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