Fix previous commit by still capping what the encoder will produce
[opus.git] / libcelt / celt.c
1 /* Copyright (c) 2007-2008 CSIRO
2    Copyright (c) 2007-2010 Xiph.Org Foundation
3    Copyright (c) 2008 Gregory Maxwell 
4    Written by Jean-Marc Valin and Gregory Maxwell */
5 /*
6    Redistribution and use in source and binary forms, with or without
7    modification, are permitted provided that the following conditions
8    are met:
9    
10    - Redistributions of source code must retain the above copyright
11    notice, this list of conditions and the following disclaimer.
12    
13    - Redistributions in binary form must reproduce the above copyright
14    notice, this list of conditions and the following disclaimer in the
15    documentation and/or other materials provided with the distribution.
16    
17    - Neither the name of the Xiph.org Foundation nor the names of its
18    contributors may be used to endorse or promote products derived from
19    this software without specific prior written permission.
20    
21    THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
22    ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
23    LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
24    A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR
25    CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
26    EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
27    PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
28    PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
29    LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
30    NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
31    SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
32 */
33
34 #ifdef HAVE_CONFIG_H
35 #include "config.h"
36 #endif
37
38 #define CELT_C
39
40 /* Always enable postfilter for Opus */
41 #if defined(OPUS_BUILD) && !defined(ENABLE_POSTFILTER)
42 #define ENABLE_POSTFILTER
43 #endif
44
45 #include "os_support.h"
46 #include "mdct.h"
47 #include <math.h>
48 #include "celt.h"
49 #include "pitch.h"
50 #include "bands.h"
51 #include "modes.h"
52 #include "entcode.h"
53 #include "quant_bands.h"
54 #include "rate.h"
55 #include "stack_alloc.h"
56 #include "mathops.h"
57 #include "float_cast.h"
58 #include <stdarg.h>
59 #include "plc.h"
60 #include "vq.h"
61
62 static const unsigned char trim_icdf[11] = {126, 124, 119, 109, 87, 41, 19, 9, 4, 2, 0};
63 /* Probs: NONE: 21.875%, LIGHT: 6.25%, NORMAL: 65.625%, AGGRESSIVE: 6.25% */
64 static const unsigned char spread_icdf[4] = {25, 23, 2, 0};
65
66 static const unsigned char tapset_icdf[3]={2,1,0};
67
68 #define COMBFILTER_MAXPERIOD 1024
69 #define COMBFILTER_MINPERIOD 15
70
71 static int resampling_factor(celt_int32 rate)
72 {
73    int ret;
74    switch (rate)
75    {
76    case 48000:
77       ret = 1;
78       break;
79    case 24000:
80       ret = 2;
81       break;
82    case 16000:
83       ret = 3;
84       break;
85    case 12000:
86       ret = 4;
87       break;
88    case 8000:
89       ret = 6;
90       break;
91    default:
92       ret = 0;
93    }
94    return ret;
95 }
96
97 /** Encoder state 
98  @brief Encoder state
99  */
100 struct CELTEncoder {
101    const CELTMode *mode;     /**< Mode used by the encoder */
102    int overlap;
103    int channels;
104    int stream_channels;
105    
106    int force_intra;
107    int clip;
108    int disable_pf;
109    int complexity;
110    int upsample;
111    int start, end;
112
113    celt_int32 bitrate;
114    int vbr;
115    int constrained_vbr;      /* If zero, VBR can do whatever it likes with the rate */
116
117    /* Everything beyond this point gets cleared on a reset */
118 #define ENCODER_RESET_START rng
119
120    ec_uint32 rng;
121    int spread_decision;
122    int delayedIntra;
123    int tonal_average;
124    int lastCodedBands;
125    int hf_average;
126    int tapset_decision;
127
128    int prefilter_period;
129    celt_word16 prefilter_gain;
130    int prefilter_tapset;
131 #ifdef RESYNTH
132    int prefilter_period_old;
133    celt_word16 prefilter_gain_old;
134    int prefilter_tapset_old;
135 #endif
136    int consec_transient;
137
138    /* VBR-related parameters */
139    celt_int32 vbr_reservoir;
140    celt_int32 vbr_drift;
141    celt_int32 vbr_offset;
142    celt_int32 vbr_count;
143
144    celt_word32 preemph_memE[2];
145    celt_word32 preemph_memD[2];
146
147 #ifdef RESYNTH
148    celt_sig syn_mem[2][2*MAX_PERIOD];
149 #endif
150
151    celt_sig in_mem[1]; /* Size = channels*mode->overlap */
152    /* celt_sig prefilter_mem[],  Size = channels*COMBFILTER_PERIOD */
153    /* celt_sig overlap_mem[],  Size = channels*mode->overlap */
154    /* celt_word16 oldEBands[], Size = 2*channels*mode->nbEBands */
155 };
156
157 int celt_encoder_get_size(int channels)
158 {
159    CELTMode *mode = celt_mode_create(48000, 960, NULL);
160    return celt_encoder_get_size_custom(mode, channels);
161 }
162
163 int celt_encoder_get_size_custom(const CELTMode *mode, int channels)
164 {
165    int size = sizeof(struct CELTEncoder)
166          + (2*channels*mode->overlap-1)*sizeof(celt_sig)
167          + channels*COMBFILTER_MAXPERIOD*sizeof(celt_sig)
168          + 3*channels*mode->nbEBands*sizeof(celt_word16);
169    return size;
170 }
171
172 CELTEncoder *celt_encoder_create(int sampling_rate, int channels, int *error)
173 {
174    CELTEncoder *st;
175    st = (CELTEncoder *)celt_alloc(celt_encoder_get_size(channels));
176    if (st!=NULL && celt_encoder_init(st, sampling_rate, channels, error)==NULL)
177    {
178       celt_encoder_destroy(st);
179       st = NULL;
180    }
181    return st;
182 }
183
184 CELTEncoder *celt_encoder_create_custom(const CELTMode *mode, int channels, int *error)
185 {
186    CELTEncoder *st = (CELTEncoder *)celt_alloc(celt_encoder_get_size_custom(mode, channels));
187    if (st!=NULL && celt_encoder_init_custom(st, mode, channels, error)==NULL)
188    {
189       celt_encoder_destroy(st);
190       st = NULL;
191    }
192    return st;
193 }
194
195 CELTEncoder *celt_encoder_init(CELTEncoder *st, int sampling_rate, int channels, int *error)
196 {
197    celt_encoder_init_custom(st, celt_mode_create(48000, 960, NULL), channels, error);
198    st->upsample = resampling_factor(sampling_rate);
199    if (st->upsample==0)
200    {
201       if (error)
202          *error = CELT_BAD_ARG;
203       return NULL;
204    }
205    return st;
206 }
207
208 CELTEncoder *celt_encoder_init_custom(CELTEncoder *st, const CELTMode *mode, int channels, int *error)
209 {
210    if (channels < 0 || channels > 2)
211    {
212       if (error)
213          *error = CELT_BAD_ARG;
214       return NULL;
215    }
216
217    if (st==NULL || mode==NULL)
218    {
219       if (error)
220          *error = CELT_ALLOC_FAIL;
221       return NULL;
222    }
223
224    CELT_MEMSET((char*)st, 0, celt_encoder_get_size_custom(mode, channels));
225    
226    st->mode = mode;
227    st->overlap = mode->overlap;
228    st->stream_channels = st->channels = channels;
229
230    st->upsample = 1;
231    st->start = 0;
232    st->end = st->mode->effEBands;
233    st->constrained_vbr = 1;
234    st->clip = 1;
235
236    st->bitrate = 255000*channels;
237    st->vbr = 0;
238    st->vbr_offset = 0;
239    st->force_intra  = 0;
240    st->delayedIntra = 1;
241    st->tonal_average = 256;
242    st->spread_decision = SPREAD_NORMAL;
243    st->hf_average = 0;
244    st->tapset_decision = 0;
245    st->complexity = 5;
246
247    if (error)
248       *error = CELT_OK;
249    return st;
250 }
251
252 void celt_encoder_destroy(CELTEncoder *st)
253 {
254    celt_free(st);
255 }
256
257 static inline celt_int16 FLOAT2INT16(float x)
258 {
259    x = x*CELT_SIG_SCALE;
260    x = MAX32(x, -32768);
261    x = MIN32(x, 32767);
262    return (celt_int16)float2int(x);
263 }
264
265 static inline celt_word16 SIG2WORD16(celt_sig x)
266 {
267 #ifdef FIXED_POINT
268    x = PSHR32(x, SIG_SHIFT);
269    x = MAX32(x, -32768);
270    x = MIN32(x, 32767);
271    return EXTRACT16(x);
272 #else
273    return (celt_word16)x;
274 #endif
275 }
276
277 static int transient_analysis(const celt_word32 * restrict in, int len, int C,
278                               int overlap)
279 {
280    int i;
281    VARDECL(celt_word16, tmp);
282    celt_word32 mem0=0,mem1=0;
283    int is_transient = 0;
284    int block;
285    int N;
286    /* FIXME: Make that smaller */
287    celt_word16 bins[50];
288    SAVE_STACK;
289    ALLOC(tmp, len, celt_word16);
290
291    block = overlap/2;
292    N=len/block;
293    if (C==1)
294    {
295       for (i=0;i<len;i++)
296          tmp[i] = SHR32(in[i],SIG_SHIFT);
297    } else {
298       for (i=0;i<len;i++)
299          tmp[i] = SHR32(ADD32(in[i],in[i+len]), SIG_SHIFT+1);
300    }
301
302    /* High-pass filter: (1 - 2*z^-1 + z^-2) / (1 - z^-1 + .5*z^-2) */
303    for (i=0;i<len;i++)
304    {
305       celt_word32 x,y;
306       x = tmp[i];
307       y = ADD32(mem0, x);
308 #ifdef FIXED_POINT
309       mem0 = mem1 + y - SHL32(x,1);
310       mem1 = x - SHR32(y,1);
311 #else
312       mem0 = mem1 + y - 2*x;
313       mem1 = x - .5f*y;
314 #endif
315       tmp[i] = EXTRACT16(SHR(y,2));
316    }
317    /* First few samples are bad because we don't propagate the memory */
318    for (i=0;i<12;i++)
319       tmp[i] = 0;
320
321    for (i=0;i<N;i++)
322    {
323       int j;
324       float max_abs=0;
325       for (j=0;j<block;j++)
326          max_abs = MAX32(max_abs, tmp[i*block+j]);
327       bins[i] = max_abs;
328    }
329    for (i=0;i<N;i++)
330    {
331       int j;
332       int conseq=0;
333       celt_word16 t1, t2, t3;
334
335       t1 = MULT16_16_Q15(QCONST16(.15f, 15), bins[i]);
336       t2 = MULT16_16_Q15(QCONST16(.4f, 15), bins[i]);
337       t3 = MULT16_16_Q15(QCONST16(.15f, 15), bins[i]);
338       for (j=0;j<i;j++)
339       {
340          if (bins[j] < t1)
341             conseq++;
342          if (bins[j] < t2)
343             conseq++;
344          else
345             conseq = 0;
346       }
347       if (conseq>=3)
348          is_transient=1;
349       conseq = 0;
350       for (j=i+1;j<N;j++)
351       {
352          if (bins[j] < t3)
353             conseq++;
354          else
355             conseq = 0;
356       }
357       if (conseq>=7)
358          is_transient=1;
359    }
360    RESTORE_STACK;
361    return is_transient;
362 }
363
364 /** Apply window and compute the MDCT for all sub-frames and 
365     all channels in a frame */
366 static void compute_mdcts(const CELTMode *mode, int shortBlocks, celt_sig * restrict in, celt_sig * restrict out, int _C, int LM)
367 {
368    const int C = CHANNELS(_C);
369    if (C==1 && !shortBlocks)
370    {
371       const int overlap = OVERLAP(mode);
372       clt_mdct_forward(&mode->mdct, in, out, mode->window, overlap, mode->maxLM-LM);
373    } else {
374       const int overlap = OVERLAP(mode);
375       int N = mode->shortMdctSize<<LM;
376       int B = 1;
377       int b, c;
378       VARDECL(celt_word32, tmp);
379       SAVE_STACK;
380       if (shortBlocks)
381       {
382          /*lookup = &mode->mdct[0];*/
383          N = mode->shortMdctSize;
384          B = shortBlocks;
385       }
386       ALLOC(tmp, N, celt_word32);
387       c=0; do {
388          for (b=0;b<B;b++)
389          {
390             int j;
391             clt_mdct_forward(&mode->mdct, in+c*(B*N+overlap)+b*N, tmp, mode->window, overlap, shortBlocks ? mode->maxLM : mode->maxLM-LM);
392             /* Interleaving the sub-frames */
393             for (j=0;j<N;j++)
394                out[(j*B+b)+c*N*B] = tmp[j];
395          }
396       } while (++c<C);
397       RESTORE_STACK;
398    }
399 }
400
401 /** Compute the IMDCT and apply window for all sub-frames and 
402     all channels in a frame */
403 static void compute_inv_mdcts(const CELTMode *mode, int shortBlocks, celt_sig *X,
404       celt_sig * restrict out_mem[],
405       celt_sig * restrict overlap_mem[], int _C, int LM)
406 {
407    int c;
408    const int C = CHANNELS(_C);
409    const int N = mode->shortMdctSize<<LM;
410    const int overlap = OVERLAP(mode);
411    c=0; do {
412       int j;
413          VARDECL(celt_word32, x);
414          VARDECL(celt_word32, tmp);
415          int b;
416          int N2 = N;
417          int B = 1;
418          SAVE_STACK;
419          
420          ALLOC(x, N+overlap, celt_word32);
421          ALLOC(tmp, N, celt_word32);
422
423          if (shortBlocks)
424          {
425             N2 = mode->shortMdctSize;
426             B = shortBlocks;
427          }
428          /* Prevents problems from the imdct doing the overlap-add */
429          CELT_MEMSET(x, 0, overlap);
430
431          for (b=0;b<B;b++)
432          {
433             /* De-interleaving the sub-frames */
434             for (j=0;j<N2;j++)
435                tmp[j] = X[(j*B+b)+c*N2*B];
436             clt_mdct_backward(&mode->mdct, tmp, x+N2*b, mode->window, overlap, shortBlocks ? mode->maxLM : mode->maxLM-LM);
437          }
438
439          for (j=0;j<overlap;j++)
440             out_mem[c][j] = x[j] + overlap_mem[c][j];
441          for (;j<N;j++)
442             out_mem[c][j] = x[j];
443          for (j=0;j<overlap;j++)
444             overlap_mem[c][j] = x[N+j];
445          RESTORE_STACK;
446    } while (++c<C);
447 }
448
449 static void deemphasis(celt_sig *in[], celt_word16 *pcm, int N, int _C, int downsample, const celt_word16 *coef, celt_sig *mem)
450 {
451    const int C = CHANNELS(_C);
452    int c;
453    int count=0;
454    c=0; do {
455       int j;
456       celt_sig * restrict x;
457       celt_word16  * restrict y;
458       celt_sig m = mem[c];
459       x =in[c];
460       y = pcm+c;
461       for (j=0;j<N;j++)
462       {
463          celt_sig tmp = *x + m;
464          m = MULT16_32_Q15(coef[0], tmp)
465            - MULT16_32_Q15(coef[1], *x);
466          tmp = SHL32(MULT16_32_Q15(coef[3], tmp), 2);
467          x++;
468          /* Technically the store could be moved outside of the if because
469             the stores we don't want will just be overwritten */
470          if (++count==downsample)
471          {
472             *y = SCALEOUT(SIG2WORD16(tmp));
473             y+=C;
474             count=0;
475          }
476       }
477       mem[c] = m;
478    } while (++c<C);
479 }
480
481 #ifdef ENABLE_POSTFILTER
482 static void comb_filter(celt_word32 *y, celt_word32 *x, int T0, int T1, int N,
483       celt_word16 g0, celt_word16 g1, int tapset0, int tapset1,
484       const celt_word16 *window, int overlap)
485 {
486    int i;
487    /* printf ("%d %d %f %f\n", T0, T1, g0, g1); */
488    celt_word16 g00, g01, g02, g10, g11, g12;
489    static const celt_word16 gains[3][3] = {
490          {QCONST16(0.3066406250f, 15), QCONST16(0.2170410156f, 15), QCONST16(0.1296386719f, 15)},
491          {QCONST16(0.4638671875f, 15), QCONST16(0.2680664062f, 15), QCONST16(0.f, 15)},
492          {QCONST16(0.7998046875f, 15), QCONST16(0.1000976562f, 15), QCONST16(0.f, 15)}};
493    g00 = MULT16_16_Q15(g0, gains[tapset0][0]);
494    g01 = MULT16_16_Q15(g0, gains[tapset0][1]);
495    g02 = MULT16_16_Q15(g0, gains[tapset0][2]);
496    g10 = MULT16_16_Q15(g1, gains[tapset1][0]);
497    g11 = MULT16_16_Q15(g1, gains[tapset1][1]);
498    g12 = MULT16_16_Q15(g1, gains[tapset1][2]);
499    for (i=0;i<overlap;i++)
500    {
501       celt_word16 f;
502       f = MULT16_16_Q15(window[i],window[i]);
503       y[i] = x[i]
504                + MULT16_32_Q15(MULT16_16_Q15((Q15ONE-f),g00),x[i-T0])
505                + MULT16_32_Q15(MULT16_16_Q15((Q15ONE-f),g01),x[i-T0-1])
506                + MULT16_32_Q15(MULT16_16_Q15((Q15ONE-f),g01),x[i-T0+1])
507                + MULT16_32_Q15(MULT16_16_Q15((Q15ONE-f),g02),x[i-T0-2])
508                + MULT16_32_Q15(MULT16_16_Q15((Q15ONE-f),g02),x[i-T0+2])
509                + MULT16_32_Q15(MULT16_16_Q15(f,g10),x[i-T1])
510                + MULT16_32_Q15(MULT16_16_Q15(f,g11),x[i-T1-1])
511                + MULT16_32_Q15(MULT16_16_Q15(f,g11),x[i-T1+1])
512                + MULT16_32_Q15(MULT16_16_Q15(f,g12),x[i-T1-2])
513                + MULT16_32_Q15(MULT16_16_Q15(f,g12),x[i-T1+2]);
514
515    }
516    for (i=overlap;i<N;i++)
517       y[i] = x[i]
518                + MULT16_32_Q15(g10,x[i-T1])
519                + MULT16_32_Q15(g11,x[i-T1-1])
520                + MULT16_32_Q15(g11,x[i-T1+1])
521                + MULT16_32_Q15(g12,x[i-T1-2])
522                + MULT16_32_Q15(g12,x[i-T1+2]);
523 }
524 #endif /* ENABLE_POSTFILTER */
525
526 static const signed char tf_select_table[4][8] = {
527       {0, -1, 0, -1,    0,-1, 0,-1},
528       {0, -1, 0, -2,    1, 0, 1,-1},
529       {0, -2, 0, -3,    2, 0, 1,-1},
530       {0, -2, 0, -3,    3, 0, 1,-1},
531 };
532
533 static celt_word32 l1_metric(const celt_norm *tmp, int N, int LM, int width)
534 {
535    int i, j;
536    static const celt_word16 sqrtM_1[4] = {Q15ONE, QCONST16(.70710678f,15), QCONST16(0.5f,15), QCONST16(0.35355339f,15)};
537    celt_word32 L1;
538    celt_word16 bias;
539    L1=0;
540    for (i=0;i<1<<LM;i++)
541    {
542       celt_word32 L2 = 0;
543       for (j=0;j<N>>LM;j++)
544          L2 = MAC16_16(L2, tmp[(j<<LM)+i], tmp[(j<<LM)+i]);
545       L1 += celt_sqrt(L2);
546    }
547    L1 = MULT16_32_Q15(sqrtM_1[LM], L1);
548    if (width==1)
549       bias = QCONST16(.12f,15)*LM;
550    else if (width==2)
551       bias = QCONST16(.05f,15)*LM;
552    else
553       bias = QCONST16(.02f,15)*LM;
554    L1 = MAC16_32_Q15(L1, bias, L1);
555    return L1;
556 }
557
558 static int tf_analysis(const CELTMode *m, celt_word16 *bandLogE, celt_word16 *oldBandE,
559       int len, int C, int isTransient, int *tf_res, int nbCompressedBytes, celt_norm *X,
560       int N0, int LM, int *tf_sum)
561 {
562    int i;
563    VARDECL(int, metric);
564    int cost0;
565    int cost1;
566    VARDECL(int, path0);
567    VARDECL(int, path1);
568    VARDECL(celt_norm, tmp);
569    int lambda;
570    int tf_select=0;
571    SAVE_STACK;
572
573    if (nbCompressedBytes<15*C)
574    {
575       *tf_sum = 0;
576       for (i=0;i<len;i++)
577          tf_res[i] = isTransient;
578       return 0;
579    }
580    if (nbCompressedBytes<40)
581       lambda = 12;
582    else if (nbCompressedBytes<60)
583       lambda = 6;
584    else if (nbCompressedBytes<100)
585       lambda = 4;
586    else
587       lambda = 3;
588
589    ALLOC(metric, len, int);
590    ALLOC(tmp, (m->eBands[len]-m->eBands[len-1])<<LM, celt_norm);
591    ALLOC(path0, len, int);
592    ALLOC(path1, len, int);
593
594    *tf_sum = 0;
595    for (i=0;i<len;i++)
596    {
597       int j, k, N;
598       celt_word32 L1, best_L1;
599       int best_level=0;
600       N = (m->eBands[i+1]-m->eBands[i])<<LM;
601       for (j=0;j<N;j++)
602          tmp[j] = X[j+(m->eBands[i]<<LM)];
603       /* Just add the right channel if we're in stereo */
604       if (C==2)
605          for (j=0;j<N;j++)
606             tmp[j] = ADD16(tmp[j],X[N0+j+(m->eBands[i]<<LM)]);
607       L1 = l1_metric(tmp, N, isTransient ? LM : 0, N>>LM);
608       best_L1 = L1;
609       /*printf ("%f ", L1);*/
610       for (k=0;k<LM;k++)
611       {
612          int B;
613
614          if (isTransient)
615             B = (LM-k-1);
616          else
617             B = k+1;
618
619          if (isTransient)
620             haar1(tmp, N>>(LM-k), 1<<(LM-k));
621          else
622             haar1(tmp, N>>k, 1<<k);
623
624          L1 = l1_metric(tmp, N, B, N>>LM);
625
626          if (L1 < best_L1)
627          {
628             best_L1 = L1;
629             best_level = k+1;
630          }
631       }
632       /*printf ("%d ", isTransient ? LM-best_level : best_level);*/
633       if (isTransient)
634          metric[i] = best_level;
635       else
636          metric[i] = -best_level;
637       *tf_sum += metric[i];
638    }
639    /*printf("\n");*/
640    /* FIXME: Figure out how to set this */
641    tf_select = 0;
642
643    cost0 = 0;
644    cost1 = isTransient ? 0 : lambda;
645    /* Viterbi forward pass */
646    for (i=1;i<len;i++)
647    {
648       int curr0, curr1;
649       int from0, from1;
650
651       from0 = cost0;
652       from1 = cost1 + lambda;
653       if (from0 < from1)
654       {
655          curr0 = from0;
656          path0[i]= 0;
657       } else {
658          curr0 = from1;
659          path0[i]= 1;
660       }
661
662       from0 = cost0 + lambda;
663       from1 = cost1;
664       if (from0 < from1)
665       {
666          curr1 = from0;
667          path1[i]= 0;
668       } else {
669          curr1 = from1;
670          path1[i]= 1;
671       }
672       cost0 = curr0 + abs(metric[i]-tf_select_table[LM][4*isTransient+2*tf_select+0]);
673       cost1 = curr1 + abs(metric[i]-tf_select_table[LM][4*isTransient+2*tf_select+1]);
674    }
675    tf_res[len-1] = cost0 < cost1 ? 0 : 1;
676    /* Viterbi backward pass to check the decisions */
677    for (i=len-2;i>=0;i--)
678    {
679       if (tf_res[i+1] == 1)
680          tf_res[i] = path1[i+1];
681       else
682          tf_res[i] = path0[i+1];
683    }
684    RESTORE_STACK;
685    return tf_select;
686 }
687
688 static void tf_encode(int start, int end, int isTransient, int *tf_res, int LM, int tf_select, ec_enc *enc)
689 {
690    int curr, i;
691    int tf_select_rsv;
692    int tf_changed;
693    int logp;
694    ec_uint32 budget;
695    ec_uint32 tell;
696    budget = enc->storage*8;
697    tell = ec_tell(enc);
698    logp = isTransient ? 2 : 4;
699    /* Reserve space to code the tf_select decision. */
700    tf_select_rsv = LM>0 && tell+logp+1 <= budget;
701    budget -= tf_select_rsv;
702    curr = tf_changed = 0;
703    for (i=start;i<end;i++)
704    {
705       if (tell+logp<=budget)
706       {
707          ec_enc_bit_logp(enc, tf_res[i] ^ curr, logp);
708          tell = ec_tell(enc);
709          curr = tf_res[i];
710          tf_changed |= curr;
711       }
712       else
713          tf_res[i] = curr;
714       logp = isTransient ? 4 : 5;
715    }
716    /* Only code tf_select if it would actually make a difference. */
717    if (tf_select_rsv &&
718          tf_select_table[LM][4*isTransient+0+tf_changed]!=
719          tf_select_table[LM][4*isTransient+2+tf_changed])
720       ec_enc_bit_logp(enc, tf_select, 1);
721    else
722       tf_select = 0;
723    for (i=start;i<end;i++)
724       tf_res[i] = tf_select_table[LM][4*isTransient+2*tf_select+tf_res[i]];
725    /*printf("%d %d ", isTransient, tf_select); for(i=0;i<end;i++)printf("%d ", tf_res[i]);printf("\n");*/
726 }
727
728 static void tf_decode(int start, int end, int isTransient, int *tf_res, int LM, ec_dec *dec)
729 {
730    int i, curr, tf_select;
731    int tf_select_rsv;
732    int tf_changed;
733    int logp;
734    ec_uint32 budget;
735    ec_uint32 tell;
736
737    budget = dec->storage*8;
738    tell = ec_tell(dec);
739    logp = isTransient ? 2 : 4;
740    tf_select_rsv = LM>0 && tell+logp+1<=budget;
741    budget -= tf_select_rsv;
742    tf_changed = curr = 0;
743    for (i=start;i<end;i++)
744    {
745       if (tell+logp<=budget)
746       {
747          curr ^= ec_dec_bit_logp(dec, logp);
748          tell = ec_tell(dec);
749          tf_changed |= curr;
750       }
751       tf_res[i] = curr;
752       logp = isTransient ? 4 : 5;
753    }
754    tf_select = 0;
755    if (tf_select_rsv &&
756      tf_select_table[LM][4*isTransient+0+tf_changed] !=
757      tf_select_table[LM][4*isTransient+2+tf_changed])
758    {
759       tf_select = ec_dec_bit_logp(dec, 1);
760    }
761    for (i=start;i<end;i++)
762    {
763       tf_res[i] = tf_select_table[LM][4*isTransient+2*tf_select+tf_res[i]];
764    }
765 }
766
767 static void init_caps(const CELTMode *m,int *cap,int LM,int C)
768 {
769    int i;
770    for (i=0;i<m->nbEBands;i++)
771    {
772       int N;
773       N=(m->eBands[i+1]-m->eBands[i])<<LM;
774       cap[i] = (m->cache.caps[m->nbEBands*(2*LM+C-1)+i]+64)*C*N>>2;
775    }
776 }
777
778 static int alloc_trim_analysis(const CELTMode *m, const celt_norm *X,
779       const celt_word16 *bandLogE, int nbEBands, int LM, int C, int N0)
780 {
781    int i;
782    celt_word32 diff=0;
783    int c;
784    int trim_index = 5;
785    if (C==2)
786    {
787       celt_word16 sum = 0; /* Q10 */
788       /* Compute inter-channel correlation for low frequencies */
789       for (i=0;i<8;i++)
790       {
791          int j;
792          celt_word32 partial = 0;
793          for (j=m->eBands[i]<<LM;j<m->eBands[i+1]<<LM;j++)
794             partial = MAC16_16(partial, X[j], X[N0+j]);
795          sum = ADD16(sum, EXTRACT16(SHR32(partial, 18)));
796       }
797       sum = MULT16_16_Q15(QCONST16(1.f/8, 15), sum);
798       /*printf ("%f\n", sum);*/
799       if (sum > QCONST16(.995f,10))
800          trim_index-=4;
801       else if (sum > QCONST16(.92f,10))
802          trim_index-=3;
803       else if (sum > QCONST16(.85f,10))
804          trim_index-=2;
805       else if (sum > QCONST16(.8f,10))
806          trim_index-=1;
807    }
808
809    /* Estimate spectral tilt */
810    c=0; do {
811       for (i=0;i<nbEBands-1;i++)
812       {
813          diff += bandLogE[i+c*nbEBands]*(celt_int32)(2+2*i-nbEBands);
814       }
815    } while (++c<0);
816    diff /= C*(nbEBands-1);
817    /*printf("%f\n", diff);*/
818    if (diff > QCONST16(2.f, DB_SHIFT))
819       trim_index--;
820    if (diff > QCONST16(8.f, DB_SHIFT))
821       trim_index--;
822    if (diff < -QCONST16(4.f, DB_SHIFT))
823       trim_index++;
824    if (diff < -QCONST16(10.f, DB_SHIFT))
825       trim_index++;
826
827    if (trim_index<0)
828       trim_index = 0;
829    if (trim_index>10)
830       trim_index = 10;
831    return trim_index;
832 }
833
834 static int stereo_analysis(const CELTMode *m, const celt_norm *X,
835       int LM, int N0)
836 {
837    int i;
838    int thetas;
839    celt_word32 sumLR = EPSILON, sumMS = EPSILON;
840
841    /* Use the L1 norm to model the entropy of the L/R signal vs the M/S signal */
842    for (i=0;i<13;i++)
843    {
844       int j;
845       for (j=m->eBands[i]<<LM;j<m->eBands[i+1]<<LM;j++)
846       {
847          celt_word16 L, R, M, S;
848          L = X[j];
849          R = X[N0+j];
850          M = L+R;
851          S = L-R;
852          sumLR += EXTEND32(ABS16(L)) + EXTEND32(ABS16(R));
853          sumMS += EXTEND32(ABS16(M)) + EXTEND32(ABS16(S));
854       }
855    }
856    sumMS = MULT16_32_Q15(QCONST16(0.707107f, 15), sumMS);
857    thetas = 13;
858    /* We don't need thetas for lower bands with LM<=1 */
859    if (LM<=1)
860       thetas -= 8;
861    return MULT16_32_Q15((m->eBands[13]<<(LM+1))+thetas, sumMS)
862          > MULT16_32_Q15(m->eBands[13]<<(LM+1), sumLR);
863 }
864
865 #ifdef FIXED_POINT
866 CELT_STATIC
867 int celt_encode_with_ec(CELTEncoder * restrict st, const celt_int16 * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes, ec_enc *enc)
868 {
869 #else
870 CELT_STATIC
871 int celt_encode_with_ec_float(CELTEncoder * restrict st, const celt_sig * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes, ec_enc *enc)
872 {
873 #endif
874    int i, c, N;
875    int bits;
876    ec_enc _enc;
877    VARDECL(celt_sig, in);
878    VARDECL(celt_sig, freq);
879    VARDECL(celt_norm, X);
880    VARDECL(celt_ener, bandE);
881    VARDECL(celt_word16, bandLogE);
882    VARDECL(int, fine_quant);
883    VARDECL(celt_word16, error);
884    VARDECL(int, pulses);
885    VARDECL(int, cap);
886    VARDECL(int, offsets);
887    VARDECL(int, fine_priority);
888    VARDECL(int, tf_res);
889    VARDECL(unsigned char, collapse_masks);
890    celt_sig *_overlap_mem;
891    celt_sig *prefilter_mem;
892    celt_word16 *oldBandE, *oldLogE, *oldLogE2;
893    int shortBlocks=0;
894    int isTransient=0;
895    int resynth;
896    const int CC = CHANNELS(st->channels);
897    const int C = CHANNELS(st->stream_channels);
898    int LM, M;
899    int tf_select;
900    int nbFilledBytes, nbAvailableBytes;
901    int effEnd;
902    int codedBands;
903    int tf_sum;
904    int alloc_trim;
905    int pitch_index=COMBFILTER_MINPERIOD;
906    celt_word16 gain1 = 0;
907    int intensity=0;
908    int dual_stereo=0;
909    int effectiveBytes;
910    celt_word16 pf_threshold;
911    int dynalloc_logp;
912    celt_int32 vbr_rate;
913    celt_int32 total_bits;
914    celt_int32 total_boost;
915    celt_int32 balance;
916    celt_int32 tell;
917    int prefilter_tapset=0;
918    int pf_on;
919    int anti_collapse_rsv;
920    int anti_collapse_on=0;
921    int silence=0;
922    SAVE_STACK;
923
924    if (nbCompressedBytes<2 || pcm==NULL)
925      return CELT_BAD_ARG;
926
927    /* Can't produce more than 1275 output bytes */
928    nbCompressedBytes = IMIN(nbCompressedBytes,1275);
929    frame_size *= st->upsample;
930    for (LM=0;LM<=st->mode->maxLM;LM++)
931       if (st->mode->shortMdctSize<<LM==frame_size)
932          break;
933    if (LM>st->mode->maxLM)
934       return CELT_BAD_ARG;
935    M=1<<LM;
936    N = M*st->mode->shortMdctSize;
937
938    prefilter_mem = st->in_mem+CC*(st->overlap);
939    _overlap_mem = prefilter_mem+CC*COMBFILTER_MAXPERIOD;
940    /*_overlap_mem = st->in_mem+C*(st->overlap);*/
941    oldBandE = (celt_word16*)(st->in_mem+CC*(2*st->overlap+COMBFILTER_MAXPERIOD));
942    oldLogE = oldBandE + CC*st->mode->nbEBands;
943    oldLogE2 = oldLogE + CC*st->mode->nbEBands;
944
945    if (enc==NULL)
946    {
947       tell=1;
948       nbFilledBytes=0;
949    } else {
950       tell=ec_tell(enc);
951       nbFilledBytes=(tell+4)>>3;
952    }
953    nbAvailableBytes = nbCompressedBytes - nbFilledBytes;
954
955    if (st->vbr)
956    {
957       celt_int32 den=st->mode->Fs>>BITRES;
958       vbr_rate=(st->bitrate*frame_size+(den>>1))/den;
959       effectiveBytes = vbr_rate>>3;
960    } else {
961       celt_int32 tmp;
962       vbr_rate = 0;
963       tmp = st->bitrate*frame_size;
964       if (tell>1)
965          tmp += tell;
966       nbCompressedBytes = IMAX(2, IMIN(nbCompressedBytes,
967             (tmp+4*st->mode->Fs)/(8*st->mode->Fs)));
968       effectiveBytes = nbCompressedBytes;
969    }
970
971    if (enc==NULL)
972    {
973       ec_enc_init(&_enc, compressed, nbCompressedBytes);
974       enc = &_enc;
975    }
976
977    if (vbr_rate>0)
978    {
979       /* Computes the max bit-rate allowed in VBR mode to avoid violating the
980           target rate and buffering.
981          We must do this up front so that bust-prevention logic triggers
982           correctly if we don't have enough bits. */
983       if (st->constrained_vbr)
984       {
985          celt_int32 vbr_bound;
986          celt_int32 max_allowed;
987          /* We could use any multiple of vbr_rate as bound (depending on the
988              delay).
989             This is clamped to ensure we use at least two bytes if the encoder
990              was entirely empty, but to allow 0 in hybrid mode. */
991          vbr_bound = vbr_rate;
992          max_allowed = IMIN(IMAX(tell==1?2:0,
993                vbr_rate+vbr_bound-st->vbr_reservoir>>(BITRES+3)),
994                nbAvailableBytes);
995          if(max_allowed < nbAvailableBytes)
996          {
997             nbCompressedBytes = nbFilledBytes+max_allowed;
998             nbAvailableBytes = max_allowed;
999             ec_enc_shrink(enc, nbCompressedBytes);
1000          }
1001       }
1002    }
1003    total_bits = nbCompressedBytes*8;
1004
1005    effEnd = st->end;
1006    if (effEnd > st->mode->effEBands)
1007       effEnd = st->mode->effEBands;
1008
1009    ALLOC(in, CC*(N+st->overlap), celt_sig);
1010
1011    /* Find pitch period and gain */
1012    {
1013       VARDECL(celt_sig, _pre);
1014       celt_sig *pre[2];
1015       SAVE_STACK;
1016       ALLOC(_pre, CC*(N+COMBFILTER_MAXPERIOD), celt_sig);
1017
1018       pre[0] = _pre;
1019       pre[1] = _pre + (N+COMBFILTER_MAXPERIOD);
1020
1021       silence = 1;
1022       c=0; do {
1023          int count = 0;
1024          const celt_word16 * restrict pcmp = pcm+c;
1025          celt_sig * restrict inp = in+c*(N+st->overlap)+st->overlap;
1026
1027          for (i=0;i<N;i++)
1028          {
1029             celt_sig x, tmp;
1030
1031             x = SCALEIN(*pcmp);
1032 #ifndef FIXED_POINT
1033             if (st->clip)
1034                x = MAX32(-65536.f, MIN32(65536.f,x));
1035 #endif
1036             if (++count==st->upsample)
1037             {
1038                count=0;
1039                pcmp+=CC;
1040             } else {
1041                x = 0;
1042             }
1043             /* Apply pre-emphasis */
1044             tmp = MULT16_16(st->mode->preemph[2], x);
1045             *inp = tmp + st->preemph_memE[c];
1046             st->preemph_memE[c] = MULT16_32_Q15(st->mode->preemph[1], *inp)
1047                                    - MULT16_32_Q15(st->mode->preemph[0], tmp);
1048             silence = silence && *inp == 0;
1049             inp++;
1050          }
1051          CELT_COPY(pre[c], prefilter_mem+c*COMBFILTER_MAXPERIOD, COMBFILTER_MAXPERIOD);
1052          CELT_COPY(pre[c]+COMBFILTER_MAXPERIOD, in+c*(N+st->overlap)+st->overlap, N);
1053       } while (++c<CC);
1054
1055       if (tell==1)
1056          ec_enc_bit_logp(enc, silence, 15);
1057       else
1058          silence=0;
1059       if (silence)
1060       {
1061          /*In VBR mode there is no need to send more than the minimum. */
1062          if (vbr_rate>0)
1063          {
1064             effectiveBytes=nbCompressedBytes=IMIN(nbCompressedBytes, nbFilledBytes+2);
1065             total_bits=nbCompressedBytes*8;
1066             nbAvailableBytes=2;
1067             ec_enc_shrink(enc, nbCompressedBytes);
1068          }
1069          /* Pretend we've filled all the remaining bits with zeros
1070             (that's what the initialiser did anyway) */
1071          tell = nbCompressedBytes*8;
1072          enc->nbits_total+=tell-ec_tell(enc);
1073       }
1074 #ifdef ENABLE_POSTFILTER
1075       if (nbAvailableBytes>12*C && st->start==0 && !silence && !st->disable_pf && st->complexity >= 5)
1076       {
1077          VARDECL(celt_word16, pitch_buf);
1078          ALLOC(pitch_buf, (COMBFILTER_MAXPERIOD+N)>>1, celt_word16);
1079
1080          pitch_downsample(pre, pitch_buf, COMBFILTER_MAXPERIOD+N, CC);
1081          pitch_search(pitch_buf+(COMBFILTER_MAXPERIOD>>1), pitch_buf, N,
1082                COMBFILTER_MAXPERIOD-COMBFILTER_MINPERIOD, &pitch_index);
1083          pitch_index = COMBFILTER_MAXPERIOD-pitch_index;
1084
1085          gain1 = remove_doubling(pitch_buf, COMBFILTER_MAXPERIOD, COMBFILTER_MINPERIOD,
1086                N, &pitch_index, st->prefilter_period, st->prefilter_gain);
1087          if (pitch_index > COMBFILTER_MAXPERIOD-2)
1088             pitch_index = COMBFILTER_MAXPERIOD-2;
1089          gain1 = MULT16_16_Q15(QCONST16(.7f,15),gain1);
1090          prefilter_tapset = st->tapset_decision;
1091       } else {
1092          gain1 = 0;
1093       }
1094
1095       /* Gain threshold for enabling the prefilter/postfilter */
1096       pf_threshold = QCONST16(.2f,15);
1097
1098       /* Adjusting the threshold based on rate and continuity */
1099       if (abs(pitch_index-st->prefilter_period)*10>pitch_index)
1100          pf_threshold += QCONST16(.2f,15);
1101       if (nbAvailableBytes<25)
1102          pf_threshold += QCONST16(.1f,15);
1103       if (nbAvailableBytes<35)
1104          pf_threshold += QCONST16(.1f,15);
1105       if (st->prefilter_gain > QCONST16(.4f,15))
1106          pf_threshold -= QCONST16(.1f,15);
1107       if (st->prefilter_gain > QCONST16(.55f,15))
1108          pf_threshold -= QCONST16(.1f,15);
1109
1110       /* Hard threshold at 0.2 */
1111       pf_threshold = MAX16(pf_threshold, QCONST16(.2f,15));
1112       if (gain1<pf_threshold)
1113       {
1114          if(st->start==0 && tell+16<=total_bits)
1115             ec_enc_bit_logp(enc, 0, 1);
1116          gain1 = 0;
1117          pf_on = 0;
1118       } else {
1119          /*This block is not gated by a total bits check only because
1120            of the nbAvailableBytes check above.*/
1121          int qg;
1122          int octave;
1123
1124          if (ABS16(gain1-st->prefilter_gain)<QCONST16(.1f,15))
1125             gain1=st->prefilter_gain;
1126
1127 #ifdef FIXED_POINT
1128          qg = ((gain1+1536)>>10)/3-1;
1129 #else
1130          qg = floor(.5+gain1*32/3)-1;
1131 #endif
1132          qg = IMAX(0, IMIN(7, qg));
1133          ec_enc_bit_logp(enc, 1, 1);
1134          pitch_index += 1;
1135          octave = EC_ILOG(pitch_index)-5;
1136          ec_enc_uint(enc, octave, 6);
1137          ec_enc_bits(enc, pitch_index-(16<<octave), 4+octave);
1138          pitch_index -= 1;
1139          ec_enc_bits(enc, qg, 3);
1140          if (ec_tell(enc)+2<=total_bits)
1141             ec_enc_icdf(enc, prefilter_tapset, tapset_icdf, 2);
1142          else
1143            prefilter_tapset = 0;
1144          gain1 = QCONST16(0.09375f,15)*(qg+1);
1145          pf_on = 1;
1146       }
1147       /*printf("%d %f\n", pitch_index, gain1);*/
1148 #else /* ENABLE_POSTFILTER */
1149       if(st->start==0 && tell+16<=total_bits)
1150          ec_enc_bit_logp(enc, 0, 1);
1151       pf_on = 0;
1152 #endif /* ENABLE_POSTFILTER */
1153
1154       c=0; do {
1155          st->prefilter_period=IMAX(st->prefilter_period, COMBFILTER_MINPERIOD);
1156          CELT_COPY(in+c*(N+st->overlap), st->in_mem+c*(st->overlap), st->overlap);
1157 #ifdef ENABLE_POSTFILTER
1158          comb_filter(in+c*(N+st->overlap)+st->overlap, pre[c]+COMBFILTER_MAXPERIOD,
1159                st->prefilter_period, pitch_index, N, -st->prefilter_gain, -gain1,
1160                st->prefilter_tapset, prefilter_tapset, st->mode->window, st->mode->overlap);
1161 #endif /* ENABLE_POSTFILTER */
1162          CELT_COPY(st->in_mem+c*(st->overlap), in+c*(N+st->overlap)+N, st->overlap);
1163
1164 #ifdef ENABLE_POSTFILTER
1165          if (N>COMBFILTER_MAXPERIOD)
1166          {
1167             CELT_MOVE(prefilter_mem+c*COMBFILTER_MAXPERIOD, pre[c]+N, COMBFILTER_MAXPERIOD);
1168          } else {
1169             CELT_MOVE(prefilter_mem+c*COMBFILTER_MAXPERIOD, prefilter_mem+c*COMBFILTER_MAXPERIOD+N, COMBFILTER_MAXPERIOD-N);
1170             CELT_MOVE(prefilter_mem+c*COMBFILTER_MAXPERIOD+COMBFILTER_MAXPERIOD-N, pre[c]+COMBFILTER_MAXPERIOD, N);
1171          }
1172 #endif /* ENABLE_POSTFILTER */
1173       } while (++c<CC);
1174
1175       RESTORE_STACK;
1176    }
1177
1178 #ifdef RESYNTH
1179    resynth = 1;
1180 #else
1181    resynth = 0;
1182 #endif
1183
1184    isTransient = 0;
1185    shortBlocks = 0;
1186    if (LM>0 && ec_tell(enc)+3<=total_bits)
1187    {
1188       if (st->complexity > 1)
1189       {
1190          isTransient = transient_analysis(in, N+st->overlap, CC,
1191                   st->overlap);
1192          if (isTransient)
1193             shortBlocks = M;
1194       }
1195       ec_enc_bit_logp(enc, isTransient, 3);
1196    }
1197
1198    ALLOC(freq, CC*N, celt_sig); /**< Interleaved signal MDCTs */
1199    ALLOC(bandE,st->mode->nbEBands*CC, celt_ener);
1200    ALLOC(bandLogE,st->mode->nbEBands*CC, celt_word16);
1201    /* Compute MDCTs */
1202    compute_mdcts(st->mode, shortBlocks, in, freq, CC, LM);
1203
1204    if (CC==2&&C==1)
1205    {
1206       for (i=0;i<N;i++)
1207          freq[i] = ADD32(HALF32(freq[i]), HALF32(freq[N+i]));
1208    }
1209    if (st->upsample != 1)
1210    {
1211       c=0; do
1212       {
1213          int bound = N/st->upsample;
1214          for (i=0;i<bound;i++)
1215             freq[c*N+i] *= st->upsample;
1216          for (;i<N;i++)
1217             freq[c*N+i] = 0;
1218       } while (++c<C);
1219    }
1220    ALLOC(X, C*N, celt_norm);         /**< Interleaved normalised MDCTs */
1221
1222    compute_band_energies(st->mode, freq, bandE, effEnd, C, M);
1223
1224    amp2Log2(st->mode, effEnd, st->end, bandE, bandLogE, C);
1225
1226    /* Band normalisation */
1227    normalise_bands(st->mode, freq, X, bandE, effEnd, C, M);
1228
1229    ALLOC(tf_res, st->mode->nbEBands, int);
1230    /* Needs to be before coarse energy quantization because otherwise the energy gets modified */
1231    tf_select = tf_analysis(st->mode, bandLogE, oldBandE, effEnd, C, isTransient, tf_res, effectiveBytes, X, N, LM, &tf_sum);
1232    for (i=effEnd;i<st->end;i++)
1233       tf_res[i] = tf_res[effEnd-1];
1234
1235    ALLOC(error, C*st->mode->nbEBands, celt_word16);
1236    quant_coarse_energy(st->mode, st->start, st->end, effEnd, bandLogE,
1237          oldBandE, total_bits, error, enc,
1238          C, LM, nbAvailableBytes, st->force_intra,
1239          &st->delayedIntra, st->complexity >= 4);
1240
1241    tf_encode(st->start, st->end, isTransient, tf_res, LM, tf_select, enc);
1242
1243    st->spread_decision = SPREAD_NORMAL;
1244    if (ec_tell(enc)+4<=total_bits)
1245    {
1246       if (shortBlocks || st->complexity < 3 || nbAvailableBytes < 10*C)
1247       {
1248          if (st->complexity == 0)
1249             st->spread_decision = SPREAD_NONE;
1250       } else {
1251          st->spread_decision = spreading_decision(st->mode, X,
1252                &st->tonal_average, st->spread_decision, &st->hf_average,
1253                &st->tapset_decision, pf_on&&!shortBlocks, effEnd, C, M);
1254       }
1255       ec_enc_icdf(enc, st->spread_decision, spread_icdf, 5);
1256    }
1257
1258    ALLOC(cap, st->mode->nbEBands, int);
1259    ALLOC(offsets, st->mode->nbEBands, int);
1260
1261    init_caps(st->mode,cap,LM,C);
1262    for (i=0;i<st->mode->nbEBands;i++)
1263       offsets[i] = 0;
1264    /* Dynamic allocation code */
1265    /* Make sure that dynamic allocation can't make us bust the budget */
1266    if (effectiveBytes > 50 && LM>=1)
1267    {
1268       int t1, t2;
1269       if (LM <= 1)
1270       {
1271          t1 = 3;
1272          t2 = 5;
1273       } else {
1274          t1 = 2;
1275          t2 = 4;
1276       }
1277       for (i=1;i<st->mode->nbEBands-1;i++)
1278       {
1279          celt_word32 d2;
1280          d2 = 2*bandLogE[i]-bandLogE[i-1]-bandLogE[i+1];
1281          if (C==2)
1282             d2 = HALF32(d2 + 2*bandLogE[i+st->mode->nbEBands]-
1283                   bandLogE[i-1+st->mode->nbEBands]-bandLogE[i+1+st->mode->nbEBands]);
1284          if (d2 > SHL16(t1,DB_SHIFT))
1285             offsets[i] += 1;
1286          if (d2 > SHL16(t2,DB_SHIFT))
1287             offsets[i] += 1;
1288       }
1289    }
1290    dynalloc_logp = 6;
1291    total_bits<<=BITRES;
1292    total_boost = 0;
1293    tell = ec_tell_frac(enc);
1294    for (i=st->start;i<st->end;i++)
1295    {
1296       int width, quanta;
1297       int dynalloc_loop_logp;
1298       int boost;
1299       int j;
1300       width = C*(st->mode->eBands[i+1]-st->mode->eBands[i])<<LM;
1301       /* quanta is 6 bits, but no more than 1 bit/sample
1302          and no less than 1/8 bit/sample */
1303       quanta = IMIN(width<<BITRES, IMAX(6<<BITRES, width));
1304       dynalloc_loop_logp = dynalloc_logp;
1305       boost = 0;
1306       for (j = 0; tell+(dynalloc_loop_logp<<BITRES) < total_bits-total_boost
1307             && boost < cap[i]; j++)
1308       {
1309          int flag;
1310          flag = j<offsets[i];
1311          ec_enc_bit_logp(enc, flag, dynalloc_loop_logp);
1312          tell = ec_tell_frac(enc);
1313          if (!flag)
1314             break;
1315          boost += quanta;
1316          total_boost += quanta;
1317          dynalloc_loop_logp = 1;
1318       }
1319       /* Making dynalloc more likely */
1320       if (j)
1321          dynalloc_logp = IMAX(2, dynalloc_logp-1);
1322       offsets[i] = boost;
1323    }
1324    alloc_trim = 5;
1325    if (tell+(6<<BITRES) <= total_bits - total_boost)
1326    {
1327       alloc_trim = alloc_trim_analysis(st->mode, X, bandLogE,
1328             st->mode->nbEBands, LM, C, N);
1329       ec_enc_icdf(enc, alloc_trim, trim_icdf, 7);
1330       tell = ec_tell_frac(enc);
1331    }
1332
1333    /* Variable bitrate */
1334    if (vbr_rate>0)
1335    {
1336      celt_word16 alpha;
1337      celt_int32 delta;
1338      /* The target rate in 8th bits per frame */
1339      celt_int32 target;
1340      celt_int32 min_allowed;
1341
1342      target = vbr_rate + st->vbr_offset - ((40*C+20)<<BITRES);
1343
1344      /* Shortblocks get a large boost in bitrate, but since they
1345         are uncommon long blocks are not greatly affected */
1346      if (shortBlocks || tf_sum < -2*(st->end-st->start))
1347         target = 7*target/4;
1348      else if (tf_sum < -(st->end-st->start))
1349         target = 3*target/2;
1350      else if (M > 1)
1351         target-=(target+14)/28;
1352
1353      /* The current offset is removed from the target and the space used
1354         so far is added*/
1355      target=target+tell;
1356
1357      /* In VBR mode the frame size must not be reduced so much that it would
1358          result in the encoder running out of bits.
1359         The margin of 2 bytes ensures that none of the bust-prevention logic
1360          in the decoder will have triggered so far. */
1361      min_allowed = (tell+total_boost+(1<<BITRES+3)-1>>(BITRES+3)) + 2 - nbFilledBytes;
1362
1363      nbAvailableBytes = target+(1<<(BITRES+2))>>(BITRES+3);
1364      nbAvailableBytes = IMAX(min_allowed,nbAvailableBytes);
1365      nbAvailableBytes = IMIN(nbCompressedBytes,nbAvailableBytes+nbFilledBytes) - nbFilledBytes;
1366
1367      if(silence)
1368      {
1369        nbAvailableBytes = 2;
1370        target = 2*8<<BITRES;
1371      }
1372
1373      /* By how much did we "miss" the target on that frame */
1374      delta = target - vbr_rate;
1375
1376      target=nbAvailableBytes<<(BITRES+3);
1377
1378      if (st->vbr_count < 970)
1379      {
1380         st->vbr_count++;
1381         alpha = celt_rcp(SHL32(EXTEND32(st->vbr_count+20),16));
1382      } else
1383         alpha = QCONST16(.001f,15);
1384      /* How many bits have we used in excess of what we're allowed */
1385      if (st->constrained_vbr)
1386         st->vbr_reservoir += target - vbr_rate;
1387      /*printf ("%d\n", st->vbr_reservoir);*/
1388
1389      /* Compute the offset we need to apply in order to reach the target */
1390      st->vbr_drift += (celt_int32)MULT16_32_Q15(alpha,delta-st->vbr_offset-st->vbr_drift);
1391      st->vbr_offset = -st->vbr_drift;
1392      /*printf ("%d\n", st->vbr_drift);*/
1393
1394      if (st->constrained_vbr && st->vbr_reservoir < 0)
1395      {
1396         /* We're under the min value -- increase rate */
1397         int adjust = (-st->vbr_reservoir)/(8<<BITRES);
1398         /* Unless we're just coding silence */
1399         nbAvailableBytes += silence?0:adjust;
1400         st->vbr_reservoir = 0;
1401         /*printf ("+%d\n", adjust);*/
1402      }
1403      nbCompressedBytes = IMIN(nbCompressedBytes,nbAvailableBytes+nbFilledBytes);
1404      /* This moves the raw bits to take into account the new compressed size */
1405      ec_enc_shrink(enc, nbCompressedBytes);
1406    }
1407    if (C==2)
1408    {
1409       int effectiveRate;
1410
1411       /* Always use MS for 2.5 ms frames until we can do a better analysis */
1412       if (LM!=0)
1413          dual_stereo = stereo_analysis(st->mode, X, LM, N);
1414
1415       /* Account for coarse energy */
1416       effectiveRate = (8*effectiveBytes - 80)>>LM;
1417
1418       /* effectiveRate in kb/s */
1419       effectiveRate = 2*effectiveRate/5;
1420       if (effectiveRate<35)
1421          intensity = 8;
1422       else if (effectiveRate<50)
1423          intensity = 12;
1424       else if (effectiveRate<68)
1425          intensity = 16;
1426       else if (effectiveRate<84)
1427          intensity = 18;
1428       else if (effectiveRate<102)
1429          intensity = 19;
1430       else if (effectiveRate<130)
1431          intensity = 20;
1432       else
1433          intensity = 100;
1434       intensity = IMIN(st->end,IMAX(st->start, intensity));
1435    }
1436
1437    /* Bit allocation */
1438    ALLOC(fine_quant, st->mode->nbEBands, int);
1439    ALLOC(pulses, st->mode->nbEBands, int);
1440    ALLOC(fine_priority, st->mode->nbEBands, int);
1441
1442    /* bits =   packet size        -       where we are         - safety*/
1443    bits = (nbCompressedBytes*8<<BITRES) - ec_tell_frac(enc) - 1;
1444    anti_collapse_rsv = isTransient&&LM>=2&&bits>=(LM+2<<BITRES) ? (1<<BITRES) : 0;
1445    bits -= anti_collapse_rsv;
1446    codedBands = compute_allocation(st->mode, st->start, st->end, offsets, cap,
1447          alloc_trim, &intensity, &dual_stereo, bits, &balance, pulses,
1448          fine_quant, fine_priority, C, LM, enc, 1, st->lastCodedBands);
1449    st->lastCodedBands = codedBands;
1450
1451    quant_fine_energy(st->mode, st->start, st->end, oldBandE, error, fine_quant, enc, C);
1452
1453 #ifdef MEASURE_NORM_MSE
1454    float X0[3000];
1455    float bandE0[60];
1456    c=0; do 
1457       for (i=0;i<N;i++)
1458          X0[i+c*N] = X[i+c*N];
1459    while (++c<C);
1460    for (i=0;i<C*st->mode->nbEBands;i++)
1461       bandE0[i] = bandE[i];
1462 #endif
1463
1464    /* Residual quantisation */
1465    ALLOC(collapse_masks, C*st->mode->nbEBands, unsigned char);
1466    quant_all_bands(1, st->mode, st->start, st->end, X, C==2 ? X+N : NULL, collapse_masks,
1467          bandE, pulses, shortBlocks, st->spread_decision, dual_stereo, intensity, tf_res, resynth,
1468          nbCompressedBytes*(8<<BITRES)-anti_collapse_rsv, balance, enc, LM, codedBands, &st->rng);
1469
1470    if (anti_collapse_rsv > 0)
1471    {
1472       anti_collapse_on = st->consec_transient<2;
1473       ec_enc_bits(enc, anti_collapse_on, 1);
1474    }
1475    quant_energy_finalise(st->mode, st->start, st->end, oldBandE, error, fine_quant, fine_priority, nbCompressedBytes*8-ec_tell(enc), enc, C);
1476
1477    if (silence)
1478    {
1479       for (i=0;i<C*st->mode->nbEBands;i++)
1480          oldBandE[i] = -QCONST16(28.f,DB_SHIFT);
1481    }
1482
1483 #ifdef RESYNTH
1484    /* Re-synthesis of the coded audio if required */
1485    if (resynth)
1486    {
1487       celt_sig *out_mem[2];
1488       celt_sig *overlap_mem[2];
1489
1490       log2Amp(st->mode, st->start, st->end, bandE, oldBandE, C);
1491       if (silence)
1492       {
1493          for (i=0;i<C*st->mode->nbEBands;i++)
1494             bandE[i] = 0;
1495       }
1496
1497 #ifdef MEASURE_NORM_MSE
1498       measure_norm_mse(st->mode, X, X0, bandE, bandE0, M, N, C);
1499 #endif
1500       if (anti_collapse_on)
1501       {
1502          anti_collapse(st->mode, X, collapse_masks, LM, C, CC, N,
1503                st->start, st->end, oldBandE, oldLogE, oldLogE2, pulses, st->rng);
1504       }
1505
1506       /* Synthesis */
1507       denormalise_bands(st->mode, X, freq, bandE, effEnd, C, M);
1508
1509       CELT_MOVE(st->syn_mem[0], st->syn_mem[0]+N, MAX_PERIOD);
1510       if (CC==2)
1511          CELT_MOVE(st->syn_mem[1], st->syn_mem[1]+N, MAX_PERIOD);
1512
1513       c=0; do
1514          for (i=0;i<M*st->mode->eBands[st->start];i++)
1515             freq[c*N+i] = 0;
1516       while (++c<C);
1517       c=0; do
1518          for (i=M*st->mode->eBands[st->end];i<N;i++)
1519             freq[c*N+i] = 0;
1520       while (++c<C);
1521
1522       if (CC==2&&C==1)
1523       {
1524          for (i=0;i<N;i++)
1525             freq[N+i] = freq[i];
1526       }
1527
1528       out_mem[0] = st->syn_mem[0]+MAX_PERIOD;
1529       if (CC==2)
1530          out_mem[1] = st->syn_mem[1]+MAX_PERIOD;
1531
1532       c=0; do
1533          overlap_mem[c] = _overlap_mem + c*st->overlap;
1534       while (++c<CC);
1535
1536       compute_inv_mdcts(st->mode, shortBlocks, freq, out_mem, overlap_mem, CC, LM);
1537
1538 #ifdef ENABLE_POSTFILTER
1539       c=0; do {
1540          st->prefilter_period=IMAX(st->prefilter_period, COMBFILTER_MINPERIOD);
1541          st->prefilter_period_old=IMAX(st->prefilter_period_old, COMBFILTER_MINPERIOD);
1542          if (LM!=0)
1543          {
1544             comb_filter(out_mem[c], out_mem[c], st->prefilter_period, st->prefilter_period, st->overlap,
1545                   st->prefilter_gain, st->prefilter_gain, st->prefilter_tapset, st->prefilter_tapset,
1546                   NULL, 0);
1547             comb_filter(out_mem[c]+st->overlap, out_mem[c]+st->overlap, st->prefilter_period, pitch_index, N-st->overlap,
1548                   st->prefilter_gain, gain1, st->prefilter_tapset, prefilter_tapset,
1549                   st->mode->window, st->mode->overlap);
1550          } else {
1551             comb_filter(out_mem[c], out_mem[c], st->prefilter_period_old, st->prefilter_period, N,
1552                   st->prefilter_gain_old, st->prefilter_gain, st->prefilter_tapset_old, st->prefilter_tapset,
1553                   st->mode->window, st->mode->overlap);
1554          }
1555       } while (++c<CC);
1556 #endif /* ENABLE_POSTFILTER */
1557
1558       deemphasis(out_mem, (celt_word16*)pcm, N, CC, st->upsample, st->mode->preemph, st->preemph_memD);
1559       st->prefilter_period_old = st->prefilter_period;
1560       st->prefilter_gain_old = st->prefilter_gain;
1561       st->prefilter_tapset_old = st->prefilter_tapset;
1562    }
1563 #endif
1564
1565    st->prefilter_period = pitch_index;
1566    st->prefilter_gain = gain1;
1567    st->prefilter_tapset = prefilter_tapset;
1568
1569    if (CC==2&&C==1) {
1570       for (i=0;i<st->mode->nbEBands;i++)
1571          oldBandE[st->mode->nbEBands+i]=oldBandE[i];
1572    }
1573
1574    /* In case start or end were to change */
1575    c=0; do
1576    {
1577       for (i=0;i<st->start;i++)
1578          oldBandE[c*st->mode->nbEBands+i]=0;
1579       for (i=st->end;i<st->mode->nbEBands;i++)
1580          oldBandE[c*st->mode->nbEBands+i]=0;
1581    } while (++c<CC);
1582    if (!isTransient)
1583    {
1584       for (i=0;i<CC*st->mode->nbEBands;i++)
1585          oldLogE2[i] = oldLogE[i];
1586       for (i=0;i<CC*st->mode->nbEBands;i++)
1587          oldLogE[i] = oldBandE[i];
1588    } else {
1589       for (i=0;i<CC*st->mode->nbEBands;i++)
1590          oldLogE[i] = MIN16(oldLogE[i], oldBandE[i]);
1591    }
1592    if (isTransient)
1593       st->consec_transient++;
1594    else
1595       st->consec_transient=0;
1596    st->rng = enc->rng;
1597
1598    /* If there's any room left (can only happen for very high rates),
1599       it's already filled with zeros */
1600    ec_enc_done(enc);
1601    
1602    RESTORE_STACK;
1603    if (ec_get_error(enc))
1604       return CELT_CORRUPTED_DATA;
1605    else
1606       return nbCompressedBytes;
1607 }
1608
1609 #ifdef FIXED_POINT
1610 #ifndef DISABLE_FLOAT_API
1611 CELT_STATIC
1612 int celt_encode_with_ec_float(CELTEncoder * restrict st, const float * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes, ec_enc *enc)
1613 {
1614    int j, ret, C, N;
1615    VARDECL(celt_int16, in);
1616    ALLOC_STACK;
1617    SAVE_STACK;
1618
1619    if (pcm==NULL)
1620       return CELT_BAD_ARG;
1621
1622    C = CHANNELS(st->channels);
1623    N = frame_size;
1624    ALLOC(in, C*N, celt_int16);
1625
1626    for (j=0;j<C*N;j++)
1627      in[j] = FLOAT2INT16(pcm[j]);
1628
1629    ret=celt_encode_with_ec(st,in,frame_size,compressed,nbCompressedBytes, enc);
1630 #ifdef RESYNTH
1631    for (j=0;j<C*N;j++)
1632       ((float*)pcm)[j]=in[j]*(1.f/32768.f);
1633 #endif
1634    RESTORE_STACK;
1635    return ret;
1636
1637 }
1638 #endif /*DISABLE_FLOAT_API*/
1639 #else
1640 CELT_STATIC
1641 int celt_encode_with_ec(CELTEncoder * restrict st, const celt_int16 * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes, ec_enc *enc)
1642 {
1643    int j, ret, C, N;
1644    VARDECL(celt_sig, in);
1645    ALLOC_STACK;
1646    SAVE_STACK;
1647
1648    if (pcm==NULL)
1649       return CELT_BAD_ARG;
1650
1651    C=CHANNELS(st->channels);
1652    N=frame_size;
1653    ALLOC(in, C*N, celt_sig);
1654    for (j=0;j<C*N;j++) {
1655      in[j] = SCALEOUT(pcm[j]);
1656    }
1657
1658    ret = celt_encode_with_ec_float(st,in,frame_size,compressed,nbCompressedBytes, enc);
1659 #ifdef RESYNTH
1660    for (j=0;j<C*N;j++)
1661       ((celt_int16*)pcm)[j] = FLOAT2INT16(in[j]);
1662 #endif
1663    RESTORE_STACK;
1664    return ret;
1665 }
1666 #endif
1667
1668 int celt_encode(CELTEncoder * restrict st, const celt_int16 * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes)
1669 {
1670    return celt_encode_with_ec(st, pcm, frame_size, compressed, nbCompressedBytes, NULL);
1671 }
1672
1673 #ifndef DISABLE_FLOAT_API
1674 int celt_encode_float(CELTEncoder * restrict st, const float * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes)
1675 {
1676    return celt_encode_with_ec_float(st, pcm, frame_size, compressed, nbCompressedBytes, NULL);
1677 }
1678 #endif /* DISABLE_FLOAT_API */
1679
1680 int celt_encoder_ctl(CELTEncoder * restrict st, int request, ...)
1681 {
1682    va_list ap;
1683    
1684    va_start(ap, request);
1685    switch (request)
1686    {
1687       case CELT_GET_MODE_REQUEST:
1688       {
1689          const CELTMode ** value = va_arg(ap, const CELTMode**);
1690          if (value==0)
1691             goto bad_arg;
1692          *value=st->mode;
1693       }
1694       break;
1695       case CELT_SET_COMPLEXITY_REQUEST:
1696       {
1697          int value = va_arg(ap, celt_int32);
1698          if (value<0 || value>10)
1699             goto bad_arg;
1700          st->complexity = value;
1701       }
1702       break;
1703       case CELT_SET_START_BAND_REQUEST:
1704       {
1705          celt_int32 value = va_arg(ap, celt_int32);
1706          if (value<0 || value>=st->mode->nbEBands)
1707             goto bad_arg;
1708          st->start = value;
1709       }
1710       break;
1711       case CELT_SET_END_BAND_REQUEST:
1712       {
1713          celt_int32 value = va_arg(ap, celt_int32);
1714          if (value<1 || value>st->mode->nbEBands)
1715             goto bad_arg;
1716          st->end = value;
1717       }
1718       break;
1719       case CELT_SET_PREDICTION_REQUEST:
1720       {
1721          int value = va_arg(ap, celt_int32);
1722          if (value<0 || value>2)
1723             goto bad_arg;
1724          st->disable_pf = value<=1;
1725          st->force_intra = value==0;
1726       }
1727       break;
1728       case CELT_SET_VBR_CONSTRAINT_REQUEST:
1729       {
1730          celt_int32 value = va_arg(ap, celt_int32);
1731          st->constrained_vbr = value;
1732       }
1733       break;
1734       case CELT_SET_VBR_REQUEST:
1735       {
1736          celt_int32 value = va_arg(ap, celt_int32);
1737          st->vbr = value;
1738       }
1739       break;
1740       case CELT_SET_BITRATE_REQUEST:
1741       {
1742          celt_int32 value = va_arg(ap, celt_int32);
1743          if (value<=500)
1744             goto bad_arg;
1745          value = IMIN(value, 260000*st->channels);
1746          st->bitrate = value;
1747       }
1748       break;
1749       case CELT_SET_CHANNELS_REQUEST:
1750       {
1751          celt_int32 value = va_arg(ap, celt_int32);
1752          if (value<1 || value>2)
1753             goto bad_arg;
1754          st->stream_channels = value;
1755       }
1756       break;
1757       case CELT_RESET_STATE:
1758       {
1759          CELT_MEMSET((char*)&st->ENCODER_RESET_START, 0,
1760                celt_encoder_get_size_custom(st->mode, st->channels)-
1761                ((char*)&st->ENCODER_RESET_START - (char*)st));
1762          st->vbr_offset = 0;
1763          st->delayedIntra = 1;
1764          st->spread_decision = SPREAD_NORMAL;
1765          st->tonal_average = QCONST16(1.f,8);
1766       }
1767       break;
1768       case CELT_SET_INPUT_CLIPPING_REQUEST:
1769       {
1770          celt_int32 value = va_arg(ap, celt_int32);
1771          st->clip = value;
1772       }
1773       break;
1774       default:
1775          goto bad_request;
1776    }
1777    va_end(ap);
1778    return CELT_OK;
1779 bad_arg:
1780    va_end(ap);
1781    return CELT_BAD_ARG;
1782 bad_request:
1783    va_end(ap);
1784    return CELT_UNIMPLEMENTED;
1785 }
1786
1787 /**********************************************************************/
1788 /*                                                                    */
1789 /*                             DECODER                                */
1790 /*                                                                    */
1791 /**********************************************************************/
1792 #define DECODE_BUFFER_SIZE 2048
1793
1794 /** Decoder state 
1795  @brief Decoder state
1796  */
1797 struct CELTDecoder {
1798    const CELTMode *mode;
1799    int overlap;
1800    int channels;
1801    int stream_channels;
1802
1803    int downsample;
1804    int start, end;
1805
1806    /* Everything beyond this point gets cleared on a reset */
1807 #define DECODER_RESET_START rng
1808
1809    ec_uint32 rng;
1810    int last_pitch_index;
1811    int loss_count;
1812    int postfilter_period;
1813    int postfilter_period_old;
1814    celt_word16 postfilter_gain;
1815    celt_word16 postfilter_gain_old;
1816    int postfilter_tapset;
1817    int postfilter_tapset_old;
1818
1819    celt_sig preemph_memD[2];
1820    
1821    celt_sig _decode_mem[1]; /* Size = channels*(DECODE_BUFFER_SIZE+mode->overlap) */
1822    /* celt_word16 lpc[],  Size = channels*LPC_ORDER */
1823    /* celt_word16 oldEBands[], Size = channels*mode->nbEBands */
1824    /* celt_word16 oldLogE[], Size = channels*mode->nbEBands */
1825    /* celt_word16 oldLogE2[], Size = channels*mode->nbEBands */
1826    /* celt_word16 backgroundLogE[], Size = channels*mode->nbEBands */
1827 };
1828
1829 int celt_decoder_get_size(int channels)
1830 {
1831    const CELTMode *mode = celt_mode_create(48000, 960, NULL);
1832    return celt_decoder_get_size_custom(mode, channels);
1833 }
1834
1835 int celt_decoder_get_size_custom(const CELTMode *mode, int channels)
1836 {
1837    int size = sizeof(struct CELTDecoder)
1838             + (channels*(DECODE_BUFFER_SIZE+mode->overlap)-1)*sizeof(celt_sig)
1839             + channels*LPC_ORDER*sizeof(celt_word16)
1840             + 4*channels*mode->nbEBands*sizeof(celt_word16);
1841    return size;
1842 }
1843
1844 CELTDecoder *celt_decoder_create(int sampling_rate, int channels, int *error)
1845 {
1846    CELTDecoder *st;
1847    st = (CELTDecoder *)celt_alloc(celt_decoder_get_size(channels));
1848    if (st!=NULL && celt_decoder_init(st, sampling_rate, channels, error)==NULL)
1849    {
1850       celt_decoder_destroy(st);
1851       st = NULL;
1852    }
1853    return st;
1854 }
1855
1856 CELTDecoder *celt_decoder_create_custom(const CELTMode *mode, int channels, int *error)
1857 {
1858    CELTDecoder *st = (CELTDecoder *)celt_alloc(celt_decoder_get_size_custom(mode, channels));
1859    if (st!=NULL && celt_decoder_init_custom(st, mode, channels, error)==NULL)
1860    {
1861       celt_decoder_destroy(st);
1862       st = NULL;
1863    }
1864    return st;
1865 }
1866
1867 CELTDecoder *celt_decoder_init(CELTDecoder *st, int sampling_rate, int channels, int *error)
1868 {
1869    celt_decoder_init_custom(st, celt_mode_create(48000, 960, NULL), channels, error);
1870    st->downsample = resampling_factor(sampling_rate);
1871    if (st->downsample==0)
1872    {
1873       if (error)
1874          *error = CELT_BAD_ARG;
1875       return NULL;
1876    }
1877    return st;
1878 }
1879
1880 CELTDecoder *celt_decoder_init_custom(CELTDecoder *st, const CELTMode *mode, int channels, int *error)
1881 {
1882    if (channels < 0 || channels > 2)
1883    {
1884       if (error)
1885          *error = CELT_BAD_ARG;
1886       return NULL;
1887    }
1888
1889    if (st==NULL)
1890    {
1891       if (error)
1892          *error = CELT_ALLOC_FAIL;
1893       return NULL;
1894    }
1895
1896    CELT_MEMSET((char*)st, 0, celt_decoder_get_size_custom(mode, channels));
1897
1898    st->mode = mode;
1899    st->overlap = mode->overlap;
1900    st->stream_channels = st->channels = channels;
1901
1902    st->downsample = 1;
1903    st->start = 0;
1904    st->end = st->mode->effEBands;
1905
1906    st->loss_count = 0;
1907
1908    if (error)
1909       *error = CELT_OK;
1910    return st;
1911 }
1912
1913 void celt_decoder_destroy(CELTDecoder *st)
1914 {
1915    celt_free(st);
1916 }
1917
1918 static void celt_decode_lost(CELTDecoder * restrict st, celt_word16 * restrict pcm, int N, int LM)
1919 {
1920    int c;
1921    int pitch_index;
1922    int overlap = st->mode->overlap;
1923    celt_word16 fade = Q15ONE;
1924    int i, len;
1925    const int C = CHANNELS(st->channels);
1926    int offset;
1927    celt_sig *out_mem[2];
1928    celt_sig *decode_mem[2];
1929    celt_sig *overlap_mem[2];
1930    celt_word16 *lpc;
1931    celt_word32 *out_syn[2];
1932    celt_word16 *oldBandE, *oldLogE2, *backgroundLogE;
1933    int plc=1;
1934    SAVE_STACK;
1935    
1936    c=0; do {
1937       decode_mem[c] = st->_decode_mem + c*(DECODE_BUFFER_SIZE+st->overlap);
1938       out_mem[c] = decode_mem[c]+DECODE_BUFFER_SIZE-MAX_PERIOD;
1939       overlap_mem[c] = decode_mem[c]+DECODE_BUFFER_SIZE;
1940    } while (++c<C);
1941    lpc = (celt_word16*)(st->_decode_mem+(DECODE_BUFFER_SIZE+st->overlap)*C);
1942    oldBandE = lpc+C*LPC_ORDER;
1943    oldLogE2 = oldBandE + C*st->mode->nbEBands;
1944    backgroundLogE = oldLogE2  + C*st->mode->nbEBands;
1945
1946    out_syn[0] = out_mem[0]+MAX_PERIOD-N;
1947    if (C==2)
1948       out_syn[1] = out_mem[1]+MAX_PERIOD-N;
1949
1950    len = N+st->mode->overlap;
1951    
1952    if (st->loss_count >= 5)
1953    {
1954       VARDECL(celt_sig, freq);
1955       VARDECL(celt_norm, X);
1956       VARDECL(celt_ener, bandE);
1957       celt_uint32 seed;
1958
1959       ALLOC(freq, C*N, celt_sig); /**< Interleaved signal MDCTs */
1960       ALLOC(X, C*N, celt_norm);   /**< Interleaved normalised MDCTs */
1961       ALLOC(bandE, st->mode->nbEBands*C, celt_ener);
1962
1963       log2Amp(st->mode, st->start, st->end, bandE, backgroundLogE, C);
1964
1965       seed = st->rng;
1966       for (i=0;i<C*N;i++)
1967       {
1968             seed = lcg_rand(seed);
1969             X[i] = (celt_int32)(seed)>>20;
1970       }
1971       st->rng = seed;
1972       for (c=0;c<C;c++)
1973          for (i=0;i<st->mode->nbEBands;i++)
1974             renormalise_vector(X+N*c+(st->mode->eBands[i]<<LM), (st->mode->eBands[i+1]-st->mode->eBands[i])<<LM, Q15ONE);
1975
1976       denormalise_bands(st->mode, X, freq, bandE, st->mode->nbEBands, C, 1<<LM);
1977
1978       compute_inv_mdcts(st->mode, 0, freq, out_syn, overlap_mem, C, LM);
1979       plc = 0;
1980    } else if (st->loss_count == 0)
1981    {
1982       celt_word16 pitch_buf[MAX_PERIOD>>1];
1983       int len2 = len;
1984       /* FIXME: This is a kludge */
1985       if (len2>MAX_PERIOD>>1)
1986          len2 = MAX_PERIOD>>1;
1987       pitch_downsample(out_mem, pitch_buf, MAX_PERIOD, C);
1988       pitch_search(pitch_buf+((MAX_PERIOD-len2)>>1), pitch_buf, len2,
1989                    MAX_PERIOD-len2-100, &pitch_index);
1990       pitch_index = MAX_PERIOD-len2-pitch_index;
1991       st->last_pitch_index = pitch_index;
1992    } else {
1993       pitch_index = st->last_pitch_index;
1994       fade = QCONST16(.8f,15);
1995    }
1996
1997    if (plc)
1998    {
1999       c=0; do {
2000          /* FIXME: This is more memory than necessary */
2001          celt_word32 e[2*MAX_PERIOD];
2002          celt_word16 exc[2*MAX_PERIOD];
2003          celt_word32 ac[LPC_ORDER+1];
2004          celt_word16 decay = 1;
2005          celt_word32 S1=0;
2006          celt_word16 mem[LPC_ORDER]={0};
2007
2008          offset = MAX_PERIOD-pitch_index;
2009          for (i=0;i<MAX_PERIOD;i++)
2010             exc[i] = ROUND16(out_mem[c][i], SIG_SHIFT);
2011
2012          if (st->loss_count == 0)
2013          {
2014             _celt_autocorr(exc, ac, st->mode->window, st->mode->overlap,
2015                   LPC_ORDER, MAX_PERIOD);
2016
2017             /* Noise floor -40 dB */
2018 #ifdef FIXED_POINT
2019             ac[0] += SHR32(ac[0],13);
2020 #else
2021             ac[0] *= 1.0001f;
2022 #endif
2023             /* Lag windowing */
2024             for (i=1;i<=LPC_ORDER;i++)
2025             {
2026                /*ac[i] *= exp(-.5*(2*M_PI*.002*i)*(2*M_PI*.002*i));*/
2027 #ifdef FIXED_POINT
2028                ac[i] -= MULT16_32_Q15(2*i*i, ac[i]);
2029 #else
2030                ac[i] -= ac[i]*(.008f*i)*(.008f*i);
2031 #endif
2032             }
2033
2034             _celt_lpc(lpc+c*LPC_ORDER, ac, LPC_ORDER);
2035          }
2036          for (i=0;i<LPC_ORDER;i++)
2037             mem[i] = ROUND16(out_mem[c][MAX_PERIOD-1-i], SIG_SHIFT);
2038          fir(exc, lpc+c*LPC_ORDER, exc, MAX_PERIOD, LPC_ORDER, mem);
2039          /*for (i=0;i<MAX_PERIOD;i++)printf("%d ", exc[i]); printf("\n");*/
2040          /* Check if the waveform is decaying (and if so how fast) */
2041          {
2042             celt_word32 E1=1, E2=1;
2043             int period;
2044             if (pitch_index <= MAX_PERIOD/2)
2045                period = pitch_index;
2046             else
2047                period = MAX_PERIOD/2;
2048             for (i=0;i<period;i++)
2049             {
2050                E1 += SHR32(MULT16_16(exc[MAX_PERIOD-period+i],exc[MAX_PERIOD-period+i]),8);
2051                E2 += SHR32(MULT16_16(exc[MAX_PERIOD-2*period+i],exc[MAX_PERIOD-2*period+i]),8);
2052             }
2053             if (E1 > E2)
2054                E1 = E2;
2055             decay = celt_sqrt(frac_div32(SHR(E1,1),E2));
2056          }
2057
2058          /* Copy excitation, taking decay into account */
2059          for (i=0;i<len+st->mode->overlap;i++)
2060          {
2061             celt_word16 tmp;
2062             if (offset+i >= MAX_PERIOD)
2063             {
2064                offset -= pitch_index;
2065                decay = MULT16_16_Q15(decay, decay);
2066             }
2067             e[i] = SHL32(EXTEND32(MULT16_16_Q15(decay, exc[offset+i])), SIG_SHIFT);
2068             tmp = ROUND16(out_mem[c][offset+i],SIG_SHIFT);
2069             S1 += SHR32(MULT16_16(tmp,tmp),8);
2070          }
2071          for (i=0;i<LPC_ORDER;i++)
2072             mem[i] = ROUND16(out_mem[c][MAX_PERIOD-1-i], SIG_SHIFT);
2073          for (i=0;i<len+st->mode->overlap;i++)
2074             e[i] = MULT16_32_Q15(fade, e[i]);
2075          iir(e, lpc+c*LPC_ORDER, e, len+st->mode->overlap, LPC_ORDER, mem);
2076
2077          {
2078             celt_word32 S2=0;
2079             for (i=0;i<len+overlap;i++)
2080             {
2081                celt_word16 tmp = ROUND16(e[i],SIG_SHIFT);
2082                S2 += SHR32(MULT16_16(tmp,tmp),8);
2083             }
2084             /* This checks for an "explosion" in the synthesis */
2085 #ifdef FIXED_POINT
2086             if (!(S1 > SHR32(S2,2)))
2087 #else
2088                /* Float test is written this way to catch NaNs at the same time */
2089                if (!(S1 > 0.2f*S2))
2090 #endif
2091                {
2092                   for (i=0;i<len+overlap;i++)
2093                      e[i] = 0;
2094                } else if (S1 < S2)
2095                {
2096                   celt_word16 ratio = celt_sqrt(frac_div32(SHR32(S1,1)+1,S2+1));
2097                   for (i=0;i<len+overlap;i++)
2098                      e[i] = MULT16_32_Q15(ratio, e[i]);
2099                }
2100          }
2101
2102 #ifdef ENABLE_POSTFILTER
2103          /* Apply post-filter to the MDCT overlap of the previous frame */
2104          comb_filter(out_mem[c]+MAX_PERIOD, out_mem[c]+MAX_PERIOD, st->postfilter_period, st->postfilter_period, st->overlap,
2105                st->postfilter_gain, st->postfilter_gain, st->postfilter_tapset, st->postfilter_tapset,
2106                NULL, 0);
2107 #endif /* ENABLE_POSTFILTER */
2108
2109          for (i=0;i<MAX_PERIOD+st->mode->overlap-N;i++)
2110             out_mem[c][i] = out_mem[c][N+i];
2111
2112          /* Apply TDAC to the concealed audio so that it blends with the
2113          previous and next frames */
2114          for (i=0;i<overlap/2;i++)
2115          {
2116             celt_word32 tmp;
2117             tmp = MULT16_32_Q15(st->mode->window[i],           e[N+overlap-1-i]) +
2118                   MULT16_32_Q15(st->mode->window[overlap-i-1], e[N+i          ]);
2119             out_mem[c][MAX_PERIOD+i] = MULT16_32_Q15(st->mode->window[overlap-i-1], tmp);
2120             out_mem[c][MAX_PERIOD+overlap-i-1] = MULT16_32_Q15(st->mode->window[i], tmp);
2121          }
2122          for (i=0;i<N;i++)
2123             out_mem[c][MAX_PERIOD-N+i] = e[i];
2124
2125 #ifdef ENABLE_POSTFILTER
2126          /* Apply pre-filter to the MDCT overlap for the next frame (post-filter will be applied then) */
2127          comb_filter(e, out_mem[c]+MAX_PERIOD, st->postfilter_period, st->postfilter_period, st->overlap,
2128                -st->postfilter_gain, -st->postfilter_gain, st->postfilter_tapset, st->postfilter_tapset,
2129                NULL, 0);
2130 #endif /* ENABLE_POSTFILTER */
2131          for (i=0;i<overlap;i++)
2132             out_mem[c][MAX_PERIOD+i] = e[i];
2133       } while (++c<C);
2134    }
2135
2136    deemphasis(out_syn, pcm, N, C, st->downsample, st->mode->preemph, st->preemph_memD);
2137    
2138    st->loss_count++;
2139
2140    RESTORE_STACK;
2141 }
2142
2143 #ifdef FIXED_POINT
2144 CELT_STATIC
2145 int celt_decode_with_ec(CELTDecoder * restrict st, const unsigned char *data, int len, celt_int16 * restrict pcm, int frame_size, ec_dec *dec)
2146 {
2147 #else
2148 CELT_STATIC
2149 int celt_decode_with_ec_float(CELTDecoder * restrict st, const unsigned char *data, int len, celt_sig * restrict pcm, int frame_size, ec_dec *dec)
2150 {
2151 #endif
2152    int c, i, N;
2153    int spread_decision;
2154    int bits;
2155    ec_dec _dec;
2156    VARDECL(celt_sig, freq);
2157    VARDECL(celt_norm, X);
2158    VARDECL(celt_ener, bandE);
2159    VARDECL(int, fine_quant);
2160    VARDECL(int, pulses);
2161    VARDECL(int, cap);
2162    VARDECL(int, offsets);
2163    VARDECL(int, fine_priority);
2164    VARDECL(int, tf_res);
2165    VARDECL(unsigned char, collapse_masks);
2166    celt_sig *out_mem[2];
2167    celt_sig *decode_mem[2];
2168    celt_sig *overlap_mem[2];
2169    celt_sig *out_syn[2];
2170    celt_word16 *lpc;
2171    celt_word16 *oldBandE, *oldLogE, *oldLogE2, *backgroundLogE;
2172
2173    int shortBlocks;
2174    int isTransient;
2175    int intra_ener;
2176    const int CC = CHANNELS(st->channels);
2177    int LM, M;
2178    int effEnd;
2179    int codedBands;
2180    int alloc_trim;
2181    int postfilter_pitch;
2182    celt_word16 postfilter_gain;
2183    int intensity=0;
2184    int dual_stereo=0;
2185    celt_int32 total_bits;
2186    celt_int32 balance;
2187    celt_int32 tell;
2188    int dynalloc_logp;
2189    int postfilter_tapset;
2190    int anti_collapse_rsv;
2191    int anti_collapse_on=0;
2192    int silence;
2193    const int C = CHANNELS(st->stream_channels);
2194
2195    SAVE_STACK;
2196
2197    if (len<0 || len>1275 || pcm==NULL)
2198       return CELT_BAD_ARG;
2199
2200    frame_size *= st->downsample;
2201    for (LM=0;LM<=st->mode->maxLM;LM++)
2202       if (st->mode->shortMdctSize<<LM==frame_size)
2203          break;
2204    if (LM>st->mode->maxLM)
2205       return CELT_BAD_ARG;
2206    M=1<<LM;
2207
2208    c=0; do {
2209       decode_mem[c] = st->_decode_mem + c*(DECODE_BUFFER_SIZE+st->overlap);
2210       out_mem[c] = decode_mem[c]+DECODE_BUFFER_SIZE-MAX_PERIOD;
2211       overlap_mem[c] = decode_mem[c]+DECODE_BUFFER_SIZE;
2212    } while (++c<CC);
2213    lpc = (celt_word16*)(st->_decode_mem+(DECODE_BUFFER_SIZE+st->overlap)*CC);
2214    oldBandE = lpc+CC*LPC_ORDER;
2215    oldLogE = oldBandE + CC*st->mode->nbEBands;
2216    oldLogE2 = oldLogE + CC*st->mode->nbEBands;
2217    backgroundLogE = oldLogE2  + CC*st->mode->nbEBands;
2218
2219    N = M*st->mode->shortMdctSize;
2220
2221    effEnd = st->end;
2222    if (effEnd > st->mode->effEBands)
2223       effEnd = st->mode->effEBands;
2224
2225    ALLOC(freq, CC*N, celt_sig); /**< Interleaved signal MDCTs */
2226    ALLOC(X, CC*N, celt_norm);   /**< Interleaved normalised MDCTs */
2227    ALLOC(bandE, st->mode->nbEBands*CC, celt_ener);
2228    c=0; do
2229       for (i=0;i<M*st->mode->eBands[st->start];i++)
2230          X[c*N+i] = 0;
2231    while (++c<CC);
2232    c=0; do   
2233       for (i=M*st->mode->eBands[effEnd];i<N;i++)
2234          X[c*N+i] = 0;
2235    while (++c<CC);
2236
2237    if (data == NULL || len<=1)
2238    {
2239       celt_decode_lost(st, pcm, N, LM);
2240       RESTORE_STACK;
2241       return CELT_OK;
2242    }
2243    if (len<0) {
2244      RESTORE_STACK;
2245      return CELT_BAD_ARG;
2246    }
2247    
2248    if (dec == NULL)
2249    {
2250       ec_dec_init(&_dec,(unsigned char*)data,len);
2251       dec = &_dec;
2252    }
2253
2254    if (C>CC)
2255    {
2256       RESTORE_STACK;
2257       return CELT_CORRUPTED_DATA;
2258    } else if (C<CC)
2259    {
2260       for (i=0;i<st->mode->nbEBands;i++)
2261          oldBandE[i]=MAX16(oldBandE[i],oldBandE[st->mode->nbEBands+i]);
2262    }
2263
2264    total_bits = len*8;
2265    tell = ec_tell(dec);
2266
2267    if (tell==1)
2268       silence = ec_dec_bit_logp(dec, 15);
2269    else
2270       silence = 0;
2271    if (silence)
2272    {
2273       /* Pretend we've read all the remaining bits */
2274       tell = len*8;
2275       dec->nbits_total+=tell-ec_tell(dec);
2276    }
2277
2278    postfilter_gain = 0;
2279    postfilter_pitch = 0;
2280    postfilter_tapset = 0;
2281    if (st->start==0 && tell+16 <= total_bits)
2282    {
2283       if(ec_dec_bit_logp(dec, 1))
2284       {
2285 #ifdef ENABLE_POSTFILTER
2286          int qg, octave;
2287          octave = ec_dec_uint(dec, 6);
2288          postfilter_pitch = (16<<octave)+ec_dec_bits(dec, 4+octave)-1;
2289          qg = ec_dec_bits(dec, 3);
2290          if (ec_tell(dec)+2<=total_bits)
2291             postfilter_tapset = ec_dec_icdf(dec, tapset_icdf, 2);
2292          postfilter_gain = QCONST16(.09375f,15)*(qg+1);
2293 #else /* ENABLE_POSTFILTER */
2294          RESTORE_STACK;
2295          return CELT_CORRUPTED_DATA;
2296 #endif /* ENABLE_POSTFILTER */
2297       }
2298       tell = ec_tell(dec);
2299    }
2300
2301    if (LM > 0 && tell+3 <= total_bits)
2302    {
2303       isTransient = ec_dec_bit_logp(dec, 3);
2304       tell = ec_tell(dec);
2305    }
2306    else
2307       isTransient = 0;
2308
2309    if (isTransient)
2310       shortBlocks = M;
2311    else
2312       shortBlocks = 0;
2313
2314    /* Decode the global flags (first symbols in the stream) */
2315    intra_ener = tell+3<=total_bits ? ec_dec_bit_logp(dec, 3) : 0;
2316    /* Get band energies */
2317    unquant_coarse_energy(st->mode, st->start, st->end, oldBandE,
2318          intra_ener, dec, C, LM);
2319
2320    ALLOC(tf_res, st->mode->nbEBands, int);
2321    tf_decode(st->start, st->end, isTransient, tf_res, LM, dec);
2322
2323    tell = ec_tell(dec);
2324    spread_decision = SPREAD_NORMAL;
2325    if (tell+4 <= total_bits)
2326       spread_decision = ec_dec_icdf(dec, spread_icdf, 5);
2327
2328    ALLOC(pulses, st->mode->nbEBands, int);
2329    ALLOC(cap, st->mode->nbEBands, int);
2330    ALLOC(offsets, st->mode->nbEBands, int);
2331    ALLOC(fine_priority, st->mode->nbEBands, int);
2332
2333    init_caps(st->mode,cap,LM,C);
2334
2335    dynalloc_logp = 6;
2336    total_bits<<=BITRES;
2337    tell = ec_tell_frac(dec);
2338    for (i=st->start;i<st->end;i++)
2339    {
2340       int width, quanta;
2341       int dynalloc_loop_logp;
2342       int boost;
2343       width = C*(st->mode->eBands[i+1]-st->mode->eBands[i])<<LM;
2344       /* quanta is 6 bits, but no more than 1 bit/sample
2345          and no less than 1/8 bit/sample */
2346       quanta = IMIN(width<<BITRES, IMAX(6<<BITRES, width));
2347       dynalloc_loop_logp = dynalloc_logp;
2348       boost = 0;
2349       while (tell+(dynalloc_loop_logp<<BITRES) < total_bits && boost < cap[i])
2350       {
2351          int flag;
2352          flag = ec_dec_bit_logp(dec, dynalloc_loop_logp);
2353          tell = ec_tell_frac(dec);
2354          if (!flag)
2355             break;
2356          boost += quanta;
2357          total_bits -= quanta;
2358          dynalloc_loop_logp = 1;
2359       }
2360       offsets[i] = boost;
2361       /* Making dynalloc more likely */
2362       if (boost>0)
2363          dynalloc_logp = IMAX(2, dynalloc_logp-1);
2364    }
2365
2366    ALLOC(fine_quant, st->mode->nbEBands, int);
2367    alloc_trim = tell+(6<<BITRES) <= total_bits ?
2368          ec_dec_icdf(dec, trim_icdf, 7) : 5;
2369
2370    bits = (len*8<<BITRES) - ec_tell_frac(dec) - 1;
2371    anti_collapse_rsv = isTransient&&LM>=2&&bits>=(LM+2<<BITRES) ? (1<<BITRES) : 0;
2372    bits -= anti_collapse_rsv;
2373    codedBands = compute_allocation(st->mode, st->start, st->end, offsets, cap,
2374          alloc_trim, &intensity, &dual_stereo, bits, &balance, pulses,
2375          fine_quant, fine_priority, C, LM, dec, 0, 0);
2376    
2377    unquant_fine_energy(st->mode, st->start, st->end, oldBandE, fine_quant, dec, C);
2378
2379    /* Decode fixed codebook */
2380    ALLOC(collapse_masks, C*st->mode->nbEBands, unsigned char);
2381    quant_all_bands(0, st->mode, st->start, st->end, X, C==2 ? X+N : NULL, collapse_masks,
2382          NULL, pulses, shortBlocks, spread_decision, dual_stereo, intensity, tf_res, 1,
2383          len*(8<<BITRES)-anti_collapse_rsv, balance, dec, LM, codedBands, &st->rng);
2384
2385    if (anti_collapse_rsv > 0)
2386    {
2387       anti_collapse_on = ec_dec_bits(dec, 1);
2388    }
2389
2390    unquant_energy_finalise(st->mode, st->start, st->end, oldBandE,
2391          fine_quant, fine_priority, len*8-ec_tell(dec), dec, C);
2392
2393    if (anti_collapse_on)
2394       anti_collapse(st->mode, X, collapse_masks, LM, C, CC, N,
2395             st->start, st->end, oldBandE, oldLogE, oldLogE2, pulses, st->rng);
2396
2397    log2Amp(st->mode, st->start, st->end, bandE, oldBandE, C);
2398
2399    if (silence)
2400    {
2401       for (i=0;i<C*st->mode->nbEBands;i++)
2402       {
2403          bandE[i] = 0;
2404          oldBandE[i] = -QCONST16(28.f,DB_SHIFT);
2405       }
2406    }
2407    /* Synthesis */
2408    denormalise_bands(st->mode, X, freq, bandE, effEnd, C, M);
2409
2410    CELT_MOVE(decode_mem[0], decode_mem[0]+N, DECODE_BUFFER_SIZE-N);
2411    if (CC==2)
2412       CELT_MOVE(decode_mem[1], decode_mem[1]+N, DECODE_BUFFER_SIZE-N);
2413
2414    c=0; do
2415       for (i=0;i<M*st->mode->eBands[st->start];i++)
2416          freq[c*N+i] = 0;
2417    while (++c<C);
2418    c=0; do {
2419       int bound = M*st->mode->eBands[effEnd];
2420       if (st->downsample!=1)
2421          bound = IMIN(bound, N/st->downsample);
2422       for (i=bound;i<N;i++)
2423          freq[c*N+i] = 0;
2424    } while (++c<C);
2425
2426    out_syn[0] = out_mem[0]+MAX_PERIOD-N;
2427    if (CC==2)
2428       out_syn[1] = out_mem[1]+MAX_PERIOD-N;
2429
2430    if (CC==2&&C==1)
2431    {
2432       for (i=0;i<N;i++)
2433          freq[N+i] = freq[i];
2434    }
2435
2436    /* Compute inverse MDCTs */
2437    compute_inv_mdcts(st->mode, shortBlocks, freq, out_syn, overlap_mem, CC, LM);
2438
2439 #ifdef ENABLE_POSTFILTER
2440    c=0; do {
2441       st->postfilter_period=IMAX(st->postfilter_period, COMBFILTER_MINPERIOD);
2442       st->postfilter_period_old=IMAX(st->postfilter_period_old, COMBFILTER_MINPERIOD);
2443       if (LM!=0)
2444       {
2445          comb_filter(out_syn[c], out_syn[c], st->postfilter_period, st->postfilter_period, st->overlap,
2446                st->postfilter_gain, st->postfilter_gain, st->postfilter_tapset, st->postfilter_tapset,
2447                NULL, 0);
2448          comb_filter(out_syn[c]+st->overlap, out_syn[c]+st->overlap, st->postfilter_period, postfilter_pitch, N-st->overlap,
2449                st->postfilter_gain, postfilter_gain, st->postfilter_tapset, postfilter_tapset,
2450                st->mode->window, st->mode->overlap);
2451       } else {
2452          comb_filter(out_syn[c], out_syn[c], st->postfilter_period_old, st->postfilter_period, N-st->overlap,
2453                st->postfilter_gain_old, st->postfilter_gain, st->postfilter_tapset_old, st->postfilter_tapset,
2454                st->mode->window, st->mode->overlap);
2455       }
2456    } while (++c<CC);
2457    st->postfilter_period_old = st->postfilter_period;
2458    st->postfilter_gain_old = st->postfilter_gain;
2459    st->postfilter_tapset_old = st->postfilter_tapset;
2460    st->postfilter_period = postfilter_pitch;
2461    st->postfilter_gain = postfilter_gain;
2462    st->postfilter_tapset = postfilter_tapset;
2463 #endif /* ENABLE_POSTFILTER */
2464
2465    if (CC==2&&C==1) {
2466       for (i=0;i<st->mode->nbEBands;i++)
2467          oldBandE[st->mode->nbEBands+i]=oldBandE[i];
2468    }
2469
2470    /* In case start or end were to change */
2471    c=0; do
2472    {
2473       for (i=0;i<st->start;i++)
2474          oldBandE[c*st->mode->nbEBands+i]=0;
2475       for (i=st->end;i<st->mode->nbEBands;i++)
2476          oldBandE[c*st->mode->nbEBands+i]=0;
2477    } while (++c<CC);
2478    if (!isTransient)
2479    {
2480       for (i=0;i<CC*st->mode->nbEBands;i++)
2481          oldLogE2[i] = oldLogE[i];
2482       for (i=0;i<CC*st->mode->nbEBands;i++)
2483          oldLogE[i] = oldBandE[i];
2484       for (i=0;i<CC*st->mode->nbEBands;i++)
2485          backgroundLogE[i] = MIN16(backgroundLogE[i] + M*QCONST16(0.001f,DB_SHIFT), oldBandE[i]);
2486    } else {
2487       for (i=0;i<CC*st->mode->nbEBands;i++)
2488          oldLogE[i] = MIN16(oldLogE[i], oldBandE[i]);
2489    }
2490    st->rng = dec->rng;
2491
2492    deemphasis(out_syn, pcm, N, CC, st->downsample, st->mode->preemph, st->preemph_memD);
2493    st->loss_count = 0;
2494    RESTORE_STACK;
2495    if (ec_tell(dec) > 8*len || ec_get_error(dec))
2496       return CELT_CORRUPTED_DATA;
2497    else
2498       return CELT_OK;
2499 }
2500
2501 #ifdef FIXED_POINT
2502 #ifndef DISABLE_FLOAT_API
2503 CELT_STATIC
2504 int celt_decode_with_ec_float(CELTDecoder * restrict st, const unsigned char *data, int len, float * restrict pcm, int frame_size, ec_dec *dec)
2505 {
2506    int j, ret, C, N;
2507    VARDECL(celt_int16, out);
2508    ALLOC_STACK;
2509    SAVE_STACK;
2510
2511    if (pcm==NULL)
2512       return CELT_BAD_ARG;
2513
2514    C = CHANNELS(st->channels);
2515    N = frame_size;
2516    
2517    ALLOC(out, C*N, celt_int16);
2518    ret=celt_decode_with_ec(st, data, len, out, frame_size, dec);
2519    if (ret==0)
2520       for (j=0;j<C*N;j++)
2521          pcm[j]=out[j]*(1.f/32768.f);
2522      
2523    RESTORE_STACK;
2524    return ret;
2525 }
2526 #endif /*DISABLE_FLOAT_API*/
2527 #else
2528 CELT_STATIC
2529 int celt_decode_with_ec(CELTDecoder * restrict st, const unsigned char *data, int len, celt_int16 * restrict pcm, int frame_size, ec_dec *dec)
2530 {
2531    int j, ret, C, N;
2532    VARDECL(celt_sig, out);
2533    ALLOC_STACK;
2534    SAVE_STACK;
2535
2536    if (pcm==NULL)
2537       return CELT_BAD_ARG;
2538
2539    C = CHANNELS(st->channels);
2540    N = frame_size;
2541    ALLOC(out, C*N, celt_sig);
2542
2543    ret=celt_decode_with_ec_float(st, data, len, out, frame_size, dec);
2544
2545    if (ret==0)
2546       for (j=0;j<C*N;j++)
2547          pcm[j] = FLOAT2INT16 (out[j]);
2548    
2549    RESTORE_STACK;
2550    return ret;
2551 }
2552 #endif
2553
2554 int celt_decode(CELTDecoder * restrict st, const unsigned char *data, int len, celt_int16 * restrict pcm, int frame_size)
2555 {
2556    return celt_decode_with_ec(st, data, len, pcm, frame_size, NULL);
2557 }
2558
2559 #ifndef DISABLE_FLOAT_API
2560 int celt_decode_float(CELTDecoder * restrict st, const unsigned char *data, int len, float * restrict pcm, int frame_size)
2561 {
2562    return celt_decode_with_ec_float(st, data, len, pcm, frame_size, NULL);
2563 }
2564 #endif /* DISABLE_FLOAT_API */
2565
2566 int celt_decoder_ctl(CELTDecoder * restrict st, int request, ...)
2567 {
2568    va_list ap;
2569
2570    va_start(ap, request);
2571    switch (request)
2572    {
2573       case CELT_GET_MODE_REQUEST:
2574       {
2575          const CELTMode ** value = va_arg(ap, const CELTMode**);
2576          if (value==0)
2577             goto bad_arg;
2578          *value=st->mode;
2579       }
2580       break;
2581       case CELT_SET_START_BAND_REQUEST:
2582       {
2583          celt_int32 value = va_arg(ap, celt_int32);
2584          if (value<0 || value>=st->mode->nbEBands)
2585             goto bad_arg;
2586          st->start = value;
2587       }
2588       break;
2589       case CELT_SET_END_BAND_REQUEST:
2590       {
2591          celt_int32 value = va_arg(ap, celt_int32);
2592          if (value<0 || value>=st->mode->nbEBands)
2593             goto bad_arg;
2594          st->end = value;
2595       }
2596       break;
2597       case CELT_SET_CHANNELS_REQUEST:
2598       {
2599          celt_int32 value = va_arg(ap, celt_int32);
2600          if (value<1 || value>2)
2601             goto bad_arg;
2602          st->stream_channels = value;
2603       }
2604       break;
2605       case CELT_RESET_STATE:
2606       {
2607          CELT_MEMSET((char*)&st->DECODER_RESET_START, 0,
2608                celt_decoder_get_size_custom(st->mode, st->channels)-
2609                ((char*)&st->DECODER_RESET_START - (char*)st));
2610       }
2611       break;
2612       default:
2613          goto bad_request;
2614    }
2615    va_end(ap);
2616    return CELT_OK;
2617 bad_arg:
2618    va_end(ap);
2619    return CELT_BAD_ARG;
2620 bad_request:
2621       va_end(ap);
2622   return CELT_UNIMPLEMENTED;
2623 }
2624
2625 const char *celt_strerror(int error)
2626 {
2627    static const char *error_strings[8] = {
2628       "success",
2629       "invalid argument",
2630       "invalid mode",
2631       "internal error",
2632       "corrupted stream",
2633       "request not implemented",
2634       "invalid state",
2635       "memory allocation failed"
2636    };
2637    if (error > 0 || error < -7)
2638       return "unknown error";
2639    else 
2640       return error_strings[-error];
2641 }
2642