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