Moves the remains of compute_inv_mdcts() to celt_synthesis()
[opus.git] / celt / celt_decoder.c
1 /* Copyright (c) 2007-2008 CSIRO
2    Copyright (c) 2007-2010 Xiph.Org Foundation
3    Copyright (c) 2008 Gregory Maxwell
4    Written by Jean-Marc Valin and Gregory Maxwell */
5 /*
6    Redistribution and use in source and binary forms, with or without
7    modification, are permitted provided that the following conditions
8    are met:
9
10    - Redistributions of source code must retain the above copyright
11    notice, this list of conditions and the following disclaimer.
12
13    - Redistributions in binary form must reproduce the above copyright
14    notice, this list of conditions and the following disclaimer in the
15    documentation and/or other materials provided with the distribution.
16
17    THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
18    ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
19    LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
20    A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
21    OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
22    EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
23    PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
24    PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
25    LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
26    NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
27    SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
28 */
29
30 #ifdef HAVE_CONFIG_H
31 #include "config.h"
32 #endif
33
34 #define CELT_DECODER_C
35
36 #include "cpu_support.h"
37 #include "os_support.h"
38 #include "mdct.h"
39 #include <math.h>
40 #include "celt.h"
41 #include "pitch.h"
42 #include "bands.h"
43 #include "modes.h"
44 #include "entcode.h"
45 #include "quant_bands.h"
46 #include "rate.h"
47 #include "stack_alloc.h"
48 #include "mathops.h"
49 #include "float_cast.h"
50 #include <stdarg.h>
51 #include "celt_lpc.h"
52 #include "vq.h"
53
54 /**********************************************************************/
55 /*                                                                    */
56 /*                             DECODER                                */
57 /*                                                                    */
58 /**********************************************************************/
59 #define DECODE_BUFFER_SIZE 2048
60
61 /** Decoder state
62  @brief Decoder state
63  */
64 struct OpusCustomDecoder {
65    const OpusCustomMode *mode;
66    int overlap;
67    int channels;
68    int stream_channels;
69
70    int downsample;
71    int start, end;
72    int signalling;
73    int arch;
74
75    /* Everything beyond this point gets cleared on a reset */
76 #define DECODER_RESET_START rng
77
78    opus_uint32 rng;
79    int error;
80    int last_pitch_index;
81    int loss_count;
82    int postfilter_period;
83    int postfilter_period_old;
84    opus_val16 postfilter_gain;
85    opus_val16 postfilter_gain_old;
86    int postfilter_tapset;
87    int postfilter_tapset_old;
88
89    celt_sig preemph_memD[2];
90
91    celt_sig _decode_mem[1]; /* Size = channels*(DECODE_BUFFER_SIZE+mode->overlap) */
92    /* opus_val16 lpc[],  Size = channels*LPC_ORDER */
93    /* opus_val16 oldEBands[], Size = 2*mode->nbEBands */
94    /* opus_val16 oldLogE[], Size = 2*mode->nbEBands */
95    /* opus_val16 oldLogE2[], Size = 2*mode->nbEBands */
96    /* opus_val16 backgroundLogE[], Size = 2*mode->nbEBands */
97 };
98
99 int celt_decoder_get_size(int channels)
100 {
101    const CELTMode *mode = opus_custom_mode_create(48000, 960, NULL);
102    return opus_custom_decoder_get_size(mode, channels);
103 }
104
105 OPUS_CUSTOM_NOSTATIC int opus_custom_decoder_get_size(const CELTMode *mode, int channels)
106 {
107    int size = sizeof(struct CELTDecoder)
108             + (channels*(DECODE_BUFFER_SIZE+mode->overlap)-1)*sizeof(celt_sig)
109             + channels*LPC_ORDER*sizeof(opus_val16)
110             + 4*2*mode->nbEBands*sizeof(opus_val16);
111    return size;
112 }
113
114 #ifdef CUSTOM_MODES
115 CELTDecoder *opus_custom_decoder_create(const CELTMode *mode, int channels, int *error)
116 {
117    int ret;
118    CELTDecoder *st = (CELTDecoder *)opus_alloc(opus_custom_decoder_get_size(mode, channels));
119    ret = opus_custom_decoder_init(st, mode, channels);
120    if (ret != OPUS_OK)
121    {
122       opus_custom_decoder_destroy(st);
123       st = NULL;
124    }
125    if (error)
126       *error = ret;
127    return st;
128 }
129 #endif /* CUSTOM_MODES */
130
131 int celt_decoder_init(CELTDecoder *st, opus_int32 sampling_rate, int channels)
132 {
133    int ret;
134    ret = opus_custom_decoder_init(st, opus_custom_mode_create(48000, 960, NULL), channels);
135    if (ret != OPUS_OK)
136       return ret;
137    st->downsample = resampling_factor(sampling_rate);
138    if (st->downsample==0)
139       return OPUS_BAD_ARG;
140    else
141       return OPUS_OK;
142 }
143
144 OPUS_CUSTOM_NOSTATIC int opus_custom_decoder_init(CELTDecoder *st, const CELTMode *mode, int channels)
145 {
146    if (channels < 0 || channels > 2)
147       return OPUS_BAD_ARG;
148
149    if (st==NULL)
150       return OPUS_ALLOC_FAIL;
151
152    OPUS_CLEAR((char*)st, opus_custom_decoder_get_size(mode, channels));
153
154    st->mode = mode;
155    st->overlap = mode->overlap;
156    st->stream_channels = st->channels = channels;
157
158    st->downsample = 1;
159    st->start = 0;
160    st->end = st->mode->effEBands;
161    st->signalling = 1;
162    st->arch = opus_select_arch();
163
164    st->loss_count = 0;
165
166    opus_custom_decoder_ctl(st, OPUS_RESET_STATE);
167
168    return OPUS_OK;
169 }
170
171 #ifdef CUSTOM_MODES
172 void opus_custom_decoder_destroy(CELTDecoder *st)
173 {
174    opus_free(st);
175 }
176 #endif /* CUSTOM_MODES */
177
178 static OPUS_INLINE opus_val16 SIG2WORD16(celt_sig x)
179 {
180 #ifdef FIXED_POINT
181    x = PSHR32(x, SIG_SHIFT);
182    x = MAX32(x, -32768);
183    x = MIN32(x, 32767);
184    return EXTRACT16(x);
185 #else
186    return (opus_val16)x;
187 #endif
188 }
189
190 #ifndef RESYNTH
191 static
192 #endif
193 void deemphasis(celt_sig *in[], opus_val16 *pcm, int N, int C, int downsample, const opus_val16 *coef,
194       celt_sig *mem)
195 {
196    int c;
197    int Nd;
198    int apply_downsampling=0;
199    opus_val16 coef0;
200    VARDECL(celt_sig, scratch);
201    SAVE_STACK;
202
203    ALLOC(scratch, N, celt_sig);
204    coef0 = coef[0];
205    Nd = N/downsample;
206    c=0; do {
207       int j;
208       celt_sig * OPUS_RESTRICT x;
209       opus_val16  * OPUS_RESTRICT y;
210       celt_sig m = mem[c];
211       x =in[c];
212       y = pcm+c;
213 #ifdef CUSTOM_MODES
214       if (coef[1] != 0)
215       {
216          opus_val16 coef1 = coef[1];
217          opus_val16 coef3 = coef[3];
218          for (j=0;j<N;j++)
219          {
220             celt_sig tmp = x[j] + m + VERY_SMALL;
221             m = MULT16_32_Q15(coef0, tmp)
222                           - MULT16_32_Q15(coef1, x[j]);
223             tmp = SHL32(MULT16_32_Q15(coef3, tmp), 2);
224             scratch[j] = tmp;
225          }
226          apply_downsampling=1;
227       } else
228 #endif
229       if (downsample>1)
230       {
231          /* Shortcut for the standard (non-custom modes) case */
232          for (j=0;j<N;j++)
233          {
234             celt_sig tmp = x[j] + m + VERY_SMALL;
235             m = MULT16_32_Q15(coef0, tmp);
236             scratch[j] = tmp;
237          }
238          apply_downsampling=1;
239       } else {
240          /* Shortcut for the standard (non-custom modes) case */
241          for (j=0;j<N;j++)
242          {
243             celt_sig tmp = x[j] + m + VERY_SMALL;
244             m = MULT16_32_Q15(coef0, tmp);
245             y[j*C] = SCALEOUT(SIG2WORD16(tmp));
246          }
247       }
248       mem[c] = m;
249
250       if (apply_downsampling)
251       {
252          /* Perform down-sampling */
253          for (j=0;j<Nd;j++)
254             y[j*C] = SCALEOUT(SIG2WORD16(scratch[j*downsample]));
255       }
256    } while (++c<C);
257    RESTORE_STACK;
258 }
259
260 #ifndef RESYNTH
261 static
262 #endif
263 void celt_synthesis(const CELTMode *mode, celt_norm *X, celt_sig * out_syn[],
264       opus_val16 *oldBandE, int start, int effEnd, int C, int CC, int isTransient,
265       int LM, int downsample, int silence)
266 {
267    int c, i;
268    int M;
269    int b;
270    int B;
271    int N, NB;
272    int shift;
273    int nbEBands;
274    int overlap;
275    VARDECL(celt_sig, freq);
276    SAVE_STACK;
277
278    overlap = mode->overlap;
279    nbEBands = mode->nbEBands;
280    N = mode->shortMdctSize<<LM;
281    ALLOC(freq, N, celt_sig); /**< Interleaved signal MDCTs */
282    M = 1<<LM;
283
284    if (isTransient)
285    {
286       B = M;
287       NB = mode->shortMdctSize;
288       shift = mode->maxLM;
289    } else {
290       B = 1;
291       NB = mode->shortMdctSize<<LM;
292       shift = mode->maxLM-LM;
293    }
294
295    if (CC==2&&C==1)
296    {
297       /* Copying a mono streams to two channels */
298       celt_sig *freq2;
299       denormalise_bands(mode, X, freq, oldBandE, start, effEnd, M,
300             downsample, silence);
301       /* Store a temporary copy in the output buffer because the IMDCT destroys its input. */
302       freq2 = out_syn[1]+overlap/2;
303       OPUS_COPY(freq2, freq, N);
304       for (b=0;b<B;b++)
305          clt_mdct_backward(&mode->mdct, &freq2[b], out_syn[0]+NB*b, mode->window, overlap, shift, B);
306       for (b=0;b<B;b++)
307          clt_mdct_backward(&mode->mdct, &freq[b], out_syn[1]+NB*b, mode->window, overlap, shift, B);
308    } else if (CC==1&&C==2)
309    {
310       /* Downmixing a stereo stream to mono */
311       celt_sig *freq2;
312       freq2 = out_syn[0]+overlap/2;
313       denormalise_bands(mode, X, freq, oldBandE, start, effEnd, M,
314             downsample, silence);
315       /* Use the output buffer as temp array before downmixing. */
316       denormalise_bands(mode, X+N, freq2, oldBandE+nbEBands, start, effEnd, M,
317             downsample, silence);
318       for (i=0;i<N;i++)
319          freq[i] = HALF32(ADD32(freq[i],freq2[i]));
320       for (b=0;b<B;b++)
321          clt_mdct_backward(&mode->mdct, &freq[b], out_syn[0]+NB*b, mode->window, overlap, shift, B);
322    } else {
323       /* Normal case (mono or stereo) */
324       c=0; do {
325          denormalise_bands(mode, X+c*N, freq, oldBandE+c*nbEBands, start, effEnd, M,
326                downsample, silence);
327          for (b=0;b<B;b++)
328             clt_mdct_backward(&mode->mdct, &freq[b], out_syn[c]+NB*b, mode->window, overlap, shift, B);
329       } while (++c<CC);
330    }
331    RESTORE_STACK;
332 }
333
334 static void tf_decode(int start, int end, int isTransient, int *tf_res, int LM, ec_dec *dec)
335 {
336    int i, curr, tf_select;
337    int tf_select_rsv;
338    int tf_changed;
339    int logp;
340    opus_uint32 budget;
341    opus_uint32 tell;
342
343    budget = dec->storage*8;
344    tell = ec_tell(dec);
345    logp = isTransient ? 2 : 4;
346    tf_select_rsv = LM>0 && tell+logp+1<=budget;
347    budget -= tf_select_rsv;
348    tf_changed = curr = 0;
349    for (i=start;i<end;i++)
350    {
351       if (tell+logp<=budget)
352       {
353          curr ^= ec_dec_bit_logp(dec, logp);
354          tell = ec_tell(dec);
355          tf_changed |= curr;
356       }
357       tf_res[i] = curr;
358       logp = isTransient ? 4 : 5;
359    }
360    tf_select = 0;
361    if (tf_select_rsv &&
362      tf_select_table[LM][4*isTransient+0+tf_changed] !=
363      tf_select_table[LM][4*isTransient+2+tf_changed])
364    {
365       tf_select = ec_dec_bit_logp(dec, 1);
366    }
367    for (i=start;i<end;i++)
368    {
369       tf_res[i] = tf_select_table[LM][4*isTransient+2*tf_select+tf_res[i]];
370    }
371 }
372
373 /* The maximum pitch lag to allow in the pitch-based PLC. It's possible to save
374    CPU time in the PLC pitch search by making this smaller than MAX_PERIOD. The
375    current value corresponds to a pitch of 66.67 Hz. */
376 #define PLC_PITCH_LAG_MAX (720)
377 /* The minimum pitch lag to allow in the pitch-based PLC. This corresponds to a
378    pitch of 480 Hz. */
379 #define PLC_PITCH_LAG_MIN (100)
380
381 static void celt_decode_lost(CELTDecoder * OPUS_RESTRICT st, opus_val16 * OPUS_RESTRICT pcm, int N, int LM)
382 {
383    int c;
384    int i;
385    const int C = st->channels;
386    celt_sig *decode_mem[2];
387    celt_sig *out_syn[2];
388    opus_val16 *lpc;
389    opus_val16 *oldBandE, *oldLogE, *oldLogE2, *backgroundLogE;
390    const OpusCustomMode *mode;
391    int nbEBands;
392    int overlap;
393    int start;
394    int downsample;
395    int loss_count;
396    int noise_based;
397    const opus_int16 *eBands;
398    SAVE_STACK;
399
400    mode = st->mode;
401    nbEBands = mode->nbEBands;
402    overlap = mode->overlap;
403    eBands = mode->eBands;
404
405    c=0; do {
406       decode_mem[c] = st->_decode_mem + c*(DECODE_BUFFER_SIZE+overlap);
407       out_syn[c] = decode_mem[c]+DECODE_BUFFER_SIZE-N;
408    } while (++c<C);
409    lpc = (opus_val16*)(st->_decode_mem+(DECODE_BUFFER_SIZE+overlap)*C);
410    oldBandE = lpc+C*LPC_ORDER;
411    oldLogE = oldBandE + 2*nbEBands;
412    oldLogE2 = oldLogE + 2*nbEBands;
413    backgroundLogE = oldLogE2  + 2*nbEBands;
414
415    loss_count = st->loss_count;
416    start = st->start;
417    downsample = st->downsample;
418    noise_based = loss_count >= 5 || start != 0;
419    if (noise_based)
420    {
421       /* Noise-based PLC/CNG */
422       VARDECL(celt_norm, X);
423       opus_uint32 seed;
424       opus_val16 *plcLogE;
425       int end;
426       int effEnd;
427
428       end = st->end;
429       effEnd = IMAX(start, IMIN(end, mode->effEBands));
430
431       /* Share the interleaved signal MDCT coefficient buffer with the
432          deemphasis scratch buffer. */
433       ALLOC(X, C*N, celt_norm);   /**< Interleaved normalised MDCTs */
434
435       if (loss_count >= 5)
436          plcLogE = backgroundLogE;
437       else {
438          /* Energy decay */
439          opus_val16 decay = loss_count==0 ?
440                QCONST16(1.5f, DB_SHIFT) : QCONST16(.5f, DB_SHIFT);
441          c=0; do
442          {
443             for (i=start;i<end;i++)
444                oldBandE[c*nbEBands+i] -= decay;
445          } while (++c<C);
446          plcLogE = oldBandE;
447       }
448       seed = st->rng;
449       for (c=0;c<C;c++)
450       {
451          for (i=start;i<effEnd;i++)
452          {
453             int j;
454             int boffs;
455             int blen;
456             boffs = N*c+(eBands[i]<<LM);
457             blen = (eBands[i+1]-eBands[i])<<LM;
458             for (j=0;j<blen;j++)
459             {
460                seed = celt_lcg_rand(seed);
461                X[boffs+j] = (celt_norm)((opus_int32)seed>>20);
462             }
463             renormalise_vector(X+boffs, blen, Q15ONE);
464          }
465       }
466       st->rng = seed;
467
468       c=0; do {
469          OPUS_MOVE(decode_mem[c], decode_mem[c]+N,
470                DECODE_BUFFER_SIZE-N+(overlap>>1));
471       } while (++c<C);
472
473       celt_synthesis(mode, X, out_syn, plcLogE, start, effEnd, C, C, 0, LM, st->downsample, 0);
474    } else {
475       /* Pitch-based PLC */
476       const opus_val16 *window;
477       opus_val16 fade = Q15ONE;
478       int pitch_index;
479       VARDECL(opus_val32, etmp);
480       VARDECL(opus_val16, exc);
481
482       if (loss_count == 0)
483       {
484          VARDECL( opus_val16, lp_pitch_buf );
485          ALLOC( lp_pitch_buf, DECODE_BUFFER_SIZE>>1, opus_val16 );
486          pitch_downsample(decode_mem, lp_pitch_buf,
487                DECODE_BUFFER_SIZE, C, st->arch);
488          pitch_search(lp_pitch_buf+(PLC_PITCH_LAG_MAX>>1), lp_pitch_buf,
489                DECODE_BUFFER_SIZE-PLC_PITCH_LAG_MAX,
490                PLC_PITCH_LAG_MAX-PLC_PITCH_LAG_MIN, &pitch_index, st->arch);
491          pitch_index = PLC_PITCH_LAG_MAX-pitch_index;
492          st->last_pitch_index = pitch_index;
493       } else {
494          pitch_index = st->last_pitch_index;
495          fade = QCONST16(.8f,15);
496       }
497
498       ALLOC(etmp, overlap, opus_val32);
499       ALLOC(exc, MAX_PERIOD, opus_val16);
500       window = mode->window;
501       c=0; do {
502          opus_val16 decay;
503          opus_val16 attenuation;
504          opus_val32 S1=0;
505          celt_sig *buf;
506          int extrapolation_offset;
507          int extrapolation_len;
508          int exc_length;
509          int j;
510
511          buf = decode_mem[c];
512          for (i=0;i<MAX_PERIOD;i++) {
513             exc[i] = ROUND16(buf[DECODE_BUFFER_SIZE-MAX_PERIOD+i], SIG_SHIFT);
514          }
515
516          if (loss_count == 0)
517          {
518             opus_val32 ac[LPC_ORDER+1];
519             /* Compute LPC coefficients for the last MAX_PERIOD samples before
520                the first loss so we can work in the excitation-filter domain. */
521             _celt_autocorr(exc, ac, window, overlap,
522                    LPC_ORDER, MAX_PERIOD, st->arch);
523             /* Add a noise floor of -40 dB. */
524 #ifdef FIXED_POINT
525             ac[0] += SHR32(ac[0],13);
526 #else
527             ac[0] *= 1.0001f;
528 #endif
529             /* Use lag windowing to stabilize the Levinson-Durbin recursion. */
530             for (i=1;i<=LPC_ORDER;i++)
531             {
532                /*ac[i] *= exp(-.5*(2*M_PI*.002*i)*(2*M_PI*.002*i));*/
533 #ifdef FIXED_POINT
534                ac[i] -= MULT16_32_Q15(2*i*i, ac[i]);
535 #else
536                ac[i] -= ac[i]*(0.008f*0.008f)*i*i;
537 #endif
538             }
539             _celt_lpc(lpc+c*LPC_ORDER, ac, LPC_ORDER);
540          }
541          /* We want the excitation for 2 pitch periods in order to look for a
542             decaying signal, but we can't get more than MAX_PERIOD. */
543          exc_length = IMIN(2*pitch_index, MAX_PERIOD);
544          /* Initialize the LPC history with the samples just before the start
545             of the region for which we're computing the excitation. */
546          {
547             opus_val16 lpc_mem[LPC_ORDER];
548             for (i=0;i<LPC_ORDER;i++)
549             {
550                lpc_mem[i] =
551                      ROUND16(buf[DECODE_BUFFER_SIZE-exc_length-1-i], SIG_SHIFT);
552             }
553             /* Compute the excitation for exc_length samples before the loss. */
554             celt_fir(exc+MAX_PERIOD-exc_length, lpc+c*LPC_ORDER,
555                   exc+MAX_PERIOD-exc_length, exc_length, LPC_ORDER, lpc_mem);
556          }
557
558          /* Check if the waveform is decaying, and if so how fast.
559             We do this to avoid adding energy when concealing in a segment
560             with decaying energy. */
561          {
562             opus_val32 E1=1, E2=1;
563             int decay_length;
564 #ifdef FIXED_POINT
565             int shift = IMAX(0,2*celt_zlog2(celt_maxabs16(&exc[MAX_PERIOD-exc_length], exc_length))-20);
566 #endif
567             decay_length = exc_length>>1;
568             for (i=0;i<decay_length;i++)
569             {
570                opus_val16 e;
571                e = exc[MAX_PERIOD-decay_length+i];
572                E1 += SHR32(MULT16_16(e, e), shift);
573                e = exc[MAX_PERIOD-2*decay_length+i];
574                E2 += SHR32(MULT16_16(e, e), shift);
575             }
576             E1 = MIN32(E1, E2);
577             decay = celt_sqrt(frac_div32(SHR32(E1, 1), E2));
578          }
579
580          /* Move the decoder memory one frame to the left to give us room to
581             add the data for the new frame. We ignore the overlap that extends
582             past the end of the buffer, because we aren't going to use it. */
583          OPUS_MOVE(buf, buf+N, DECODE_BUFFER_SIZE-N);
584
585          /* Extrapolate from the end of the excitation with a period of
586             "pitch_index", scaling down each period by an additional factor of
587             "decay". */
588          extrapolation_offset = MAX_PERIOD-pitch_index;
589          /* We need to extrapolate enough samples to cover a complete MDCT
590             window (including overlap/2 samples on both sides). */
591          extrapolation_len = N+overlap;
592          /* We also apply fading if this is not the first loss. */
593          attenuation = MULT16_16_Q15(fade, decay);
594          for (i=j=0;i<extrapolation_len;i++,j++)
595          {
596             opus_val16 tmp;
597             if (j >= pitch_index) {
598                j -= pitch_index;
599                attenuation = MULT16_16_Q15(attenuation, decay);
600             }
601             buf[DECODE_BUFFER_SIZE-N+i] =
602                   SHL32(EXTEND32(MULT16_16_Q15(attenuation,
603                         exc[extrapolation_offset+j])), SIG_SHIFT);
604             /* Compute the energy of the previously decoded signal whose
605                excitation we're copying. */
606             tmp = ROUND16(
607                   buf[DECODE_BUFFER_SIZE-MAX_PERIOD-N+extrapolation_offset+j],
608                   SIG_SHIFT);
609             S1 += SHR32(MULT16_16(tmp, tmp), 8);
610          }
611
612          {
613             opus_val16 lpc_mem[LPC_ORDER];
614             /* Copy the last decoded samples (prior to the overlap region) to
615                synthesis filter memory so we can have a continuous signal. */
616             for (i=0;i<LPC_ORDER;i++)
617                lpc_mem[i] = ROUND16(buf[DECODE_BUFFER_SIZE-N-1-i], SIG_SHIFT);
618             /* Apply the synthesis filter to convert the excitation back into
619                the signal domain. */
620             celt_iir(buf+DECODE_BUFFER_SIZE-N, lpc+c*LPC_ORDER,
621                   buf+DECODE_BUFFER_SIZE-N, extrapolation_len, LPC_ORDER,
622                   lpc_mem);
623          }
624
625          /* Check if the synthesis energy is higher than expected, which can
626             happen with the signal changes during our window. If so,
627             attenuate. */
628          {
629             opus_val32 S2=0;
630             for (i=0;i<extrapolation_len;i++)
631             {
632                opus_val16 tmp = ROUND16(buf[DECODE_BUFFER_SIZE-N+i], SIG_SHIFT);
633                S2 += SHR32(MULT16_16(tmp, tmp), 8);
634             }
635             /* This checks for an "explosion" in the synthesis. */
636 #ifdef FIXED_POINT
637             if (!(S1 > SHR32(S2,2)))
638 #else
639             /* The float test is written this way to catch NaNs in the output
640                of the IIR filter at the same time. */
641             if (!(S1 > 0.2f*S2))
642 #endif
643             {
644                for (i=0;i<extrapolation_len;i++)
645                   buf[DECODE_BUFFER_SIZE-N+i] = 0;
646             } else if (S1 < S2)
647             {
648                opus_val16 ratio = celt_sqrt(frac_div32(SHR32(S1,1)+1,S2+1));
649                for (i=0;i<overlap;i++)
650                {
651                   opus_val16 tmp_g = Q15ONE
652                         - MULT16_16_Q15(window[i], Q15ONE-ratio);
653                   buf[DECODE_BUFFER_SIZE-N+i] =
654                         MULT16_32_Q15(tmp_g, buf[DECODE_BUFFER_SIZE-N+i]);
655                }
656                for (i=overlap;i<extrapolation_len;i++)
657                {
658                   buf[DECODE_BUFFER_SIZE-N+i] =
659                         MULT16_32_Q15(ratio, buf[DECODE_BUFFER_SIZE-N+i]);
660                }
661             }
662          }
663
664          /* Apply the pre-filter to the MDCT overlap for the next frame because
665             the post-filter will be re-applied in the decoder after the MDCT
666             overlap. */
667          comb_filter(etmp, buf+DECODE_BUFFER_SIZE,
668               st->postfilter_period, st->postfilter_period, overlap,
669               -st->postfilter_gain, -st->postfilter_gain,
670               st->postfilter_tapset, st->postfilter_tapset, NULL, 0);
671
672          /* Simulate TDAC on the concealed audio so that it blends with the
673             MDCT of the next frame. */
674          for (i=0;i<overlap/2;i++)
675          {
676             buf[DECODE_BUFFER_SIZE+i] =
677                MULT16_32_Q15(window[i], etmp[overlap-1-i])
678                + MULT16_32_Q15(window[overlap-i-1], etmp[i]);
679          }
680       } while (++c<C);
681    }
682
683    deemphasis(out_syn, pcm, N, C, downsample,
684          mode->preemph, st->preemph_memD);
685
686    st->loss_count = loss_count+1;
687
688    RESTORE_STACK;
689 }
690
691 int celt_decode_with_ec(CELTDecoder * OPUS_RESTRICT st, const unsigned char *data, int len, opus_val16 * OPUS_RESTRICT pcm, int frame_size, ec_dec *dec)
692 {
693    int c, i, N;
694    int spread_decision;
695    opus_int32 bits;
696    ec_dec _dec;
697    VARDECL(celt_sig, freq);
698 #ifdef SMALL_FOOTPRINT
699    celt_norm *X;
700 #else
701    VARDECL(celt_norm, X);
702 #endif
703    VARDECL(int, fine_quant);
704    VARDECL(int, pulses);
705    VARDECL(int, cap);
706    VARDECL(int, offsets);
707    VARDECL(int, fine_priority);
708    VARDECL(int, tf_res);
709    VARDECL(unsigned char, collapse_masks);
710    celt_sig *decode_mem[2];
711    celt_sig *out_syn[2];
712    opus_val16 *lpc;
713    opus_val16 *oldBandE, *oldLogE, *oldLogE2, *backgroundLogE;
714
715    int shortBlocks;
716    int isTransient;
717    int intra_ener;
718    const int CC = st->channels;
719    int LM, M;
720    int start;
721    int end;
722    int effEnd;
723    int codedBands;
724    int alloc_trim;
725    int postfilter_pitch;
726    opus_val16 postfilter_gain;
727    int intensity=0;
728    int dual_stereo=0;
729    opus_int32 total_bits;
730    opus_int32 balance;
731    opus_int32 tell;
732    int dynalloc_logp;
733    int postfilter_tapset;
734    int anti_collapse_rsv;
735    int anti_collapse_on=0;
736    int silence;
737    int C = st->stream_channels;
738    const OpusCustomMode *mode;
739    int nbEBands;
740    int overlap;
741    const opus_int16 *eBands;
742    ALLOC_STACK;
743
744    mode = st->mode;
745    nbEBands = mode->nbEBands;
746    overlap = mode->overlap;
747    eBands = mode->eBands;
748    start = st->start;
749    end = st->end;
750    frame_size *= st->downsample;
751
752    c=0; do {
753       decode_mem[c] = st->_decode_mem + c*(DECODE_BUFFER_SIZE+overlap);
754    } while (++c<CC);
755    lpc = (opus_val16*)(st->_decode_mem+(DECODE_BUFFER_SIZE+overlap)*CC);
756    oldBandE = lpc+CC*LPC_ORDER;
757    oldLogE = oldBandE + 2*nbEBands;
758    oldLogE2 = oldLogE + 2*nbEBands;
759    backgroundLogE = oldLogE2  + 2*nbEBands;
760
761 #ifdef CUSTOM_MODES
762    if (st->signalling && data!=NULL)
763    {
764       int data0=data[0];
765       /* Convert "standard mode" to Opus header */
766       if (mode->Fs==48000 && mode->shortMdctSize==120)
767       {
768          data0 = fromOpus(data0);
769          if (data0<0)
770             return OPUS_INVALID_PACKET;
771       }
772       st->end = end = IMAX(1, mode->effEBands-2*(data0>>5));
773       LM = (data0>>3)&0x3;
774       C = 1 + ((data0>>2)&0x1);
775       data++;
776       len--;
777       if (LM>mode->maxLM)
778          return OPUS_INVALID_PACKET;
779       if (frame_size < mode->shortMdctSize<<LM)
780          return OPUS_BUFFER_TOO_SMALL;
781       else
782          frame_size = mode->shortMdctSize<<LM;
783    } else {
784 #else
785    {
786 #endif
787       for (LM=0;LM<=mode->maxLM;LM++)
788          if (mode->shortMdctSize<<LM==frame_size)
789             break;
790       if (LM>mode->maxLM)
791          return OPUS_BAD_ARG;
792    }
793    M=1<<LM;
794
795    if (len<0 || len>1275 || pcm==NULL)
796       return OPUS_BAD_ARG;
797
798    N = M*mode->shortMdctSize;
799
800    effEnd = end;
801    if (effEnd > mode->effEBands)
802       effEnd = mode->effEBands;
803
804    if (data == NULL || len<=1)
805    {
806       celt_decode_lost(st, pcm, N, LM);
807       RESTORE_STACK;
808       return frame_size/st->downsample;
809    }
810
811    if (dec == NULL)
812    {
813       ec_dec_init(&_dec,(unsigned char*)data,len);
814       dec = &_dec;
815    }
816
817    if (C==1)
818    {
819       for (i=0;i<nbEBands;i++)
820          oldBandE[i]=MAX16(oldBandE[i],oldBandE[nbEBands+i]);
821    }
822
823    total_bits = len*8;
824    tell = ec_tell(dec);
825
826    if (tell >= total_bits)
827       silence = 1;
828    else if (tell==1)
829       silence = ec_dec_bit_logp(dec, 15);
830    else
831       silence = 0;
832    if (silence)
833    {
834       /* Pretend we've read all the remaining bits */
835       tell = len*8;
836       dec->nbits_total+=tell-ec_tell(dec);
837    }
838
839    postfilter_gain = 0;
840    postfilter_pitch = 0;
841    postfilter_tapset = 0;
842    if (start==0 && tell+16 <= total_bits)
843    {
844       if(ec_dec_bit_logp(dec, 1))
845       {
846          int qg, octave;
847          octave = ec_dec_uint(dec, 6);
848          postfilter_pitch = (16<<octave)+ec_dec_bits(dec, 4+octave)-1;
849          qg = ec_dec_bits(dec, 3);
850          if (ec_tell(dec)+2<=total_bits)
851             postfilter_tapset = ec_dec_icdf(dec, tapset_icdf, 2);
852          postfilter_gain = QCONST16(.09375f,15)*(qg+1);
853       }
854       tell = ec_tell(dec);
855    }
856
857    if (LM > 0 && tell+3 <= total_bits)
858    {
859       isTransient = ec_dec_bit_logp(dec, 3);
860       tell = ec_tell(dec);
861    }
862    else
863       isTransient = 0;
864
865    if (isTransient)
866       shortBlocks = M;
867    else
868       shortBlocks = 0;
869
870    /* Decode the global flags (first symbols in the stream) */
871    intra_ener = tell+3<=total_bits ? ec_dec_bit_logp(dec, 3) : 0;
872    /* Get band energies */
873    unquant_coarse_energy(mode, start, end, oldBandE,
874          intra_ener, dec, C, LM);
875
876    ALLOC(tf_res, nbEBands, int);
877    tf_decode(start, end, isTransient, tf_res, LM, dec);
878
879    tell = ec_tell(dec);
880    spread_decision = SPREAD_NORMAL;
881    if (tell+4 <= total_bits)
882       spread_decision = ec_dec_icdf(dec, spread_icdf, 5);
883
884    ALLOC(cap, nbEBands, int);
885
886    init_caps(mode,cap,LM,C);
887
888    ALLOC(offsets, nbEBands, int);
889
890    dynalloc_logp = 6;
891    total_bits<<=BITRES;
892    tell = ec_tell_frac(dec);
893    for (i=start;i<end;i++)
894    {
895       int width, quanta;
896       int dynalloc_loop_logp;
897       int boost;
898       width = C*(eBands[i+1]-eBands[i])<<LM;
899       /* quanta is 6 bits, but no more than 1 bit/sample
900          and no less than 1/8 bit/sample */
901       quanta = IMIN(width<<BITRES, IMAX(6<<BITRES, width));
902       dynalloc_loop_logp = dynalloc_logp;
903       boost = 0;
904       while (tell+(dynalloc_loop_logp<<BITRES) < total_bits && boost < cap[i])
905       {
906          int flag;
907          flag = ec_dec_bit_logp(dec, dynalloc_loop_logp);
908          tell = ec_tell_frac(dec);
909          if (!flag)
910             break;
911          boost += quanta;
912          total_bits -= quanta;
913          dynalloc_loop_logp = 1;
914       }
915       offsets[i] = boost;
916       /* Making dynalloc more likely */
917       if (boost>0)
918          dynalloc_logp = IMAX(2, dynalloc_logp-1);
919    }
920
921    ALLOC(fine_quant, nbEBands, int);
922    alloc_trim = tell+(6<<BITRES) <= total_bits ?
923          ec_dec_icdf(dec, trim_icdf, 7) : 5;
924
925    bits = (((opus_int32)len*8)<<BITRES) - ec_tell_frac(dec) - 1;
926    anti_collapse_rsv = isTransient&&LM>=2&&bits>=((LM+2)<<BITRES) ? (1<<BITRES) : 0;
927    bits -= anti_collapse_rsv;
928
929    ALLOC(pulses, nbEBands, int);
930    ALLOC(fine_priority, nbEBands, int);
931
932    codedBands = compute_allocation(mode, start, end, offsets, cap,
933          alloc_trim, &intensity, &dual_stereo, bits, &balance, pulses,
934          fine_quant, fine_priority, C, LM, dec, 0, 0, 0);
935
936    unquant_fine_energy(mode, start, end, oldBandE, fine_quant, dec, C);
937
938    c=0; do {
939       OPUS_MOVE(decode_mem[c], decode_mem[c]+N, DECODE_BUFFER_SIZE-N+overlap/2);
940       out_syn[c] = decode_mem[c]+DECODE_BUFFER_SIZE-N;
941    } while (++c<CC);
942
943    /* Decode fixed codebook */
944    ALLOC(collapse_masks, C*nbEBands, unsigned char);
945
946 #ifdef SMALL_FOOTPRINT
947    /* This is an ugly hack that breaks aliasing rules and would be easily broken,
948       but it saves almost 4kB of stack. */
949    X = (celt_norm*)(out_syn[CC-1]+overlap/2);
950 #else
951    ALLOC(X, C*N, celt_norm);   /**< Interleaved normalised MDCTs */
952 #endif
953
954    quant_all_bands(0, mode, start, end, X, C==2 ? X+N : NULL, collapse_masks,
955          NULL, pulses, shortBlocks, spread_decision, dual_stereo, intensity, tf_res,
956          len*(8<<BITRES)-anti_collapse_rsv, balance, dec, LM, codedBands, &st->rng);
957
958    if (anti_collapse_rsv > 0)
959    {
960       anti_collapse_on = ec_dec_bits(dec, 1);
961    }
962
963    unquant_energy_finalise(mode, start, end, oldBandE,
964          fine_quant, fine_priority, len*8-ec_tell(dec), dec, C);
965
966    if (anti_collapse_on)
967       anti_collapse(mode, X, collapse_masks, LM, C, N,
968             start, end, oldBandE, oldLogE, oldLogE2, pulses, st->rng);
969
970    if (silence)
971    {
972       for (i=0;i<C*nbEBands;i++)
973          oldBandE[i] = -QCONST16(28.f,DB_SHIFT);
974    }
975
976    celt_synthesis(mode, X, out_syn, oldBandE, start, effEnd, C, CC, isTransient, LM, st->downsample, silence);
977
978    c=0; do {
979       st->postfilter_period=IMAX(st->postfilter_period, COMBFILTER_MINPERIOD);
980       st->postfilter_period_old=IMAX(st->postfilter_period_old, COMBFILTER_MINPERIOD);
981       comb_filter(out_syn[c], out_syn[c], st->postfilter_period_old, st->postfilter_period, mode->shortMdctSize,
982             st->postfilter_gain_old, st->postfilter_gain, st->postfilter_tapset_old, st->postfilter_tapset,
983             mode->window, overlap);
984       if (LM!=0)
985          comb_filter(out_syn[c]+mode->shortMdctSize, out_syn[c]+mode->shortMdctSize, st->postfilter_period, postfilter_pitch, N-mode->shortMdctSize,
986                st->postfilter_gain, postfilter_gain, st->postfilter_tapset, postfilter_tapset,
987                mode->window, overlap);
988
989    } while (++c<CC);
990    st->postfilter_period_old = st->postfilter_period;
991    st->postfilter_gain_old = st->postfilter_gain;
992    st->postfilter_tapset_old = st->postfilter_tapset;
993    st->postfilter_period = postfilter_pitch;
994    st->postfilter_gain = postfilter_gain;
995    st->postfilter_tapset = postfilter_tapset;
996    if (LM!=0)
997    {
998       st->postfilter_period_old = st->postfilter_period;
999       st->postfilter_gain_old = st->postfilter_gain;
1000       st->postfilter_tapset_old = st->postfilter_tapset;
1001    }
1002
1003    if (C==1)
1004       OPUS_COPY(&oldBandE[nbEBands], oldBandE, nbEBands);
1005
1006    /* In case start or end were to change */
1007    if (!isTransient)
1008    {
1009       OPUS_COPY(oldLogE2, oldLogE, 2*nbEBands);
1010       OPUS_COPY(oldLogE, oldBandE, 2*nbEBands);
1011       for (i=0;i<2*nbEBands;i++)
1012          backgroundLogE[i] = MIN16(backgroundLogE[i] + M*QCONST16(0.001f,DB_SHIFT), oldBandE[i]);
1013    } else {
1014       for (i=0;i<2*nbEBands;i++)
1015          oldLogE[i] = MIN16(oldLogE[i], oldBandE[i]);
1016    }
1017    c=0; do
1018    {
1019       for (i=0;i<start;i++)
1020       {
1021          oldBandE[c*nbEBands+i]=0;
1022          oldLogE[c*nbEBands+i]=oldLogE2[c*nbEBands+i]=-QCONST16(28.f,DB_SHIFT);
1023       }
1024       for (i=end;i<nbEBands;i++)
1025       {
1026          oldBandE[c*nbEBands+i]=0;
1027          oldLogE[c*nbEBands+i]=oldLogE2[c*nbEBands+i]=-QCONST16(28.f,DB_SHIFT);
1028       }
1029    } while (++c<2);
1030    st->rng = dec->rng;
1031
1032    /* We reuse freq[] as scratch space for the de-emphasis */
1033    deemphasis(out_syn, pcm, N, CC, st->downsample, mode->preemph, st->preemph_memD);
1034    st->loss_count = 0;
1035    RESTORE_STACK;
1036    if (ec_tell(dec) > 8*len)
1037       return OPUS_INTERNAL_ERROR;
1038    if(ec_get_error(dec))
1039       st->error = 1;
1040    return frame_size/st->downsample;
1041 }
1042
1043
1044 #ifdef CUSTOM_MODES
1045
1046 #ifdef FIXED_POINT
1047 int opus_custom_decode(CELTDecoder * OPUS_RESTRICT st, const unsigned char *data, int len, opus_int16 * OPUS_RESTRICT pcm, int frame_size)
1048 {
1049    return celt_decode_with_ec(st, data, len, pcm, frame_size, NULL);
1050 }
1051
1052 #ifndef DISABLE_FLOAT_API
1053 int opus_custom_decode_float(CELTDecoder * OPUS_RESTRICT st, const unsigned char *data, int len, float * OPUS_RESTRICT pcm, int frame_size)
1054 {
1055    int j, ret, C, N;
1056    VARDECL(opus_int16, out);
1057    ALLOC_STACK;
1058
1059    if (pcm==NULL)
1060       return OPUS_BAD_ARG;
1061
1062    C = st->channels;
1063    N = frame_size;
1064
1065    ALLOC(out, C*N, opus_int16);
1066    ret=celt_decode_with_ec(st, data, len, out, frame_size, NULL);
1067    if (ret>0)
1068       for (j=0;j<C*ret;j++)
1069          pcm[j]=out[j]*(1.f/32768.f);
1070
1071    RESTORE_STACK;
1072    return ret;
1073 }
1074 #endif /* DISABLE_FLOAT_API */
1075
1076 #else
1077
1078 int opus_custom_decode_float(CELTDecoder * OPUS_RESTRICT st, const unsigned char *data, int len, float * OPUS_RESTRICT pcm, int frame_size)
1079 {
1080    return celt_decode_with_ec(st, data, len, pcm, frame_size, NULL);
1081 }
1082
1083 int opus_custom_decode(CELTDecoder * OPUS_RESTRICT st, const unsigned char *data, int len, opus_int16 * OPUS_RESTRICT pcm, int frame_size)
1084 {
1085    int j, ret, C, N;
1086    VARDECL(celt_sig, out);
1087    ALLOC_STACK;
1088
1089    if (pcm==NULL)
1090       return OPUS_BAD_ARG;
1091
1092    C = st->channels;
1093    N = frame_size;
1094    ALLOC(out, C*N, celt_sig);
1095
1096    ret=celt_decode_with_ec(st, data, len, out, frame_size, NULL);
1097
1098    if (ret>0)
1099       for (j=0;j<C*ret;j++)
1100          pcm[j] = FLOAT2INT16 (out[j]);
1101
1102    RESTORE_STACK;
1103    return ret;
1104 }
1105
1106 #endif
1107 #endif /* CUSTOM_MODES */
1108
1109 int opus_custom_decoder_ctl(CELTDecoder * OPUS_RESTRICT st, int request, ...)
1110 {
1111    va_list ap;
1112
1113    va_start(ap, request);
1114    switch (request)
1115    {
1116       case CELT_SET_START_BAND_REQUEST:
1117       {
1118          opus_int32 value = va_arg(ap, opus_int32);
1119          if (value<0 || value>=st->mode->nbEBands)
1120             goto bad_arg;
1121          st->start = value;
1122       }
1123       break;
1124       case CELT_SET_END_BAND_REQUEST:
1125       {
1126          opus_int32 value = va_arg(ap, opus_int32);
1127          if (value<1 || value>st->mode->nbEBands)
1128             goto bad_arg;
1129          st->end = value;
1130       }
1131       break;
1132       case CELT_SET_CHANNELS_REQUEST:
1133       {
1134          opus_int32 value = va_arg(ap, opus_int32);
1135          if (value<1 || value>2)
1136             goto bad_arg;
1137          st->stream_channels = value;
1138       }
1139       break;
1140       case CELT_GET_AND_CLEAR_ERROR_REQUEST:
1141       {
1142          opus_int32 *value = va_arg(ap, opus_int32*);
1143          if (value==NULL)
1144             goto bad_arg;
1145          *value=st->error;
1146          st->error = 0;
1147       }
1148       break;
1149       case OPUS_GET_LOOKAHEAD_REQUEST:
1150       {
1151          opus_int32 *value = va_arg(ap, opus_int32*);
1152          if (value==NULL)
1153             goto bad_arg;
1154          *value = st->overlap/st->downsample;
1155       }
1156       break;
1157       case OPUS_RESET_STATE:
1158       {
1159          int i;
1160          opus_val16 *lpc, *oldBandE, *oldLogE, *oldLogE2;
1161          lpc = (opus_val16*)(st->_decode_mem+(DECODE_BUFFER_SIZE+st->overlap)*st->channels);
1162          oldBandE = lpc+st->channels*LPC_ORDER;
1163          oldLogE = oldBandE + 2*st->mode->nbEBands;
1164          oldLogE2 = oldLogE + 2*st->mode->nbEBands;
1165          OPUS_CLEAR((char*)&st->DECODER_RESET_START,
1166                opus_custom_decoder_get_size(st->mode, st->channels)-
1167                ((char*)&st->DECODER_RESET_START - (char*)st));
1168          for (i=0;i<2*st->mode->nbEBands;i++)
1169             oldLogE[i]=oldLogE2[i]=-QCONST16(28.f,DB_SHIFT);
1170       }
1171       break;
1172       case OPUS_GET_PITCH_REQUEST:
1173       {
1174          opus_int32 *value = va_arg(ap, opus_int32*);
1175          if (value==NULL)
1176             goto bad_arg;
1177          *value = st->postfilter_period;
1178       }
1179       break;
1180       case CELT_GET_MODE_REQUEST:
1181       {
1182          const CELTMode ** value = va_arg(ap, const CELTMode**);
1183          if (value==0)
1184             goto bad_arg;
1185          *value=st->mode;
1186       }
1187       break;
1188       case CELT_SET_SIGNALLING_REQUEST:
1189       {
1190          opus_int32 value = va_arg(ap, opus_int32);
1191          st->signalling = value;
1192       }
1193       break;
1194       case OPUS_GET_FINAL_RANGE_REQUEST:
1195       {
1196          opus_uint32 * value = va_arg(ap, opus_uint32 *);
1197          if (value==0)
1198             goto bad_arg;
1199          *value=st->rng;
1200       }
1201       break;
1202       default:
1203          goto bad_request;
1204    }
1205    va_end(ap);
1206    return OPUS_OK;
1207 bad_arg:
1208    va_end(ap);
1209    return OPUS_BAD_ARG;
1210 bad_request:
1211       va_end(ap);
1212   return OPUS_UNIMPLEMENTED;
1213 }