535703cc988c14b4ad2944188fad352bab0fada4
[opus.git] / libcelt / celt.c
1 /* (C) 2007-2008 Jean-Marc Valin, CSIRO
2 */
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_pitch.h"
48 #include "quant_bands.h"
49 #include "psy.h"
50 #include "rate.h"
51 #include "stack_alloc.h"
52
53 static const celt_word16_t preemph = QCONST16(0.8f,15);
54
55
56 /** Encoder state 
57  @brief Encoder state
58  */
59 struct CELTEncoder {
60    const CELTMode *mode;     /**< Mode used by the encoder */
61    int frame_size;
62    int block_size;
63    int overlap;
64    int channels;
65    
66    ec_byte_buffer buf;
67    ec_enc         enc;
68
69    celt_word16_t * restrict preemph_memE; /* Input is 16-bit, so why bother with 32 */
70    celt_sig_t    * restrict preemph_memD;
71
72    celt_sig_t *in_mem;
73    celt_sig_t *out_mem;
74
75    celt_word16_t *oldBandE;
76 #ifdef EXP_PSY
77    celt_word16_t *psy_mem;
78    struct PsyDecay psy;
79 #endif
80 };
81
82 CELTEncoder EXPORT *celt_encoder_create(const CELTMode *mode)
83 {
84    int N, C;
85    CELTEncoder *st;
86
87    if (check_mode(mode) != CELT_OK)
88       return NULL;
89
90    N = mode->mdctSize;
91    C = mode->nbChannels;
92    st = celt_alloc(sizeof(CELTEncoder));
93    
94    st->mode = mode;
95    st->frame_size = N;
96    st->block_size = N;
97    st->overlap = mode->overlap;
98
99    ec_byte_writeinit(&st->buf);
100    ec_enc_init(&st->enc,&st->buf);
101
102    st->in_mem = celt_alloc(st->overlap*C*sizeof(celt_sig_t));
103    st->out_mem = celt_alloc((MAX_PERIOD+st->overlap)*C*sizeof(celt_sig_t));
104
105    st->oldBandE = (celt_word16_t*)celt_alloc(C*mode->nbEBands*sizeof(celt_word16_t));
106
107    st->preemph_memE = (celt_word16_t*)celt_alloc(C*sizeof(celt_word16_t));;
108    st->preemph_memD = (celt_sig_t*)celt_alloc(C*sizeof(celt_sig_t));;
109
110 #ifdef EXP_PSY
111    st->psy_mem = celt_alloc(MAX_PERIOD*sizeof(celt_word16_t));
112    psydecay_init(&st->psy, MAX_PERIOD/2, st->mode->Fs);
113 #endif
114
115    return st;
116 }
117
118 void EXPORT celt_encoder_destroy(CELTEncoder *st)
119 {
120    if (st == NULL)
121    {
122       celt_warning("NULL passed to celt_encoder_destroy");
123       return;
124    }
125    if (check_mode(st->mode) != CELT_OK)
126       return;
127
128    ec_byte_writeclear(&st->buf);
129
130    celt_free(st->in_mem);
131    celt_free(st->out_mem);
132    
133    celt_free(st->oldBandE);
134    
135    celt_free(st->preemph_memE);
136    celt_free(st->preemph_memD);
137    
138 #ifdef EXP_PSY
139    celt_free (st->psy_mem);
140    psydecay_clear(&st->psy);
141 #endif
142    
143    celt_free(st);
144 }
145
146 static inline celt_int16_t SIG2INT16(celt_sig_t x)
147 {
148    x = PSHR32(x, SIG_SHIFT);
149    x = MAX32(x, -32768);
150    x = MIN32(x, 32767);
151 #ifdef FIXED_POINT
152    return EXTRACT16(x);
153 #else
154    return (celt_int16_t)floor(.5+x);
155 #endif
156 }
157
158 /** Apply window and compute the MDCT for all sub-frames and all channels in a frame */
159 static void compute_mdcts(const CELTMode *mode, const celt_word16_t * restrict window, celt_sig_t * restrict in, celt_sig_t * restrict out)
160 {
161    const mdct_lookup *lookup = MDCT(mode);
162    const int N = FRAMESIZE(mode);
163    const int C = CHANNELS(mode);
164    const int overlap = OVERLAP(mode);
165    if (C==1)
166    {
167       mdct_forward(lookup, in, out, window, overlap);
168    } else {
169       int c;
170       VARDECL(celt_word32_t, x);
171       VARDECL(celt_word32_t, tmp);
172       SAVE_STACK;
173       ALLOC(x, N+overlap, celt_word32_t);
174       ALLOC(tmp, N, celt_word32_t);
175       for (c=0;c<C;c++)
176       {
177          int j;
178          for (j=0;j<N+overlap;j++)
179             x[j] = in[C*j+c];
180          mdct_forward(lookup, x, tmp, window, overlap);
181          /* Interleaving the sub-frames */
182          for (j=0;j<N;j++)
183             out[C*j+c] = tmp[j];
184       }
185       RESTORE_STACK;
186    }
187 }
188
189 /** Compute the IMDCT and apply window for all sub-frames and all channels in a frame */
190 static void compute_inv_mdcts(const CELTMode *mode, const celt_word16_t * restrict window, celt_sig_t *X, celt_sig_t * restrict out_mem)
191 {
192    int c, N4;
193    const int C = CHANNELS(mode);
194    const mdct_lookup *lookup = MDCT(mode);
195    const int N = FRAMESIZE(mode);
196    const int overlap = OVERLAP(mode);
197    N4 = (N-overlap)>>1;
198    for (c=0;c<C;c++)
199    {
200       int j;
201       if (C==1) {
202          mdct_backward(lookup, X, out_mem+C*(MAX_PERIOD-N-N4), window, overlap);
203       } else {
204          VARDECL(celt_word32_t, x);
205          VARDECL(celt_word32_t, tmp);
206          SAVE_STACK;
207          ALLOC(x, 2*N, celt_word32_t);
208          ALLOC(tmp, N, celt_word32_t);
209          /* De-interleaving the sub-frames */
210          for (j=0;j<N;j++)
211             tmp[j] = X[C*j+c];
212          /* Prevents problems from the imdct doing the overlap-add */
213          CELT_MEMSET(x+N4, 0, overlap);
214          mdct_backward(lookup, tmp, x, window, overlap);
215          /* The first and last part would need to be set to zero if we actually
216          wanted to use them. */
217          for (j=0;j<overlap;j++)
218             out_mem[C*(MAX_PERIOD-N)+C*j+c] += x[j+N4];
219          for (j=0;j<overlap;j++)
220             out_mem[C*(MAX_PERIOD)+C*(overlap-j-1)+c] = x[2*N-j-N4-1];
221          for (j=0;j<2*N4;j++)
222             out_mem[C*(MAX_PERIOD-N)+C*(j+overlap)+c] = x[j+N4+overlap];
223          RESTORE_STACK;
224       }
225    }
226 }
227
228 int EXPORT celt_encode(CELTEncoder * restrict st, celt_int16_t * restrict pcm, unsigned char *compressed, int nbCompressedBytes)
229 {
230    int i, c, N, N4;
231    int has_pitch;
232    int pitch_index;
233    celt_word32_t curr_power, pitch_power;
234    VARDECL(celt_sig_t, in);
235    VARDECL(celt_sig_t, freq);
236    VARDECL(celt_norm_t, X);
237    VARDECL(celt_norm_t, P);
238    VARDECL(celt_ener_t, bandE);
239    VARDECL(celt_pgain_t, gains);
240    VARDECL(int, stereo_mode);
241 #ifdef EXP_PSY
242    VARDECL(celt_word32_t, mask);
243 #endif
244    const int C = CHANNELS(st->mode);
245    SAVE_STACK;
246
247    if (check_mode(st->mode) != CELT_OK)
248       return CELT_INVALID_MODE;
249
250    N = st->block_size;
251    N4 = (N-st->overlap)>>1;
252    ALLOC(in, 2*C*N-2*C*N4, celt_sig_t);
253
254    CELT_COPY(in, st->in_mem, C*st->overlap);
255    for (c=0;c<C;c++)
256    {
257       const celt_int16_t * restrict pcmp = pcm+c;
258       celt_sig_t * restrict inp = in+C*st->overlap+c;
259       for (i=0;i<N;i++)
260       {
261          /* Apply pre-emphasis */
262          celt_sig_t tmp = SHL32(EXTEND32(*pcmp), SIG_SHIFT);
263          *inp = SUB32(tmp, SHR32(MULT16_16(preemph,st->preemph_memE[c]),1));
264          st->preemph_memE[c] = *pcmp;
265          inp += C;
266          pcmp += C;
267       }
268    }
269    CELT_COPY(st->in_mem, in+C*(2*N-2*N4-st->overlap), C*st->overlap);
270
271    /* Pitch analysis: we do it early to save on the peak stack space */
272    find_spectral_pitch(st->mode, st->mode->fft, &st->mode->psy, in, st->out_mem, st->mode->window, 2*N-2*N4, MAX_PERIOD-(2*N-2*N4), &pitch_index);
273
274    ALLOC(freq, C*N, celt_sig_t); /**< Interleaved signal MDCTs */
275    
276    /*for (i=0;i<(B+1)*C*N;i++) printf ("%f(%d) ", in[i], i); printf ("\n");*/
277    /* Compute MDCTs */
278    compute_mdcts(st->mode, st->mode->window, in, freq);
279
280 #ifdef EXP_PSY
281    CELT_MOVE(st->psy_mem, st->out_mem+N, MAX_PERIOD+st->overlap-N);
282    for (i=0;i<N;i++)
283       st->psy_mem[MAX_PERIOD+st->overlap-N+i] = in[C*(st->overlap+i)];
284    for (c=1;c<C;c++)
285       for (i=0;i<N;i++)
286          st->psy_mem[MAX_PERIOD+st->overlap-N+i] += in[C*(st->overlap+i)+c];
287
288    ALLOC(mask, N, celt_sig_t);
289    compute_mdct_masking(&st->psy, freq, st->psy_mem, mask, C*N);
290
291    /* Invert and stretch the mask to length of X 
292       For some reason, I get better results by using the sqrt instead,
293       although there's no valid reason to. Must investigate further */
294    for (i=0;i<C*N;i++)
295       mask[i] = 1/(.1+mask[i]);
296 #endif
297    
298    /* Deferred allocation after find_spectral_pitch() to reduce the peak memory usage */
299    ALLOC(X, C*N, celt_norm_t);         /**< Interleaved normalised MDCTs */
300    ALLOC(P, C*N, celt_norm_t);         /**< Interleaved normalised pitch MDCTs*/
301    ALLOC(bandE,st->mode->nbEBands*C, celt_ener_t);
302    ALLOC(gains,st->mode->nbPBands, celt_pgain_t);
303
304    /*printf ("%f %f\n", curr_power, pitch_power);*/
305    /*int j;
306    for (j=0;j<B*N;j++)
307       printf ("%f ", X[j]);
308    for (j=0;j<B*N;j++)
309       printf ("%f ", P[j]);
310    printf ("\n");*/
311
312    /* Band normalisation */
313    compute_band_energies(st->mode, freq, bandE);
314    normalise_bands(st->mode, freq, X, bandE);
315    /*for (i=0;i<st->mode->nbEBands;i++)printf("%f ", bandE[i]);printf("\n");*/
316    /*for (i=0;i<N*B*C;i++)printf("%f ", X[i]);printf("\n");*/
317
318    /* Compute MDCTs of the pitch part */
319    compute_mdcts(st->mode, st->mode->window, st->out_mem+pitch_index*C, freq);
320
321    {
322       /* Normalise the pitch vector as well (discard the energies) */
323       VARDECL(celt_ener_t, bandEp);
324       ALLOC(bandEp, st->mode->nbEBands*st->mode->nbChannels, celt_ener_t);
325       compute_band_energies(st->mode, freq, bandEp);
326       normalise_bands(st->mode, freq, P, bandEp);
327       pitch_power = bandEp[0]+bandEp[1]+bandEp[2];
328    }
329    curr_power = bandE[0]+bandE[1]+bandE[2];
330    /* Check if we can safely use the pitch (i.e. effective gain isn't too high) */
331    if (MULT16_32_Q15(QCONST16(.1f, 15),curr_power) + QCONST32(10.f,ENER_SHIFT) < pitch_power)
332    {
333       /* Simulates intensity stereo */
334       /*for (i=30;i<N*B;i++)
335          X[i*C+1] = P[i*C+1] = 0;*/
336
337       /* Pitch prediction */
338       compute_pitch_gain(st->mode, X, P, gains);
339       has_pitch = quant_pitch(gains, st->mode->nbPBands, &st->enc);
340       if (has_pitch)
341          ec_enc_uint(&st->enc, pitch_index, MAX_PERIOD-(2*N-2*N4));
342    } else {
343       /* No pitch, so we just pretend we found a gain of zero */
344       for (i=0;i<st->mode->nbPBands;i++)
345          gains[i] = 0;
346       ec_enc_bits(&st->enc, 0, 7);
347       for (i=0;i<C*N;i++)
348          P[i] = 0;
349    }
350    quant_energy(st->mode, bandE, st->oldBandE, 20*C+nbCompressedBytes*8/5, st->mode->prob, &st->enc);
351
352    ALLOC(stereo_mode, st->mode->nbEBands, int);
353    stereo_decision(st->mode, X, stereo_mode, st->mode->nbEBands);
354
355    pitch_quant_bands(st->mode, P, gains);
356
357    /*for (i=0;i<B*N;i++) printf("%f ",P[i]);printf("\n");*/
358
359    /* Residual quantisation */
360    quant_bands(st->mode, X, P, NULL, bandE, stereo_mode, nbCompressedBytes*8, &st->enc);
361    
362    if (C==2)
363    {
364       renormalise_bands(st->mode, X);
365    }
366    /* Synthesis */
367    denormalise_bands(st->mode, X, freq, bandE);
368
369
370    CELT_MOVE(st->out_mem, st->out_mem+C*N, C*(MAX_PERIOD+st->overlap-N));
371
372    compute_inv_mdcts(st->mode, st->mode->window, freq, st->out_mem);
373    /* De-emphasis and put everything back at the right place in the synthesis history */
374 #ifndef SHORTCUTS
375    for (c=0;c<C;c++)
376    {
377       int j;
378       const celt_sig_t * restrict outp=st->out_mem+C*(MAX_PERIOD-N)+c;
379       celt_int16_t * restrict pcmp = pcm+c;
380       for (j=0;j<N;j++)
381       {
382          celt_sig_t tmp = ADD32(*outp, MULT16_32_Q15(preemph,st->preemph_memD[c]));
383          st->preemph_memD[c] = tmp;
384          *pcmp = SIG2INT16(tmp);
385          pcmp += C;
386          outp += C;
387       }
388    }
389 #endif
390    if (ec_enc_tell(&st->enc, 0) < nbCompressedBytes*8 - 7)
391       celt_warning_int ("many unused bits: ", nbCompressedBytes*8-ec_enc_tell(&st->enc, 0));
392    /*printf ("%d\n", ec_enc_tell(&st->enc, 0)-8*nbCompressedBytes);*/
393    /* Finishing the stream with a 0101... pattern so that the decoder can check is everything's right */
394    {
395       int val = 0;
396       while (ec_enc_tell(&st->enc, 0) < nbCompressedBytes*8)
397       {
398          ec_enc_uint(&st->enc, val, 2);
399          val = 1-val;
400       }
401    }
402    ec_enc_done(&st->enc);
403    {
404       unsigned char *data;
405       int nbBytes = ec_byte_bytes(&st->buf);
406       if (nbBytes > nbCompressedBytes)
407       {
408          celt_warning_int ("got too many bytes:", nbBytes);
409          RESTORE_STACK;
410          return CELT_INTERNAL_ERROR;
411       }
412       /*printf ("%d\n", *nbBytes);*/
413       data = ec_byte_get_buffer(&st->buf);
414       for (i=0;i<nbBytes;i++)
415          compressed[i] = data[i];
416       for (;i<nbCompressedBytes;i++)
417          compressed[i] = 0;
418    }
419    /* Reset the packing for the next encoding */
420    ec_byte_reset(&st->buf);
421    ec_enc_init(&st->enc,&st->buf);
422
423    RESTORE_STACK;
424    return nbCompressedBytes;
425 }
426
427
428 /****************************************************************************/
429 /*                                                                          */
430 /*                                DECODER                                   */
431 /*                                                                          */
432 /****************************************************************************/
433
434
435 /** Decoder state 
436  @brief Decoder state
437  */
438 struct CELTDecoder {
439    const CELTMode *mode;
440    int frame_size;
441    int block_size;
442    int overlap;
443
444    ec_byte_buffer buf;
445    ec_enc         enc;
446
447    celt_sig_t * restrict preemph_memD;
448
449    celt_sig_t *out_mem;
450
451    celt_word16_t *oldBandE;
452    
453    int last_pitch_index;
454 };
455
456 CELTDecoder EXPORT *celt_decoder_create(const CELTMode *mode)
457 {
458    int N, C;
459    CELTDecoder *st;
460
461    if (check_mode(mode) != CELT_OK)
462       return NULL;
463
464    N = mode->mdctSize;
465    C = CHANNELS(mode);
466    st = celt_alloc(sizeof(CELTDecoder));
467    
468    st->mode = mode;
469    st->frame_size = N;
470    st->block_size = N;
471    st->overlap = mode->overlap;
472
473    st->out_mem = celt_alloc((MAX_PERIOD+st->overlap)*C*sizeof(celt_sig_t));
474    
475    st->oldBandE = (celt_word16_t*)celt_alloc(C*mode->nbEBands*sizeof(celt_word16_t));
476
477    st->preemph_memD = (celt_sig_t*)celt_alloc(C*sizeof(celt_sig_t));;
478
479    st->last_pitch_index = 0;
480    return st;
481 }
482
483 void EXPORT celt_decoder_destroy(CELTDecoder *st)
484 {
485    if (st == NULL)
486    {
487       celt_warning("NULL passed to celt_encoder_destroy");
488       return;
489    }
490    if (check_mode(st->mode) != CELT_OK)
491       return;
492
493
494    celt_free(st->out_mem);
495    
496    celt_free(st->oldBandE);
497    
498    celt_free(st->preemph_memD);
499
500    celt_free(st);
501 }
502
503 /** Handles lost packets by just copying past data with the same offset as the last
504     pitch period */
505 static void celt_decode_lost(CELTDecoder * restrict st, short * restrict pcm)
506 {
507    int c, N;
508    int pitch_index;
509    int i, len;
510    VARDECL(celt_sig_t, freq);
511    const int C = CHANNELS(st->mode);
512    int offset;
513    SAVE_STACK;
514    N = st->block_size;
515    ALLOC(freq,C*N, celt_sig_t);         /**< Interleaved signal MDCTs */
516    
517    len = N+st->mode->overlap;
518 #if 0
519    pitch_index = st->last_pitch_index;
520    
521    /* Use the pitch MDCT as the "guessed" signal */
522    compute_mdcts(st->mode, st->mode->window, st->out_mem+pitch_index*C, freq);
523
524 #else
525    find_spectral_pitch(st->mode, st->mode->fft, &st->mode->psy, st->out_mem+MAX_PERIOD-len, st->out_mem, st->mode->window, len, MAX_PERIOD-len-100, &pitch_index);
526    pitch_index = MAX_PERIOD-len-pitch_index;
527    offset = MAX_PERIOD-pitch_index;
528    while (offset+len >= MAX_PERIOD)
529       offset -= pitch_index;
530    compute_mdcts(st->mode, st->mode->window, st->out_mem+offset*C, freq);
531    for (i=0;i<N;i++)
532       freq[i] = MULT16_32_Q15(QCONST16(.9f,15),freq[i]);
533 #endif
534    
535    
536    
537    CELT_MOVE(st->out_mem, st->out_mem+C*N, C*(MAX_PERIOD+st->mode->overlap-N));
538    /* Compute inverse MDCTs */
539    compute_inv_mdcts(st->mode, st->mode->window, freq, st->out_mem);
540
541    for (c=0;c<C;c++)
542    {
543       int j;
544       for (j=0;j<N;j++)
545       {
546          celt_sig_t tmp = ADD32(st->out_mem[C*(MAX_PERIOD-N)+C*j+c],
547                                 MULT16_32_Q15(preemph,st->preemph_memD[c]));
548          st->preemph_memD[c] = tmp;
549          pcm[C*j+c] = SIG2INT16(tmp);
550       }
551    }
552    RESTORE_STACK;
553 }
554
555 int EXPORT celt_decode(CELTDecoder * restrict st, unsigned char *data, int len, celt_int16_t * restrict pcm)
556 {
557    int c, N, N4;
558    int has_pitch;
559    int pitch_index;
560    ec_dec dec;
561    ec_byte_buffer buf;
562    VARDECL(celt_sig_t, freq);
563    VARDECL(celt_norm_t, X);
564    VARDECL(celt_norm_t, P);
565    VARDECL(celt_ener_t, bandE);
566    VARDECL(celt_pgain_t, gains);
567    VARDECL(int, stereo_mode);
568    const int C = CHANNELS(st->mode);
569    SAVE_STACK;
570
571    if (check_mode(st->mode) != CELT_OK)
572       return CELT_INVALID_MODE;
573
574    N = st->block_size;
575    N4 = (N-st->overlap)>>1;
576
577    ALLOC(freq, C*N, celt_sig_t); /**< Interleaved signal MDCTs */
578    ALLOC(X, C*N, celt_norm_t);         /**< Interleaved normalised MDCTs */
579    ALLOC(P, C*N, celt_norm_t);         /**< Interleaved normalised pitch MDCTs*/
580    ALLOC(bandE, st->mode->nbEBands*C, celt_ener_t);
581    ALLOC(gains, st->mode->nbPBands, celt_pgain_t);
582    
583    if (check_mode(st->mode) != CELT_OK)
584    {
585       RESTORE_STACK;
586       return CELT_INVALID_MODE;
587    }
588    if (data == NULL)
589    {
590       celt_decode_lost(st, pcm);
591       RESTORE_STACK;
592       return 0;
593    }
594    
595    ec_byte_readinit(&buf,data,len);
596    ec_dec_init(&dec,&buf);
597    
598    /* Get the pitch gains */
599    has_pitch = unquant_pitch(gains, st->mode->nbPBands, &dec);
600    
601    /* Get the pitch index */
602    if (has_pitch)
603    {
604       pitch_index = ec_dec_uint(&dec, MAX_PERIOD-(2*N-2*N4));
605       st->last_pitch_index = pitch_index;
606    } else {
607       /* FIXME: We could be more intelligent here and just not compute the MDCT */
608       pitch_index = 0;
609    }
610
611    /* Get band energies */
612    unquant_energy(st->mode, bandE, st->oldBandE, 20*C+len*8/5, st->mode->prob, &dec);
613
614    /* Pitch MDCT */
615    compute_mdcts(st->mode, st->mode->window, st->out_mem+pitch_index*C, freq);
616
617    {
618       VARDECL(celt_ener_t, bandEp);
619       ALLOC(bandEp, st->mode->nbEBands*C, celt_ener_t);
620       compute_band_energies(st->mode, freq, bandEp);
621       normalise_bands(st->mode, freq, P, bandEp);
622    }
623
624    ALLOC(stereo_mode, st->mode->nbEBands, int);
625    stereo_decision(st->mode, X, stereo_mode, st->mode->nbEBands);
626    /* Apply pitch gains */
627    pitch_quant_bands(st->mode, P, gains);
628
629    /* Decode fixed codebook and merge with pitch */
630    unquant_bands(st->mode, X, P, bandE, stereo_mode, len*8, &dec);
631
632    if (C==2)
633    {
634       renormalise_bands(st->mode, X);
635    }
636    /* Synthesis */
637    denormalise_bands(st->mode, X, freq, bandE);
638
639
640    CELT_MOVE(st->out_mem, st->out_mem+C*N, C*(MAX_PERIOD+st->overlap-N));
641    /* Compute inverse MDCTs */
642    compute_inv_mdcts(st->mode, st->mode->window, freq, st->out_mem);
643
644    for (c=0;c<C;c++)
645    {
646       int j;
647       const celt_sig_t * restrict outp=st->out_mem+C*(MAX_PERIOD-N)+c;
648       celt_int16_t * restrict pcmp = pcm+c;
649       for (j=0;j<N;j++)
650       {
651          celt_sig_t tmp = ADD32(*outp, MULT16_32_Q15(preemph,st->preemph_memD[c]));
652          st->preemph_memD[c] = tmp;
653          *pcmp = SIG2INT16(tmp);
654          pcmp += C;
655          outp += C;
656       }
657    }
658
659    {
660       unsigned int val = 0;
661       while (ec_dec_tell(&dec, 0) < len*8)
662       {
663          if (ec_dec_uint(&dec, 2) != val)
664          {
665             celt_warning("decode error");
666             RESTORE_STACK;
667             return CELT_CORRUPTED_DATA;
668          }
669          val = 1-val;
670       }
671    }
672
673    RESTORE_STACK;
674    return 0;
675    /*printf ("\n");*/
676 }
677