Miscellaneous comment, copyright notice, readme updates.
[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 static const celt_word16_t base_resolution = QCONST16(6.f,8);
88 static const celt_word16_t base_resolution_1 = QCONST16(0.1666667f,15);
89
90 int *quant_prob_alloc(const CELTMode *m)
91 {
92    int i;
93    int *prob;
94    prob = celt_alloc(2*m->nbEBands*sizeof(int));
95    for (i=0;i<m->nbEBands;i++)
96    {
97       prob[2*i] = 6000-i*200;
98       prob[2*i+1] = ec_laplace_get_start_freq(prob[2*i]);
99    }
100    return prob;
101 }
102
103 void quant_prob_free(int *freq)
104 {
105    celt_free(freq);
106 }
107
108 static void quant_coarse_energy_mono(const CELTMode *m, celt_ener_t *eBands, celt_word16_t *oldEBands, unsigned budget, int *prob, celt_word16_t *error, ec_enc *enc)
109 {
110    int i;
111    unsigned bits;
112    celt_word16_t prev = 0;
113    celt_word16_t coef = m->ePredCoef;
114    celt_word16_t beta;
115    /* The .8 is a heuristic */
116    beta = MULT16_16_Q15(QCONST16(.8f,15),coef);
117    
118    bits = ec_enc_tell(enc, 0);
119    /* Encode at a fixed coarse resolution */
120    for (i=0;i<m->nbEBands;i++)
121    {
122       int qi;
123       celt_word16_t q;   /* dB */
124       celt_word16_t x;   /* dB */
125       celt_word16_t f;   /* Q8 */
126       celt_word16_t mean = MULT16_16_Q15(Q15ONE-coef,eMeans[i]);
127       x = amp2dB(eBands[i]);
128 #ifdef FIXED_POINT
129       f = MULT16_16_Q15(x-mean-MULT16_16_Q15(coef,oldEBands[i])-prev,base_resolution_1);
130       /* Rounding to nearest integer here is really important! */
131       qi = (f+128)>>8;
132 #else
133       f = (x-mean-coef*oldEBands[i]-prev)*base_resolution_1;
134       /* Rounding to nearest integer here is really important! */
135       qi = (int)floor(.5+f);
136 #endif
137       /* If we don't have enough bits to encode all the energy, just assume something safe.
138          We allow slightly busting the budget here */
139       if (ec_enc_tell(enc, 0) - bits > budget)
140       {
141          qi = -1;
142          error[i] = 128;
143       } else {
144          ec_laplace_encode_start(enc, &qi, prob[2*i], prob[2*i+1]);
145          error[i] = f - SHL16(qi,8);
146       }
147       q = qi*base_resolution;
148       
149       oldEBands[i] = mean+MULT16_16_Q15(coef,oldEBands[i])+prev+q;
150       if (oldEBands[i] < -QCONST16(12.f,8))
151          oldEBands[i] = -QCONST16(12.f,8);
152       prev = mean+prev+MULT16_16_Q15(Q15ONE-beta,q);
153    }
154 }
155
156 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)
157 {
158    int i;
159    /* Encode finer resolution */
160    for (i=0;i<m->nbEBands;i++)
161    {
162       int q2;
163       celt_int16_t frac = 1<<fine_quant[i];
164       celt_word16_t offset = (error[i]+QCONST16(.5f,8))*frac;
165       if (fine_quant[i] <= 0)
166          continue;
167 #ifdef FIXED_POINT
168       /* Has to be without rounding */
169       q2 = offset>>8;
170 #else
171       q2 = (int)floor(offset);
172 #endif
173       if (q2 > frac-1)
174          q2 = frac-1;
175       ec_enc_bits(enc, q2, fine_quant[i]);
176 #ifdef FIXED_POINT
177       offset = SUB16(SHR16(SHL16(q2,8)+QCONST16(.5,8),fine_quant[i]),QCONST16(.5f,8));
178 #else
179       offset = (q2+.5f)*(1<<(14-fine_quant[i]))*(1.f/16384) - .5f;
180 #endif
181       oldEBands[i] += PSHR32(MULT16_16(DB_SCALING*6,offset),8);
182       /*printf ("%f ", error[i] - offset);*/
183    }
184    for (i=0;i<m->nbEBands;i++)
185    {
186       eBands[i] = dB2Amp(oldEBands[i]);
187    }
188    /*printf ("%d\n", ec_enc_tell(enc, 0)-9);*/
189
190    /*printf ("\n");*/
191 }
192
193 static void unquant_coarse_energy_mono(const CELTMode *m, celt_ener_t *eBands, celt_word16_t *oldEBands, unsigned budget, int *prob, ec_dec *dec)
194 {
195    int i;
196    unsigned bits;
197    celt_word16_t prev = 0;
198    celt_word16_t coef = m->ePredCoef;
199    /* The .8 is a heuristic */
200    celt_word16_t beta = MULT16_16_Q15(QCONST16(.8f,15),coef);
201    
202    bits = ec_dec_tell(dec, 0);
203    /* Decode at a fixed coarse resolution */
204    for (i=0;i<m->nbEBands;i++)
205    {
206       int qi;
207       celt_word16_t q;
208       celt_word16_t mean = MULT16_16_Q15(Q15ONE-coef,eMeans[i]);
209       /* If we didn't have enough bits to encode all the energy, just assume something safe.
210          We allow slightly busting the budget here */
211       if (ec_dec_tell(dec, 0) - bits > budget)
212          qi = -1;
213       else
214          qi = ec_laplace_decode_start(dec, prob[2*i], prob[2*i+1]);
215       q = qi*base_resolution;
216       
217       oldEBands[i] = mean+MULT16_16_Q15(coef,oldEBands[i])+prev+q;
218       if (oldEBands[i] < -QCONST16(12.f,8))
219          oldEBands[i] = -QCONST16(12.f,8);
220       
221       prev = mean+prev+MULT16_16_Q15(Q15ONE-beta,q);
222    }
223 }
224
225 static void unquant_fine_energy_mono(const CELTMode *m, celt_ener_t *eBands, celt_word16_t *oldEBands, int *fine_quant, ec_dec *dec)
226 {
227    int i;
228    /* Decode finer resolution */
229    for (i=0;i<m->nbEBands;i++)
230    {
231       int q2;
232       celt_word16_t offset;
233       if (fine_quant[i] <= 0)
234          continue;
235       q2 = ec_dec_bits(dec, fine_quant[i]);
236 #ifdef FIXED_POINT
237       offset = SUB16(SHR16(SHL16(q2,8)+QCONST16(.5,8),fine_quant[i]),QCONST16(.5f,8));
238 #else
239       offset = (q2+.5f)*(1<<(14-fine_quant[i]))*(1.f/16384) - .5f;
240 #endif
241       oldEBands[i] += PSHR32(MULT16_16(DB_SCALING*6,offset),8);
242    }
243    for (i=0;i<m->nbEBands;i++)
244    {
245       eBands[i] = dB2Amp(oldEBands[i]);
246    }
247    /*printf ("\n");*/
248 }
249
250
251
252 void quant_coarse_energy(const CELTMode *m, celt_ener_t *eBands, celt_word16_t *oldEBands, int budget, int *prob, celt_word16_t *error, ec_enc *enc)
253 {
254    int C;
255    C = m->nbChannels;
256
257    if (C==1)
258    {
259       quant_coarse_energy_mono(m, eBands, oldEBands, budget, prob, error, enc);
260
261    } else {
262       int c;
263       for (c=0;c<C;c++)
264       {
265          int i;
266          VARDECL(celt_ener_t, E);
267          SAVE_STACK;
268          ALLOC(E, m->nbEBands, celt_ener_t);
269          for (i=0;i<m->nbEBands;i++)
270             E[i] = eBands[C*i+c];
271          quant_coarse_energy_mono(m, E, oldEBands+c*m->nbEBands, budget/C, prob, error+c*m->nbEBands, enc);
272          RESTORE_STACK;
273       }
274    }
275 }
276
277 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)
278 {
279    int C;
280    C = m->nbChannels;
281
282    if (C==1)
283    {
284       quant_fine_energy_mono(m, eBands, oldEBands, error, fine_quant, enc);
285
286    } else {
287       int c;
288       VARDECL(celt_ener_t, E);
289       ALLOC(E, m->nbEBands, celt_ener_t);
290       for (c=0;c<C;c++)
291       {
292          int i;
293          SAVE_STACK;
294          quant_fine_energy_mono(m, E, oldEBands+c*m->nbEBands, error+c*m->nbEBands, fine_quant, enc);
295          for (i=0;i<m->nbEBands;i++)
296             eBands[C*i+c] = E[i];
297          RESTORE_STACK;
298       }
299    }
300 }
301
302
303 void unquant_coarse_energy(const CELTMode *m, celt_ener_t *eBands, celt_word16_t *oldEBands, int budget, int *prob, ec_dec *dec)
304 {
305    int C;   
306
307    C = m->nbChannels;
308    if (C==1)
309    {
310       unquant_coarse_energy_mono(m, eBands, oldEBands, budget, prob, dec);
311    }
312    else {
313       int c;
314       VARDECL(celt_ener_t, E);
315       SAVE_STACK;
316       ALLOC(E, m->nbEBands, celt_ener_t);
317       for (c=0;c<C;c++)
318       {
319          unquant_coarse_energy_mono(m, E, oldEBands+c*m->nbEBands, budget/C, prob, dec);
320       }
321       RESTORE_STACK;
322    }
323 }
324
325 void unquant_fine_energy(const CELTMode *m, celt_ener_t *eBands, celt_word16_t *oldEBands, int *fine_quant, ec_dec *dec)
326 {
327    int C;   
328
329    C = m->nbChannels;
330
331    if (C==1)
332    {
333       unquant_fine_energy_mono(m, eBands, oldEBands, fine_quant, dec);
334    }
335    else {
336       int c;
337       VARDECL(celt_ener_t, E);
338       SAVE_STACK;
339       ALLOC(E, m->nbEBands, celt_ener_t);
340       for (c=0;c<C;c++)
341       {
342          int i;
343          unquant_fine_energy_mono(m, E, oldEBands+c*m->nbEBands, fine_quant, dec);
344          for (i=0;i<m->nbEBands;i++)
345             eBands[C*i+c] = E[i];
346       }
347       RESTORE_STACK;
348    }
349 }