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