Making the encoder fail less often at stupidly low bit-rates
[opus.git] / libcelt / quant_bands.c
1 /* Copyright (c) 2007-2008 CSIRO
2    Copyright (c) 2007-2009 Xiph.Org Foundation
3    Written by Jean-Marc Valin */
4 /*
5    Redistribution and use in source and binary forms, with or without
6    modification, are permitted provided that the following conditions
7    are met:
8    
9    - Redistributions of source code must retain the above copyright
10    notice, this list of conditions and the following disclaimer.
11    
12    - Redistributions in binary form must reproduce the above copyright
13    notice, this list of conditions and the following disclaimer in the
14    documentation and/or other materials provided with the distribution.
15    
16    - Neither the name of the Xiph.org Foundation nor the names of its
17    contributors may be used to endorse or promote products derived from
18    this software without specific prior written permission.
19    
20    THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
21    ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
22    LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
23    A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR
24    CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
25    EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
26    PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
27    PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
28    LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
29    NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
30    SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
31 */
32
33 #ifdef HAVE_CONFIG_H
34 #include "config.h"
35 #endif
36
37 #include "quant_bands.h"
38 #include "laplace.h"
39 #include <math.h>
40 #include "os_support.h"
41 #include "arch.h"
42 #include "mathops.h"
43 #include "stack_alloc.h"
44
45 #ifdef FIXED_POINT
46 /* Mean energy in each band quantized in Q6 */
47 static const signed char eMeans[25] = {
48       103,100, 92, 85, 81,
49        77, 72, 70, 78, 75,
50        73, 71, 78, 74, 69,
51        72, 70, 74, 76, 71,
52        60, 60, 60, 60, 60
53 };
54 #else
55 /* Mean energy in each band quantized in Q6 and converted back to float */
56 static const celt_word16 eMeans[25] = {
57       6.437500f, 6.250000f, 5.750000f, 5.312500f, 5.062500f,
58       4.812500f, 4.500000f, 4.375000f, 4.875000f, 4.687500f,
59       4.562500f, 4.437500f, 4.875000f, 4.625000f, 4.312500f,
60       4.500000f, 4.375000f, 4.625000f, 4.750000f, 4.437500f,
61       3.750000f, 3.750000f, 3.750000f, 3.750000f, 3.750000f
62 };
63 #endif
64 /* prediction coefficients: 0.9, 0.8, 0.65, 0.5 */
65 #ifdef FIXED_POINT
66 static const celt_word16 pred_coef[4] = {29440, 26112, 21248, 16384};
67 static const celt_word16 beta_coef[4] = {30147, 22282, 12124, 6554};
68 #else
69 static const celt_word16 pred_coef[4] = {29440/32768., 26112/32768., 21248/32768., 16384/32768.};
70 static const celt_word16 beta_coef[4] = {30147/32768., 22282/32768., 12124/32768., 6554/32768.};
71 #endif
72
73 /*Parameters of the Laplace-like probability models used for the coarse energy.
74   There is one pair of parameters for each frame size, prediction type
75    (inter/intra), and band number.
76   The first number of each pair is the probability of 0, and the second is the
77    decay rate, both in Q8 precision.*/
78 static const unsigned char e_prob_model[4][2][42] = {
79    /*120 sample frames.*/
80    {
81       /*Inter*/
82       {
83           72, 127,  65, 129,  66, 128,  65, 128,  64, 128,  62, 128,  64, 128,
84           64, 128,  92,  78,  92,  79,  92,  78,  90,  79, 116,  41, 115,  40,
85          114,  40, 132,  26, 132,  26, 145,  17, 161,  12, 176,  10, 177,  11
86       },
87       /*Intra*/
88       {
89           24, 179,  48, 138,  54, 135,  54, 132,  53, 134,  56, 133,  55, 132,
90           55, 132,  61, 114,  70,  96,  74,  88,  75,  88,  87,  74,  89,  66,
91           91,  67, 100,  59, 108,  50, 120,  40, 122,  37,  97,  43,  78,  50
92       }
93    },
94    /*240 sample frames.*/
95    {
96       /*Inter*/
97       {
98           83,  78,  84,  81,  88,  75,  86,  74,  87,  71,  90,  73,  93,  74,
99           93,  74, 109,  40, 114,  36, 117,  34, 117,  34, 143,  17, 145,  18,
100          146,  19, 162,  12, 165,  10, 178,   7, 189,   6, 190,   8, 177,   9
101       },
102       /*Intra*/
103       {
104           23, 178,  54, 115,  63, 102,  66,  98,  69,  99,  74,  89,  71,  91,
105           73,  91,  78,  89,  86,  80,  92,  66,  93,  64, 102,  59, 103,  60,
106          104,  60, 117,  52, 123,  44, 138,  35, 133,  31,  97,  38,  77,  45
107       }
108    },
109    /*480 sample frames.*/
110    {
111       /*Inter*/
112       {
113           61,  90,  93,  60, 105,  42, 107,  41, 110,  45, 116,  38, 113,  38,
114          112,  38, 124,  26, 132,  27, 136,  19, 140,  20, 155,  14, 159,  16,
115          158,  18, 170,  13, 177,  10, 187,   8, 192,   6, 175,   9, 159,  10
116       },
117       /*Intra*/
118       {
119           21, 178,  59, 110,  71,  86,  75,  85,  84,  83,  91,  66,  88,  73,
120           87,  72,  92,  75,  98,  72, 105,  58, 107,  54, 115,  52, 114,  55,
121          112,  56, 129,  51, 132,  40, 150,  33, 140,  29,  98,  35,  77,  42
122       }
123    },
124    /*960 sample frames.*/
125    {
126       /*Inter*/
127       {
128           42, 121,  96,  66, 108,  43, 111,  40, 117,  44, 123,  32, 120,  36,
129          119,  33, 127,  33, 134,  34, 139,  21, 147,  23, 152,  20, 158,  25,
130          154,  26, 166,  21, 173,  16, 184,  13, 184,  10, 150,  13, 139,  15
131       },
132       /*Intra*/
133       {
134           22, 178,  63, 114,  74,  82,  84,  83,  92,  82, 103,  62,  96,  72,
135           96,  67, 101,  73, 107,  72, 113,  55, 118,  52, 125,  52, 118,  52,
136          117,  55, 135,  49, 137,  39, 157,  32, 145,  29,  97,  33,  77,  40
137       }
138    }
139 };
140
141 static int intra_decision(const celt_word16 *eBands, celt_word16 *oldEBands, int start, int end, int len, int C)
142 {
143    int c, i;
144    celt_word32 dist = 0;
145    c=0; do {
146       for (i=start;i<end;i++)
147       {
148          celt_word16 d = SHR16(SUB16(eBands[i+c*len], oldEBands[i+c*len]),2);
149          dist = MAC16_16(dist, d,d);
150       }
151    } while (++c<C);
152    return SHR32(dist,2*DB_SHIFT-4) > 2*C*(end-start);
153 }
154
155 static void quant_coarse_energy_impl(const CELTMode *m, int start, int end,
156       const celt_word16 *eBands, celt_word16 *oldEBands, int budget,
157       const unsigned char *prob_model, celt_word16 *error, ec_enc *enc,
158       int _C, int LM, int intra, celt_word16 max_decay)
159 {
160    const int C = CHANNELS(_C);
161    int i, c;
162    celt_word32 prev[2] = {0,0};
163    celt_word16 coef;
164    celt_word16 beta;
165
166    ec_enc_bit_prob(enc, intra, 8192);
167    if (intra)
168    {
169       coef = 0;
170       beta = QCONST16(.15f,15);
171    } else {
172       beta = beta_coef[LM];
173       coef = pred_coef[LM];
174    }
175
176    /* Encode at a fixed coarse resolution */
177    for (i=start;i<end;i++)
178    {
179       c=0;
180       do {
181          int bits_left;
182          int qi;
183          int pi;
184          celt_word16 q;
185          celt_word16 x;
186          celt_word32 f;
187          x = eBands[i+c*m->nbEBands];
188 #ifdef FIXED_POINT
189          f = SHL32(EXTEND32(x),15) -MULT16_16(coef,oldEBands[i+c*m->nbEBands])-prev[c];
190          /* Rounding to nearest integer here is really important! */
191          qi = (f+QCONST32(.5,DB_SHIFT+15))>>(DB_SHIFT+15);
192 #else
193          f = x-coef*oldEBands[i+c*m->nbEBands]-prev[c];
194          /* Rounding to nearest integer here is really important! */
195          qi = (int)floor(.5f+f);
196 #endif
197          /* Prevent the energy from going down too quickly (e.g. for bands
198             that have just one bin) */
199          if (qi < 0 && x < oldEBands[i+c*m->nbEBands]-max_decay)
200          {
201             qi += (int)SHR16(oldEBands[i+c*m->nbEBands]-max_decay-x, DB_SHIFT);
202             if (qi > 0)
203                qi = 0;
204          }
205          /* If we don't have enough bits to encode all the energy, just assume something safe.
206             We allow slightly busting the budget here */
207          bits_left = budget-(int)ec_enc_tell(enc, 0)-3*C*(end-i);
208          if (i!=start && bits_left < 30)
209          {
210             qi = IMAX(-1,qi);
211             if (bits_left < 24)
212                qi = IMIN(1, qi);
213             if (bits_left < 16)
214                qi = IMAX(-1, qi);
215             if (bits_left<8)
216                qi = 0;
217          }
218          pi = 2*IMIN(i,20);
219          ec_laplace_encode(enc, &qi,
220                prob_model[pi]<<7, prob_model[pi+1]<<6);
221          error[i+c*m->nbEBands] = PSHR32(f,15) - SHL16(qi,DB_SHIFT);
222          q = SHL16(qi,DB_SHIFT);
223          
224          oldEBands[i+c*m->nbEBands] = PSHR32(MULT16_16(coef,oldEBands[i+c*m->nbEBands]) + prev[c] + SHL32(EXTEND32(q),15), 15);
225          prev[c] = prev[c] + SHL32(EXTEND32(q),15) - MULT16_16(beta,q);
226       } while (++c < C);
227    }
228 }
229
230 void quant_coarse_energy(const CELTMode *m, int start, int end, int effEnd,
231       const celt_word16 *eBands, celt_word16 *oldEBands, int budget,
232       celt_word16 *error, ec_enc *enc, int _C, int LM, int nbAvailableBytes,
233       int force_intra, int *delayedIntra, int two_pass)
234 {
235    const int C = CHANNELS(_C);
236    int intra;
237    celt_word16 max_decay;
238    VARDECL(celt_word16, oldEBands_intra);
239    VARDECL(celt_word16, error_intra);
240    ec_enc enc_start_state;
241    ec_byte_buffer buf_start_state;
242    SAVE_STACK;
243
244    intra = force_intra || (*delayedIntra && nbAvailableBytes > end*C);
245    if (/*shortBlocks || */intra_decision(eBands, oldEBands, start, effEnd, m->nbEBands, C))
246       *delayedIntra = 1;
247    else
248       *delayedIntra = 0;
249
250    /* Encode the global flags using a simple probability model
251       (first symbols in the stream) */
252
253 #ifdef FIXED_POINT
254       max_decay = MIN32(QCONST16(16,DB_SHIFT), SHL32(EXTEND32(nbAvailableBytes),DB_SHIFT-3));
255 #else
256    max_decay = MIN32(16.f, .125f*nbAvailableBytes);
257 #endif
258
259    enc_start_state = *enc;
260    buf_start_state = *(enc->buf);
261
262    ALLOC(oldEBands_intra, C*m->nbEBands, celt_word16);
263    ALLOC(error_intra, C*m->nbEBands, celt_word16);
264    CELT_COPY(oldEBands_intra, oldEBands, C*end);
265
266    if (two_pass || intra)
267    {
268       quant_coarse_energy_impl(m, start, end, eBands, oldEBands_intra, budget,
269             e_prob_model[LM][1], error_intra, enc, C, LM, 1, max_decay);
270    }
271
272    if (!intra)
273    {
274       ec_enc enc_intra_state;
275       ec_byte_buffer buf_intra_state;
276       int tell_intra;
277       VARDECL(unsigned char, intra_bits);
278
279       tell_intra = ec_enc_tell(enc, 3);
280
281       enc_intra_state = *enc;
282       buf_intra_state = *(enc->buf);
283
284       ALLOC(intra_bits, buf_intra_state.ptr-buf_start_state.ptr, unsigned char);
285       /* Copy bits from intra bit-stream */
286       CELT_COPY(intra_bits, buf_start_state.ptr, buf_intra_state.ptr-buf_start_state.ptr);
287
288       *enc = enc_start_state;
289       *(enc->buf) = buf_start_state;
290
291       quant_coarse_energy_impl(m, start, end, eBands, oldEBands, budget,
292             e_prob_model[LM][intra], error, enc, C, LM, 0, max_decay);
293
294       if (two_pass && ec_enc_tell(enc, 3) > tell_intra)
295       {
296          *enc = enc_intra_state;
297          *(enc->buf) = buf_intra_state;
298          /* Copy bits from to bit-stream */
299          CELT_COPY(buf_start_state.ptr, intra_bits, buf_intra_state.ptr-buf_start_state.ptr);
300          CELT_COPY(oldEBands, oldEBands_intra, C*end);
301          CELT_COPY(error, error_intra, C*end);
302       }
303    } else {
304       CELT_COPY(oldEBands, oldEBands_intra, C*end);
305       CELT_COPY(error, error_intra, C*end);
306    }
307    RESTORE_STACK;
308 }
309
310 void quant_fine_energy(const CELTMode *m, int start, int end, celt_ener *eBands, celt_word16 *oldEBands, celt_word16 *error, int *fine_quant, ec_enc *enc, int _C)
311 {
312    int i, c;
313    const int C = CHANNELS(_C);
314
315    /* Encode finer resolution */
316    for (i=start;i<end;i++)
317    {
318       celt_int16 frac = 1<<fine_quant[i];
319       if (fine_quant[i] <= 0)
320          continue;
321       c=0;
322       do {
323          int q2;
324          celt_word16 offset;
325 #ifdef FIXED_POINT
326          /* Has to be without rounding */
327          q2 = (error[i+c*m->nbEBands]+QCONST16(.5f,DB_SHIFT))>>(DB_SHIFT-fine_quant[i]);
328 #else
329          q2 = (int)floor((error[i+c*m->nbEBands]+.5f)*frac);
330 #endif
331          if (q2 > frac-1)
332             q2 = frac-1;
333          if (q2<0)
334             q2 = 0;
335          ec_enc_bits(enc, q2, fine_quant[i]);
336 #ifdef FIXED_POINT
337          offset = SUB16(SHR32(SHL32(EXTEND32(q2),DB_SHIFT)+QCONST16(.5,DB_SHIFT),fine_quant[i]),QCONST16(.5f,DB_SHIFT));
338 #else
339          offset = (q2+.5f)*(1<<(14-fine_quant[i]))*(1.f/16384) - .5f;
340 #endif
341          oldEBands[i+c*m->nbEBands] += offset;
342          error[i+c*m->nbEBands] -= offset;
343          /*printf ("%f ", error[i] - offset);*/
344       } while (++c < C);
345    }
346 }
347
348 void quant_energy_finalise(const CELTMode *m, int start, int end, celt_ener *eBands, celt_word16 *oldEBands, celt_word16 *error, int *fine_quant, int *fine_priority, int bits_left, ec_enc *enc, int _C)
349 {
350    int i, prio, c;
351    const int C = CHANNELS(_C);
352
353    /* Use up the remaining bits */
354    for (prio=0;prio<2;prio++)
355    {
356       for (i=start;i<end && bits_left>=C ;i++)
357       {
358          if (fine_quant[i] >= 7 || fine_priority[i]!=prio)
359             continue;
360          c=0;
361          do {
362             int q2;
363             celt_word16 offset;
364             q2 = error[i+c*m->nbEBands]<0 ? 0 : 1;
365             ec_enc_bits(enc, q2, 1);
366 #ifdef FIXED_POINT
367             offset = SHR16(SHL16(q2,DB_SHIFT)-QCONST16(.5,DB_SHIFT),fine_quant[i]+1);
368 #else
369             offset = (q2-.5f)*(1<<(14-fine_quant[i]-1))*(1.f/16384);
370 #endif
371             oldEBands[i+c*m->nbEBands] += offset;
372             bits_left--;
373          } while (++c < C);
374       }
375    }
376 }
377
378 void unquant_coarse_energy(const CELTMode *m, int start, int end, celt_ener *eBands, celt_word16 *oldEBands, int intra, ec_dec *dec, int _C, int LM)
379 {
380    const unsigned char *prob_model = e_prob_model[LM][intra];
381    int i, c;
382    celt_word32 prev[2] = {0, 0};
383    celt_word16 coef;
384    celt_word16 beta;
385    const int C = CHANNELS(_C);
386
387
388    if (intra)
389    {
390       coef = 0;
391       beta = QCONST16(.15f,15);
392    } else {
393       beta = beta_coef[LM];
394       coef = pred_coef[LM];
395    }
396
397    /* Decode at a fixed coarse resolution */
398    for (i=start;i<end;i++)
399    {
400       c=0;
401       do {
402          int qi;
403          int pi;
404          celt_word16 q;
405          pi = 2*IMIN(i,20);
406          qi = ec_laplace_decode(dec,
407                prob_model[pi]<<7, prob_model[pi+1]<<6);
408          q = SHL16(qi,DB_SHIFT);
409
410          oldEBands[i+c*m->nbEBands] = PSHR32(MULT16_16(coef,oldEBands[i+c*m->nbEBands]) + prev[c] + SHL32(EXTEND32(q),15), 15);
411          prev[c] = prev[c] + SHL32(EXTEND32(q),15) - MULT16_16(beta,q);
412       } while (++c < C);
413    }
414 }
415
416 void unquant_fine_energy(const CELTMode *m, int start, int end, celt_ener *eBands, celt_word16 *oldEBands, int *fine_quant, ec_dec *dec, int _C)
417 {
418    int i, c;
419    const int C = CHANNELS(_C);
420    /* Decode finer resolution */
421    for (i=start;i<end;i++)
422    {
423       if (fine_quant[i] <= 0)
424          continue;
425       c=0; 
426       do {
427          int q2;
428          celt_word16 offset;
429          q2 = ec_dec_bits(dec, fine_quant[i]);
430 #ifdef FIXED_POINT
431          offset = SUB16(SHR32(SHL32(EXTEND32(q2),DB_SHIFT)+QCONST16(.5,DB_SHIFT),fine_quant[i]),QCONST16(.5f,DB_SHIFT));
432 #else
433          offset = (q2+.5f)*(1<<(14-fine_quant[i]))*(1.f/16384) - .5f;
434 #endif
435          oldEBands[i+c*m->nbEBands] += offset;
436       } while (++c < C);
437    }
438 }
439
440 void unquant_energy_finalise(const CELTMode *m, int start, int end, celt_ener *eBands, celt_word16 *oldEBands, int *fine_quant,  int *fine_priority, int bits_left, ec_dec *dec, int _C)
441 {
442    int i, prio, c;
443    const int C = CHANNELS(_C);
444
445    /* Use up the remaining bits */
446    for (prio=0;prio<2;prio++)
447    {
448       for (i=start;i<end && bits_left>=C ;i++)
449       {
450          if (fine_quant[i] >= 7 || fine_priority[i]!=prio)
451             continue;
452          c=0;
453          do {
454             int q2;
455             celt_word16 offset;
456             q2 = ec_dec_bits(dec, 1);
457 #ifdef FIXED_POINT
458             offset = SHR16(SHL16(q2,DB_SHIFT)-QCONST16(.5,DB_SHIFT),fine_quant[i]+1);
459 #else
460             offset = (q2-.5f)*(1<<(14-fine_quant[i]-1))*(1.f/16384);
461 #endif
462             oldEBands[i+c*m->nbEBands] += offset;
463             bits_left--;
464          } while (++c < C);
465       }
466    }
467 }
468
469 void log2Amp(const CELTMode *m, int start, int end,
470       celt_ener *eBands, celt_word16 *oldEBands, int _C)
471 {
472    int c, i;
473    const int C = CHANNELS(_C);
474    c=0;
475    do {
476       for (i=start;i<m->nbEBands;i++)
477       {
478          celt_word16 lg = oldEBands[i+c*m->nbEBands]
479                         + SHL16((celt_word16)eMeans[i],6);
480          eBands[i+c*m->nbEBands] = PSHR32(celt_exp2(SHL16(lg,11-DB_SHIFT)),4);
481          if (oldEBands[i+c*m->nbEBands] < -QCONST16(14.f,DB_SHIFT))
482             oldEBands[i+c*m->nbEBands] = -QCONST16(14.f,DB_SHIFT);
483       }
484    } while (++c < C);
485 }
486
487 void amp2Log2(const CELTMode *m, int effEnd, int end,
488       celt_ener *bandE, celt_word16 *bandLogE, int _C)
489 {
490    int c, i;
491    const int C = CHANNELS(_C);
492    c=0;
493    do {
494       for (i=0;i<effEnd;i++)
495          bandLogE[i+c*m->nbEBands] =
496                celt_log2(MAX32(QCONST32(.001f,14),SHL32(bandE[i+c*m->nbEBands],2)))
497                - SHL16((celt_word16)eMeans[i],6);
498       for (i=effEnd;i<end;i++)
499          bandLogE[c*m->nbEBands+i] = -QCONST16(14.f,DB_SHIFT);
500    } while (++c < C);
501 }