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