[xiph-commits] r12036 - trunk/speex/libspeex

[jm at svn.xiph.org]

Author: jm
Date: 2006-11-06 00:37:10 -0800 (Mon, 06 Nov 2006)
New Revision: 12036

Modified:
   trunk/speex/libspeex/math_approx.c
Log:
Re-wrote spx_sqrt in terms of spx_ilog4, a few MACs and 2 shifts.


Modified: trunk/speex/libspeex/math_approx.c
===================================================================
--- trunk/speex/libspeex/math_approx.c	2006-11-05 22:30:30 UTC (rev 12035)
+++ trunk/speex/libspeex/math_approx.c	2006-11-06 08:37:10 UTC (rev 12036)
@@ -95,64 +95,31 @@
 #ifdef FIXED_POINT
 
 /* sqrt(x) ~= 0.22178 + 1.29227*x - 0.77070*x^2 + 0.25723*x^3 (for .25 < x < 1) */
+/*#define C0 3634
+#define C1 21173
+#define C2 -12627
+#define C3 4215*/
+
+/* sqrt(x) ~= 0.22178 + 1.29227*x - 0.77070*x^2 + 0.25659*x^3 (for .25 < x < 1) */
 #define C0 3634
 #define C1 21173
 #define C2 -12627
-#define C3 4215
+#define C3 4204
 
 spx_word16_t spx_sqrt(spx_word32_t x)
 {
-   int k=0;
+   int k;
    spx_word32_t rt;
-
-   if (x<=0)
-      return 0;
-#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
-   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)))))));
-   if (rt > 16383)
-      rt = 16383;
+   k = spx_ilog4(x)-6;
    if (k>0)
-      rt <<= k;
+      x = SHR32(x, (k<<1));
    else
-      rt >>= -k;
-   rt >>=7;
+      x = SHL32(x, (-k<<1));
+   rt = ADD16(C0, MULT16_16_Q14(x, ADD16(C1, MULT16_16_Q14(x, ADD16(C2, MULT16_16_Q14(x, (C3)))))));
+   if (k>7)
+      rt <<= k-7;
+   else
+      rt >>= 7-k;
    return rt;
 }