2 * Copyright (c) 2012 Andrew D'Addesio
3 * Copyright (c) 2013-2014 Mozilla Corporation
5 * This file is part of FFmpeg.
7 * FFmpeg is free software; you can redistribute it and/or
8 * modify it under the terms of the GNU Lesser General Public
9 * License as published by the Free Software Foundation; either
10 * version 2.1 of the License, or (at your option) any later version.
12 * FFmpeg is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
15 * Lesser General Public License for more details.
17 * You should have received a copy of the GNU Lesser General Public
18 * License along with FFmpeg; if not, write to the Free Software
19 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
29 #include "libavutil/float_dsp.h"
30 #include "libavutil/libm.h"
39 CELT_SPREAD_AGGRESSIVE
42 typedef struct CeltFrame {
43 float energy[CELT_MAX_BANDS];
44 float prev_energy[2][CELT_MAX_BANDS];
46 uint8_t collapse_masks[CELT_MAX_BANDS];
48 /* buffer for mdct output + postfilter */
49 DECLARE_ALIGNED(32, float, buf)[2048];
51 /* postfilter parameters */
53 float pf_gains_new[3];
57 float pf_gains_old[3];
63 // constant values that do not change during context lifetime
64 AVCodecContext *avctx;
65 IMDCT15Context *imdct[4];
66 AVFloatDSPContext *dsp;
69 // values that have inter-frame effect and must be reset on flush
74 // values that only affect a single frame
79 /* number of iMDCT blocks in the frame */
81 /* size of each block */
92 enum CeltSpread spread;
96 int fine_bits [CELT_MAX_BANDS];
97 int fine_priority[CELT_MAX_BANDS];
98 int pulses [CELT_MAX_BANDS];
99 int tf_change [CELT_MAX_BANDS];
101 DECLARE_ALIGNED(32, float, coeffs)[2][CELT_MAX_FRAME_SIZE];
102 DECLARE_ALIGNED(32, float, scratch)[22 * 8]; // MAX(celt_freq_range) * 1<<CELT_MAX_LOG_BLOCKS
105 static const uint16_t celt_model_tapset[] = { 4, 2, 3, 4 };
107 static const uint16_t celt_model_spread[] = { 32, 7, 9, 30, 32 };
109 static const uint16_t celt_model_alloc_trim[] = {
110 128, 2, 4, 9, 19, 41, 87, 109, 119, 124, 126, 128
113 static const uint16_t celt_model_energy_small[] = { 4, 2, 3, 4 };
115 static const uint8_t celt_freq_bands[] = { /* in steps of 200Hz */
116 0, 1, 2, 3, 4, 5, 6, 7, 8, 10, 12, 14, 16, 20, 24, 28, 34, 40, 48, 60, 78, 100
119 static const uint8_t celt_freq_range[] = {
120 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 4, 4, 4, 6, 6, 8, 12, 18, 22
123 static const uint8_t celt_log_freq_range[] = {
124 0, 0, 0, 0, 0, 0, 0, 0, 8, 8, 8, 8, 16, 16, 16, 21, 21, 24, 29, 34, 36
127 static const int8_t celt_tf_select[4][2][2][2] = {
128 { { { 0, -1 }, { 0, -1 } }, { { 0, -1 }, { 0, -1 } } },
129 { { { 0, -1 }, { 0, -2 } }, { { 1, 0 }, { 1, -1 } } },
130 { { { 0, -2 }, { 0, -3 } }, { { 2, 0 }, { 1, -1 } } },
131 { { { 0, -2 }, { 0, -3 } }, { { 3, 0 }, { 1, -1 } } }
134 static const float celt_mean_energy[] = {
135 6.437500f, 6.250000f, 5.750000f, 5.312500f, 5.062500f,
136 4.812500f, 4.500000f, 4.375000f, 4.875000f, 4.687500f,
137 4.562500f, 4.437500f, 4.875000f, 4.625000f, 4.312500f,
138 4.500000f, 4.375000f, 4.625000f, 4.750000f, 4.437500f,
139 3.750000f, 3.750000f, 3.750000f, 3.750000f, 3.750000f
142 static const float celt_alpha_coef[] = {
143 29440.0f/32768.0f, 26112.0f/32768.0f, 21248.0f/32768.0f, 16384.0f/32768.0f
146 static const float celt_beta_coef[] = { /* TODO: precompute 1 minus this if the code ends up neater */
147 30147.0f/32768.0f, 22282.0f/32768.0f, 12124.0f/32768.0f, 6554.0f/32768.0f
150 static const uint8_t celt_coarse_energy_dist[4][2][42] = {
152 { // 120-sample inter
153 72, 127, 65, 129, 66, 128, 65, 128, 64, 128, 62, 128, 64, 128,
154 64, 128, 92, 78, 92, 79, 92, 78, 90, 79, 116, 41, 115, 40,
155 114, 40, 132, 26, 132, 26, 145, 17, 161, 12, 176, 10, 177, 11
156 }, { // 120-sample intra
157 24, 179, 48, 138, 54, 135, 54, 132, 53, 134, 56, 133, 55, 132,
158 55, 132, 61, 114, 70, 96, 74, 88, 75, 88, 87, 74, 89, 66,
159 91, 67, 100, 59, 108, 50, 120, 40, 122, 37, 97, 43, 78, 50
162 { // 240-sample inter
163 83, 78, 84, 81, 88, 75, 86, 74, 87, 71, 90, 73, 93, 74,
164 93, 74, 109, 40, 114, 36, 117, 34, 117, 34, 143, 17, 145, 18,
165 146, 19, 162, 12, 165, 10, 178, 7, 189, 6, 190, 8, 177, 9
166 }, { // 240-sample intra
167 23, 178, 54, 115, 63, 102, 66, 98, 69, 99, 74, 89, 71, 91,
168 73, 91, 78, 89, 86, 80, 92, 66, 93, 64, 102, 59, 103, 60,
169 104, 60, 117, 52, 123, 44, 138, 35, 133, 31, 97, 38, 77, 45
172 { // 480-sample inter
173 61, 90, 93, 60, 105, 42, 107, 41, 110, 45, 116, 38, 113, 38,
174 112, 38, 124, 26, 132, 27, 136, 19, 140, 20, 155, 14, 159, 16,
175 158, 18, 170, 13, 177, 10, 187, 8, 192, 6, 175, 9, 159, 10
176 }, { // 480-sample intra
177 21, 178, 59, 110, 71, 86, 75, 85, 84, 83, 91, 66, 88, 73,
178 87, 72, 92, 75, 98, 72, 105, 58, 107, 54, 115, 52, 114, 55,
179 112, 56, 129, 51, 132, 40, 150, 33, 140, 29, 98, 35, 77, 42
182 { // 960-sample inter
183 42, 121, 96, 66, 108, 43, 111, 40, 117, 44, 123, 32, 120, 36,
184 119, 33, 127, 33, 134, 34, 139, 21, 147, 23, 152, 20, 158, 25,
185 154, 26, 166, 21, 173, 16, 184, 13, 184, 10, 150, 13, 139, 15
186 }, { // 960-sample intra
187 22, 178, 63, 114, 74, 82, 84, 83, 92, 82, 103, 62, 96, 72,
188 96, 67, 101, 73, 107, 72, 113, 55, 118, 52, 125, 52, 118, 52,
189 117, 55, 135, 49, 137, 39, 157, 32, 145, 29, 97, 33, 77, 40
194 static const uint8_t celt_static_alloc[11][21] = { /* 1/32 bit/sample */
195 { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
196 { 90, 80, 75, 69, 63, 56, 49, 40, 34, 29, 20, 18, 10, 0, 0, 0, 0, 0, 0, 0, 0 },
197 { 110, 100, 90, 84, 78, 71, 65, 58, 51, 45, 39, 32, 26, 20, 12, 0, 0, 0, 0, 0, 0 },
198 { 118, 110, 103, 93, 86, 80, 75, 70, 65, 59, 53, 47, 40, 31, 23, 15, 4, 0, 0, 0, 0 },
199 { 126, 119, 112, 104, 95, 89, 83, 78, 72, 66, 60, 54, 47, 39, 32, 25, 17, 12, 1, 0, 0 },
200 { 134, 127, 120, 114, 103, 97, 91, 85, 78, 72, 66, 60, 54, 47, 41, 35, 29, 23, 16, 10, 1 },
201 { 144, 137, 130, 124, 113, 107, 101, 95, 88, 82, 76, 70, 64, 57, 51, 45, 39, 33, 26, 15, 1 },
202 { 152, 145, 138, 132, 123, 117, 111, 105, 98, 92, 86, 80, 74, 67, 61, 55, 49, 43, 36, 20, 1 },
203 { 162, 155, 148, 142, 133, 127, 121, 115, 108, 102, 96, 90, 84, 77, 71, 65, 59, 53, 46, 30, 1 },
204 { 172, 165, 158, 152, 143, 137, 131, 125, 118, 112, 106, 100, 94, 87, 81, 75, 69, 63, 56, 45, 20 },
205 { 200, 200, 200, 200, 200, 200, 200, 200, 198, 193, 188, 183, 178, 173, 168, 163, 158, 153, 148, 129, 104 }
208 static const uint8_t celt_static_caps[4][2][21] = {
210 {224, 224, 224, 224, 224, 224, 224, 224, 160, 160,
211 160, 160, 185, 185, 185, 178, 178, 168, 134, 61, 37},
212 {224, 224, 224, 224, 224, 224, 224, 224, 240, 240,
213 240, 240, 207, 207, 207, 198, 198, 183, 144, 66, 40},
215 {160, 160, 160, 160, 160, 160, 160, 160, 185, 185,
216 185, 185, 193, 193, 193, 183, 183, 172, 138, 64, 38},
217 {240, 240, 240, 240, 240, 240, 240, 240, 207, 207,
218 207, 207, 204, 204, 204, 193, 193, 180, 143, 66, 40},
220 {185, 185, 185, 185, 185, 185, 185, 185, 193, 193,
221 193, 193, 193, 193, 193, 183, 183, 172, 138, 65, 39},
222 {207, 207, 207, 207, 207, 207, 207, 207, 204, 204,
223 204, 204, 201, 201, 201, 188, 188, 176, 141, 66, 40},
225 {193, 193, 193, 193, 193, 193, 193, 193, 193, 193,
226 193, 193, 194, 194, 194, 184, 184, 173, 139, 65, 39},
227 {204, 204, 204, 204, 204, 204, 204, 204, 201, 201,
228 201, 201, 198, 198, 198, 187, 187, 175, 140, 66, 40}
232 static const uint8_t celt_cache_bits[392] = {
233 40, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7,
234 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7,
235 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 40, 15, 23, 28,
236 31, 34, 36, 38, 39, 41, 42, 43, 44, 45, 46, 47, 47, 49, 50,
237 51, 52, 53, 54, 55, 55, 57, 58, 59, 60, 61, 62, 63, 63, 65,
238 66, 67, 68, 69, 70, 71, 71, 40, 20, 33, 41, 48, 53, 57, 61,
239 64, 66, 69, 71, 73, 75, 76, 78, 80, 82, 85, 87, 89, 91, 92,
240 94, 96, 98, 101, 103, 105, 107, 108, 110, 112, 114, 117, 119, 121, 123,
241 124, 126, 128, 40, 23, 39, 51, 60, 67, 73, 79, 83, 87, 91, 94,
242 97, 100, 102, 105, 107, 111, 115, 118, 121, 124, 126, 129, 131, 135, 139,
243 142, 145, 148, 150, 153, 155, 159, 163, 166, 169, 172, 174, 177, 179, 35,
244 28, 49, 65, 78, 89, 99, 107, 114, 120, 126, 132, 136, 141, 145, 149,
245 153, 159, 165, 171, 176, 180, 185, 189, 192, 199, 205, 211, 216, 220, 225,
246 229, 232, 239, 245, 251, 21, 33, 58, 79, 97, 112, 125, 137, 148, 157,
247 166, 174, 182, 189, 195, 201, 207, 217, 227, 235, 243, 251, 17, 35, 63,
248 86, 106, 123, 139, 152, 165, 177, 187, 197, 206, 214, 222, 230, 237, 250,
249 25, 31, 55, 75, 91, 105, 117, 128, 138, 146, 154, 161, 168, 174, 180,
250 185, 190, 200, 208, 215, 222, 229, 235, 240, 245, 255, 16, 36, 65, 89,
251 110, 128, 144, 159, 173, 185, 196, 207, 217, 226, 234, 242, 250, 11, 41,
252 74, 103, 128, 151, 172, 191, 209, 225, 241, 255, 9, 43, 79, 110, 138,
253 163, 186, 207, 227, 246, 12, 39, 71, 99, 123, 144, 164, 182, 198, 214,
254 228, 241, 253, 9, 44, 81, 113, 142, 168, 192, 214, 235, 255, 7, 49,
255 90, 127, 160, 191, 220, 247, 6, 51, 95, 134, 170, 203, 234, 7, 47,
256 87, 123, 155, 184, 212, 237, 6, 52, 97, 137, 174, 208, 240, 5, 57,
257 106, 151, 192, 231, 5, 59, 111, 158, 202, 243, 5, 55, 103, 147, 187,
258 224, 5, 60, 113, 161, 206, 248, 4, 65, 122, 175, 224, 4, 67, 127,
262 static const int16_t celt_cache_index[105] = {
263 -1, -1, -1, -1, -1, -1, -1, -1, 0, 0, 0, 0, 41, 41, 41,
264 82, 82, 123, 164, 200, 222, 0, 0, 0, 0, 0, 0, 0, 0, 41,
265 41, 41, 41, 123, 123, 123, 164, 164, 240, 266, 283, 295, 41, 41, 41,
266 41, 41, 41, 41, 41, 123, 123, 123, 123, 240, 240, 240, 266, 266, 305,
267 318, 328, 336, 123, 123, 123, 123, 123, 123, 123, 123, 240, 240, 240, 240,
268 305, 305, 305, 318, 318, 343, 351, 358, 364, 240, 240, 240, 240, 240, 240,
269 240, 240, 305, 305, 305, 305, 343, 343, 343, 351, 351, 370, 376, 382, 387,
272 static const uint8_t celt_log2_frac[] = {
273 0, 8, 13, 16, 19, 21, 23, 24, 26, 27, 28, 29, 30, 31, 32, 32, 33, 34, 34, 35, 36, 36, 37, 37
276 static const uint8_t celt_bit_interleave[] = {
277 0, 1, 1, 1, 2, 3, 3, 3, 2, 3, 3, 3, 2, 3, 3, 3
280 static const uint8_t celt_bit_deinterleave[] = {
281 0x00, 0x03, 0x0C, 0x0F, 0x30, 0x33, 0x3C, 0x3F,
282 0xC0, 0xC3, 0xCC, 0xCF, 0xF0, 0xF3, 0xFC, 0xFF
285 static const uint8_t celt_hadamard_ordery[] = {
288 7, 0, 4, 3, 6, 1, 5, 2,
289 15, 0, 8, 7, 12, 3, 11, 4, 14, 1, 9, 6, 13, 2, 10, 5
292 static const uint16_t celt_qn_exp2[] = {
293 16384, 17866, 19483, 21247, 23170, 25267, 27554, 30048
296 static const uint32_t celt_pvq_u[1272] = {
297 /* N = 0, K = 0...176 */
298 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
299 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
300 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
301 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
302 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
303 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
304 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
305 /* N = 1, K = 1...176 */
306 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
307 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
308 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
309 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
310 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
311 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
312 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
313 /* N = 2, K = 2...176 */
314 3, 5, 7, 9, 11, 13, 15, 17, 19, 21, 23, 25, 27, 29, 31, 33, 35, 37, 39, 41,
315 43, 45, 47, 49, 51, 53, 55, 57, 59, 61, 63, 65, 67, 69, 71, 73, 75, 77, 79,
316 81, 83, 85, 87, 89, 91, 93, 95, 97, 99, 101, 103, 105, 107, 109, 111, 113,
317 115, 117, 119, 121, 123, 125, 127, 129, 131, 133, 135, 137, 139, 141, 143,
318 145, 147, 149, 151, 153, 155, 157, 159, 161, 163, 165, 167, 169, 171, 173,
319 175, 177, 179, 181, 183, 185, 187, 189, 191, 193, 195, 197, 199, 201, 203,
320 205, 207, 209, 211, 213, 215, 217, 219, 221, 223, 225, 227, 229, 231, 233,
321 235, 237, 239, 241, 243, 245, 247, 249, 251, 253, 255, 257, 259, 261, 263,
322 265, 267, 269, 271, 273, 275, 277, 279, 281, 283, 285, 287, 289, 291, 293,
323 295, 297, 299, 301, 303, 305, 307, 309, 311, 313, 315, 317, 319, 321, 323,
324 325, 327, 329, 331, 333, 335, 337, 339, 341, 343, 345, 347, 349, 351,
325 /* N = 3, K = 3...176 */
326 13, 25, 41, 61, 85, 113, 145, 181, 221, 265, 313, 365, 421, 481, 545, 613,
327 685, 761, 841, 925, 1013, 1105, 1201, 1301, 1405, 1513, 1625, 1741, 1861,
328 1985, 2113, 2245, 2381, 2521, 2665, 2813, 2965, 3121, 3281, 3445, 3613, 3785,
329 3961, 4141, 4325, 4513, 4705, 4901, 5101, 5305, 5513, 5725, 5941, 6161, 6385,
330 6613, 6845, 7081, 7321, 7565, 7813, 8065, 8321, 8581, 8845, 9113, 9385, 9661,
331 9941, 10225, 10513, 10805, 11101, 11401, 11705, 12013, 12325, 12641, 12961,
332 13285, 13613, 13945, 14281, 14621, 14965, 15313, 15665, 16021, 16381, 16745,
333 17113, 17485, 17861, 18241, 18625, 19013, 19405, 19801, 20201, 20605, 21013,
334 21425, 21841, 22261, 22685, 23113, 23545, 23981, 24421, 24865, 25313, 25765,
335 26221, 26681, 27145, 27613, 28085, 28561, 29041, 29525, 30013, 30505, 31001,
336 31501, 32005, 32513, 33025, 33541, 34061, 34585, 35113, 35645, 36181, 36721,
337 37265, 37813, 38365, 38921, 39481, 40045, 40613, 41185, 41761, 42341, 42925,
338 43513, 44105, 44701, 45301, 45905, 46513, 47125, 47741, 48361, 48985, 49613,
339 50245, 50881, 51521, 52165, 52813, 53465, 54121, 54781, 55445, 56113, 56785,
340 57461, 58141, 58825, 59513, 60205, 60901, 61601,
341 /* N = 4, K = 4...176 */
342 63, 129, 231, 377, 575, 833, 1159, 1561, 2047, 2625, 3303, 4089, 4991, 6017,
343 7175, 8473, 9919, 11521, 13287, 15225, 17343, 19649, 22151, 24857, 27775,
344 30913, 34279, 37881, 41727, 45825, 50183, 54809, 59711, 64897, 70375, 76153,
345 82239, 88641, 95367, 102425, 109823, 117569, 125671, 134137, 142975, 152193,
346 161799, 171801, 182207, 193025, 204263, 215929, 228031, 240577, 253575,
347 267033, 280959, 295361, 310247, 325625, 341503, 357889, 374791, 392217,
348 410175, 428673, 447719, 467321, 487487, 508225, 529543, 551449, 573951,
349 597057, 620775, 645113, 670079, 695681, 721927, 748825, 776383, 804609,
350 833511, 863097, 893375, 924353, 956039, 988441, 1021567, 1055425, 1090023,
351 1125369, 1161471, 1198337, 1235975, 1274393, 1313599, 1353601, 1394407,
352 1436025, 1478463, 1521729, 1565831, 1610777, 1656575, 1703233, 1750759,
353 1799161, 1848447, 1898625, 1949703, 2001689, 2054591, 2108417, 2163175,
354 2218873, 2275519, 2333121, 2391687, 2451225, 2511743, 2573249, 2635751,
355 2699257, 2763775, 2829313, 2895879, 2963481, 3032127, 3101825, 3172583,
356 3244409, 3317311, 3391297, 3466375, 3542553, 3619839, 3698241, 3777767,
357 3858425, 3940223, 4023169, 4107271, 4192537, 4278975, 4366593, 4455399,
358 4545401, 4636607, 4729025, 4822663, 4917529, 5013631, 5110977, 5209575,
359 5309433, 5410559, 5512961, 5616647, 5721625, 5827903, 5935489, 6044391,
360 6154617, 6266175, 6379073, 6493319, 6608921, 6725887, 6844225, 6963943,
362 /* N = 5, K = 5...176 */
363 321, 681, 1289, 2241, 3649, 5641, 8361, 11969, 16641, 22569, 29961, 39041,
364 50049, 63241, 78889, 97281, 118721, 143529, 172041, 204609, 241601, 283401,
365 330409, 383041, 441729, 506921, 579081, 658689, 746241, 842249, 947241,
366 1061761, 1186369, 1321641, 1468169, 1626561, 1797441, 1981449, 2179241,
367 2391489, 2618881, 2862121, 3121929, 3399041, 3694209, 4008201, 4341801,
368 4695809, 5071041, 5468329, 5888521, 6332481, 6801089, 7295241, 7815849,
369 8363841, 8940161, 9545769, 10181641, 10848769, 11548161, 12280841, 13047849,
370 13850241, 14689089, 15565481, 16480521, 17435329, 18431041, 19468809,
371 20549801, 21675201, 22846209, 24064041, 25329929, 26645121, 28010881,
372 29428489, 30899241, 32424449, 34005441, 35643561, 37340169, 39096641,
373 40914369, 42794761, 44739241, 46749249, 48826241, 50971689, 53187081,
374 55473921, 57833729, 60268041, 62778409, 65366401, 68033601, 70781609,
375 73612041, 76526529, 79526721, 82614281, 85790889, 89058241, 92418049,
376 95872041, 99421961, 103069569, 106816641, 110664969, 114616361, 118672641,
377 122835649, 127107241, 131489289, 135983681, 140592321, 145317129, 150160041,
378 155123009, 160208001, 165417001, 170752009, 176215041, 181808129, 187533321,
379 193392681, 199388289, 205522241, 211796649, 218213641, 224775361, 231483969,
380 238341641, 245350569, 252512961, 259831041, 267307049, 274943241, 282741889,
381 290705281, 298835721, 307135529, 315607041, 324252609, 333074601, 342075401,
382 351257409, 360623041, 370174729, 379914921, 389846081, 399970689, 410291241,
383 420810249, 431530241, 442453761, 453583369, 464921641, 476471169, 488234561,
384 500214441, 512413449, 524834241, 537479489, 550351881, 563454121, 576788929,
385 590359041, 604167209, 618216201, 632508801,
386 /* N = 6, K = 6...96 (technically V(109,5) fits in 32 bits, but that can't be
387 achieved by splitting an Opus band) */
388 1683, 3653, 7183, 13073, 22363, 36365, 56695, 85305, 124515, 177045, 246047,
389 335137, 448427, 590557, 766727, 982729, 1244979, 1560549, 1937199, 2383409,
390 2908411, 3522221, 4235671, 5060441, 6009091, 7095093, 8332863, 9737793,
391 11326283, 13115773, 15124775, 17372905, 19880915, 22670725, 25765455,
392 29189457, 32968347, 37129037, 41699767, 46710137, 52191139, 58175189,
393 64696159, 71789409, 79491819, 87841821, 96879431, 106646281, 117185651,
394 128542501, 140763503, 153897073, 167993403, 183104493, 199284183, 216588185,
395 235074115, 254801525, 275831935, 298228865, 322057867, 347386557, 374284647,
396 402823977, 433078547, 465124549, 499040399, 534906769, 572806619, 612825229,
397 655050231, 699571641, 746481891, 795875861, 847850911, 902506913, 959946283,
398 1020274013, 1083597703, 1150027593, 1219676595, 1292660325, 1369097135,
399 1449108145, 1532817275, 1620351277, 1711839767, 1807415257, 1907213187,
400 2011371957, 2120032959,
401 /* N = 7, K = 7...54 (technically V(60,6) fits in 32 bits, but that can't be
402 achieved by splitting an Opus band) */
403 8989, 19825, 40081, 75517, 134245, 227305, 369305, 579125, 880685, 1303777,
404 1884961, 2668525, 3707509, 5064793, 6814249, 9041957, 11847485, 15345233,
405 19665841, 24957661, 31388293, 39146185, 48442297, 59511829, 72616013,
406 88043969, 106114625, 127178701, 151620757, 179861305, 212358985, 249612805,
407 292164445, 340600625, 395555537, 457713341, 527810725, 606639529, 695049433,
408 793950709, 904317037, 1027188385, 1163673953, 1314955181, 1482288821,
409 1667010073, 1870535785, 2094367717,
410 /* N = 8, K = 8...37 (technically V(40,7) fits in 32 bits, but that can't be
411 achieved by splitting an Opus band) */
412 48639, 108545, 224143, 433905, 795455, 1392065, 2340495, 3800305, 5984767,
413 9173505, 13726991, 20103025, 28875327, 40754369, 56610575, 77500017,
414 104692735, 139703809, 184327311, 240673265, 311207743, 398796225, 506750351,
415 638878193, 799538175, 993696769, 1226990095, 1505789553, 1837271615,
417 /* N = 9, K = 9...28 (technically V(29,8) fits in 32 bits, but that can't be
418 achieved by splitting an Opus band) */
419 265729, 598417, 1256465, 2485825, 4673345, 8405905, 14546705, 24331777,
420 39490049, 62390545, 96220561, 145198913, 214828609, 312193553, 446304145,
421 628496897, 872893441, 1196924561, 1621925137, 2173806145,
422 /* N = 10, K = 10...24 */
423 1462563, 3317445, 7059735, 14218905, 27298155, 50250765, 89129247, 152951073,
424 254831667, 413442773, 654862247, 1014889769, 1541911931, 2300409629,
426 /* N = 11, K = 11...19 (technically V(20,10) fits in 32 bits, but that can't be
427 achieved by splitting an Opus band) */
428 8097453, 18474633, 39753273, 81270333, 158819253, 298199265, 540279585,
429 948062325, 1616336765,
430 /* N = 12, K = 12...18 */
431 45046719, 103274625, 224298231, 464387817, 921406335, 1759885185,
433 /* N = 13, K = 13...16 */
434 251595969, 579168825, 1267854873, 2653649025,
439 DECLARE_ALIGNED(32, static const float, celt_window)[120] = {
440 6.7286966e-05f, 0.00060551348f, 0.0016815970f, 0.0032947962f, 0.0054439943f,
441 0.0081276923f, 0.011344001f, 0.015090633f, 0.019364886f, 0.024163635f,
442 0.029483315f, 0.035319905f, 0.041668911f, 0.048525347f, 0.055883718f,
443 0.063737999f, 0.072081616f, 0.080907428f, 0.090207705f, 0.099974111f,
444 0.11019769f, 0.12086883f, 0.13197729f, 0.14351214f, 0.15546177f,
445 0.16781389f, 0.18055550f, 0.19367290f, 0.20715171f, 0.22097682f,
446 0.23513243f, 0.24960208f, 0.26436860f, 0.27941419f, 0.29472040f,
447 0.31026818f, 0.32603788f, 0.34200931f, 0.35816177f, 0.37447407f,
448 0.39092462f, 0.40749142f, 0.42415215f, 0.44088423f, 0.45766484f,
449 0.47447104f, 0.49127978f, 0.50806798f, 0.52481261f, 0.54149077f,
450 0.55807973f, 0.57455701f, 0.59090049f, 0.60708841f, 0.62309951f,
451 0.63891306f, 0.65450896f, 0.66986776f, 0.68497077f, 0.69980010f,
452 0.71433873f, 0.72857055f, 0.74248043f, 0.75605424f, 0.76927895f,
453 0.78214257f, 0.79463430f, 0.80674445f, 0.81846456f, 0.82978733f,
454 0.84070669f, 0.85121779f, 0.86131698f, 0.87100183f, 0.88027111f,
455 0.88912479f, 0.89756398f, 0.90559094f, 0.91320904f, 0.92042270f,
456 0.92723738f, 0.93365955f, 0.93969656f, 0.94535671f, 0.95064907f,
457 0.95558353f, 0.96017067f, 0.96442171f, 0.96834849f, 0.97196334f,
458 0.97527906f, 0.97830883f, 0.98106616f, 0.98356480f, 0.98581869f,
459 0.98784191f, 0.98964856f, 0.99125274f, 0.99266849f, 0.99390969f,
460 0.99499004f, 0.99592297f, 0.99672162f, 0.99739874f, 0.99796667f,
461 0.99843728f, 0.99882195f, 0.99913147f, 0.99937606f, 0.99956527f,
462 0.99970802f, 0.99981248f, 0.99988613f, 0.99993565f, 0.99996697f,
463 0.99998518f, 0.99999457f, 0.99999859f, 0.99999982f, 1.0000000f,
466 /* square of the window, used for the postfilter */
467 const float ff_celt_window2[120] = {
468 4.5275357e-09f, 3.66647e-07f, 2.82777e-06f, 1.08557e-05f, 2.96371e-05f, 6.60594e-05f,
469 0.000128686f, 0.000227727f, 0.000374999f, 0.000583881f, 0.000869266f, 0.0012475f,
470 0.0017363f, 0.00235471f, 0.00312299f, 0.00406253f, 0.00519576f, 0.00654601f,
471 0.00813743f, 0.00999482f, 0.0121435f, 0.0146093f, 0.017418f, 0.0205957f, 0.0241684f,
472 0.0281615f, 0.0326003f, 0.0375092f, 0.0429118f, 0.0488308f, 0.0552873f, 0.0623012f,
473 0.0698908f, 0.0780723f, 0.0868601f, 0.0962664f, 0.106301f, 0.11697f, 0.12828f,
474 0.140231f, 0.152822f, 0.166049f, 0.179905f, 0.194379f, 0.209457f, 0.225123f, 0.241356f,
475 0.258133f, 0.275428f, 0.293212f, 0.311453f, 0.330116f, 0.349163f, 0.368556f, 0.388253f,
476 0.40821f, 0.428382f, 0.448723f, 0.469185f, 0.48972f, 0.51028f, 0.530815f, 0.551277f,
477 0.571618f, 0.59179f, 0.611747f, 0.631444f, 0.650837f, 0.669884f, 0.688547f, 0.706788f,
478 0.724572f, 0.741867f, 0.758644f, 0.774877f, 0.790543f, 0.805621f, 0.820095f, 0.833951f,
479 0.847178f, 0.859769f, 0.87172f, 0.88303f, 0.893699f, 0.903734f, 0.91314f, 0.921928f,
480 0.930109f, 0.937699f, 0.944713f, 0.951169f, 0.957088f, 0.962491f, 0.9674f, 0.971838f,
481 0.975832f, 0.979404f, 0.982582f, 0.985391f, 0.987857f, 0.990005f, 0.991863f, 0.993454f,
482 0.994804f, 0.995937f, 0.996877f, 0.997645f, 0.998264f, 0.998753f, 0.999131f, 0.999416f,
483 0.999625f, 0.999772f, 0.999871f, 0.999934f, 0.99997f, 0.999989f, 0.999997f, 0.99999964f, 1.0f,
486 static const uint32_t * const celt_pvq_u_row[15] = {
487 celt_pvq_u + 0, celt_pvq_u + 176, celt_pvq_u + 351,
488 celt_pvq_u + 525, celt_pvq_u + 698, celt_pvq_u + 870,
489 celt_pvq_u + 1041, celt_pvq_u + 1131, celt_pvq_u + 1178,
490 celt_pvq_u + 1207, celt_pvq_u + 1226, celt_pvq_u + 1240,
491 celt_pvq_u + 1248, celt_pvq_u + 1254, celt_pvq_u + 1257
494 static inline int16_t celt_cos(int16_t x)
496 x = (MUL16(x, x) + 4096) >> 13;
497 x = (32767-x) + ROUND_MUL16(x, (-7651 + ROUND_MUL16(x, (8277 + ROUND_MUL16(-626, x)))));
501 static inline int celt_log2tan(int isin, int icos)
504 lc = opus_ilog(icos);
505 ls = opus_ilog(isin);
508 return (ls << 11) - (lc << 11) +
509 ROUND_MUL16(isin, ROUND_MUL16(isin, -2597) + 7932) -
510 ROUND_MUL16(icos, ROUND_MUL16(icos, -2597) + 7932);
513 static inline uint32_t celt_rng(CeltContext *s)
515 s->seed = 1664525 * s->seed + 1013904223;
519 static void celt_decode_coarse_energy(CeltContext *s, OpusRangeCoder *rc)
524 const uint8_t *model;
526 /* use the 2D z-transform to apply prediction in both */
527 /* the time domain (alpha) and the frequency domain (beta) */
529 if (opus_rc_tell(rc)+3 <= s->framebits && opus_rc_p2model(rc, 3)) {
532 beta = 1.0f - 4915.0f/32768.0f;
533 model = celt_coarse_energy_dist[s->duration][1];
535 alpha = celt_alpha_coef[s->duration];
536 beta = 1.0f - celt_beta_coef[s->duration];
537 model = celt_coarse_energy_dist[s->duration][0];
540 for (i = 0; i < CELT_MAX_BANDS; i++) {
541 for (j = 0; j < s->coded_channels; j++) {
542 CeltFrame *frame = &s->frame[j];
546 if (i < s->startband || i >= s->endband) {
547 frame->energy[i] = 0.0;
551 available = s->framebits - opus_rc_tell(rc);
552 if (available >= 15) {
553 /* decode using a Laplace distribution */
554 int k = FFMIN(i, 20) << 1;
555 value = opus_rc_laplace(rc, model[k] << 7, model[k+1] << 6);
556 } else if (available >= 2) {
557 int x = opus_rc_getsymbol(rc, celt_model_energy_small);
558 value = (x>>1) ^ -(x&1);
559 } else if (available >= 1) {
560 value = -(float)opus_rc_p2model(rc, 1);
563 frame->energy[i] = FFMAX(-9.0f, frame->energy[i]) * alpha + prev[j] + value;
564 prev[j] += beta * value;
569 static void celt_decode_fine_energy(CeltContext *s, OpusRangeCoder *rc)
572 for (i = s->startband; i < s->endband; i++) {
574 if (!s->fine_bits[i])
577 for (j = 0; j < s->coded_channels; j++) {
578 CeltFrame *frame = &s->frame[j];
581 q2 = opus_getrawbits(rc, s->fine_bits[i]);
582 offset = (q2 + 0.5f) * (1 << (14 - s->fine_bits[i])) / 16384.0f - 0.5f;
583 frame->energy[i] += offset;
588 static void celt_decode_final_energy(CeltContext *s, OpusRangeCoder *rc,
593 for (priority = 0; priority < 2; priority++) {
594 for (i = s->startband; i < s->endband && bits_left >= s->coded_channels; i++) {
595 if (s->fine_priority[i] != priority || s->fine_bits[i] >= CELT_MAX_FINE_BITS)
598 for (j = 0; j < s->coded_channels; j++) {
601 q2 = opus_getrawbits(rc, 1);
602 offset = (q2 - 0.5f) * (1 << (14 - s->fine_bits[i] - 1)) / 16384.0f;
603 s->frame[j].energy[i] += offset;
610 static void celt_decode_tf_changes(CeltContext *s, OpusRangeCoder *rc,
613 int i, diff = 0, tf_select = 0, tf_changed = 0, tf_select_bit;
614 int consumed, bits = transient ? 2 : 4;
616 consumed = opus_rc_tell(rc);
617 tf_select_bit = (s->duration != 0 && consumed+bits+1 <= s->framebits);
619 for (i = s->startband; i < s->endband; i++) {
620 if (consumed+bits+tf_select_bit <= s->framebits) {
621 diff ^= opus_rc_p2model(rc, bits);
622 consumed = opus_rc_tell(rc);
625 s->tf_change[i] = diff;
626 bits = transient ? 4 : 5;
629 if (tf_select_bit && celt_tf_select[s->duration][transient][0][tf_changed] !=
630 celt_tf_select[s->duration][transient][1][tf_changed])
631 tf_select = opus_rc_p2model(rc, 1);
633 for (i = s->startband; i < s->endband; i++) {
634 s->tf_change[i] = celt_tf_select[s->duration][transient][tf_select][s->tf_change[i]];
638 static void celt_decode_allocation(CeltContext *s, OpusRangeCoder *rc)
640 // approx. maximum bit allocation for each band before boost/trim
641 int cap[CELT_MAX_BANDS];
642 int boost[CELT_MAX_BANDS];
643 int threshold[CELT_MAX_BANDS];
644 int bits1[CELT_MAX_BANDS];
645 int bits2[CELT_MAX_BANDS];
646 int trim_offset[CELT_MAX_BANDS];
648 int skip_startband = s->startband;
654 int intensitystereo_bit = 0;
655 int dualstereo_bit = 0;
657 int remaining, bandbits;
658 int low, high, total, done;
663 consumed = opus_rc_tell(rc);
665 /* obtain spread flag */
666 s->spread = CELT_SPREAD_NORMAL;
667 if (consumed + 4 <= s->framebits)
668 s->spread = opus_rc_getsymbol(rc, celt_model_spread);
670 /* generate static allocation caps */
671 for (i = 0; i < CELT_MAX_BANDS; i++) {
672 cap[i] = (celt_static_caps[s->duration][s->coded_channels - 1][i] + 64)
673 * celt_freq_range[i] << (s->coded_channels - 1) << s->duration >> 2;
676 /* obtain band boost */
677 totalbits = s->framebits << 3; // convert to 1/8 bits
678 consumed = opus_rc_tell_frac(rc);
679 for (i = s->startband; i < s->endband; i++) {
680 int quanta, band_dynalloc;
684 quanta = celt_freq_range[i] << (s->coded_channels - 1) << s->duration;
685 quanta = FFMIN(quanta << 3, FFMAX(6 << 3, quanta));
686 band_dynalloc = dynalloc;
687 while (consumed + (band_dynalloc<<3) < totalbits && boost[i] < cap[i]) {
688 int add = opus_rc_p2model(rc, band_dynalloc);
689 consumed = opus_rc_tell_frac(rc);
697 /* dynalloc is more likely to occur if it's already been used for earlier bands */
699 dynalloc = FFMAX(2, dynalloc - 1);
702 /* obtain allocation trim */
703 if (consumed + (6 << 3) <= totalbits)
704 alloctrim = opus_rc_getsymbol(rc, celt_model_alloc_trim);
706 /* anti-collapse bit reservation */
707 totalbits = (s->framebits << 3) - opus_rc_tell_frac(rc) - 1;
708 s->anticollapse_bit = 0;
709 if (s->blocks > 1 && s->duration >= 2 &&
710 totalbits >= ((s->duration + 2) << 3))
711 s->anticollapse_bit = 1 << 3;
712 totalbits -= s->anticollapse_bit;
714 /* band skip bit reservation */
715 if (totalbits >= 1 << 3)
717 totalbits -= skip_bit;
719 /* intensity/dual stereo bit reservation */
720 if (s->coded_channels == 2) {
721 intensitystereo_bit = celt_log2_frac[s->endband - s->startband];
722 if (intensitystereo_bit <= totalbits) {
723 totalbits -= intensitystereo_bit;
724 if (totalbits >= 1 << 3) {
725 dualstereo_bit = 1 << 3;
729 intensitystereo_bit = 0;
732 for (i = s->startband; i < s->endband; i++) {
733 int trim = alloctrim - 5 - s->duration;
734 int band = celt_freq_range[i] * (s->endband - i - 1);
735 int duration = s->duration + 3;
736 int scale = duration + s->coded_channels - 1;
738 /* PVQ minimum allocation threshold, below this value the band is
740 threshold[i] = FFMAX(3 * celt_freq_range[i] << duration >> 4,
741 s->coded_channels << 3);
743 trim_offset[i] = trim * (band << scale) >> 6;
745 if (celt_freq_range[i] << s->duration == 1)
746 trim_offset[i] -= s->coded_channels << 3;
751 high = CELT_VECTORS - 1;
752 while (low <= high) {
753 int center = (low + high) >> 1;
756 for (i = s->endband - 1; i >= s->startband; i--) {
757 bandbits = celt_freq_range[i] * celt_static_alloc[center][i]
758 << (s->coded_channels - 1) << s->duration >> 2;
761 bandbits = FFMAX(0, bandbits + trim_offset[i]);
762 bandbits += boost[i];
764 if (bandbits >= threshold[i] || done) {
766 total += FFMIN(bandbits, cap[i]);
767 } else if (bandbits >= s->coded_channels << 3)
768 total += s->coded_channels << 3;
771 if (total > totalbits)
778 for (i = s->startband; i < s->endband; i++) {
779 bits1[i] = celt_freq_range[i] * celt_static_alloc[low][i]
780 << (s->coded_channels - 1) << s->duration >> 2;
781 bits2[i] = high >= CELT_VECTORS ? cap[i] :
782 celt_freq_range[i] * celt_static_alloc[high][i]
783 << (s->coded_channels - 1) << s->duration >> 2;
786 bits1[i] = FFMAX(0, bits1[i] + trim_offset[i]);
788 bits2[i] = FFMAX(0, bits2[i] + trim_offset[i]);
790 bits1[i] += boost[i];
791 bits2[i] += boost[i];
795 bits2[i] = FFMAX(0, bits2[i] - bits1[i]);
800 high = 1 << CELT_ALLOC_STEPS;
801 for (i = 0; i < CELT_ALLOC_STEPS; i++) {
802 int center = (low + high) >> 1;
805 for (j = s->endband - 1; j >= s->startband; j--) {
806 bandbits = bits1[j] + (center * bits2[j] >> CELT_ALLOC_STEPS);
808 if (bandbits >= threshold[j] || done) {
810 total += FFMIN(bandbits, cap[j]);
811 } else if (bandbits >= s->coded_channels << 3)
812 total += s->coded_channels << 3;
814 if (total > totalbits)
821 for (i = s->endband - 1; i >= s->startband; i--) {
822 bandbits = bits1[i] + (low * bits2[i] >> CELT_ALLOC_STEPS);
824 if (bandbits >= threshold[i] || done)
827 bandbits = (bandbits >= s->coded_channels << 3) ?
828 s->coded_channels << 3 : 0;
830 bandbits = FFMIN(bandbits, cap[i]);
831 s->pulses[i] = bandbits;
836 for (s->codedbands = s->endband; ; s->codedbands--) {
838 j = s->codedbands - 1;
840 if (j == skip_startband) {
841 /* all remaining bands are not skipped */
842 totalbits += skip_bit;
846 /* determine the number of bits available for coding "do not skip" markers */
847 remaining = totalbits - total;
848 bandbits = remaining / (celt_freq_bands[j+1] - celt_freq_bands[s->startband]);
849 remaining -= bandbits * (celt_freq_bands[j+1] - celt_freq_bands[s->startband]);
850 allocation = s->pulses[j] + bandbits * celt_freq_range[j]
851 + FFMAX(0, remaining - (celt_freq_bands[j] - celt_freq_bands[s->startband]));
853 /* a "do not skip" marker is only coded if the allocation is
854 above the chosen threshold */
855 if (allocation >= FFMAX(threshold[j], (s->coded_channels + 1) <<3 )) {
856 if (opus_rc_p2model(rc, 1))
860 allocation -= 1 << 3;
863 /* the band is skipped, so reclaim its bits */
864 total -= s->pulses[j];
865 if (intensitystereo_bit) {
866 total -= intensitystereo_bit;
867 intensitystereo_bit = celt_log2_frac[j - s->startband];
868 total += intensitystereo_bit;
871 total += s->pulses[j] = (allocation >= s->coded_channels << 3) ?
872 s->coded_channels << 3 : 0;
875 /* obtain stereo flags */
876 s->intensitystereo = 0;
878 if (intensitystereo_bit)
879 s->intensitystereo = s->startband +
880 opus_rc_unimodel(rc, s->codedbands + 1 - s->startband);
881 if (s->intensitystereo <= s->startband)
882 totalbits += dualstereo_bit; /* no intensity stereo means no dual stereo */
883 else if (dualstereo_bit)
884 s->dualstereo = opus_rc_p2model(rc, 1);
886 /* supply the remaining bits in this frame to lower bands */
887 remaining = totalbits - total;
888 bandbits = remaining / (celt_freq_bands[s->codedbands] - celt_freq_bands[s->startband]);
889 remaining -= bandbits * (celt_freq_bands[s->codedbands] - celt_freq_bands[s->startband]);
890 for (i = s->startband; i < s->codedbands; i++) {
891 int bits = FFMIN(remaining, celt_freq_range[i]);
893 s->pulses[i] += bits + bandbits * celt_freq_range[i];
897 for (i = s->startband; i < s->codedbands; i++) {
898 int N = celt_freq_range[i] << s->duration;
899 int prev_extra = extrabits;
900 s->pulses[i] += extrabits;
903 int dof; // degrees of freedom
904 int temp; // dof * channels * log(dof)
905 int offset; // fine energy quantization offset, i.e.
906 // extra bits assigned over the standard
908 int fine_bits, max_bits;
910 extrabits = FFMAX(0, s->pulses[i] - cap[i]);
911 s->pulses[i] -= extrabits;
913 /* intensity stereo makes use of an extra degree of freedom */
914 dof = N * s->coded_channels
915 + (s->coded_channels == 2 && N > 2 && !s->dualstereo && i < s->intensitystereo);
916 temp = dof * (celt_log_freq_range[i] + (s->duration<<3));
917 offset = (temp >> 1) - dof * CELT_FINE_OFFSET;
918 if (N == 2) /* dof=2 is the only case that doesn't fit the model */
921 /* grant an additional bias for the first and second pulses */
922 if (s->pulses[i] + offset < 2 * (dof << 3))
924 else if (s->pulses[i] + offset < 3 * (dof << 3))
927 fine_bits = (s->pulses[i] + offset + (dof << 2)) / (dof << 3);
928 max_bits = FFMIN((s->pulses[i]>>3) >> (s->coded_channels - 1),
931 max_bits = FFMAX(max_bits, 0);
933 s->fine_bits[i] = av_clip(fine_bits, 0, max_bits);
935 /* if fine_bits was rounded down or capped,
936 give priority for the final fine energy pass */
937 s->fine_priority[i] = (s->fine_bits[i] * (dof<<3) >= s->pulses[i] + offset);
939 /* the remaining bits are assigned to PVQ */
940 s->pulses[i] -= s->fine_bits[i] << (s->coded_channels - 1) << 3;
942 /* all bits go to fine energy except for the sign bit */
943 extrabits = FFMAX(0, s->pulses[i] - (s->coded_channels << 3));
944 s->pulses[i] -= extrabits;
946 s->fine_priority[i] = 1;
949 /* hand back a limited number of extra fine energy bits to this band */
951 int fineextra = FFMIN(extrabits >> (s->coded_channels + 2),
952 CELT_MAX_FINE_BITS - s->fine_bits[i]);
953 s->fine_bits[i] += fineextra;
955 fineextra <<= s->coded_channels + 2;
956 s->fine_priority[i] = (fineextra >= extrabits - prev_extra);
957 extrabits -= fineextra;
960 s->remaining = extrabits;
962 /* skipped bands dedicate all of their bits for fine energy */
963 for (; i < s->endband; i++) {
964 s->fine_bits[i] = s->pulses[i] >> (s->coded_channels - 1) >> 3;
966 s->fine_priority[i] = s->fine_bits[i] < 1;
970 static inline int celt_bits2pulses(const uint8_t *cache, int bits)
972 // TODO: Find the size of cache and make it into an array in the parameters list
973 int i, low = 0, high;
978 for (i = 0; i < 6; i++) {
979 int center = (low + high + 1) >> 1;
980 if (cache[center] >= bits)
986 return (bits - (low == 0 ? -1 : cache[low]) <= cache[high] - bits) ? low : high;
989 static inline int celt_pulses2bits(const uint8_t *cache, int pulses)
991 // TODO: Find the size of cache and make it into an array in the parameters list
992 return (pulses == 0) ? 0 : cache[pulses] + 1;
995 static inline void celt_normalize_residual(const int * av_restrict iy, float * av_restrict X,
999 for (i = 0; i < N; i++)
1003 static void celt_exp_rotation1(float *X, unsigned int len, unsigned int stride,
1010 for (i = 0; i < len - stride; i++) {
1014 Xptr[stride] = c * x2 + s * x1;
1015 *Xptr++ = c * x1 - s * x2;
1018 Xptr = &X[len - 2 * stride - 1];
1019 for (i = len - 2 * stride - 1; i >= 0; i--) {
1023 Xptr[stride] = c * x2 + s * x1;
1024 *Xptr-- = c * x1 - s * x2;
1028 static inline void celt_exp_rotation(float *X, unsigned int len,
1029 unsigned int stride, unsigned int K,
1030 enum CeltSpread spread)
1032 unsigned int stride2 = 0;
1037 if (2*K >= len || spread == CELT_SPREAD_NONE)
1040 gain = (float)len / (len + (20 - 5*spread) * K);
1041 theta = M_PI * gain * gain / 4;
1046 if (len >= stride << 3) {
1048 /* This is just a simple (equivalent) way of computing sqrt(len/stride) with rounding.
1049 It's basically incrementing long as (stride2+0.5)^2 < len/stride. */
1050 while ((stride2 * stride2 + stride2) * stride + (stride >> 2) < len)
1054 /*NOTE: As a minor optimization, we could be passing around log2(B), not B, for both this and for
1055 extract_collapse_mask().*/
1057 for (i = 0; i < stride; i++) {
1059 celt_exp_rotation1(X + i * len, len, stride2, s, c);
1060 celt_exp_rotation1(X + i * len, len, 1, c, s);
1064 static inline unsigned int celt_extract_collapse_mask(const int *iy,
1068 unsigned int collapse_mask;
1075 /*NOTE: As a minor optimization, we could be passing around log2(B), not B, for both this and for
1079 for (i = 0; i < B; i++)
1080 for (j = 0; j < N0; j++)
1081 collapse_mask |= (iy[i*N0+j]!=0)<<i;
1082 return collapse_mask;
1085 static inline void celt_renormalize_vector(float *X, int N, float gain)
1089 for (i = 0; i < N; i++)
1091 g = gain / sqrtf(g);
1093 for (i = 0; i < N; i++)
1097 static inline void celt_stereo_merge(float *X, float *Y, float mid, int N)
1100 float xp = 0, side = 0;
1105 /* Compute the norm of X+Y and X-Y as |X|^2 + |Y|^2 +/- sum(xy) */
1106 for (i = 0; i < N; i++) {
1108 side += Y[i] * Y[i];
1111 /* Compensating for the mid normalization */
1114 E[0] = mid2 * mid2 + side - 2 * xp;
1115 E[1] = mid2 * mid2 + side + 2 * xp;
1116 if (E[0] < 6e-4f || E[1] < 6e-4f) {
1117 for (i = 0; i < N; i++)
1123 gain[0] = 1.0f / sqrtf(t);
1125 gain[1] = 1.0f / sqrtf(t);
1127 for (i = 0; i < N; i++) {
1129 /* Apply mid scaling (side is already scaled) */
1130 value[0] = mid * X[i];
1132 X[i] = gain[0] * (value[0] - value[1]);
1133 Y[i] = gain[1] * (value[0] + value[1]);
1137 static void celt_interleave_hadamard(float *tmp, float *X, int N0,
1138 int stride, int hadamard)
1144 const uint8_t *ordery = celt_hadamard_ordery + stride - 2;
1145 for (i = 0; i < stride; i++)
1146 for (j = 0; j < N0; j++)
1147 tmp[j*stride+i] = X[ordery[i]*N0+j];
1149 for (i = 0; i < stride; i++)
1150 for (j = 0; j < N0; j++)
1151 tmp[j*stride+i] = X[i*N0+j];
1154 for (i = 0; i < N; i++)
1158 static void celt_deinterleave_hadamard(float *tmp, float *X, int N0,
1159 int stride, int hadamard)
1165 const uint8_t *ordery = celt_hadamard_ordery + stride - 2;
1166 for (i = 0; i < stride; i++)
1167 for (j = 0; j < N0; j++)
1168 tmp[ordery[i]*N0+j] = X[j*stride+i];
1170 for (i = 0; i < stride; i++)
1171 for (j = 0; j < N0; j++)
1172 tmp[i*N0+j] = X[j*stride+i];
1175 for (i = 0; i < N; i++)
1179 static void celt_haar1(float *X, int N0, int stride)
1183 for (i = 0; i < stride; i++) {
1184 for (j = 0; j < N0; j++) {
1185 float x0 = X[stride * (2 * j + 0) + i];
1186 float x1 = X[stride * (2 * j + 1) + i];
1187 X[stride * (2 * j + 0) + i] = (x0 + x1) * M_SQRT1_2;
1188 X[stride * (2 * j + 1) + i] = (x0 - x1) * M_SQRT1_2;
1193 static inline int celt_compute_qn(int N, int b, int offset, int pulse_cap,
1198 if (dualstereo && N == 2)
1201 /* The upper limit ensures that in a stereo split with itheta==16384, we'll
1202 * always have enough bits left over to code at least one pulse in the
1203 * side; otherwise it would collapse, since it doesn't get folded. */
1204 qb = FFMIN3(b - pulse_cap - (4 << 3), (b + N2 * offset) / N2, 8 << 3);
1205 qn = (qb < (1 << 3 >> 1)) ? 1 : ((celt_qn_exp2[qb & 0x7] >> (14 - (qb >> 3))) + 1) >> 1 << 1;
1209 // this code was adapted from libopus
1210 static inline uint64_t celt_cwrsi(unsigned int N, unsigned int K, unsigned int i, int *y)
1220 /*Lots of pulses case:*/
1222 const uint32_t *row = celt_pvq_u_row[N];
1224 /* Are the pulses in this dimension negative? */
1229 /*Count how many pulses were placed in this dimension.*/
1235 p = celt_pvq_u_row[--K][N];
1238 for (p = row[K]; p > i; p = row[K])
1242 val = (k0 - K + s) ^ s;
1245 } else { /*Lots of dimensions case:*/
1246 /*Are there any pulses in this dimension at all?*/
1247 p = celt_pvq_u_row[K ][N];
1248 q = celt_pvq_u_row[K + 1][N];
1250 if (p <= i && i < q) {
1254 /*Are the pulses in this dimension negative?*/
1258 /*Count how many pulses were placed in this dimension.*/
1260 do p = celt_pvq_u_row[--K][N];
1264 val = (k0 - K + s) ^ s;
1282 val = (k0 - K + s) ^ s;
1295 static inline float celt_decode_pulses(OpusRangeCoder *rc, int *y, unsigned int N, unsigned int K)
1298 #define CELT_PVQ_U(n, k) (celt_pvq_u_row[FFMIN(n, k)][FFMAX(n, k)])
1299 #define CELT_PVQ_V(n, k) (CELT_PVQ_U(n, k) + CELT_PVQ_U(n, (k) + 1))
1300 idx = opus_rc_unimodel(rc, CELT_PVQ_V(N, K));
1301 return celt_cwrsi(N, K, idx, y);
1304 /** Decode pulse vector and combine the result with the pitch vector to produce
1305 the final normalised signal in the current band. */
1306 static inline unsigned int celt_alg_unquant(OpusRangeCoder *rc, float *X,
1307 unsigned int N, unsigned int K,
1308 enum CeltSpread spread,
1309 unsigned int blocks, float gain)
1313 gain /= sqrtf(celt_decode_pulses(rc, y, N, K));
1314 celt_normalize_residual(y, X, N, gain);
1315 celt_exp_rotation(X, N, blocks, K, spread);
1316 return celt_extract_collapse_mask(y, N, blocks);
1319 static unsigned int celt_decode_band(CeltContext *s, OpusRangeCoder *rc,
1320 const int band, float *X, float *Y,
1321 int N, int b, unsigned int blocks,
1322 float *lowband, int duration,
1323 float *lowband_out, int level,
1324 float gain, float *lowband_scratch,
1327 const uint8_t *cache;
1328 int dualstereo, split;
1329 int imid = 0, iside = 0;
1330 unsigned int N0 = N;
1334 int time_divide = 0;
1337 float mid = 0, side = 0;
1338 int longblocks = (B0 == 1);
1339 unsigned int cm = 0;
1341 N_B0 = N_B = N / blocks;
1342 split = dualstereo = (Y != NULL);
1345 /* special case for one sample */
1348 for (i = 0; i <= dualstereo; i++) {
1350 if (s->remaining2 >= 1<<3) {
1351 sign = opus_getrawbits(rc, 1);
1352 s->remaining2 -= 1 << 3;
1355 x[0] = sign ? -1.0f : 1.0f;
1359 lowband_out[0] = X[0];
1363 if (!dualstereo && level == 0) {
1364 int tf_change = s->tf_change[band];
1367 recombine = tf_change;
1368 /* Band recombining to increase frequency resolution */
1371 (recombine || ((N_B & 1) == 0 && tf_change < 0) || B0 > 1)) {
1373 for (j = 0; j < N; j++)
1374 lowband_scratch[j] = lowband[j];
1375 lowband = lowband_scratch;
1378 for (k = 0; k < recombine; k++) {
1380 celt_haar1(lowband, N >> k, 1 << k);
1381 fill = celt_bit_interleave[fill & 0xF] | celt_bit_interleave[fill >> 4] << 2;
1383 blocks >>= recombine;
1386 /* Increasing the time resolution */
1387 while ((N_B & 1) == 0 && tf_change < 0) {
1389 celt_haar1(lowband, N_B, blocks);
1390 fill |= fill << blocks;
1399 /* Reorganize the samples in time order instead of frequency order */
1400 if (B0 > 1 && lowband)
1401 celt_deinterleave_hadamard(s->scratch, lowband, N_B >> recombine,
1402 B0 << recombine, longblocks);
1405 /* If we need 1.5 more bit than we can produce, split the band in two. */
1406 cache = celt_cache_bits +
1407 celt_cache_index[(duration + 1) * CELT_MAX_BANDS + band];
1408 if (!dualstereo && duration >= 0 && b > cache[cache[0]] + 12 && N > 2) {
1414 fill = (fill & 1) | (fill << 1);
1415 blocks = (blocks + 1) >> 1;
1421 int mbits, sbits, delta;
1428 /* Decide on the resolution to give to the split parameter theta */
1429 pulse_cap = celt_log_freq_range[band] + duration * 8;
1430 offset = (pulse_cap >> 1) - (dualstereo && N == 2 ? CELT_QTHETA_OFFSET_TWOPHASE :
1431 CELT_QTHETA_OFFSET);
1432 qn = (dualstereo && band >= s->intensitystereo) ? 1 :
1433 celt_compute_qn(N, b, offset, pulse_cap, dualstereo);
1434 tell = opus_rc_tell_frac(rc);
1436 /* Entropy coding of the angle. We use a uniform pdf for the
1437 time split, a step for stereo, and a triangular one for the rest. */
1438 if (dualstereo && N > 2)
1439 itheta = opus_rc_stepmodel(rc, qn/2);
1440 else if (dualstereo || B0 > 1)
1441 itheta = opus_rc_unimodel(rc, qn+1);
1443 itheta = opus_rc_trimodel(rc, qn);
1444 itheta = itheta * 16384 / qn;
1445 /* NOTE: Renormalising X and Y *may* help fixed-point a bit at very high rate.
1446 Let's do that at higher complexity */
1447 } else if (dualstereo) {
1448 inv = (b > 2 << 3 && s->remaining2 > 2 << 3) ? opus_rc_p2model(rc, 2) : 0;
1451 qalloc = opus_rc_tell_frac(rc) - tell;
1458 fill = av_mod_uintp2(fill, blocks);
1460 } else if (itheta == 16384) {
1463 fill &= ((1 << blocks) - 1) << blocks;
1466 imid = celt_cos(itheta);
1467 iside = celt_cos(16384-itheta);
1468 /* This is the mid vs side allocation that minimizes squared error
1470 delta = ROUND_MUL16((N - 1) << 7, celt_log2tan(iside, imid));
1473 mid = imid / 32768.0f;
1474 side = iside / 32768.0f;
1476 /* This is a special case for N=2 that only works for stereo and takes
1477 advantage of the fact that mid and side are orthogonal to encode
1478 the side with just one bit. */
1479 if (N == 2 && dualstereo) {
1485 /* Only need one bit for the side */
1486 sbits = (itheta != 0 && itheta != 16384) ? 1 << 3 : 0;
1488 c = (itheta > 8192);
1489 s->remaining2 -= qalloc+sbits;
1494 sign = opus_getrawbits(rc, 1);
1495 sign = 1 - 2 * sign;
1496 /* We use orig_fill here because we want to fold the side, but if
1497 itheta==16384, we'll have cleared the low bits of fill. */
1498 cm = celt_decode_band(s, rc, band, x2, NULL, N, mbits, blocks,
1499 lowband, duration, lowband_out, level, gain,
1500 lowband_scratch, orig_fill);
1501 /* We don't split N=2 bands, so cm is either 1 or 0 (for a fold-collapse),
1502 and there's no need to worry about mixing with the other channel. */
1503 y2[0] = -sign * x2[1];
1504 y2[1] = sign * x2[0];
1516 /* "Normal" split code */
1517 float *next_lowband2 = NULL;
1518 float *next_lowband_out1 = NULL;
1522 /* Give more bits to low-energy MDCTs than they would
1523 * otherwise deserve */
1524 if (B0 > 1 && !dualstereo && (itheta & 0x3fff)) {
1526 /* Rough approximation for pre-echo masking */
1527 delta -= delta >> (4 - duration);
1529 /* Corresponds to a forward-masking slope of
1530 * 1.5 dB per 10 ms */
1531 delta = FFMIN(0, delta + (N << 3 >> (5 - duration)));
1533 mbits = av_clip((b - delta) / 2, 0, b);
1535 s->remaining2 -= qalloc;
1537 if (lowband && !dualstereo)
1538 next_lowband2 = lowband + N; /* >32-bit split case */
1540 /* Only stereo needs to pass on lowband_out.
1541 * Otherwise, it's handled at the end */
1543 next_lowband_out1 = lowband_out;
1545 next_level = level + 1;
1547 rebalance = s->remaining2;
1548 if (mbits >= sbits) {
1549 /* In stereo mode, we do not apply a scaling to the mid
1550 * because we need the normalized mid for folding later */
1551 cm = celt_decode_band(s, rc, band, X, NULL, N, mbits, blocks,
1552 lowband, duration, next_lowband_out1,
1553 next_level, dualstereo ? 1.0f : (gain * mid),
1554 lowband_scratch, fill);
1556 rebalance = mbits - (rebalance - s->remaining2);
1557 if (rebalance > 3 << 3 && itheta != 0)
1558 sbits += rebalance - (3 << 3);
1560 /* For a stereo split, the high bits of fill are always zero,
1561 * so no folding will be done to the side. */
1562 cm |= celt_decode_band(s, rc, band, Y, NULL, N, sbits, blocks,
1563 next_lowband2, duration, NULL,
1564 next_level, gain * side, NULL,
1565 fill >> blocks) << ((B0 >> 1) & (dualstereo - 1));
1567 /* For a stereo split, the high bits of fill are always zero,
1568 * so no folding will be done to the side. */
1569 cm = celt_decode_band(s, rc, band, Y, NULL, N, sbits, blocks,
1570 next_lowband2, duration, NULL,
1571 next_level, gain * side, NULL,
1572 fill >> blocks) << ((B0 >> 1) & (dualstereo - 1));
1574 rebalance = sbits - (rebalance - s->remaining2);
1575 if (rebalance > 3 << 3 && itheta != 16384)
1576 mbits += rebalance - (3 << 3);
1578 /* In stereo mode, we do not apply a scaling to the mid because
1579 * we need the normalized mid for folding later */
1580 cm |= celt_decode_band(s, rc, band, X, NULL, N, mbits, blocks,
1581 lowband, duration, next_lowband_out1,
1582 next_level, dualstereo ? 1.0f : (gain * mid),
1583 lowband_scratch, fill);
1587 /* This is the basic no-split case */
1588 unsigned int q = celt_bits2pulses(cache, b);
1589 unsigned int curr_bits = celt_pulses2bits(cache, q);
1590 s->remaining2 -= curr_bits;
1592 /* Ensures we can never bust the budget */
1593 while (s->remaining2 < 0 && q > 0) {
1594 s->remaining2 += curr_bits;
1595 curr_bits = celt_pulses2bits(cache, --q);
1596 s->remaining2 -= curr_bits;
1600 /* Finally do the actual quantization */
1601 cm = celt_alg_unquant(rc, X, N, (q < 8) ? q : (8 + (q & 7)) << ((q >> 3) - 1),
1602 s->spread, blocks, gain);
1604 /* If there's no pulse, fill the band anyway */
1606 unsigned int cm_mask = (1 << blocks) - 1;
1609 for (j = 0; j < N; j++)
1614 for (j = 0; j < N; j++)
1615 X[j] = (((int32_t)celt_rng(s)) >> 20);
1618 /* Folded spectrum */
1619 for (j = 0; j < N; j++) {
1620 /* About 48 dB below the "normal" folding level */
1621 X[j] = lowband[j] + (((celt_rng(s)) & 0x8000) ? 1.0f / 256 : -1.0f / 256);
1625 celt_renormalize_vector(X, N, gain);
1630 /* This code is used by the decoder and by the resynthesis-enabled encoder */
1634 celt_stereo_merge(X, Y, mid, N);
1636 for (j = 0; j < N; j++)
1639 } else if (level == 0) {
1642 /* Undo the sample reorganization going from time order to frequency order */
1644 celt_interleave_hadamard(s->scratch, X, N_B>>recombine,
1645 B0<<recombine, longblocks);
1647 /* Undo time-freq changes that we did earlier */
1650 for (k = 0; k < time_divide; k++) {
1654 celt_haar1(X, N_B, blocks);
1657 for (k = 0; k < recombine; k++) {
1658 cm = celt_bit_deinterleave[cm];
1659 celt_haar1(X, N0>>k, 1<<k);
1661 blocks <<= recombine;
1663 /* Scale output for later folding */
1666 float n = sqrtf(N0);
1667 for (j = 0; j < N0; j++)
1668 lowband_out[j] = n * X[j];
1670 cm = av_mod_uintp2(cm, blocks);
1675 static void celt_denormalize(CeltContext *s, CeltFrame *frame, float *data)
1679 for (i = s->startband; i < s->endband; i++) {
1680 float *dst = data + (celt_freq_bands[i] << s->duration);
1681 float norm = exp2(frame->energy[i] + celt_mean_energy[i]);
1683 for (j = 0; j < celt_freq_range[i] << s->duration; j++)
1688 static void celt_postfilter_apply_transition(CeltFrame *frame, float *data)
1690 const int T0 = frame->pf_period_old;
1691 const int T1 = frame->pf_period;
1693 float g00, g01, g02;
1694 float g10, g11, g12;
1696 float x0, x1, x2, x3, x4;
1700 if (frame->pf_gains[0] == 0.0 &&
1701 frame->pf_gains_old[0] == 0.0)
1704 g00 = frame->pf_gains_old[0];
1705 g01 = frame->pf_gains_old[1];
1706 g02 = frame->pf_gains_old[2];
1707 g10 = frame->pf_gains[0];
1708 g11 = frame->pf_gains[1];
1709 g12 = frame->pf_gains[2];
1716 for (i = 0; i < CELT_OVERLAP; i++) {
1717 float w = ff_celt_window2[i];
1718 x0 = data[i - T1 + 2];
1720 data[i] += (1.0 - w) * g00 * data[i - T0] +
1721 (1.0 - w) * g01 * (data[i - T0 - 1] + data[i - T0 + 1]) +
1722 (1.0 - w) * g02 * (data[i - T0 - 2] + data[i - T0 + 2]) +
1724 w * g11 * (x1 + x3) +
1725 w * g12 * (x0 + x4);
1733 static void celt_postfilter_apply(CeltFrame *frame,
1734 float *data, int len)
1736 const int T = frame->pf_period;
1738 float x0, x1, x2, x3, x4;
1741 if (frame->pf_gains[0] == 0.0 || len <= 0)
1744 g0 = frame->pf_gains[0];
1745 g1 = frame->pf_gains[1];
1746 g2 = frame->pf_gains[2];
1753 for (i = 0; i < len; i++) {
1754 x0 = data[i - T + 2];
1755 data[i] += g0 * x2 +
1765 static void celt_postfilter(CeltContext *s, CeltFrame *frame)
1767 int len = s->blocksize * s->blocks;
1769 celt_postfilter_apply_transition(frame, frame->buf + 1024);
1771 frame->pf_period_old = frame->pf_period;
1772 memcpy(frame->pf_gains_old, frame->pf_gains, sizeof(frame->pf_gains));
1774 frame->pf_period = frame->pf_period_new;
1775 memcpy(frame->pf_gains, frame->pf_gains_new, sizeof(frame->pf_gains));
1777 if (len > CELT_OVERLAP) {
1778 celt_postfilter_apply_transition(frame, frame->buf + 1024 + CELT_OVERLAP);
1779 celt_postfilter_apply(frame, frame->buf + 1024 + 2 * CELT_OVERLAP,
1780 len - 2 * CELT_OVERLAP);
1782 frame->pf_period_old = frame->pf_period;
1783 memcpy(frame->pf_gains_old, frame->pf_gains, sizeof(frame->pf_gains));
1786 memmove(frame->buf, frame->buf + len, (1024 + CELT_OVERLAP / 2) * sizeof(float));
1789 static int parse_postfilter(CeltContext *s, OpusRangeCoder *rc, int consumed)
1791 static const float postfilter_taps[3][3] = {
1792 { 0.3066406250f, 0.2170410156f, 0.1296386719f },
1793 { 0.4638671875f, 0.2680664062f, 0.0 },
1794 { 0.7998046875f, 0.1000976562f, 0.0 }
1798 memset(s->frame[0].pf_gains_new, 0, sizeof(s->frame[0].pf_gains_new));
1799 memset(s->frame[1].pf_gains_new, 0, sizeof(s->frame[1].pf_gains_new));
1801 if (s->startband == 0 && consumed + 16 <= s->framebits) {
1802 int has_postfilter = opus_rc_p2model(rc, 1);
1803 if (has_postfilter) {
1805 int tapset, octave, period;
1807 octave = opus_rc_unimodel(rc, 6);
1808 period = (16 << octave) + opus_getrawbits(rc, 4 + octave) - 1;
1809 gain = 0.09375f * (opus_getrawbits(rc, 3) + 1);
1810 tapset = (opus_rc_tell(rc) + 2 <= s->framebits) ?
1811 opus_rc_getsymbol(rc, celt_model_tapset) : 0;
1813 for (i = 0; i < 2; i++) {
1814 CeltFrame *frame = &s->frame[i];
1816 frame->pf_period_new = FFMAX(period, CELT_POSTFILTER_MINPERIOD);
1817 frame->pf_gains_new[0] = gain * postfilter_taps[tapset][0];
1818 frame->pf_gains_new[1] = gain * postfilter_taps[tapset][1];
1819 frame->pf_gains_new[2] = gain * postfilter_taps[tapset][2];
1823 consumed = opus_rc_tell(rc);
1829 static void process_anticollapse(CeltContext *s, CeltFrame *frame, float *X)
1833 for (i = s->startband; i < s->endband; i++) {
1834 int renormalize = 0;
1838 float thresh, sqrt_1;
1841 /* depth in 1/8 bits */
1842 depth = (1 + s->pulses[i]) / (celt_freq_range[i] << s->duration);
1843 thresh = exp2f(-1.0 - 0.125f * depth);
1844 sqrt_1 = 1.0f / sqrtf(celt_freq_range[i] << s->duration);
1846 xptr = X + (celt_freq_bands[i] << s->duration);
1848 prev[0] = frame->prev_energy[0][i];
1849 prev[1] = frame->prev_energy[1][i];
1850 if (s->coded_channels == 1) {
1851 CeltFrame *frame1 = &s->frame[1];
1853 prev[0] = FFMAX(prev[0], frame1->prev_energy[0][i]);
1854 prev[1] = FFMAX(prev[1], frame1->prev_energy[1][i]);
1856 Ediff = frame->energy[i] - FFMIN(prev[0], prev[1]);
1857 Ediff = FFMAX(0, Ediff);
1859 /* r needs to be multiplied by 2 or 2*sqrt(2) depending on LM because
1860 short blocks don't have the same energy as long */
1861 r = exp2(1 - Ediff);
1862 if (s->duration == 3)
1864 r = FFMIN(thresh, r) * sqrt_1;
1865 for (k = 0; k < 1 << s->duration; k++) {
1866 /* Detect collapse */
1867 if (!(frame->collapse_masks[i] & 1 << k)) {
1868 /* Fill with noise */
1869 for (j = 0; j < celt_freq_range[i]; j++)
1870 xptr[(j << s->duration) + k] = (celt_rng(s) & 0x8000) ? r : -r;
1875 /* We just added some energy, so we need to renormalize */
1877 celt_renormalize_vector(xptr, celt_freq_range[i] << s->duration, 1.0f);
1881 static void celt_decode_bands(CeltContext *s, OpusRangeCoder *rc)
1883 float lowband_scratch[8 * 22];
1884 float norm[2 * 8 * 100];
1886 int totalbits = (s->framebits << 3) - s->anticollapse_bit;
1888 int update_lowband = 1;
1889 int lowband_offset = 0;
1893 memset(s->coeffs, 0, sizeof(s->coeffs));
1895 for (i = s->startband; i < s->endband; i++) {
1896 int band_offset = celt_freq_bands[i] << s->duration;
1897 int band_size = celt_freq_range[i] << s->duration;
1898 float *X = s->coeffs[0] + band_offset;
1899 float *Y = (s->coded_channels == 2) ? s->coeffs[1] + band_offset : NULL;
1901 int consumed = opus_rc_tell_frac(rc);
1902 float *norm2 = norm + 8 * 100;
1903 int effective_lowband = -1;
1907 /* Compute how many bits we want to allocate to this band */
1908 if (i != s->startband)
1909 s->remaining -= consumed;
1910 s->remaining2 = totalbits - consumed - 1;
1911 if (i <= s->codedbands - 1) {
1912 int curr_balance = s->remaining / FFMIN(3, s->codedbands-i);
1913 b = av_clip_uintp2(FFMIN(s->remaining2 + 1, s->pulses[i] + curr_balance), 14);
1917 if (celt_freq_bands[i] - celt_freq_range[i] >= celt_freq_bands[s->startband] &&
1918 (update_lowband || lowband_offset == 0))
1921 /* Get a conservative estimate of the collapse_mask's for the bands we're
1922 going to be folding from. */
1923 if (lowband_offset != 0 && (s->spread != CELT_SPREAD_AGGRESSIVE ||
1924 s->blocks > 1 || s->tf_change[i] < 0)) {
1925 int foldstart, foldend;
1927 /* This ensures we never repeat spectral content within one band */
1928 effective_lowband = FFMAX(celt_freq_bands[s->startband],
1929 celt_freq_bands[lowband_offset] - celt_freq_range[i]);
1930 foldstart = lowband_offset;
1931 while (celt_freq_bands[--foldstart] > effective_lowband);
1932 foldend = lowband_offset - 1;
1933 while (celt_freq_bands[++foldend] < effective_lowband + celt_freq_range[i]);
1936 for (j = foldstart; j < foldend; j++) {
1937 cm[0] |= s->frame[0].collapse_masks[j];
1938 cm[1] |= s->frame[s->coded_channels - 1].collapse_masks[j];
1941 /* Otherwise, we'll be using the LCG to fold, so all blocks will (almost
1942 always) be non-zero.*/
1943 cm[0] = cm[1] = (1 << s->blocks) - 1;
1945 if (s->dualstereo && i == s->intensitystereo) {
1946 /* Switch off dual stereo to do intensity */
1948 for (j = celt_freq_bands[s->startband] << s->duration; j < band_offset; j++)
1949 norm[j] = (norm[j] + norm2[j]) / 2;
1952 if (s->dualstereo) {
1953 cm[0] = celt_decode_band(s, rc, i, X, NULL, band_size, b / 2, s->blocks,
1954 effective_lowband != -1 ? norm + (effective_lowband << s->duration) : NULL, s->duration,
1955 norm + band_offset, 0, 1.0f, lowband_scratch, cm[0]);
1957 cm[1] = celt_decode_band(s, rc, i, Y, NULL, band_size, b/2, s->blocks,
1958 effective_lowband != -1 ? norm2 + (effective_lowband << s->duration) : NULL, s->duration,
1959 norm2 + band_offset, 0, 1.0f, lowband_scratch, cm[1]);
1961 cm[0] = celt_decode_band(s, rc, i, X, Y, band_size, b, s->blocks,
1962 effective_lowband != -1 ? norm + (effective_lowband << s->duration) : NULL, s->duration,
1963 norm + band_offset, 0, 1.0f, lowband_scratch, cm[0]|cm[1]);
1968 s->frame[0].collapse_masks[i] = (uint8_t)cm[0];
1969 s->frame[s->coded_channels - 1].collapse_masks[i] = (uint8_t)cm[1];
1970 s->remaining += s->pulses[i] + consumed;
1972 /* Update the folding position only as long as we have 1 bit/sample depth */
1973 update_lowband = (b > band_size << 3);
1977 int ff_celt_decode_frame(CeltContext *s, OpusRangeCoder *rc,
1978 float **output, int coded_channels, int frame_size,
1979 int startband, int endband)
1983 int consumed; // bits of entropy consumed thus far for this frame
1986 int anticollapse = 0;
1987 IMDCT15Context *imdct;
1988 float imdct_scale = 1.0;
1990 if (coded_channels != 1 && coded_channels != 2) {
1991 av_log(s->avctx, AV_LOG_ERROR, "Invalid number of coded channels: %d\n",
1993 return AVERROR_INVALIDDATA;
1995 if (startband < 0 || startband > endband || endband > CELT_MAX_BANDS) {
1996 av_log(s->avctx, AV_LOG_ERROR, "Invalid start/end band: %d %d\n",
1997 startband, endband);
1998 return AVERROR_INVALIDDATA;
2002 s->coded_channels = coded_channels;
2003 s->startband = startband;
2004 s->endband = endband;
2005 s->framebits = rc->rb.bytes * 8;
2007 s->duration = av_log2(frame_size / CELT_SHORT_BLOCKSIZE);
2008 if (s->duration > CELT_MAX_LOG_BLOCKS ||
2009 frame_size != CELT_SHORT_BLOCKSIZE * (1 << s->duration)) {
2010 av_log(s->avctx, AV_LOG_ERROR, "Invalid CELT frame size: %d\n",
2012 return AVERROR_INVALIDDATA;
2015 if (!s->output_channels)
2016 s->output_channels = coded_channels;
2018 memset(s->frame[0].collapse_masks, 0, sizeof(s->frame[0].collapse_masks));
2019 memset(s->frame[1].collapse_masks, 0, sizeof(s->frame[1].collapse_masks));
2021 consumed = opus_rc_tell(rc);
2023 /* obtain silence flag */
2024 if (consumed >= s->framebits)
2026 else if (consumed == 1)
2027 silence = opus_rc_p2model(rc, 15);
2031 consumed = s->framebits;
2032 rc->total_read_bits += s->framebits - opus_rc_tell(rc);
2035 /* obtain post-filter options */
2036 consumed = parse_postfilter(s, rc, consumed);
2038 /* obtain transient flag */
2039 if (s->duration != 0 && consumed+3 <= s->framebits)
2040 transient = opus_rc_p2model(rc, 3);
2042 s->blocks = transient ? 1 << s->duration : 1;
2043 s->blocksize = frame_size / s->blocks;
2045 imdct = s->imdct[transient ? 0 : s->duration];
2047 if (coded_channels == 1) {
2048 for (i = 0; i < CELT_MAX_BANDS; i++)
2049 s->frame[0].energy[i] = FFMAX(s->frame[0].energy[i], s->frame[1].energy[i]);
2052 celt_decode_coarse_energy(s, rc);
2053 celt_decode_tf_changes (s, rc, transient);
2054 celt_decode_allocation (s, rc);
2055 celt_decode_fine_energy (s, rc);
2056 celt_decode_bands (s, rc);
2058 if (s->anticollapse_bit)
2059 anticollapse = opus_getrawbits(rc, 1);
2061 celt_decode_final_energy(s, rc, s->framebits - opus_rc_tell(rc));
2063 /* apply anti-collapse processing and denormalization to
2064 * each coded channel */
2065 for (i = 0; i < s->coded_channels; i++) {
2066 CeltFrame *frame = &s->frame[i];
2069 process_anticollapse(s, frame, s->coeffs[i]);
2071 celt_denormalize(s, frame, s->coeffs[i]);
2074 /* stereo -> mono downmix */
2075 if (s->output_channels < s->coded_channels) {
2076 s->dsp->vector_fmac_scalar(s->coeffs[0], s->coeffs[1], 1.0, FFALIGN(frame_size, 16));
2078 } else if (s->output_channels > s->coded_channels)
2079 memcpy(s->coeffs[1], s->coeffs[0], frame_size * sizeof(float));
2082 for (i = 0; i < 2; i++) {
2083 CeltFrame *frame = &s->frame[i];
2085 for (j = 0; j < FF_ARRAY_ELEMS(frame->energy); j++)
2086 frame->energy[j] = CELT_ENERGY_SILENCE;
2088 memset(s->coeffs, 0, sizeof(s->coeffs));
2091 /* transform and output for each output channel */
2092 for (i = 0; i < s->output_channels; i++) {
2093 CeltFrame *frame = &s->frame[i];
2094 float m = frame->deemph_coeff;
2096 /* iMDCT and overlap-add */
2097 for (j = 0; j < s->blocks; j++) {
2098 float *dst = frame->buf + 1024 + j * s->blocksize;
2100 imdct->imdct_half(imdct, dst + CELT_OVERLAP / 2, s->coeffs[i] + j,
2101 s->blocks, imdct_scale);
2102 s->dsp->vector_fmul_window(dst, dst, dst + CELT_OVERLAP / 2,
2103 celt_window, CELT_OVERLAP / 2);
2107 celt_postfilter(s, frame);
2109 /* deemphasis and output scaling */
2110 for (j = 0; j < frame_size; j++) {
2111 float tmp = frame->buf[1024 - frame_size + j] + m;
2112 m = tmp * CELT_DEEMPH_COEFF;
2113 output[i][j] = tmp / 32768.;
2115 frame->deemph_coeff = m;
2118 if (coded_channels == 1)
2119 memcpy(s->frame[1].energy, s->frame[0].energy, sizeof(s->frame[0].energy));
2121 for (i = 0; i < 2; i++ ) {
2122 CeltFrame *frame = &s->frame[i];
2125 memcpy(frame->prev_energy[1], frame->prev_energy[0], sizeof(frame->prev_energy[0]));
2126 memcpy(frame->prev_energy[0], frame->energy, sizeof(frame->prev_energy[0]));
2128 for (j = 0; j < CELT_MAX_BANDS; j++)
2129 frame->prev_energy[0][j] = FFMIN(frame->prev_energy[0][j], frame->energy[j]);
2132 for (j = 0; j < s->startband; j++) {
2133 frame->prev_energy[0][j] = CELT_ENERGY_SILENCE;
2134 frame->energy[j] = 0.0;
2136 for (j = s->endband; j < CELT_MAX_BANDS; j++) {
2137 frame->prev_energy[0][j] = CELT_ENERGY_SILENCE;
2138 frame->energy[j] = 0.0;
2142 s->seed = rc->range;
2147 void ff_celt_flush(CeltContext *s)
2154 for (i = 0; i < 2; i++) {
2155 CeltFrame *frame = &s->frame[i];
2157 for (j = 0; j < CELT_MAX_BANDS; j++)
2158 frame->prev_energy[0][j] = frame->prev_energy[1][j] = CELT_ENERGY_SILENCE;
2160 memset(frame->energy, 0, sizeof(frame->energy));
2161 memset(frame->buf, 0, sizeof(frame->buf));
2163 memset(frame->pf_gains, 0, sizeof(frame->pf_gains));
2164 memset(frame->pf_gains_old, 0, sizeof(frame->pf_gains_old));
2165 memset(frame->pf_gains_new, 0, sizeof(frame->pf_gains_new));
2167 frame->deemph_coeff = 0.0;
2174 void ff_celt_free(CeltContext **ps)
2176 CeltContext *s = *ps;
2182 for (i = 0; i < FF_ARRAY_ELEMS(s->imdct); i++)
2183 ff_imdct15_uninit(&s->imdct[i]);
2189 int ff_celt_init(AVCodecContext *avctx, CeltContext **ps, int output_channels)
2194 if (output_channels != 1 && output_channels != 2) {
2195 av_log(avctx, AV_LOG_ERROR, "Invalid number of output channels: %d\n",
2197 return AVERROR(EINVAL);
2200 s = av_mallocz(sizeof(*s));
2202 return AVERROR(ENOMEM);
2205 s->output_channels = output_channels;
2207 for (i = 0; i < FF_ARRAY_ELEMS(s->imdct); i++) {
2208 ret = ff_imdct15_init(&s->imdct[i], i + 3);
2213 s->dsp = avpriv_float_dsp_alloc(avctx->flags & AV_CODEC_FLAG_BITEXACT);
2215 ret = AVERROR(ENOMEM);