Mode doc
[opus.git] / tests / dft-test.c
1 #ifdef HAVE_CONFIG_H
2 #include "config.h"
3 #endif
4
5 #include <stdio.h>
6 #include "kiss_fft.h"
7
8
9 void check(kiss_fft_cpx  * in,kiss_fft_cpx  * out,int nfft,int isinverse)
10 {
11     int bin,k;
12     double errpow=0,sigpow=0;
13     
14     for (bin=0;bin<nfft;++bin) {
15         double ansr = 0;
16         double ansi = 0;
17         double difr;
18         double difi;
19
20         for (k=0;k<nfft;++k) {
21             double phase = -2*M_PI*bin*k/nfft;
22             double re = cos(phase);
23             double im = sin(phase);
24             if (isinverse)
25                 im = -im;
26
27 #ifdef FIXED_POINT
28             re /= nfft;
29             im /= nfft;
30 #endif            
31
32             ansr += in[k].r * re - in[k].i * im;
33             ansi += in[k].r * im + in[k].i * re;
34         }
35         difr = ansr - out[bin].r;
36         difi = ansi - out[bin].i;
37         errpow += difr*difr + difi*difi;
38         sigpow += ansr*ansr+ansi*ansi;
39     }
40     printf("nfft=%d inverse=%d,snr = %f\n",nfft,isinverse,10*log10(sigpow/errpow) );
41 }
42
43 void test1d(int nfft,int isinverse)
44 {
45     size_t buflen = sizeof(kiss_fft_cpx)*nfft;
46
47     kiss_fft_cpx  * in = (kiss_fft_cpx*)malloc(buflen);
48     kiss_fft_cpx  * out= (kiss_fft_cpx*)malloc(buflen);
49     kiss_fft_cfg  cfg = kiss_fft_alloc(nfft,0,0);
50     int k;
51
52     for (k=0;k<nfft;++k) {
53         in[k].r = (rand() % 65536) - 32768;
54         in[k].i = (rand() % 65536) - 32768;
55     }
56
57     if (isinverse)
58        kiss_ifft(cfg,in,out);
59     else
60        kiss_fft(cfg,in,out);
61
62     check(in,out,nfft,isinverse);
63
64     free(in);
65     free(out);
66     free(cfg);
67 }
68
69 int main(int argc,char ** argv)
70 {
71     if (argc>1) {
72         int k;
73         for (k=1;k<argc;++k) {
74             test1d(atoi(argv[k]),0);
75             test1d(atoi(argv[k]),1);
76         }
77     }else{
78         test1d(32,0);
79         test1d(32,1);
80         test1d(36,0);
81         test1d(36,1);
82         test1d(50,0);
83         test1d(50,1);
84         test1d(120,0);
85         test1d(120,1);
86         test1d(105,0);
87         test1d(105,1);
88     }
89     return 0;
90 }