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