Fixes for the tests in pseudo-stack mode.
[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 #define CELT_C 
9 #include "../libcelt/stack_alloc.h"
10 #include "../libcelt/kiss_fft.c"
11
12
13 #ifndef M_PI
14 #define M_PI 3.141592653
15 #endif
16
17 int ret = 0;
18
19 void check(kiss_fft_cpx  * in,kiss_fft_cpx  * out,int nfft,int isinverse)
20 {
21     int bin,k;
22     double errpow=0,sigpow=0, snr;
23     
24     for (bin=0;bin<nfft;++bin) {
25         double ansr = 0;
26         double ansi = 0;
27         double difr;
28         double difi;
29
30         for (k=0;k<nfft;++k) {
31             double phase = -2*M_PI*bin*k/nfft;
32             double re = cos(phase);
33             double im = sin(phase);
34             if (isinverse)
35                 im = -im;
36
37             if (!isinverse)
38             {
39                re /= nfft;
40                im /= nfft;
41             }
42
43             ansr += in[k].r * re - in[k].i * im;
44             ansi += in[k].r * im + in[k].i * re;
45         }
46         /*printf ("%d %d ", (int)ansr, (int)ansi);*/
47         difr = ansr - out[bin].r;
48         difi = ansi - out[bin].i;
49         errpow += difr*difr + difi*difi;
50         sigpow += ansr*ansr+ansi*ansi;
51     }
52     snr = 10*log10(sigpow/errpow);
53     printf("nfft=%d inverse=%d,snr = %f\n",nfft,isinverse,snr );
54     if (snr<60) {
55        printf( "** poor snr: %f ** \n", snr);
56        ret = 1;
57     }
58 }
59
60 void test1d(int nfft,int isinverse)
61 {
62     size_t buflen = sizeof(kiss_fft_cpx)*nfft;
63
64     kiss_fft_cpx  * in = (kiss_fft_cpx*)malloc(buflen);
65     kiss_fft_cpx  * out= (kiss_fft_cpx*)malloc(buflen);
66     kiss_fft_cfg  cfg = kiss_fft_alloc(nfft,0,0);
67     int k;
68
69     for (k=0;k<nfft;++k) {
70         in[k].r = (rand() % 32767) - 16384;
71         in[k].i = (rand() % 32767) - 16384;
72     }
73
74 #ifdef DOUBLE_PRECISION
75     for (k=0;k<nfft;++k) {
76        in[k].r *= 32768;
77        in[k].i *= 32768;
78     }
79 #endif
80     
81     if (isinverse)
82     {
83        for (k=0;k<nfft;++k) {
84           in[k].r /= nfft;
85           in[k].i /= nfft;
86        }
87     }
88     
89     /*for (k=0;k<nfft;++k) printf("%d %d ", in[k].r, in[k].i);printf("\n");*/
90        
91     if (isinverse)
92        kiss_ifft(cfg,in,out);
93     else
94        kiss_fft(cfg,in,out);
95
96     /*for (k=0;k<nfft;++k) printf("%d %d ", out[k].r, out[k].i);printf("\n");*/
97
98     check(in,out,nfft,isinverse);
99
100     free(in);
101     free(out);
102     free(cfg);
103 }
104
105 int main(int argc,char ** argv)
106 {
107     ALLOC_STACK;
108     if (argc>1) {
109         int k;
110         for (k=1;k<argc;++k) {
111             test1d(atoi(argv[k]),0);
112             test1d(atoi(argv[k]),1);
113         }
114     }else{
115         test1d(32,0);
116         test1d(32,1);
117         test1d(128,0);
118         test1d(128,1);
119         test1d(256,0);
120         test1d(256,1);
121 #ifndef RADIX_TWO_ONLY
122         test1d(36,0);
123         test1d(36,1);
124         test1d(50,0);
125         test1d(50,1);
126         test1d(120,0);
127         test1d(120,1);
128         test1d(105,0);
129         test1d(105,1);
130 #endif
131     }
132     return ret;
133 }