Mercurial > mplayer.hg
view libfaad2/mdct.c @ 22930:8a38e7270fe7
documented -mpegopts :interleaving2
author | nicodvb |
---|---|
date | Sun, 08 Apr 2007 16:21:32 +0000 |
parents | 59b6fa5b4201 |
children | e83eef58b30a |
line wrap: on
line source
/* ** FAAD2 - Freeware Advanced Audio (AAC) Decoder including SBR decoding ** Copyright (C) 2003-2004 M. Bakker, Ahead Software AG, http://www.nero.com ** ** This program is free software; you can redistribute it and/or modify ** it under the terms of the GNU General Public License as published by ** the Free Software Foundation; either version 2 of the License, or ** (at your option) any later version. ** ** This program is distributed in the hope that it will be useful, ** but WITHOUT ANY WARRANTY; without even the implied warranty of ** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ** GNU General Public License for more details. ** ** You should have received a copy of the GNU General Public License ** along with this program; if not, write to the Free Software ** Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. ** ** Any non-GPL usage of this software or parts of this software is strictly ** forbidden. ** ** Commercial non-GPL licensing of this software is possible. ** For more info contact Ahead Software through Mpeg4AAClicense@nero.com. ** ** $Id: mdct.c,v 1.43 2004/09/04 14:56:28 menno Exp $ **/ /* * Fast (I)MDCT Implementation using (I)FFT ((Inverse) Fast Fourier Transform) * and consists of three steps: pre-(I)FFT complex multiplication, complex * (I)FFT, post-(I)FFT complex multiplication, * * As described in: * P. Duhamel, Y. Mahieux, and J.P. Petit, "A Fast Algorithm for the * Implementation of Filter Banks Based on 'Time Domain Aliasing * Cancellation’," IEEE Proc. on ICASSP‘91, 1991, pp. 2209-2212. * * * As of April 6th 2002 completely rewritten. * This (I)MDCT can now be used for any data size n, where n is divisible by 8. * */ #include "common.h" #include "structs.h" #include <stdlib.h> #ifdef _WIN32_WCE #define assert(x) #else #include <assert.h> #endif #include "cfft.h" #include "mdct.h" #include "mdct_tab.h" mdct_info *faad_mdct_init(uint16_t N) { mdct_info *mdct = (mdct_info*)faad_malloc(sizeof(mdct_info)); assert(N % 8 == 0); mdct->N = N; /* NOTE: For "small framelengths" in FIXED_POINT the coefficients need to be * scaled by sqrt("(nearest power of 2) > N" / N) */ /* RE(mdct->sincos[k]) = scale*(real_t)(cos(2.0*M_PI*(k+1./8.) / (real_t)N)); * IM(mdct->sincos[k]) = scale*(real_t)(sin(2.0*M_PI*(k+1./8.) / (real_t)N)); */ /* scale is 1 for fixed point, sqrt(N) for floating point */ switch (N) { case 2048: mdct->sincos = (complex_t*)mdct_tab_2048; break; case 256: mdct->sincos = (complex_t*)mdct_tab_256; break; #ifdef LD_DEC case 1024: mdct->sincos = (complex_t*)mdct_tab_1024; break; #endif #ifdef ALLOW_SMALL_FRAMELENGTH case 1920: mdct->sincos = (complex_t*)mdct_tab_1920; break; case 240: mdct->sincos = (complex_t*)mdct_tab_240; break; #ifdef LD_DEC case 960: mdct->sincos = (complex_t*)mdct_tab_960; break; #endif #endif #ifdef SSR_DEC case 512: mdct->sincos = (complex_t*)mdct_tab_512; break; case 64: mdct->sincos = (complex_t*)mdct_tab_64; break; #endif } /* initialise fft */ mdct->cfft = cffti(N/4); #ifdef PROFILE mdct->cycles = 0; mdct->fft_cycles = 0; #endif return mdct; } void faad_mdct_end(mdct_info *mdct) { if (mdct != NULL) { #ifdef PROFILE printf("MDCT[%.4d]: %I64d cycles\n", mdct->N, mdct->cycles); printf("CFFT[%.4d]: %I64d cycles\n", mdct->N/4, mdct->fft_cycles); #endif cfftu(mdct->cfft); faad_free(mdct); } } void faad_imdct(mdct_info *mdct, real_t *X_in, real_t *X_out) { uint16_t k; complex_t x; #ifdef ALLOW_SMALL_FRAMELENGTH #ifdef FIXED_POINT real_t scale, b_scale = 0; #endif #endif ALIGN complex_t Z1[512]; complex_t *sincos = mdct->sincos; uint16_t N = mdct->N; uint16_t N2 = N >> 1; uint16_t N4 = N >> 2; uint16_t N8 = N >> 3; #ifdef PROFILE int64_t count1, count2 = faad_get_ts(); #endif #ifdef ALLOW_SMALL_FRAMELENGTH #ifdef FIXED_POINT /* detect non-power of 2 */ if (N & (N-1)) { /* adjust scale for non-power of 2 MDCT */ /* 2048/1920 */ b_scale = 1; scale = COEF_CONST(1.0666666666666667); } #endif #endif /* pre-IFFT complex multiplication */ for (k = 0; k < N4; k++) { ComplexMult(&IM(Z1[k]), &RE(Z1[k]), X_in[2*k], X_in[N2 - 1 - 2*k], RE(sincos[k]), IM(sincos[k])); } #ifdef PROFILE count1 = faad_get_ts(); #endif /* complex IFFT, any non-scaling FFT can be used here */ cfftb(mdct->cfft, Z1); #ifdef PROFILE count1 = faad_get_ts() - count1; #endif /* post-IFFT complex multiplication */ for (k = 0; k < N4; k++) { RE(x) = RE(Z1[k]); IM(x) = IM(Z1[k]); ComplexMult(&IM(Z1[k]), &RE(Z1[k]), IM(x), RE(x), RE(sincos[k]), IM(sincos[k])); #ifdef ALLOW_SMALL_FRAMELENGTH #ifdef FIXED_POINT /* non-power of 2 MDCT scaling */ if (b_scale) { RE(Z1[k]) = MUL_C(RE(Z1[k]), scale); IM(Z1[k]) = MUL_C(IM(Z1[k]), scale); } #endif #endif } /* reordering */ for (k = 0; k < N8; k+=2) { X_out[ 2*k] = IM(Z1[N8 + k]); X_out[ 2 + 2*k] = IM(Z1[N8 + 1 + k]); X_out[ 1 + 2*k] = -RE(Z1[N8 - 1 - k]); X_out[ 3 + 2*k] = -RE(Z1[N8 - 2 - k]); X_out[N4 + 2*k] = RE(Z1[ k]); X_out[N4 + + 2 + 2*k] = RE(Z1[ 1 + k]); X_out[N4 + 1 + 2*k] = -IM(Z1[N4 - 1 - k]); X_out[N4 + 3 + 2*k] = -IM(Z1[N4 - 2 - k]); X_out[N2 + 2*k] = RE(Z1[N8 + k]); X_out[N2 + + 2 + 2*k] = RE(Z1[N8 + 1 + k]); X_out[N2 + 1 + 2*k] = -IM(Z1[N8 - 1 - k]); X_out[N2 + 3 + 2*k] = -IM(Z1[N8 - 2 - k]); X_out[N2 + N4 + 2*k] = -IM(Z1[ k]); X_out[N2 + N4 + 2 + 2*k] = -IM(Z1[ 1 + k]); X_out[N2 + N4 + 1 + 2*k] = RE(Z1[N4 - 1 - k]); X_out[N2 + N4 + 3 + 2*k] = RE(Z1[N4 - 2 - k]); } #ifdef PROFILE count2 = faad_get_ts() - count2; mdct->fft_cycles += count1; mdct->cycles += (count2 - count1); #endif } #ifdef LTP_DEC void faad_mdct(mdct_info *mdct, real_t *X_in, real_t *X_out) { uint16_t k; complex_t x; ALIGN complex_t Z1[512]; complex_t *sincos = mdct->sincos; uint16_t N = mdct->N; uint16_t N2 = N >> 1; uint16_t N4 = N >> 2; uint16_t N8 = N >> 3; #ifndef FIXED_POINT real_t scale = REAL_CONST(N); #else real_t scale = REAL_CONST(4.0/N); #endif #ifdef ALLOW_SMALL_FRAMELENGTH #ifdef FIXED_POINT /* detect non-power of 2 */ if (N & (N-1)) { /* adjust scale for non-power of 2 MDCT */ /* *= sqrt(2048/1920) */ scale = MUL_C(scale, COEF_CONST(1.0327955589886444)); } #endif #endif /* pre-FFT complex multiplication */ for (k = 0; k < N8; k++) { uint16_t n = k << 1; RE(x) = X_in[N - N4 - 1 - n] + X_in[N - N4 + n]; IM(x) = X_in[ N4 + n] - X_in[ N4 - 1 - n]; ComplexMult(&RE(Z1[k]), &IM(Z1[k]), RE(x), IM(x), RE(sincos[k]), IM(sincos[k])); RE(Z1[k]) = MUL_R(RE(Z1[k]), scale); IM(Z1[k]) = MUL_R(IM(Z1[k]), scale); RE(x) = X_in[N2 - 1 - n] - X_in[ n]; IM(x) = X_in[N2 + n] + X_in[N - 1 - n]; ComplexMult(&RE(Z1[k + N8]), &IM(Z1[k + N8]), RE(x), IM(x), RE(sincos[k + N8]), IM(sincos[k + N8])); RE(Z1[k + N8]) = MUL_R(RE(Z1[k + N8]), scale); IM(Z1[k + N8]) = MUL_R(IM(Z1[k + N8]), scale); } /* complex FFT, any non-scaling FFT can be used here */ cfftf(mdct->cfft, Z1); /* post-FFT complex multiplication */ for (k = 0; k < N4; k++) { uint16_t n = k << 1; ComplexMult(&RE(x), &IM(x), RE(Z1[k]), IM(Z1[k]), RE(sincos[k]), IM(sincos[k])); X_out[ n] = -RE(x); X_out[N2 - 1 - n] = IM(x); X_out[N2 + n] = -IM(x); X_out[N - 1 - n] = RE(x); } } #endif