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