annotate src/ext/kissfft/test/test_vs_dft.c @ 196:da283326bcd3 tip master

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