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