Increase caps/allocation accuracy.
[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 void init_caps(const CELTMode *m,int *cap,int LM,int C)
760 {
761    int i;
762    for (i=0;i<m->nbEBands;i++)
763    {
764       int N;
765       N=(m->eBands[i+1]-m->eBands[i])<<LM;
766       cap[i] = (m->cache.caps[m->nbEBands*(2*LM+C-1)+i]+64)*C*N>>2;
767    }
768 }
769
770 static int alloc_trim_analysis(const CELTMode *m, const celt_norm *X,
771       const celt_word16 *bandLogE, int nbEBands, int LM, int C, int N0)
772 {
773    int i;
774    celt_word32 diff=0;
775    int c;
776    int trim_index = 5;
777    if (C==2)
778    {
779       celt_word16 sum = 0; /* Q10 */
780       /* Compute inter-channel correlation for low frequencies */
781       for (i=0;i<8;i++)
782       {
783          int j;
784          celt_word32 partial = 0;
785          for (j=m->eBands[i]<<LM;j<m->eBands[i+1]<<LM;j++)
786             partial = MAC16_16(partial, X[j], X[N0+j]);
787          sum = ADD16(sum, EXTRACT16(SHR32(partial, 18)));
788       }
789       sum = MULT16_16_Q15(QCONST16(1.f/8, 15), sum);
790       /*printf ("%f\n", sum);*/
791       if (sum > QCONST16(.995f,10))
792          trim_index-=4;
793       else if (sum > QCONST16(.92f,10))
794          trim_index-=3;
795       else if (sum > QCONST16(.85f,10))
796          trim_index-=2;
797       else if (sum > QCONST16(.8f,10))
798          trim_index-=1;
799    }
800
801    /* Estimate spectral tilt */
802    c=0; do {
803       for (i=0;i<nbEBands-1;i++)
804       {
805          diff += bandLogE[i+c*nbEBands]*(celt_int32)(2+2*i-nbEBands);
806       }
807    } while (++c<0);
808    diff /= C*(nbEBands-1);
809    /*printf("%f\n", diff);*/
810    if (diff > QCONST16(2.f, DB_SHIFT))
811       trim_index--;
812    if (diff > QCONST16(8.f, DB_SHIFT))
813       trim_index--;
814    if (diff < -QCONST16(4.f, DB_SHIFT))
815       trim_index++;
816    if (diff < -QCONST16(10.f, DB_SHIFT))
817       trim_index++;
818
819    if (trim_index<0)
820       trim_index = 0;
821    if (trim_index>10)
822       trim_index = 10;
823    return trim_index;
824 }
825
826 static int stereo_analysis(const CELTMode *m, const celt_norm *X,
827       int LM, int N0)
828 {
829    int i;
830    int thetas;
831    celt_word32 sumLR = EPSILON, sumMS = EPSILON;
832
833    /* Use the L1 norm to model the entropy of the L/R signal vs the M/S signal */
834    for (i=0;i<13;i++)
835    {
836       int j;
837       for (j=m->eBands[i]<<LM;j<m->eBands[i+1]<<LM;j++)
838       {
839          celt_word16 L, R, M, S;
840          L = X[j];
841          R = X[N0+j];
842          M = L+R;
843          S = L-R;
844          sumLR += EXTEND32(ABS16(L)) + EXTEND32(ABS16(R));
845          sumMS += EXTEND32(ABS16(M)) + EXTEND32(ABS16(S));
846       }
847    }
848    sumMS = MULT16_32_Q15(QCONST16(0.707107f, 15), sumMS);
849    thetas = 13;
850    /* We don't need thetas for lower bands with LM<=1 */
851    if (LM<=1)
852       thetas -= 8;
853    return MULT16_32_Q15((m->eBands[13]<<(LM+1))+thetas, sumMS)
854          > MULT16_32_Q15(m->eBands[13]<<(LM+1), sumLR);
855 }
856
857 #ifdef FIXED_POINT
858 int celt_encode_with_ec(CELTEncoder * restrict st, const celt_int16 * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes, ec_enc *enc)
859 {
860 #else
861 int celt_encode_with_ec_float(CELTEncoder * restrict st, const celt_sig * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes, ec_enc *enc)
862 {
863 #endif
864    int i, c, N;
865    int bits;
866    ec_byte_buffer buf;
867    ec_enc         _enc;
868    VARDECL(celt_sig, in);
869    VARDECL(celt_sig, freq);
870    VARDECL(celt_norm, X);
871    VARDECL(celt_ener, bandE);
872    VARDECL(celt_word16, bandLogE);
873    VARDECL(int, fine_quant);
874    VARDECL(celt_word16, error);
875    VARDECL(int, pulses);
876    VARDECL(int, cap);
877    VARDECL(int, offsets);
878    VARDECL(int, fine_priority);
879    VARDECL(int, tf_res);
880    VARDECL(unsigned char, collapse_masks);
881    celt_sig *_overlap_mem;
882    celt_sig *prefilter_mem;
883    celt_word16 *oldBandE, *oldLogE, *oldLogE2;
884    int shortBlocks=0;
885    int isTransient=0;
886    int resynth;
887    const int CC = CHANNELS(st->channels);
888    const int C = CHANNELS(st->stream_channels);
889    int LM, M;
890    int tf_select;
891    int nbFilledBytes, nbAvailableBytes;
892    int effEnd;
893    int codedBands;
894    int tf_sum;
895    int alloc_trim;
896    int pitch_index=COMBFILTER_MINPERIOD;
897    celt_word16 gain1 = 0;
898    int intensity=0;
899    int dual_stereo=0;
900    int effectiveBytes;
901    celt_word16 pf_threshold;
902    int dynalloc_logp;
903    celt_int32 vbr_rate;
904    celt_int32 total_bits;
905    celt_int32 total_boost;
906    celt_int32 balance;
907    celt_int32 tell;
908    int prefilter_tapset=0;
909    int pf_on;
910    int anti_collapse_rsv;
911    int anti_collapse_on=0;
912    int silence=0;
913    SAVE_STACK;
914
915    if (nbCompressedBytes<2 || pcm==NULL)
916      return CELT_BAD_ARG;
917
918    frame_size *= st->upsample;
919    for (LM=0;LM<4;LM++)
920       if (st->mode->shortMdctSize<<LM==frame_size)
921          break;
922    if (LM>=MAX_CONFIG_SIZES)
923       return CELT_BAD_ARG;
924    M=1<<LM;
925    N = M*st->mode->shortMdctSize;
926
927    prefilter_mem = st->in_mem+CC*(st->overlap);
928    _overlap_mem = prefilter_mem+CC*COMBFILTER_MAXPERIOD;
929    /*_overlap_mem = st->in_mem+C*(st->overlap);*/
930    oldBandE = (celt_word16*)(st->in_mem+CC*(2*st->overlap+COMBFILTER_MAXPERIOD));
931    oldLogE = oldBandE + CC*st->mode->nbEBands;
932    oldLogE2 = oldLogE + CC*st->mode->nbEBands;
933
934    if (enc==NULL)
935    {
936       tell=1;
937       nbFilledBytes=0;
938    } else {
939       tell=ec_enc_tell(enc, 0);
940       nbFilledBytes=(tell+4)>>3;
941    }
942    nbAvailableBytes = nbCompressedBytes - nbFilledBytes;
943
944    if (st->vbr)
945    {
946       vbr_rate = ((2*st->bitrate*frame_size<<BITRES)+st->mode->Fs)/(2*st->mode->Fs);
947       effectiveBytes = vbr_rate>>3;
948    } else {
949       celt_int32 tmp;
950       vbr_rate = 0;
951       tmp = st->bitrate*frame_size;
952       if (tell>1)
953          tmp += tell;
954       nbCompressedBytes = IMAX(2, IMIN(nbCompressedBytes,
955             (tmp+4*st->mode->Fs)/(8*st->mode->Fs)));
956       effectiveBytes = nbCompressedBytes;
957    }
958
959    if (enc==NULL)
960    {
961       ec_byte_writeinit_buffer(&buf, compressed, nbCompressedBytes);
962       ec_enc_init(&_enc,&buf);
963       enc = &_enc;
964    }
965
966    if (vbr_rate>0)
967    {
968       /* Computes the max bit-rate allowed in VBR mode to avoid violating the
969           target rate and buffering.
970          We must do this up front so that bust-prevention logic triggers
971           correctly if we don't have enough bits. */
972       if (st->constrained_vbr)
973       {
974          celt_int32 vbr_bound;
975          celt_int32 max_allowed;
976          /* We could use any multiple of vbr_rate as bound (depending on the
977              delay).
978             This is clamped to ensure we use at least two bytes if the encoder
979              was entirely empty, but to allow 0 in hybrid mode. */
980          vbr_bound = vbr_rate;
981          max_allowed = IMIN(IMAX(tell==1?2:0,
982                vbr_rate+vbr_bound-st->vbr_reservoir>>(BITRES+3)),
983                nbAvailableBytes);
984          if(max_allowed < nbAvailableBytes)
985          {
986             nbCompressedBytes = nbFilledBytes+max_allowed;
987             nbAvailableBytes = max_allowed;
988             ec_byte_shrink(enc->buf, nbCompressedBytes);
989          }
990       }
991    }
992    total_bits = nbCompressedBytes*8;
993
994    effEnd = st->end;
995    if (effEnd > st->mode->effEBands)
996       effEnd = st->mode->effEBands;
997
998    ALLOC(in, CC*(N+st->overlap), celt_sig);
999
1000    /* Find pitch period and gain */
1001    {
1002       VARDECL(celt_sig, _pre);
1003       celt_sig *pre[2];
1004       SAVE_STACK;
1005       c = 0;
1006       ALLOC(_pre, CC*(N+COMBFILTER_MAXPERIOD), celt_sig);
1007
1008       pre[0] = _pre;
1009       pre[1] = _pre + (N+COMBFILTER_MAXPERIOD);
1010
1011       silence = 1;
1012       c=0; do {
1013          int count = 0;
1014          const celt_word16 * restrict pcmp = pcm+c;
1015          celt_sig * restrict inp = in+c*(N+st->overlap)+st->overlap;
1016
1017          for (i=0;i<N;i++)
1018          {
1019             celt_sig x, tmp;
1020
1021             x = SCALEIN(*pcmp);
1022             if (++count==st->upsample)
1023             {
1024                count=0;
1025                pcmp+=CC;
1026             } else {
1027                x = 0;
1028             }
1029             /* Apply pre-emphasis */
1030             tmp = MULT16_16(st->mode->preemph[2], x);
1031             *inp = tmp + st->preemph_memE[c];
1032             st->preemph_memE[c] = MULT16_32_Q15(st->mode->preemph[1], *inp)
1033                                    - MULT16_32_Q15(st->mode->preemph[0], tmp);
1034             silence = silence && *inp == 0;
1035             inp++;
1036          }
1037          CELT_COPY(pre[c], prefilter_mem+c*COMBFILTER_MAXPERIOD, COMBFILTER_MAXPERIOD);
1038          CELT_COPY(pre[c]+COMBFILTER_MAXPERIOD, in+c*(N+st->overlap)+st->overlap, N);
1039       } while (++c<CC);
1040
1041       if (tell==1)
1042          ec_enc_bit_logp(enc, silence, 15);
1043       else
1044          silence=0;
1045       if (silence)
1046       {
1047          /*In VBR mode there is no need to send more than the minimum. */
1048          if (vbr_rate>0)
1049          {
1050             effectiveBytes=nbCompressedBytes=IMIN(nbCompressedBytes, nbFilledBytes+2);
1051             total_bits=nbCompressedBytes*8;
1052             nbAvailableBytes=2;
1053             ec_byte_shrink(enc->buf, nbCompressedBytes);
1054          }
1055          /* Pretend we've filled all the remaining bits with zeros
1056             (that's what the initialiser did anyway) */
1057          tell = nbCompressedBytes*8;
1058          enc->nbits_total+=tell-ec_enc_tell(enc,0);
1059       }
1060 #ifdef ENABLE_POSTFILTER
1061       if (nbAvailableBytes>12*C && st->start==0 && !silence)
1062       {
1063          VARDECL(celt_word16, pitch_buf);
1064          ALLOC(pitch_buf, (COMBFILTER_MAXPERIOD+N)>>1, celt_word16);
1065
1066          pitch_downsample(pre, pitch_buf, COMBFILTER_MAXPERIOD+N, CC);
1067          pitch_search(pitch_buf+(COMBFILTER_MAXPERIOD>>1), pitch_buf, N,
1068                COMBFILTER_MAXPERIOD-COMBFILTER_MINPERIOD, &pitch_index);
1069          pitch_index = COMBFILTER_MAXPERIOD-pitch_index;
1070
1071          gain1 = remove_doubling(pitch_buf, COMBFILTER_MAXPERIOD, COMBFILTER_MINPERIOD,
1072                N, &pitch_index, st->prefilter_period, st->prefilter_gain);
1073          if (pitch_index > COMBFILTER_MAXPERIOD-2)
1074             pitch_index = COMBFILTER_MAXPERIOD-2;
1075          gain1 = MULT16_16_Q15(QCONST16(.7f,15),gain1);
1076          prefilter_tapset = st->tapset_decision;
1077       } else {
1078          gain1 = 0;
1079       }
1080
1081       /* Gain threshold for enabling the prefilter/postfilter */
1082       pf_threshold = QCONST16(.2f,15);
1083
1084       /* Adjusting the threshold based on rate and continuity */
1085       if (abs(pitch_index-st->prefilter_period)*10>pitch_index)
1086          pf_threshold += QCONST16(.2f,15);
1087       if (nbAvailableBytes<25)
1088          pf_threshold += QCONST16(.1f,15);
1089       if (nbAvailableBytes<35)
1090          pf_threshold += QCONST16(.1f,15);
1091       if (st->prefilter_gain > QCONST16(.4f,15))
1092          pf_threshold -= QCONST16(.1f,15);
1093       if (st->prefilter_gain > QCONST16(.55f,15))
1094          pf_threshold -= QCONST16(.1f,15);
1095
1096       /* Hard threshold at 0.2 */
1097       pf_threshold = MAX16(pf_threshold, QCONST16(.2f,15));
1098       if (gain1<pf_threshold)
1099       {
1100          if(st->start==0 && tell+17<=total_bits)
1101             ec_enc_bit_logp(enc, 0, 1);
1102          gain1 = 0;
1103          pf_on = 0;
1104       } else {
1105          int qg;
1106          int octave;
1107
1108          if (gain1 > QCONST16(.6f,15))
1109             gain1 = QCONST16(.6f,15);
1110          if (ABS16(gain1-st->prefilter_gain)<QCONST16(.1f,15))
1111             gain1=st->prefilter_gain;
1112
1113 #ifdef FIXED_POINT
1114          qg = ((gain1+2048)>>12)-2;
1115 #else
1116          qg = floor(.5+gain1*8)-2;
1117 #endif
1118          ec_enc_bit_logp(enc, 1, 1);
1119          pitch_index += 1;
1120          octave = EC_ILOG(pitch_index)-5;
1121          ec_enc_uint(enc, octave, 6);
1122          ec_enc_bits(enc, pitch_index-(16<<octave), 4+octave);
1123          pitch_index -= 1;
1124          ec_enc_bits(enc, qg, 2);
1125          gain1 = QCONST16(.125f,15)*(qg+2);
1126          ec_enc_icdf(enc, prefilter_tapset, tapset_icdf, 2);
1127          pf_on = 1;
1128       }
1129       /*printf("%d %f\n", pitch_index, gain1);*/
1130 #else /* ENABLE_POSTFILTER */
1131       if(st->start==0 && tell+17<=total_bits)
1132          ec_enc_bit_logp(enc, 0, 1);
1133       pf_on = 0;
1134 #endif /* ENABLE_POSTFILTER */
1135
1136       c=0; do {
1137          st->prefilter_period=IMAX(st->prefilter_period, COMBFILTER_MINPERIOD);
1138          CELT_COPY(in+c*(N+st->overlap), st->in_mem+c*(st->overlap), st->overlap);
1139 #ifdef ENABLE_POSTFILTER
1140          comb_filter(in+c*(N+st->overlap)+st->overlap, pre[c]+COMBFILTER_MAXPERIOD,
1141                st->prefilter_period, pitch_index, N, -st->prefilter_gain, -gain1,
1142                st->prefilter_tapset, prefilter_tapset, st->mode->window, st->mode->overlap);
1143 #endif /* ENABLE_POSTFILTER */
1144          CELT_COPY(st->in_mem+c*(st->overlap), in+c*(N+st->overlap)+N, st->overlap);
1145
1146 #ifdef ENABLE_POSTFILTER
1147          if (N>COMBFILTER_MAXPERIOD)
1148          {
1149             CELT_MOVE(prefilter_mem+c*COMBFILTER_MAXPERIOD, pre[c]+N, COMBFILTER_MAXPERIOD);
1150          } else {
1151             CELT_MOVE(prefilter_mem+c*COMBFILTER_MAXPERIOD, prefilter_mem+c*COMBFILTER_MAXPERIOD+N, COMBFILTER_MAXPERIOD-N);
1152             CELT_MOVE(prefilter_mem+c*COMBFILTER_MAXPERIOD+COMBFILTER_MAXPERIOD-N, pre[c]+COMBFILTER_MAXPERIOD, N);
1153          }
1154 #endif /* ENABLE_POSTFILTER */
1155       } while (++c<CC);
1156
1157       RESTORE_STACK;
1158    }
1159
1160 #ifdef RESYNTH
1161    resynth = 1;
1162 #else
1163    resynth = 0;
1164 #endif
1165
1166    isTransient = 0;
1167    shortBlocks = 0;
1168    if (LM>0 && ec_enc_tell(enc, 0)+3<=total_bits)
1169    {
1170       if (st->complexity > 1)
1171       {
1172          isTransient = transient_analysis(in, N+st->overlap, CC,
1173                   st->overlap);
1174          if (isTransient)
1175             shortBlocks = M;
1176       }
1177       ec_enc_bit_logp(enc, isTransient, 3);
1178    }
1179
1180    ALLOC(freq, CC*N, celt_sig); /**< Interleaved signal MDCTs */
1181    ALLOC(bandE,st->mode->nbEBands*CC, celt_ener);
1182    ALLOC(bandLogE,st->mode->nbEBands*CC, celt_word16);
1183    /* Compute MDCTs */
1184    compute_mdcts(st->mode, shortBlocks, in, freq, CC, LM);
1185
1186    if (CC==2&&C==1)
1187    {
1188       for (i=0;i<N;i++)
1189          freq[i] = ADD32(HALF32(freq[i]), HALF32(freq[N+i]));
1190    }
1191    if (st->upsample != 1)
1192    {
1193       c=0; do
1194       {
1195          int bound = N/st->upsample;
1196          for (i=0;i<bound;i++)
1197             freq[c*N+i] *= st->upsample;
1198          for (;i<N;i++)
1199             freq[c*N+i] = 0;
1200       } while (++c<C);
1201    }
1202    ALLOC(X, C*N, celt_norm);         /**< Interleaved normalised MDCTs */
1203
1204    compute_band_energies(st->mode, freq, bandE, effEnd, C, M);
1205
1206    amp2Log2(st->mode, effEnd, st->end, bandE, bandLogE, C);
1207
1208    /* Band normalisation */
1209    normalise_bands(st->mode, freq, X, bandE, effEnd, C, M);
1210
1211    ALLOC(tf_res, st->mode->nbEBands, int);
1212    /* Needs to be before coarse energy quantization because otherwise the energy gets modified */
1213    tf_select = tf_analysis(st->mode, bandLogE, oldBandE, effEnd, C, isTransient, tf_res, effectiveBytes, X, N, LM, &tf_sum);
1214    for (i=effEnd;i<st->end;i++)
1215       tf_res[i] = tf_res[effEnd-1];
1216
1217    ALLOC(error, C*st->mode->nbEBands, celt_word16);
1218    quant_coarse_energy(st->mode, st->start, st->end, effEnd, bandLogE,
1219          oldBandE, total_bits, error, enc,
1220          C, LM, nbAvailableBytes, st->force_intra,
1221          &st->delayedIntra, st->complexity >= 4);
1222
1223    tf_encode(st->start, st->end, isTransient, tf_res, LM, tf_select, enc);
1224
1225    st->spread_decision = SPREAD_NORMAL;
1226    if (ec_enc_tell(enc, 0)+4<=total_bits)
1227    {
1228       if (shortBlocks || st->complexity < 3 || nbAvailableBytes < 10*C)
1229       {
1230          if (st->complexity == 0)
1231             st->spread_decision = SPREAD_NONE;
1232       } else {
1233          st->spread_decision = spreading_decision(st->mode, X,
1234                &st->tonal_average, st->spread_decision, &st->hf_average,
1235                &st->tapset_decision, pf_on&&!shortBlocks, effEnd, C, M);
1236       }
1237       ec_enc_icdf(enc, st->spread_decision, spread_icdf, 5);
1238    }
1239
1240    ALLOC(cap, st->mode->nbEBands, int);
1241    ALLOC(offsets, st->mode->nbEBands, int);
1242
1243    init_caps(st->mode,cap,LM,C);
1244    for (i=0;i<st->mode->nbEBands;i++)
1245       offsets[i] = 0;
1246    /* Dynamic allocation code */
1247    /* Make sure that dynamic allocation can't make us bust the budget */
1248    if (effectiveBytes > 50 && LM>=1)
1249    {
1250       int t1, t2;
1251       if (LM <= 1)
1252       {
1253          t1 = 3;
1254          t2 = 5;
1255       } else {
1256          t1 = 2;
1257          t2 = 4;
1258       }
1259       for (i=1;i<st->mode->nbEBands-1;i++)
1260       {
1261          celt_word32 d2;
1262          d2 = 2*bandLogE[i]-bandLogE[i-1]-bandLogE[i+1];
1263          if (C==2)
1264             d2 = HALF32(d2 + 2*bandLogE[i+st->mode->nbEBands]-
1265                   bandLogE[i-1+st->mode->nbEBands]-bandLogE[i+1+st->mode->nbEBands]);
1266          if (d2 > SHL16(t1,DB_SHIFT))
1267             offsets[i] += 1;
1268          if (d2 > SHL16(t2,DB_SHIFT))
1269             offsets[i] += 1;
1270       }
1271    }
1272    dynalloc_logp = 6;
1273    total_bits<<=BITRES;
1274    total_boost = 0;
1275    tell = ec_enc_tell(enc, BITRES);
1276    for (i=st->start;i<st->end;i++)
1277    {
1278       int width, quanta;
1279       int dynalloc_loop_logp;
1280       int boost;
1281       int j;
1282       width = C*(st->mode->eBands[i+1]-st->mode->eBands[i])<<LM;
1283       /* quanta is 6 bits, but no more than 1 bit/sample
1284          and no less than 1/8 bit/sample */
1285       quanta = IMIN(width<<BITRES, IMAX(6<<BITRES, width));
1286       dynalloc_loop_logp = dynalloc_logp;
1287       boost = 0;
1288       for (j = 0; tell+(dynalloc_loop_logp<<BITRES) < total_bits-total_boost
1289             && boost < cap[i]; j++)
1290       {
1291          int flag;
1292          flag = j<offsets[i];
1293          ec_enc_bit_logp(enc, flag, dynalloc_loop_logp);
1294          tell = ec_enc_tell(enc, BITRES);
1295          if (!flag)
1296             break;
1297          boost += quanta;
1298          total_boost += quanta;
1299          dynalloc_loop_logp = 1;
1300       }
1301       /* Making dynalloc more likely */
1302       if (j)
1303          dynalloc_logp = IMAX(2, dynalloc_logp-1);
1304       offsets[i] = boost;
1305    }
1306    alloc_trim = 5;
1307    if (tell+(6<<BITRES) <= total_bits - total_boost)
1308    {
1309       alloc_trim = alloc_trim_analysis(st->mode, X, bandLogE,
1310             st->mode->nbEBands, LM, C, N);
1311       ec_enc_icdf(enc, alloc_trim, trim_icdf, 7);
1312       tell = ec_enc_tell(enc, BITRES);
1313    }
1314
1315    /* Variable bitrate */
1316    if (vbr_rate>0)
1317    {
1318      celt_word16 alpha;
1319      celt_int32 delta;
1320      /* The target rate in 8th bits per frame */
1321      celt_int32 target;
1322      celt_int32 min_allowed;
1323
1324      target = vbr_rate + st->vbr_offset - ((40*C+20)<<BITRES);
1325
1326      /* Shortblocks get a large boost in bitrate, but since they
1327         are uncommon long blocks are not greatly affected */
1328      if (shortBlocks || tf_sum < -2*(st->end-st->start))
1329         target = 7*target/4;
1330      else if (tf_sum < -(st->end-st->start))
1331         target = 3*target/2;
1332      else if (M > 1)
1333         target-=(target+14)/28;
1334
1335      /* The current offset is removed from the target and the space used
1336         so far is added*/
1337      target=target+tell;
1338
1339      /* In VBR mode the frame size must not be reduced so much that it would
1340          result in the encoder running out of bits.
1341         The margin of 2 bytes ensures that none of the bust-prevention logic
1342          in the decoder will have triggered so far. */
1343      min_allowed = (tell+total_boost+(1<<BITRES+3)-1>>(BITRES+3)) + 2 - nbFilledBytes;
1344
1345      nbAvailableBytes = target+(1<<(BITRES+2))>>(BITRES+3);
1346      nbAvailableBytes = IMAX(min_allowed,nbAvailableBytes);
1347      nbAvailableBytes = IMIN(nbCompressedBytes,nbAvailableBytes+nbFilledBytes) - nbFilledBytes;
1348
1349      if(silence)
1350      {
1351        nbAvailableBytes = 2;
1352        target = 2*8<<BITRES;
1353      }
1354
1355      /* By how much did we "miss" the target on that frame */
1356      delta = target - vbr_rate;
1357
1358      target=nbAvailableBytes<<(BITRES+3);
1359
1360      if (st->vbr_count < 970)
1361      {
1362         st->vbr_count++;
1363         alpha = celt_rcp(SHL32(EXTEND32(st->vbr_count+20),16));
1364      } else
1365         alpha = QCONST16(.001f,15);
1366      /* How many bits have we used in excess of what we're allowed */
1367      if (st->constrained_vbr)
1368         st->vbr_reservoir += target - vbr_rate;
1369      /*printf ("%d\n", st->vbr_reservoir);*/
1370
1371      /* Compute the offset we need to apply in order to reach the target */
1372      st->vbr_drift += (celt_int32)MULT16_32_Q15(alpha,delta-st->vbr_offset-st->vbr_drift);
1373      st->vbr_offset = -st->vbr_drift;
1374      /*printf ("%d\n", st->vbr_drift);*/
1375
1376      if (st->constrained_vbr && st->vbr_reservoir < 0)
1377      {
1378         /* We're under the min value -- increase rate */
1379         int adjust = (-st->vbr_reservoir)/(8<<BITRES);
1380         /* Unless we're just coding silence */
1381         nbAvailableBytes += silence?0:adjust;
1382         st->vbr_reservoir = 0;
1383         /*printf ("+%d\n", adjust);*/
1384      }
1385      nbCompressedBytes = IMIN(nbCompressedBytes,nbAvailableBytes+nbFilledBytes);
1386      /* This moves the raw bits to take into account the new compressed size */
1387      ec_byte_shrink(enc->buf, nbCompressedBytes);
1388    }
1389    if (C==2)
1390    {
1391       int effectiveRate;
1392
1393       /* Always use MS for 2.5 ms frames until we can do a better analysis */
1394       if (LM!=0)
1395          dual_stereo = stereo_analysis(st->mode, X, LM, N);
1396
1397       /* Account for coarse energy */
1398       effectiveRate = (8*effectiveBytes - 80)>>LM;
1399
1400       /* effectiveRate in kb/s */
1401       effectiveRate = 2*effectiveRate/5;
1402       if (effectiveRate<35)
1403          intensity = 8;
1404       else if (effectiveRate<50)
1405          intensity = 12;
1406       else if (effectiveRate<68)
1407          intensity = 16;
1408       else if (effectiveRate<84)
1409          intensity = 18;
1410       else if (effectiveRate<102)
1411          intensity = 19;
1412       else if (effectiveRate<130)
1413          intensity = 20;
1414       else
1415          intensity = 100;
1416       intensity = IMIN(st->end,IMAX(st->start, intensity));
1417    }
1418
1419    /* Bit allocation */
1420    ALLOC(fine_quant, st->mode->nbEBands, int);
1421    ALLOC(pulses, st->mode->nbEBands, int);
1422    ALLOC(fine_priority, st->mode->nbEBands, int);
1423
1424    /* bits =   packet size        -       where we are         - safety*/
1425    bits = (nbCompressedBytes*8<<BITRES) - ec_enc_tell(enc, BITRES) - 1;
1426    anti_collapse_rsv = isTransient&&LM>=2&&bits>=(LM+2<<BITRES) ? (1<<BITRES) : 0;
1427    bits -= anti_collapse_rsv;
1428    codedBands = compute_allocation(st->mode, st->start, st->end, offsets, cap,
1429          alloc_trim, &intensity, &dual_stereo, bits, &balance, pulses,
1430          fine_quant, fine_priority, C, LM, enc, 1, st->lastCodedBands);
1431    st->lastCodedBands = codedBands;
1432
1433    quant_fine_energy(st->mode, st->start, st->end, oldBandE, error, fine_quant, enc, C);
1434
1435 #ifdef MEASURE_NORM_MSE
1436    float X0[3000];
1437    float bandE0[60];
1438    c=0; do 
1439       for (i=0;i<N;i++)
1440          X0[i+c*N] = X[i+c*N];
1441    while (++c<C);
1442    for (i=0;i<C*st->mode->nbEBands;i++)
1443       bandE0[i] = bandE[i];
1444 #endif
1445
1446    /* Residual quantisation */
1447    ALLOC(collapse_masks, C*st->mode->nbEBands, unsigned char);
1448    quant_all_bands(1, st->mode, st->start, st->end, X, C==2 ? X+N : NULL, collapse_masks,
1449          bandE, pulses, shortBlocks, st->spread_decision, dual_stereo, intensity, tf_res, resynth,
1450          nbCompressedBytes*(8<<BITRES)-anti_collapse_rsv, balance, enc, LM, codedBands, &st->rng);
1451
1452    if (anti_collapse_rsv > 0)
1453    {
1454       anti_collapse_on = st->consec_transient<2;
1455       ec_enc_bits(enc, anti_collapse_on, 1);
1456    }
1457    quant_energy_finalise(st->mode, st->start, st->end, oldBandE, error, fine_quant, fine_priority, nbCompressedBytes*8-ec_enc_tell(enc, 0), enc, C);
1458
1459    if (silence)
1460    {
1461       for (i=0;i<C*st->mode->nbEBands;i++)
1462          oldBandE[i] = -QCONST16(28.f,DB_SHIFT);
1463    }
1464
1465 #ifdef RESYNTH
1466    /* Re-synthesis of the coded audio if required */
1467    if (resynth)
1468    {
1469       celt_sig *out_mem[2];
1470       celt_sig *overlap_mem[2];
1471
1472       log2Amp(st->mode, st->start, st->end, bandE, oldBandE, C);
1473       if (silence)
1474       {
1475          for (i=0;i<C*st->mode->nbEBands;i++)
1476             bandE[i] = 0;
1477       }
1478
1479 #ifdef MEASURE_NORM_MSE
1480       measure_norm_mse(st->mode, X, X0, bandE, bandE0, M, N, C);
1481 #endif
1482       if (anti_collapse_on)
1483       {
1484          anti_collapse(st->mode, X, collapse_masks, LM, C, CC, N,
1485                st->start, st->end, oldBandE, oldLogE, oldLogE2, pulses, st->rng);
1486       }
1487
1488       /* Synthesis */
1489       denormalise_bands(st->mode, X, freq, bandE, effEnd, C, M);
1490
1491       CELT_MOVE(st->syn_mem[0], st->syn_mem[0]+N, MAX_PERIOD);
1492       if (CC==2)
1493          CELT_MOVE(st->syn_mem[1], st->syn_mem[1]+N, MAX_PERIOD);
1494
1495       c=0; do
1496          for (i=0;i<M*st->mode->eBands[st->start];i++)
1497             freq[c*N+i] = 0;
1498       while (++c<C);
1499       c=0; do
1500          for (i=M*st->mode->eBands[st->end];i<N;i++)
1501             freq[c*N+i] = 0;
1502       while (++c<C);
1503
1504       if (CC==2&&C==1)
1505       {
1506          for (i=0;i<N;i++)
1507             freq[N+i] = freq[i];
1508       }
1509
1510       out_mem[0] = st->syn_mem[0]+MAX_PERIOD;
1511       if (CC==2)
1512          out_mem[1] = st->syn_mem[1]+MAX_PERIOD;
1513
1514       c=0; do
1515          overlap_mem[c] = _overlap_mem + c*st->overlap;
1516       while (++c<CC);
1517
1518       compute_inv_mdcts(st->mode, shortBlocks, freq, out_mem, overlap_mem, CC, LM);
1519
1520 #ifdef ENABLE_POSTFILTER
1521       c=0; do {
1522          st->prefilter_period=IMAX(st->prefilter_period, COMBFILTER_MINPERIOD);
1523          st->prefilter_period_old=IMAX(st->prefilter_period_old, COMBFILTER_MINPERIOD);
1524          if (LM!=0)
1525          {
1526             comb_filter(out_mem[c], out_mem[c], st->prefilter_period, st->prefilter_period, st->overlap,
1527                   st->prefilter_gain, st->prefilter_gain, st->prefilter_tapset, st->prefilter_tapset,
1528                   NULL, 0);
1529             comb_filter(out_mem[c]+st->overlap, out_mem[c]+st->overlap, st->prefilter_period, pitch_index, N-st->overlap,
1530                   st->prefilter_gain, gain1, st->prefilter_tapset, prefilter_tapset,
1531                   st->mode->window, st->mode->overlap);
1532          } else {
1533             comb_filter(out_mem[c], out_mem[c], st->prefilter_period_old, st->prefilter_period, N,
1534                   st->prefilter_gain_old, st->prefilter_gain, st->prefilter_tapset_old, st->prefilter_tapset,
1535                   st->mode->window, st->mode->overlap);
1536          }
1537       } while (++c<CC);
1538 #endif /* ENABLE_POSTFILTER */
1539
1540       deemphasis(out_mem, (celt_word16*)pcm, N, CC, st->upsample, st->mode->preemph, st->preemph_memD);
1541       st->prefilter_period_old = st->prefilter_period;
1542       st->prefilter_gain_old = st->prefilter_gain;
1543       st->prefilter_tapset_old = st->prefilter_tapset;
1544    }
1545 #endif
1546
1547    st->prefilter_period = pitch_index;
1548    st->prefilter_gain = gain1;
1549    st->prefilter_tapset = prefilter_tapset;
1550
1551    if (CC==2&&C==1) {
1552       for (i=0;i<st->mode->nbEBands;i++)
1553          oldBandE[st->mode->nbEBands+i]=oldBandE[i];
1554    }
1555
1556    /* In case start or end were to change */
1557    c=0; do
1558    {
1559       for (i=0;i<st->start;i++)
1560          oldBandE[c*st->mode->nbEBands+i]=0;
1561       for (i=st->end;i<st->mode->nbEBands;i++)
1562          oldBandE[c*st->mode->nbEBands+i]=0;
1563    } while (++c<CC);
1564    if (!isTransient)
1565    {
1566       for (i=0;i<CC*st->mode->nbEBands;i++)
1567          oldLogE2[i] = oldLogE[i];
1568       for (i=0;i<CC*st->mode->nbEBands;i++)
1569          oldLogE[i] = oldBandE[i];
1570    } else {
1571       for (i=0;i<CC*st->mode->nbEBands;i++)
1572          oldLogE[i] = MIN16(oldLogE[i], oldBandE[i]);
1573    }
1574    if (isTransient)
1575       st->consec_transient++;
1576    else
1577       st->consec_transient=0;
1578    st->rng = enc->rng;
1579
1580    /* If there's any room left (can only happen for very high rates),
1581       it's already filled with zeros */
1582    ec_enc_done(enc);
1583    
1584    RESTORE_STACK;
1585    if (ec_enc_get_error(enc))
1586       return CELT_CORRUPTED_DATA;
1587    else
1588       return nbCompressedBytes;
1589 }
1590
1591 #ifdef FIXED_POINT
1592 #ifndef DISABLE_FLOAT_API
1593 int celt_encode_with_ec_float(CELTEncoder * restrict st, const float * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes, ec_enc *enc)
1594 {
1595    int j, ret, C, N;
1596    VARDECL(celt_int16, in);
1597    ALLOC_STACK;
1598    SAVE_STACK;
1599
1600    if (pcm==NULL)
1601       return CELT_BAD_ARG;
1602
1603    C = CHANNELS(st->channels);
1604    N = frame_size;
1605    ALLOC(in, C*N, celt_int16);
1606
1607    for (j=0;j<C*N;j++)
1608      in[j] = FLOAT2INT16(pcm[j]);
1609
1610    ret=celt_encode_with_ec(st,in,frame_size,compressed,nbCompressedBytes, enc);
1611 #ifdef RESYNTH
1612    for (j=0;j<C*N;j++)
1613       ((float*)pcm)[j]=in[j]*(1.f/32768.f);
1614 #endif
1615    RESTORE_STACK;
1616    return ret;
1617
1618 }
1619 #endif /*DISABLE_FLOAT_API*/
1620 #else
1621 int celt_encode_with_ec(CELTEncoder * restrict st, const celt_int16 * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes, ec_enc *enc)
1622 {
1623    int j, ret, C, N;
1624    VARDECL(celt_sig, in);
1625    ALLOC_STACK;
1626    SAVE_STACK;
1627
1628    if (pcm==NULL)
1629       return CELT_BAD_ARG;
1630
1631    C=CHANNELS(st->channels);
1632    N=frame_size;
1633    ALLOC(in, C*N, celt_sig);
1634    for (j=0;j<C*N;j++) {
1635      in[j] = SCALEOUT(pcm[j]);
1636    }
1637
1638    ret = celt_encode_with_ec_float(st,in,frame_size,compressed,nbCompressedBytes, enc);
1639 #ifdef RESYNTH
1640    for (j=0;j<C*N;j++)
1641       ((celt_int16*)pcm)[j] = FLOAT2INT16(in[j]);
1642 #endif
1643    RESTORE_STACK;
1644    return ret;
1645 }
1646 #endif
1647
1648 int celt_encode(CELTEncoder * restrict st, const celt_int16 * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes)
1649 {
1650    return celt_encode_with_ec(st, pcm, frame_size, compressed, nbCompressedBytes, NULL);
1651 }
1652
1653 #ifndef DISABLE_FLOAT_API
1654 int celt_encode_float(CELTEncoder * restrict st, const float * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes)
1655 {
1656    return celt_encode_with_ec_float(st, pcm, frame_size, compressed, nbCompressedBytes, NULL);
1657 }
1658 #endif /* DISABLE_FLOAT_API */
1659
1660 int celt_encoder_ctl(CELTEncoder * restrict st, int request, ...)
1661 {
1662    va_list ap;
1663    
1664    va_start(ap, request);
1665    switch (request)
1666    {
1667       case CELT_GET_MODE_REQUEST:
1668       {
1669          const CELTMode ** value = va_arg(ap, const CELTMode**);
1670          if (value==0)
1671             goto bad_arg;
1672          *value=st->mode;
1673       }
1674       break;
1675       case CELT_SET_COMPLEXITY_REQUEST:
1676       {
1677          int value = va_arg(ap, celt_int32);
1678          if (value<0 || value>10)
1679             goto bad_arg;
1680          st->complexity = value;
1681       }
1682       break;
1683       case CELT_SET_START_BAND_REQUEST:
1684       {
1685          celt_int32 value = va_arg(ap, celt_int32);
1686          if (value<0 || value>=st->mode->nbEBands)
1687             goto bad_arg;
1688          st->start = value;
1689       }
1690       break;
1691       case CELT_SET_END_BAND_REQUEST:
1692       {
1693          celt_int32 value = va_arg(ap, celt_int32);
1694          if (value<1 || value>st->mode->nbEBands)
1695             goto bad_arg;
1696          st->end = value;
1697       }
1698       break;
1699       case CELT_SET_PREDICTION_REQUEST:
1700       {
1701          int value = va_arg(ap, celt_int32);
1702          if (value<0 || value>2)
1703             goto bad_arg;
1704          if (value==0)
1705          {
1706             st->force_intra   = 1;
1707          } else if (value==1) {
1708             st->force_intra   = 0;
1709          } else {
1710             st->force_intra   = 0;
1711          }   
1712       }
1713       break;
1714       case CELT_SET_VBR_CONSTRAINT_REQUEST:
1715       {
1716          celt_int32 value = va_arg(ap, celt_int32);
1717          st->constrained_vbr = value;
1718       }
1719       break;
1720       case CELT_SET_VBR_REQUEST:
1721       {
1722          celt_int32 value = va_arg(ap, celt_int32);
1723          st->vbr = value;
1724       }
1725       break;
1726       case CELT_SET_BITRATE_REQUEST:
1727       {
1728          celt_int32 value = va_arg(ap, celt_int32);
1729          if (value<=500)
1730             goto bad_arg;
1731          st->bitrate = value;
1732       }
1733       break;
1734       case CELT_SET_CHANNELS_REQUEST:
1735       {
1736          celt_int32 value = va_arg(ap, celt_int32);
1737          if (value<1 || value>2)
1738             goto bad_arg;
1739          st->stream_channels = value;
1740       }
1741       break;
1742       case CELT_RESET_STATE:
1743       {
1744          CELT_MEMSET((char*)&st->ENCODER_RESET_START, 0,
1745                celt_encoder_get_size_custom(st->mode, st->channels)-
1746                ((char*)&st->ENCODER_RESET_START - (char*)st));
1747          st->vbr_offset = 0;
1748          st->delayedIntra = 1;
1749          st->spread_decision = SPREAD_NORMAL;
1750          st->tonal_average = QCONST16(1.f,8);
1751       }
1752       break;
1753       default:
1754          goto bad_request;
1755    }
1756    va_end(ap);
1757    return CELT_OK;
1758 bad_arg:
1759    va_end(ap);
1760    return CELT_BAD_ARG;
1761 bad_request:
1762    va_end(ap);
1763    return CELT_UNIMPLEMENTED;
1764 }
1765
1766 /**********************************************************************/
1767 /*                                                                    */
1768 /*                             DECODER                                */
1769 /*                                                                    */
1770 /**********************************************************************/
1771 #define DECODE_BUFFER_SIZE 2048
1772
1773 /** Decoder state 
1774  @brief Decoder state
1775  */
1776 struct CELTDecoder {
1777    const CELTMode *mode;
1778    int overlap;
1779    int channels;
1780    int stream_channels;
1781
1782    int downsample;
1783    int start, end;
1784
1785    /* Everything beyond this point gets cleared on a reset */
1786 #define DECODER_RESET_START rng
1787
1788    ec_uint32 rng;
1789    int last_pitch_index;
1790    int loss_count;
1791    int postfilter_period;
1792    int postfilter_period_old;
1793    celt_word16 postfilter_gain;
1794    celt_word16 postfilter_gain_old;
1795    int postfilter_tapset;
1796    int postfilter_tapset_old;
1797
1798    celt_sig preemph_memD[2];
1799    
1800    celt_sig _decode_mem[1]; /* Size = channels*(DECODE_BUFFER_SIZE+mode->overlap) */
1801    /* celt_word16 lpc[],  Size = channels*LPC_ORDER */
1802    /* celt_word16 oldEBands[], Size = channels*mode->nbEBands */
1803    /* celt_word16 oldLogE[], Size = channels*mode->nbEBands */
1804    /* celt_word16 oldLogE2[], Size = channels*mode->nbEBands */
1805    /* celt_word16 backgroundLogE[], Size = channels*mode->nbEBands */
1806 };
1807
1808 int celt_decoder_get_size(int channels)
1809 {
1810    const CELTMode *mode = celt_mode_create(48000, 960, NULL);
1811    return celt_decoder_get_size_custom(mode, channels);
1812 }
1813
1814 int celt_decoder_get_size_custom(const CELTMode *mode, int channels)
1815 {
1816    int size = sizeof(struct CELTDecoder)
1817             + (channels*(DECODE_BUFFER_SIZE+mode->overlap)-1)*sizeof(celt_sig)
1818             + channels*LPC_ORDER*sizeof(celt_word16)
1819             + 4*channels*mode->nbEBands*sizeof(celt_word16);
1820    return size;
1821 }
1822
1823 CELTDecoder *celt_decoder_create(int sampling_rate, int channels, int *error)
1824 {
1825    CELTDecoder *st;
1826    st = (CELTDecoder *)celt_alloc(celt_decoder_get_size(channels));
1827    if (st!=NULL && celt_decoder_init(st, sampling_rate, channels, error)==NULL)
1828    {
1829       celt_decoder_destroy(st);
1830       st = NULL;
1831    }
1832    return st;
1833 }
1834
1835 CELTDecoder *celt_decoder_create_custom(const CELTMode *mode, int channels, int *error)
1836 {
1837    CELTDecoder *st = (CELTDecoder *)celt_alloc(celt_decoder_get_size_custom(mode, channels));
1838    if (st!=NULL && celt_decoder_init_custom(st, mode, channels, error)==NULL)
1839    {
1840       celt_decoder_destroy(st);
1841       st = NULL;
1842    }
1843    return st;
1844 }
1845
1846 CELTDecoder *celt_decoder_init(CELTDecoder *st, int sampling_rate, int channels, int *error)
1847 {
1848    celt_decoder_init_custom(st, celt_mode_create(48000, 960, NULL), channels, error);
1849    st->downsample = resampling_factor(sampling_rate);
1850    if (st->downsample==0)
1851    {
1852       if (error)
1853          *error = CELT_BAD_ARG;
1854       return NULL;
1855    }
1856    return st;
1857 }
1858
1859 CELTDecoder *celt_decoder_init_custom(CELTDecoder *st, const CELTMode *mode, int channels, int *error)
1860 {
1861    if (channels < 0 || channels > 2)
1862    {
1863       if (error)
1864          *error = CELT_BAD_ARG;
1865       return NULL;
1866    }
1867
1868    if (st==NULL)
1869    {
1870       if (error)
1871          *error = CELT_ALLOC_FAIL;
1872       return NULL;
1873    }
1874
1875    CELT_MEMSET((char*)st, 0, celt_decoder_get_size_custom(mode, channels));
1876
1877    st->mode = mode;
1878    st->overlap = mode->overlap;
1879    st->stream_channels = st->channels = channels;
1880
1881    st->downsample = 1;
1882    st->start = 0;
1883    st->end = st->mode->effEBands;
1884
1885    st->loss_count = 0;
1886
1887    if (error)
1888       *error = CELT_OK;
1889    return st;
1890 }
1891
1892 void celt_decoder_destroy(CELTDecoder *st)
1893 {
1894    celt_free(st);
1895 }
1896
1897 static void celt_decode_lost(CELTDecoder * restrict st, celt_word16 * restrict pcm, int N, int LM)
1898 {
1899    int c;
1900    int pitch_index;
1901    int overlap = st->mode->overlap;
1902    celt_word16 fade = Q15ONE;
1903    int i, len;
1904    const int C = CHANNELS(st->channels);
1905    int offset;
1906    celt_sig *out_mem[2];
1907    celt_sig *decode_mem[2];
1908    celt_sig *overlap_mem[2];
1909    celt_word16 *lpc;
1910    celt_word32 *out_syn[2];
1911    celt_word16 *oldBandE, *oldLogE2, *backgroundLogE;
1912    SAVE_STACK;
1913    
1914    c=0; do {
1915       decode_mem[c] = st->_decode_mem + c*(DECODE_BUFFER_SIZE+st->overlap);
1916       out_mem[c] = decode_mem[c]+DECODE_BUFFER_SIZE-MAX_PERIOD;
1917       overlap_mem[c] = decode_mem[c]+DECODE_BUFFER_SIZE;
1918    } while (++c<C);
1919    lpc = (celt_word16*)(st->_decode_mem+(DECODE_BUFFER_SIZE+st->overlap)*C);
1920    oldBandE = lpc+C*LPC_ORDER;
1921    oldLogE2 = oldBandE + C*st->mode->nbEBands;
1922    backgroundLogE = oldLogE2  + C*st->mode->nbEBands;
1923
1924    out_syn[0] = out_mem[0]+MAX_PERIOD-N;
1925    if (C==2)
1926       out_syn[1] = out_mem[1]+MAX_PERIOD-N;
1927
1928    len = N+st->mode->overlap;
1929    
1930    if (st->loss_count >= 5)
1931    {
1932       VARDECL(celt_sig, freq);
1933       VARDECL(celt_norm, X);
1934       VARDECL(celt_ener, bandE);
1935       celt_uint32 seed;
1936
1937       ALLOC(freq, C*N, celt_sig); /**< Interleaved signal MDCTs */
1938       ALLOC(X, C*N, celt_norm);   /**< Interleaved normalised MDCTs */
1939       ALLOC(bandE, st->mode->nbEBands*C, celt_ener);
1940
1941       log2Amp(st->mode, st->start, st->end, bandE, backgroundLogE, C);
1942
1943       seed = st->rng;
1944       for (i=0;i<C*N;i++)
1945       {
1946             seed = lcg_rand(seed);
1947             X[i] = (celt_int32)(seed)>>20;
1948       }
1949       st->rng = seed;
1950       for (c=0;c<C;c++)
1951          for (i=0;i<st->mode->nbEBands;i++)
1952             renormalise_vector(X+N*c+(st->mode->eBands[i]<<LM), (st->mode->eBands[i+1]-st->mode->eBands[i])<<LM, Q15ONE);
1953
1954       denormalise_bands(st->mode, X, freq, bandE, st->mode->nbEBands, C, 1<<LM);
1955
1956       compute_inv_mdcts(st->mode, 0, freq, out_syn, overlap_mem, C, LM);
1957    } else if (st->loss_count == 0)
1958    {
1959       celt_word16 pitch_buf[MAX_PERIOD>>1];
1960       int len2 = len;
1961       /* FIXME: This is a kludge */
1962       if (len2>MAX_PERIOD>>1)
1963          len2 = MAX_PERIOD>>1;
1964       pitch_downsample(out_mem, pitch_buf, MAX_PERIOD, C);
1965       pitch_search(pitch_buf+((MAX_PERIOD-len2)>>1), pitch_buf, len2,
1966                    MAX_PERIOD-len2-100, &pitch_index);
1967       pitch_index = MAX_PERIOD-len2-pitch_index;
1968       st->last_pitch_index = pitch_index;
1969    } else {
1970       pitch_index = st->last_pitch_index;
1971       fade = QCONST16(.8f,15);
1972    }
1973
1974    c=0; do {
1975       /* FIXME: This is more memory than necessary */
1976       celt_word32 e[2*MAX_PERIOD];
1977       celt_word16 exc[2*MAX_PERIOD];
1978       celt_word32 ac[LPC_ORDER+1];
1979       celt_word16 decay = 1;
1980       celt_word32 S1=0;
1981       celt_word16 mem[LPC_ORDER]={0};
1982
1983       offset = MAX_PERIOD-pitch_index;
1984       for (i=0;i<MAX_PERIOD;i++)
1985          exc[i] = ROUND16(out_mem[c][i], SIG_SHIFT);
1986
1987       if (st->loss_count == 0)
1988       {
1989          _celt_autocorr(exc, ac, st->mode->window, st->mode->overlap,
1990                         LPC_ORDER, MAX_PERIOD);
1991
1992          /* Noise floor -40 dB */
1993 #ifdef FIXED_POINT
1994          ac[0] += SHR32(ac[0],13);
1995 #else
1996          ac[0] *= 1.0001f;
1997 #endif
1998          /* Lag windowing */
1999          for (i=1;i<=LPC_ORDER;i++)
2000          {
2001             /*ac[i] *= exp(-.5*(2*M_PI*.002*i)*(2*M_PI*.002*i));*/
2002 #ifdef FIXED_POINT
2003             ac[i] -= MULT16_32_Q15(2*i*i, ac[i]);
2004 #else
2005             ac[i] -= ac[i]*(.008f*i)*(.008f*i);
2006 #endif
2007          }
2008
2009          _celt_lpc(lpc+c*LPC_ORDER, ac, LPC_ORDER);
2010       }
2011       for (i=0;i<LPC_ORDER;i++)
2012          mem[i] = ROUND16(out_mem[c][MAX_PERIOD-1-i], SIG_SHIFT);
2013       fir(exc, lpc+c*LPC_ORDER, exc, MAX_PERIOD, LPC_ORDER, mem);
2014       /*for (i=0;i<MAX_PERIOD;i++)printf("%d ", exc[i]); printf("\n");*/
2015       /* Check if the waveform is decaying (and if so how fast) */
2016       {
2017          celt_word32 E1=1, E2=1;
2018          int period;
2019          if (pitch_index <= MAX_PERIOD/2)
2020             period = pitch_index;
2021          else
2022             period = MAX_PERIOD/2;
2023          for (i=0;i<period;i++)
2024          {
2025             E1 += SHR32(MULT16_16(exc[MAX_PERIOD-period+i],exc[MAX_PERIOD-period+i]),8);
2026             E2 += SHR32(MULT16_16(exc[MAX_PERIOD-2*period+i],exc[MAX_PERIOD-2*period+i]),8);
2027          }
2028          if (E1 > E2)
2029             E1 = E2;
2030          decay = celt_sqrt(frac_div32(SHR(E1,1),E2));
2031       }
2032
2033       /* Copy excitation, taking decay into account */
2034       for (i=0;i<len+st->mode->overlap;i++)
2035       {
2036          celt_word16 tmp;
2037          if (offset+i >= MAX_PERIOD)
2038          {
2039             offset -= pitch_index;
2040             decay = MULT16_16_Q15(decay, decay);
2041          }
2042          e[i] = SHL32(EXTEND32(MULT16_16_Q15(decay, exc[offset+i])), SIG_SHIFT);
2043          tmp = ROUND16(out_mem[c][offset+i],SIG_SHIFT);
2044          S1 += SHR32(MULT16_16(tmp,tmp),8);
2045       }
2046       for (i=0;i<LPC_ORDER;i++)
2047          mem[i] = ROUND16(out_mem[c][MAX_PERIOD-1-i], SIG_SHIFT);
2048       for (i=0;i<len+st->mode->overlap;i++)
2049          e[i] = MULT16_32_Q15(fade, e[i]);
2050       iir(e, lpc+c*LPC_ORDER, e, len+st->mode->overlap, LPC_ORDER, mem);
2051
2052       {
2053          celt_word32 S2=0;
2054          for (i=0;i<len+overlap;i++)
2055          {
2056             celt_word16 tmp = ROUND16(e[i],SIG_SHIFT);
2057             S2 += SHR32(MULT16_16(tmp,tmp),8);
2058          }
2059          /* This checks for an "explosion" in the synthesis */
2060 #ifdef FIXED_POINT
2061          if (!(S1 > SHR32(S2,2)))
2062 #else
2063          /* Float test is written this way to catch NaNs at the same time */
2064          if (!(S1 > 0.2f*S2))
2065 #endif
2066          {
2067             for (i=0;i<len+overlap;i++)
2068                e[i] = 0;
2069          } else if (S1 < S2)
2070          {
2071             celt_word16 ratio = celt_sqrt(frac_div32(SHR32(S1,1)+1,S2+1));
2072             for (i=0;i<len+overlap;i++)
2073                e[i] = MULT16_32_Q15(ratio, e[i]);
2074          }
2075       }
2076
2077 #ifdef ENABLE_POSTFILTER
2078       /* Apply post-filter to the MDCT overlap of the previous frame */
2079       comb_filter(out_mem[c]+MAX_PERIOD, out_mem[c]+MAX_PERIOD, st->postfilter_period, st->postfilter_period, st->overlap,
2080                   st->postfilter_gain, st->postfilter_gain, st->postfilter_tapset, st->postfilter_tapset,
2081                   NULL, 0);
2082 #endif /* ENABLE_POSTFILTER */
2083
2084       for (i=0;i<MAX_PERIOD+st->mode->overlap-N;i++)
2085          out_mem[c][i] = out_mem[c][N+i];
2086
2087       /* Apply TDAC to the concealed audio so that it blends with the
2088          previous and next frames */
2089       for (i=0;i<overlap/2;i++)
2090       {
2091          celt_word32 tmp;
2092          tmp = MULT16_32_Q15(st->mode->window[i],           e[N+overlap-1-i]) +
2093                MULT16_32_Q15(st->mode->window[overlap-i-1], e[N+i          ]);
2094          out_mem[c][MAX_PERIOD+i] = MULT16_32_Q15(st->mode->window[overlap-i-1], tmp);
2095          out_mem[c][MAX_PERIOD+overlap-i-1] = MULT16_32_Q15(st->mode->window[i], tmp);
2096       }
2097       for (i=0;i<N;i++)
2098          out_mem[c][MAX_PERIOD-N+i] = e[i];
2099
2100 #ifdef ENABLE_POSTFILTER
2101       /* Apply pre-filter to the MDCT overlap for the next frame (post-filter will be applied then) */
2102       comb_filter(e, out_mem[c]+MAX_PERIOD, st->postfilter_period, st->postfilter_period, st->overlap,
2103                   -st->postfilter_gain, -st->postfilter_gain, st->postfilter_tapset, st->postfilter_tapset,
2104                   NULL, 0);
2105 #endif /* ENABLE_POSTFILTER */
2106       for (i=0;i<overlap;i++)
2107          out_mem[c][MAX_PERIOD+i] = e[i];
2108    } while (++c<C);
2109
2110    deemphasis(out_syn, pcm, N, C, st->downsample, st->mode->preemph, st->preemph_memD);
2111    
2112    st->loss_count++;
2113
2114    RESTORE_STACK;
2115 }
2116
2117 #ifdef FIXED_POINT
2118 int celt_decode_with_ec(CELTDecoder * restrict st, const unsigned char *data, int len, celt_int16 * restrict pcm, int frame_size, ec_dec *dec)
2119 {
2120 #else
2121 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)
2122 {
2123 #endif
2124    int c, i, N;
2125    int spread_decision;
2126    int bits;
2127    ec_dec _dec;
2128    ec_byte_buffer buf;
2129    VARDECL(celt_sig, freq);
2130    VARDECL(celt_norm, X);
2131    VARDECL(celt_ener, bandE);
2132    VARDECL(int, fine_quant);
2133    VARDECL(int, pulses);
2134    VARDECL(int, cap);
2135    VARDECL(int, offsets);
2136    VARDECL(int, fine_priority);
2137    VARDECL(int, tf_res);
2138    VARDECL(unsigned char, collapse_masks);
2139    celt_sig *out_mem[2];
2140    celt_sig *decode_mem[2];
2141    celt_sig *overlap_mem[2];
2142    celt_sig *out_syn[2];
2143    celt_word16 *lpc;
2144    celt_word16 *oldBandE, *oldLogE, *oldLogE2, *backgroundLogE;
2145
2146    int shortBlocks;
2147    int isTransient;
2148    int intra_ener;
2149    const int CC = CHANNELS(st->channels);
2150    int LM, M;
2151    int effEnd;
2152    int codedBands;
2153    int alloc_trim;
2154    int postfilter_pitch;
2155    celt_word16 postfilter_gain;
2156    int intensity=0;
2157    int dual_stereo=0;
2158    celt_int32 total_bits;
2159    celt_int32 balance;
2160    celt_int32 tell;
2161    int dynalloc_logp;
2162    int postfilter_tapset;
2163    int anti_collapse_rsv;
2164    int anti_collapse_on=0;
2165    int silence;
2166    const int C = CHANNELS(st->stream_channels);
2167
2168    SAVE_STACK;
2169
2170    if (pcm==NULL)
2171       return CELT_BAD_ARG;
2172
2173    frame_size *= st->downsample;
2174    for (LM=0;LM<4;LM++)
2175       if (st->mode->shortMdctSize<<LM==frame_size)
2176          break;
2177    if (LM>=MAX_CONFIG_SIZES)
2178       return CELT_BAD_ARG;
2179    M=1<<LM;
2180
2181    c=0; do {
2182       decode_mem[c] = st->_decode_mem + c*(DECODE_BUFFER_SIZE+st->overlap);
2183       out_mem[c] = decode_mem[c]+DECODE_BUFFER_SIZE-MAX_PERIOD;
2184       overlap_mem[c] = decode_mem[c]+DECODE_BUFFER_SIZE;
2185    } while (++c<CC);
2186    lpc = (celt_word16*)(st->_decode_mem+(DECODE_BUFFER_SIZE+st->overlap)*CC);
2187    oldBandE = lpc+CC*LPC_ORDER;
2188    oldLogE = oldBandE + CC*st->mode->nbEBands;
2189    oldLogE2 = oldLogE + CC*st->mode->nbEBands;
2190    backgroundLogE = oldLogE2  + CC*st->mode->nbEBands;
2191
2192    N = M*st->mode->shortMdctSize;
2193
2194    effEnd = st->end;
2195    if (effEnd > st->mode->effEBands)
2196       effEnd = st->mode->effEBands;
2197
2198    ALLOC(freq, CC*N, celt_sig); /**< Interleaved signal MDCTs */
2199    ALLOC(X, CC*N, celt_norm);   /**< Interleaved normalised MDCTs */
2200    ALLOC(bandE, st->mode->nbEBands*CC, celt_ener);
2201    c=0; do
2202       for (i=0;i<M*st->mode->eBands[st->start];i++)
2203          X[c*N+i] = 0;
2204    while (++c<CC);
2205    c=0; do   
2206       for (i=M*st->mode->eBands[effEnd];i<N;i++)
2207          X[c*N+i] = 0;
2208    while (++c<CC);
2209
2210    if (data == NULL || len<=1)
2211    {
2212       celt_decode_lost(st, pcm, N, LM);
2213       RESTORE_STACK;
2214       return CELT_OK;
2215    }
2216    if (len<0) {
2217      RESTORE_STACK;
2218      return CELT_BAD_ARG;
2219    }
2220    
2221    if (dec == NULL)
2222    {
2223       ec_byte_readinit(&buf,(unsigned char*)data,len);
2224       ec_dec_init(&_dec,&buf);
2225       dec = &_dec;
2226    }
2227
2228    if (C>CC)
2229    {
2230       RESTORE_STACK;
2231       return CELT_CORRUPTED_DATA;
2232    } else if (C<CC)
2233    {
2234       for (i=0;i<st->mode->nbEBands;i++)
2235          oldBandE[i]=MAX16(oldBandE[i],oldBandE[st->mode->nbEBands+i]);
2236    }
2237
2238    total_bits = len*8;
2239    tell = ec_dec_tell(dec, 0);
2240
2241    if (tell==1)
2242       silence = ec_dec_bit_logp(dec, 15);
2243    else
2244       silence = 0;
2245    if (silence)
2246    {
2247       /* Pretend we've read all the remaining bits */
2248       tell = len*8;
2249       dec->nbits_total+=tell-ec_dec_tell(dec,0);
2250    }
2251
2252    postfilter_gain = 0;
2253    postfilter_pitch = 0;
2254    postfilter_tapset = 0;
2255    if (st->start==0 && tell+17 <= total_bits)
2256    {
2257       if(ec_dec_bit_logp(dec, 1))
2258       {
2259 #ifdef ENABLE_POSTFILTER
2260          int qg, octave;
2261          octave = ec_dec_uint(dec, 6);
2262          postfilter_pitch = (16<<octave)+ec_dec_bits(dec, 4+octave)-1;
2263          qg = ec_dec_bits(dec, 2);
2264          postfilter_tapset = ec_dec_icdf(dec, tapset_icdf, 2);
2265          postfilter_gain = QCONST16(.125f,15)*(qg+2);
2266 #else /* ENABLE_POSTFILTER */
2267          RESTORE_STACK;
2268          return CELT_CORRUPTED_DATA;
2269 #endif /* ENABLE_POSTFILTER */
2270       }
2271       tell = ec_dec_tell(dec, 0);
2272    }
2273
2274    if (LM > 0 && tell+3 <= total_bits)
2275    {
2276       isTransient = ec_dec_bit_logp(dec, 3);
2277       tell = ec_dec_tell(dec, 0);
2278    }
2279    else
2280       isTransient = 0;
2281
2282    if (isTransient)
2283       shortBlocks = M;
2284    else
2285       shortBlocks = 0;
2286
2287    /* Decode the global flags (first symbols in the stream) */
2288    intra_ener = tell+3<=total_bits ? ec_dec_bit_logp(dec, 3) : 0;
2289    /* Get band energies */
2290    unquant_coarse_energy(st->mode, st->start, st->end, oldBandE,
2291          intra_ener, dec, C, LM);
2292
2293    ALLOC(tf_res, st->mode->nbEBands, int);
2294    tf_decode(st->start, st->end, isTransient, tf_res, LM, dec);
2295
2296    tell = ec_dec_tell(dec, 0);
2297    spread_decision = SPREAD_NORMAL;
2298    if (tell+4 <= total_bits)
2299       spread_decision = ec_dec_icdf(dec, spread_icdf, 5);
2300
2301    ALLOC(pulses, st->mode->nbEBands, int);
2302    ALLOC(cap, st->mode->nbEBands, int);
2303    ALLOC(offsets, st->mode->nbEBands, int);
2304    ALLOC(fine_priority, st->mode->nbEBands, int);
2305
2306    init_caps(st->mode,cap,LM,C);
2307
2308    dynalloc_logp = 6;
2309    total_bits<<=BITRES;
2310    tell = ec_dec_tell(dec, BITRES);
2311    for (i=st->start;i<st->end;i++)
2312    {
2313       int width, quanta;
2314       int dynalloc_loop_logp;
2315       int boost;
2316       width = C*(st->mode->eBands[i+1]-st->mode->eBands[i])<<LM;
2317       /* quanta is 6 bits, but no more than 1 bit/sample
2318          and no less than 1/8 bit/sample */
2319       quanta = IMIN(width<<BITRES, IMAX(6<<BITRES, width));
2320       dynalloc_loop_logp = dynalloc_logp;
2321       boost = 0;
2322       while (tell+(dynalloc_loop_logp<<BITRES) < total_bits && boost < cap[i])
2323       {
2324          int flag;
2325          flag = ec_dec_bit_logp(dec, dynalloc_loop_logp);
2326          tell = ec_dec_tell(dec, BITRES);
2327          if (!flag)
2328             break;
2329          boost += quanta;
2330          total_bits -= quanta;
2331          dynalloc_loop_logp = 1;
2332       }
2333       offsets[i] = boost;
2334       /* Making dynalloc more likely */
2335       if (boost>0)
2336          dynalloc_logp = IMAX(2, dynalloc_logp-1);
2337    }
2338
2339    ALLOC(fine_quant, st->mode->nbEBands, int);
2340    alloc_trim = tell+(6<<BITRES) <= total_bits ?
2341          ec_dec_icdf(dec, trim_icdf, 7) : 5;
2342
2343    bits = (len*8<<BITRES) - ec_dec_tell(dec, BITRES) - 1;
2344    anti_collapse_rsv = isTransient&&LM>=2&&bits>=(LM+2<<BITRES) ? (1<<BITRES) : 0;
2345    bits -= anti_collapse_rsv;
2346    codedBands = compute_allocation(st->mode, st->start, st->end, offsets, cap,
2347          alloc_trim, &intensity, &dual_stereo, bits, &balance, pulses,
2348          fine_quant, fine_priority, C, LM, dec, 0, 0);
2349    
2350    unquant_fine_energy(st->mode, st->start, st->end, oldBandE, fine_quant, dec, C);
2351
2352    /* Decode fixed codebook */
2353    ALLOC(collapse_masks, C*st->mode->nbEBands, unsigned char);
2354    quant_all_bands(0, st->mode, st->start, st->end, X, C==2 ? X+N : NULL, collapse_masks,
2355          NULL, pulses, shortBlocks, spread_decision, dual_stereo, intensity, tf_res, 1,
2356          len*(8<<BITRES)-anti_collapse_rsv, balance, dec, LM, codedBands, &st->rng);
2357
2358    if (anti_collapse_rsv > 0)
2359    {
2360       anti_collapse_on = ec_dec_bits(dec, 1);
2361    }
2362
2363    unquant_energy_finalise(st->mode, st->start, st->end, oldBandE,
2364          fine_quant, fine_priority, len*8-ec_dec_tell(dec, 0), dec, C);
2365
2366    if (anti_collapse_on)
2367       anti_collapse(st->mode, X, collapse_masks, LM, C, CC, N,
2368             st->start, st->end, oldBandE, oldLogE, oldLogE2, pulses, st->rng);
2369
2370    log2Amp(st->mode, st->start, st->end, bandE, oldBandE, C);
2371
2372    if (silence)
2373    {
2374       for (i=0;i<C*st->mode->nbEBands;i++)
2375       {
2376          bandE[i] = 0;
2377          oldBandE[i] = -QCONST16(28.f,DB_SHIFT);
2378       }
2379    }
2380    /* Synthesis */
2381    denormalise_bands(st->mode, X, freq, bandE, effEnd, C, M);
2382
2383    CELT_MOVE(decode_mem[0], decode_mem[0]+N, DECODE_BUFFER_SIZE-N);
2384    if (CC==2)
2385       CELT_MOVE(decode_mem[1], decode_mem[1]+N, DECODE_BUFFER_SIZE-N);
2386
2387    c=0; do
2388       for (i=0;i<M*st->mode->eBands[st->start];i++)
2389          freq[c*N+i] = 0;
2390    while (++c<C);
2391    c=0; do {
2392       int bound = M*st->mode->eBands[effEnd];
2393       if (st->downsample!=1)
2394          bound = IMIN(bound, N/st->downsample);
2395       for (i=M*st->mode->eBands[effEnd];i<N;i++)
2396          freq[c*N+i] = 0;
2397    } while (++c<C);
2398
2399    out_syn[0] = out_mem[0]+MAX_PERIOD-N;
2400    if (CC==2)
2401       out_syn[1] = out_mem[1]+MAX_PERIOD-N;
2402
2403    if (CC==2&&C==1)
2404    {
2405       for (i=0;i<N;i++)
2406          freq[N+i] = freq[i];
2407    }
2408
2409    /* Compute inverse MDCTs */
2410    compute_inv_mdcts(st->mode, shortBlocks, freq, out_syn, overlap_mem, CC, LM);
2411
2412 #ifdef ENABLE_POSTFILTER
2413    c=0; do {
2414       st->postfilter_period=IMAX(st->postfilter_period, COMBFILTER_MINPERIOD);
2415       st->postfilter_period_old=IMAX(st->postfilter_period_old, COMBFILTER_MINPERIOD);
2416       if (LM!=0)
2417       {
2418          comb_filter(out_syn[c], out_syn[c], st->postfilter_period, st->postfilter_period, st->overlap,
2419                st->postfilter_gain, st->postfilter_gain, st->postfilter_tapset, st->postfilter_tapset,
2420                NULL, 0);
2421          comb_filter(out_syn[c]+st->overlap, out_syn[c]+st->overlap, st->postfilter_period, postfilter_pitch, N-st->overlap,
2422                st->postfilter_gain, postfilter_gain, st->postfilter_tapset, postfilter_tapset,
2423                st->mode->window, st->mode->overlap);
2424       } else {
2425          comb_filter(out_syn[c], out_syn[c], st->postfilter_period_old, st->postfilter_period, N-st->overlap,
2426                st->postfilter_gain_old, st->postfilter_gain, st->postfilter_tapset_old, st->postfilter_tapset,
2427                st->mode->window, st->mode->overlap);
2428       }
2429    } while (++c<CC);
2430    st->postfilter_period_old = st->postfilter_period;
2431    st->postfilter_gain_old = st->postfilter_gain;
2432    st->postfilter_tapset_old = st->postfilter_tapset;
2433    st->postfilter_period = postfilter_pitch;
2434    st->postfilter_gain = postfilter_gain;
2435    st->postfilter_tapset = postfilter_tapset;
2436 #endif /* ENABLE_POSTFILTER */
2437
2438    if (CC==2&&C==1) {
2439       for (i=0;i<st->mode->nbEBands;i++)
2440          oldBandE[st->mode->nbEBands+i]=oldBandE[i];
2441    }
2442
2443    /* In case start or end were to change */
2444    c=0; do
2445    {
2446       for (i=0;i<st->start;i++)
2447          oldBandE[c*st->mode->nbEBands+i]=0;
2448       for (i=st->end;i<st->mode->nbEBands;i++)
2449          oldBandE[c*st->mode->nbEBands+i]=0;
2450    } while (++c<CC);
2451    if (!isTransient)
2452    {
2453       for (i=0;i<CC*st->mode->nbEBands;i++)
2454          oldLogE2[i] = oldLogE[i];
2455       for (i=0;i<CC*st->mode->nbEBands;i++)
2456          oldLogE[i] = oldBandE[i];
2457       for (i=0;i<CC*st->mode->nbEBands;i++)
2458          backgroundLogE[i] = MIN16(backgroundLogE[i] + M*QCONST16(0.001f,DB_SHIFT), oldBandE[i]);
2459    } else {
2460       for (i=0;i<CC*st->mode->nbEBands;i++)
2461          oldLogE[i] = MIN16(oldLogE[i], oldBandE[i]);
2462    }
2463    st->rng = dec->rng;
2464
2465    deemphasis(out_syn, pcm, N, CC, st->downsample, st->mode->preemph, st->preemph_memD);
2466    st->loss_count = 0;
2467    RESTORE_STACK;
2468    if (ec_dec_tell(dec,0) > 8*len || ec_dec_get_error(dec))
2469       return CELT_CORRUPTED_DATA;
2470    else
2471       return CELT_OK;
2472 }
2473
2474 #ifdef FIXED_POINT
2475 #ifndef DISABLE_FLOAT_API
2476 int celt_decode_with_ec_float(CELTDecoder * restrict st, const unsigned char *data, int len, float * restrict pcm, int frame_size, ec_dec *dec)
2477 {
2478    int j, ret, C, N;
2479    VARDECL(celt_int16, out);
2480    ALLOC_STACK;
2481    SAVE_STACK;
2482
2483    if (pcm==NULL)
2484       return CELT_BAD_ARG;
2485
2486    C = CHANNELS(st->channels);
2487    N = frame_size;
2488    
2489    ALLOC(out, C*N, celt_int16);
2490    ret=celt_decode_with_ec(st, data, len, out, frame_size, dec);
2491    if (ret==0)
2492       for (j=0;j<C*N;j++)
2493          pcm[j]=out[j]*(1.f/32768.f);
2494      
2495    RESTORE_STACK;
2496    return ret;
2497 }
2498 #endif /*DISABLE_FLOAT_API*/
2499 #else
2500 int celt_decode_with_ec(CELTDecoder * restrict st, const unsigned char *data, int len, celt_int16 * restrict pcm, int frame_size, ec_dec *dec)
2501 {
2502    int j, ret, C, N;
2503    VARDECL(celt_sig, out);
2504    ALLOC_STACK;
2505    SAVE_STACK;
2506
2507    if (pcm==NULL)
2508       return CELT_BAD_ARG;
2509
2510    C = CHANNELS(st->channels);
2511    N = frame_size;
2512    ALLOC(out, C*N, celt_sig);
2513
2514    ret=celt_decode_with_ec_float(st, data, len, out, frame_size, dec);
2515
2516    if (ret==0)
2517       for (j=0;j<C*N;j++)
2518          pcm[j] = FLOAT2INT16 (out[j]);
2519    
2520    RESTORE_STACK;
2521    return ret;
2522 }
2523 #endif
2524
2525 int celt_decode(CELTDecoder * restrict st, const unsigned char *data, int len, celt_int16 * restrict pcm, int frame_size)
2526 {
2527    return celt_decode_with_ec(st, data, len, pcm, frame_size, NULL);
2528 }
2529
2530 #ifndef DISABLE_FLOAT_API
2531 int celt_decode_float(CELTDecoder * restrict st, const unsigned char *data, int len, float * restrict pcm, int frame_size)
2532 {
2533    return celt_decode_with_ec_float(st, data, len, pcm, frame_size, NULL);
2534 }
2535 #endif /* DISABLE_FLOAT_API */
2536
2537 int celt_decoder_ctl(CELTDecoder * restrict st, int request, ...)
2538 {
2539    va_list ap;
2540
2541    va_start(ap, request);
2542    switch (request)
2543    {
2544       case CELT_GET_MODE_REQUEST:
2545       {
2546          const CELTMode ** value = va_arg(ap, const CELTMode**);
2547          if (value==0)
2548             goto bad_arg;
2549          *value=st->mode;
2550       }
2551       break;
2552       case CELT_SET_START_BAND_REQUEST:
2553       {
2554          celt_int32 value = va_arg(ap, celt_int32);
2555          if (value<0 || value>=st->mode->nbEBands)
2556             goto bad_arg;
2557          st->start = value;
2558       }
2559       break;
2560       case CELT_SET_END_BAND_REQUEST:
2561       {
2562          celt_int32 value = va_arg(ap, celt_int32);
2563          if (value<0 || value>=st->mode->nbEBands)
2564             goto bad_arg;
2565          st->end = value;
2566       }
2567       break;
2568       case CELT_SET_CHANNELS_REQUEST:
2569       {
2570          celt_int32 value = va_arg(ap, celt_int32);
2571          if (value<1 || value>2)
2572             goto bad_arg;
2573          st->stream_channels = value;
2574       }
2575       break;
2576       case CELT_RESET_STATE:
2577       {
2578          CELT_MEMSET((char*)&st->DECODER_RESET_START, 0,
2579                celt_decoder_get_size_custom(st->mode, st->channels)-
2580                ((char*)&st->DECODER_RESET_START - (char*)st));
2581       }
2582       break;
2583       default:
2584          goto bad_request;
2585    }
2586    va_end(ap);
2587    return CELT_OK;
2588 bad_arg:
2589    va_end(ap);
2590    return CELT_BAD_ARG;
2591 bad_request:
2592       va_end(ap);
2593   return CELT_UNIMPLEMENTED;
2594 }
2595
2596 const char *celt_strerror(int error)
2597 {
2598    static const char *error_strings[8] = {
2599       "success",
2600       "invalid argument",
2601       "invalid mode",
2602       "internal error",
2603       "corrupted stream",
2604       "request not implemented",
2605       "invalid state",
2606       "memory allocation failed"
2607    };
2608    if (error > 0 || error < -7)
2609       return "unknown error";
2610    else 
2611       return error_strings[-error];
2612 }
2613