diff options
Diffstat (limited to 'lib/rbcodec/codecs/libopus/celt/bands.c')
-rw-r--r-- | lib/rbcodec/codecs/libopus/celt/bands.c | 1302 |
1 files changed, 1302 insertions, 0 deletions
diff --git a/lib/rbcodec/codecs/libopus/celt/bands.c b/lib/rbcodec/codecs/libopus/celt/bands.c new file mode 100644 index 0000000000..6e612980b6 --- /dev/null +++ b/lib/rbcodec/codecs/libopus/celt/bands.c | |||
@@ -0,0 +1,1302 @@ | |||
1 | /* Copyright (c) 2007-2008 CSIRO | ||
2 | Copyright (c) 2007-2009 Xiph.Org Foundation | ||
3 | Copyright (c) 2008-2009 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 "opus_config.h" | ||
32 | #endif | ||
33 | |||
34 | #include <math.h> | ||
35 | #include "bands.h" | ||
36 | #include "modes.h" | ||
37 | #include "vq.h" | ||
38 | #include "cwrs.h" | ||
39 | #include "stack_alloc.h" | ||
40 | #include "os_support.h" | ||
41 | #include "mathops.h" | ||
42 | #include "rate.h" | ||
43 | |||
44 | opus_uint32 celt_lcg_rand(opus_uint32 seed) | ||
45 | { | ||
46 | return 1664525 * seed + 1013904223; | ||
47 | } | ||
48 | |||
49 | /* This is a cos() approximation designed to be bit-exact on any platform. Bit exactness | ||
50 | with this approximation is important because it has an impact on the bit allocation */ | ||
51 | static opus_int16 bitexact_cos(opus_int16 x) | ||
52 | { | ||
53 | opus_int32 tmp; | ||
54 | opus_int16 x2; | ||
55 | tmp = (4096+((opus_int32)(x)*(x)))>>13; | ||
56 | celt_assert(tmp<=32767); | ||
57 | x2 = tmp; | ||
58 | x2 = (32767-x2) + FRAC_MUL16(x2, (-7651 + FRAC_MUL16(x2, (8277 + FRAC_MUL16(-626, x2))))); | ||
59 | celt_assert(x2<=32766); | ||
60 | return 1+x2; | ||
61 | } | ||
62 | |||
63 | static int bitexact_log2tan(int isin,int icos) | ||
64 | { | ||
65 | int lc; | ||
66 | int ls; | ||
67 | lc=EC_ILOG(icos); | ||
68 | ls=EC_ILOG(isin); | ||
69 | icos<<=15-lc; | ||
70 | isin<<=15-ls; | ||
71 | return (ls-lc)*(1<<11) | ||
72 | +FRAC_MUL16(isin, FRAC_MUL16(isin, -2597) + 7932) | ||
73 | -FRAC_MUL16(icos, FRAC_MUL16(icos, -2597) + 7932); | ||
74 | } | ||
75 | |||
76 | #ifdef FIXED_POINT | ||
77 | /* Compute the amplitude (sqrt energy) in each of the bands */ | ||
78 | void compute_band_energies(const CELTMode *m, const celt_sig *X, celt_ener *bandE, int end, int C, int M) | ||
79 | { | ||
80 | int i, c, N; | ||
81 | const opus_int16 *eBands = m->eBands; | ||
82 | N = M*m->shortMdctSize; | ||
83 | c=0; do { | ||
84 | for (i=0;i<end;i++) | ||
85 | { | ||
86 | int j; | ||
87 | opus_val32 maxval=0; | ||
88 | opus_val32 sum = 0; | ||
89 | |||
90 | j=M*eBands[i]; do { | ||
91 | maxval = MAX32(maxval, X[j+c*N]); | ||
92 | maxval = MAX32(maxval, -X[j+c*N]); | ||
93 | } while (++j<M*eBands[i+1]); | ||
94 | |||
95 | if (maxval > 0) | ||
96 | { | ||
97 | int shift = celt_ilog2(maxval)-10; | ||
98 | j=M*eBands[i]; do { | ||
99 | sum = MAC16_16(sum, EXTRACT16(VSHR32(X[j+c*N],shift)), | ||
100 | EXTRACT16(VSHR32(X[j+c*N],shift))); | ||
101 | } while (++j<M*eBands[i+1]); | ||
102 | /* We're adding one here to ensure the normalized band isn't larger than unity norm */ | ||
103 | bandE[i+c*m->nbEBands] = EPSILON+VSHR32(EXTEND32(celt_sqrt(sum)),-shift); | ||
104 | } else { | ||
105 | bandE[i+c*m->nbEBands] = EPSILON; | ||
106 | } | ||
107 | /*printf ("%f ", bandE[i+c*m->nbEBands]);*/ | ||
108 | } | ||
109 | } while (++c<C); | ||
110 | /*printf ("\n");*/ | ||
111 | } | ||
112 | |||
113 | /* Normalise each band such that the energy is one. */ | ||
114 | void normalise_bands(const CELTMode *m, const celt_sig * OPUS_RESTRICT freq, celt_norm * OPUS_RESTRICT X, const celt_ener *bandE, int end, int C, int M) | ||
115 | { | ||
116 | int i, c, N; | ||
117 | const opus_int16 *eBands = m->eBands; | ||
118 | N = M*m->shortMdctSize; | ||
119 | c=0; do { | ||
120 | i=0; do { | ||
121 | opus_val16 g; | ||
122 | int j,shift; | ||
123 | opus_val16 E; | ||
124 | shift = celt_zlog2(bandE[i+c*m->nbEBands])-13; | ||
125 | E = VSHR32(bandE[i+c*m->nbEBands], shift); | ||
126 | g = EXTRACT16(celt_rcp(SHL32(E,3))); | ||
127 | j=M*eBands[i]; do { | ||
128 | X[j+c*N] = MULT16_16_Q15(VSHR32(freq[j+c*N],shift-1),g); | ||
129 | } while (++j<M*eBands[i+1]); | ||
130 | } while (++i<end); | ||
131 | } while (++c<C); | ||
132 | } | ||
133 | |||
134 | #else /* FIXED_POINT */ | ||
135 | /* Compute the amplitude (sqrt energy) in each of the bands */ | ||
136 | void compute_band_energies(const CELTMode *m, const celt_sig *X, celt_ener *bandE, int end, int C, int M) | ||
137 | { | ||
138 | int i, c, N; | ||
139 | const opus_int16 *eBands = m->eBands; | ||
140 | N = M*m->shortMdctSize; | ||
141 | c=0; do { | ||
142 | for (i=0;i<end;i++) | ||
143 | { | ||
144 | int j; | ||
145 | opus_val32 sum = 1e-27f; | ||
146 | for (j=M*eBands[i];j<M*eBands[i+1];j++) | ||
147 | sum += X[j+c*N]*X[j+c*N]; | ||
148 | bandE[i+c*m->nbEBands] = celt_sqrt(sum); | ||
149 | /*printf ("%f ", bandE[i+c*m->nbEBands]);*/ | ||
150 | } | ||
151 | } while (++c<C); | ||
152 | /*printf ("\n");*/ | ||
153 | } | ||
154 | |||
155 | /* Normalise each band such that the energy is one. */ | ||
156 | void normalise_bands(const CELTMode *m, const celt_sig * OPUS_RESTRICT freq, celt_norm * OPUS_RESTRICT X, const celt_ener *bandE, int end, int C, int M) | ||
157 | { | ||
158 | int i, c, N; | ||
159 | const opus_int16 *eBands = m->eBands; | ||
160 | N = M*m->shortMdctSize; | ||
161 | c=0; do { | ||
162 | for (i=0;i<end;i++) | ||
163 | { | ||
164 | int j; | ||
165 | opus_val16 g = 1.f/(1e-27f+bandE[i+c*m->nbEBands]); | ||
166 | for (j=M*eBands[i];j<M*eBands[i+1];j++) | ||
167 | X[j+c*N] = freq[j+c*N]*g; | ||
168 | } | ||
169 | } while (++c<C); | ||
170 | } | ||
171 | |||
172 | #endif /* FIXED_POINT */ | ||
173 | |||
174 | /* De-normalise the energy to produce the synthesis from the unit-energy bands */ | ||
175 | void denormalise_bands(const CELTMode *m, const celt_norm * OPUS_RESTRICT X, celt_sig * OPUS_RESTRICT freq, const celt_ener *bandE, int end, int C, int M) | ||
176 | { | ||
177 | int i, c, N; | ||
178 | const opus_int16 *eBands = m->eBands; | ||
179 | N = M*m->shortMdctSize; | ||
180 | celt_assert2(C<=2, "denormalise_bands() not implemented for >2 channels"); | ||
181 | c=0; do { | ||
182 | celt_sig * OPUS_RESTRICT f; | ||
183 | const celt_norm * OPUS_RESTRICT x; | ||
184 | f = freq+c*N; | ||
185 | x = X+c*N; | ||
186 | for (i=0;i<end;i++) | ||
187 | { | ||
188 | int j, band_end; | ||
189 | opus_val32 g = SHR32(bandE[i+c*m->nbEBands],1); | ||
190 | j=M*eBands[i]; | ||
191 | band_end = M*eBands[i+1]; | ||
192 | do { | ||
193 | *f++ = SHL32(MULT16_32_Q15(*x, g),2); | ||
194 | x++; | ||
195 | } while (++j<band_end); | ||
196 | } | ||
197 | for (i=M*eBands[end];i<N;i++) | ||
198 | *f++ = 0; | ||
199 | } while (++c<C); | ||
200 | } | ||
201 | |||
202 | /* This prevents energy collapse for transients with multiple short MDCTs */ | ||
203 | void anti_collapse(const CELTMode *m, celt_norm *X_, unsigned char *collapse_masks, int LM, int C, int size, | ||
204 | int start, int end, opus_val16 *logE, opus_val16 *prev1logE, | ||
205 | opus_val16 *prev2logE, int *pulses, opus_uint32 seed) | ||
206 | { | ||
207 | int c, i, j, k; | ||
208 | for (i=start;i<end;i++) | ||
209 | { | ||
210 | int N0; | ||
211 | opus_val16 thresh, sqrt_1; | ||
212 | int depth; | ||
213 | #ifdef FIXED_POINT | ||
214 | int shift; | ||
215 | opus_val32 thresh32; | ||
216 | #endif | ||
217 | |||
218 | N0 = m->eBands[i+1]-m->eBands[i]; | ||
219 | /* depth in 1/8 bits */ | ||
220 | depth = (1+pulses[i])/((m->eBands[i+1]-m->eBands[i])<<LM); | ||
221 | |||
222 | #ifdef FIXED_POINT | ||
223 | thresh32 = SHR32(celt_exp2(-SHL16(depth, 10-BITRES)),1); | ||
224 | thresh = MULT16_32_Q15(QCONST16(0.5f, 15), MIN32(32767,thresh32)); | ||
225 | { | ||
226 | opus_val32 t; | ||
227 | t = N0<<LM; | ||
228 | shift = celt_ilog2(t)>>1; | ||
229 | t = SHL32(t, (7-shift)<<1); | ||
230 | sqrt_1 = celt_rsqrt_norm(t); | ||
231 | } | ||
232 | #else | ||
233 | thresh = .5f*celt_exp2(-.125f*depth); | ||
234 | sqrt_1 = celt_rsqrt(N0<<LM); | ||
235 | #endif | ||
236 | |||
237 | c=0; do | ||
238 | { | ||
239 | celt_norm *X; | ||
240 | opus_val16 prev1; | ||
241 | opus_val16 prev2; | ||
242 | opus_val32 Ediff; | ||
243 | opus_val16 r; | ||
244 | int renormalize=0; | ||
245 | prev1 = prev1logE[c*m->nbEBands+i]; | ||
246 | prev2 = prev2logE[c*m->nbEBands+i]; | ||
247 | if (C==1) | ||
248 | { | ||
249 | prev1 = MAX16(prev1,prev1logE[m->nbEBands+i]); | ||
250 | prev2 = MAX16(prev2,prev2logE[m->nbEBands+i]); | ||
251 | } | ||
252 | Ediff = EXTEND32(logE[c*m->nbEBands+i])-EXTEND32(MIN16(prev1,prev2)); | ||
253 | Ediff = MAX32(0, Ediff); | ||
254 | |||
255 | #ifdef FIXED_POINT | ||
256 | if (Ediff < 16384) | ||
257 | { | ||
258 | opus_val32 r32 = SHR32(celt_exp2(-EXTRACT16(Ediff)),1); | ||
259 | r = 2*MIN16(16383,r32); | ||
260 | } else { | ||
261 | r = 0; | ||
262 | } | ||
263 | if (LM==3) | ||
264 | r = MULT16_16_Q14(23170, MIN32(23169, r)); | ||
265 | r = SHR16(MIN16(thresh, r),1); | ||
266 | r = SHR32(MULT16_16_Q15(sqrt_1, r),shift); | ||
267 | #else | ||
268 | /* r needs to be multiplied by 2 or 2*sqrt(2) depending on LM because | ||
269 | short blocks don't have the same energy as long */ | ||
270 | r = 2.f*celt_exp2(-Ediff); | ||
271 | if (LM==3) | ||
272 | r *= 1.41421356f; | ||
273 | r = MIN16(thresh, r); | ||
274 | r = r*sqrt_1; | ||
275 | #endif | ||
276 | X = X_+c*size+(m->eBands[i]<<LM); | ||
277 | for (k=0;k<1<<LM;k++) | ||
278 | { | ||
279 | /* Detect collapse */ | ||
280 | if (!(collapse_masks[i*C+c]&1<<k)) | ||
281 | { | ||
282 | /* Fill with noise */ | ||
283 | for (j=0;j<N0;j++) | ||
284 | { | ||
285 | seed = celt_lcg_rand(seed); | ||
286 | X[(j<<LM)+k] = (seed&0x8000 ? r : -r); | ||
287 | } | ||
288 | renormalize = 1; | ||
289 | } | ||
290 | } | ||
291 | /* We just added some energy, so we need to renormalise */ | ||
292 | if (renormalize) | ||
293 | renormalise_vector(X, N0<<LM, Q15ONE); | ||
294 | } while (++c<C); | ||
295 | } | ||
296 | } | ||
297 | |||
298 | static void intensity_stereo(const CELTMode *m, celt_norm *X, celt_norm *Y, const celt_ener *bandE, int bandID, int N) | ||
299 | { | ||
300 | int i = bandID; | ||
301 | int j; | ||
302 | opus_val16 a1, a2; | ||
303 | opus_val16 left, right; | ||
304 | opus_val16 norm; | ||
305 | #ifdef FIXED_POINT | ||
306 | int shift = celt_zlog2(MAX32(bandE[i], bandE[i+m->nbEBands]))-13; | ||
307 | #endif | ||
308 | left = VSHR32(bandE[i],shift); | ||
309 | right = VSHR32(bandE[i+m->nbEBands],shift); | ||
310 | norm = EPSILON + celt_sqrt(EPSILON+MULT16_16(left,left)+MULT16_16(right,right)); | ||
311 | a1 = DIV32_16(SHL32(EXTEND32(left),14),norm); | ||
312 | a2 = DIV32_16(SHL32(EXTEND32(right),14),norm); | ||
313 | for (j=0;j<N;j++) | ||
314 | { | ||
315 | celt_norm r, l; | ||
316 | l = X[j]; | ||
317 | r = Y[j]; | ||
318 | X[j] = MULT16_16_Q14(a1,l) + MULT16_16_Q14(a2,r); | ||
319 | /* Side is not encoded, no need to calculate */ | ||
320 | } | ||
321 | } | ||
322 | |||
323 | static void stereo_split(celt_norm *X, celt_norm *Y, int N) | ||
324 | { | ||
325 | int j; | ||
326 | for (j=0;j<N;j++) | ||
327 | { | ||
328 | celt_norm r, l; | ||
329 | l = MULT16_16_Q15(QCONST16(.70710678f,15), X[j]); | ||
330 | r = MULT16_16_Q15(QCONST16(.70710678f,15), Y[j]); | ||
331 | X[j] = l+r; | ||
332 | Y[j] = r-l; | ||
333 | } | ||
334 | } | ||
335 | |||
336 | static void stereo_merge(celt_norm *X, celt_norm *Y, opus_val16 mid, int N) | ||
337 | { | ||
338 | int j; | ||
339 | opus_val32 xp=0, side=0; | ||
340 | opus_val32 El, Er; | ||
341 | opus_val16 mid2; | ||
342 | #ifdef FIXED_POINT | ||
343 | int kl, kr; | ||
344 | #endif | ||
345 | opus_val32 t, lgain, rgain; | ||
346 | |||
347 | /* Compute the norm of X+Y and X-Y as |X|^2 + |Y|^2 +/- sum(xy) */ | ||
348 | for (j=0;j<N;j++) | ||
349 | { | ||
350 | xp = MAC16_16(xp, X[j], Y[j]); | ||
351 | side = MAC16_16(side, Y[j], Y[j]); | ||
352 | } | ||
353 | /* Compensating for the mid normalization */ | ||
354 | xp = MULT16_32_Q15(mid, xp); | ||
355 | /* mid and side are in Q15, not Q14 like X and Y */ | ||
356 | mid2 = SHR32(mid, 1); | ||
357 | El = MULT16_16(mid2, mid2) + side - 2*xp; | ||
358 | Er = MULT16_16(mid2, mid2) + side + 2*xp; | ||
359 | if (Er < QCONST32(6e-4f, 28) || El < QCONST32(6e-4f, 28)) | ||
360 | { | ||
361 | for (j=0;j<N;j++) | ||
362 | Y[j] = X[j]; | ||
363 | return; | ||
364 | } | ||
365 | |||
366 | #ifdef FIXED_POINT | ||
367 | kl = celt_ilog2(El)>>1; | ||
368 | kr = celt_ilog2(Er)>>1; | ||
369 | #endif | ||
370 | t = VSHR32(El, (kl-7)<<1); | ||
371 | lgain = celt_rsqrt_norm(t); | ||
372 | t = VSHR32(Er, (kr-7)<<1); | ||
373 | rgain = celt_rsqrt_norm(t); | ||
374 | |||
375 | #ifdef FIXED_POINT | ||
376 | if (kl < 7) | ||
377 | kl = 7; | ||
378 | if (kr < 7) | ||
379 | kr = 7; | ||
380 | #endif | ||
381 | |||
382 | for (j=0;j<N;j++) | ||
383 | { | ||
384 | celt_norm r, l; | ||
385 | /* Apply mid scaling (side is already scaled) */ | ||
386 | l = MULT16_16_Q15(mid, X[j]); | ||
387 | r = Y[j]; | ||
388 | X[j] = EXTRACT16(PSHR32(MULT16_16(lgain, SUB16(l,r)), kl+1)); | ||
389 | Y[j] = EXTRACT16(PSHR32(MULT16_16(rgain, ADD16(l,r)), kr+1)); | ||
390 | } | ||
391 | } | ||
392 | |||
393 | /* Decide whether we should spread the pulses in the current frame */ | ||
394 | int spreading_decision(const CELTMode *m, celt_norm *X, int *average, | ||
395 | int last_decision, int *hf_average, int *tapset_decision, int update_hf, | ||
396 | int end, int C, int M) | ||
397 | { | ||
398 | int i, c, N0; | ||
399 | int sum = 0, nbBands=0; | ||
400 | const opus_int16 * OPUS_RESTRICT eBands = m->eBands; | ||
401 | int decision; | ||
402 | int hf_sum=0; | ||
403 | |||
404 | celt_assert(end>0); | ||
405 | |||
406 | N0 = M*m->shortMdctSize; | ||
407 | |||
408 | if (M*(eBands[end]-eBands[end-1]) <= 8) | ||
409 | return SPREAD_NONE; | ||
410 | c=0; do { | ||
411 | for (i=0;i<end;i++) | ||
412 | { | ||
413 | int j, N, tmp=0; | ||
414 | int tcount[3] = {0,0,0}; | ||
415 | celt_norm * OPUS_RESTRICT x = X+M*eBands[i]+c*N0; | ||
416 | N = M*(eBands[i+1]-eBands[i]); | ||
417 | if (N<=8) | ||
418 | continue; | ||
419 | /* Compute rough CDF of |x[j]| */ | ||
420 | for (j=0;j<N;j++) | ||
421 | { | ||
422 | opus_val32 x2N; /* Q13 */ | ||
423 | |||
424 | x2N = MULT16_16(MULT16_16_Q15(x[j], x[j]), N); | ||
425 | if (x2N < QCONST16(0.25f,13)) | ||
426 | tcount[0]++; | ||
427 | if (x2N < QCONST16(0.0625f,13)) | ||
428 | tcount[1]++; | ||
429 | if (x2N < QCONST16(0.015625f,13)) | ||
430 | tcount[2]++; | ||
431 | } | ||
432 | |||
433 | /* Only include four last bands (8 kHz and up) */ | ||
434 | if (i>m->nbEBands-4) | ||
435 | hf_sum += 32*(tcount[1]+tcount[0])/N; | ||
436 | tmp = (2*tcount[2] >= N) + (2*tcount[1] >= N) + (2*tcount[0] >= N); | ||
437 | sum += tmp*256; | ||
438 | nbBands++; | ||
439 | } | ||
440 | } while (++c<C); | ||
441 | |||
442 | if (update_hf) | ||
443 | { | ||
444 | if (hf_sum) | ||
445 | hf_sum /= C*(4-m->nbEBands+end); | ||
446 | *hf_average = (*hf_average+hf_sum)>>1; | ||
447 | hf_sum = *hf_average; | ||
448 | if (*tapset_decision==2) | ||
449 | hf_sum += 4; | ||
450 | else if (*tapset_decision==0) | ||
451 | hf_sum -= 4; | ||
452 | if (hf_sum > 22) | ||
453 | *tapset_decision=2; | ||
454 | else if (hf_sum > 18) | ||
455 | *tapset_decision=1; | ||
456 | else | ||
457 | *tapset_decision=0; | ||
458 | } | ||
459 | /*printf("%d %d %d\n", hf_sum, *hf_average, *tapset_decision);*/ | ||
460 | celt_assert(nbBands>0); /*M*(eBands[end]-eBands[end-1]) <= 8 assures this*/ | ||
461 | sum /= nbBands; | ||
462 | /* Recursive averaging */ | ||
463 | sum = (sum+*average)>>1; | ||
464 | *average = sum; | ||
465 | /* Hysteresis */ | ||
466 | sum = (3*sum + (((3-last_decision)<<7) + 64) + 2)>>2; | ||
467 | if (sum < 80) | ||
468 | { | ||
469 | decision = SPREAD_AGGRESSIVE; | ||
470 | } else if (sum < 256) | ||
471 | { | ||
472 | decision = SPREAD_NORMAL; | ||
473 | } else if (sum < 384) | ||
474 | { | ||
475 | decision = SPREAD_LIGHT; | ||
476 | } else { | ||
477 | decision = SPREAD_NONE; | ||
478 | } | ||
479 | #ifdef FUZZING | ||
480 | decision = rand()&0x3; | ||
481 | *tapset_decision=rand()%3; | ||
482 | #endif | ||
483 | return decision; | ||
484 | } | ||
485 | |||
486 | #ifdef MEASURE_NORM_MSE | ||
487 | |||
488 | float MSE[30] = {0}; | ||
489 | int nbMSEBands = 0; | ||
490 | int MSECount[30] = {0}; | ||
491 | |||
492 | void dump_norm_mse(void) | ||
493 | { | ||
494 | int i; | ||
495 | for (i=0;i<nbMSEBands;i++) | ||
496 | { | ||
497 | printf ("%g ", MSE[i]/MSECount[i]); | ||
498 | } | ||
499 | printf ("\n"); | ||
500 | } | ||
501 | |||
502 | void measure_norm_mse(const CELTMode *m, float *X, float *X0, float *bandE, float *bandE0, int M, int N, int C) | ||
503 | { | ||
504 | static int init = 0; | ||
505 | int i; | ||
506 | if (!init) | ||
507 | { | ||
508 | atexit(dump_norm_mse); | ||
509 | init = 1; | ||
510 | } | ||
511 | for (i=0;i<m->nbEBands;i++) | ||
512 | { | ||
513 | int j; | ||
514 | int c; | ||
515 | float g; | ||
516 | if (bandE0[i]<10 || (C==2 && bandE0[i+m->nbEBands]<1)) | ||
517 | continue; | ||
518 | c=0; do { | ||
519 | g = bandE[i+c*m->nbEBands]/(1e-15+bandE0[i+c*m->nbEBands]); | ||
520 | for (j=M*m->eBands[i];j<M*m->eBands[i+1];j++) | ||
521 | MSE[i] += (g*X[j+c*N]-X0[j+c*N])*(g*X[j+c*N]-X0[j+c*N]); | ||
522 | } while (++c<C); | ||
523 | MSECount[i]+=C; | ||
524 | } | ||
525 | nbMSEBands = m->nbEBands; | ||
526 | } | ||
527 | |||
528 | #endif | ||
529 | |||
530 | /* Indexing table for converting from natural Hadamard to ordery Hadamard | ||
531 | This is essentially a bit-reversed Gray, on top of which we've added | ||
532 | an inversion of the order because we want the DC at the end rather than | ||
533 | the beginning. The lines are for N=2, 4, 8, 16 */ | ||
534 | static const int ordery_table[] = { | ||
535 | 1, 0, | ||
536 | 3, 0, 2, 1, | ||
537 | 7, 0, 4, 3, 6, 1, 5, 2, | ||
538 | 15, 0, 8, 7, 12, 3, 11, 4, 14, 1, 9, 6, 13, 2, 10, 5, | ||
539 | }; | ||
540 | |||
541 | static void deinterleave_hadamard(celt_norm *X, int N0, int stride, int hadamard) | ||
542 | { | ||
543 | int i,j; | ||
544 | VARDECL(celt_norm, tmp); | ||
545 | int N; | ||
546 | SAVE_STACK; | ||
547 | N = N0*stride; | ||
548 | ALLOC(tmp, N, celt_norm); | ||
549 | celt_assert(stride>0); | ||
550 | if (hadamard) | ||
551 | { | ||
552 | const int *ordery = ordery_table+stride-2; | ||
553 | for (i=0;i<stride;i++) | ||
554 | { | ||
555 | for (j=0;j<N0;j++) | ||
556 | tmp[ordery[i]*N0+j] = X[j*stride+i]; | ||
557 | } | ||
558 | } else { | ||
559 | for (i=0;i<stride;i++) | ||
560 | for (j=0;j<N0;j++) | ||
561 | tmp[i*N0+j] = X[j*stride+i]; | ||
562 | } | ||
563 | for (j=0;j<N;j++) | ||
564 | X[j] = tmp[j]; | ||
565 | RESTORE_STACK; | ||
566 | } | ||
567 | |||
568 | static void interleave_hadamard(celt_norm *X, int N0, int stride, int hadamard) | ||
569 | { | ||
570 | int i,j; | ||
571 | VARDECL(celt_norm, tmp); | ||
572 | int N; | ||
573 | SAVE_STACK; | ||
574 | N = N0*stride; | ||
575 | ALLOC(tmp, N, celt_norm); | ||
576 | if (hadamard) | ||
577 | { | ||
578 | const int *ordery = ordery_table+stride-2; | ||
579 | for (i=0;i<stride;i++) | ||
580 | for (j=0;j<N0;j++) | ||
581 | tmp[j*stride+i] = X[ordery[i]*N0+j]; | ||
582 | } else { | ||
583 | for (i=0;i<stride;i++) | ||
584 | for (j=0;j<N0;j++) | ||
585 | tmp[j*stride+i] = X[i*N0+j]; | ||
586 | } | ||
587 | for (j=0;j<N;j++) | ||
588 | X[j] = tmp[j]; | ||
589 | RESTORE_STACK; | ||
590 | } | ||
591 | |||
592 | void haar1(celt_norm *X, int N0, int stride) | ||
593 | { | ||
594 | int i, j; | ||
595 | N0 >>= 1; | ||
596 | for (i=0;i<stride;i++) | ||
597 | for (j=0;j<N0;j++) | ||
598 | { | ||
599 | celt_norm tmp1, tmp2; | ||
600 | tmp1 = MULT16_16_Q15(QCONST16(.70710678f,15), X[stride*2*j+i]); | ||
601 | tmp2 = MULT16_16_Q15(QCONST16(.70710678f,15), X[stride*(2*j+1)+i]); | ||
602 | X[stride*2*j+i] = tmp1 + tmp2; | ||
603 | X[stride*(2*j+1)+i] = tmp1 - tmp2; | ||
604 | } | ||
605 | } | ||
606 | |||
607 | static int compute_qn(int N, int b, int offset, int pulse_cap, int stereo) | ||
608 | { | ||
609 | static const opus_int16 exp2_table8[8] = | ||
610 | {16384, 17866, 19483, 21247, 23170, 25267, 27554, 30048}; | ||
611 | int qn, qb; | ||
612 | int N2 = 2*N-1; | ||
613 | if (stereo && N==2) | ||
614 | N2--; | ||
615 | /* The upper limit ensures that in a stereo split with itheta==16384, we'll | ||
616 | always have enough bits left over to code at least one pulse in the | ||
617 | side; otherwise it would collapse, since it doesn't get folded. */ | ||
618 | qb = IMIN(b-pulse_cap-(4<<BITRES), (b+N2*offset)/N2); | ||
619 | |||
620 | qb = IMIN(8<<BITRES, qb); | ||
621 | |||
622 | if (qb<(1<<BITRES>>1)) { | ||
623 | qn = 1; | ||
624 | } else { | ||
625 | qn = exp2_table8[qb&0x7]>>(14-(qb>>BITRES)); | ||
626 | qn = (qn+1)>>1<<1; | ||
627 | } | ||
628 | celt_assert(qn <= 256); | ||
629 | return qn; | ||
630 | } | ||
631 | |||
632 | /* This function is responsible for encoding and decoding a band for both | ||
633 | the mono and stereo case. Even in the mono case, it can split the band | ||
634 | in two and transmit the energy difference with the two half-bands. It | ||
635 | can be called recursively so bands can end up being split in 8 parts. */ | ||
636 | static unsigned quant_band(int encode, const CELTMode *m, int i, celt_norm *X, celt_norm *Y, | ||
637 | int N, int b, int spread, int B, int intensity, int tf_change, celt_norm *lowband, ec_ctx *ec, | ||
638 | opus_int32 *remaining_bits, int LM, celt_norm *lowband_out, const celt_ener *bandE, int level, | ||
639 | opus_uint32 *seed, opus_val16 gain, celt_norm *lowband_scratch, int fill) | ||
640 | { | ||
641 | const unsigned char *cache; | ||
642 | int q; | ||
643 | int curr_bits; | ||
644 | int stereo, split; | ||
645 | int imid=0, iside=0; | ||
646 | int N0=N; | ||
647 | int N_B=N; | ||
648 | int N_B0; | ||
649 | int B0=B; | ||
650 | int time_divide=0; | ||
651 | int recombine=0; | ||
652 | int inv = 0; | ||
653 | opus_val16 mid=0, side=0; | ||
654 | int longBlocks; | ||
655 | unsigned cm=0; | ||
656 | #ifdef RESYNTH | ||
657 | int resynth = 1; | ||
658 | #else | ||
659 | int resynth = !encode; | ||
660 | #endif | ||
661 | |||
662 | longBlocks = B0==1; | ||
663 | |||
664 | N_B /= B; | ||
665 | N_B0 = N_B; | ||
666 | |||
667 | split = stereo = Y != NULL; | ||
668 | |||
669 | /* Special case for one sample */ | ||
670 | if (N==1) | ||
671 | { | ||
672 | int c; | ||
673 | celt_norm *x = X; | ||
674 | c=0; do { | ||
675 | int sign=0; | ||
676 | if (*remaining_bits>=1<<BITRES) | ||
677 | { | ||
678 | if (encode) | ||
679 | { | ||
680 | sign = x[0]<0; | ||
681 | ec_enc_bits(ec, sign, 1); | ||
682 | } else { | ||
683 | sign = ec_dec_bits(ec, 1); | ||
684 | } | ||
685 | *remaining_bits -= 1<<BITRES; | ||
686 | b-=1<<BITRES; | ||
687 | } | ||
688 | if (resynth) | ||
689 | x[0] = sign ? -NORM_SCALING : NORM_SCALING; | ||
690 | x = Y; | ||
691 | } while (++c<1+stereo); | ||
692 | if (lowband_out) | ||
693 | lowband_out[0] = SHR16(X[0],4); | ||
694 | return 1; | ||
695 | } | ||
696 | |||
697 | if (!stereo && level == 0) | ||
698 | { | ||
699 | int k; | ||
700 | if (tf_change>0) | ||
701 | recombine = tf_change; | ||
702 | /* Band recombining to increase frequency resolution */ | ||
703 | |||
704 | if (lowband && (recombine || ((N_B&1) == 0 && tf_change<0) || B0>1)) | ||
705 | { | ||
706 | int j; | ||
707 | for (j=0;j<N;j++) | ||
708 | lowband_scratch[j] = lowband[j]; | ||
709 | lowband = lowband_scratch; | ||
710 | } | ||
711 | |||
712 | for (k=0;k<recombine;k++) | ||
713 | { | ||
714 | static const unsigned char bit_interleave_table[16]={ | ||
715 | 0,1,1,1,2,3,3,3,2,3,3,3,2,3,3,3 | ||
716 | }; | ||
717 | if (encode) | ||
718 | haar1(X, N>>k, 1<<k); | ||
719 | if (lowband) | ||
720 | haar1(lowband, N>>k, 1<<k); | ||
721 | fill = bit_interleave_table[fill&0xF]|bit_interleave_table[fill>>4]<<2; | ||
722 | } | ||
723 | B>>=recombine; | ||
724 | N_B<<=recombine; | ||
725 | |||
726 | /* Increasing the time resolution */ | ||
727 | while ((N_B&1) == 0 && tf_change<0) | ||
728 | { | ||
729 | if (encode) | ||
730 | haar1(X, N_B, B); | ||
731 | if (lowband) | ||
732 | haar1(lowband, N_B, B); | ||
733 | fill |= fill<<B; | ||
734 | B <<= 1; | ||
735 | N_B >>= 1; | ||
736 | time_divide++; | ||
737 | tf_change++; | ||
738 | } | ||
739 | B0=B; | ||
740 | N_B0 = N_B; | ||
741 | |||
742 | /* Reorganize the samples in time order instead of frequency order */ | ||
743 | if (B0>1) | ||
744 | { | ||
745 | if (encode) | ||
746 | deinterleave_hadamard(X, N_B>>recombine, B0<<recombine, longBlocks); | ||
747 | if (lowband) | ||
748 | deinterleave_hadamard(lowband, N_B>>recombine, B0<<recombine, longBlocks); | ||
749 | } | ||
750 | } | ||
751 | |||
752 | /* If we need 1.5 more bit than we can produce, split the band in two. */ | ||
753 | cache = m->cache.bits + m->cache.index[(LM+1)*m->nbEBands+i]; | ||
754 | if (!stereo && LM != -1 && b > cache[cache[0]]+12 && N>2) | ||
755 | { | ||
756 | N >>= 1; | ||
757 | Y = X+N; | ||
758 | split = 1; | ||
759 | LM -= 1; | ||
760 | if (B==1) | ||
761 | fill = (fill&1)|(fill<<1); | ||
762 | B = (B+1)>>1; | ||
763 | } | ||
764 | |||
765 | if (split) | ||
766 | { | ||
767 | int qn; | ||
768 | int itheta=0; | ||
769 | int mbits, sbits, delta; | ||
770 | int qalloc; | ||
771 | int pulse_cap; | ||
772 | int offset; | ||
773 | int orig_fill; | ||
774 | opus_int32 tell; | ||
775 | |||
776 | /* Decide on the resolution to give to the split parameter theta */ | ||
777 | pulse_cap = m->logN[i]+LM*(1<<BITRES); | ||
778 | offset = (pulse_cap>>1) - (stereo&&N==2 ? QTHETA_OFFSET_TWOPHASE : QTHETA_OFFSET); | ||
779 | qn = compute_qn(N, b, offset, pulse_cap, stereo); | ||
780 | if (stereo && i>=intensity) | ||
781 | qn = 1; | ||
782 | if (encode) | ||
783 | { | ||
784 | /* theta is the atan() of the ratio between the (normalized) | ||
785 | side and mid. With just that parameter, we can re-scale both | ||
786 | mid and side because we know that 1) they have unit norm and | ||
787 | 2) they are orthogonal. */ | ||
788 | itheta = stereo_itheta(X, Y, stereo, N); | ||
789 | } | ||
790 | tell = ec_tell_frac(ec); | ||
791 | if (qn!=1) | ||
792 | { | ||
793 | if (encode) | ||
794 | itheta = (itheta*qn+8192)>>14; | ||
795 | |||
796 | /* Entropy coding of the angle. We use a uniform pdf for the | ||
797 | time split, a step for stereo, and a triangular one for the rest. */ | ||
798 | if (stereo && N>2) | ||
799 | { | ||
800 | int p0 = 3; | ||
801 | int x = itheta; | ||
802 | int x0 = qn/2; | ||
803 | int ft = p0*(x0+1) + x0; | ||
804 | /* Use a probability of p0 up to itheta=8192 and then use 1 after */ | ||
805 | if (encode) | ||
806 | { | ||
807 | ec_encode(ec,x<=x0?p0*x:(x-1-x0)+(x0+1)*p0,x<=x0?p0*(x+1):(x-x0)+(x0+1)*p0,ft); | ||
808 | } else { | ||
809 | int fs; | ||
810 | fs=ec_decode(ec,ft); | ||
811 | if (fs<(x0+1)*p0) | ||
812 | x=fs/p0; | ||
813 | else | ||
814 | x=x0+1+(fs-(x0+1)*p0); | ||
815 | ec_dec_update(ec,x<=x0?p0*x:(x-1-x0)+(x0+1)*p0,x<=x0?p0*(x+1):(x-x0)+(x0+1)*p0,ft); | ||
816 | itheta = x; | ||
817 | } | ||
818 | } else if (B0>1 || stereo) { | ||
819 | /* Uniform pdf */ | ||
820 | if (encode) | ||
821 | ec_enc_uint(ec, itheta, qn+1); | ||
822 | else | ||
823 | itheta = ec_dec_uint(ec, qn+1); | ||
824 | } else { | ||
825 | int fs=1, ft; | ||
826 | ft = ((qn>>1)+1)*((qn>>1)+1); | ||
827 | if (encode) | ||
828 | { | ||
829 | int fl; | ||
830 | |||
831 | fs = itheta <= (qn>>1) ? itheta + 1 : qn + 1 - itheta; | ||
832 | fl = itheta <= (qn>>1) ? itheta*(itheta + 1)>>1 : | ||
833 | ft - ((qn + 1 - itheta)*(qn + 2 - itheta)>>1); | ||
834 | |||
835 | ec_encode(ec, fl, fl+fs, ft); | ||
836 | } else { | ||
837 | /* Triangular pdf */ | ||
838 | int fl=0; | ||
839 | int fm; | ||
840 | fm = ec_decode(ec, ft); | ||
841 | |||
842 | if (fm < ((qn>>1)*((qn>>1) + 1)>>1)) | ||
843 | { | ||
844 | itheta = (isqrt32(8*(opus_uint32)fm + 1) - 1)>>1; | ||
845 | fs = itheta + 1; | ||
846 | fl = itheta*(itheta + 1)>>1; | ||
847 | } | ||
848 | else | ||
849 | { | ||
850 | itheta = (2*(qn + 1) | ||
851 | - isqrt32(8*(opus_uint32)(ft - fm - 1) + 1))>>1; | ||
852 | fs = qn + 1 - itheta; | ||
853 | fl = ft - ((qn + 1 - itheta)*(qn + 2 - itheta)>>1); | ||
854 | } | ||
855 | |||
856 | ec_dec_update(ec, fl, fl+fs, ft); | ||
857 | } | ||
858 | } | ||
859 | itheta = (opus_int32)itheta*16384/qn; | ||
860 | if (encode && stereo) | ||
861 | { | ||
862 | if (itheta==0) | ||
863 | intensity_stereo(m, X, Y, bandE, i, N); | ||
864 | else | ||
865 | stereo_split(X, Y, N); | ||
866 | } | ||
867 | /* NOTE: Renormalising X and Y *may* help fixed-point a bit at very high rate. | ||
868 | Let's do that at higher complexity */ | ||
869 | } else if (stereo) { | ||
870 | if (encode) | ||
871 | { | ||
872 | inv = itheta > 8192; | ||
873 | if (inv) | ||
874 | { | ||
875 | int j; | ||
876 | for (j=0;j<N;j++) | ||
877 | Y[j] = -Y[j]; | ||
878 | } | ||
879 | intensity_stereo(m, X, Y, bandE, i, N); | ||
880 | } | ||
881 | if (b>2<<BITRES && *remaining_bits > 2<<BITRES) | ||
882 | { | ||
883 | if (encode) | ||
884 | ec_enc_bit_logp(ec, inv, 2); | ||
885 | else | ||
886 | inv = ec_dec_bit_logp(ec, 2); | ||
887 | } else | ||
888 | inv = 0; | ||
889 | itheta = 0; | ||
890 | } | ||
891 | qalloc = ec_tell_frac(ec) - tell; | ||
892 | b -= qalloc; | ||
893 | |||
894 | orig_fill = fill; | ||
895 | if (itheta == 0) | ||
896 | { | ||
897 | imid = 32767; | ||
898 | iside = 0; | ||
899 | fill &= (1<<B)-1; | ||
900 | delta = -16384; | ||
901 | } else if (itheta == 16384) | ||
902 | { | ||
903 | imid = 0; | ||
904 | iside = 32767; | ||
905 | fill &= ((1<<B)-1)<<B; | ||
906 | delta = 16384; | ||
907 | } else { | ||
908 | imid = bitexact_cos(itheta); | ||
909 | iside = bitexact_cos(16384-itheta); | ||
910 | /* This is the mid vs side allocation that minimizes squared error | ||
911 | in that band. */ | ||
912 | delta = FRAC_MUL16((N-1)<<7,bitexact_log2tan(iside,imid)); | ||
913 | } | ||
914 | |||
915 | #ifdef FIXED_POINT | ||
916 | mid = imid; | ||
917 | side = iside; | ||
918 | #else | ||
919 | mid = (1.f/32768)*imid; | ||
920 | side = (1.f/32768)*iside; | ||
921 | #endif | ||
922 | |||
923 | /* This is a special case for N=2 that only works for stereo and takes | ||
924 | advantage of the fact that mid and side are orthogonal to encode | ||
925 | the side with just one bit. */ | ||
926 | if (N==2 && stereo) | ||
927 | { | ||
928 | int c; | ||
929 | int sign=0; | ||
930 | celt_norm *x2, *y2; | ||
931 | mbits = b; | ||
932 | sbits = 0; | ||
933 | /* Only need one bit for the side */ | ||
934 | if (itheta != 0 && itheta != 16384) | ||
935 | sbits = 1<<BITRES; | ||
936 | mbits -= sbits; | ||
937 | c = itheta > 8192; | ||
938 | *remaining_bits -= qalloc+sbits; | ||
939 | |||
940 | x2 = c ? Y : X; | ||
941 | y2 = c ? X : Y; | ||
942 | if (sbits) | ||
943 | { | ||
944 | if (encode) | ||
945 | { | ||
946 | /* Here we only need to encode a sign for the side */ | ||
947 | sign = x2[0]*y2[1] - x2[1]*y2[0] < 0; | ||
948 | ec_enc_bits(ec, sign, 1); | ||
949 | } else { | ||
950 | sign = ec_dec_bits(ec, 1); | ||
951 | } | ||
952 | } | ||
953 | sign = 1-2*sign; | ||
954 | /* We use orig_fill here because we want to fold the side, but if | ||
955 | itheta==16384, we'll have cleared the low bits of fill. */ | ||
956 | cm = quant_band(encode, m, i, x2, NULL, N, mbits, spread, B, intensity, tf_change, lowband, ec, remaining_bits, LM, lowband_out, NULL, level, seed, gain, lowband_scratch, orig_fill); | ||
957 | /* We don't split N=2 bands, so cm is either 1 or 0 (for a fold-collapse), | ||
958 | and there's no need to worry about mixing with the other channel. */ | ||
959 | y2[0] = -sign*x2[1]; | ||
960 | y2[1] = sign*x2[0]; | ||
961 | if (resynth) | ||
962 | { | ||
963 | celt_norm tmp; | ||
964 | X[0] = MULT16_16_Q15(mid, X[0]); | ||
965 | X[1] = MULT16_16_Q15(mid, X[1]); | ||
966 | Y[0] = MULT16_16_Q15(side, Y[0]); | ||
967 | Y[1] = MULT16_16_Q15(side, Y[1]); | ||
968 | tmp = X[0]; | ||
969 | X[0] = SUB16(tmp,Y[0]); | ||
970 | Y[0] = ADD16(tmp,Y[0]); | ||
971 | tmp = X[1]; | ||
972 | X[1] = SUB16(tmp,Y[1]); | ||
973 | Y[1] = ADD16(tmp,Y[1]); | ||
974 | } | ||
975 | } else { | ||
976 | /* "Normal" split code */ | ||
977 | celt_norm *next_lowband2=NULL; | ||
978 | celt_norm *next_lowband_out1=NULL; | ||
979 | int next_level=0; | ||
980 | opus_int32 rebalance; | ||
981 | |||
982 | /* Give more bits to low-energy MDCTs than they would otherwise deserve */ | ||
983 | if (B0>1 && !stereo && (itheta&0x3fff)) | ||
984 | { | ||
985 | if (itheta > 8192) | ||
986 | /* Rough approximation for pre-echo masking */ | ||
987 | delta -= delta>>(4-LM); | ||
988 | else | ||
989 | /* Corresponds to a forward-masking slope of 1.5 dB per 10 ms */ | ||
990 | delta = IMIN(0, delta + (N<<BITRES>>(5-LM))); | ||
991 | } | ||
992 | mbits = IMAX(0, IMIN(b, (b-delta)/2)); | ||
993 | sbits = b-mbits; | ||
994 | *remaining_bits -= qalloc; | ||
995 | |||
996 | if (lowband && !stereo) | ||
997 | next_lowband2 = lowband+N; /* >32-bit split case */ | ||
998 | |||
999 | /* Only stereo needs to pass on lowband_out. Otherwise, it's | ||
1000 | handled at the end */ | ||
1001 | if (stereo) | ||
1002 | next_lowband_out1 = lowband_out; | ||
1003 | else | ||
1004 | next_level = level+1; | ||
1005 | |||
1006 | rebalance = *remaining_bits; | ||
1007 | if (mbits >= sbits) | ||
1008 | { | ||
1009 | /* In stereo mode, we do not apply a scaling to the mid because we need the normalized | ||
1010 | mid for folding later */ | ||
1011 | cm = quant_band(encode, m, i, X, NULL, N, mbits, spread, B, intensity, tf_change, | ||
1012 | lowband, ec, remaining_bits, LM, next_lowband_out1, | ||
1013 | NULL, next_level, seed, stereo ? Q15ONE : MULT16_16_P15(gain,mid), lowband_scratch, fill); | ||
1014 | rebalance = mbits - (rebalance-*remaining_bits); | ||
1015 | if (rebalance > 3<<BITRES && itheta!=0) | ||
1016 | sbits += rebalance - (3<<BITRES); | ||
1017 | |||
1018 | /* For a stereo split, the high bits of fill are always zero, so no | ||
1019 | folding will be done to the side. */ | ||
1020 | cm |= quant_band(encode, m, i, Y, NULL, N, sbits, spread, B, intensity, tf_change, | ||
1021 | next_lowband2, ec, remaining_bits, LM, NULL, | ||
1022 | NULL, next_level, seed, MULT16_16_P15(gain,side), NULL, fill>>B)<<((B0>>1)&(stereo-1)); | ||
1023 | } else { | ||
1024 | /* For a stereo split, the high bits of fill are always zero, so no | ||
1025 | folding will be done to the side. */ | ||
1026 | cm = quant_band(encode, m, i, Y, NULL, N, sbits, spread, B, intensity, tf_change, | ||
1027 | next_lowband2, ec, remaining_bits, LM, NULL, | ||
1028 | NULL, next_level, seed, MULT16_16_P15(gain,side), NULL, fill>>B)<<((B0>>1)&(stereo-1)); | ||
1029 | rebalance = sbits - (rebalance-*remaining_bits); | ||
1030 | if (rebalance > 3<<BITRES && itheta!=16384) | ||
1031 | mbits += rebalance - (3<<BITRES); | ||
1032 | /* In stereo mode, we do not apply a scaling to the mid because we need the normalized | ||
1033 | mid for folding later */ | ||
1034 | cm |= quant_band(encode, m, i, X, NULL, N, mbits, spread, B, intensity, tf_change, | ||
1035 | lowband, ec, remaining_bits, LM, next_lowband_out1, | ||
1036 | NULL, next_level, seed, stereo ? Q15ONE : MULT16_16_P15(gain,mid), lowband_scratch, fill); | ||
1037 | } | ||
1038 | } | ||
1039 | |||
1040 | } else { | ||
1041 | /* This is the basic no-split case */ | ||
1042 | q = bits2pulses(m, i, LM, b); | ||
1043 | curr_bits = pulses2bits(m, i, LM, q); | ||
1044 | *remaining_bits -= curr_bits; | ||
1045 | |||
1046 | /* Ensures we can never bust the budget */ | ||
1047 | while (*remaining_bits < 0 && q > 0) | ||
1048 | { | ||
1049 | *remaining_bits += curr_bits; | ||
1050 | q--; | ||
1051 | curr_bits = pulses2bits(m, i, LM, q); | ||
1052 | *remaining_bits -= curr_bits; | ||
1053 | } | ||
1054 | |||
1055 | if (q!=0) | ||
1056 | { | ||
1057 | int K = get_pulses(q); | ||
1058 | |||
1059 | /* Finally do the actual quantization */ | ||
1060 | if (encode) | ||
1061 | { | ||
1062 | cm = alg_quant(X, N, K, spread, B, ec | ||
1063 | #ifdef RESYNTH | ||
1064 | , gain | ||
1065 | #endif | ||
1066 | ); | ||
1067 | } else { | ||
1068 | cm = alg_unquant(X, N, K, spread, B, ec, gain); | ||
1069 | } | ||
1070 | } else { | ||
1071 | /* If there's no pulse, fill the band anyway */ | ||
1072 | int j; | ||
1073 | if (resynth) | ||
1074 | { | ||
1075 | unsigned cm_mask; | ||
1076 | /*B can be as large as 16, so this shift might overflow an int on a | ||
1077 | 16-bit platform; use a long to get defined behavior.*/ | ||
1078 | cm_mask = (unsigned)(1UL<<B)-1; | ||
1079 | fill &= cm_mask; | ||
1080 | if (!fill) | ||
1081 | { | ||
1082 | for (j=0;j<N;j++) | ||
1083 | X[j] = 0; | ||
1084 | } else { | ||
1085 | if (lowband == NULL) | ||
1086 | { | ||
1087 | /* Noise */ | ||
1088 | for (j=0;j<N;j++) | ||
1089 | { | ||
1090 | *seed = celt_lcg_rand(*seed); | ||
1091 | X[j] = (celt_norm)((opus_int32)*seed>>20); | ||
1092 | } | ||
1093 | cm = cm_mask; | ||
1094 | } else { | ||
1095 | /* Folded spectrum */ | ||
1096 | for (j=0;j<N;j++) | ||
1097 | { | ||
1098 | opus_val16 tmp; | ||
1099 | *seed = celt_lcg_rand(*seed); | ||
1100 | /* About 48 dB below the "normal" folding level */ | ||
1101 | tmp = QCONST16(1.0f/256, 10); | ||
1102 | tmp = (*seed)&0x8000 ? tmp : -tmp; | ||
1103 | X[j] = lowband[j]+tmp; | ||
1104 | } | ||
1105 | cm = fill; | ||
1106 | } | ||
1107 | renormalise_vector(X, N, gain); | ||
1108 | } | ||
1109 | } | ||
1110 | } | ||
1111 | } | ||
1112 | |||
1113 | /* This code is used by the decoder and by the resynthesis-enabled encoder */ | ||
1114 | if (resynth) | ||
1115 | { | ||
1116 | if (stereo) | ||
1117 | { | ||
1118 | if (N!=2) | ||
1119 | stereo_merge(X, Y, mid, N); | ||
1120 | if (inv) | ||
1121 | { | ||
1122 | int j; | ||
1123 | for (j=0;j<N;j++) | ||
1124 | Y[j] = -Y[j]; | ||
1125 | } | ||
1126 | } else if (level == 0) | ||
1127 | { | ||
1128 | int k; | ||
1129 | |||
1130 | /* Undo the sample reorganization going from time order to frequency order */ | ||
1131 | if (B0>1) | ||
1132 | interleave_hadamard(X, N_B>>recombine, B0<<recombine, longBlocks); | ||
1133 | |||
1134 | /* Undo time-freq changes that we did earlier */ | ||
1135 | N_B = N_B0; | ||
1136 | B = B0; | ||
1137 | for (k=0;k<time_divide;k++) | ||
1138 | { | ||
1139 | B >>= 1; | ||
1140 | N_B <<= 1; | ||
1141 | cm |= cm>>B; | ||
1142 | haar1(X, N_B, B); | ||
1143 | } | ||
1144 | |||
1145 | for (k=0;k<recombine;k++) | ||
1146 | { | ||
1147 | static const unsigned char bit_deinterleave_table[16]={ | ||
1148 | 0x00,0x03,0x0C,0x0F,0x30,0x33,0x3C,0x3F, | ||
1149 | 0xC0,0xC3,0xCC,0xCF,0xF0,0xF3,0xFC,0xFF | ||
1150 | }; | ||
1151 | cm = bit_deinterleave_table[cm]; | ||
1152 | haar1(X, N0>>k, 1<<k); | ||
1153 | } | ||
1154 | B<<=recombine; | ||
1155 | |||
1156 | /* Scale output for later folding */ | ||
1157 | if (lowband_out) | ||
1158 | { | ||
1159 | int j; | ||
1160 | opus_val16 n; | ||
1161 | n = celt_sqrt(SHL32(EXTEND32(N0),22)); | ||
1162 | for (j=0;j<N0;j++) | ||
1163 | lowband_out[j] = MULT16_16_Q15(n,X[j]); | ||
1164 | } | ||
1165 | cm &= (1<<B)-1; | ||
1166 | } | ||
1167 | } | ||
1168 | return cm; | ||
1169 | } | ||
1170 | |||
1171 | void quant_all_bands(int encode, const CELTMode *m, int start, int end, | ||
1172 | celt_norm *X_, celt_norm *Y_, unsigned char *collapse_masks, const celt_ener *bandE, int *pulses, | ||
1173 | int shortBlocks, int spread, int dual_stereo, int intensity, int *tf_res, | ||
1174 | opus_int32 total_bits, opus_int32 balance, ec_ctx *ec, int LM, int codedBands, opus_uint32 *seed) | ||
1175 | { | ||
1176 | int i; | ||
1177 | opus_int32 remaining_bits; | ||
1178 | const opus_int16 * OPUS_RESTRICT eBands = m->eBands; | ||
1179 | celt_norm * OPUS_RESTRICT norm, * OPUS_RESTRICT norm2; | ||
1180 | VARDECL(celt_norm, _norm); | ||
1181 | VARDECL(celt_norm, lowband_scratch); | ||
1182 | int B; | ||
1183 | int M; | ||
1184 | int lowband_offset; | ||
1185 | int update_lowband = 1; | ||
1186 | int C = Y_ != NULL ? 2 : 1; | ||
1187 | #ifdef RESYNTH | ||
1188 | int resynth = 1; | ||
1189 | #else | ||
1190 | int resynth = !encode; | ||
1191 | #endif | ||
1192 | SAVE_STACK; | ||
1193 | |||
1194 | M = 1<<LM; | ||
1195 | B = shortBlocks ? M : 1; | ||
1196 | ALLOC(_norm, C*M*eBands[m->nbEBands], celt_norm); | ||
1197 | ALLOC(lowband_scratch, M*(eBands[m->nbEBands]-eBands[m->nbEBands-1]), celt_norm); | ||
1198 | norm = _norm; | ||
1199 | norm2 = norm + M*eBands[m->nbEBands]; | ||
1200 | |||
1201 | lowband_offset = 0; | ||
1202 | for (i=start;i<end;i++) | ||
1203 | { | ||
1204 | opus_int32 tell; | ||
1205 | int b; | ||
1206 | int N; | ||
1207 | opus_int32 curr_balance; | ||
1208 | int effective_lowband=-1; | ||
1209 | celt_norm * OPUS_RESTRICT X, * OPUS_RESTRICT Y; | ||
1210 | int tf_change=0; | ||
1211 | unsigned x_cm; | ||
1212 | unsigned y_cm; | ||
1213 | |||
1214 | X = X_+M*eBands[i]; | ||
1215 | if (Y_!=NULL) | ||
1216 | Y = Y_+M*eBands[i]; | ||
1217 | else | ||
1218 | Y = NULL; | ||
1219 | N = M*eBands[i+1]-M*eBands[i]; | ||
1220 | tell = ec_tell_frac(ec); | ||
1221 | |||
1222 | /* Compute how many bits we want to allocate to this band */ | ||
1223 | if (i != start) | ||
1224 | balance -= tell; | ||
1225 | remaining_bits = total_bits-tell-1; | ||
1226 | if (i <= codedBands-1) | ||
1227 | { | ||
1228 | curr_balance = balance / IMIN(3, codedBands-i); | ||
1229 | b = IMAX(0, IMIN(16383, IMIN(remaining_bits+1,pulses[i]+curr_balance))); | ||
1230 | } else { | ||
1231 | b = 0; | ||
1232 | } | ||
1233 | |||
1234 | if (resynth && M*eBands[i]-N >= M*eBands[start] && (update_lowband || lowband_offset==0)) | ||
1235 | lowband_offset = i; | ||
1236 | |||
1237 | tf_change = tf_res[i]; | ||
1238 | if (i>=m->effEBands) | ||
1239 | { | ||
1240 | X=norm; | ||
1241 | if (Y_!=NULL) | ||
1242 | Y = norm; | ||
1243 | } | ||
1244 | |||
1245 | /* Get a conservative estimate of the collapse_mask's for the bands we're | ||
1246 | going to be folding from. */ | ||
1247 | if (lowband_offset != 0 && (spread!=SPREAD_AGGRESSIVE || B>1 || tf_change<0)) | ||
1248 | { | ||
1249 | int fold_start; | ||
1250 | int fold_end; | ||
1251 | int fold_i; | ||
1252 | /* This ensures we never repeat spectral content within one band */ | ||
1253 | effective_lowband = IMAX(M*eBands[start], M*eBands[lowband_offset]-N); | ||
1254 | fold_start = lowband_offset; | ||
1255 | while(M*eBands[--fold_start] > effective_lowband); | ||
1256 | fold_end = lowband_offset-1; | ||
1257 | while(M*eBands[++fold_end] < effective_lowband+N); | ||
1258 | x_cm = y_cm = 0; | ||
1259 | fold_i = fold_start; do { | ||
1260 | x_cm |= collapse_masks[fold_i*C+0]; | ||
1261 | y_cm |= collapse_masks[fold_i*C+C-1]; | ||
1262 | } while (++fold_i<fold_end); | ||
1263 | } | ||
1264 | /* Otherwise, we'll be using the LCG to fold, so all blocks will (almost | ||
1265 | always) be non-zero.*/ | ||
1266 | else | ||
1267 | x_cm = y_cm = (1<<B)-1; | ||
1268 | |||
1269 | if (dual_stereo && i==intensity) | ||
1270 | { | ||
1271 | int j; | ||
1272 | |||
1273 | /* Switch off dual stereo to do intensity */ | ||
1274 | dual_stereo = 0; | ||
1275 | if (resynth) | ||
1276 | for (j=M*eBands[start];j<M*eBands[i];j++) | ||
1277 | norm[j] = HALF32(norm[j]+norm2[j]); | ||
1278 | } | ||
1279 | if (dual_stereo) | ||
1280 | { | ||
1281 | x_cm = quant_band(encode, m, i, X, NULL, N, b/2, spread, B, intensity, tf_change, | ||
1282 | effective_lowband != -1 ? norm+effective_lowband : NULL, ec, &remaining_bits, LM, | ||
1283 | norm+M*eBands[i], bandE, 0, seed, Q15ONE, lowband_scratch, x_cm); | ||
1284 | y_cm = quant_band(encode, m, i, Y, NULL, N, b/2, spread, B, intensity, tf_change, | ||
1285 | effective_lowband != -1 ? norm2+effective_lowband : NULL, ec, &remaining_bits, LM, | ||
1286 | norm2+M*eBands[i], bandE, 0, seed, Q15ONE, lowband_scratch, y_cm); | ||
1287 | } else { | ||
1288 | x_cm = quant_band(encode, m, i, X, Y, N, b, spread, B, intensity, tf_change, | ||
1289 | effective_lowband != -1 ? norm+effective_lowband : NULL, ec, &remaining_bits, LM, | ||
1290 | norm+M*eBands[i], bandE, 0, seed, Q15ONE, lowband_scratch, x_cm|y_cm); | ||
1291 | y_cm = x_cm; | ||
1292 | } | ||
1293 | collapse_masks[i*C+0] = (unsigned char)x_cm; | ||
1294 | collapse_masks[i*C+C-1] = (unsigned char)y_cm; | ||
1295 | balance += pulses[i] + tell; | ||
1296 | |||
1297 | /* Update the folding position only as long as we have 1 bit/sample depth */ | ||
1298 | update_lowband = b>(N<<BITRES); | ||
1299 | } | ||
1300 | RESTORE_STACK; | ||
1301 | } | ||
1302 | |||