c@409: #include c@409: #include c@409: #include c@409: #include "kiss_fft.h" c@409: #include "kiss_fftr.h" c@409: #include c@409: c@409: c@409: static c@409: double two_tone_test( int nfft, int bin1,int bin2) c@409: { c@409: kiss_fftr_cfg cfg = NULL; c@409: kiss_fft_cpx *kout = NULL; c@409: kiss_fft_scalar *tbuf = NULL; c@409: c@409: int i; c@409: double f1 = bin1*2*M_PI/nfft; c@409: double f2 = bin2*2*M_PI/nfft; c@409: double sigpow=0; c@409: double noisepow=0; c@409: #if FIXED_POINT==32 c@409: long maxrange = LONG_MAX; c@409: #else c@409: long maxrange = SHRT_MAX;/* works fine for float too*/ c@409: #endif c@409: c@409: cfg = kiss_fftr_alloc(nfft , 0, NULL, NULL); c@409: tbuf = KISS_FFT_MALLOC(nfft * sizeof(kiss_fft_scalar)); c@409: kout = KISS_FFT_MALLOC(nfft * sizeof(kiss_fft_cpx)); c@409: c@409: /* generate a signal with two tones*/ c@409: for (i = 0; i < nfft; i++) { c@409: #ifdef USE_SIMD c@409: tbuf[i] = _mm_set1_ps( (maxrange>>1)*cos(f1*i) c@409: + (maxrange>>1)*cos(f2*i) ); c@409: #else c@409: tbuf[i] = (maxrange>>1)*cos(f1*i) c@409: + (maxrange>>1)*cos(f2*i); c@409: #endif c@409: } c@409: c@409: kiss_fftr(cfg, tbuf, kout); c@409: c@409: for (i=0;i < (nfft/2+1);++i) { c@409: #ifdef USE_SIMD c@409: double tmpr = (double)*(float*)&kout[i].r / (double)maxrange; c@409: double tmpi = (double)*(float*)&kout[i].i / (double)maxrange; c@409: #else c@409: double tmpr = (double)kout[i].r / (double)maxrange; c@409: double tmpi = (double)kout[i].i / (double)maxrange; c@409: #endif c@409: double mag2 = tmpr*tmpr + tmpi*tmpi; c@409: if (i!=0 && i!= nfft/2) c@409: mag2 *= 2; /* all bins except DC and Nyquist have symmetric counterparts implied*/ c@409: c@409: /* if there is power in one of the expected bins, it is signal, otherwise noise*/ c@409: if ( i!=bin1 && i != bin2 ) c@409: noisepow += mag2; c@409: else c@409: sigpow += mag2; c@409: } c@409: kiss_fft_cleanup(); c@409: /*printf("TEST %d,%d,%d noise @ %fdB\n",nfft,bin1,bin2,10*log10(noisepow/sigpow +1e-30) );*/ c@409: return 10*log10(sigpow/(noisepow+1e-50) ); c@409: } c@409: c@409: int main(int argc,char ** argv) c@409: { c@409: int nfft = 4*2*2*3*5; c@409: if (argc>1) nfft = atoi(argv[1]); c@409: c@409: int i,j; c@409: double minsnr = 500; c@409: double maxsnr = -500; c@409: double snr; c@409: for (i=0;i>4)+1) { c@409: for (j=i;j>4)+7) { c@409: snr = two_tone_test(nfft,i,j); c@409: if (snrmaxsnr) { c@409: maxsnr=snr; c@409: } c@409: } c@409: } c@409: snr = two_tone_test(nfft,nfft/2,nfft/2); c@409: if (snrmaxsnr) maxsnr=snr; c@409: c@409: printf("TwoToneTest: snr ranges from %ddB to %ddB\n",(int)minsnr,(int)maxsnr); c@409: printf("sizeof(kiss_fft_scalar) = %d\n",(int)sizeof(kiss_fft_scalar) ); c@409: return 0; c@409: }