/*

 * DCA encoder

 * Copyright (C) 2008-2012 Alexander E. Patrakov

 *               2010 Benjamin Larsson

 *               2011 Xiang Wang

 *

 * This file is part of FFmpeg.

 *

 * FFmpeg 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.

 *

 * FFmpeg 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 FFmpeg; if not, write to the Free Software

 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA

 */

#include "libavutil/avassert.h"

#include "libavutil/channel_layout.h"

#include "libavutil/common.h"

#include "libavutil/ffmath.h"

#include "libavutil/mem.h"

#include "libavutil/mem_internal.h"

#include "libavutil/opt.h"

#include "libavutil/thread.h"

#include "libavutil/tx.h"

#include "avcodec.h"

#include "codec_internal.h"

#include "dcaadpcm.h"

#include "dcamath.h"

#include "dca_core.h"

#include "dcadata.h"

#include "dcaenc.h"

#include "encode.h"

#include "put_bits.h"

#define MAX_CHANNELS 6

#define DCA_MAX_FRAME_SIZE 16384

#define DCA_HEADER_SIZE 13

#define DCA_LFE_SAMPLES 8

#define DCAENC_SUBBANDS 32

#define SUBFRAMES 1

#define SUBSUBFRAMES 2

#define SUBBAND_SAMPLES (SUBFRAMES * SUBSUBFRAMES * 8)

#define AUBANDS 25

#define COS_T(x) (c->cos_table[(x) & 2047])

typedef struct CompressionOptions {

    int adpcm_mode;

} CompressionOptions;

typedef struct DCAEncContext {

    AVClass *class;

    PutBitContext pb;

    DCAADPCMEncContext adpcm_ctx;

    AVTXContext *mdct;

    av_tx_fn mdct_fn;

    CompressionOptions options;

    int frame_size;

    int frame_bits;

    int fullband_channels;

    int channels;

    int lfe_channel;

    int samplerate_index;

    int bitrate_index;

    int channel_config;

    const int32_t *band_interpolation;

    const int32_t *band_spectrum;

    int lfe_scale_factor;

    softfloat lfe_quant;

    int32_t lfe_peak_cb;

    const int8_t *channel_order_tab;  ///< channel reordering table, lfe and non lfe

    int32_t prediction_mode[MAX_CHANNELS][DCAENC_SUBBANDS];

    int32_t adpcm_history[MAX_CHANNELS][DCAENC_SUBBANDS][DCA_ADPCM_COEFFS * 2];

    int32_t history[MAX_CHANNELS][512]; /* This is a circular buffer */

    int32_t *subband[MAX_CHANNELS][DCAENC_SUBBANDS];

    int32_t quantized[MAX_CHANNELS][DCAENC_SUBBANDS][SUBBAND_SAMPLES];

    int32_t peak_cb[MAX_CHANNELS][DCAENC_SUBBANDS];

    int32_t diff_peak_cb[MAX_CHANNELS][DCAENC_SUBBANDS]; ///< expected peak of residual signal

    int32_t downsampled_lfe[DCA_LFE_SAMPLES];

    int32_t masking_curve_cb[SUBSUBFRAMES][256];

    int32_t bit_allocation_sel[MAX_CHANNELS];

    int abits[MAX_CHANNELS][DCAENC_SUBBANDS];

    int scale_factor[MAX_CHANNELS][DCAENC_SUBBANDS];

    softfloat quant[MAX_CHANNELS][DCAENC_SUBBANDS];

    int32_t quant_index_sel[MAX_CHANNELS][DCA_CODE_BOOKS];

    int32_t eff_masking_curve_cb[256];

    int32_t band_masking_cb[32];

    int32_t worst_quantization_noise;

    int32_t worst_noise_ever;

    int consumed_bits;

    int consumed_adpcm_bits; ///< Number of bits to transmit ADPCM related info

    int32_t cos_table[2048];

    int32_t band_interpolation_tab[2][512];

    int32_t band_spectrum_tab[2][8];

    int32_t auf[9][AUBANDS][256];

    int32_t cb_to_add[256];

    int32_t cb_to_level[2048];

    int32_t lfe_fir_64i[512];

} DCAEncContext;

/* Transfer function of outer and middle ear, Hz -> dB */

static double hom(double f)

{

    double f1 = f / 1000;

    return -3.64 * pow(f1, -0.8)

           + 6.8 * exp(-0.6 * (f1 - 3.4) * (f1 - 3.4))

           - 6.0 * exp(-0.15 * (f1 - 8.7) * (f1 - 8.7))

           - 0.0006 * (f1 * f1) * (f1 * f1);

}

static double gammafilter(int i, double f)

{

    double h = (f - fc[i]) / erb[i];

    h = 1 + h * h;

    h = 1 / (h * h);

    return 20 * log10(h);

}

static int subband_bufer_alloc(DCAEncContext *c)

{

    int ch, band;

    int32_t *bufer = av_calloc(MAX_CHANNELS * DCAENC_SUBBANDS *

                               (SUBBAND_SAMPLES + DCA_ADPCM_COEFFS),

                               sizeof(int32_t));

    if (!bufer)

        return AVERROR(ENOMEM);

    /* we need a place for DCA_ADPCM_COEFF samples from previous frame

     * to calc prediction coefficients for each subband */

    for (ch = 0; ch < MAX_CHANNELS; ch++) {

        for (band = 0; band < DCAENC_SUBBANDS; band++) {

            c->subband[ch][band] = bufer +

                                   ch * DCAENC_SUBBANDS * (SUBBAND_SAMPLES + DCA_ADPCM_COEFFS) +

                                   band * (SUBBAND_SAMPLES + DCA_ADPCM_COEFFS) + DCA_ADPCM_COEFFS;

        }

    }

    return 0;

}

static void subband_bufer_free(DCAEncContext *c)

{

    if (c->subband[0][0]) {

        int32_t *bufer = c->subband[0][0] - DCA_ADPCM_COEFFS;

        av_free(bufer);

        c->subband[0][0] = NULL;

    }

}

static uint16_t bitalloc_12_table[DCA_BITALLOC_12_COUNT][12 + 1][2];

static uint16_t bitalloc_table[DCA_NUM_BITALLOC_CODES][2];

static const uint16_t (*bitalloc_tables[DCA_CODE_BOOKS][8])[2];

static av_cold void create_enc_table(uint16_t dst[][2], unsigned count,

                                     const uint8_t (**src_tablep)[2])

{

    const uint8_t (*src_table)[2] = *src_tablep;

    uint16_t code = 0;

    for (unsigned i = 0; i < count; i++) {

        unsigned dst_idx = src_table[i][0];

        dst[dst_idx][0] = code >> (16 - src_table[i][1]);

        dst[dst_idx][1] = src_table[i][1];

        code += 1 << (16 - src_table[i][1]);

    }

    *src_tablep += count;

}

static av_cold void dcaenc_init_static_tables(void)

{

    uint16_t (*bitalloc_dst)[2] = bitalloc_table;

    const uint8_t (*src_table)[2] = ff_dca_vlc_src_tables;

    for (unsigned i = 0; i < DCA_CODE_BOOKS; i++) {

        for (unsigned j = 0; j < ff_dca_quant_index_group_size[i]; j++) {

            create_enc_table(bitalloc_dst, ff_dca_bitalloc_sizes[i],

                             &src_table);

            bitalloc_tables[i][j] = bitalloc_dst - ff_dca_bitalloc_offsets[i];

            bitalloc_dst += ff_dca_bitalloc_sizes[i];

        }

    }

    for (unsigned i = 0; i < DCA_BITALLOC_12_COUNT; i++)

        create_enc_table(&bitalloc_12_table[i][1], 12, &src_table);

}

static av_cold int encode_init(AVCodecContext *avctx)

{

    static AVOnce init_static_once = AV_ONCE_INIT;

    DCAEncContext *c = avctx->priv_data;

    AVChannelLayout layout = avctx->ch_layout;

    int i, j, k, min_frame_bits;

    float scale = 1.0f;

    int ret;

    if ((ret = subband_bufer_alloc(c)) < 0)

        return ret;

    c->fullband_channels = c->channels = layout.nb_channels;

    c->lfe_channel = (c->channels == 3 || c->channels == 6);

    c->band_interpolation = c->band_interpolation_tab[1];

    c->band_spectrum = c->band_spectrum_tab[1];

    c->worst_quantization_noise = -2047;

    c->worst_noise_ever = -2047;

    c->consumed_adpcm_bits = 0;

    if (ff_dcaadpcm_init(&c->adpcm_ctx))

        return AVERROR(ENOMEM);

    switch (layout.nb_channels) {

    case 1: /* mono */

        c->channel_config = 0;

        break;

    case 2: /* stereo */

        c->channel_config = 2;

        break;

    case 4: /* 2.2 */

        c->channel_config = 8;

        break;

    case 5: /* 5.0 */

        c->channel_config = 9;

        break;

    case 6: /* 5.1 */

        c->channel_config = 9;

        break;

    default:

        av_assert1(!"impossible channel layout");

    }

    if (c->lfe_channel) {

        c->fullband_channels--;

        c->channel_order_tab = channel_reorder_lfe[c->channel_config];

    } else {

        c->channel_order_tab = channel_reorder_nolfe[c->channel_config];

    }

    for (i = 0; i < MAX_CHANNELS; i++) {

        for (j = 0; j < DCA_CODE_BOOKS; j++) {

            c->quant_index_sel[i][j] = ff_dca_quant_index_group_size[j];

        }

        /* 6 - no Huffman */

        c->bit_allocation_sel[i] = 6;

        for (j = 0; j < DCAENC_SUBBANDS; j++) {

            /* -1 - no ADPCM */

            c->prediction_mode[i][j] = -1;

            memset(c->adpcm_history[i][j], 0, sizeof(int32_t)*DCA_ADPCM_COEFFS);

        }

    }

    for (i = 0; i < 9; i++) {

        if (sample_rates[i] == avctx->sample_rate)

            break;

    }

    if (i == 9)

        return AVERROR(EINVAL);

    c->samplerate_index = i;

    if (avctx->bit_rate < 32000 || avctx->bit_rate > 3840000) {

        av_log(avctx, AV_LOG_ERROR, "Bit rate %"PRId64" not supported.", avctx->bit_rate);

        return AVERROR(EINVAL);

    }

    for (i = 0; ff_dca_bit_rates[i] < avctx->bit_rate; i++)

        ;

    c->bitrate_index = i;

    c->frame_bits = FFALIGN((avctx->bit_rate * 512 + avctx->sample_rate - 1) / avctx->sample_rate, 32);

    min_frame_bits = 132 + (493 + 28 * 32) * c->fullband_channels + c->lfe_channel * 72;

    if (c->frame_bits < min_frame_bits || c->frame_bits > (DCA_MAX_FRAME_SIZE << 3))

        return AVERROR(EINVAL);

    c->frame_size = (c->frame_bits + 7) / 8;

    avctx->frame_size = 32 * SUBBAND_SAMPLES;

    if ((ret = av_tx_init(&c->mdct, &c->mdct_fn, AV_TX_INT32_MDCT, 0, 256, &scale, 0)) < 0)

        return ret;

    /* Init all tables */

    c->cos_table[0] = 0x7fffffff;

    c->cos_table[512] = 0;

    c->cos_table[1024] = -c->cos_table[0];

    for (i = 1; i < 512; i++) {

        c->cos_table[i]   = (int32_t)(0x7fffffff * cos(M_PI * i / 1024));

        c->cos_table[1024-i] = -c->cos_table[i];

        c->cos_table[1024+i] = -c->cos_table[i];

        c->cos_table[2048-i] = +c->cos_table[i];

    }

    for (i = 0; i < 2048; i++)

        c->cb_to_level[i] = (int32_t)(0x7fffffff * ff_exp10(-0.005 * i));

    for (k = 0; k < 32; k++) {

        for (j = 0; j < 8; j++) {

            c->lfe_fir_64i[64 * j + k] = (int32_t)(0xffffff800000ULL * ff_dca_lfe_fir_64[8 * k + j]);

            c->lfe_fir_64i[64 * (7-j) + (63 - k)] = (int32_t)(0xffffff800000ULL * ff_dca_lfe_fir_64[8 * k + j]);

        }

    }

    for (i = 0; i < 512; i++) {

        c->band_interpolation_tab[0][i] = (int32_t)(0x1000000000ULL * ff_dca_fir_32bands_perfect[i]);

        c->band_interpolation_tab[1][i] = (int32_t)(0x1000000000ULL * ff_dca_fir_32bands_nonperfect[i]);

    }

    for (i = 0; i < 9; i++) {

        for (j = 0; j < AUBANDS; j++) {

            for (k = 0; k < 256; k++) {

                double freq = sample_rates[i] * (k + 0.5) / 512;

                c->auf[i][j][k] = (int32_t)(10 * (hom(freq) + gammafilter(j, freq)));

            }

        }

    }

    for (i = 0; i < 256; i++) {

        double add = 1 + ff_exp10(-0.01 * i);

        c->cb_to_add[i] = (int32_t)(100 * log10(add));

    }

    for (j = 0; j < 8; j++) {

        double accum = 0;

        for (i = 0; i < 512; i++) {

            double reconst = ff_dca_fir_32bands_perfect[i] * ((i & 64) ? (-1) : 1);

            accum += reconst * cos(2 * M_PI * (i + 0.5 - 256) * (j + 0.5) / 512);

        }

        c->band_spectrum_tab[0][j] = (int32_t)(200 * log10(accum));

    }

    for (j = 0; j < 8; j++) {

        double accum = 0;

        for (i = 0; i < 512; i++) {

            double reconst = ff_dca_fir_32bands_nonperfect[i] * ((i & 64) ? (-1) : 1);

            accum += reconst * cos(2 * M_PI * (i + 0.5 - 256) * (j + 0.5) / 512);

        }

        c->band_spectrum_tab[1][j] = (int32_t)(200 * log10(accum));

    }

    ff_thread_once(&init_static_once, dcaenc_init_static_tables);

    return 0;

}

static av_cold int encode_close(AVCodecContext *avctx)

{

    DCAEncContext *c = avctx->priv_data;

    av_tx_uninit(&c->mdct);

    subband_bufer_free(c);

    ff_dcaadpcm_free(&c->adpcm_ctx);

    return 0;

}

static void subband_transform(DCAEncContext *c, const int32_t *input)

{

    int ch, subs, i, k, j;

    for (ch = 0; ch < c->fullband_channels; ch++) {

        /* History is copied because it is also needed for PSY */

        int32_t hist[512];

        int hist_start = 0;

        const int chi = c->channel_order_tab[ch];

        memcpy(hist, &c->history[ch][0], 512 * sizeof(int32_t));

        for (subs = 0; subs < SUBBAND_SAMPLES; subs++) {

            int32_t accum[64];

            int32_t resp;

            int band;

            /* Calculate the convolutions at once */

            memset(accum, 0, 64 * sizeof(int32_t));

            for (k = 0, i = hist_start, j = 0;

                    i < 512; k = (k + 1) & 63, i++, j++)

                accum[k] += mul32(hist[i], c->band_interpolation[j]);

            for (i = 0; i < hist_start; k = (k + 1) & 63, i++, j++)

                accum[k] += mul32(hist[i], c->band_interpolation[j]);

            for (k = 16; k < 32; k++)

                accum[k] = accum[k] - accum[31 - k];

            for (k = 32; k < 48; k++)

                accum[k] = accum[k] + accum[95 - k];

            for (band = 0; band < 32; band++) {

                resp = 0;

                for (i = 16; i < 48; i++) {

                    int s = (2 * band + 1) * (2 * (i + 16) + 1);

                    resp += mul32(accum[i], COS_T(s << 3)) >> 3;

                }

                c->subband[ch][band][subs] = ((band + 1) & 2) ? -resp : resp;

            }

            /* Copy in 32 new samples from input */

            for (i = 0; i < 32; i++)

                hist[i + hist_start] = input[(subs * 32 + i) * c->channels + chi];

            hist_start = (hist_start + 32) & 511;

        }

    }

}

static void lfe_downsample(DCAEncContext *c, const int32_t *input)

{

    /* FIXME: make 128x LFE downsampling possible */

    const int lfech = lfe_index[c->channel_config];

    int i, j, lfes;

    int32_t hist[512];

    int32_t accum;

    int hist_start = 0;

    memcpy(hist, &c->history[c->channels - 1][0], 512 * sizeof(int32_t));

    for (lfes = 0; lfes < DCA_LFE_SAMPLES; lfes++) {

        /* Calculate the convolution */

        accum = 0;

        for (i = hist_start, j = 0; i < 512; i++, j++)

            accum += mul32(hist[i], c->lfe_fir_64i[j]);

        for (i = 0; i < hist_start; i++, j++)

            accum += mul32(hist[i], c->lfe_fir_64i[j]);

        c->downsampled_lfe[lfes] = accum;

        /* Copy in 64 new samples from input */

        for (i = 0; i < 64; i++)

            hist[i + hist_start] = input[(lfes * 64 + i) * c->channels + lfech];

        hist_start = (hist_start + 64) & 511;

    }

}

static uint32_t dca_vlc_calc_alloc_bits(const int values[], uint8_t n, uint8_t sel)

{

    uint32_t sum = 0;

    for (unsigned i = 0; i < n; i++)

        sum += bitalloc_12_table[sel][values[i]][1];

    return sum;

}

static void dca_vlc_enc_alloc(PutBitContext *pb, const int values[],

                              uint8_t n, uint8_t sel)

{

    for (unsigned i = 0; i < n; i++)

        put_bits(pb, bitalloc_12_table[sel][values[i]][1],

                     bitalloc_12_table[sel][values[i]][0]);

}

static uint32_t dca_vlc_calc_quant_bits(const int values[], uint8_t n,

                                        uint8_t sel, uint8_t table)

{

    uint32_t sum = 0;

    for (unsigned i = 0; i < n; i++)

        sum += bitalloc_tables[table][sel][values[i]][1];

    return sum;

}

static void dca_vlc_enc_quant(PutBitContext *pb, const int values[],

                              uint8_t n, uint8_t sel, uint8_t table)

{

    for (unsigned i = 0; i < n; i++)

        put_bits(pb, bitalloc_tables[table][sel][values[i]][1],

                     bitalloc_tables[table][sel][values[i]][0]);

}

static int32_t get_cb(DCAEncContext *c, int32_t in)

{

    int i, res = 0;

    in = FFABS(in);

    for (i = 1024; i > 0; i >>= 1) {

        if (c->cb_to_level[i + res] >= in)

            res += i;

    }

    return -res;

}

static int32_t add_cb(DCAEncContext *c, int32_t a, int32_t b)

{

    if (a < b)

        FFSWAP(int32_t, a, b);

    if (a - b >= 256)

        return a;

    return a + c->cb_to_add[a - b];

}

static void calc_power(DCAEncContext *c,

                       const int32_t in[2 * 256], int32_t power[256])

{

    int i;

    LOCAL_ALIGNED_32(int32_t, data,  [512]);

    LOCAL_ALIGNED_32(int32_t, coeff, [256]);

    for (i = 0; i < 512; i++)

        data[i] = norm__(mul32(in[i], 0x3fffffff - (COS_T(4 * i + 2) >> 1)), 4);

    c->mdct_fn(c->mdct, coeff, data, sizeof(int32_t));

    for (i = 0; i < 256; i++) {

        const int32_t cb = get_cb(c, coeff[i]);

        power[i] = add_cb(c, cb, cb);

    }

}

static void adjust_jnd(DCAEncContext *c,

                       const int32_t in[512], int32_t out_cb[256])

{

    int32_t power[256];

    int32_t out_cb_unnorm[256];

    int32_t denom;

    const int32_t ca_cb = -1114;

    const int32_t cs_cb = 928;

    const int samplerate_index = c->samplerate_index;

    int i, j;

    calc_power(c, in, power);

    for (j = 0; j < 256; j++)

        out_cb_unnorm[j] = -2047; /* and can only grow */

    for (i = 0; i < AUBANDS; i++) {

        denom = ca_cb; /* and can only grow */

        for (j = 0; j < 256; j++)

            denom = add_cb(c, denom, power[j] + c->auf[samplerate_index][i][j]);

        for (j = 0; j < 256; j++)

            out_cb_unnorm[j] = add_cb(c, out_cb_unnorm[j],

                                      -denom + c->auf[samplerate_index][i][j]);

    }

    for (j = 0; j < 256; j++)

        out_cb[j] = add_cb(c, out_cb[j], -out_cb_unnorm[j] - ca_cb - cs_cb);

}

typedef void (*walk_band_t)(DCAEncContext *c, int band1, int band2, int f,

                            int32_t spectrum1, int32_t spectrum2, int channel,

                            int32_t * arg);

static void walk_band_low(DCAEncContext *c, int band, int channel,

                          walk_band_t walk, int32_t *arg)

{

    int f;

    if (band == 0) {

        for (f = 0; f < 4; f++)

            walk(c, 0, 0, f, 0, -2047, channel, arg);

    } else {

        for (f = 0; f < 8; f++)

            walk(c, band, band - 1, 8 * band - 4 + f,

                    c->band_spectrum[7 - f], c->band_spectrum[f], channel, arg);

    }

}

static void walk_band_high(DCAEncContext *c, int band, int channel,

                           walk_band_t walk, int32_t *arg)

{

    int f;

    if (band == 31) {

        for (f = 0; f < 4; f++)

            walk(c, 31, 31, 256 - 4 + f, 0, -2047, channel, arg);

    } else {

        for (f = 0; f < 8; f++)

            walk(c, band, band + 1, 8 * band + 4 + f,

                    c->band_spectrum[f], c->band_spectrum[7 - f], channel, arg);

    }

}

static void update_band_masking(DCAEncContext *c, int band1, int band2,

                                int f, int32_t spectrum1, int32_t spectrum2,

                                int channel, int32_t * arg)

{

    int32_t value = c->eff_masking_curve_cb[f] - spectrum1;

    if (value < c->band_masking_cb[band1])

        c->band_masking_cb[band1] = value;

}

static void calc_masking(DCAEncContext *c, const int32_t *input)

{

    int i, k, band, ch, ssf;

    int32_t data[512];

    for (i = 0; i < 256; i++)

        for (ssf = 0; ssf < SUBSUBFRAMES; ssf++)

            c->masking_curve_cb[ssf][i] = -2047;

    for (ssf = 0; ssf < SUBSUBFRAMES; ssf++)

        for (ch = 0; ch < c->fullband_channels; ch++) {

            const int chi = c->channel_order_tab[ch];

            for (i = 0, k = 128 + 256 * ssf; k < 512; i++, k++)

                data[i] = c->history[ch][k];

            for (k -= 512; i < 512; i++, k++)

                data[i] = input[k * c->channels + chi];

            adjust_jnd(c, data, c->masking_curve_cb[ssf]);

        }

    for (i = 0; i < 256; i++) {

        int32_t m = 2048;

        for (ssf = 0; ssf < SUBSUBFRAMES; ssf++)

            if (c->masking_curve_cb[ssf][i] < m)

                m = c->masking_curve_cb[ssf][i];

        c->eff_masking_curve_cb[i] = m;

    }

    for (band = 0; band < 32; band++) {

        c->band_masking_cb[band] = 2048;

        walk_band_low(c, band, 0, update_band_masking, NULL);

        walk_band_high(c, band, 0, update_band_masking, NULL);

    }

}

static inline int32_t find_peak(DCAEncContext *c, const int32_t *in, int len)

{

    int sample;

    int32_t m = 0;

    for (sample = 0; sample < len; sample++) {

        int32_t s = abs(in[sample]);

        if (m < s)

            m = s;

    }

    return get_cb(c, m);

}

static void find_peaks(DCAEncContext *c)

{

    int band, ch;

    for (ch = 0; ch < c->fullband_channels; ch++) {

        for (band = 0; band < 32; band++)

            c->peak_cb[ch][band] = find_peak(c, c->subband[ch][band],

                                             SUBBAND_SAMPLES);

    }

    if (c->lfe_channel)

        c->lfe_peak_cb = find_peak(c, c->downsampled_lfe, DCA_LFE_SAMPLES);

}

static void adpcm_analysis(DCAEncContext *c)

{

    int ch, band;

    int pred_vq_id;

    int32_t *samples;

    int32_t estimated_diff[SUBBAND_SAMPLES];

    c->consumed_adpcm_bits = 0;

    for (ch = 0; ch < c->fullband_channels; ch++) {

        for (band = 0; band < 32; band++) {

            samples = c->subband[ch][band] - DCA_ADPCM_COEFFS;

            pred_vq_id = ff_dcaadpcm_subband_analysis(&c->adpcm_ctx, samples,

                                                      SUBBAND_SAMPLES, estimated_diff);

            if (pred_vq_id >= 0) {

                c->prediction_mode[ch][band] = pred_vq_id;

                c->consumed_adpcm_bits += 12; //12 bits to transmit prediction vq index

                c->diff_peak_cb[ch][band] = find_peak(c, estimated_diff, 16);

            } else {

                c->prediction_mode[ch][band] = -1;

            }

        }

    }

}

static const int snr_fudge = 128;

#define USED_1ABITS 1

#define USED_26ABITS 4

static inline int32_t get_step_size(DCAEncContext *c, int ch, int band)

{

    int32_t step_size;

    if (c->bitrate_index == 3)

        step_size = ff_dca_lossless_quant[c->abits[ch][band]];

    else

        step_size = ff_dca_lossy_quant[c->abits[ch][band]];

    return step_size;

}

static int calc_one_scale(DCAEncContext *c, int32_t peak_cb, int abits,

                          softfloat *quant)

{

    int32_t peak;

    int our_nscale, try_remove;

    softfloat our_quant;

    av_assert0(peak_cb <= 0);

    av_assert0(peak_cb >= -2047);

    our_nscale = 127;

    peak = c->cb_to_level[-peak_cb];

    for (try_remove = 64; try_remove > 0; try_remove >>= 1) {

        if (scalefactor_inv[our_nscale - try_remove].e + stepsize_inv[abits].e <= 17)

            continue;

        our_quant.m = mul32(scalefactor_inv[our_nscale - try_remove].m, stepsize_inv[abits].m);

        our_quant.e = scalefactor_inv[our_nscale - try_remove].e + stepsize_inv[abits].e - 17;

        if ((ff_dca_quant_levels[abits] - 1) / 2 < quantize_value(peak, our_quant))

            continue;

        our_nscale -= try_remove;

    }

    if (our_nscale >= 125)

        our_nscale = 124;

    quant->m = mul32(scalefactor_inv[our_nscale].m, stepsize_inv[abits].m);

    quant->e = scalefactor_inv[our_nscale].e + stepsize_inv[abits].e - 17;

    av_assert0((ff_dca_quant_levels[abits] - 1) / 2 >= quantize_value(peak, *quant));

    return our_nscale;

}

static inline void quantize_adpcm_subband(DCAEncContext *c, int ch, int band)

{

    int32_t step_size;

    int32_t diff_peak_cb = c->diff_peak_cb[ch][band];

    c->scale_factor[ch][band] = calc_one_scale(c, diff_peak_cb,

                                               c->abits[ch][band],

                                               &c->quant[ch][band]);

    step_size = get_step_size(c, ch, band);

    ff_dcaadpcm_do_real(c->prediction_mode[ch][band],

                        c->quant[ch][band],

                        ff_dca_scale_factor_quant7[c->scale_factor[ch][band]],

                        step_size, c->adpcm_history[ch][band], c->subband[ch][band],

                        c->adpcm_history[ch][band] + 4, c->quantized[ch][band],

                        SUBBAND_SAMPLES, c->cb_to_level[-diff_peak_cb]);

}

static void quantize_adpcm(DCAEncContext *c)

{

    int band, ch;

    for (ch = 0; ch < c->fullband_channels; ch++)

        for (band = 0; band < 32; band++)

            if (c->prediction_mode[ch][band] >= 0)

                quantize_adpcm_subband(c, ch, band);

}

static void quantize_pcm(DCAEncContext *c)

{

    int sample, band, ch;

    for (ch = 0; ch < c->fullband_channels; ch++) {

        for (band = 0; band < 32; band++) {

            if (c->prediction_mode[ch][band] == -1) {

                for (sample = 0; sample < SUBBAND_SAMPLES; sample++) {

                    int32_t val = quantize_value(c->subband[ch][band][sample],

                                                 c->quant[ch][band]);

                    c->quantized[ch][band][sample] = val;

                }

            }

        }

    }

}

static void accumulate_huff_bit_consumption(int abits, int32_t *quantized,

                                            uint32_t *result)

{

    uint8_t sel, id = abits - 1;

    for (sel = 0; sel < ff_dca_quant_index_group_size[id]; sel++)

        result[sel] += dca_vlc_calc_quant_bits(quantized, SUBBAND_SAMPLES,

                                               sel, id);

}

static uint32_t set_best_code(uint32_t vlc_bits[DCA_CODE_BOOKS][7],

                              uint32_t clc_bits[DCA_CODE_BOOKS],

                              int32_t res[DCA_CODE_BOOKS])

{

    uint8_t i, sel;

    uint32_t best_sel_bits[DCA_CODE_BOOKS];

    int32_t best_sel_id[DCA_CODE_BOOKS];

    uint32_t t, bits = 0;

    for (i = 0; i < DCA_CODE_BOOKS; i++) {

        av_assert0(!((!!vlc_bits[i][0]) ^ (!!clc_bits[i])));

        if (vlc_bits[i][0] == 0) {

            /* do not transmit adjustment index for empty codebooks */

            res[i] = ff_dca_quant_index_group_size[i];

            /* and skip it */

            continue;

        }

        best_sel_bits[i] = vlc_bits[i][0];

        best_sel_id[i] = 0;

        for (sel = 0; sel < ff_dca_quant_index_group_size[i]; sel++) {

            if (best_sel_bits[i] > vlc_bits[i][sel] && vlc_bits[i][sel]) {

                best_sel_bits[i] = vlc_bits[i][sel];

                best_sel_id[i] = sel;

            }

        }

        /* 2 bits to transmit scale factor adjustment index */

        t = best_sel_bits[i] + 2;

        if (t < clc_bits[i]) {

            res[i] = best_sel_id[i];

            bits += t;

        } else {

            res[i] = ff_dca_quant_index_group_size[i];

            bits += clc_bits[i];

        }

    }

    return bits;

}

static uint32_t set_best_abits_code(int abits[DCAENC_SUBBANDS], int bands,

                                    int32_t *res)

{

    uint8_t i;

    uint32_t t;

    int32_t best_sel = 6;

    int32_t best_bits = bands * 5;

    /* Check do we have subband which cannot be encoded by Huffman tables */

    for (i = 0; i < bands; i++) {

        if (abits[i] > 12 || abits[i] == 0) {

            *res = best_sel;

            return best_bits;

        }

    }

    for (i = 0; i < DCA_BITALLOC_12_COUNT; i++) {

        t = dca_vlc_calc_alloc_bits(abits, bands, i);

        if (t < best_bits) {

            best_bits = t;

            best_sel = i;

        }

    }

    *res = best_sel;

    return best_bits;

}

static int init_quantization_noise(DCAEncContext *c, int noise, int forbid_zero)

{

    int ch, band, ret = USED_26ABITS | USED_1ABITS;

    uint32_t huff_bit_count_accum[MAX_CHANNELS][DCA_CODE_BOOKS][7];

    uint32_t clc_bit_count_accum[MAX_CHANNELS][DCA_CODE_BOOKS];

    uint32_t bits_counter = 0;

    c->consumed_bits = 132 + 333 * c->fullband_channels;

    c->consumed_bits += c->consumed_adpcm_bits;

    if (c->lfe_channel)

        c->consumed_bits += 72;

    /* attempt to guess the bit distribution based on the previous frame */

    for (ch = 0; ch < c->fullband_channels; ch++) {

        for (band = 0; band < 32; band++) {

            int snr_cb = c->peak_cb[ch][band] - c->band_masking_cb[band] - noise;

            if (snr_cb >= 1312) {

                c->abits[ch][band] = 26;