X-Git-Url: http://git.tuebingen.mpg.de/?p=paraslash.git;a=blobdiff_plain;f=imdct.c;h=702d3d21327eb508703e1b43636d67ea088721da;hp=e2cfe03426eee48092b6973924e7e458a7a4fab3;hb=df6606e438a34c7bfc325e5f320843d1e97d3f7a;hpb=44aa9a29ecc951d74a28bd79462cb0f9b628523e diff --git a/imdct.c b/imdct.c index e2cfe034..702d3d21 100644 --- a/imdct.c +++ b/imdct.c @@ -29,21 +29,13 @@ typedef float fftsample_t; -#define DECLARE_ALIGNED(n,t,v) t v __attribute__ ((aligned (n))) -#define DECLARE_ALIGNED_16(t, v) DECLARE_ALIGNED(16, t, v) -#define M_SQRT1_2 0.70710678118654752440 /* 1/sqrt(2) */ - struct fft_complex { fftsample_t re, im; }; struct fft_context { int nbits; - int inverse; uint16_t *revtab; - struct fft_complex *exptab; - struct fft_complex *exptab1; /* only used by SSE code */ - struct fft_complex *tmp_buf; }; struct mdct_context { @@ -78,22 +70,22 @@ static fftsample_t *ff_cos_tabs[] = { ff_cos_16384, ff_cos_32768, ff_cos_65536, }; -static int split_radix_permutation(int i, int n, int inverse) +static int split_radix_permutation(int i, int n) { int m; if (n <= 2) return i & 1; m = n >> 1; - if (!(i & m)) - return split_radix_permutation(i, m, inverse) * 2; + if ((i & m) == 0) + return split_radix_permutation(i, m) * 2; m >>= 1; - if (inverse == !(i & m)) - return split_radix_permutation(i, m, inverse) * 4 + 1; + if ((i & m) == 0) + return split_radix_permutation(i, m) * 4 + 1; else - return split_radix_permutation(i, m, inverse) * 4 - 1; + return split_radix_permutation(i, m) * 4 - 1; } -#define sqrthalf (float)M_SQRT1_2 +#define SQRTHALF (float)0.70710678118654752440 /* 1/sqrt(2) */ #define BF(x,y,a,b) {\ x = a - b;\ @@ -203,7 +195,7 @@ static void fft8(struct fft_complex *z) BF(z[6].re, z[2].re, z[2].re, t7); BF(z[6].im, z[2].im, z[2].im, t8); - TRANSFORM(z[1], z[3], z[5], z[7], sqrthalf, sqrthalf); + TRANSFORM(z[1], z[3], z[5], z[7], SQRTHALF, SQRTHALF); } static void fft16(struct fft_complex *z) @@ -215,7 +207,7 @@ static void fft16(struct fft_complex *z) fft4(z + 12); TRANSFORM_ZERO(z[0], z[4], z[8], z[12]); - TRANSFORM(z[2], z[6], z[10], z[14], sqrthalf, sqrthalf); + TRANSFORM(z[2], z[6], z[10], z[14], SQRTHALF, SQRTHALF); TRANSFORM(z[1], z[5], z[9], z[13], ff_cos_16[1], ff_cos_16[3]); TRANSFORM(z[3], z[7], z[11], z[15], ff_cos_16[3], ff_cos_16[1]); } @@ -325,116 +317,34 @@ void imdct(struct mdct_context *s, float *output, const float *input) } } -static int fft_init(struct fft_context *s, int nbits, int inverse) +static int fft_init(struct fft_context *s, int nbits) { - int i, j, m, n; - float alpha, c1, s1, s2; - int split_radix = 1; + int i, j, n; if (nbits < 2 || nbits > 16) return -E_FFT_BAD_PARAMS; s->nbits = nbits; n = 1 << nbits; - s->tmp_buf = NULL; - s->exptab = para_malloc((n / 2) * sizeof(struct fft_complex)); s->revtab = para_malloc(n * sizeof(uint16_t)); - s->inverse = inverse; - - s2 = inverse ? 1.0 : -1.0; - - s->exptab1 = NULL; - - if (split_radix) { - for (j = 4; j <= nbits; j++) { - int k = 1 << j; - double freq = 2 * M_PI / k; - fftsample_t *tab = ff_cos_tabs[j - 4]; - for (i = 0; i <= k / 4; i++) - tab[i] = cos(i * freq); - for (i = 1; i < k / 4; i++) - tab[k / 2 - i] = tab[i]; - } - for (i = 0; i < n; i++) - s->revtab[-split_radix_permutation( - i, n, s->inverse) & (n - 1)] = i; - s->tmp_buf = para_malloc(n * sizeof(struct fft_complex)); - } else { - int np, nblocks, np2, l; - struct fft_complex *q; - - for (i = 0; i < (n / 2); i++) { - alpha = 2 * M_PI * (float) i / (float) n; - c1 = cos(alpha); - s1 = sin(alpha) * s2; - s->exptab[i].re = c1; - s->exptab[i].im = s1; - } - - np = 1 << nbits; - nblocks = np >> 3; - np2 = np >> 1; - s->exptab1 = para_malloc(np * 2 * sizeof(struct fft_complex)); - q = s->exptab1; - do { - for (l = 0; l < np2; l += 2 * nblocks) { - *q++ = s->exptab[l]; - *q++ = s->exptab[l + nblocks]; - - q->re = -s->exptab[l].im; - q->im = s->exptab[l].re; - q++; - q->re = -s->exptab[l + nblocks].im; - q->im = s->exptab[l + nblocks].re; - q++; - } - nblocks = nblocks >> 1; - } while (nblocks != 0); - freep(&s->exptab); - /* compute bit reverse table */ - for (i = 0; i < n; i++) { - m = 0; - for (j = 0; j < nbits; j++) - m |= ((i >> j) & 1) << (nbits - j - 1); - s->revtab[i] = m; - } + for (j = 4; j <= nbits; j++) { + int k = 1 << j; + double freq = 2 * M_PI / k; + fftsample_t *tab = ff_cos_tabs[j - 4]; + for (i = 0; i <= k / 4; i++) + tab[i] = cos(i * freq); + for (i = 1; i < k / 4; i++) + tab[k / 2 - i] = tab[i]; } - return 0; -} - -static void fft_end(struct fft_context *ctx) -{ - freep(&ctx->revtab); - freep(&ctx->exptab); - freep(&ctx->exptab1); - freep(&ctx->tmp_buf); -} - -DECLARE_ALIGNED(16, float, ff_sine_128[128]); -DECLARE_ALIGNED(16, float, ff_sine_256[256]); -DECLARE_ALIGNED(16, float, ff_sine_512[512]); -DECLARE_ALIGNED(16, float, ff_sine_1024[1024]); -DECLARE_ALIGNED(16, float, ff_sine_2048[2048]); -DECLARE_ALIGNED(16, float, ff_sine_4096[4096]); - -float *ff_sine_windows[6] = { - ff_sine_128, ff_sine_256, ff_sine_512, ff_sine_1024, - ff_sine_2048, ff_sine_4096 -}; - -// Generate a sine window. -void sine_window_init(float *window, int n) -{ - int i; - for (i = 0; i < n; i++) - window[i] = sinf((i + 0.5) * (M_PI / (2.0 * n))); + s->revtab[-split_radix_permutation(i, n) & (n - 1)] = i; + return 0; } /** - * Init MDCT or IMDCT computation. + * Initialize the inverse modified cosine transform. */ -int imdct_init(int nbits, int inverse, struct mdct_context **result) +int imdct_init(int nbits, struct mdct_context **result) { int ret, n, n4, i; double alpha; @@ -453,7 +363,7 @@ int imdct_init(int nbits, int inverse, struct mdct_context **result) s->tcos[i] = -cos(alpha); s->tsin[i] = -sin(alpha); } - ret = fft_init(&s->fft, s->nbits - 2, inverse); + ret = fft_init(&s->fft, s->nbits - 2); if (ret < 0) goto fail; *result = s; @@ -467,8 +377,8 @@ fail: void imdct_end(struct mdct_context *ctx) { - freep(&ctx->tcos); - freep(&ctx->tsin); - fft_end(&ctx->fft); + free(ctx->tcos); + free(ctx->tsin); + free(ctx->fft.revtab); free(ctx); }