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