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