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