git.sesse.net Git - ffmpeg/blob - libavcodec/opus_celt.c

   1 /*
   2  * Copyright (c) 2012 Andrew D'Addesio
   3  * Copyright (c) 2013-2014 Mozilla Corporation
   4  *
   5  * This file is part of FFmpeg.
   6  *
   7  * FFmpeg is free software; you can redistribute it and/or
   8  * modify it under the terms of the GNU Lesser General Public
   9  * License as published by the Free Software Foundation; either
  10  * version 2.1 of the License, or (at your option) any later version.
  11  *
  12  * FFmpeg is distributed in the hope that it will be useful,
  13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  15  * Lesser General Public License for more details.
  16  *
  17  * You should have received a copy of the GNU Lesser General Public
  18  * License along with FFmpeg; if not, write to the Free Software
  19  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
  20  */
  21
  22 /**
  23  * @file
  24  * Opus CELT decoder
  25  */
  26
  27 #include <stdint.h>
  28
  29 #include "libavutil/float_dsp.h"
  30
  31 #include "opus.h"
  32
  33 enum CeltSpread {
  34     CELT_SPREAD_NONE,
  35     CELT_SPREAD_LIGHT,
  36     CELT_SPREAD_NORMAL,
  37     CELT_SPREAD_AGGRESSIVE
  38 };
  39
  40 typedef struct CeltFrame {
  41     float energy[CELT_MAX_BANDS];
  42     float prev_energy[2][CELT_MAX_BANDS];
  43
  44     uint8_t collapse_masks[CELT_MAX_BANDS];
  45
  46     /* buffer for mdct output + postfilter */
  47     DECLARE_ALIGNED(32, float, buf)[2048];
  48
  49     /* postfilter parameters */
  50     int pf_period_new;
  51     float pf_gains_new[3];
  52     int pf_period;
  53     float pf_gains[3];
  54     int pf_period_old;
  55     float pf_gains_old[3];
  56
  57     float deemph_coeff;
  58 } CeltFrame;
  59
  60 struct CeltContext {
  61     // constant values that do not change during context lifetime
  62     AVCodecContext    *avctx;
  63     CeltIMDCTContext  *imdct[4];
  64     AVFloatDSPContext  dsp;
  65     int output_channels;
  66
  67     // values that have inter-frame effect and must be reset on flush
  68     CeltFrame frame[2];
  69     uint32_t seed;
  70     int flushed;
  71
  72     // values that only affect a single frame
  73     int coded_channels;
  74     int framebits;
  75     int duration;
  76
  77     /* number of iMDCT blocks in the frame */
  78     int blocks;
  79     /* size of each block */
  80     int blocksize;
  81
  82     int startband;
  83     int endband;
  84     int codedbands;
  85
  86     int anticollapse_bit;
  87
  88     int intensitystereo;
  89     int dualstereo;
  90     enum CeltSpread spread;
  91
  92     int remaining;
  93     int remaining2;
  94     int fine_bits    [CELT_MAX_BANDS];
  95     int fine_priority[CELT_MAX_BANDS];
  96     int pulses       [CELT_MAX_BANDS];
  97     int tf_change    [CELT_MAX_BANDS];
  98
  99     DECLARE_ALIGNED(32, float, coeffs)[2][CELT_MAX_FRAME_SIZE];
 100     DECLARE_ALIGNED(32, float, scratch)[22 * 8]; // MAX(celt_freq_range) * 1<<CELT_MAX_LOG_BLOCKS
 101 };
 102
 103 static const uint16_t celt_model_tapset[] = { 4, 2, 3, 4 };
 104
 105 static const uint16_t celt_model_spread[] = { 32, 7, 9, 30, 32 };
 106
 107 static const uint16_t celt_model_alloc_trim[] = {
 108     128,   2,   4,   9,  19,  41,  87, 109, 119, 124, 126, 128
 109 };
 110
 111 static const uint16_t celt_model_energy_small[] = { 4, 2, 3, 4 };
 112
 113 static const uint8_t celt_freq_bands[] = { /* in steps of 200Hz */
 114     0,  1,  2,  3,  4,  5,  6,  7,  8, 10, 12, 14, 16, 20, 24, 28, 34, 40, 48, 60, 78, 100
 115 };
 116
 117 static const uint8_t celt_freq_range[] = {
 118     1,  1,  1,  1,  1,  1,  1,  1,  2,  2,  2,  2,  4,  4,  4,  6,  6,  8, 12, 18, 22
 119 };
 120
 121 static const uint8_t celt_log_freq_range[] = {
 122     0,  0,  0,  0,  0,  0,  0,  0,  8,  8,  8,  8, 16, 16, 16, 21, 21, 24, 29, 34, 36
 123 };
 124
 125 static const int8_t celt_tf_select[4][2][2][2] = {
 126     { { { 0, -1 }, { 0, -1 } }, { { 0, -1 }, { 0, -1 } } },
 127     { { { 0, -1 }, { 0, -2 } }, { { 1,  0 }, { 1, -1 } } },
 128     { { { 0, -2 }, { 0, -3 } }, { { 2,  0 }, { 1, -1 } } },
 129     { { { 0, -2 }, { 0, -3 } }, { { 3,  0 }, { 1, -1 } } }
 130 };
 131
 132 static const float celt_mean_energy[] = {
 133     6.437500f, 6.250000f, 5.750000f, 5.312500f, 5.062500f,
 134     4.812500f, 4.500000f, 4.375000f, 4.875000f, 4.687500f,
 135     4.562500f, 4.437500f, 4.875000f, 4.625000f, 4.312500f,
 136     4.500000f, 4.375000f, 4.625000f, 4.750000f, 4.437500f,
 137     3.750000f, 3.750000f, 3.750000f, 3.750000f, 3.750000f
 138 };
 139
 140 static const float celt_alpha_coef[] = {
 141     29440.0f/32768.0f,    26112.0f/32768.0f,    21248.0f/32768.0f,    16384.0f/32768.0f
 142 };
 143
 144 static const float celt_beta_coef[] = { /* TODO: precompute 1 minus this if the code ends up neater */
 145     30147.0f/32768.0f,    22282.0f/32768.0f,    12124.0f/32768.0f,     6554.0f/32768.0f
 146 };
 147
 148 static const uint8_t celt_coarse_energy_dist[4][2][42] = {
 149     {
 150         {       // 120-sample inter
 151              72, 127,  65, 129,  66, 128,  65, 128,  64, 128,  62, 128,  64, 128,
 152              64, 128,  92,  78,  92,  79,  92,  78,  90,  79, 116,  41, 115,  40,
 153             114,  40, 132,  26, 132,  26, 145,  17, 161,  12, 176,  10, 177,  11
 154         }, {    // 120-sample intra
 155              24, 179,  48, 138,  54, 135,  54, 132,  53, 134,  56, 133,  55, 132,
 156              55, 132,  61, 114,  70,  96,  74,  88,  75,  88,  87,  74,  89,  66,
 157              91,  67, 100,  59, 108,  50, 120,  40, 122,  37,  97,  43,  78,  50
 158         }
 159     }, {
 160         {       // 240-sample inter
 161              83,  78,  84,  81,  88,  75,  86,  74,  87,  71,  90,  73,  93,  74,
 162              93,  74, 109,  40, 114,  36, 117,  34, 117,  34, 143,  17, 145,  18,
 163             146,  19, 162,  12, 165,  10, 178,   7, 189,   6, 190,   8, 177,   9
 164         }, {    // 240-sample intra
 165              23, 178,  54, 115,  63, 102,  66,  98,  69,  99,  74,  89,  71,  91,
 166              73,  91,  78,  89,  86,  80,  92,  66,  93,  64, 102,  59, 103,  60,
 167             104,  60, 117,  52, 123,  44, 138,  35, 133,  31,  97,  38,  77,  45
 168         }
 169     }, {
 170         {       // 480-sample inter
 171              61,  90,  93,  60, 105,  42, 107,  41, 110,  45, 116,  38, 113,  38,
 172             112,  38, 124,  26, 132,  27, 136,  19, 140,  20, 155,  14, 159,  16,
 173             158,  18, 170,  13, 177,  10, 187,   8, 192,   6, 175,   9, 159,  10
 174         }, {    // 480-sample intra
 175              21, 178,  59, 110,  71,  86,  75,  85,  84,  83,  91,  66,  88,  73,
 176              87,  72,  92,  75,  98,  72, 105,  58, 107,  54, 115,  52, 114,  55,
 177             112,  56, 129,  51, 132,  40, 150,  33, 140,  29,  98,  35,  77,  42
 178         }
 179     }, {
 180         {       // 960-sample inter
 181              42, 121,  96,  66, 108,  43, 111,  40, 117,  44, 123,  32, 120,  36,
 182             119,  33, 127,  33, 134,  34, 139,  21, 147,  23, 152,  20, 158,  25,
 183             154,  26, 166,  21, 173,  16, 184,  13, 184,  10, 150,  13, 139,  15
 184         }, {    // 960-sample intra
 185              22, 178,  63, 114,  74,  82,  84,  83,  92,  82, 103,  62,  96,  72,
 186              96,  67, 101,  73, 107,  72, 113,  55, 118,  52, 125,  52, 118,  52,
 187             117,  55, 135,  49, 137,  39, 157,  32, 145,  29,  97,  33,  77,  40
 188         }
 189     }
 190 };
 191
 192 static const uint8_t celt_static_alloc[11][21] = {  /* 1/32 bit/sample */
 193     {   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0 },
 194     {  90,  80,  75,  69,  63,  56,  49,  40,  34,  29,  20,  18,  10,   0,   0,   0,   0,   0,   0,   0,   0 },
 195     { 110, 100,  90,  84,  78,  71,  65,  58,  51,  45,  39,  32,  26,  20,  12,   0,   0,   0,   0,   0,   0 },
 196     { 118, 110, 103,  93,  86,  80,  75,  70,  65,  59,  53,  47,  40,  31,  23,  15,   4,   0,   0,   0,   0 },
 197     { 126, 119, 112, 104,  95,  89,  83,  78,  72,  66,  60,  54,  47,  39,  32,  25,  17,  12,   1,   0,   0 },
 198     { 134, 127, 120, 114, 103,  97,  91,  85,  78,  72,  66,  60,  54,  47,  41,  35,  29,  23,  16,  10,   1 },
 199     { 144, 137, 130, 124, 113, 107, 101,  95,  88,  82,  76,  70,  64,  57,  51,  45,  39,  33,  26,  15,   1 },
 200     { 152, 145, 138, 132, 123, 117, 111, 105,  98,  92,  86,  80,  74,  67,  61,  55,  49,  43,  36,  20,   1 },
 201     { 162, 155, 148, 142, 133, 127, 121, 115, 108, 102,  96,  90,  84,  77,  71,  65,  59,  53,  46,  30,   1 },
 202     { 172, 165, 158, 152, 143, 137, 131, 125, 118, 112, 106, 100,  94,  87,  81,  75,  69,  63,  56,  45,  20 },
 203     { 200, 200, 200, 200, 200, 200, 200, 200, 198, 193, 188, 183, 178, 173, 168, 163, 158, 153, 148, 129, 104 }
 204 };
 205
 206 static const uint8_t celt_static_caps[4][2][21] = {
 207     {       // 120-sample
 208         {224, 224, 224, 224, 224, 224, 224, 224, 160, 160,
 209          160, 160, 185, 185, 185, 178, 178, 168, 134,  61,  37},
 210         {224, 224, 224, 224, 224, 224, 224, 224, 240, 240,
 211          240, 240, 207, 207, 207, 198, 198, 183, 144,  66,  40},
 212     }, {    // 240-sample
 213         {160, 160, 160, 160, 160, 160, 160, 160, 185, 185,
 214          185, 185, 193, 193, 193, 183, 183, 172, 138,  64,  38},
 215         {240, 240, 240, 240, 240, 240, 240, 240, 207, 207,
 216          207, 207, 204, 204, 204, 193, 193, 180, 143,  66,  40},
 217     }, {    // 480-sample
 218         {185, 185, 185, 185, 185, 185, 185, 185, 193, 193,
 219          193, 193, 193, 193, 193, 183, 183, 172, 138,  65,  39},
 220         {207, 207, 207, 207, 207, 207, 207, 207, 204, 204,
 221          204, 204, 201, 201, 201, 188, 188, 176, 141,  66,  40},
 222     }, {    // 960-sample
 223         {193, 193, 193, 193, 193, 193, 193, 193, 193, 193,
 224          193, 193, 194, 194, 194, 184, 184, 173, 139,  65,  39},
 225         {204, 204, 204, 204, 204, 204, 204, 204, 201, 201,
 226          201, 201, 198, 198, 198, 187, 187, 175, 140,  66,  40}
 227     }
 228 };
 229
 230 static const uint8_t celt_cache_bits[392] = {
 231     40, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7,
 232     7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7,
 233     7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 40, 15, 23, 28,
 234     31, 34, 36, 38, 39, 41, 42, 43, 44, 45, 46, 47, 47, 49, 50,
 235     51, 52, 53, 54, 55, 55, 57, 58, 59, 60, 61, 62, 63, 63, 65,
 236     66, 67, 68, 69, 70, 71, 71, 40, 20, 33, 41, 48, 53, 57, 61,
 237     64, 66, 69, 71, 73, 75, 76, 78, 80, 82, 85, 87, 89, 91, 92,
 238     94, 96, 98, 101, 103, 105, 107, 108, 110, 112, 114, 117, 119, 121, 123,
 239     124, 126, 128, 40, 23, 39, 51, 60, 67, 73, 79, 83, 87, 91, 94,
 240     97, 100, 102, 105, 107, 111, 115, 118, 121, 124, 126, 129, 131, 135, 139,
 241     142, 145, 148, 150, 153, 155, 159, 163, 166, 169, 172, 174, 177, 179, 35,
 242     28, 49, 65, 78, 89, 99, 107, 114, 120, 126, 132, 136, 141, 145, 149,
 243     153, 159, 165, 171, 176, 180, 185, 189, 192, 199, 205, 211, 216, 220, 225,
 244     229, 232, 239, 245, 251, 21, 33, 58, 79, 97, 112, 125, 137, 148, 157,
 245     166, 174, 182, 189, 195, 201, 207, 217, 227, 235, 243, 251, 17, 35, 63,
 246     86, 106, 123, 139, 152, 165, 177, 187, 197, 206, 214, 222, 230, 237, 250,
 247     25, 31, 55, 75, 91, 105, 117, 128, 138, 146, 154, 161, 168, 174, 180,
 248     185, 190, 200, 208, 215, 222, 229, 235, 240, 245, 255, 16, 36, 65, 89,
 249     110, 128, 144, 159, 173, 185, 196, 207, 217, 226, 234, 242, 250, 11, 41,
 250     74, 103, 128, 151, 172, 191, 209, 225, 241, 255, 9, 43, 79, 110, 138,
 251     163, 186, 207, 227, 246, 12, 39, 71, 99, 123, 144, 164, 182, 198, 214,
 252     228, 241, 253, 9, 44, 81, 113, 142, 168, 192, 214, 235, 255, 7, 49,
 253     90, 127, 160, 191, 220, 247, 6, 51, 95, 134, 170, 203, 234, 7, 47,
 254     87, 123, 155, 184, 212, 237, 6, 52, 97, 137, 174, 208, 240, 5, 57,
 255     106, 151, 192, 231, 5, 59, 111, 158, 202, 243, 5, 55, 103, 147, 187,
 256     224, 5, 60, 113, 161, 206, 248, 4, 65, 122, 175, 224, 4, 67, 127,
 257     182, 234
 258 };
 259
 260 static const int16_t celt_cache_index[105] = {
 261     -1, -1, -1, -1, -1, -1, -1, -1, 0, 0, 0, 0, 41, 41, 41,
 262     82, 82, 123, 164, 200, 222, 0, 0, 0, 0, 0, 0, 0, 0, 41,
 263     41, 41, 41, 123, 123, 123, 164, 164, 240, 266, 283, 295, 41, 41, 41,
 264     41, 41, 41, 41, 41, 123, 123, 123, 123, 240, 240, 240, 266, 266, 305,
 265     318, 328, 336, 123, 123, 123, 123, 123, 123, 123, 123, 240, 240, 240, 240,
 266     305, 305, 305, 318, 318, 343, 351, 358, 364, 240, 240, 240, 240, 240, 240,
 267     240, 240, 305, 305, 305, 305, 343, 343, 343, 351, 351, 370, 376, 382, 387,
 268 };
 269
 270 static const uint8_t celt_log2_frac[] = {
 271     0, 8, 13, 16, 19, 21, 23, 24, 26, 27, 28, 29, 30, 31, 32, 32, 33, 34, 34, 35, 36, 36, 37, 37
 272 };
 273
 274 static const uint8_t celt_bit_interleave[] = {
 275     0, 1, 1, 1, 2, 3, 3, 3, 2, 3, 3, 3, 2, 3, 3, 3
 276 };
 277
 278 static const uint8_t celt_bit_deinterleave[] = {
 279     0x00, 0x03, 0x0C, 0x0F, 0x30, 0x33, 0x3C, 0x3F,
 280     0xC0, 0xC3, 0xCC, 0xCF, 0xF0, 0xF3, 0xFC, 0xFF
 281 };
 282
 283 static const uint8_t celt_hadamard_ordery[] = {
 284     1,   0,
 285     3,   0,  2,  1,
 286     7,   0,  4,  3,  6,  1,  5,  2,
 287     15,  0,  8,  7, 12,  3, 11,  4, 14,  1,  9,  6, 13,  2, 10,  5
 288 };
 289
 290 static const uint16_t celt_qn_exp2[] = {
 291     16384, 17866, 19483, 21247, 23170, 25267, 27554, 30048
 292 };
 293
 294 static const uint32_t celt_pvq_u[1272] = {
 295     /* N = 0, K = 0...176 */
 296     1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
 297     0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
 298     0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
 299     0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
 300     0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
 301     0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
 302     0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
 303     /* N = 1, K = 1...176 */
 304     1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
 305     1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
 306     1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
 307     1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
 308     1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
 309     1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
 310     1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
 311     /* N = 2, K = 2...176 */
 312     3, 5, 7, 9, 11, 13, 15, 17, 19, 21, 23, 25, 27, 29, 31, 33, 35, 37, 39, 41,
 313     43, 45, 47, 49, 51, 53, 55, 57, 59, 61, 63, 65, 67, 69, 71, 73, 75, 77, 79,
 314     81, 83, 85, 87, 89, 91, 93, 95, 97, 99, 101, 103, 105, 107, 109, 111, 113,
 315     115, 117, 119, 121, 123, 125, 127, 129, 131, 133, 135, 137, 139, 141, 143,
 316     145, 147, 149, 151, 153, 155, 157, 159, 161, 163, 165, 167, 169, 171, 173,
 317     175, 177, 179, 181, 183, 185, 187, 189, 191, 193, 195, 197, 199, 201, 203,
 318     205, 207, 209, 211, 213, 215, 217, 219, 221, 223, 225, 227, 229, 231, 233,
 319     235, 237, 239, 241, 243, 245, 247, 249, 251, 253, 255, 257, 259, 261, 263,
 320     265, 267, 269, 271, 273, 275, 277, 279, 281, 283, 285, 287, 289, 291, 293,
 321     295, 297, 299, 301, 303, 305, 307, 309, 311, 313, 315, 317, 319, 321, 323,
 322     325, 327, 329, 331, 333, 335, 337, 339, 341, 343, 345, 347, 349, 351,
 323     /* N = 3, K = 3...176 */
 324     13, 25, 41, 61, 85, 113, 145, 181, 221, 265, 313, 365, 421, 481, 545, 613,
 325     685, 761, 841, 925, 1013, 1105, 1201, 1301, 1405, 1513, 1625, 1741, 1861,
 326     1985, 2113, 2245, 2381, 2521, 2665, 2813, 2965, 3121, 3281, 3445, 3613, 3785,
 327     3961, 4141, 4325, 4513, 4705, 4901, 5101, 5305, 5513, 5725, 5941, 6161, 6385,
 328     6613, 6845, 7081, 7321, 7565, 7813, 8065, 8321, 8581, 8845, 9113, 9385, 9661,
 329     9941, 10225, 10513, 10805, 11101, 11401, 11705, 12013, 12325, 12641, 12961,
 330     13285, 13613, 13945, 14281, 14621, 14965, 15313, 15665, 16021, 16381, 16745,
 331     17113, 17485, 17861, 18241, 18625, 19013, 19405, 19801, 20201, 20605, 21013,
 332     21425, 21841, 22261, 22685, 23113, 23545, 23981, 24421, 24865, 25313, 25765,
 333     26221, 26681, 27145, 27613, 28085, 28561, 29041, 29525, 30013, 30505, 31001,
 334     31501, 32005, 32513, 33025, 33541, 34061, 34585, 35113, 35645, 36181, 36721,
 335     37265, 37813, 38365, 38921, 39481, 40045, 40613, 41185, 41761, 42341, 42925,
 336     43513, 44105, 44701, 45301, 45905, 46513, 47125, 47741, 48361, 48985, 49613,
 337     50245, 50881, 51521, 52165, 52813, 53465, 54121, 54781, 55445, 56113, 56785,
 338     57461, 58141, 58825, 59513, 60205, 60901, 61601,
 339     /* N = 4, K = 4...176 */
 340     63, 129, 231, 377, 575, 833, 1159, 1561, 2047, 2625, 3303, 4089, 4991, 6017,
 341     7175, 8473, 9919, 11521, 13287, 15225, 17343, 19649, 22151, 24857, 27775,
 342     30913, 34279, 37881, 41727, 45825, 50183, 54809, 59711, 64897, 70375, 76153,
 343     82239, 88641, 95367, 102425, 109823, 117569, 125671, 134137, 142975, 152193,
 344     161799, 171801, 182207, 193025, 204263, 215929, 228031, 240577, 253575,
 345     267033, 280959, 295361, 310247, 325625, 341503, 357889, 374791, 392217,
 346     410175, 428673, 447719, 467321, 487487, 508225, 529543, 551449, 573951,
 347     597057, 620775, 645113, 670079, 695681, 721927, 748825, 776383, 804609,
 348     833511, 863097, 893375, 924353, 956039, 988441, 1021567, 1055425, 1090023,
 349     1125369, 1161471, 1198337, 1235975, 1274393, 1313599, 1353601, 1394407,
 350     1436025, 1478463, 1521729, 1565831, 1610777, 1656575, 1703233, 1750759,
 351     1799161, 1848447, 1898625, 1949703, 2001689, 2054591, 2108417, 2163175,
 352     2218873, 2275519, 2333121, 2391687, 2451225, 2511743, 2573249, 2635751,
 353     2699257, 2763775, 2829313, 2895879, 2963481, 3032127, 3101825, 3172583,
 354     3244409, 3317311, 3391297, 3466375, 3542553, 3619839, 3698241, 3777767,
 355     3858425, 3940223, 4023169, 4107271, 4192537, 4278975, 4366593, 4455399,
 356     4545401, 4636607, 4729025, 4822663, 4917529, 5013631, 5110977, 5209575,
 357     5309433, 5410559, 5512961, 5616647, 5721625, 5827903, 5935489, 6044391,
 358     6154617, 6266175, 6379073, 6493319, 6608921, 6725887, 6844225, 6963943,
 359     7085049, 7207551,
 360     /* N = 5, K = 5...176 */
 361     321, 681, 1289, 2241, 3649, 5641, 8361, 11969, 16641, 22569, 29961, 39041,
 362     50049, 63241, 78889, 97281, 118721, 143529, 172041, 204609, 241601, 283401,
 363     330409, 383041, 441729, 506921, 579081, 658689, 746241, 842249, 947241,
 364     1061761, 1186369, 1321641, 1468169, 1626561, 1797441, 1981449, 2179241,
 365     2391489, 2618881, 2862121, 3121929, 3399041, 3694209, 4008201, 4341801,
 366     4695809, 5071041, 5468329, 5888521, 6332481, 6801089, 7295241, 7815849,
 367     8363841, 8940161, 9545769, 10181641, 10848769, 11548161, 12280841, 13047849,
 368     13850241, 14689089, 15565481, 16480521, 17435329, 18431041, 19468809,
 369     20549801, 21675201, 22846209, 24064041, 25329929, 26645121, 28010881,
 370     29428489, 30899241, 32424449, 34005441, 35643561, 37340169, 39096641,
 371     40914369, 42794761, 44739241, 46749249, 48826241, 50971689, 53187081,
 372     55473921, 57833729, 60268041, 62778409, 65366401, 68033601, 70781609,
 373     73612041, 76526529, 79526721, 82614281, 85790889, 89058241, 92418049,
 374     95872041, 99421961, 103069569, 106816641, 110664969, 114616361, 118672641,
 375     122835649, 127107241, 131489289, 135983681, 140592321, 145317129, 150160041,
 376     155123009, 160208001, 165417001, 170752009, 176215041, 181808129, 187533321,
 377     193392681, 199388289, 205522241, 211796649, 218213641, 224775361, 231483969,
 378     238341641, 245350569, 252512961, 259831041, 267307049, 274943241, 282741889,
 379     290705281, 298835721, 307135529, 315607041, 324252609, 333074601, 342075401,
 380     351257409, 360623041, 370174729, 379914921, 389846081, 399970689, 410291241,
 381     420810249, 431530241, 442453761, 453583369, 464921641, 476471169, 488234561,
 382     500214441, 512413449, 524834241, 537479489, 550351881, 563454121, 576788929,
 383     590359041, 604167209, 618216201, 632508801,
 384     /* N = 6, K = 6...96 (technically V(109,5) fits in 32 bits, but that can't be
 385      achieved by splitting an Opus band) */
 386     1683, 3653, 7183, 13073, 22363, 36365, 56695, 85305, 124515, 177045, 246047,
 387     335137, 448427, 590557, 766727, 982729, 1244979, 1560549, 1937199, 2383409,
 388     2908411, 3522221, 4235671, 5060441, 6009091, 7095093, 8332863, 9737793,
 389     11326283, 13115773, 15124775, 17372905, 19880915, 22670725, 25765455,
 390     29189457, 32968347, 37129037, 41699767, 46710137, 52191139, 58175189,
 391     64696159, 71789409, 79491819, 87841821, 96879431, 106646281, 117185651,
 392     128542501, 140763503, 153897073, 167993403, 183104493, 199284183, 216588185,
 393     235074115, 254801525, 275831935, 298228865, 322057867, 347386557, 374284647,
 394     402823977, 433078547, 465124549, 499040399, 534906769, 572806619, 612825229,
 395     655050231, 699571641, 746481891, 795875861, 847850911, 902506913, 959946283,
 396     1020274013, 1083597703, 1150027593, 1219676595, 1292660325, 1369097135,
 397     1449108145, 1532817275, 1620351277, 1711839767, 1807415257, 1907213187,
 398     2011371957, 2120032959,
 399     /* N = 7, K = 7...54 (technically V(60,6) fits in 32 bits, but that can't be
 400      achieved by splitting an Opus band) */
 401     8989, 19825, 40081, 75517, 134245, 227305, 369305, 579125, 880685, 1303777,
 402     1884961, 2668525, 3707509, 5064793, 6814249, 9041957, 11847485, 15345233,
 403     19665841, 24957661, 31388293, 39146185, 48442297, 59511829, 72616013,
 404     88043969, 106114625, 127178701, 151620757, 179861305, 212358985, 249612805,
 405     292164445, 340600625, 395555537, 457713341, 527810725, 606639529, 695049433,
 406     793950709, 904317037, 1027188385, 1163673953, 1314955181, 1482288821,
 407     1667010073, 1870535785, 2094367717,
 408     /* N = 8, K = 8...37 (technically V(40,7) fits in 32 bits, but that can't be
 409      achieved by splitting an Opus band) */
 410     48639, 108545, 224143, 433905, 795455, 1392065, 2340495, 3800305, 5984767,
 411     9173505, 13726991, 20103025, 28875327, 40754369, 56610575, 77500017,
 412     104692735, 139703809, 184327311, 240673265, 311207743, 398796225, 506750351,
 413     638878193, 799538175, 993696769, 1226990095, 1505789553, 1837271615,
 414     2229491905,
 415     /* N = 9, K = 9...28 (technically V(29,8) fits in 32 bits, but that can't be
 416      achieved by splitting an Opus band) */
 417     265729, 598417, 1256465, 2485825, 4673345, 8405905, 14546705, 24331777,
 418     39490049, 62390545, 96220561, 145198913, 214828609, 312193553, 446304145,
 419     628496897, 872893441, 1196924561, 1621925137, 2173806145,
 420     /* N = 10, K = 10...24 */
 421     1462563, 3317445, 7059735, 14218905, 27298155, 50250765, 89129247, 152951073,
 422     254831667, 413442773, 654862247, 1014889769, 1541911931, 2300409629,
 423     3375210671,
 424     /* N = 11, K = 11...19 (technically V(20,10) fits in 32 bits, but that can't be
 425      achieved by splitting an Opus band) */
 426     8097453, 18474633, 39753273, 81270333, 158819253, 298199265, 540279585,
 427     948062325, 1616336765,
 428     /* N = 12, K = 12...18 */
 429     45046719, 103274625, 224298231, 464387817, 921406335, 1759885185,
 430     3248227095,
 431     /* N = 13, K = 13...16 */
 432     251595969, 579168825, 1267854873, 2653649025,
 433     /* N = 14, K = 14 */
 434     1409933619
 435 };
 436
 437 DECLARE_ALIGNED(32, static const float, celt_window)[120] = {
 438     6.7286966e-05f, 0.00060551348f, 0.0016815970f, 0.0032947962f, 0.0054439943f,
 439     0.0081276923f, 0.011344001f, 0.015090633f, 0.019364886f, 0.024163635f,
 440     0.029483315f, 0.035319905f, 0.041668911f, 0.048525347f, 0.055883718f,
 441     0.063737999f, 0.072081616f, 0.080907428f, 0.090207705f, 0.099974111f,
 442     0.11019769f, 0.12086883f, 0.13197729f, 0.14351214f, 0.15546177f,
 443     0.16781389f, 0.18055550f, 0.19367290f, 0.20715171f, 0.22097682f,
 444     0.23513243f, 0.24960208f, 0.26436860f, 0.27941419f, 0.29472040f,
 445     0.31026818f, 0.32603788f, 0.34200931f, 0.35816177f, 0.37447407f,
 446     0.39092462f, 0.40749142f, 0.42415215f, 0.44088423f, 0.45766484f,
 447     0.47447104f, 0.49127978f, 0.50806798f, 0.52481261f, 0.54149077f,
 448     0.55807973f, 0.57455701f, 0.59090049f, 0.60708841f, 0.62309951f,
 449     0.63891306f, 0.65450896f, 0.66986776f, 0.68497077f, 0.69980010f,
 450     0.71433873f, 0.72857055f, 0.74248043f, 0.75605424f, 0.76927895f,
 451     0.78214257f, 0.79463430f, 0.80674445f, 0.81846456f, 0.82978733f,
 452     0.84070669f, 0.85121779f, 0.86131698f, 0.87100183f, 0.88027111f,
 453     0.88912479f, 0.89756398f, 0.90559094f, 0.91320904f, 0.92042270f,
 454     0.92723738f, 0.93365955f, 0.93969656f, 0.94535671f, 0.95064907f,
 455     0.95558353f, 0.96017067f, 0.96442171f, 0.96834849f, 0.97196334f,
 456     0.97527906f, 0.97830883f, 0.98106616f, 0.98356480f, 0.98581869f,
 457     0.98784191f, 0.98964856f, 0.99125274f, 0.99266849f, 0.99390969f,
 458     0.99499004f, 0.99592297f, 0.99672162f, 0.99739874f, 0.99796667f,
 459     0.99843728f, 0.99882195f, 0.99913147f, 0.99937606f, 0.99956527f,
 460     0.99970802f, 0.99981248f, 0.99988613f, 0.99993565f, 0.99996697f,
 461     0.99998518f, 0.99999457f, 0.99999859f, 0.99999982f, 1.0000000f,
 462 };
 463
 464 /* square of the window, used for the postfilter */
 465 const float ff_celt_window2[120] = {
 466     4.5275357e-09f, 3.66647e-07f, 2.82777e-06f, 1.08557e-05f, 2.96371e-05f, 6.60594e-05f,
 467     0.000128686f, 0.000227727f, 0.000374999f, 0.000583881f, 0.000869266f, 0.0012475f,
 468     0.0017363f, 0.00235471f, 0.00312299f, 0.00406253f, 0.00519576f, 0.00654601f,
 469     0.00813743f, 0.00999482f, 0.0121435f, 0.0146093f, 0.017418f, 0.0205957f, 0.0241684f,
 470     0.0281615f, 0.0326003f, 0.0375092f, 0.0429118f, 0.0488308f, 0.0552873f, 0.0623012f,
 471     0.0698908f, 0.0780723f, 0.0868601f, 0.0962664f, 0.106301f, 0.11697f, 0.12828f,
 472     0.140231f, 0.152822f, 0.166049f, 0.179905f, 0.194379f, 0.209457f, 0.225123f, 0.241356f,
 473     0.258133f, 0.275428f, 0.293212f, 0.311453f, 0.330116f, 0.349163f, 0.368556f, 0.388253f,
 474     0.40821f, 0.428382f, 0.448723f, 0.469185f, 0.48972f, 0.51028f, 0.530815f, 0.551277f,
 475     0.571618f, 0.59179f, 0.611747f, 0.631444f, 0.650837f, 0.669884f, 0.688547f, 0.706788f,
 476     0.724572f, 0.741867f, 0.758644f, 0.774877f, 0.790543f, 0.805621f, 0.820095f, 0.833951f,
 477     0.847178f, 0.859769f, 0.87172f, 0.88303f, 0.893699f, 0.903734f, 0.91314f, 0.921928f,
 478     0.930109f, 0.937699f, 0.944713f, 0.951169f, 0.957088f, 0.962491f, 0.9674f, 0.971838f,
 479     0.975832f, 0.979404f, 0.982582f, 0.985391f, 0.987857f, 0.990005f, 0.991863f, 0.993454f,
 480     0.994804f, 0.995937f, 0.996877f, 0.997645f, 0.998264f, 0.998753f, 0.999131f, 0.999416f,
 481     0.999625f, 0.999772f, 0.999871f, 0.999934f, 0.99997f, 0.999989f, 0.999997f, 0.99999964f, 1.0f,
 482 };
 483
 484 static const uint32_t * const celt_pvq_u_row[15] = {
 485     celt_pvq_u +    0, celt_pvq_u +  176, celt_pvq_u +  351,
 486     celt_pvq_u +  525, celt_pvq_u +  698, celt_pvq_u +  870,
 487     celt_pvq_u + 1041, celt_pvq_u + 1131, celt_pvq_u + 1178,
 488     celt_pvq_u + 1207, celt_pvq_u + 1226, celt_pvq_u + 1240,
 489     celt_pvq_u + 1248, celt_pvq_u + 1254, celt_pvq_u + 1257
 490 };
 491
 492 static inline int16_t celt_cos(int16_t x)
 493 {
 494     x = (MUL16(x, x) + 4096) >> 13;
 495     x = (32767-x) + ROUND_MUL16(x, (-7651 + ROUND_MUL16(x, (8277 + ROUND_MUL16(-626, x)))));
 496     return 1+x;
 497 }
 498
 499 static inline int celt_log2tan(int isin, int icos)
 500 {
 501     int lc, ls;
 502     lc = opus_ilog(icos);
 503     ls = opus_ilog(isin);
 504     icos <<= 15 - lc;
 505     isin <<= 15 - ls;
 506     return (ls << 11) - (lc << 11) +
 507            ROUND_MUL16(isin, ROUND_MUL16(isin, -2597) + 7932) -
 508            ROUND_MUL16(icos, ROUND_MUL16(icos, -2597) + 7932);
 509 }
 510
 511 static inline uint32_t celt_rng(CeltContext *s)
 512 {
 513     s->seed = 1664525 * s->seed + 1013904223;
 514     return s->seed;
 515 }
 516
 517 static void celt_decode_coarse_energy(CeltContext *s, OpusRangeCoder *rc)
 518 {
 519     int i, j;
 520     float prev[2] = {0};
 521     float alpha, beta;
 522     const uint8_t *model;
 523
 524     /* use the 2D z-transform to apply prediction in both */
 525     /* the time domain (alpha) and the frequency domain (beta) */
 526
 527     if (opus_rc_tell(rc)+3 <= s->framebits && opus_rc_p2model(rc, 3)) {
 528         /* intra frame */
 529         alpha = 0;
 530         beta  = 1.0f - 4915.0f/32768.0f;
 531         model = celt_coarse_energy_dist[s->duration][1];
 532     } else {
 533         alpha = celt_alpha_coef[s->duration];
 534         beta  = 1.0f - celt_beta_coef[s->duration];
 535         model = celt_coarse_energy_dist[s->duration][0];
 536     }
 537
 538     for (i = 0; i < CELT_MAX_BANDS; i++) {
 539         for (j = 0; j < s->coded_channels; j++) {
 540             CeltFrame *frame = &s->frame[j];
 541             float value;
 542             int available;
 543
 544             if (i < s->startband || i >= s->endband) {
 545                 frame->energy[i] = 0.0;
 546                 continue;
 547             }
 548
 549             available = s->framebits - opus_rc_tell(rc);
 550             if (available >= 15) {
 551                 /* decode using a Laplace distribution */
 552                 int k = FFMIN(i, 20) << 1;
 553                 value = opus_rc_laplace(rc, model[k] << 7, model[k+1] << 6);
 554             } else if (available >= 2) {
 555                 int x = opus_rc_getsymbol(rc, celt_model_energy_small);
 556                 value = (x>>1) ^ -(x&1);
 557             } else if (available >= 1) {
 558                 value = -(float)opus_rc_p2model(rc, 1);
 559             } else value = -1;
 560
 561             frame->energy[i] = FFMAX(-9.0f, frame->energy[i]) * alpha + prev[j] + value;
 562             prev[j] += beta * value;
 563         }
 564     }
 565 }
 566
 567 static void celt_decode_fine_energy(CeltContext *s, OpusRangeCoder *rc)
 568 {
 569     int i;
 570     for (i = s->startband; i < s->endband; i++) {
 571         int j;
 572         if (!s->fine_bits[i])
 573             continue;
 574
 575         for (j = 0; j < s->coded_channels; j++) {
 576             CeltFrame *frame = &s->frame[j];
 577             int q2;
 578             float offset;
 579             q2 = opus_getrawbits(rc, s->fine_bits[i]);
 580             offset = (q2 + 0.5f) * (1 << (14 - s->fine_bits[i])) / 16384.0f - 0.5f;
 581             frame->energy[i] += offset;
 582         }
 583     }
 584 }
 585
 586 static void celt_decode_final_energy(CeltContext *s, OpusRangeCoder *rc,
 587                                      int bits_left)
 588 {
 589     int priority, i, j;
 590
 591     for (priority = 0; priority < 2; priority++) {
 592         for (i = s->startband; i < s->endband && bits_left >= s->coded_channels; i++) {
 593             if (s->fine_priority[i] != priority || s->fine_bits[i] >= CELT_MAX_FINE_BITS)
 594                 continue;
 595
 596             for (j = 0; j < s->coded_channels; j++) {
 597                 int q2;
 598                 float offset;
 599                 q2 = opus_getrawbits(rc, 1);
 600                 offset = (q2 - 0.5f) * (1 << (14 - s->fine_bits[i] - 1)) / 16384.0f;
 601                 s->frame[j].energy[i] += offset;
 602                 bits_left--;
 603             }
 604         }
 605     }
 606 }
 607
 608 static void celt_decode_tf_changes(CeltContext *s, OpusRangeCoder *rc,
 609                                    int transient)
 610 {
 611     int i, diff = 0, tf_select = 0, tf_changed = 0, tf_select_bit;
 612     int consumed, bits = transient ? 2 : 4;
 613
 614     consumed = opus_rc_tell(rc);
 615     tf_select_bit = (s->duration != 0 && consumed+bits+1 <= s->framebits);
 616
 617     for (i = s->startband; i < s->endband; i++) {
 618         if (consumed+bits+tf_select_bit <= s->framebits) {
 619             diff ^= opus_rc_p2model(rc, bits);
 620             consumed = opus_rc_tell(rc);
 621             tf_changed |= diff;
 622         }
 623         s->tf_change[i] = diff;
 624         bits = transient ? 4 : 5;
 625     }
 626
 627     if (tf_select_bit && celt_tf_select[s->duration][transient][0][tf_changed] !=
 628                          celt_tf_select[s->duration][transient][1][tf_changed])
 629         tf_select = opus_rc_p2model(rc, 1);
 630
 631     for (i = s->startband; i < s->endband; i++) {
 632         s->tf_change[i] = celt_tf_select[s->duration][transient][tf_select][s->tf_change[i]];
 633     }
 634 }
 635
 636 static void celt_decode_allocation(CeltContext *s, OpusRangeCoder *rc)
 637 {
 638     // approx. maximum bit allocation for each band before boost/trim
 639     int cap[CELT_MAX_BANDS];
 640     int boost[CELT_MAX_BANDS];
 641     int threshold[CELT_MAX_BANDS];
 642     int bits1[CELT_MAX_BANDS];
 643     int bits2[CELT_MAX_BANDS];
 644     int trim_offset[CELT_MAX_BANDS];
 645
 646     int skip_startband = s->startband;
 647     int dynalloc       = 6;
 648     int alloctrim      = 5;
 649     int extrabits      = 0;
 650
 651     int skip_bit            = 0;
 652     int intensitystereo_bit = 0;
 653     int dualstereo_bit      = 0;
 654
 655     int remaining, bandbits;
 656     int low, high, total, done;
 657     int totalbits;
 658     int consumed;
 659     int i, j;
 660
 661     consumed = opus_rc_tell(rc);
 662
 663     /* obtain spread flag */
 664     s->spread = CELT_SPREAD_NORMAL;
 665     if (consumed + 4 <= s->framebits)
 666         s->spread = opus_rc_getsymbol(rc, celt_model_spread);
 667
 668     /* generate static allocation caps */
 669     for (i = 0; i < CELT_MAX_BANDS; i++) {
 670         cap[i] = (celt_static_caps[s->duration][s->coded_channels - 1][i] + 64)
 671                  * celt_freq_range[i] << (s->coded_channels - 1) << s->duration >> 2;
 672     }
 673
 674     /* obtain band boost */
 675     totalbits = s->framebits << 3; // convert to 1/8 bits
 676     consumed = opus_rc_tell_frac(rc);
 677     for (i = s->startband; i < s->endband; i++) {
 678         int quanta, band_dynalloc;
 679
 680         boost[i] = 0;
 681
 682         quanta = celt_freq_range[i] << (s->coded_channels - 1) << s->duration;
 683         quanta = FFMIN(quanta << 3, FFMAX(6 << 3, quanta));
 684         band_dynalloc = dynalloc;
 685         while (consumed + (band_dynalloc<<3) < totalbits && boost[i] < cap[i]) {
 686             int add = opus_rc_p2model(rc, band_dynalloc);
 687             consumed = opus_rc_tell_frac(rc);
 688             if (!add)
 689                 break;
 690
 691             boost[i]     += quanta;
 692             totalbits    -= quanta;
 693             band_dynalloc = 1;
 694         }
 695         /* dynalloc is more likely to occur if it's already been used for earlier bands */
 696         if (boost[i])
 697             dynalloc = FFMAX(2, dynalloc - 1);
 698     }
 699
 700     /* obtain allocation trim */
 701     if (consumed + (6 << 3) <= totalbits)
 702         alloctrim = opus_rc_getsymbol(rc, celt_model_alloc_trim);
 703
 704     /* anti-collapse bit reservation */
 705     totalbits = (s->framebits << 3) - opus_rc_tell_frac(rc) - 1;
 706     s->anticollapse_bit = 0;
 707     if (s->blocks > 1 && s->duration >= 2 &&
 708         totalbits >= ((s->duration + 2) << 3))
 709         s->anticollapse_bit = 1 << 3;
 710     totalbits -= s->anticollapse_bit;
 711
 712     /* band skip bit reservation */
 713     if (totalbits >= 1 << 3)
 714         skip_bit = 1 << 3;
 715     totalbits -= skip_bit;
 716
 717     /* intensity/dual stereo bit reservation */
 718     if (s->coded_channels == 2) {
 719         intensitystereo_bit = celt_log2_frac[s->endband - s->startband];
 720         if (intensitystereo_bit <= totalbits) {
 721             totalbits -= intensitystereo_bit;
 722             if (totalbits >= 1 << 3) {
 723                 dualstereo_bit = 1 << 3;
 724                 totalbits -= 1 << 3;
 725             }
 726         } else
 727             intensitystereo_bit = 0;
 728     }
 729
 730     for (i = s->startband; i < s->endband; i++) {
 731         int trim     = alloctrim - 5 - s->duration;
 732         int band     = celt_freq_range[i] * (s->endband - i - 1);
 733         int duration = s->duration + 3;
 734         int scale    = duration + s->coded_channels - 1;
 735
 736         /* PVQ minimum allocation threshold, below this value the band is
 737          * skipped */
 738         threshold[i] = FFMAX(3 * celt_freq_range[i] << duration >> 4,
 739                              s->coded_channels << 3);
 740
 741         trim_offset[i] = trim * (band << scale) >> 6;
 742
 743         if (celt_freq_range[i] << s->duration == 1)
 744             trim_offset[i] -= s->coded_channels << 3;
 745     }
 746
 747     /* bisection */
 748     low  = 1;
 749     high = CELT_VECTORS - 1;
 750     while (low <= high) {
 751         int center = (low + high) >> 1;
 752         done = total = 0;
 753
 754         for (i = s->endband - 1; i >= s->startband; i--) {
 755             bandbits = celt_freq_range[i] * celt_static_alloc[center][i]
 756                        << (s->coded_channels - 1) << s->duration >> 2;
 757
 758             if (bandbits)
 759                 bandbits = FFMAX(0, bandbits + trim_offset[i]);
 760             bandbits += boost[i];
 761
 762             if (bandbits >= threshold[i] || done) {
 763                 done = 1;
 764                 total += FFMIN(bandbits, cap[i]);
 765             } else if (bandbits >= s->coded_channels << 3)
 766                 total += s->coded_channels << 3;
 767         }
 768
 769         if (total > totalbits)
 770             high = center - 1;
 771         else
 772             low = center + 1;
 773     }
 774     high = low--;
 775
 776     for (i = s->startband; i < s->endband; i++) {
 777         bits1[i] = celt_freq_range[i] * celt_static_alloc[low][i]
 778                    << (s->coded_channels - 1) << s->duration >> 2;
 779         bits2[i] = high >= CELT_VECTORS ? cap[i] :
 780                    celt_freq_range[i] * celt_static_alloc[high][i]
 781                    << (s->coded_channels - 1) << s->duration >> 2;
 782
 783         if (bits1[i])
 784             bits1[i] = FFMAX(0, bits1[i] + trim_offset[i]);
 785         if (bits2[i])
 786             bits2[i] = FFMAX(0, bits2[i] + trim_offset[i]);
 787         if (low)
 788             bits1[i] += boost[i];
 789         bits2[i] += boost[i];
 790
 791         if (boost[i])
 792             skip_startband = i;
 793         bits2[i] = FFMAX(0, bits2[i] - bits1[i]);
 794     }
 795
 796     /* bisection */
 797     low  = 0;
 798     high = 1 << CELT_ALLOC_STEPS;
 799     for (i = 0; i < CELT_ALLOC_STEPS; i++) {
 800         int center = (low + high) >> 1;
 801         done = total = 0;
 802
 803         for (j = s->endband - 1; j >= s->startband; j--) {
 804             bandbits = bits1[j] + (center * bits2[j] >> CELT_ALLOC_STEPS);
 805
 806             if (bandbits >= threshold[j] || done) {
 807                 done = 1;
 808                 total += FFMIN(bandbits, cap[j]);
 809             } else if (bandbits >= s->coded_channels << 3)
 810                 total += s->coded_channels << 3;
 811         }
 812         if (total > totalbits)
 813             high = center;
 814         else
 815             low = center;
 816     }
 817
 818     done = total = 0;
 819     for (i = s->endband - 1; i >= s->startband; i--) {
 820         bandbits = bits1[i] + (low * bits2[i] >> CELT_ALLOC_STEPS);
 821
 822         if (bandbits >= threshold[i] || done)
 823             done = 1;
 824         else
 825             bandbits = (bandbits >= s->coded_channels << 3) ?
 826                        s->coded_channels << 3 : 0;
 827
 828         bandbits     = FFMIN(bandbits, cap[i]);
 829         s->pulses[i] = bandbits;
 830         total      += bandbits;
 831     }
 832
 833     /* band skipping */
 834     for (s->codedbands = s->endband; ; s->codedbands--) {
 835         int allocation;
 836         j = s->codedbands - 1;
 837
 838         if (j == skip_startband) {
 839             /* all remaining bands are not skipped */
 840             totalbits += skip_bit;
 841             break;
 842         }
 843
 844         /* determine the number of bits available for coding "do not skip" markers */
 845         remaining   = totalbits - total;
 846         bandbits    = remaining / (celt_freq_bands[j+1] - celt_freq_bands[s->startband]);
 847         remaining  -= bandbits  * (celt_freq_bands[j+1] - celt_freq_bands[s->startband]);
 848         allocation  = s->pulses[j] + bandbits * celt_freq_range[j]
 849                       + FFMAX(0, remaining - (celt_freq_bands[j] - celt_freq_bands[s->startband]));
 850
 851         /* a "do not skip" marker is only coded if the allocation is
 852            above the chosen threshold */
 853         if (allocation >= FFMAX(threshold[j], (s->coded_channels + 1) <<3 )) {
 854             if (opus_rc_p2model(rc, 1))
 855                 break;
 856
 857             total      += 1 << 3;
 858             allocation -= 1 << 3;
 859         }
 860
 861         /* the band is skipped, so reclaim its bits */
 862         total -= s->pulses[j];
 863         if (intensitystereo_bit) {
 864             total -= intensitystereo_bit;
 865             intensitystereo_bit = celt_log2_frac[j - s->startband];
 866             total += intensitystereo_bit;
 867         }
 868
 869         total += s->pulses[j] = (allocation >= s->coded_channels << 3) ?
 870                               s->coded_channels << 3 : 0;
 871     }
 872
 873     /* obtain stereo flags */
 874     s->intensitystereo = 0;
 875     s->dualstereo      = 0;
 876     if (intensitystereo_bit)
 877         s->intensitystereo = s->startband +
 878                           opus_rc_unimodel(rc, s->codedbands + 1 - s->startband);
 879     if (s->intensitystereo <= s->startband)
 880         totalbits += dualstereo_bit; /* no intensity stereo means no dual stereo */
 881     else if (dualstereo_bit)
 882         s->dualstereo = opus_rc_p2model(rc, 1);
 883
 884     /* supply the remaining bits in this frame to lower bands */
 885     remaining = totalbits - total;
 886     bandbits  = remaining / (celt_freq_bands[s->codedbands] - celt_freq_bands[s->startband]);
 887     remaining -= bandbits * (celt_freq_bands[s->codedbands] - celt_freq_bands[s->startband]);
 888     for (i = s->startband; i < s->codedbands; i++) {
 889         int bits = FFMIN(remaining, celt_freq_range[i]);
 890
 891         s->pulses[i] += bits + bandbits * celt_freq_range[i];
 892         remaining    -= bits;
 893     }
 894
 895     for (i = s->startband; i < s->codedbands; i++) {
 896         int N = celt_freq_range[i] << s->duration;
 897         int prev_extra = extrabits;
 898         s->pulses[i] += extrabits;
 899
 900         if (N > 1) {
 901             int dof;        // degrees of freedom
 902             int temp;       // dof * channels * log(dof)
 903             int offset;     // fine energy quantization offset, i.e.
 904                             // extra bits assigned over the standard
 905                             // totalbits/dof
 906             int fine_bits, max_bits;
 907
 908             extrabits = FFMAX(0, s->pulses[i] - cap[i]);
 909             s->pulses[i] -= extrabits;
 910
 911             /* intensity stereo makes use of an extra degree of freedom */
 912             dof = N * s->coded_channels
 913                   + (s->coded_channels == 2 && N > 2 && !s->dualstereo && i < s->intensitystereo);
 914             temp = dof * (celt_log_freq_range[i] + (s->duration<<3));
 915             offset = (temp >> 1) - dof * CELT_FINE_OFFSET;
 916             if (N == 2) /* dof=2 is the only case that doesn't fit the model */
 917                 offset += dof<<1;
 918
 919             /* grant an additional bias for the first and second pulses */
 920             if (s->pulses[i] + offset < 2 * (dof << 3))
 921                 offset += temp >> 2;
 922             else if (s->pulses[i] + offset < 3 * (dof << 3))
 923                 offset += temp >> 3;
 924
 925             fine_bits = (s->pulses[i] + offset + (dof << 2)) / (dof << 3);
 926             max_bits  = FFMIN((s->pulses[i]>>3) >> (s->coded_channels - 1),
 927                               CELT_MAX_FINE_BITS);
 928
 929             max_bits  = FFMAX(max_bits, 0);
 930
 931             s->fine_bits[i] = av_clip(fine_bits, 0, max_bits);
 932
 933             /* if fine_bits was rounded down or capped,
 934                give priority for the final fine energy pass */
 935             s->fine_priority[i] = (s->fine_bits[i] * (dof<<3) >= s->pulses[i] + offset);
 936
 937             /* the remaining bits are assigned to PVQ */
 938             s->pulses[i] -= s->fine_bits[i] << (s->coded_channels - 1) << 3;
 939         } else {
 940             /* all bits go to fine energy except for the sign bit */
 941             extrabits = FFMAX(0, s->pulses[i] - (s->coded_channels << 3));
 942             s->pulses[i] -= extrabits;
 943             s->fine_bits[i] = 0;
 944             s->fine_priority[i] = 1;
 945         }
 946
 947         /* hand back a limited number of extra fine energy bits to this band */
 948         if (extrabits > 0) {
 949             int fineextra = FFMIN(extrabits >> (s->coded_channels + 2),
 950                                   CELT_MAX_FINE_BITS - s->fine_bits[i]);
 951             s->fine_bits[i] += fineextra;
 952
 953             fineextra <<= s->coded_channels + 2;
 954             s->fine_priority[i] = (fineextra >= extrabits - prev_extra);
 955             extrabits -= fineextra;
 956         }
 957     }
 958     s->remaining = extrabits;
 959
 960     /* skipped bands dedicate all of their bits for fine energy */
 961     for (; i < s->endband; i++) {
 962         s->fine_bits[i]     = s->pulses[i] >> (s->coded_channels - 1) >> 3;
 963         s->pulses[i]        = 0;
 964         s->fine_priority[i] = s->fine_bits[i] < 1;
 965     }
 966 }
 967
 968 static inline int celt_bits2pulses(const uint8_t *cache, int bits)
 969 {
 970     // TODO: Find the size of cache and make it into an array in the parameters list
 971     int i, low = 0, high;
 972
 973     high = cache[0];
 974     bits--;
 975
 976     for (i = 0; i < 6; i++) {
 977         int center = (low + high + 1) >> 1;
 978         if (cache[center] >= bits)
 979             high = center;
 980         else
 981             low = center;
 982     }
 983
 984     return (bits - (low == 0 ? -1 : cache[low]) <= cache[high] - bits) ? low : high;
 985 }
 986
 987 static inline int celt_pulses2bits(const uint8_t *cache, int pulses)
 988 {
 989     // TODO: Find the size of cache and make it into an array in the parameters list
 990    return (pulses == 0) ? 0 : cache[pulses] + 1;
 991 }
 992
 993 static inline void celt_normalize_residual(const int * restrict iy, float * restrict X,
 994                                            int N, float g)
 995 {
 996     int i;
 997     for (i = 0; i < N; i++)
 998         X[i] = g * iy[i];
 999 }
1000
1001 static void celt_exp_rotation1(float *X, unsigned int len, unsigned int stride,
1002                                float c, float s)
1003 {
1004     float *Xptr;
1005     int i;
1006
1007     Xptr = X;
1008     for (i = 0; i < len - stride; i++) {
1009         float x1, x2;
1010         x1           = Xptr[0];
1011         x2           = Xptr[stride];
1012         Xptr[stride] = c * x2 + s * x1;
1013         *Xptr++      = c * x1 - s * x2;
1014     }
1015
1016     Xptr = &X[len - 2 * stride - 1];
1017     for (i = len - 2 * stride - 1; i >= 0; i--) {
1018         float x1, x2;
1019         x1           = Xptr[0];
1020         x2           = Xptr[stride];
1021         Xptr[stride] = c * x2 + s * x1;
1022         *Xptr--      = c * x1 - s * x2;
1023     }
1024 }
1025
1026 static inline void celt_exp_rotation(float *X, unsigned int len,
1027                                      unsigned int stride, unsigned int K,
1028                                      enum CeltSpread spread)
1029 {
1030     unsigned int stride2 = 0;
1031     float c, s;
1032     float gain, theta;
1033     int i;
1034
1035     if (2*K >= len || spread == CELT_SPREAD_NONE)
1036         return;
1037
1038     gain = (float)len / (len + (20 - 5*spread) * K);
1039     theta = M_PI * gain * gain / 4;
1040
1041     c = cos(theta);
1042     s = sin(theta);
1043
1044     if (len >= stride << 3) {
1045         stride2 = 1;
1046         /* This is just a simple (equivalent) way of computing sqrt(len/stride) with rounding.
1047         It's basically incrementing long as (stride2+0.5)^2 < len/stride. */
1048         while ((stride2 * stride2 + stride2) * stride + (stride >> 2) < len)
1049             stride2++;
1050     }
1051
1052     /*NOTE: As a minor optimization, we could be passing around log2(B), not B, for both this and for
1053     extract_collapse_mask().*/
1054     len /= stride;
1055     for (i = 0; i < stride; i++) {
1056         if (stride2)
1057             celt_exp_rotation1(X + i * len, len, stride2, s, c);
1058         celt_exp_rotation1(X + i * len, len, 1, c, s);
1059     }
1060 }
1061
1062 static inline unsigned int celt_extract_collapse_mask(const int *iy,
1063                                                       unsigned int N,
1064                                                       unsigned int B)
1065 {
1066     unsigned int collapse_mask;
1067     int N0;
1068     int i, j;
1069
1070     if (B <= 1)
1071         return 1;
1072
1073     /*NOTE: As a minor optimization, we could be passing around log2(B), not B, for both this and for
1074     exp_rotation().*/
1075     N0 = N/B;
1076     collapse_mask = 0;
1077     for (i = 0; i < B; i++)
1078         for (j = 0; j < N0; j++)
1079             collapse_mask |= (iy[i*N0+j]!=0)<<i;
1080     return collapse_mask;
1081 }
1082
1083 static inline void celt_renormalize_vector(float *X, int N, float gain)
1084 {
1085     int i;
1086     float g = 1e-15f;
1087     for (i = 0; i < N; i++)
1088         g += X[i] * X[i];
1089     g = gain / sqrtf(g);
1090
1091     for (i = 0; i < N; i++)
1092         X[i] *= g;
1093 }
1094
1095 static inline void celt_stereo_merge(float *X, float *Y, float mid, int N)
1096 {
1097     int i;
1098     float xp = 0, side = 0;
1099     float E[2];
1100     float mid2;
1101     float t, gain[2];
1102
1103     /* Compute the norm of X+Y and X-Y as |X|^2 + |Y|^2 +/- sum(xy) */
1104     for (i = 0; i < N; i++) {
1105         xp   += X[i] * Y[i];
1106         side += Y[i] * Y[i];
1107     }
1108
1109     /* Compensating for the mid normalization */
1110     xp *= mid;
1111     mid2 = mid;
1112     E[0] = mid2 * mid2 + side - 2 * xp;
1113     E[1] = mid2 * mid2 + side + 2 * xp;
1114     if (E[0] < 6e-4f || E[1] < 6e-4f) {
1115         for (i = 0; i < N; i++)
1116             Y[i] = X[i];
1117         return;
1118     }
1119
1120     t = E[0];
1121     gain[0] = 1.0f / sqrtf(t);
1122     t = E[1];
1123     gain[1] = 1.0f / sqrtf(t);
1124
1125     for (i = 0; i < N; i++) {
1126         float value[2];
1127         /* Apply mid scaling (side is already scaled) */
1128         value[0] = mid * X[i];
1129         value[1] = Y[i];
1130         X[i] = gain[0] * (value[0] - value[1]);
1131         Y[i] = gain[1] * (value[0] + value[1]);
1132     }
1133 }
1134
1135 static void celt_interleave_hadamard(float *tmp, float *X, int N0,
1136                                      int stride, int hadamard)
1137 {
1138     int i, j;
1139     int N = N0*stride;
1140
1141     if (hadamard) {
1142         const uint8_t *ordery = celt_hadamard_ordery + stride - 2;
1143         for (i = 0; i < stride; i++)
1144             for (j = 0; j < N0; j++)
1145                 tmp[j*stride+i] = X[ordery[i]*N0+j];
1146     } else {
1147         for (i = 0; i < stride; i++)
1148             for (j = 0; j < N0; j++)
1149                 tmp[j*stride+i] = X[i*N0+j];
1150     }
1151
1152     for (i = 0; i < N; i++)
1153         X[i] = tmp[i];
1154 }
1155
1156 static void celt_deinterleave_hadamard(float *tmp, float *X, int N0,
1157                                        int stride, int hadamard)
1158 {
1159     int i, j;
1160     int N = N0*stride;
1161
1162     if (hadamard) {
1163         const uint8_t *ordery = celt_hadamard_ordery + stride - 2;
1164         for (i = 0; i < stride; i++)
1165             for (j = 0; j < N0; j++)
1166                 tmp[ordery[i]*N0+j] = X[j*stride+i];
1167     } else {
1168         for (i = 0; i < stride; i++)
1169             for (j = 0; j < N0; j++)
1170                 tmp[i*N0+j] = X[j*stride+i];
1171     }
1172
1173     for (i = 0; i < N; i++)
1174         X[i] = tmp[i];
1175 }
1176
1177 static void celt_haar1(float *X, int N0, int stride)
1178 {
1179     int i, j;
1180     N0 >>= 1;
1181     for (i = 0; i < stride; i++) {
1182         for (j = 0; j < N0; j++) {
1183             float x0 = X[stride * (2 * j + 0) + i];
1184             float x1 = X[stride * (2 * j + 1) + i];
1185             X[stride * (2 * j + 0) + i] = (x0 + x1) * M_SQRT1_2;
1186             X[stride * (2 * j + 1) + i] = (x0 - x1) * M_SQRT1_2;
1187         }
1188     }
1189 }
1190
1191 static inline int celt_compute_qn(int N, int b, int offset, int pulse_cap,
1192                                   int dualstereo)
1193 {
1194     int qn, qb;
1195     int N2 = 2 * N - 1;
1196     if (dualstereo && N == 2)
1197         N2--;
1198
1199     /* The upper limit ensures that in a stereo split with itheta==16384, we'll
1200      * always have enough bits left over to code at least one pulse in the
1201      * side; otherwise it would collapse, since it doesn't get folded. */
1202     qb = FFMIN3(b - pulse_cap - (4 << 3), (b + N2 * offset) / N2, 8 << 3);
1203     qn = (qb < (1 << 3 >> 1)) ? 1 : ((celt_qn_exp2[qb & 0x7] >> (14 - (qb >> 3))) + 1) >> 1 << 1;
1204     return qn;
1205 }
1206
1207 // this code was adapted from libopus
1208 static inline uint64_t celt_cwrsi(unsigned int N, unsigned int K, unsigned int i, int *y)
1209 {
1210     uint64_t norm = 0;
1211     uint32_t p;
1212     int s, val;
1213     int k0;
1214
1215     while (N > 2) {
1216         uint32_t q;
1217
1218         /*Lots of pulses case:*/
1219         if (K >= N) {
1220             const uint32_t *row = celt_pvq_u_row[N];
1221
1222             /* Are the pulses in this dimension negative? */
1223             p  = row[K + 1];
1224             s  = -(i >= p);
1225             i -= p & s;
1226
1227             /*Count how many pulses were placed in this dimension.*/
1228             k0 = K;
1229             q = row[N];
1230             if (q > i) {
1231                 K = N;
1232                 do {
1233                     p = celt_pvq_u_row[--K][N];
1234                 } while (p > i);
1235             } else
1236                 for (p = row[K]; p > i; p = row[K])
1237                     K--;
1238
1239             i    -= p;
1240             val   = (k0 - K + s) ^ s;
1241             norm += val * val;
1242             *y++  = val;
1243         } else { /*Lots of dimensions case:*/
1244             /*Are there any pulses in this dimension at all?*/
1245             p = celt_pvq_u_row[K    ][N];
1246             q = celt_pvq_u_row[K + 1][N];
1247
1248             if (p <= i && i < q) {
1249                 i -= p;
1250                 *y++ = 0;
1251             } else {
1252                 /*Are the pulses in this dimension negative?*/
1253                 s  = -(i >= q);
1254                 i -= q & s;
1255
1256                 /*Count how many pulses were placed in this dimension.*/
1257                 k0 = K;
1258                 do p = celt_pvq_u_row[--K][N];
1259                 while (p > i);
1260
1261                 i    -= p;
1262                 val   = (k0 - K + s) ^ s;
1263                 norm += val * val;
1264                 *y++  = val;
1265             }
1266         }
1267         N--;
1268     }
1269
1270     /* N == 2 */
1271     p  = 2 * K + 1;
1272     s  = -(i >= p);
1273     i -= p & s;
1274     k0 = K;
1275     K  = (i + 1) / 2;
1276
1277     if (K)
1278         i -= 2 * K - 1;
1279
1280     val   = (k0 - K + s) ^ s;
1281     norm += val * val;
1282     *y++  = val;
1283
1284     /* N==1 */
1285     s     = -i;
1286     val   = (K + s) ^ s;
1287     norm += val * val;
1288     *y    = val;
1289
1290     return norm;
1291 }
1292
1293 static inline float celt_decode_pulses(OpusRangeCoder *rc, int *y, unsigned int N, unsigned int K)
1294 {
1295     unsigned int idx;
1296 #define CELT_PVQ_U(n, k) (celt_pvq_u_row[FFMIN(n, k)][FFMAX(n, k)])
1297 #define CELT_PVQ_V(n, k) (CELT_PVQ_U(n, k) + CELT_PVQ_U(n, k + 1))
1298     idx = opus_rc_unimodel(rc, CELT_PVQ_V(N, K));
1299     return celt_cwrsi(N, K, idx, y);
1300 }
1301
1302 /** Decode pulse vector and combine the result with the pitch vector to produce
1303     the final normalised signal in the current band. */
1304 static inline unsigned int celt_alg_unquant(OpusRangeCoder *rc, float *X,
1305                                             unsigned int N, unsigned int K,
1306                                             enum CeltSpread spread,
1307                                             unsigned int blocks, float gain)
1308 {
1309     int y[176];
1310
1311     gain /= sqrtf(celt_decode_pulses(rc, y, N, K));
1312     celt_normalize_residual(y, X, N, gain);
1313     celt_exp_rotation(X, N, blocks, K, spread);
1314     return celt_extract_collapse_mask(y, N, blocks);
1315 }
1316
1317 static unsigned int celt_decode_band(CeltContext *s, OpusRangeCoder *rc,
1318                                      const int band, float *X, float *Y,
1319                                      int N, int b, unsigned int blocks,
1320                                      float *lowband, int duration,
1321                                      float *lowband_out, int level,
1322                                      float gain, float *lowband_scratch,
1323                                      int fill)
1324 {
1325     const uint8_t *cache;
1326     int dualstereo, split;
1327     int imid = 0, iside = 0;
1328     unsigned int N0 = N;
1329     int N_B;
1330     int N_B0;
1331     int B0 = blocks;
1332     int time_divide = 0;
1333     int recombine = 0;
1334     int inv = 0;
1335     float mid = 0, side = 0;
1336     int longblocks = (B0 == 1);
1337     unsigned int cm = 0;
1338
1339     N_B0 = N_B = N / blocks;
1340     split = dualstereo = (Y != NULL);
1341
1342     if (N == 1) {
1343         /* special case for one sample */
1344         int i;
1345         float *x = X;
1346         for (i = 0; i <= dualstereo; i++) {
1347             int sign = 0;
1348             if (s->remaining2 >= 1<<3) {
1349                 sign           = opus_getrawbits(rc, 1);
1350                 s->remaining2 -= 1 << 3;
1351                 b             -= 1 << 3;
1352             }
1353             x[0] = sign ? -1.0f : 1.0f;
1354             x = Y;
1355         }
1356         if (lowband_out)
1357             lowband_out[0] = X[0];
1358         return 1;
1359     }
1360
1361     if (!dualstereo && level == 0) {
1362         int tf_change = s->tf_change[band];
1363         int k;
1364         if (tf_change > 0)
1365             recombine = tf_change;
1366         /* Band recombining to increase frequency resolution */
1367
1368         if (lowband &&
1369             (recombine || ((N_B & 1) == 0 && tf_change < 0) || B0 > 1)) {
1370             int j;
1371             for (j = 0; j < N; j++)
1372                 lowband_scratch[j] = lowband[j];
1373             lowband = lowband_scratch;
1374         }
1375
1376         for (k = 0; k < recombine; k++) {
1377             if (lowband)
1378                 celt_haar1(lowband, N >> k, 1 << k);
1379             fill = celt_bit_interleave[fill & 0xF] | celt_bit_interleave[fill >> 4] << 2;
1380         }
1381         blocks >>= recombine;
1382         N_B <<= recombine;
1383
1384         /* Increasing the time resolution */
1385         while ((N_B & 1) == 0 && tf_change < 0) {
1386             if (lowband)
1387                 celt_haar1(lowband, N_B, blocks);
1388             fill |= fill << blocks;
1389             blocks <<= 1;
1390             N_B >>= 1;
1391             time_divide++;
1392             tf_change++;
1393         }
1394         B0 = blocks;
1395         N_B0 = N_B;
1396
1397         /* Reorganize the samples in time order instead of frequency order */
1398         if (B0 > 1 && lowband)
1399             celt_deinterleave_hadamard(s->scratch, lowband, N_B >> recombine,
1400                                        B0 << recombine, longblocks);
1401     }
1402
1403     /* If we need 1.5 more bit than we can produce, split the band in two. */
1404     cache = celt_cache_bits +
1405             celt_cache_index[(duration + 1) * CELT_MAX_BANDS + band];
1406     if (!dualstereo && duration >= 0 && b > cache[cache[0]] + 12 && N > 2) {
1407         N >>= 1;
1408         Y = X + N;
1409         split = 1;
1410         duration -= 1;
1411         if (blocks == 1)
1412             fill = (fill & 1) | (fill << 1);
1413         blocks = (blocks + 1) >> 1;
1414     }
1415
1416     if (split) {
1417         int qn;
1418         int itheta = 0;
1419         int mbits, sbits, delta;
1420         int qalloc;
1421         int pulse_cap;
1422         int offset;
1423         int orig_fill;
1424         int tell;
1425
1426         /* Decide on the resolution to give to the split parameter theta */
1427         pulse_cap = celt_log_freq_range[band] + duration * 8;
1428         offset = (pulse_cap >> 1) - (dualstereo && N == 2 ? CELT_QTHETA_OFFSET_TWOPHASE :
1429                                                           CELT_QTHETA_OFFSET);
1430         qn = (dualstereo && band >= s->intensitystereo) ? 1 :
1431              celt_compute_qn(N, b, offset, pulse_cap, dualstereo);
1432         tell = opus_rc_tell_frac(rc);
1433         if (qn != 1) {
1434             /* Entropy coding of the angle. We use a uniform pdf for the
1435             time split, a step for stereo, and a triangular one for the rest. */
1436             if (dualstereo && N > 2)
1437                 itheta = opus_rc_stepmodel(rc, qn/2);
1438             else if (dualstereo || B0 > 1)
1439                 itheta = opus_rc_unimodel(rc, qn+1);
1440             else
1441                 itheta = opus_rc_trimodel(rc, qn);
1442             itheta = itheta * 16384 / qn;
1443             /* NOTE: Renormalising X and Y *may* help fixed-point a bit at very high rate.
1444             Let's do that at higher complexity */
1445         } else if (dualstereo) {
1446             inv = (b > 2 << 3 && s->remaining2 > 2 << 3) ? opus_rc_p2model(rc, 2) : 0;
1447             itheta = 0;
1448         }
1449         qalloc = opus_rc_tell_frac(rc) - tell;
1450         b -= qalloc;
1451
1452         orig_fill = fill;
1453         if (itheta == 0) {
1454             imid = 32767;
1455             iside = 0;
1456             fill &= (1 << blocks) - 1;
1457             delta = -16384;
1458         } else if (itheta == 16384) {
1459             imid = 0;
1460             iside = 32767;
1461             fill &= ((1 << blocks) - 1) << blocks;
1462             delta = 16384;
1463         } else {
1464             imid = celt_cos(itheta);
1465             iside = celt_cos(16384-itheta);
1466             /* This is the mid vs side allocation that minimizes squared error
1467             in that band. */
1468             delta = ROUND_MUL16((N - 1) << 7, celt_log2tan(iside, imid));
1469         }
1470
1471         mid  = imid  / 32768.0f;
1472         side = iside / 32768.0f;
1473
1474         /* This is a special case for N=2 that only works for stereo and takes
1475         advantage of the fact that mid and side are orthogonal to encode
1476         the side with just one bit. */
1477         if (N == 2 && dualstereo) {
1478             int c;
1479             int sign = 0;
1480             float tmp;
1481             float *x2, *y2;
1482             mbits = b;
1483             /* Only need one bit for the side */
1484             sbits = (itheta != 0 && itheta != 16384) ? 1 << 3 : 0;
1485             mbits -= sbits;
1486             c = (itheta > 8192);
1487             s->remaining2 -= qalloc+sbits;
1488
1489             x2 = c ? Y : X;
1490             y2 = c ? X : Y;
1491             if (sbits)
1492                 sign = opus_getrawbits(rc, 1);
1493             sign = 1 - 2 * sign;
1494             /* We use orig_fill here because we want to fold the side, but if
1495             itheta==16384, we'll have cleared the low bits of fill. */
1496             cm = celt_decode_band(s, rc, band, x2, NULL, N, mbits, blocks,
1497                                   lowband, duration, lowband_out, level, gain,
1498                                   lowband_scratch, orig_fill);
1499             /* We don't split N=2 bands, so cm is either 1 or 0 (for a fold-collapse),
1500             and there's no need to worry about mixing with the other channel. */
1501             y2[0] = -sign * x2[1];
1502             y2[1] =  sign * x2[0];
1503             X[0] *= mid;
1504             X[1] *= mid;
1505             Y[0] *= side;
1506             Y[1] *= side;
1507             tmp = X[0];
1508             X[0] = tmp - Y[0];
1509             Y[0] = tmp + Y[0];
1510             tmp = X[1];
1511             X[1] = tmp - Y[1];
1512             Y[1] = tmp + Y[1];
1513         } else {
1514             /* "Normal" split code */
1515             float *next_lowband2     = NULL;
1516             float *next_lowband_out1 = NULL;
1517             int next_level = 0;
1518             int rebalance;
1519
1520             /* Give more bits to low-energy MDCTs than they would
1521              * otherwise deserve */
1522             if (B0 > 1 && !dualstereo && (itheta & 0x3fff)) {
1523                 if (itheta > 8192)
1524                     /* Rough approximation for pre-echo masking */
1525                     delta -= delta >> (4 - duration);
1526                 else
1527                     /* Corresponds to a forward-masking slope of
1528                      * 1.5 dB per 10 ms */
1529                     delta = FFMIN(0, delta + (N << 3 >> (5 - duration)));
1530             }
1531             mbits = av_clip((b - delta) / 2, 0, b);
1532             sbits = b - mbits;
1533             s->remaining2 -= qalloc;
1534
1535             if (lowband && !dualstereo)
1536                 next_lowband2 = lowband + N; /* >32-bit split case */
1537
1538             /* Only stereo needs to pass on lowband_out.
1539              * Otherwise, it's handled at the end */
1540             if (dualstereo)
1541                 next_lowband_out1 = lowband_out;
1542             else
1543                 next_level = level + 1;
1544
1545             rebalance = s->remaining2;
1546             if (mbits >= sbits) {
1547                 /* In stereo mode, we do not apply a scaling to the mid
1548                  * because we need the normalized mid for folding later */
1549                 cm = celt_decode_band(s, rc, band, X, NULL, N, mbits, blocks,
1550                                       lowband, duration, next_lowband_out1,
1551                                       next_level, dualstereo ? 1.0f : (gain * mid),
1552                                       lowband_scratch, fill);
1553
1554                 rebalance = mbits - (rebalance - s->remaining2);
1555                 if (rebalance > 3 << 3 && itheta != 0)
1556                     sbits += rebalance - (3 << 3);
1557
1558                 /* For a stereo split, the high bits of fill are always zero,
1559                  * so no folding will be done to the side. */
1560                 cm |= celt_decode_band(s, rc, band, Y, NULL, N, sbits, blocks,
1561                                        next_lowband2, duration, NULL,
1562                                        next_level, gain * side, NULL,
1563                                        fill >> blocks) << ((B0 >> 1) & (dualstereo - 1));
1564             } else {
1565                 /* For a stereo split, the high bits of fill are always zero,
1566                  * so no folding will be done to the side. */
1567                 cm = celt_decode_band(s, rc, band, Y, NULL, N, sbits, blocks,
1568                                       next_lowband2, duration, NULL,
1569                                       next_level, gain * side, NULL,
1570                                       fill >> blocks) << ((B0 >> 1) & (dualstereo - 1));
1571
1572                 rebalance = sbits - (rebalance - s->remaining2);
1573                 if (rebalance > 3 << 3 && itheta != 16384)
1574                     mbits += rebalance - (3 << 3);
1575
1576                 /* In stereo mode, we do not apply a scaling to the mid because
1577                  * we need the normalized mid for folding later */
1578                 cm |= celt_decode_band(s, rc, band, X, NULL, N, mbits, blocks,
1579                                        lowband, duration, next_lowband_out1,
1580                                        next_level, dualstereo ? 1.0f : (gain * mid),
1581                                        lowband_scratch, fill);
1582             }
1583         }
1584     } else {
1585         /* This is the basic no-split case */
1586         unsigned int q         = celt_bits2pulses(cache, b);
1587         unsigned int curr_bits = celt_pulses2bits(cache, q);
1588         s->remaining2 -= curr_bits;
1589
1590         /* Ensures we can never bust the budget */
1591         while (s->remaining2 < 0 && q > 0) {
1592             s->remaining2 += curr_bits;
1593             curr_bits      = celt_pulses2bits(cache, --q);
1594             s->remaining2 -= curr_bits;
1595         }
1596
1597         if (q != 0) {
1598             /* Finally do the actual quantization */
1599             cm = celt_alg_unquant(rc, X, N, (q < 8) ? q : (8 + (q & 7)) << ((q >> 3) - 1),
1600                                   s->spread, blocks, gain);
1601         } else {
1602             /* If there's no pulse, fill the band anyway */
1603             int j;
1604             unsigned int cm_mask = (1 << blocks) - 1;
1605             fill &= cm_mask;
1606             if (!fill) {
1607                 for (j = 0; j < N; j++)
1608                     X[j] = 0.0f;
1609             } else {
1610                 if (lowband == NULL) {
1611                     /* Noise */
1612                     for (j = 0; j < N; j++)
1613                         X[j] = (((int32_t)celt_rng(s)) >> 20);
1614                     cm = cm_mask;
1615                 } else {
1616                     /* Folded spectrum */
1617                     for (j = 0; j < N; j++) {
1618                         /* About 48 dB below the "normal" folding level */
1619                         X[j] = lowband[j] + (((celt_rng(s)) & 0x8000) ? 1.0f / 256 : -1.0f / 256);
1620                     }
1621                     cm = fill;
1622                 }
1623                 celt_renormalize_vector(X, N, gain);
1624             }
1625         }
1626     }
1627
1628     /* This code is used by the decoder and by the resynthesis-enabled encoder */
1629     if (dualstereo) {
1630         int j;
1631         if (N != 2)
1632             celt_stereo_merge(X, Y, mid, N);
1633         if (inv) {
1634             for (j = 0; j < N; j++)
1635                 Y[j] *= -1;
1636         }
1637     } else if (level == 0) {
1638         int k;
1639
1640         /* Undo the sample reorganization going from time order to frequency order */
1641         if (B0 > 1)
1642             celt_interleave_hadamard(s->scratch, X, N_B>>recombine,
1643                                      B0<<recombine, longblocks);
1644
1645         /* Undo time-freq changes that we did earlier */
1646         N_B = N_B0;
1647         blocks = B0;
1648         for (k = 0; k < time_divide; k++) {
1649             blocks >>= 1;
1650             N_B <<= 1;
1651             cm |= cm >> blocks;
1652             celt_haar1(X, N_B, blocks);
1653         }
1654
1655         for (k = 0; k < recombine; k++) {
1656             cm = celt_bit_deinterleave[cm];
1657             celt_haar1(X, N0>>k, 1<<k);
1658         }
1659         blocks <<= recombine;
1660
1661         /* Scale output for later folding */
1662         if (lowband_out) {
1663             int j;
1664             float n = sqrtf(N0);
1665             for (j = 0; j < N0; j++)
1666                 lowband_out[j] = n * X[j];
1667         }
1668         cm &= (1 << blocks) - 1;
1669     }
1670     return cm;
1671 }
1672
1673 static void celt_denormalize(CeltContext *s, CeltFrame *frame, float *data)
1674 {
1675     int i, j;
1676
1677     for (i = s->startband; i < s->endband; i++) {
1678         float *dst = data + (celt_freq_bands[i] << s->duration);
1679         float norm = pow(2, frame->energy[i] + celt_mean_energy[i]);
1680
1681         for (j = 0; j < celt_freq_range[i] << s->duration; j++)
1682             dst[j] *= norm;
1683     }
1684 }
1685
1686 static void celt_postfilter_apply_transition(CeltFrame *frame, float *data)
1687 {
1688     const int T0 = frame->pf_period_old;
1689     const int T1 = frame->pf_period;
1690
1691     float g00, g01, g02;
1692     float g10, g11, g12;
1693
1694     float x0, x1, x2, x3, x4;
1695
1696     int i;
1697
1698     if (frame->pf_gains[0]     == 0.0 &&
1699         frame->pf_gains_old[0] == 0.0)
1700         return;
1701
1702     g00 = frame->pf_gains_old[0];
1703     g01 = frame->pf_gains_old[1];
1704     g02 = frame->pf_gains_old[2];
1705     g10 = frame->pf_gains[0];
1706     g11 = frame->pf_gains[1];
1707     g12 = frame->pf_gains[2];
1708
1709     x1 = data[-T1 + 1];
1710     x2 = data[-T1];
1711     x3 = data[-T1 - 1];
1712     x4 = data[-T1 - 2];
1713
1714     for (i = 0; i < CELT_OVERLAP; i++) {
1715         float w = ff_celt_window2[i];
1716         x0 = data[i - T1 + 2];
1717
1718         data[i] +=  (1.0 - w) * g00 * data[i - T0]                          +
1719                     (1.0 - w) * g01 * (data[i - T0 - 1] + data[i - T0 + 1]) +
1720                     (1.0 - w) * g02 * (data[i - T0 - 2] + data[i - T0 + 2]) +
1721                     w         * g10 * x2                                    +
1722                     w         * g11 * (x1 + x3)                             +
1723                     w         * g12 * (x0 + x4);
1724         x4 = x3;
1725         x3 = x2;
1726         x2 = x1;
1727         x1 = x0;
1728     }
1729 }
1730
1731 static void celt_postfilter_apply(CeltFrame *frame,
1732                                   float *data, int len)
1733 {
1734     const int T = frame->pf_period;
1735     float g0, g1, g2;
1736     float x0, x1, x2, x3, x4;
1737     int i;
1738
1739     if (frame->pf_gains[0] == 0.0 || len <= 0)
1740         return;
1741
1742     g0 = frame->pf_gains[0];
1743     g1 = frame->pf_gains[1];
1744     g2 = frame->pf_gains[2];
1745
1746     x4 = data[-T - 2];
1747     x3 = data[-T - 1];
1748     x2 = data[-T];
1749     x1 = data[-T + 1];
1750
1751     for (i = 0; i < len; i++) {
1752         x0 = data[i - T + 2];
1753         data[i] += g0 * x2        +
1754                    g1 * (x1 + x3) +
1755                    g2 * (x0 + x4);
1756         x4 = x3;
1757         x3 = x2;
1758         x2 = x1;
1759         x1 = x0;
1760     }
1761 }
1762
1763 static void celt_postfilter(CeltContext *s, CeltFrame *frame)
1764 {
1765     int len = s->blocksize * s->blocks;
1766
1767     celt_postfilter_apply_transition(frame, frame->buf + 1024);
1768
1769     frame->pf_period_old = frame->pf_period;
1770     memcpy(frame->pf_gains_old, frame->pf_gains, sizeof(frame->pf_gains));
1771
1772     frame->pf_period = frame->pf_period_new;
1773     memcpy(frame->pf_gains, frame->pf_gains_new, sizeof(frame->pf_gains));
1774
1775     if (len > CELT_OVERLAP) {
1776         celt_postfilter_apply_transition(frame, frame->buf + 1024 + CELT_OVERLAP);
1777         celt_postfilter_apply(frame, frame->buf + 1024 + 2 * CELT_OVERLAP,
1778                               len - 2 * CELT_OVERLAP);
1779
1780         frame->pf_period_old = frame->pf_period;
1781         memcpy(frame->pf_gains_old, frame->pf_gains, sizeof(frame->pf_gains));
1782     }
1783
1784     memmove(frame->buf, frame->buf + len, (1024 + CELT_OVERLAP / 2) * sizeof(float));
1785 }
1786
1787 static int parse_postfilter(CeltContext *s, OpusRangeCoder *rc, int consumed)
1788 {
1789     static const float postfilter_taps[3][3] = {
1790         { 0.3066406250f, 0.2170410156f, 0.1296386719f },
1791         { 0.4638671875f, 0.2680664062f, 0.0           },
1792         { 0.7998046875f, 0.1000976562f, 0.0           }
1793     };
1794     int i;
1795
1796     memset(s->frame[0].pf_gains_new, 0, sizeof(s->frame[0].pf_gains_new));
1797     memset(s->frame[1].pf_gains_new, 0, sizeof(s->frame[1].pf_gains_new));
1798
1799     if (s->startband == 0 && consumed + 16 <= s->framebits) {
1800         int has_postfilter = opus_rc_p2model(rc, 1);
1801         if (has_postfilter) {
1802             float gain;
1803             int tapset, octave, period;
1804
1805             octave = opus_rc_unimodel(rc, 6);
1806             period = (16 << octave) + opus_getrawbits(rc, 4 + octave) - 1;
1807             gain   = 0.09375f * (opus_getrawbits(rc, 3) + 1);
1808             tapset = (opus_rc_tell(rc) + 2 <= s->framebits) ?
1809                      opus_rc_getsymbol(rc, celt_model_tapset) : 0;
1810
1811             for (i = 0; i < 2; i++) {
1812                 CeltFrame *frame = &s->frame[i];
1813
1814                 frame->pf_period_new = FFMAX(period, CELT_POSTFILTER_MINPERIOD);
1815                 frame->pf_gains_new[0] = gain * postfilter_taps[tapset][0];
1816                 frame->pf_gains_new[1] = gain * postfilter_taps[tapset][1];
1817                 frame->pf_gains_new[2] = gain * postfilter_taps[tapset][2];
1818             }
1819         }
1820
1821         consumed = opus_rc_tell(rc);
1822     }
1823
1824     return consumed;
1825 }
1826
1827 static void process_anticollapse(CeltContext *s, CeltFrame *frame, float *X)
1828 {
1829     int i, j, k;
1830
1831     for (i = s->startband; i < s->endband; i++) {
1832         int renormalize = 0;
1833         float *xptr;
1834         float prev[2];
1835         float Ediff, r;
1836         float thresh, sqrt_1;
1837         int depth;
1838
1839         /* depth in 1/8 bits */
1840         depth = (1 + s->pulses[i]) / (celt_freq_range[i] << s->duration);
1841         thresh = pow(2, -1.0 - 0.125f * depth);
1842         sqrt_1 = 1.0f / sqrtf(celt_freq_range[i] << s->duration);
1843
1844         xptr = X + (celt_freq_bands[i] << s->duration);
1845
1846         prev[0] = frame->prev_energy[0][i];
1847         prev[1] = frame->prev_energy[1][i];
1848         if (s->coded_channels == 1) {
1849             CeltFrame *frame1 = &s->frame[1];
1850
1851             prev[0] = FFMAX(prev[0], frame1->prev_energy[0][i]);
1852             prev[1] = FFMAX(prev[1], frame1->prev_energy[1][i]);
1853         }
1854         Ediff = frame->energy[i] - FFMIN(prev[0], prev[1]);
1855         Ediff = FFMAX(0, Ediff);
1856
1857         /* r needs to be multiplied by 2 or 2*sqrt(2) depending on LM because
1858         short blocks don't have the same energy as long */
1859         r = pow(2, 1 - Ediff);
1860         if (s->duration == 3)
1861             r *= M_SQRT2;
1862         r = FFMIN(thresh, r) * sqrt_1;
1863         for (k = 0; k < 1 << s->duration; k++) {
1864             /* Detect collapse */
1865             if (!(frame->collapse_masks[i] & 1 << k)) {
1866                 /* Fill with noise */
1867                 for (j = 0; j < celt_freq_range[i]; j++)
1868                     xptr[(j << s->duration) + k] = (celt_rng(s) & 0x8000) ? r : -r;
1869                 renormalize = 1;
1870             }
1871         }
1872
1873         /* We just added some energy, so we need to renormalize */
1874         if (renormalize)
1875             celt_renormalize_vector(xptr, celt_freq_range[i] << s->duration, 1.0f);
1876     }
1877 }
1878
1879 static void celt_decode_bands(CeltContext *s, OpusRangeCoder *rc)
1880 {
1881     float lowband_scratch[8 * 22];
1882     float norm[2 * 8 * 100];
1883
1884     int totalbits = (s->framebits << 3) - s->anticollapse_bit;
1885
1886     int update_lowband = 1;
1887     int lowband_offset = 0;
1888
1889     int i, j;
1890
1891     memset(s->coeffs, 0, sizeof(s->coeffs));
1892
1893     for (i = s->startband; i < s->endband; i++) {
1894         int band_offset = celt_freq_bands[i] << s->duration;
1895         int band_size   = celt_freq_range[i] << s->duration;
1896         float *X = s->coeffs[0] + band_offset;
1897         float *Y = (s->coded_channels == 2) ? s->coeffs[1] + band_offset : NULL;
1898
1899         int consumed = opus_rc_tell_frac(rc);
1900         float *norm2 = norm + 8 * 100;
1901         int effective_lowband = -1;
1902         unsigned int cm[2];
1903         int b;
1904
1905         /* Compute how many bits we want to allocate to this band */
1906         if (i != s->startband)
1907             s->remaining -= consumed;
1908         s->remaining2 = totalbits - consumed - 1;
1909         if (i <= s->codedbands - 1) {
1910             int curr_balance = s->remaining / FFMIN(3, s->codedbands-i);
1911             b = av_clip(FFMIN(s->remaining2 + 1, s->pulses[i] + curr_balance), 0, 16383);
1912         } else
1913             b = 0;
1914
1915         if (celt_freq_bands[i] - celt_freq_range[i] >= celt_freq_bands[s->startband] &&
1916             (update_lowband || lowband_offset == 0))
1917             lowband_offset = i;
1918
1919         /* Get a conservative estimate of the collapse_mask's for the bands we're
1920         going to be folding from. */
1921         if (lowband_offset != 0 && (s->spread != CELT_SPREAD_AGGRESSIVE ||
1922                                     s->blocks > 1 || s->tf_change[i] < 0)) {
1923             int foldstart, foldend;
1924
1925             /* This ensures we never repeat spectral content within one band */
1926             effective_lowband = FFMAX(celt_freq_bands[s->startband],
1927                                       celt_freq_bands[lowband_offset] - celt_freq_range[i]);
1928             foldstart = lowband_offset;
1929             while (celt_freq_bands[--foldstart] > effective_lowband);
1930             foldend = lowband_offset - 1;
1931             while (celt_freq_bands[++foldend] < effective_lowband + celt_freq_range[i]);
1932
1933             cm[0] = cm[1] = 0;
1934             for (j = foldstart; j < foldend; j++) {
1935                 cm[0] |= s->frame[0].collapse_masks[j];
1936                 cm[1] |= s->frame[s->coded_channels - 1].collapse_masks[j];
1937             }
1938         } else
1939             /* Otherwise, we'll be using the LCG to fold, so all blocks will (almost
1940             always) be non-zero.*/
1941             cm[0] = cm[1] = (1 << s->blocks) - 1;
1942
1943         if (s->dualstereo && i == s->intensitystereo) {
1944             /* Switch off dual stereo to do intensity */
1945             s->dualstereo = 0;
1946             for (j = celt_freq_bands[s->startband] << s->duration; j < band_offset; j++)
1947                 norm[j] = (norm[j] + norm2[j]) / 2;
1948         }
1949
1950         if (s->dualstereo) {
1951             cm[0] = celt_decode_band(s, rc, i, X, NULL, band_size, b / 2, s->blocks,
1952                                      effective_lowband != -1 ? norm + (effective_lowband << s->duration) : NULL, s->duration,
1953             norm + band_offset, 0, 1.0f, lowband_scratch, cm[0]);
1954
1955             cm[1] = celt_decode_band(s, rc, i, Y, NULL, band_size, b/2, s->blocks,
1956                                      effective_lowband != -1 ? norm2 + (effective_lowband << s->duration) : NULL, s->duration,
1957             norm2 + band_offset, 0, 1.0f, lowband_scratch, cm[1]);
1958         } else {
1959             cm[0] = celt_decode_band(s, rc, i, X, Y, band_size, b, s->blocks,
1960             effective_lowband != -1 ? norm + (effective_lowband << s->duration) : NULL, s->duration,
1961             norm + band_offset, 0, 1.0f, lowband_scratch, cm[0]|cm[1]);
1962
1963             cm[1] = cm[0];
1964         }
1965
1966         s->frame[0].collapse_masks[i]                     = (uint8_t)cm[0];
1967         s->frame[s->coded_channels - 1].collapse_masks[i] = (uint8_t)cm[1];
1968         s->remaining += s->pulses[i] + consumed;
1969
1970         /* Update the folding position only as long as we have 1 bit/sample depth */
1971         update_lowband = (b > band_size << 3);
1972     }
1973 }
1974
1975 int ff_celt_decode_frame(CeltContext *s, OpusRangeCoder *rc,
1976                          float **output, int coded_channels, int frame_size,
1977                          int startband,  int endband)
1978 {
1979     int i, j;
1980
1981     int consumed;           // bits of entropy consumed thus far for this frame
1982     int silence = 0;
1983     int transient = 0;
1984     int anticollapse = 0;
1985     CeltIMDCTContext *imdct;
1986     float imdct_scale = 1.0;
1987
1988     if (coded_channels != 1 && coded_channels != 2) {
1989         av_log(s->avctx, AV_LOG_ERROR, "Invalid number of coded channels: %d\n",
1990                coded_channels);
1991         return AVERROR_INVALIDDATA;
1992     }
1993     if (startband < 0 || startband > endband || endband > CELT_MAX_BANDS) {
1994         av_log(s->avctx, AV_LOG_ERROR, "Invalid start/end band: %d %d\n",
1995                startband, endband);
1996         return AVERROR_INVALIDDATA;
1997     }
1998
1999     s->flushed        = 0;
2000     s->coded_channels = coded_channels;
2001     s->startband      = startband;
2002     s->endband        = endband;
2003     s->framebits      = rc->rb.bytes * 8;
2004
2005     s->duration = av_log2(frame_size / CELT_SHORT_BLOCKSIZE);
2006     if (s->duration > CELT_MAX_LOG_BLOCKS ||
2007         frame_size != CELT_SHORT_BLOCKSIZE * (1 << s->duration)) {
2008         av_log(s->avctx, AV_LOG_ERROR, "Invalid CELT frame size: %d\n",
2009                frame_size);
2010         return AVERROR_INVALIDDATA;
2011     }
2012
2013     if (!s->output_channels)
2014         s->output_channels = coded_channels;
2015
2016     memset(s->frame[0].collapse_masks, 0, sizeof(s->frame[0].collapse_masks));
2017     memset(s->frame[1].collapse_masks, 0, sizeof(s->frame[1].collapse_masks));
2018
2019     consumed = opus_rc_tell(rc);
2020
2021     /* obtain silence flag */
2022     if (consumed >= s->framebits)
2023         silence = 1;
2024     else if (consumed == 1)
2025         silence = opus_rc_p2model(rc, 15);
2026
2027
2028     if (silence) {
2029         consumed = s->framebits;
2030         rc->total_read_bits += s->framebits - opus_rc_tell(rc);
2031     }
2032
2033     /* obtain post-filter options */
2034     consumed = parse_postfilter(s, rc, consumed);
2035
2036     /* obtain transient flag */
2037     if (s->duration != 0 && consumed+3 <= s->framebits)
2038         transient = opus_rc_p2model(rc, 3);
2039
2040     s->blocks    = transient ? 1 << s->duration : 1;
2041     s->blocksize = frame_size / s->blocks;
2042
2043     imdct = s->imdct[transient ? 0 : s->duration];
2044
2045     if (coded_channels == 1) {
2046         for (i = 0; i < CELT_MAX_BANDS; i++)
2047             s->frame[0].energy[i] = FFMAX(s->frame[0].energy[i], s->frame[1].energy[i]);
2048     }
2049
2050     celt_decode_coarse_energy(s, rc);
2051     celt_decode_tf_changes   (s, rc, transient);
2052     celt_decode_allocation   (s, rc);
2053     celt_decode_fine_energy  (s, rc);
2054     celt_decode_bands        (s, rc);
2055
2056     if (s->anticollapse_bit)
2057         anticollapse = opus_getrawbits(rc, 1);
2058
2059     celt_decode_final_energy(s, rc, s->framebits - opus_rc_tell(rc));
2060
2061     /* apply anti-collapse processing and denormalization to
2062      * each coded channel */
2063     for (i = 0; i < s->coded_channels; i++) {
2064         CeltFrame *frame = &s->frame[i];
2065
2066         if (anticollapse)
2067             process_anticollapse(s, frame, s->coeffs[i]);
2068
2069         celt_denormalize(s, frame, s->coeffs[i]);
2070     }
2071
2072     /* stereo -> mono downmix */
2073     if (s->output_channels < s->coded_channels) {
2074         s->dsp.vector_fmac_scalar(s->coeffs[0], s->coeffs[1], 1.0, FFALIGN(frame_size, 16));
2075         imdct_scale = 0.5;
2076     } else if (s->output_channels > s->coded_channels)
2077         memcpy(s->coeffs[1], s->coeffs[0], frame_size * sizeof(float));
2078
2079     if (silence) {
2080         for (i = 0; i < 2; i++) {
2081             CeltFrame *frame = &s->frame[i];
2082
2083             for (j = 0; j < FF_ARRAY_ELEMS(frame->energy); j++)
2084                 frame->energy[j] = CELT_ENERGY_SILENCE;
2085         }
2086         memset(s->coeffs, 0, sizeof(s->coeffs));
2087     }
2088
2089     /* transform and output for each output channel */
2090     for (i = 0; i < s->output_channels; i++) {
2091         CeltFrame *frame = &s->frame[i];
2092         float m = frame->deemph_coeff;
2093
2094         /* iMDCT and overlap-add */
2095         for (j = 0; j < s->blocks; j++) {
2096             float *dst  = frame->buf + 1024 + j * s->blocksize;
2097
2098             ff_celt_imdct_half(imdct, dst + CELT_OVERLAP / 2, s->coeffs[i] + j,
2099                                s->blocks, imdct_scale);
2100             s->dsp.vector_fmul_window(dst, dst, dst + CELT_OVERLAP / 2,
2101                                       celt_window, CELT_OVERLAP / 2);
2102         }
2103
2104         /* postfilter */
2105         celt_postfilter(s, frame);
2106
2107         /* deemphasis and output scaling */
2108         for (j = 0; j < frame_size; j++) {
2109             float tmp = frame->buf[1024 - frame_size + j] + m;
2110             m = tmp * CELT_DEEMPH_COEFF;
2111             output[i][j] = tmp / 32768.;
2112         }
2113         frame->deemph_coeff = m;
2114     }
2115
2116     if (coded_channels == 1)
2117         memcpy(s->frame[1].energy, s->frame[0].energy, sizeof(s->frame[0].energy));
2118
2119     for (i = 0; i < 2; i++ ) {
2120         CeltFrame *frame = &s->frame[i];
2121
2122         if (!transient) {
2123             memcpy(frame->prev_energy[1], frame->prev_energy[0], sizeof(frame->prev_energy[0]));
2124             memcpy(frame->prev_energy[0], frame->energy,         sizeof(frame->prev_energy[0]));
2125         } else {
2126             for (j = 0; j < CELT_MAX_BANDS; j++)
2127                 frame->prev_energy[0][j] = FFMIN(frame->prev_energy[0][j], frame->energy[j]);
2128         }
2129
2130         for (j = 0; j < s->startband; j++) {
2131             frame->prev_energy[0][j] = CELT_ENERGY_SILENCE;
2132             frame->energy[j]         = 0.0;
2133         }
2134         for (j = s->endband; j < CELT_MAX_BANDS; j++) {
2135             frame->prev_energy[0][j] = CELT_ENERGY_SILENCE;
2136             frame->energy[j]         = 0.0;
2137         }
2138     }
2139
2140     s->seed = rc->range;
2141
2142     return 0;
2143 }
2144
2145 void ff_celt_flush(CeltContext *s)
2146 {
2147     int i, j;
2148
2149     if (s->flushed)
2150         return;
2151
2152     for (i = 0; i < 2; i++) {
2153         CeltFrame *frame = &s->frame[i];
2154
2155         for (j = 0; j < CELT_MAX_BANDS; j++)
2156             frame->prev_energy[0][j] = frame->prev_energy[1][j] = CELT_ENERGY_SILENCE;
2157
2158         memset(frame->energy, 0, sizeof(frame->energy));
2159         memset(frame->buf,    0, sizeof(frame->buf));
2160
2161         memset(frame->pf_gains,     0, sizeof(frame->pf_gains));
2162         memset(frame->pf_gains_old, 0, sizeof(frame->pf_gains_old));
2163         memset(frame->pf_gains_new, 0, sizeof(frame->pf_gains_new));
2164
2165         frame->deemph_coeff = 0.0;
2166     }
2167     s->seed = 0;
2168
2169     s->flushed = 1;
2170 }
2171
2172 void ff_celt_free(CeltContext **ps)
2173 {
2174     CeltContext *s = *ps;
2175     int i;
2176
2177     if (!s)
2178         return;
2179
2180     for (i = 0; i < FF_ARRAY_ELEMS(s->imdct); i++)
2181         ff_celt_imdct_uninit(&s->imdct[i]);
2182
2183     av_freep(ps);
2184 }
2185
2186 int ff_celt_init(AVCodecContext *avctx, CeltContext **ps, int output_channels)
2187 {
2188     CeltContext *s;
2189     int i, ret;
2190
2191     if (output_channels != 1 && output_channels != 2) {
2192         av_log(avctx, AV_LOG_ERROR, "Invalid number of output channels: %d\n",
2193                output_channels);
2194         return AVERROR(EINVAL);
2195     }
2196
2197     s = av_mallocz(sizeof(*s));
2198     if (!s)
2199         return AVERROR(ENOMEM);
2200
2201     s->avctx           = avctx;
2202     s->output_channels = output_channels;
2203
2204     for (i = 0; i < FF_ARRAY_ELEMS(s->imdct); i++) {
2205         ret = ff_celt_imdct_init(&s->imdct[i], i + 3);
2206         if (ret < 0)
2207             goto fail;
2208     }
2209
2210     avpriv_float_dsp_init(&s->dsp, avctx->flags & CODEC_FLAG_BITEXACT);
2211
2212     ff_celt_flush(s);
2213
2214     *ps = s;
2215
2216     return 0;
2217 fail:
2218     ff_celt_free(&s);
2219     return ret;
2220 }