00001
00002
00003
00004
00005
00006
00007 #ifndef FMATH_H
00008 #define FMATH_H
00009
00010 #include <math.h>
00011
00012 #ifdef HAVE_FLOAT_H
00013 #include <float.h>
00014 #endif
00015
00016 namespace FD {
00017 #ifdef WIN32
00018 #define for if(0);else for
00019 #endif
00020
00021 union FloatManip {
00022 float f;
00023 unsigned int i;
00024 };
00025
00026 #ifndef M_LN2
00027 #define M_LN2 0.69314718055994530942
00028 #endif
00029
00030
00031
00032
00033
00034 #define FLOGLOOKUP2SIZE 256
00035 #define FLOGLOOKUP2SHIFT 15
00036
00037
00038
00039
00040
00041
00042
00043 extern float logtable2[FLOGLOOKUP2SIZE];
00044
00045 inline void build_flog_table()
00046 {
00047 static bool init=false;
00048 if (!init)
00049 {
00050 FloatManip m;
00051
00052 for (int i=0;i<FLOGLOOKUP2SIZE;i++)
00053 {
00054 m.i = (i<<FLOGLOOKUP2SHIFT) | 0x3f800000;
00055 logtable2[i]=log(m.f);
00056 }
00057 init=true;
00058 }
00059 }
00060
00061
00062
00063 inline float flog(float f)
00064 {
00065 build_flog_table();
00066 FloatManip m;
00067 m.f = f;
00068
00069 unsigned int id1 = m.i>>23;
00070
00071 unsigned int id2 = (m.i & 0x007fffff)>>FLOGLOOKUP2SHIFT;
00072 float f2=m.f;
00073
00074 m.i &= 0xffff8000;
00075 return (int(id1)-127)*M_LN2 + logtable2[id2] + (f2-m.f)/f2;
00076 }
00077
00078
00079 inline void fflogv(const float *fin, float *fout, int len)
00080 {
00081 build_flog_table();
00082 FloatManip m;
00083 for (int i=0;i<len;i++)
00084 {
00085 m.f = fin[i];
00086
00087 float tb = logtable2[(m.i & 0x007fffff)>>FLOGLOOKUP2SHIFT];
00088
00089 int id = (m.i>>23)-127;
00090
00091 fout[i] = id*M_LN2 + tb;
00092 }
00093 }
00094
00095
00096
00097
00098
00099
00100 #define FEXPSHIFT 22
00101 #define FEXPSIZE 1024
00102 #define FEXPMASK 0xffc00000
00103
00104 extern float exptable[FEXPSIZE];
00105
00106 inline void build_fexp_table()
00107 {
00108 static bool init=false;
00109 if (!init)
00110 {
00111 FloatManip m;
00112
00113 for (int i=0;i<FEXPSIZE;i++)
00114 {
00115 m.i = i<<FEXPSHIFT;
00116 exptable[i]=exp(m.f);
00117 }
00118 init=true;
00119 }
00120 }
00121
00122
00123
00124
00125 inline float fexp(float f)
00126 {
00127 build_fexp_table();
00128 FloatManip m;
00129 m.f = f;
00130
00131 unsigned int id = m.i>>FEXPSHIFT;
00132 float val = exptable[id];
00133 float f2;
00134
00135
00136
00137 f2=m.f;
00138 m.i &= FEXPMASK;
00139 m.f=f2-m.f;
00140 id = m.i>>FEXPSHIFT;
00141 val *= exptable[id];
00142
00143
00144 f2=m.f;
00145 m.i &= FEXPMASK;
00146 m.f=f2-m.f;
00147
00148
00149 return val + ( val*m.f * (1+.5*m.f) );
00150 }
00151
00152 }
00153 #endif