0e39d37b62224a270ae55c9f998d74d03cd3b557
[opus.git] / libcelt / quant_bands.c
1 /* (C) 2007-2008 Jean-Marc Valin, CSIRO
2 */
3 /*
4    Redistribution and use in source and binary forms, with or without
5    modification, are permitted provided that the following conditions
6    are met:
7    
8    - Redistributions of source code must retain the above copyright
9    notice, this list of conditions and the following disclaimer.
10    
11    - Redistributions in binary form must reproduce the above copyright
12    notice, this list of conditions and the following disclaimer in the
13    documentation and/or other materials provided with the distribution.
14    
15    - Neither the name of the Xiph.org Foundation nor the names of its
16    contributors may be used to endorse or promote products derived from
17    this software without specific prior written permission.
18    
19    THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
20    ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
21    LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
22    A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR
23    CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
24    EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
25    PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
26    PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
27    LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
28    NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
29    SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
30 */
31
32 #ifdef HAVE_CONFIG_H
33 #include "config.h"
34 #endif
35
36 #include "quant_bands.h"
37 #include "laplace.h"
38 #include <math.h>
39 #include "os_support.h"
40 #include "arch.h"
41 #include "mathops.h"
42 #include "stack_alloc.h"
43
44 #ifdef FIXED_POINT
45 const celt_word16_t eMeans[24] = {11520, -2048, -3072, -640, 256, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
46 #else
47 const celt_word16_t eMeans[24] = {45.f, -8.f, -12.f, -2.5f, 1.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f};
48 #endif
49
50
51 #ifdef FIXED_POINT
52 static inline celt_ener_t dB2Amp(celt_ener_t dB)
53 {
54    celt_ener_t amp;
55    if (dB>24659)
56       dB=24659;
57    amp = PSHR32(celt_exp2(MULT16_16_Q14(21771,dB)),2)-QCONST16(.3f, 14);
58    if (amp < 0)
59       amp = 0;
60    return PSHR32(amp,2);
61 }
62
63 #define DBofTWO 24661
64 static inline celt_word16_t amp2dB(celt_ener_t amp)
65 {
66    /* equivalent to return 6.0207*log2(.3+amp) */
67    return ROUND16(MULT16_16(24661,celt_log2(ADD32(QCONST32(.3f,14),SHL32(amp,2)))),12);
68    /* return DB_SCALING*20*log10(.3+ENER_SCALING_1*amp); */
69 }
70 #else
71 static inline celt_ener_t dB2Amp(celt_ener_t dB)
72 {
73    celt_ener_t amp;
74    /*amp = pow(10, .05*dB)-.3;*/
75    amp = exp(0.115129f*dB)-.3f;
76    if (amp < 0)
77       amp = 0;
78    return amp;
79 }
80 static inline celt_word16_t amp2dB(celt_ener_t amp)
81 {
82    /*return 20*log10(.3+amp);*/
83    return 8.68589f*log(.3f+amp);
84 }
85 #endif
86
87 int intra_decision(celt_ener_t *eBands, celt_word16_t *oldEBands, int len)
88 {
89    int i;
90    celt_word32_t dist = 0;
91    for (i=0;i<len;i++)
92    {
93       celt_word16_t d = SUB16(amp2dB(eBands[i]), oldEBands[i]);
94       dist = MAC16_16(dist, d,d);
95    }
96    return SHR32(dist,16) > 64*len;
97 }
98
99 static const celt_word16_t base_resolution = QCONST16(6.f,8);
100 static const celt_word16_t base_resolution_1 = QCONST16(0.1666667f,15);
101
102 int *quant_prob_alloc(const CELTMode *m)
103 {
104    int i;
105    int *prob;
106    prob = celt_alloc(4*m->nbEBands*sizeof(int));
107    for (i=0;i<m->nbEBands;i++)
108    {
109       prob[2*i] = 6000-i*200;
110       prob[2*i+1] = ec_laplace_get_start_freq(prob[2*i]);
111    }
112    for (i=0;i<m->nbEBands;i++)
113    {
114       prob[2*m->nbEBands+2*i] = 9000-i*240;
115       prob[2*m->nbEBands+2*i+1] = ec_laplace_get_start_freq(prob[2*m->nbEBands+2*i]);
116    }
117    return prob;
118 }
119
120 void quant_prob_free(int *freq)
121 {
122    celt_free(freq);
123 }
124
125 static void quant_coarse_energy_mono(const CELTMode *m, celt_ener_t *eBands, celt_word16_t *oldEBands, unsigned budget, int intra, int *prob, celt_word16_t *error, ec_enc *enc)
126 {
127    int i;
128    unsigned bits;
129    celt_word16_t prev = 0;
130    celt_word16_t coef = m->ePredCoef;
131    celt_word16_t beta;
132    
133    if (intra)
134    {
135       coef = 0;
136       prob += 2*m->nbEBands;
137    }
138    /* The .8 is a heuristic */
139    beta = MULT16_16_Q15(QCONST16(.8f,15),coef);
140    
141    bits = ec_enc_tell(enc, 0);
142    /* Encode at a fixed coarse resolution */
143    for (i=0;i<m->nbEBands;i++)
144    {
145       int qi;
146       celt_word16_t q;   /* dB */
147       celt_word16_t x;   /* dB */
148       celt_word16_t f;   /* Q8 */
149       celt_word16_t mean = MULT16_16_Q15(Q15ONE-coef,eMeans[i]);
150       x = amp2dB(eBands[i]);
151 #ifdef FIXED_POINT
152       f = MULT16_16_Q15(x-mean-MULT16_16_Q15(coef,oldEBands[i])-prev,base_resolution_1);
153       /* Rounding to nearest integer here is really important! */
154       qi = (f+128)>>8;
155 #else
156       f = (x-mean-coef*oldEBands[i]-prev)*base_resolution_1;
157       /* Rounding to nearest integer here is really important! */
158       qi = (int)floor(.5+f);
159 #endif
160       /* If we don't have enough bits to encode all the energy, just assume something safe.
161          We allow slightly busting the budget here */
162       if (ec_enc_tell(enc, 0) - bits > budget)
163       {
164          qi = -1;
165          error[i] = 128;
166       } else {
167          ec_laplace_encode_start(enc, &qi, prob[2*i], prob[2*i+1]);
168          error[i] = f - SHL16(qi,8);
169       }
170       q = qi*base_resolution;
171       
172       oldEBands[i] = mean+MULT16_16_Q15(coef,oldEBands[i])+prev+q;
173       if (oldEBands[i] < -QCONST16(12.f,8))
174          oldEBands[i] = -QCONST16(12.f,8);
175       prev = mean+prev+MULT16_16_Q15(Q15ONE-beta,q);
176    }
177 }
178
179 static void quant_fine_energy_mono(const CELTMode *m, celt_ener_t *eBands, celt_word16_t *oldEBands, celt_word16_t *error, int *fine_quant, ec_enc *enc)
180 {
181    int i;
182    /* Encode finer resolution */
183    for (i=0;i<m->nbEBands;i++)
184    {
185       int q2;
186       celt_int16_t frac = 1<<fine_quant[i];
187       celt_word16_t offset = (error[i]+QCONST16(.5f,8))*frac;
188       if (fine_quant[i] <= 0)
189          continue;
190 #ifdef FIXED_POINT
191       /* Has to be without rounding */
192       q2 = offset>>8;
193 #else
194       q2 = (int)floor(offset);
195 #endif
196       if (q2 > frac-1)
197          q2 = frac-1;
198       ec_enc_bits(enc, q2, fine_quant[i]);
199 #ifdef FIXED_POINT
200       offset = SUB16(SHR16(SHL16(q2,8)+QCONST16(.5,8),fine_quant[i]),QCONST16(.5f,8));
201 #else
202       offset = (q2+.5f)*(1<<(14-fine_quant[i]))*(1.f/16384) - .5f;
203 #endif
204       oldEBands[i] += PSHR32(MULT16_16(DB_SCALING*6,offset),8);
205       /*printf ("%f ", error[i] - offset);*/
206    }
207    for (i=0;i<m->nbEBands;i++)
208    {
209       eBands[i] = dB2Amp(oldEBands[i]);
210    }
211    /*printf ("%d\n", ec_enc_tell(enc, 0)-9);*/
212
213    /*printf ("\n");*/
214 }
215
216 static void unquant_coarse_energy_mono(const CELTMode *m, celt_ener_t *eBands, celt_word16_t *oldEBands, unsigned budget, int intra, int *prob, ec_dec *dec)
217 {
218    int i;
219    unsigned bits;
220    celt_word16_t prev = 0;
221    celt_word16_t coef = m->ePredCoef;
222    celt_word16_t beta;
223    
224    if (intra)
225    {
226       coef = 0;
227       prob += 2*m->nbEBands;
228    }
229    /* The .8 is a heuristic */
230    beta = MULT16_16_Q15(QCONST16(.8f,15),coef);
231    
232    bits = ec_dec_tell(dec, 0);
233    /* Decode at a fixed coarse resolution */
234    for (i=0;i<m->nbEBands;i++)
235    {
236       int qi;
237       celt_word16_t q;
238       celt_word16_t mean = MULT16_16_Q15(Q15ONE-coef,eMeans[i]);
239       /* If we didn't have enough bits to encode all the energy, just assume something safe.
240          We allow slightly busting the budget here */
241       if (ec_dec_tell(dec, 0) - bits > budget)
242          qi = -1;
243       else
244          qi = ec_laplace_decode_start(dec, prob[2*i], prob[2*i+1]);
245       q = qi*base_resolution;
246       
247       oldEBands[i] = mean+MULT16_16_Q15(coef,oldEBands[i])+prev+q;
248       if (oldEBands[i] < -QCONST16(12.f,8))
249          oldEBands[i] = -QCONST16(12.f,8);
250       
251       prev = mean+prev+MULT16_16_Q15(Q15ONE-beta,q);
252    }
253 }
254
255 static void unquant_fine_energy_mono(const CELTMode *m, celt_ener_t *eBands, celt_word16_t *oldEBands, int *fine_quant, ec_dec *dec)
256 {
257    int i;
258    /* Decode finer resolution */
259    for (i=0;i<m->nbEBands;i++)
260    {
261       int q2;
262       celt_word16_t offset;
263       if (fine_quant[i] <= 0)
264          continue;
265       q2 = ec_dec_bits(dec, fine_quant[i]);
266 #ifdef FIXED_POINT
267       offset = SUB16(SHR16(SHL16(q2,8)+QCONST16(.5,8),fine_quant[i]),QCONST16(.5f,8));
268 #else
269       offset = (q2+.5f)*(1<<(14-fine_quant[i]))*(1.f/16384) - .5f;
270 #endif
271       oldEBands[i] += PSHR32(MULT16_16(DB_SCALING*6,offset),8);
272    }
273    for (i=0;i<m->nbEBands;i++)
274    {
275       eBands[i] = dB2Amp(oldEBands[i]);
276    }
277    /*printf ("\n");*/
278 }
279
280
281
282 void quant_coarse_energy(const CELTMode *m, celt_ener_t *eBands, celt_word16_t *oldEBands, int budget, int intra, int *prob, celt_word16_t *error, ec_enc *enc)
283 {
284    int C;
285    C = m->nbChannels;
286
287    if (C==1)
288    {
289       quant_coarse_energy_mono(m, eBands, oldEBands, budget, intra, prob, error, enc);
290    } else {
291       int c;
292       for (c=0;c<C;c++)
293       {
294          int i;
295          VARDECL(celt_ener_t, E);
296          SAVE_STACK;
297          ALLOC(E, m->nbEBands, celt_ener_t);
298          for (i=0;i<m->nbEBands;i++)
299             E[i] = eBands[C*i+c];
300          quant_coarse_energy_mono(m, E, oldEBands+c*m->nbEBands, budget/C, intra, prob, error+c*m->nbEBands, enc);
301          RESTORE_STACK;
302       }
303    }
304 }
305
306 void quant_fine_energy(const CELTMode *m, celt_ener_t *eBands, celt_word16_t *oldEBands, celt_word16_t *error, int *fine_quant, ec_enc *enc)
307 {
308    int C;
309    C = m->nbChannels;
310
311    if (C==1)
312    {
313       quant_fine_energy_mono(m, eBands, oldEBands, error, fine_quant, enc);
314
315    } else {
316       int c;
317       VARDECL(celt_ener_t, E);
318       ALLOC(E, m->nbEBands, celt_ener_t);
319       for (c=0;c<C;c++)
320       {
321          int i;
322          SAVE_STACK;
323          quant_fine_energy_mono(m, E, oldEBands+c*m->nbEBands, error+c*m->nbEBands, fine_quant, enc);
324          for (i=0;i<m->nbEBands;i++)
325             eBands[C*i+c] = E[i];
326          RESTORE_STACK;
327       }
328    }
329 }
330
331
332 void unquant_coarse_energy(const CELTMode *m, celt_ener_t *eBands, celt_word16_t *oldEBands, int budget, int intra, int *prob, ec_dec *dec)
333 {
334    int C;   
335
336    C = m->nbChannels;
337    if (C==1)
338    {
339       unquant_coarse_energy_mono(m, eBands, oldEBands, budget, intra, prob, dec);
340    }
341    else {
342       int c;
343       VARDECL(celt_ener_t, E);
344       SAVE_STACK;
345       ALLOC(E, m->nbEBands, celt_ener_t);
346       for (c=0;c<C;c++)
347       {
348          unquant_coarse_energy_mono(m, E, oldEBands+c*m->nbEBands, budget/C, intra, prob, dec);
349       }
350       RESTORE_STACK;
351    }
352 }
353
354 void unquant_fine_energy(const CELTMode *m, celt_ener_t *eBands, celt_word16_t *oldEBands, int *fine_quant, ec_dec *dec)
355 {
356    int C;   
357
358    C = m->nbChannels;
359
360    if (C==1)
361    {
362       unquant_fine_energy_mono(m, eBands, oldEBands, fine_quant, dec);
363    }
364    else {
365       int c;
366       VARDECL(celt_ener_t, E);
367       SAVE_STACK;
368       ALLOC(E, m->nbEBands, celt_ener_t);
369       for (c=0;c<C;c++)
370       {
371          int i;
372          unquant_fine_energy_mono(m, E, oldEBands+c*m->nbEBands, fine_quant, dec);
373          for (i=0;i<m->nbEBands;i++)
374             eBands[C*i+c] = E[i];
375       }
376       RESTORE_STACK;
377    }
378 }