testcase for the MDCT and IMDCT
authorJean-Marc Valin <Jean-Marc.Valin@csiro.au>
Fri, 22 Feb 2008 07:02:56 +0000 (18:02 +1100)
committerJean-Marc Valin <Jean-Marc.Valin@csiro.au>
Fri, 22 Feb 2008 07:02:56 +0000 (18:02 +1100)
tests/Makefile.am
tests/mdct-test.c [new file with mode: 0644]

index 817a7d7..79795ca 100644 (file)
@@ -1,9 +1,9 @@
 INCLUDES = -I$(top_srcdir)/libcelt
 METASOURCES = AUTO
 
-TESTS = type-test ectest cwrs32-test cwrs64-test real-fft-test dft-test laplace-test
+TESTS = type-test ectest cwrs32-test cwrs64-test real-fft-test dft-test laplace-test mdct-test
 
-noinst_PROGRAMS = type-test ectest cwrs32-test cwrs64-test real-fft-test dft-test laplace-test
+noinst_PROGRAMS = type-test ectest cwrs32-test cwrs64-test real-fft-test dft-test laplace-test mdct-test
 
 type_test_SOURCES = type-test.c
 ectest_SOURCES = ectest.c
@@ -12,5 +12,6 @@ cwrs64_test_SOURCES = cwrs64-test.c
 real_fft_test_SOURCES = real-fft-test.c
 dft_test_SOURCES = dft-test.c
 laplace_test_SOURCES = laplace-test.c
+mdct_test_SOURCES = mdct-test.c
 
 LDADD = $(top_builddir)/libcelt/libcelt.la
diff --git a/tests/mdct-test.c b/tests/mdct-test.c
new file mode 100644 (file)
index 0000000..909e58e
--- /dev/null
@@ -0,0 +1,124 @@
+#ifdef HAVE_CONFIG_H
+#include "config.h"
+#endif
+
+#include <stdio.h>
+#include "mdct.h"
+
+
+void check(float  * in,float  * out,int nfft,int isinverse)
+{
+    int bin,k;
+    double errpow=0,sigpow=0;
+    
+    for (bin=0;bin<nfft/2;++bin) {
+        double ansr = 0;
+        double difr;
+
+        for (k=0;k<nfft;++k) {
+           double phase = 2*M_PI*(k+.5+.25*nfft)*(bin+.5)/nfft;
+           double re = cos(phase);
+            
+           re /= nfft/4;
+
+           ansr += in[k] * re;
+        }
+        /*printf ("%f %f\n", ansr, out[bin]);*/
+        difr = ansr - out[bin];
+        errpow += difr*difr;
+        sigpow += ansr*ansr;
+    }
+    printf("nfft=%d inverse=%d,snr = %f\n",nfft,isinverse,10*log10(sigpow/errpow) );
+}
+
+void check_inv(float  * in,float  * out,int nfft,int isinverse)
+{
+   int bin,k;
+   double errpow=0,sigpow=0;
+    
+   for (bin=0;bin<nfft;++bin) {
+      double ansr = 0;
+      double difr;
+
+      for (k=0;k<nfft/2;++k) {
+         double phase = 2*M_PI*(bin+.5+.25*nfft)*(k+.5)/nfft;
+         double re = cos(phase);
+            
+         ansr += in[k] * re;
+      }
+      /*printf ("%f %f\n", ansr, out[bin]);*/
+      difr = ansr - out[bin];
+      errpow += difr*difr;
+      sigpow += ansr*ansr;
+   }
+   printf("nfft=%d inverse=%d,snr = %f\n",nfft,isinverse,10*log10(sigpow/errpow) );
+}
+
+
+void test1d(int nfft,int isinverse)
+{
+    mdct_lookup cfg;
+    size_t buflen = sizeof(float)*nfft;
+
+    float  * in = (float*)malloc(buflen);
+    float  * out= (float*)malloc(buflen);
+    int k;
+
+    mdct_init(&cfg, nfft);
+    for (k=0;k<nfft;++k) {
+        in[k] = (rand() % 65536) - 32768;
+    }
+
+#ifdef DOUBLE_PRECISION
+    for (k=0;k<nfft;++k) {
+       in[k].r *= 65536;
+       in[k].i *= 65536;
+    }
+#endif
+    
+    if (isinverse)
+    {
+       for (k=0;k<nfft;++k) {
+          in[k] /= nfft;
+       }
+    }
+    
+    /*for (k=0;k<nfft;++k) printf("%d %d ", in[k].r, in[k].i);printf("\n");*/
+       
+    if (isinverse)
+    {
+       mdct_backward(&cfg,in,out);
+       check_inv(in,out,nfft,isinverse);
+    } else {
+       mdct_forward(&cfg,in,out);
+       check(in,out,nfft,isinverse);
+    }
+    /*for (k=0;k<nfft;++k) printf("%d %d ", out[k].r, out[k].i);printf("\n");*/
+
+
+    free(in);
+    free(out);
+    mdct_clear(&cfg);
+}
+
+int main(int argc,char ** argv)
+{
+    if (argc>1) {
+        int k;
+        for (k=1;k<argc;++k) {
+            test1d(atoi(argv[k]),0);
+            test1d(atoi(argv[k]),1);
+        }
+    }else{
+        test1d(32,0);
+        test1d(32,1);
+        exit(0);
+        test1d(36,0);
+        test1d(36,1);
+        test1d(50,0);
+        test1d(50,1);
+        test1d(120,0);
+        test1d(120,1);
+    }
+    return 0;
+}