WMAVoice postfilter.

Originally committed as revision 22938 to svn://svn.ffmpeg.org/ffmpeg/trunk

2 files changed, 706 insertions, 7 deletions

diff --git a/libavcodec/wmavoice.c b/libavcodec/wmavoice.c
index 97dabd2526..04bb6b3420 100644
--- a/libavcodec/wmavoice.c
+++ b/libavcodec/wmavoice.c
@@ -36,9 +36,13 @@
 #include "acelp_filters.h"
 #include "lsp.h"
 #include "libavutil/lzo.h"
+#include "avfft.h"
+#include "fft.h"
 
 #define MAX_BLOCKS           8   ///< maximum number of blocks per frame
 #define MAX_LSPS             16  ///< maximum filter order
+#define MAX_LSPS_ALIGN16     16  ///< same as #MAX_LSPS; needs to be multiple
+                                 ///< of 16 for ASM input buffer alignment
 #define MAX_FRAMES           3   ///< maximum number of frames per superframe
 #define MAX_FRAMESIZE        160 ///< maximum number of samples per frame
 #define MAX_SIGNAL_HISTORY   416 ///< maximum excitation signal history
@@ -140,8 +144,15 @@ typedef struct {
     int history_nsamples;         ///< number of samples in history for signal
                                   ///< prediction (through ACB)
 
+    /* postfilter specific values */
     int do_apf;                   ///< whether to apply the averaged
                                   ///< projection filter (APF)
+    int denoise_strength;         ///< strength of denoising in Wiener filter
+                                  ///< [0-11]
+    int denoise_tilt_corr;        ///< Whether to apply tilt correction to the
+                                  ///< Wiener filter coefficients (postfilter)
+    int dc_level;                 ///< Predicted amount of DC noise, based
+                                  ///< on which a DC removal filter is used
 
     int lsps;                     ///< number of LSPs per frame [10 or 16]
     int lsp_q_mode;               ///< defines quantizer defaults [0, 1]
@@ -244,6 +255,34 @@ typedef struct {
     float synth_history[MAX_LSPS]; ///< see #excitation_history
     /**
      * @}
+     * @defgroup post_filter Postfilter values
+     * Varibales used for postfilter implementation, mostly history for
+     * smoothing and so on, and context variables for FFT/iFFT.
+     * @{
+     */
+    RDFTContext rdft, irdft;      ///< contexts for FFT-calculation in the
+                                  ///< postfilter (for denoise filter)
+    DCTContext dct, dst;          ///< contexts for phase shift (in Hilbert
+                                  ///< transform, part of postfilter)
+    float sin[511], cos[511];     ///< 8-bit cosine/sine windows over [-pi,pi]
+                                  ///< range
+    float postfilter_agc;         ///< gain control memory, used in
+                                  ///< #adaptive_gain_control()
+    float dcf_mem[2];             ///< DC filter history
+    float zero_exc_pf[MAX_SIGNAL_HISTORY + MAX_SFRAMESIZE];
+                                  ///< zero filter output (i.e. excitation)
+                                  ///< by postfilter
+    float denoise_filter_cache[MAX_FRAMESIZE];
+    int   denoise_filter_cache_size; ///< samples in #denoise_filter_cache
+    DECLARE_ALIGNED(16, float, tilted_lpcs_pf)[0x80];
+                                  ///< aligned buffer for LPC tilting
+    DECLARE_ALIGNED(16, float, denoise_coeffs_pf)[0x80];
+                                  ///< aligned buffer for denoise coefficients
+    DECLARE_ALIGNED(16, float, synth_filter_out_buf)[80 + MAX_LSPS_ALIGN16];
+                                  ///< aligned buffer for postfilter speech
+                                  ///< synthesis
+    /**
+     * @}
      */
 } WMAVoiceContext;
 
@@ -313,6 +352,28 @@ static av_cold int wmavoice_decode_init(AVCodecContext *ctx)
     flags                = AV_RL32(ctx->extradata + 18);
     s->spillover_bitsize = 3 + av_ceil_log2(ctx->block_align);
     s->do_apf            =    flags & 0x1;
+    if (s->do_apf) {
+        ff_rdft_init(&s->rdft,  7, DFT_R2C);
+        ff_rdft_init(&s->irdft, 7, IDFT_C2R);
+        ff_dct_init(&s->dct,  6, DCT_I);
+        ff_dct_init(&s->dst,  6, DST_I);
+
+        ff_sine_window_init(s->cos, 256);
+        memcpy(&s->sin[255], s->cos, 256 * sizeof(s->cos[0]));
+        for (n = 0; n < 255; n++) {
+            s->sin[n]       = -s->sin[510 - n];
+            s->cos[510 - n] =  s->cos[n];
+        }
+    }
+    s->denoise_strength  =   (flags >> 2) & 0xF;
+    if (s->denoise_strength >= 12) {
+        av_log(ctx, AV_LOG_ERROR,
+               "Invalid denoise filter strength %d (max=11)\n",
+               s->denoise_strength);
+        return -1;
+    }
+    s->denoise_tilt_corr = !!(flags & 0x40);
+    s->dc_level          =   (flags >> 7) & 0xF;
     s->lsp_q_mode        = !!(flags & 0x2000);
     s->lsp_def_mode      = !!(flags & 0x4000);
     lsp16_flag           =    flags & 0x1000;
@@ -370,6 +431,366 @@ static av_cold int wmavoice_decode_init(AVCodecContext *ctx)
 }
 
 /**
+ * @defgroup postfilter Postfilter functions
+ * Postfilter functions (gain control, wiener denoise filter, DC filter,
+ * kalman smoothening, plus surrounding code to wrap it)
+ * @{
+ */
+/**
+ * Adaptive gain control (as used in postfilter).
+ *
+ * Identical to #ff_adaptive_gain_control() in acelp_vectors.c, except
+ * that the energy here is calculated using sum(abs(...)), whereas the
+ * other codecs (e.g. AMR-NB, SIPRO) use sqrt(dotproduct(...)).
+ *
+ * @param out output buffer for filtered samples
+ * @param in input buffer containing the samples as they are after the
+ *           postfilter steps so far
+ * @param speech_synth input buffer containing speech synth before postfilter
+ * @param size input buffer size
+ * @param alpha exponential filter factor
+ * @param gain_mem pointer to filter memory (single float)
+ */
+static void adaptive_gain_control(float *out, const float *in,
+                                  const float *speech_synth,
+                                  int size, float alpha, float *gain_mem)
+{
+    int i;
+    float speech_energy = 0.0, postfilter_energy = 0.0, gain_scale_factor;
+    float mem = *gain_mem;
+
+    for (i = 0; i < size; i++) {
+        speech_energy     += fabsf(speech_synth[i]);
+        postfilter_energy += fabsf(in[i]);
+    }
+    gain_scale_factor = (1.0 - alpha) * speech_energy / postfilter_energy;
+
+    for (i = 0; i < size; i++) {
+        mem = alpha * mem + gain_scale_factor;
+        out[i] = in[i] * mem;
+    }
+
+    *gain_mem = mem;
+}
+
+/**
+ * Kalman smoothing function.
+ *
+ * This function looks back pitch +/- 3 samples back into history to find
+ * the best fitting curve (that one giving the optimal gain of the two
+ * signals, i.e. the highest dot product between the two), and then
+ * uses that signal history to smoothen the output of the speech synthesis
+ * filter.
+ *
+ * @param s WMA Voice decoding context
+ * @param pitch pitch of the speech signal
+ * @param in input speech signal
+ * @param out output pointer for smoothened signal
+ * @param size input/output buffer size
+ *
+ * @returns -1 if no smoothening took place, e.g. because no optimal
+ *          fit could be found, or 0 on success.
+ */
+static int kalman_smoothen(WMAVoiceContext *s, int pitch,
+                           const float *in, float *out, int size)
+{
+    int n;
+    float optimal_gain = 0, dot;
+    const float *ptr = &in[-FFMAX(s->min_pitch_val, pitch - 3)],
+                *end = &in[-FFMIN(s->max_pitch_val, pitch + 3)],
+                *best_hist_ptr;
+
+    /* find best fitting point in history */
+    do {
+        dot = ff_dot_productf(in, ptr, size);
+        if (dot > optimal_gain) {
+            optimal_gain  = dot;
+            best_hist_ptr = ptr;
+        }
+    } while (--ptr >= end);
+
+    if (optimal_gain <= 0)
+        return -1;
+    dot = ff_dot_productf(best_hist_ptr, best_hist_ptr, size);
+    if (dot <= 0) // would be 1.0
+        return -1;
+
+    if (optimal_gain <= dot) {
+        dot = dot / (dot + 0.6 * optimal_gain); // 0.625-1.000
+    } else
+        dot = 0.625;
+
+    /* actual smoothing */
+    for (n = 0; n < size; n++)
+        out[n] = best_hist_ptr[n] + dot * (in[n] - best_hist_ptr[n]);
+
+    return 0;
+}
+
+/**
+ * Get the tilt factor of a formant filter from its transfer function
+ * @see #tilt_factor() in amrnbdec.c, which does essentially the same,
+ *      but somehow (??) it does a speech synthesis filter in the
+ *      middle, which is missing here
+ *
+ * @param lpcs LPC coefficients
+ * @param n_lpcs Size of LPC buffer
+ * @returns the tilt factor
+ */
+static float tilt_factor(const float *lpcs, int n_lpcs)
+{
+    float rh0, rh1;
+
+    rh0 = 1.0     + ff_dot_productf(lpcs,  lpcs,    n_lpcs);
+    rh1 = lpcs[0] + ff_dot_productf(lpcs, &lpcs[1], n_lpcs - 1);
+
+    return rh1 / rh0;
+}
+
+/**
+ * Derive denoise filter coefficients (in real domain) from the LPCs.
+ */
+static void calc_input_response(WMAVoiceContext *s, float *lpcs,
+                                int fcb_type, float *coeffs, int remainder)
+{
+    float last_coeff, min = 15.0, max = -15.0;
+    float irange, angle_mul, gain_mul, range, sq;
+    int n, idx;
+
+    /* Create frequency power spectrum of speech input (i.e. RDFT of LPCs) */
+    ff_rdft_calc(&s->rdft, lpcs);
+#define log_range(var, assign) do { \
+        float tmp = log10f(assign);  var = tmp; \
+        max       = FFMAX(max, tmp); min = FFMIN(min, tmp); \
+    } while (0)
+    log_range(last_coeff,  lpcs[1]         * lpcs[1]);
+    for (n = 1; n < 64; n++)
+        log_range(lpcs[n], lpcs[n * 2]     * lpcs[n * 2] +
+                           lpcs[n * 2 + 1] * lpcs[n * 2 + 1]);
+    log_range(lpcs[0],     lpcs[0]         * lpcs[0]);
+#undef log_range
+    range    = max - min;
+    lpcs[64] = last_coeff;
+
+    /* Now, use this spectrum to pick out these frequencies with higher
+     * (relative) power/energy (which we then take to be "not noise"),
+     * and set up a table (still in lpc[]) of (relative) gains per frequency.
+     * These frequencies will be maintained, while others ("noise") will be
+     * decreased in the filter output. */
+    irange    = 64.0 / range; // so irange*(max-value) is in the range [0, 63]
+    gain_mul  = range * (fcb_type == FCB_TYPE_HARDCODED ? (5.0 / 13.0) :
+                                                          (5.0 / 14.7));
+    angle_mul = gain_mul * (8.0 * M_LN10 / M_PI);
+    for (n = 0; n <= 64; n++) {
+        float pow;
+
+        idx = FFMAX(0, lrint((max - lpcs[n]) * irange) - 1);
+        pow = wmavoice_denoise_power_table[s->denoise_strength][idx];
+        lpcs[n] = angle_mul * pow;
+
+        /* 70.57 =~ 1/log10(1.0331663) */
+        idx = (pow * gain_mul - 0.0295) * 70.570526123;
+        if (idx > 127) { // fallback if index falls outside table range
+            coeffs[n] = wmavoice_energy_table[127] *
+                        powf(1.0331663, idx - 127);
+        } else
+            coeffs[n] = wmavoice_energy_table[FFMAX(0, idx)];
+    }
+
+    /* calculate the Hilbert transform of the gains, which we do (since this
+     * is a sinus input) by doing a phase shift (in theory, H(sin())=cos()).
+     * Hilbert_Transform(RDFT(x)) = Laplace_Transform(x), which calculates the
+     * "moment" of the LPCs in this filter. */
+    ff_dct_calc(&s->dct, lpcs);
+    ff_dct_calc(&s->dst, lpcs);
+
+    /* Split out the coefficient indexes into phase/magnitude pairs */
+    idx = 255 + av_clip(lpcs[64],               -255, 255);
+    coeffs[0]  = coeffs[0]  * s->cos[idx];
+    idx = 255 + av_clip(lpcs[64] - 2 * lpcs[63], -255, 255);
+    last_coeff = coeffs[64] * s->cos[idx];
+    for (n = 63;; n--) {
+        idx = 255 + av_clip(-lpcs[64] - 2 * lpcs[n - 1], -255, 255);
+        coeffs[n * 2 + 1] = coeffs[n] * s->sin[idx];
+        coeffs[n * 2]     = coeffs[n] * s->cos[idx];
+
+        if (!--n) break;
+
+        idx = 255 + av_clip( lpcs[64] - 2 * lpcs[n - 1], -255, 255);
+        coeffs[n * 2 + 1] = coeffs[n] * s->sin[idx];
+        coeffs[n * 2]     = coeffs[n] * s->cos[idx];
+    }
+    coeffs[1] = last_coeff;
+
+    /* move into real domain */
+    ff_rdft_calc(&s->irdft, coeffs);
+
+    /* tilt correction and normalize scale */
+    memset(&coeffs[remainder], 0, sizeof(coeffs[0]) * (128 - remainder));
+    if (s->denoise_tilt_corr) {
+        float tilt_mem = 0;
+
+        coeffs[remainder - 1] = 0;
+        ff_tilt_compensation(&tilt_mem,
+                             -1.8 * tilt_factor(coeffs, remainder - 1),
+                             coeffs, remainder);
+    }
+    sq = (1.0 / 64.0) * sqrtf(1 / ff_dot_productf(coeffs, coeffs, remainder));
+    for (n = 0; n < remainder; n++)
+        coeffs[n] *= sq;
+}
+
+/**
+ * This function applies a Wiener filter on the (noisy) speech signal as
+ * a means to denoise it.
+ *
+ * - take RDFT of LPCs to get the power spectrum of the noise + speech;
+ * - using this power spectrum, calculate (for each frequency) the Wiener
+ *    filter gain, which depends on the frequency power and desired level
+ *    of noise subtraction (when set too high, this leads to artifacts)
+ *    We can do this symmetrically over the X-axis (so 0-4kHz is the inverse
+ *    of 4-8kHz);
+ * - by doing a phase shift, calculate the Hilbert transform of this array
+ *    of per-frequency filter-gains to get the filtering coefficients;
+ * - smoothen/normalize/de-tilt these filter coefficients as desired;
+ * - take RDFT of noisy sound, apply the coefficients and take its IRDFT
+ *    to get the denoised speech signal;
+ * - the leftover (i.e. output of the IRDFT on denoised speech data beyond
+ *    the frame boundary) are saved and applied to subsequent frames by an
+ *    overlap-add method (otherwise you get clicking-artifacts).
+ *
+ * @param s WMA Voice decoding context
+ * @param s fcb_type Frame (codebook) type
+ * @param synth_pf input: the noisy speech signal, output: denoised speech
+ *                 data; should be 16-byte aligned (for ASM purposes)
+ * @param size size of the speech data
+ * @param lpcs LPCs used to synthesize this frame's speech data
+ */
+static void wiener_denoise(WMAVoiceContext *s, int fcb_type,
+                           float *synth_pf, int size,
+                           const float *lpcs)
+{
+    int remainder, lim, n;
+
+    if (fcb_type != FCB_TYPE_SILENCE) {
+        float *tilted_lpcs = s->tilted_lpcs_pf,
+              *coeffs = s->denoise_coeffs_pf, tilt_mem = 0;
+
+        tilted_lpcs[0]           = 1.0;
+        memcpy(&tilted_lpcs[1], lpcs, sizeof(lpcs[0]) * s->lsps);
+        memset(&tilted_lpcs[s->lsps + 1], 0,
+               sizeof(tilted_lpcs[0]) * (128 - s->lsps - 1));
+        ff_tilt_compensation(&tilt_mem, 0.7 * tilt_factor(lpcs, s->lsps),
+                             tilted_lpcs, s->lsps + 2);
+
+        /* The IRDFT output (127 samples for 7-bit filter) beyond the frame
+         * size is applied to the next frame. All input beyond this is zero,
+         * and thus all output beyond this will go towards zero, hence we can
+         * limit to min(size-1, 127-size) as a performance consideration. */
+        remainder = FFMIN(127 - size, size - 1);
+        calc_input_response(s, tilted_lpcs, fcb_type, coeffs, remainder);
+
+        /* apply coefficients (in frequency spectrum domain), i.e. complex
+         * number multiplication */
+        memset(&synth_pf[size], 0, sizeof(synth_pf[0]) * (128 - size));
+        ff_rdft_calc(&s->rdft, synth_pf);
+        ff_rdft_calc(&s->rdft, coeffs);
+        synth_pf[0] *= coeffs[0];
+        synth_pf[1] *= coeffs[1];
+        for (n = 1; n < 128; n++) {
+            float v1 = synth_pf[n * 2], v2 = synth_pf[n * 2 + 1];
+            synth_pf[n * 2]     = v1 * coeffs[n * 2] - v2 * coeffs[n * 2 + 1];
+            synth_pf[n * 2 + 1] = v2 * coeffs[n * 2] + v1 * coeffs[n * 2 + 1];
+        }
+        ff_rdft_calc(&s->irdft, synth_pf);
+    }
+
+    /* merge filter output with the history of previous runs */
+    if (s->denoise_filter_cache_size) {
+        lim = FFMIN(s->denoise_filter_cache_size, size);
+        for (n = 0; n < lim; n++)
+            synth_pf[n] += s->denoise_filter_cache[n];
+        s->denoise_filter_cache_size -= lim;
+        memmove(s->denoise_filter_cache, &s->denoise_filter_cache[size],
+                sizeof(s->denoise_filter_cache[0]) * s->denoise_filter_cache_size);
+    }
+
+    /* move remainder of filter output into a cache for future runs */
+    if (fcb_type != FCB_TYPE_SILENCE) {
+        lim = FFMIN(remainder, s->denoise_filter_cache_size);
+        for (n = 0; n < lim; n++)
+            s->denoise_filter_cache[n] += synth_pf[size + n];
+        if (lim < remainder) {
+            memcpy(&s->denoise_filter_cache[lim], &synth_pf[size + lim],
+                   sizeof(s->denoise_filter_cache[0]) * (remainder - lim));
+            s->denoise_filter_cache_size = remainder;
+        }
+    }
+}
+
+/**
+ * Averaging projection filter, the postfilter used in WMAVoice.
+ *
+ * This uses the following steps:
+ * - A zero-synthesis filter (generate excitation from synth signal)
+ * - Kalman smoothing on excitation, based on pitch
+ * - Re-synthesized smoothened output
+ * - Iterative Wiener denoise filter
+ * - Adaptive gain filter
+ * - DC filter
+ *
+ * @param s WMAVoice decoding context
+ * @param synth Speech synthesis output (before postfilter)
+ * @param samples Output buffer for filtered samples
+ * @param size Buffer size of synth & samples
+ * @param lpcs Generated LPCs used for speech synthesis
+ * @param fcb_type Frame type (silence, hardcoded, AW-pulses or FCB-pulses)
+ * @param pitch Pitch of the input signal
+ */
+static void postfilter(WMAVoiceContext *s, const float *synth,
+                       float *samples,    int size,
+                       const float *lpcs, float *zero_exc_pf,
+                       int fcb_type,      int pitch)
+{
+    float synth_filter_in_buf[MAX_FRAMESIZE / 2],
+          *synth_pf = &s->synth_filter_out_buf[MAX_LSPS_ALIGN16],
+          *synth_filter_in = zero_exc_pf;
+
+    assert(size <= MAX_FRAMESIZE / 2);
+
+    /* generate excitation from input signal */
+    ff_celp_lp_zero_synthesis_filterf(zero_exc_pf, lpcs, synth, size, s->lsps);
+
+    if (fcb_type >= FCB_TYPE_AW_PULSES &&
+        !kalman_smoothen(s, pitch, zero_exc_pf, synth_filter_in_buf, size))
+        synth_filter_in = synth_filter_in_buf;
+
+    /* re-synthesize speech after smoothening, and keep history */
+    ff_celp_lp_synthesis_filterf(synth_pf, lpcs,
+                                 synth_filter_in, size, s->lsps);
+    memcpy(&synth_pf[-s->lsps], &synth_pf[size - s->lsps],
+           sizeof(synth_pf[0]) * s->lsps);
+
+    wiener_denoise(s, fcb_type, synth_pf, size, lpcs);
+
+    adaptive_gain_control(samples, synth_pf, synth, size, 0.99,
+                          &s->postfilter_agc);
+
+    if (s->dc_level > 8) {
+        /* remove ultra-low frequency DC noise / highpass filter;
+         * coefficients are identical to those used in SIPR decoding,
+         * and very closely resemble those used in AMR-NB decoding. */
+        ff_acelp_apply_order_2_transfer_function(samples, samples,
+            (const float[2]) { -1.99997,      1.0 },
+            (const float[2]) { -1.9330735188, 0.93589198496 },
+            0.93980580475, s->dcf_mem, size);
+    }
+}
+/**
+ * @}
+ */
+
+/**
  * Dequantize LSPs
  * @param lsps output pointer to the array that will hold the LSPs
  * @param num number of LSPs to be dequantized
@@ -980,6 +1401,7 @@ static void synth_block(WMAVoiceContext *s, GetBitContext *gb,
  *
  * @param ctx WMA Voice decoder context
  * @param gb bit I/O context (s->gb or one for cross-packet superframes)
+ * @param frame_idx Frame number within superframe [0-2]
  * @param samples pointer to output sample buffer, has space for at least 160
  *                samples
  * @param lsps LSP array
@@ -988,7 +1410,7 @@ static void synth_block(WMAVoiceContext *s, GetBitContext *gb,
  * @param synth target buffer for synthesized speech data
  * @return 0 on success, <0 on error.
  */
-static int synth_frame(AVCodecContext *ctx, GetBitContext *gb,
+static int synth_frame(AVCodecContext *ctx, GetBitContext *gb, int frame_idx,
                        float *samples,
                        const double *lsps, const double *prev_lsps,
                        float *excitation, float *synth)
@@ -1113,10 +1535,23 @@ static int synth_frame(AVCodecContext *ctx, GetBitContext *gb,
     /* Averaging projection filter, if applicable. Else, just copy samples
      * from synthesis buffer */
     if (s->do_apf) {
-        // FIXME this is where APF would take place, currently not implemented
-        av_log_missing_feature(ctx, "APF", 0);
-        s->do_apf = 0;
-    } //else
+        double i_lsps[MAX_LSPS];
+        float lpcs[MAX_LSPS];
+
+        for (n = 0; n < s->lsps; n++) // LSF -> LSP
+            i_lsps[n] = cos(0.5 * (prev_lsps[n] + lsps[n]));
+        ff_acelp_lspd2lpc(i_lsps, lpcs, s->lsps >> 1);
+        postfilter(s, synth, samples, 80, lpcs,
+                   &s->zero_exc_pf[s->history_nsamples + MAX_FRAMESIZE * frame_idx],
+                   frame_descs[bd_idx].fcb_type, pitch[0]);
+
+        for (n = 0; n < s->lsps; n++) // LSF -> LSP
+            i_lsps[n] = cos(lsps[n]);
+        ff_acelp_lspd2lpc(i_lsps, lpcs, s->lsps >> 1);
+        postfilter(s, &synth[80], &samples[80], 80, lpcs,
+                   &s->zero_exc_pf[s->history_nsamples + MAX_FRAMESIZE * frame_idx + 80],
+                   frame_descs[bd_idx].fcb_type, pitch[0]);
+    } else
         memcpy(samples, synth, 160 * sizeof(synth[0]));
 
     /* Cache values for next frame */
@@ -1355,7 +1790,7 @@ static int synth_superframe(AVCodecContext *ctx,
             stabilize_lsps(lsps[n], s->lsps);
         }
 
-        if ((res = synth_frame(ctx, gb,
+        if ((res = synth_frame(ctx, gb, n,
                                &samples[n * MAX_FRAMESIZE],
                                lsps[n], n == 0 ? s->prev_lsps : lsps[n - 1],
                                &excitation[s->history_nsamples + n * MAX_FRAMESIZE],
@@ -1381,6 +1816,9 @@ static int synth_superframe(AVCodecContext *ctx,
            s->lsps             * sizeof(*synth));
     memcpy(s->excitation_history, &excitation[MAX_SFRAMESIZE],
            s->history_nsamples * sizeof(*excitation));
+    if (s->do_apf)
+        memmove(s->zero_exc_pf,       &s->zero_exc_pf[MAX_SFRAMESIZE],
+                s->history_nsamples * sizeof(*s->zero_exc_pf));
 
     return 0;
 }
@@ -1535,11 +1973,26 @@ static int wmavoice_decode_packet(AVCodecContext *ctx, void *data,
     return size;
 }
 
+static av_cold int wmavoice_decode_end(AVCodecContext *ctx)
+{
+    WMAVoiceContext *s = ctx->priv_data;
+
+    if (s->do_apf) {
+        ff_rdft_end(&s->rdft);
+        ff_rdft_end(&s->irdft);
+        ff_dct_end(&s->dct);
+        ff_dct_end(&s->dst);
+    }
+
+    return 0;
+}
+
 static av_cold void wmavoice_flush(AVCodecContext *ctx)
 {
     WMAVoiceContext *s = ctx->priv_data;
     int n;
 
+    s->postfilter_agc    = 0;
     s->sframe_cache_size = 0;
     s->skip_bits_next    = 0;
     for (n = 0; n < s->lsps; n++)
@@ -1550,6 +2003,16 @@ static av_cold void wmavoice_flush(AVCodecContext *ctx)
            sizeof(*s->synth_history)      * MAX_LSPS);
     memset(s->gain_pred_err,      0,
            sizeof(s->gain_pred_err));
+
+    if (s->do_apf) {
+        memset(&s->synth_filter_out_buf[MAX_LSPS_ALIGN16 - s->lsps], 0,
+               sizeof(*s->synth_filter_out_buf) * s->lsps);
+        memset(s->dcf_mem,              0,
+               sizeof(*s->dcf_mem)              * 2);
+        memset(s->zero_exc_pf,          0,
+               sizeof(*s->zero_exc_pf)          * s->history_nsamples);
+        memset(s->denoise_filter_cache, 0, sizeof(s->denoise_filter_cache));
+    }
 }
 
 AVCodec wmavoice_decoder = {
@@ -1559,7 +2022,7 @@ AVCodec wmavoice_decoder = {
     sizeof(WMAVoiceContext),
     wmavoice_decode_init,
     NULL,
-    NULL,
+    wmavoice_decode_end,
     wmavoice_decode_packet,
     CODEC_CAP_SUBFRAMES,
     .flush     = wmavoice_flush,
diff --git a/libavcodec/wmavoice_data.h b/libavcodec/wmavoice_data.h
index 5d74f384e5..cbf65b043e 100644
--- a/libavcodec/wmavoice_data.h
+++ b/libavcodec/wmavoice_data.h
@@ -3020,4 +3020,240 @@ static const float wmavoice_ipol2_coeffs[32] = {
     0, -0.0273968070, -0.0392575669, -0.0276240534
 };
 
+/**
+ * LUT for 1.071575641632 * pow(1.0331663, n - 127)
+ */
+static const float wmavoice_energy_table[128] = {
+    0.0169982178, 0.0175619858, 0.0181444519, 0.0187462362,
+    0.0193679795, 0.0200103437, 0.0206740128, 0.0213596933,
+    0.0220681153, 0.0228000330, 0.0235562258, 0.0243374986,
+    0.0251446834, 0.0259786395, 0.0268402549, 0.0277304468,
+    0.0286501631, 0.0296003830, 0.0305821182, 0.0315964139,
+    0.0326443501, 0.0337270424, 0.0348456436, 0.0360013446,
+    0.0371953760, 0.0384290090, 0.0397035571, 0.0410203772,
+    0.0423808713, 0.0437864880, 0.0452387238, 0.0467391249,
+    0.0482892887, 0.0498908657, 0.0515455612, 0.0532551367,
+    0.0550214125, 0.0568462692, 0.0587316496, 0.0606795611,
+    0.0626920777, 0.0647713419, 0.0669195677, 0.0691390421,
+    0.0714321284, 0.0738012678, 0.0762489827, 0.0787778794,
+    0.0813906502, 0.0840900769, 0.0868790336, 0.0897604897,
+    0.0927375130, 0.0958132732, 0.0989910450, 0.1022742117,
+    0.1056662688, 0.1091708280, 0.1127916204, 0.1165325012,
+    0.1203974531, 0.1243905911, 0.1285161668, 0.1327785725,
+    0.1371823465, 0.1417321773, 0.1464329093, 0.1512895470,
+    0.1563072616, 0.1614913951, 0.1668474671, 0.1723811803,
+    0.1780984262, 0.1840052921, 0.1901080668, 0.1964132480,
+    0.2029275487, 0.2096579046, 0.2166114816, 0.2237956830,
+    0.2312181577, 0.2388868085, 0.2468098001, 0.2549955679,
+    0.2634528274, 0.2721905830, 0.2812181375, 0.2905451026,
+    0.3001814086, 0.3101373153, 0.3204234225, 0.3310506819,
+    0.3420304081, 0.3533742912, 0.3650944090, 0.3772032397,
+    0.3897136755, 0.4026390362, 0.4159930832, 0.4297900346,
+    0.4440445799, 0.4587718956, 0.4739876619, 0.4897080789,
+    0.5059498840, 0.5227303696, 0.5400674019, 0.5579794393,
+    0.5764855528, 0.5956054456, 0.6153594745, 0.6357686714,
+    0.6568547659, 0.6786402082, 0.7011481929, 0.7244026842,
+    0.7484284410, 0.7732510432, 0.7988969192, 0.8253933741,
+    0.8527686184, 0.8810517982, 0.9102730265, 0.9404634147,
+    0.9716551065, 1.0038813113, 1.0371763400, 1.0715756416
+};
+
+/**
+ * LUT for f(x,y) = pow((y + 6.9) / 64, 0.025 * (x + 1)).
+ */
+static const float wmavoice_denoise_power_table[12][64] = {
+    { 0.9458379339, 0.9490436287, 0.9518757236, 0.9544130754,
+      0.9567118717, 0.9588135761, 0.9607496688, 0.9625446194,
+      0.9642178285, 0.9657849396, 0.9672587526, 0.9686498743,
+      0.9699671937, 0.9712182343, 0.9724094211, 0.9735462842,
+      0.9746336187, 0.9756756090, 0.9766759291, 0.9776378218,
+      0.9785641645, 0.9794575217, 0.9803201890, 0.9811542296,
+      0.9819615045, 0.9827436985, 0.9835023412, 0.9842388263,
+      0.9849544265, 0.9856503078, 0.9863275406, 0.9869871101,
+      0.9876299254, 0.9882568267, 0.9888685922, 0.9894659445,
+      0.9900495551, 0.9906200497, 0.9911780119, 0.9917239872,
+      0.9922584859, 0.9927819864, 0.9932949377, 0.9937977618,
+      0.9942908555, 0.9947745929, 0.9952493267, 0.9957153901,
+      0.9961730980, 0.9966227482, 0.9970646231, 0.9974989903,
+      0.9979261037, 0.9983462046, 0.9987595223, 0.9991662752,
+      0.9995666709, 0.9999609077, 1.0003491745, 1.0007316515,
+      1.0011085110, 1.0014799178, 1.0018460292, 1.0022069960 },
+    { 0.8946093973, 0.9006838092, 0.9060673931, 0.9109043185,
+      0.9152976055, 0.9193234737, 0.9230399260, 0.9264921443,
+      0.9297160207, 0.9327405496, 0.9355894944, 0.9382825789,
+      0.9408363568, 0.9432648587, 0.9455800822, 0.9477923675,
+      0.9499106907, 0.9519428941, 0.9538958704, 0.9557757107,
+      0.9575878241, 0.9593370368, 0.9610276730, 0.9626636222,
+      0.9642483964, 0.9657851769, 0.9672768552, 0.9687260672,
+      0.9701352224, 0.9715065293, 0.9728420173, 0.9741435556,
+      0.9754128696, 0.9766515555, 0.9778610927, 0.9790428553,
+      0.9801981216, 0.9813280829, 0.9824338513, 0.9835164667,
+      0.9845769028, 0.9856160726, 0.9866348334, 0.9876339913,
+      0.9886143053, 0.9895764906, 0.9905212223, 0.9914491381,
+      0.9923608411, 0.9932569022, 0.9941378627, 0.9950042356,
+      0.9958565084, 0.9966951442, 0.9975205834, 0.9983332454,
+      0.9991335296, 0.9999218170, 1.0006984708, 1.0014638383,
+      1.0022182509, 1.0029620257, 1.0036954662, 1.0044188628 },
+    { 0.8461555040, 0.8547882305, 0.8624635555, 0.8693789920,
+      0.8756760853, 0.8814598273, 0.8868103032, 0.8917900284,
+      0.8964487626, 0.9008267754, 0.9049571273, 0.9088673021,
+      0.9125804007, 0.9161160306, 0.9194909803, 0.9227197376,
+      0.9258148939, 0.9287874629, 0.9316471355, 0.9344024839,
+      0.9370611291, 0.9396298766, 0.9421148300, 0.9445214846,
+      0.9468548060, 0.9491192967, 0.9513190517, 0.9534578074,
+      0.9555389816, 0.9575657096, 0.9595408742, 0.9614671327,
+      0.9633469396, 0.9651825670, 0.9669761222, 0.9687295635,
+      0.9704447142, 0.9721232742, 0.9737668316, 0.9753768718,
+      0.9769547868, 0.9785018824, 0.9800193854, 0.9815084500,
+      0.9829701633, 0.9844055505, 0.9858155796, 0.9872011653,
+      0.9885631734, 0.9899024236, 0.9912196934, 0.9925157203,
+      0.9937912053, 0.9950468143, 0.9962831814, 0.9975009102,
+      0.9987005760, 0.9998827277, 1.0010478892, 1.0021965608,
+      1.0033292209, 1.0044463270, 1.0055483173, 1.0066356112 },
+    { 0.8003259737, 0.8112313241, 0.8209581209, 0.8297466775,
+      0.8377697066, 0.8451556492, 0.8520027051, 0.8583876935,
+      0.8643718792, 0.8700049328, 0.8753277020, 0.8803741979,
+      0.8851730502, 0.8897485937, 0.8941216918, 0.8983103719,
+      0.9023303202, 0.9061952736, 0.9099173316, 0.9135072091,
+      0.9169744409, 0.9203275502, 0.9235741882, 0.9267212496,
+      0.9297749699, 0.9327410079, 0.9356245146, 0.9384301933,
+      0.9411623497, 0.9438249364, 0.9464215906, 0.9489556668,
+      0.9514302661, 0.9538482608, 0.9562123167, 0.9585249126,
+      0.9607883576, 0.9630048062, 0.9651762722, 0.9673046403,
+      0.9693916775, 0.9714390425, 0.9734482944, 0.9754209007,
+      0.9773582446, 0.9792616307, 0.9811322918, 0.9829713934,
+      0.9847800389, 0.9865592739, 0.9883100900, 0.9900334289,
+      0.9917301853, 0.9934012104, 0.9950473143, 0.9966692689,
+      0.9982678100, 0.9998436400, 1.0013974295, 1.0029298194,
+      1.0044414224, 1.0059328250, 1.0074045889, 1.0088572520 },
+    { 0.7569786654, 0.7698939195, 0.7814501054, 0.7919210783,
+      0.8015042240, 0.8103467104, 0.8185613167, 0.8262364557,
+      0.8334427763, 0.8402376615, 0.8466683811, 0.8527743561,
+      0.8585888194, 0.8641400582, 0.8694523567, 0.8745467247,
+      0.8794414652, 0.8841526254, 0.8886943552, 0.8930791981,
+      0.8973183276, 0.9014217415, 0.9053984227, 0.9092564737,
+      0.9130032283, 0.9166453478, 0.9201889007, 0.9236394320,
+      0.9270020224, 0.9302813390, 0.9334816797, 0.9366070112,
+      0.9396610028, 0.9426470554, 0.9455683275, 0.9484277579,
+      0.9512280860, 0.9539718690, 0.9566614986, 0.9592992147,
+      0.9618871182, 0.9644271823, 0.9669212630, 0.9693711079,
+      0.9717783651, 0.9741445900, 0.9764712529, 0.9787597445,
+      0.9810113822, 0.9832274148, 0.9854090274, 0.9875573457,
+      0.9896734398, 0.9917583281, 0.9938129803, 0.9958383209,
+      0.9978352315, 0.9998045539, 1.0017470919, 1.0036636145,
+      1.0055548568, 1.0074215229, 1.0092642871, 1.0110837959 },
+    { 0.7159791370, 0.7306629191, 0.7438433845, 0.7558198318,
+      0.7668086064, 0.7769714272, 0.7864325139, 0.7952894548,
+      0.8036203840, 0.8114888792, 0.8189474022, 0.8260397728,
+      0.8328029877, 0.8392685815, 0.8454636629, 0.8514117142,
+      0.8571332177, 0.8626461513, 0.8679663850, 0.8731080020,
+      0.8780835596, 0.8829043049, 0.8875803529, 0.8921208349,
+      0.8965340237, 0.9008274393, 0.9050079382, 0.9090817905,
+      0.9130547454, 0.9169320882, 0.9207186893, 0.9244190474,
+      0.9280373261, 0.9315773876, 0.9350428208, 0.9384369673,
+      0.9417629433, 0.9450236603, 0.9482218422, 0.9513600421,
+      0.9544406555, 0.9574659338, 0.9604379957, 0.9633588374,
+      0.9662303420, 0.9690542879, 0.9718323569, 0.9745661408,
+      0.9772571477, 0.9799068082, 0.9825164805, 0.9850874551,
+      0.9876209597, 0.9901181627, 0.9925801775, 0.9950080658,
+      0.9974028405, 0.9997654692, 1.0020968764, 1.0043979464,
+      1.0066695255, 1.0089124239, 1.0111274185, 1.0133152537 },
+    { 0.6772002277, 0.6934309881, 0.7080464599, 0.7213643301,
+      0.7336148970, 0.7449707526, 0.7555647772, 0.7655015856,
+      0.7748651015, 0.7837237382, 0.7921340426, 0.8001433220,
+      0.8077915768, 0.8151129499, 0.8221368310, 0.8288887107,
+      0.8353908496, 0.8416628090, 0.8477218755, 0.8535834053,
+      0.8592611049, 0.8647672624, 0.8701129393, 0.8753081305,
+      0.8803618988, 0.8852824894, 0.8900774261, 0.8947535945,
+      0.8993173131, 0.9037743949, 0.9081302004, 0.9123896841,
+      0.9165574352, 0.9206377129, 0.9246344779, 0.9285514202,
+      0.9323919830, 0.9361593853, 0.9398566405, 0.9434865742,
+      0.9470518396, 0.9505549317, 0.9539981992, 0.9573838564,
+      0.9607139933, 0.9639905847, 0.9672154989, 0.9703905051,
+      0.9735172803, 0.9765974162, 0.9796324243, 0.9826237418,
+      0.9855727362, 0.9884807098, 0.9913489039, 0.9941785028,
+      0.9969706369, 0.9997263861, 1.0024467831, 1.0051328157,
+      1.0077854297, 1.0104055314, 1.0129939892, 1.0155516364 },
+    { 0.6405216642, 0.6580962612, 0.6739722363, 0.6884795488,
+      0.7018580813, 0.7142880714, 0.7259086094, 0.7368294324,
+      0.7471387455, 0.7569085832, 0.7661985859, 0.7750587283,
+      0.7835313288, 0.7916525600, 0.7994535998, 0.8069615243,
+      0.8142000068, 0.8211898738, 0.8279495504, 0.8344954211,
+      0.8408421252, 0.8470027997, 0.8529892811, 0.8588122744,
+      0.8644814947, 0.8700057878, 0.8753932324, 0.8806512276,
+      0.8857865684, 0.8908055105, 0.8957138271, 0.9005168576,
+      0.9052195513, 0.9098265046, 0.9143419945, 0.9187700080,
+      0.9231142680, 0.9273782568, 0.9315652364, 0.9356782672,
+      0.9397202245, 0.9436938133, 0.9476015819, 0.9514459336,
+      0.9552291382, 0.9589533414, 0.9626205741, 0.9662327603,
+      0.9697917251, 0.9732992008, 0.9767568340, 0.9801661903,
+      0.9835287605, 0.9868459649, 0.9901191578, 0.9933496315,
+      0.9965386205, 0.9996873045, 1.0027968119, 1.0058682226,
+      1.0089025710, 1.0119008485, 1.0148640056, 1.0177929548 },
+    { 0.6058296875, 0.6245620637, 0.6415378101, 0.6570938835,
+      0.6714759586, 0.6848691001, 0.6974164561, 0.7092312055,
+      0.7204044988, 0.7310109103, 0.7411122884, 0.7507605397,
+      0.7599996842, 0.7688674015, 0.7773962122, 0.7856143935,
+      0.7935466990, 0.8012149303, 0.8086383963, 0.8158342858,
+      0.8228179717, 0.8296032631, 0.8362026133, 0.8426272954,
+      0.8488875492, 0.8549927056, 0.8609512936, 0.8667711307,
+      0.8724594015, 0.8780227256, 0.8834672161, 0.8887985309,
+      0.8940219180, 0.8991422543, 0.9041640810, 0.9090916337,
+      0.9139288704, 0.9186794948, 0.9233469789, 0.9279345818,
+      0.9324453671, 0.9368822185, 0.9412478543, 0.9455448393,
+      0.9497755970, 0.9539424198, 0.9580474782, 0.9620928299,
+      0.9660804271, 0.9700121244, 0.9738896845, 0.9777147851,
+      0.9814890239, 0.9852139236, 0.9888909370, 0.9925214512,
+      0.9961067913, 0.9996482244, 1.0031469629, 1.0066041676,
+      1.0100209506, 1.0133983785, 1.0167374742, 1.0200392198 },
+    { 0.5730166999, 0.5927366473, 0.6106642672, 0.6271389942,
+      0.6424090212, 0.6566617910, 0.6700426292, 0.6826666808,
+      0.6946268614, 0.7059993279, 0.7168473476, 0.7272241023,
+      0.7371747608, 0.7467380401, 0.7559474006, 0.7648319736,
+      0.7734172908, 0.7817258650, 0.7897776570, 0.7975904541,
+      0.8051801811, 0.8125611560, 0.8197463039, 0.8267473349,
+      0.8335748949, 0.8402386937, 0.8467476129, 0.8531098003,
+      0.8593327495, 0.8654233698, 0.8713880464, 0.8772326935,
+      0.8829628002, 0.8885834710, 0.8940994619, 0.8995152120,
+      0.9048348715, 0.9100623268, 0.9152012229, 0.9202549833,
+      0.9252268281, 0.9301197899, 0.9349367288, 0.9396803449,
+      0.9443531909, 0.9489576823, 0.9534961076, 0.9579706374,
+      0.9623833320, 0.9667361492, 0.9710309512, 0.9752695109,
+      0.9794535174, 0.9835845813, 0.9876642399, 0.9916939614,
+      0.9956751493, 0.9996091459, 1.0034972362, 1.0073406510,
+      1.0111405700, 1.0148981248, 1.0186144013, 1.0222904422 },
+    { 0.5419809316, 0.5625329386, 0.5812764912, 0.5985496562,
+      0.6146003370, 0.6296162401, 0.6437432340, 0.6570971404,
+      0.6697716039, 0.6818435182, 0.6933768712, 0.7044255353,
+      0.7150353340, 0.7252456009, 0.7350903742, 0.7445993259,
+      0.7537984929, 0.7627108595, 0.7713568269, 0.7797545943,
+      0.7879204712, 0.7958691361, 0.8036138516, 0.8111666444,
+      0.8185384580, 0.8257392814, 0.8327782597, 0.8396637886,
+      0.8464035955, 0.8530048108, 0.8594740287, 0.8658173611,
+      0.8720404845, 0.8781486812, 0.8841468762, 0.8900396688,
+      0.8958313620, 0.9015259874, 0.9071273286, 0.9126389413,
+      0.9180641715, 0.9234061727, 0.9286679198, 0.9338522236,
+      0.9389617420, 0.9439989920, 0.9489663591, 0.9538661069,
+      0.9587003852, 0.9634712378, 0.9681806094, 0.9728303524,
+      0.9774222323, 0.9819579336, 0.9864390644, 0.9908671615,
+      0.9952436943, 0.9995700689, 1.0038476318, 1.0080776733,
+      1.0122614305, 1.0164000906, 1.0204947932, 1.0245466331 },
+    { 0.5126261246, 0.5338683013, 0.5533029807, 0.5712636181,
+      0.5879954388, 0.6036845987, 0.6184760989, 0.6324853169,
+      0.6458057215, 0.6585142011, 0.6706748475, 0.6823417062,
+      0.6935608163, 0.7043717519, 0.7148088052, 0.7249019070,
+      0.7346773529, 0.7441583823, 0.7533656456, 0.7623175831,
+      0.7710307376, 0.7795200117, 0.7877988829, 0.7958795841,
+      0.8037732557, 0.8114900754, 0.8190393682, 0.8264297018,
+      0.8336689680, 0.8407644543, 0.8477229049, 0.8545505751,
+      0.8612532786, 0.8678364291, 0.8743050768, 0.8806639416,
+      0.8869174414, 0.8930697184, 0.8991246621, 0.9050859297,
+      0.9109569648, 0.9167410144, 0.9224411436, 0.9280602496,
+      0.9336010737, 0.9390662129, 0.9444581300, 0.9497791628,
+      0.9550315328, 0.9602173528, 0.9653386345, 0.9703972943,
+      0.9753951600, 0.9803339761, 0.9852154088, 0.9900410510,
+      0.9948124263, 0.9995309934, 1.0041981497, 1.0088152348,
+      1.0133835335, 1.0179042791, 1.0223786564, 1.0268078035 },
+};
+
 #endif /* AVCODEC_WMAVOICE_DATA_H */