[xiph-cvs] cvs commit: speex/libspeex filters.c lsp.c math_approx.c math_approx.h nb_celp.c

[Jean-Marc Valin]

jm          03/11/07 03:34:15

  Modified:    libspeex filters.c lsp.c math_approx.c math_approx.h
                        nb_celp.c
  Log:
  fixed-point: acos function approximated with fixed-point arithmetic

Revision  Changes    Path
1.53      +1 -1      speex/libspeex/filters.c

Index: filters.c
===================================================================
RCS file: /usr/local/cvsroot/speex/libspeex/filters.c,v
retrieving revision 1.52
retrieving revision 1.53
diff -u -r1.52 -r1.53
--- filters.c	6 Nov 2003 08:41:56 -0000	1.52
+++ filters.c	7 Nov 2003 08:34:14 -0000	1.53
@@ -118,7 +118,7 @@
    }
    
    /*FIXME: remove division*/
-   return SHR(SHL((spx_word32_t)sqroot(1+sum/len),(sig_shift+3)),SIG_SHIFT);
+   return SHR(SHL((spx_word32_t)spx_sqrt(1+sum/len),(sig_shift+3)),SIG_SHIFT);
 }
 
 

<p><p>1.39      +3 -3      speex/libspeex/lsp.c

Index: lsp.c
===================================================================
RCS file: /usr/local/cvsroot/speex/libspeex/lsp.c,v
retrieving revision 1.38
retrieving revision 1.39
diff -u -r1.38 -r1.39
--- lsp.c	2 Nov 2003 07:44:35 -0000	1.38
+++ lsp.c	7 Nov 2003 08:34:14 -0000	1.39
@@ -46,7 +46,7 @@
 #include <math.h>
 #include "lsp.h"
 #include "stack_alloc.h"
-
+#include "math_approx.h"
 
 #ifndef M_PI
 #define M_PI           3.14159265358979323846  /* pi */
@@ -78,8 +78,8 @@
 
 #define FREQ_SCALE 16384
 #define ANGLE2X(a) (SHL(cos_32(a),2))
-#define X2ANGLE(x) (acos(.00006103515625*(x))*LSP_SCALING)
-
+//#define X2ANGLE(x) (acos(.00006103515625*(x))*LSP_SCALING)
+#define X2ANGLE(x) (spx_acos(x))
 /* maybe we could approximate acos as 
    sqrt(6.7349563814-5.213449731*sqrt(0.6688572663+x))
 */

<p><p>1.7       +37 -2     speex/libspeex/math_approx.c

Index: math_approx.c
===================================================================
RCS file: /usr/local/cvsroot/speex/libspeex/math_approx.c,v
retrieving revision 1.6
retrieving revision 1.7
diff -u -r1.6 -r1.7
--- math_approx.c	6 Nov 2003 09:14:53 -0000	1.6
+++ math_approx.c	7 Nov 2003 08:34:14 -0000	1.7
@@ -41,7 +41,7 @@
 #define C2 -12627
 #define C3 4215
 
-spx_word16_t sqroot(spx_word32_t x)
+spx_word16_t spx_sqrt(spx_word32_t x)
 {
    int k=0;
    spx_word32_t rt;
@@ -79,12 +79,47 @@
       k++;
       }
 #endif
+   while (x<4096)
+   {
+      x<<=2;
+      k--;
+   }
    rt = ADD16(C0, MULT16_16_Q14(x, ADD16(C1, MULT16_16_Q14(x, ADD16(C2, MULT16_16_Q14(x, (C3)))))));
-   rt <<= k;
+   if (k>0)
+      rt <<= k;
+   else
+      rt >>= -k;
    rt >>=7;
    return rt;
 }
 
 /* log(x) ~= -2.18151 + 4.20592*x - 2.88938*x^2 + 0.86535*x^3 (for .5 < x < 1) */
 
+
+#define A1 16469
+#define A2 2242
+#define A3 1486
+
+spx_word16_t spx_acos(spx_word16_t x)
+{
+   int s=0;
+   spx_word16_t ret;
+   spx_word32_t sq;
+   if (x<0)
+   {
+      s=1;
+      x = -x;
+   }
+   x = 16384-x;
+   
+   x = x >> 1;
+   sq = MULT16_16_Q13(x, ADD16(A1, MULT16_16_Q13(x, ADD16(A2, MULT16_16_Q13(x, (A3))))));
+   ret = spx_sqrt(SHL(sq,13));
+   
+   /*ret = spx_sqrt(67108864*(-1.6129e-04 + 2.0104e+00*f + 2.7373e-01*f*f + 1.8136e-01*f*f*f));*/
+   if (s)
+      ret = 25736-ret;
+   return ret;
+}
+
 #endif

<p><p>1.4       +4 -4      speex/libspeex/math_approx.h

Index: math_approx.h
===================================================================
RCS file: /usr/local/cvsroot/speex/libspeex/math_approx.h,v
retrieving revision 1.3
retrieving revision 1.4
diff -u -r1.3 -r1.4
--- math_approx.h	6 Nov 2003 08:41:56 -0000	1.3
+++ math_approx.h	7 Nov 2003 08:34:14 -0000	1.4
@@ -35,12 +35,12 @@
 
 #include "misc.h"
 
-float speex_cos(float x);
-
 #ifdef FIXED_POINT
-spx_word16_t sqroot(spx_word32_t x);
+spx_word16_t spx_sqrt(spx_word32_t x);
+spx_word16_t spx_acos(spx_word16_t x);
 #else
-#define sqroot sqrt
+#define spx_sqrt sqrt
+#define spx_acos acos
 #endif
 
 #endif

<p><p>1.141     +1 -1      speex/libspeex/nb_celp.c

Index: nb_celp.c
===================================================================
RCS file: /usr/local/cvsroot/speex/libspeex/nb_celp.c,v
retrieving revision 1.140
retrieving revision 1.141
diff -u -r1.140 -r1.141
--- nb_celp.c	9 Oct 2003 20:53:41 -0000	1.140
+++ nb_celp.c	7 Nov 2003 08:34:14 -0000	1.141
@@ -877,7 +877,7 @@
    /* The next frame will not be the first (Duh!) */
    st->first = 0;
 
-   {
+   if (0) {
       float ener=0, err=0;
       float snr;
       for (i=0;i<st->frameSize;i++)

<p><p>--- >8 ----
List archives:  http://www.xiph.org/archives/
Ogg project homepage: http://www.xiph.org/ogg/
To unsubscribe from this list, send a message to 'cvs-request at xiph.org'
containing only the word 'unsubscribe' in the body.  No subject is needed.
Unsubscribe messages sent to the list will be ignored/filtered.