New VQ search is now enabled by default after fixing the last remaining issues:
[opus.git] / tests / mdct-test.c
1 #ifdef HAVE_CONFIG_H
2 #include "config.h"
3 #endif
4
5 #include <stdio.h>
6 #include "mdct.h"
7 #include "stack_alloc.h"
8
9 #ifndef M_PI
10 #define M_PI 3.141592653
11 #endif
12
13 int ret = 0;
14 void check(kiss_fft_scalar  * in,kiss_fft_scalar  * out,int nfft,int isinverse)
15 {
16     int bin,k;
17     double errpow=0,sigpow=0;
18     double snr;
19     for (bin=0;bin<nfft/2;++bin) {
20         double ansr = 0;
21         double difr;
22
23         for (k=0;k<nfft;++k) {
24            double phase = 2*M_PI*(k+.5+.25*nfft)*(bin+.5)/nfft;
25            double re = cos(phase);
26             
27            re /= nfft/4;
28
29            ansr += in[k] * re;
30         }
31         /*printf ("%f %f\n", ansr, out[bin]);*/
32         difr = ansr - out[bin];
33         errpow += difr*difr;
34         sigpow += ansr*ansr;
35     }
36     snr = 10*log10(sigpow/errpow);
37     printf("nfft=%d inverse=%d,snr = %f\n",nfft,isinverse,snr );
38     if (snr<60) {
39        printf( "** poor snr: %f **\n", snr);
40        ret = 1;
41     }
42 }
43
44 void check_inv(kiss_fft_scalar  * in,kiss_fft_scalar  * out,int nfft,int isinverse)
45 {
46    int bin,k;
47    double errpow=0,sigpow=0;
48    double snr;
49    for (bin=0;bin<nfft;++bin) {
50       double ansr = 0;
51       double difr;
52
53       for (k=0;k<nfft/2;++k) {
54          double phase = 2*M_PI*(bin+.5+.25*nfft)*(k+.5)/nfft;
55          double re = cos(phase);
56
57          //re *= 2;
58
59          ansr += in[k] * re;
60       }
61       /*printf ("%f %f\n", ansr, out[bin]);*/
62       difr = ansr - out[bin];
63       errpow += difr*difr;
64       sigpow += ansr*ansr;
65    }
66    snr = 10*log10(sigpow/errpow);
67    printf("nfft=%d inverse=%d,snr = %f\n",nfft,isinverse,snr );
68    if (snr<60) {
69       printf( "** poor snr: %f **\n", snr);
70       ret = 1;
71    }
72 }
73
74
75 void test1d(int nfft,int isinverse)
76 {
77     mdct_lookup cfg;
78     size_t buflen = sizeof(kiss_fft_scalar)*nfft;
79
80     kiss_fft_scalar  * in = (kiss_fft_scalar*)malloc(buflen);
81     kiss_fft_scalar  * out= (kiss_fft_scalar*)malloc(buflen);
82     celt_word16_t  * window= (celt_word16_t*)malloc(sizeof(celt_word16_t)*nfft/2);
83     int k;
84
85     mdct_init(&cfg, nfft);
86     for (k=0;k<nfft;++k) {
87         in[k] = (rand() % 32768) - 16384;
88     }
89
90     for (k=0;k<nfft/2;++k) {
91        window[k] = Q15ONE;
92     }
93 #ifdef DOUBLE_PRECISION
94     for (k=0;k<nfft;++k) {
95        in[k] *= 32768;
96     }
97 #endif
98     
99     if (isinverse)
100     {
101        for (k=0;k<nfft;++k) {
102           in[k] /= nfft;
103        }
104     }
105     
106     /*for (k=0;k<nfft;++k) printf("%d %d ", in[k].r, in[k].i);printf("\n");*/
107        
108     if (isinverse)
109     {
110        for (k=0;k<nfft;++k)
111           out[k] = 0;
112        mdct_backward(&cfg,in,out, window, nfft/2);
113        check_inv(in,out,nfft,isinverse);
114     } else {
115        mdct_forward(&cfg,in,out,window, nfft/2);
116        check(in,out,nfft,isinverse);
117     }
118     /*for (k=0;k<nfft;++k) printf("%d %d ", out[k].r, out[k].i);printf("\n");*/
119
120
121     free(in);
122     free(out);
123     mdct_clear(&cfg);
124 }
125
126 int main(int argc,char ** argv)
127 {
128     ALLOC_STACK;
129     if (argc>1) {
130         int k;
131         for (k=1;k<argc;++k) {
132             test1d(atoi(argv[k]),0);
133             test1d(atoi(argv[k]),1);
134         }
135     }else{
136         test1d(32,0);
137         test1d(32,1);
138         test1d(256,0);
139         test1d(256,1);
140         test1d(512,0);
141         test1d(512,1);
142 #ifndef RADIX_TWO_ONLY
143         test1d(40,0);
144         test1d(40,1);
145         test1d(56,0);
146         test1d(56,1);
147         test1d(120,0);
148         test1d(120,1);
149         test1d(240,0);
150         test1d(240,1);
151         test1d(480,0);
152         test1d(480,1);
153 #endif
154     }
155     return ret;
156 }