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