Mercurial > hg > silvet
annotate constant-q-cpp/src/ext/kissfft/test/test_vs_dft.c @ 366:5d0a2ebb4d17
Bring dependent libraries in to repo
author | Chris Cannam |
---|---|
date | Fri, 24 Jun 2016 14:47:45 +0100 |
parents | |
children |
rev | line source |
---|---|
Chris@366 | 1 #include "kiss_fft.h" |
Chris@366 | 2 |
Chris@366 | 3 |
Chris@366 | 4 void check(kiss_fft_cpx * in,kiss_fft_cpx * out,int nfft,int isinverse) |
Chris@366 | 5 { |
Chris@366 | 6 int bin,k; |
Chris@366 | 7 double errpow=0,sigpow=0; |
Chris@366 | 8 |
Chris@366 | 9 for (bin=0;bin<nfft;++bin) { |
Chris@366 | 10 double ansr = 0; |
Chris@366 | 11 double ansi = 0; |
Chris@366 | 12 double difr; |
Chris@366 | 13 double difi; |
Chris@366 | 14 |
Chris@366 | 15 for (k=0;k<nfft;++k) { |
Chris@366 | 16 double phase = -2*M_PI*bin*k/nfft; |
Chris@366 | 17 double re = cos(phase); |
Chris@366 | 18 double im = sin(phase); |
Chris@366 | 19 if (isinverse) |
Chris@366 | 20 im = -im; |
Chris@366 | 21 |
Chris@366 | 22 #ifdef FIXED_POINT |
Chris@366 | 23 re /= nfft; |
Chris@366 | 24 im /= nfft; |
Chris@366 | 25 #endif |
Chris@366 | 26 |
Chris@366 | 27 ansr += in[k].r * re - in[k].i * im; |
Chris@366 | 28 ansi += in[k].r * im + in[k].i * re; |
Chris@366 | 29 } |
Chris@366 | 30 difr = ansr - out[bin].r; |
Chris@366 | 31 difi = ansi - out[bin].i; |
Chris@366 | 32 errpow += difr*difr + difi*difi; |
Chris@366 | 33 sigpow += ansr*ansr+ansi*ansi; |
Chris@366 | 34 } |
Chris@366 | 35 printf("nfft=%d inverse=%d,snr = %f\n",nfft,isinverse,10*log10(sigpow/errpow) ); |
Chris@366 | 36 } |
Chris@366 | 37 |
Chris@366 | 38 void test1d(int nfft,int isinverse) |
Chris@366 | 39 { |
Chris@366 | 40 size_t buflen = sizeof(kiss_fft_cpx)*nfft; |
Chris@366 | 41 |
Chris@366 | 42 kiss_fft_cpx * in = (kiss_fft_cpx*)malloc(buflen); |
Chris@366 | 43 kiss_fft_cpx * out= (kiss_fft_cpx*)malloc(buflen); |
Chris@366 | 44 kiss_fft_cfg cfg = kiss_fft_alloc(nfft,isinverse,0,0); |
Chris@366 | 45 int k; |
Chris@366 | 46 |
Chris@366 | 47 for (k=0;k<nfft;++k) { |
Chris@366 | 48 in[k].r = (rand() % 65536) - 32768; |
Chris@366 | 49 in[k].i = (rand() % 65536) - 32768; |
Chris@366 | 50 } |
Chris@366 | 51 |
Chris@366 | 52 kiss_fft(cfg,in,out); |
Chris@366 | 53 |
Chris@366 | 54 check(in,out,nfft,isinverse); |
Chris@366 | 55 |
Chris@366 | 56 free(in); |
Chris@366 | 57 free(out); |
Chris@366 | 58 free(cfg); |
Chris@366 | 59 } |
Chris@366 | 60 |
Chris@366 | 61 int main(int argc,char ** argv) |
Chris@366 | 62 { |
Chris@366 | 63 if (argc>1) { |
Chris@366 | 64 int k; |
Chris@366 | 65 for (k=1;k<argc;++k) { |
Chris@366 | 66 test1d(atoi(argv[k]),0); |
Chris@366 | 67 test1d(atoi(argv[k]),1); |
Chris@366 | 68 } |
Chris@366 | 69 }else{ |
Chris@366 | 70 test1d(32,0); |
Chris@366 | 71 test1d(32,1); |
Chris@366 | 72 } |
Chris@366 | 73 return 0; |
Chris@366 | 74 } |