Making the modified transient code work with stereo as well
[opus.git] / libcelt / celt.c
1 /* (C) 2007-2008 Jean-Marc Valin, CSIRO
2    (C) 2008 Gregory Maxwell */
3 /*
4    Redistribution and use in source and binary forms, with or without
5    modification, are permitted provided that the following conditions
6    are met:
7    
8    - Redistributions of source code must retain the above copyright
9    notice, this list of conditions and the following disclaimer.
10    
11    - Redistributions in binary form must reproduce the above copyright
12    notice, this list of conditions and the following disclaimer in the
13    documentation and/or other materials provided with the distribution.
14    
15    - Neither the name of the Xiph.org Foundation nor the names of its
16    contributors may be used to endorse or promote products derived from
17    this software without specific prior written permission.
18    
19    THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
20    ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
21    LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
22    A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR
23    CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
24    EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
25    PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
26    PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
27    LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
28    NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
29    SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
30 */
31
32 #ifdef HAVE_CONFIG_H
33 #include "config.h"
34 #endif
35
36 #define CELT_C
37
38 #include "os_support.h"
39 #include "mdct.h"
40 #include <math.h>
41 #include "celt.h"
42 #include "pitch.h"
43 #include "kiss_fftr.h"
44 #include "bands.h"
45 #include "modes.h"
46 #include "entcode.h"
47 #include "quant_bands.h"
48 #include "psy.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 static const celt_word16_t preemph = QCONST16(0.8f,15);
56
57 #ifdef FIXED_POINT
58 static const celt_word16_t transientWindow[16] = {
59      279,  1106,  2454,  4276,  6510,  9081, 11900, 14872,
60    17896, 20868, 23687, 26258, 28492, 30314, 31662, 32489};
61 #else
62 static const float transientWindow[16] = {
63    0.0085135, 0.0337639, 0.0748914, 0.1304955, 0.1986827, 0.2771308, 0.3631685, 0.4538658,
64    0.5461342, 0.6368315, 0.7228692, 0.8013173, 0.8695045, 0.9251086, 0.9662361, 0.9914865};
65 #endif
66
67    
68 /** Encoder state 
69  @brief Encoder state
70  */
71 struct CELTEncoder {
72    const CELTMode *mode;     /**< Mode used by the encoder */
73    int frame_size;
74    int block_size;
75    int overlap;
76    int channels;
77    
78    int pitch_enabled;
79    int pitch_available;
80    int delayedIntra;
81    int VBR_rate; /* Target number of 16th bits per frame */
82    celt_word16_t * restrict preemph_memE; /* Input is 16-bit, so why bother with 32 */
83    celt_sig_t    * restrict preemph_memD;
84
85    celt_sig_t *in_mem;
86    celt_sig_t *out_mem;
87
88    celt_word16_t *oldBandE;
89 #ifdef EXP_PSY
90    celt_word16_t *psy_mem;
91    struct PsyDecay psy;
92 #endif
93 };
94
95 CELTEncoder *celt_encoder_create(const CELTMode *mode)
96 {
97    int N, C;
98    CELTEncoder *st;
99
100    if (check_mode(mode) != CELT_OK)
101       return NULL;
102
103    N = mode->mdctSize;
104    C = mode->nbChannels;
105    st = celt_alloc(sizeof(CELTEncoder));
106    
107    st->mode = mode;
108    st->frame_size = N;
109    st->block_size = N;
110    st->overlap = mode->overlap;
111
112    st->VBR_rate = 0;
113    st->pitch_enabled = 1;
114    st->pitch_available = 1;
115    st->delayedIntra = 1;
116
117    st->in_mem = celt_alloc(st->overlap*C*sizeof(celt_sig_t));
118    st->out_mem = celt_alloc((MAX_PERIOD+st->overlap)*C*sizeof(celt_sig_t));
119
120    st->oldBandE = (celt_word16_t*)celt_alloc(C*mode->nbEBands*sizeof(celt_word16_t));
121
122    st->preemph_memE = (celt_word16_t*)celt_alloc(C*sizeof(celt_word16_t));
123    st->preemph_memD = (celt_sig_t*)celt_alloc(C*sizeof(celt_sig_t));
124
125 #ifdef EXP_PSY
126    st->psy_mem = celt_alloc(MAX_PERIOD*sizeof(celt_word16_t));
127    psydecay_init(&st->psy, MAX_PERIOD/2, st->mode->Fs);
128 #endif
129
130    return st;
131 }
132
133 void celt_encoder_destroy(CELTEncoder *st)
134 {
135    if (st == NULL)
136    {
137       celt_warning("NULL passed to celt_encoder_destroy");
138       return;
139    }
140    if (check_mode(st->mode) != CELT_OK)
141       return;
142
143    celt_free(st->in_mem);
144    celt_free(st->out_mem);
145    
146    celt_free(st->oldBandE);
147    
148    celt_free(st->preemph_memE);
149    celt_free(st->preemph_memD);
150    
151 #ifdef EXP_PSY
152    celt_free (st->psy_mem);
153    psydecay_clear(&st->psy);
154 #endif
155    
156    celt_free(st);
157 }
158
159 static inline celt_int16_t FLOAT2INT16(float x)
160 {
161    x = x*CELT_SIG_SCALE;
162    x = MAX32(x, -32768);
163    x = MIN32(x, 32767);
164    return (celt_int16_t)float2int(x);
165 }
166
167 static inline celt_word16_t SIG2WORD16(celt_sig_t x)
168 {
169 #ifdef FIXED_POINT
170    x = PSHR32(x, SIG_SHIFT);
171    x = MAX32(x, -32768);
172    x = MIN32(x, 32767);
173    return EXTRACT16(x);
174 #else
175    return (celt_word16_t)x;
176 #endif
177 }
178
179 static int transient_analysis(celt_word32_t *in, int len, int C, int *transient_time, int *transient_shift)
180 {
181    int c, i, n;
182    celt_word32_t ratio;
183    /* FIXME: Remove the floats here */
184    VARDECL(celt_word32_t, begin);
185    SAVE_STACK;
186    ALLOC(begin, len, celt_word32_t);
187    for (i=0;i<len;i++)
188       begin[i] = ABS32(SHR32(in[C*i],SIG_SHIFT));
189    for (c=1;c<C;c++)
190    {
191       for (i=0;i<len;i++)
192          begin[i] = MAX32(begin[i], ABS32(SHR32(in[C*i+c],SIG_SHIFT)));
193    }
194    for (i=1;i<len;i++)
195       begin[i] = MAX32(begin[i-1],begin[i]);
196    n = -1;
197    for (i=8;i<len-8;i++)
198    {
199       if (begin[i] < MULT16_32_Q15(QCONST16(.2f,15),begin[len-1]))
200          n=i;
201    }
202    if (n<32)
203    {
204       n = -1;
205       ratio = 0;
206    } else {
207       ratio = DIV32(begin[len-1],1+begin[n-16]);
208    }
209    /*printf ("%d %f\n", n, ratio*ratio);*/
210    if (ratio < 0)
211       ratio = 0;
212    if (ratio > 1000)
213       ratio = 1000;
214    ratio *= ratio;
215    
216    if (ratio > 2048)
217       *transient_shift = 3;
218    else
219       *transient_shift = 0;
220    
221    *transient_time = n;
222    
223    RESTORE_STACK;
224    return ratio > 20;
225 }
226
227 /** Apply window and compute the MDCT for all sub-frames and all channels in a frame */
228 static void compute_mdcts(const CELTMode *mode, int shortBlocks, celt_sig_t * restrict in, celt_sig_t * restrict out)
229 {
230    const int C = CHANNELS(mode);
231    if (C==1 && !shortBlocks)
232    {
233       const mdct_lookup *lookup = MDCT(mode);
234       const int overlap = OVERLAP(mode);
235       mdct_forward(lookup, in, out, mode->window, overlap);
236    } else if (!shortBlocks) {
237       const mdct_lookup *lookup = MDCT(mode);
238       const int overlap = OVERLAP(mode);
239       const int N = FRAMESIZE(mode);
240       int c;
241       VARDECL(celt_word32_t, x);
242       VARDECL(celt_word32_t, tmp);
243       SAVE_STACK;
244       ALLOC(x, N+overlap, celt_word32_t);
245       ALLOC(tmp, N, celt_word32_t);
246       for (c=0;c<C;c++)
247       {
248          int j;
249          for (j=0;j<N+overlap;j++)
250             x[j] = in[C*j+c];
251          mdct_forward(lookup, x, tmp, mode->window, overlap);
252          /* Interleaving the sub-frames */
253          for (j=0;j<N;j++)
254             out[C*j+c] = tmp[j];
255       }
256       RESTORE_STACK;
257    } else {
258       const mdct_lookup *lookup = &mode->shortMdct;
259       const int overlap = mode->overlap;
260       const int N = mode->shortMdctSize;
261       int b, c;
262       VARDECL(celt_word32_t, x);
263       VARDECL(celt_word32_t, tmp);
264       SAVE_STACK;
265       ALLOC(x, N+overlap, celt_word32_t);
266       ALLOC(tmp, N, celt_word32_t);
267       for (c=0;c<C;c++)
268       {
269          int B = mode->nbShortMdcts;
270          for (b=0;b<B;b++)
271          {
272             int j;
273             for (j=0;j<N+overlap;j++)
274                x[j] = in[C*(b*N+j)+c];
275             mdct_forward(lookup, x, tmp, mode->window, overlap);
276             /* Interleaving the sub-frames */
277             for (j=0;j<N;j++)
278                out[C*(j*B+b)+c] = tmp[j];
279          }
280       }
281       RESTORE_STACK;
282    }
283 }
284
285 /** Compute the IMDCT and apply window for all sub-frames and all channels in a frame */
286 static void compute_inv_mdcts(const CELTMode *mode, int shortBlocks, celt_sig_t *X, int transient_time, int transient_shift, celt_sig_t * restrict out_mem)
287 {
288    int c, N4;
289    const int C = CHANNELS(mode);
290    const int N = FRAMESIZE(mode);
291    const int overlap = OVERLAP(mode);
292    N4 = (N-overlap)>>1;
293    for (c=0;c<C;c++)
294    {
295       int j;
296       if (transient_shift==0 && C==1 && !shortBlocks) {
297          const mdct_lookup *lookup = MDCT(mode);
298          mdct_backward(lookup, X, out_mem+C*(MAX_PERIOD-N-N4), mode->window, overlap);
299       } else if (!shortBlocks) {
300          const mdct_lookup *lookup = MDCT(mode);
301          VARDECL(celt_word32_t, x);
302          VARDECL(celt_word32_t, tmp);
303          SAVE_STACK;
304          ALLOC(x, 2*N, celt_word32_t);
305          ALLOC(tmp, N, celt_word32_t);
306          /* De-interleaving the sub-frames */
307          for (j=0;j<N;j++)
308             tmp[j] = X[C*j+c];
309          /* Prevents problems from the imdct doing the overlap-add */
310          CELT_MEMSET(x+N4, 0, N);
311          mdct_backward(lookup, tmp, x, mode->window, overlap);
312          celt_assert(transient_shift == 0);
313          /* The first and last part would need to be set to zero if we actually
314             wanted to use them. */
315          for (j=0;j<overlap;j++)
316             out_mem[C*(MAX_PERIOD-N)+C*j+c] += x[j+N4];
317          for (j=0;j<overlap;j++)
318             out_mem[C*(MAX_PERIOD)+C*(overlap-j-1)+c] = x[2*N-j-N4-1];
319          for (j=0;j<2*N4;j++)
320             out_mem[C*(MAX_PERIOD-N)+C*(j+overlap)+c] = x[j+N4+overlap];
321          RESTORE_STACK;
322       } else {
323          int b;
324          const int N2 = mode->shortMdctSize;
325          const int B = mode->nbShortMdcts;
326          const mdct_lookup *lookup = &mode->shortMdct;
327          VARDECL(celt_word32_t, x);
328          VARDECL(celt_word32_t, tmp);
329          SAVE_STACK;
330          ALLOC(x, 2*N, celt_word32_t);
331          ALLOC(tmp, N, celt_word32_t);
332          /* Prevents problems from the imdct doing the overlap-add */
333          CELT_MEMSET(x+N4, 0, N2);
334          for (b=0;b<B;b++)
335          {
336             /* De-interleaving the sub-frames */
337             for (j=0;j<N2;j++)
338                tmp[j] = X[C*(j*B+b)+c];
339             mdct_backward(lookup, tmp, x+N4+N2*b, mode->window, overlap);
340          }
341          if (transient_shift > 0)
342          {
343 #ifdef FIXED_POINT
344             for (j=0;j<16;j++)
345                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));
346             for (j=transient_time;j<N+overlap;j++)
347                x[N4+j] = SHL32(x[N4+j], transient_shift);
348 #else
349             for (j=0;j<16;j++)
350                x[N4+transient_time+j-16] *= 1+transientWindow[j]*((1<<transient_shift)-1);
351             for (j=transient_time;j<N+overlap;j++)
352                x[N4+j] *= 1<<transient_shift;
353 #endif
354          }
355          /* The first and last part would need to be set to zero if we actually
356          wanted to use them. */
357          for (j=0;j<overlap;j++)
358             out_mem[C*(MAX_PERIOD-N)+C*j+c] += x[j+N4];
359          for (j=0;j<overlap;j++)
360             out_mem[C*(MAX_PERIOD)+C*(overlap-j-1)+c] = x[2*N-j-N4-1];
361          for (j=0;j<2*N4;j++)
362             out_mem[C*(MAX_PERIOD-N)+C*(j+overlap)+c] = x[j+N4+overlap];
363          RESTORE_STACK;
364       }
365    }
366 }
367
368 #define FLAG_NONE        0
369 #define FLAG_INTRA       1U<<16
370 #define FLAG_PITCH       1U<<15
371 #define FLAG_SHORT       1U<<14
372 #define FLAG_FOLD        1U<<13
373 #define FLAG_MASK        (FLAG_INTRA|FLAG_PITCH|FLAG_SHORT|FLAG_FOLD)
374
375 celt_int32_t flaglist[8] = {
376       0 /*00  */ | FLAG_FOLD,
377       1 /*01  */ | FLAG_PITCH|FLAG_FOLD,
378       8 /*1000*/ | FLAG_NONE,
379       9 /*1001*/ | FLAG_SHORT|FLAG_FOLD,
380      10 /*1010*/ | FLAG_PITCH,
381      11 /*1011*/ | FLAG_INTRA,
382       6 /*110 */ | FLAG_INTRA|FLAG_FOLD,
383       7 /*111 */ | FLAG_INTRA|FLAG_SHORT|FLAG_FOLD
384 };
385
386 void encode_flags(ec_enc *enc, int intra_ener, int has_pitch, int shortBlocks, int has_fold)
387 {
388    int i;
389    int flags=FLAG_NONE;
390    int flag_bits;
391    flags |= intra_ener   ? FLAG_INTRA : 0;
392    flags |= has_pitch    ? FLAG_PITCH : 0;
393    flags |= shortBlocks  ? FLAG_SHORT : 0;
394    flags |= has_fold     ? FLAG_FOLD  : 0;
395    for (i=0;i<8;i++)
396       if (flags == (flaglist[i]&FLAG_MASK))
397          break;
398    celt_assert(i<8);
399    flag_bits = flaglist[i]&0xf;
400    /*printf ("enc %d: %d %d %d %d\n", flag_bits, intra_ener, has_pitch, shortBlocks, has_fold);*/
401    if (i<2)
402       ec_enc_bits(enc, flag_bits, 2);
403    else if (i<6)
404       ec_enc_bits(enc, flag_bits, 4);
405    else
406       ec_enc_bits(enc, flag_bits, 3);
407 }
408
409 void decode_flags(ec_dec *dec, int *intra_ener, int *has_pitch, int *shortBlocks, int *has_fold)
410 {
411    int i;
412    int flag_bits;
413    flag_bits = ec_dec_bits(dec, 2);
414    /*printf ("(%d) ", flag_bits);*/
415    if (flag_bits==2)
416       flag_bits = (flag_bits<<2) | ec_dec_bits(dec, 2);
417    else if (flag_bits==3)
418       flag_bits = (flag_bits<<1) | ec_dec_bits(dec, 1);
419    for (i=0;i<8;i++)
420       if (flag_bits == (flaglist[i]&0xf))
421          break;
422    celt_assert(i<8);
423    *intra_ener  = (flaglist[i]&FLAG_INTRA) != 0;
424    *has_pitch   = (flaglist[i]&FLAG_PITCH) != 0;
425    *shortBlocks = (flaglist[i]&FLAG_SHORT) != 0;
426    *has_fold    = (flaglist[i]&FLAG_FOLD ) != 0;
427    /*printf ("dec %d: %d %d %d %d\n", flag_bits, *intra_ener, *has_pitch, *shortBlocks, *has_fold);*/
428 }
429
430 #ifdef FIXED_POINT
431 int celt_encode(CELTEncoder * restrict st, const celt_int16_t * pcm, celt_int16_t * optional_synthesis, unsigned char *compressed, int nbCompressedBytes)
432 {
433 #else
434 int celt_encode_float(CELTEncoder * restrict st, const celt_sig_t * pcm, celt_sig_t * optional_synthesis, unsigned char *compressed, int nbCompressedBytes)
435 {
436 #endif
437    int i, c, N, N4;
438    int has_pitch;
439    int pitch_index;
440    int bits;
441    int has_fold=1;
442    unsigned coarse_needed;
443    ec_byte_buffer buf;
444    ec_enc         enc;
445    VARDECL(celt_sig_t, in);
446    VARDECL(celt_sig_t, freq);
447    VARDECL(celt_norm_t, X);
448    VARDECL(celt_norm_t, P);
449    VARDECL(celt_ener_t, bandE);
450    VARDECL(celt_pgain_t, gains);
451    VARDECL(int, stereo_mode);
452    VARDECL(int, fine_quant);
453    VARDECL(celt_word16_t, error);
454    VARDECL(int, pulses);
455    VARDECL(int, offsets);
456 #ifdef EXP_PSY
457    VARDECL(celt_word32_t, mask);
458    VARDECL(celt_word32_t, tonality);
459    VARDECL(celt_word32_t, bandM);
460    VARDECL(celt_ener_t, bandN);
461 #endif
462    int intra_ener = 0;
463    int shortBlocks=0;
464    int transient_time;
465    int transient_shift;
466    const int C = CHANNELS(st->mode);
467    int mdct_weight_shift = 0;
468    int mdct_weight_pos=0;
469    SAVE_STACK;
470
471    if (check_mode(st->mode) != CELT_OK)
472       return CELT_INVALID_MODE;
473
474    if (nbCompressedBytes<0)
475      return CELT_BAD_ARG; 
476
477    /* The memset is important for now in case the encoder doesn't fill up all the bytes */
478    CELT_MEMSET(compressed, 0, nbCompressedBytes);
479    ec_byte_writeinit_buffer(&buf, compressed, nbCompressedBytes);
480    ec_enc_init(&enc,&buf);
481
482    N = st->block_size;
483    N4 = (N-st->overlap)>>1;
484    ALLOC(in, 2*C*N-2*C*N4, celt_sig_t);
485
486    CELT_COPY(in, st->in_mem, C*st->overlap);
487    for (c=0;c<C;c++)
488    {
489       const celt_word16_t * restrict pcmp = pcm+c;
490       celt_sig_t * restrict inp = in+C*st->overlap+c;
491       for (i=0;i<N;i++)
492       {
493          /* Apply pre-emphasis */
494          celt_sig_t tmp = SCALEIN(SHL32(EXTEND32(*pcmp), SIG_SHIFT));
495          *inp = SUB32(tmp, SHR32(MULT16_16(preemph,st->preemph_memE[c]),3));
496          st->preemph_memE[c] = SCALEIN(*pcmp);
497          inp += C;
498          pcmp += C;
499       }
500    }
501    CELT_COPY(st->in_mem, in+C*(2*N-2*N4-st->overlap), C*st->overlap);
502    
503    /* Transient handling */
504    if (st->mode->nbShortMdcts > 1)
505    {
506       if (transient_analysis(in, N+st->overlap, C, &transient_time, &transient_shift))
507       {
508 #ifndef FIXED_POINT
509          float gain_1;
510 #endif
511          /* Apply the inverse shaping window */
512          if (transient_shift)
513          {
514 #ifdef FIXED_POINT
515             for (c=0;c<C;c++)
516                for (i=0;i<16;i++)
517                   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]);
518             for (c=0;c<C;c++)
519                for (i=transient_time;i<N+st->overlap;i++)
520                   in[C*i+c] = SHR32(in[C*i+c], transient_shift);
521 #else
522             for (c=0;c<C;c++)
523                for (i=0;i<16;i++)
524                   in[C*(transient_time+i-16)+c] /= 1+transientWindow[i]*((1<<transient_shift)-1);
525             gain_1 = 1./(1<<transient_shift);
526             for (c=0;c<C;c++)
527                for (i=transient_time;i<N+st->overlap;i++)
528                   in[C*i+c] *= gain_1;
529 #endif
530          }
531          shortBlocks = 1;
532          has_fold = 1;
533       } else {
534          transient_time = -1;
535          transient_shift = 0;
536          shortBlocks = 0;
537       }
538    } else {
539       transient_time = -1;
540       transient_shift = 0;
541       shortBlocks = 0;
542    }
543
544    ALLOC(freq, C*N, celt_sig_t); /**< Interleaved signal MDCTs */
545    ALLOC(bandE,st->mode->nbEBands*C, celt_ener_t);
546    /* Compute MDCTs */
547    compute_mdcts(st->mode, shortBlocks, in, freq);
548    if (shortBlocks && !transient_shift) 
549    {
550       celt_word32_t sum[4]={1,1,1,1};
551       int m;
552       for (c=0;c<C;c++)
553       {
554          m=0;
555          do {
556             celt_word32_t tmp=0;
557             for (i=m*C+c;i<N;i+=C*st->mode->nbShortMdcts)
558                tmp += ABS32(freq[i]);
559             sum[m++] += tmp;
560          } while (m<st->mode->nbShortMdcts);
561       }
562       m=0;
563 #ifdef FIXED_POINT
564       do {
565          if (SHR32(sum[m+1],3) > sum[m])
566          {
567             mdct_weight_shift=2;
568             mdct_weight_pos = m;
569          } else if (SHR32(sum[m+1],1) > sum[m] && mdct_weight_shift < 2)
570          {
571             mdct_weight_shift=1;
572             mdct_weight_pos = m;
573          }
574          m++;
575       } while (m<st->mode->nbShortMdcts-1);
576       if (mdct_weight_shift)
577       {
578          for (c=0;c<C;c++)
579             for (m=mdct_weight_pos+1;m<st->mode->nbShortMdcts;m++)
580                for (i=m*C+c;i<N;i+=C*st->mode->nbShortMdcts)
581                   freq[i] = SHR32(freq[i],mdct_weight_shift);
582       }
583 #else
584       do {
585          if (sum[m+1] > 8*sum[m])
586          {
587             mdct_weight_shift=2;
588             mdct_weight_pos = m;
589          } else if (sum[m+1] > 2*sum[m] && mdct_weight_shift < 2)
590          {
591             mdct_weight_shift=1;
592             mdct_weight_pos = m;
593          }
594          m++;
595       } while (m<st->mode->nbShortMdcts-1);
596       if (mdct_weight_shift)
597       {
598          for (c=0;c<C;c++)
599             for (m=mdct_weight_pos+1;m<st->mode->nbShortMdcts;m++)
600                for (i=m*C+c;i<N;i+=C*st->mode->nbShortMdcts)
601                   freq[i] = (1./(1<<mdct_weight_shift))*freq[i];
602       }
603 #endif
604       /*printf ("%f\n", short_ratio);*/
605       /*if (short_ratio < 1)
606          short_ratio = 1;
607       short_ratio = 1<<(int)floor(.5+log2(short_ratio));
608       if (short_ratio>4)
609          short_ratio = 4;*/
610    }/* else if (transient_shift)
611       printf ("8\n");
612       else printf ("1\n");*/
613
614    compute_band_energies(st->mode, freq, bandE);
615
616    intra_ener = st->delayedIntra;
617    if (intra_decision(bandE, st->oldBandE, st->mode->nbEBands) || shortBlocks)
618       st->delayedIntra = 1;
619    else
620       st->delayedIntra = 0;
621    /* Pitch analysis: we do it early to save on the peak stack space */
622    /* Don't use pitch if there isn't enough data available yet, or if we're using shortBlocks */
623    has_pitch = st->pitch_enabled && (st->pitch_available >= MAX_PERIOD) && (!shortBlocks) && !intra_ener;
624 #ifdef EXP_PSY
625    ALLOC(tonality, MAX_PERIOD/4, celt_word16_t);
626    {
627       VARDECL(celt_word16_t, X);
628       ALLOC(X, MAX_PERIOD/2, celt_word16_t);
629       find_spectral_pitch(st->mode, st->mode->fft, &st->mode->psy, in, st->out_mem, st->mode->window, X, 2*N-2*N4, MAX_PERIOD-(2*N-2*N4), &pitch_index);
630       compute_tonality(st->mode, X, st->psy_mem, MAX_PERIOD, tonality, MAX_PERIOD/4);
631    }
632 #else
633    if (has_pitch)
634    {
635       find_spectral_pitch(st->mode, st->mode->fft, &st->mode->psy, in, st->out_mem, st->mode->window, NULL, 2*N-2*N4, MAX_PERIOD-(2*N-2*N4), &pitch_index);
636    }
637 #endif
638
639 #ifdef EXP_PSY
640    ALLOC(mask, N, celt_sig_t);
641    compute_mdct_masking(&st->psy, freq, tonality, st->psy_mem, mask, C*N);
642    /*for (i=0;i<256;i++)
643       printf ("%f %f %f ", freq[i], tonality[i], mask[i]);
644    printf ("\n");*/
645 #endif
646
647    /* Deferred allocation after find_spectral_pitch() to reduce the peak memory usage */
648    ALLOC(X, C*N, celt_norm_t);         /**< Interleaved normalised MDCTs */
649    ALLOC(P, C*N, celt_norm_t);         /**< Interleaved normalised pitch MDCTs*/
650    ALLOC(gains,st->mode->nbPBands, celt_pgain_t);
651
652
653    /* Band normalisation */
654    normalise_bands(st->mode, freq, X, bandE);
655
656 #ifdef EXP_PSY
657    ALLOC(bandN,C*st->mode->nbEBands, celt_ener_t);
658    ALLOC(bandM,st->mode->nbEBands, celt_ener_t);
659    compute_noise_energies(st->mode, freq, tonality, bandN);
660
661    /*for (i=0;i<st->mode->nbEBands;i++)
662       printf ("%f ", (.1+bandN[i])/(.1+bandE[i]));
663    printf ("\n");*/
664    has_fold = 0;
665    for (i=st->mode->nbPBands;i<st->mode->nbEBands;i++)
666       if (bandN[i] < .4*bandE[i])
667          has_fold++;
668    /*printf ("%d\n", has_fold);*/
669    if (has_fold>=2)
670       has_fold = 0;
671    else
672       has_fold = 1;
673    for (i=0;i<N;i++)
674       mask[i] = sqrt(mask[i]);
675    compute_band_energies(st->mode, mask, bandM);
676    /*for (i=0;i<st->mode->nbEBands;i++)
677       printf ("%f %f ", bandE[i], bandM[i]);
678    printf ("\n");*/
679 #endif
680
681    /* Compute MDCTs of the pitch part */
682    if (has_pitch)
683    {
684       celt_word32_t curr_power, pitch_power=0;
685       /* Normalise the pitch vector as well (discard the energies) */
686       VARDECL(celt_ener_t, bandEp);
687       
688       compute_mdcts(st->mode, 0, st->out_mem+pitch_index*C, freq);
689       ALLOC(bandEp, st->mode->nbEBands*st->mode->nbChannels, celt_ener_t);
690       compute_band_energies(st->mode, freq, bandEp);
691       normalise_bands(st->mode, freq, P, bandEp);
692       pitch_power = bandEp[0]+bandEp[1]+bandEp[2];
693       /* Check if we can safely use the pitch (i.e. effective gain isn't too high) */
694       curr_power = bandE[0]+bandE[1]+bandE[2];
695       if ((MULT16_32_Q15(QCONST16(.1f, 15),curr_power) + QCONST32(10.f,ENER_SHIFT) < pitch_power))
696       {
697          /* Pitch prediction */
698          has_pitch = compute_pitch_gain(st->mode, X, P, gains);
699       } else {
700          has_pitch = 0;
701       }
702    }
703    
704    encode_flags(&enc, intra_ener, has_pitch, shortBlocks, has_fold);
705    if (has_pitch)
706    {
707       ec_enc_uint(&enc, pitch_index, MAX_PERIOD-(2*N-2*N4));
708    } else {
709       for (i=0;i<st->mode->nbPBands;i++)
710          gains[i] = 0;
711       for (i=0;i<C*N;i++)
712          P[i] = 0;
713    }
714    if (shortBlocks)
715    {
716       if (transient_shift)
717       {
718          ec_enc_bits(&enc, transient_shift, 2);
719          ec_enc_uint(&enc, transient_time, N+st->overlap);
720       } else {
721          ec_enc_bits(&enc, mdct_weight_shift, 2);
722          if (mdct_weight_shift && st->mode->nbShortMdcts!=2)
723             ec_enc_uint(&enc, mdct_weight_pos, st->mode->nbShortMdcts-1);
724       }
725    }
726
727 #ifdef STDIN_TUNING2
728    static int fine_quant[30];
729    static int pulses[30];
730    static int init=0;
731    if (!init)
732    {
733       for (i=0;i<st->mode->nbEBands;i++)
734          scanf("%d ", &fine_quant[i]);
735       for (i=0;i<st->mode->nbEBands;i++)
736          scanf("%d ", &pulses[i]);
737       init = 1;
738    }
739 #else
740    ALLOC(fine_quant, st->mode->nbEBands, int);
741    ALLOC(pulses, st->mode->nbEBands, int);
742 #endif
743
744    /* Bit allocation */
745    ALLOC(error, C*st->mode->nbEBands, celt_word16_t);
746    coarse_needed = quant_coarse_energy(st->mode, bandE, st->oldBandE, nbCompressedBytes*8/3, intra_ener, st->mode->prob, error, &enc);
747    coarse_needed = ((coarse_needed*3-1)>>3)+1;
748
749    /* Variable bitrate */
750    if (st->VBR_rate>0)
751    {
752      /* The target rate in 16th bits per frame */
753      int target=st->VBR_rate;
754    
755      /* Shortblocks get a large boost in bitrate, but since they are uncommon long blocks are not greatly effected */
756      if (shortBlocks)
757        target*=2;
758      else if (st->mode->nbShortMdcts > 1)
759        target-=(target+14)/28;     
760
761      /*The average energy is removed from the target and the actual energy added*/
762      target=target-588+ec_enc_tell(&enc, 4);
763
764      /* In VBR mode the frame size must not be reduced so much that it would result in the coarse energy busting its budget */
765      target=IMAX(coarse_needed,(target+64)/128);
766      nbCompressedBytes=IMIN(nbCompressedBytes,target);
767    }
768
769    ALLOC(offsets, st->mode->nbEBands, int);
770    ALLOC(stereo_mode, st->mode->nbEBands, int);
771    stereo_decision(st->mode, X, stereo_mode, st->mode->nbEBands);
772
773    for (i=0;i<st->mode->nbEBands;i++)
774       offsets[i] = 0;
775    bits = nbCompressedBytes*8 - ec_enc_tell(&enc, 0) - 1;
776    if (has_pitch)
777       bits -= st->mode->nbPBands;
778 #ifndef STDIN_TUNING
779    compute_allocation(st->mode, offsets, stereo_mode, bits, pulses, fine_quant);
780 #endif
781
782    quant_fine_energy(st->mode, bandE, st->oldBandE, error, fine_quant, &enc);
783
784    /* Residual quantisation */
785    if (C==1)
786       quant_bands(st->mode, X, P, NULL, has_pitch, gains, bandE, pulses, shortBlocks, has_fold, nbCompressedBytes*8, &enc);
787    else
788       quant_bands_stereo(st->mode, X, P, NULL, has_pitch, gains, bandE, pulses, shortBlocks, has_fold, nbCompressedBytes*8, &enc);
789
790    /* Re-synthesis of the coded audio if required */
791    if (st->pitch_available>0 || optional_synthesis!=NULL)
792    {
793       if (st->pitch_available>0 && st->pitch_available<MAX_PERIOD)
794         st->pitch_available+=st->frame_size;
795
796       /* Synthesis */
797       denormalise_bands(st->mode, X, freq, bandE);
798       
799       
800       CELT_MOVE(st->out_mem, st->out_mem+C*N, C*(MAX_PERIOD+st->overlap-N));
801       
802       if (mdct_weight_shift)
803       {
804          int m;
805          for (c=0;c<C;c++)
806             for (m=mdct_weight_pos+1;m<st->mode->nbShortMdcts;m++)
807                for (i=m*C+c;i<N;i+=C*st->mode->nbShortMdcts)
808 #ifdef FIXED_POINT
809                   freq[i] = SHL32(freq[i], mdct_weight_shift);
810 #else
811                   freq[i] = (1<<mdct_weight_shift)*freq[i];
812 #endif
813       }
814       compute_inv_mdcts(st->mode, shortBlocks, freq, transient_time, transient_shift, st->out_mem);
815       /* De-emphasis and put everything back at the right place in the synthesis history */
816       if (optional_synthesis != NULL) {
817          for (c=0;c<C;c++)
818          {
819             int j;
820             for (j=0;j<N;j++)
821             {
822                celt_sig_t tmp = MAC16_32_Q15(st->out_mem[C*(MAX_PERIOD-N)+C*j+c],
823                                    preemph,st->preemph_memD[c]);
824                st->preemph_memD[c] = tmp;
825                optional_synthesis[C*j+c] = SCALEOUT(SIG2WORD16(tmp));
826             }
827          }
828       }
829    }
830
831    /*fprintf (stderr, "remaining bits after encode = %d\n", nbCompressedBytes*8-ec_enc_tell(&enc, 0));*/
832    /*if (ec_enc_tell(&enc, 0) < nbCompressedBytes*8 - 7)
833       celt_warning_int ("many unused bits: ", nbCompressedBytes*8-ec_enc_tell(&enc, 0));*/
834
835    /* Finishing the stream with a 0101... pattern so that the decoder can check is everything's right */
836    {
837       int val = 0;
838       while (ec_enc_tell(&enc, 0) < nbCompressedBytes*8)
839       {
840          ec_enc_uint(&enc, val, 2);
841          val = 1-val;
842       }
843    }
844    ec_enc_done(&enc);
845    {
846       /*unsigned char *data;*/
847       int nbBytes = ec_byte_bytes(&buf);
848       if (nbBytes > nbCompressedBytes)
849       {
850          celt_warning_int ("got too many bytes:", nbBytes);
851          RESTORE_STACK;
852          return CELT_INTERNAL_ERROR;
853       }
854    }
855
856    RESTORE_STACK;
857    return nbCompressedBytes;
858 }
859
860 #ifdef FIXED_POINT
861 #ifndef DISABLE_FLOAT_API
862 int celt_encode_float(CELTEncoder * restrict st, const float * pcm, float * optional_synthesis, unsigned char *compressed, int nbCompressedBytes)
863 {
864    int j, ret;
865    const int C = CHANNELS(st->mode);
866    const int N = st->block_size;
867    VARDECL(celt_int16_t, in);
868    SAVE_STACK;
869    ALLOC(in, C*N, celt_int16_t);
870
871    for (j=0;j<C*N;j++)
872      in[j] = FLOAT2INT16(pcm[j]);
873
874    if (optional_synthesis != NULL) {
875      ret=celt_encode(st,in,in,compressed,nbCompressedBytes);
876       for (j=0;j<C*N;j++)
877          optional_synthesis[j]=in[j]*(1/32768.);
878    } else {
879      ret=celt_encode(st,in,NULL,compressed,nbCompressedBytes);
880    }
881    RESTORE_STACK;
882    return ret;
883
884 }
885 #endif /*DISABLE_FLOAT_API*/
886 #else
887 int celt_encode(CELTEncoder * restrict st, const celt_int16_t * pcm, celt_int16_t * optional_synthesis, unsigned char *compressed, int nbCompressedBytes)
888 {
889    int j, ret;
890    VARDECL(celt_sig_t, in);
891    const int C = CHANNELS(st->mode);
892    const int N = st->block_size;
893    SAVE_STACK;
894    ALLOC(in, C*N, celt_sig_t);
895    for (j=0;j<C*N;j++) {
896      in[j] = SCALEOUT(pcm[j]);
897    }
898
899    if (optional_synthesis != NULL) {
900       ret = celt_encode_float(st,in,in,compressed,nbCompressedBytes);
901       for (j=0;j<C*N;j++)
902          optional_synthesis[j] = FLOAT2INT16(in[j]);
903    } else {
904       ret = celt_encode_float(st,in,NULL,compressed,nbCompressedBytes);
905    }
906    RESTORE_STACK;
907    return ret;
908 }
909 #endif
910
911 int celt_encoder_ctl(CELTEncoder * restrict st, int request, ...)
912 {
913    va_list ap;
914    va_start(ap, request);
915    switch (request)
916    {
917       case CELT_SET_COMPLEXITY_REQUEST:
918       {
919          int value = va_arg(ap, celt_int32_t);
920          if (value<0 || value>10)
921             goto bad_arg;
922          if (value<=2) {
923             st->pitch_enabled = 0; 
924             st->pitch_available = 0;
925          } else {
926               st->pitch_enabled = 1;
927               if (st->pitch_available<1)
928                 st->pitch_available = 1;
929          }   
930       }
931       break;
932       case CELT_SET_LTP_REQUEST:
933       {
934          int value = va_arg(ap, celt_int32_t);
935          if (value<0 || value>1 || (value==1 && st->pitch_available==0))
936             goto bad_arg;
937          if (value==0)
938             st->pitch_enabled = 0;
939          else
940             st->pitch_enabled = 1;
941       }
942       break;
943       case CELT_SET_VBR_RATE_REQUEST:
944       {
945          int value = va_arg(ap, celt_int32_t);
946          if (value<0)
947             goto bad_arg;
948          if (value>3072000)
949             value = 3072000;
950          st->VBR_rate = ((st->mode->Fs<<3)+(st->block_size>>1))/st->block_size;
951          st->VBR_rate = ((value<<7)+(st->VBR_rate>>1))/st->VBR_rate;
952       }
953       break;
954       case CELT_RESET_STATE:
955       {
956          const CELTMode *mode = st->mode;
957          int C = mode->nbChannels;
958
959          if (st->pitch_available > 0) st->pitch_available = 1;
960
961          CELT_MEMSET(st->in_mem, 0, st->overlap*C);
962          CELT_MEMSET(st->out_mem, 0, (MAX_PERIOD+st->overlap)*C);
963
964          CELT_MEMSET(st->oldBandE, 0, C*mode->nbEBands);
965
966          CELT_MEMSET(st->preemph_memE, 0, C);
967          CELT_MEMSET(st->preemph_memD, 0, C);
968          st->delayedIntra = 1;
969       }
970       break;
971       default:
972          goto bad_request;
973    }
974    va_end(ap);
975    return CELT_OK;
976 bad_arg:
977    va_end(ap);
978    return CELT_BAD_ARG;
979 bad_request:
980    va_end(ap);
981    return CELT_UNIMPLEMENTED;
982 }
983
984 /****************************************************************************/
985 /*                                                                          */
986 /*                                DECODER                                   */
987 /*                                                                          */
988 /****************************************************************************/
989 #ifdef NEW_PLC
990 #define DECODE_BUFFER_SIZE 2048
991 #else
992 #define DECODE_BUFFER_SIZE MAX_PERIOD
993 #endif
994
995 /** Decoder state 
996  @brief Decoder state
997  */
998 struct CELTDecoder {
999    const CELTMode *mode;
1000    int frame_size;
1001    int block_size;
1002    int overlap;
1003
1004    ec_byte_buffer buf;
1005    ec_enc         enc;
1006
1007    celt_sig_t * restrict preemph_memD;
1008
1009    celt_sig_t *out_mem;
1010    celt_sig_t *decode_mem;
1011
1012    celt_word16_t *oldBandE;
1013    
1014    int last_pitch_index;
1015 };
1016
1017 CELTDecoder *celt_decoder_create(const CELTMode *mode)
1018 {
1019    int N, C;
1020    CELTDecoder *st;
1021
1022    if (check_mode(mode) != CELT_OK)
1023       return NULL;
1024
1025    N = mode->mdctSize;
1026    C = CHANNELS(mode);
1027    st = celt_alloc(sizeof(CELTDecoder));
1028    
1029    st->mode = mode;
1030    st->frame_size = N;
1031    st->block_size = N;
1032    st->overlap = mode->overlap;
1033
1034    st->decode_mem = celt_alloc((DECODE_BUFFER_SIZE+st->overlap)*C*sizeof(celt_sig_t));
1035    st->out_mem = st->decode_mem+DECODE_BUFFER_SIZE-MAX_PERIOD;
1036    
1037    st->oldBandE = (celt_word16_t*)celt_alloc(C*mode->nbEBands*sizeof(celt_word16_t));
1038
1039    st->preemph_memD = (celt_sig_t*)celt_alloc(C*sizeof(celt_sig_t));
1040
1041    st->last_pitch_index = 0;
1042    return st;
1043 }
1044
1045 void celt_decoder_destroy(CELTDecoder *st)
1046 {
1047    if (st == NULL)
1048    {
1049       celt_warning("NULL passed to celt_encoder_destroy");
1050       return;
1051    }
1052    if (check_mode(st->mode) != CELT_OK)
1053       return;
1054
1055
1056    celt_free(st->decode_mem);
1057    
1058    celt_free(st->oldBandE);
1059    
1060    celt_free(st->preemph_memD);
1061
1062    celt_free(st);
1063 }
1064
1065 /** Handles lost packets by just copying past data with the same offset as the last
1066     pitch period */
1067 #ifdef NEW_PLC
1068 #include "plc.c"
1069 #else
1070 static void celt_decode_lost(CELTDecoder * restrict st, celt_word16_t * restrict pcm)
1071 {
1072    int c, N;
1073    int pitch_index;
1074    int i, len;
1075    VARDECL(celt_sig_t, freq);
1076    const int C = CHANNELS(st->mode);
1077    int offset;
1078    SAVE_STACK;
1079    N = st->block_size;
1080    ALLOC(freq,C*N, celt_sig_t);         /**< Interleaved signal MDCTs */
1081    
1082    len = N+st->mode->overlap;
1083 #if 0
1084    pitch_index = st->last_pitch_index;
1085    
1086    /* Use the pitch MDCT as the "guessed" signal */
1087    compute_mdcts(st->mode, st->mode->window, st->out_mem+pitch_index*C, freq);
1088
1089 #else
1090    find_spectral_pitch(st->mode, st->mode->fft, &st->mode->psy, st->out_mem+MAX_PERIOD-len, st->out_mem, st->mode->window, NULL, len, MAX_PERIOD-len-100, &pitch_index);
1091    pitch_index = MAX_PERIOD-len-pitch_index;
1092    offset = MAX_PERIOD-pitch_index;
1093    while (offset+len >= MAX_PERIOD)
1094       offset -= pitch_index;
1095    compute_mdcts(st->mode, 0, st->out_mem+offset*C, freq);
1096    for (i=0;i<N;i++)
1097       freq[i] = ADD32(EPSILON, MULT16_32_Q15(QCONST16(.9f,15),freq[i]));
1098 #endif
1099    
1100    
1101    
1102    CELT_MOVE(st->out_mem, st->out_mem+C*N, C*(MAX_PERIOD+st->mode->overlap-N));
1103    /* Compute inverse MDCTs */
1104    compute_inv_mdcts(st->mode, 0, freq, -1, 0, st->out_mem);
1105
1106    for (c=0;c<C;c++)
1107    {
1108       int j;
1109       for (j=0;j<N;j++)
1110       {
1111          celt_sig_t tmp = MAC16_32_Q15(st->out_mem[C*(MAX_PERIOD-N)+C*j+c],
1112                                 preemph,st->preemph_memD[c]);
1113          st->preemph_memD[c] = tmp;
1114          pcm[C*j+c] = SCALEOUT(SIG2WORD16(tmp));
1115       }
1116    }
1117    RESTORE_STACK;
1118 }
1119 #endif
1120
1121 #ifdef FIXED_POINT
1122 int celt_decode(CELTDecoder * restrict st, const unsigned char *data, int len, celt_int16_t * restrict pcm)
1123 {
1124 #else
1125 int celt_decode_float(CELTDecoder * restrict st, const unsigned char *data, int len, celt_sig_t * restrict pcm)
1126 {
1127 #endif
1128    int i, c, N, N4;
1129    int has_pitch, has_fold;
1130    int pitch_index;
1131    int bits;
1132    ec_dec dec;
1133    ec_byte_buffer buf;
1134    VARDECL(celt_sig_t, freq);
1135    VARDECL(celt_norm_t, X);
1136    VARDECL(celt_norm_t, P);
1137    VARDECL(celt_ener_t, bandE);
1138    VARDECL(celt_pgain_t, gains);
1139    VARDECL(int, stereo_mode);
1140    VARDECL(int, fine_quant);
1141    VARDECL(int, pulses);
1142    VARDECL(int, offsets);
1143
1144    int shortBlocks;
1145    int intra_ener;
1146    int transient_time;
1147    int transient_shift;
1148    int mdct_weight_shift=0;
1149    const int C = CHANNELS(st->mode);
1150    int mdct_weight_pos=0;
1151    SAVE_STACK;
1152
1153    if (check_mode(st->mode) != CELT_OK)
1154       return CELT_INVALID_MODE;
1155
1156    N = st->block_size;
1157    N4 = (N-st->overlap)>>1;
1158
1159    ALLOC(freq, C*N, celt_sig_t); /**< Interleaved signal MDCTs */
1160    ALLOC(X, C*N, celt_norm_t);         /**< Interleaved normalised MDCTs */
1161    ALLOC(P, C*N, celt_norm_t);         /**< Interleaved normalised pitch MDCTs*/
1162    ALLOC(bandE, st->mode->nbEBands*C, celt_ener_t);
1163    ALLOC(gains, st->mode->nbPBands, celt_pgain_t);
1164    
1165    if (check_mode(st->mode) != CELT_OK)
1166    {
1167       RESTORE_STACK;
1168       return CELT_INVALID_MODE;
1169    }
1170    if (data == NULL)
1171    {
1172       celt_decode_lost(st, pcm);
1173       RESTORE_STACK;
1174       return 0;
1175    }
1176    if (len<0) {
1177      RESTORE_STACK;
1178      return CELT_BAD_ARG;
1179    }
1180    
1181    ec_byte_readinit(&buf,(unsigned char*)data,len);
1182    ec_dec_init(&dec,&buf);
1183    
1184    decode_flags(&dec, &intra_ener, &has_pitch, &shortBlocks, &has_fold);
1185    if (shortBlocks)
1186    {
1187       transient_shift = ec_dec_bits(&dec, 2);
1188       if (transient_shift == 3)
1189       {
1190          transient_time = ec_dec_uint(&dec, N+st->mode->overlap);
1191       } else {
1192          mdct_weight_shift = transient_shift;
1193          if (mdct_weight_shift && st->mode->nbShortMdcts>2)
1194             mdct_weight_pos = ec_dec_uint(&dec, st->mode->nbShortMdcts-1);
1195          transient_shift = 0;
1196          transient_time = 0;
1197       }
1198    } else {
1199       transient_time = -1;
1200       transient_shift = 0;
1201    }
1202    
1203    if (has_pitch)
1204    {
1205       pitch_index = ec_dec_uint(&dec, MAX_PERIOD-(2*N-2*N4));
1206       st->last_pitch_index = pitch_index;
1207    } else {
1208       pitch_index = 0;
1209       for (i=0;i<st->mode->nbPBands;i++)
1210          gains[i] = 0;
1211    }
1212
1213    ALLOC(fine_quant, st->mode->nbEBands, int);
1214    /* Get band energies */
1215    unquant_coarse_energy(st->mode, bandE, st->oldBandE, len*8/3, intra_ener, st->mode->prob, &dec);
1216    
1217    ALLOC(pulses, st->mode->nbEBands, int);
1218    ALLOC(offsets, st->mode->nbEBands, int);
1219    ALLOC(stereo_mode, st->mode->nbEBands, int);
1220    stereo_decision(st->mode, X, stereo_mode, st->mode->nbEBands);
1221
1222    for (i=0;i<st->mode->nbEBands;i++)
1223       offsets[i] = 0;
1224
1225    bits = len*8 - ec_dec_tell(&dec, 0) - 1;
1226    if (has_pitch)
1227       bits -= st->mode->nbPBands;
1228    compute_allocation(st->mode, offsets, stereo_mode, bits, pulses, fine_quant);
1229    /*bits = ec_dec_tell(&dec, 0);
1230    compute_fine_allocation(st->mode, fine_quant, (20*C+len*8/5-(ec_dec_tell(&dec, 0)-bits))/C);*/
1231    
1232    unquant_fine_energy(st->mode, bandE, st->oldBandE, fine_quant, &dec);
1233
1234
1235    if (has_pitch) 
1236    {
1237       VARDECL(celt_ener_t, bandEp);
1238       
1239       /* Pitch MDCT */
1240       compute_mdcts(st->mode, 0, st->out_mem+pitch_index*C, freq);
1241       ALLOC(bandEp, st->mode->nbEBands*C, celt_ener_t);
1242       compute_band_energies(st->mode, freq, bandEp);
1243       normalise_bands(st->mode, freq, P, bandEp);
1244       /* Apply pitch gains */
1245    } else {
1246       for (i=0;i<C*N;i++)
1247          P[i] = 0;
1248    }
1249
1250    /* Decode fixed codebook and merge with pitch */
1251    if (C==1)
1252       unquant_bands(st->mode, X, P, has_pitch, gains, bandE, pulses, shortBlocks, has_fold, len*8, &dec);
1253    else
1254       unquant_bands_stereo(st->mode, X, P, has_pitch, gains, bandE, pulses, shortBlocks, has_fold, len*8, &dec);
1255
1256    /* Synthesis */
1257    denormalise_bands(st->mode, X, freq, bandE);
1258
1259
1260    CELT_MOVE(st->decode_mem, st->decode_mem+C*N, C*(DECODE_BUFFER_SIZE+st->overlap-N));
1261    if (mdct_weight_shift)
1262    {
1263       int m;
1264       for (c=0;c<C;c++)
1265          for (m=mdct_weight_pos+1;m<st->mode->nbShortMdcts;m++)
1266             for (i=m*C+c;i<N;i+=C*st->mode->nbShortMdcts)
1267 #ifdef FIXED_POINT
1268                freq[i] = SHL32(freq[i], mdct_weight_shift);
1269 #else
1270                freq[i] = (1<<mdct_weight_shift)*freq[i];
1271 #endif
1272    }
1273    /* Compute inverse MDCTs */
1274    compute_inv_mdcts(st->mode, shortBlocks, freq, transient_time, transient_shift, st->out_mem);
1275
1276    for (c=0;c<C;c++)
1277    {
1278       int j;
1279       for (j=0;j<N;j++)
1280       {
1281          celt_sig_t tmp = MAC16_32_Q15(st->out_mem[C*(MAX_PERIOD-N)+C*j+c],
1282                                 preemph,st->preemph_memD[c]);
1283          st->preemph_memD[c] = tmp;
1284          pcm[C*j+c] = SCALEOUT(SIG2WORD16(tmp));
1285       }
1286    }
1287
1288    {
1289       unsigned int val = 0;
1290       while (ec_dec_tell(&dec, 0) < len*8)
1291       {
1292          if (ec_dec_uint(&dec, 2) != val)
1293          {
1294             celt_warning("decode error");
1295             RESTORE_STACK;
1296             return CELT_CORRUPTED_DATA;
1297          }
1298          val = 1-val;
1299       }
1300    }
1301
1302    RESTORE_STACK;
1303    return 0;
1304    /*printf ("\n");*/
1305 }
1306
1307 #ifdef FIXED_POINT
1308 #ifndef DISABLE_FLOAT_API
1309 int celt_decode_float(CELTDecoder * restrict st, const unsigned char *data, int len, float * restrict pcm)
1310 {
1311    int j, ret;
1312    const int C = CHANNELS(st->mode);
1313    const int N = st->block_size;
1314    VARDECL(celt_int16_t, out);
1315    SAVE_STACK;
1316    ALLOC(out, C*N, celt_int16_t);
1317
1318    ret=celt_decode(st, data, len, out);
1319
1320    for (j=0;j<C*N;j++)
1321      pcm[j]=out[j]*(1/32768.);
1322    RESTORE_STACK;
1323    return ret;
1324 }
1325 #endif /*DISABLE_FLOAT_API*/
1326 #else
1327 int celt_decode(CELTDecoder * restrict st, const unsigned char *data, int len, celt_int16_t * restrict pcm)
1328 {
1329    int j, ret;
1330    VARDECL(celt_sig_t, out);
1331    const int C = CHANNELS(st->mode);
1332    const int N = st->block_size;
1333    SAVE_STACK;
1334    ALLOC(out, C*N, celt_sig_t);
1335
1336    ret=celt_decode_float(st, data, len, out);
1337
1338    for (j=0;j<C*N;j++)
1339      pcm[j] = FLOAT2INT16 (out[j]);
1340
1341    RESTORE_STACK;
1342    return ret;
1343 }
1344 #endif
1345
1346 int celt_decoder_ctl(CELTDecoder * restrict st, int request, ...)
1347 {
1348    va_list ap;
1349    va_start(ap, request);
1350    switch (request)
1351    {
1352       case CELT_RESET_STATE:
1353       {
1354          const CELTMode *mode = st->mode;
1355          int C = mode->nbChannels;
1356
1357          CELT_MEMSET(st->decode_mem, 0, (DECODE_BUFFER_SIZE+st->overlap)*C);
1358          CELT_MEMSET(st->oldBandE, 0, C*mode->nbEBands);
1359
1360          CELT_MEMSET(st->preemph_memD, 0, C);
1361
1362          st->last_pitch_index = 0;
1363       }
1364       break;
1365       default:
1366          goto bad_request;
1367    }
1368    va_end(ap);
1369    return CELT_OK;
1370 #if 0    /* Put this back in if you ever need "bad_arg" */
1371 bad_arg:
1372    va_end(ap);
1373    return CELT_BAD_ARG;
1374 #endif
1375 bad_request:
1376       va_end(ap);
1377   return CELT_UNIMPLEMENTED;
1378 }