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