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