fixed-point: interger version of sqrt function
authorjm <jm@0101bb08-14d6-0310-b084-bc0e0c8e3800>
Thu, 6 Nov 2003 08:41:56 +0000 (08:41 +0000)
committerjm <jm@0101bb08-14d6-0310-b084-bc0e0c8e3800>
Thu, 6 Nov 2003 08:41:56 +0000 (08:41 +0000)
git-svn-id: http://svn.xiph.org/trunk/speex@5537 0101bb08-14d6-0310-b084-bc0e0c8e3800

libspeex/filters.c
libspeex/math_approx.c
libspeex/math_approx.h

index 6e1c2a1..d4e8903 100644 (file)
@@ -34,6 +34,7 @@
 #include "stack_alloc.h"
 #include <math.h>
 #include "misc.h"
+#include "math_approx.h"
 
 void bw_lpc(float gamma, spx_coef_t *lpc_in, spx_coef_t *lpc_out, int order)
 {
@@ -47,6 +48,7 @@ void bw_lpc(float gamma, spx_coef_t *lpc_in, spx_coef_t *lpc_out, int order)
 }
 
 
+
 #ifdef FIXED_POINT
 
 int normalize16(spx_sig_t *x, spx_word16_t *y, int max_scale, int len)
@@ -116,7 +118,7 @@ spx_word16_t compute_rms(spx_sig_t *x, int len)
    }
    
    /*FIXME: remove division*/
-   return (1<<(sig_shift+3))*sqrt(1+sum/len)/(float)SIG_SCALING;
+   return SHR(SHL((spx_word32_t)sqroot(1+sum/len),(sig_shift+3)),SIG_SHIFT);
 }
 
 
index 748d619..6cef31c 100644 (file)
    SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 */
 
-#include <math.h>
 #include "math_approx.h"
+#include "misc.h"
 
-#ifdef SLOW_TRIG
+#ifdef FIXED_POINT
 
-float cos_sin[102] = {
-   1.00000000, 0.00000000,
-   0.99804751, 0.06245932,
-   0.99219767, 0.12467473,
-   0.98247331, 0.18640330,
-   0.96891242, 0.24740396,
-   0.95156795, 0.30743851,
-   0.93050762, 0.36627253,
-   0.90581368, 0.42367626,
-   0.87758256, 0.47942554,
-   0.84592450, 0.53330267,
-   0.81096312, 0.58509727,
-   0.77283495, 0.63460708,
-   0.73168887, 0.68163876,
-   0.68768556, 0.72600866,
-   0.64099686, 0.76754350,
-   0.59180508, 0.80608111,
-   0.54030231, 0.84147098,
-   0.48668967, 0.87357494,
-   0.43117652, 0.90226759,
-   0.37397963, 0.92743692,
-   0.31532236, 0.94898462,
-   0.25543377, 0.96682656,
-   0.19454771, 0.98089306,
-   0.13290194, 0.99112919,
-   0.07073720, 0.99749499,
-   0.00829623, 0.99996559,
-   -0.05417714, 0.99853134,
-   -0.11643894, 0.99319785,
-   -0.17824606, 0.98398595,
-   -0.23935712, 0.97093160,
-   -0.29953351, 0.95408578,
-   -0.35854022, 0.93351428,
-   -0.41614684, 0.90929743,
-   -0.47212841, 0.88152979,
-   -0.52626633, 0.85031979,
-   -0.57834920, 0.81578931,
-   -0.62817362, 0.77807320,
-   -0.67554504, 0.73731872,
-   -0.72027847, 0.69368503,
-   -0.76219923, 0.64734252,
-   -0.80114362, 0.59847214,
-   -0.83695955, 0.54726475,
-   -0.86950718, 0.49392030,
-   -0.89865940, 0.43864710,
-   -0.92430238, 0.38166099,
-   -0.94633597, 0.32318451,
-   -0.96467415, 0.26344599,
-   -0.97924529, 0.20267873,
-   -0.98999250, 0.14112001,
-   -0.99687381, 0.07901022,
-   -0.99986235, 0.01659189
-};
+#define C0 3634
+#define C1 21173
+#define C2 -12627
+#define C3 4215
 
-float speex_cos(float x)
+spx_word16_t sqroot(spx_word32_t x)
 {
-   int ind;
-   float delta;
-   ind = (int)floor(x*16+.5);
-   delta = x-0.062500*ind;
-   ind <<= 1;
-   return cos_sin[ind] - delta*(cos_sin[ind+1] + 
-                                  .5*delta*(cos_sin[ind] - 
-                                            .3333333*delta*cos_sin[ind+1]));
-}
+   int k=0;
+   spx_word32_t rt;
 
+#if 1
+   if (x>16777216)
+   {
+      x>>=10;
+      k+=5;
+   }
+   if (x>1048576)
+   {
+      x>>=6;
+      k+=3;
+   }
+   if (x>262144)
+   {
+      x>>=4;
+      k+=2;
+   }
+   if (x>32768)
+   {
+      x>>=2;
+      k+=1;
+   }
+   if (x>16384)
+   {
+      x>>=2;
+      k+=1;
+   }
+#else
+   while (x>16384)
+   {
+      x>>=2;
+      k++;
+      }
 #endif
+   rt = ADD16(C0, MULT16_16_Q14(x, ADD16(C1, MULT16_16_Q14(x, ADD16(C2, MULT16_16_Q14(x, (C3)))))));
+   rt <<= k;
+   rt >>=7;
+   return rt;
+}
 
+#endif
index c917ffd..bf48900 100644 (file)
 #ifndef MATH_APPROX_H
 #define MATH_APPROX_H
 
+#include "misc.h"
+
 float speex_cos(float x);
 
+#ifdef FIXED_POINT
+spx_word16_t sqroot(spx_word32_t x);
+#else
+#define sqroot sqrt
+#endif
 
 #endif