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