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