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