Adds code to revise the transient decision using the band energy
[opus.git] / celt / celt_encoder.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_ENCODER_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
54 /** Encoder state
55  @brief Encoder state
56  */
57 struct OpusCustomEncoder {
58    const OpusCustomMode *mode;     /**< Mode used by the encoder */
59    int overlap;
60    int channels;
61    int stream_channels;
62
63    int force_intra;
64    int clip;
65    int disable_pf;
66    int complexity;
67    int upsample;
68    int start, end;
69
70    opus_int32 bitrate;
71    int vbr;
72    int signalling;
73    int constrained_vbr;      /* If zero, VBR can do whatever it likes with the rate */
74    int loss_rate;
75    int lsb_depth;
76
77    /* Everything beyond this point gets cleared on a reset */
78 #define ENCODER_RESET_START rng
79
80    opus_uint32 rng;
81    int spread_decision;
82    opus_val32 delayedIntra;
83    int tonal_average;
84    int lastCodedBands;
85    int hf_average;
86    int tapset_decision;
87
88    int prefilter_period;
89    opus_val16 prefilter_gain;
90    int prefilter_tapset;
91 #ifdef RESYNTH
92    int prefilter_period_old;
93    opus_val16 prefilter_gain_old;
94    int prefilter_tapset_old;
95 #endif
96    int consec_transient;
97    AnalysisInfo analysis;
98
99    opus_val32 preemph_memE[2];
100    opus_val32 preemph_memD[2];
101
102    /* VBR-related parameters */
103    opus_int32 vbr_reservoir;
104    opus_int32 vbr_drift;
105    opus_int32 vbr_offset;
106    opus_int32 vbr_count;
107    opus_val16 overlap_max;
108    opus_val16 stereo_saving;
109    int intensity;
110
111 #ifdef RESYNTH
112    /* +MAX_PERIOD/2 to make space for overlap */
113    celt_sig syn_mem[2][2*MAX_PERIOD+MAX_PERIOD/2];
114 #endif
115
116    celt_sig in_mem[1]; /* Size = channels*mode->overlap */
117    /* celt_sig prefilter_mem[],  Size = channels*COMBFILTER_MAXPERIOD */
118    /* opus_val16 oldBandE[],     Size = channels*mode->nbEBands */
119    /* opus_val16 oldLogE[],      Size = channels*mode->nbEBands */
120    /* opus_val16 oldLogE2[],     Size = channels*mode->nbEBands */
121 };
122
123 int celt_encoder_get_size(int channels)
124 {
125    CELTMode *mode = opus_custom_mode_create(48000, 960, NULL);
126    return opus_custom_encoder_get_size(mode, channels);
127 }
128
129 OPUS_CUSTOM_NOSTATIC int opus_custom_encoder_get_size(const CELTMode *mode, int channels)
130 {
131    int size = sizeof(struct CELTEncoder)
132          + (channels*mode->overlap-1)*sizeof(celt_sig)    /* celt_sig in_mem[channels*mode->overlap]; */
133          + channels*COMBFILTER_MAXPERIOD*sizeof(celt_sig) /* celt_sig prefilter_mem[channels*COMBFILTER_MAXPERIOD]; */
134          + 3*channels*mode->nbEBands*sizeof(opus_val16);  /* opus_val16 oldBandE[channels*mode->nbEBands]; */
135                                                           /* opus_val16 oldLogE[channels*mode->nbEBands]; */
136                                                           /* opus_val16 oldLogE2[channels*mode->nbEBands]; */
137    return size;
138 }
139
140 #ifdef CUSTOM_MODES
141 CELTEncoder *opus_custom_encoder_create(const CELTMode *mode, int channels, int *error)
142 {
143    int ret;
144    CELTEncoder *st = (CELTEncoder *)opus_alloc(opus_custom_encoder_get_size(mode, channels));
145    /* init will handle the NULL case */
146    ret = opus_custom_encoder_init(st, mode, channels);
147    if (ret != OPUS_OK)
148    {
149       opus_custom_encoder_destroy(st);
150       st = NULL;
151    }
152    if (error)
153       *error = ret;
154    return st;
155 }
156 #endif /* CUSTOM_MODES */
157
158 int celt_encoder_init(CELTEncoder *st, opus_int32 sampling_rate, int channels)
159 {
160    int ret;
161    ret = opus_custom_encoder_init(st, opus_custom_mode_create(48000, 960, NULL), channels);
162    if (ret != OPUS_OK)
163       return ret;
164    st->upsample = resampling_factor(sampling_rate);
165    return OPUS_OK;
166 }
167
168 OPUS_CUSTOM_NOSTATIC int opus_custom_encoder_init(CELTEncoder *st, const CELTMode *mode, int channels)
169 {
170    if (channels < 0 || channels > 2)
171       return OPUS_BAD_ARG;
172
173    if (st==NULL || mode==NULL)
174       return OPUS_ALLOC_FAIL;
175
176    OPUS_CLEAR((char*)st, opus_custom_encoder_get_size(mode, channels));
177
178    st->mode = mode;
179    st->overlap = mode->overlap;
180    st->stream_channels = st->channels = channels;
181
182    st->upsample = 1;
183    st->start = 0;
184    st->end = st->mode->effEBands;
185    st->signalling = 1;
186
187    st->constrained_vbr = 1;
188    st->clip = 1;
189
190    st->bitrate = OPUS_BITRATE_MAX;
191    st->vbr = 0;
192    st->force_intra  = 0;
193    st->complexity = 5;
194    st->lsb_depth=24;
195
196    opus_custom_encoder_ctl(st, OPUS_RESET_STATE);
197
198    return OPUS_OK;
199 }
200
201 #ifdef CUSTOM_MODES
202 void opus_custom_encoder_destroy(CELTEncoder *st)
203 {
204    opus_free(st);
205 }
206 #endif /* CUSTOM_MODES */
207
208
209 static int transient_analysis(const opus_val32 * OPUS_RESTRICT in, int len, int C,
210                               opus_val16 *tf_estimate, int *tf_chan)
211 {
212    int i;
213    VARDECL(opus_val16, tmp);
214    opus_val32 mem0,mem1;
215    int is_transient = 0;
216    opus_int32 mask_metric = 0;
217    int c;
218    opus_val16 tf_max;
219    int len2;
220    /* Table of 6*64/x, trained on real data to minimize the average error */
221    static const unsigned char inv_table[128] = {
222          255,255,156,110, 86, 70, 59, 51, 45, 40, 37, 33, 31, 28, 26, 25,
223           23, 22, 21, 20, 19, 18, 17, 16, 16, 15, 15, 14, 13, 13, 12, 12,
224           12, 12, 11, 11, 11, 10, 10, 10,  9,  9,  9,  9,  9,  9,  8,  8,
225            8,  8,  8,  7,  7,  7,  7,  7,  7,  6,  6,  6,  6,  6,  6,  6,
226            6,  6,  6,  6,  6,  6,  6,  6,  6,  5,  5,  5,  5,  5,  5,  5,
227            5,  5,  5,  5,  5,  4,  4,  4,  4,  4,  4,  4,  4,  4,  4,  4,
228            4,  4,  4,  4,  4,  4,  4,  4,  4,  4,  4,  4,  4,  4,  3,  3,
229            3,  3,  3,  3,  3,  3,  3,  3,  3,  3,  3,  3,  3,  3,  3,  2,
230    };
231    SAVE_STACK;
232    ALLOC(tmp, len, opus_val16);
233
234    len2=len/2;
235    tf_max = 0;
236    for (c=0;c<C;c++)
237    {
238       opus_val32 mean;
239       opus_int32 unmask=0;
240       opus_val32 norm;
241       opus_val16 maxE;
242       mem0=0;
243       mem1=0;
244       /* High-pass filter: (1 - 2*z^-1 + z^-2) / (1 - z^-1 + .5*z^-2) */
245       for (i=0;i<len;i++)
246       {
247          opus_val32 x,y;
248          x = SHR32(in[i+c*len],SIG_SHIFT);
249          y = ADD32(mem0, x);
250 #ifdef FIXED_POINT
251          mem0 = mem1 + y - SHL32(x,1);
252          mem1 = x - SHR32(y,1);
253 #else
254          mem0 = mem1 + y - 2*x;
255          mem1 = x - .5f*y;
256 #endif
257          tmp[i] = EXTRACT16(SHR32(y,2));
258          /*printf("%f ", tmp[i]);*/
259       }
260       /*printf("\n");*/
261       /* First few samples are bad because we don't propagate the memory */
262       for (i=0;i<12;i++)
263          tmp[i] = 0;
264
265 #ifdef FIXED_POINT
266       /* Normalize tmp to max range */
267       {
268          int shift=0;
269          shift = 14-celt_ilog2(1+celt_maxabs16(tmp, len));
270          if (shift!=0)
271          {
272             for (i=0;i<len;i++)
273                tmp[i] = SHL16(tmp[i], shift);
274          }
275       }
276 #endif
277
278       mean=0;
279       mem0=0;
280       /* Grouping by two to reduce complexity */
281       /* Forward pass to compute the post-echo threshold*/
282       for (i=0;i<len2;i++)
283       {
284          opus_val16 x2 = PSHR32(MULT16_16(tmp[2*i],tmp[2*i]) + MULT16_16(tmp[2*i+1],tmp[2*i+1]),16);
285          mean += x2;
286 #ifdef FIXED_POINT
287          /* FIXME: Use PSHR16() instead */
288          tmp[i] = mem0 + PSHR32(x2-mem0,4);
289 #else
290          tmp[i] = mem0 + MULT16_16_P15(QCONST16(.0625f,15),x2-mem0);
291 #endif
292          mem0 = tmp[i];
293       }
294
295       mem0=0;
296       maxE=0;
297       /* Backward pass to compute the pre-echo threshold */
298       for (i=len2-1;i>=0;i--)
299       {
300 #ifdef FIXED_POINT
301          /* FIXME: Use PSHR16() instead */
302          tmp[i] = mem0 + PSHR32(tmp[i]-mem0,3);
303 #else
304          tmp[i] = mem0 + MULT16_16_P15(QCONST16(0.125f,15),tmp[i]-mem0);
305 #endif
306          mem0 = tmp[i];
307          maxE = MAX16(maxE, mem0);
308       }
309       /*for (i=0;i<len2;i++)printf("%f ", tmp[i]/mean);printf("\n");*/
310
311       /* Compute the ratio of the "frame energy" over the harmonic mean of the energy.
312          This essentially corresponds to a bitrate-normalized temporal noise-to-mask
313          ratio */
314
315       /* As a compromise with the old transient detector, frame energy is the
316          geometric mean of the energy and half the max */
317 #ifdef FIXED_POINT
318       /* Costs two sqrt() to avoid overflows */
319       mean = MULT16_16(celt_sqrt(mean), celt_sqrt(MULT16_16(maxE,len2>>1)));
320 #else
321       mean = celt_sqrt(mean * maxE*.5*len2);
322 #endif
323       /* Inverse of the mean energy in Q15+6 */
324       norm = SHL32(EXTEND32(len2),6+14)/ADD32(EPSILON,SHR32(mean,1));
325       /* Compute harmonic mean discarding the unreliable boundaries
326          The data is smooth, so we only take 1/4th of the samples */
327       unmask=0;
328       for (i=12;i<len2-5;i+=4)
329       {
330          int id;
331 #ifdef FIXED_POINT
332          id = IMAX(0,IMIN(127,MULT16_32_Q15(tmp[i],norm))); /* Do not round to nearest */
333 #else
334          id = IMAX(0,IMIN(127,(int)floor(64*norm*tmp[i]))); /* Do not round to nearest */
335 #endif
336          unmask += inv_table[id];
337       }
338       /*printf("%d\n", unmask);*/
339       /* Normalize, compensate for the 1/4th of the sample and the factor of 6 in the inverse table */
340       unmask = 64*unmask*4/(6*(len2-17));
341       if (unmask>mask_metric)
342       {
343          *tf_chan = c;
344          mask_metric = unmask;
345       }
346    }
347    is_transient = mask_metric>200;
348
349    /* Arbitrary metric for VBR boost */
350    tf_max = MAX16(0,celt_sqrt(27*mask_metric)-42);
351    /* *tf_estimate = 1 + MIN16(1, sqrt(MAX16(0, tf_max-30))/20); */
352    *tf_estimate = celt_sqrt(MAX16(0, SHL32(MULT16_16(QCONST16(0.0069,14),MIN16(163,tf_max)),14)-QCONST32(0.139,28)));
353    /*printf("%d %f\n", tf_max, mask_metric);*/
354    RESTORE_STACK;
355 #ifdef FUZZING
356    is_transient = rand()&0x1;
357 #endif
358    /*printf("%d %f %d\n", is_transient, (float)*tf_estimate, tf_max);*/
359    return is_transient;
360 }
361
362 /* Looks for sudden increases of energy to decide whether we need to patch
363    the transient decision */
364 int patch_transient_decision(opus_val16 *new, opus_val16 *old, int nbEBands,
365       int end, int C)
366 {
367    int i, c;
368    opus_val32 mean_diff=0;
369    opus_val16 spread_old[26];
370    /* Apply an aggressive (-6 dB/Bark) spreading function to the old frame to
371       avoid false detection caused by irrelevant bands */
372    if (C==1)
373    {
374       spread_old[0] = old[0];
375       for (i=1;i<end;i++)
376          spread_old[i] = MAX16(spread_old[i-1]-QCONST16(1.0f, DB_SHIFT), old[i]);
377    } else {
378       spread_old[0] = MAX16(old[0],old[nbEBands]);
379       for (i=1;i<end;i++)
380          spread_old[i] = MAX16(spread_old[i-1]-QCONST16(1.0f, DB_SHIFT),
381                                MAX16(old[i],old[i+nbEBands]));
382    }
383    for (i=end-2;i>=0;i--)
384       spread_old[i] = MAX16(spread_old[i], spread_old[i+1]-QCONST16(1.0f, DB_SHIFT));
385    /* Compute mean increase */
386    c=0; do {
387       for (i=2;i<end-1;i++)
388       {
389          opus_val16 x1, x2;
390          x1 = MAX16(0, new[i]);
391          x2 = MAX16(0, spread_old[i]);
392          mean_diff = ADD32(mean_diff, EXTEND32(MAX16(0, SUB16(x1, x2))));
393       }
394    } while (++c<C);
395    mean_diff = DIV32(mean_diff, C*(end-3));
396    /*printf("%f %f %d\n", mean_diff, max_diff, count);*/
397    return mean_diff > QCONST16(1.f, DB_SHIFT);
398 }
399
400 /** Apply window and compute the MDCT for all sub-frames and
401     all channels in a frame */
402 static void compute_mdcts(const CELTMode *mode, int shortBlocks, celt_sig * OPUS_RESTRICT in, celt_sig * OPUS_RESTRICT out, int C, int LM)
403 {
404    const int overlap = OVERLAP(mode);
405    int N;
406    int B;
407    int shift;
408    int b, c;
409    if (shortBlocks)
410    {
411       B = shortBlocks;
412       N = mode->shortMdctSize;
413       shift = mode->maxLM;
414    } else {
415       B = 1;
416       N = mode->shortMdctSize<<LM;
417       shift = mode->maxLM-LM;
418    }
419    c=0; do {
420       for (b=0;b<B;b++)
421       {
422          /* Interleaving the sub-frames while doing the MDCTs */
423          clt_mdct_forward(&mode->mdct, in+c*(B*N+overlap)+b*N, &out[b+c*N*B], mode->window, overlap, shift, B);
424       }
425    } while (++c<C);
426 }
427
428
429 static void preemphasis(const opus_val16 * OPUS_RESTRICT pcmp, celt_sig * OPUS_RESTRICT inp,
430                         int N, int CC, int upsample, const opus_val16 *coef, celt_sig *mem, int clip)
431 {
432    int i;
433    opus_val16 coef0;
434    celt_sig m;
435    int Nu;
436
437    coef0 = coef[0];
438
439
440    Nu = N/upsample;
441    if (upsample!=1)
442    {
443       for (i=0;i<N;i++)
444          inp[i] = 0;
445    }
446    for (i=0;i<Nu;i++)
447    {
448       celt_sig x;
449
450       x = SCALEIN(pcmp[CC*i]);
451 #ifndef FIXED_POINT
452       /* Replace NaNs with zeros */
453       if (!(x==x))
454          x = 0;
455 #endif
456       inp[i*upsample] = x;
457    }
458
459 #ifndef FIXED_POINT
460    if (clip)
461    {
462       /* Clip input to avoid encoding non-portable files */
463       for (i=0;i<Nu;i++)
464          inp[i*upsample] = MAX32(-65536.f, MIN32(65536.f,inp[i*upsample]));
465    }
466 #endif
467    m = *mem;
468 #ifdef CUSTOM_MODES
469    if (coef[1] != 0)
470    {
471       opus_val16 coef1 = coef[1];
472       opus_val16 coef2 = coef[2];
473       for (i=0;i<N;i++)
474       {
475          opus_val16 x, tmp;
476          x = inp[i];
477          /* Apply pre-emphasis */
478          tmp = MULT16_16(coef2, x);
479          inp[i] = tmp + m;
480          m = MULT16_32_Q15(coef1, inp[i]) - MULT16_32_Q15(coef0, tmp);
481       }
482    } else
483 #endif
484    {
485       for (i=0;i<N;i++)
486       {
487          celt_sig x;
488          x = SHL32(inp[i], SIG_SHIFT);
489          /* Apply pre-emphasis */
490          inp[i] = x + m;
491          m = - MULT16_32_Q15(coef0, x);
492       }
493    }
494    *mem = m;
495 }
496
497
498
499 static opus_val32 l1_metric(const celt_norm *tmp, int N, int LM, opus_val16 bias)
500 {
501    int i;
502    opus_val32 L1;
503    L1 = 0;
504    for (i=0;i<N;i++)
505       L1 += EXTEND32(ABS16(tmp[i]));
506    /* When in doubt, prefer good freq resolution */
507    L1 = MAC16_32_Q15(L1, LM*bias, L1);
508    return L1;
509
510 }
511
512 static int tf_analysis(const CELTMode *m, int len, int isTransient,
513       int *tf_res, int lambda, celt_norm *X, int N0, int LM,
514       int *tf_sum, opus_val16 tf_estimate, int tf_chan)
515 {
516    int i;
517    VARDECL(int, metric);
518    int cost0;
519    int cost1;
520    VARDECL(int, path0);
521    VARDECL(int, path1);
522    VARDECL(celt_norm, tmp);
523    VARDECL(celt_norm, tmp_1);
524    int sel;
525    int selcost[2];
526    int tf_select=0;
527    opus_val16 bias;
528
529    SAVE_STACK;
530    bias = MULT16_16_Q14(QCONST16(.04f,15), MAX16(-QCONST16(.25f,14), QCONST16(.5f,14)-tf_estimate));
531    /*printf("%f ", bias);*/
532
533    ALLOC(metric, len, int);
534    ALLOC(tmp, (m->eBands[len]-m->eBands[len-1])<<LM, celt_norm);
535    ALLOC(tmp_1, (m->eBands[len]-m->eBands[len-1])<<LM, celt_norm);
536    ALLOC(path0, len, int);
537    ALLOC(path1, len, int);
538
539    *tf_sum = 0;
540    for (i=0;i<len;i++)
541    {
542       int j, k, N;
543       int narrow;
544       opus_val32 L1, best_L1;
545       int best_level=0;
546       N = (m->eBands[i+1]-m->eBands[i])<<LM;
547       /* band is too narrow to be split down to LM=-1 */
548       narrow = (m->eBands[i+1]-m->eBands[i])==1;
549       for (j=0;j<N;j++)
550          tmp[j] = X[tf_chan*N0 + j+(m->eBands[i]<<LM)];
551       /* Just add the right channel if we're in stereo */
552       /*if (C==2)
553          for (j=0;j<N;j++)
554             tmp[j] = ADD16(SHR16(tmp[j], 1),SHR16(X[N0+j+(m->eBands[i]<<LM)], 1));*/
555       L1 = l1_metric(tmp, N, isTransient ? LM : 0, bias);
556       best_L1 = L1;
557       /* Check the -1 case for transients */
558       if (isTransient && !narrow)
559       {
560          for (j=0;j<N;j++)
561             tmp_1[j] = tmp[j];
562          haar1(tmp_1, N>>LM, 1<<LM);
563          L1 = l1_metric(tmp_1, N, LM+1, bias);
564          if (L1<best_L1)
565          {
566             best_L1 = L1;
567             best_level = -1;
568          }
569       }
570       /*printf ("%f ", L1);*/
571       for (k=0;k<LM+!(isTransient||narrow);k++)
572       {
573          int B;
574
575          if (isTransient)
576             B = (LM-k-1);
577          else
578             B = k+1;
579
580          haar1(tmp, N>>k, 1<<k);
581
582          L1 = l1_metric(tmp, N, B, bias);
583
584          if (L1 < best_L1)
585          {
586             best_L1 = L1;
587             best_level = k+1;
588          }
589       }
590       /*printf ("%d ", isTransient ? LM-best_level : best_level);*/
591       /* metric is in Q1 to be able to select the mid-point (-0.5) for narrower bands */
592       if (isTransient)
593          metric[i] = 2*best_level;
594       else
595          metric[i] = -2*best_level;
596       *tf_sum += (isTransient ? LM : 0) - metric[i]/2;
597       /* For bands that can't be split to -1, set the metric to the half-way point to avoid
598          biasing the decision */
599       if (narrow && (metric[i]==0 || metric[i]==-2*LM))
600          metric[i]-=1;
601       /*printf("%d ", metric[i]);*/
602    }
603    /*printf("\n");*/
604    /* Search for the optimal tf resolution, including tf_select */
605    tf_select = 0;
606    for (sel=0;sel<2;sel++)
607    {
608       cost0 = 0;
609       cost1 = isTransient ? 0 : lambda;
610       for (i=1;i<len;i++)
611       {
612          int curr0, curr1;
613          curr0 = IMIN(cost0, cost1 + lambda);
614          curr1 = IMIN(cost0 + lambda, cost1);
615          cost0 = curr0 + abs(metric[i]-2*tf_select_table[LM][4*isTransient+2*sel+0]);
616          cost1 = curr1 + abs(metric[i]-2*tf_select_table[LM][4*isTransient+2*sel+1]);
617       }
618       cost0 = IMIN(cost0, cost1);
619       selcost[sel]=cost0;
620    }
621    /* For now, we're conservative and only allow tf_select=1 for transients.
622     * If tests confirm it's useful for non-transients, we could allow it. */
623    if (selcost[1]<selcost[0] && isTransient)
624       tf_select=1;
625    cost0 = 0;
626    cost1 = isTransient ? 0 : lambda;
627    /* Viterbi forward pass */
628    for (i=1;i<len;i++)
629    {
630       int curr0, curr1;
631       int from0, from1;
632
633       from0 = cost0;
634       from1 = cost1 + lambda;
635       if (from0 < from1)
636       {
637          curr0 = from0;
638          path0[i]= 0;
639       } else {
640          curr0 = from1;
641          path0[i]= 1;
642       }
643
644       from0 = cost0 + lambda;
645       from1 = cost1;
646       if (from0 < from1)
647       {
648          curr1 = from0;
649          path1[i]= 0;
650       } else {
651          curr1 = from1;
652          path1[i]= 1;
653       }
654       cost0 = curr0 + abs(metric[i]-2*tf_select_table[LM][4*isTransient+2*tf_select+0]);
655       cost1 = curr1 + abs(metric[i]-2*tf_select_table[LM][4*isTransient+2*tf_select+1]);
656    }
657    tf_res[len-1] = cost0 < cost1 ? 0 : 1;
658    /* Viterbi backward pass to check the decisions */
659    for (i=len-2;i>=0;i--)
660    {
661       if (tf_res[i+1] == 1)
662          tf_res[i] = path1[i+1];
663       else
664          tf_res[i] = path0[i+1];
665    }
666    /*printf("%d %f\n", *tf_sum, tf_estimate);*/
667    RESTORE_STACK;
668 #ifdef FUZZING
669    tf_select = rand()&0x1;
670    tf_res[0] = rand()&0x1;
671    for (i=1;i<len;i++)
672       tf_res[i] = tf_res[i-1] ^ ((rand()&0xF) == 0);
673 #endif
674    return tf_select;
675 }
676
677 static void tf_encode(int start, int end, int isTransient, int *tf_res, int LM, int tf_select, ec_enc *enc)
678 {
679    int curr, i;
680    int tf_select_rsv;
681    int tf_changed;
682    int logp;
683    opus_uint32 budget;
684    opus_uint32 tell;
685    budget = enc->storage*8;
686    tell = ec_tell(enc);
687    logp = isTransient ? 2 : 4;
688    /* Reserve space to code the tf_select decision. */
689    tf_select_rsv = LM>0 && tell+logp+1 <= budget;
690    budget -= tf_select_rsv;
691    curr = tf_changed = 0;
692    for (i=start;i<end;i++)
693    {
694       if (tell+logp<=budget)
695       {
696          ec_enc_bit_logp(enc, tf_res[i] ^ curr, logp);
697          tell = ec_tell(enc);
698          curr = tf_res[i];
699          tf_changed |= curr;
700       }
701       else
702          tf_res[i] = curr;
703       logp = isTransient ? 4 : 5;
704    }
705    /* Only code tf_select if it would actually make a difference. */
706    if (tf_select_rsv &&
707          tf_select_table[LM][4*isTransient+0+tf_changed]!=
708          tf_select_table[LM][4*isTransient+2+tf_changed])
709       ec_enc_bit_logp(enc, tf_select, 1);
710    else
711       tf_select = 0;
712    for (i=start;i<end;i++)
713       tf_res[i] = tf_select_table[LM][4*isTransient+2*tf_select+tf_res[i]];
714    /*for(i=0;i<end;i++)printf("%d ", isTransient ? tf_res[i] : LM+tf_res[i]);printf("\n");*/
715 }
716
717
718 static int alloc_trim_analysis(const CELTMode *m, const celt_norm *X,
719       const opus_val16 *bandLogE, int end, int LM, int C, int N0,
720       AnalysisInfo *analysis, opus_val16 *stereo_saving, opus_val16 tf_estimate,
721       int intensity)
722 {
723    int i;
724    opus_val32 diff=0;
725    int c;
726    int trim_index = 5;
727    opus_val16 trim = QCONST16(5.f, 8);
728    opus_val16 logXC, logXC2;
729    if (C==2)
730    {
731       opus_val16 sum = 0; /* Q10 */
732       opus_val16 minXC; /* Q10 */
733       /* Compute inter-channel correlation for low frequencies */
734       for (i=0;i<8;i++)
735       {
736          int j;
737          opus_val32 partial = 0;
738          for (j=m->eBands[i]<<LM;j<m->eBands[i+1]<<LM;j++)
739             partial = MAC16_16(partial, X[j], X[N0+j]);
740          sum = ADD16(sum, EXTRACT16(SHR32(partial, 18)));
741       }
742       sum = MULT16_16_Q15(QCONST16(1.f/8, 15), sum);
743       sum = MIN16(QCONST16(1.f, 10), ABS16(sum));
744       minXC = sum;
745       for (i=8;i<intensity;i++)
746       {
747          int j;
748          opus_val32 partial = 0;
749          for (j=m->eBands[i]<<LM;j<m->eBands[i+1]<<LM;j++)
750             partial = MAC16_16(partial, X[j], X[N0+j]);
751          minXC = MIN16(minXC, ABS16(EXTRACT16(SHR32(partial, 18))));
752       }
753       minXC = MIN16(QCONST16(1.f, 10), ABS16(minXC));
754       /*printf ("%f\n", sum);*/
755       if (sum > QCONST16(.995f,10))
756          trim_index-=4;
757       else if (sum > QCONST16(.92f,10))
758          trim_index-=3;
759       else if (sum > QCONST16(.85f,10))
760          trim_index-=2;
761       else if (sum > QCONST16(.8f,10))
762          trim_index-=1;
763       /* mid-side savings estimations based on the LF average*/
764       logXC = celt_log2(QCONST32(1.001f, 20)-MULT16_16(sum, sum));
765       /* mid-side savings estimations based on min correlation */
766       logXC2 = MAX16(HALF16(logXC), celt_log2(QCONST32(1.001f, 20)-MULT16_16(minXC, minXC)));
767 #ifdef FIXED_POINT
768       /* Compensate for Q20 vs Q14 input and convert output to Q8 */
769       logXC = PSHR32(logXC-QCONST16(6.f, DB_SHIFT),DB_SHIFT-8);
770       logXC2 = PSHR32(logXC2-QCONST16(6.f, DB_SHIFT),DB_SHIFT-8);
771 #endif
772
773       trim += MAX16(-QCONST16(4.f, 8), MULT16_16_Q15(QCONST16(.75f,15),logXC));
774       *stereo_saving = MIN16(*stereo_saving + QCONST16(0.25f, 8), -HALF16(logXC2));
775    }
776
777    /* Estimate spectral tilt */
778    c=0; do {
779       for (i=0;i<end-1;i++)
780       {
781          diff += bandLogE[i+c*m->nbEBands]*(opus_int32)(2+2*i-end);
782       }
783    } while (++c<C);
784    diff /= C*(end-1);
785    /*printf("%f\n", diff);*/
786    if (diff > QCONST16(2.f, DB_SHIFT))
787       trim_index--;
788    if (diff > QCONST16(8.f, DB_SHIFT))
789       trim_index--;
790    if (diff < -QCONST16(4.f, DB_SHIFT))
791       trim_index++;
792    if (diff < -QCONST16(10.f, DB_SHIFT))
793       trim_index++;
794    trim -= MAX16(-QCONST16(2.f, 8), MIN16(QCONST16(2.f, 8), SHR16(diff+QCONST16(1.f, DB_SHIFT),DB_SHIFT-8)/6 ));
795    trim -= 2*SHR16(tf_estimate, 14-8);
796 #ifndef FIXED_POINT
797    if (analysis->valid)
798    {
799       trim -= MAX16(-QCONST16(2.f, 8), MIN16(QCONST16(2.f, 8), 2*(analysis->tonality_slope+.05f)));
800    }
801 #endif
802
803 #ifdef FIXED_POINT
804    trim_index = PSHR32(trim, 8);
805 #else
806    trim_index = (int)floor(.5f+trim);
807 #endif
808    if (trim_index<0)
809       trim_index = 0;
810    if (trim_index>10)
811       trim_index = 10;
812    /*printf("%d\n", trim_index);*/
813 #ifdef FUZZING
814    trim_index = rand()%11;
815 #endif
816    return trim_index;
817 }
818
819 static int stereo_analysis(const CELTMode *m, const celt_norm *X,
820       int LM, int N0)
821 {
822    int i;
823    int thetas;
824    opus_val32 sumLR = EPSILON, sumMS = EPSILON;
825
826    /* Use the L1 norm to model the entropy of the L/R signal vs the M/S signal */
827    for (i=0;i<13;i++)
828    {
829       int j;
830       for (j=m->eBands[i]<<LM;j<m->eBands[i+1]<<LM;j++)
831       {
832          opus_val32 L, R, M, S;
833          /* We cast to 32-bit first because of the -32768 case */
834          L = EXTEND32(X[j]);
835          R = EXTEND32(X[N0+j]);
836          M = ADD32(L, R);
837          S = SUB32(L, R);
838          sumLR = ADD32(sumLR, ADD32(ABS32(L), ABS32(R)));
839          sumMS = ADD32(sumMS, ADD32(ABS32(M), ABS32(S)));
840       }
841    }
842    sumMS = MULT16_32_Q15(QCONST16(0.707107f, 15), sumMS);
843    thetas = 13;
844    /* We don't need thetas for lower bands with LM<=1 */
845    if (LM<=1)
846       thetas -= 8;
847    return MULT16_32_Q15((m->eBands[13]<<(LM+1))+thetas, sumMS)
848          > MULT16_32_Q15(m->eBands[13]<<(LM+1), sumLR);
849 }
850
851 static opus_val16 dynalloc_analysis(const opus_val16 *bandLogE, const opus_val16 *bandLogE2,
852       int nbEBands, int start, int end, int C, int *offsets, int lsb_depth, const opus_int16 *logN,
853       int isTransient, int vbr, int constrained_vbr, const opus_int16 *eBands, int LM,
854       int effectiveBytes, opus_int32 *tot_boost_)
855 {
856    int i, c;
857    opus_int32 tot_boost=0;
858    opus_val16 maxDepth;
859    VARDECL(opus_val16, follower);
860    VARDECL(opus_val16, noise_floor);
861    SAVE_STACK;
862    ALLOC(follower, C*nbEBands, opus_val16);
863    ALLOC(noise_floor, C*nbEBands, opus_val16);
864    for (i=0;i<nbEBands;i++)
865       offsets[i] = 0;
866    /* Dynamic allocation code */
867    maxDepth=-QCONST16(32.f, DB_SHIFT);
868    for (i=0;i<end;i++)
869    {
870       /* Noise floor must take into account eMeans, the depth, the width of the bands
871          and the preemphasis filter (approx. square of bark band ID) */
872       noise_floor[i] = MULT16_16(QCONST16(0.0625f, DB_SHIFT),logN[i])
873             +QCONST16(.5f,DB_SHIFT)+SHL16(9-lsb_depth,DB_SHIFT)-SHL16(eMeans[i],6)
874             +MULT16_16(QCONST16(.0062,DB_SHIFT),(i+5)*(i+5));
875    }
876    c=0;do
877    {
878       for (i=0;i<end;i++)
879          maxDepth = MAX16(maxDepth, bandLogE[c*nbEBands+i]-noise_floor[i]);
880    } while (++c<C);
881    /* Make sure that dynamic allocation can't make us bust the budget */
882    if (effectiveBytes > 50 && LM>=1)
883    {
884       int last=0;
885       c=0;do
886       {
887          follower[c*nbEBands] = bandLogE2[c*nbEBands];
888          for (i=1;i<end;i++)
889          {
890             /* The last band to be at least 3 dB higher than the previous one
891                is the last we'll consider. Otherwise, we run into problems on
892                bandlimited signals. */
893             if (bandLogE2[c*nbEBands+i] > bandLogE2[c*nbEBands+i-1]+QCONST16(.5f,DB_SHIFT))
894                last=i;
895             follower[c*nbEBands+i] = MIN16(follower[c*nbEBands+i-1]+QCONST16(1.5f,DB_SHIFT), bandLogE2[c*nbEBands+i]);
896          }
897          for (i=last-1;i>=0;i--)
898             follower[c*nbEBands+i] = MIN16(follower[c*nbEBands+i], MIN16(follower[c*nbEBands+i+1]+QCONST16(2.f,DB_SHIFT), bandLogE2[c*nbEBands+i]));
899          for (i=0;i<end;i++)
900             follower[c*nbEBands+i] = MAX16(follower[c*nbEBands+i], noise_floor[i]);
901       } while (++c<C);
902       if (C==2)
903       {
904          for (i=start;i<end;i++)
905          {
906             /* Consider 24 dB "cross-talk" */
907             follower[nbEBands+i] = MAX16(follower[nbEBands+i], follower[         i]-QCONST16(4.f,DB_SHIFT));
908             follower[         i] = MAX16(follower[         i], follower[nbEBands+i]-QCONST16(4.f,DB_SHIFT));
909             follower[i] = HALF16(MAX16(0, bandLogE[i]-follower[i]) + MAX16(0, bandLogE[nbEBands+i]-follower[nbEBands+i]));
910          }
911       } else {
912          for (i=start;i<end;i++)
913          {
914             follower[i] = MAX16(0, bandLogE[i]-follower[i]);
915          }
916       }
917       /* For non-transient CBR/CVBR frames, halve the dynalloc contribution */
918       if ((!vbr || constrained_vbr)&&!isTransient)
919       {
920          for (i=start;i<end;i++)
921             follower[i] = HALF16(follower[i]);
922       }
923       for (i=start;i<end;i++)
924       {
925          int width;
926          int boost;
927          int boost_bits;
928
929          if (i<8)
930             follower[i] *= 2;
931          if (i>=12)
932             follower[i] = HALF16(follower[i]);
933          follower[i] = MIN16(follower[i], QCONST16(4, DB_SHIFT));
934
935          width = C*(eBands[i+1]-eBands[i])<<LM;
936          if (width<6)
937          {
938             boost = (int)SHR32(EXTEND32(follower[i]),DB_SHIFT);
939             boost_bits = boost*width<<BITRES;
940          } else if (width > 48) {
941             boost = (int)SHR32(EXTEND32(follower[i])*8,DB_SHIFT);
942             boost_bits = (boost*width<<BITRES)/8;
943          } else {
944             boost = (int)SHR32(EXTEND32(follower[i])*width/6,DB_SHIFT);
945             boost_bits = boost*6<<BITRES;
946          }
947          /* For CBR and non-transient CVBR frames, limit dynalloc to 1/4 of the bits */
948          if ((!vbr || (constrained_vbr&&!isTransient))
949                && (tot_boost+boost_bits)>>BITRES>>3 > effectiveBytes/4)
950          {
951             offsets[i] = 0;
952             break;
953          } else {
954             offsets[i] = boost;
955             tot_boost += boost_bits;
956          }
957       }
958    }
959    *tot_boost_ = tot_boost;
960    RESTORE_STACK;
961    return maxDepth;
962 }
963
964
965 static int run_prefilter(CELTEncoder *st, celt_sig *in, celt_sig *prefilter_mem, int CC, int N,
966       int prefilter_tapset, int *pitch, opus_val16 *gain, int *qgain, int enabled, int nbAvailableBytes)
967 {
968    int c;
969    VARDECL(celt_sig, _pre);
970    celt_sig *pre[2];
971    const CELTMode *mode;
972    int pitch_index;
973    opus_val16 gain1;
974    opus_val16 pf_threshold;
975    int pf_on;
976    int qg;
977    SAVE_STACK;
978
979    mode = st->mode;
980    ALLOC(_pre, CC*(N+COMBFILTER_MAXPERIOD), celt_sig);
981
982    pre[0] = _pre;
983    pre[1] = _pre + (N+COMBFILTER_MAXPERIOD);
984
985
986    c=0; do {
987       OPUS_COPY(pre[c], prefilter_mem+c*COMBFILTER_MAXPERIOD, COMBFILTER_MAXPERIOD);
988       OPUS_COPY(pre[c]+COMBFILTER_MAXPERIOD, in+c*(N+st->overlap)+st->overlap, N);
989    } while (++c<CC);
990
991    if (enabled)
992    {
993       VARDECL(opus_val16, pitch_buf);
994       ALLOC(pitch_buf, (COMBFILTER_MAXPERIOD+N)>>1, opus_val16);
995
996       pitch_downsample(pre, pitch_buf, COMBFILTER_MAXPERIOD+N, CC);
997       /* Don't search for the fir last 1.5 octave of the range because
998          there's too many false-positives due to short-term correlation */
999       pitch_search(pitch_buf+(COMBFILTER_MAXPERIOD>>1), pitch_buf, N,
1000             COMBFILTER_MAXPERIOD-3*COMBFILTER_MINPERIOD, &pitch_index);
1001       pitch_index = COMBFILTER_MAXPERIOD-pitch_index;
1002
1003       gain1 = remove_doubling(pitch_buf, COMBFILTER_MAXPERIOD, COMBFILTER_MINPERIOD,
1004             N, &pitch_index, st->prefilter_period, st->prefilter_gain);
1005       if (pitch_index > COMBFILTER_MAXPERIOD-2)
1006          pitch_index = COMBFILTER_MAXPERIOD-2;
1007       gain1 = MULT16_16_Q15(QCONST16(.7f,15),gain1);
1008       /*printf("%d %d %f %f\n", pitch_change, pitch_index, gain1, st->analysis.tonality);*/
1009       if (st->loss_rate>2)
1010          gain1 = HALF32(gain1);
1011       if (st->loss_rate>4)
1012          gain1 = HALF32(gain1);
1013       if (st->loss_rate>8)
1014          gain1 = 0;
1015    } else {
1016       gain1 = 0;
1017       pitch_index = COMBFILTER_MINPERIOD;
1018    }
1019
1020    /* Gain threshold for enabling the prefilter/postfilter */
1021    pf_threshold = QCONST16(.2f,15);
1022
1023    /* Adjusting the threshold based on rate and continuity */
1024    if (abs(pitch_index-st->prefilter_period)*10>pitch_index)
1025       pf_threshold += QCONST16(.2f,15);
1026    if (nbAvailableBytes<25)
1027       pf_threshold += QCONST16(.1f,15);
1028    if (nbAvailableBytes<35)
1029       pf_threshold += QCONST16(.1f,15);
1030    if (st->prefilter_gain > QCONST16(.4f,15))
1031       pf_threshold -= QCONST16(.1f,15);
1032    if (st->prefilter_gain > QCONST16(.55f,15))
1033       pf_threshold -= QCONST16(.1f,15);
1034
1035    /* Hard threshold at 0.2 */
1036    pf_threshold = MAX16(pf_threshold, QCONST16(.2f,15));
1037    if (gain1<pf_threshold)
1038    {
1039       gain1 = 0;
1040       pf_on = 0;
1041       qg = 0;
1042    } else {
1043       /*This block is not gated by a total bits check only because
1044         of the nbAvailableBytes check above.*/
1045       if (ABS16(gain1-st->prefilter_gain)<QCONST16(.1f,15))
1046          gain1=st->prefilter_gain;
1047
1048 #ifdef FIXED_POINT
1049       qg = ((gain1+1536)>>10)/3-1;
1050 #else
1051       qg = (int)floor(.5f+gain1*32/3)-1;
1052 #endif
1053       qg = IMAX(0, IMIN(7, qg));
1054       gain1 = QCONST16(0.09375f,15)*(qg+1);
1055       pf_on = 1;
1056    }
1057    /*printf("%d %f\n", pitch_index, gain1);*/
1058
1059    c=0; do {
1060       int offset = mode->shortMdctSize-st->overlap;
1061       st->prefilter_period=IMAX(st->prefilter_period, COMBFILTER_MINPERIOD);
1062       OPUS_COPY(in+c*(N+st->overlap), st->in_mem+c*(st->overlap), st->overlap);
1063       if (offset)
1064          comb_filter(in+c*(N+st->overlap)+st->overlap, pre[c]+COMBFILTER_MAXPERIOD,
1065                st->prefilter_period, st->prefilter_period, offset, -st->prefilter_gain, -st->prefilter_gain,
1066                st->prefilter_tapset, st->prefilter_tapset, NULL, 0);
1067
1068       comb_filter(in+c*(N+st->overlap)+st->overlap+offset, pre[c]+COMBFILTER_MAXPERIOD+offset,
1069             st->prefilter_period, pitch_index, N-offset, -st->prefilter_gain, -gain1,
1070             st->prefilter_tapset, prefilter_tapset, mode->window, st->overlap);
1071       OPUS_COPY(st->in_mem+c*(st->overlap), in+c*(N+st->overlap)+N, st->overlap);
1072
1073       if (N>COMBFILTER_MAXPERIOD)
1074       {
1075          OPUS_MOVE(prefilter_mem+c*COMBFILTER_MAXPERIOD, pre[c]+N, COMBFILTER_MAXPERIOD);
1076       } else {
1077          OPUS_MOVE(prefilter_mem+c*COMBFILTER_MAXPERIOD, prefilter_mem+c*COMBFILTER_MAXPERIOD+N, COMBFILTER_MAXPERIOD-N);
1078          OPUS_MOVE(prefilter_mem+c*COMBFILTER_MAXPERIOD+COMBFILTER_MAXPERIOD-N, pre[c]+COMBFILTER_MAXPERIOD, N);
1079       }
1080    } while (++c<CC);
1081
1082    RESTORE_STACK;
1083    *gain = gain1;
1084    *pitch = pitch_index;
1085    *qgain = qg;
1086    return pf_on;
1087 }
1088
1089 int celt_encode_with_ec(CELTEncoder * OPUS_RESTRICT st, const opus_val16 * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes, ec_enc *enc)
1090 {
1091    int i, c, N;
1092    opus_int32 bits;
1093    ec_enc _enc;
1094    VARDECL(celt_sig, in);
1095    VARDECL(celt_sig, freq);
1096    VARDECL(celt_norm, X);
1097    VARDECL(celt_ener, bandE);
1098    VARDECL(opus_val16, bandLogE);
1099    VARDECL(opus_val16, bandLogE2);
1100    VARDECL(int, fine_quant);
1101    VARDECL(opus_val16, error);
1102    VARDECL(int, pulses);
1103    VARDECL(int, cap);
1104    VARDECL(int, offsets);
1105    VARDECL(int, fine_priority);
1106    VARDECL(int, tf_res);
1107    VARDECL(unsigned char, collapse_masks);
1108    celt_sig *prefilter_mem;
1109    opus_val16 *oldBandE, *oldLogE, *oldLogE2;
1110    int shortBlocks=0;
1111    int isTransient=0;
1112    const int CC = st->channels;
1113    const int C = st->stream_channels;
1114    int LM, M;
1115    int tf_select;
1116    int nbFilledBytes, nbAvailableBytes;
1117    int effEnd;
1118    int codedBands;
1119    int tf_sum;
1120    int alloc_trim;
1121    int pitch_index=COMBFILTER_MINPERIOD;
1122    opus_val16 gain1 = 0;
1123    int dual_stereo=0;
1124    int effectiveBytes;
1125    int dynalloc_logp;
1126    opus_int32 vbr_rate;
1127    opus_int32 total_bits;
1128    opus_int32 total_boost;
1129    opus_int32 balance;
1130    opus_int32 tell;
1131    int prefilter_tapset=0;
1132    int pf_on;
1133    int anti_collapse_rsv;
1134    int anti_collapse_on=0;
1135    int silence=0;
1136    int tf_chan = 0;
1137    opus_val16 tf_estimate;
1138    int pitch_change=0;
1139    opus_int32 tot_boost;
1140    opus_val16 sample_max;
1141    opus_val16 maxDepth;
1142    const OpusCustomMode *mode;
1143    int nbEBands;
1144    int overlap;
1145    const opus_int16 *eBands;
1146    int secondMdct;
1147    int signalBandwidth;
1148    ALLOC_STACK;
1149
1150    mode = st->mode;
1151    nbEBands = mode->nbEBands;
1152    overlap = mode->overlap;
1153    eBands = mode->eBands;
1154    tf_estimate = 0;
1155    if (nbCompressedBytes<2 || pcm==NULL)
1156      return OPUS_BAD_ARG;
1157
1158    frame_size *= st->upsample;
1159    for (LM=0;LM<=mode->maxLM;LM++)
1160       if (mode->shortMdctSize<<LM==frame_size)
1161          break;
1162    if (LM>mode->maxLM)
1163       return OPUS_BAD_ARG;
1164    M=1<<LM;
1165    N = M*mode->shortMdctSize;
1166
1167    prefilter_mem = st->in_mem+CC*(st->overlap);
1168    oldBandE = (opus_val16*)(st->in_mem+CC*(st->overlap+COMBFILTER_MAXPERIOD));
1169    oldLogE = oldBandE + CC*nbEBands;
1170    oldLogE2 = oldLogE + CC*nbEBands;
1171
1172    if (enc==NULL)
1173    {
1174       tell=1;
1175       nbFilledBytes=0;
1176    } else {
1177       tell=ec_tell(enc);
1178       nbFilledBytes=(tell+4)>>3;
1179    }
1180
1181 #ifdef CUSTOM_MODES
1182    if (st->signalling && enc==NULL)
1183    {
1184       int tmp = (mode->effEBands-st->end)>>1;
1185       st->end = IMAX(1, mode->effEBands-tmp);
1186       compressed[0] = tmp<<5;
1187       compressed[0] |= LM<<3;
1188       compressed[0] |= (C==2)<<2;
1189       /* Convert "standard mode" to Opus header */
1190       if (mode->Fs==48000 && mode->shortMdctSize==120)
1191       {
1192          int c0 = toOpus(compressed[0]);
1193          if (c0<0)
1194             return OPUS_BAD_ARG;
1195          compressed[0] = c0;
1196       }
1197       compressed++;
1198       nbCompressedBytes--;
1199    }
1200 #else
1201    celt_assert(st->signalling==0);
1202 #endif
1203
1204    /* Can't produce more than 1275 output bytes */
1205    nbCompressedBytes = IMIN(nbCompressedBytes,1275);
1206    nbAvailableBytes = nbCompressedBytes - nbFilledBytes;
1207
1208    if (st->vbr && st->bitrate!=OPUS_BITRATE_MAX)
1209    {
1210       opus_int32 den=mode->Fs>>BITRES;
1211       vbr_rate=(st->bitrate*frame_size+(den>>1))/den;
1212 #ifdef CUSTOM_MODES
1213       if (st->signalling)
1214          vbr_rate -= 8<<BITRES;
1215 #endif
1216       effectiveBytes = vbr_rate>>(3+BITRES);
1217    } else {
1218       opus_int32 tmp;
1219       vbr_rate = 0;
1220       tmp = st->bitrate*frame_size;
1221       if (tell>1)
1222          tmp += tell;
1223       if (st->bitrate!=OPUS_BITRATE_MAX)
1224          nbCompressedBytes = IMAX(2, IMIN(nbCompressedBytes,
1225                (tmp+4*mode->Fs)/(8*mode->Fs)-!!st->signalling));
1226       effectiveBytes = nbCompressedBytes;
1227    }
1228
1229    if (enc==NULL)
1230    {
1231       ec_enc_init(&_enc, compressed, nbCompressedBytes);
1232       enc = &_enc;
1233    }
1234
1235    if (vbr_rate>0)
1236    {
1237       /* Computes the max bit-rate allowed in VBR mode to avoid violating the
1238           target rate and buffering.
1239          We must do this up front so that bust-prevention logic triggers
1240           correctly if we don't have enough bits. */
1241       if (st->constrained_vbr)
1242       {
1243          opus_int32 vbr_bound;
1244          opus_int32 max_allowed;
1245          /* We could use any multiple of vbr_rate as bound (depending on the
1246              delay).
1247             This is clamped to ensure we use at least two bytes if the encoder
1248              was entirely empty, but to allow 0 in hybrid mode. */
1249          vbr_bound = vbr_rate;
1250          max_allowed = IMIN(IMAX(tell==1?2:0,
1251                (vbr_rate+vbr_bound-st->vbr_reservoir)>>(BITRES+3)),
1252                nbAvailableBytes);
1253          if(max_allowed < nbAvailableBytes)
1254          {
1255             nbCompressedBytes = nbFilledBytes+max_allowed;
1256             nbAvailableBytes = max_allowed;
1257             ec_enc_shrink(enc, nbCompressedBytes);
1258          }
1259       }
1260    }
1261    total_bits = nbCompressedBytes*8;
1262
1263    effEnd = st->end;
1264    if (effEnd > mode->effEBands)
1265       effEnd = mode->effEBands;
1266
1267    ALLOC(in, CC*(N+st->overlap), celt_sig);
1268
1269    sample_max=MAX16(st->overlap_max, celt_maxabs16(pcm, C*(N-overlap)/st->upsample));
1270    st->overlap_max=celt_maxabs16(pcm+C*(N-overlap)/st->upsample, C*overlap/st->upsample);
1271    sample_max=MAX16(sample_max, st->overlap_max);
1272 #ifdef FIXED_POINT
1273    silence = (sample_max==0);
1274 #else
1275    silence = (sample_max <= (opus_val16)1/(1<<st->lsb_depth));
1276 #endif
1277 #ifdef FUZZING
1278    if ((rand()&0x3F)==0)
1279       silence = 1;
1280 #endif
1281    if (tell==1)
1282       ec_enc_bit_logp(enc, silence, 15);
1283    else
1284       silence=0;
1285    if (silence)
1286    {
1287       /*In VBR mode there is no need to send more than the minimum. */
1288       if (vbr_rate>0)
1289       {
1290          effectiveBytes=nbCompressedBytes=IMIN(nbCompressedBytes, nbFilledBytes+2);
1291          total_bits=nbCompressedBytes*8;
1292          nbAvailableBytes=2;
1293          ec_enc_shrink(enc, nbCompressedBytes);
1294       }
1295       /* Pretend we've filled all the remaining bits with zeros
1296             (that's what the initialiser did anyway) */
1297       tell = nbCompressedBytes*8;
1298       enc->nbits_total+=tell-ec_tell(enc);
1299    }
1300    c=0; do {
1301       preemphasis(pcm+c, in+c*(N+st->overlap)+st->overlap, N, CC, st->upsample,
1302                   mode->preemph, st->preemph_memE+c, st->clip);
1303    } while (++c<CC);
1304
1305
1306
1307    /* Find pitch period and gain */
1308    {
1309       int enabled;
1310       int qg;
1311       enabled = nbAvailableBytes>12*C && st->start==0 && !silence && !st->disable_pf && st->complexity >= 5;
1312
1313       prefilter_tapset = st->tapset_decision;
1314       pf_on = run_prefilter(st, in, prefilter_mem, CC, N, prefilter_tapset, &pitch_index, &gain1, &qg, enabled, nbAvailableBytes);
1315       if ((gain1 > QCONST16(.4f,15) || st->prefilter_gain > QCONST16(.4f,15)) && st->analysis.tonality > .3
1316             && (pitch_index > 1.26*st->prefilter_period || pitch_index < .79*st->prefilter_period))
1317          pitch_change = 1;
1318       if (pf_on==0)
1319       {
1320          if(st->start==0 && tell+16<=total_bits)
1321             ec_enc_bit_logp(enc, 0, 1);
1322       } else {
1323          /*This block is not gated by a total bits check only because
1324            of the nbAvailableBytes check above.*/
1325          int octave;
1326          ec_enc_bit_logp(enc, 1, 1);
1327          pitch_index += 1;
1328          octave = EC_ILOG(pitch_index)-5;
1329          ec_enc_uint(enc, octave, 6);
1330          ec_enc_bits(enc, pitch_index-(16<<octave), 4+octave);
1331          pitch_index -= 1;
1332          ec_enc_bits(enc, qg, 3);
1333          ec_enc_icdf(enc, prefilter_tapset, tapset_icdf, 2);
1334       }
1335    }
1336
1337    isTransient = 0;
1338    shortBlocks = 0;
1339    if (LM>0 && ec_tell(enc)+3<=total_bits)
1340    {
1341       if (st->complexity >= 1)
1342       {
1343          isTransient = transient_analysis(in, N+st->overlap, CC,
1344                   &tf_estimate, &tf_chan);
1345          if (isTransient)
1346             shortBlocks = M;
1347       }
1348    }
1349
1350    ALLOC(freq, CC*N, celt_sig); /**< Interleaved signal MDCTs */
1351    ALLOC(bandE,nbEBands*CC, celt_ener);
1352    ALLOC(bandLogE,nbEBands*CC, opus_val16);
1353
1354    secondMdct = shortBlocks && st->complexity>=8;
1355    ALLOC(bandLogE2, C*nbEBands, opus_val16);
1356    if (secondMdct)
1357    {
1358       compute_mdcts(mode, 0, in, freq, CC, LM);
1359       if (CC==2&&C==1)
1360       {
1361          for (i=0;i<N;i++)
1362             freq[i] = ADD32(HALF32(freq[i]), HALF32(freq[N+i]));
1363       }
1364       if (st->upsample != 1)
1365       {
1366          c=0; do
1367          {
1368             int bound = N/st->upsample;
1369             for (i=0;i<bound;i++)
1370                freq[c*N+i] *= st->upsample;
1371             for (;i<N;i++)
1372                freq[c*N+i] = 0;
1373          } while (++c<C);
1374       }
1375       compute_band_energies(mode, freq, bandE, effEnd, C, M);
1376       amp2Log2(mode, effEnd, st->end, bandE, bandLogE2, C);
1377       for (i=0;i<C*nbEBands;i++)
1378          bandLogE2[i] += HALF16(SHL16(LM, DB_SHIFT));
1379    }
1380
1381    compute_mdcts(mode, shortBlocks, in, freq, CC, LM);
1382
1383    if (CC==2&&C==1)
1384    {
1385       for (i=0;i<N;i++)
1386          freq[i] = ADD32(HALF32(freq[i]), HALF32(freq[N+i]));
1387       tf_chan = 0;
1388    }
1389    if (st->upsample != 1)
1390    {
1391       c=0; do
1392       {
1393          int bound = N/st->upsample;
1394          for (i=0;i<bound;i++)
1395             freq[c*N+i] *= st->upsample;
1396          for (;i<N;i++)
1397             freq[c*N+i] = 0;
1398       } while (++c<C);
1399    }
1400    compute_band_energies(mode, freq, bandE, effEnd, C, M);
1401
1402    amp2Log2(mode, effEnd, st->end, bandE, bandLogE, C);
1403    /*for (i=0;i<21;i++)
1404       printf("%f ", bandLogE[i]);
1405    printf("\n");*/
1406
1407    if (!secondMdct)
1408    {
1409       for (i=0;i<C*nbEBands;i++)
1410          bandLogE2[i] = bandLogE[i];
1411    }
1412
1413    /* Last chance to catch any transient we might have missed in the
1414       time-domain analysis */
1415    if (LM>0 && ec_tell(enc)+3<=total_bits && !isTransient && st->complexity>=5)
1416    {
1417       if (patch_transient_decision(bandLogE, oldBandE, nbEBands, st->end, C))
1418       {
1419          isTransient = 1;
1420          shortBlocks = M;
1421          compute_mdcts(mode, shortBlocks, in, freq, CC, LM);
1422          if (CC==2&&C==1)
1423          {
1424             for (i=0;i<N;i++)
1425                freq[i] = ADD32(HALF32(freq[i]), HALF32(freq[N+i]));
1426          }
1427          if (st->upsample != 1)
1428          {
1429             c=0; do
1430             {
1431                int bound = N/st->upsample;
1432                for (i=0;i<bound;i++)
1433                   freq[c*N+i] *= st->upsample;
1434                for (;i<N;i++)
1435                   freq[c*N+i] = 0;
1436             } while (++c<C);
1437          }
1438          compute_band_energies(mode, freq, bandE, effEnd, C, M);
1439          amp2Log2(mode, effEnd, st->end, bandE, bandLogE, C);
1440          /* Compensate for the scaling of short vs long mdcts */
1441          for (i=0;i<C*nbEBands;i++)
1442             bandLogE2[i] += HALF16(SHL16(LM, DB_SHIFT));
1443          tf_estimate = QCONST16(.2,14);
1444       }
1445    }
1446
1447    if (LM>0 && ec_tell(enc)+3<=total_bits)
1448       ec_enc_bit_logp(enc, isTransient, 3);
1449
1450    ALLOC(X, C*N, celt_norm);         /**< Interleaved normalised MDCTs */
1451
1452    /* Band normalisation */
1453    normalise_bands(mode, freq, X, bandE, effEnd, C, M);
1454
1455    ALLOC(tf_res, nbEBands, int);
1456    /* Disable variable tf resolution for hybrid and at very low bitrate */
1457    if (effectiveBytes>=15*C && st->start==0 && st->complexity>=2)
1458    {
1459       int lambda;
1460       if (effectiveBytes<40)
1461          lambda = 12;
1462       else if (effectiveBytes<60)
1463          lambda = 6;
1464       else if (effectiveBytes<100)
1465          lambda = 4;
1466       else
1467          lambda = 3;
1468       lambda*=2;
1469       tf_select = tf_analysis(mode, effEnd, isTransient, tf_res, lambda, X, N, LM, &tf_sum, tf_estimate, tf_chan);
1470       for (i=effEnd;i<st->end;i++)
1471          tf_res[i] = tf_res[effEnd-1];
1472    } else {
1473       tf_sum = 0;
1474       for (i=0;i<st->end;i++)
1475          tf_res[i] = isTransient;
1476       tf_select=0;
1477    }
1478
1479    ALLOC(error, C*nbEBands, opus_val16);
1480    quant_coarse_energy(mode, st->start, st->end, effEnd, bandLogE,
1481          oldBandE, total_bits, error, enc,
1482          C, LM, nbAvailableBytes, st->force_intra,
1483          &st->delayedIntra, st->complexity >= 4, st->loss_rate);
1484
1485    tf_encode(st->start, st->end, isTransient, tf_res, LM, tf_select, enc);
1486
1487    if (ec_tell(enc)+4<=total_bits)
1488    {
1489       if (shortBlocks || st->complexity < 3 || nbAvailableBytes < 10*C || st->start != 0)
1490       {
1491          if (st->complexity == 0)
1492             st->spread_decision = SPREAD_NONE;
1493          else
1494             st->spread_decision = SPREAD_NORMAL;
1495       } else {
1496          /* Disable new spreading+tapset estimator until we can show it works
1497             better than the old one. So far it seems like spreading_decision()
1498             works best. */
1499          if (0&&st->analysis.valid)
1500          {
1501             static const opus_val16 spread_thresholds[3] = {-QCONST16(.6f, 15), -QCONST16(.2f, 15), -QCONST16(.07f, 15)};
1502             static const opus_val16 spread_histeresis[3] = {QCONST16(.15f, 15), QCONST16(.07f, 15), QCONST16(.02f, 15)};
1503             static const opus_val16 tapset_thresholds[2] = {QCONST16(.0f, 15), QCONST16(.15f, 15)};
1504             static const opus_val16 tapset_histeresis[2] = {QCONST16(.1f, 15), QCONST16(.05f, 15)};
1505             st->spread_decision = hysteresis_decision(-st->analysis.tonality, spread_thresholds, spread_histeresis, 3, st->spread_decision);
1506             st->tapset_decision = hysteresis_decision(st->analysis.tonality_slope, tapset_thresholds, tapset_histeresis, 2, st->tapset_decision);
1507          } else {
1508             st->spread_decision = spreading_decision(mode, X,
1509                   &st->tonal_average, st->spread_decision, &st->hf_average,
1510                   &st->tapset_decision, pf_on&&!shortBlocks, effEnd, C, M);
1511          }
1512          /*printf("%d %d\n", st->tapset_decision, st->spread_decision);*/
1513          /*printf("%f %d %f %d\n\n", st->analysis.tonality, st->spread_decision, st->analysis.tonality_slope, st->tapset_decision);*/
1514       }
1515       ec_enc_icdf(enc, st->spread_decision, spread_icdf, 5);
1516    }
1517
1518    ALLOC(offsets, nbEBands, int);
1519
1520    maxDepth = dynalloc_analysis(bandLogE, bandLogE2, nbEBands, st->start, st->end, C, offsets,
1521          st->lsb_depth, mode->logN, isTransient, st->vbr, st->constrained_vbr,
1522          eBands, LM, effectiveBytes, &tot_boost);
1523    ALLOC(cap, nbEBands, int);
1524    init_caps(mode,cap,LM,C);
1525
1526    dynalloc_logp = 6;
1527    total_bits<<=BITRES;
1528    total_boost = 0;
1529    tell = ec_tell_frac(enc);
1530    for (i=st->start;i<st->end;i++)
1531    {
1532       int width, quanta;
1533       int dynalloc_loop_logp;
1534       int boost;
1535       int j;
1536       width = C*(eBands[i+1]-eBands[i])<<LM;
1537       /* quanta is 6 bits, but no more than 1 bit/sample
1538          and no less than 1/8 bit/sample */
1539       quanta = IMIN(width<<BITRES, IMAX(6<<BITRES, width));
1540       dynalloc_loop_logp = dynalloc_logp;
1541       boost = 0;
1542       for (j = 0; tell+(dynalloc_loop_logp<<BITRES) < total_bits-total_boost
1543             && boost < cap[i]; j++)
1544       {
1545          int flag;
1546          flag = j<offsets[i];
1547          ec_enc_bit_logp(enc, flag, dynalloc_loop_logp);
1548          tell = ec_tell_frac(enc);
1549          if (!flag)
1550             break;
1551          boost += quanta;
1552          total_boost += quanta;
1553          dynalloc_loop_logp = 1;
1554       }
1555       /* Making dynalloc more likely */
1556       if (j)
1557          dynalloc_logp = IMAX(2, dynalloc_logp-1);
1558       offsets[i] = boost;
1559    }
1560
1561    if (C==2)
1562    {
1563       int effectiveRate;
1564
1565       static const opus_val16 intensity_thresholds[21]=
1566       /* 0  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19  20  off*/
1567         { 16,21,23,25,27,29,31,33,35,38,42,46,50,54,58,63,68,75,84,102,130};
1568       static const opus_val16 intensity_histeresis[21]=
1569         {  2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3, 4, 5, 6,  8, 12};
1570
1571       /* Always use MS for 2.5 ms frames until we can do a better analysis */
1572       if (LM!=0)
1573          dual_stereo = stereo_analysis(mode, X, LM, N);
1574
1575       /* Account for coarse energy */
1576       effectiveRate = (8*effectiveBytes - 80)>>LM;
1577
1578       /* effectiveRate in kb/s */
1579       effectiveRate = 2*effectiveRate/5;
1580
1581       st->intensity = hysteresis_decision((opus_val16)effectiveRate, intensity_thresholds, intensity_histeresis, 21, st->intensity);
1582       st->intensity = IMIN(st->end,IMAX(st->start, st->intensity));
1583    }
1584
1585    alloc_trim = 5;
1586    if (tell+(6<<BITRES) <= total_bits - total_boost)
1587    {
1588       alloc_trim = alloc_trim_analysis(mode, X, bandLogE,
1589             st->end, LM, C, N, &st->analysis, &st->stereo_saving, tf_estimate, st->intensity);
1590       ec_enc_icdf(enc, alloc_trim, trim_icdf, 7);
1591       tell = ec_tell_frac(enc);
1592    }
1593
1594    /* Variable bitrate */
1595    if (vbr_rate>0)
1596    {
1597      opus_val16 alpha;
1598      opus_int32 delta;
1599      /* The target rate in 8th bits per frame */
1600      opus_int32 target, base_target;
1601      opus_int32 min_allowed;
1602      int coded_bins;
1603      int coded_bands;
1604      int lm_diff = mode->maxLM - LM;
1605      coded_bands = st->lastCodedBands ? st->lastCodedBands : nbEBands;
1606      coded_bins = eBands[coded_bands]<<LM;
1607      if (C==2)
1608         coded_bins += eBands[IMIN(st->intensity, coded_bands)]<<LM;
1609
1610      /* Don't attempt to use more than 510 kb/s, even for frames smaller than 20 ms.
1611         The CELT allocator will just not be able to use more than that anyway. */
1612      nbCompressedBytes = IMIN(nbCompressedBytes,1275>>(3-LM));
1613      target = vbr_rate - ((40*C+20)<<BITRES);
1614      base_target = target;
1615
1616      if (st->constrained_vbr)
1617         target += (st->vbr_offset>>lm_diff);
1618
1619      /*printf("%f %f %f %f %d %d ", st->analysis.activity, st->analysis.tonality, tf_estimate, st->stereo_saving, tot_boost, coded_bands);*/
1620 #ifndef FIXED_POINT
1621      if (st->analysis.valid && st->analysis.activity<.4)
1622         target -= (opus_int32)((coded_bins<<BITRES)*(.4f-st->analysis.activity));
1623 #endif
1624      /* Stereo savings */
1625      if (C==2)
1626      {
1627         int coded_stereo_bands;
1628         int coded_stereo_dof;
1629         opus_val16 max_frac;
1630         coded_stereo_bands = IMIN(st->intensity, coded_bands);
1631         coded_stereo_dof = (eBands[coded_stereo_bands]<<LM)-coded_stereo_bands;
1632         /* Maximum fraction of the bits we can save if the signal is mono. */
1633         max_frac = DIV32_16(MULT16_16(QCONST16(0.8f, 15), coded_stereo_dof), coded_bins);
1634         /*printf("%d %d %d ", coded_stereo_dof, coded_bins, tot_boost);*/
1635         target -= (opus_int32)MIN32(MULT16_32_Q15(max_frac,target),
1636                         SHR16(MULT16_16(st->stereo_saving-QCONST16(0.1f,8),(coded_stereo_dof<<BITRES)),8));
1637      }
1638      /* Boost the rate according to dynalloc (minus the dynalloc average for calibration). */
1639      target += tot_boost-(16<<LM);
1640      /* Apply transient boost, compensating for average boost. */
1641      target += (opus_int32)SHL32(MULT16_32_Q15(tf_estimate-QCONST16(0.04f,14), target),1);
1642
1643 #ifndef FIXED_POINT
1644      /* Apply tonality boost */
1645      if (st->analysis.valid) {
1646         opus_int32 tonal_target;
1647         float tonal;
1648
1649         /* Tonality boost (compensating for the average). */
1650         tonal = MAX16(0.f,st->analysis.tonality-.15f)-0.09f;
1651         tonal_target = target + (opus_int32)((coded_bins<<BITRES)*1.2f*tonal);
1652         if (pitch_change)
1653            tonal_target +=  (opus_int32)((coded_bins<<BITRES)*.8f);
1654         /*printf("%f %f ", st->analysis.tonality, tonal);*/
1655         target = tonal_target;
1656      }
1657 #endif
1658
1659      {
1660         opus_int32 floor_depth;
1661         int bins;
1662         bins = eBands[nbEBands-2]<<LM;
1663         /*floor_depth = SHR32(MULT16_16((C*bins<<BITRES),celt_log2(SHL32(MAX16(1,sample_max),13))), DB_SHIFT);*/
1664         floor_depth = (opus_int32)SHR32(MULT16_16((C*bins<<BITRES),maxDepth), DB_SHIFT);
1665         floor_depth = IMAX(floor_depth, target>>2);
1666         target = IMIN(target, floor_depth);
1667         /*printf("%f %d\n", maxDepth, floor_depth);*/
1668      }
1669
1670      if (st->constrained_vbr || st->bitrate<64000)
1671      {
1672         opus_val16 rate_factor;
1673 #ifdef FIXED_POINT
1674         rate_factor = MAX16(0,(st->bitrate-32000));
1675 #else
1676         rate_factor = MAX16(0,(1.f/32768)*(st->bitrate-32000));
1677 #endif
1678         if (st->constrained_vbr)
1679            rate_factor = MIN16(rate_factor, QCONST16(0.67f, 15));
1680         target = base_target + (opus_int32)MULT16_32_Q15(rate_factor, target-base_target);
1681
1682      }
1683      /* Don't allow more than doubling the rate */
1684      target = IMIN(2*base_target, target);
1685
1686      /* The current offset is removed from the target and the space used
1687         so far is added*/
1688      target=target+tell;
1689      /* In VBR mode the frame size must not be reduced so much that it would
1690          result in the encoder running out of bits.
1691         The margin of 2 bytes ensures that none of the bust-prevention logic
1692          in the decoder will have triggered so far. */
1693      min_allowed = ((tell+total_boost+(1<<(BITRES+3))-1)>>(BITRES+3)) + 2 - nbFilledBytes;
1694
1695      nbAvailableBytes = (target+(1<<(BITRES+2)))>>(BITRES+3);
1696      nbAvailableBytes = IMAX(min_allowed,nbAvailableBytes);
1697      nbAvailableBytes = IMIN(nbCompressedBytes,nbAvailableBytes+nbFilledBytes) - nbFilledBytes;
1698
1699      /* By how much did we "miss" the target on that frame */
1700      delta = target - vbr_rate;
1701
1702      target=nbAvailableBytes<<(BITRES+3);
1703
1704      /*If the frame is silent we don't adjust our drift, otherwise
1705        the encoder will shoot to very high rates after hitting a
1706        span of silence, but we do allow the bitres to refill.
1707        This means that we'll undershoot our target in CVBR/VBR modes
1708        on files with lots of silence. */
1709      if(silence)
1710      {
1711        nbAvailableBytes = 2;
1712        target = 2*8<<BITRES;
1713        delta = 0;
1714      }
1715
1716      if (st->vbr_count < 970)
1717      {
1718         st->vbr_count++;
1719         alpha = celt_rcp(SHL32(EXTEND32(st->vbr_count+20),16));
1720      } else
1721         alpha = QCONST16(.001f,15);
1722      /* How many bits have we used in excess of what we're allowed */
1723      if (st->constrained_vbr)
1724         st->vbr_reservoir += target - vbr_rate;
1725      /*printf ("%d\n", st->vbr_reservoir);*/
1726
1727      /* Compute the offset we need to apply in order to reach the target */
1728      if (st->constrained_vbr)
1729      {
1730         st->vbr_drift += (opus_int32)MULT16_32_Q15(alpha,(delta*(1<<lm_diff))-st->vbr_offset-st->vbr_drift);
1731         st->vbr_offset = -st->vbr_drift;
1732      }
1733      /*printf ("%d\n", st->vbr_drift);*/
1734
1735      if (st->constrained_vbr && st->vbr_reservoir < 0)
1736      {
1737         /* We're under the min value -- increase rate */
1738         int adjust = (-st->vbr_reservoir)/(8<<BITRES);
1739         /* Unless we're just coding silence */
1740         nbAvailableBytes += silence?0:adjust;
1741         st->vbr_reservoir = 0;
1742         /*printf ("+%d\n", adjust);*/
1743      }
1744      nbCompressedBytes = IMIN(nbCompressedBytes,nbAvailableBytes+nbFilledBytes);
1745      /*printf("%d\n", nbCompressedBytes*50*8);*/
1746      /* This moves the raw bits to take into account the new compressed size */
1747      ec_enc_shrink(enc, nbCompressedBytes);
1748    }
1749
1750    /* Bit allocation */
1751    ALLOC(fine_quant, nbEBands, int);
1752    ALLOC(pulses, nbEBands, int);
1753    ALLOC(fine_priority, nbEBands, int);
1754
1755    /* bits =           packet size                    - where we are - safety*/
1756    bits = (((opus_int32)nbCompressedBytes*8)<<BITRES) - ec_tell_frac(enc) - 1;
1757    anti_collapse_rsv = isTransient&&LM>=2&&bits>=((LM+2)<<BITRES) ? (1<<BITRES) : 0;
1758    bits -= anti_collapse_rsv;
1759    signalBandwidth = st->end-1;
1760 #ifndef FIXED_POINT
1761    if (st->analysis.valid)
1762       signalBandwidth = st->analysis.bandwidth;
1763 #endif
1764    codedBands = compute_allocation(mode, st->start, st->end, offsets, cap,
1765          alloc_trim, &st->intensity, &dual_stereo, bits, &balance, pulses,
1766          fine_quant, fine_priority, C, LM, enc, 1, st->lastCodedBands, signalBandwidth);
1767    st->lastCodedBands = codedBands;
1768
1769    quant_fine_energy(mode, st->start, st->end, oldBandE, error, fine_quant, enc, C);
1770
1771 #ifdef MEASURE_NORM_MSE
1772    float X0[3000];
1773    float bandE0[60];
1774    c=0; do
1775       for (i=0;i<N;i++)
1776          X0[i+c*N] = X[i+c*N];
1777    while (++c<C);
1778    for (i=0;i<C*nbEBands;i++)
1779       bandE0[i] = bandE[i];
1780 #endif
1781
1782    /* Residual quantisation */
1783    ALLOC(collapse_masks, C*nbEBands, unsigned char);
1784    quant_all_bands(1, mode, st->start, st->end, X, C==2 ? X+N : NULL, collapse_masks,
1785          bandE, pulses, shortBlocks, st->spread_decision, dual_stereo, st->intensity, tf_res,
1786          nbCompressedBytes*(8<<BITRES)-anti_collapse_rsv, balance, enc, LM, codedBands, &st->rng);
1787
1788    if (anti_collapse_rsv > 0)
1789    {
1790       anti_collapse_on = st->consec_transient<2;
1791 #ifdef FUZZING
1792       anti_collapse_on = rand()&0x1;
1793 #endif
1794       ec_enc_bits(enc, anti_collapse_on, 1);
1795    }
1796    quant_energy_finalise(mode, st->start, st->end, oldBandE, error, fine_quant, fine_priority, nbCompressedBytes*8-ec_tell(enc), enc, C);
1797
1798    if (silence)
1799    {
1800       for (i=0;i<C*nbEBands;i++)
1801          oldBandE[i] = -QCONST16(28.f,DB_SHIFT);
1802    }
1803
1804 #ifdef RESYNTH
1805    /* Re-synthesis of the coded audio if required */
1806    {
1807       celt_sig *out_mem[2];
1808
1809       log2Amp(mode, st->start, st->end, bandE, oldBandE, C);
1810       if (silence)
1811       {
1812          for (i=0;i<C*nbEBands;i++)
1813             bandE[i] = 0;
1814       }
1815
1816 #ifdef MEASURE_NORM_MSE
1817       measure_norm_mse(mode, X, X0, bandE, bandE0, M, N, C);
1818 #endif
1819       if (anti_collapse_on)
1820       {
1821          anti_collapse(mode, X, collapse_masks, LM, C, N,
1822                st->start, st->end, oldBandE, oldLogE, oldLogE2, pulses, st->rng);
1823       }
1824
1825       /* Synthesis */
1826       denormalise_bands(mode, X, freq, bandE, st->start, effEnd, C, M);
1827
1828       c=0; do {
1829          OPUS_MOVE(st->syn_mem[c], st->syn_mem[c]+N, 2*MAX_PERIOD-N+overlap/2);
1830       } while (++c<CC);
1831
1832       if (CC==2&&C==1)
1833       {
1834          for (i=0;i<N;i++)
1835             freq[N+i] = freq[i];
1836       }
1837
1838       c=0; do {
1839          out_mem[c] = st->syn_mem[c]+2*MAX_PERIOD-N;
1840       } while (++c<CC);
1841
1842       compute_inv_mdcts(mode, shortBlocks, freq, out_mem, CC, LM);
1843
1844       c=0; do {
1845          st->prefilter_period=IMAX(st->prefilter_period, COMBFILTER_MINPERIOD);
1846          st->prefilter_period_old=IMAX(st->prefilter_period_old, COMBFILTER_MINPERIOD);
1847          comb_filter(out_mem[c], out_mem[c], st->prefilter_period_old, st->prefilter_period, mode->shortMdctSize,
1848                st->prefilter_gain_old, st->prefilter_gain, st->prefilter_tapset_old, st->prefilter_tapset,
1849                mode->window, st->overlap);
1850          if (LM!=0)
1851             comb_filter(out_mem[c]+mode->shortMdctSize, out_mem[c]+mode->shortMdctSize, st->prefilter_period, pitch_index, N-mode->shortMdctSize,
1852                   st->prefilter_gain, gain1, st->prefilter_tapset, prefilter_tapset,
1853                   mode->window, overlap);
1854       } while (++c<CC);
1855
1856       /* We reuse freq[] as scratch space for the de-emphasis */
1857       deemphasis(out_mem, (opus_val16*)pcm, N, CC, st->upsample, mode->preemph, st->preemph_memD, freq);
1858       st->prefilter_period_old = st->prefilter_period;
1859       st->prefilter_gain_old = st->prefilter_gain;
1860       st->prefilter_tapset_old = st->prefilter_tapset;
1861    }
1862 #endif
1863
1864    st->prefilter_period = pitch_index;
1865    st->prefilter_gain = gain1;
1866    st->prefilter_tapset = prefilter_tapset;
1867 #ifdef RESYNTH
1868    if (LM!=0)
1869    {
1870       st->prefilter_period_old = st->prefilter_period;
1871       st->prefilter_gain_old = st->prefilter_gain;
1872       st->prefilter_tapset_old = st->prefilter_tapset;
1873    }
1874 #endif
1875
1876    if (CC==2&&C==1) {
1877       for (i=0;i<nbEBands;i++)
1878          oldBandE[nbEBands+i]=oldBandE[i];
1879    }
1880
1881    if (!isTransient)
1882    {
1883       for (i=0;i<CC*nbEBands;i++)
1884          oldLogE2[i] = oldLogE[i];
1885       for (i=0;i<CC*nbEBands;i++)
1886          oldLogE[i] = oldBandE[i];
1887    } else {
1888       for (i=0;i<CC*nbEBands;i++)
1889          oldLogE[i] = MIN16(oldLogE[i], oldBandE[i]);
1890    }
1891    /* In case start or end were to change */
1892    c=0; do
1893    {
1894       for (i=0;i<st->start;i++)
1895       {
1896          oldBandE[c*nbEBands+i]=0;
1897          oldLogE[c*nbEBands+i]=oldLogE2[c*nbEBands+i]=-QCONST16(28.f,DB_SHIFT);
1898       }
1899       for (i=st->end;i<nbEBands;i++)
1900       {
1901          oldBandE[c*nbEBands+i]=0;
1902          oldLogE[c*nbEBands+i]=oldLogE2[c*nbEBands+i]=-QCONST16(28.f,DB_SHIFT);
1903       }
1904    } while (++c<CC);
1905
1906    if (isTransient)
1907       st->consec_transient++;
1908    else
1909       st->consec_transient=0;
1910    st->rng = enc->rng;
1911
1912    /* If there's any room left (can only happen for very high rates),
1913       it's already filled with zeros */
1914    ec_enc_done(enc);
1915
1916 #ifdef CUSTOM_MODES
1917    if (st->signalling)
1918       nbCompressedBytes++;
1919 #endif
1920
1921    RESTORE_STACK;
1922    if (ec_get_error(enc))
1923       return OPUS_INTERNAL_ERROR;
1924    else
1925       return nbCompressedBytes;
1926 }
1927
1928
1929 #ifdef CUSTOM_MODES
1930
1931 #ifdef FIXED_POINT
1932 int opus_custom_encode(CELTEncoder * OPUS_RESTRICT st, const opus_int16 * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes)
1933 {
1934    return celt_encode_with_ec(st, pcm, frame_size, compressed, nbCompressedBytes, NULL);
1935 }
1936
1937 #ifndef DISABLE_FLOAT_API
1938 int opus_custom_encode_float(CELTEncoder * OPUS_RESTRICT st, const float * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes)
1939 {
1940    int j, ret, C, N;
1941    VARDECL(opus_int16, in);
1942    ALLOC_STACK;
1943
1944    if (pcm==NULL)
1945       return OPUS_BAD_ARG;
1946
1947    C = st->channels;
1948    N = frame_size;
1949    ALLOC(in, C*N, opus_int16);
1950
1951    for (j=0;j<C*N;j++)
1952      in[j] = FLOAT2INT16(pcm[j]);
1953
1954    ret=celt_encode_with_ec(st,in,frame_size,compressed,nbCompressedBytes, NULL);
1955 #ifdef RESYNTH
1956    for (j=0;j<C*N;j++)
1957       ((float*)pcm)[j]=in[j]*(1.f/32768.f);
1958 #endif
1959    RESTORE_STACK;
1960    return ret;
1961 }
1962 #endif /* DISABLE_FLOAT_API */
1963 #else
1964
1965 int opus_custom_encode(CELTEncoder * OPUS_RESTRICT st, const opus_int16 * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes)
1966 {
1967    int j, ret, C, N;
1968    VARDECL(celt_sig, in);
1969    ALLOC_STACK;
1970
1971    if (pcm==NULL)
1972       return OPUS_BAD_ARG;
1973
1974    C=st->channels;
1975    N=frame_size;
1976    ALLOC(in, C*N, celt_sig);
1977    for (j=0;j<C*N;j++) {
1978      in[j] = SCALEOUT(pcm[j]);
1979    }
1980
1981    ret = celt_encode_with_ec(st,in,frame_size,compressed,nbCompressedBytes, NULL);
1982 #ifdef RESYNTH
1983    for (j=0;j<C*N;j++)
1984       ((opus_int16*)pcm)[j] = FLOAT2INT16(in[j]);
1985 #endif
1986    RESTORE_STACK;
1987    return ret;
1988 }
1989
1990 int opus_custom_encode_float(CELTEncoder * OPUS_RESTRICT st, const float * pcm, int frame_size, unsigned char *compressed, int nbCompressedBytes)
1991 {
1992    return celt_encode_with_ec(st, pcm, frame_size, compressed, nbCompressedBytes, NULL);
1993 }
1994
1995 #endif
1996
1997 #endif /* CUSTOM_MODES */
1998
1999 int opus_custom_encoder_ctl(CELTEncoder * OPUS_RESTRICT st, int request, ...)
2000 {
2001    va_list ap;
2002
2003    va_start(ap, request);
2004    switch (request)
2005    {
2006       case OPUS_SET_COMPLEXITY_REQUEST:
2007       {
2008          int value = va_arg(ap, opus_int32);
2009          if (value<0 || value>10)
2010             goto bad_arg;
2011          st->complexity = value;
2012       }
2013       break;
2014       case CELT_SET_START_BAND_REQUEST:
2015       {
2016          opus_int32 value = va_arg(ap, opus_int32);
2017          if (value<0 || value>=st->mode->nbEBands)
2018             goto bad_arg;
2019          st->start = value;
2020       }
2021       break;
2022       case CELT_SET_END_BAND_REQUEST:
2023       {
2024          opus_int32 value = va_arg(ap, opus_int32);
2025          if (value<1 || value>st->mode->nbEBands)
2026             goto bad_arg;
2027          st->end = value;
2028       }
2029       break;
2030       case CELT_SET_PREDICTION_REQUEST:
2031       {
2032          int value = va_arg(ap, opus_int32);
2033          if (value<0 || value>2)
2034             goto bad_arg;
2035          st->disable_pf = value<=1;
2036          st->force_intra = value==0;
2037       }
2038       break;
2039       case OPUS_SET_PACKET_LOSS_PERC_REQUEST:
2040       {
2041          int value = va_arg(ap, opus_int32);
2042          if (value<0 || value>100)
2043             goto bad_arg;
2044          st->loss_rate = value;
2045       }
2046       break;
2047       case OPUS_SET_VBR_CONSTRAINT_REQUEST:
2048       {
2049          opus_int32 value = va_arg(ap, opus_int32);
2050          st->constrained_vbr = value;
2051       }
2052       break;
2053       case OPUS_SET_VBR_REQUEST:
2054       {
2055          opus_int32 value = va_arg(ap, opus_int32);
2056          st->vbr = value;
2057       }
2058       break;
2059       case OPUS_SET_BITRATE_REQUEST:
2060       {
2061          opus_int32 value = va_arg(ap, opus_int32);
2062          if (value<=500 && value!=OPUS_BITRATE_MAX)
2063             goto bad_arg;
2064          value = IMIN(value, 260000*st->channels);
2065          st->bitrate = value;
2066       }
2067       break;
2068       case CELT_SET_CHANNELS_REQUEST:
2069       {
2070          opus_int32 value = va_arg(ap, opus_int32);
2071          if (value<1 || value>2)
2072             goto bad_arg;
2073          st->stream_channels = value;
2074       }
2075       break;
2076       case OPUS_SET_LSB_DEPTH_REQUEST:
2077       {
2078           opus_int32 value = va_arg(ap, opus_int32);
2079           if (value<8 || value>24)
2080              goto bad_arg;
2081           st->lsb_depth=value;
2082       }
2083       break;
2084       case OPUS_GET_LSB_DEPTH_REQUEST:
2085       {
2086           opus_int32 *value = va_arg(ap, opus_int32*);
2087           *value=st->lsb_depth;
2088       }
2089       break;
2090       case OPUS_RESET_STATE:
2091       {
2092          int i;
2093          opus_val16 *oldBandE, *oldLogE, *oldLogE2;
2094          oldBandE = (opus_val16*)(st->in_mem+st->channels*(st->overlap+COMBFILTER_MAXPERIOD));
2095          oldLogE = oldBandE + st->channels*st->mode->nbEBands;
2096          oldLogE2 = oldLogE + st->channels*st->mode->nbEBands;
2097          OPUS_CLEAR((char*)&st->ENCODER_RESET_START,
2098                opus_custom_encoder_get_size(st->mode, st->channels)-
2099                ((char*)&st->ENCODER_RESET_START - (char*)st));
2100          for (i=0;i<st->channels*st->mode->nbEBands;i++)
2101             oldLogE[i]=oldLogE2[i]=-QCONST16(28.f,DB_SHIFT);
2102          st->vbr_offset = 0;
2103          st->delayedIntra = 1;
2104          st->spread_decision = SPREAD_NORMAL;
2105          st->tonal_average = 256;
2106          st->hf_average = 0;
2107          st->tapset_decision = 0;
2108       }
2109       break;
2110 #ifdef CUSTOM_MODES
2111       case CELT_SET_INPUT_CLIPPING_REQUEST:
2112       {
2113          opus_int32 value = va_arg(ap, opus_int32);
2114          st->clip = value;
2115       }
2116       break;
2117 #endif
2118       case CELT_SET_SIGNALLING_REQUEST:
2119       {
2120          opus_int32 value = va_arg(ap, opus_int32);
2121          st->signalling = value;
2122       }
2123       break;
2124       case CELT_SET_ANALYSIS_REQUEST:
2125       {
2126          AnalysisInfo *info = va_arg(ap, AnalysisInfo *);
2127          if (info)
2128             OPUS_COPY(&st->analysis, info, 1);
2129       }
2130       break;
2131       case CELT_GET_MODE_REQUEST:
2132       {
2133          const CELTMode ** value = va_arg(ap, const CELTMode**);
2134          if (value==0)
2135             goto bad_arg;
2136          *value=st->mode;
2137       }
2138       break;
2139       case OPUS_GET_FINAL_RANGE_REQUEST:
2140       {
2141          opus_uint32 * value = va_arg(ap, opus_uint32 *);
2142          if (value==0)
2143             goto bad_arg;
2144          *value=st->rng;
2145       }
2146       break;
2147       default:
2148          goto bad_request;
2149    }
2150    va_end(ap);
2151    return OPUS_OK;
2152 bad_arg:
2153    va_end(ap);
2154    return OPUS_BAD_ARG;
2155 bad_request:
2156    va_end(ap);
2157    return OPUS_UNIMPLEMENTED;
2158 }