Further cleanup of the MDCT code, fixes PLC bug
[opus.git] / celt / celt.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_C
35
36 #include "os_support.h"
37 #include "mdct.h"
38 #include <math.h>
39 #include "celt.h"
40 #include "pitch.h"
41 #include "bands.h"
42 #include "modes.h"
43 #include "entcode.h"
44 #include "quant_bands.h"
45 #include "rate.h"
46 #include "stack_alloc.h"
47 #include "mathops.h"
48 #include "float_cast.h"
49 #include <stdarg.h>
50 #include "celt_lpc.h"
51 #include "vq.h"
52
53 #ifndef OPUS_VERSION
54 #define OPUS_VERSION "unknown"
55 #endif
56
57 #ifdef CUSTOM_MODES
58 #define OPUS_CUSTOM_NOSTATIC
59 #else
60 #define OPUS_CUSTOM_NOSTATIC static inline
61 #endif
62
63 static const unsigned char trim_icdf[11] = {126, 124, 119, 109, 87, 41, 19, 9, 4, 2, 0};
64 /* Probs: NONE: 21.875%, LIGHT: 6.25%, NORMAL: 65.625%, AGGRESSIVE: 6.25% */
65 static const unsigned char spread_icdf[4] = {25, 23, 2, 0};
66
67 static const unsigned char tapset_icdf[3]={2,1,0};
68
69 #ifdef CUSTOM_MODES
70 static const unsigned char toOpusTable[20] = {
71       0xE0, 0xE8, 0xF0, 0xF8,
72       0xC0, 0xC8, 0xD0, 0xD8,
73       0xA0, 0xA8, 0xB0, 0xB8,
74       0x00, 0x00, 0x00, 0x00,
75       0x80, 0x88, 0x90, 0x98,
76 };
77
78 static const unsigned char fromOpusTable[16] = {
79       0x80, 0x88, 0x90, 0x98,
80       0x40, 0x48, 0x50, 0x58,
81       0x20, 0x28, 0x30, 0x38,
82       0x00, 0x08, 0x10, 0x18
83 };
84
85 static inline int toOpus(unsigned char c)
86 {
87    int ret=0;
88    if (c<0xA0)
89       ret = toOpusTable[c>>3];
90    if (ret == 0)
91       return -1;
92    else
93       return ret|(c&0x7);
94 }
95
96 static inline int fromOpus(unsigned char c)
97 {
98    if (c<0x80)
99       return -1;
100    else
101       return fromOpusTable[(c>>3)-16] | (c&0x7);
102 }
103 #endif /* CUSTOM_MODES */
104
105 #define COMBFILTER_MAXPERIOD 1024
106 #define COMBFILTER_MINPERIOD 15
107
108 static int resampling_factor(opus_int32 rate)
109 {
110    int ret;
111    switch (rate)
112    {
113    case 48000:
114       ret = 1;
115       break;
116    case 24000:
117       ret = 2;
118       break;
119    case 16000:
120       ret = 3;
121       break;
122    case 12000:
123       ret = 4;
124       break;
125    case 8000:
126       ret = 6;
127       break;
128    default:
129 #ifndef CUSTOM_MODES
130       celt_assert(0);
131 #endif
132       ret = 0;
133       break;
134    }
135    return ret;
136 }
137
138 /** Encoder state
139  @brief Encoder state
140  */
141 struct OpusCustomEncoder {
142    const OpusCustomMode *mode;     /**< Mode used by the encoder */
143    int overlap;
144    int channels;
145    int stream_channels;
146
147    int force_intra;
148    int clip;
149    int disable_pf;
150    int complexity;
151    int upsample;
152    int start, end;
153
154    opus_int32 bitrate;
155    int vbr;
156    int signalling;
157    int constrained_vbr;      /* If zero, VBR can do whatever it likes with the rate */
158    int loss_rate;
159    int lsb_depth;
160
161    /* Everything beyond this point gets cleared on a reset */
162 #define ENCODER_RESET_START rng
163
164    opus_uint32 rng;
165    int spread_decision;
166    opus_val32 delayedIntra;
167    int tonal_average;
168    int lastCodedBands;
169    int hf_average;
170    int tapset_decision;
171
172    int prefilter_period;
173    opus_val16 prefilter_gain;
174    int prefilter_tapset;
175 #ifdef RESYNTH
176    int prefilter_period_old;
177    opus_val16 prefilter_gain_old;
178    int prefilter_tapset_old;
179 #endif
180    int consec_transient;
181    AnalysisInfo analysis;
182
183    opus_val32 preemph_memE[2];
184    opus_val32 preemph_memD[2];
185
186    /* VBR-related parameters */
187    opus_int32 vbr_reservoir;
188    opus_int32 vbr_drift;
189    opus_int32 vbr_offset;
190    opus_int32 vbr_count;
191    opus_val16 overlap_max;
192    opus_val16 stereo_saving;
193    int intensity;
194
195 #ifdef RESYNTH
196    /* +MAX_PERIOD/2 to make space for overlap */
197    celt_sig syn_mem[2][2*MAX_PERIOD+MAX_PERIOD/2];
198 #endif
199
200    celt_sig in_mem[1]; /* Size = channels*mode->overlap */
201    /* celt_sig prefilter_mem[],  Size = channels*COMBFILTER_MAXPERIOD */
202    /* opus_val16 oldBandE[],     Size = channels*mode->nbEBands */
203    /* opus_val16 oldLogE[],      Size = channels*mode->nbEBands */
204    /* opus_val16 oldLogE2[],     Size = channels*mode->nbEBands */
205 };
206
207 int celt_encoder_get_size(int channels)
208 {
209    CELTMode *mode = opus_custom_mode_create(48000, 960, NULL);
210    return opus_custom_encoder_get_size(mode, channels);
211 }
212
213 OPUS_CUSTOM_NOSTATIC int opus_custom_encoder_get_size(const CELTMode *mode, int channels)
214 {
215    int size = sizeof(struct CELTEncoder)
216          + (channels*mode->overlap-1)*sizeof(celt_sig)    /* celt_sig in_mem[channels*mode->overlap]; */
217          + channels*COMBFILTER_MAXPERIOD*sizeof(celt_sig) /* celt_sig prefilter_mem[channels*COMBFILTER_MAXPERIOD]; */
218          + 3*channels*mode->nbEBands*sizeof(opus_val16);  /* opus_val16 oldBandE[channels*mode->nbEBands]; */
219                                                           /* opus_val16 oldLogE[channels*mode->nbEBands]; */
220                                                           /* opus_val16 oldLogE2[channels*mode->nbEBands]; */
221    return size;
222 }
223
224 #ifdef CUSTOM_MODES
225 CELTEncoder *opus_custom_encoder_create(const CELTMode *mode, int channels, int *error)
226 {
227    int ret;
228    CELTEncoder *st = (CELTEncoder *)opus_alloc(opus_custom_encoder_get_size(mode, channels));
229    /* init will handle the NULL case */
230    ret = opus_custom_encoder_init(st, mode, channels);
231    if (ret != OPUS_OK)
232    {
233       opus_custom_encoder_destroy(st);
234       st = NULL;
235    }
236    if (error)
237       *error = ret;
238    return st;
239 }
240 #endif /* CUSTOM_MODES */
241
242 int celt_encoder_init(CELTEncoder *st, opus_int32 sampling_rate, int channels)
243 {
244    int ret;
245    ret = opus_custom_encoder_init(st, opus_custom_mode_create(48000, 960, NULL), channels);
246    if (ret != OPUS_OK)
247       return ret;
248    st->upsample = resampling_factor(sampling_rate);
249    return OPUS_OK;
250 }
251
252 OPUS_CUSTOM_NOSTATIC int opus_custom_encoder_init(CELTEncoder *st, const CELTMode *mode, int channels)
253 {
254    if (channels < 0 || channels > 2)
255       return OPUS_BAD_ARG;
256
257    if (st==NULL || mode==NULL)
258       return OPUS_ALLOC_FAIL;
259
260    OPUS_CLEAR((char*)st, opus_custom_encoder_get_size(mode, channels));
261
262    st->mode = mode;
263    st->overlap = mode->overlap;
264    st->stream_channels = st->channels = channels;
265
266    st->upsample = 1;
267    st->start = 0;
268    st->end = st->mode->effEBands;
269    st->signalling = 1;
270
271    st->constrained_vbr = 1;
272    st->clip = 1;
273
274    st->bitrate = OPUS_BITRATE_MAX;
275    st->vbr = 0;
276    st->force_intra  = 0;
277    st->complexity = 5;
278    st->lsb_depth=24;
279
280    opus_custom_encoder_ctl(st, OPUS_RESET_STATE);
281
282    return OPUS_OK;
283 }
284
285 #ifdef CUSTOM_MODES
286 void opus_custom_encoder_destroy(CELTEncoder *st)
287 {
288    opus_free(st);
289 }
290 #endif /* CUSTOM_MODES */
291
292 static inline opus_val16 SIG2WORD16(celt_sig x)
293 {
294 #ifdef FIXED_POINT
295    x = PSHR32(x, SIG_SHIFT);
296    x = MAX32(x, -32768);
297    x = MIN32(x, 32767);
298    return EXTRACT16(x);
299 #else
300    return (opus_val16)x;
301 #endif
302 }
303
304 static int transient_analysis(const opus_val32 * OPUS_RESTRICT in, int len, int C,
305                               opus_val16 *tf_estimate, int *tf_chan)
306 {
307    int i;
308    VARDECL(opus_val16, tmp);
309    opus_val32 mem0,mem1;
310    int is_transient = 0;
311    opus_int32 mask_metric = 0;
312    int c;
313    int tf_max;
314    /* Table of 6*64/x, trained on real data to minimize the average error */
315    static const unsigned char inv_table[128] = {
316          255,255,156,110, 86, 70, 59, 51, 45, 40, 37, 33, 31, 28, 26, 25,
317           23, 22, 21, 20, 19, 18, 17, 16, 16, 15, 15, 14, 13, 13, 12, 12,
318           12, 12, 11, 11, 11, 10, 10, 10,  9,  9,  9,  9,  9,  9,  8,  8,
319            8,  8,  8,  7,  7,  7,  7,  7,  7,  6,  6,  6,  6,  6,  6,  6,
320            6,  6,  6,  6,  6,  6,  6,  6,  6,  5,  5,  5,  5,  5,  5,  5,
321            5,  5,  5,  5,  5,  4,  4,  4,  4,  4,  4,  4,  4,  4,  4,  4,
322            4,  4,  4,  4,  4,  4,  4,  4,  4,  4,  4,  4,  4,  4,  3,  3,
323            3,  3,  3,  3,  3,  3,  3,  3,  3,  3,  3,  3,  3,  3,  3,  2,
324    };
325    SAVE_STACK;
326    ALLOC(tmp, len, opus_val16);
327
328    tf_max = 0;
329    for (c=0;c<C;c++)
330    {
331       opus_val32 mean;
332       opus_int32 unmask=0;
333       opus_val32 norm;
334       mem0=0;
335       mem1=0;
336       /* High-pass filter: (1 - 2*z^-1 + z^-2) / (1 - z^-1 + .5*z^-2) */
337       for (i=0;i<len;i++)
338       {
339          opus_val32 x,y;
340          x = SHR32(in[i+c*len],SIG_SHIFT);
341          y = ADD32(mem0, x);
342 #ifdef FIXED_POINT
343          mem0 = mem1 + y - SHL32(x,1);
344          mem1 = x - SHR32(y,1);
345 #else
346          mem0 = mem1 + y - 2*x;
347          mem1 = x - .5f*y;
348 #endif
349          tmp[i] = EXTRACT16(SHR32(y,2));
350          /*printf("%f ", tmp[i]);*/
351       }
352       /*printf("\n");*/
353       /* First few samples are bad because we don't propagate the memory */
354       for (i=0;i<12;i++)
355          tmp[i] = 0;
356
357 #ifdef FIXED_POINT
358       /* Normalize tmp to max range */
359       {
360          int shift=0;
361          shift = 14-celt_ilog2(1+celt_maxabs16(tmp, len));
362          if (shift!=0)
363          {
364             for (i=0;i<len;i++)
365                tmp[i] = SHL16(tmp[i], shift);
366          }
367       }
368 #endif
369
370       mean=0;
371       mem0=0;
372       /*  Grouping by two to reduce complexity */
373       len/=2;
374       /* Forward pass to compute the post-echo threshold*/
375       for (i=0;i<len;i++)
376       {
377          opus_val16 x2 = PSHR32(MULT16_16(tmp[2*i],tmp[2*i]) + MULT16_16(tmp[2*i+1],tmp[2*i+1]),16);
378          mean += x2;
379 #ifdef FIXED_POINT
380          /* FIXME: Use PSHR16() instead */
381          tmp[i] = mem0 + PSHR32(x2-mem0,4);
382 #else
383          tmp[i] = mem0 + MULT16_16_P15(QCONST16(.0625f,15),x2-mem0);
384 #endif
385          mem0 = tmp[i];
386       }
387
388       mem0=0;
389       /* Backward pass to compute the pre-echo threshold */
390       for (i=len-1;i>=0;i--)
391       {
392 #ifdef FIXED_POINT
393          /* FIXME: Use PSHR16() instead */
394          tmp[i] = mem0 + PSHR32(tmp[i]-mem0,3);
395 #else
396          tmp[i] = mem0 + MULT16_16_P15(QCONST16(0.125f,15),tmp[i]-mem0);
397 #endif
398          mem0 = tmp[i];
399       }
400       /*for (i=0;i<len;i++)printf("%f ", tmp[i]/mean);printf("\n");*/
401
402       /* Compute the ratio of the mean energy over the harmonic mean of the energy.
403          This essentially corresponds to a bitrate-normalized temporal noise-to-mask
404          ratio */
405
406       /* Inverse of the mean energy in Q15+6 */
407       norm = SHL32(EXTEND32(len),6+14)/ADD32(EPSILON,SHR32(mean,1));
408       /* Compute harmonic mean discarding the unreliable boundaries
409          The data is smooth, so we only take 1/4th of the samples */
410       unmask=0;
411       for (i=12;i<len-5;i+=4)
412       {
413          int id;
414 #ifdef FIXED_POINT
415          id = IMAX(0,IMIN(127,MULT16_32_Q15(tmp[i],norm))); /* Do not round to nearest */
416 #else
417          id = IMAX(0,IMIN(127,floor(64*norm*tmp[i]))); /* Do not round to nearest */
418 #endif
419          unmask += inv_table[id];
420       }
421       /*printf("%d\n", unmask);*/
422       /* Normalize, compensate for the 1/4th of the sample and the factor of 6 in the inverse table */
423       unmask = 64*unmask*4/(6*(len-17));
424       if (unmask>mask_metric)
425       {
426          *tf_chan = c;
427          mask_metric = unmask;
428       }
429    }
430    is_transient = mask_metric>141;
431
432    /* Arbitrary metric for VBR boost */
433    tf_max = MAX16(0,celt_sqrt(64*mask_metric)-64);
434    /* *tf_estimate = 1 + MIN16(1, sqrt(MAX16(0, tf_max-30))/20); */
435    *tf_estimate = QCONST16(1.f, 14) + celt_sqrt(MAX16(0, SHL32(MULT16_16(QCONST16(0.0069,14),IMIN(163,tf_max)),14)-QCONST32(0.139,28)));
436    /*printf("%d %f\n", tf_max, mask_metric);*/
437    RESTORE_STACK;
438 #ifdef FUZZING
439    is_transient = rand()&0x1;
440 #endif
441    /*printf("%d %f %d\n", is_transient, (float)*tf_estimate, tf_max);*/
442    return is_transient;
443 }
444
445 /** Apply window and compute the MDCT for all sub-frames and
446     all channels in a frame */
447 static void compute_mdcts(const CELTMode *mode, int shortBlocks, celt_sig * OPUS_RESTRICT in, celt_sig * OPUS_RESTRICT out, int C, int LM)
448 {
449    const int overlap = OVERLAP(mode);
450    int N;
451    int B;
452    int shift;
453    int b, c;
454    if (shortBlocks)
455    {
456       B = shortBlocks;
457       N = mode->shortMdctSize;
458       shift = mode->maxLM;
459    } else {
460       B = 1;
461       N = mode->shortMdctSize<<LM;
462       shift = mode->maxLM-LM;
463    }
464    c=0; do {
465       for (b=0;b<B;b++)
466       {
467          /* Interleaving the sub-frames while doing the MDCTs */
468          clt_mdct_forward(&mode->mdct, in+c*(B*N+overlap)+b*N, &out[b+c*N*B], mode->window, overlap, shift, B);
469       }
470    } while (++c<C);
471 }
472
473 /** Compute the IMDCT and apply window for all sub-frames and
474     all channels in a frame */
475 static void compute_inv_mdcts(const CELTMode *mode, int shortBlocks, celt_sig *X,
476       celt_sig * OPUS_RESTRICT out_mem[], int C, int LM)
477 {
478    int b, c;
479    int B;
480    int N;
481    int shift;
482    const int overlap = OVERLAP(mode);
483
484    if (shortBlocks)
485    {
486       B = shortBlocks;
487       N = mode->shortMdctSize;
488       shift = mode->maxLM;
489    } else {
490       B = 1;
491       N = mode->shortMdctSize<<LM;
492       shift = mode->maxLM-LM;
493    }
494    c=0; do {
495       /* IMDCT on the interleaved the sub-frames, overlap-add is performed by the IMDCT */
496       for (b=0;b<B;b++)
497          clt_mdct_backward(&mode->mdct, &X[b+c*N*B], out_mem[c]+N*b, mode->window, overlap, shift, B);
498    } while (++c<C);
499 }
500
501 static void preemphasis(const opus_val16 * OPUS_RESTRICT pcmp, celt_sig * OPUS_RESTRICT inp,
502                         int N, int CC, int upsample, const opus_val16 *coef, celt_sig *mem, int clip)
503 {
504    int i;
505    opus_val16 coef0, coef1;
506    celt_sig m;
507    int Nu;
508
509    coef0 = coef[0];
510    coef1 = coef[1];
511
512
513    Nu = N/upsample;
514    if (upsample!=1)
515    {
516       for (i=0;i<N;i++)
517          inp[i] = 0;
518    }
519    for (i=0;i<Nu;i++)
520    {
521       celt_sig x;
522
523       x = SCALEIN(pcmp[CC*i]);
524 #ifndef FIXED_POINT
525       /* Replace NaNs with zeros */
526       if (!(x==x))
527          x = 0;
528 #endif
529       inp[i*upsample] = x;
530    }
531
532 #ifndef FIXED_POINT
533    if (clip)
534    {
535       /* Clip input to avoid encoding non-portable files */
536       for (i=0;i<Nu;i++)
537          inp[i*upsample] = MAX32(-65536.f, MIN32(65536.f,inp[i*upsample]));
538    }
539 #endif
540    m = *mem;
541    if (coef1 == 0)
542    {
543       for (i=0;i<N;i++)
544       {
545          celt_sig x;
546          x = SHL32(inp[i], SIG_SHIFT);
547          /* Apply pre-emphasis */
548          inp[i] = x + m;
549          m = - MULT16_32_Q15(coef0, x);
550       }
551    } else {
552       opus_val16 coef2 = coef[2];
553       for (i=0;i<N;i++)
554       {
555          opus_val16 x, tmp;
556          x = inp[i];
557          /* Apply pre-emphasis */
558          tmp = MULT16_16(coef2, x);
559          inp[i] = tmp + m;
560          m = MULT16_32_Q15(coef1, inp[i]) - MULT16_32_Q15(coef0, tmp);
561       }
562    }
563    *mem = m;
564 }
565
566 static void deemphasis(celt_sig *in[], opus_val16 *pcm, int N, int C, int downsample, const opus_val16 *coef, celt_sig *mem, celt_sig * OPUS_RESTRICT scratch)
567 {
568    int c;
569    int Nd;
570    opus_val16 coef0, coef1;
571
572    coef0 = coef[0];
573    coef1 = coef[1];
574    Nd = N/downsample;
575    c=0; do {
576       int j;
577       celt_sig * OPUS_RESTRICT x;
578       opus_val16  * OPUS_RESTRICT y;
579       celt_sig m = mem[c];
580       x =in[c];
581       y = pcm+c;
582       /* Shortcut for the standard (non-custom modes) case */
583       if (coef1 == 0)
584       {
585          for (j=0;j<N;j++)
586          {
587             celt_sig tmp = x[j] + m;
588             m = MULT16_32_Q15(coef0, tmp);
589             scratch[j] = tmp;
590          }
591       } else {
592          opus_val16 coef3 = coef[3];
593          for (j=0;j<N;j++)
594          {
595             celt_sig tmp = x[j] + m;
596             m = MULT16_32_Q15(coef0, tmp)
597               - MULT16_32_Q15(coef1, x[j]);
598             tmp = SHL32(MULT16_32_Q15(coef3, tmp), 2);
599             scratch[j] = tmp;
600          }
601       }
602       mem[c] = m;
603
604       /* Perform down-sampling */
605       for (j=0;j<Nd;j++)
606          y[j*C] = SCALEOUT(SIG2WORD16(scratch[j*downsample]));
607    } while (++c<C);
608 }
609
610 static void comb_filter(opus_val32 *y, opus_val32 *x, int T0, int T1, int N,
611       opus_val16 g0, opus_val16 g1, int tapset0, int tapset1,
612       const opus_val16 *window, int overlap)
613 {
614    int i;
615    /* printf ("%d %d %f %f\n", T0, T1, g0, g1); */
616    opus_val16 g00, g01, g02, g10, g11, g12;
617    opus_val32 x0, x1, x2, x3, x4;
618    static const opus_val16 gains[3][3] = {
619          {QCONST16(0.3066406250f, 15), QCONST16(0.2170410156f, 15), QCONST16(0.1296386719f, 15)},
620          {QCONST16(0.4638671875f, 15), QCONST16(0.2680664062f, 15), QCONST16(0.f, 15)},
621          {QCONST16(0.7998046875f, 15), QCONST16(0.1000976562f, 15), QCONST16(0.f, 15)}};
622
623    if (g0==0 && g1==0)
624    {
625       /* OPT: Happens to work without the OPUS_MOVE(), but only because the current encoder already copies x to y */
626       if (x!=y)
627          OPUS_MOVE(y, x, N);
628       return;
629    }
630    g00 = MULT16_16_Q15(g0, gains[tapset0][0]);
631    g01 = MULT16_16_Q15(g0, gains[tapset0][1]);
632    g02 = MULT16_16_Q15(g0, gains[tapset0][2]);
633    g10 = MULT16_16_Q15(g1, gains[tapset1][0]);
634    g11 = MULT16_16_Q15(g1, gains[tapset1][1]);
635    g12 = MULT16_16_Q15(g1, gains[tapset1][2]);
636    x1 = x[-T1+1];
637    x2 = x[-T1  ];
638    x3 = x[-T1-1];
639    x4 = x[-T1-2];
640    for (i=0;i<overlap;i++)
641    {
642       opus_val16 f;
643       x0=x[i-T1+2];
644       f = MULT16_16_Q15(window[i],window[i]);
645       y[i] = x[i]
646                + MULT16_32_Q15(MULT16_16_Q15((Q15ONE-f),g00),x[i-T0])
647                + MULT16_32_Q15(MULT16_16_Q15((Q15ONE-f),g01),ADD32(x[i-T0+1],x[i-T0-1]))
648                + MULT16_32_Q15(MULT16_16_Q15((Q15ONE-f),g02),ADD32(x[i-T0+2],x[i-T0-2]))
649                + MULT16_32_Q15(MULT16_16_Q15(f,g10),x2)
650                + MULT16_32_Q15(MULT16_16_Q15(f,g11),ADD32(x1,x3))
651                + MULT16_32_Q15(MULT16_16_Q15(f,g12),ADD32(x0,x4));
652       x4=x3;
653       x3=x2;
654       x2=x1;
655       x1=x0;
656
657    }
658    if (g1==0)
659    {
660       /* OPT: Happens to work without the OPUS_MOVE(), but only because the current encoder already copies x to y */
661       if (x!=y)
662          OPUS_MOVE(y+overlap, x+overlap, N-overlap);
663       return;
664    }
665    /* OPT: For machines where the movs are costly, unroll by 5 */
666    for (;i<N;i++)
667    {
668       x0=x[i-T1+2];
669       y[i] = x[i]
670                + MULT16_32_Q15(g10,x2)
671                + MULT16_32_Q15(g11,ADD32(x1,x3))
672                + MULT16_32_Q15(g12,ADD32(x0,x4));
673       x4=x3;
674       x3=x2;
675       x2=x1;
676       x1=x0;
677    }
678 }
679
680 static const signed char tf_select_table[4][8] = {
681       {0, -1, 0, -1,    0,-1, 0,-1},
682       {0, -1, 0, -2,    1, 0, 1,-1},
683       {0, -2, 0, -3,    2, 0, 1,-1},
684       {0, -2, 0, -3,    3, 0, 1,-1},
685 };
686
687 static opus_val32 l1_metric(const celt_norm *tmp, int N, int LM, opus_val16 bias)
688 {
689    int i;
690    opus_val32 L1;
691    L1 = 0;
692    for (i=0;i<N;i++)
693       L1 += EXTEND32(ABS16(tmp[i]));
694    /* When in doubt, prefer good freq resolution */
695    L1 = MAC16_32_Q15(L1, LM*bias, L1);
696    return L1;
697
698 }
699
700 static int tf_analysis(const CELTMode *m, int len, int C, int isTransient,
701       int *tf_res, int nbCompressedBytes, celt_norm *X, int N0, int LM,
702       int *tf_sum, opus_val16 tf_estimate, int tf_chan)
703 {
704    int i;
705    VARDECL(int, metric);
706    int cost0;
707    int cost1;
708    VARDECL(int, path0);
709    VARDECL(int, path1);
710    VARDECL(celt_norm, tmp);
711    VARDECL(celt_norm, tmp_1);
712    int lambda;
713    int sel;
714    int selcost[2];
715    int tf_select=0;
716    opus_val16 bias;
717
718    SAVE_STACK;
719    bias = MULT16_16_Q14(QCONST16(.04f,15), MAX16(-QCONST16(.25f,14), QCONST16(1.5f,14)-tf_estimate));
720    /*printf("%f ", bias);*/
721
722    if (nbCompressedBytes<15*C)
723    {
724       *tf_sum = 0;
725       for (i=0;i<len;i++)
726          tf_res[i] = isTransient;
727       return 0;
728    }
729    if (nbCompressedBytes<40)
730       lambda = 12;
731    else if (nbCompressedBytes<60)
732       lambda = 6;
733    else if (nbCompressedBytes<100)
734       lambda = 4;
735    else
736       lambda = 3;
737    lambda*=2;
738    ALLOC(metric, len, int);
739    ALLOC(tmp, (m->eBands[len]-m->eBands[len-1])<<LM, celt_norm);
740    ALLOC(tmp_1, (m->eBands[len]-m->eBands[len-1])<<LM, celt_norm);
741    ALLOC(path0, len, int);
742    ALLOC(path1, len, int);
743
744    *tf_sum = 0;
745    for (i=0;i<len;i++)
746    {
747       int j, k, N;
748       int narrow;
749       opus_val32 L1, best_L1;
750       int best_level=0;
751       N = (m->eBands[i+1]-m->eBands[i])<<LM;
752       /* band is too narrow to be split down to LM=-1 */
753       narrow = (m->eBands[i+1]-m->eBands[i])==1;
754       for (j=0;j<N;j++)
755          tmp[j] = X[tf_chan*N0 + j+(m->eBands[i]<<LM)];
756       /* Just add the right channel if we're in stereo */
757       /*if (C==2)
758          for (j=0;j<N;j++)
759             tmp[j] = ADD16(SHR16(tmp[j], 1),SHR16(X[N0+j+(m->eBands[i]<<LM)], 1));*/
760       L1 = l1_metric(tmp, N, isTransient ? LM : 0, bias);
761       best_L1 = L1;
762       /* Check the -1 case for transients */
763       if (isTransient && !narrow)
764       {
765          for (j=0;j<N;j++)
766             tmp_1[j] = tmp[j];
767          haar1(tmp_1, N>>LM, 1<<LM);
768          L1 = l1_metric(tmp_1, N, LM+1, bias);
769          if (L1<best_L1)
770          {
771             best_L1 = L1;
772             best_level = -1;
773          }
774       }
775       /*printf ("%f ", L1);*/
776       for (k=0;k<LM+!(isTransient||narrow);k++)
777       {
778          int B;
779
780          if (isTransient)
781             B = (LM-k-1);
782          else
783             B = k+1;
784
785          haar1(tmp, N>>k, 1<<k);
786
787          L1 = l1_metric(tmp, N, B, bias);
788
789          if (L1 < best_L1)
790          {
791             best_L1 = L1;
792             best_level = k+1;
793          }
794       }
795       /*printf ("%d ", isTransient ? LM-best_level : best_level);*/
796       /* metric is in Q1 to be able to select the mid-point (-0.5) for narrower bands */
797       if (isTransient)
798          metric[i] = 2*best_level;
799       else
800          metric[i] = -2*best_level;
801       *tf_sum += (isTransient ? LM : 0) - metric[i]/2;
802       /* For bands that can't be split to -1, set the metric to the half-way point to avoid
803          biasing the decision */
804       if (narrow && (metric[i]==0 || metric[i]==-2*LM))
805          metric[i]-=1;
806       /*printf("%d ", metric[i]);*/
807    }
808    /*printf("\n");*/
809    /* Search for the optimal tf resolution, including tf_select */
810    tf_select = 0;
811    for (sel=0;sel<2;sel++)
812    {
813       cost0 = 0;
814       cost1 = isTransient ? 0 : lambda;
815       for (i=1;i<len;i++)
816       {
817          int curr0, curr1;
818          curr0 = IMIN(cost0, cost1 + lambda);
819          curr1 = IMIN(cost0 + lambda, cost1);
820          cost0 = curr0 + abs(metric[i]-2*tf_select_table[LM][4*isTransient+2*sel+0]);
821          cost1 = curr1 + abs(metric[i]-2*tf_select_table[LM][4*isTransient+2*sel+1]);
822       }
823       cost0 = IMIN(cost0, cost1);
824       selcost[sel]=cost0;
825    }
826    /* For now, we're conservative and only allow tf_select=1 for transients.
827     * If tests confirm it's useful for non-transients, we could allow it. */
828    if (selcost[1]<selcost[0] && isTransient)
829       tf_select=1;
830    cost0 = 0;
831    cost1 = isTransient ? 0 : lambda;
832    /* Viterbi forward pass */
833    for (i=1;i<len;i++)
834    {
835       int curr0, curr1;
836       int from0, from1;
837
838       from0 = cost0;
839       from1 = cost1 + lambda;
840       if (from0 < from1)
841       {
842          curr0 = from0;
843          path0[i]= 0;
844       } else {
845          curr0 = from1;
846          path0[i]= 1;
847       }
848
849       from0 = cost0 + lambda;
850       from1 = cost1;
851       if (from0 < from1)
852       {
853          curr1 = from0;
854          path1[i]= 0;
855       } else {
856          curr1 = from1;
857          path1[i]= 1;
858       }
859       cost0 = curr0 + abs(metric[i]-2*tf_select_table[LM][4*isTransient+2*tf_select+0]);
860       cost1 = curr1 + abs(metric[i]-2*tf_select_table[LM][4*isTransient+2*tf_select+1]);
861    }
862    tf_res[len-1] = cost0 < cost1 ? 0 : 1;
863    /* Viterbi backward pass to check the decisions */
864    for (i=len-2;i>=0;i--)
865    {
866       if (tf_res[i+1] == 1)
867          tf_res[i] = path1[i+1];
868       else
869          tf_res[i] = path0[i+1];
870    }
871    /*printf("%d %f\n", *tf_sum, tf_estimate);*/
872    RESTORE_STACK;
873 #ifdef FUZZING
874    tf_select = rand()&0x1;
875    tf_res[0] = rand()&0x1;
876    for (i=1;i<len;i++)
877       tf_res[i] = tf_res[i-1] ^ ((rand()&0xF) == 0);
878 #endif
879    return tf_select;
880 }
881
882 static void tf_encode(int start, int end, int isTransient, int *tf_res, int LM, int tf_select, ec_enc *enc)
883 {
884    int curr, i;
885    int tf_select_rsv;
886    int tf_changed;
887    int logp;
888    opus_uint32 budget;
889    opus_uint32 tell;
890    budget = enc->storage*8;
891    tell = ec_tell(enc);
892    logp = isTransient ? 2 : 4;
893    /* Reserve space to code the tf_select decision. */
894    tf_select_rsv = LM>0 && tell+logp+1 <= budget;
895    budget -= tf_select_rsv;
896    curr = tf_changed = 0;
897    for (i=start;i<end;i++)
898    {
899       if (tell+logp<=budget)
900       {
901          ec_enc_bit_logp(enc, tf_res[i] ^ curr, logp);
902          tell = ec_tell(enc);
903          curr = tf_res[i];
904          tf_changed |= curr;
905       }
906       else
907          tf_res[i] = curr;
908       logp = isTransient ? 4 : 5;
909    }
910    /* Only code tf_select if it would actually make a difference. */
911    if (tf_select_rsv &&
912          tf_select_table[LM][4*isTransient+0+tf_changed]!=
913          tf_select_table[LM][4*isTransient+2+tf_changed])
914       ec_enc_bit_logp(enc, tf_select, 1);
915    else
916       tf_select = 0;
917    for (i=start;i<end;i++)
918       tf_res[i] = tf_select_table[LM][4*isTransient+2*tf_select+tf_res[i]];
919    /*for(i=0;i<end;i++)printf("%d ", isTransient ? tf_res[i] : LM+tf_res[i]);printf("\n");*/
920 }
921
922 static void tf_decode(int start, int end, int isTransient, int *tf_res, int LM, ec_dec *dec)
923 {
924    int i, curr, tf_select;
925    int tf_select_rsv;
926    int tf_changed;
927    int logp;
928    opus_uint32 budget;
929    opus_uint32 tell;
930
931    budget = dec->storage*8;
932    tell = ec_tell(dec);
933    logp = isTransient ? 2 : 4;
934    tf_select_rsv = LM>0 && tell+logp+1<=budget;
935    budget -= tf_select_rsv;
936    tf_changed = curr = 0;
937    for (i=start;i<end;i++)
938    {
939       if (tell+logp<=budget)
940       {
941          curr ^= ec_dec_bit_logp(dec, logp);
942          tell = ec_tell(dec);
943          tf_changed |= curr;
944       }
945       tf_res[i] = curr;
946       logp = isTransient ? 4 : 5;
947    }
948    tf_select = 0;
949    if (tf_select_rsv &&
950      tf_select_table[LM][4*isTransient+0+tf_changed] !=
951      tf_select_table[LM][4*isTransient+2+tf_changed])
952    {
953       tf_select = ec_dec_bit_logp(dec, 1);
954    }
955    for (i=start;i<end;i++)
956    {
957       tf_res[i] = tf_select_table[LM][4*isTransient+2*tf_select+tf_res[i]];
958    }
959 }
960
961 static void init_caps(const CELTMode *m,int *cap,int LM,int C)
962 {
963    int i;
964    for (i=0;i<m->nbEBands;i++)
965    {
966       int N;
967       N=(m->eBands[i+1]-m->eBands[i])<<LM;
968       cap[i] = (m->cache.caps[m->nbEBands*(2*LM+C-1)+i]+64)*C*N>>2;
969    }
970 }
971
972 static int alloc_trim_analysis(const CELTMode *m, const celt_norm *X,
973       const opus_val16 *bandLogE, int end, int LM, int C, int N0,
974       AnalysisInfo *analysis, opus_val16 *stereo_saving, opus_val16 tf_estimate,
975       int intensity)
976 {
977    int i;
978    opus_val32 diff=0;
979    int c;
980    int trim_index = 5;
981    opus_val16 trim = QCONST16(5.f, 8);
982    opus_val16 logXC, logXC2;
983    if (C==2)
984    {
985       opus_val16 sum = 0; /* Q10 */
986       opus_val16 minXC; /* Q10 */
987       /* Compute inter-channel correlation for low frequencies */
988       for (i=0;i<8;i++)
989       {
990          int j;
991          opus_val32 partial = 0;
992          for (j=m->eBands[i]<<LM;j<m->eBands[i+1]<<LM;j++)
993             partial = MAC16_16(partial, X[j], X[N0+j]);
994          sum = ADD16(sum, EXTRACT16(SHR32(partial, 18)));
995       }
996       sum = MULT16_16_Q15(QCONST16(1.f/8, 15), sum);
997       sum = MIN16(QCONST16(1.f, 10), ABS16(sum));
998       minXC = sum;
999       for (i=8;i<intensity;i++)
1000       {
1001          int j;
1002          opus_val32 partial = 0;
1003          for (j=m->eBands[i]<<LM;j<m->eBands[i+1]<<LM;j++)
1004             partial = MAC16_16(partial, X[j], X[N0+j]);
1005          minXC = MIN16(minXC, ABS16(EXTRACT16(SHR32(partial, 18))));
1006       }
1007       minXC = MIN16(QCONST16(1.f, 10), ABS16(minXC));
1008       /*printf ("%f\n", sum);*/
1009       if (sum > QCONST16(.995f,10))
1010          trim_index-=4;
1011       else if (sum > QCONST16(.92f,10))
1012          trim_index-=3;
1013       else if (sum > QCONST16(.85f,10))
1014          trim_index-=2;
1015       else if (sum > QCONST16(.8f,10))
1016          trim_index-=1;
1017       /* mid-side savings estimations based on the LF average*/
1018       logXC = celt_log2(QCONST32(1.001f, 20)-MULT16_16(sum, sum));
1019       /* mid-side savings estimations based on min correlation */
1020       logXC2 = MAX16(HALF16(logXC), celt_log2(QCONST32(1.001f, 20)-MULT16_16(minXC, minXC)));
1021 #ifdef FIXED_POINT
1022       /* Compensate for Q20 vs Q14 input and convert output to Q8 */
1023       logXC = PSHR32(logXC-QCONST16(6.f, DB_SHIFT),DB_SHIFT-8);
1024       logXC2 = PSHR32(logXC2-QCONST16(6.f, DB_SHIFT),DB_SHIFT-8);
1025 #endif
1026
1027       trim += MAX16(-QCONST16(4.f, 8), MULT16_16_Q15(QCONST16(.75f,15),logXC));
1028       *stereo_saving = MIN16(*stereo_saving + QCONST16(0.25f, 8), -HALF16(logXC2));
1029    }
1030
1031    /* Estimate spectral tilt */
1032    c=0; do {
1033       for (i=0;i<end-1;i++)
1034       {
1035          diff += bandLogE[i+c*m->nbEBands]*(opus_int32)(2+2*i-end);
1036       }
1037    } while (++c<C);
1038    diff /= C*(end-1);
1039    /*printf("%f\n", diff);*/
1040    if (diff > QCONST16(2.f, DB_SHIFT))
1041       trim_index--;
1042    if (diff > QCONST16(8.f, DB_SHIFT))
1043       trim_index--;
1044    if (diff < -QCONST16(4.f, DB_SHIFT))
1045       trim_index++;
1046    if (diff < -QCONST16(10.f, DB_SHIFT))
1047       trim_index++;
1048    trim -= MAX16(-QCONST16(2.f, 8), MIN16(QCONST16(2.f, 8), SHR16(diff+QCONST16(1.f, DB_SHIFT),DB_SHIFT-8)/6 ));
1049    trim -= 2*SHR16(tf_estimate-QCONST16(1.f,14), 14-8);
1050 #ifndef FIXED_POINT
1051    if (analysis->valid)
1052    {
1053       trim -= MAX16(-QCONST16(2.f, 8), MIN16(QCONST16(2.f, 8), 2*(analysis->tonality_slope+.05)));
1054    }
1055 #endif
1056
1057 #ifdef FIXED_POINT
1058    trim_index = PSHR32(trim, 8);
1059 #else
1060    trim_index = floor(.5+trim);
1061 #endif
1062    if (trim_index<0)
1063       trim_index = 0;
1064    if (trim_index>10)
1065       trim_index = 10;
1066    /*printf("%d\n", trim_index);*/
1067 #ifdef FUZZING
1068    trim_index = rand()%11;
1069 #endif
1070    return trim_index;
1071 }
1072
1073 static int stereo_analysis(const CELTMode *m, const celt_norm *X,
1074       int LM, int N0)
1075 {
1076    int i;
1077    int thetas;
1078    opus_val32 sumLR = EPSILON, sumMS = EPSILON;
1079
1080    /* Use the L1 norm to model the entropy of the L/R signal vs the M/S signal */
1081    for (i=0;i<13;i++)
1082    {
1083       int j;
1084       for (j=m->eBands[i]<<LM;j<m->eBands[i+1]<<LM;j++)
1085       {
1086          opus_val32 L, R, M, S;
1087          /* We cast to 32-bit first because of the -32768 case */
1088          L = EXTEND32(X[j]);
1089          R = EXTEND32(X[N0+j]);
1090          M = ADD32(L, R);
1091          S = SUB32(L, R);
1092          sumLR = ADD32(sumLR, ADD32(ABS32(L), ABS32(R)));
1093          sumMS = ADD32(sumMS, ADD32(ABS32(M), ABS32(S)));
1094       }
1095    }
1096    sumMS = MULT16_32_Q15(QCONST16(0.707107f, 15), sumMS);
1097    thetas = 13;
1098    /* We don't need thetas for lower bands with LM<=1 */
1099    if (LM<=1)
1100       thetas -= 8;
1101    return MULT16_32_Q15((m->eBands[13]<<(LM+1))+thetas, sumMS)
1102          > MULT16_32_Q15(m->eBands[13]<<(LM+1), sumLR);
1103 }
1104
1105 static int run_prefilter(CELTEncoder *st, celt_sig *in, celt_sig *prefilter_mem, int CC, int N,
1106       int prefilter_tapset, int *pitch, opus_val16 *gain, int *qgain, int enabled, int nbAvailableBytes)
1107 {
1108    int c;
1109    VARDECL(celt_sig, _pre);
1110    celt_sig *pre[2];
1111    const CELTMode *mode;
1112    int pitch_index;
1113    opus_val16 gain1;
1114    opus_val16 pf_threshold;
1115    int pf_on;
1116    int qg;
1117    SAVE_STACK;
1118
1119    mode = st->mode;
1120    ALLOC(_pre, CC*(N+COMBFILTER_MAXPERIOD), celt_sig);
1121
1122    pre[0] = _pre;
1123    pre[1] = _pre + (N+COMBFILTER_MAXPERIOD);
1124
1125
1126    c=0; do {
1127       OPUS_COPY(pre[c], prefilter_mem+c*COMBFILTER_MAXPERIOD, COMBFILTER_MAXPERIOD);
1128       OPUS_COPY(pre[c]+COMBFILTER_MAXPERIOD, in+c*(N+st->overlap)+st->overlap, N);
1129    } while (++c<CC);
1130
1131    if (enabled)
1132    {
1133       VARDECL(opus_val16, pitch_buf);
1134       ALLOC(pitch_buf, (COMBFILTER_MAXPERIOD+N)>>1, opus_val16);
1135
1136       pitch_downsample(pre, pitch_buf, COMBFILTER_MAXPERIOD+N, CC);
1137       /* Don't search for the fir last 1.5 octave of the range because
1138          there's too many false-positives due to short-term correlation */
1139       pitch_search(pitch_buf+(COMBFILTER_MAXPERIOD>>1), pitch_buf, N,
1140             COMBFILTER_MAXPERIOD-3*COMBFILTER_MINPERIOD, &pitch_index);
1141       pitch_index = COMBFILTER_MAXPERIOD-pitch_index;
1142
1143       gain1 = remove_doubling(pitch_buf, COMBFILTER_MAXPERIOD, COMBFILTER_MINPERIOD,
1144             N, &pitch_index, st->prefilter_period, st->prefilter_gain);
1145       if (pitch_index > COMBFILTER_MAXPERIOD-2)
1146          pitch_index = COMBFILTER_MAXPERIOD-2;
1147       gain1 = MULT16_16_Q15(QCONST16(.7f,15),gain1);
1148       /*printf("%d %d %f %f\n", pitch_change, pitch_index, gain1, st->analysis.tonality);*/
1149       if (st->loss_rate>2)
1150          gain1 = HALF32(gain1);
1151       if (st->loss_rate>4)
1152          gain1 = HALF32(gain1);
1153       if (st->loss_rate>8)
1154          gain1 = 0;
1155    } else {
1156       gain1 = 0;
1157       pitch_index = COMBFILTER_MINPERIOD;
1158    }
1159
1160    /* Gain threshold for enabling the prefilter/postfilter */
1161    pf_threshold = QCONST16(.2f,15);
1162
1163    /* Adjusting the threshold based on rate and continuity */
1164    if (abs(pitch_index-st->prefilter_period)*10>pitch_index)
1165       pf_threshold += QCONST16(.2f,15);
1166    if (nbAvailableBytes<25)
1167       pf_threshold += QCONST16(.1f,15);
1168    if (nbAvailableBytes<35)
1169       pf_threshold += QCONST16(.1f,15);
1170    if (st->prefilter_gain > QCONST16(.4f,15))
1171       pf_threshold -= QCONST16(.1f,15);
1172    if (st->prefilter_gain > QCONST16(.55f,15))
1173       pf_threshold -= QCONST16(.1f,15);
1174
1175    /* Hard threshold at 0.2 */
1176    pf_threshold = MAX16(pf_threshold, QCONST16(.2f,15));
1177    if (gain1<pf_threshold)
1178    {
1179       gain1 = 0;
1180       pf_on = 0;
1181       qg = 0;
1182    } else {
1183       /*This block is not gated by a total bits check only because
1184         of the nbAvailableBytes check above.*/
1185       if (ABS16(gain1-st->prefilter_gain)<QCONST16(.1f,15))
1186          gain1=st->prefilter_gain;
1187
1188 #ifdef FIXED_POINT
1189       qg = ((gain1+1536)>>10)/3-1;
1190 #else
1191       qg = (int)floor(.5f+gain1*32/3)-1;
1192 #endif
1193       qg = IMAX(0, IMIN(7, qg));
1194       gain1 = QCONST16(0.09375f,15)*(qg+1);
1195       pf_on = 1;
1196    }
1197    /*printf("%d %f\n", pitch_index, gain1);*/
1198
1199    c=0; do {
1200       int offset = mode->shortMdctSize-st->overlap;
1201       st->prefilter_period=IMAX(st->prefilter_period, COMBFILTER_MINPERIOD);
1202       OPUS_COPY(in+c*(N+st->overlap), st->in_mem+c*(st->overlap), st->overlap);
1203       if (offset)
1204          comb_filter(in+c*(N+st->overlap)+st->overlap, pre[c]+COMBFILTER_MAXPERIOD,
1205                st->prefilter_period, st->prefilter_period, offset, -st->prefilter_gain, -st->prefilter_gain,
1206                st->prefilter_tapset, st->prefilter_tapset, NULL, 0);
1207
1208       comb_filter(in+c*(N+st->overlap)+st->overlap+offset, pre[c]+COMBFILTER_MAXPERIOD+offset,
1209             st->prefilter_period, pitch_index, N-offset, -st->prefilter_gain, -gain1,
1210             st->prefilter_tapset, prefilter_tapset, mode->window, st->overlap);
1211       OPUS_COPY(st->in_mem+c*(st->overlap), in+c*(N+st->overlap)+N, st->overlap);
1212
1213       if (N>COMBFILTER_MAXPERIOD)
1214       {
1215          OPUS_MOVE(prefilter_mem+c*COMBFILTER_MAXPERIOD, pre[c]+N, COMBFILTER_MAXPERIOD);
1216       } else {
1217          OPUS_MOVE(prefilter_mem+c*COMBFILTER_MAXPERIOD, prefilter_mem+c*COMBFILTER_MAXPERIOD+N, COMBFILTER_MAXPERIOD-N);
1218          OPUS_MOVE(prefilter_mem+c*COMBFILTER_MAXPERIOD+COMBFILTER_MAXPERIOD-N, pre[c]+COMBFILTER_MAXPERIOD, N);
1219       }
1220    } while (++c<CC);
1221
1222    RESTORE_STACK;
1223    *gain = gain1;
1224    *pitch = pitch_index;
1225    *qgain = qg;
1226    return pf_on;
1227 }
1228
1229 int celt_encode_with_ec(CELTEncoder * OPUS_RESTRICT st, const opus_val16 * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes, ec_enc *enc)
1230 {
1231    int i, c, N;
1232    opus_int32 bits;
1233    ec_enc _enc;
1234    VARDECL(celt_sig, in);
1235    VARDECL(celt_sig, freq);
1236    VARDECL(celt_norm, X);
1237    VARDECL(celt_ener, bandE);
1238    VARDECL(opus_val16, bandLogE);
1239    VARDECL(opus_val16, bandLogE2);
1240    VARDECL(int, fine_quant);
1241    VARDECL(opus_val16, error);
1242    VARDECL(int, pulses);
1243    VARDECL(int, cap);
1244    VARDECL(int, offsets);
1245    VARDECL(int, fine_priority);
1246    VARDECL(int, tf_res);
1247    VARDECL(unsigned char, collapse_masks);
1248    celt_sig *prefilter_mem;
1249    opus_val16 *oldBandE, *oldLogE, *oldLogE2;
1250    int shortBlocks=0;
1251    int isTransient=0;
1252    const int CC = st->channels;
1253    const int C = st->stream_channels;
1254    int LM, M;
1255    int tf_select;
1256    int nbFilledBytes, nbAvailableBytes;
1257    int effEnd;
1258    int codedBands;
1259    int tf_sum;
1260    int alloc_trim;
1261    int pitch_index=COMBFILTER_MINPERIOD;
1262    opus_val16 gain1 = 0;
1263    int dual_stereo=0;
1264    int effectiveBytes;
1265    int dynalloc_logp;
1266    opus_int32 vbr_rate;
1267    opus_int32 total_bits;
1268    opus_int32 total_boost;
1269    opus_int32 balance;
1270    opus_int32 tell;
1271    int prefilter_tapset=0;
1272    int pf_on;
1273    int anti_collapse_rsv;
1274    int anti_collapse_on=0;
1275    int silence=0;
1276    int tf_chan = 0;
1277    opus_val16 tf_estimate;
1278    int pitch_change=0;
1279    opus_int32 tot_boost=0;
1280    opus_val16 sample_max;
1281    opus_val16 maxDepth;
1282    const OpusCustomMode *mode;
1283    int nbEBands;
1284    int overlap;
1285    const opus_int16 *eBands;
1286    ALLOC_STACK;
1287
1288    mode = st->mode;
1289    nbEBands = mode->nbEBands;
1290    overlap = mode->overlap;
1291    eBands = mode->eBands;
1292    tf_estimate = QCONST16(1.0f,14);
1293    if (nbCompressedBytes<2 || pcm==NULL)
1294      return OPUS_BAD_ARG;
1295
1296    frame_size *= st->upsample;
1297    for (LM=0;LM<=mode->maxLM;LM++)
1298       if (mode->shortMdctSize<<LM==frame_size)
1299          break;
1300    if (LM>mode->maxLM)
1301       return OPUS_BAD_ARG;
1302    M=1<<LM;
1303    N = M*mode->shortMdctSize;
1304
1305    prefilter_mem = st->in_mem+CC*(st->overlap);
1306    oldBandE = (opus_val16*)(st->in_mem+CC*(st->overlap+COMBFILTER_MAXPERIOD));
1307    oldLogE = oldBandE + CC*nbEBands;
1308    oldLogE2 = oldLogE + CC*nbEBands;
1309
1310    if (enc==NULL)
1311    {
1312       tell=1;
1313       nbFilledBytes=0;
1314    } else {
1315       tell=ec_tell(enc);
1316       nbFilledBytes=(tell+4)>>3;
1317    }
1318
1319 #ifdef CUSTOM_MODES
1320    if (st->signalling && enc==NULL)
1321    {
1322       int tmp = (mode->effEBands-st->end)>>1;
1323       st->end = IMAX(1, mode->effEBands-tmp);
1324       compressed[0] = tmp<<5;
1325       compressed[0] |= LM<<3;
1326       compressed[0] |= (C==2)<<2;
1327       /* Convert "standard mode" to Opus header */
1328       if (mode->Fs==48000 && mode->shortMdctSize==120)
1329       {
1330          int c0 = toOpus(compressed[0]);
1331          if (c0<0)
1332             return OPUS_BAD_ARG;
1333          compressed[0] = c0;
1334       }
1335       compressed++;
1336       nbCompressedBytes--;
1337    }
1338 #else
1339    celt_assert(st->signalling==0);
1340 #endif
1341
1342    /* Can't produce more than 1275 output bytes */
1343    nbCompressedBytes = IMIN(nbCompressedBytes,1275);
1344    nbAvailableBytes = nbCompressedBytes - nbFilledBytes;
1345
1346    if (st->vbr && st->bitrate!=OPUS_BITRATE_MAX)
1347    {
1348       opus_int32 den=mode->Fs>>BITRES;
1349       vbr_rate=(st->bitrate*frame_size+(den>>1))/den;
1350 #ifdef CUSTOM_MODES
1351       if (st->signalling)
1352          vbr_rate -= 8<<BITRES;
1353 #endif
1354       effectiveBytes = vbr_rate>>(3+BITRES);
1355    } else {
1356       opus_int32 tmp;
1357       vbr_rate = 0;
1358       tmp = st->bitrate*frame_size;
1359       if (tell>1)
1360          tmp += tell;
1361       if (st->bitrate!=OPUS_BITRATE_MAX)
1362          nbCompressedBytes = IMAX(2, IMIN(nbCompressedBytes,
1363                (tmp+4*mode->Fs)/(8*mode->Fs)-!!st->signalling));
1364       effectiveBytes = nbCompressedBytes;
1365    }
1366
1367    if (enc==NULL)
1368    {
1369       ec_enc_init(&_enc, compressed, nbCompressedBytes);
1370       enc = &_enc;
1371    }
1372
1373    if (vbr_rate>0)
1374    {
1375       /* Computes the max bit-rate allowed in VBR mode to avoid violating the
1376           target rate and buffering.
1377          We must do this up front so that bust-prevention logic triggers
1378           correctly if we don't have enough bits. */
1379       if (st->constrained_vbr)
1380       {
1381          opus_int32 vbr_bound;
1382          opus_int32 max_allowed;
1383          /* We could use any multiple of vbr_rate as bound (depending on the
1384              delay).
1385             This is clamped to ensure we use at least two bytes if the encoder
1386              was entirely empty, but to allow 0 in hybrid mode. */
1387          vbr_bound = vbr_rate;
1388          max_allowed = IMIN(IMAX(tell==1?2:0,
1389                (vbr_rate+vbr_bound-st->vbr_reservoir)>>(BITRES+3)),
1390                nbAvailableBytes);
1391          if(max_allowed < nbAvailableBytes)
1392          {
1393             nbCompressedBytes = nbFilledBytes+max_allowed;
1394             nbAvailableBytes = max_allowed;
1395             ec_enc_shrink(enc, nbCompressedBytes);
1396          }
1397       }
1398    }
1399    total_bits = nbCompressedBytes*8;
1400
1401    effEnd = st->end;
1402    if (effEnd > mode->effEBands)
1403       effEnd = mode->effEBands;
1404
1405    ALLOC(in, CC*(N+st->overlap), celt_sig);
1406
1407    sample_max=MAX16(st->overlap_max, celt_maxabs16(pcm, C*(N-overlap)/st->upsample));
1408    st->overlap_max=celt_maxabs16(pcm+C*(N-overlap)/st->upsample, C*overlap/st->upsample);
1409    sample_max=MAX16(sample_max, st->overlap_max);
1410 #ifdef FIXED_POINT
1411    silence = (sample_max==0);
1412 #else
1413    silence = (sample_max <= (opus_val16)1/(1<<st->lsb_depth));
1414 #endif
1415 #ifdef FUZZING
1416    if ((rand()&0x3F)==0)
1417       silence = 1;
1418 #endif
1419    if (tell==1)
1420       ec_enc_bit_logp(enc, silence, 15);
1421    else
1422       silence=0;
1423    if (silence)
1424    {
1425       /*In VBR mode there is no need to send more than the minimum. */
1426       if (vbr_rate>0)
1427       {
1428          effectiveBytes=nbCompressedBytes=IMIN(nbCompressedBytes, nbFilledBytes+2);
1429          total_bits=nbCompressedBytes*8;
1430          nbAvailableBytes=2;
1431          ec_enc_shrink(enc, nbCompressedBytes);
1432       }
1433       /* Pretend we've filled all the remaining bits with zeros
1434             (that's what the initialiser did anyway) */
1435       tell = nbCompressedBytes*8;
1436       enc->nbits_total+=tell-ec_tell(enc);
1437    }
1438    c=0; do {
1439       preemphasis(pcm+c, in+c*(N+st->overlap)+st->overlap, N, CC, st->upsample,
1440                   mode->preemph, st->preemph_memE+c, st->clip);
1441    } while (++c<CC);
1442
1443
1444
1445    /* Find pitch period and gain */
1446    {
1447       int enabled;
1448       int qg;
1449       enabled = nbAvailableBytes>12*C && st->start==0 && !silence && !st->disable_pf && st->complexity >= 5;
1450
1451       prefilter_tapset = st->tapset_decision;
1452       pf_on = run_prefilter(st, in, prefilter_mem, CC, N, prefilter_tapset, &pitch_index, &gain1, &qg, enabled, nbAvailableBytes);
1453       if ((gain1 > QCONST16(.4f,15) || st->prefilter_gain > QCONST16(.4f,15)) && st->analysis.tonality > .3
1454             && (pitch_index > 1.26*st->prefilter_period || pitch_index < .79*st->prefilter_period))
1455          pitch_change = 1;
1456       if (pf_on==0)
1457       {
1458          if(st->start==0 && tell+16<=total_bits)
1459             ec_enc_bit_logp(enc, 0, 1);
1460       } else {
1461          /*This block is not gated by a total bits check only because
1462            of the nbAvailableBytes check above.*/
1463          int octave;
1464          ec_enc_bit_logp(enc, 1, 1);
1465          pitch_index += 1;
1466          octave = EC_ILOG(pitch_index)-5;
1467          ec_enc_uint(enc, octave, 6);
1468          ec_enc_bits(enc, pitch_index-(16<<octave), 4+octave);
1469          pitch_index -= 1;
1470          ec_enc_bits(enc, qg, 3);
1471          ec_enc_icdf(enc, prefilter_tapset, tapset_icdf, 2);
1472       }
1473    }
1474
1475    isTransient = 0;
1476    shortBlocks = 0;
1477    if (LM>0 && ec_tell(enc)+3<=total_bits)
1478    {
1479       if (st->complexity > 1)
1480       {
1481          isTransient = transient_analysis(in, N+st->overlap, CC,
1482                   &tf_estimate, &tf_chan);
1483          if (isTransient)
1484             shortBlocks = M;
1485       }
1486       ec_enc_bit_logp(enc, isTransient, 3);
1487    }
1488
1489    ALLOC(freq, CC*N, celt_sig); /**< Interleaved signal MDCTs */
1490    ALLOC(bandE,nbEBands*CC, celt_ener);
1491    ALLOC(bandLogE,nbEBands*CC, opus_val16);
1492    /* Compute MDCTs */
1493    compute_mdcts(mode, shortBlocks, in, freq, CC, LM);
1494
1495    if (CC==2&&C==1)
1496    {
1497       for (i=0;i<N;i++)
1498          freq[i] = ADD32(HALF32(freq[i]), HALF32(freq[N+i]));
1499       tf_chan = 0;
1500    }
1501    if (st->upsample != 1)
1502    {
1503       c=0; do
1504       {
1505          int bound = N/st->upsample;
1506          for (i=0;i<bound;i++)
1507             freq[c*N+i] *= st->upsample;
1508          for (;i<N;i++)
1509             freq[c*N+i] = 0;
1510       } while (++c<C);
1511    }
1512    compute_band_energies(mode, freq, bandE, effEnd, C, M);
1513
1514    amp2Log2(mode, effEnd, st->end, bandE, bandLogE, C);
1515    /*for (i=0;i<21;i++)
1516       printf("%f ", bandLogE[i]);
1517    printf("\n");*/
1518
1519    ALLOC(bandLogE2, C*nbEBands, opus_val16);
1520    if (shortBlocks && st->complexity>=8)
1521    {
1522       VARDECL(celt_sig, freq2);
1523       VARDECL(opus_val32, bandE2);
1524       ALLOC(freq2, CC*N, celt_sig);
1525       compute_mdcts(mode, 0, in, freq2, CC, LM);
1526       if (CC==2&&C==1)
1527       {
1528          for (i=0;i<N;i++)
1529             freq2[i] = ADD32(HALF32(freq2[i]), HALF32(freq2[N+i]));
1530       }
1531       if (st->upsample != 1)
1532       {
1533          c=0; do
1534          {
1535             int bound = N/st->upsample;
1536             for (i=0;i<bound;i++)
1537                freq2[c*N+i] *= st->upsample;
1538             for (;i<N;i++)
1539                freq2[c*N+i] = 0;
1540          } while (++c<C);
1541       }
1542       ALLOC(bandE2, C*nbEBands, opus_val32);
1543       compute_band_energies(mode, freq2, bandE2, effEnd, C, M);
1544       amp2Log2(mode, effEnd, st->end, bandE2, bandLogE2, C);
1545       for (i=0;i<C*nbEBands;i++)
1546          bandLogE2[i] += HALF16(SHL16(LM, DB_SHIFT));
1547    } else {
1548       for (i=0;i<C*nbEBands;i++)
1549          bandLogE2[i] = bandLogE[i];
1550    }
1551
1552    ALLOC(X, C*N, celt_norm);         /**< Interleaved normalised MDCTs */
1553
1554    /* Band normalisation */
1555    normalise_bands(mode, freq, X, bandE, effEnd, C, M);
1556
1557    ALLOC(tf_res, nbEBands, int);
1558    tf_select = tf_analysis(mode, effEnd, C, isTransient, tf_res, effectiveBytes, X, N, LM, &tf_sum, tf_estimate, tf_chan);
1559    for (i=effEnd;i<st->end;i++)
1560       tf_res[i] = tf_res[effEnd-1];
1561
1562    ALLOC(error, C*nbEBands, opus_val16);
1563    quant_coarse_energy(mode, st->start, st->end, effEnd, bandLogE,
1564          oldBandE, total_bits, error, enc,
1565          C, LM, nbAvailableBytes, st->force_intra,
1566          &st->delayedIntra, st->complexity >= 4, st->loss_rate);
1567
1568    tf_encode(st->start, st->end, isTransient, tf_res, LM, tf_select, enc);
1569
1570    if (ec_tell(enc)+4<=total_bits)
1571    {
1572       if (shortBlocks || st->complexity < 3 || nbAvailableBytes < 10*C)
1573       {
1574          if (st->complexity == 0)
1575             st->spread_decision = SPREAD_NONE;
1576       } else {
1577          if (st->analysis.valid)
1578          {
1579             static const opus_val16 spread_thresholds[3] = {-QCONST16(.6f, 15), -QCONST16(.2f, 15), -QCONST16(.07f, 15)};
1580             static const opus_val16 spread_histeresis[3] = {QCONST16(.15f, 15), QCONST16(.07f, 15), QCONST16(.02f, 15)};
1581             static const opus_val16 tapset_thresholds[2] = {QCONST16(.0f, 15), QCONST16(.15f, 15)};
1582             static const opus_val16 tapset_histeresis[2] = {QCONST16(.1f, 15), QCONST16(.05f, 15)};
1583             st->spread_decision = hysteresis_decision(-st->analysis.tonality, spread_thresholds, spread_histeresis, 3, st->spread_decision);
1584             st->tapset_decision = hysteresis_decision(st->analysis.tonality_slope, tapset_thresholds, tapset_histeresis, 2, st->tapset_decision);
1585          } else {
1586             st->spread_decision = spreading_decision(mode, X,
1587                   &st->tonal_average, st->spread_decision, &st->hf_average,
1588                   &st->tapset_decision, pf_on&&!shortBlocks, effEnd, C, M);
1589          }
1590          /*printf("%d %d\n", st->tapset_decision, st->spread_decision);*/
1591          /*printf("%f %d %f %d\n\n", st->analysis.tonality, st->spread_decision, st->analysis.tonality_slope, st->tapset_decision);*/
1592       }
1593       ec_enc_icdf(enc, st->spread_decision, spread_icdf, 5);
1594    }
1595
1596    ALLOC(cap, nbEBands, int);
1597    ALLOC(offsets, nbEBands, int);
1598
1599    init_caps(mode,cap,LM,C);
1600    for (i=0;i<nbEBands;i++)
1601       offsets[i] = 0;
1602    /* Dynamic allocation code */
1603    maxDepth=-QCONST16(32.f, DB_SHIFT);
1604    /* Make sure that dynamic allocation can't make us bust the budget */
1605    if (effectiveBytes > 50 && LM>=1)
1606    {
1607       int last=0;
1608       VARDECL(opus_val16, follower);
1609       ALLOC(follower, C*nbEBands, opus_val16);
1610       c=0;do
1611       {
1612          follower[c*nbEBands] = bandLogE2[c*nbEBands];
1613          for (i=1;i<st->end;i++)
1614          {
1615             /* The last band to be at least 3 dB higher than the previous one
1616                is the last we'll consider. Otherwise, we run into problems on
1617                bandlimited signals. */
1618             if (bandLogE2[c*nbEBands+i] > bandLogE2[c*nbEBands+i-1]+QCONST16(.5f,DB_SHIFT))
1619                last=i;
1620             follower[c*nbEBands+i] = MIN16(follower[c*nbEBands+i-1]+QCONST16(1.5f,DB_SHIFT), bandLogE2[c*nbEBands+i]);
1621          }
1622          for (i=last-1;i>=0;i--)
1623             follower[c*nbEBands+i] = MIN16(follower[c*nbEBands+i], MIN16(follower[c*nbEBands+i+1]+QCONST16(2.f,DB_SHIFT), bandLogE2[c*nbEBands+i]));
1624          for (i=0;i<st->end;i++)
1625          {
1626             opus_val16 noise_floor;
1627             /* Noise floor must take into account eMeans, the depth, the width of the bands
1628                and the preemphasis filter (approx. square of bark band ID) */
1629             noise_floor = MULT16_16(QCONST16(0.0625f, DB_SHIFT),mode->logN[i])
1630                   +QCONST16(.5f,DB_SHIFT)+SHL16(9-st->lsb_depth,DB_SHIFT)-SHL16(eMeans[i],6)
1631                   +MULT16_16(QCONST16(.0062,DB_SHIFT),(i+5)*(i+5));
1632             follower[c*nbEBands+i] = MAX16(follower[c*nbEBands+i], noise_floor);
1633             maxDepth = MAX16(maxDepth, bandLogE[c*nbEBands+i]-noise_floor);
1634          }
1635       } while (++c<C);
1636       if (C==2)
1637       {
1638          for (i=st->start;i<st->end;i++)
1639          {
1640             /* Consider 24 dB "cross-talk" */
1641             follower[nbEBands+i] = MAX16(follower[nbEBands+i], follower[                   i]-QCONST16(4.f,DB_SHIFT));
1642             follower[                   i] = MAX16(follower[                   i], follower[nbEBands+i]-QCONST16(4.f,DB_SHIFT));
1643             follower[i] = HALF16(MAX16(0, bandLogE[i]-follower[i]) + MAX16(0, bandLogE[nbEBands+i]-follower[nbEBands+i]));
1644          }
1645       } else {
1646          for (i=st->start;i<st->end;i++)
1647          {
1648             follower[i] = MAX16(0, bandLogE[i]-follower[i]);
1649          }
1650       }
1651       /* For non-transient CBR/CVBR frames, halve the dynalloc contribution */
1652       if ((!st->vbr || st->constrained_vbr)&&!isTransient)
1653       {
1654          for (i=st->start;i<st->end;i++)
1655             follower[i] = HALF16(follower[i]);
1656       }
1657       for (i=st->start;i<st->end;i++)
1658       {
1659          int width;
1660          int boost;
1661          int boost_bits;
1662
1663          if (i<8)
1664             follower[i] *= 2;
1665          if (i>=12)
1666             follower[i] = HALF16(follower[i]);
1667          follower[i] = MIN16(follower[i], QCONST16(4, DB_SHIFT));
1668
1669          width = C*(eBands[i+1]-eBands[i])<<LM;
1670          if (width<6)
1671          {
1672             boost = SHR32(EXTEND32(follower[i]),DB_SHIFT);
1673             boost_bits = boost*width<<BITRES;
1674          } else if (width > 48) {
1675             boost = SHR32(EXTEND32(follower[i])*8,DB_SHIFT);
1676             boost_bits = (boost*width<<BITRES)/8;
1677          } else {
1678             boost = SHR32(EXTEND32(follower[i])*width/6,DB_SHIFT);
1679             boost_bits = boost*6<<BITRES;
1680          }
1681          /* For CBR and non-transient CVBR frames, limit dynalloc to 1/4 of the bits */
1682          if ((!st->vbr || (st->constrained_vbr&&!isTransient))
1683                && (tot_boost+boost_bits)>>BITRES>>3 > effectiveBytes/4)
1684          {
1685             offsets[i] = 0;
1686             break;
1687          } else {
1688             offsets[i] = boost;
1689             tot_boost += boost_bits;
1690          }
1691       }
1692    }
1693    dynalloc_logp = 6;
1694    total_bits<<=BITRES;
1695    total_boost = 0;
1696    tell = ec_tell_frac(enc);
1697    for (i=st->start;i<st->end;i++)
1698    {
1699       int width, quanta;
1700       int dynalloc_loop_logp;
1701       int boost;
1702       int j;
1703       width = C*(eBands[i+1]-eBands[i])<<LM;
1704       /* quanta is 6 bits, but no more than 1 bit/sample
1705          and no less than 1/8 bit/sample */
1706       quanta = IMIN(width<<BITRES, IMAX(6<<BITRES, width));
1707       dynalloc_loop_logp = dynalloc_logp;
1708       boost = 0;
1709       for (j = 0; tell+(dynalloc_loop_logp<<BITRES) < total_bits-total_boost
1710             && boost < cap[i]; j++)
1711       {
1712          int flag;
1713          flag = j<offsets[i];
1714          ec_enc_bit_logp(enc, flag, dynalloc_loop_logp);
1715          tell = ec_tell_frac(enc);
1716          if (!flag)
1717             break;
1718          boost += quanta;
1719          total_boost += quanta;
1720          dynalloc_loop_logp = 1;
1721       }
1722       /* Making dynalloc more likely */
1723       if (j)
1724          dynalloc_logp = IMAX(2, dynalloc_logp-1);
1725       offsets[i] = boost;
1726    }
1727
1728    if (C==2)
1729    {
1730       int effectiveRate;
1731
1732       static const opus_val16 intensity_thresholds[21]=
1733       /* 0  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19  20  off*/
1734         { 16,21,23,25,27,29,31,33,35,38,42,46,50,54,58,63,68,75,84,102,130};
1735       static const opus_val16 intensity_histeresis[21]=
1736         {  2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3, 4, 5, 6,  8, 12};
1737
1738       /* Always use MS for 2.5 ms frames until we can do a better analysis */
1739       if (LM!=0)
1740          dual_stereo = stereo_analysis(mode, X, LM, N);
1741
1742       /* Account for coarse energy */
1743       effectiveRate = (8*effectiveBytes - 80)>>LM;
1744
1745       /* effectiveRate in kb/s */
1746       effectiveRate = 2*effectiveRate/5;
1747
1748       st->intensity = hysteresis_decision(effectiveRate, intensity_thresholds, intensity_histeresis, 21, st->intensity);
1749       st->intensity = IMIN(st->end,IMAX(st->start, st->intensity));
1750    }
1751
1752    alloc_trim = 5;
1753    if (tell+(6<<BITRES) <= total_bits - total_boost)
1754    {
1755       alloc_trim = alloc_trim_analysis(mode, X, bandLogE,
1756             st->end, LM, C, N, &st->analysis, &st->stereo_saving, tf_estimate, st->intensity);
1757       ec_enc_icdf(enc, alloc_trim, trim_icdf, 7);
1758       tell = ec_tell_frac(enc);
1759    }
1760
1761    /* Variable bitrate */
1762    if (vbr_rate>0)
1763    {
1764      opus_val16 alpha;
1765      opus_int32 delta;
1766      /* The target rate in 8th bits per frame */
1767      opus_int32 target, base_target;
1768      opus_int32 min_allowed;
1769      int coded_bins;
1770      int coded_bands;
1771      int lm_diff = mode->maxLM - LM;
1772      coded_bands = st->lastCodedBands ? st->lastCodedBands : nbEBands;
1773      coded_bins = eBands[coded_bands]<<LM;
1774      if (C==2)
1775         coded_bins += eBands[IMIN(st->intensity, coded_bands)]<<LM;
1776
1777      /* Don't attempt to use more than 510 kb/s, even for frames smaller than 20 ms.
1778         The CELT allocator will just not be able to use more than that anyway. */
1779      nbCompressedBytes = IMIN(nbCompressedBytes,1275>>(3-LM));
1780      target = vbr_rate - ((40*C+20)<<BITRES);
1781      base_target = target;
1782
1783      if (st->constrained_vbr)
1784         target += (st->vbr_offset>>lm_diff);
1785
1786      /*printf("%f %f %f %f %d %d ", st->analysis.activity, st->analysis.tonality, tf_estimate, st->stereo_saving, tot_boost, coded_bands);*/
1787 #ifndef FIXED_POINT
1788      if (st->analysis.valid && st->analysis.activity<.4)
1789         target -= (coded_bins<<BITRES)*1*(.4-st->analysis.activity);
1790 #endif
1791      /* Stereo savings */
1792      if (C==2)
1793      {
1794         int coded_stereo_bands;
1795         int coded_stereo_dof;
1796         coded_stereo_bands = IMIN(st->intensity, coded_bands);
1797         coded_stereo_dof = (eBands[coded_stereo_bands]<<LM)-coded_stereo_bands;
1798         /*printf("%d %d %d ", coded_stereo_dof, coded_bins, tot_boost);*/
1799         target -= MIN32(target/3, SHR16(MULT16_16(st->stereo_saving,(coded_stereo_dof<<BITRES)),8));
1800         target += MULT16_16_Q15(QCONST16(0.035,15),coded_stereo_dof<<BITRES);
1801      }
1802      /* Limits starving of other bands when using dynalloc */
1803      target += tot_boost;
1804      /* Compensates for the average transient boost */
1805      target = MULT16_32_Q15(QCONST16(0.96f,15),target);
1806      /* Apply transient boost */
1807      target = SHL32(MULT16_32_Q15(tf_estimate, target),1);
1808
1809 #ifndef FIXED_POINT
1810      /* Apply tonality boost */
1811      if (st->analysis.valid) {
1812         int tonal_target;
1813         float tonal;
1814
1815         /* Compensates for the average tonality boost */
1816         target -= MULT16_16_Q15(QCONST16(0.13f,15),coded_bins<<BITRES);
1817
1818         tonal = MAX16(0,st->analysis.tonality-.2);
1819         tonal_target = target + (coded_bins<<BITRES)*2.0f*tonal;
1820         if (pitch_change)
1821            tonal_target +=  (coded_bins<<BITRES)*.8;
1822         /*printf("%f %f ", st->analysis.tonality, tonal);*/
1823         target = IMAX(tonal_target,target);
1824      }
1825 #endif
1826
1827      {
1828         opus_int32 floor_depth;
1829         int bins;
1830         bins = eBands[nbEBands-2]<<LM;
1831         /*floor_depth = SHR32(MULT16_16((C*bins<<BITRES),celt_log2(SHL32(MAX16(1,sample_max),13))), DB_SHIFT);*/
1832         floor_depth = SHR32(MULT16_16((C*bins<<BITRES),maxDepth), DB_SHIFT);
1833         floor_depth = IMAX(floor_depth, target>>2);
1834         target = IMIN(target, floor_depth);
1835         /*printf("%f %d\n", maxDepth, floor_depth);*/
1836      }
1837
1838      if (st->constrained_vbr || st->bitrate<64000)
1839      {
1840         opus_val16 rate_factor;
1841 #ifdef FIXED_POINT
1842         rate_factor = MAX16(0,(st->bitrate-32000));
1843 #else
1844         rate_factor = MAX16(0,(1.f/32768)*(st->bitrate-32000));
1845 #endif
1846         if (st->constrained_vbr)
1847            rate_factor = MIN16(rate_factor, QCONST16(0.67f, 15));
1848         target = base_target + MULT16_32_Q15(rate_factor, target-base_target);
1849
1850      }
1851      /* Don't allow more than doubling the rate */
1852      target = IMIN(2*base_target, target);
1853
1854      /* The current offset is removed from the target and the space used
1855         so far is added*/
1856      target=target+tell;
1857      /* In VBR mode the frame size must not be reduced so much that it would
1858          result in the encoder running out of bits.
1859         The margin of 2 bytes ensures that none of the bust-prevention logic
1860          in the decoder will have triggered so far. */
1861      min_allowed = ((tell+total_boost+(1<<(BITRES+3))-1)>>(BITRES+3)) + 2 - nbFilledBytes;
1862
1863      nbAvailableBytes = (target+(1<<(BITRES+2)))>>(BITRES+3);
1864      nbAvailableBytes = IMAX(min_allowed,nbAvailableBytes);
1865      nbAvailableBytes = IMIN(nbCompressedBytes,nbAvailableBytes+nbFilledBytes) - nbFilledBytes;
1866
1867      /* By how much did we "miss" the target on that frame */
1868      delta = target - vbr_rate;
1869
1870      target=nbAvailableBytes<<(BITRES+3);
1871
1872      /*If the frame is silent we don't adjust our drift, otherwise
1873        the encoder will shoot to very high rates after hitting a
1874        span of silence, but we do allow the bitres to refill.
1875        This means that we'll undershoot our target in CVBR/VBR modes
1876        on files with lots of silence. */
1877      if(silence)
1878      {
1879        nbAvailableBytes = 2;
1880        target = 2*8<<BITRES;
1881        delta = 0;
1882      }
1883
1884      if (st->vbr_count < 970)
1885      {
1886         st->vbr_count++;
1887         alpha = celt_rcp(SHL32(EXTEND32(st->vbr_count+20),16));
1888      } else
1889         alpha = QCONST16(.001f,15);
1890      /* How many bits have we used in excess of what we're allowed */
1891      if (st->constrained_vbr)
1892         st->vbr_reservoir += target - vbr_rate;
1893      /*printf ("%d\n", st->vbr_reservoir);*/
1894
1895      /* Compute the offset we need to apply in order to reach the target */
1896      if (st->constrained_vbr)
1897      {
1898         st->vbr_drift += (opus_int32)MULT16_32_Q15(alpha,(delta*(1<<lm_diff))-st->vbr_offset-st->vbr_drift);
1899         st->vbr_offset = -st->vbr_drift;
1900      }
1901      /*printf ("%d\n", st->vbr_drift);*/
1902
1903      if (st->constrained_vbr && st->vbr_reservoir < 0)
1904      {
1905         /* We're under the min value -- increase rate */
1906         int adjust = (-st->vbr_reservoir)/(8<<BITRES);
1907         /* Unless we're just coding silence */
1908         nbAvailableBytes += silence?0:adjust;
1909         st->vbr_reservoir = 0;
1910         /*printf ("+%d\n", adjust);*/
1911      }
1912      nbCompressedBytes = IMIN(nbCompressedBytes,nbAvailableBytes+nbFilledBytes);
1913      /*printf("%d\n", nbCompressedBytes*50*8);*/
1914      /* This moves the raw bits to take into account the new compressed size */
1915      ec_enc_shrink(enc, nbCompressedBytes);
1916    }
1917
1918    /* Bit allocation */
1919    ALLOC(fine_quant, nbEBands, int);
1920    ALLOC(pulses, nbEBands, int);
1921    ALLOC(fine_priority, nbEBands, int);
1922
1923    /* bits =           packet size                    - where we are - safety*/
1924    bits = (((opus_int32)nbCompressedBytes*8)<<BITRES) - ec_tell_frac(enc) - 1;
1925    anti_collapse_rsv = isTransient&&LM>=2&&bits>=((LM+2)<<BITRES) ? (1<<BITRES) : 0;
1926    bits -= anti_collapse_rsv;
1927    codedBands = compute_allocation(mode, st->start, st->end, offsets, cap,
1928          alloc_trim, &st->intensity, &dual_stereo, bits, &balance, pulses,
1929          fine_quant, fine_priority, C, LM, enc, 1, st->lastCodedBands);
1930    st->lastCodedBands = codedBands;
1931
1932    quant_fine_energy(mode, st->start, st->end, oldBandE, error, fine_quant, enc, C);
1933
1934 #ifdef MEASURE_NORM_MSE
1935    float X0[3000];
1936    float bandE0[60];
1937    c=0; do
1938       for (i=0;i<N;i++)
1939          X0[i+c*N] = X[i+c*N];
1940    while (++c<C);
1941    for (i=0;i<C*nbEBands;i++)
1942       bandE0[i] = bandE[i];
1943 #endif
1944
1945    /* Residual quantisation */
1946    ALLOC(collapse_masks, C*nbEBands, unsigned char);
1947    quant_all_bands(1, mode, st->start, st->end, X, C==2 ? X+N : NULL, collapse_masks,
1948          bandE, pulses, shortBlocks, st->spread_decision, dual_stereo, st->intensity, tf_res,
1949          nbCompressedBytes*(8<<BITRES)-anti_collapse_rsv, balance, enc, LM, codedBands, &st->rng);
1950
1951    if (anti_collapse_rsv > 0)
1952    {
1953       anti_collapse_on = st->consec_transient<2;
1954 #ifdef FUZZING
1955       anti_collapse_on = rand()&0x1;
1956 #endif
1957       ec_enc_bits(enc, anti_collapse_on, 1);
1958    }
1959    quant_energy_finalise(mode, st->start, st->end, oldBandE, error, fine_quant, fine_priority, nbCompressedBytes*8-ec_tell(enc), enc, C);
1960
1961    if (silence)
1962    {
1963       for (i=0;i<C*nbEBands;i++)
1964          oldBandE[i] = -QCONST16(28.f,DB_SHIFT);
1965    }
1966
1967 #ifdef RESYNTH
1968    /* Re-synthesis of the coded audio if required */
1969    {
1970       celt_sig *out_mem[2];
1971
1972       log2Amp(mode, st->start, st->end, bandE, oldBandE, C);
1973       if (silence)
1974       {
1975          for (i=0;i<C*nbEBands;i++)
1976             bandE[i] = 0;
1977       }
1978
1979 #ifdef MEASURE_NORM_MSE
1980       measure_norm_mse(mode, X, X0, bandE, bandE0, M, N, C);
1981 #endif
1982       if (anti_collapse_on)
1983       {
1984          anti_collapse(mode, X, collapse_masks, LM, C, N,
1985                st->start, st->end, oldBandE, oldLogE, oldLogE2, pulses, st->rng);
1986       }
1987
1988       /* Synthesis */
1989       denormalise_bands(mode, X, freq, bandE, effEnd, C, M);
1990
1991       c=0; do {
1992          OPUS_MOVE(st->syn_mem[c], st->syn_mem[c]+N, 2*MAX_PERIOD-N+overlap);
1993       } while (++c<CC);
1994
1995       c=0; do
1996          for (i=0;i<M*eBands[st->start];i++)
1997             freq[c*N+i] = 0;
1998       while (++c<C);
1999       c=0; do
2000          for (i=M*eBands[st->end];i<N;i++)
2001             freq[c*N+i] = 0;
2002       while (++c<C);
2003
2004       if (CC==2&&C==1)
2005       {
2006          for (i=0;i<N;i++)
2007             freq[N+i] = freq[i];
2008       }
2009
2010       c=0; do {
2011          out_mem[c] = st->syn_mem[c]+2*MAX_PERIOD-N;
2012       } while (++c<CC);
2013
2014       compute_inv_mdcts(mode, shortBlocks, freq, out_mem, CC, LM);
2015
2016       c=0; do {
2017          st->prefilter_period=IMAX(st->prefilter_period, COMBFILTER_MINPERIOD);
2018          st->prefilter_period_old=IMAX(st->prefilter_period_old, COMBFILTER_MINPERIOD);
2019          comb_filter(out_mem[c], out_mem[c], st->prefilter_period_old, st->prefilter_period, mode->shortMdctSize,
2020                st->prefilter_gain_old, st->prefilter_gain, st->prefilter_tapset_old, st->prefilter_tapset,
2021                mode->window, st->overlap);
2022          if (LM!=0)
2023             comb_filter(out_mem[c]+mode->shortMdctSize, out_mem[c]+mode->shortMdctSize, st->prefilter_period, pitch_index, N-mode->shortMdctSize,
2024                   st->prefilter_gain, gain1, st->prefilter_tapset, prefilter_tapset,
2025                   mode->window, overlap);
2026       } while (++c<CC);
2027
2028       /* We reuse freq[] as scratch space for the de-emphasis */
2029       deemphasis(out_mem, (opus_val16*)pcm, N, CC, st->upsample, mode->preemph, st->preemph_memD, freq);
2030       st->prefilter_period_old = st->prefilter_period;
2031       st->prefilter_gain_old = st->prefilter_gain;
2032       st->prefilter_tapset_old = st->prefilter_tapset;
2033    }
2034 #endif
2035
2036    st->prefilter_period = pitch_index;
2037    st->prefilter_gain = gain1;
2038    st->prefilter_tapset = prefilter_tapset;
2039 #ifdef RESYNTH
2040    if (LM!=0)
2041    {
2042       st->prefilter_period_old = st->prefilter_period;
2043       st->prefilter_gain_old = st->prefilter_gain;
2044       st->prefilter_tapset_old = st->prefilter_tapset;
2045    }
2046 #endif
2047
2048    if (CC==2&&C==1) {
2049       for (i=0;i<nbEBands;i++)
2050          oldBandE[nbEBands+i]=oldBandE[i];
2051    }
2052
2053    if (!isTransient)
2054    {
2055       for (i=0;i<CC*nbEBands;i++)
2056          oldLogE2[i] = oldLogE[i];
2057       for (i=0;i<CC*nbEBands;i++)
2058          oldLogE[i] = oldBandE[i];
2059    } else {
2060       for (i=0;i<CC*nbEBands;i++)
2061          oldLogE[i] = MIN16(oldLogE[i], oldBandE[i]);
2062    }
2063    /* In case start or end were to change */
2064    c=0; do
2065    {
2066       for (i=0;i<st->start;i++)
2067       {
2068          oldBandE[c*nbEBands+i]=0;
2069          oldLogE[c*nbEBands+i]=oldLogE2[c*nbEBands+i]=-QCONST16(28.f,DB_SHIFT);
2070       }
2071       for (i=st->end;i<nbEBands;i++)
2072       {
2073          oldBandE[c*nbEBands+i]=0;
2074          oldLogE[c*nbEBands+i]=oldLogE2[c*nbEBands+i]=-QCONST16(28.f,DB_SHIFT);
2075       }
2076    } while (++c<CC);
2077
2078    if (isTransient)
2079       st->consec_transient++;
2080    else
2081       st->consec_transient=0;
2082    st->rng = enc->rng;
2083
2084    /* If there's any room left (can only happen for very high rates),
2085       it's already filled with zeros */
2086    ec_enc_done(enc);
2087
2088 #ifdef CUSTOM_MODES
2089    if (st->signalling)
2090       nbCompressedBytes++;
2091 #endif
2092
2093    RESTORE_STACK;
2094    if (ec_get_error(enc))
2095       return OPUS_INTERNAL_ERROR;
2096    else
2097       return nbCompressedBytes;
2098 }
2099
2100
2101 #ifdef CUSTOM_MODES
2102
2103 #ifdef FIXED_POINT
2104 int opus_custom_encode(CELTEncoder * OPUS_RESTRICT st, const opus_int16 * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes)
2105 {
2106    return celt_encode_with_ec(st, pcm, frame_size, compressed, nbCompressedBytes, NULL);
2107 }
2108
2109 #ifndef DISABLE_FLOAT_API
2110 int opus_custom_encode_float(CELTEncoder * OPUS_RESTRICT st, const float * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes)
2111 {
2112    int j, ret, C, N;
2113    VARDECL(opus_int16, in);
2114    ALLOC_STACK;
2115
2116    if (pcm==NULL)
2117       return OPUS_BAD_ARG;
2118
2119    C = st->channels;
2120    N = frame_size;
2121    ALLOC(in, C*N, opus_int16);
2122
2123    for (j=0;j<C*N;j++)
2124      in[j] = FLOAT2INT16(pcm[j]);
2125
2126    ret=celt_encode_with_ec(st,in,frame_size,compressed,nbCompressedBytes, NULL);
2127 #ifdef RESYNTH
2128    for (j=0;j<C*N;j++)
2129       ((float*)pcm)[j]=in[j]*(1.f/32768.f);
2130 #endif
2131    RESTORE_STACK;
2132    return ret;
2133 }
2134 #endif /* DISABLE_FLOAT_API */
2135 #else
2136
2137 int opus_custom_encode(CELTEncoder * OPUS_RESTRICT st, const opus_int16 * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes)
2138 {
2139    int j, ret, C, N;
2140    VARDECL(celt_sig, in);
2141    ALLOC_STACK;
2142
2143    if (pcm==NULL)
2144       return OPUS_BAD_ARG;
2145
2146    C=st->channels;
2147    N=frame_size;
2148    ALLOC(in, C*N, celt_sig);
2149    for (j=0;j<C*N;j++) {
2150      in[j] = SCALEOUT(pcm[j]);
2151    }
2152
2153    ret = celt_encode_with_ec(st,in,frame_size,compressed,nbCompressedBytes, NULL);
2154 #ifdef RESYNTH
2155    for (j=0;j<C*N;j++)
2156       ((opus_int16*)pcm)[j] = FLOAT2INT16(in[j]);
2157 #endif
2158    RESTORE_STACK;
2159    return ret;
2160 }
2161
2162 int opus_custom_encode_float(CELTEncoder * OPUS_RESTRICT st, const float * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes)
2163 {
2164    return celt_encode_with_ec(st, pcm, frame_size, compressed, nbCompressedBytes, NULL);
2165 }
2166
2167 #endif
2168
2169 #endif /* CUSTOM_MODES */
2170
2171 int opus_custom_encoder_ctl(CELTEncoder * OPUS_RESTRICT st, int request, ...)
2172 {
2173    va_list ap;
2174
2175    va_start(ap, request);
2176    switch (request)
2177    {
2178       case OPUS_SET_COMPLEXITY_REQUEST:
2179       {
2180          int value = va_arg(ap, opus_int32);
2181          if (value<0 || value>10)
2182             goto bad_arg;
2183          st->complexity = value;
2184       }
2185       break;
2186       case CELT_SET_START_BAND_REQUEST:
2187       {
2188          opus_int32 value = va_arg(ap, opus_int32);
2189          if (value<0 || value>=st->mode->nbEBands)
2190             goto bad_arg;
2191          st->start = value;
2192       }
2193       break;
2194       case CELT_SET_END_BAND_REQUEST:
2195       {
2196          opus_int32 value = va_arg(ap, opus_int32);
2197          if (value<1 || value>st->mode->nbEBands)
2198             goto bad_arg;
2199          st->end = value;
2200       }
2201       break;
2202       case CELT_SET_PREDICTION_REQUEST:
2203       {
2204          int value = va_arg(ap, opus_int32);
2205          if (value<0 || value>2)
2206             goto bad_arg;
2207          st->disable_pf = value<=1;
2208          st->force_intra = value==0;
2209       }
2210       break;
2211       case OPUS_SET_PACKET_LOSS_PERC_REQUEST:
2212       {
2213          int value = va_arg(ap, opus_int32);
2214          if (value<0 || value>100)
2215             goto bad_arg;
2216          st->loss_rate = value;
2217       }
2218       break;
2219       case OPUS_SET_VBR_CONSTRAINT_REQUEST:
2220       {
2221          opus_int32 value = va_arg(ap, opus_int32);
2222          st->constrained_vbr = value;
2223       }
2224       break;
2225       case OPUS_SET_VBR_REQUEST:
2226       {
2227          opus_int32 value = va_arg(ap, opus_int32);
2228          st->vbr = value;
2229       }
2230       break;
2231       case OPUS_SET_BITRATE_REQUEST:
2232       {
2233          opus_int32 value = va_arg(ap, opus_int32);
2234          if (value<=500 && value!=OPUS_BITRATE_MAX)
2235             goto bad_arg;
2236          value = IMIN(value, 260000*st->channels);
2237          st->bitrate = value;
2238       }
2239       break;
2240       case CELT_SET_CHANNELS_REQUEST:
2241       {
2242          opus_int32 value = va_arg(ap, opus_int32);
2243          if (value<1 || value>2)
2244             goto bad_arg;
2245          st->stream_channels = value;
2246       }
2247       break;
2248       case OPUS_SET_LSB_DEPTH_REQUEST:
2249       {
2250           opus_int32 value = va_arg(ap, opus_int32);
2251           if (value<8 || value>24)
2252              goto bad_arg;
2253           st->lsb_depth=value;
2254       }
2255       break;
2256       case OPUS_GET_LSB_DEPTH_REQUEST:
2257       {
2258           opus_int32 *value = va_arg(ap, opus_int32*);
2259           *value=st->lsb_depth;
2260       }
2261       break;
2262       case OPUS_RESET_STATE:
2263       {
2264          int i;
2265          opus_val16 *oldBandE, *oldLogE, *oldLogE2;
2266          oldBandE = (opus_val16*)(st->in_mem+st->channels*(st->overlap+COMBFILTER_MAXPERIOD));
2267          oldLogE = oldBandE + st->channels*st->mode->nbEBands;
2268          oldLogE2 = oldLogE + st->channels*st->mode->nbEBands;
2269          OPUS_CLEAR((char*)&st->ENCODER_RESET_START,
2270                opus_custom_encoder_get_size(st->mode, st->channels)-
2271                ((char*)&st->ENCODER_RESET_START - (char*)st));
2272          for (i=0;i<st->channels*st->mode->nbEBands;i++)
2273             oldLogE[i]=oldLogE2[i]=-QCONST16(28.f,DB_SHIFT);
2274          st->vbr_offset = 0;
2275          st->delayedIntra = 1;
2276          st->spread_decision = SPREAD_NORMAL;
2277          st->tonal_average = 256;
2278          st->hf_average = 0;
2279          st->tapset_decision = 0;
2280       }
2281       break;
2282 #ifdef CUSTOM_MODES
2283       case CELT_SET_INPUT_CLIPPING_REQUEST:
2284       {
2285          opus_int32 value = va_arg(ap, opus_int32);
2286          st->clip = value;
2287       }
2288       break;
2289 #endif
2290       case CELT_SET_SIGNALLING_REQUEST:
2291       {
2292          opus_int32 value = va_arg(ap, opus_int32);
2293          st->signalling = value;
2294       }
2295       break;
2296       case CELT_SET_ANALYSIS_REQUEST:
2297       {
2298          AnalysisInfo *info = va_arg(ap, AnalysisInfo *);
2299          if (info)
2300             OPUS_COPY(&st->analysis, info, 1);
2301       }
2302       break;
2303       case CELT_GET_MODE_REQUEST:
2304       {
2305          const CELTMode ** value = va_arg(ap, const CELTMode**);
2306          if (value==0)
2307             goto bad_arg;
2308          *value=st->mode;
2309       }
2310       break;
2311       case OPUS_GET_FINAL_RANGE_REQUEST:
2312       {
2313          opus_uint32 * value = va_arg(ap, opus_uint32 *);
2314          if (value==0)
2315             goto bad_arg;
2316          *value=st->rng;
2317       }
2318       break;
2319       default:
2320          goto bad_request;
2321    }
2322    va_end(ap);
2323    return OPUS_OK;
2324 bad_arg:
2325    va_end(ap);
2326    return OPUS_BAD_ARG;
2327 bad_request:
2328    va_end(ap);
2329    return OPUS_UNIMPLEMENTED;
2330 }
2331
2332 /**********************************************************************/
2333 /*                                                                    */
2334 /*                             DECODER                                */
2335 /*                                                                    */
2336 /**********************************************************************/
2337 #define DECODE_BUFFER_SIZE 2048
2338
2339 /** Decoder state
2340  @brief Decoder state
2341  */
2342 struct OpusCustomDecoder {
2343    const OpusCustomMode *mode;
2344    int overlap;
2345    int channels;
2346    int stream_channels;
2347
2348    int downsample;
2349    int start, end;
2350    int signalling;
2351
2352    /* Everything beyond this point gets cleared on a reset */
2353 #define DECODER_RESET_START rng
2354
2355    opus_uint32 rng;
2356    int error;
2357    int last_pitch_index;
2358    int loss_count;
2359    int postfilter_period;
2360    int postfilter_period_old;
2361    opus_val16 postfilter_gain;
2362    opus_val16 postfilter_gain_old;
2363    int postfilter_tapset;
2364    int postfilter_tapset_old;
2365
2366    celt_sig preemph_memD[2];
2367
2368    celt_sig _decode_mem[1]; /* Size = channels*(DECODE_BUFFER_SIZE+mode->overlap) */
2369    /* opus_val16 lpc[],  Size = channels*LPC_ORDER */
2370    /* opus_val16 oldEBands[], Size = 2*mode->nbEBands */
2371    /* opus_val16 oldLogE[], Size = 2*mode->nbEBands */
2372    /* opus_val16 oldLogE2[], Size = 2*mode->nbEBands */
2373    /* opus_val16 backgroundLogE[], Size = 2*mode->nbEBands */
2374 };
2375
2376 int celt_decoder_get_size(int channels)
2377 {
2378    const CELTMode *mode = opus_custom_mode_create(48000, 960, NULL);
2379    return opus_custom_decoder_get_size(mode, channels);
2380 }
2381
2382 OPUS_CUSTOM_NOSTATIC int opus_custom_decoder_get_size(const CELTMode *mode, int channels)
2383 {
2384    int size = sizeof(struct CELTDecoder)
2385             + (channels*(DECODE_BUFFER_SIZE+mode->overlap)-1)*sizeof(celt_sig)
2386             + channels*LPC_ORDER*sizeof(opus_val16)
2387             + 4*2*mode->nbEBands*sizeof(opus_val16);
2388    return size;
2389 }
2390
2391 #ifdef CUSTOM_MODES
2392 CELTDecoder *opus_custom_decoder_create(const CELTMode *mode, int channels, int *error)
2393 {
2394    int ret;
2395    CELTDecoder *st = (CELTDecoder *)opus_alloc(opus_custom_decoder_get_size(mode, channels));
2396    ret = opus_custom_decoder_init(st, mode, channels);
2397    if (ret != OPUS_OK)
2398    {
2399       opus_custom_decoder_destroy(st);
2400       st = NULL;
2401    }
2402    if (error)
2403       *error = ret;
2404    return st;
2405 }
2406 #endif /* CUSTOM_MODES */
2407
2408 int celt_decoder_init(CELTDecoder *st, opus_int32 sampling_rate, int channels)
2409 {
2410    int ret;
2411    ret = opus_custom_decoder_init(st, opus_custom_mode_create(48000, 960, NULL), channels);
2412    if (ret != OPUS_OK)
2413       return ret;
2414    st->downsample = resampling_factor(sampling_rate);
2415    if (st->downsample==0)
2416       return OPUS_BAD_ARG;
2417    else
2418       return OPUS_OK;
2419 }
2420
2421 OPUS_CUSTOM_NOSTATIC int opus_custom_decoder_init(CELTDecoder *st, const CELTMode *mode, int channels)
2422 {
2423    if (channels < 0 || channels > 2)
2424       return OPUS_BAD_ARG;
2425
2426    if (st==NULL)
2427       return OPUS_ALLOC_FAIL;
2428
2429    OPUS_CLEAR((char*)st, opus_custom_decoder_get_size(mode, channels));
2430
2431    st->mode = mode;
2432    st->overlap = mode->overlap;
2433    st->stream_channels = st->channels = channels;
2434
2435    st->downsample = 1;
2436    st->start = 0;
2437    st->end = st->mode->effEBands;
2438    st->signalling = 1;
2439
2440    st->loss_count = 0;
2441
2442    opus_custom_decoder_ctl(st, OPUS_RESET_STATE);
2443
2444    return OPUS_OK;
2445 }
2446
2447 #ifdef CUSTOM_MODES
2448 void opus_custom_decoder_destroy(CELTDecoder *st)
2449 {
2450    opus_free(st);
2451 }
2452 #endif /* CUSTOM_MODES */
2453
2454 static void celt_decode_lost(CELTDecoder * OPUS_RESTRICT st, opus_val16 * OPUS_RESTRICT pcm, int N, int LM)
2455 {
2456    int c;
2457    int pitch_index;
2458    opus_val16 fade = Q15ONE;
2459    int i, len;
2460    const int C = st->channels;
2461    int offset;
2462    celt_sig *out_mem[2];
2463    celt_sig *decode_mem[2];
2464    opus_val16 *lpc;
2465    opus_val32 *out_syn[2];
2466    opus_val16 *oldBandE, *oldLogE, *oldLogE2, *backgroundLogE;
2467    const OpusCustomMode *mode;
2468    int nbEBands;
2469    int overlap;
2470    const opus_int16 *eBands;
2471    VARDECL(celt_sig, scratch);
2472    SAVE_STACK;
2473
2474    mode = st->mode;
2475    nbEBands = mode->nbEBands;
2476    overlap = mode->overlap;
2477    eBands = mode->eBands;
2478
2479    c=0; do {
2480       decode_mem[c] = st->_decode_mem + c*(DECODE_BUFFER_SIZE+st->overlap);
2481       out_mem[c] = decode_mem[c]+DECODE_BUFFER_SIZE-MAX_PERIOD;
2482    } while (++c<C);
2483    lpc = (opus_val16*)(st->_decode_mem+(DECODE_BUFFER_SIZE+st->overlap)*C);
2484    oldBandE = lpc+C*LPC_ORDER;
2485    oldLogE = oldBandE + 2*nbEBands;
2486    oldLogE2 = oldLogE + 2*nbEBands;
2487    backgroundLogE = oldLogE2  + 2*nbEBands;
2488
2489    c=0; do {
2490       out_syn[c] = out_mem[c]+MAX_PERIOD-N;
2491    } while (++c<C);
2492
2493    len = N+overlap;
2494
2495    if (st->loss_count >= 5 || st->start!=0)
2496    {
2497       /* Noise-based PLC/CNG */
2498       VARDECL(celt_sig, freq);
2499       VARDECL(celt_norm, X);
2500       VARDECL(celt_ener, bandE);
2501       opus_uint32 seed;
2502       int effEnd;
2503
2504       effEnd = st->end;
2505       if (effEnd > mode->effEBands)
2506          effEnd = mode->effEBands;
2507
2508       ALLOC(freq, C*N, celt_sig); /**< Interleaved signal MDCTs */
2509       ALLOC(X, C*N, celt_norm);   /**< Interleaved normalised MDCTs */
2510       ALLOC(bandE, nbEBands*C, celt_ener);
2511
2512       if (st->loss_count >= 5)
2513          log2Amp(mode, st->start, st->end, bandE, backgroundLogE, C);
2514       else {
2515          /* Energy decay */
2516          opus_val16 decay = st->loss_count==0 ? QCONST16(1.5f, DB_SHIFT) : QCONST16(.5f, DB_SHIFT);
2517          c=0; do
2518          {
2519             for (i=st->start;i<st->end;i++)
2520                oldBandE[c*nbEBands+i] -= decay;
2521          } while (++c<C);
2522          log2Amp(mode, st->start, st->end, bandE, oldBandE, C);
2523       }
2524       seed = st->rng;
2525       for (c=0;c<C;c++)
2526       {
2527          for (i=0;i<(eBands[st->start]<<LM);i++)
2528             X[c*N+i] = 0;
2529          for (i=st->start;i<mode->effEBands;i++)
2530          {
2531             int j;
2532             int boffs;
2533             int blen;
2534             boffs = N*c+(eBands[i]<<LM);
2535             blen = (eBands[i+1]-eBands[i])<<LM;
2536             for (j=0;j<blen;j++)
2537             {
2538                seed = celt_lcg_rand(seed);
2539                X[boffs+j] = (celt_norm)((opus_int32)seed>>20);
2540             }
2541             renormalise_vector(X+boffs, blen, Q15ONE);
2542          }
2543          for (i=(eBands[st->end]<<LM);i<N;i++)
2544             X[c*N+i] = 0;
2545       }
2546       st->rng = seed;
2547
2548       denormalise_bands(mode, X, freq, bandE, mode->effEBands, C, 1<<LM);
2549
2550       c=0; do
2551          for (i=0;i<eBands[st->start]<<LM;i++)
2552             freq[c*N+i] = 0;
2553       while (++c<C);
2554       c=0; do {
2555          int bound = eBands[effEnd]<<LM;
2556          if (st->downsample!=1)
2557             bound = IMIN(bound, N/st->downsample);
2558          for (i=bound;i<N;i++)
2559             freq[c*N+i] = 0;
2560       } while (++c<C);
2561       c=0; do {
2562          OPUS_MOVE(decode_mem[c], decode_mem[c]+N, DECODE_BUFFER_SIZE-N+overlap);
2563       } while (++c<C);
2564       compute_inv_mdcts(mode, 0, freq, out_syn, C, LM);
2565    } else {
2566       /* Pitch-based PLC */
2567       VARDECL(opus_val32, e);
2568
2569       if (st->loss_count == 0)
2570       {
2571          opus_val16 pitch_buf[DECODE_BUFFER_SIZE>>1];
2572          /* Corresponds to a min pitch of 67 Hz. It's possible to save CPU in this
2573          search by using only part of the decode buffer */
2574          int poffset = 720;
2575          pitch_downsample(decode_mem, pitch_buf, DECODE_BUFFER_SIZE, C);
2576          /* Max pitch is 100 samples (480 Hz) */
2577          pitch_search(pitch_buf+((poffset)>>1), pitch_buf, DECODE_BUFFER_SIZE-poffset,
2578                poffset-100, &pitch_index);
2579          pitch_index = poffset-pitch_index;
2580          st->last_pitch_index = pitch_index;
2581       } else {
2582          pitch_index = st->last_pitch_index;
2583          fade = QCONST16(.8f,15);
2584       }
2585
2586       ALLOC(e, MAX_PERIOD+2*overlap, opus_val32);
2587       c=0; do {
2588          opus_val16 exc[MAX_PERIOD];
2589          opus_val32 ac[LPC_ORDER+1];
2590          opus_val16 decay = 1;
2591          opus_val32 S1=0;
2592          opus_val16 mem[LPC_ORDER]={0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0};
2593
2594          offset = MAX_PERIOD-pitch_index;
2595          for (i=0;i<MAX_PERIOD;i++)
2596             exc[i] = ROUND16(out_mem[c][i], SIG_SHIFT);
2597
2598          if (st->loss_count == 0)
2599          {
2600             _celt_autocorr(exc, ac, mode->window, overlap,
2601                   LPC_ORDER, MAX_PERIOD);
2602
2603             /* Noise floor -40 dB */
2604 #ifdef FIXED_POINT
2605             ac[0] += SHR32(ac[0],13);
2606 #else
2607             ac[0] *= 1.0001f;
2608 #endif
2609             /* Lag windowing */
2610             for (i=1;i<=LPC_ORDER;i++)
2611             {
2612                /*ac[i] *= exp(-.5*(2*M_PI*.002*i)*(2*M_PI*.002*i));*/
2613 #ifdef FIXED_POINT
2614                ac[i] -= MULT16_32_Q15(2*i*i, ac[i]);
2615 #else
2616                ac[i] -= ac[i]*(.008f*i)*(.008f*i);
2617 #endif
2618             }
2619
2620             _celt_lpc(lpc+c*LPC_ORDER, ac, LPC_ORDER);
2621          }
2622          for (i=0;i<LPC_ORDER;i++)
2623             mem[i] = ROUND16(out_mem[c][MAX_PERIOD-1-i], SIG_SHIFT);
2624          celt_fir(exc, lpc+c*LPC_ORDER, exc, MAX_PERIOD, LPC_ORDER, mem);
2625          /*for (i=0;i<MAX_PERIOD;i++)printf("%d ", exc[i]); printf("\n");*/
2626          /* Check if the waveform is decaying (and if so how fast) */
2627          {
2628             opus_val32 E1=1, E2=1;
2629             int period;
2630             if (pitch_index <= MAX_PERIOD/2)
2631                period = pitch_index;
2632             else
2633                period = MAX_PERIOD/2;
2634             for (i=0;i<period;i++)
2635             {
2636                E1 += SHR32(MULT16_16(exc[MAX_PERIOD-period+i],exc[MAX_PERIOD-period+i]),8);
2637                E2 += SHR32(MULT16_16(exc[MAX_PERIOD-2*period+i],exc[MAX_PERIOD-2*period+i]),8);
2638             }
2639             if (E1 > E2)
2640                E1 = E2;
2641             decay = celt_sqrt(frac_div32(SHR32(E1,1),E2));
2642          }
2643
2644          /* Copy excitation, taking decay into account */
2645          for (i=0;i<len+overlap;i++)
2646          {
2647             opus_val16 tmp;
2648             if (offset+i >= MAX_PERIOD)
2649             {
2650                offset -= pitch_index;
2651                decay = MULT16_16_Q15(decay, decay);
2652             }
2653             e[i] = SHL32(EXTEND32(MULT16_16_Q15(decay, exc[offset+i])), SIG_SHIFT);
2654             tmp = ROUND16(out_mem[c][offset+i],SIG_SHIFT);
2655             S1 += SHR32(MULT16_16(tmp,tmp),8);
2656          }
2657          for (i=0;i<LPC_ORDER;i++)
2658             mem[i] = ROUND16(out_mem[c][MAX_PERIOD-1-i], SIG_SHIFT);
2659          for (i=0;i<len+overlap;i++)
2660             e[i] = MULT16_32_Q15(fade, e[i]);
2661          celt_iir(e, lpc+c*LPC_ORDER, e, len+overlap, LPC_ORDER, mem);
2662
2663          {
2664             opus_val32 S2=0;
2665             for (i=0;i<len+overlap;i++)
2666             {
2667                opus_val16 tmp = ROUND16(e[i],SIG_SHIFT);
2668                S2 += SHR32(MULT16_16(tmp,tmp),8);
2669             }
2670             /* This checks for an "explosion" in the synthesis */
2671 #ifdef FIXED_POINT
2672             if (!(S1 > SHR32(S2,2)))
2673 #else
2674                /* Float test is written this way to catch NaNs at the same time */
2675                if (!(S1 > 0.2f*S2))
2676 #endif
2677                {
2678                   for (i=0;i<len+overlap;i++)
2679                      e[i] = 0;
2680                } else if (S1 < S2)
2681                {
2682                   opus_val16 ratio = celt_sqrt(frac_div32(SHR32(S1,1)+1,S2+1));
2683                   for (i=0;i<len+overlap;i++)
2684                      e[i] = MULT16_32_Q15(ratio, e[i]);
2685                }
2686          }
2687
2688          /* Apply post-filter to the MDCT overlap of the previous frame */
2689          comb_filter(out_mem[c]+MAX_PERIOD, out_mem[c]+MAX_PERIOD, st->postfilter_period, st->postfilter_period, st->overlap,
2690                st->postfilter_gain, st->postfilter_gain, st->postfilter_tapset, st->postfilter_tapset,
2691                NULL, 0);
2692
2693          for (i=0;i<MAX_PERIOD+overlap-N;i++)
2694             out_mem[c][i] = out_mem[c][N+i];
2695
2696          /* Apply TDAC to the concealed audio so that it blends with the
2697          previous and next frames */
2698          for (i=0;i<overlap/2;i++)
2699          {
2700             opus_val32 tmp;
2701             tmp = MULT16_32_Q15(mode->window[i],           e[N+overlap-1-i]) +
2702                   MULT16_32_Q15(mode->window[overlap-i-1], e[N+i          ]);
2703             out_mem[c][MAX_PERIOD+i] = MULT16_32_Q15(mode->window[overlap-i-1], tmp);
2704             out_mem[c][MAX_PERIOD+overlap-i-1] = MULT16_32_Q15(mode->window[i], tmp);
2705          }
2706          for (i=0;i<N;i++)
2707             out_mem[c][MAX_PERIOD-N+i] = e[i];
2708
2709          /* Apply pre-filter to the MDCT overlap for the next frame (post-filter will be applied then) */
2710          comb_filter(e, out_mem[c]+MAX_PERIOD, st->postfilter_period, st->postfilter_period, st->overlap,
2711                -st->postfilter_gain, -st->postfilter_gain, st->postfilter_tapset, st->postfilter_tapset,
2712                NULL, 0);
2713          for (i=0;i<overlap;i++)
2714             out_mem[c][MAX_PERIOD+i] = e[i];
2715       } while (++c<C);
2716    }
2717
2718    ALLOC(scratch, N, celt_sig);
2719    deemphasis(out_syn, pcm, N, C, st->downsample, mode->preemph, st->preemph_memD, scratch);
2720
2721    st->loss_count++;
2722
2723    RESTORE_STACK;
2724 }
2725
2726 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)
2727 {
2728    int c, i, N;
2729    int spread_decision;
2730    opus_int32 bits;
2731    ec_dec _dec;
2732    VARDECL(celt_sig, freq);
2733    VARDECL(celt_norm, X);
2734    VARDECL(celt_ener, bandE);
2735    VARDECL(int, fine_quant);
2736    VARDECL(int, pulses);
2737    VARDECL(int, cap);
2738    VARDECL(int, offsets);
2739    VARDECL(int, fine_priority);
2740    VARDECL(int, tf_res);
2741    VARDECL(unsigned char, collapse_masks);
2742    celt_sig *out_mem[2];
2743    celt_sig *decode_mem[2];
2744    celt_sig *out_syn[2];
2745    opus_val16 *lpc;
2746    opus_val16 *oldBandE, *oldLogE, *oldLogE2, *backgroundLogE;
2747
2748    int shortBlocks;
2749    int isTransient;
2750    int intra_ener;
2751    const int CC = st->channels;
2752    int LM, M;
2753    int effEnd;
2754    int codedBands;
2755    int alloc_trim;
2756    int postfilter_pitch;
2757    opus_val16 postfilter_gain;
2758    int intensity=0;
2759    int dual_stereo=0;
2760    opus_int32 total_bits;
2761    opus_int32 balance;
2762    opus_int32 tell;
2763    int dynalloc_logp;
2764    int postfilter_tapset;
2765    int anti_collapse_rsv;
2766    int anti_collapse_on=0;
2767    int silence;
2768    int C = st->stream_channels;
2769    const OpusCustomMode *mode;
2770    int nbEBands;
2771    int overlap;
2772    const opus_int16 *eBands;
2773    ALLOC_STACK;
2774
2775    mode = st->mode;
2776    nbEBands = mode->nbEBands;
2777    overlap = mode->overlap;
2778    eBands = mode->eBands;
2779    frame_size *= st->downsample;
2780
2781    c=0; do {
2782       decode_mem[c] = st->_decode_mem + c*(DECODE_BUFFER_SIZE+overlap);
2783       out_mem[c] = decode_mem[c]+DECODE_BUFFER_SIZE-MAX_PERIOD;
2784    } while (++c<CC);
2785    lpc = (opus_val16*)(st->_decode_mem+(DECODE_BUFFER_SIZE+overlap)*CC);
2786    oldBandE = lpc+CC*LPC_ORDER;
2787    oldLogE = oldBandE + 2*nbEBands;
2788    oldLogE2 = oldLogE + 2*nbEBands;
2789    backgroundLogE = oldLogE2  + 2*nbEBands;
2790
2791 #ifdef CUSTOM_MODES
2792    if (st->signalling && data!=NULL)
2793    {
2794       int data0=data[0];
2795       /* Convert "standard mode" to Opus header */
2796       if (mode->Fs==48000 && mode->shortMdctSize==120)
2797       {
2798          data0 = fromOpus(data0);
2799          if (data0<0)
2800             return OPUS_INVALID_PACKET;
2801       }
2802       st->end = IMAX(1, mode->effEBands-2*(data0>>5));
2803       LM = (data0>>3)&0x3;
2804       C = 1 + ((data0>>2)&0x1);
2805       data++;
2806       len--;
2807       if (LM>mode->maxLM)
2808          return OPUS_INVALID_PACKET;
2809       if (frame_size < mode->shortMdctSize<<LM)
2810          return OPUS_BUFFER_TOO_SMALL;
2811       else
2812          frame_size = mode->shortMdctSize<<LM;
2813    } else {
2814 #else
2815    {
2816 #endif
2817       for (LM=0;LM<=mode->maxLM;LM++)
2818          if (mode->shortMdctSize<<LM==frame_size)
2819             break;
2820       if (LM>mode->maxLM)
2821          return OPUS_BAD_ARG;
2822    }
2823    M=1<<LM;
2824
2825    if (len<0 || len>1275 || pcm==NULL)
2826       return OPUS_BAD_ARG;
2827
2828    N = M*mode->shortMdctSize;
2829
2830    effEnd = st->end;
2831    if (effEnd > mode->effEBands)
2832       effEnd = mode->effEBands;
2833
2834    if (data == NULL || len<=1)
2835    {
2836       celt_decode_lost(st, pcm, N, LM);
2837       RESTORE_STACK;
2838       return frame_size/st->downsample;
2839    }
2840
2841    ALLOC(freq, IMAX(CC,C)*N, celt_sig); /**< Interleaved signal MDCTs */
2842    ALLOC(X, C*N, celt_norm);   /**< Interleaved normalised MDCTs */
2843    ALLOC(bandE, nbEBands*C, celt_ener);
2844    c=0; do
2845       for (i=0;i<M*eBands[st->start];i++)
2846          X[c*N+i] = 0;
2847    while (++c<C);
2848    c=0; do
2849       for (i=M*eBands[effEnd];i<N;i++)
2850          X[c*N+i] = 0;
2851    while (++c<C);
2852
2853    if (dec == NULL)
2854    {
2855       ec_dec_init(&_dec,(unsigned char*)data,len);
2856       dec = &_dec;
2857    }
2858
2859    if (C==1)
2860    {
2861       for (i=0;i<nbEBands;i++)
2862          oldBandE[i]=MAX16(oldBandE[i],oldBandE[nbEBands+i]);
2863    }
2864
2865    total_bits = len*8;
2866    tell = ec_tell(dec);
2867
2868    if (tell >= total_bits)
2869       silence = 1;
2870    else if (tell==1)
2871       silence = ec_dec_bit_logp(dec, 15);
2872    else
2873       silence = 0;
2874    if (silence)
2875    {
2876       /* Pretend we've read all the remaining bits */
2877       tell = len*8;
2878       dec->nbits_total+=tell-ec_tell(dec);
2879    }
2880
2881    postfilter_gain = 0;
2882    postfilter_pitch = 0;
2883    postfilter_tapset = 0;
2884    if (st->start==0 && tell+16 <= total_bits)
2885    {
2886       if(ec_dec_bit_logp(dec, 1))
2887       {
2888          int qg, octave;
2889          octave = ec_dec_uint(dec, 6);
2890          postfilter_pitch = (16<<octave)+ec_dec_bits(dec, 4+octave)-1;
2891          qg = ec_dec_bits(dec, 3);
2892          if (ec_tell(dec)+2<=total_bits)
2893             postfilter_tapset = ec_dec_icdf(dec, tapset_icdf, 2);
2894          postfilter_gain = QCONST16(.09375f,15)*(qg+1);
2895       }
2896       tell = ec_tell(dec);
2897    }
2898
2899    if (LM > 0 && tell+3 <= total_bits)
2900    {
2901       isTransient = ec_dec_bit_logp(dec, 3);
2902       tell = ec_tell(dec);
2903    }
2904    else
2905       isTransient = 0;
2906
2907    if (isTransient)
2908       shortBlocks = M;
2909    else
2910       shortBlocks = 0;
2911
2912    /* Decode the global flags (first symbols in the stream) */
2913    intra_ener = tell+3<=total_bits ? ec_dec_bit_logp(dec, 3) : 0;
2914    /* Get band energies */
2915    unquant_coarse_energy(mode, st->start, st->end, oldBandE,
2916          intra_ener, dec, C, LM);
2917
2918    ALLOC(tf_res, nbEBands, int);
2919    tf_decode(st->start, st->end, isTransient, tf_res, LM, dec);
2920
2921    tell = ec_tell(dec);
2922    spread_decision = SPREAD_NORMAL;
2923    if (tell+4 <= total_bits)
2924       spread_decision = ec_dec_icdf(dec, spread_icdf, 5);
2925
2926    ALLOC(pulses, nbEBands, int);
2927    ALLOC(cap, nbEBands, int);
2928    ALLOC(offsets, nbEBands, int);
2929    ALLOC(fine_priority, nbEBands, int);
2930
2931    init_caps(mode,cap,LM,C);
2932
2933    dynalloc_logp = 6;
2934    total_bits<<=BITRES;
2935    tell = ec_tell_frac(dec);
2936    for (i=st->start;i<st->end;i++)
2937    {
2938       int width, quanta;
2939       int dynalloc_loop_logp;
2940       int boost;
2941       width = C*(eBands[i+1]-eBands[i])<<LM;
2942       /* quanta is 6 bits, but no more than 1 bit/sample
2943          and no less than 1/8 bit/sample */
2944       quanta = IMIN(width<<BITRES, IMAX(6<<BITRES, width));
2945       dynalloc_loop_logp = dynalloc_logp;
2946       boost = 0;
2947       while (tell+(dynalloc_loop_logp<<BITRES) < total_bits && boost < cap[i])
2948       {
2949          int flag;
2950          flag = ec_dec_bit_logp(dec, dynalloc_loop_logp);
2951          tell = ec_tell_frac(dec);
2952          if (!flag)
2953             break;
2954          boost += quanta;
2955          total_bits -= quanta;
2956          dynalloc_loop_logp = 1;
2957       }
2958       offsets[i] = boost;
2959       /* Making dynalloc more likely */
2960       if (boost>0)
2961          dynalloc_logp = IMAX(2, dynalloc_logp-1);
2962    }
2963
2964    ALLOC(fine_quant, nbEBands, int);
2965    alloc_trim = tell+(6<<BITRES) <= total_bits ?
2966          ec_dec_icdf(dec, trim_icdf, 7) : 5;
2967
2968    bits = (((opus_int32)len*8)<<BITRES) - ec_tell_frac(dec) - 1;
2969    anti_collapse_rsv = isTransient&&LM>=2&&bits>=((LM+2)<<BITRES) ? (1<<BITRES) : 0;
2970    bits -= anti_collapse_rsv;
2971    codedBands = compute_allocation(mode, st->start, st->end, offsets, cap,
2972          alloc_trim, &intensity, &dual_stereo, bits, &balance, pulses,
2973          fine_quant, fine_priority, C, LM, dec, 0, 0);
2974
2975    unquant_fine_energy(mode, st->start, st->end, oldBandE, fine_quant, dec, C);
2976
2977    /* Decode fixed codebook */
2978    ALLOC(collapse_masks, C*nbEBands, unsigned char);
2979    quant_all_bands(0, mode, st->start, st->end, X, C==2 ? X+N : NULL, collapse_masks,
2980          NULL, pulses, shortBlocks, spread_decision, dual_stereo, intensity, tf_res,
2981          len*(8<<BITRES)-anti_collapse_rsv, balance, dec, LM, codedBands, &st->rng);
2982
2983    if (anti_collapse_rsv > 0)
2984    {
2985       anti_collapse_on = ec_dec_bits(dec, 1);
2986    }
2987
2988    unquant_energy_finalise(mode, st->start, st->end, oldBandE,
2989          fine_quant, fine_priority, len*8-ec_tell(dec), dec, C);
2990
2991    if (anti_collapse_on)
2992       anti_collapse(mode, X, collapse_masks, LM, C, N,
2993             st->start, st->end, oldBandE, oldLogE, oldLogE2, pulses, st->rng);
2994
2995    log2Amp(mode, st->start, st->end, bandE, oldBandE, C);
2996
2997    if (silence)
2998    {
2999       for (i=0;i<C*nbEBands;i++)
3000       {
3001          bandE[i] = 0;
3002          oldBandE[i] = -QCONST16(28.f,DB_SHIFT);
3003       }
3004    }
3005    /* Synthesis */
3006    denormalise_bands(mode, X, freq, bandE, effEnd, C, M);
3007
3008    c=0; do {
3009       OPUS_MOVE(decode_mem[c], decode_mem[c]+N, DECODE_BUFFER_SIZE-N+overlap);
3010    } while (++c<CC);
3011
3012    c=0; do
3013       for (i=0;i<M*eBands[st->start];i++)
3014          freq[c*N+i] = 0;
3015    while (++c<C);
3016    c=0; do {
3017       int bound = M*eBands[effEnd];
3018       if (st->downsample!=1)
3019          bound = IMIN(bound, N/st->downsample);
3020       for (i=bound;i<N;i++)
3021          freq[c*N+i] = 0;
3022    } while (++c<C);
3023
3024    c=0; do {
3025       out_syn[c] = out_mem[c]+MAX_PERIOD-N;
3026    } while (++c<CC);
3027
3028    if (CC==2&&C==1)
3029    {
3030       for (i=0;i<N;i++)
3031          freq[N+i] = freq[i];
3032    }
3033    if (CC==1&&C==2)
3034    {
3035       for (i=0;i<N;i++)
3036          freq[i] = HALF32(ADD32(freq[i],freq[N+i]));
3037    }
3038
3039    /* Compute inverse MDCTs */
3040    compute_inv_mdcts(mode, shortBlocks, freq, out_syn, CC, LM);
3041
3042    c=0; do {
3043       st->postfilter_period=IMAX(st->postfilter_period, COMBFILTER_MINPERIOD);
3044       st->postfilter_period_old=IMAX(st->postfilter_period_old, COMBFILTER_MINPERIOD);
3045       comb_filter(out_syn[c], out_syn[c], st->postfilter_period_old, st->postfilter_period, mode->shortMdctSize,
3046             st->postfilter_gain_old, st->postfilter_gain, st->postfilter_tapset_old, st->postfilter_tapset,
3047             mode->window, overlap);
3048       if (LM!=0)
3049          comb_filter(out_syn[c]+mode->shortMdctSize, out_syn[c]+mode->shortMdctSize, st->postfilter_period, postfilter_pitch, N-mode->shortMdctSize,
3050                st->postfilter_gain, postfilter_gain, st->postfilter_tapset, postfilter_tapset,
3051                mode->window, overlap);
3052
3053    } while (++c<CC);
3054    st->postfilter_period_old = st->postfilter_period;
3055    st->postfilter_gain_old = st->postfilter_gain;
3056    st->postfilter_tapset_old = st->postfilter_tapset;
3057    st->postfilter_period = postfilter_pitch;
3058    st->postfilter_gain = postfilter_gain;
3059    st->postfilter_tapset = postfilter_tapset;
3060    if (LM!=0)
3061    {
3062       st->postfilter_period_old = st->postfilter_period;
3063       st->postfilter_gain_old = st->postfilter_gain;
3064       st->postfilter_tapset_old = st->postfilter_tapset;
3065    }
3066
3067    if (C==1) {
3068       for (i=0;i<nbEBands;i++)
3069          oldBandE[nbEBands+i]=oldBandE[i];
3070    }
3071
3072    /* In case start or end were to change */
3073    if (!isTransient)
3074    {
3075       for (i=0;i<2*nbEBands;i++)
3076          oldLogE2[i] = oldLogE[i];
3077       for (i=0;i<2*nbEBands;i++)
3078          oldLogE[i] = oldBandE[i];
3079       for (i=0;i<2*nbEBands;i++)
3080          backgroundLogE[i] = MIN16(backgroundLogE[i] + M*QCONST16(0.001f,DB_SHIFT), oldBandE[i]);
3081    } else {
3082       for (i=0;i<2*nbEBands;i++)
3083          oldLogE[i] = MIN16(oldLogE[i], oldBandE[i]);
3084    }
3085    c=0; do
3086    {
3087       for (i=0;i<st->start;i++)
3088       {
3089          oldBandE[c*nbEBands+i]=0;
3090          oldLogE[c*nbEBands+i]=oldLogE2[c*nbEBands+i]=-QCONST16(28.f,DB_SHIFT);
3091       }
3092       for (i=st->end;i<nbEBands;i++)
3093       {
3094          oldBandE[c*nbEBands+i]=0;
3095          oldLogE[c*nbEBands+i]=oldLogE2[c*nbEBands+i]=-QCONST16(28.f,DB_SHIFT);
3096       }
3097    } while (++c<2);
3098    st->rng = dec->rng;
3099
3100    /* We reuse freq[] as scratch space for the de-emphasis */
3101    deemphasis(out_syn, pcm, N, CC, st->downsample, mode->preemph, st->preemph_memD, freq);
3102    st->loss_count = 0;
3103    RESTORE_STACK;
3104    if (ec_tell(dec) > 8*len)
3105       return OPUS_INTERNAL_ERROR;
3106    if(ec_get_error(dec))
3107       st->error = 1;
3108    return frame_size/st->downsample;
3109 }
3110
3111
3112 #ifdef CUSTOM_MODES
3113
3114 #ifdef FIXED_POINT
3115 int opus_custom_decode(CELTDecoder * OPUS_RESTRICT st, const unsigned char *data, int len, opus_int16 * OPUS_RESTRICT pcm, int frame_size)
3116 {
3117    return celt_decode_with_ec(st, data, len, pcm, frame_size, NULL);
3118 }
3119
3120 #ifndef DISABLE_FLOAT_API
3121 int opus_custom_decode_float(CELTDecoder * OPUS_RESTRICT st, const unsigned char *data, int len, float * OPUS_RESTRICT pcm, int frame_size)
3122 {
3123    int j, ret, C, N;
3124    VARDECL(opus_int16, out);
3125    ALLOC_STACK;
3126
3127    if (pcm==NULL)
3128       return OPUS_BAD_ARG;
3129
3130    C = st->channels;
3131    N = frame_size;
3132
3133    ALLOC(out, C*N, opus_int16);
3134    ret=celt_decode_with_ec(st, data, len, out, frame_size, NULL);
3135    if (ret>0)
3136       for (j=0;j<C*ret;j++)
3137          pcm[j]=out[j]*(1.f/32768.f);
3138
3139    RESTORE_STACK;
3140    return ret;
3141 }
3142 #endif /* DISABLE_FLOAT_API */
3143
3144 #else
3145
3146 int opus_custom_decode_float(CELTDecoder * OPUS_RESTRICT st, const unsigned char *data, int len, float * OPUS_RESTRICT pcm, int frame_size)
3147 {
3148    return celt_decode_with_ec(st, data, len, pcm, frame_size, NULL);
3149 }
3150
3151 int opus_custom_decode(CELTDecoder * OPUS_RESTRICT st, const unsigned char *data, int len, opus_int16 * OPUS_RESTRICT pcm, int frame_size)
3152 {
3153    int j, ret, C, N;
3154    VARDECL(celt_sig, out);
3155    ALLOC_STACK;
3156
3157    if (pcm==NULL)
3158       return OPUS_BAD_ARG;
3159
3160    C = st->channels;
3161    N = frame_size;
3162    ALLOC(out, C*N, celt_sig);
3163
3164    ret=celt_decode_with_ec(st, data, len, out, frame_size, NULL);
3165
3166    if (ret>0)
3167       for (j=0;j<C*ret;j++)
3168          pcm[j] = FLOAT2INT16 (out[j]);
3169
3170    RESTORE_STACK;
3171    return ret;
3172 }
3173
3174 #endif
3175 #endif /* CUSTOM_MODES */
3176
3177 int opus_custom_decoder_ctl(CELTDecoder * OPUS_RESTRICT st, int request, ...)
3178 {
3179    va_list ap;
3180
3181    va_start(ap, request);
3182    switch (request)
3183    {
3184       case CELT_SET_START_BAND_REQUEST:
3185       {
3186          opus_int32 value = va_arg(ap, opus_int32);
3187          if (value<0 || value>=st->mode->nbEBands)
3188             goto bad_arg;
3189          st->start = value;
3190       }
3191       break;
3192       case CELT_SET_END_BAND_REQUEST:
3193       {
3194          opus_int32 value = va_arg(ap, opus_int32);
3195          if (value<1 || value>st->mode->nbEBands)
3196             goto bad_arg;
3197          st->end = value;
3198       }
3199       break;
3200       case CELT_SET_CHANNELS_REQUEST:
3201       {
3202          opus_int32 value = va_arg(ap, opus_int32);
3203          if (value<1 || value>2)
3204             goto bad_arg;
3205          st->stream_channels = value;
3206       }
3207       break;
3208       case CELT_GET_AND_CLEAR_ERROR_REQUEST:
3209       {
3210          opus_int32 *value = va_arg(ap, opus_int32*);
3211          if (value==NULL)
3212             goto bad_arg;
3213          *value=st->error;
3214          st->error = 0;
3215       }
3216       break;
3217       case OPUS_GET_LOOKAHEAD_REQUEST:
3218       {
3219          opus_int32 *value = va_arg(ap, opus_int32*);
3220          if (value==NULL)
3221             goto bad_arg;
3222          *value = st->overlap/st->downsample;
3223       }
3224       break;
3225       case OPUS_RESET_STATE:
3226       {
3227          int i;
3228          opus_val16 *lpc, *oldBandE, *oldLogE, *oldLogE2;
3229          lpc = (opus_val16*)(st->_decode_mem+(DECODE_BUFFER_SIZE+st->overlap)*st->channels);
3230          oldBandE = lpc+st->channels*LPC_ORDER;
3231          oldLogE = oldBandE + 2*st->mode->nbEBands;
3232          oldLogE2 = oldLogE + 2*st->mode->nbEBands;
3233          OPUS_CLEAR((char*)&st->DECODER_RESET_START,
3234                opus_custom_decoder_get_size(st->mode, st->channels)-
3235                ((char*)&st->DECODER_RESET_START - (char*)st));
3236          for (i=0;i<2*st->mode->nbEBands;i++)
3237             oldLogE[i]=oldLogE2[i]=-QCONST16(28.f,DB_SHIFT);
3238       }
3239       break;
3240       case OPUS_GET_PITCH_REQUEST:
3241       {
3242          opus_int32 *value = va_arg(ap, opus_int32*);
3243          if (value==NULL)
3244             goto bad_arg;
3245          *value = st->postfilter_period;
3246       }
3247       break;
3248 #ifdef OPUS_BUILD
3249       case CELT_GET_MODE_REQUEST:
3250       {
3251          const CELTMode ** value = va_arg(ap, const CELTMode**);
3252          if (value==0)
3253             goto bad_arg;
3254          *value=st->mode;
3255       }
3256       break;
3257       case CELT_SET_SIGNALLING_REQUEST:
3258       {
3259          opus_int32 value = va_arg(ap, opus_int32);
3260          st->signalling = value;
3261       }
3262       break;
3263       case OPUS_GET_FINAL_RANGE_REQUEST:
3264       {
3265          opus_uint32 * value = va_arg(ap, opus_uint32 *);
3266          if (value==0)
3267             goto bad_arg;
3268          *value=st->rng;
3269       }
3270       break;
3271 #endif
3272       default:
3273          goto bad_request;
3274    }
3275    va_end(ap);
3276    return OPUS_OK;
3277 bad_arg:
3278    va_end(ap);
3279    return OPUS_BAD_ARG;
3280 bad_request:
3281       va_end(ap);
3282   return OPUS_UNIMPLEMENTED;
3283 }
3284
3285
3286
3287 const char *opus_strerror(int error)
3288 {
3289    static const char * const error_strings[8] = {
3290       "success",
3291       "invalid argument",
3292       "buffer too small",
3293       "internal error",
3294       "corrupted stream",
3295       "request not implemented",
3296       "invalid state",
3297       "memory allocation failed"
3298    };
3299    if (error > 0 || error < -7)
3300       return "unknown error";
3301    else
3302       return error_strings[-error];
3303 }
3304
3305 const char *opus_get_version_string(void)
3306 {
3307     return "libopus " OPUS_VERSION
3308 #ifdef FIXED_POINT
3309           "-fixed"
3310 #endif
3311 #ifdef FUZZING
3312           "-fuzzing"
3313 #endif
3314           ;
3315 }