|
c@174
|
1 #include <stdlib.h>
|
|
c@174
|
2 #include <string.h>
|
|
c@174
|
3 #include <stdio.h>
|
|
c@174
|
4 #include "kiss_fft.h"
|
|
c@174
|
5 #include "kiss_fftr.h"
|
|
c@174
|
6 #include <limits.h>
|
|
c@174
|
7
|
|
c@174
|
8
|
|
c@174
|
9 static
|
|
c@174
|
10 double two_tone_test( int nfft, int bin1,int bin2)
|
|
c@174
|
11 {
|
|
c@174
|
12 kiss_fftr_cfg cfg = NULL;
|
|
c@174
|
13 kiss_fft_cpx *kout = NULL;
|
|
c@174
|
14 kiss_fft_scalar *tbuf = NULL;
|
|
c@174
|
15
|
|
c@174
|
16 int i;
|
|
c@174
|
17 double f1 = bin1*2*M_PI/nfft;
|
|
c@174
|
18 double f2 = bin2*2*M_PI/nfft;
|
|
c@174
|
19 double sigpow=0;
|
|
c@174
|
20 double noisepow=0;
|
|
c@174
|
21 #if FIXED_POINT==32
|
|
c@174
|
22 long maxrange = LONG_MAX;
|
|
c@174
|
23 #else
|
|
c@174
|
24 long maxrange = SHRT_MAX;/* works fine for float too*/
|
|
c@174
|
25 #endif
|
|
c@174
|
26
|
|
c@174
|
27 cfg = kiss_fftr_alloc(nfft , 0, NULL, NULL);
|
|
c@174
|
28 tbuf = KISS_FFT_MALLOC(nfft * sizeof(kiss_fft_scalar));
|
|
c@174
|
29 kout = KISS_FFT_MALLOC(nfft * sizeof(kiss_fft_cpx));
|
|
c@174
|
30
|
|
c@174
|
31 /* generate a signal with two tones*/
|
|
c@174
|
32 for (i = 0; i < nfft; i++) {
|
|
c@174
|
33 #ifdef USE_SIMD
|
|
c@174
|
34 tbuf[i] = _mm_set1_ps( (maxrange>>1)*cos(f1*i)
|
|
c@174
|
35 + (maxrange>>1)*cos(f2*i) );
|
|
c@174
|
36 #else
|
|
c@174
|
37 tbuf[i] = (maxrange>>1)*cos(f1*i)
|
|
c@174
|
38 + (maxrange>>1)*cos(f2*i);
|
|
c@174
|
39 #endif
|
|
c@174
|
40 }
|
|
c@174
|
41
|
|
c@174
|
42 kiss_fftr(cfg, tbuf, kout);
|
|
c@174
|
43
|
|
c@174
|
44 for (i=0;i < (nfft/2+1);++i) {
|
|
c@174
|
45 #ifdef USE_SIMD
|
|
c@174
|
46 double tmpr = (double)*(float*)&kout[i].r / (double)maxrange;
|
|
c@174
|
47 double tmpi = (double)*(float*)&kout[i].i / (double)maxrange;
|
|
c@174
|
48 #else
|
|
c@174
|
49 double tmpr = (double)kout[i].r / (double)maxrange;
|
|
c@174
|
50 double tmpi = (double)kout[i].i / (double)maxrange;
|
|
c@174
|
51 #endif
|
|
c@174
|
52 double mag2 = tmpr*tmpr + tmpi*tmpi;
|
|
c@174
|
53 if (i!=0 && i!= nfft/2)
|
|
c@174
|
54 mag2 *= 2; /* all bins except DC and Nyquist have symmetric counterparts implied*/
|
|
c@174
|
55
|
|
c@174
|
56 /* if there is power in one of the expected bins, it is signal, otherwise noise*/
|
|
c@174
|
57 if ( i!=bin1 && i != bin2 )
|
|
c@174
|
58 noisepow += mag2;
|
|
c@174
|
59 else
|
|
c@174
|
60 sigpow += mag2;
|
|
c@174
|
61 }
|
|
c@174
|
62 kiss_fft_cleanup();
|
|
c@174
|
63 /*printf("TEST %d,%d,%d noise @ %fdB\n",nfft,bin1,bin2,10*log10(noisepow/sigpow +1e-30) );*/
|
|
c@174
|
64 return 10*log10(sigpow/(noisepow+1e-50) );
|
|
c@174
|
65 }
|
|
c@174
|
66
|
|
c@174
|
67 int main(int argc,char ** argv)
|
|
c@174
|
68 {
|
|
c@174
|
69 int nfft = 4*2*2*3*5;
|
|
c@174
|
70 if (argc>1) nfft = atoi(argv[1]);
|
|
c@174
|
71
|
|
c@174
|
72 int i,j;
|
|
c@174
|
73 double minsnr = 500;
|
|
c@174
|
74 double maxsnr = -500;
|
|
c@174
|
75 double snr;
|
|
c@174
|
76 for (i=0;i<nfft/2;i+= (nfft>>4)+1) {
|
|
c@174
|
77 for (j=i;j<nfft/2;j+=(nfft>>4)+7) {
|
|
c@174
|
78 snr = two_tone_test(nfft,i,j);
|
|
c@174
|
79 if (snr<minsnr) {
|
|
c@174
|
80 minsnr=snr;
|
|
c@174
|
81 }
|
|
c@174
|
82 if (snr>maxsnr) {
|
|
c@174
|
83 maxsnr=snr;
|
|
c@174
|
84 }
|
|
c@174
|
85 }
|
|
c@174
|
86 }
|
|
c@174
|
87 snr = two_tone_test(nfft,nfft/2,nfft/2);
|
|
c@174
|
88 if (snr<minsnr) minsnr=snr;
|
|
c@174
|
89 if (snr>maxsnr) maxsnr=snr;
|
|
c@174
|
90
|
|
c@174
|
91 printf("TwoToneTest: snr ranges from %ddB to %ddB\n",(int)minsnr,(int)maxsnr);
|
|
c@174
|
92 printf("sizeof(kiss_fft_scalar) = %d\n",(int)sizeof(kiss_fft_scalar) );
|
|
c@174
|
93 return 0;
|
|
c@174
|
94 }
|
