Fixed a long-standing rare mismatch
[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 #define E_MEANS_SIZE (5)
46
47 const celt_word16 eMeans[E_MEANS_SIZE] = {QCONST16(7.5f,DB_SHIFT), -QCONST16(1.33f,DB_SHIFT), -QCONST16(2.f,DB_SHIFT), -QCONST16(0.42f,DB_SHIFT), QCONST16(0.17f,DB_SHIFT)};
48
49 /* FIXME: Implement for stereo */
50 int intra_decision(celt_word16 *eBands, celt_word16 *oldEBands, int len)
51 {
52    int i;
53    celt_word32 dist = 0;
54    for (i=0;i<len;i++)
55    {
56       celt_word16 d = SUB16(eBands[i], oldEBands[i]);
57       dist = MAC16_16(dist, d,d);
58    }
59    return SHR32(dist,2*DB_SHIFT) > 2*len;
60 }
61
62 int *quant_prob_alloc(const CELTMode *m)
63 {
64    int i;
65    int *prob;
66    prob = celt_alloc(4*m->nbEBands*sizeof(int));
67    if (prob==NULL)
68      return NULL;
69    for (i=0;i<m->nbEBands;i++)
70    {
71       prob[2*i] = 6000-i*200;
72       prob[2*i+1] = ec_laplace_get_start_freq(prob[2*i]);
73    }
74    for (i=0;i<m->nbEBands;i++)
75    {
76       prob[2*m->nbEBands+2*i] = 9000-i*240;
77       prob[2*m->nbEBands+2*i+1] = ec_laplace_get_start_freq(prob[2*m->nbEBands+2*i]);
78    }
79    return prob;
80 }
81
82 void quant_prob_free(int *freq)
83 {
84    celt_free(freq);
85 }
86
87 unsigned quant_coarse_energy(const CELTMode *m, int start, celt_word16 *eBands, celt_word16 *oldEBands, int budget, int intra, int *prob, celt_word16 *error, ec_enc *enc, int _C)
88 {
89    int i, c;
90    unsigned bits_used = 0;
91    celt_word32 prev[2] = {0,0};
92    celt_word16 coef = m->ePredCoef;
93    celt_word16 beta;
94    const int C = CHANNELS(_C);
95
96    if (intra)
97    {
98       coef = 0;
99       prob += 2*m->nbEBands;
100    }
101    /* The .8 is a heuristic */
102    beta = MULT16_16_P15(QCONST16(.8f,15),coef);
103
104    /* Encode at a fixed coarse resolution */
105    for (i=start;i<m->nbEBands;i++)
106    {
107       c=0;
108       do {
109          int qi;
110          celt_word16 q;
111          celt_word16 x;
112          celt_word32 f;
113          celt_word32 mean =  (i-start < E_MEANS_SIZE) ? SUB32(SHL32(EXTEND32(eMeans[i-start]),15), MULT16_16(coef,eMeans[i-start])) : 0;
114          x = eBands[i+c*m->nbEBands];
115 #ifdef FIXED_POINT
116          f = SHL32(EXTEND32(x),15)-mean -MULT16_16(coef,oldEBands[i+c*m->nbEBands])-prev[c];
117          /* Rounding to nearest integer here is really important! */
118          qi = (f+QCONST32(.5,DB_SHIFT+15))>>(DB_SHIFT+15);
119 #else
120          f = x-mean-coef*oldEBands[i+c*m->nbEBands]-prev[c];
121          /* Rounding to nearest integer here is really important! */
122          qi = (int)floor(.5f+f);
123 #endif
124          /* If we don't have enough bits to encode all the energy, just assume something safe.
125             We allow slightly busting the budget here */
126          bits_used=ec_enc_tell(enc, 0);
127          if (bits_used > budget)
128          {
129             qi = -1;
130             error[i+c*m->nbEBands] = QCONST16(.5f,DB_SHIFT);
131          } else {
132             ec_laplace_encode_start(enc, &qi, prob[2*i], prob[2*i+1]);
133             error[i+c*m->nbEBands] = PSHR32(f,15) - SHL16(qi,DB_SHIFT);
134          }
135          q = SHL16(qi,DB_SHIFT);
136          
137          oldEBands[i+c*m->nbEBands] = PSHR32(MULT16_16(coef,oldEBands[i+c*m->nbEBands]) + mean + prev[c] + SHL32(EXTEND32(q),15), 15);
138          prev[c] = mean + prev[c] + SHL32(EXTEND32(q),15) - MULT16_16(beta,q);
139       } while (++c < C);
140    }
141    return bits_used;
142 }
143
144 void quant_fine_energy(const CELTMode *m, int start, celt_ener *eBands, celt_word16 *oldEBands, celt_word16 *error, int *fine_quant, ec_enc *enc, int _C)
145 {
146    int i, c;
147    const int C = CHANNELS(_C);
148
149    /* Encode finer resolution */
150    for (i=start;i<m->nbEBands;i++)
151    {
152       celt_int16 frac = 1<<fine_quant[i];
153       if (fine_quant[i] <= 0)
154          continue;
155       c=0;
156       do {
157          int q2;
158          celt_word16 offset;
159 #ifdef FIXED_POINT
160          /* Has to be without rounding */
161          q2 = (error[i+c*m->nbEBands]+QCONST16(.5f,DB_SHIFT))>>(DB_SHIFT-fine_quant[i]);
162 #else
163          q2 = (int)floor((error[i+c*m->nbEBands]+.5f)*frac);
164 #endif
165          if (q2 > frac-1)
166             q2 = frac-1;
167          if (q2<0)
168             q2 = 0;
169          ec_enc_bits(enc, q2, fine_quant[i]);
170 #ifdef FIXED_POINT
171          offset = SUB16(SHR16(SHL16(q2,DB_SHIFT)+QCONST16(.5,DB_SHIFT),fine_quant[i]),QCONST16(.5f,DB_SHIFT));
172 #else
173          offset = (q2+.5f)*(1<<(14-fine_quant[i]))*(1.f/16384) - .5f;
174 #endif
175          oldEBands[i+c*m->nbEBands] += offset;
176          error[i+c*m->nbEBands] -= offset;
177          /*printf ("%f ", error[i] - offset);*/
178       } while (++c < C);
179    }
180    for (i=start;i<C*m->nbEBands;i++)
181       eBands[i] = log2Amp(oldEBands[i]);
182 }
183
184 void quant_energy_finalise(const CELTMode *m, int start, celt_ener *eBands, celt_word16 *oldEBands, celt_word16 *error, int *fine_quant, int *fine_priority, int bits_left, ec_enc *enc, int _C)
185 {
186    int i, prio, c;
187    const int C = CHANNELS(_C);
188
189    /* Use up the remaining bits */
190    for (prio=0;prio<2;prio++)
191    {
192       for (i=start;i<m->nbEBands && bits_left>=C ;i++)
193       {
194          if (fine_quant[i] >= 7 || fine_priority[i]!=prio)
195             continue;
196          c=0;
197          do {
198             int q2;
199             celt_word16 offset;
200             q2 = error[i+c*m->nbEBands]<0 ? 0 : 1;
201             ec_enc_bits(enc, q2, 1);
202 #ifdef FIXED_POINT
203             offset = SHR16(SHL16(q2,DB_SHIFT)-QCONST16(.5,DB_SHIFT),fine_quant[i]+1);
204 #else
205             offset = (q2-.5f)*(1<<(14-fine_quant[i]-1))*(1.f/16384);
206 #endif
207             oldEBands[i+c*m->nbEBands] += offset;
208             eBands[i+c*m->nbEBands] = log2Amp(oldEBands[i+c*m->nbEBands]);
209             bits_left--;
210          } while (++c < C);
211       }
212    }
213    for (i=start;i<C*m->nbEBands;i++)
214    {
215       if (oldEBands[i] < -QCONST16(7.f,DB_SHIFT))
216          oldEBands[i] = -QCONST16(7.f,DB_SHIFT);
217    }
218 }
219
220 void unquant_coarse_energy(const CELTMode *m, int start, celt_ener *eBands, celt_word16 *oldEBands, int budget, int intra, int *prob, ec_dec *dec, int _C)
221 {
222    int i, c;
223    celt_word32 prev[2] = {0, 0};
224    celt_word16 coef = m->ePredCoef;
225    celt_word16 beta;
226    const int C = CHANNELS(_C);
227
228    if (intra)
229    {
230       coef = 0;
231       prob += 2*m->nbEBands;
232    }
233    /* The .8 is a heuristic */
234    beta = MULT16_16_P15(QCONST16(.8f,15),coef);
235
236    /* Decode at a fixed coarse resolution */
237    for (i=start;i<m->nbEBands;i++)
238    {
239       c=0;
240       do {
241          int qi;
242          celt_word16 q;
243          celt_word32 mean =  (i-start < E_MEANS_SIZE) ? SUB32(SHL32(EXTEND32(eMeans[i-start]),15), MULT16_16(coef,eMeans[i-start])) : 0;
244          /* If we didn't have enough bits to encode all the energy, just assume something safe.
245             We allow slightly busting the budget here */
246          if (ec_dec_tell(dec, 0) > budget)
247             qi = -1;
248          else
249             qi = ec_laplace_decode_start(dec, prob[2*i], prob[2*i+1]);
250          q = SHL16(qi,DB_SHIFT);
251
252          oldEBands[i+c*m->nbEBands] = PSHR32(MULT16_16(coef,oldEBands[i+c*m->nbEBands]) + mean + prev[c] + SHL32(EXTEND32(q),15), 15);
253          prev[c] = mean + prev[c] + SHL32(EXTEND32(q),15) - MULT16_16(beta,q);
254       } while (++c < C);
255    }
256 }
257
258 void unquant_fine_energy(const CELTMode *m, int start, celt_ener *eBands, celt_word16 *oldEBands, int *fine_quant, ec_dec *dec, int _C)
259 {
260    int i, c;
261    const int C = CHANNELS(_C);
262    /* Decode finer resolution */
263    for (i=start;i<m->nbEBands;i++)
264    {
265       if (fine_quant[i] <= 0)
266          continue;
267       c=0; 
268       do {
269          int q2;
270          celt_word16 offset;
271          q2 = ec_dec_bits(dec, fine_quant[i]);
272 #ifdef FIXED_POINT
273          offset = SUB16(SHR16(SHL16(q2,DB_SHIFT)+QCONST16(.5,DB_SHIFT),fine_quant[i]),QCONST16(.5f,DB_SHIFT));
274 #else
275          offset = (q2+.5f)*(1<<(14-fine_quant[i]))*(1.f/16384) - .5f;
276 #endif
277          oldEBands[i+c*m->nbEBands] += offset;
278       } while (++c < C);
279    }
280    for (i=start;i<C*m->nbEBands;i++)
281       eBands[i] = log2Amp(oldEBands[i]);
282 }
283
284 void unquant_energy_finalise(const CELTMode *m, int start, celt_ener *eBands, celt_word16 *oldEBands, int *fine_quant,  int *fine_priority, int bits_left, ec_dec *dec, int _C)
285 {
286    int i, prio, c;
287    const int C = CHANNELS(_C);
288
289    /* Use up the remaining bits */
290    for (prio=0;prio<2;prio++)
291    {
292       for (i=start;i<m->nbEBands && bits_left>=C ;i++)
293       {
294          if (fine_quant[i] >= 7 || fine_priority[i]!=prio)
295             continue;
296          c=0;
297          do {
298             int q2;
299             celt_word16 offset;
300             q2 = ec_dec_bits(dec, 1);
301 #ifdef FIXED_POINT
302             offset = SHR16(SHL16(q2,DB_SHIFT)-QCONST16(.5,DB_SHIFT),fine_quant[i]+1);
303 #else
304             offset = (q2-.5f)*(1<<(14-fine_quant[i]-1))*(1.f/16384);
305 #endif
306             oldEBands[i+c*m->nbEBands] += offset;
307             eBands[i+c*m->nbEBands] = log2Amp(oldEBands[i+c*m->nbEBands]);
308             bits_left--;
309          } while (++c < C);
310       }
311    }
312    for (i=start;i<C*m->nbEBands;i++)
313    {
314       if (oldEBands[i] < -QCONST16(7.f,DB_SHIFT))
315          oldEBands[i] = -QCONST16(7.f,DB_SHIFT);
316    }
317 }