2e206e77272f43bd9110f18bb6e5efb6c9adb5e7
[opus.git] / libcelt / celt.c
1 /* Copyright (c) 2007-2008 CSIRO
2    Copyright (c) 2007-2009 Xiph.Org Foundation
3    Copyright (c) 2008 Gregory Maxwell 
4    Written by Jean-Marc Valin and Gregory Maxwell */
5 /*
6    Redistribution and use in source and binary forms, with or without
7    modification, are permitted provided that the following conditions
8    are met:
9    
10    - Redistributions of source code must retain the above copyright
11    notice, this list of conditions and the following disclaimer.
12    
13    - Redistributions in binary form must reproduce the above copyright
14    notice, this list of conditions and the following disclaimer in the
15    documentation and/or other materials provided with the distribution.
16    
17    - Neither the name of the Xiph.org Foundation nor the names of its
18    contributors may be used to endorse or promote products derived from
19    this software without specific prior written permission.
20    
21    THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
22    ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
23    LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
24    A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR
25    CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
26    EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
27    PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
28    PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
29    LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
30    NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
31    SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
32 */
33
34 #ifdef HAVE_CONFIG_H
35 #include "config.h"
36 #endif
37
38 #define CELT_C
39
40 #include "os_support.h"
41 #include "mdct.h"
42 #include <math.h>
43 #include "celt.h"
44 #include "pitch.h"
45 #include "kiss_fftr.h"
46 #include "bands.h"
47 #include "modes.h"
48 #include "entcode.h"
49 #include "quant_bands.h"
50 #include "psy.h"
51 #include "rate.h"
52 #include "stack_alloc.h"
53 #include "mathops.h"
54 #include "float_cast.h"
55 #include <stdarg.h>
56
57 #define LPC_ORDER 24
58 #if !defined(FIXED_POINT) || defined(NEW_PLC)
59 #include "plc.c"
60 #endif
61
62 static const celt_word16 preemph = QCONST16(0.8f,15);
63
64 #ifdef FIXED_POINT
65 static const celt_word16 transientWindow[16] = {
66      279,  1106,  2454,  4276,  6510,  9081, 11900, 14872,
67    17896, 20868, 23687, 26258, 28492, 30314, 31662, 32489};
68 #else
69 static const float transientWindow[16] = {
70    0.0085135, 0.0337639, 0.0748914, 0.1304955, 
71    0.1986827, 0.2771308, 0.3631685, 0.4538658,
72    0.5461342, 0.6368315, 0.7228692, 0.8013173, 
73    0.8695045, 0.9251086, 0.9662361, 0.9914865};
74 #endif
75
76 #define ENCODERVALID   0x4c434554
77 #define ENCODERPARTIAL 0x5445434c
78 #define ENCODERFREED   0x4c004500
79    
80 /** Encoder state 
81  @brief Encoder state
82  */
83 struct CELTEncoder {
84    celt_uint32 marker;
85    const CELTMode *mode;     /**< Mode used by the encoder */
86    int frame_size;
87    int block_size;
88    int overlap;
89    int channels;
90    
91    int pitch_enabled;       /* Complexity level is allowed to use pitch */
92    int pitch_permitted;     /*  Use of the LTP is permitted by the user */
93    int pitch_available;     /*  Amount of pitch buffer available */
94    int force_intra;
95    int delayedIntra;
96    celt_word16 tonal_average;
97    int fold_decision;
98    celt_word16 gain_prod;
99
100    /* VBR-related parameters */
101    celt_int32 vbr_reservoir;
102    celt_int32 vbr_drift;
103    celt_int32 vbr_offset;
104    celt_int32 vbr_count;
105
106    celt_int32 vbr_rate; /* Target number of 16th bits per frame */
107    celt_word16 * restrict preemph_memE; 
108    celt_sig    * restrict preemph_memD;
109
110    celt_sig *in_mem;
111    celt_sig *out_mem;
112    celt_word16 *pitch_buf;
113    celt_sig xmem;
114
115    celt_word16 *oldBandE;
116 };
117
118 static int check_encoder(const CELTEncoder *st) 
119 {
120    if (st==NULL)
121    {
122       celt_warning("NULL passed as an encoder structure");  
123       return CELT_INVALID_STATE;
124    }
125    if (st->marker == ENCODERVALID)
126       return CELT_OK;
127    if (st->marker == ENCODERFREED)
128       celt_warning("Referencing an encoder that has already been freed");
129    else
130       celt_warning("This is not a valid CELT encoder structure");
131    return CELT_INVALID_STATE;
132 }
133
134 CELTEncoder *celt_encoder_create(const CELTMode *mode, int channels, int *error)
135 {
136    int N, C;
137    CELTEncoder *st;
138
139    if (check_mode(mode) != CELT_OK)
140    {
141       if (error)
142          *error = CELT_INVALID_MODE;
143       return NULL;
144    }
145 #ifdef DISABLE_STEREO
146    if (channels > 1)
147    {
148       celt_warning("Stereo support was disable from this build");
149       if (error)
150          *error = CELT_BAD_ARG;
151       return NULL;
152    }
153 #endif
154
155    if (channels < 0 || channels > 2)
156    {
157       celt_warning("Only mono and stereo supported");
158       if (error)
159          *error = CELT_BAD_ARG;
160       return NULL;
161    }
162
163    N = mode->mdctSize;
164    C = channels;
165    st = celt_alloc(sizeof(CELTEncoder));
166    
167    if (st==NULL)
168    {
169       if (error)
170          *error = CELT_ALLOC_FAIL;
171       return NULL;
172    }
173    st->marker = ENCODERPARTIAL;
174    st->mode = mode;
175    st->frame_size = N;
176    st->block_size = N;
177    st->overlap = mode->overlap;
178    st->channels = channels;
179
180    st->vbr_rate = 0;
181    st->pitch_enabled = 1;
182    st->pitch_permitted = 1;
183    st->pitch_available = 1;
184    st->force_intra  = 0;
185    st->delayedIntra = 1;
186    st->tonal_average = QCONST16(1.,8);
187    st->fold_decision = 1;
188
189    st->in_mem = celt_alloc(st->overlap*C*sizeof(celt_sig));
190    st->out_mem = celt_alloc((MAX_PERIOD+st->overlap)*C*sizeof(celt_sig));
191    st->pitch_buf = celt_alloc(((MAX_PERIOD>>1)+2)*sizeof(celt_word16));
192
193    st->oldBandE = (celt_word16*)celt_alloc(C*mode->nbEBands*sizeof(celt_word16));
194
195    st->preemph_memE = (celt_word16*)celt_alloc(C*sizeof(celt_word16));
196    st->preemph_memD = (celt_sig*)celt_alloc(C*sizeof(celt_sig));
197
198    if ((st->in_mem!=NULL) && (st->out_mem!=NULL) && (st->oldBandE!=NULL) 
199        && (st->preemph_memE!=NULL) && (st->preemph_memD!=NULL))
200    {
201       if (error)
202          *error = CELT_OK;
203       st->marker   = ENCODERVALID;
204       return st;
205    }
206    /* If the setup fails for some reason deallocate it. */
207    celt_encoder_destroy(st);  
208    if (error)
209       *error = CELT_ALLOC_FAIL;
210    return NULL;
211 }
212
213 void celt_encoder_destroy(CELTEncoder *st)
214 {
215    if (st == NULL)
216    {
217       celt_warning("NULL passed to celt_encoder_destroy");
218       return;
219    }
220
221    if (st->marker == ENCODERFREED)
222    {
223       celt_warning("Freeing an encoder which has already been freed"); 
224       return;
225    }
226
227    if (st->marker != ENCODERVALID && st->marker != ENCODERPARTIAL)
228    {
229       celt_warning("This is not a valid CELT encoder structure");
230       return;
231    }
232    /*Check_mode is non-fatal here because we can still free
233     the encoder memory even if the mode is bad, although calling
234     the free functions in this order is a violation of the API.*/
235    check_mode(st->mode);
236    
237    celt_free(st->in_mem);
238    celt_free(st->out_mem);
239    celt_free(st->pitch_buf);
240    celt_free(st->oldBandE);
241    
242    celt_free(st->preemph_memE);
243    celt_free(st->preemph_memD);
244
245    st->marker = ENCODERFREED;
246    
247    celt_free(st);
248 }
249
250 static inline celt_int16 FLOAT2INT16(float x)
251 {
252    x = x*CELT_SIG_SCALE;
253    x = MAX32(x, -32768);
254    x = MIN32(x, 32767);
255    return (celt_int16)float2int(x);
256 }
257
258 static inline celt_word16 SIG2WORD16(celt_sig x)
259 {
260 #ifdef FIXED_POINT
261    x = PSHR32(x, SIG_SHIFT);
262    x = MAX32(x, -32768);
263    x = MIN32(x, 32767);
264    return EXTRACT16(x);
265 #else
266    return (celt_word16)x;
267 #endif
268 }
269
270 static int transient_analysis(celt_word32 *in, int len, int C, int *transient_time, int *transient_shift)
271 {
272    int c, i, n;
273    celt_word32 ratio;
274    VARDECL(celt_word32, begin);
275    SAVE_STACK;
276    ALLOC(begin, len, celt_word32);
277    for (i=0;i<len;i++)
278       begin[i] = ABS32(SHR32(in[C*i],SIG_SHIFT));
279    for (c=1;c<C;c++)
280    {
281       for (i=0;i<len;i++)
282          begin[i] = MAX32(begin[i], ABS32(SHR32(in[C*i+c],SIG_SHIFT)));
283    }
284    for (i=1;i<len;i++)
285       begin[i] = MAX32(begin[i-1],begin[i]);
286    n = -1;
287    for (i=8;i<len-8;i++)
288    {
289       if (begin[i] < MULT16_32_Q15(QCONST16(.2f,15),begin[len-1]))
290          n=i;
291    }
292    if (n<32)
293    {
294       n = -1;
295       ratio = 0;
296    } else {
297       ratio = DIV32(begin[len-1],1+begin[n-16]);
298    }
299    if (ratio < 0)
300       ratio = 0;
301    if (ratio > 1000)
302       ratio = 1000;
303    ratio *= ratio;
304    
305    if (ratio > 2048)
306       *transient_shift = 3;
307    else
308       *transient_shift = 0;
309    
310    *transient_time = n;
311    
312    RESTORE_STACK;
313    return ratio > 20;
314 }
315
316 /** Apply window and compute the MDCT for all sub-frames and 
317     all channels in a frame */
318 static void compute_mdcts(const CELTMode *mode, int shortBlocks, celt_sig * restrict in, celt_sig * restrict out, int _C)
319 {
320    const int C = CHANNELS(_C);
321    if (C==1 && !shortBlocks)
322    {
323       const mdct_lookup *lookup = MDCT(mode);
324       const int overlap = OVERLAP(mode);
325       clt_mdct_forward(lookup, in, out, mode->window, overlap);
326    } else {
327       const mdct_lookup *lookup = MDCT(mode);
328       const int overlap = OVERLAP(mode);
329       int N = FRAMESIZE(mode);
330       int B = 1;
331       int b, c;
332       VARDECL(celt_word32, x);
333       VARDECL(celt_word32, tmp);
334       SAVE_STACK;
335       if (shortBlocks)
336       {
337          lookup = &mode->shortMdct;
338          N = mode->shortMdctSize;
339          B = mode->nbShortMdcts;
340       }
341       ALLOC(x, N+overlap, celt_word32);
342       ALLOC(tmp, N, celt_word32);
343       for (c=0;c<C;c++)
344       {
345          for (b=0;b<B;b++)
346          {
347             int j;
348             for (j=0;j<N+overlap;j++)
349                x[j] = in[C*(b*N+j)+c];
350             clt_mdct_forward(lookup, x, tmp, mode->window, overlap);
351             /* Interleaving the sub-frames */
352             for (j=0;j<N;j++)
353                out[(j*B+b)+c*N*B] = tmp[j];
354          }
355       }
356       RESTORE_STACK;
357    }
358 }
359
360 /** Compute the IMDCT and apply window for all sub-frames and 
361     all channels in a frame */
362 static void compute_inv_mdcts(const CELTMode *mode, int shortBlocks, celt_sig *X, int transient_time, int transient_shift, celt_sig * restrict out_mem, int _C)
363 {
364    int c, N4;
365    const int C = CHANNELS(_C);
366    const int N = FRAMESIZE(mode);
367    const int overlap = OVERLAP(mode);
368    N4 = (N-overlap)>>1;
369    for (c=0;c<C;c++)
370    {
371       int j;
372       if (transient_shift==0 && C==1 && !shortBlocks) {
373          const mdct_lookup *lookup = MDCT(mode);
374          clt_mdct_backward(lookup, X, out_mem+C*(MAX_PERIOD-N-N4), mode->window, overlap);
375       } else {
376          VARDECL(celt_word32, x);
377          VARDECL(celt_word32, tmp);
378          int b;
379          int N2 = N;
380          int B = 1;
381          int n4offset=0;
382          const mdct_lookup *lookup = MDCT(mode);
383          SAVE_STACK;
384          
385          ALLOC(x, 2*N, celt_word32);
386          ALLOC(tmp, N, celt_word32);
387
388          if (shortBlocks)
389          {
390             lookup = &mode->shortMdct;
391             N2 = mode->shortMdctSize;
392             B = mode->nbShortMdcts;
393             n4offset = N4;
394          }
395          /* Prevents problems from the imdct doing the overlap-add */
396          CELT_MEMSET(x+N4, 0, N2);
397
398          for (b=0;b<B;b++)
399          {
400             /* De-interleaving the sub-frames */
401             for (j=0;j<N2;j++)
402                tmp[j] = X[(j*B+b)+c*N2*B];
403             clt_mdct_backward(lookup, tmp, x+n4offset+N2*b, mode->window, overlap);
404          }
405
406          if (transient_shift > 0)
407          {
408 #ifdef FIXED_POINT
409             for (j=0;j<16;j++)
410                x[N4+transient_time+j-16] = MULT16_32_Q15(SHR16(Q15_ONE-transientWindow[j],transient_shift)+transientWindow[j], SHL32(x[N4+transient_time+j-16],transient_shift));
411             for (j=transient_time;j<N+overlap;j++)
412                x[N4+j] = SHL32(x[N4+j], transient_shift);
413 #else
414             for (j=0;j<16;j++)
415                x[N4+transient_time+j-16] *= 1+transientWindow[j]*((1<<transient_shift)-1);
416             for (j=transient_time;j<N+overlap;j++)
417                x[N4+j] *= 1<<transient_shift;
418 #endif
419          }
420          /* The first and last part would need to be set to zero 
421             if we actually wanted to use them. */
422          for (j=0;j<overlap;j++)
423             out_mem[C*(MAX_PERIOD-N)+C*j+c] += x[j+N4];
424          for (j=0;j<overlap;j++)
425             out_mem[C*(MAX_PERIOD)+C*(overlap-j-1)+c] = x[2*N-j-N4-1];
426          for (j=0;j<2*N4;j++)
427             out_mem[C*(MAX_PERIOD-N)+C*(j+overlap)+c] = x[j+N4+overlap];
428          RESTORE_STACK;
429       }
430    }
431 }
432
433 #define FLAG_NONE        0
434 #define FLAG_INTRA       (1U<<13)
435 #define FLAG_PITCH       (1U<<12)
436 #define FLAG_SHORT       (1U<<11)
437 #define FLAG_FOLD        (1U<<10)
438 #define FLAG_MASK        (FLAG_INTRA|FLAG_PITCH|FLAG_SHORT|FLAG_FOLD)
439
440 static const int flaglist[8] = {
441       0 /*00  */ | FLAG_FOLD,
442       1 /*01  */ | FLAG_PITCH|FLAG_FOLD,
443       8 /*1000*/ | FLAG_NONE,
444       9 /*1001*/ | FLAG_SHORT|FLAG_FOLD,
445      10 /*1010*/ | FLAG_PITCH,
446      11 /*1011*/ | FLAG_INTRA,
447       6 /*110 */ | FLAG_INTRA|FLAG_FOLD,
448       7 /*111 */ | FLAG_INTRA|FLAG_SHORT|FLAG_FOLD
449 };
450
451 static void encode_flags(ec_enc *enc, int intra_ener, int has_pitch, int shortBlocks, int has_fold)
452 {
453    int i;
454    int flags=FLAG_NONE;
455    int flag_bits;
456    flags |= intra_ener   ? FLAG_INTRA : 0;
457    flags |= has_pitch    ? FLAG_PITCH : 0;
458    flags |= shortBlocks  ? FLAG_SHORT : 0;
459    flags |= has_fold     ? FLAG_FOLD  : 0;
460    for (i=0;i<8;i++)
461       if (flags == (flaglist[i]&FLAG_MASK))
462          break;
463    celt_assert(i<8);
464    flag_bits = flaglist[i]&0xf;
465    /*printf ("enc %d: %d %d %d %d\n", flag_bits, intra_ener, has_pitch, shortBlocks, has_fold);*/
466    if (i<2)
467       ec_enc_uint(enc, flag_bits, 4);
468    else if (i<6)
469       ec_enc_uint(enc, flag_bits, 16);
470    else
471       ec_enc_uint(enc, flag_bits, 8);
472 }
473
474 static void decode_flags(ec_dec *dec, int *intra_ener, int *has_pitch, int *shortBlocks, int *has_fold)
475 {
476    int i;
477    int flag_bits;
478    flag_bits = ec_dec_uint(dec, 4);
479    /*printf ("(%d) ", flag_bits);*/
480    if (flag_bits==2)
481       flag_bits = (flag_bits<<2) | ec_dec_uint(dec, 4);
482    else if (flag_bits==3)
483       flag_bits = (flag_bits<<1) | ec_dec_uint(dec, 2);
484    for (i=0;i<8;i++)
485       if (flag_bits == (flaglist[i]&0xf))
486          break;
487    celt_assert(i<8);
488    *intra_ener  = (flaglist[i]&FLAG_INTRA) != 0;
489    *has_pitch   = (flaglist[i]&FLAG_PITCH) != 0;
490    *shortBlocks = (flaglist[i]&FLAG_SHORT) != 0;
491    *has_fold    = (flaglist[i]&FLAG_FOLD ) != 0;
492    /*printf ("dec %d: %d %d %d %d\n", flag_bits, *intra_ener, *has_pitch, *shortBlocks, *has_fold);*/
493 }
494
495 static void deemphasis(celt_sig *in, celt_word16 *pcm, int N, int _C, celt_word16 coef, celt_sig *mem)
496 {
497    const int C = CHANNELS(_C);
498    int c;
499    for (c=0;c<C;c++)
500    {
501       int j;
502       for (j=0;j<N;j++)
503       {
504          celt_sig tmp = MAC16_32_Q15(in[C*(MAX_PERIOD-N)+C*j+c],
505                                        coef,mem[c]);
506          mem[c] = tmp;
507          pcm[C*j+c] = SCALEOUT(SIG2WORD16(tmp));
508       }
509    }
510
511 }
512
513 static void mdct_shape(const CELTMode *mode, celt_norm *X, int start, int end, int N, int nbShortMdcts, int mdct_weight_shift, int _C)
514 {
515    int m, i, c;
516    const int C = CHANNELS(_C);
517    for (c=0;c<C;c++)
518       for (m=start;m<end;m++)
519          for (i=m+c*N;i<(c+1)*N;i+=nbShortMdcts)
520 #ifdef FIXED_POINT
521             X[i] = SHR16(X[i], mdct_weight_shift);
522 #else
523             X[i] = (1.f/(1<<mdct_weight_shift))*X[i];
524 #endif
525    renormalise_bands(mode, X, C);
526 }
527
528
529 #ifdef FIXED_POINT
530 int celt_encode(CELTEncoder * restrict st, const celt_int16 * pcm, celt_int16 * optional_synthesis, unsigned char *compressed, int nbCompressedBytes)
531 {
532 #else
533 int celt_encode_float(CELTEncoder * restrict st, const celt_sig * pcm, celt_sig * optional_synthesis, unsigned char *compressed, int nbCompressedBytes)
534 {
535 #endif
536    int i, c, N, NN, N4;
537    int has_pitch;
538    int pitch_index;
539    int bits;
540    int has_fold=1;
541    int coarse_needed;
542    ec_byte_buffer buf;
543    ec_enc         enc;
544    VARDECL(celt_sig, in);
545    VARDECL(celt_sig, freq);
546    VARDECL(celt_sig, pitch_freq);
547    VARDECL(celt_norm, X);
548    VARDECL(celt_ener, bandE);
549    VARDECL(celt_word16, bandLogE);
550    VARDECL(int, fine_quant);
551    VARDECL(celt_word16, error);
552    VARDECL(int, pulses);
553    VARDECL(int, offsets);
554    VARDECL(int, fine_priority);
555    int intra_ener = 0;
556    int shortBlocks=0;
557    int transient_time;
558    int transient_shift;
559    const int C = CHANNELS(st->channels);
560    int mdct_weight_shift = 0;
561    int mdct_weight_pos=0;
562    int gain_id=0;
563    int norm_rate;
564    SAVE_STACK;
565
566    if (check_encoder(st) != CELT_OK)
567       return CELT_INVALID_STATE;
568
569    if (check_mode(st->mode) != CELT_OK)
570       return CELT_INVALID_MODE;
571
572    if (nbCompressedBytes<0 || pcm==NULL)
573      return CELT_BAD_ARG; 
574
575    /* The memset is important for now in case the encoder doesn't 
576       fill up all the bytes */
577    CELT_MEMSET(compressed, 0, nbCompressedBytes);
578    ec_byte_writeinit_buffer(&buf, compressed, nbCompressedBytes);
579    ec_enc_init(&enc,&buf);
580
581    N = st->block_size;
582    N4 = (N-st->overlap)>>1;
583    ALLOC(in, 2*C*N-2*C*N4, celt_sig);
584
585    CELT_COPY(in, st->in_mem, C*st->overlap);
586    for (c=0;c<C;c++)
587    {
588       const celt_word16 * restrict pcmp = pcm+c;
589       celt_sig * restrict inp = in+C*st->overlap+c;
590       for (i=0;i<N;i++)
591       {
592          /* Apply pre-emphasis */
593          celt_sig tmp = SCALEIN(SHL32(EXTEND32(*pcmp), SIG_SHIFT));
594          *inp = SUB32(tmp, SHR32(MULT16_16(preemph,st->preemph_memE[c]),3));
595          st->preemph_memE[c] = SCALEIN(*pcmp);
596          inp += C;
597          pcmp += C;
598       }
599    }
600    CELT_COPY(st->in_mem, in+C*(2*N-2*N4-st->overlap), C*st->overlap);
601
602    /* Transient handling */
603    transient_time = -1;
604    transient_shift = 0;
605    shortBlocks = 0;
606
607    if (st->mode->nbShortMdcts > 1 && transient_analysis(in, N+st->overlap, C, &transient_time, &transient_shift))
608    {
609 #ifndef FIXED_POINT
610       float gain_1;
611 #endif
612       /* Apply the inverse shaping window */
613       if (transient_shift)
614       {
615 #ifdef FIXED_POINT
616          for (c=0;c<C;c++)
617             for (i=0;i<16;i++)
618                in[C*(transient_time+i-16)+c] = MULT16_32_Q15(EXTRACT16(SHR32(celt_rcp(Q15ONE+MULT16_16(transientWindow[i],((1<<transient_shift)-1))),1)), in[C*(transient_time+i-16)+c]);
619          for (c=0;c<C;c++)
620             for (i=transient_time;i<N+st->overlap;i++)
621                in[C*i+c] = SHR32(in[C*i+c], transient_shift);
622 #else
623          for (c=0;c<C;c++)
624             for (i=0;i<16;i++)
625                in[C*(transient_time+i-16)+c] /= 1+transientWindow[i]*((1<<transient_shift)-1);
626          gain_1 = 1./(1<<transient_shift);
627          for (c=0;c<C;c++)
628             for (i=transient_time;i<N+st->overlap;i++)
629                in[C*i+c] *= gain_1;
630 #endif
631       }
632       shortBlocks = 1;
633       has_fold = 1;
634    }
635
636    ALLOC(freq, C*N, celt_sig); /**< Interleaved signal MDCTs */
637    ALLOC(bandE,st->mode->nbEBands*C, celt_ener);
638    ALLOC(bandLogE,st->mode->nbEBands*C, celt_word16);
639    /* Compute MDCTs */
640    compute_mdcts(st->mode, shortBlocks, in, freq, C);
641
642
643    norm_rate = (nbCompressedBytes-5)*8*(celt_uint32)st->mode->Fs/(C*N)>>10;
644    /* Pitch analysis: we do it early to save on the peak stack space */
645    /* Don't use pitch if there isn't enough data available yet, 
646       or if we're using shortBlocks */
647    has_pitch = st->pitch_enabled && st->pitch_permitted && (N <= 512) 
648             && (st->pitch_available >= MAX_PERIOD) && (!shortBlocks)
649             && norm_rate < 50;
650    if (has_pitch)
651    {
652       VARDECL(celt_word16, x_lp);
653       SAVE_STACK;
654       ALLOC(x_lp, (2*N-2*N4)>>1, celt_word16);
655       pitch_downsample(in, x_lp, 2*N-2*N4, N, C, &st->xmem, &st->pitch_buf[MAX_PERIOD>>1]);
656       pitch_search(st->mode, x_lp, st->pitch_buf, 2*N-2*N4, MAX_PERIOD-(2*N-2*N4), &pitch_index, &st->xmem);
657       RESTORE_STACK;
658    }
659
660    /* Deferred allocation after find_spectral_pitch() to reduce 
661       the peak memory usage */
662    ALLOC(X, C*N, celt_norm);         /**< Interleaved normalised MDCTs */
663
664    ALLOC(pitch_freq, C*N, celt_sig); /**< Interleaved signal MDCTs */
665    if (has_pitch)
666    {
667       compute_mdcts(st->mode, 0, st->out_mem+pitch_index*C, pitch_freq, C);
668       has_pitch = compute_pitch_gain(st->mode, freq, pitch_freq, norm_rate, &gain_id, C, &st->gain_prod);
669    }
670    
671    if (has_pitch)
672       apply_pitch(st->mode, freq, pitch_freq, gain_id, 1, C);
673
674    compute_band_energies(st->mode, freq, bandE, C);
675    for (i=0;i<st->mode->nbEBands*C;i++)
676       bandLogE[i] = amp2Log(bandE[i]);
677
678    /* Band normalisation */
679    normalise_bands(st->mode, freq, X, bandE, C);
680    if (!shortBlocks && !folding_decision(st->mode, X, &st->tonal_average, &st->fold_decision, C))
681       has_fold = 0;
682
683    /* Don't use intra energy when we're operating at low bit-rate */
684    intra_ener = st->force_intra || (!has_pitch && st->delayedIntra && nbCompressedBytes > st->mode->nbEBands);
685    if (shortBlocks || intra_decision(bandLogE, st->oldBandE, st->mode->nbEBands))
686       st->delayedIntra = 1;
687    else
688       st->delayedIntra = 0;
689
690    NN = st->mode->eBands[st->mode->nbEBands];
691    if (shortBlocks && !transient_shift) 
692    {
693       celt_word32 sum[8]={1,1,1,1,1,1,1,1};
694       int m;
695       for (c=0;c<C;c++)
696       {
697          m=0;
698          do {
699             celt_word32 tmp=0;
700             for (i=m+c*N;i<c*N+NN;i+=st->mode->nbShortMdcts)
701                tmp += ABS32(X[i]);
702             sum[m++] += tmp;
703          } while (m<st->mode->nbShortMdcts);
704       }
705       m=0;
706 #ifdef FIXED_POINT
707       do {
708          if (SHR32(sum[m+1],3) > sum[m])
709          {
710             mdct_weight_shift=2;
711             mdct_weight_pos = m;
712          } else if (SHR32(sum[m+1],1) > sum[m] && mdct_weight_shift < 2)
713          {
714             mdct_weight_shift=1;
715             mdct_weight_pos = m;
716          }
717          m++;
718       } while (m<st->mode->nbShortMdcts-1);
719 #else
720       do {
721          if (sum[m+1] > 8*sum[m])
722          {
723             mdct_weight_shift=2;
724             mdct_weight_pos = m;
725          } else if (sum[m+1] > 2*sum[m] && mdct_weight_shift < 2)
726          {
727             mdct_weight_shift=1;
728             mdct_weight_pos = m;
729          }
730          m++;
731       } while (m<st->mode->nbShortMdcts-1);
732 #endif
733       if (mdct_weight_shift)
734       {
735          mdct_shape(st->mode, X, mdct_weight_pos+1, st->mode->nbShortMdcts, N, st->mode->nbShortMdcts, mdct_weight_shift, C);
736          renormalise_bands(st->mode, X, C);
737       }
738    }
739
740
741    encode_flags(&enc, intra_ener, has_pitch, shortBlocks, has_fold);
742    if (has_pitch)
743    {
744       ec_enc_uint(&enc, pitch_index, MAX_PERIOD-(2*N-2*N4));
745       ec_enc_uint(&enc, gain_id, 16);
746    }
747    if (shortBlocks)
748    {
749       if (transient_shift)
750       {
751          ec_enc_uint(&enc, transient_shift, 4);
752          ec_enc_uint(&enc, transient_time, N+st->overlap);
753       } else {
754          ec_enc_uint(&enc, mdct_weight_shift, 4);
755          if (mdct_weight_shift && st->mode->nbShortMdcts!=2)
756             ec_enc_uint(&enc, mdct_weight_pos, st->mode->nbShortMdcts-1);
757       }
758    }
759
760    ALLOC(fine_quant, st->mode->nbEBands, int);
761    ALLOC(pulses, st->mode->nbEBands, int);
762
763    /* Computes the max bit-rate allowed in VBR more to avoid busting the budget */
764    if (st->vbr_rate>0)
765    {
766       celt_int32 vbr_bound, max_allowed;
767
768       vbr_bound = st->vbr_rate;
769       max_allowed = (st->vbr_rate + vbr_bound - st->vbr_reservoir)>>(BITRES+3);
770       if (max_allowed < 4)
771          max_allowed = 4;
772       if (max_allowed < nbCompressedBytes)
773          nbCompressedBytes = max_allowed;
774    }
775
776    /* Bit allocation */
777    ALLOC(error, C*st->mode->nbEBands, celt_word16);
778    coarse_needed = quant_coarse_energy(st->mode, bandLogE, st->oldBandE, nbCompressedBytes*4-8, intra_ener, st->mode->prob, error, &enc, C);
779    coarse_needed = ((coarse_needed*3-1)>>3)+1;
780    if (coarse_needed > nbCompressedBytes)
781       coarse_needed = nbCompressedBytes;
782    /* Variable bitrate */
783    if (st->vbr_rate>0)
784    {
785      celt_word16 alpha;
786      celt_int32 delta;
787      /* The target rate in 16th bits per frame */
788      celt_int32 target=st->vbr_rate;
789    
790      /* Shortblocks get a large boost in bitrate, but since they 
791         are uncommon long blocks are not greatly effected */
792      if (shortBlocks)
793        target*=2;
794      else if (st->mode->nbShortMdcts > 1)
795        target-=(target+14)/28;
796
797      /* The average energy is removed from the target and the actual 
798         energy added*/
799      target=target+st->vbr_offset-588+ec_enc_tell(&enc, BITRES);
800
801      /* In VBR mode the frame size must not be reduced so much that it would result in the coarse energy busting its budget */
802      target=IMAX(coarse_needed,(target+64)/128);
803      target=IMIN(nbCompressedBytes,target);
804      /* Make the adaptation coef (alpha) higher at the beginning */
805      if (st->vbr_count < 990)
806      {
807         st->vbr_count++;
808         alpha = celt_rcp(SHL32(EXTEND32(st->vbr_count+10),16));
809         /*printf ("%d %d\n", st->vbr_count+10, alpha);*/
810      } else
811         alpha = QCONST16(.001f,15);
812
813      /* By how much did we "miss" the target on that frame */
814      delta = (8<<BITRES)*(celt_int32)target - st->vbr_rate;
815      /* How many bits have we used in excess of what we're allowed */
816      st->vbr_reservoir += delta;
817      /*printf ("%d\n", st->vbr_reservoir);*/
818
819      /* Compute the offset we need to apply in order to reach the target */
820      st->vbr_drift += MULT16_32_Q15(alpha,delta-st->vbr_offset-st->vbr_drift);
821      st->vbr_offset = -st->vbr_drift;
822      /*printf ("%d\n", st->vbr_drift);*/
823
824      /* We could use any multiple of vbr_rate as bound (depending on the delay) */
825      if (st->vbr_reservoir < 0)
826      {
827         /* We're under the min value -- increase rate */
828         int adjust = 1-(st->vbr_reservoir-1)/(8<<BITRES);
829         st->vbr_reservoir += adjust*(8<<BITRES);
830         target += adjust;
831         /*printf ("+%d\n", adjust);*/
832      }
833      if (target < nbCompressedBytes)
834         nbCompressedBytes = target;
835      /* This moves the raw bits to take into account the new compressed size */
836      ec_byte_shrink(&buf, nbCompressedBytes);
837    }
838
839    ALLOC(offsets, st->mode->nbEBands, int);
840    ALLOC(fine_priority, st->mode->nbEBands, int);
841
842    for (i=0;i<st->mode->nbEBands;i++)
843       offsets[i] = 0;
844    bits = nbCompressedBytes*8 - ec_enc_tell(&enc, 0) - 1;
845    compute_allocation(st->mode, offsets, bits, pulses, fine_quant, fine_priority, C);
846
847    quant_fine_energy(st->mode, bandE, st->oldBandE, error, fine_quant, &enc, C);
848
849    /* Residual quantisation */
850    if (C==1)
851       quant_bands(st->mode, X, bandE, pulses, shortBlocks, has_fold, nbCompressedBytes*8, 1, &enc);
852 #ifndef DISABLE_STEREO
853    else
854       quant_bands_stereo(st->mode, X, bandE, pulses, shortBlocks, has_fold, nbCompressedBytes*8, &enc);
855 #endif
856
857    quant_energy_finalise(st->mode, bandE, st->oldBandE, error, fine_quant, fine_priority, nbCompressedBytes*8-ec_enc_tell(&enc, 0), &enc, C);
858
859    /* Re-synthesis of the coded audio if required */
860    if (st->pitch_available>0 || optional_synthesis!=NULL)
861    {
862       if (st->pitch_available>0 && st->pitch_available<MAX_PERIOD)
863         st->pitch_available+=st->frame_size;
864
865       if (mdct_weight_shift)
866       {
867          mdct_shape(st->mode, X, 0, mdct_weight_pos+1, N, st->mode->nbShortMdcts, mdct_weight_shift, C);
868       }
869
870       /* Synthesis */
871       denormalise_bands(st->mode, X, freq, bandE, C);
872
873       CELT_MOVE(st->out_mem, st->out_mem+C*N, C*(MAX_PERIOD+st->overlap-N));
874
875       if (has_pitch)
876          apply_pitch(st->mode, freq, pitch_freq, gain_id, 0, C);
877       
878       compute_inv_mdcts(st->mode, shortBlocks, freq, transient_time, transient_shift, st->out_mem, C);
879
880       /* De-emphasis and put everything back at the right place 
881          in the synthesis history */
882       if (optional_synthesis != NULL) {
883          deemphasis(st->out_mem, optional_synthesis, N, C, preemph, st->preemph_memD);
884
885       }
886    }
887
888    ec_enc_done(&enc);
889    
890    RESTORE_STACK;
891    return nbCompressedBytes;
892 }
893
894 #ifdef FIXED_POINT
895 #ifndef DISABLE_FLOAT_API
896 int celt_encode_float(CELTEncoder * restrict st, const float * pcm, float * optional_synthesis, unsigned char *compressed, int nbCompressedBytes)
897 {
898    int j, ret, C, N;
899    VARDECL(celt_int16, in);
900    SAVE_STACK;
901
902    if (check_encoder(st) != CELT_OK)
903       return CELT_INVALID_STATE;
904
905    if (check_mode(st->mode) != CELT_OK)
906       return CELT_INVALID_MODE;
907
908    if (pcm==NULL)
909       return CELT_BAD_ARG;
910
911    C = CHANNELS(st->channels);
912    N = st->block_size;
913    ALLOC(in, C*N, celt_int16);
914
915    for (j=0;j<C*N;j++)
916      in[j] = FLOAT2INT16(pcm[j]);
917
918    if (optional_synthesis != NULL) {
919      ret=celt_encode(st,in,in,compressed,nbCompressedBytes);
920       for (j=0;j<C*N;j++)
921          optional_synthesis[j]=in[j]*(1/32768.);
922    } else {
923      ret=celt_encode(st,in,NULL,compressed,nbCompressedBytes);
924    }
925    RESTORE_STACK;
926    return ret;
927
928 }
929 #endif /*DISABLE_FLOAT_API*/
930 #else
931 int celt_encode(CELTEncoder * restrict st, const celt_int16 * pcm, celt_int16 * optional_synthesis, unsigned char *compressed, int nbCompressedBytes)
932 {
933    int j, ret, C, N;
934    VARDECL(celt_sig, in);
935    SAVE_STACK;
936
937    if (check_encoder(st) != CELT_OK)
938       return CELT_INVALID_STATE;
939
940    if (check_mode(st->mode) != CELT_OK)
941       return CELT_INVALID_MODE;
942
943    if (pcm==NULL)
944       return CELT_BAD_ARG;
945
946    C=CHANNELS(st->channels);
947    N=st->block_size;
948    ALLOC(in, C*N, celt_sig);
949    for (j=0;j<C*N;j++) {
950      in[j] = SCALEOUT(pcm[j]);
951    }
952
953    if (optional_synthesis != NULL) {
954       ret = celt_encode_float(st,in,in,compressed,nbCompressedBytes);
955       for (j=0;j<C*N;j++)
956          optional_synthesis[j] = FLOAT2INT16(in[j]);
957    } else {
958       ret = celt_encode_float(st,in,NULL,compressed,nbCompressedBytes);
959    }
960    RESTORE_STACK;
961    return ret;
962 }
963 #endif
964
965 int celt_encoder_ctl(CELTEncoder * restrict st, int request, ...)
966 {
967    va_list ap;
968    
969    if (check_encoder(st) != CELT_OK)
970       return CELT_INVALID_STATE;
971
972    va_start(ap, request);
973    if ((request!=CELT_GET_MODE_REQUEST) && (check_mode(st->mode) != CELT_OK))
974      goto bad_mode;
975    switch (request)
976    {
977       case CELT_GET_MODE_REQUEST:
978       {
979          const CELTMode ** value = va_arg(ap, const CELTMode**);
980          if (value==0)
981             goto bad_arg;
982          *value=st->mode;
983       }
984       break;
985       case CELT_SET_COMPLEXITY_REQUEST:
986       {
987          int value = va_arg(ap, celt_int32);
988          if (value<0 || value>10)
989             goto bad_arg;
990          if (value<=2) {
991             st->pitch_enabled = 0; 
992             st->pitch_available = 0;
993          } else {
994               st->pitch_enabled = 1;
995               if (st->pitch_available<1)
996                 st->pitch_available = 1;
997          }   
998       }
999       break;
1000       case CELT_SET_PREDICTION_REQUEST:
1001       {
1002          int value = va_arg(ap, celt_int32);
1003          if (value<0 || value>2)
1004             goto bad_arg;
1005          if (value==0)
1006          {
1007             st->force_intra   = 1;
1008             st->pitch_permitted = 0;
1009          } else if (value==1) {
1010             st->force_intra   = 0;
1011             st->pitch_permitted = 0;
1012          } else {
1013             st->force_intra   = 0;
1014             st->pitch_permitted = 1;
1015          }   
1016       }
1017       break;
1018       case CELT_SET_VBR_RATE_REQUEST:
1019       {
1020          celt_int32 value = va_arg(ap, celt_int32);
1021          if (value<0)
1022             goto bad_arg;
1023          if (value>3072000)
1024             value = 3072000;
1025          st->vbr_rate = ((st->mode->Fs<<3)+(st->block_size>>1))/st->block_size;
1026          st->vbr_rate = ((value<<7)+(st->vbr_rate>>1))/st->vbr_rate;
1027       }
1028       break;
1029       case CELT_RESET_STATE:
1030       {
1031          const CELTMode *mode = st->mode;
1032          int C = st->channels;
1033
1034          if (st->pitch_available > 0) st->pitch_available = 1;
1035
1036          CELT_MEMSET(st->in_mem, 0, st->overlap*C);
1037          CELT_MEMSET(st->out_mem, 0, (MAX_PERIOD+st->overlap)*C);
1038
1039          CELT_MEMSET(st->oldBandE, 0, C*mode->nbEBands);
1040
1041          CELT_MEMSET(st->preemph_memE, 0, C);
1042          CELT_MEMSET(st->preemph_memD, 0, C);
1043          st->delayedIntra = 1;
1044       }
1045       break;
1046       default:
1047          goto bad_request;
1048    }
1049    va_end(ap);
1050    return CELT_OK;
1051 bad_mode:
1052   va_end(ap);
1053   return CELT_INVALID_MODE;
1054 bad_arg:
1055    va_end(ap);
1056    return CELT_BAD_ARG;
1057 bad_request:
1058    va_end(ap);
1059    return CELT_UNIMPLEMENTED;
1060 }
1061
1062 /**********************************************************************/
1063 /*                                                                    */
1064 /*                             DECODER                                */
1065 /*                                                                    */
1066 /**********************************************************************/
1067 #ifdef NEW_PLC
1068 #define DECODE_BUFFER_SIZE 2048
1069 #else
1070 #define DECODE_BUFFER_SIZE MAX_PERIOD
1071 #endif
1072
1073 #define DECODERVALID   0x4c434454
1074 #define DECODERPARTIAL 0x5444434c
1075 #define DECODERFREED   0x4c004400
1076
1077 /** Decoder state 
1078  @brief Decoder state
1079  */
1080 struct CELTDecoder {
1081    celt_uint32 marker;
1082    const CELTMode *mode;
1083    int frame_size;
1084    int block_size;
1085    int overlap;
1086    int channels;
1087
1088    ec_byte_buffer buf;
1089    ec_enc         enc;
1090
1091    celt_sig * restrict preemph_memD;
1092
1093    celt_sig *out_mem;
1094    celt_sig *decode_mem;
1095
1096    celt_word16 *oldBandE;
1097    
1098 #ifdef NEW_PLC
1099    float *lpc;
1100 #endif
1101
1102    int last_pitch_index;
1103    int loss_count;
1104 };
1105
1106 int check_decoder(const CELTDecoder *st) 
1107 {
1108    if (st==NULL)
1109    {
1110       celt_warning("NULL passed a decoder structure");  
1111       return CELT_INVALID_STATE;
1112    }
1113    if (st->marker == DECODERVALID)
1114       return CELT_OK;
1115    if (st->marker == DECODERFREED)
1116       celt_warning("Referencing a decoder that has already been freed");
1117    else
1118       celt_warning("This is not a valid CELT decoder structure");
1119    return CELT_INVALID_STATE;
1120 }
1121
1122 CELTDecoder *celt_decoder_create(const CELTMode *mode, int channels, int *error)
1123 {
1124    int N, C;
1125    CELTDecoder *st;
1126
1127    if (check_mode(mode) != CELT_OK)
1128    {
1129       if (error)
1130          *error = CELT_INVALID_MODE;
1131       return NULL;
1132    }
1133 #ifdef DISABLE_STEREO
1134    if (channels > 1)
1135    {
1136       celt_warning("Stereo support was disable from this build");
1137       if (error)
1138          *error = CELT_BAD_ARG;
1139       return NULL;
1140    }
1141 #endif
1142
1143    if (channels < 0 || channels > 2)
1144    {
1145       celt_warning("Only mono and stereo supported");
1146       if (error)
1147          *error = CELT_BAD_ARG;
1148       return NULL;
1149    }
1150
1151    N = mode->mdctSize;
1152    C = CHANNELS(channels);
1153    st = celt_alloc(sizeof(CELTDecoder));
1154
1155    if (st==NULL)
1156    {
1157       if (error)
1158          *error = CELT_ALLOC_FAIL;
1159       return NULL;
1160    }
1161
1162    st->marker = DECODERPARTIAL;
1163    st->mode = mode;
1164    st->frame_size = N;
1165    st->block_size = N;
1166    st->overlap = mode->overlap;
1167    st->channels = channels;
1168
1169    st->decode_mem = celt_alloc((DECODE_BUFFER_SIZE+st->overlap)*C*sizeof(celt_sig));
1170    st->out_mem = st->decode_mem+DECODE_BUFFER_SIZE-MAX_PERIOD;
1171    
1172    st->oldBandE = (celt_word16*)celt_alloc(C*mode->nbEBands*sizeof(celt_word16));
1173    
1174    st->preemph_memD = (celt_sig*)celt_alloc(C*sizeof(celt_sig));
1175
1176 #ifdef NEW_PLC
1177    st->lpc = (float*)celt_alloc(C*LPC_ORDER*sizeof(float));
1178 #endif
1179
1180    st->loss_count = 0;
1181
1182    if ((st->decode_mem!=NULL) && (st->out_mem!=NULL) && (st->oldBandE!=NULL) &&
1183 #ifdef NEW_PLC
1184          (st->lpc!=NULL) &&
1185 #endif
1186        (st->preemph_memD!=NULL))
1187    {
1188       if (error)
1189          *error = CELT_OK;
1190       st->marker = DECODERVALID;
1191       return st;
1192    }
1193    /* If the setup fails for some reason deallocate it. */
1194    celt_decoder_destroy(st);
1195    if (error)
1196       *error = CELT_ALLOC_FAIL;
1197    return NULL;
1198 }
1199
1200 void celt_decoder_destroy(CELTDecoder *st)
1201 {
1202    if (st == NULL)
1203    {
1204       celt_warning("NULL passed to celt_decoder_destroy");
1205       return;
1206    }
1207
1208    if (st->marker == DECODERFREED) 
1209    {
1210       celt_warning("Freeing a decoder which has already been freed"); 
1211       return;
1212    }
1213    
1214    if (st->marker != DECODERVALID && st->marker != DECODERPARTIAL)
1215    {
1216       celt_warning("This is not a valid CELT decoder structure");
1217       return;
1218    }
1219    
1220    /*Check_mode is non-fatal here because we can still free
1221      the encoder memory even if the mode is bad, although calling
1222      the free functions in this order is a violation of the API.*/
1223    check_mode(st->mode);
1224    
1225    celt_free(st->decode_mem);
1226    celt_free(st->oldBandE);
1227    celt_free(st->preemph_memD);
1228
1229 #ifdef NEW_PLC
1230    celt_free(st->lpc);
1231 #endif
1232    
1233    st->marker = DECODERFREED;
1234    
1235    celt_free(st);
1236 }
1237
1238 static void celt_decode_lost(CELTDecoder * restrict st, celt_word16 * restrict pcm)
1239 {
1240    int c, N;
1241    int pitch_index;
1242    int overlap = st->mode->overlap;
1243    celt_word16 fade = Q15ONE;
1244    int i, len;
1245    VARDECL(celt_sig, freq);
1246    const int C = CHANNELS(st->channels);
1247    int offset;
1248    SAVE_STACK;
1249    N = st->block_size;
1250    
1251    len = N+st->mode->overlap;
1252    
1253    if (st->loss_count == 0)
1254    {
1255       celt_word16 pitch_buf[MAX_PERIOD>>1];
1256       celt_word32 tmp=0;
1257       celt_word32 mem0[2]={0,0};
1258       celt_word16 mem1[2]={0,0};
1259       pitch_downsample(st->out_mem, pitch_buf, MAX_PERIOD, MAX_PERIOD,
1260                        C, mem0, mem1);
1261       pitch_search(st->mode, pitch_buf+((MAX_PERIOD-len)>>1), pitch_buf, len,
1262                    MAX_PERIOD-len-100, &pitch_index, &tmp);
1263       pitch_index = MAX_PERIOD-len-pitch_index;
1264       st->last_pitch_index = pitch_index;
1265    } else {
1266       pitch_index = st->last_pitch_index;
1267       if (st->loss_count < 5)
1268          fade = QCONST16(.8f,15);
1269       else
1270          fade = 0;
1271    }
1272
1273 #ifndef NEW_PLC
1274    offset = MAX_PERIOD-pitch_index;
1275    ALLOC(freq,C*N, celt_sig); /**< Interleaved signal MDCTs */
1276    while (offset+len >= MAX_PERIOD)
1277       offset -= pitch_index;
1278    compute_mdcts(st->mode, 0, st->out_mem+offset*C, freq, C);
1279    for (i=0;i<C*N;i++)
1280       freq[i] = ADD32(VERY_SMALL, MULT16_32_Q15(fade,freq[i]));
1281
1282    CELT_MOVE(st->out_mem, st->out_mem+C*N, C*(MAX_PERIOD+st->mode->overlap-N));
1283    /* Compute inverse MDCTs */
1284    compute_inv_mdcts(st->mode, 0, freq, -1, 0, st->out_mem, C);
1285 #else
1286    for (c=0;c<C;c++)
1287    {
1288       celt_word32 e[MAX_PERIOD];
1289       float exc[MAX_PERIOD];
1290       float ac[LPC_ORDER+1];
1291       float decay = 1;
1292       float S1=0;
1293       float mem[LPC_ORDER]={0};
1294
1295       offset = MAX_PERIOD-pitch_index;
1296       for (i=0;i<MAX_PERIOD;i++)
1297          exc[i] = st->out_mem[i*C+c];
1298
1299       if (st->loss_count == 0)
1300       {
1301          _celt_autocorr(exc, ac, st->mode->window, st->mode->overlap,
1302                         LPC_ORDER, MAX_PERIOD);
1303
1304          /* Noise floor -50 dB */
1305          ac[0] *= 1.00001;
1306          /* Lag windowing */
1307          for (i=1;i<=LPC_ORDER;i++)
1308          {
1309             /*ac[i] *= exp(-.5*(2*M_PI*.002*i)*(2*M_PI*.002*i));*/
1310             ac[i] -= ac[i]*(.008*i)*(.008*i);
1311          }
1312
1313          _celt_lpc(st->lpc, ac, LPC_ORDER);
1314       }
1315       fir(exc, st->lpc, exc, MAX_PERIOD, LPC_ORDER, mem);
1316
1317       /* Check if the waveform is decaying (and if so how fast) */
1318       {
1319          float E1=0, E2=0;
1320          int period;
1321          if (pitch_index <= MAX_PERIOD/2)
1322             period = pitch_index;
1323          else
1324             period = MAX_PERIOD/2;
1325          for (i=0;i<period;i++)
1326          {
1327             E1 += exc[MAX_PERIOD-period+i]*exc[MAX_PERIOD-period+i];
1328             E2 += exc[MAX_PERIOD-2*period+i]*exc[MAX_PERIOD-2*period+i];
1329          }
1330          decay = sqrt((E1+1)/(E2+1));
1331          if (decay > 1)
1332             decay = 1;
1333       }
1334
1335       /* Copy excitation, taking decay into account */
1336       for (i=0;i<len+st->mode->overlap;i++)
1337       {
1338          if (offset+i >= MAX_PERIOD)
1339          {
1340             offset -= pitch_index;
1341             decay *= decay;
1342          }
1343          e[i] = decay*exc[offset+i];
1344          S1 += st->out_mem[offset+i]*1.*st->out_mem[offset+i];
1345       }
1346
1347       iir(e, st->lpc, e, len+st->mode->overlap, LPC_ORDER, mem);
1348
1349       {
1350          float S2=0;
1351          for (i=0;i<len+overlap;i++)
1352             S2 += e[i]*1.*e[i];
1353          /* This checks for an "explosion" in the synthesis (including NaNs) */
1354          if (!(S1 > 0.2f*S2))
1355          {
1356             for (i=0;i<len+overlap;i++)
1357                e[i] = 0;
1358          } else if (S1 < S2)
1359          {
1360             float ratio = sqrt((S1+1)/(S2+1));
1361             for (i=0;i<len+overlap;i++)
1362                e[i] *= ratio;
1363          }
1364       }
1365
1366       for (i=0;i<MAX_PERIOD+st->mode->overlap-N;i++)
1367          st->out_mem[C*i+c] = st->out_mem[C*(N+i)+c];
1368
1369       /* Apply TDAC to the concealed audio so that it blends with the
1370          previous and next frames */
1371       for (i=0;i<overlap/2;i++)
1372       {
1373          celt_word32 tmp1, tmp2;
1374          tmp1 = MULT16_32_Q15(st->mode->window[i          ], e[i          ]) -
1375                 MULT16_32_Q15(st->mode->window[overlap-i-1], e[overlap-i-1]);
1376          tmp2 = MULT16_32_Q15(st->mode->window[i],           e[N+overlap-1-i]) +
1377                 MULT16_32_Q15(st->mode->window[overlap-i-1], e[N+i          ]);
1378          tmp1 = MULT16_32_Q15(fade, tmp1);
1379          tmp2 = MULT16_32_Q15(fade, tmp2);
1380          st->out_mem[C*(MAX_PERIOD+i)+c] = MULT16_32_Q15(st->mode->window[overlap-i-1], tmp2);
1381          st->out_mem[C*(MAX_PERIOD+overlap-i-1)+c] = MULT16_32_Q15(st->mode->window[i], tmp2);
1382          st->out_mem[C*(MAX_PERIOD-N+i)+c] += MULT16_32_Q15(st->mode->window[i], tmp1);
1383          st->out_mem[C*(MAX_PERIOD-N+overlap-i-1)+c] -= MULT16_32_Q15(st->mode->window[overlap-i-1], tmp1);
1384       }
1385       for (i=0;i<N-overlap;i++)
1386          st->out_mem[C*(MAX_PERIOD-N+overlap+i)+c] = MULT16_32_Q15(fade, e[overlap+i]);
1387    }
1388 #endif
1389
1390    deemphasis(st->out_mem, pcm, N, C, preemph, st->preemph_memD);
1391    
1392    st->loss_count++;
1393
1394    RESTORE_STACK;
1395 }
1396
1397 #ifdef FIXED_POINT
1398 int celt_decode(CELTDecoder * restrict st, const unsigned char *data, int len, celt_int16 * restrict pcm)
1399 {
1400 #else
1401 int celt_decode_float(CELTDecoder * restrict st, const unsigned char *data, int len, celt_sig * restrict pcm)
1402 {
1403 #endif
1404    int i, N, N4;
1405    int has_pitch, has_fold;
1406    int pitch_index;
1407    int bits;
1408    ec_dec dec;
1409    ec_byte_buffer buf;
1410    VARDECL(celt_sig, freq);
1411    VARDECL(celt_sig, pitch_freq);
1412    VARDECL(celt_norm, X);
1413    VARDECL(celt_ener, bandE);
1414    VARDECL(int, fine_quant);
1415    VARDECL(int, pulses);
1416    VARDECL(int, offsets);
1417    VARDECL(int, fine_priority);
1418
1419    int shortBlocks;
1420    int intra_ener;
1421    int transient_time;
1422    int transient_shift;
1423    int mdct_weight_shift=0;
1424    const int C = CHANNELS(st->channels);
1425    int mdct_weight_pos=0;
1426    int gain_id=0;
1427    SAVE_STACK;
1428
1429    if (check_decoder(st) != CELT_OK)
1430       return CELT_INVALID_STATE;
1431
1432    if (check_mode(st->mode) != CELT_OK)
1433       return CELT_INVALID_MODE;
1434
1435    if (pcm==NULL)
1436       return CELT_BAD_ARG;
1437
1438    N = st->block_size;
1439    N4 = (N-st->overlap)>>1;
1440
1441    ALLOC(freq, C*N, celt_sig); /**< Interleaved signal MDCTs */
1442    ALLOC(X, C*N, celt_norm);   /**< Interleaved normalised MDCTs */
1443    ALLOC(bandE, st->mode->nbEBands*C, celt_ener);
1444    
1445    if (data == NULL)
1446    {
1447       celt_decode_lost(st, pcm);
1448       RESTORE_STACK;
1449       return 0;
1450    }
1451    if (len<0) {
1452      RESTORE_STACK;
1453      return CELT_BAD_ARG;
1454    }
1455    
1456    ec_byte_readinit(&buf,(unsigned char*)data,len);
1457    ec_dec_init(&dec,&buf);
1458    
1459    decode_flags(&dec, &intra_ener, &has_pitch, &shortBlocks, &has_fold);
1460    if (shortBlocks)
1461    {
1462       transient_shift = ec_dec_uint(&dec, 4);
1463       if (transient_shift == 3)
1464       {
1465          transient_time = ec_dec_uint(&dec, N+st->mode->overlap);
1466       } else {
1467          mdct_weight_shift = transient_shift;
1468          if (mdct_weight_shift && st->mode->nbShortMdcts>2)
1469             mdct_weight_pos = ec_dec_uint(&dec, st->mode->nbShortMdcts-1);
1470          transient_shift = 0;
1471          transient_time = 0;
1472       }
1473    } else {
1474       transient_time = -1;
1475       transient_shift = 0;
1476    }
1477    
1478    if (has_pitch)
1479    {
1480       pitch_index = ec_dec_uint(&dec, MAX_PERIOD-(2*N-2*N4));
1481       gain_id = ec_dec_uint(&dec, 16);
1482    } else {
1483       pitch_index = 0;
1484    }
1485
1486    ALLOC(fine_quant, st->mode->nbEBands, int);
1487    /* Get band energies */
1488    unquant_coarse_energy(st->mode, bandE, st->oldBandE, len*4-8, intra_ener, st->mode->prob, &dec, C);
1489    
1490    ALLOC(pulses, st->mode->nbEBands, int);
1491    ALLOC(offsets, st->mode->nbEBands, int);
1492    ALLOC(fine_priority, st->mode->nbEBands, int);
1493
1494    for (i=0;i<st->mode->nbEBands;i++)
1495       offsets[i] = 0;
1496
1497    bits = len*8 - ec_dec_tell(&dec, 0) - 1;
1498    compute_allocation(st->mode, offsets, bits, pulses, fine_quant, fine_priority, C);
1499    /*bits = ec_dec_tell(&dec, 0);
1500    compute_fine_allocation(st->mode, fine_quant, (20*C+len*8/5-(ec_dec_tell(&dec, 0)-bits))/C);*/
1501    
1502    unquant_fine_energy(st->mode, bandE, st->oldBandE, fine_quant, &dec, C);
1503
1504    ALLOC(pitch_freq, C*N, celt_sig); /**< Interleaved signal MDCTs */
1505    if (has_pitch) 
1506    {
1507       /* Pitch MDCT */
1508       compute_mdcts(st->mode, 0, st->out_mem+pitch_index*C, pitch_freq, C);
1509    }
1510
1511    /* Decode fixed codebook and merge with pitch */
1512    if (C==1)
1513       quant_bands(st->mode, X, bandE, pulses, shortBlocks, has_fold, len*8, 0, &dec);
1514 #ifndef DISABLE_STEREO
1515    else
1516       unquant_bands_stereo(st->mode, X, bandE, pulses, shortBlocks, has_fold, len*8, &dec);
1517 #endif
1518    unquant_energy_finalise(st->mode, bandE, st->oldBandE, fine_quant, fine_priority, len*8-ec_dec_tell(&dec, 0), &dec, C);
1519    
1520    if (mdct_weight_shift)
1521    {
1522       mdct_shape(st->mode, X, 0, mdct_weight_pos+1, N, st->mode->nbShortMdcts, mdct_weight_shift, C);
1523    }
1524
1525    /* Synthesis */
1526    denormalise_bands(st->mode, X, freq, bandE, C);
1527
1528
1529    CELT_MOVE(st->decode_mem, st->decode_mem+C*N, C*(DECODE_BUFFER_SIZE+st->overlap-N));
1530
1531    if (has_pitch)
1532       apply_pitch(st->mode, freq, pitch_freq, gain_id, 0, C);
1533
1534    /* Compute inverse MDCTs */
1535    compute_inv_mdcts(st->mode, shortBlocks, freq, transient_time, transient_shift, st->out_mem, C);
1536
1537    deemphasis(st->out_mem, pcm, N, C, preemph, st->preemph_memD);
1538    st->loss_count = 0;
1539    RESTORE_STACK;
1540    return 0;
1541 }
1542
1543 #ifdef FIXED_POINT
1544 #ifndef DISABLE_FLOAT_API
1545 int celt_decode_float(CELTDecoder * restrict st, const unsigned char *data, int len, float * restrict pcm)
1546 {
1547    int j, ret, C, N;
1548    VARDECL(celt_int16, out);
1549    SAVE_STACK;
1550
1551    if (check_decoder(st) != CELT_OK)
1552       return CELT_INVALID_STATE;
1553
1554    if (check_mode(st->mode) != CELT_OK)
1555       return CELT_INVALID_MODE;
1556
1557    if (pcm==NULL)
1558       return CELT_BAD_ARG;
1559
1560    C = CHANNELS(st->channels);
1561    N = st->block_size;
1562    
1563    ALLOC(out, C*N, celt_int16);
1564    ret=celt_decode(st, data, len, out);
1565    for (j=0;j<C*N;j++)
1566       pcm[j]=out[j]*(1/32768.);
1567      
1568    RESTORE_STACK;
1569    return ret;
1570 }
1571 #endif /*DISABLE_FLOAT_API*/
1572 #else
1573 int celt_decode(CELTDecoder * restrict st, const unsigned char *data, int len, celt_int16 * restrict pcm)
1574 {
1575    int j, ret, C, N;
1576    VARDECL(celt_sig, out);
1577    SAVE_STACK;
1578
1579    if (check_decoder(st) != CELT_OK)
1580       return CELT_INVALID_STATE;
1581
1582    if (check_mode(st->mode) != CELT_OK)
1583       return CELT_INVALID_MODE;
1584
1585    if (pcm==NULL)
1586       return CELT_BAD_ARG;
1587
1588    C = CHANNELS(st->channels);
1589    N = st->block_size;
1590    ALLOC(out, C*N, celt_sig);
1591
1592    ret=celt_decode_float(st, data, len, out);
1593
1594    for (j=0;j<C*N;j++)
1595       pcm[j] = FLOAT2INT16 (out[j]);
1596    
1597    RESTORE_STACK;
1598    return ret;
1599 }
1600 #endif
1601
1602 int celt_decoder_ctl(CELTDecoder * restrict st, int request, ...)
1603 {
1604    va_list ap;
1605
1606    if (check_decoder(st) != CELT_OK)
1607       return CELT_INVALID_STATE;
1608
1609    va_start(ap, request);
1610    if ((request!=CELT_GET_MODE_REQUEST) && (check_mode(st->mode) != CELT_OK))
1611      goto bad_mode;
1612    switch (request)
1613    {
1614       case CELT_GET_MODE_REQUEST:
1615       {
1616          const CELTMode ** value = va_arg(ap, const CELTMode**);
1617          if (value==0)
1618             goto bad_arg;
1619          *value=st->mode;
1620       }
1621       break;
1622       case CELT_RESET_STATE:
1623       {
1624          const CELTMode *mode = st->mode;
1625          int C = st->channels;
1626
1627          CELT_MEMSET(st->decode_mem, 0, (DECODE_BUFFER_SIZE+st->overlap)*C);
1628          CELT_MEMSET(st->oldBandE, 0, C*mode->nbEBands);
1629
1630          CELT_MEMSET(st->preemph_memD, 0, C);
1631
1632          st->loss_count = 0;
1633       }
1634       break;
1635       default:
1636          goto bad_request;
1637    }
1638    va_end(ap);
1639    return CELT_OK;
1640 bad_mode:
1641   va_end(ap);
1642   return CELT_INVALID_MODE;
1643 bad_arg:
1644    va_end(ap);
1645    return CELT_BAD_ARG;
1646 bad_request:
1647       va_end(ap);
1648   return CELT_UNIMPLEMENTED;
1649 }
1650
1651 const char *celt_strerror(int error)
1652 {
1653    static const char *error_strings[8] = {
1654       "success",
1655       "invalid argument",
1656       "invalid mode",
1657       "internal error",
1658       "corrupted stream",
1659       "request not implemented",
1660       "invalid state",
1661       "memory allocation failed"
1662    };
1663    if (error > 0 || error < -7)
1664       return "unknown error";
1665    else 
1666       return error_strings[-error];
1667 }
1668