Adds a generic CELT_SET_BITRATE() ctl() API for CBR and VBR
[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, 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    SAVE_STACK;
1589
1590    if (pcm==NULL)
1591       return CELT_BAD_ARG;
1592
1593    C = CHANNELS(st->channels);
1594    N = frame_size;
1595    ALLOC(in, C*N, celt_int16);
1596
1597    for (j=0;j<C*N;j++)
1598      in[j] = FLOAT2INT16(pcm[j]);
1599
1600    ret=celt_encode_with_ec(st,in,frame_size,compressed,nbCompressedBytes, enc);
1601 #ifdef RESYNTH
1602    for (j=0;j<C*N;j++)
1603       ((float*)pcm)[j]=in[j]*(1.f/32768.f);
1604 #endif
1605    RESTORE_STACK;
1606    return ret;
1607
1608 }
1609 #endif /*DISABLE_FLOAT_API*/
1610 #else
1611 int celt_encode_with_ec(CELTEncoder * restrict st, const celt_int16 * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes, ec_enc *enc)
1612 {
1613    int j, ret, C, N;
1614    VARDECL(celt_sig, in);
1615    SAVE_STACK;
1616
1617    if (pcm==NULL)
1618       return CELT_BAD_ARG;
1619
1620    C=CHANNELS(st->channels);
1621    N=frame_size;
1622    ALLOC(in, C*N, celt_sig);
1623    for (j=0;j<C*N;j++) {
1624      in[j] = SCALEOUT(pcm[j]);
1625    }
1626
1627    ret = celt_encode_with_ec_float(st,in,frame_size,compressed,nbCompressedBytes, enc);
1628 #ifdef RESYNTH
1629    for (j=0;j<C*N;j++)
1630       ((celt_int16*)pcm)[j] = FLOAT2INT16(in[j]);
1631 #endif
1632    RESTORE_STACK;
1633    return ret;
1634 }
1635 #endif
1636
1637 int celt_encode(CELTEncoder * restrict st, const celt_int16 * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes)
1638 {
1639    return celt_encode_with_ec(st, pcm, frame_size, compressed, nbCompressedBytes, NULL);
1640 }
1641
1642 #ifndef DISABLE_FLOAT_API
1643 int celt_encode_float(CELTEncoder * restrict st, const float * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes)
1644 {
1645    return celt_encode_with_ec_float(st, pcm, frame_size, compressed, nbCompressedBytes, NULL);
1646 }
1647 #endif /* DISABLE_FLOAT_API */
1648
1649 int celt_encoder_ctl(CELTEncoder * restrict st, int request, ...)
1650 {
1651    va_list ap;
1652    
1653    va_start(ap, request);
1654    switch (request)
1655    {
1656       case CELT_GET_MODE_REQUEST:
1657       {
1658          const CELTMode ** value = va_arg(ap, const CELTMode**);
1659          if (value==0)
1660             goto bad_arg;
1661          *value=st->mode;
1662       }
1663       break;
1664       case CELT_SET_COMPLEXITY_REQUEST:
1665       {
1666          int value = va_arg(ap, celt_int32);
1667          if (value<0 || value>10)
1668             goto bad_arg;
1669          st->complexity = value;
1670       }
1671       break;
1672       case CELT_SET_START_BAND_REQUEST:
1673       {
1674          celt_int32 value = va_arg(ap, celt_int32);
1675          if (value<0 || value>=st->mode->nbEBands)
1676             goto bad_arg;
1677          st->start = value;
1678       }
1679       break;
1680       case CELT_SET_END_BAND_REQUEST:
1681       {
1682          celt_int32 value = va_arg(ap, celt_int32);
1683          if (value<1 || value>st->mode->nbEBands)
1684             goto bad_arg;
1685          st->end = value;
1686       }
1687       break;
1688       case CELT_SET_PREDICTION_REQUEST:
1689       {
1690          int value = va_arg(ap, celt_int32);
1691          if (value<0 || value>2)
1692             goto bad_arg;
1693          if (value==0)
1694          {
1695             st->force_intra   = 1;
1696          } else if (value==1) {
1697             st->force_intra   = 0;
1698          } else {
1699             st->force_intra   = 0;
1700          }   
1701       }
1702       break;
1703       case CELT_SET_VBR_CONSTRAINT_REQUEST:
1704       {
1705          celt_int32 value = va_arg(ap, celt_int32);
1706          st->constrained_vbr = value;
1707       }
1708       break;
1709       case CELT_SET_VBR_REQUEST:
1710       {
1711          celt_int32 value = va_arg(ap, celt_int32);
1712          st->vbr = value;
1713       }
1714       break;
1715       case CELT_SET_BITRATE_REQUEST:
1716       {
1717          celt_int32 value = va_arg(ap, celt_int32);
1718          if (value<=500)
1719             goto bad_arg;
1720          st->bitrate = value;
1721       }
1722       break;
1723       case CELT_SET_CHANNELS_REQUEST:
1724       {
1725          celt_int32 value = va_arg(ap, celt_int32);
1726          if (value<1 || value>2)
1727             goto bad_arg;
1728          st->stream_channels = value;
1729       }
1730       break;
1731       case CELT_RESET_STATE:
1732       {
1733          CELT_MEMSET((char*)&st->ENCODER_RESET_START, 0,
1734                celt_encoder_get_size_custom(st->mode, st->channels)-
1735                ((char*)&st->ENCODER_RESET_START - (char*)st));
1736          st->vbr_offset = 0;
1737          st->delayedIntra = 1;
1738          st->spread_decision = SPREAD_NORMAL;
1739          st->tonal_average = QCONST16(1.f,8);
1740       }
1741       break;
1742       default:
1743          goto bad_request;
1744    }
1745    va_end(ap);
1746    return CELT_OK;
1747 bad_arg:
1748    va_end(ap);
1749    return CELT_BAD_ARG;
1750 bad_request:
1751    va_end(ap);
1752    return CELT_UNIMPLEMENTED;
1753 }
1754
1755 /**********************************************************************/
1756 /*                                                                    */
1757 /*                             DECODER                                */
1758 /*                                                                    */
1759 /**********************************************************************/
1760 #define DECODE_BUFFER_SIZE 2048
1761
1762 /** Decoder state 
1763  @brief Decoder state
1764  */
1765 struct CELTDecoder {
1766    const CELTMode *mode;
1767    int overlap;
1768    int channels;
1769    int stream_channels;
1770
1771    int downsample;
1772    int start, end;
1773
1774    /* Everything beyond this point gets cleared on a reset */
1775 #define DECODER_RESET_START rng
1776
1777    ec_uint32 rng;
1778    int last_pitch_index;
1779    int loss_count;
1780    int postfilter_period;
1781    int postfilter_period_old;
1782    celt_word16 postfilter_gain;
1783    celt_word16 postfilter_gain_old;
1784    int postfilter_tapset;
1785    int postfilter_tapset_old;
1786
1787    celt_sig preemph_memD[2];
1788    
1789    celt_sig _decode_mem[1]; /* Size = channels*(DECODE_BUFFER_SIZE+mode->overlap) */
1790    /* celt_word16 lpc[],  Size = channels*LPC_ORDER */
1791    /* celt_word16 oldEBands[], Size = channels*mode->nbEBands */
1792    /* celt_word16 oldLogE[], Size = channels*mode->nbEBands */
1793    /* celt_word16 oldLogE2[], Size = channels*mode->nbEBands */
1794    /* celt_word16 backgroundLogE[], Size = channels*mode->nbEBands */
1795 };
1796
1797 int celt_decoder_get_size(int channels)
1798 {
1799    const CELTMode *mode = celt_mode_create(48000, 960, NULL);
1800    return celt_decoder_get_size_custom(mode, channels);
1801 }
1802
1803 int celt_decoder_get_size_custom(const CELTMode *mode, int channels)
1804 {
1805    int size = sizeof(struct CELTDecoder)
1806             + (channels*(DECODE_BUFFER_SIZE+mode->overlap)-1)*sizeof(celt_sig)
1807             + channels*LPC_ORDER*sizeof(celt_word16)
1808             + 4*channels*mode->nbEBands*sizeof(celt_word16);
1809    return size;
1810 }
1811
1812 CELTDecoder *celt_decoder_create(int sampling_rate, int channels, int *error)
1813 {
1814    CELTDecoder *st;
1815    st = (CELTDecoder *)celt_alloc(celt_decoder_get_size(channels));
1816    if (st!=NULL && celt_decoder_init(st, sampling_rate, channels, error)==NULL)
1817    {
1818       celt_decoder_destroy(st);
1819       st = NULL;
1820    }
1821    return st;
1822 }
1823
1824 CELTDecoder *celt_decoder_create_custom(const CELTMode *mode, int channels, int *error)
1825 {
1826    CELTDecoder *st = (CELTDecoder *)celt_alloc(celt_decoder_get_size_custom(mode, channels));
1827    if (st!=NULL && celt_decoder_init_custom(st, mode, channels, error)==NULL)
1828    {
1829       celt_decoder_destroy(st);
1830       st = NULL;
1831    }
1832    return st;
1833 }
1834
1835 CELTDecoder *celt_decoder_init(CELTDecoder *st, int sampling_rate, int channels, int *error)
1836 {
1837    celt_decoder_init_custom(st, celt_mode_create(48000, 960, NULL), channels, error);
1838    st->downsample = resampling_factor(sampling_rate);
1839    if (st->downsample==0)
1840    {
1841       if (error)
1842          *error = CELT_BAD_ARG;
1843       return NULL;
1844    }
1845    return st;
1846 }
1847
1848 CELTDecoder *celt_decoder_init_custom(CELTDecoder *st, const CELTMode *mode, int channels, int *error)
1849 {
1850    if (channels < 0 || channels > 2)
1851    {
1852       if (error)
1853          *error = CELT_BAD_ARG;
1854       return NULL;
1855    }
1856
1857    if (st==NULL)
1858    {
1859       if (error)
1860          *error = CELT_ALLOC_FAIL;
1861       return NULL;
1862    }
1863
1864    CELT_MEMSET((char*)st, 0, celt_decoder_get_size_custom(mode, channels));
1865
1866    st->mode = mode;
1867    st->overlap = mode->overlap;
1868    st->stream_channels = st->channels = channels;
1869
1870    st->downsample = 1;
1871    st->start = 0;
1872    st->end = st->mode->effEBands;
1873
1874    st->loss_count = 0;
1875
1876    if (error)
1877       *error = CELT_OK;
1878    return st;
1879 }
1880
1881 void celt_decoder_destroy(CELTDecoder *st)
1882 {
1883    celt_free(st);
1884 }
1885
1886 static void celt_decode_lost(CELTDecoder * restrict st, celt_word16 * restrict pcm, int N, int LM)
1887 {
1888    int c;
1889    int pitch_index;
1890    int overlap = st->mode->overlap;
1891    celt_word16 fade = Q15ONE;
1892    int i, len;
1893    const int C = CHANNELS(st->channels);
1894    int offset;
1895    celt_sig *out_mem[2];
1896    celt_sig *decode_mem[2];
1897    celt_sig *overlap_mem[2];
1898    celt_word16 *lpc;
1899    celt_word32 *out_syn[2];
1900    celt_word16 *oldBandE, *oldLogE2, *backgroundLogE;
1901    SAVE_STACK;
1902    
1903    c=0; do {
1904       decode_mem[c] = st->_decode_mem + c*(DECODE_BUFFER_SIZE+st->overlap);
1905       out_mem[c] = decode_mem[c]+DECODE_BUFFER_SIZE-MAX_PERIOD;
1906       overlap_mem[c] = decode_mem[c]+DECODE_BUFFER_SIZE;
1907    } while (++c<C);
1908    lpc = (celt_word16*)(st->_decode_mem+(DECODE_BUFFER_SIZE+st->overlap)*C);
1909    oldBandE = lpc+C*LPC_ORDER;
1910    oldLogE2 = oldBandE + C*st->mode->nbEBands;
1911    backgroundLogE = oldLogE2  + C*st->mode->nbEBands;
1912
1913    out_syn[0] = out_mem[0]+MAX_PERIOD-N;
1914    if (C==2)
1915       out_syn[1] = out_mem[1]+MAX_PERIOD-N;
1916
1917    len = N+st->mode->overlap;
1918    
1919    if (st->loss_count >= 5)
1920    {
1921       VARDECL(celt_sig, freq);
1922       VARDECL(celt_norm, X);
1923       VARDECL(celt_ener, bandE);
1924       celt_uint32 seed;
1925
1926       ALLOC(freq, C*N, celt_sig); /**< Interleaved signal MDCTs */
1927       ALLOC(X, C*N, celt_norm);   /**< Interleaved normalised MDCTs */
1928       ALLOC(bandE, st->mode->nbEBands*C, celt_ener);
1929
1930       log2Amp(st->mode, st->start, st->end, bandE, backgroundLogE, C);
1931
1932       seed = st->rng;
1933       for (i=0;i<C*N;i++)
1934       {
1935             seed = lcg_rand(seed);
1936             X[i] = (celt_int32)(seed)>>20;
1937       }
1938       st->rng = seed;
1939       for (c=0;c<C;c++)
1940          for (i=0;i<st->mode->nbEBands;i++)
1941             renormalise_vector(X+N*c+(st->mode->eBands[i]<<LM), (st->mode->eBands[i+1]-st->mode->eBands[i])<<LM, Q15ONE);
1942
1943       denormalise_bands(st->mode, X, freq, bandE, st->mode->nbEBands, C, 1<<LM);
1944
1945       compute_inv_mdcts(st->mode, 0, freq, out_syn, overlap_mem, C, LM);
1946    } else if (st->loss_count == 0)
1947    {
1948       celt_word16 pitch_buf[MAX_PERIOD>>1];
1949       int len2 = len;
1950       /* FIXME: This is a kludge */
1951       if (len2>MAX_PERIOD>>1)
1952          len2 = MAX_PERIOD>>1;
1953       pitch_downsample(out_mem, pitch_buf, MAX_PERIOD, C);
1954       pitch_search(pitch_buf+((MAX_PERIOD-len2)>>1), pitch_buf, len2,
1955                    MAX_PERIOD-len2-100, &pitch_index);
1956       pitch_index = MAX_PERIOD-len2-pitch_index;
1957       st->last_pitch_index = pitch_index;
1958    } else {
1959       pitch_index = st->last_pitch_index;
1960       fade = QCONST16(.8f,15);
1961    }
1962
1963    c=0; do {
1964       /* FIXME: This is more memory than necessary */
1965       celt_word32 e[2*MAX_PERIOD];
1966       celt_word16 exc[2*MAX_PERIOD];
1967       celt_word32 ac[LPC_ORDER+1];
1968       celt_word16 decay = 1;
1969       celt_word32 S1=0;
1970       celt_word16 mem[LPC_ORDER]={0};
1971
1972       offset = MAX_PERIOD-pitch_index;
1973       for (i=0;i<MAX_PERIOD;i++)
1974          exc[i] = ROUND16(out_mem[c][i], SIG_SHIFT);
1975
1976       if (st->loss_count == 0)
1977       {
1978          _celt_autocorr(exc, ac, st->mode->window, st->mode->overlap,
1979                         LPC_ORDER, MAX_PERIOD);
1980
1981          /* Noise floor -40 dB */
1982 #ifdef FIXED_POINT
1983          ac[0] += SHR32(ac[0],13);
1984 #else
1985          ac[0] *= 1.0001f;
1986 #endif
1987          /* Lag windowing */
1988          for (i=1;i<=LPC_ORDER;i++)
1989          {
1990             /*ac[i] *= exp(-.5*(2*M_PI*.002*i)*(2*M_PI*.002*i));*/
1991 #ifdef FIXED_POINT
1992             ac[i] -= MULT16_32_Q15(2*i*i, ac[i]);
1993 #else
1994             ac[i] -= ac[i]*(.008f*i)*(.008f*i);
1995 #endif
1996          }
1997
1998          _celt_lpc(lpc+c*LPC_ORDER, ac, LPC_ORDER);
1999       }
2000       for (i=0;i<LPC_ORDER;i++)
2001          mem[i] = ROUND16(out_mem[c][MAX_PERIOD-1-i], SIG_SHIFT);
2002       fir(exc, lpc+c*LPC_ORDER, exc, MAX_PERIOD, LPC_ORDER, mem);
2003       /*for (i=0;i<MAX_PERIOD;i++)printf("%d ", exc[i]); printf("\n");*/
2004       /* Check if the waveform is decaying (and if so how fast) */
2005       {
2006          celt_word32 E1=1, E2=1;
2007          int period;
2008          if (pitch_index <= MAX_PERIOD/2)
2009             period = pitch_index;
2010          else
2011             period = MAX_PERIOD/2;
2012          for (i=0;i<period;i++)
2013          {
2014             E1 += SHR32(MULT16_16(exc[MAX_PERIOD-period+i],exc[MAX_PERIOD-period+i]),8);
2015             E2 += SHR32(MULT16_16(exc[MAX_PERIOD-2*period+i],exc[MAX_PERIOD-2*period+i]),8);
2016          }
2017          if (E1 > E2)
2018             E1 = E2;
2019          decay = celt_sqrt(frac_div32(SHR(E1,1),E2));
2020       }
2021
2022       /* Copy excitation, taking decay into account */
2023       for (i=0;i<len+st->mode->overlap;i++)
2024       {
2025          celt_word16 tmp;
2026          if (offset+i >= MAX_PERIOD)
2027          {
2028             offset -= pitch_index;
2029             decay = MULT16_16_Q15(decay, decay);
2030          }
2031          e[i] = SHL32(EXTEND32(MULT16_16_Q15(decay, exc[offset+i])), SIG_SHIFT);
2032          tmp = ROUND16(out_mem[c][offset+i],SIG_SHIFT);
2033          S1 += SHR32(MULT16_16(tmp,tmp),8);
2034       }
2035       for (i=0;i<LPC_ORDER;i++)
2036          mem[i] = ROUND16(out_mem[c][MAX_PERIOD-1-i], SIG_SHIFT);
2037       for (i=0;i<len+st->mode->overlap;i++)
2038          e[i] = MULT16_32_Q15(fade, e[i]);
2039       iir(e, lpc+c*LPC_ORDER, e, len+st->mode->overlap, LPC_ORDER, mem);
2040
2041       {
2042          celt_word32 S2=0;
2043          for (i=0;i<len+overlap;i++)
2044          {
2045             celt_word16 tmp = ROUND16(e[i],SIG_SHIFT);
2046             S2 += SHR32(MULT16_16(tmp,tmp),8);
2047          }
2048          /* This checks for an "explosion" in the synthesis */
2049 #ifdef FIXED_POINT
2050          if (!(S1 > SHR32(S2,2)))
2051 #else
2052          /* Float test is written this way to catch NaNs at the same time */
2053          if (!(S1 > 0.2f*S2))
2054 #endif
2055          {
2056             for (i=0;i<len+overlap;i++)
2057                e[i] = 0;
2058          } else if (S1 < S2)
2059          {
2060             celt_word16 ratio = celt_sqrt(frac_div32(SHR32(S1,1)+1,S2+1));
2061             for (i=0;i<len+overlap;i++)
2062                e[i] = MULT16_32_Q15(ratio, e[i]);
2063          }
2064       }
2065
2066 #ifdef ENABLE_POSTFILTER
2067       /* Apply post-filter to the MDCT overlap of the previous frame */
2068       comb_filter(out_mem[c]+MAX_PERIOD, out_mem[c]+MAX_PERIOD, st->postfilter_period, st->postfilter_period, st->overlap,
2069                   st->postfilter_gain, st->postfilter_gain, st->postfilter_tapset, st->postfilter_tapset,
2070                   NULL, 0);
2071 #endif /* ENABLE_POSTFILTER */
2072
2073       for (i=0;i<MAX_PERIOD+st->mode->overlap-N;i++)
2074          out_mem[c][i] = out_mem[c][N+i];
2075
2076       /* Apply TDAC to the concealed audio so that it blends with the
2077          previous and next frames */
2078       for (i=0;i<overlap/2;i++)
2079       {
2080          celt_word32 tmp;
2081          tmp = MULT16_32_Q15(st->mode->window[i],           e[N+overlap-1-i]) +
2082                MULT16_32_Q15(st->mode->window[overlap-i-1], e[N+i          ]);
2083          out_mem[c][MAX_PERIOD+i] = MULT16_32_Q15(st->mode->window[overlap-i-1], tmp);
2084          out_mem[c][MAX_PERIOD+overlap-i-1] = MULT16_32_Q15(st->mode->window[i], tmp);
2085       }
2086       for (i=0;i<N;i++)
2087          out_mem[c][MAX_PERIOD-N+i] = e[i];
2088
2089 #ifdef ENABLE_POSTFILTER
2090       /* Apply pre-filter to the MDCT overlap for the next frame (post-filter will be applied then) */
2091       comb_filter(e, out_mem[c]+MAX_PERIOD, st->postfilter_period, st->postfilter_period, st->overlap,
2092                   -st->postfilter_gain, -st->postfilter_gain, st->postfilter_tapset, st->postfilter_tapset,
2093                   NULL, 0);
2094 #endif /* ENABLE_POSTFILTER */
2095       for (i=0;i<overlap;i++)
2096          out_mem[c][MAX_PERIOD+i] = e[i];
2097    } while (++c<C);
2098
2099    deemphasis(out_syn, pcm, N, C, st->downsample, st->mode->preemph, st->preemph_memD);
2100    
2101    st->loss_count++;
2102
2103    RESTORE_STACK;
2104 }
2105
2106 #ifdef FIXED_POINT
2107 int celt_decode_with_ec(CELTDecoder * restrict st, const unsigned char *data, int len, celt_int16 * restrict pcm, int frame_size, ec_dec *dec)
2108 {
2109 #else
2110 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)
2111 {
2112 #endif
2113    int c, i, N;
2114    int spread_decision;
2115    int bits;
2116    ec_dec _dec;
2117    ec_byte_buffer buf;
2118    VARDECL(celt_sig, freq);
2119    VARDECL(celt_norm, X);
2120    VARDECL(celt_ener, bandE);
2121    VARDECL(int, fine_quant);
2122    VARDECL(int, pulses);
2123    VARDECL(int, cap);
2124    VARDECL(int, offsets);
2125    VARDECL(int, fine_priority);
2126    VARDECL(int, tf_res);
2127    VARDECL(unsigned char, collapse_masks);
2128    celt_sig *out_mem[2];
2129    celt_sig *decode_mem[2];
2130    celt_sig *overlap_mem[2];
2131    celt_sig *out_syn[2];
2132    celt_word16 *lpc;
2133    celt_word16 *oldBandE, *oldLogE, *oldLogE2, *backgroundLogE;
2134
2135    int shortBlocks;
2136    int isTransient;
2137    int intra_ener;
2138    const int CC = CHANNELS(st->channels);
2139    int LM, M;
2140    int effEnd;
2141    int codedBands;
2142    int alloc_trim;
2143    int postfilter_pitch;
2144    celt_word16 postfilter_gain;
2145    int intensity=0;
2146    int dual_stereo=0;
2147    celt_int32 total_bits;
2148    celt_int32 balance;
2149    celt_int32 tell;
2150    int dynalloc_logp;
2151    int postfilter_tapset;
2152    int anti_collapse_rsv;
2153    int anti_collapse_on=0;
2154    int silence;
2155    const int C = CHANNELS(st->stream_channels);
2156
2157    SAVE_STACK;
2158
2159    if (pcm==NULL)
2160       return CELT_BAD_ARG;
2161
2162    frame_size *= st->downsample;
2163    for (LM=0;LM<4;LM++)
2164       if (st->mode->shortMdctSize<<LM==frame_size)
2165          break;
2166    if (LM>=MAX_CONFIG_SIZES)
2167       return CELT_BAD_ARG;
2168    M=1<<LM;
2169
2170    c=0; do {
2171       decode_mem[c] = st->_decode_mem + c*(DECODE_BUFFER_SIZE+st->overlap);
2172       out_mem[c] = decode_mem[c]+DECODE_BUFFER_SIZE-MAX_PERIOD;
2173       overlap_mem[c] = decode_mem[c]+DECODE_BUFFER_SIZE;
2174    } while (++c<CC);
2175    lpc = (celt_word16*)(st->_decode_mem+(DECODE_BUFFER_SIZE+st->overlap)*CC);
2176    oldBandE = lpc+CC*LPC_ORDER;
2177    oldLogE = oldBandE + CC*st->mode->nbEBands;
2178    oldLogE2 = oldLogE + CC*st->mode->nbEBands;
2179    backgroundLogE = oldLogE2  + CC*st->mode->nbEBands;
2180
2181    N = M*st->mode->shortMdctSize;
2182
2183    effEnd = st->end;
2184    if (effEnd > st->mode->effEBands)
2185       effEnd = st->mode->effEBands;
2186
2187    ALLOC(freq, CC*N, celt_sig); /**< Interleaved signal MDCTs */
2188    ALLOC(X, CC*N, celt_norm);   /**< Interleaved normalised MDCTs */
2189    ALLOC(bandE, st->mode->nbEBands*CC, celt_ener);
2190    c=0; do
2191       for (i=0;i<M*st->mode->eBands[st->start];i++)
2192          X[c*N+i] = 0;
2193    while (++c<CC);
2194    c=0; do   
2195       for (i=M*st->mode->eBands[effEnd];i<N;i++)
2196          X[c*N+i] = 0;
2197    while (++c<CC);
2198
2199    if (data == NULL || len<=1)
2200    {
2201       celt_decode_lost(st, pcm, N, LM);
2202       RESTORE_STACK;
2203       return CELT_OK;
2204    }
2205    if (len<0) {
2206      RESTORE_STACK;
2207      return CELT_BAD_ARG;
2208    }
2209    
2210    if (dec == NULL)
2211    {
2212       ec_byte_readinit(&buf,(unsigned char*)data,len);
2213       ec_dec_init(&_dec,&buf);
2214       dec = &_dec;
2215    }
2216
2217    if (C>CC)
2218    {
2219       RESTORE_STACK;
2220       return CELT_CORRUPTED_DATA;
2221    } else if (C<CC)
2222    {
2223       for (i=0;i<st->mode->nbEBands;i++)
2224          oldBandE[i]=MAX16(oldBandE[i],oldBandE[st->mode->nbEBands+i]);
2225    }
2226
2227    total_bits = len*8;
2228    tell = ec_dec_tell(dec, 0);
2229
2230    if (tell==1)
2231       silence = ec_dec_bit_logp(dec, 15);
2232    else
2233       silence = 0;
2234    if (silence)
2235    {
2236       /* Pretend we've read all the remaining bits */
2237       tell = len*8;
2238       dec->nbits_total+=tell-ec_dec_tell(dec,0);
2239    }
2240
2241    postfilter_gain = 0;
2242    postfilter_pitch = 0;
2243    postfilter_tapset = 0;
2244    if (st->start==0 && tell+17 <= total_bits)
2245    {
2246       if(ec_dec_bit_logp(dec, 1))
2247       {
2248 #ifdef ENABLE_POSTFILTER
2249          int qg, octave;
2250          octave = ec_dec_uint(dec, 6);
2251          postfilter_pitch = (16<<octave)+ec_dec_bits(dec, 4+octave)-1;
2252          qg = ec_dec_bits(dec, 2);
2253          postfilter_tapset = ec_dec_icdf(dec, tapset_icdf, 2);
2254          postfilter_gain = QCONST16(.125f,15)*(qg+2);
2255 #else /* ENABLE_POSTFILTER */
2256          RESTORE_STACK;
2257          return CELT_CORRUPTED_DATA;
2258 #endif /* ENABLE_POSTFILTER */
2259       }
2260       tell = ec_dec_tell(dec, 0);
2261    }
2262
2263    if (LM > 0 && tell+3 <= total_bits)
2264    {
2265       isTransient = ec_dec_bit_logp(dec, 3);
2266       tell = ec_dec_tell(dec, 0);
2267    }
2268    else
2269       isTransient = 0;
2270
2271    if (isTransient)
2272       shortBlocks = M;
2273    else
2274       shortBlocks = 0;
2275
2276    /* Decode the global flags (first symbols in the stream) */
2277    intra_ener = tell+3<=total_bits ? ec_dec_bit_logp(dec, 3) : 0;
2278    /* Get band energies */
2279    unquant_coarse_energy(st->mode, st->start, st->end, oldBandE,
2280          intra_ener, dec, C, LM);
2281
2282    ALLOC(tf_res, st->mode->nbEBands, int);
2283    tf_decode(st->start, st->end, isTransient, tf_res, LM, dec);
2284
2285    tell = ec_dec_tell(dec, 0);
2286    spread_decision = SPREAD_NORMAL;
2287    if (tell+4 <= total_bits)
2288       spread_decision = ec_dec_icdf(dec, spread_icdf, 5);
2289
2290    ALLOC(pulses, st->mode->nbEBands, int);
2291    ALLOC(cap, st->mode->nbEBands, int);
2292    ALLOC(offsets, st->mode->nbEBands, int);
2293    ALLOC(fine_priority, st->mode->nbEBands, int);
2294
2295    for (i=0;i<st->mode->nbEBands;i++)
2296       cap[i] = st->mode->cache.caps[st->mode->nbEBands*(2*LM+C-1)+i]
2297             << C+LM+BITRES-2;
2298
2299    dynalloc_logp = 6;
2300    total_bits<<=BITRES;
2301    tell = ec_dec_tell(dec, BITRES);
2302    for (i=st->start;i<st->end;i++)
2303    {
2304       int width, quanta;
2305       int dynalloc_loop_logp;
2306       int boost;
2307       width = C*(st->mode->eBands[i+1]-st->mode->eBands[i])<<LM;
2308       /* quanta is 6 bits, but no more than 1 bit/sample
2309          and no less than 1/8 bit/sample */
2310       quanta = IMIN(width<<BITRES, IMAX(6<<BITRES, width));
2311       dynalloc_loop_logp = dynalloc_logp;
2312       boost = 0;
2313       while (tell+(dynalloc_loop_logp<<BITRES) < total_bits && boost < cap[i])
2314       {
2315          int flag;
2316          flag = ec_dec_bit_logp(dec, dynalloc_loop_logp);
2317          tell = ec_dec_tell(dec, BITRES);
2318          if (!flag)
2319             break;
2320          boost += quanta;
2321          total_bits -= quanta;
2322          dynalloc_loop_logp = 1;
2323       }
2324       offsets[i] = boost;
2325       /* Making dynalloc more likely */
2326       if (boost>0)
2327          dynalloc_logp = IMAX(2, dynalloc_logp-1);
2328    }
2329
2330    ALLOC(fine_quant, st->mode->nbEBands, int);
2331    alloc_trim = tell+(6<<BITRES) <= total_bits ?
2332          ec_dec_icdf(dec, trim_icdf, 7) : 5;
2333
2334    bits = (len*8<<BITRES) - ec_dec_tell(dec, BITRES) - 1;
2335    anti_collapse_rsv = isTransient&&LM>=2&&bits>=(LM+2<<BITRES) ? (1<<BITRES) : 0;
2336    bits -= anti_collapse_rsv;
2337    codedBands = compute_allocation(st->mode, st->start, st->end, offsets, cap,
2338          alloc_trim, &intensity, &dual_stereo, bits, &balance, pulses,
2339          fine_quant, fine_priority, C, LM, dec, 0, 0);
2340    
2341    unquant_fine_energy(st->mode, st->start, st->end, oldBandE, fine_quant, dec, C);
2342
2343    /* Decode fixed codebook */
2344    ALLOC(collapse_masks, st->mode->nbEBands, unsigned char);
2345    quant_all_bands(0, st->mode, st->start, st->end, X, C==2 ? X+N : NULL, collapse_masks,
2346          NULL, pulses, shortBlocks, spread_decision, dual_stereo, intensity, tf_res, 1,
2347          len*(8<<BITRES)-anti_collapse_rsv, balance, dec, LM, codedBands, &st->rng);
2348
2349    if (anti_collapse_rsv > 0)
2350    {
2351       anti_collapse_on = ec_dec_bits(dec, 1);
2352    }
2353
2354    unquant_energy_finalise(st->mode, st->start, st->end, oldBandE,
2355          fine_quant, fine_priority, len*8-ec_dec_tell(dec, 0), dec, C);
2356
2357    if (anti_collapse_on)
2358       anti_collapse(st->mode, X, collapse_masks, LM, C, CC, N,
2359             st->start, st->end, oldBandE, oldLogE, oldLogE2, pulses, st->rng);
2360
2361    log2Amp(st->mode, st->start, st->end, bandE, oldBandE, C);
2362
2363    if (silence)
2364    {
2365       for (i=0;i<C*st->mode->nbEBands;i++)
2366       {
2367          bandE[i] = 0;
2368          oldBandE[i] = -QCONST16(28.f,DB_SHIFT);
2369       }
2370    }
2371    /* Synthesis */
2372    denormalise_bands(st->mode, X, freq, bandE, effEnd, C, M);
2373
2374    CELT_MOVE(decode_mem[0], decode_mem[0]+N, DECODE_BUFFER_SIZE-N);
2375    if (CC==2)
2376       CELT_MOVE(decode_mem[1], decode_mem[1]+N, DECODE_BUFFER_SIZE-N);
2377
2378    c=0; do
2379       for (i=0;i<M*st->mode->eBands[st->start];i++)
2380          freq[c*N+i] = 0;
2381    while (++c<C);
2382    c=0; do {
2383       int bound = M*st->mode->eBands[effEnd];
2384       if (st->downsample!=1)
2385          bound = IMIN(bound, N/st->downsample);
2386       for (i=M*st->mode->eBands[effEnd];i<N;i++)
2387          freq[c*N+i] = 0;
2388    } while (++c<C);
2389
2390    out_syn[0] = out_mem[0]+MAX_PERIOD-N;
2391    if (CC==2)
2392       out_syn[1] = out_mem[1]+MAX_PERIOD-N;
2393
2394    if (CC==2&&C==1)
2395    {
2396       for (i=0;i<N;i++)
2397          freq[N+i] = freq[i];
2398    }
2399
2400    /* Compute inverse MDCTs */
2401    compute_inv_mdcts(st->mode, shortBlocks, freq, out_syn, overlap_mem, CC, LM);
2402
2403 #ifdef ENABLE_POSTFILTER
2404    c=0; do {
2405       st->postfilter_period=IMAX(st->postfilter_period, COMBFILTER_MINPERIOD);
2406       st->postfilter_period_old=IMAX(st->postfilter_period_old, COMBFILTER_MINPERIOD);
2407       if (LM!=0)
2408       {
2409          comb_filter(out_syn[c], out_syn[c], st->postfilter_period, st->postfilter_period, st->overlap,
2410                st->postfilter_gain, st->postfilter_gain, st->postfilter_tapset, st->postfilter_tapset,
2411                NULL, 0);
2412          comb_filter(out_syn[c]+st->overlap, out_syn[c]+st->overlap, st->postfilter_period, postfilter_pitch, N-st->overlap,
2413                st->postfilter_gain, postfilter_gain, st->postfilter_tapset, postfilter_tapset,
2414                st->mode->window, st->mode->overlap);
2415       } else {
2416          comb_filter(out_syn[c], out_syn[c], st->postfilter_period_old, st->postfilter_period, N-st->overlap,
2417                st->postfilter_gain_old, st->postfilter_gain, st->postfilter_tapset_old, st->postfilter_tapset,
2418                st->mode->window, st->mode->overlap);
2419       }
2420    } while (++c<CC);
2421    st->postfilter_period_old = st->postfilter_period;
2422    st->postfilter_gain_old = st->postfilter_gain;
2423    st->postfilter_tapset_old = st->postfilter_tapset;
2424    st->postfilter_period = postfilter_pitch;
2425    st->postfilter_gain = postfilter_gain;
2426    st->postfilter_tapset = postfilter_tapset;
2427 #endif /* ENABLE_POSTFILTER */
2428
2429    if (CC==2&&C==1) {
2430       for (i=0;i<st->mode->nbEBands;i++)
2431          oldBandE[st->mode->nbEBands+i]=oldBandE[i];
2432    }
2433
2434    /* In case start or end were to change */
2435    c=0; do
2436    {
2437       for (i=0;i<st->start;i++)
2438          oldBandE[c*st->mode->nbEBands+i]=0;
2439       for (i=st->end;i<st->mode->nbEBands;i++)
2440          oldBandE[c*st->mode->nbEBands+i]=0;
2441    } while (++c<CC);
2442    if (!isTransient)
2443    {
2444       for (i=0;i<CC*st->mode->nbEBands;i++)
2445          oldLogE2[i] = oldLogE[i];
2446       for (i=0;i<CC*st->mode->nbEBands;i++)
2447          oldLogE[i] = oldBandE[i];
2448       for (i=0;i<CC*st->mode->nbEBands;i++)
2449          backgroundLogE[i] = MIN16(backgroundLogE[i] + M*QCONST16(0.001f,DB_SHIFT), oldBandE[i]);
2450    } else {
2451       for (i=0;i<CC*st->mode->nbEBands;i++)
2452          oldLogE[i] = MIN16(oldLogE[i], oldBandE[i]);
2453    }
2454    st->rng = dec->rng;
2455
2456    deemphasis(out_syn, pcm, N, CC, st->downsample, st->mode->preemph, st->preemph_memD);
2457    st->loss_count = 0;
2458    RESTORE_STACK;
2459    if (ec_dec_tell(dec,0) > 8*len || ec_dec_get_error(dec))
2460       return CELT_CORRUPTED_DATA;
2461    else
2462       return CELT_OK;
2463 }
2464
2465 #ifdef FIXED_POINT
2466 #ifndef DISABLE_FLOAT_API
2467 int celt_decode_with_ec_float(CELTDecoder * restrict st, const unsigned char *data, int len, float * restrict pcm, int frame_size, ec_dec *dec)
2468 {
2469    int j, ret, C, N;
2470    VARDECL(celt_int16, out);
2471    SAVE_STACK;
2472
2473    if (pcm==NULL)
2474       return CELT_BAD_ARG;
2475
2476    C = CHANNELS(st->channels);
2477    N = frame_size;
2478    
2479    ALLOC(out, C*N, celt_int16);
2480    ret=celt_decode_with_ec(st, data, len, out, frame_size, dec);
2481    if (ret==0)
2482       for (j=0;j<C*N;j++)
2483          pcm[j]=out[j]*(1.f/32768.f);
2484      
2485    RESTORE_STACK;
2486    return ret;
2487 }
2488 #endif /*DISABLE_FLOAT_API*/
2489 #else
2490 int celt_decode_with_ec(CELTDecoder * restrict st, const unsigned char *data, int len, celt_int16 * restrict pcm, int frame_size, ec_dec *dec)
2491 {
2492    int j, ret, C, N;
2493    VARDECL(celt_sig, out);
2494    SAVE_STACK;
2495
2496    if (pcm==NULL)
2497       return CELT_BAD_ARG;
2498
2499    C = CHANNELS(st->channels);
2500    N = frame_size;
2501    ALLOC(out, C*N, celt_sig);
2502
2503    ret=celt_decode_with_ec_float(st, data, len, out, frame_size, dec);
2504
2505    if (ret==0)
2506       for (j=0;j<C*N;j++)
2507          pcm[j] = FLOAT2INT16 (out[j]);
2508    
2509    RESTORE_STACK;
2510    return ret;
2511 }
2512 #endif
2513
2514 int celt_decode(CELTDecoder * restrict st, const unsigned char *data, int len, celt_int16 * restrict pcm, int frame_size)
2515 {
2516    return celt_decode_with_ec(st, data, len, pcm, frame_size, NULL);
2517 }
2518
2519 #ifndef DISABLE_FLOAT_API
2520 int celt_decode_float(CELTDecoder * restrict st, const unsigned char *data, int len, float * restrict pcm, int frame_size)
2521 {
2522    return celt_decode_with_ec_float(st, data, len, pcm, frame_size, NULL);
2523 }
2524 #endif /* DISABLE_FLOAT_API */
2525
2526 int celt_decoder_ctl(CELTDecoder * restrict st, int request, ...)
2527 {
2528    va_list ap;
2529
2530    va_start(ap, request);
2531    switch (request)
2532    {
2533       case CELT_GET_MODE_REQUEST:
2534       {
2535          const CELTMode ** value = va_arg(ap, const CELTMode**);
2536          if (value==0)
2537             goto bad_arg;
2538          *value=st->mode;
2539       }
2540       break;
2541       case CELT_SET_START_BAND_REQUEST:
2542       {
2543          celt_int32 value = va_arg(ap, celt_int32);
2544          if (value<0 || value>=st->mode->nbEBands)
2545             goto bad_arg;
2546          st->start = value;
2547       }
2548       break;
2549       case CELT_SET_END_BAND_REQUEST:
2550       {
2551          celt_int32 value = va_arg(ap, celt_int32);
2552          if (value<0 || value>=st->mode->nbEBands)
2553             goto bad_arg;
2554          st->end = value;
2555       }
2556       break;
2557       case CELT_SET_CHANNELS_REQUEST:
2558       {
2559          celt_int32 value = va_arg(ap, celt_int32);
2560          if (value<1 || value>2)
2561             goto bad_arg;
2562          st->stream_channels = value;
2563       }
2564       break;
2565       case CELT_RESET_STATE:
2566       {
2567          CELT_MEMSET((char*)&st->DECODER_RESET_START, 0,
2568                celt_decoder_get_size_custom(st->mode, st->channels)-
2569                ((char*)&st->DECODER_RESET_START - (char*)st));
2570       }
2571       break;
2572       default:
2573          goto bad_request;
2574    }
2575    va_end(ap);
2576    return CELT_OK;
2577 bad_arg:
2578    va_end(ap);
2579    return CELT_BAD_ARG;
2580 bad_request:
2581       va_end(ap);
2582   return CELT_UNIMPLEMENTED;
2583 }
2584
2585 const char *celt_strerror(int error)
2586 {
2587    static const char *error_strings[8] = {
2588       "success",
2589       "invalid argument",
2590       "invalid mode",
2591       "internal error",
2592       "corrupted stream",
2593       "request not implemented",
2594       "invalid state",
2595       "memory allocation failed"
2596    };
2597    if (error > 0 || error < -7)
2598       return "unknown error";
2599    else 
2600       return error_strings[-error];
2601 }
2602