fixed-point: exp_rotation() mostly converted (still need to convert the cos/sin)
authorJean-Marc Valin <Jean-Marc.Valin@csiro.au>
Wed, 27 Feb 2008 02:49:54 +0000 (13:49 +1100)
committerJean-Marc Valin <Jean-Marc.Valin@csiro.au>
Wed, 27 Feb 2008 02:49:54 +0000 (13:49 +1100)
libcelt/bands.c
libcelt/bands.h
tests/Makefile.am
tests/rotation-test.c [new file with mode: 0644]

index 916dfc4..200af28 100644 (file)
 #include "vq.h"
 #include "cwrs.h"
 
-/** Applies a series of rotations so that pulses are spread like a two-sided
-exponential. The effect of this is to reduce the tonal noise created by the
-sparse spectrum resulting from the pulse codebook */
-static void exp_rotation(celt_norm_t *X, int len, float theta, int dir, int stride, int iter)
+
+void exp_rotation(celt_norm_t *X, int len, float theta, int dir, int stride, int iter)
 {
    int i, k;
-   float c, s;
-   c = cos(theta);
-   s = dir*sin(theta);
+   celt_word16_t c, s;
+   c = Q15ONE*cos(theta);
+   s = dir*Q15ONE*sin(theta);
    for (k=0;k<iter;k++)
    {
+      /* We could use MULT16_16_P15 instead of MULT16_16_Q15 for more accuracy, 
+         but at this point, I really don't think it's necessary */
       for (i=0;i<len-stride;i++)
       {
-         float x1, x2;
+         celt_norm_t x1, x2;
          x1 = X[i];
          x2 = X[i+stride];
-         X[i] = c*x1 - s*x2;
-         X[i+stride] = c*x2 + s*x1;
+         X[i] = MULT16_16_Q15(c,x1) - MULT16_16_Q15(s,x2);
+         X[i+stride] = MULT16_16_Q15(c,x2) + MULT16_16_Q15(s,x1);
       }
       for (i=len-2*stride-1;i>=0;i--)
       {
-         float x1, x2;
+         celt_norm_t x1, x2;
          x1 = X[i];
          x2 = X[i+stride];
-         X[i] = c*x1 - s*x2;
-         X[i+stride] = c*x2 + s*x1;
+         X[i] = MULT16_16_Q15(c,x1) - MULT16_16_Q15(s,x2);
+         X[i+stride] = MULT16_16_Q15(c,x2) + MULT16_16_Q15(s,x1);
       }
    }
 }
index 6890f77..8e394d3 100644 (file)
 #include "entdec.h"
 #include "rate.h"
 
+/** Applies a series of rotations so that pulses are spread like a two-sided
+exponential. The effect of this is to reduce the tonal noise created by the
+sparse spectrum resulting from the pulse codebook */
+void exp_rotation(celt_norm_t *X, int len, float theta, int dir, int stride, int iter);
+
 /** Compute the amplitude (sqrt energy) in each of the bands 
  * @param m Mode data 
  * @param X Spectrum
index 79795ca..e19c6c6 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 mdct-test
+TESTS = type-test ectest cwrs32-test cwrs64-test real-fft-test dft-test laplace-test mdct-test rotation-test
 
-noinst_PROGRAMS = 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 mdct-test rotation-test
 
 type_test_SOURCES = type-test.c
 ectest_SOURCES = ectest.c
@@ -13,5 +13,6 @@ 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
+rotation_test_SOURCES = rotation-test.c
 
 LDADD = $(top_builddir)/libcelt/libcelt.la
diff --git a/tests/rotation-test.c b/tests/rotation-test.c
new file mode 100644 (file)
index 0000000..b73762c
--- /dev/null
@@ -0,0 +1,43 @@
+#ifdef HAVE_CONFIG_H
+#include "config.h"
+#endif
+
+#include <stdio.h>
+#include <stdlib.h>
+#include "celt_types.h"
+#include "bands.h"
+#include <math.h>
+#define MAX_SIZE 100
+
+int ret=0;
+void test_rotation(int N, int K)
+{
+   int i;
+   double err = 0, ener = 0, snr;
+   float theta = .007*N/K;
+   celt_word16_t x0[MAX_SIZE];
+   celt_word16_t x1[MAX_SIZE];
+   for (i=0;i<N;i++)
+      x1[i] = x0[i] = rand()%32767-16384;
+   exp_rotation(x1, N, theta, 1, 1, 8);
+   exp_rotation(x1, N, theta, -1, 1, 8);
+   for (i=0;i<N;i++)
+   {
+      err += (x0[i]-(double)x1[i])*(x0[i]-(double)x1[i]);
+      ener += x0[i]*(double)x0[i];
+   }
+   snr = 20*log10(ener/err);
+   printf ("SNR for size %d (%d pulses) is %f\n", N, K, snr);
+   if (snr < 60)
+      ret = 1;
+}
+
+int main()
+{
+   test_rotation(4, 40);
+   test_rotation(7, 20);
+   test_rotation(10, 10);
+   test_rotation(23, 5);
+   test_rotation(50, 3);
+   return ret;
+}