fix.cpp

Go to the documentation of this file.
00001 /*
00002 This program is distributed under the terms of the 'MIT license'. The text
00003 of this licence follows...
00004 
00005 Copyright (c) 2004-2005 J.D.Medhurst (a.k.a. Tixy)
00006 
00007 Permission is hereby granted, free of charge, to any person obtaining a copy
00008 of this software and associated documentation files (the "Software"), to deal
00009 in the Software without restriction, including without limitation the rights
00010 to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
00011 copies of the Software, and to permit persons to whom the Software is
00012 furnished to do so, subject to the following conditions:
00013 
00014 The above copyright notice and this permission notice shall be included in
00015 all copies or substantial portions of the Software.
00016 
00017 THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
00018 IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
00019 FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.  IN NO EVENT SHALL THE
00020 AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
00021 LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
00022 OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
00023 THE SOFTWARE.
00024 */
00025 
00032 #include "common.h"
00033 #include "fix.h"
00034 
00035 
00041 #define FIX_USE_64BIT_MUL
00042 
00047 #define FIX_UNROLL_DIV_LOOP
00048 
00049 
00074 static inline int Interpolate(const int* table, int value, int shift)
00075     {
00076     // get values from table (required value lies between b and c)
00077     int i = value>>shift;
00078     int a = table[i+0];
00079     int b = table[i+1];
00080     int c = table[i+2];
00081     int d = table[i+3];
00082 
00083     // interpolate
00084     int f = value&((1<<shift)-1); // f = factional part of index
00085     int cadb = ( (c-a) - (d-b) ) >> 2;
00086     int r = (c-b) + cadb - (cadb*f>>shift);
00087     r = (unsigned)(r*f)>>shift;  // this cast to unsigned assumes table only contains increasing values (and gains us 1 more bit headroom)
00088     r += b;
00089 
00090     return r;
00091     }
00092 
00093 
00094 /*
00095 Members of class Fix
00096 */
00097 
00098 
00099 fix Fix::Add(fix a,fix b)
00100     {
00101     fix r = a+b;
00102     if((~(a^b) & (a^r)) < 0)
00103         r = (r>>31)^0x80000000;  // produce saturated result if overflow
00104     return r;
00105     }
00106 
00107 
00108 fix Fix::Sub(fix a,fix b)
00109     {
00110     fix r = a-b;
00111     if(((a^b) & (a^r)) < 0)
00112         r = (r>>31)^0x80000000;  // produce saturated result if overflow
00113     return r;
00114     }
00115 
00116 
00117 fix Fix::Mul(fix a,fix b)
00118     {
00119 #ifdef FIX_USE_64BIT_MUL
00120 
00121     int64_t r = ((int64_t)a*(int64_t)b);
00122     r += 0x8000;
00123     int32_t hi = r>>32;
00124     r >>= 16;
00125     if(hi>=0x8000)
00126         return 0x7fffffff;
00127     if(hi<-0x8000)
00128         return 0x80000000;
00129     return (int32_t)r;
00130 
00131 #else
00132 
00133     // calculate sign result
00134     int sign = a^b;
00135 
00136     // calculate absolute values
00137     if(a<0) a=-a;
00138     if(b<0) b=-b;
00139 
00140     // do the multiply in four parts
00141     uint32_t al = a&0xFFFF;
00142     uint32_t bl = b&0xFFFF;
00143     uint32_t c = al*bl;
00144     c += 0x8000;
00145     c >>= 16;
00146     uint32_t c1 = bl*((uint32_t)a>>16);
00147     c += c1;
00148     if(c>=c1) // No carry from addition
00149         {
00150         uint32_t c2 = al*((uint32_t)b>>16);
00151         c += c2;
00152         if(c>=c2) // No carry from addition
00153             {
00154             uint32_t d = ((uint32_t)a>>16)*((uint32_t)b>>16);
00155             if(d<0x10000) // No overflow from multiply
00156                 {
00157                 uint32_t dl = (uint32_t)d<<16;
00158                 c += dl;
00159                 if(c>=dl) // No carry from addition
00160                     {
00161                     if(sign<0)
00162                         {
00163                         if(c<=0x80000000)
00164                             return -(int32_t)c;
00165                         }
00166                     else
00167                         {
00168                         if(c<=0x7fffffff)
00169                             return c;
00170                         }
00171                     }
00172                 }
00173             }
00174         }
00175 
00176     // produce saturated result
00177     return (sign<0) ? 0x80000000 : 0x7fffffff;
00178 
00179 #endif
00180     }
00181 
00182 
00183 fix Fix::MulNS(fix a,fix b)
00184     {
00185 #ifdef FIX_USE_64BIT_MUL
00186     int64_t r = ((int64_t)a*(int64_t)b);
00187     r += 0x8000;
00188     r >>= 16;
00189     return (int32_t)r;
00190 #else
00191     uint32_t al = a&0xFFFF;
00192     uint32_t bl = b&0xFFFF;
00193     uint32_t r = al*bl;
00194     r += 0x8000;
00195     r >>= 16;
00196     r += bl*(a>>16);
00197     r += al*(b>>16);
00198     r += ((a>>16)*(b>>16))<<16;
00199     return r;
00200 #endif
00201     }
00202 
00203 
00204 fix Fix::Div(fix a,fix b)
00205     {
00206     // calculate sign bit of result
00207     int32_t r = a^b;
00208     r &= (1<<31);
00209 
00210     // calculate absolute values
00211     if(a<0) a = -a;
00212     if(b<0) b = -b;
00213 
00214     // calculate the number of integer bits is the result
00215     int32_t intBits = 0;
00216     uint32_t q = a;
00217     if((uint32_t)b<=((uint32_t)q>>16))
00218         intBits += 16, q >>= 16;
00219     if((uint32_t)b<=((uint32_t)q>>8))
00220         intBits += 8,  q >>= 8;
00221     if((uint32_t)b<=((uint32_t)q>>4))
00222         intBits += 4,  q >>= 4;
00223     if((uint32_t)b<=((uint32_t)q>>2))
00224         intBits += 2,  q >>= 2;
00225     if((uint32_t)b<=((uint32_t)q>>1))
00226         intBits += 1,  q >>= 1;
00227 
00228 #ifdef FIX_UNROLL_DIV_LOOP
00229 
00230     // calculate the integer part of result (bits 31 to 16)
00231     switch(intBits)
00232         {
00233         case 14:
00234             if((uint32_t)a>=((uint32_t)b<<14))
00235                 a -= (b<<14), r += 1<<(14+16);
00236         case 13:
00237             if((uint32_t)a>=((uint32_t)b<<13))
00238                 a -= (b<<13), r += 1<<(13+16);
00239         case 12:
00240             if((uint32_t)a>=((uint32_t)b<<12))
00241                 a -= (b<<12), r += 1<<(12+16);
00242         case 11:
00243             if((uint32_t)a>=((uint32_t)b<<11))
00244                 a -= (b<<11), r += 1<<(11+16);
00245         case 10:
00246             if((uint32_t)a>=((uint32_t)b<<10))
00247                 a -= (b<<10), r += 1<<(10+16);
00248         case 9:
00249             if((uint32_t)a>=((uint32_t)b<<9))
00250                 a -= (b<<9), r += 1<<(9+16);
00251         case 8:
00252             if((uint32_t)a>=((uint32_t)b<<8))
00253                 a -= (b<<8), r += 1<<(8+16);
00254         case 7:
00255             if((uint32_t)a>=((uint32_t)b<<7))
00256                 a -= (b<<7), r += 1<<(7+16);
00257         case 6:
00258             if((uint32_t)a>=((uint32_t)b<<6))
00259                 a -= (b<<6), r += 1<<(6+16);
00260         case 5:
00261             if((uint32_t)a>=((uint32_t)b<<5))
00262                 a -= (b<<5), r += 1<<(5+16);
00263         case 4:
00264             if((uint32_t)a>=((uint32_t)b<<4))
00265                 a -= (b<<4), r += 1<<(4+16);
00266         case 3:
00267             if((uint32_t)a>=((uint32_t)b<<3))
00268                 a -= (b<<3), r += 1<<(3+16);
00269         case 2:
00270             if((uint32_t)a>=((uint32_t)b<<2))
00271                 a -= (b<<2), r += 1<<(2+16);
00272         case 1:
00273             if((uint32_t)a>=((uint32_t)b<<1))
00274                 a -= (b<<1), r += 1<<(1+16);
00275         case 0:
00276             if((uint32_t)a>=((uint32_t)b<<0))
00277                 a -= (b<<0), r += 1<<(0+16);
00278             break;
00279     default:
00280         // produce saturated result
00281         return (r<0) ? 0x80000000 : 0x7fffffff; // saturated result
00282         }
00283 
00284     // calculate the fractional part of result (bits 15 to 0)
00285     a <<= 1;
00286     if((uint32_t)a>=(uint32_t)b)
00287         a -= b, r += 1<<15;
00288     a <<= 1;
00289     if((uint32_t)a>=(uint32_t)b)
00290         a -= b, r += 1<<14;
00291     a <<= 1;
00292     if((uint32_t)a>=(uint32_t)b)
00293         a -= b, r += 1<<13;
00294     a <<= 1;
00295     if((uint32_t)a>=(uint32_t)b)
00296         a -= b, r += 1<<12;
00297     a <<= 1;
00298     if((uint32_t)a>=(uint32_t)b)
00299         a -= b, r += 1<<11;
00300     a <<= 1;
00301     if((uint32_t)a>=(uint32_t)b)
00302         a -= b, r += 1<<10;
00303     a <<= 1;
00304     if((uint32_t)a>=(uint32_t)b)
00305         a -= b, r += 1<<9;
00306     a <<= 1;
00307     if((uint32_t)a>=(uint32_t)b)
00308         a -= b, r += 1<<8;
00309     a <<= 1;
00310     if((uint32_t)a>=(uint32_t)b)
00311         a -= b, r += 1<<7;
00312     a <<= 1;
00313     if((uint32_t)a>=(uint32_t)b)
00314         a -= b, r += 1<<6;
00315     a <<= 1;
00316     if((uint32_t)a>=(uint32_t)b)
00317         a -= b, r += 1<<5;
00318     a <<= 1;
00319     if((uint32_t)a>=(uint32_t)b)
00320         a -= b, r += 1<<4;
00321     a <<= 1;
00322     if((uint32_t)a>=(uint32_t)b)
00323         a -= b, r += 1<<3;
00324     a <<= 1;
00325     if((uint32_t)a>=(uint32_t)b)
00326         a -= b, r += 1<<2;
00327     a <<= 1;
00328     if((uint32_t)a>=(uint32_t)b)
00329         a -= b, r += 1<<1;
00330     a <<= 1;
00331     if((uint32_t)a>=(uint32_t)b)
00332         a -= b, r += 1<<0;
00333 
00334 #else
00335 
00336     if(intBits>14)
00337         return (r<0) ? 0x80000000 : 0x7fffffff; // saturated result
00338 
00339     // calculate the integer part of result (bits 31 to 16)
00340     b <<= intBits;
00341     int32_t bit = 0x10000<<intBits;
00342     while(bit>0x10000)
00343         {
00344         if((uint32_t)a>=(uint32_t)b)
00345             {
00346             a -= b;
00347             r += bit;
00348             }
00349         b >>= 1;
00350         bit >>= 1;
00351         }
00352     if((uint32_t)a>=(uint32_t)b)
00353         {
00354         a -= b;
00355         r += bit;
00356         }
00357     bit >>= 1;
00358 
00359     // calculate the fractional part of result (bits 15 to 0)
00360     do
00361         {
00362         a <<= 1;
00363         if((uint32_t)a>=(uint32_t)b)
00364             {
00365             a -= b;
00366             r += bit;
00367             }
00368         bit >>= 1;
00369         }
00370     while(bit);
00371 
00372 #endif
00373 
00374     return (r<0) ? 0x80000000-r : r;
00375 
00376     }
00377 
00378 
00379 fix Fix::Sqrt(ufix a)
00380     {
00381     ufix r,nr,m;
00382 
00383     // calculate integer part (bits 31 to 16)
00384     r = 0;
00385     m = 0x40000000;
00386     do
00387         {
00388         nr = r+m;
00389         if(nr<=a)
00390             {
00391             a -= nr;
00392             r = nr+m;
00393             }
00394         r >>= 1;
00395         m >>= 2;
00396         }
00397     while(m);
00398 
00399     // calculate bits 15 to 8
00400     r <<= 8;
00401     a <<= 8;
00402     m = 0x40;
00403     do
00404         {
00405         nr = r+m;
00406         if(nr<=a)
00407             {
00408             a -= nr;
00409             r = nr+m;
00410             }
00411         r >>= 1;
00412         m >>= 2;
00413         }
00414     while(m);
00415 
00416     // calculate bits 7 to 0
00417     r <<= 8;
00418     a <<= 8;
00419     m = 0x40;
00420     do
00421         {
00422         nr = r+m;
00423         if(nr<=a)
00424             {
00425             a -= nr;
00426             r = nr+m;
00427             }
00428         r >>= 1;
00429         m >>= 2;
00430         }
00431     while(m);
00432 
00433     // round result
00434     if(r<a)
00435         r++;
00436 
00437     return r;
00438     }
00439 
00440 
00441 fix Fix::Log2(ufix a)
00442     {
00443     // trap 0
00444     if(a==0)
00445         return (fix)-0x80000000;
00446 
00447     // calculate integer part of result in i
00448     // and set n = normalised value of a
00449     int32_t i=15*0x10000;
00450     uint32_t n=a;
00451     if(n<(uint32_t)(1<<(32-16)))
00452         n <<= 16, i -= 16*0x10000;
00453     if(n<(uint32_t)(1<<(32-8)))
00454         n <<= 8, i -= 8*0x10000;
00455     if(n<(uint32_t)(1<<(32-4)))
00456         n <<= 4, i -= 4*0x10000;
00457     if(n<(uint32_t)(1<<(32-2)))
00458         n <<= 2, i -= 2*0x10000;
00459     if(n<(uint32_t)(1<<(32-1)))
00460         n <<= 1, i -= 1*0x10000;
00461 
00462     // reduce n to the 23 most significant bits and clear
00463     // the most significant bit, leaving a 22 bit value in n
00464     n = (n-0x80000000+(1<<8))>>9;
00465 
00466     // calculate fractional part of result (in f) by interpolation of lookup table values
00467     static const int32_t LogTable[] =
00468         {
00469         0xffff45e1,
00470         0x00000000,0x0000b73d,0x00016bad,0x00021d67,0x0002cc7f,0x00037908,0x00042316,0x0004caba,
00471         0x00057007,0x0006130b,0x0006b3d8,0x0007527c,0x0007ef06,0x00088984,0x00092204,0x0009b892,
00472         0x000a4d3c,0x000ae00d,0x000b7111,0x000c0053,0x000c8ddd,0x000d19bb,0x000da3f6,0x000e2c98,
00473         0x000eb3aa,0x000f3935,0x000fbd43,0x00103fdb,0x0010c105,0x001140ca,0x0011bf31,0x00123c42,
00474         0x0012b803,0x0013327c,0x0013abb4,0x001423b0,0x00149a78,0x00151012,0x00158482,0x0015f7d0,
00475         0x00166a01,0x0016db19,0x00174b20,0x0017ba19,0x0018280a,0x001894f7,0x001900e6,0x00196bdb,
00476         0x0019d5da,0x001a3ee8,0x001aa709,0x001b0e41,0x001b7495,0x001bda07,0x001c3e9d,0x001ca259,
00477         0x001d053f,0x001d6754,0x001dc89a,0x001e2914,0x001e88c7,0x001ee7b4,0x001f45e1,0x001fa34e,
00478         0x00200000,0x00205bf9,0x0020b73d
00479         };
00480     int32_t f = Interpolate(LogTable,n,16);
00481 
00482     // scale result to 16 bits (from the 22 bit precision used in lookup table)
00483     f = (f+(1<<(5-1)))>>5;
00484 
00485     // return sum of integer and fractional part
00486     return i+f;
00487     }
00488 
00489 
00490 ufix Fix::Exp2(fix a)
00491     {
00492     if(a>=0x00100000)
00493         return 0xffffffffu; // result will be too big so result max value
00494 
00495     // special cases for small values
00496     if(a<-0x000ead96)
00497         {
00498         if(a<-0x00110000)
00499             return 0;
00500         if(a<-0x000f6a3f)
00501             return 1;
00502         return 2;
00503         }
00504 
00505     // table of values for ((2^n)-1)<<32 for n in range 0x0000.0000 to 0x0000.ff00
00506     static const int32_t Exp2TableFF00[] =
00507         {
00508         0x00000000,0x00b1afa5,0x0163da9f,0x02168143,0x02c9a3e7,0x037d42e1,0x04315e86,0x04e5f72f,
00509         0x059b0d31,0x0650a0e3,0x0706b29d,0x07bd42b7,0x08745187,0x092bdf66,0x09e3ecac,0x0a9c79b1,
00510         0x0b5586cf,0x0c0f145e,0x0cc922b7,0x0d83b233,0x0e3ec32d,0x0efa55fd,0x0fb66aff,0x1073028d,
00511         0x11301d01,0x11edbab5,0x12abdc06,0x136a814f,0x1429aaea,0x14e95934,0x15a98c8a,0x166a4547,
00512         0x172b83c7,0x17ed4869,0x18af9388,0x19726583,0x1a35beb6,0x1af99f81,0x1bbe0840,0x1c82f952,
00513         0x1d487316,0x1e0e75eb,0x1ed5022f,0x1f9c1843,0x2063b886,0x212be357,0x21f49917,0x22bdda27,
00514         0x2387a6e7,0x2451ffb8,0x251ce4fb,0x25e85711,0x26b4565e,0x2780e341,0x284dfe1f,0x291ba759,
00515         0x29e9df51,0x2ab8a66d,0x2b87fd0d,0x2c57e397,0x2d285a6e,0x2df961f6,0x2ecafa93,0x2f9d24ab,
00516         0x306fe0a3,0x31432ede,0x32170fc4,0x32eb83ba,0x33c08b26,0x3496266e,0x356c55f9,0x36431a2d,
00517         0x371a7373,0x37f26231,0x38cae6d0,0x39a401b7,0x3a7db34e,0x3b57fbfe,0x3c32dc31,0x3d0e544e,
00518         0x3dea64c1,0x3ec70df1,0x3fa4504a,0x40822c36,0x4160a21f,0x423fb270,0x431f5d95,0x43ffa3f8,
00519         0x44e08606,0x45c2042a,0x46a41ed1,0x4786d668,0x486a2b5c,0x494e1e19,0x4a32af0d,0x4b17dea6,
00520         0x4bfdad53,0x4ce41b81,0x4dcb299f,0x4eb2d81d,0x4f9b2769,0x508417f4,0x516daa2c,0x5257de83,
00521         0x5342b569,0x542e2f4f,0x551a4ca5,0x56070dde,0x56f4736b,0x57e27dbe,0x58d12d49,0x59c0827f,
00522         0x5ab07dd4,0x5ba11fba,0x5c9268a5,0x5d845909,0x5e76f15a,0x5f6a320d,0x605e1b97,0x6152ae6c,
00523         0x6247eb03,0x633dd1d1,0x6434634c,0x652b9feb,0x66238825,0x671c1c70,0x68155d44,0x690f4b19,
00524         0x6a09e667,0x6b052fa7,0x6c012750,0x6cfdcddd,0x6dfb23c6,0x6ef92985,0x6ff7df95,0x70f7466f,
00525         0x71f75e8e,0x72f8286e,0x73f9a48a,0x74fbd35d,0x75feb564,0x77024b1a,0x780694fd,0x790b938a,
00526         0x7a11473e,0x7b17b097,0x7c1ed013,0x7d26a62f,0x7e2f336c,0x7f387849,0x80427543,0x814d2add,
00527         0x82589994,0x8364c1eb,0x8471a462,0x857f4179,0x868d99b4,0x879cad93,0x88ac7d98,0x89bd0a47,
00528         0x8ace5422,0x8be05bad,0x8cf3216b,0x8e06a5e0,0x8f1ae991,0x902fed02,0x9145b0b9,0x925c353a,
00529         0x93737b0c,0x948b82b5,0x95a44cbc,0x96bdd9a7,0x97d829fd,0x98f33e47,0x9a0f170c,0x9b2bb4d5,
00530         0x9c49182a,0x9d674194,0x9e86319e,0x9fa5e8d0,0xa0c667b5,0xa1e7aed8,0xa309bec4,0xa42c9804,
00531         0xa5503b23,0xa674a8af,0xa799e133,0xa8bfe53c,0xa9e6b557,0xab0e5213,0xac36bbfd,0xad5ff3a3,
00532         0xae89f995,0xafb4ce62,0xb0e07298,0xb20ce6c9,0xb33a2b84,0xb468415b,0xb59728de,0xb6c6e29f,
00533         0xb7f76f2f,0xb928cf22,0xba5b030a,0xbb8e0b79,0xbcc1e904,0xbdf69c3f,0xbf2c25bd,0xc0628614,
00534         0xc199bdd8,0xc2d1cd9f,0xc40ab5ff,0xc544778f,0xc67f12e5,0xc7ba8898,0xc8f6d940,0xca340575,
00535         0xcb720dce,0xccb0f2e6,0xcdf0b555,0xcf3155b5,0xd072d4a0,0xd1b532b0,0xd2f87080,0xd43c8eac,
00536         0xd5818dcf,0xd6c76e86,0xd80e316c,0xd955d71f,0xda9e603d,0xdbe7cd63,0xdd321f30,0xde7d5641,
00537         0xdfc97337,0xe11676b1,0xe264614f,0xe3b333b1,0xe502ee78,0xe6539246,0xe7a51fbc,0xe8f7977c,
00538         0xea4afa2a,0xeb9f4867,0xecf482d8,0xee4aaa21,0xefa1bee6,0xf0f9c1cb,0xf252b376,0xf3ac948d,
00539         0xf50765b6,0xf6632798,0xf7bfdad9,0xf91d8022,0xfa7c1819,0xfbdba369,0xfd3c22b8,0xfe9d96b2
00540         };
00541 
00542     // table of values for ((2^n)-1)<<40 for n in range 0x0000.0000 to 0x0000.00ff
00543     static const int32_t Exp2Table00FF[] =
00544         {
00545         0x00000000,0x00b17255,0x0162e525,0x02145871,0x02c5cc37,0x03774079,0x0428b535,0x04da2a6d,
00546         0x058ba01f,0x063d164d,0x06ee8cf5,0x07a00419,0x08517bb7,0x0902f3d1,0x09b46c65,0x0a65e575,
00547         0x0b175eff,0x0bc8d905,0x0c7a5386,0x0d2bce81,0x0ddd49f8,0x0e8ec5e9,0x0f404256,0x0ff1bf3e,
00548         0x10a33ca1,0x1154ba7e,0x120638d7,0x12b7b7ab,0x136936fa,0x141ab6c4,0x14cc3709,0x157db7c9,
00549         0x162f3904,0x16e0baba,0x17923ceb,0x1843bf97,0x18f542be,0x19a6c660,0x1a584a7d,0x1b09cf16,
00550         0x1bbb5429,0x1c6cd9b7,0x1d1e5fc1,0x1dcfe645,0x1e816d45,0x1f32f4bf,0x1fe47cb5,0x20960526,
00551         0x21478e11,0x21f91778,0x22aaa15a,0x235c2bb7,0x240db68f,0x24bf41e2,0x2570cdb0,0x262259f9,
00552         0x26d3e6bd,0x278573fd,0x283701b7,0x28e88fed,0x299a1e9d,0x2a4badc9,0x2afd3d6f,0x2baecd91,
00553         0x2c605e2e,0x2d11ef46,0x2dc380d9,0x2e7512e7,0x2f26a570,0x2fd83874,0x3089cbf4,0x313b5fee,
00554         0x31ecf464,0x329e8954,0x33501ec0,0x3401b4a7,0x34b34b09,0x3564e1e6,0x3616793e,0x36c81111,
00555         0x3779a95f,0x382b4228,0x38dcdb6d,0x398e752c,0x3a400f67,0x3af1aa1d,0x3ba3454e,0x3c54e0fa,
00556         0x3d067d21,0x3db819c3,0x3e69b6e0,0x3f1b5479,0x3fccf28c,0x407e911b,0x41303025,0x41e1cfaa,
00557         0x42936faa,0x43451025,0x43f6b11b,0x44a8528d,0x4559f479,0x460b96e1,0x46bd39c3,0x476edd21,
00558         0x482080fa,0x48d2254e,0x4983ca1e,0x4a356f68,0x4ae7152e,0x4b98bb6e,0x4c4a622a,0x4cfc0961,
00559         0x4dadb113,0x4e5f5941,0x4f1101e9,0x4fc2ab0d,0x507454ab,0x5125fec5,0x51d7a95a,0x5289546a,
00560         0x533afff5,0x53ecabfc,0x549e587d,0x5550057a,0x5601b2f2,0x56b360e5,0x57650f53,0x5816be3d,
00561         0x58c86da1,0x597a1d81,0x5a2bcddc,0x5add7eb2,0x5b8f3003,0x5c40e1cf,0x5cf29417,0x5da446da,
00562         0x5e55fa17,0x5f07add1,0x5fb96205,0x606b16b4,0x611ccbdf,0x61ce8184,0x628037a5,0x6331ee41,
00563         0x63e3a559,0x64955ceb,0x654714f9,0x65f8cd82,0x66aa8686,0x675c4005,0x680df9ff,0x68bfb475,
00564         0x69716f66,0x6a232ad2,0x6ad4e6b9,0x6b86a31b,0x6c385ff9,0x6cea1d52,0x6d9bdb26,0x6e4d9975,
00565         0x6eff583f,0x6fb11785,0x7062d746,0x71149782,0x71c65839,0x7278196b,0x7329db19,0x73db9d42,
00566         0x748d5fe6,0x753f2305,0x75f0e6a0,0x76a2aab5,0x77546f46,0x78063452,0x78b7f9da,0x7969bfdc,
00567         0x7a1b865a,0x7acd4d53,0x7b7f14c7,0x7c30dcb7,0x7ce2a522,0x7d946e07,0x7e463769,0x7ef80145,
00568         0x7fa9cb9d,0x805b9670,0x810d61be,0x81bf2d87,0x8270f9cc,0x8322c68c,0x83d493c7,0x8486617d,
00569         0x85382fae,0x85e9fe5b,0x869bcd83,0x874d9d27,0x87ff6d45,0x88b13ddf,0x89630ef4,0x8a14e084,
00570         0x8ac6b290,0x8b788517,0x8c2a5819,0x8cdc2b96,0x8d8dff8f,0x8e3fd403,0x8ef1a8f2,0x8fa37e5c,
00571         0x90555442,0x91072aa3,0x91b9017f,0x926ad8d6,0x931cb0a9,0x93ce88f7,0x948061c0,0x95323b05,
00572         0x95e414c5,0x9695ef00,0x9747c9b6,0x97f9a4e8,0x98ab8095,0x995d5cbd,0x9a0f3961,0x9ac1167f,
00573         0x9b72f41a,0x9c24d22f,0x9cd6b0c0,0x9d888fcc,0x9e3a6f53,0x9eec4f55,0x9f9e2fd3,0xa05010cc,
00574         0xa101f241,0xa1b3d430,0xa265b69b,0xa3179982,0xa3c97ce3,0xa47b60c0,0xa52d4519,0xa5df29ec,
00575         0xa6910f3b,0xa742f505,0xa7f4db4b,0xa8a6c20b,0xa958a947,0xaa0a90ff,0xaabc7932,0xab6e61e0,
00576         0xac204b09,0xacd234ae,0xad841ece,0xae360969,0xaee7f480,0xaf99e011,0xb04bcc1f,0xb0fdb8a7
00577         };
00578 
00579     // Treat 'a' as i+j+k where i is integer part of value, j is bits 15 to 8 of
00580     // the fractional part and k is bits 7 to 0 of the fractional part
00581     // We calculate 2^a as 2^(i+j+k) which is (2^i)*(2^j)*(2^k)
00582 
00583     int i = a>>16; // i = integer part of value
00584     uint32_t x=Exp2TableFF00[(a>>8)&0xff];  // x = (2^j-1) << 32
00585     uint32_t y=Exp2Table00FF[a&0xff];       // y = (2^k-1) << 40
00586 
00587     // calculate p = (x*y)>>32 = (2^j-1)(2^k-1)<<40
00588     uint32_t xl = x&0xffff;
00589     uint32_t yl = y&0xffff;
00590     uint32_t xh = x>>16;
00591     uint32_t yh = y>>16;
00592     uint32_t pl = (xl*yl+(1<<15))>>16;
00593     uint32_t pm1 = xh*yl+pl;
00594     uint32_t pm2 = xl*yh;
00595     unsigned p = xh*yh;
00596     uint32_t pm = pm1+pm2;
00597     if(pm<pm2) p += (1<<16); // check for carry
00598     if(pm&(1<<15)) p += 1; // round up if required
00599     p += pm>>16;
00600 
00601     // calculate n = p+x+y = x*y+x+y = (x+1)(y+1)-1 = ((2^j)(2^k)-1)<<32
00602     unsigned round = ((p&0xff)+(y&0xff)+0x80)>>8;
00603     unsigned n = x+(y>>8)+round;
00604     n += p>>8;
00605 
00606     // n = n>>(16-i) = 2^(i-16)*n = ((2^i)*((2^j)*(2^k)-1)<<16
00607     n >>= 15-i;
00608     n = (n+1)>>1;
00609 
00610     // finally add 2^i to get result of 2^a
00611     n += 0x80000000>>(15-i);
00612 
00613     return (ufix)n;
00614     }
00615 
00616 
00617 fix Fix::Cos(fixangle angle)
00618     {
00619     return Fix::Sin(angle+0x4000);
00620     }
00621 
00622 
00623 fix Fix::Sin(fixangle angle)
00624     {
00625     // reduce angle to first quadrant
00626     unsigned n = angle&0x3FFF;
00627     if(angle&(1<<14))
00628         n = 0x4000-n;
00629 
00630     // calc Sin through table lookup
00631     static const int32_t SinTable[] =
00632         {
00633         -0x000c8fb3,
00634         0x00000000,0x000c8fb3,0x001917a7,0x00259021,0x0031f170,0x003e33f3,0x004a5019,0x00563e6a,
00635         0x0061f78b,0x006d7440,0x0078ad75,0x00839c3d,0x008e39da,0x00987fc0,0x00a26799,0x00abeb4a,
00636         0x00b504f3,0x00bdaef9,0x00c5e403,0x00cd9f02,0x00d4db31,0x00db941a,0x00e1c598,0x00e76bd8,
00637         0x00ec835e,0x00f10908,0x00f4fa0b,0x00f853f8,0x00fb14be,0x00fd3aac,0x00fec46d,0x00ffb10f,
00638         0x01000000,0x00ffb10f,0x00ffb10f
00639         };
00640     int r = Interpolate(SinTable,n,9);
00641 
00642     // scale result to 16 bits (from the 24 bit precision used in lookup table)
00643     r = (r+(1<<(8-1)))>>8;
00644 
00645     // produce correct sign for result
00646     if(angle&(1<<15))
00647         r = -r;
00648 
00649     return r;
00650     }
00651 
00652 
00653 fix Fix::Tan(fixangle angle)
00654     {
00655     // n = angle in first quadrant
00656     int n = angle&0x7fff;
00657     bool neg = n>0x4000;
00658     if(neg) n = 0x8000-n;
00659 
00660     // calculate tangent by interpolation of lookup table values
00661     unsigned r;
00662     if(n<=0x2000)
00663         {
00664         static const int32_t TanTable0000[] =
00665             {
00666             0xffe6dcba,
00667             0x00000000,0x00192346,0x00324e4f,0x004b88e7,0x0064daee,0x007e4c61,0x0097e564,0x00b1ae4d,
00668             0x00cbafaf,0x00e5f267,0x01007fa7,0x011b6104,0x0136a083,0x015248ae,0x016e649f,0x018b0019,
00669             0x01a8279a,0x01c5e872,0x01e450e1,0x02037033,0x022356e3,0x024416be,0x0265c311,0x028870db,
00670             0x02ac3705,0x02d12ea3,0x02f7733e,0x031f232f,0x03486005,0x03734ef8,0x03a01979,0x03ceedd6,
00671             0x04000000,0x04338a74
00672             };
00673         r = Interpolate(TanTable0000,n-0x0000,8);
00674         // scale result to 16 bits (from the 26 bit precision used in lookup table)
00675         r = (r+(1<<(10-1)))>>10;
00676         }
00677     else if(n<=0x2d80)
00678         {
00679         static const int32_t TanTable2000[] =
00680             {
00681             0x01f395da,
00682             0x02000000,0x020cb920,0x0219c53a,0x02272890,0x0234e7ab,0x02430762,0x02518ce0,0x02607dac,
00683             0x026fdfb2,0x027fb948,0x0290113f,0x02a0eee8,0x02b25a22,0x02c45b6c,0x02d6fbf1,0x02ea4598,
00684             0x02fe431c,0x03130022,0x0328894f,0x033eec66,0x0356386c,0x036e7dc8,0x0387ce70,0x03a23e1b,
00685             0x03bde277,0x03dad368,0x03f92b57,0x04190786,0x043a8876
00686             };
00687         r = Interpolate(TanTable2000,n-0x2000,7);
00688         // scale result to 16 bits (from the 25 bit precision used in lookup table)
00689         r = (r+(1<<(9-1)))>>9;
00690         }
00691     else if(n<=0x35c0)
00692         {
00693         static const int32_t TanTable2D80[] =
00694             {
00695             0x0204737a,0x020c83c3,0x0214c89d,
00696             0x021d443b,0x0225f8ef,0x022ee92e,0x02381791,0x024186d9,0x024b39ef,0x025533ea,0x025f7813,
00697             0x026a09e6,0x0274ed19,0x0280259c,0x028bb7a6,0x0297a7b1,0x02a3fa88,0x02b0b54a,0x02bddd70,
00698             0x02cb78da,0x02d98dd2,0x02e8231d,0x02f73fff,0x0306ec4e,0x0317307b,0x032815a6,0x0339a5ac,
00699             0x034beb3d,0x035ef1f2,0x0372c665,0x03877650,0x039d10ad,0x03b3a5d7,0x03cb47bd,0x03e40a0a,
00700             0x03fe0261
00701             };
00702         r = Interpolate(TanTable2D80,n-0x2d80,6);
00703         // scale result to 16 bits (from the 24 bit precision used in lookup table)
00704         r = (r+(1<<(8-1)))>>8;
00705         }
00706     else if(n<=0x39e0)
00707         {
00708         static const int32_t TanTable35C0[] =
00709             {
00710             0x01ebc1c4,0x01f20505,0x01f86f02,
00711             0x01ff0130,0x0205bd16,0x020ca44f,0x0213b88a,0x021afb90,0x02226f3e,0x022a158f,0x0231f097,
00712             0x023a0288,0x02424db5,0x024ad491,0x025399b6,0x025c9fe3,0x0265ea01,0x026f7b26,0x0279569c,
00713             0x02837fdc,0x028dfa9d,0x0298cace,0x02a3f4a5,0x02af7c9d,0x02bb6780,0x02c7ba6a,0x02d47ad7,
00714             0x02e1aea3,0x02ef5c1b,0x02fd8a00,0x030c3f99,0x031b84b8,0x032b61d0,0x033be000,0x034d0923,
00715             0x035ee7ea
00716             };
00717         r = Interpolate(TanTable35C0,n-0x35c0,5);
00718         // scale result to 16 bits (from the 23 bit precision used in lookup table)
00719         r = (r+(1<<(7-1)))>>7;
00720         }
00721     else if(n<=0x3c70)
00722         {
00723         static const int32_t TanTable39E0[] =
00724             {
00725             0x01a22f46,0x01a68491,0x01aaf08f,
00726             0x01af73f5,0x01b40f80,0x01b8c3f6,0x01bd9224,0x01c27ae2,0x01c77f0f,0x01cc9f96,0x01d1dd6a,
00727             0x01d7398d,0x01dcb509,0x01e250f7,0x01e80e7b,0x01edeec8,0x01f3f321,0x01fa1cd8,0x02006d4d,
00728             0x0206e5f6,0x020d885a,0x02145612,0x021b50d0,0x02227a5a,0x0229d48f,0x02316169,0x023922fc,
00729             0x02411b7b,0x02494d38,0x0251baa7,0x025a6660,0x02635323,0x026c83da,0x0275fb9a,0x027fbdac,
00730             0x0289cd8b,0x02942eec,0x029ee5c0,0x02a9f63c,0x02b564da,0x02c13664,0x02cd6ff9,0x02da1711,
00731             0x02e7318b
00732             };
00733         r = Interpolate(TanTable39E0,n-0x39e0,4);
00734         // scale result to 16 bits (from the 22 bit precision used in lookup table)
00735         r = (r+(1<<(6-1)))>>6;
00736         }
00737     else if(n<=0x3df0)
00738         {
00739         static const int32_t TanTable3C70[] =
00740             {
00741             0x0169dabd,0x016d0b89,0x01704ac0,
00742             0x017398c5,0x0176f600,0x017a62d9,0x017ddfbf,0x01816d24,0x01850b7f,0x0188bb49,0x018c7d04,
00743             0x01905133,0x01943860,0x0198331a,0x019c41f6,0x01a0658d,0x01a49e82,0x01a8ed7c,0x01ad5328,
00744             0x01b1d03c,0x01b66576,0x01bb139b,0x01bfdb79,0x01c4bde6,0x01c9bbc2,0x01ced5f9,0x01d40d7d,
00745             0x01d96350,0x01ded87c,0x01e46e1a,0x01ea254f,0x01efff4d,0x01f5fd57,0x01fc20be,0x02026ae5,
00746             0x0208dd40,0x020f7955,0x021640bf,0x021d352e,0x0224586a,0x022bac51,0x023332de,0x023aee23,
00747             0x0242e055,0x024b0bc4,0x025372e6,0x025c1852,0x0264fec8,0x026e2932,0x02779aa6,0x0281566b
00748             };
00749         r = Interpolate(TanTable3C70,n-0x3c70,3);
00750         // scale result to 16 bits (from the 21 bit precision used in lookup table)
00751         r = (r+(1<<(5-1)))>>5;
00752         }
00753     else if(n<=0x3ebc)
00754         {
00755         static const int32_t TanTable3DF0[] =
00756             {
00757             0x01396c6c,0x013bcd53,0x013e3784,0x0140ab35,0x014328a0,
00758             0x0145afff,0x0148418d,0x014add89,0x014d8433,0x015035cd,0x0152f29c,0x0155bae5,0x01588ef2,
00759             0x015b6f0e,0x015e5b87,0x016154ad,0x01645ad4,0x01676e52,0x016a8f7f,0x016dbeb8,0x0170fc5c,
00760             0x017448cf,0x0177a477,0x017b0fbd,0x017e8b10,0x018216e2,0x0185b3aa,0x018961e2,0x018d220a,
00761             0x0190f4a6,0x0194da41,0x0198d368,0x019ce0b1,0x01a102b6,0x01a53a18,0x01a9877e,0x01adeb98,
00762             0x01b26719,0x01b6fac1,0x01bba753,0x01c06d9e,0x01c54e79,0x01ca4ac3,0x01cf6366,0x01d49958,
00763             0x01d9ed98,0x01df6132,0x01e4f53d,0x01eaaadf,0x01f0834b,0x01f67fc3,0x01fca197,0x0202ea2c,
00764             0x02095af4
00765             };
00766         r = Interpolate(TanTable3DF0,n-0x3df0,2);
00767         // scale result to 16 bits (from the 20 bit precision used in lookup table)
00768         r = (r+(1<<(4-1)))>>4;
00769         }
00770     else if(n<=0x3f90)
00771         {
00772         static const int32_t TanTable3EBC[] =
00773             {
00774             0x00ffe079,0x01017516,0x01030eb9,
00775             0x0104ad7a,0x01065172,0x0107fabb,0x0109a96e,0x010b5da7,0x010d1780,0x010ed715,0x01109c84,
00776             0x011267ea,0x01143965,0x01161115,0x0117ef18,0x0119d391,0x011bbea1,0x011db06b,0x011fa912,
00777             0x0121a8bb,0x0123af8b,0x0125bda9,0x0127d33e,0x0129f071,0x012c156d,0x012e425e,0x0130776f,
00778             0x0132b4cf,0x0134faae,0x0137493b,0x0139a0a9,0x013c012c,0x013e6af8,0x0140de45,0x01435b4b,
00779             0x0145e245,0x0148736f,0x014b0f06,0x014db54c,0x01506681,0x015322eb,0x0155ead0,0x0158be78,
00780             0x015b9e30,0x015e8a44,0x01618306,0x016488c8,0x01679be1,0x016abcaa,0x016deb7e,0x017128be,
00781             0x017474cc,0x0177d00f,0x017b3af1,0x017eb5e0,0x0182414d,0x0185ddb0,0x01898b84,0x018d4b47,
00782             0x01911d7f,0x019502b5,0x0198fb77,0x019d085b,0x01a129fc,0x01a560f9,0x01a9adfb,0x01ae11b0,
00783             0x01b28ccd,0x01b72010,0x01bbcc3e,0x01c09224,0x01c5729a,0x01ca6e80,0x01cf86bf,0x01d4bc4c,
00784             0x01da1028,0x01df835d,0x01e51704,0x01eacc41,0x01f0a448,0x01f6a05b,0x01fcc1cc,0x020309fc,
00785             0x02097a5f,0x0210147d,0x0216d9f0,0x021dcc68,0x0224edad,0x022c3f9d,0x0233c433,0x023b7d81,
00786             0x02436dbc,0x024b9734,0x0253fc5f,0x025c9fd4,0x02658453,0x026eacc6,0x02781c43,0x0281d611,
00787             0x028bddad,0x029636cb,0x02a0e55c,0x02abed94,0x02b753f0,0x02c31d37,0x02cf4e89,0x02dbed5f,
00788             0x02e8ff97,0x02f68b7b
00789             };
00790         r = Interpolate(TanTable3EBC,n-0x3ebc,1);
00791         // scale result to 16 bits (from the 19 bit precision used in lookup table)
00792         r = (r+(1<<(3-1)))>>3;
00793         }
00794     else if(n<0x4000)
00795         {
00796         static const int32_t TanTable3F90[] =
00797             {
00798             0x005d1ff3,0x005df6bd,0x005ed16f,0x005fb025,0x006092fa,0x00617a0c,0x0062657b,0x00635566,
00799             0x006449ed,0x00654334,0x0066415f,0x00674491,0x00684cf3,0x00695aac,0x006a6de7,0x006b86ce,
00800             0x006ca58f,0x006dca59,0x006ef55e,0x007026d1,0x00715ee9,0x00729ddc,0x0073e3e5,0x00753142,
00801             0x00768633,0x0077e2fa,0x007947dd,0x007ab526,0x007c2b22,0x007daa20,0x007f3276,0x0080c47c,
00802             0x0082608e,0x00840710,0x0085b866,0x008774fe,0x00893d49,0x008b11bf,0x008cf2de,0x008ee12b,
00803             0x0090dd33,0x0092e78b,0x009500d0,0x009729a7,0x009962c0,0x009bacd6,0x009e08af,0x00a0771d,
00804             0x00a2f8fd,0x00a58f3f,0x00a83add,0x00aafce4,0x00add675,0x00b0c8c1,0x00b3d50f,0x00b6fcbe,
00805             0x00ba4145,0x00bda438,0x00c12747,0x00c4cc43,0x00c89521,0x00cc83fe,0x00d09b21,0x00d4dd01,
00806             0x00d94c4b,0x00ddebe4,0x00e2bef2,0x00e7c8e5,0x00ed0d7a,0x00f290c9,0x00f8574c,0x00fe65ee,
00807             0x0104c218,0x010b71c1,0x01127b80,0x0119e6a4,0x0121bb49,0x012a027b,0x0132c656,0x013c122e,
00808             0x0145f2c4,0x0150767c,0x015bada6,0x0167aad4,0x0174833b,0x01824f38,0x01912ae6,0x01a136df,
00809             0x01b2992c,0x01c57e75,0x01da1b7f,0x01f0af1b,0x020984ae,0x0224f778,0x02437703,0x02658d16,
00810             0x028be5ec,0x02b75bab,0x02e906ce,0x0322561d,0x036532a4,0x03b43743,0x0413099b,0x0486ee3e,
00811             0x0517cc0b,0x05d20dc8,0x06ca656d,0x08261355,0x0a2f982f,0x0d94caee,0x145f306a,0x28be60d9,
00812             0xffffffff
00813             };
00814         r = TanTable3F90[n-0x3f90];
00815         }
00816     else
00817         {
00818         // result is infinity so return saturated result
00819         return angle<0 ? (int32_t)-0x80000000 : (int32_t)0x7fffffff;
00820         }
00821 
00822     return neg ? -(int32_t)r: (int32_t)r;
00823     }
00824 
00825 
00826 fixangle Fix::ACos(fix value)
00827     {
00828     if(value<-0x10000)
00829         value = -0x10000;
00830     else if(value>0x10000)
00831         value = 0x10000;
00832     return 0x4000-Fix::ASin(value);
00833     }
00834 
00835 
00836 fixangle Fix::ASin(fix value)
00837     {
00838     unsigned n = value;
00839     if(value<0)
00840         n = (unsigned)-(int)n;
00841 
00842     // calculare arc-sine by interpolation of lookup table values
00843     int r;
00844     if(n<0xc000)
00845         {
00846         static const int32_t ASinTable0000[] =
00847             {
00848             0xfffeb9ff,
00849             0x00000000,0x00014601,0x00028c53,0x0003d349,0x00051b37,0x00066474,0x0007af57,0x0008fc3f,
00850             0x000a4b8d,0x000b9daa,0x000cf305,0x000e4c19,0x000fa967,0x00110b83,0x0012730c,0x0013e0b9,
00851             0x00155555,0x0016d1cc,0x0018572d,0x0019e6b6,0x001b81e1,0x001d2a76,0x001ee2a7,0x0020ad31,
00852             0x00228d9c,0x00248890
00853             };
00854         r = Interpolate(ASinTable0000,n,11);
00855         }
00856     else if(n<0xf200)
00857         {
00858         static const int32_t ASinTableC000[] =
00859             {
00860                                                                                          0x00221339,
00861             0x00228d9c,0x002309a5,0x0023876b,0x00240706,0x00248890,0x00250c26,0x002591e7,0x002619f5,
00862             0x0026a476,0x00273194,0x0027c17d,0x00285465,0x0028ea85,0x00298420,0x002a217e,0x002ac2f5,
00863             0x002b68e6,0x002c13c1,0x002cc40b,0x002d7a5e,0x002e3777,0x002efc36,0x002fc9b5,0x0030a152,
00864             0x003184d3,0x0032768f,0x003379c1
00865             };
00866         r = Interpolate(ASinTableC000,n-0xc000,9);
00867         }
00868     else if(n<0xfe00)
00869         {
00870         static const int32_t ASinTableF200[] =
00871             {
00872                                              0x003238a2,0x0032768f,0x0032b592,0x0032f5ba,0x00333719,
00873             0x003379c1,0x0033bdc8,0x00340345,0x00344a52,0x0034930c,0x0034dd94,0x00352a10,0x003578a9,
00874             0x0035c991,0x00361d01,0x0036733a,0x0036cc8b,0x00372953,0x00378a02,0x0037ef25,0x0038596e,
00875             0x0038c9be,0x00394144,0x0039c19e,0x003a4d1f,0x003ae75a,0x003b9654
00876             };
00877         r = Interpolate(ASinTableF200,n-0xF200,7);
00878         }
00879     else if(n<=0xffe0)
00880         {
00881         static const int32_t ASinTableFE00[] =
00882             {
00883                                                                                          0x003ad319,
00884             0x003ae75a,0x003afbed,0x003b10d5,0x003b2617,0x003b3bb7,0x003b51ba,0x003b6827,0x003b7f02,
00885             0x003b9654,0x003bae22,0x003bc677,0x003bdf5a,0x003bf8d6,0x003c12f8,0x003c2dcb,0x003c4960,
00886             0x003c65c7,0x003c8314,0x003ca160,0x003cc0c5,0x003ce165,0x003d0369,0x003d2702,0x003d4c6e,
00887             0x003d73ff,0x003d9e1e,0x003dcb5f,0x003dfc94,0x003e3300,0x003e70c5,0x003eba0a,0x003f1984
00888             };
00889         r = Interpolate(ASinTableFE00,n-0xFE00,4);
00890         }
00891     else if(n<=0x10000)
00892         {
00893         static const int32_t ASinTableFFE1[] =
00894             {
00895                        0x003ebf2c,0x003ec464,0x003ec9b2,0x003ecf17,0x003ed496,0x003eda2f,0x003edfe4,
00896             0x003ee5b6,0x003eeba8,0x003ef1bb,0x003ef7f2,0x003efe4f,0x003f04d5,0x003f0b88,0x003f126c,
00897             0x003f1984,0x003f20d5,0x003f2867,0x003f303e,0x003f3865,0x003f40e5,0x003f49c9,0x003f5323,
00898             0x003f5d06,0x003f678d,0x003f72dc,0x003f7f28,0x003f8cc2,0x003f9c33,0x003fae83,0x003fc661,
00899             0x00400000,
00900             };
00901         r = ASinTableFFE1[n-0xffe1];
00902         }
00903     else
00904         {
00905         // value out of range, return PI/2
00906         r = 0x00400000;
00907         }
00908 
00909     // scale result to 16 bits (from the 24 bit precision used in lookup table)
00910     r = (r+(1<<(8-1)))>>8;
00911 
00912     // produce correct sign for result
00913     if(value<0)
00914         r = -r;
00915 
00916     return r;
00917     }
00918 
00919 fixangle Fix::ATan(fix value)
00920     {
00921     unsigned n = value;
00922     if(value<0)
00923         n = (unsigned)-(int)n;
00924 
00925     // calculare arc-tangent by interpolation of lookup table values
00926     static const int32_t ATanTable[] =
00927         {
00928         0xfffeba28,
00929         0x00000000,0x000145d8,0x00028b0d,0x0003cf00,0x00051112,0x000650ad,0x00078d40,0x0008c643,
00930         0x0009fb38,0x000b2bac,0x000c5734,0x000d7d75,0x000e9e1d,0x000fb8e7,0x0010cd99,0x0011dc04,
00931         0x0012e405,0x0013e582,0x0014e06a,0x0015d4b7,0x0016c267,0x0017a982,0x00188a16,0x00196434,
00932         0x001a37f6,0x001b0575,0x001bccd2,0x001c8e2d,0x001d49ab,0x001dff72,0x001eafa7,0x001f5a74,
00933         0x00200000,0x0020a074,
00934         };
00935     int r;
00936     if(n<=0x10000)
00937         {
00938         r = Interpolate(ATanTable,n,11);
00939         }
00940     else
00941         {
00942         // For n>1 use the identity ATAN(n) = PI/2-ATAN(1/n)
00943         n = (unsigned)-(int)(n>>1)/(unsigned)n+1;   // n = 1/n
00944         r = Interpolate(ATanTable,n,11);
00945         r = 0x400000-r; // r = PI/2-r
00946         }
00947 
00948     // scale result to 16 bits (from the 24 bit precision used in lookup table)
00949     r = (r+(1<<(8-1)))>>8;
00950 
00951     if(value<0)
00952         r = -r;
00953 
00954     return r;
00955     }
00956 
00957 
00958 fix Fix::Random(uint32_t& seed)
00959     {
00960     return seed=(seed*69069+1);
00961     }
00962 
00963 
00964 ufix Fix::Random(uint32_t& seed,ufix range)
00965     {
00966     // generate next random number
00967     ufix n = seed=(seed*69069+1);
00968 
00969     // multiply range by random number (we won't bother with carry from low order half
00970     // of 64bit result because we don't really need acuracy for random numbers.)
00971     ufix lo = range&0xffff;
00972     ufix hi = range>>16;
00973     ufix r = hi*(n&0xffff)>>16;
00974     n >>= 16;
00975     r += lo*n>>16;
00976     r += hi*n;
00977 
00978     return r;
00979     }
00980 
00981 

Generated by  doxygen 1.6.1