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