Changeset 304
- Timestamp:
- 05/02/07 22:55:04 (18 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
libmpc/trunk/libmpcpsy/fft_routines.c
r137 r304 198 198 PowSpec256 ( const float* x, float* erg ) 199 199 { 200 const float* win = Hann_256; 201 float* aix = a; 202 int i; 203 200 int i; 204 201 // windowing 205 i = 256; 206 while (i--) 207 *aix++ = *x++ * *win++; 208 209 // perform FFT 210 rdft ( 256, a, ip, w ); 202 for( i = 0; i < 256; ++i ) 203 a[i] = x[i] * Hann_256[i]; 204 205 rdft(256, a, ip, w); // perform FFT 211 206 212 207 // calculate power 213 aix = a; // reset pointer 214 i = 128; 215 while (i--) { 216 *erg++ = aix[0]*aix[0] + aix[1]*aix[1]; 217 aix += 2; 218 } 208 for( i = 0; i < 128; ++i) 209 erg[i] = a[i*2] * a[i*2] + a[i*2+1] * a[i*2+1]; 219 210 } 220 211 … … 224 215 PowSpec1024 ( const float* x, float* erg ) 225 216 { 226 const float* win = Hann_1024; 227 float* aix = a; 228 int i; 229 230 i = 1024; // windowing 231 while (i--) 232 *aix++ = *x++ * *win++; 233 234 // for (i=0; i<1024; i++) 235 // a[i] = Hann_1024[i] * ((i==0 ? 0 : i-512) + 1000); 236 237 rdft ( 1024, a, ip, w ); // perform FFT 238 239 aix = a; // calculate power 240 i = 512; 241 242 243 DECONV; 244 // for (i = 0; i <= 512; i++ ) 245 // printf ("%3u %12.6f %12.6f\n", i, a[i+i], a[i+i+1]); 246 // exit(1); 247 while (i--) { 248 *erg++ = aix[0]*aix[0] + aix[1]*aix[1]; 249 aix += 2; 250 } 217 int i; 218 // windowing 219 for( i = 0; i < 1024; ++i ) 220 a[i] = x[i] * Hann_1024[i]; 221 222 rdft(1024, a, ip, w); // perform FFT 223 224 // calculate power 225 for( i = 0; i < 512; ++i) 226 erg[i] = a[i*2] * a[i*2] + a[i*2+1] * a[i*2+1]; 251 227 } 252 228 … … 256 232 PowSpec2048 ( const float* x, float* erg ) 257 233 { 258 const float* win = Hann_1600; 259 float* aix = a; 260 int i; 261 234 int i; 262 235 // windowing (only 1600 samples available -> centered in 2048!) 263 memset ( a , 0, 224*sizeof(*a) ); 264 aix = a + 224; 265 i = 1600; 266 while (i--) 267 *aix++ = *x++ * *win++; 268 memset ( a+1824, 0, 224*sizeof(*a) ); 269 270 rdft ( 2048, a, ip, w ); // perform FFT 271 272 aix = a; // calculate power 273 i = 1024; 274 while (i--) { 275 *erg++ = aix[0]*aix[0] + aix[1]*aix[1]; 276 aix += 2; 277 } 236 memset(a, 0, 224 * sizeof *a); 237 for( i = 0; i < 1600; ++i ) 238 a[i+224] = x[i] * Hann_1600[i]; 239 memset(a + 1824, 0, 224 * sizeof *a); 240 241 rdft(2048, a, ip, w); // perform FFT 242 243 // calculate power 244 for( i = 0; i < 1024; ++i) 245 erg[i] = a[i*2] * a[i*2] + a[i*2+1] * a[i*2+1]; 278 246 } 279 247 … … 283 251 PolarSpec1024 ( const float* x, float* erg, float* phs ) 284 252 { 285 const float* win = Hann_1024; 286 float* aix = a; 287 int i; 288 289 i = 1024; // windowing 290 while (i--) 291 *aix++ = *x++ * *win++; 292 293 rdft ( 1024, a, ip, w ); // perform FFT 253 int i; 254 for( i = 0; i < 1024; i++ ) 255 a[i] = x[i] * Hann_1024[i]; 256 257 rdft( 1024, a, ip, w); // perform FFT 294 258 295 259 // calculate power and phase 296 aix = a; // reset pointer 297 i = 512; 298 while (i--) { 299 *erg++ = aix[0]*aix[0] + aix[1]*aix[1]; 300 *phs++ = ATAN2F (aix[1], aix[0]); 301 aix += 2; 260 for( i = 0; i < 512; ++i ) 261 { 262 erg[i] = a[i*2] * a[i*2] + a[i*2+1] * a[i*2+1]; 263 phs[i] = ATAN2F( a[i*2+1], a[i*2] ); 302 264 } 303 265 } … … 308 270 Cepstrum2048 ( float* cep, const int MaxLine ) 309 271 { 310 float* aix = cep; 311 float* bix = cep + 2048; 312 int i; 313 272 int i, j; 314 273 // generate real, even spectrum (symmetric around 1024, cep[2048-i] = cep[i]) 315 for ( i = 0; i < 1024; i++ ) 316 *bix-- = *aix++; 317 318 // perform IFFT 319 rdft ( 2048, cep, ip, w ); 320 274 for( i = 0, j = 1024; i < 1024; ++i, --j ) 275 cep[1024 + j] = cep[i]; 276 277 rdft(2048, cep, ip, w); 278 321 279 // only real part as outcome (all even indexes of cep[]) 322 aix = cep; 323 bix = cep; 324 i = MaxLine + 1; 325 while (i--) { 326 *aix = *bix * (float) (0.9888 / 2048.); 327 // *aix = *bix * 0.0004828125f; 328 aix ++; 329 bix += 2; 330 } 331 } 280 for( i = 0; i < MaxLine + 1; ++i ) 281 cep[i] = cep[i*2] * (float) (0.9888 / 2048); 282 }
Note: See TracChangeset
for help on using the changeset viewer.