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