Changeset 109 for mppenc/trunk
- Timestamp:
- 11/10/06 22:04:46 (18 years ago)
- Location:
- mppenc/trunk
- Files:
-
- 5 edited
-
CMakeLists.txt (modified) (1 diff)
-
src/cvd.c (modified) (1 diff)
-
src/fastmath.c (modified) (3 diffs)
-
src/fastmath.h (modified) (5 diffs)
-
src/quant.c (modified) (6 diffs)
Legend:
- Unmodified
- Added
- Removed
-
mppenc/trunk/CMakeLists.txt
r105 r109 4 4 5 5 if(NOT MSVC) 6 set(CMAKE_C_FLAGS "- fno-strict-aliasing -Os-fomit-frame-pointer -pipe")6 set(CMAKE_C_FLAGS "-Wstrict-aliasing -fstrict-aliasing -ffast-math -O3 -fomit-frame-pointer -pipe") 7 7 endif(NOT MSVC) 8 8 -
mppenc/trunk/src/cvd.c
r97 r109 219 219 logfast ( float x ) 220 220 { 221 double y = x * x; 222 y *= y; 223 y *= y; 224 return (((int*)(&y))[1] + (45127.5 - 1072693248.)) * ( M_LN2 / (1L<<23) ); 221 union { double d; Int32_t n[2]; } tmp; 222 tmp.d = x * x; 223 tmp.d *= tmp.d; 224 tmp.d *= tmp.d; 225 #if ENDIAN == HAVE_LITTLE_ENDIAN 226 return (tmp.n[1] + (45127.5 - 1072693248.)) * ( M_LN2 / (1L<<23) ); 227 #else 228 return (tmp.n[0] + (45127.5 - 1072693248.)) * ( M_LN2 / (1L<<23) ); 229 #endif 225 230 } 226 231 -
mppenc/trunk/src/fastmath.c
r97 r109 19 19 20 20 #include "mppenc.h" 21 #include "fastmath.h" 21 22 22 23 #ifdef FAST_MATH … … 30 31 void Init_FastMath ( void ) 31 32 { 32 int i; 33 float X; 34 float Y; 35 double xm; 36 double x0; 37 double xp; 38 double x; 39 double y; 40 float* p; 33 int i; mpc_floatint X, Y; double xm, x0, xp, x, y; float* p; 41 34 42 35 p = (float*) tabatan2; … … 64 57 p = (float*) tabsqrt_ex; 65 58 for ( i = 0; i < 255; i++ ) { 66 *(int*)&X= (i << 23);67 *(int*)&Y= (i << 23) + (1<<23) - 1;68 *p++ = sqrt(X );59 X.n = (i << 23); 60 Y.n = (i << 23) + (1<<23) - 1; 61 *p++ = sqrt(X.f); 69 62 } 70 *(int*)&X= (255 << 23) - 1;71 *p++ = sqrt(X );63 X.n = (255 << 23) - 1; 64 *p++ = sqrt(X.f); 72 65 73 66 p = (float*) tabsqrt_m; -
mppenc/trunk/src/fastmath.h
r97 r109 18 18 */ 19 19 20 #if 1 21 # define ROUND32(x) ( floattmp = (x) + (int)0x00FD8000L, *(int*)(&floattmp) - (int)0x4B7D8000L ) 22 #else 23 # define ROUND32(x) ( (int) floor ((x) + 0.5) ) 24 #endif 20 typedef union mpc_floatint 21 { 22 float f; 23 Int32_t n; 24 } mpc_floatint; 25 26 static __inline Int32_t mpc_nearbyintf(float fVal) 27 { 28 mpc_floatint tmp; 29 tmp.f = fVal + 0x00FF8000; 30 tmp.n = tmp.n - 0x4B7F8000; 31 return tmp.n; 32 } 33 #define mpc_lrintf mpc_nearbyintf 34 static __inline Int32_t mpc_round32(float fVal) 35 { 36 mpc_floatint tmp; 37 tmp.f = fVal + 0x00FD8000; 38 tmp.n = tmp.n - 0x4B7D8000; 39 return tmp.n; 40 } 25 41 26 42 #ifdef FAST_MATH … … 29 45 my_atan2 ( float x, float y ) 30 46 { 31 float t; 32 int i; 33 float ret; 34 float floattmp; 47 float t, ret; int i; mpc_floatint mx, my; 35 48 36 if ( (*(int*)&x & 0x7FFFFFFF) < (*(int*)&y & 0x7FFFFFFF) ) { 37 i = ROUND32 (t = TABSTEP * (x / y)); 49 mx.f = x; 50 my.f = y; 51 if ( (mx.n & 0x7FFFFFFF) < (my.n & 0x7FFFFFFF) ) { 52 i = mpc_round32 (t = TABSTEP * (mx.f / my.f)); 38 53 ret = tabatan2 [1*TABSTEP+i][0] + tabatan2 [1*TABSTEP+i][1] * (t-i); 39 if ( *(int*)&y< 0 )54 if ( my.n < 0 ) 40 55 ret = (float)(ret - M_PI); 41 56 } 42 else if ( *(int*)&x < 0) {43 i = ROUND32 (t = TABSTEP * (y / x));57 else if ( mx.n < 0 ) { 58 i = mpc_round32 (t = TABSTEP * (my.f / mx.f)); 44 59 ret = - M_PI/2 - tabatan2 [1*TABSTEP+i][0] + tabatan2 [1*TABSTEP+i][1] * (i-t); 45 60 } 46 else if ( *(int*)&x > 0) {47 i = ROUND32 (t = TABSTEP * (y / x));61 else if ( mx.n > 0 ) { 62 i = mpc_round32 (t = TABSTEP * (my.f / mx.f)); 48 63 ret = + M_PI/2 - tabatan2 [1*TABSTEP+i][0] + tabatan2 [1*TABSTEP+i][1] * (i-t); 49 64 } … … 58 73 my_cos ( float x ) 59 74 { 60 float t; 61 int i; 62 float ret; 63 float floattmp; 64 65 i = ROUND32 (t = TABSTEP * x); 75 float t, ret; int i; 76 i = mpc_round32 (t = TABSTEP * x); 66 77 ret = tabcos [13*TABSTEP+i][0] + tabcos [13*TABSTEP+i][1] * (t-i); 67 78 return ret; … … 72 83 my_ifloor ( float x ) 73 84 { 74 x = x + (0x0C00000L + 0.500000001); 75 return *(int*)&x - 1262485505; 85 mpc_floatint mx; 86 mx.f = (float) (x + (0x0C00000L + 0.500000001)); 87 return mx.n - 1262485505; 76 88 } 77 89 … … 80 92 my_sqrt ( float x ) 81 93 { 82 float ret; 83 int i; 84 int ex = *(int*)&x >> 23; // get the exponent 85 float floattmp; 86 87 *(int*)&x = (*(int*)&x & 0x7FFFFF) | 0x42800000; // delete the exponent 88 i = ROUND32 (x); // Integer-part of the mantissa (round ????????????) 89 ret = tabsqrt_m [i-TABSTEP][0] + tabsqrt_m [i-TABSTEP][1] * (x-i); // calculate value 94 float ret; int i, ex; mpc_floatint mx; 95 mx.f = x; 96 ex = mx.n >> 23; // get the exponent 97 mx.n = (mx.n & 0x7FFFFF) | 0x42800000; // delete the exponent 98 i = mpc_round32 (mx.f); // Integer-part of the mantissa (round ????????????) 99 ret = tabsqrt_m [i-TABSTEP][0] + tabsqrt_m [i-TABSTEP][1] * (mx.f-i); // calculate value 90 100 ret *= tabsqrt_ex [ex]; 91 101 return ret; -
mppenc/trunk/src/quant.c
r97 r109 19 19 20 20 #include "mppenc.h" 21 #include "fastmath.h" 21 22 22 23 /* V A R I A B L E S */ … … 149 150 ISNR_Schaetzer ( const float* input, const float SNRcomp, const int res ) 150 151 { 151 int k; 152 float fac = A [res]; 153 float invfac = C [res]; 154 float Signal = 1.e-30f; 155 float Fehler = 1.e-30f; 156 float tmp ; 157 float tmp2; 158 float tmp3; 152 int k; float signal, error, sig, err; 153 const float fac = A [res]; 154 const float invfac = C [res]; 155 const float noiseComp = NoiseInjectionCompensation1D [res]; 159 156 160 157 // Summation of the absolute power and the quadratic error 161 for ( k = 0; k < 36; k++ ) { 162 tmp2 = input[k] * NoiseInjectionCompensation1D [res]; 163 // q = ftol(in), correct rounding 164 tmp = tmp2 * fac + 0xFF8000; 165 tmp3 = (*(int*) & tmp - 0x4B7F8000) * invfac; 166 tmp = tmp3 - tmp2; 167 168 Fehler += tmp * tmp; 169 Signal += tmp2 * tmp2; 158 signal = error = 1.e-30f; 159 for ( k = 0; k < 36; k++ ) 160 { 161 sig = input[k] * noiseComp; 162 err = mpc_nearbyintf(sig * fac) * invfac - sig; 163 signal += sig * sig; 164 error += err * err; 170 165 } 171 166 172 167 // Utilization of SNRcomp only if SNR > 1 !!! 173 return Signal > Fehler ? Fehler / (SNRcomp * Signal) : Fehler / Signal;168 return signal > error ? error / (signal * SNRcomp) : error / signal; 174 169 } 175 170 … … 178 173 ISNR_Schaetzer_Trans ( const float* input, const float SNRcomp, const int res ) 179 174 { 180 int k; 181 float fac = A [res]; 182 float invfac = C [res]; 183 float Signal; 184 float Fehler; 185 float ret ; 186 float tmp ; 187 float tmp2; 188 float tmp3; 189 175 int k; float signal, error, result, sig, err; 176 const float fac = A [res]; 177 const float invfac = C [res]; 178 const float noiseComp = NoiseInjectionCompensation1D [res]; 179 190 180 // Summation of the absolute power and the quadratic error 191 k = 0; 192 Signal = Fehler = 1.e-30f; 193 for ( ; k < 12; k++ ) { 194 tmp2 = input[k] * NoiseInjectionCompensation1D [res]; 195 // q = ftol(in), correct rounding 196 tmp = tmp2 * fac + 0xFF8000; 197 tmp3 = (*(int*) & tmp - 0x4B7F8000) * invfac; 198 tmp = tmp3 - tmp2; 199 200 Fehler += tmp * tmp; 201 Signal += tmp2 * tmp2; 202 } 203 tmp = Signal > Fehler ? Fehler / (SNRcomp * Signal) : Fehler / Signal; 204 ret = tmp; 205 Signal = Fehler = 1.e-30f; 206 for ( ; k < 24; k++ ) { 207 tmp2 = input[k] * NoiseInjectionCompensation1D [res]; 208 // q = ftol(in), correct rounding 209 tmp = tmp2 * fac + 0xFF8000; 210 tmp3 = (*(int*) & tmp - 0x4B7F8000) * invfac; 211 tmp = tmp3 - tmp2; 212 213 Fehler += tmp * tmp; 214 Signal += tmp2 * tmp2; 215 } 216 tmp = Signal > Fehler ? Fehler / (SNRcomp * Signal) : Fehler / Signal; 217 if ( tmp > ret ) ret = tmp; 218 //ret += tmp; 219 Signal = Fehler = 1.e-30f; 220 for ( ; k < 36; k++ ) { 221 tmp2 = input[k] * NoiseInjectionCompensation1D [res]; 222 // q = ftol(in), correct rounding 223 tmp = tmp2 * fac + 0xFF8000; 224 tmp3 = (*(int*) & tmp - 0x4B7F8000) * invfac; 225 tmp = tmp3 - tmp2; 226 227 Fehler += tmp * tmp; 228 Signal += tmp2 * tmp2; 229 } 230 tmp = Signal > Fehler ? Fehler / (SNRcomp * Signal) : Fehler / Signal; 231 if ( tmp > ret ) ret = tmp; 232 //ret += tmp; 233 //ret *= 0.33333333333f; 234 235 return ret; 181 signal = error = 1.e-30f; 182 for ( k = 0 ; k < 12; k++ ) 183 { 184 sig = input[k] * noiseComp; 185 err = mpc_nearbyintf(sig * fac) * invfac - sig; 186 signal += sig * sig; 187 error += err * err; 188 } 189 result = signal > error ? error / (signal * SNRcomp) : error / signal; 190 191 signal = error = 1.e-30f; 192 for ( ; k < 24; k++ ) 193 { 194 sig = input[k] * noiseComp; 195 err = mpc_nearbyintf(sig * fac) * invfac - sig; 196 signal += sig * sig; 197 error += err * err; 198 } 199 sig = signal > error ? error / (signal * SNRcomp) : error / signal; 200 if ( sig > result ) result = sig; 201 202 signal = error = 1.e-30f; 203 for ( ; k < 36; k++ ) 204 { 205 sig = input[k] * noiseComp; 206 err = mpc_nearbyintf(sig * fac) * invfac - sig; 207 signal += sig * sig; 208 error += err * err; 209 } 210 sig = signal > error ? error / (signal * SNRcomp) : error / signal; 211 if ( sig > result ) result = sig; 212 213 return result; 236 214 } 237 215 … … 241 219 QuantizeSubband ( unsigned int* qu_output, const float* input, const int res, float* errors ) 242 220 { 243 int n; 244 int offset = D [res]; 245 float mult = A [res] * NoiseInjectionCompensation1D [res]; 246 float invmult = C [res]; 247 float tmp; 248 int quant; 249 float signal; 250 251 for ( n = 0; n < 36 - MAX_NS_ORDER; n++, input++, qu_output++ ) { 252 // q = ftol(in), correct rounding 253 tmp = *input * mult + 0xFF8000; 254 quant = (unsigned int)(*(int*) & tmp - 0x4B7F8000 + offset); 255 221 int n, quant; float signal; 222 const int offset = D [res]; 223 const float noiseComp = NoiseInjectionCompensation1D [res]; 224 const float mult = A [res] * noiseComp; 225 const float invmult = C [res]; 226 227 for ( n = 0; n < 36 - MAX_NS_ORDER; n++) 228 { 229 quant = (unsigned int) ( mpc_lrintf( input[n] * mult ) + offset ); 256 230 // limitation to 0...2D 257 if ((unsigned int)quant > (unsigned int)2*offset ) { 258 quant = mini ( quant, 2*offset ); 259 quant = maxi ( quant, 0 ); 231 if ( quant > offset * 2 ) 232 { 233 quant = mini ( quant, offset * 2 ); 234 quant = maxi ( quant, 0 ); 260 235 } 261 *qu_output = quant; 262 } 263 264 for ( ; n < 36; n++, input++, qu_output++ ) { 265 // q = ftol(in), correct rounding 266 signal = *input * mult; 267 tmp = signal + 0xFF8000; 268 quant = (unsigned int)(*(int*) & tmp - 0x4B7F8000 + offset); 236 qu_output[n] = quant; 237 } 238 239 for ( ; n < 36; n++ ) 240 { 241 signal = input[n] * mult; 242 quant = (unsigned int) ( mpc_lrintf( signal ) + offset ); 269 243 270 244 // calculate the current error and save it for error refeeding 271 errors [n + 6] = invmult * (quant - offset) - signal * NoiseInjectionCompensation1D [res];245 errors [n + 6] = invmult * (quant - offset) - signal * noiseComp; 272 246 273 247 // limitation to 0...2D 274 if ((unsigned int)quant > (unsigned int)2*offset ) { 275 quant = mini ( quant, 2*offset ); 276 quant = maxi ( quant, 0 ); 248 if ( quant > offset * 2 ) 249 { 250 quant = mini ( quant, offset * 2); 251 quant = maxi ( quant, 0 ); 277 252 } 278 *qu_output= quant;253 qu_output[n] = quant; 279 254 } 280 255 } … … 285 260 QuantizeSubbandWithNoiseShaping ( unsigned int* qu_output, const float* input, const int res, float* errors, const float* FIR ) 286 261 { 287 #define E(x) *((int*)errors+(x)) 288 289 float signal; 290 float tmp; 291 float mult = A [res]; 292 float invmult = C [res]; 293 int offset = D [res]; 294 int n; 295 int quant; 296 297 E(0) = E(1) = E(2) = E(3) = E(4) = E(5) = 0; // arghh, it produces pops on each frame boundary! 298 299 for ( n = 0; n < 36; n++, input++, qu_output++ ) { 300 signal = *input * NoiseInjectionCompensation1D [res] - (FIR[5]*errors[n+0] + FIR[4]*errors[n+1] + FIR[3]*errors[n+2] + FIR[2]*errors[n+3] + FIR[1]*errors[n+4] + FIR[0]*errors[n+5]); 301 302 // quant = ftol(signal), correct rounding 303 tmp = signal * mult + 0xFF8000; 304 quant = *(int*) & tmp - 0x4B7F8000; 262 int n, quant; float signal; 263 const float mult = A [res]; 264 const float invmult = C [res]; 265 const int offset = D [res]; 266 const float noiseComp = NoiseInjectionCompensation1D [res]; 267 268 memset(errors, 0, 6 * sizeof *errors); // arghh, it produces pops on each frame boundary! 269 270 for ( n = 0; n < 36; n++ ) 271 { 272 signal = input[n] * noiseComp - 273 (FIR[5]*errors[n+0] + FIR[4]*errors[n+1] + FIR[3]*errors[n+2] + 274 FIR[2]*errors[n+3] + FIR[1]*errors[n+4] + FIR[0]*errors[n+5]); 275 quant = mpc_lrintf(signal * mult); 305 276 306 277 // calculate the current error and save it for error refeeding 307 errors [n + 6] = invmult * quant - signal * NoiseInjectionCompensation1D [res];278 errors [n + 6] = invmult * quant - signal * noiseComp; 308 279 309 280 // limitation to +/-D … … 311 282 quant = maxf ( quant, -offset ); 312 283 313 *qu_output = (unsigned int)(quant + offset);284 qu_output[n] = (unsigned int) ( quant + offset ); 314 285 } 315 286 }
Note: See TracChangeset
for help on using the changeset viewer.
