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 {
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;\
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)
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]);
}
}
}
-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.
+ *
+ * \param nbits The number of bits to use (4 <= \a nbits <= 18).
+ *
+ * \param result Opaque structure that must be passed to \ref imdct().
+ *
+ * \return Standard.
*/
-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;
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;
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);
}