Reduces decoder stack usage by only storing one channel of denormalized MDCT
[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 /** Compute the IMDCT and apply window for all sub-frames and
261     all channels in a frame */
262 #ifndef RESYNTH
263 static
264 #endif
265 void compute_inv_mdcts(const CELTMode *mode, int shortBlocks, celt_sig *X,
266       celt_sig * OPUS_RESTRICT out_mem, int LM)
267 {
268    int b;
269    int B;
270    int N;
271    int shift;
272    const int overlap = mode->overlap;
273
274    if (shortBlocks)
275    {
276       B = shortBlocks;
277       N = mode->shortMdctSize;
278       shift = mode->maxLM;
279    } else {
280       B = 1;
281       N = mode->shortMdctSize<<LM;
282       shift = mode->maxLM-LM;
283    }
284    /* IMDCT on the interleaved the sub-frames, overlap-add is performed by the IMDCT */
285    for (b=0;b<B;b++)
286       clt_mdct_backward(&mode->mdct, &X[b], out_mem+N*b, mode->window, overlap, shift, B);
287 }
288
289 #ifndef RESYNTH
290 static
291 #endif
292 void celt_synthesis(const CELTMode *mode, celt_norm *X, celt_sig * out_syn[],
293       opus_val16 *oldBandE, int start, int effEnd, int C, int CC, int isTransient,
294       int LM, int downsample, int silence)
295 {
296    int c, i;
297    int M, N;
298    int nbEBands;
299    int shortBlocks;
300    int overlap;
301    VARDECL(celt_sig, freq);
302    SAVE_STACK;
303
304    overlap = mode->overlap;
305    nbEBands = mode->nbEBands;
306    N = mode->shortMdctSize<<LM;
307    ALLOC(freq, N, celt_sig); /**< Interleaved signal MDCTs */
308    M = 1<<LM;
309    shortBlocks = isTransient ? M : 0;
310
311    if (CC==2&&C==1)
312    {
313       /* Copying a mono streams to two channels */
314       celt_sig *freq2;
315       denormalise_bands(mode, X, freq, oldBandE, start, effEnd, M,
316             downsample, silence);
317       /* Store a temporary copy in the output buffer because the IMDCT destroys its input. */
318       freq2 = out_syn[1]+overlap/2;
319       OPUS_COPY(freq2, freq, N);
320       compute_inv_mdcts(mode, shortBlocks, freq2, out_syn[0], LM);
321       compute_inv_mdcts(mode, shortBlocks, freq, out_syn[1], LM);
322    } else if (CC==1&&C==2)
323    {
324       /* Downmixing a stereo stream to mono */
325       celt_sig *freq2;
326       freq2 = out_syn[0]+overlap/2;
327       denormalise_bands(mode, X, freq, oldBandE, start, effEnd, M,
328             downsample, silence);
329       /* Use the output buffer as temp array before downmixing. */
330       denormalise_bands(mode, X+N, freq2, oldBandE+nbEBands, start, effEnd, M,
331             downsample, silence);
332       for (i=0;i<N;i++)
333          freq[i] = HALF32(ADD32(freq[i],freq2[i]));
334       /* Compute inverse MDCTs */
335       compute_inv_mdcts(mode, shortBlocks, freq, out_syn[0], LM);
336    } else {
337       /* Normal case (mono or stereo) */
338       c=0; do {
339          denormalise_bands(mode, X+c*N, freq, oldBandE+c*nbEBands, start, effEnd, M,
340                downsample, silence);
341          compute_inv_mdcts(mode, shortBlocks, freq, out_syn[c], LM);
342       } while (++c<CC);
343    }
344    RESTORE_STACK;
345 }
346
347 static void tf_decode(int start, int end, int isTransient, int *tf_res, int LM, ec_dec *dec)
348 {
349    int i, curr, tf_select;
350    int tf_select_rsv;
351    int tf_changed;
352    int logp;
353    opus_uint32 budget;
354    opus_uint32 tell;
355
356    budget = dec->storage*8;
357    tell = ec_tell(dec);
358    logp = isTransient ? 2 : 4;
359    tf_select_rsv = LM>0 && tell+logp+1<=budget;
360    budget -= tf_select_rsv;
361    tf_changed = curr = 0;
362    for (i=start;i<end;i++)
363    {
364       if (tell+logp<=budget)
365       {
366          curr ^= ec_dec_bit_logp(dec, logp);
367          tell = ec_tell(dec);
368          tf_changed |= curr;
369       }
370       tf_res[i] = curr;
371       logp = isTransient ? 4 : 5;
372    }
373    tf_select = 0;
374    if (tf_select_rsv &&
375      tf_select_table[LM][4*isTransient+0+tf_changed] !=
376      tf_select_table[LM][4*isTransient+2+tf_changed])
377    {
378       tf_select = ec_dec_bit_logp(dec, 1);
379    }
380    for (i=start;i<end;i++)
381    {
382       tf_res[i] = tf_select_table[LM][4*isTransient+2*tf_select+tf_res[i]];
383    }
384 }
385
386 /* The maximum pitch lag to allow in the pitch-based PLC. It's possible to save
387    CPU time in the PLC pitch search by making this smaller than MAX_PERIOD. The
388    current value corresponds to a pitch of 66.67 Hz. */
389 #define PLC_PITCH_LAG_MAX (720)
390 /* The minimum pitch lag to allow in the pitch-based PLC. This corresponds to a
391    pitch of 480 Hz. */
392 #define PLC_PITCH_LAG_MIN (100)
393
394 static void celt_decode_lost(CELTDecoder * OPUS_RESTRICT st, opus_val16 * OPUS_RESTRICT pcm, int N, int LM)
395 {
396    int c;
397    int i;
398    const int C = st->channels;
399    celt_sig *decode_mem[2];
400    celt_sig *out_syn[2];
401    opus_val16 *lpc;
402    opus_val16 *oldBandE, *oldLogE, *oldLogE2, *backgroundLogE;
403    const OpusCustomMode *mode;
404    int nbEBands;
405    int overlap;
406    int start;
407    int downsample;
408    int loss_count;
409    int noise_based;
410    const opus_int16 *eBands;
411    SAVE_STACK;
412
413    mode = st->mode;
414    nbEBands = mode->nbEBands;
415    overlap = mode->overlap;
416    eBands = mode->eBands;
417
418    c=0; do {
419       decode_mem[c] = st->_decode_mem + c*(DECODE_BUFFER_SIZE+overlap);
420       out_syn[c] = decode_mem[c]+DECODE_BUFFER_SIZE-N;
421    } while (++c<C);
422    lpc = (opus_val16*)(st->_decode_mem+(DECODE_BUFFER_SIZE+overlap)*C);
423    oldBandE = lpc+C*LPC_ORDER;
424    oldLogE = oldBandE + 2*nbEBands;
425    oldLogE2 = oldLogE + 2*nbEBands;
426    backgroundLogE = oldLogE2  + 2*nbEBands;
427
428    loss_count = st->loss_count;
429    start = st->start;
430    downsample = st->downsample;
431    noise_based = loss_count >= 5 || start != 0;
432    if (noise_based)
433    {
434       /* Noise-based PLC/CNG */
435       VARDECL(celt_norm, X);
436       opus_uint32 seed;
437       opus_val16 *plcLogE;
438       int end;
439       int effEnd;
440
441       end = st->end;
442       effEnd = IMAX(start, IMIN(end, mode->effEBands));
443
444       /* Share the interleaved signal MDCT coefficient buffer with the
445          deemphasis scratch buffer. */
446       ALLOC(X, C*N, celt_norm);   /**< Interleaved normalised MDCTs */
447
448       if (loss_count >= 5)
449          plcLogE = backgroundLogE;
450       else {
451          /* Energy decay */
452          opus_val16 decay = loss_count==0 ?
453                QCONST16(1.5f, DB_SHIFT) : QCONST16(.5f, DB_SHIFT);
454          c=0; do
455          {
456             for (i=start;i<end;i++)
457                oldBandE[c*nbEBands+i] -= decay;
458          } while (++c<C);
459          plcLogE = oldBandE;
460       }
461       seed = st->rng;
462       for (c=0;c<C;c++)
463       {
464          for (i=start;i<effEnd;i++)
465          {
466             int j;
467             int boffs;
468             int blen;
469             boffs = N*c+(eBands[i]<<LM);
470             blen = (eBands[i+1]-eBands[i])<<LM;
471             for (j=0;j<blen;j++)
472             {
473                seed = celt_lcg_rand(seed);
474                X[boffs+j] = (celt_norm)((opus_int32)seed>>20);
475             }
476             renormalise_vector(X+boffs, blen, Q15ONE);
477          }
478       }
479       st->rng = seed;
480
481       c=0; do {
482          OPUS_MOVE(decode_mem[c], decode_mem[c]+N,
483                DECODE_BUFFER_SIZE-N+(overlap>>1));
484       } while (++c<C);
485
486       celt_synthesis(mode, X, out_syn, plcLogE, start, effEnd, C, C, 0, LM, st->downsample, 0);
487    } else {
488       /* Pitch-based PLC */
489       const opus_val16 *window;
490       opus_val16 fade = Q15ONE;
491       int pitch_index;
492       VARDECL(opus_val32, etmp);
493       VARDECL(opus_val16, exc);
494
495       if (loss_count == 0)
496       {
497          VARDECL( opus_val16, lp_pitch_buf );
498          ALLOC( lp_pitch_buf, DECODE_BUFFER_SIZE>>1, opus_val16 );
499          pitch_downsample(decode_mem, lp_pitch_buf,
500                DECODE_BUFFER_SIZE, C, st->arch);
501          pitch_search(lp_pitch_buf+(PLC_PITCH_LAG_MAX>>1), lp_pitch_buf,
502                DECODE_BUFFER_SIZE-PLC_PITCH_LAG_MAX,
503                PLC_PITCH_LAG_MAX-PLC_PITCH_LAG_MIN, &pitch_index, st->arch);
504          pitch_index = PLC_PITCH_LAG_MAX-pitch_index;
505          st->last_pitch_index = pitch_index;
506       } else {
507          pitch_index = st->last_pitch_index;
508          fade = QCONST16(.8f,15);
509       }
510
511       ALLOC(etmp, overlap, opus_val32);
512       ALLOC(exc, MAX_PERIOD, opus_val16);
513       window = mode->window;
514       c=0; do {
515          opus_val16 decay;
516          opus_val16 attenuation;
517          opus_val32 S1=0;
518          celt_sig *buf;
519          int extrapolation_offset;
520          int extrapolation_len;
521          int exc_length;
522          int j;
523
524          buf = decode_mem[c];
525          for (i=0;i<MAX_PERIOD;i++) {
526             exc[i] = ROUND16(buf[DECODE_BUFFER_SIZE-MAX_PERIOD+i], SIG_SHIFT);
527          }
528
529          if (loss_count == 0)
530          {
531             opus_val32 ac[LPC_ORDER+1];
532             /* Compute LPC coefficients for the last MAX_PERIOD samples before
533                the first loss so we can work in the excitation-filter domain. */
534             _celt_autocorr(exc, ac, window, overlap,
535                    LPC_ORDER, MAX_PERIOD, st->arch);
536             /* Add a noise floor of -40 dB. */
537 #ifdef FIXED_POINT
538             ac[0] += SHR32(ac[0],13);
539 #else
540             ac[0] *= 1.0001f;
541 #endif
542             /* Use lag windowing to stabilize the Levinson-Durbin recursion. */
543             for (i=1;i<=LPC_ORDER;i++)
544             {
545                /*ac[i] *= exp(-.5*(2*M_PI*.002*i)*(2*M_PI*.002*i));*/
546 #ifdef FIXED_POINT
547                ac[i] -= MULT16_32_Q15(2*i*i, ac[i]);
548 #else
549                ac[i] -= ac[i]*(0.008f*0.008f)*i*i;
550 #endif
551             }
552             _celt_lpc(lpc+c*LPC_ORDER, ac, LPC_ORDER);
553          }
554          /* We want the excitation for 2 pitch periods in order to look for a
555             decaying signal, but we can't get more than MAX_PERIOD. */
556          exc_length = IMIN(2*pitch_index, MAX_PERIOD);
557          /* Initialize the LPC history with the samples just before the start
558             of the region for which we're computing the excitation. */
559          {
560             opus_val16 lpc_mem[LPC_ORDER];
561             for (i=0;i<LPC_ORDER;i++)
562             {
563                lpc_mem[i] =
564                      ROUND16(buf[DECODE_BUFFER_SIZE-exc_length-1-i], SIG_SHIFT);
565             }
566             /* Compute the excitation for exc_length samples before the loss. */
567             celt_fir(exc+MAX_PERIOD-exc_length, lpc+c*LPC_ORDER,
568                   exc+MAX_PERIOD-exc_length, exc_length, LPC_ORDER, lpc_mem);
569          }
570
571          /* Check if the waveform is decaying, and if so how fast.
572             We do this to avoid adding energy when concealing in a segment
573             with decaying energy. */
574          {
575             opus_val32 E1=1, E2=1;
576             int decay_length;
577 #ifdef FIXED_POINT
578             int shift = IMAX(0,2*celt_zlog2(celt_maxabs16(&exc[MAX_PERIOD-exc_length], exc_length))-20);
579 #endif
580             decay_length = exc_length>>1;
581             for (i=0;i<decay_length;i++)
582             {
583                opus_val16 e;
584                e = exc[MAX_PERIOD-decay_length+i];
585                E1 += SHR32(MULT16_16(e, e), shift);
586                e = exc[MAX_PERIOD-2*decay_length+i];
587                E2 += SHR32(MULT16_16(e, e), shift);
588             }
589             E1 = MIN32(E1, E2);
590             decay = celt_sqrt(frac_div32(SHR32(E1, 1), E2));
591          }
592
593          /* Move the decoder memory one frame to the left to give us room to
594             add the data for the new frame. We ignore the overlap that extends
595             past the end of the buffer, because we aren't going to use it. */
596          OPUS_MOVE(buf, buf+N, DECODE_BUFFER_SIZE-N);
597
598          /* Extrapolate from the end of the excitation with a period of
599             "pitch_index", scaling down each period by an additional factor of
600             "decay". */
601          extrapolation_offset = MAX_PERIOD-pitch_index;
602          /* We need to extrapolate enough samples to cover a complete MDCT
603             window (including overlap/2 samples on both sides). */
604          extrapolation_len = N+overlap;
605          /* We also apply fading if this is not the first loss. */
606          attenuation = MULT16_16_Q15(fade, decay);
607          for (i=j=0;i<extrapolation_len;i++,j++)
608          {
609             opus_val16 tmp;
610             if (j >= pitch_index) {
611                j -= pitch_index;
612                attenuation = MULT16_16_Q15(attenuation, decay);
613             }
614             buf[DECODE_BUFFER_SIZE-N+i] =
615                   SHL32(EXTEND32(MULT16_16_Q15(attenuation,
616                         exc[extrapolation_offset+j])), SIG_SHIFT);
617             /* Compute the energy of the previously decoded signal whose
618                excitation we're copying. */
619             tmp = ROUND16(
620                   buf[DECODE_BUFFER_SIZE-MAX_PERIOD-N+extrapolation_offset+j],
621                   SIG_SHIFT);
622             S1 += SHR32(MULT16_16(tmp, tmp), 8);
623          }
624
625          {
626             opus_val16 lpc_mem[LPC_ORDER];
627             /* Copy the last decoded samples (prior to the overlap region) to
628                synthesis filter memory so we can have a continuous signal. */
629             for (i=0;i<LPC_ORDER;i++)
630                lpc_mem[i] = ROUND16(buf[DECODE_BUFFER_SIZE-N-1-i], SIG_SHIFT);
631             /* Apply the synthesis filter to convert the excitation back into
632                the signal domain. */
633             celt_iir(buf+DECODE_BUFFER_SIZE-N, lpc+c*LPC_ORDER,
634                   buf+DECODE_BUFFER_SIZE-N, extrapolation_len, LPC_ORDER,
635                   lpc_mem);
636          }
637
638          /* Check if the synthesis energy is higher than expected, which can
639             happen with the signal changes during our window. If so,
640             attenuate. */
641          {
642             opus_val32 S2=0;
643             for (i=0;i<extrapolation_len;i++)
644             {
645                opus_val16 tmp = ROUND16(buf[DECODE_BUFFER_SIZE-N+i], SIG_SHIFT);
646                S2 += SHR32(MULT16_16(tmp, tmp), 8);
647             }
648             /* This checks for an "explosion" in the synthesis. */
649 #ifdef FIXED_POINT
650             if (!(S1 > SHR32(S2,2)))
651 #else
652             /* The float test is written this way to catch NaNs in the output
653                of the IIR filter at the same time. */
654             if (!(S1 > 0.2f*S2))
655 #endif
656             {
657                for (i=0;i<extrapolation_len;i++)
658                   buf[DECODE_BUFFER_SIZE-N+i] = 0;
659             } else if (S1 < S2)
660             {
661                opus_val16 ratio = celt_sqrt(frac_div32(SHR32(S1,1)+1,S2+1));
662                for (i=0;i<overlap;i++)
663                {
664                   opus_val16 tmp_g = Q15ONE
665                         - MULT16_16_Q15(window[i], Q15ONE-ratio);
666                   buf[DECODE_BUFFER_SIZE-N+i] =
667                         MULT16_32_Q15(tmp_g, buf[DECODE_BUFFER_SIZE-N+i]);
668                }
669                for (i=overlap;i<extrapolation_len;i++)
670                {
671                   buf[DECODE_BUFFER_SIZE-N+i] =
672                         MULT16_32_Q15(ratio, buf[DECODE_BUFFER_SIZE-N+i]);
673                }
674             }
675          }
676
677          /* Apply the pre-filter to the MDCT overlap for the next frame because
678             the post-filter will be re-applied in the decoder after the MDCT
679             overlap. */
680          comb_filter(etmp, buf+DECODE_BUFFER_SIZE,
681               st->postfilter_period, st->postfilter_period, overlap,
682               -st->postfilter_gain, -st->postfilter_gain,
683               st->postfilter_tapset, st->postfilter_tapset, NULL, 0);
684
685          /* Simulate TDAC on the concealed audio so that it blends with the
686             MDCT of the next frame. */
687          for (i=0;i<overlap/2;i++)
688          {
689             buf[DECODE_BUFFER_SIZE+i] =
690                MULT16_32_Q15(window[i], etmp[overlap-1-i])
691                + MULT16_32_Q15(window[overlap-i-1], etmp[i]);
692          }
693       } while (++c<C);
694    }
695
696    deemphasis(out_syn, pcm, N, C, downsample,
697          mode->preemph, st->preemph_memD);
698
699    st->loss_count = loss_count+1;
700
701    RESTORE_STACK;
702 }
703
704 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)
705 {
706    int c, i, N;
707    int spread_decision;
708    opus_int32 bits;
709    ec_dec _dec;
710    VARDECL(celt_sig, freq);
711    VARDECL(celt_norm, X);
712    VARDECL(int, fine_quant);
713    VARDECL(int, pulses);
714    VARDECL(int, cap);
715    VARDECL(int, offsets);
716    VARDECL(int, fine_priority);
717    VARDECL(int, tf_res);
718    VARDECL(unsigned char, collapse_masks);
719    celt_sig *decode_mem[2];
720    celt_sig *out_syn[2];
721    opus_val16 *lpc;
722    opus_val16 *oldBandE, *oldLogE, *oldLogE2, *backgroundLogE;
723
724    int shortBlocks;
725    int isTransient;
726    int intra_ener;
727    const int CC = st->channels;
728    int LM, M;
729    int start;
730    int end;
731    int effEnd;
732    int codedBands;
733    int alloc_trim;
734    int postfilter_pitch;
735    opus_val16 postfilter_gain;
736    int intensity=0;
737    int dual_stereo=0;
738    opus_int32 total_bits;
739    opus_int32 balance;
740    opus_int32 tell;
741    int dynalloc_logp;
742    int postfilter_tapset;
743    int anti_collapse_rsv;
744    int anti_collapse_on=0;
745    int silence;
746    int C = st->stream_channels;
747    const OpusCustomMode *mode;
748    int nbEBands;
749    int overlap;
750    const opus_int16 *eBands;
751    ALLOC_STACK;
752
753    mode = st->mode;
754    nbEBands = mode->nbEBands;
755    overlap = mode->overlap;
756    eBands = mode->eBands;
757    start = st->start;
758    end = st->end;
759    frame_size *= st->downsample;
760
761    c=0; do {
762       decode_mem[c] = st->_decode_mem + c*(DECODE_BUFFER_SIZE+overlap);
763    } while (++c<CC);
764    lpc = (opus_val16*)(st->_decode_mem+(DECODE_BUFFER_SIZE+overlap)*CC);
765    oldBandE = lpc+CC*LPC_ORDER;
766    oldLogE = oldBandE + 2*nbEBands;
767    oldLogE2 = oldLogE + 2*nbEBands;
768    backgroundLogE = oldLogE2  + 2*nbEBands;
769
770 #ifdef CUSTOM_MODES
771    if (st->signalling && data!=NULL)
772    {
773       int data0=data[0];
774       /* Convert "standard mode" to Opus header */
775       if (mode->Fs==48000 && mode->shortMdctSize==120)
776       {
777          data0 = fromOpus(data0);
778          if (data0<0)
779             return OPUS_INVALID_PACKET;
780       }
781       st->end = end = IMAX(1, mode->effEBands-2*(data0>>5));
782       LM = (data0>>3)&0x3;
783       C = 1 + ((data0>>2)&0x1);
784       data++;
785       len--;
786       if (LM>mode->maxLM)
787          return OPUS_INVALID_PACKET;
788       if (frame_size < mode->shortMdctSize<<LM)
789          return OPUS_BUFFER_TOO_SMALL;
790       else
791          frame_size = mode->shortMdctSize<<LM;
792    } else {
793 #else
794    {
795 #endif
796       for (LM=0;LM<=mode->maxLM;LM++)
797          if (mode->shortMdctSize<<LM==frame_size)
798             break;
799       if (LM>mode->maxLM)
800          return OPUS_BAD_ARG;
801    }
802    M=1<<LM;
803
804    if (len<0 || len>1275 || pcm==NULL)
805       return OPUS_BAD_ARG;
806
807    N = M*mode->shortMdctSize;
808
809    effEnd = end;
810    if (effEnd > mode->effEBands)
811       effEnd = mode->effEBands;
812
813    if (data == NULL || len<=1)
814    {
815       celt_decode_lost(st, pcm, N, LM);
816       RESTORE_STACK;
817       return frame_size/st->downsample;
818    }
819
820    if (dec == NULL)
821    {
822       ec_dec_init(&_dec,(unsigned char*)data,len);
823       dec = &_dec;
824    }
825
826    if (C==1)
827    {
828       for (i=0;i<nbEBands;i++)
829          oldBandE[i]=MAX16(oldBandE[i],oldBandE[nbEBands+i]);
830    }
831
832    total_bits = len*8;
833    tell = ec_tell(dec);
834
835    if (tell >= total_bits)
836       silence = 1;
837    else if (tell==1)
838       silence = ec_dec_bit_logp(dec, 15);
839    else
840       silence = 0;
841    if (silence)
842    {
843       /* Pretend we've read all the remaining bits */
844       tell = len*8;
845       dec->nbits_total+=tell-ec_tell(dec);
846    }
847
848    postfilter_gain = 0;
849    postfilter_pitch = 0;
850    postfilter_tapset = 0;
851    if (start==0 && tell+16 <= total_bits)
852    {
853       if(ec_dec_bit_logp(dec, 1))
854       {
855          int qg, octave;
856          octave = ec_dec_uint(dec, 6);
857          postfilter_pitch = (16<<octave)+ec_dec_bits(dec, 4+octave)-1;
858          qg = ec_dec_bits(dec, 3);
859          if (ec_tell(dec)+2<=total_bits)
860             postfilter_tapset = ec_dec_icdf(dec, tapset_icdf, 2);
861          postfilter_gain = QCONST16(.09375f,15)*(qg+1);
862       }
863       tell = ec_tell(dec);
864    }
865
866    if (LM > 0 && tell+3 <= total_bits)
867    {
868       isTransient = ec_dec_bit_logp(dec, 3);
869       tell = ec_tell(dec);
870    }
871    else
872       isTransient = 0;
873
874    if (isTransient)
875       shortBlocks = M;
876    else
877       shortBlocks = 0;
878
879    /* Decode the global flags (first symbols in the stream) */
880    intra_ener = tell+3<=total_bits ? ec_dec_bit_logp(dec, 3) : 0;
881    /* Get band energies */
882    unquant_coarse_energy(mode, start, end, oldBandE,
883          intra_ener, dec, C, LM);
884
885    ALLOC(tf_res, nbEBands, int);
886    tf_decode(start, end, isTransient, tf_res, LM, dec);
887
888    tell = ec_tell(dec);
889    spread_decision = SPREAD_NORMAL;
890    if (tell+4 <= total_bits)
891       spread_decision = ec_dec_icdf(dec, spread_icdf, 5);
892
893    ALLOC(cap, nbEBands, int);
894
895    init_caps(mode,cap,LM,C);
896
897    ALLOC(offsets, nbEBands, int);
898
899    dynalloc_logp = 6;
900    total_bits<<=BITRES;
901    tell = ec_tell_frac(dec);
902    for (i=start;i<end;i++)
903    {
904       int width, quanta;
905       int dynalloc_loop_logp;
906       int boost;
907       width = C*(eBands[i+1]-eBands[i])<<LM;
908       /* quanta is 6 bits, but no more than 1 bit/sample
909          and no less than 1/8 bit/sample */
910       quanta = IMIN(width<<BITRES, IMAX(6<<BITRES, width));
911       dynalloc_loop_logp = dynalloc_logp;
912       boost = 0;
913       while (tell+(dynalloc_loop_logp<<BITRES) < total_bits && boost < cap[i])
914       {
915          int flag;
916          flag = ec_dec_bit_logp(dec, dynalloc_loop_logp);
917          tell = ec_tell_frac(dec);
918          if (!flag)
919             break;
920          boost += quanta;
921          total_bits -= quanta;
922          dynalloc_loop_logp = 1;
923       }
924       offsets[i] = boost;
925       /* Making dynalloc more likely */
926       if (boost>0)
927          dynalloc_logp = IMAX(2, dynalloc_logp-1);
928    }
929
930    ALLOC(fine_quant, nbEBands, int);
931    alloc_trim = tell+(6<<BITRES) <= total_bits ?
932          ec_dec_icdf(dec, trim_icdf, 7) : 5;
933
934    bits = (((opus_int32)len*8)<<BITRES) - ec_tell_frac(dec) - 1;
935    anti_collapse_rsv = isTransient&&LM>=2&&bits>=((LM+2)<<BITRES) ? (1<<BITRES) : 0;
936    bits -= anti_collapse_rsv;
937
938    ALLOC(pulses, nbEBands, int);
939    ALLOC(fine_priority, nbEBands, int);
940
941    codedBands = compute_allocation(mode, start, end, offsets, cap,
942          alloc_trim, &intensity, &dual_stereo, bits, &balance, pulses,
943          fine_quant, fine_priority, C, LM, dec, 0, 0, 0);
944
945    unquant_fine_energy(mode, start, end, oldBandE, fine_quant, dec, C);
946
947    /* Decode fixed codebook */
948    ALLOC(collapse_masks, C*nbEBands, unsigned char);
949    ALLOC(X, C*N, celt_norm);   /**< Interleaved normalised MDCTs */
950
951    quant_all_bands(0, mode, start, end, X, C==2 ? X+N : NULL, collapse_masks,
952          NULL, pulses, shortBlocks, spread_decision, dual_stereo, intensity, tf_res,
953          len*(8<<BITRES)-anti_collapse_rsv, balance, dec, LM, codedBands, &st->rng);
954
955    if (anti_collapse_rsv > 0)
956    {
957       anti_collapse_on = ec_dec_bits(dec, 1);
958    }
959
960    unquant_energy_finalise(mode, start, end, oldBandE,
961          fine_quant, fine_priority, len*8-ec_tell(dec), dec, C);
962
963    if (anti_collapse_on)
964       anti_collapse(mode, X, collapse_masks, LM, C, N,
965             start, end, oldBandE, oldLogE, oldLogE2, pulses, st->rng);
966
967    if (silence)
968    {
969       for (i=0;i<C*nbEBands;i++)
970          oldBandE[i] = -QCONST16(28.f,DB_SHIFT);
971    }
972
973    c=0; do {
974       OPUS_MOVE(decode_mem[c], decode_mem[c]+N, DECODE_BUFFER_SIZE-N+overlap/2);
975    } while (++c<CC);
976
977    c=0; do {
978       out_syn[c] = decode_mem[c]+DECODE_BUFFER_SIZE-N;
979    } while (++c<CC);
980
981    celt_synthesis(mode, X, out_syn, oldBandE, start, effEnd, C, CC, isTransient, LM, st->downsample, silence);
982
983    c=0; do {
984       st->postfilter_period=IMAX(st->postfilter_period, COMBFILTER_MINPERIOD);
985       st->postfilter_period_old=IMAX(st->postfilter_period_old, COMBFILTER_MINPERIOD);
986       comb_filter(out_syn[c], out_syn[c], st->postfilter_period_old, st->postfilter_period, mode->shortMdctSize,
987             st->postfilter_gain_old, st->postfilter_gain, st->postfilter_tapset_old, st->postfilter_tapset,
988             mode->window, overlap);
989       if (LM!=0)
990          comb_filter(out_syn[c]+mode->shortMdctSize, out_syn[c]+mode->shortMdctSize, st->postfilter_period, postfilter_pitch, N-mode->shortMdctSize,
991                st->postfilter_gain, postfilter_gain, st->postfilter_tapset, postfilter_tapset,
992                mode->window, overlap);
993
994    } while (++c<CC);
995    st->postfilter_period_old = st->postfilter_period;
996    st->postfilter_gain_old = st->postfilter_gain;
997    st->postfilter_tapset_old = st->postfilter_tapset;
998    st->postfilter_period = postfilter_pitch;
999    st->postfilter_gain = postfilter_gain;
1000    st->postfilter_tapset = postfilter_tapset;
1001    if (LM!=0)
1002    {
1003       st->postfilter_period_old = st->postfilter_period;
1004       st->postfilter_gain_old = st->postfilter_gain;
1005       st->postfilter_tapset_old = st->postfilter_tapset;
1006    }
1007
1008    if (C==1)
1009       OPUS_COPY(&oldBandE[nbEBands], oldBandE, nbEBands);
1010
1011    /* In case start or end were to change */
1012    if (!isTransient)
1013    {
1014       OPUS_COPY(oldLogE2, oldLogE, 2*nbEBands);
1015       OPUS_COPY(oldLogE, oldBandE, 2*nbEBands);
1016       for (i=0;i<2*nbEBands;i++)
1017          backgroundLogE[i] = MIN16(backgroundLogE[i] + M*QCONST16(0.001f,DB_SHIFT), oldBandE[i]);
1018    } else {
1019       for (i=0;i<2*nbEBands;i++)
1020          oldLogE[i] = MIN16(oldLogE[i], oldBandE[i]);
1021    }
1022    c=0; do
1023    {
1024       for (i=0;i<start;i++)
1025       {
1026          oldBandE[c*nbEBands+i]=0;
1027          oldLogE[c*nbEBands+i]=oldLogE2[c*nbEBands+i]=-QCONST16(28.f,DB_SHIFT);
1028       }
1029       for (i=end;i<nbEBands;i++)
1030       {
1031          oldBandE[c*nbEBands+i]=0;
1032          oldLogE[c*nbEBands+i]=oldLogE2[c*nbEBands+i]=-QCONST16(28.f,DB_SHIFT);
1033       }
1034    } while (++c<2);
1035    st->rng = dec->rng;
1036
1037    /* We reuse freq[] as scratch space for the de-emphasis */
1038    deemphasis(out_syn, pcm, N, CC, st->downsample, mode->preemph, st->preemph_memD);
1039    st->loss_count = 0;
1040    RESTORE_STACK;
1041    if (ec_tell(dec) > 8*len)
1042       return OPUS_INTERNAL_ERROR;
1043    if(ec_get_error(dec))
1044       st->error = 1;
1045    return frame_size/st->downsample;
1046 }
1047
1048
1049 #ifdef CUSTOM_MODES
1050
1051 #ifdef FIXED_POINT
1052 int opus_custom_decode(CELTDecoder * OPUS_RESTRICT st, const unsigned char *data, int len, opus_int16 * OPUS_RESTRICT pcm, int frame_size)
1053 {
1054    return celt_decode_with_ec(st, data, len, pcm, frame_size, NULL);
1055 }
1056
1057 #ifndef DISABLE_FLOAT_API
1058 int opus_custom_decode_float(CELTDecoder * OPUS_RESTRICT st, const unsigned char *data, int len, float * OPUS_RESTRICT pcm, int frame_size)
1059 {
1060    int j, ret, C, N;
1061    VARDECL(opus_int16, out);
1062    ALLOC_STACK;
1063
1064    if (pcm==NULL)
1065       return OPUS_BAD_ARG;
1066
1067    C = st->channels;
1068    N = frame_size;
1069
1070    ALLOC(out, C*N, opus_int16);
1071    ret=celt_decode_with_ec(st, data, len, out, frame_size, NULL);
1072    if (ret>0)
1073       for (j=0;j<C*ret;j++)
1074          pcm[j]=out[j]*(1.f/32768.f);
1075
1076    RESTORE_STACK;
1077    return ret;
1078 }
1079 #endif /* DISABLE_FLOAT_API */
1080
1081 #else
1082
1083 int opus_custom_decode_float(CELTDecoder * OPUS_RESTRICT st, const unsigned char *data, int len, float * OPUS_RESTRICT pcm, int frame_size)
1084 {
1085    return celt_decode_with_ec(st, data, len, pcm, frame_size, NULL);
1086 }
1087
1088 int opus_custom_decode(CELTDecoder * OPUS_RESTRICT st, const unsigned char *data, int len, opus_int16 * OPUS_RESTRICT pcm, int frame_size)
1089 {
1090    int j, ret, C, N;
1091    VARDECL(celt_sig, out);
1092    ALLOC_STACK;
1093
1094    if (pcm==NULL)
1095       return OPUS_BAD_ARG;
1096
1097    C = st->channels;
1098    N = frame_size;
1099    ALLOC(out, C*N, celt_sig);
1100
1101    ret=celt_decode_with_ec(st, data, len, out, frame_size, NULL);
1102
1103    if (ret>0)
1104       for (j=0;j<C*ret;j++)
1105          pcm[j] = FLOAT2INT16 (out[j]);
1106
1107    RESTORE_STACK;
1108    return ret;
1109 }
1110
1111 #endif
1112 #endif /* CUSTOM_MODES */
1113
1114 int opus_custom_decoder_ctl(CELTDecoder * OPUS_RESTRICT st, int request, ...)
1115 {
1116    va_list ap;
1117
1118    va_start(ap, request);
1119    switch (request)
1120    {
1121       case CELT_SET_START_BAND_REQUEST:
1122       {
1123          opus_int32 value = va_arg(ap, opus_int32);
1124          if (value<0 || value>=st->mode->nbEBands)
1125             goto bad_arg;
1126          st->start = value;
1127       }
1128       break;
1129       case CELT_SET_END_BAND_REQUEST:
1130       {
1131          opus_int32 value = va_arg(ap, opus_int32);
1132          if (value<1 || value>st->mode->nbEBands)
1133             goto bad_arg;
1134          st->end = value;
1135       }
1136       break;
1137       case CELT_SET_CHANNELS_REQUEST:
1138       {
1139          opus_int32 value = va_arg(ap, opus_int32);
1140          if (value<1 || value>2)
1141             goto bad_arg;
1142          st->stream_channels = value;
1143       }
1144       break;
1145       case CELT_GET_AND_CLEAR_ERROR_REQUEST:
1146       {
1147          opus_int32 *value = va_arg(ap, opus_int32*);
1148          if (value==NULL)
1149             goto bad_arg;
1150          *value=st->error;
1151          st->error = 0;
1152       }
1153       break;
1154       case OPUS_GET_LOOKAHEAD_REQUEST:
1155       {
1156          opus_int32 *value = va_arg(ap, opus_int32*);
1157          if (value==NULL)
1158             goto bad_arg;
1159          *value = st->overlap/st->downsample;
1160       }
1161       break;
1162       case OPUS_RESET_STATE:
1163       {
1164          int i;
1165          opus_val16 *lpc, *oldBandE, *oldLogE, *oldLogE2;
1166          lpc = (opus_val16*)(st->_decode_mem+(DECODE_BUFFER_SIZE+st->overlap)*st->channels);
1167          oldBandE = lpc+st->channels*LPC_ORDER;
1168          oldLogE = oldBandE + 2*st->mode->nbEBands;
1169          oldLogE2 = oldLogE + 2*st->mode->nbEBands;
1170          OPUS_CLEAR((char*)&st->DECODER_RESET_START,
1171                opus_custom_decoder_get_size(st->mode, st->channels)-
1172                ((char*)&st->DECODER_RESET_START - (char*)st));
1173          for (i=0;i<2*st->mode->nbEBands;i++)
1174             oldLogE[i]=oldLogE2[i]=-QCONST16(28.f,DB_SHIFT);
1175       }
1176       break;
1177       case OPUS_GET_PITCH_REQUEST:
1178       {
1179          opus_int32 *value = va_arg(ap, opus_int32*);
1180          if (value==NULL)
1181             goto bad_arg;
1182          *value = st->postfilter_period;
1183       }
1184       break;
1185       case CELT_GET_MODE_REQUEST:
1186       {
1187          const CELTMode ** value = va_arg(ap, const CELTMode**);
1188          if (value==0)
1189             goto bad_arg;
1190          *value=st->mode;
1191       }
1192       break;
1193       case CELT_SET_SIGNALLING_REQUEST:
1194       {
1195          opus_int32 value = va_arg(ap, opus_int32);
1196          st->signalling = value;
1197       }
1198       break;
1199       case OPUS_GET_FINAL_RANGE_REQUEST:
1200       {
1201          opus_uint32 * value = va_arg(ap, opus_uint32 *);
1202          if (value==0)
1203             goto bad_arg;
1204          *value=st->rng;
1205       }
1206       break;
1207       default:
1208          goto bad_request;
1209    }
1210    va_end(ap);
1211    return OPUS_OK;
1212 bad_arg:
1213    va_end(ap);
1214    return OPUS_BAD_ARG;
1215 bad_request:
1216       va_end(ap);
1217   return OPUS_UNIMPLEMENTED;
1218 }