Fixes some side-information rate control issues in VBR mode
[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
56 static const int trim_cdf[12] = {0, 2, 4, 9, 19, 41, 87, 109, 119, 124, 126, 128};
57
58 #define COMBFILTER_MAXPERIOD 1024
59 #define COMBFILTER_MINPERIOD 16
60
61 /** Encoder state 
62  @brief Encoder state
63  */
64 struct CELTEncoder {
65    const CELTMode *mode;     /**< Mode used by the encoder */
66    int overlap;
67    int channels;
68    
69    int force_intra;
70    int complexity;
71    int start, end;
72
73    celt_int32 vbr_rate_norm; /* Target number of 8th bits per frame */
74
75    /* Everything beyond this point gets cleared on a reset */
76 #define ENCODER_RESET_START frame_max
77
78    celt_word32 frame_max;
79    int fold_decision;
80    int delayedIntra;
81    int tonal_average;
82
83    int prefilter_period;
84    celt_word16 prefilter_gain;
85
86    /* VBR-related parameters */
87    celt_int32 vbr_reservoir;
88    celt_int32 vbr_drift;
89    celt_int32 vbr_offset;
90    celt_int32 vbr_count;
91
92    celt_word32 preemph_memE[2];
93    celt_word32 preemph_memD[2];
94
95 #ifdef RESYNTH
96    celt_sig syn_mem[2][2*MAX_PERIOD];
97 #endif
98
99    celt_sig in_mem[1]; /* Size = channels*mode->overlap */
100    /* celt_sig prefilter_mem[],  Size = channels*COMBFILTER_PERIOD */
101    /* celt_sig overlap_mem[],  Size = channels*mode->overlap */
102    /* celt_word16 oldEBands[], Size = channels*mode->nbEBands */
103 };
104
105 int celt_encoder_get_size(const CELTMode *mode, int channels)
106 {
107    int size = sizeof(struct CELTEncoder)
108          + (2*channels*mode->overlap-1)*sizeof(celt_sig)
109          + channels*COMBFILTER_MAXPERIOD*sizeof(celt_sig)
110          + channels*mode->nbEBands*sizeof(celt_word16);
111    return size;
112 }
113
114 CELTEncoder *celt_encoder_create(const CELTMode *mode, int channels, int *error)
115 {
116    return celt_encoder_init(
117          (CELTEncoder *)celt_alloc(celt_encoder_get_size(mode, channels)),
118          mode, channels, error);
119 }
120
121 CELTEncoder *celt_encoder_init(CELTEncoder *st, const CELTMode *mode, int channels, int *error)
122 {
123    if (channels < 0 || channels > 2)
124    {
125       if (error)
126          *error = CELT_BAD_ARG;
127       return NULL;
128    }
129
130    if (st==NULL)
131    {
132       if (error)
133          *error = CELT_ALLOC_FAIL;
134       return NULL;
135    }
136
137    CELT_MEMSET((char*)st, 0, celt_encoder_get_size(mode, channels));
138    
139    st->mode = mode;
140    st->overlap = mode->overlap;
141    st->channels = channels;
142
143    st->start = 0;
144    st->end = st->mode->effEBands;
145
146    st->vbr_rate_norm = 0;
147    st->vbr_offset = -(64<<BITRES);
148    st->force_intra  = 0;
149    st->delayedIntra = 1;
150    st->tonal_average = 256;
151    st->fold_decision = 1;
152    st->complexity = 5;
153
154    if (error)
155       *error = CELT_OK;
156    return st;
157 }
158
159 void celt_encoder_destroy(CELTEncoder *st)
160 {
161    celt_free(st);
162 }
163
164 static inline celt_int16 FLOAT2INT16(float x)
165 {
166    x = x*CELT_SIG_SCALE;
167    x = MAX32(x, -32768);
168    x = MIN32(x, 32767);
169    return (celt_int16)float2int(x);
170 }
171
172 static inline celt_word16 SIG2WORD16(celt_sig x)
173 {
174 #ifdef FIXED_POINT
175    x = PSHR32(x, SIG_SHIFT);
176    x = MAX32(x, -32768);
177    x = MIN32(x, 32767);
178    return EXTRACT16(x);
179 #else
180    return (celt_word16)x;
181 #endif
182 }
183
184 static int transient_analysis(const celt_word32 * restrict in, int len, int C,
185                               celt_word32 *frame_max, int overlap)
186 {
187    int i, n;
188    celt_word32 threshold;
189    VARDECL(celt_word32, begin);
190    VARDECL(celt_word16, tmp);
191    celt_word32 mem0=0,mem1=0;
192    SAVE_STACK;
193    ALLOC(tmp, len, celt_word16);
194    ALLOC(begin, len+1, celt_word32);
195
196    if (C==1)
197    {
198       for (i=0;i<len;i++)
199          tmp[i] = SHR32(in[i],SIG_SHIFT);
200    } else {
201       for (i=0;i<len;i++)
202          tmp[i] = SHR32(ADD32(in[i],in[i+len]), SIG_SHIFT+1);
203    }
204
205    /* High-pass filter: (1 - 2*z^-1 + z^-2) / (1 - z^-1 + .5*z^-2) */
206    for (i=0;i<len;i++)
207    {
208       celt_word32 x,y;
209       x = tmp[i];
210       y = ADD32(mem0, x);
211 #ifdef FIXED_POINT
212       mem0 = mem1 + y - SHL32(x,1);
213       mem1 = x - SHR32(y,1);
214 #else
215       mem0 = mem1 + y - 2*x;
216       mem1 = x - .5f*y;
217 #endif
218       tmp[i] = EXTRACT16(SHR(y,2));
219    }
220    /* First few samples are bad because we don't propagate the memory */
221    for (i=0;i<24;i++)
222       tmp[i] = 0;
223
224    begin[0] = 0;
225    for (i=0;i<len;i++)
226       begin[i+1] = MAX32(begin[i], ABS32(tmp[i]));
227
228    n = -1;
229
230    threshold = MULT16_32_Q15(QCONST16(.4f,15),begin[len]);
231    /* If the following condition isn't met, there's just no way
232       we'll have a transient*/
233    if (*frame_max < threshold)
234    {
235       /* It's likely we have a transient, now find it */
236       for (i=8;i<len-8;i++)
237       {
238          if (begin[i+1] < threshold)
239             n=i;
240       }
241    }
242
243    *frame_max = begin[len-overlap];
244    /* Only consider the last 7.5 ms for the next transient */
245    if (len>360+overlap)
246    {
247       *frame_max = 0;
248       for (i=len-360-overlap;i<len-overlap;i++)
249          *frame_max = MAX32(*frame_max, ABS32(tmp[i]));
250    }
251    RESTORE_STACK;
252    return n>=32;
253 }
254
255 /** Apply window and compute the MDCT for all sub-frames and 
256     all channels in a frame */
257 static void compute_mdcts(const CELTMode *mode, int shortBlocks, celt_sig * restrict in, celt_sig * restrict out, int _C, int LM)
258 {
259    const int C = CHANNELS(_C);
260    if (C==1 && !shortBlocks)
261    {
262       const int overlap = OVERLAP(mode);
263       clt_mdct_forward(&mode->mdct, in, out, mode->window, overlap, mode->maxLM-LM);
264    } else {
265       const int overlap = OVERLAP(mode);
266       int N = mode->shortMdctSize<<LM;
267       int B = 1;
268       int b, c;
269       VARDECL(celt_word32, tmp);
270       SAVE_STACK;
271       if (shortBlocks)
272       {
273          /*lookup = &mode->mdct[0];*/
274          N = mode->shortMdctSize;
275          B = shortBlocks;
276       }
277       ALLOC(tmp, N, celt_word32);
278       c=0; do {
279          for (b=0;b<B;b++)
280          {
281             int j;
282             clt_mdct_forward(&mode->mdct, in+c*(B*N+overlap)+b*N, tmp, mode->window, overlap, shortBlocks ? mode->maxLM : mode->maxLM-LM);
283             /* Interleaving the sub-frames */
284             for (j=0;j<N;j++)
285                out[(j*B+b)+c*N*B] = tmp[j];
286          }
287       } while (++c<C);
288       RESTORE_STACK;
289    }
290 }
291
292 /** Compute the IMDCT and apply window for all sub-frames and 
293     all channels in a frame */
294 static void compute_inv_mdcts(const CELTMode *mode, int shortBlocks, celt_sig *X,
295       celt_sig * restrict out_mem[],
296       celt_sig * restrict overlap_mem[], int _C, int LM)
297 {
298    int c;
299    const int C = CHANNELS(_C);
300    const int N = mode->shortMdctSize<<LM;
301    const int overlap = OVERLAP(mode);
302    c=0; do {
303       int j;
304          VARDECL(celt_word32, x);
305          VARDECL(celt_word32, tmp);
306          int b;
307          int N2 = N;
308          int B = 1;
309          SAVE_STACK;
310          
311          ALLOC(x, N+overlap, celt_word32);
312          ALLOC(tmp, N, celt_word32);
313
314          if (shortBlocks)
315          {
316             N2 = mode->shortMdctSize;
317             B = shortBlocks;
318          }
319          /* Prevents problems from the imdct doing the overlap-add */
320          CELT_MEMSET(x, 0, overlap);
321
322          for (b=0;b<B;b++)
323          {
324             /* De-interleaving the sub-frames */
325             for (j=0;j<N2;j++)
326                tmp[j] = X[(j*B+b)+c*N2*B];
327             clt_mdct_backward(&mode->mdct, tmp, x+N2*b, mode->window, overlap, shortBlocks ? mode->maxLM : mode->maxLM-LM);
328          }
329
330          for (j=0;j<overlap;j++)
331             out_mem[c][j] = x[j] + overlap_mem[c][j];
332          for (;j<N;j++)
333             out_mem[c][j] = x[j];
334          for (j=0;j<overlap;j++)
335             overlap_mem[c][j] = x[N+j];
336          RESTORE_STACK;
337    } while (++c<C);
338 }
339
340 static void deemphasis(celt_sig *in[], celt_word16 *pcm, int N, int _C, const celt_word16 *coef, celt_sig *mem)
341 {
342    const int C = CHANNELS(_C);
343    int c;
344    c=0; do {
345       int j;
346       celt_sig * restrict x;
347       celt_word16  * restrict y;
348       celt_sig m = mem[c];
349       x =in[c];
350       y = pcm+c;
351       for (j=0;j<N;j++)
352       {
353          celt_sig tmp = *x + m;
354          m = MULT16_32_Q15(coef[0], tmp)
355            - MULT16_32_Q15(coef[1], *x);
356          tmp = SHL32(MULT16_32_Q15(coef[3], tmp), 2);
357          *y = SCALEOUT(SIG2WORD16(tmp));
358          x++;
359          y+=C;
360       }
361       mem[c] = m;
362    } while (++c<C);
363 }
364
365 #ifdef ENABLE_POSTFILTER
366 /* FIXME: Handle the case where T = maxperiod */
367 static void comb_filter(celt_word32 *y, celt_word32 *x, int T0, int T1, int N,
368       int C, celt_word16 g0, celt_word16 g1, const celt_word16 *window, int overlap)
369 {
370    int i;
371    /* printf ("%d %d %f %f\n", T0, T1, g0, g1); */
372    celt_word16 g00, g01, g02, g10, g11, g12;
373    celt_word16 t0, t1, t2;
374    /* zeros at theta = +/- 5*pi/6 */
375    t0 = QCONST16(.26795f, 15);
376    t1 = QCONST16(.46410f, 15);
377    t2 = QCONST16(.26795f, 15);
378    g00 = MULT16_16_Q15(g0, t0);
379    g01 = MULT16_16_Q15(g0, t1);
380    g02 = MULT16_16_Q15(g0, t2);
381    g10 = MULT16_16_Q15(g1, t0);
382    g11 = MULT16_16_Q15(g1, t1);
383    g12 = MULT16_16_Q15(g1, t2);
384    for (i=0;i<overlap;i++)
385    {
386       celt_word16 f;
387       f = MULT16_16_Q15(window[i],window[i]);
388       y[i] = x[i]
389                + MULT16_32_Q15(MULT16_16_Q15((Q15ONE-f),g01),x[i-T0])
390                + MULT16_32_Q15(MULT16_16_Q15((Q15ONE-f),g00),x[i-T0-1])
391                + MULT16_32_Q15(MULT16_16_Q15((Q15ONE-f),g02),x[i-T0+1])
392                + MULT16_32_Q15(MULT16_16_Q15(f,g11),x[i-T1])
393                + MULT16_32_Q15(MULT16_16_Q15(f,g10),x[i-T1-1])
394                + MULT16_32_Q15(MULT16_16_Q15(f,g12),x[i-T1+1]);
395
396    }
397    for (i=overlap;i<N;i++)
398       y[i] = x[i]
399                + MULT16_32_Q15(g11,x[i-T1])
400                + MULT16_32_Q15(g10,x[i-T1-1])
401                + MULT16_32_Q15(g12,x[i-T1+1]);
402 }
403 #endif /* ENABLE_POSTFILTER */
404
405 static const signed char tf_select_table[4][8] = {
406       {0, -1, 0, -1,    0,-1, 0,-1},
407       {0, -1, 0, -2,    1, 0, 1 -1},
408       {0, -2, 0, -3,    2, 0, 1 -1},
409       {0, -2, 0, -3,    2, 0, 1 -1},
410 };
411
412 static celt_word32 l1_metric(const celt_norm *tmp, int N, int LM, int width)
413 {
414    int i, j;
415    static const celt_word16 sqrtM_1[4] = {Q15ONE, QCONST16(0.70711f,15), QCONST16(0.5f,15), QCONST16(0.35355f,15)};
416    celt_word32 L1;
417    celt_word16 bias;
418    L1=0;
419    for (i=0;i<1<<LM;i++)
420    {
421       celt_word32 L2 = 0;
422       for (j=0;j<N>>LM;j++)
423          L2 = MAC16_16(L2, tmp[(j<<LM)+i], tmp[(j<<LM)+i]);
424       L1 += celt_sqrt(L2);
425    }
426    L1 = MULT16_32_Q15(sqrtM_1[LM], L1);
427    if (width==1)
428       bias = QCONST16(.12f,15)*LM;
429    else if (width==2)
430       bias = QCONST16(.05f,15)*LM;
431    else
432       bias = QCONST16(.02f,15)*LM;
433    L1 = MAC16_32_Q15(L1, bias, L1);
434    return L1;
435 }
436
437 static int tf_analysis(const CELTMode *m, celt_word16 *bandLogE, celt_word16 *oldBandE,
438       int len, int C, int isTransient, int *tf_res, int nbCompressedBytes, celt_norm *X,
439       int N0, int LM, int *tf_sum)
440 {
441    int i;
442    VARDECL(int, metric);
443    int cost0;
444    int cost1;
445    VARDECL(int, path0);
446    VARDECL(int, path1);
447    VARDECL(celt_norm, tmp);
448    int lambda;
449    int tf_select=0;
450    SAVE_STACK;
451
452    /* FIXME: Should check number of bytes *left* */
453    if (nbCompressedBytes<15*C)
454    {
455       *tf_sum = 0;
456       for (i=0;i<len;i++)
457          tf_res[i] = isTransient;
458       return 0;
459    }
460    if (nbCompressedBytes<40)
461       lambda = 12;
462    else if (nbCompressedBytes<60)
463       lambda = 6;
464    else if (nbCompressedBytes<100)
465       lambda = 4;
466    else
467       lambda = 3;
468
469    ALLOC(metric, len, int);
470    ALLOC(tmp, (m->eBands[len]-m->eBands[len-1])<<LM, celt_norm);
471    ALLOC(path0, len, int);
472    ALLOC(path1, len, int);
473
474    *tf_sum = 0;
475    for (i=0;i<len;i++)
476    {
477       int j, k, N;
478       celt_word32 L1, best_L1;
479       int best_level=0;
480       N = (m->eBands[i+1]-m->eBands[i])<<LM;
481       for (j=0;j<N;j++)
482          tmp[j] = X[j+(m->eBands[i]<<LM)];
483       /* FIXME: Do something with the right channel */
484       if (C==2)
485          for (j=0;j<N;j++)
486             tmp[j] = ADD16(tmp[j],X[N0+j+(m->eBands[i]<<LM)]);
487       L1 = l1_metric(tmp, N, isTransient ? LM : 0, N>>LM);
488       best_L1 = L1;
489       /*printf ("%f ", L1);*/
490       for (k=0;k<LM;k++)
491       {
492          int B;
493
494          if (isTransient)
495             B = (LM-k-1);
496          else
497             B = k+1;
498
499          if (isTransient)
500             haar1(tmp, N>>(LM-k), 1<<(LM-k));
501          else
502             haar1(tmp, N>>k, 1<<k);
503
504          L1 = l1_metric(tmp, N, B, N>>LM);
505
506          if (L1 < best_L1)
507          {
508             best_L1 = L1;
509             best_level = k+1;
510          }
511       }
512       /*printf ("%d ", isTransient ? LM-best_level : best_level);*/
513       if (isTransient)
514          metric[i] = best_level;
515       else
516          metric[i] = -best_level;
517       *tf_sum += metric[i];
518    }
519    /*printf("\n");*/
520    /* FIXME: Figure out how to set this */
521    tf_select = 0;
522
523    cost0 = 0;
524    cost1 = isTransient ? 0 : lambda;
525    /* Viterbi forward pass */
526    for (i=1;i<len;i++)
527    {
528       int curr0, curr1;
529       int from0, from1;
530
531       from0 = cost0;
532       from1 = cost1 + lambda;
533       if (from0 < from1)
534       {
535          curr0 = from0;
536          path0[i]= 0;
537       } else {
538          curr0 = from1;
539          path0[i]= 1;
540       }
541
542       from0 = cost0 + lambda;
543       from1 = cost1;
544       if (from0 < from1)
545       {
546          curr1 = from0;
547          path1[i]= 0;
548       } else {
549          curr1 = from1;
550          path1[i]= 1;
551       }
552       cost0 = curr0 + abs(metric[i]-tf_select_table[LM][4*isTransient+2*tf_select+0]);
553       cost1 = curr1 + abs(metric[i]-tf_select_table[LM][4*isTransient+2*tf_select+1]);
554    }
555    tf_res[len-1] = cost0 < cost1 ? 0 : 1;
556    /* Viterbi backward pass to check the decisions */
557    for (i=len-2;i>=0;i--)
558    {
559       if (tf_res[i+1] == 1)
560          tf_res[i] = path1[i+1];
561       else
562          tf_res[i] = path0[i+1];
563    }
564    RESTORE_STACK;
565    return tf_select;
566 }
567
568 static void tf_encode(int start, int end, int isTransient, int *tf_res, int LM, int tf_select, ec_enc *enc)
569 {
570    int curr, i;
571    ec_enc_bit_prob(enc, tf_res[start], isTransient ? 16384 : 4096);
572    curr = tf_res[start];
573    for (i=start+1;i<end;i++)
574    {
575       ec_enc_bit_prob(enc, tf_res[i] ^ curr, isTransient ? 4096 : 2048);
576       curr = tf_res[i];
577    }
578    ec_enc_bits(enc, tf_select, 1);
579    for (i=start;i<end;i++)
580       tf_res[i] = tf_select_table[LM][4*isTransient+2*tf_select+tf_res[i]];
581    /*printf("%d %d ", isTransient, tf_select); for(i=0;i<end;i++)printf("%d ", tf_res[i]);printf("\n");*/
582 }
583
584 static void tf_decode(int start, int end, int C, int isTransient, int *tf_res, int LM, ec_dec *dec)
585 {
586    int i, curr, tf_select;
587    tf_res[start] = ec_dec_bit_prob(dec, isTransient ? 16384 : 4096);
588    curr = tf_res[start];
589    for (i=start+1;i<end;i++)
590    {
591       tf_res[i] = ec_dec_bit_prob(dec, isTransient ? 4096 : 2048) ^ curr;
592       curr = tf_res[i];
593    }
594    tf_select = ec_dec_bits(dec, 1);
595    for (i=start;i<end;i++)
596       tf_res[i] = tf_select_table[LM][4*isTransient+2*tf_select+tf_res[i]];
597 }
598
599 static int alloc_trim_analysis(const CELTMode *m, const celt_norm *X,
600       const celt_word16 *bandLogE, int nbEBands, int LM, int C, int N0)
601 {
602    int i;
603    int trim_index = 5;
604    if (C==2)
605    {
606       celt_word16 sum = 0; /* Q10 */
607       /* Compute inter-channel correlation for low frequencies */
608       for (i=0;i<8;i++)
609       {
610          int j;
611          celt_word32 partial = 0;
612          for (j=m->eBands[i]<<LM;j<m->eBands[i+1]<<LM;j++)
613             partial = MAC16_16(partial, X[j], X[N0+j]);
614          sum = ADD16(sum, EXTRACT16(SHR32(partial, 18)));
615       }
616       sum = MULT16_16_Q15(QCONST16(1.f/8, 15), sum);
617       /*printf ("%f\n", sum);*/
618       if (sum > QCONST16(.995f,10))
619          trim_index-=4;
620       else if (sum > QCONST16(.92f,10))
621          trim_index-=3;
622       else if (sum > QCONST16(.85f,10))
623          trim_index-=2;
624       else if (sum > QCONST16(.8f,10))
625          trim_index-=1;
626    }
627 #if 0
628    float diff=0;
629    int c;
630    c=0; do {
631       for (i=0;i<nbEBands-1;i++)
632       {
633          diff += bandLogE[i+c*nbEBands]*(i-.5*nbEBands);
634       }
635    } while (++c<0);
636    diff /= C*(nbEBands-1);
637    /*printf("%f\n", diff);*/
638    if (diff > 4)
639       trim_index--;
640    if (diff > 8)
641       trim_index--;
642    if (diff < -4)
643       trim_index++;
644 #endif
645    if (trim_index<0)
646       trim_index = 0;
647    if (trim_index>5)
648       trim_index = 5;
649    return trim_index;
650 }
651
652 static int stereo_analysis(const CELTMode *m, const celt_norm *X,
653       int nbEBands, int LM, int C, int N0)
654 {
655    int i;
656    int thetas;
657    celt_word32 sumLR = EPSILON, sumMS = EPSILON;
658
659    /* Use the L1 norm to model the entropy of the L/R signal vs the M/S signal */
660    for (i=0;i<13;i++)
661    {
662       int j;
663       for (j=m->eBands[i]<<LM;j<m->eBands[i+1]<<LM;j++)
664       {
665          celt_word16 L, R, M, S;
666          L = X[j];
667          R = X[N0+j];
668          M = L+R;
669          S = L-R;
670          sumLR += EXTEND32(ABS16(L)) + EXTEND32(ABS16(R));
671          sumMS += EXTEND32(ABS16(M)) + EXTEND32(ABS16(S));
672       }
673    }
674    sumMS = MULT16_32_Q15(QCONST16(0.707107f, 15), sumMS);
675    thetas = 13;
676    /* We don't need thetas for lower bands with LM<=1 */
677    if (LM<=1)
678       thetas -= 8;
679    return MULT16_32_Q15((m->eBands[13]<<(LM+1))+thetas, sumMS)
680          > MULT16_32_Q15(m->eBands[13]<<(LM+1), sumLR);
681 }
682
683 #ifdef FIXED_POINT
684 int celt_encode_with_ec(CELTEncoder * restrict st, const celt_int16 * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes, ec_enc *enc)
685 {
686 #else
687 int celt_encode_with_ec_float(CELTEncoder * restrict st, const celt_sig * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes, ec_enc *enc)
688 {
689 #endif
690    int i, c, N;
691    int bits;
692    int has_fold=1;
693    ec_byte_buffer buf;
694    ec_enc         _enc;
695    VARDECL(celt_sig, in);
696    VARDECL(celt_sig, freq);
697    VARDECL(celt_norm, X);
698    VARDECL(celt_ener, bandE);
699    VARDECL(celt_word16, bandLogE);
700    VARDECL(int, fine_quant);
701    VARDECL(celt_word16, error);
702    VARDECL(int, pulses);
703    VARDECL(int, offsets);
704    VARDECL(int, fine_priority);
705    VARDECL(int, tf_res);
706    celt_sig *_overlap_mem;
707    celt_sig *prefilter_mem;
708    celt_word16 *oldBandE;
709    int shortBlocks=0;
710    int isTransient=0;
711    int resynth;
712    const int C = CHANNELS(st->channels);
713    int LM, M;
714    int tf_select;
715    int nbFilledBytes, nbAvailableBytes;
716    int effEnd;
717    int codedBands;
718    int tf_sum;
719    int alloc_trim;
720    int pitch_index=0;
721    celt_word16 gain1 = 0;
722    int intensity=0;
723    int dual_stereo=0;
724    int effectiveBytes;
725    SAVE_STACK;
726
727    if (nbCompressedBytes<0 || pcm==NULL)
728      return CELT_BAD_ARG;
729
730    for (LM=0;LM<4;LM++)
731       if (st->mode->shortMdctSize<<LM==frame_size)
732          break;
733    if (LM>=MAX_CONFIG_SIZES)
734       return CELT_BAD_ARG;
735    M=1<<LM;
736
737    prefilter_mem = st->in_mem+C*(st->overlap);
738    _overlap_mem = prefilter_mem+C*COMBFILTER_MAXPERIOD;
739    /*_overlap_mem = st->in_mem+C*(st->overlap);*/
740    oldBandE = (celt_word16*)(st->in_mem+C*(2*st->overlap+COMBFILTER_MAXPERIOD));
741
742    if (enc==NULL)
743    {
744       ec_byte_writeinit_buffer(&buf, compressed, nbCompressedBytes);
745       ec_enc_init(&_enc,&buf);
746       enc = &_enc;
747       nbFilledBytes=0;
748    } else {
749       nbFilledBytes=(ec_enc_tell(enc, 0)+4)>>3;
750    }
751    nbAvailableBytes = nbCompressedBytes - nbFilledBytes;
752
753    if (st->vbr_rate_norm>0)
754       effectiveBytes = st->vbr_rate_norm>>BITRES<<LM>>3;
755    else
756       effectiveBytes = nbCompressedBytes;
757
758    effEnd = st->end;
759    if (effEnd > st->mode->effEBands)
760       effEnd = st->mode->effEBands;
761
762    N = M*st->mode->shortMdctSize;
763    ALLOC(in, C*(N+st->overlap), celt_sig);
764
765    /* Find pitch period and gain */
766    {
767       VARDECL(celt_sig, _pre);
768       celt_sig *pre[2];
769       SAVE_STACK;
770       c = 0;
771       ALLOC(_pre, C*(N+COMBFILTER_MAXPERIOD), celt_sig);
772
773       pre[0] = _pre;
774       pre[1] = _pre + (N+COMBFILTER_MAXPERIOD);
775
776       c=0; do {
777          const celt_word16 * restrict pcmp = pcm+c;
778          celt_sig * restrict inp = in+c*(N+st->overlap)+st->overlap;
779
780          for (i=0;i<N;i++)
781          {
782             /* Apply pre-emphasis */
783             celt_sig tmp = MULT16_16(st->mode->preemph[2], SCALEIN(*pcmp));
784             *inp = tmp + st->preemph_memE[c];
785             st->preemph_memE[c] = MULT16_32_Q15(st->mode->preemph[1], *inp)
786                                    - MULT16_32_Q15(st->mode->preemph[0], tmp);
787             inp++;
788             pcmp+=C;
789          }
790          CELT_COPY(pre[c], prefilter_mem+c*COMBFILTER_MAXPERIOD, COMBFILTER_MAXPERIOD);
791          CELT_COPY(pre[c]+COMBFILTER_MAXPERIOD, in+c*(N+st->overlap)+st->overlap, N);
792       } while (++c<C);
793
794 #ifdef ENABLE_POSTFILTER
795       {
796          VARDECL(celt_word16, pitch_buf);
797          ALLOC(pitch_buf, (COMBFILTER_MAXPERIOD+N)>>1, celt_word16);
798          celt_word32 tmp=0;
799          celt_word32 mem0[2]={0,0};
800          celt_word16 mem1[2]={0,0};
801
802          pitch_downsample(pre, pitch_buf, COMBFILTER_MAXPERIOD+N, COMBFILTER_MAXPERIOD+N,
803                           C, mem0, mem1);
804          pitch_search(st->mode, pitch_buf+(COMBFILTER_MAXPERIOD>>1), pitch_buf, N,
805                COMBFILTER_MAXPERIOD-COMBFILTER_MINPERIOD, &pitch_index, &tmp, 1<<LM);
806          pitch_index = COMBFILTER_MAXPERIOD-pitch_index;
807
808          gain1 = remove_doubling(pitch_buf, COMBFILTER_MAXPERIOD, COMBFILTER_MINPERIOD,
809                N, &pitch_index, st->prefilter_period, st->prefilter_gain);
810       }
811       if (pitch_index > COMBFILTER_MAXPERIOD)
812          pitch_index = COMBFILTER_MAXPERIOD;
813       gain1 = MULT16_16_Q15(QCONST16(.7f,15),gain1);
814       if (gain1 > QCONST16(.6f,15))
815          gain1 = QCONST16(.6f,15);
816       if (ABS16(gain1-st->prefilter_gain)<QCONST16(.1,15))
817          gain1=st->prefilter_gain;
818       if (gain1<QCONST16(.2f,15))
819       {
820          ec_enc_bit_prob(enc, 0, 32768);
821          gain1 = 0;
822       } else {
823          int qg;
824          int octave;
825 #ifdef FIXED_POINT
826          qg = ((gain1+2048)>>12)-2;
827 #else
828          qg = floor(.5+gain1*8)-2;
829 #endif
830          ec_enc_bit_prob(enc, 1, 32768);
831          octave = EC_ILOG(pitch_index)-5;
832          ec_enc_uint(enc, octave, 6);
833          ec_enc_bits(enc, pitch_index-(16<<octave), 4+octave);
834          ec_enc_bits(enc, qg, 2);
835          gain1 = QCONST16(.125f,15)*(qg+2);
836       }
837       /*printf("%d %f\n", pitch_index, gain1);*/
838 #else /* ENABLE_POSTFILTER */
839       ec_enc_bit_prob(enc, 0, 32768);
840 #endif /* ENABLE_POSTFILTER */
841
842       c=0; do {
843          CELT_COPY(in+c*(N+st->overlap), st->in_mem+c*(st->overlap), st->overlap);
844 #ifdef ENABLE_POSTFILTER
845          comb_filter(in+c*(N+st->overlap)+st->overlap, pre[c]+COMBFILTER_MAXPERIOD,
846                st->prefilter_period, pitch_index, N, C, -st->prefilter_gain, -gain1, st->mode->window, st->mode->overlap);
847 #endif /* ENABLE_POSTFILTER */
848          CELT_COPY(st->in_mem+c*(st->overlap), in+c*(N+st->overlap)+N, st->overlap);
849
850 #ifdef ENABLE_POSTFILTER
851          if (N>COMBFILTER_MAXPERIOD)
852          {
853             CELT_MOVE(prefilter_mem+c*COMBFILTER_MAXPERIOD, pre[c]+N, COMBFILTER_MAXPERIOD);
854          } else {
855             CELT_MOVE(prefilter_mem+c*COMBFILTER_MAXPERIOD, prefilter_mem+c*COMBFILTER_MAXPERIOD+N, COMBFILTER_MAXPERIOD-N);
856             CELT_MOVE(prefilter_mem+c*COMBFILTER_MAXPERIOD+COMBFILTER_MAXPERIOD-N, pre[c]+COMBFILTER_MAXPERIOD, N);
857          }
858 #endif /* ENABLE_POSTFILTER */
859       } while (++c<C);
860
861       RESTORE_STACK;
862    }
863
864 #ifdef RESYNTH
865    resynth = 1;
866 #else
867    resynth = 0;
868 #endif
869
870    if (st->complexity > 1 && LM>0)
871    {
872       isTransient = M > 1 &&
873          transient_analysis(in, N+st->overlap, C, &st->frame_max, st->overlap);
874    } else {
875       isTransient = 0;
876    }
877
878    if (isTransient)
879       shortBlocks = M;
880    else
881       shortBlocks = 0;
882
883    ALLOC(freq, C*N, celt_sig); /**< Interleaved signal MDCTs */
884    ALLOC(bandE,st->mode->nbEBands*C, celt_ener);
885    ALLOC(bandLogE,st->mode->nbEBands*C, celt_word16);
886    /* Compute MDCTs */
887    compute_mdcts(st->mode, shortBlocks, in, freq, C, LM);
888
889    ALLOC(X, C*N, celt_norm);         /**< Interleaved normalised MDCTs */
890
891    compute_band_energies(st->mode, freq, bandE, effEnd, C, M);
892
893    amp2Log2(st->mode, effEnd, st->end, bandE, bandLogE, C);
894
895    /* Band normalisation */
896    normalise_bands(st->mode, freq, X, bandE, effEnd, C, M);
897
898    ALLOC(tf_res, st->mode->nbEBands, int);
899    /* Needs to be before coarse energy quantization because otherwise the energy gets modified */
900    tf_select = tf_analysis(st->mode, bandLogE, oldBandE, effEnd, C, isTransient, tf_res, effectiveBytes, X, N, LM, &tf_sum);
901    for (i=effEnd;i<st->end;i++)
902       tf_res[i] = tf_res[effEnd-1];
903
904    ALLOC(error, C*st->mode->nbEBands, celt_word16);
905    quant_coarse_energy(st->mode, st->start, st->end, effEnd, bandLogE,
906          oldBandE, nbCompressedBytes*8, error, enc,
907          C, LM, nbAvailableBytes, st->force_intra,
908          &st->delayedIntra, st->complexity >= 4);
909
910    if (LM > 0)
911       ec_enc_bit_prob(enc, shortBlocks!=0, 8192);
912
913    tf_encode(st->start, st->end, isTransient, tf_res, LM, tf_select, enc);
914
915    if (shortBlocks || st->complexity < 3)
916    {
917       if (st->complexity == 0)
918       {
919          has_fold = 0;
920          st->fold_decision = 3;
921       } else {
922          has_fold = 1;
923          st->fold_decision = 1;
924       }
925    } else {
926       has_fold = folding_decision(st->mode, X, &st->tonal_average, &st->fold_decision, effEnd, C, M);
927    }
928    ec_enc_bit_prob(enc, has_fold>>1, 8192);
929    ec_enc_bit_prob(enc, has_fold&1, (has_fold>>1) ? 32768 : 49152);
930
931    ALLOC(offsets, st->mode->nbEBands, int);
932
933    for (i=0;i<st->mode->nbEBands;i++)
934       offsets[i] = 0;
935    /* Dynamic allocation code */
936    /* Make sure that dynamic allocation can't make us bust the budget */
937    if (effectiveBytes > 50 && LM>=1)
938    {
939       int t1, t2;
940       if (LM <= 1)
941       {
942          t1 = 3;
943          t2 = 5;
944       } else {
945          t1 = 2;
946          t2 = 4;
947       }
948       for (i=1;i<st->mode->nbEBands-1;i++)
949       {
950          celt_word32 d2;
951          d2 = 2*bandLogE[i]-bandLogE[i-1]-bandLogE[i+1];
952          if (C==2)
953             d2 = HALF32(d2 + 2*bandLogE[i+st->mode->nbEBands]-
954                   bandLogE[i-1+st->mode->nbEBands]-bandLogE[i+1+st->mode->nbEBands]);
955          if (d2 > SHL16(t1,DB_SHIFT))
956             offsets[i] += 1;
957          if (d2 > SHL16(t2,DB_SHIFT))
958             offsets[i] += 1;
959       }
960    }
961    for (i=0;i<st->mode->nbEBands;i++)
962    {
963       int j;
964       ec_enc_bit_prob(enc, offsets[i]!=0, 1024);
965       if (offsets[i]!=0)
966       {
967          for (j=0;j<offsets[i]-1;j++)
968             ec_enc_bit_prob(enc, 1, 32768);
969          ec_enc_bit_prob(enc, 0, 32768);
970       }
971       offsets[i] *= (6<<BITRES);
972    }
973    alloc_trim = alloc_trim_analysis(st->mode, X, bandLogE, st->mode->nbEBands, LM, C, N);
974    ec_encode_bin(enc, trim_cdf[alloc_trim], trim_cdf[alloc_trim+1], 7);
975
976    /* Variable bitrate */
977    if (st->vbr_rate_norm>0)
978    {
979      celt_word16 alpha;
980      celt_int32 delta;
981      /* The target rate in 8th bits per frame */
982      celt_int32 vbr_rate;
983      celt_int32 target;
984      celt_int32 vbr_bound, max_allowed;
985
986      target = vbr_rate = M*st->vbr_rate_norm;
987
988      /* Shortblocks get a large boost in bitrate, but since they
989         are uncommon long blocks are not greatly affected */
990      if (shortBlocks || tf_sum < -2*(st->end-st->start))
991         target*=2;
992      else if (tf_sum < -(st->end-st->start))
993         target = 3*target/2;
994      else if (M > 1)
995         target-=(target+14)/28;
996
997      /* The current offset is removed from the target and the space used
998         so far is added*/
999      target=target+st->vbr_offset+ec_enc_tell(enc, BITRES);
1000
1001      /* Computes the max bit-rate allowed in VBR more to avoid violating the target rate and buffering */
1002      vbr_bound = vbr_rate;
1003      max_allowed = IMIN(vbr_rate+vbr_bound-st->vbr_reservoir>>(BITRES+3),nbAvailableBytes);
1004
1005      /* In VBR mode the frame size must not be reduced so much that it would result in the encoder running out of bits */
1006      nbAvailableBytes = target+(1<<(BITRES+2))>>(BITRES+3);
1007      nbAvailableBytes=IMAX(16,IMIN(max_allowed,nbAvailableBytes));
1008      target=nbAvailableBytes<<(BITRES+3);
1009
1010      if (st->vbr_count < 970)
1011      {
1012         st->vbr_count++;
1013         alpha = celt_rcp(SHL32(EXTEND32(st->vbr_count+20),16));
1014      } else
1015         alpha = QCONST16(.001f,15);
1016      /* By how much did we "miss" the target on that frame */
1017      delta = (celt_int32)target - vbr_rate;
1018      /* How many bits have we used in excess of what we're allowed */
1019      st->vbr_reservoir += delta;
1020      /*printf ("%d\n", st->vbr_reservoir);*/
1021
1022      /* Compute the offset we need to apply in order to reach the target */
1023      st->vbr_drift += (celt_int32)MULT16_32_Q15(alpha,delta-st->vbr_offset-st->vbr_drift);
1024      st->vbr_offset = -st->vbr_drift;
1025      /*printf ("%d\n", st->vbr_drift);*/
1026
1027      /* We could use any multiple of vbr_rate as bound (depending on the delay) */
1028      if (st->vbr_reservoir < 0)
1029      {
1030         /* We're under the min value -- increase rate */
1031         int adjust = (-st->vbr_reservoir)/(8<<BITRES);
1032         nbAvailableBytes += adjust;
1033         st->vbr_reservoir = 0;
1034         /*printf ("+%d\n", adjust);*/
1035      }
1036      nbCompressedBytes = IMIN(nbCompressedBytes,nbAvailableBytes+nbFilledBytes);
1037
1038      /* This moves the raw bits to take into account the new compressed size */
1039      ec_byte_shrink(&buf, nbCompressedBytes);
1040    }
1041
1042    if (C==2)
1043    {
1044       dual_stereo = stereo_analysis(st->mode, X, st->mode->nbEBands, LM, C, N);
1045       ec_enc_bit_prob(enc, dual_stereo, 32768);
1046    }
1047    if (C==2)
1048    {
1049       int effectiveRate;
1050
1051       /* Account for coarse energy */
1052       effectiveRate = (8*effectiveBytes - 80)>>LM;
1053
1054       /* effectiveRate in kb/s */
1055       effectiveRate = 2*effectiveRate/5;
1056       if (effectiveRate<35)
1057          intensity = 6;
1058       else if (effectiveRate<50)
1059          intensity = 12;
1060       else if (effectiveRate<68)
1061          intensity = 16;
1062       else if (effectiveRate<84)
1063          intensity = 18;
1064       else if (effectiveRate<102)
1065          intensity = 19;
1066       else if (effectiveRate<130)
1067          intensity = 20;
1068       else
1069          intensity = 100;
1070       intensity = IMIN(st->end,IMAX(st->start, intensity));
1071       ec_enc_uint(enc, intensity, 1+st->end-st->start);
1072    }
1073
1074    /* Bit allocation */
1075    ALLOC(fine_quant, st->mode->nbEBands, int);
1076    ALLOC(pulses, st->mode->nbEBands, int);
1077    ALLOC(fine_priority, st->mode->nbEBands, int);
1078
1079    bits = nbCompressedBytes*8 - ec_enc_tell(enc, 0) - 1;
1080    codedBands = compute_allocation(st->mode, st->start, st->end, offsets, alloc_trim, bits, pulses, fine_quant, fine_priority, C, LM);
1081
1082    quant_fine_energy(st->mode, st->start, st->end, bandE, oldBandE, error, fine_quant, enc, C);
1083
1084 #ifdef MEASURE_NORM_MSE
1085    float X0[3000];
1086    float bandE0[60];
1087    c=0; do 
1088       for (i=0;i<N;i++)
1089          X0[i+c*N] = X[i+c*N];
1090    while (++c<C);
1091    for (i=0;i<C*st->mode->nbEBands;i++)
1092       bandE0[i] = bandE[i];
1093 #endif
1094
1095    /* Residual quantisation */
1096    quant_all_bands(1, st->mode, st->start, st->end, X, C==2 ? X+N : NULL,
1097          bandE, pulses, shortBlocks, has_fold, dual_stereo, intensity, tf_res, resynth,
1098          nbCompressedBytes*8, enc, LM, codedBands);
1099
1100    quant_energy_finalise(st->mode, st->start, st->end, bandE, oldBandE, error, fine_quant, fine_priority, nbCompressedBytes*8-ec_enc_tell(enc, 0), enc, C);
1101
1102 #ifdef RESYNTH
1103    /* Re-synthesis of the coded audio if required */
1104    if (resynth)
1105    {
1106       celt_sig *out_mem[2];
1107       celt_sig *overlap_mem[2];
1108
1109       log2Amp(st->mode, st->start, st->end, bandE, oldBandE, C);
1110
1111 #ifdef MEASURE_NORM_MSE
1112       measure_norm_mse(st->mode, X, X0, bandE, bandE0, M, N, C);
1113 #endif
1114
1115       /* Synthesis */
1116       denormalise_bands(st->mode, X, freq, bandE, effEnd, C, M);
1117
1118       CELT_MOVE(st->syn_mem[0], st->syn_mem[0]+N, MAX_PERIOD);
1119       if (C==2)
1120          CELT_MOVE(st->syn_mem[1], st->syn_mem[1]+N, MAX_PERIOD);
1121
1122       c=0; do
1123          for (i=0;i<M*st->mode->eBands[st->start];i++)
1124             freq[c*N+i] = 0;
1125       while (++c<C);
1126       c=0; do
1127          for (i=M*st->mode->eBands[st->end];i<N;i++)
1128             freq[c*N+i] = 0;
1129       while (++c<C);
1130
1131       out_mem[0] = st->syn_mem[0]+MAX_PERIOD;
1132       if (C==2)
1133          out_mem[1] = st->syn_mem[1]+MAX_PERIOD;
1134
1135       c=0; do
1136          overlap_mem[c] = _overlap_mem + c*st->overlap;
1137       while (++c<C);
1138
1139       compute_inv_mdcts(st->mode, shortBlocks, freq, out_mem, overlap_mem, C, LM);
1140
1141 #ifdef ENABLE_POSTFILTER
1142       c=0; do {
1143          comb_filter(out_mem[c], out_mem[c], st->prefilter_period, st->prefilter_period, st->overlap, C,
1144                st->prefilter_gain, st->prefilter_gain, NULL, 0);
1145          comb_filter(out_mem[c]+st->overlap, out_mem[c]+st->overlap, st->prefilter_period, pitch_index, N-st->overlap, C,
1146                st->prefilter_gain, gain1, st->mode->window, st->mode->overlap);
1147       } while (++c<C);
1148 #endif /* ENABLE_POSTFILTER */
1149
1150       deemphasis(out_mem, (celt_word16*)pcm, N, C, st->mode->preemph, st->preemph_memD);
1151    }
1152 #endif
1153
1154    st->prefilter_period = pitch_index;
1155    st->prefilter_gain = gain1;
1156
1157    /* If there's any room left (can only happen for very high rates),
1158       fill it with zeros */
1159    while (ec_enc_tell(enc,0) + 8 <= nbCompressedBytes*8)
1160       ec_enc_bits(enc, 0, 8);
1161    ec_enc_done(enc);
1162    
1163    RESTORE_STACK;
1164    if (ec_enc_get_error(enc))
1165       return CELT_CORRUPTED_DATA;
1166    else
1167       return nbCompressedBytes;
1168 }
1169
1170 #ifdef FIXED_POINT
1171 #ifndef DISABLE_FLOAT_API
1172 int celt_encode_with_ec_float(CELTEncoder * restrict st, const float * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes, ec_enc *enc)
1173 {
1174    int j, ret, C, N, LM, M;
1175    VARDECL(celt_int16, in);
1176    SAVE_STACK;
1177
1178    if (pcm==NULL)
1179       return CELT_BAD_ARG;
1180
1181    for (LM=0;LM<4;LM++)
1182       if (st->mode->shortMdctSize<<LM==frame_size)
1183          break;
1184    if (LM>=MAX_CONFIG_SIZES)
1185       return CELT_BAD_ARG;
1186    M=1<<LM;
1187
1188    C = CHANNELS(st->channels);
1189    N = M*st->mode->shortMdctSize;
1190    ALLOC(in, C*N, celt_int16);
1191
1192    for (j=0;j<C*N;j++)
1193      in[j] = FLOAT2INT16(pcm[j]);
1194
1195    ret=celt_encode_with_ec(st,in,frame_size,compressed,nbCompressedBytes, enc);
1196 #ifdef RESYNTH
1197    for (j=0;j<C*N;j++)
1198       ((float*)pcm)[j]=in[j]*(1.f/32768.f);
1199 #endif
1200    RESTORE_STACK;
1201    return ret;
1202
1203 }
1204 #endif /*DISABLE_FLOAT_API*/
1205 #else
1206 int celt_encode_with_ec(CELTEncoder * restrict st, const celt_int16 * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes, ec_enc *enc)
1207 {
1208    int j, ret, C, N, LM, M;
1209    VARDECL(celt_sig, in);
1210    SAVE_STACK;
1211
1212    if (pcm==NULL)
1213       return CELT_BAD_ARG;
1214
1215    for (LM=0;LM<4;LM++)
1216       if (st->mode->shortMdctSize<<LM==frame_size)
1217          break;
1218    if (LM>=MAX_CONFIG_SIZES)
1219       return CELT_BAD_ARG;
1220    M=1<<LM;
1221
1222    C=CHANNELS(st->channels);
1223    N=M*st->mode->shortMdctSize;
1224    ALLOC(in, C*N, celt_sig);
1225    for (j=0;j<C*N;j++) {
1226      in[j] = SCALEOUT(pcm[j]);
1227    }
1228
1229    ret = celt_encode_with_ec_float(st,in,frame_size,compressed,nbCompressedBytes, enc);
1230 #ifdef RESYNTH
1231    for (j=0;j<C*N;j++)
1232       ((celt_int16*)pcm)[j] = FLOAT2INT16(in[j]);
1233 #endif
1234    RESTORE_STACK;
1235    return ret;
1236 }
1237 #endif
1238
1239 int celt_encode(CELTEncoder * restrict st, const celt_int16 * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes)
1240 {
1241    return celt_encode_with_ec(st, pcm, frame_size, compressed, nbCompressedBytes, NULL);
1242 }
1243
1244 #ifndef DISABLE_FLOAT_API
1245 int celt_encode_float(CELTEncoder * restrict st, const float * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes)
1246 {
1247    return celt_encode_with_ec_float(st, pcm, frame_size, compressed, nbCompressedBytes, NULL);
1248 }
1249 #endif /* DISABLE_FLOAT_API */
1250
1251 int celt_encoder_ctl(CELTEncoder * restrict st, int request, ...)
1252 {
1253    va_list ap;
1254    
1255    va_start(ap, request);
1256    switch (request)
1257    {
1258       case CELT_GET_MODE_REQUEST:
1259       {
1260          const CELTMode ** value = va_arg(ap, const CELTMode**);
1261          if (value==0)
1262             goto bad_arg;
1263          *value=st->mode;
1264       }
1265       break;
1266       case CELT_SET_COMPLEXITY_REQUEST:
1267       {
1268          int value = va_arg(ap, celt_int32);
1269          if (value<0 || value>10)
1270             goto bad_arg;
1271          st->complexity = value;
1272       }
1273       break;
1274       case CELT_SET_START_BAND_REQUEST:
1275       {
1276          celt_int32 value = va_arg(ap, celt_int32);
1277          if (value<0 || value>=st->mode->nbEBands)
1278             goto bad_arg;
1279          st->start = value;
1280       }
1281       break;
1282       case CELT_SET_END_BAND_REQUEST:
1283       {
1284          celt_int32 value = va_arg(ap, celt_int32);
1285          if (value<0 || value>=st->mode->nbEBands)
1286             goto bad_arg;
1287          st->end = value;
1288       }
1289       break;
1290       case CELT_SET_PREDICTION_REQUEST:
1291       {
1292          int value = va_arg(ap, celt_int32);
1293          if (value<0 || value>2)
1294             goto bad_arg;
1295          if (value==0)
1296          {
1297             st->force_intra   = 1;
1298          } else if (value==1) {
1299             st->force_intra   = 0;
1300          } else {
1301             st->force_intra   = 0;
1302          }   
1303       }
1304       break;
1305       case CELT_SET_VBR_RATE_REQUEST:
1306       {
1307          celt_int32 value = va_arg(ap, celt_int32);
1308          int frame_rate;
1309          int N = st->mode->shortMdctSize;
1310          if (value<0)
1311             goto bad_arg;
1312          if (value>3072000)
1313             value = 3072000;
1314          frame_rate = ((st->mode->Fs<<3)+(N>>1))/N;
1315          st->vbr_rate_norm = ((value<<(BITRES+3))+(frame_rate>>1))/frame_rate;
1316       }
1317       break;
1318       case CELT_RESET_STATE:
1319       {
1320          CELT_MEMSET((char*)&st->ENCODER_RESET_START, 0,
1321                celt_encoder_get_size(st->mode, st->channels)-
1322                ((char*)&st->ENCODER_RESET_START - (char*)st));
1323          st->vbr_offset = -(64<<BITRES);
1324          st->delayedIntra = 1;
1325          st->fold_decision = 1;
1326          st->tonal_average = QCONST16(1.f,8);
1327       }
1328       break;
1329       default:
1330          goto bad_request;
1331    }
1332    va_end(ap);
1333    return CELT_OK;
1334 bad_arg:
1335    va_end(ap);
1336    return CELT_BAD_ARG;
1337 bad_request:
1338    va_end(ap);
1339    return CELT_UNIMPLEMENTED;
1340 }
1341
1342 /**********************************************************************/
1343 /*                                                                    */
1344 /*                             DECODER                                */
1345 /*                                                                    */
1346 /**********************************************************************/
1347 #define DECODE_BUFFER_SIZE 2048
1348
1349 /** Decoder state 
1350  @brief Decoder state
1351  */
1352 struct CELTDecoder {
1353    const CELTMode *mode;
1354    int overlap;
1355    int channels;
1356
1357    int start, end;
1358
1359    /* Everything beyond this point gets cleared on a reset */
1360 #define DECODER_RESET_START last_pitch_index
1361
1362    int last_pitch_index;
1363    int loss_count;
1364    int postfilter_period;
1365    celt_word16 postfilter_gain;
1366
1367    celt_sig preemph_memD[2];
1368    
1369    celt_sig _decode_mem[1]; /* Size = channels*(DECODE_BUFFER_SIZE+mode->overlap) */
1370    /* celt_word16 lpc[],  Size = channels*LPC_ORDER */
1371    /* celt_word16 oldEBands[], Size = channels*mode->nbEBands */
1372 };
1373
1374 int celt_decoder_get_size(const CELTMode *mode, int channels)
1375 {
1376    int size = sizeof(struct CELTDecoder)
1377             + (channels*(DECODE_BUFFER_SIZE+mode->overlap)-1)*sizeof(celt_sig)
1378             + channels*LPC_ORDER*sizeof(celt_word16)
1379             + channels*mode->nbEBands*sizeof(celt_word16);
1380    return size;
1381 }
1382
1383 CELTDecoder *celt_decoder_create(const CELTMode *mode, int channels, int *error)
1384 {
1385    return celt_decoder_init(
1386          (CELTDecoder *)celt_alloc(celt_decoder_get_size(mode, channels)),
1387          mode, channels, error);
1388 }
1389
1390 CELTDecoder *celt_decoder_init(CELTDecoder *st, const CELTMode *mode, int channels, int *error)
1391 {
1392    if (channels < 0 || channels > 2)
1393    {
1394       if (error)
1395          *error = CELT_BAD_ARG;
1396       return NULL;
1397    }
1398
1399    if (st==NULL)
1400    {
1401       if (error)
1402          *error = CELT_ALLOC_FAIL;
1403       return NULL;
1404    }
1405
1406    CELT_MEMSET((char*)st, 0, celt_decoder_get_size(mode, channels));
1407
1408    st->mode = mode;
1409    st->overlap = mode->overlap;
1410    st->channels = channels;
1411
1412    st->start = 0;
1413    st->end = st->mode->effEBands;
1414
1415    st->loss_count = 0;
1416
1417    if (error)
1418       *error = CELT_OK;
1419    return st;
1420 }
1421
1422 void celt_decoder_destroy(CELTDecoder *st)
1423 {
1424    celt_free(st);
1425 }
1426
1427 static void celt_decode_lost(CELTDecoder * restrict st, celt_word16 * restrict pcm, int N, int LM)
1428 {
1429    int c;
1430    int pitch_index;
1431    int overlap = st->mode->overlap;
1432    celt_word16 fade = Q15ONE;
1433    int i, len;
1434    const int C = CHANNELS(st->channels);
1435    int offset;
1436    celt_sig *out_mem[2];
1437    celt_sig *decode_mem[2];
1438    celt_sig *overlap_mem[2];
1439    celt_word16 *lpc;
1440    SAVE_STACK;
1441    
1442    c=0; do {
1443       decode_mem[c] = st->_decode_mem + c*(DECODE_BUFFER_SIZE+st->overlap);
1444       out_mem[c] = decode_mem[c]+DECODE_BUFFER_SIZE-MAX_PERIOD;
1445       overlap_mem[c] = decode_mem[c]+DECODE_BUFFER_SIZE;
1446    } while (++c<C);
1447    lpc = (celt_word16*)(st->_decode_mem+(DECODE_BUFFER_SIZE+st->overlap)*C);
1448
1449    len = N+st->mode->overlap;
1450    
1451    if (st->loss_count == 0)
1452    {
1453       celt_word16 pitch_buf[MAX_PERIOD>>1];
1454       celt_word32 tmp=0;
1455       celt_word32 mem0[2]={0,0};
1456       celt_word16 mem1[2]={0,0};
1457       int len2 = len;
1458       /* FIXME: This is a kludge */
1459       if (len2>MAX_PERIOD>>1)
1460          len2 = MAX_PERIOD>>1;
1461       pitch_downsample(out_mem, pitch_buf, MAX_PERIOD, MAX_PERIOD,
1462                        C, mem0, mem1);
1463       pitch_search(st->mode, pitch_buf+((MAX_PERIOD-len2)>>1), pitch_buf, len2,
1464                    MAX_PERIOD-len2-100, &pitch_index, &tmp, 1<<LM);
1465       pitch_index = MAX_PERIOD-len2-pitch_index;
1466       st->last_pitch_index = pitch_index;
1467    } else {
1468       pitch_index = st->last_pitch_index;
1469       if (st->loss_count < 5)
1470          fade = QCONST16(.8f,15);
1471       else
1472          fade = 0;
1473    }
1474
1475    c=0; do {
1476       /* FIXME: This is more memory than necessary */
1477       celt_word32 e[2*MAX_PERIOD];
1478       celt_word16 exc[2*MAX_PERIOD];
1479       celt_word32 ac[LPC_ORDER+1];
1480       celt_word16 decay = 1;
1481       celt_word32 S1=0;
1482       celt_word16 mem[LPC_ORDER]={0};
1483
1484       offset = MAX_PERIOD-pitch_index;
1485       for (i=0;i<MAX_PERIOD;i++)
1486          exc[i] = ROUND16(out_mem[c][i], SIG_SHIFT);
1487
1488       if (st->loss_count == 0)
1489       {
1490          _celt_autocorr(exc, ac, st->mode->window, st->mode->overlap,
1491                         LPC_ORDER, MAX_PERIOD);
1492
1493          /* Noise floor -40 dB */
1494 #ifdef FIXED_POINT
1495          ac[0] += SHR32(ac[0],13);
1496 #else
1497          ac[0] *= 1.0001f;
1498 #endif
1499          /* Lag windowing */
1500          for (i=1;i<=LPC_ORDER;i++)
1501          {
1502             /*ac[i] *= exp(-.5*(2*M_PI*.002*i)*(2*M_PI*.002*i));*/
1503 #ifdef FIXED_POINT
1504             ac[i] -= MULT16_32_Q15(2*i*i, ac[i]);
1505 #else
1506             ac[i] -= ac[i]*(.008f*i)*(.008f*i);
1507 #endif
1508          }
1509
1510          _celt_lpc(lpc+c*LPC_ORDER, ac, LPC_ORDER);
1511       }
1512       for (i=0;i<LPC_ORDER;i++)
1513          mem[i] = ROUND16(out_mem[c][MAX_PERIOD-1-i], SIG_SHIFT);
1514       fir(exc, lpc+c*LPC_ORDER, exc, MAX_PERIOD, LPC_ORDER, mem);
1515       /*for (i=0;i<MAX_PERIOD;i++)printf("%d ", exc[i]); printf("\n");*/
1516       /* Check if the waveform is decaying (and if so how fast) */
1517       {
1518          celt_word32 E1=1, E2=1;
1519          int period;
1520          if (pitch_index <= MAX_PERIOD/2)
1521             period = pitch_index;
1522          else
1523             period = MAX_PERIOD/2;
1524          for (i=0;i<period;i++)
1525          {
1526             E1 += SHR32(MULT16_16(exc[MAX_PERIOD-period+i],exc[MAX_PERIOD-period+i]),8);
1527             E2 += SHR32(MULT16_16(exc[MAX_PERIOD-2*period+i],exc[MAX_PERIOD-2*period+i]),8);
1528          }
1529          if (E1 > E2)
1530             E1 = E2;
1531          decay = celt_sqrt(frac_div32(SHR(E1,1),E2));
1532       }
1533
1534       /* Copy excitation, taking decay into account */
1535       for (i=0;i<len+st->mode->overlap;i++)
1536       {
1537          celt_word16 tmp;
1538          if (offset+i >= MAX_PERIOD)
1539          {
1540             offset -= pitch_index;
1541             decay = MULT16_16_Q15(decay, decay);
1542          }
1543          e[i] = SHL32(EXTEND32(MULT16_16_Q15(decay, exc[offset+i])), SIG_SHIFT);
1544          tmp = ROUND16(out_mem[c][offset+i],SIG_SHIFT);
1545          S1 += SHR32(MULT16_16(tmp,tmp),8);
1546       }
1547       for (i=0;i<LPC_ORDER;i++)
1548          mem[i] = ROUND16(out_mem[c][MAX_PERIOD-1-i], SIG_SHIFT);
1549       for (i=0;i<len+st->mode->overlap;i++)
1550          e[i] = MULT16_32_Q15(fade, e[i]);
1551       iir(e, lpc+c*LPC_ORDER, e, len+st->mode->overlap, LPC_ORDER, mem);
1552
1553       {
1554          celt_word32 S2=0;
1555          for (i=0;i<len+overlap;i++)
1556          {
1557             celt_word16 tmp = ROUND16(e[i],SIG_SHIFT);
1558             S2 += SHR32(MULT16_16(tmp,tmp),8);
1559          }
1560          /* This checks for an "explosion" in the synthesis */
1561 #ifdef FIXED_POINT
1562          if (!(S1 > SHR32(S2,2)))
1563 #else
1564          /* Float test is written this way to catch NaNs at the same time */
1565          if (!(S1 > 0.2f*S2))
1566 #endif
1567          {
1568             for (i=0;i<len+overlap;i++)
1569                e[i] = 0;
1570          } else if (S1 < S2)
1571          {
1572             celt_word16 ratio = celt_sqrt(frac_div32(SHR32(S1,1)+1,S2+1));
1573             for (i=0;i<len+overlap;i++)
1574                e[i] = MULT16_32_Q15(ratio, e[i]);
1575          }
1576       }
1577
1578 #ifdef ENABLE_POSTFILTER
1579       /* Apply post-filter to the MDCT overlap of the previous frame */
1580       comb_filter(out_mem[c]+MAX_PERIOD, out_mem[c]+MAX_PERIOD, st->postfilter_period, st->postfilter_period, st->overlap, C,
1581                   st->postfilter_gain, st->postfilter_gain, NULL, 0);
1582 #endif /* ENABLE_POSTFILTER */
1583
1584       for (i=0;i<MAX_PERIOD+st->mode->overlap-N;i++)
1585          out_mem[c][i] = out_mem[c][N+i];
1586
1587       /* Apply TDAC to the concealed audio so that it blends with the
1588          previous and next frames */
1589       for (i=0;i<overlap/2;i++)
1590       {
1591          celt_word32 tmp;
1592          tmp = MULT16_32_Q15(st->mode->window[i],           e[N+overlap-1-i]) +
1593                MULT16_32_Q15(st->mode->window[overlap-i-1], e[N+i          ]);
1594          out_mem[c][MAX_PERIOD+i] = MULT16_32_Q15(st->mode->window[overlap-i-1], tmp);
1595          out_mem[c][MAX_PERIOD+overlap-i-1] = MULT16_32_Q15(st->mode->window[i], tmp);
1596       }
1597       for (i=0;i<N;i++)
1598          out_mem[c][MAX_PERIOD-N+i] = e[i];
1599
1600 #ifdef ENABLE_POSTFILTER
1601       /* Apply pre-filter to the MDCT overlap for the next frame (post-filter will be applied then) */
1602       comb_filter(e, out_mem[c]+MAX_PERIOD, st->postfilter_period, st->postfilter_period, st->overlap, C,
1603                   -st->postfilter_gain, -st->postfilter_gain, NULL, 0);
1604 #endif /* ENABLE_POSTFILTER */
1605       for (i=0;i<overlap;i++)
1606          out_mem[c][MAX_PERIOD+i] = e[i];
1607    } while (++c<C);
1608
1609    {
1610       celt_word32 *out_syn[2];
1611       out_syn[0] = out_mem[0]+MAX_PERIOD-N;
1612       if (C==2)
1613          out_syn[1] = out_mem[1]+MAX_PERIOD-N;
1614       deemphasis(out_syn, pcm, N, C, st->mode->preemph, st->preemph_memD);
1615    }
1616    
1617    st->loss_count++;
1618
1619    RESTORE_STACK;
1620 }
1621
1622 #ifdef FIXED_POINT
1623 int celt_decode_with_ec(CELTDecoder * restrict st, const unsigned char *data, int len, celt_int16 * restrict pcm, int frame_size, ec_dec *dec)
1624 {
1625 #else
1626 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)
1627 {
1628 #endif
1629    int c, i, N;
1630    int has_fold;
1631    int bits;
1632    ec_dec _dec;
1633    ec_byte_buffer buf;
1634    VARDECL(celt_sig, freq);
1635    VARDECL(celt_norm, X);
1636    VARDECL(celt_ener, bandE);
1637    VARDECL(int, fine_quant);
1638    VARDECL(int, pulses);
1639    VARDECL(int, offsets);
1640    VARDECL(int, fine_priority);
1641    VARDECL(int, tf_res);
1642    celt_sig *out_mem[2];
1643    celt_sig *decode_mem[2];
1644    celt_sig *overlap_mem[2];
1645    celt_sig *out_syn[2];
1646    celt_word16 *lpc;
1647    celt_word16 *oldBandE;
1648
1649    int shortBlocks;
1650    int isTransient;
1651    int intra_ener;
1652    const int C = CHANNELS(st->channels);
1653    int LM, M;
1654    int nbFilledBytes, nbAvailableBytes;
1655    int effEnd;
1656    int codedBands;
1657    int alloc_trim;
1658    int postfilter_pitch;
1659    celt_word16 postfilter_gain;
1660    int intensity=0;
1661    int dual_stereo=0;
1662    SAVE_STACK;
1663
1664    if (pcm==NULL)
1665       return CELT_BAD_ARG;
1666
1667    for (LM=0;LM<4;LM++)
1668       if (st->mode->shortMdctSize<<LM==frame_size)
1669          break;
1670    if (LM>=MAX_CONFIG_SIZES)
1671       return CELT_BAD_ARG;
1672    M=1<<LM;
1673
1674    c=0; do {
1675       decode_mem[c] = st->_decode_mem + c*(DECODE_BUFFER_SIZE+st->overlap);
1676       out_mem[c] = decode_mem[c]+DECODE_BUFFER_SIZE-MAX_PERIOD;
1677       overlap_mem[c] = decode_mem[c]+DECODE_BUFFER_SIZE;
1678    } while (++c<C);
1679    lpc = (celt_word16*)(st->_decode_mem+(DECODE_BUFFER_SIZE+st->overlap)*C);
1680    oldBandE = lpc+C*LPC_ORDER;
1681
1682    N = M*st->mode->shortMdctSize;
1683
1684    effEnd = st->end;
1685    if (effEnd > st->mode->effEBands)
1686       effEnd = st->mode->effEBands;
1687
1688    ALLOC(freq, C*N, celt_sig); /**< Interleaved signal MDCTs */
1689    ALLOC(X, C*N, celt_norm);   /**< Interleaved normalised MDCTs */
1690    ALLOC(bandE, st->mode->nbEBands*C, celt_ener);
1691    c=0; do
1692       for (i=0;i<M*st->mode->eBands[st->start];i++)
1693          X[c*N+i] = 0;
1694    while (++c<C);
1695    c=0; do   
1696       for (i=M*st->mode->eBands[effEnd];i<N;i++)
1697          X[c*N+i] = 0;
1698    while (++c<C);
1699
1700    if (data == NULL)
1701    {
1702       celt_decode_lost(st, pcm, N, LM);
1703       RESTORE_STACK;
1704       return CELT_OK;
1705    }
1706    if (len<0) {
1707      RESTORE_STACK;
1708      return CELT_BAD_ARG;
1709    }
1710    
1711    if (dec == NULL)
1712    {
1713       ec_byte_readinit(&buf,(unsigned char*)data,len);
1714       ec_dec_init(&_dec,&buf);
1715       dec = &_dec;
1716       nbFilledBytes = 0;
1717    } else {
1718       nbFilledBytes = (ec_dec_tell(dec, 0)+4)>>3;
1719    }
1720    nbAvailableBytes = len-nbFilledBytes;
1721
1722    if (ec_dec_bit_prob(dec, 32768))
1723    {
1724 #ifdef ENABLE_POSTFILTER
1725       int qg, octave;
1726       octave = ec_dec_uint(dec, 6);
1727       postfilter_pitch = (16<<octave)+ec_dec_bits(dec, 4+octave);
1728       qg = ec_dec_bits(dec, 2);
1729       postfilter_gain = QCONST16(.125f,15)*(qg+2);
1730 #else /* ENABLE_POSTFILTER */
1731       RESTORE_STACK;
1732       return CELT_CORRUPTED_DATA;
1733 #endif /* ENABLE_POSTFILTER */
1734
1735    } else {
1736       postfilter_gain = 0;
1737       postfilter_pitch = 0;
1738    }
1739
1740    /* Decode the global flags (first symbols in the stream) */
1741    intra_ener = ec_dec_bit_prob(dec, 8192);
1742    /* Get band energies */
1743    unquant_coarse_energy(st->mode, st->start, st->end, bandE, oldBandE,
1744          intra_ener, dec, C, LM);
1745
1746    if (LM > 0)
1747       isTransient = ec_dec_bit_prob(dec, 8192);
1748    else
1749       isTransient = 0;
1750
1751    if (isTransient)
1752       shortBlocks = M;
1753    else
1754       shortBlocks = 0;
1755
1756    ALLOC(tf_res, st->mode->nbEBands, int);
1757    tf_decode(st->start, st->end, C, isTransient, tf_res, LM, dec);
1758
1759    has_fold = ec_dec_bit_prob(dec, 8192)<<1;
1760    has_fold |= ec_dec_bit_prob(dec, (has_fold>>1) ? 32768 : 49152);
1761
1762    ALLOC(pulses, st->mode->nbEBands, int);
1763    ALLOC(offsets, st->mode->nbEBands, int);
1764    ALLOC(fine_priority, st->mode->nbEBands, int);
1765
1766    for (i=0;i<st->mode->nbEBands;i++)
1767       offsets[i] = 0;
1768    for (i=0;i<st->mode->nbEBands;i++)
1769    {
1770       if (ec_dec_bit_prob(dec, 1024))
1771       {
1772          while (ec_dec_bit_prob(dec, 32768))
1773             offsets[i]++;
1774          offsets[i]++;
1775          offsets[i] *= (6<<BITRES);
1776       }
1777    }
1778
1779    ALLOC(fine_quant, st->mode->nbEBands, int);
1780    {
1781       int fl;
1782       alloc_trim = 0;
1783       fl = ec_decode_bin(dec, 7);
1784       while (trim_cdf[alloc_trim+1] <= fl)
1785          alloc_trim++;
1786       ec_dec_update(dec, trim_cdf[alloc_trim], trim_cdf[alloc_trim+1], 128);
1787    }
1788
1789    if (C==2)
1790    {
1791       dual_stereo = ec_dec_bit_prob(dec, 32768);
1792       intensity = ec_dec_uint(dec, 1+st->end-st->start);
1793    }
1794
1795    bits = len*8 - ec_dec_tell(dec, 0) - 1;
1796    codedBands = compute_allocation(st->mode, st->start, st->end, offsets, alloc_trim, bits, pulses, fine_quant, fine_priority, C, LM);
1797    
1798    unquant_fine_energy(st->mode, st->start, st->end, bandE, oldBandE, fine_quant, dec, C);
1799
1800    /* Decode fixed codebook */
1801    quant_all_bands(0, st->mode, st->start, st->end, X, C==2 ? X+N : NULL,
1802          NULL, pulses, shortBlocks, has_fold, dual_stereo, intensity, tf_res, 1,
1803          len*8, dec, LM, codedBands);
1804
1805    unquant_energy_finalise(st->mode, st->start, st->end, bandE, oldBandE,
1806          fine_quant, fine_priority, len*8-ec_dec_tell(dec, 0), dec, C);
1807
1808    log2Amp(st->mode, st->start, st->end, bandE, oldBandE, C);
1809
1810    /* Synthesis */
1811    denormalise_bands(st->mode, X, freq, bandE, effEnd, C, M);
1812
1813    CELT_MOVE(decode_mem[0], decode_mem[0]+N, DECODE_BUFFER_SIZE-N);
1814    if (C==2)
1815       CELT_MOVE(decode_mem[1], decode_mem[1]+N, DECODE_BUFFER_SIZE-N);
1816
1817    c=0; do
1818       for (i=0;i<M*st->mode->eBands[st->start];i++)
1819          freq[c*N+i] = 0;
1820    while (++c<C);
1821    c=0; do
1822       for (i=M*st->mode->eBands[effEnd];i<N;i++)
1823          freq[c*N+i] = 0;
1824    while (++c<C);
1825
1826    out_syn[0] = out_mem[0]+MAX_PERIOD-N;
1827    if (C==2)
1828       out_syn[1] = out_mem[1]+MAX_PERIOD-N;
1829
1830    /* Compute inverse MDCTs */
1831    compute_inv_mdcts(st->mode, shortBlocks, freq, out_syn, overlap_mem, C, LM);
1832
1833 #ifdef ENABLE_POSTFILTER
1834    c=0; do {
1835       comb_filter(out_syn[c], out_syn[c], st->postfilter_period, st->postfilter_period, st->overlap, C,
1836             st->postfilter_gain, st->postfilter_gain, NULL, 0);
1837       comb_filter(out_syn[c]+st->overlap, out_syn[c]+st->overlap, st->postfilter_period, postfilter_pitch, N-st->overlap, C,
1838             st->postfilter_gain, postfilter_gain, st->mode->window, st->mode->overlap);
1839    } while (++c<C);
1840    st->postfilter_period = postfilter_pitch;
1841    st->postfilter_gain = postfilter_gain;
1842 #endif /* ENABLE_POSTFILTER */
1843
1844    deemphasis(out_syn, pcm, N, C, st->mode->preemph, st->preemph_memD);
1845    st->loss_count = 0;
1846    RESTORE_STACK;
1847    if (ec_dec_get_error(dec))
1848       return CELT_CORRUPTED_DATA;
1849    else
1850       return CELT_OK;
1851 }
1852
1853 #ifdef FIXED_POINT
1854 #ifndef DISABLE_FLOAT_API
1855 int celt_decode_with_ec_float(CELTDecoder * restrict st, const unsigned char *data, int len, float * restrict pcm, int frame_size, ec_dec *dec)
1856 {
1857    int j, ret, C, N, LM, M;
1858    VARDECL(celt_int16, out);
1859    SAVE_STACK;
1860
1861    if (pcm==NULL)
1862       return CELT_BAD_ARG;
1863
1864    for (LM=0;LM<4;LM++)
1865       if (st->mode->shortMdctSize<<LM==frame_size)
1866          break;
1867    if (LM>=MAX_CONFIG_SIZES)
1868       return CELT_BAD_ARG;
1869    M=1<<LM;
1870
1871    C = CHANNELS(st->channels);
1872    N = M*st->mode->shortMdctSize;
1873    
1874    ALLOC(out, C*N, celt_int16);
1875    ret=celt_decode_with_ec(st, data, len, out, frame_size, dec);
1876    if (ret==0)
1877       for (j=0;j<C*N;j++)
1878          pcm[j]=out[j]*(1.f/32768.f);
1879      
1880    RESTORE_STACK;
1881    return ret;
1882 }
1883 #endif /*DISABLE_FLOAT_API*/
1884 #else
1885 int celt_decode_with_ec(CELTDecoder * restrict st, const unsigned char *data, int len, celt_int16 * restrict pcm, int frame_size, ec_dec *dec)
1886 {
1887    int j, ret, C, N, LM, M;
1888    VARDECL(celt_sig, out);
1889    SAVE_STACK;
1890
1891    if (pcm==NULL)
1892       return CELT_BAD_ARG;
1893
1894    for (LM=0;LM<4;LM++)
1895       if (st->mode->shortMdctSize<<LM==frame_size)
1896          break;
1897    if (LM>=MAX_CONFIG_SIZES)
1898       return CELT_BAD_ARG;
1899    M=1<<LM;
1900
1901    C = CHANNELS(st->channels);
1902    N = M*st->mode->shortMdctSize;
1903    ALLOC(out, C*N, celt_sig);
1904
1905    ret=celt_decode_with_ec_float(st, data, len, out, frame_size, dec);
1906
1907    if (ret==0)
1908       for (j=0;j<C*N;j++)
1909          pcm[j] = FLOAT2INT16 (out[j]);
1910    
1911    RESTORE_STACK;
1912    return ret;
1913 }
1914 #endif
1915
1916 int celt_decode(CELTDecoder * restrict st, const unsigned char *data, int len, celt_int16 * restrict pcm, int frame_size)
1917 {
1918    return celt_decode_with_ec(st, data, len, pcm, frame_size, NULL);
1919 }
1920
1921 #ifndef DISABLE_FLOAT_API
1922 int celt_decode_float(CELTDecoder * restrict st, const unsigned char *data, int len, float * restrict pcm, int frame_size)
1923 {
1924    return celt_decode_with_ec_float(st, data, len, pcm, frame_size, NULL);
1925 }
1926 #endif /* DISABLE_FLOAT_API */
1927
1928 int celt_decoder_ctl(CELTDecoder * restrict st, int request, ...)
1929 {
1930    va_list ap;
1931
1932    va_start(ap, request);
1933    switch (request)
1934    {
1935       case CELT_GET_MODE_REQUEST:
1936       {
1937          const CELTMode ** value = va_arg(ap, const CELTMode**);
1938          if (value==0)
1939             goto bad_arg;
1940          *value=st->mode;
1941       }
1942       break;
1943       case CELT_SET_START_BAND_REQUEST:
1944       {
1945          celt_int32 value = va_arg(ap, celt_int32);
1946          if (value<0 || value>=st->mode->nbEBands)
1947             goto bad_arg;
1948          st->start = value;
1949       }
1950       break;
1951       case CELT_SET_END_BAND_REQUEST:
1952       {
1953          celt_int32 value = va_arg(ap, celt_int32);
1954          if (value<0 || value>=st->mode->nbEBands)
1955             goto bad_arg;
1956          st->end = value;
1957       }
1958       break;
1959       case CELT_RESET_STATE:
1960       {
1961          CELT_MEMSET((char*)&st->DECODER_RESET_START, 0,
1962                celt_decoder_get_size(st->mode, st->channels)-
1963                ((char*)&st->DECODER_RESET_START - (char*)st));
1964       }
1965       break;
1966       default:
1967          goto bad_request;
1968    }
1969    va_end(ap);
1970    return CELT_OK;
1971 bad_arg:
1972    va_end(ap);
1973    return CELT_BAD_ARG;
1974 bad_request:
1975       va_end(ap);
1976   return CELT_UNIMPLEMENTED;
1977 }
1978
1979 const char *celt_strerror(int error)
1980 {
1981    static const char *error_strings[8] = {
1982       "success",
1983       "invalid argument",
1984       "invalid mode",
1985       "internal error",
1986       "corrupted stream",
1987       "request not implemented",
1988       "invalid state",
1989       "memory allocation failed"
1990    };
1991    if (error > 0 || error < -7)
1992       return "unknown error";
1993    else 
1994       return error_strings[-error];
1995 }
1996