ac3enc: use generic fixed-point mdct

This makes the AC3 encoder use the shared fixed-point MDCT rather
than its own implementation.  The checksum changes are due to
different rounding in the MDCT.

Signed-off-by: Mans Rullgard <mans@mansr.com>

5 files changed, 21 insertions, 422 deletions

diff --git a/libavcodec/ac3enc.c b/libavcodec/ac3enc.c
index a869d7de9d..372f019f64 100644
--- a/libavcodec/ac3enc.c
+++ b/libavcodec/ac3enc.c
@@ -29,6 +29,8 @@
 //#define DEBUG
 //#define ASSERT_LEVEL 2
 
+#include <stdint.h>
+
 #include "libavutil/audioconvert.h"
 #include "libavutil/avassert.h"
 #include "libavutil/crc.h"
@@ -39,6 +41,7 @@
 #include "ac3dsp.h"
 #include "ac3.h"
 #include "audioconvert.h"
+#include "fft.h"
 
 
 #ifndef CONFIG_AC3ENC_FLOAT
@@ -55,16 +58,22 @@
 #define AC3_REMATRIXING_NONE    1
 #define AC3_REMATRIXING_ALWAYS  3
 
-/** Scale a float value by 2^bits and convert to an integer. */
-#define SCALE_FLOAT(a, bits) lrintf((a) * (float)(1 << (bits)))
-
-
 #if CONFIG_AC3ENC_FLOAT
-#include "ac3enc_float.h"
+#define MAC_COEF(d,a,b) ((d)+=(a)*(b))
+typedef float SampleType;
+typedef float CoefType;
+typedef float CoefSumType;
 #else
-#include "ac3enc_fixed.h"
+#define MAC_COEF(d,a,b) MAC64(d,a,b)
+typedef int16_t SampleType;
+typedef int32_t CoefType;
+typedef int64_t CoefSumType;
 #endif
 
+typedef struct AC3MDCTContext {
+    const SampleType *window;           ///< MDCT window function
+    FFTContext fft;                     ///< FFT context for MDCT calculation
+} AC3MDCTContext;
 
 /**
  * Encoding Options used by AVOption.
@@ -279,8 +288,6 @@ static av_cold void mdct_end(AC3MDCTContext *mdct);
 static av_cold int mdct_init(AVCodecContext *avctx, AC3MDCTContext *mdct,
                              int nbits);
 
-static void mdct512(AC3MDCTContext *mdct, CoefType *out, SampleType *in);
-
 static void apply_window(DSPContext *dsp, SampleType *output, const SampleType *input,
                          const SampleType *window, unsigned int len);
 
@@ -386,7 +393,8 @@ static void apply_mdct(AC3EncodeContext *s)
 
             block->coeff_shift[ch] = normalize_samples(s);
 
-            mdct512(&s->mdct, block->mdct_coef[ch], s->windowed_samples);
+            s->mdct.fft.mdct_calcw(&s->mdct.fft, block->mdct_coef[ch],
+                                   s->windowed_samples);
         }
     }
 }
diff --git a/libavcodec/ac3enc_fixed.c b/libavcodec/ac3enc_fixed.c
index 967c7f0805..720d87d5d7 100644
--- a/libavcodec/ac3enc_fixed.c
+++ b/libavcodec/ac3enc_fixed.c
@@ -26,54 +26,17 @@
  * fixed-point AC-3 encoder.
  */
 
+#define CONFIG_FFT_FLOAT 0
 #undef CONFIG_AC3ENC_FLOAT
 #include "ac3enc.c"
 
 
-/** Scale a float value by 2^15, convert to an integer, and clip to range -32767..32767. */
-#define FIX15(a) av_clip(SCALE_FLOAT(a, 15), -32767, 32767)
-
-
 /**
  * Finalize MDCT and free allocated memory.
  */
 static av_cold void mdct_end(AC3MDCTContext *mdct)
 {
-    mdct->nbits = 0;
-    av_freep(&mdct->costab);
-    av_freep(&mdct->sintab);
-    av_freep(&mdct->xcos1);
-    av_freep(&mdct->xsin1);
-    av_freep(&mdct->rot_tmp);
-    av_freep(&mdct->cplx_tmp);
-}
-
-
-/**
- * Initialize FFT tables.
- * @param ln log2(FFT size)
- */
-static av_cold int fft_init(AVCodecContext *avctx, AC3MDCTContext *mdct, int ln)
-{
-    int i, n, n2;
-    float alpha;
-
-    n  = 1 << ln;
-    n2 = n >> 1;
-
-    FF_ALLOC_OR_GOTO(avctx, mdct->costab, n2 * sizeof(*mdct->costab), fft_alloc_fail);
-    FF_ALLOC_OR_GOTO(avctx, mdct->sintab, n2 * sizeof(*mdct->sintab), fft_alloc_fail);
-
-    for (i = 0; i < n2; i++) {
-        alpha     = 2.0 * M_PI * i / n;
-        mdct->costab[i] = FIX15(cos(alpha));
-        mdct->sintab[i] = FIX15(sin(alpha));
-    }
-
-    return 0;
-fft_alloc_fail:
-    mdct_end(mdct);
-    return AVERROR(ENOMEM);
+    ff_fft_end(&mdct->fft);
 }
 
 
@@ -84,167 +47,9 @@ fft_alloc_fail:
 static av_cold int mdct_init(AVCodecContext *avctx, AC3MDCTContext *mdct,
                              int nbits)
 {
-    int i, n, n4, ret;
-
-    n  = 1 << nbits;
-    n4 = n >> 2;
-
-    mdct->nbits = nbits;
-
-    ret = fft_init(avctx, mdct, nbits - 2);
-    if (ret)
-        return ret;
-
+    int ret = ff_mdct_init(&mdct->fft, nbits, 0, 1.0);
     mdct->window = ff_ac3_window;
-
-    FF_ALLOC_OR_GOTO(avctx, mdct->xcos1,    n4 * sizeof(*mdct->xcos1),    mdct_alloc_fail);
-    FF_ALLOC_OR_GOTO(avctx, mdct->xsin1,    n4 * sizeof(*mdct->xsin1),    mdct_alloc_fail);
-    FF_ALLOC_OR_GOTO(avctx, mdct->rot_tmp,  n  * sizeof(*mdct->rot_tmp),  mdct_alloc_fail);
-    FF_ALLOC_OR_GOTO(avctx, mdct->cplx_tmp, n4 * sizeof(*mdct->cplx_tmp), mdct_alloc_fail);
-
-    for (i = 0; i < n4; i++) {
-        float alpha = 2.0 * M_PI * (i + 1.0 / 8.0) / n;
-        mdct->xcos1[i] = FIX15(-cos(alpha));
-        mdct->xsin1[i] = FIX15(-sin(alpha));
-    }
-
-    return 0;
-mdct_alloc_fail:
-    mdct_end(mdct);
-    return AVERROR(ENOMEM);
-}
-
-
-/** Butterfly op */
-#define BF(pre, pim, qre, qim, pre1, pim1, qre1, qim1)  \
-{                                                       \
-  int ax, ay, bx, by;                                   \
-  bx  = pre1;                                           \
-  by  = pim1;                                           \
-  ax  = qre1;                                           \
-  ay  = qim1;                                           \
-  pre = (bx + ax) >> 1;                                 \
-  pim = (by + ay) >> 1;                                 \
-  qre = (bx - ax) >> 1;                                 \
-  qim = (by - ay) >> 1;                                 \
-}
-
-
-/** Complex multiply */
-#define CMUL(pre, pim, are, aim, bre, bim, rshift)      \
-{                                                       \
-   pre = (MUL16(are, bre) - MUL16(aim, bim)) >> rshift; \
-   pim = (MUL16(are, bim) + MUL16(bre, aim)) >> rshift; \
-}
-
-
-/**
- * Calculate a 2^n point complex FFT on 2^ln points.
- * @param z  complex input/output samples
- * @param ln log2(FFT size)
- */
-static void fft(AC3MDCTContext *mdct, IComplex *z, int ln)
-{
-    int j, l, np, np2;
-    int nblocks, nloops;
-    register IComplex *p,*q;
-    int tmp_re, tmp_im;
-
-    np = 1 << ln;
-
-    /* reverse */
-    for (j = 0; j < np; j++) {
-        int k = av_reverse[j] >> (8 - ln);
-        if (k < j)
-            FFSWAP(IComplex, z[k], z[j]);
-    }
-
-    /* pass 0 */
-
-    p = &z[0];
-    j = np >> 1;
-    do {
-        BF(p[0].re, p[0].im, p[1].re, p[1].im,
-           p[0].re, p[0].im, p[1].re, p[1].im);
-        p += 2;
-    } while (--j);
-
-    /* pass 1 */
-
-    p = &z[0];
-    j = np >> 2;
-    do {
-        BF(p[0].re, p[0].im, p[2].re,  p[2].im,
-           p[0].re, p[0].im, p[2].re,  p[2].im);
-        BF(p[1].re, p[1].im, p[3].re,  p[3].im,
-           p[1].re, p[1].im, p[3].im, -p[3].re);
-        p+=4;
-    } while (--j);
-
-    /* pass 2 .. ln-1 */
-
-    nblocks = np >> 3;
-    nloops  =  1 << 2;
-    np2     = np >> 1;
-    do {
-        p = z;
-        q = z + nloops;
-        for (j = 0; j < nblocks; j++) {
-            BF(p->re, p->im, q->re, q->im,
-               p->re, p->im, q->re, q->im);
-            p++;
-            q++;
-            for(l = nblocks; l < np2; l += nblocks) {
-                CMUL(tmp_re, tmp_im, mdct->costab[l], -mdct->sintab[l], q->re, q->im, 15);
-                BF(p->re, p->im, q->re,  q->im,
-                   p->re, p->im, tmp_re, tmp_im);
-                p++;
-                q++;
-            }
-            p += nloops;
-            q += nloops;
-        }
-        nblocks = nblocks >> 1;
-        nloops  = nloops  << 1;
-    } while (nblocks);
-}
-
-
-/**
- * Calculate a 512-point MDCT
- * @param out 256 output frequency coefficients
- * @param in  512 windowed input audio samples
- */
-static void mdct512(AC3MDCTContext *mdct, int32_t *out, int16_t *in)
-{
-    int i, re, im, n, n2, n4;
-    int16_t *rot = mdct->rot_tmp;
-    IComplex *x  = mdct->cplx_tmp;
-
-    n  = 1 << mdct->nbits;
-    n2 = n >> 1;
-    n4 = n >> 2;
-
-    /* shift to simplify computations */
-    for (i = 0; i <n4; i++)
-        rot[i] = -in[i + 3*n4];
-    memcpy(&rot[n4], &in[0], 3*n4*sizeof(*in));
-
-    /* pre rotation */
-    for (i = 0; i < n4; i++) {
-        re =  ((int)rot[   2*i] - (int)rot[ n-1-2*i]) >> 1;
-        im = -((int)rot[n2+2*i] - (int)rot[n2-1-2*i]) >> 1;
-        CMUL(x[i].re, x[i].im, re, im, -mdct->xcos1[i], mdct->xsin1[i], 15);
-    }
-
-    fft(mdct, x, mdct->nbits - 2);
-
-    /* post rotation */
-    for (i = 0; i < n4; i++) {
-        re = x[i].re;
-        im = x[i].im;
-        CMUL(out[n2-1-2*i], out[2*i], re, im, mdct->xsin1[i], mdct->xcos1[i], 0);
-    }
+    return ret;
 }
 
 
@@ -304,101 +109,6 @@ static void scale_coefficients(AC3EncodeContext *s)
 }
 
 
-#ifdef TEST
-/*************************************************************************/
-/* TEST */
-
-#include "libavutil/lfg.h"
-
-#define MDCT_NBITS 9
-#define MDCT_SAMPLES (1 << MDCT_NBITS)
-#define FN (MDCT_SAMPLES/4)
-
-
-static void fft_test(AC3MDCTContext *mdct, AVLFG *lfg)
-{
-    IComplex in[FN], in1[FN];
-    int k, n, i;
-    float sum_re, sum_im, a;
-
-    for (i = 0; i < FN; i++) {
-        in[i].re = av_lfg_get(lfg) % 65535 - 32767;
-        in[i].im = av_lfg_get(lfg) % 65535 - 32767;
-        in1[i]   = in[i];
-    }
-    fft(mdct, in, 7);
-
-    /* do it by hand */
-    for (k = 0; k < FN; k++) {
-        sum_re = 0;
-        sum_im = 0;
-        for (n = 0; n < FN; n++) {
-            a = -2 * M_PI * (n * k) / FN;
-            sum_re += in1[n].re * cos(a) - in1[n].im * sin(a);
-            sum_im += in1[n].re * sin(a) + in1[n].im * cos(a);
-        }
-        av_log(NULL, AV_LOG_DEBUG, "%3d: %6d,%6d %6.0f,%6.0f\n",
-               k, in[k].re, in[k].im, sum_re / FN, sum_im / FN);
-    }
-}
-
-
-static void mdct_test(AC3MDCTContext *mdct, AVLFG *lfg)
-{
-    int16_t input[MDCT_SAMPLES];
-    int32_t output[AC3_MAX_COEFS];
-    float input1[MDCT_SAMPLES];
-    float output1[AC3_MAX_COEFS];
-    float s, a, err, e, emax;
-    int i, k, n;
-
-    for (i = 0; i < MDCT_SAMPLES; i++) {
-        input[i]  = (av_lfg_get(lfg) % 65535 - 32767) * 9 / 10;
-        input1[i] = input[i];
-    }
-
-    mdct512(mdct, output, input);
-
-    /* do it by hand */
-    for (k = 0; k < AC3_MAX_COEFS; k++) {
-        s = 0;
-        for (n = 0; n < MDCT_SAMPLES; n++) {
-            a = (2*M_PI*(2*n+1+MDCT_SAMPLES/2)*(2*k+1) / (4 * MDCT_SAMPLES));
-            s += input1[n] * cos(a);
-        }
-        output1[k] = -2 * s / MDCT_SAMPLES;
-    }
-
-    err  = 0;
-    emax = 0;
-    for (i = 0; i < AC3_MAX_COEFS; i++) {
-        av_log(NULL, AV_LOG_DEBUG, "%3d: %7d %7.0f\n", i, output[i], output1[i]);
-        e = output[i] - output1[i];
-        if (e > emax)
-            emax = e;
-        err += e * e;
-    }
-    av_log(NULL, AV_LOG_DEBUG, "err2=%f emax=%f\n", err / AC3_MAX_COEFS, emax);
-}
-
-
-int main(void)
-{
-    AVLFG lfg;
-    AC3MDCTContext mdct;
-
-    mdct.avctx = NULL;
-    av_log_set_level(AV_LOG_DEBUG);
-    mdct_init(&mdct, 9);
-
-    fft_test(&mdct, &lfg);
-    mdct_test(&mdct, &lfg);
-
-    return 0;
-}
-#endif /* TEST */
-
-
 AVCodec ff_ac3_fixed_encoder = {
     "ac3_fixed",
     AVMEDIA_TYPE_AUDIO,
diff --git a/libavcodec/ac3enc_fixed.h b/libavcodec/ac3enc_fixed.h
deleted file mode 100644
index 0c9e6a16c6..0000000000
--- a/libavcodec/ac3enc_fixed.h
+++ /dev/null
@@ -1,61 +0,0 @@
-/*
- * The simplest AC-3 encoder
- * Copyright (c) 2000 Fabrice Bellard
- * Copyright (c) 2006-2010 Justin Ruggles <justin.ruggles@gmail.com>
- * Copyright (c) 2006-2010 Prakash Punnoor <prakash@punnoor.de>
- *
- * This file is part of Libav.
- *
- * Libav is free software; you can redistribute it and/or
- * modify it under the terms of the GNU Lesser General Public
- * License as published by the Free Software Foundation; either
- * version 2.1 of the License, or (at your option) any later version.
- *
- * Libav 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
- * Lesser General Public License for more details.
- *
- * You should have received a copy of the GNU Lesser General Public
- * License along with Libav; if not, write to the Free Software
- * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
- */
-
-/**
- * @file
- * fixed-point AC-3 encoder header.
- */
-
-#ifndef AVCODEC_AC3ENC_FIXED_H
-#define AVCODEC_AC3ENC_FIXED_H
-
-#include <stdint.h>
-
-
-typedef int16_t SampleType;
-typedef int32_t CoefType;
-typedef int64_t CoefSumType;
-
-#define MAC_COEF(d,a,b) MAC64(d,a,b)
-
-
-/**
- * Compex number.
- * Used in fixed-point MDCT calculation.
- */
-typedef struct IComplex {
-    int16_t re,im;
-} IComplex;
-
-typedef struct AC3MDCTContext {
-    const int16_t *window;                  ///< MDCT window function
-    int nbits;                              ///< log2(transform size)
-    int16_t *costab;                        ///< FFT cos table
-    int16_t *sintab;                        ///< FFT sin table
-    int16_t *xcos1;                         ///< MDCT cos table
-    int16_t *xsin1;                         ///< MDCT sin table
-    int16_t *rot_tmp;                       ///< temp buffer for pre-rotated samples
-    IComplex *cplx_tmp;                     ///< temp buffer for complex pre-rotated samples
-} AC3MDCTContext;
-
-#endif /* AVCODEC_AC3ENC_FIXED_H */
diff --git a/libavcodec/ac3enc_float.c b/libavcodec/ac3enc_float.c
index 19afd27bb1..d2435dee15 100644
--- a/libavcodec/ac3enc_float.c
+++ b/libavcodec/ac3enc_float.c
@@ -69,17 +69,6 @@ static av_cold int mdct_init(AVCodecContext *avctx, AC3MDCTContext *mdct,
 
 
 /**
- * Calculate a 512-point MDCT
- * @param out 256 output frequency coefficients
- * @param in  512 windowed input audio samples
- */
-static void mdct512(AC3MDCTContext *mdct, float *out, float *in)
-{
-    mdct->fft.mdct_calc(&mdct->fft, out, in);
-}
-
-
-/**
  * Apply KBD window to input samples prior to MDCT.
  */
 static void apply_window(DSPContext *dsp, float *output, const float *input,
diff --git a/libavcodec/ac3enc_float.h b/libavcodec/ac3enc_float.h
deleted file mode 100644
index ecdc421e92..0000000000
--- a/libavcodec/ac3enc_float.h
+++ /dev/null
@@ -1,47 +0,0 @@
-/*
- * The simplest AC-3 encoder
- * Copyright (c) 2000 Fabrice Bellard
- * Copyright (c) 2006-2010 Justin Ruggles <justin.ruggles@gmail.com>
- * Copyright (c) 2006-2010 Prakash Punnoor <prakash@punnoor.de>
- *
- * This file is part of Libav.
- *
- * Libav is free software; you can redistribute it and/or
- * modify it under the terms of the GNU Lesser General Public
- * License as published by the Free Software Foundation; either
- * version 2.1 of the License, or (at your option) any later version.
- *
- * Libav 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
- * Lesser General Public License for more details.
- *
- * You should have received a copy of the GNU Lesser General Public
- * License along with Libav; if not, write to the Free Software
- * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
- */
-
-/**
- * @file
- * floating-point AC-3 encoder header.
- */
-
-#ifndef AVCODEC_AC3ENC_FLOAT_H
-#define AVCODEC_AC3ENC_FLOAT_H
-
-#include "fft.h"
-
-
-typedef float SampleType;
-typedef float CoefType;
-typedef float CoefSumType;
-
-#define MAC_COEF(d,a,b) ((d)+=(a)*(b))
-
-
-typedef struct AC3MDCTContext {
-    const float *window;    ///< MDCT window function
-    FFTContext fft;         ///< FFT context for MDCT calculation
-} AC3MDCTContext;
-
-#endif /* AVCODEC_AC3ENC_FLOAT_H */