Git Browser for ODROID
Code Review Sign In
git.odroid.com / yocto / linux / multimedia / 585b7207c4befb130b5762e563b0f127c9b9ecd4 / . / libfaad / ps_dec.c
blob: 142a5ae07a28529cbbf588698b26dec1b9535899 [file] [log] [blame]
xingri.gaoc18d4472023-02-28 02:51:02 +0000[diff] [blame]1/*
2** FAAD2 - Freeware Advanced Audio (AAC) Decoder including SBR decoding
3** Copyright (C) 2003-2005 M. Bakker, Nero AG, http://www.nero.com
4**
5** This program is free software; you can redistribute it and/or modify
6** it under the terms of the GNU General Public License as published by
7** the Free Software Foundation; either version 2 of the License, or
8** (at your option) any later version.
9**
10** This program is distributed in the hope that it will be useful,
11** but WITHOUT ANY WARRANTY; without even the implied warranty of
12** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13** GNU General Public License for more details.
14**
15** You should have received a copy of the GNU General Public License
16** along with this program; if not, write to the Free Software
17** Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
18**
19** Any non-GPL usage of this software or parts of this software is strictly
20** forbidden.
21**
22** The "appropriate copyright message" mentioned in section 2c of the GPLv2
23** must read: "Code from FAAD2 is copyright (c) Nero AG, www.nero.com"
24**
25** Commercial non-GPL licensing of this software is possible.
26** For more info contact Nero AG through Mpeg4AAClicense@nero.com.
27**
28** $Id: ps_dec.c,v 1.16 2009/01/26 22:32:31 menno Exp $
29**/
30
31#include <stdio.h>
32#include <stdlib.h>
33#include <string.h>
34#include <fcntl.h>
35#include "common.h"
36
37#ifdef PS_DEC
38
39#include "ps_dec.h"
40#include "ps_tables.h"
41
42/* constants */
43#define NEGATE_IPD_MASK (0x1000)
44#define DECAY_SLOPE FRAC_CONST(0.05)
45#define COEF_SQRT2 COEF_CONST(1.4142135623731)
46
47#ifndef __unused
48#define __unused __attribute__((unused))
49#endif
50
51/* tables */
52/* filters are mirrored in coef 6, second half left out */
53static const real_t p8_13_20[7] = {
54 FRAC_CONST(0.00746082949812),
55 FRAC_CONST(0.02270420949825),
56 FRAC_CONST(0.04546865930473),
57 FRAC_CONST(0.07266113929591),
58 FRAC_CONST(0.09885108575264),
59 FRAC_CONST(0.11793710567217),
60 FRAC_CONST(0.125)
61};
62
63static const real_t p2_13_20[7] = {
64 FRAC_CONST(0.0),
65 FRAC_CONST(0.01899487526049),
66 FRAC_CONST(0.0),
67 FRAC_CONST(-0.07293139167538),
68 FRAC_CONST(0.0),
69 FRAC_CONST(0.30596630545168),
70 FRAC_CONST(0.5)
71};
72
73static const real_t p12_13_34[7] = {
74 FRAC_CONST(0.04081179924692),
75 FRAC_CONST(0.03812810994926),
76 FRAC_CONST(0.05144908135699),
77 FRAC_CONST(0.06399831151592),
78 FRAC_CONST(0.07428313801106),
79 FRAC_CONST(0.08100347892914),
80 FRAC_CONST(0.08333333333333)
81};
82
83static const real_t p8_13_34[7] = {
84 FRAC_CONST(0.01565675600122),
85 FRAC_CONST(0.03752716391991),
86 FRAC_CONST(0.05417891378782),
87 FRAC_CONST(0.08417044116767),
88 FRAC_CONST(0.10307344158036),
89 FRAC_CONST(0.12222452249753),
90 FRAC_CONST(0.125)
91};
92
93static const real_t p4_13_34[7] = {
94 FRAC_CONST(-0.05908211155639),
95 FRAC_CONST(-0.04871498374946),
96 FRAC_CONST(0.0),
97 FRAC_CONST(0.07778723915851),
98 FRAC_CONST(0.16486303567403),
99 FRAC_CONST(0.23279856662996),
100 FRAC_CONST(0.25)
101};
102
103#ifdef PARAM_32KHZ
104static const uint8_t delay_length_d[2][NO_ALLPASS_LINKS] = {
105 { 1, 2, 3 } /* d_24kHz */,
106 { 3, 4, 5 } /* d_48kHz */
107};
108#else
109static const uint8_t delay_length_d[NO_ALLPASS_LINKS] = {
110 3, 4, 5 /* d_48kHz */
111};
112#endif
113static const real_t filter_a[NO_ALLPASS_LINKS] = { /* a(m) = exp(-d_48kHz(m)/7) */
114 FRAC_CONST(0.65143905753106),
115 FRAC_CONST(0.56471812200776),
116 FRAC_CONST(0.48954165955695)
117};
118
119static const uint8_t group_border20[10 + 12 + 1] = {
120 6, 7, 0, 1, 2, 3, /* 6 subqmf subbands */
121 9, 8, /* 2 subqmf subbands */
122 10, 11, /* 2 subqmf subbands */
123 3, 4, 5, 6, 7, 8, 9, 11, 14, 18, 23, 35, 64
124};
125
126static const uint8_t group_border34[32 + 18 + 1] = {
127 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, /* 12 subqmf subbands */
128 12, 13, 14, 15, 16, 17, 18, 19, /* 8 subqmf subbands */
129 20, 21, 22, 23, /* 4 subqmf subbands */
130 24, 25, 26, 27, /* 4 subqmf subbands */
131 28, 29, 30, 31, /* 4 subqmf subbands */
132 32 - 27, 33 - 27, 34 - 27, 35 - 27, 36 - 27, 37 - 27, 38 - 27, 40 - 27, 42 - 27, 44 - 27, 46 - 27, 48 - 27, 51 - 27, 54 - 27, 57 - 27, 60 - 27, 64 - 27, 68 - 27, 91 - 27
133};
134
135static const uint16_t map_group2bk20[10 + 12] = {
136 (NEGATE_IPD_MASK | 1), (NEGATE_IPD_MASK | 0),
137 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19
138};
139
140static const uint16_t map_group2bk34[32 + 18] = {
141 0, 1, 2, 3, 4, 5, 6, 6, 7, (NEGATE_IPD_MASK | 2), (NEGATE_IPD_MASK | 1), (NEGATE_IPD_MASK | 0),
142 10, 10, 4, 5, 6, 7, 8, 9,
143 10, 11, 12, 9,
144 14, 11, 12, 13,
145 14, 15, 16, 13,
146 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33
147};
148
149/* type definitions */
150typedef struct {
151 uint8_t frame_len;
152 uint8_t resolution20[3];
153 uint8_t resolution34[5];
154
155 qmf_t *work;
156 qmf_t **buffer;
157 qmf_t **temp;
158} hyb_info;
159
160/* static function declarations */
161static void ps_data_decode(ps_info *ps);
162static hyb_info *hybrid_init(uint8_t numTimeSlotsRate);
163static void channel_filter2(hyb_info *hyb, uint8_t frame_len, const real_t *filter,
164 qmf_t *buffer, qmf_t **X_hybrid);
165static void INLINE DCT3_4_unscaled(real_t *y, real_t *x);
166static void channel_filter8(hyb_info *hyb, uint8_t frame_len, const real_t *filter,
167 qmf_t *buffer, qmf_t **X_hybrid);
168static void hybrid_analysis(hyb_info *hyb, qmf_t X[32][64], qmf_t X_hybrid[32][32],
169 uint8_t use34, uint8_t numTimeSlotsRate);
170static void hybrid_synthesis(hyb_info *hyb, qmf_t X[32][64], qmf_t X_hybrid[32][32],
171 uint8_t use34, uint8_t numTimeSlotsRate);
172static int8_t delta_clip(int8_t i, int8_t min, int8_t max);
173static void delta_decode(uint8_t enable, int8_t *index, int8_t *index_prev,
174 uint8_t dt_flag, uint8_t nr_par, uint8_t stride,
175 int8_t min_index, int8_t max_index);
176static void delta_modulo_decode(uint8_t enable, int8_t *index, int8_t *index_prev,
177 uint8_t dt_flag, uint8_t nr_par, uint8_t stride,
178 int8_t and_modulo);
179static void map20indexto34(int8_t *index, uint8_t bins);
180#ifdef PS_LOW_POWER
181static void map34indexto20(int8_t *index, uint8_t bins);
182#endif
183static void ps_data_decode(ps_info *ps);
184static void ps_decorrelate(ps_info *ps, qmf_t X_left[38][64], qmf_t X_right[38][64],
185 qmf_t X_hybrid_left[32][32], qmf_t X_hybrid_right[32][32]);
186static void ps_mix_phase(ps_info *ps, qmf_t X_left[38][64], qmf_t X_right[38][64],
187 qmf_t X_hybrid_left[32][32], qmf_t X_hybrid_right[32][32]);
188
189/* */
190
191
192static hyb_info *hybrid_init(uint8_t numTimeSlotsRate)
193{
194 uint8_t i;
195
196 hyb_info *hyb = (hyb_info*)faad_malloc(sizeof(hyb_info));
197
198 hyb->resolution34[0] = 12;
199 hyb->resolution34[1] = 8;
200 hyb->resolution34[2] = 4;
201 hyb->resolution34[3] = 4;
202 hyb->resolution34[4] = 4;
203
204 hyb->resolution20[0] = 8;
205 hyb->resolution20[1] = 2;
206 hyb->resolution20[2] = 2;
207
208 hyb->frame_len = numTimeSlotsRate;
209
210 hyb->work = (qmf_t*)faad_malloc((hyb->frame_len + 12) * sizeof(qmf_t));
211 memset(hyb->work, 0, (hyb->frame_len + 12) * sizeof(qmf_t));
212
213 hyb->buffer = (qmf_t**)faad_malloc(5 * sizeof(qmf_t*));
214 for (i = 0; i < 5; i++) {
215 hyb->buffer[i] = (qmf_t*)faad_malloc(hyb->frame_len * sizeof(qmf_t));
216 memset(hyb->buffer[i], 0, hyb->frame_len * sizeof(qmf_t));
217 }
218
219 hyb->temp = (qmf_t**)faad_malloc(hyb->frame_len * sizeof(qmf_t*));
220 for (i = 0; i < hyb->frame_len; i++) {
221 hyb->temp[i] = (qmf_t*)faad_malloc(12 /*max*/ * sizeof(qmf_t));
222 }
223
224 return hyb;
225}
226
227static void hybrid_free(hyb_info *hyb)
228{
229 uint8_t i;
230
231 if (!hyb) {
232 return;
233 }
234
235 if (hyb->work) {
236 faad_free(hyb->work);
237 }
238
239 for (i = 0; i < 5; i++) {
240 if (hyb->buffer[i]) {
241 faad_free(hyb->buffer[i]);
242 }
243 }
244 if (hyb->buffer) {
245 faad_free(hyb->buffer);
246 }
247
248 for (i = 0; i < hyb->frame_len; i++) {
249 if (hyb->temp[i]) {
250 faad_free(hyb->temp[i]);
251 }
252 }
253 if (hyb->temp) {
254 faad_free(hyb->temp);
255 }
256
257 faad_free(hyb);
258}
259
260/* real filter, size 2 */
261static void channel_filter2(hyb_info *hyb __unused, uint8_t frame_len, const real_t *filter,
262 qmf_t *buffer, qmf_t **X_hybrid)
263{
264 uint8_t i;
265
266 for (i = 0; i < frame_len; i++) {
267 real_t r0 = MUL_F(filter[0], (QMF_RE(buffer[0 + i]) + QMF_RE(buffer[12 + i])));
268 real_t r1 = MUL_F(filter[1], (QMF_RE(buffer[1 + i]) + QMF_RE(buffer[11 + i])));
269 real_t r2 = MUL_F(filter[2], (QMF_RE(buffer[2 + i]) + QMF_RE(buffer[10 + i])));
270 real_t r3 = MUL_F(filter[3], (QMF_RE(buffer[3 + i]) + QMF_RE(buffer[9 + i])));
271 real_t r4 = MUL_F(filter[4], (QMF_RE(buffer[4 + i]) + QMF_RE(buffer[8 + i])));
272 real_t r5 = MUL_F(filter[5], (QMF_RE(buffer[5 + i]) + QMF_RE(buffer[7 + i])));
273 real_t r6 = MUL_F(filter[6], QMF_RE(buffer[6 + i]));
274 real_t i0 = MUL_F(filter[0], (QMF_IM(buffer[0 + i]) + QMF_IM(buffer[12 + i])));
275 real_t i1 = MUL_F(filter[1], (QMF_IM(buffer[1 + i]) + QMF_IM(buffer[11 + i])));
276 real_t i2 = MUL_F(filter[2], (QMF_IM(buffer[2 + i]) + QMF_IM(buffer[10 + i])));
277 real_t i3 = MUL_F(filter[3], (QMF_IM(buffer[3 + i]) + QMF_IM(buffer[9 + i])));
278 real_t i4 = MUL_F(filter[4], (QMF_IM(buffer[4 + i]) + QMF_IM(buffer[8 + i])));
279 real_t i5 = MUL_F(filter[5], (QMF_IM(buffer[5 + i]) + QMF_IM(buffer[7 + i])));
280 real_t i6 = MUL_F(filter[6], QMF_IM(buffer[6 + i]));
281
282 /* q = 0 */
283 QMF_RE(X_hybrid[i][0]) = r0 + r1 + r2 + r3 + r4 + r5 + r6;
284 QMF_IM(X_hybrid[i][0]) = i0 + i1 + i2 + i3 + i4 + i5 + i6;
285
286 /* q = 1 */
287 QMF_RE(X_hybrid[i][1]) = r0 - r1 + r2 - r3 + r4 - r5 + r6;
288 QMF_IM(X_hybrid[i][1]) = i0 - i1 + i2 - i3 + i4 - i5 + i6;
289 }
290}
291
292/* complex filter, size 4 */
293static void channel_filter4(hyb_info *hyb __unused, uint8_t frame_len, const real_t *filter,
294 qmf_t *buffer, qmf_t **X_hybrid)
295{
296 uint8_t i;
297 real_t input_re1[2], input_re2[2], input_im1[2], input_im2[2];
298
299 for (i = 0; i < frame_len; i++) {
300 input_re1[0] = -MUL_F(filter[2], (QMF_RE(buffer[i + 2]) + QMF_RE(buffer[i + 10]))) +
301 MUL_F(filter[6], QMF_RE(buffer[i + 6]));
302 input_re1[1] = MUL_F(FRAC_CONST(-0.70710678118655),
303 (MUL_F(filter[1], (QMF_RE(buffer[i + 1]) + QMF_RE(buffer[i + 11]))) +
304 MUL_F(filter[3], (QMF_RE(buffer[i + 3]) + QMF_RE(buffer[i + 9]))) -
305 MUL_F(filter[5], (QMF_RE(buffer[i + 5]) + QMF_RE(buffer[i + 7])))));
306
307 input_im1[0] = MUL_F(filter[0], (QMF_IM(buffer[i + 0]) - QMF_IM(buffer[i + 12]))) -
308 MUL_F(filter[4], (QMF_IM(buffer[i + 4]) - QMF_IM(buffer[i + 8])));
309 input_im1[1] = MUL_F(FRAC_CONST(0.70710678118655),
310 (MUL_F(filter[1], (QMF_IM(buffer[i + 1]) - QMF_IM(buffer[i + 11]))) -
311 MUL_F(filter[3], (QMF_IM(buffer[i + 3]) - QMF_IM(buffer[i + 9]))) -
312 MUL_F(filter[5], (QMF_IM(buffer[i + 5]) - QMF_IM(buffer[i + 7])))));
313
314 input_re2[0] = MUL_F(filter[0], (QMF_RE(buffer[i + 0]) - QMF_RE(buffer[i + 12]))) -
315 MUL_F(filter[4], (QMF_RE(buffer[i + 4]) - QMF_RE(buffer[i + 8])));
316 input_re2[1] = MUL_F(FRAC_CONST(0.70710678118655),
317 (MUL_F(filter[1], (QMF_RE(buffer[i + 1]) - QMF_RE(buffer[i + 11]))) -
318 MUL_F(filter[3], (QMF_RE(buffer[i + 3]) - QMF_RE(buffer[i + 9]))) -
319 MUL_F(filter[5], (QMF_RE(buffer[i + 5]) - QMF_RE(buffer[i + 7])))));
320
321 input_im2[0] = -MUL_F(filter[2], (QMF_IM(buffer[i + 2]) + QMF_IM(buffer[i + 10]))) +
322 MUL_F(filter[6], QMF_IM(buffer[i + 6]));
323 input_im2[1] = MUL_F(FRAC_CONST(-0.70710678118655),
324 (MUL_F(filter[1], (QMF_IM(buffer[i + 1]) + QMF_IM(buffer[i + 11]))) +
325 MUL_F(filter[3], (QMF_IM(buffer[i + 3]) + QMF_IM(buffer[i + 9]))) -
326 MUL_F(filter[5], (QMF_IM(buffer[i + 5]) + QMF_IM(buffer[i + 7])))));
327
328 /* q == 0 */
329 QMF_RE(X_hybrid[i][0]) = input_re1[0] + input_re1[1] + input_im1[0] + input_im1[1];
330 QMF_IM(X_hybrid[i][0]) = -input_re2[0] - input_re2[1] + input_im2[0] + input_im2[1];
331
332 /* q == 1 */
333 QMF_RE(X_hybrid[i][1]) = input_re1[0] - input_re1[1] - input_im1[0] + input_im1[1];
334 QMF_IM(X_hybrid[i][1]) = input_re2[0] - input_re2[1] + input_im2[0] - input_im2[1];
335
336 /* q == 2 */
337 QMF_RE(X_hybrid[i][2]) = input_re1[0] - input_re1[1] + input_im1[0] - input_im1[1];
338 QMF_IM(X_hybrid[i][2]) = -input_re2[0] + input_re2[1] + input_im2[0] - input_im2[1];
339
340 /* q == 3 */
341 QMF_RE(X_hybrid[i][3]) = input_re1[0] + input_re1[1] - input_im1[0] - input_im1[1];
342 QMF_IM(X_hybrid[i][3]) = input_re2[0] + input_re2[1] + input_im2[0] + input_im2[1];
343 }
344}
345
346static void INLINE DCT3_4_unscaled(real_t *y, real_t *x)
347{
348 real_t f0, f1, f2, f3, f4, f5, f6, f7, f8;
349
350 f0 = MUL_F(x[2], FRAC_CONST(0.7071067811865476));
351 f1 = x[0] - f0;
352 f2 = x[0] + f0;
353 f3 = x[1] + x[3];
354 f4 = MUL_C(x[1], COEF_CONST(1.3065629648763766));
355 f5 = MUL_F(f3, FRAC_CONST(-0.9238795325112866));
356 f6 = MUL_F(x[3], FRAC_CONST(-0.5411961001461967));
357 f7 = f4 + f5;
358 f8 = f6 - f5;
359 y[3] = f2 - f8;
360 y[0] = f2 + f8;
361 y[2] = f1 - f7;
362 y[1] = f1 + f7;
363}
364
365/* complex filter, size 8 */
366static void channel_filter8(hyb_info *hyb __unused, uint8_t frame_len, const real_t *filter,
367 qmf_t *buffer, qmf_t **X_hybrid)
368{
369 uint8_t i, n;
370 real_t input_re1[4], input_re2[4], input_im1[4], input_im2[4];
371 real_t x[4];
372
373 for (i = 0; i < frame_len; i++) {
374 input_re1[0] = MUL_F(filter[6], QMF_RE(buffer[6 + i]));
375 input_re1[1] = MUL_F(filter[5], (QMF_RE(buffer[5 + i]) + QMF_RE(buffer[7 + i])));
376 input_re1[2] = -MUL_F(filter[0], (QMF_RE(buffer[0 + i]) + QMF_RE(buffer[12 + i]))) + MUL_F(filter[4], (QMF_RE(buffer[4 + i]) + QMF_RE(buffer[8 + i])));
377 input_re1[3] = -MUL_F(filter[1], (QMF_RE(buffer[1 + i]) + QMF_RE(buffer[11 + i]))) + MUL_F(filter[3], (QMF_RE(buffer[3 + i]) + QMF_RE(buffer[9 + i])));
378
379 input_im1[0] = MUL_F(filter[5], (QMF_IM(buffer[7 + i]) - QMF_IM(buffer[5 + i])));
380 input_im1[1] = MUL_F(filter[0], (QMF_IM(buffer[12 + i]) - QMF_IM(buffer[0 + i]))) + MUL_F(filter[4], (QMF_IM(buffer[8 + i]) - QMF_IM(buffer[4 + i])));
381 input_im1[2] = MUL_F(filter[1], (QMF_IM(buffer[11 + i]) - QMF_IM(buffer[1 + i]))) + MUL_F(filter[3], (QMF_IM(buffer[9 + i]) - QMF_IM(buffer[3 + i])));
382 input_im1[3] = MUL_F(filter[2], (QMF_IM(buffer[10 + i]) - QMF_IM(buffer[2 + i])));
383
384 for (n = 0; n < 4; n++) {
385 x[n] = input_re1[n] - input_im1[3 - n];
386 }
387 DCT3_4_unscaled(x, x);
388 QMF_RE(X_hybrid[i][7]) = x[0];
389 QMF_RE(X_hybrid[i][5]) = x[2];
390 QMF_RE(X_hybrid[i][3]) = x[3];
391 QMF_RE(X_hybrid[i][1]) = x[1];
392
393 for (n = 0; n < 4; n++) {
394 x[n] = input_re1[n] + input_im1[3 - n];
395 }
396 DCT3_4_unscaled(x, x);
397 QMF_RE(X_hybrid[i][6]) = x[1];
398 QMF_RE(X_hybrid[i][4]) = x[3];
399 QMF_RE(X_hybrid[i][2]) = x[2];
400 QMF_RE(X_hybrid[i][0]) = x[0];
401
402 input_im2[0] = MUL_F(filter[6], QMF_IM(buffer[6 + i]));
403 input_im2[1] = MUL_F(filter[5], (QMF_IM(buffer[5 + i]) + QMF_IM(buffer[7 + i])));
404 input_im2[2] = -MUL_F(filter[0], (QMF_IM(buffer[0 + i]) + QMF_IM(buffer[12 + i]))) + MUL_F(filter[4], (QMF_IM(buffer[4 + i]) + QMF_IM(buffer[8 + i])));
405 input_im2[3] = -MUL_F(filter[1], (QMF_IM(buffer[1 + i]) + QMF_IM(buffer[11 + i]))) + MUL_F(filter[3], (QMF_IM(buffer[3 + i]) + QMF_IM(buffer[9 + i])));
406
407 input_re2[0] = MUL_F(filter[5], (QMF_RE(buffer[7 + i]) - QMF_RE(buffer[5 + i])));
408 input_re2[1] = MUL_F(filter[0], (QMF_RE(buffer[12 + i]) - QMF_RE(buffer[0 + i]))) + MUL_F(filter[4], (QMF_RE(buffer[8 + i]) - QMF_RE(buffer[4 + i])));
409 input_re2[2] = MUL_F(filter[1], (QMF_RE(buffer[11 + i]) - QMF_RE(buffer[1 + i]))) + MUL_F(filter[3], (QMF_RE(buffer[9 + i]) - QMF_RE(buffer[3 + i])));
410 input_re2[3] = MUL_F(filter[2], (QMF_RE(buffer[10 + i]) - QMF_RE(buffer[2 + i])));
411
412 for (n = 0; n < 4; n++) {
413 x[n] = input_im2[n] + input_re2[3 - n];
414 }
415 DCT3_4_unscaled(x, x);
416 QMF_IM(X_hybrid[i][7]) = x[0];
417 QMF_IM(X_hybrid[i][5]) = x[2];
418 QMF_IM(X_hybrid[i][3]) = x[3];
419 QMF_IM(X_hybrid[i][1]) = x[1];
420
421 for (n = 0; n < 4; n++) {
422 x[n] = input_im2[n] - input_re2[3 - n];
423 }
424 DCT3_4_unscaled(x, x);
425 QMF_IM(X_hybrid[i][6]) = x[1];
426 QMF_IM(X_hybrid[i][4]) = x[3];
427 QMF_IM(X_hybrid[i][2]) = x[2];
428 QMF_IM(X_hybrid[i][0]) = x[0];
429 }
430}
431
432static void INLINE DCT3_6_unscaled(real_t *y, real_t *x)
433{
434 real_t f0, f1, f2, f3, f4, f5, f6, f7;
435
436 f0 = MUL_F(x[3], FRAC_CONST(0.70710678118655));
437 f1 = x[0] + f0;
438 f2 = x[0] - f0;
439 f3 = MUL_F((x[1] - x[5]), FRAC_CONST(0.70710678118655));
440 f4 = MUL_F(x[2], FRAC_CONST(0.86602540378444)) + MUL_F(x[4], FRAC_CONST(0.5));
441 f5 = f4 - x[4];
442 f6 = MUL_F(x[1], FRAC_CONST(0.96592582628907)) + MUL_F(x[5], FRAC_CONST(0.25881904510252));
443 f7 = f6 - f3;
444 y[0] = f1 + f6 + f4;
445 y[1] = f2 + f3 - x[4];
446 y[2] = f7 + f2 - f5;
447 y[3] = f1 - f7 - f5;
448 y[4] = f1 - f3 - x[4];
449 y[5] = f2 - f6 + f4;
450}
451
452/* complex filter, size 12 */
453static void channel_filter12(hyb_info *hyb __unused, uint8_t frame_len, const real_t *filter,
454 qmf_t *buffer, qmf_t **X_hybrid)
455{
456 uint8_t i, n;
457 real_t input_re1[6], input_re2[6], input_im1[6], input_im2[6];
458 real_t out_re1[6], out_re2[6], out_im1[6], out_im2[6];
459
460 for (i = 0; i < frame_len; i++) {
461 for (n = 0; n < 6; n++) {
462 if (n == 0) {
463 input_re1[0] = MUL_F(QMF_RE(buffer[6 + i]), filter[6]);
464 input_re2[0] = MUL_F(QMF_IM(buffer[6 + i]), filter[6]);
465 } else {
466 input_re1[6 - n] = MUL_F((QMF_RE(buffer[n + i]) + QMF_RE(buffer[12 - n + i])), filter[n]);
467 input_re2[6 - n] = MUL_F((QMF_IM(buffer[n + i]) + QMF_IM(buffer[12 - n + i])), filter[n]);
468 }
469 input_im2[n] = MUL_F((QMF_RE(buffer[n + i]) - QMF_RE(buffer[12 - n + i])), filter[n]);
470 input_im1[n] = MUL_F((QMF_IM(buffer[n + i]) - QMF_IM(buffer[12 - n + i])), filter[n]);
471 }
472
473 DCT3_6_unscaled(out_re1, input_re1);
474 DCT3_6_unscaled(out_re2, input_re2);
475
476 DCT3_6_unscaled(out_im1, input_im1);
477 DCT3_6_unscaled(out_im2, input_im2);
478
479 for (n = 0; n < 6; n += 2) {
480 QMF_RE(X_hybrid[i][n]) = out_re1[n] - out_im1[n];
481 QMF_IM(X_hybrid[i][n]) = out_re2[n] + out_im2[n];
482 QMF_RE(X_hybrid[i][n + 1]) = out_re1[n + 1] + out_im1[n + 1];
483 QMF_IM(X_hybrid[i][n + 1]) = out_re2[n + 1] - out_im2[n + 1];
484
485 QMF_RE(X_hybrid[i][10 - n]) = out_re1[n + 1] - out_im1[n + 1];
486 QMF_IM(X_hybrid[i][10 - n]) = out_re2[n + 1] + out_im2[n + 1];
487 QMF_RE(X_hybrid[i][11 - n]) = out_re1[n] + out_im1[n];
488 QMF_IM(X_hybrid[i][11 - n]) = out_re2[n] - out_im2[n];
489 }
490 }
491}
492
493/* Hybrid analysis: further split up QMF subbands
494 * to improve frequency resolution
495 */
496static void hybrid_analysis(hyb_info *hyb, qmf_t X[32][64], qmf_t X_hybrid[32][32],
497 uint8_t use34, uint8_t numTimeSlotsRate)
498{
499 uint8_t k, n, band;
500 uint8_t offset = 0;
501 uint8_t qmf_bands = (use34) ? 5 : 3;
502 uint8_t *resolution = (use34) ? hyb->resolution34 : hyb->resolution20;
503
504 for (band = 0; band < qmf_bands; band++) {
505 /* build working buffer */
506 memcpy(hyb->work, hyb->buffer[band], 12 * sizeof(qmf_t));
507
508 /* add new samples */
509 for (n = 0; n < hyb->frame_len; n++) {
510 QMF_RE(hyb->work[12 + n]) = QMF_RE(X[n + 6 /*delay*/][band]);
511 QMF_IM(hyb->work[12 + n]) = QMF_IM(X[n + 6 /*delay*/][band]);
512 }
513
514 /* store samples */
515 memcpy(hyb->buffer[band], hyb->work + hyb->frame_len, 12 * sizeof(qmf_t));
516
517
518 switch (resolution[band]) {
519 case 2:
520 /* Type B real filter, Q[p] = 2 */
521 channel_filter2(hyb, hyb->frame_len, p2_13_20, hyb->work, hyb->temp);
522 break;
523 case 4:
524 /* Type A complex filter, Q[p] = 4 */
525 channel_filter4(hyb, hyb->frame_len, p4_13_34, hyb->work, hyb->temp);
526 break;
527 case 8:
528 /* Type A complex filter, Q[p] = 8 */
529 channel_filter8(hyb, hyb->frame_len, (use34) ? p8_13_34 : p8_13_20,
530 hyb->work, hyb->temp);
531 break;
532 case 12:
533 /* Type A complex filter, Q[p] = 12 */
534 channel_filter12(hyb, hyb->frame_len, p12_13_34, hyb->work, hyb->temp);
535 break;
536 }
537
538 for (n = 0; n < hyb->frame_len; n++) {
539 for (k = 0; k < resolution[band]; k++) {
540 QMF_RE(X_hybrid[n][offset + k]) = QMF_RE(hyb->temp[n][k]);
541 QMF_IM(X_hybrid[n][offset + k]) = QMF_IM(hyb->temp[n][k]);
542 }
543 }
544 offset += resolution[band];
545 }
546
547 /* group hybrid channels */
548 if (!use34) {
549 for (n = 0; n < numTimeSlotsRate; n++) {
550 QMF_RE(X_hybrid[n][3]) += QMF_RE(X_hybrid[n][4]);
551 QMF_IM(X_hybrid[n][3]) += QMF_IM(X_hybrid[n][4]);
552 QMF_RE(X_hybrid[n][4]) = 0;
553 QMF_IM(X_hybrid[n][4]) = 0;
554
555 QMF_RE(X_hybrid[n][2]) += QMF_RE(X_hybrid[n][5]);
556 QMF_IM(X_hybrid[n][2]) += QMF_IM(X_hybrid[n][5]);
557 QMF_RE(X_hybrid[n][5]) = 0;
558 QMF_IM(X_hybrid[n][5]) = 0;
559 }
560 }
561}
562
563static void hybrid_synthesis(hyb_info *hyb, qmf_t X[32][64], qmf_t X_hybrid[32][32],
564 uint8_t use34, uint8_t numTimeSlotsRate __unused)
565{
566 uint8_t k, n, band;
567 uint8_t offset = 0;
568 uint8_t qmf_bands = (use34) ? 5 : 3;
569 uint8_t *resolution = (use34) ? hyb->resolution34 : hyb->resolution20;
570
571 for (band = 0; band < qmf_bands; band++) {
572 for (n = 0; n < hyb->frame_len; n++) {
573 QMF_RE(X[n][band]) = 0;
574 QMF_IM(X[n][band]) = 0;
575
576 for (k = 0; k < resolution[band]; k++) {
577 QMF_RE(X[n][band]) += QMF_RE(X_hybrid[n][offset + k]);
578 QMF_IM(X[n][band]) += QMF_IM(X_hybrid[n][offset + k]);
579 }
580 }
581 offset += resolution[band];
582 }
583}
584
585/* limits the value i to the range [min,max] */
586static int8_t delta_clip(int8_t i, int8_t min, int8_t max)
587{
588 if (i < min) {
589 return min;
590 } else if (i > max) {
591 return max;
592 } else {
593 return i;
594 }
595}
596
597//int iid = 0;
598
599/* delta decode array */
600static void delta_decode(uint8_t enable, int8_t *index, int8_t *index_prev,
601 uint8_t dt_flag, uint8_t nr_par, uint8_t stride,
602 int8_t min_index, int8_t max_index)
603{
604 int8_t i;
605
606 if (enable == 1) {
607 if (dt_flag == 0) {
608 /* delta coded in frequency direction */
609 index[0] = 0 + index[0];
610 index[0] = delta_clip(index[0], min_index, max_index);
611
612 for (i = 1; i < nr_par; i++) {
613 index[i] = index[i - 1] + index[i];
614 index[i] = delta_clip(index[i], min_index, max_index);
615 }
616 } else {
617 /* delta coded in time direction */
618 for (i = 0; i < nr_par; i++) {
619 //int8_t tmp2;
620 //int8_t tmp = index[i];
621
622 //printf("%d %d\n", index_prev[i*stride], index[i]);
623 //printf("%d\n", index[i]);
624
625 index[i] = index_prev[i * stride] + index[i];
626 //tmp2 = index[i];
627 index[i] = delta_clip(index[i], min_index, max_index);
628
629 //if (iid)
630 //{
631 // if (index[i] == 7)
632 // {
633 // printf("%d %d %d\n", index_prev[i*stride], tmp, tmp2);
634 // }
635 //}
636 }
637 }
638 } else {
639 /* set indices to zero */
640 for (i = 0; i < nr_par; i++) {
641 index[i] = 0;
642 }
643 }
644
645 /* coarse */
646 if (stride == 2) {
647 for (i = (nr_par << 1) - 1; i > 0; i--) {
648 index[i] = index[i >> 1];
649 }
650 }
651}
652
653/* delta modulo decode array */
654/* in: log2 value of the modulo value to allow using AND instead of MOD */
655static void delta_modulo_decode(uint8_t enable, int8_t *index, int8_t *index_prev,
656 uint8_t dt_flag, uint8_t nr_par, uint8_t stride,
657 int8_t and_modulo)
658{
659 int8_t i;
660
661 if (enable == 1) {
662 if (dt_flag == 0) {
663 /* delta coded in frequency direction */
664 index[0] = 0 + index[0];
665 index[0] &= and_modulo;
666
667 for (i = 1; i < nr_par; i++) {
668 index[i] = index[i - 1] + index[i];
669 index[i] &= and_modulo;
670 }
671 } else {
672 /* delta coded in time direction */
673 for (i = 0; i < nr_par; i++) {
674 index[i] = index_prev[i * stride] + index[i];
675 index[i] &= and_modulo;
676 }
677 }
678 } else {
679 /* set indices to zero */
680 for (i = 0; i < nr_par; i++) {
681 index[i] = 0;
682 }
683 }
684
685 /* coarse */
686 if (stride == 2) {
687 index[0] = 0;
688 for (i = (nr_par << 1) - 1; i > 0; i--) {
689 index[i] = index[i >> 1];
690 }
691 }
692}
693
694#ifdef PS_LOW_POWER
695static void map34indexto20(int8_t *index, uint8_t bins)
696{
697 index[0] = (2 * index[0] + index[1]) / 3;
698 index[1] = (index[1] + 2 * index[2]) / 3;
699 index[2] = (2 * index[3] + index[4]) / 3;
700 index[3] = (index[4] + 2 * index[5]) / 3;
701 index[4] = (index[6] + index[7]) / 2;
702 index[5] = (index[8] + index[9]) / 2;
703 index[6] = index[10];
704 index[7] = index[11];
705 index[8] = (index[12] + index[13]) / 2;
706 index[9] = (index[14] + index[15]) / 2;
707 index[10] = index[16];
708
709 if (bins == 34) {
710 index[11] = index[17];
711 index[12] = index[18];
712 index[13] = index[19];
713 index[14] = (index[20] + index[21]) / 2;
714 index[15] = (index[22] + index[23]) / 2;
715 index[16] = (index[24] + index[25]) / 2;
716 index[17] = (index[26] + index[27]) / 2;
717 index[18] = (index[28] + index[29] + index[30] + index[31]) / 4;
718 index[19] = (index[32] + index[33]) / 2;
719 }
720}
721#endif
722
723static void map20indexto34(int8_t *index, uint8_t bins)
724{
yuliang.hu77fab5b2024-07-19 18:32:11 +0800[diff] [blame]725 /*coverity[no_effect]*/
xingri.gaoc18d4472023-02-28 02:51:02 +0000[diff] [blame]726 index[0] = index[0];
727 index[1] = (index[0] + index[1]) / 2;
728 index[2] = index[1];
729 index[3] = index[2];
730 index[4] = (index[2] + index[3]) / 2;
731 index[5] = index[3];
732 index[6] = index[4];
733 index[7] = index[4];
734 index[8] = index[5];
735 index[9] = index[5];
736 index[10] = index[6];
737 index[11] = index[7];
738 index[12] = index[8];
739 index[13] = index[8];
740 index[14] = index[9];
741 index[15] = index[9];
742 index[16] = index[10];
743
744 if (bins == 34) {
745 index[17] = index[11];
746 index[18] = index[12];
747 index[19] = index[13];
748 index[20] = index[14];
749 index[21] = index[14];
750 index[22] = index[15];
751 index[23] = index[15];
752 index[24] = index[16];
753 index[25] = index[16];
754 index[26] = index[17];
755 index[27] = index[17];
756 index[28] = index[18];
757 index[29] = index[18];
758 index[30] = index[18];
759 index[31] = index[18];
760 index[32] = index[19];
761 index[33] = index[19];
762 }
763}
764
765/* parse the bitstream data decoded in ps_data() */
766static void ps_data_decode(ps_info *ps)
767{
768 uint8_t env, bin;
769
770 /* ps data not available, use data from previous frame */
771 if (ps->ps_data_available == 0) {
772 ps->num_env = 0;
773 }
774
775 for (env = 0; env < ps->num_env; env++) {
776 int8_t *iid_index_prev;
777 int8_t *icc_index_prev;
778 int8_t *ipd_index_prev;
779 int8_t *opd_index_prev;
780
781 int8_t num_iid_steps = (ps->iid_mode < 3) ? 7 : 15 /*fine quant*/;
782
783 if (env == 0) {
784 /* take last envelope from previous frame */
785 iid_index_prev = ps->iid_index_prev;
786 icc_index_prev = ps->icc_index_prev;
787 ipd_index_prev = ps->ipd_index_prev;
788 opd_index_prev = ps->opd_index_prev;
789 } else {
790 /* take index values from previous envelope */
791 iid_index_prev = ps->iid_index[env - 1];
792 icc_index_prev = ps->icc_index[env - 1];
793 ipd_index_prev = ps->ipd_index[env - 1];
794 opd_index_prev = ps->opd_index[env - 1];
795 }
796
797 // iid = 1;
798 /* delta decode iid parameters */
799 delta_decode(ps->enable_iid, ps->iid_index[env], iid_index_prev,
800 ps->iid_dt[env], ps->nr_iid_par,
801 (ps->iid_mode == 0 || ps->iid_mode == 3) ? 2 : 1,
802 -num_iid_steps, num_iid_steps);
803 // iid = 0;
804
805 /* delta decode icc parameters */
806 delta_decode(ps->enable_icc, ps->icc_index[env], icc_index_prev,
807 ps->icc_dt[env], ps->nr_icc_par,
808 (ps->icc_mode == 0 || ps->icc_mode == 3) ? 2 : 1,
809 0, 7);
810
811 /* delta modulo decode ipd parameters */
812 delta_modulo_decode(ps->enable_ipdopd, ps->ipd_index[env], ipd_index_prev,
813 ps->ipd_dt[env], ps->nr_ipdopd_par, 1, 7);
814
815 /* delta modulo decode opd parameters */
816 delta_modulo_decode(ps->enable_ipdopd, ps->opd_index[env], opd_index_prev,
817 ps->opd_dt[env], ps->nr_ipdopd_par, 1, 7);
818 }
819
820 /* handle error case */
821 if (ps->num_env == 0) {
822 /* force to 1 */
823 ps->num_env = 1;
824
825 if (ps->enable_iid) {
826 for (bin = 0; bin < 34; bin++) {
827 ps->iid_index[0][bin] = ps->iid_index_prev[bin];
828 }
829 } else {
830 for (bin = 0; bin < 34; bin++) {
831 ps->iid_index[0][bin] = 0;
832 }
833 }
834
835 if (ps->enable_icc) {
836 for (bin = 0; bin < 34; bin++) {
837 ps->icc_index[0][bin] = ps->icc_index_prev[bin];
838 }
839 } else {
840 for (bin = 0; bin < 34; bin++) {
841 ps->icc_index[0][bin] = 0;
842 }
843 }
844
845 if (ps->enable_ipdopd) {
846 for (bin = 0; bin < 17; bin++) {
847 ps->ipd_index[0][bin] = ps->ipd_index_prev[bin];
848 ps->opd_index[0][bin] = ps->opd_index_prev[bin];
849 }
850 } else {
851 for (bin = 0; bin < 17; bin++) {
852 ps->ipd_index[0][bin] = 0;
853 ps->opd_index[0][bin] = 0;
854 }
855 }
856 }
857
858 /* update previous indices */
859 for (bin = 0; bin < 34; bin++) {
860 ps->iid_index_prev[bin] = ps->iid_index[ps->num_env - 1][bin];
861 }
862 for (bin = 0; bin < 34; bin++) {
863 ps->icc_index_prev[bin] = ps->icc_index[ps->num_env - 1][bin];
864 }
865 for (bin = 0; bin < 17; bin++) {
866 ps->ipd_index_prev[bin] = ps->ipd_index[ps->num_env - 1][bin];
867 ps->opd_index_prev[bin] = ps->opd_index[ps->num_env - 1][bin];
868 }
869
870 ps->ps_data_available = 0;
871
872 if (ps->frame_class == 0) {
873 ps->border_position[0] = 0;
874 for (env = 1; env < ps->num_env; env++) {
875 ps->border_position[env] = (env * ps->numTimeSlotsRate) / ps->num_env;
876 }
877 ps->border_position[ps->num_env] = ps->numTimeSlotsRate;
878 } else {
879 ps->border_position[0] = 0;
880
881 if (ps->border_position[ps->num_env] < ps->numTimeSlotsRate) {
882 for (bin = 0; bin < 34; bin++) {
883 ps->iid_index[ps->num_env][bin] = ps->iid_index[ps->num_env - 1][bin];
884 ps->icc_index[ps->num_env][bin] = ps->icc_index[ps->num_env - 1][bin];
885 }
886 for (bin = 0; bin < 17; bin++) {
887 ps->ipd_index[ps->num_env][bin] = ps->ipd_index[ps->num_env - 1][bin];
888 ps->opd_index[ps->num_env][bin] = ps->opd_index[ps->num_env - 1][bin];
889 }
890 ps->num_env++;
891 ps->border_position[ps->num_env] = ps->numTimeSlotsRate;
892 }
893
894 for (env = 1; env < ps->num_env; env++) {
895 int8_t thr = ps->numTimeSlotsRate - (ps->num_env - env);
896
897 if (ps->border_position[env] > thr) {
898 ps->border_position[env] = thr;
899 } else {
900 thr = ps->border_position[env - 1] + 1;
901 if (ps->border_position[env] < thr) {
902 ps->border_position[env] = thr;
903 }
904 }
905 }
906 }
907
908 /* make sure that the indices of all parameters can be mapped
909 * to the same hybrid synthesis filterbank
910 */
911#ifdef PS_LOW_POWER
912 for (env = 0; env < ps->num_env; env++) {
913 if (ps->iid_mode == 2 || ps->iid_mode == 5) {
914 map34indexto20(ps->iid_index[env], 34);
915 }
916 if (ps->icc_mode == 2 || ps->icc_mode == 5) {
917 map34indexto20(ps->icc_index[env], 34);
918 }
919
920 /* disable ipd/opd */
921 for (bin = 0; bin < 17; bin++) {
922 ps->aaIpdIndex[env][bin] = 0;
923 ps->aaOpdIndex[env][bin] = 0;
924 }
925 }
926#else
927 if (ps->use34hybrid_bands) {
928 for (env = 0; env < ps->num_env; env++) {
929 if (ps->iid_mode != 2 && ps->iid_mode != 5) {
930 map20indexto34(ps->iid_index[env], 34);
931 }
932 if (ps->icc_mode != 2 && ps->icc_mode != 5) {
933 map20indexto34(ps->icc_index[env], 34);
934 }
935 if (ps->ipd_mode != 2 && ps->ipd_mode != 5) {
936 map20indexto34(ps->ipd_index[env], 17);
937 map20indexto34(ps->opd_index[env], 17);
938 }
939 }
940 }
941#endif
942
943#if 0
944 for (env = 0; env < ps->num_env; env++) {
945 printf("iid[env:%d]:", env);
946 for (bin = 0; bin < 34; bin++) {
947 printf(" %d", ps->iid_index[env][bin]);
948 }
949 printf("\n");
950 }
951 for (env = 0; env < ps->num_env; env++) {
952 printf("icc[env:%d]:", env);
953 for (bin = 0; bin < 34; bin++) {
954 printf(" %d", ps->icc_index[env][bin]);
955 }
956 printf("\n");
957 }
958 for (env = 0; env < ps->num_env; env++) {
959 printf("ipd[env:%d]:", env);
960 for (bin = 0; bin < 17; bin++) {
961 printf(" %d", ps->ipd_index[env][bin]);
962 }
963 printf("\n");
964 }
965 for (env = 0; env < ps->num_env; env++) {
966 printf("opd[env:%d]:", env);
967 for (bin = 0; bin < 17; bin++) {
968 printf(" %d", ps->opd_index[env][bin]);
969 }
970 printf("\n");
971 }
972 printf("\n");
973#endif
974}
975
976/* decorrelate the mono signal using an allpass filter */
977static void ps_decorrelate(ps_info *ps, qmf_t X_left[38][64], qmf_t X_right[38][64],
978 qmf_t X_hybrid_left[32][32], qmf_t X_hybrid_right[32][32])
979{
980 uint8_t gr, n, m, bk;
981 uint8_t temp_delay;
982 uint8_t sb, maxsb;
983 const complex_t *Phi_Fract_SubQmf;
984 uint8_t temp_delay_ser[NO_ALLPASS_LINKS];
985 real_t P_SmoothPeakDecayDiffNrg, nrg;
986 real_t P[32][34];
987 real_t G_TransientRatio[32][34] = {{0}};
988 complex_t inputLeft;
989
990
991 /* chose hybrid filterbank: 20 or 34 band case */
992 if (ps->use34hybrid_bands) {
993 Phi_Fract_SubQmf = Phi_Fract_SubQmf34;
994 } else {
995 Phi_Fract_SubQmf = Phi_Fract_SubQmf20;
996 }
997
998 /* clear the energy values */
999 for (n = 0; n < 32; n++) {
1000 for (bk = 0; bk < 34; bk++) {
1001 P[n][bk] = 0;
1002 }
1003 }
1004
1005 /* calculate the energy in each parameter band b(k) */
1006 for (gr = 0; gr < ps->num_groups; gr++) {
1007 /* select the parameter index b(k) to which this group belongs */
1008 bk = (~NEGATE_IPD_MASK) & ps->map_group2bk[gr];
1009
1010 /* select the upper subband border for this group */
1011 maxsb = (gr < ps->num_hybrid_groups) ? ps->group_border[gr] + 1 : ps->group_border[gr + 1];
1012
1013 for (sb = ps->group_border[gr]; sb < maxsb; sb++) {
1014 for (n = ps->border_position[0]; n < ps->border_position[ps->num_env]; n++) {
1015#ifdef FIXED_POINT
1016 uint32_t in_re, in_im;
1017#endif
1018
1019 /* input from hybrid subbands or QMF subbands */
1020 if (gr < ps->num_hybrid_groups) {
1021 RE(inputLeft) = QMF_RE(X_hybrid_left[n][sb]);
1022 IM(inputLeft) = QMF_IM(X_hybrid_left[n][sb]);
1023 } else {
1024 RE(inputLeft) = QMF_RE(X_left[n][sb]);
1025 IM(inputLeft) = QMF_IM(X_left[n][sb]);
1026 }
1027
1028 /* accumulate energy */
1029#ifdef FIXED_POINT
1030 /* NOTE: all input is scaled by 2^(-5) because of fixed point QMF
1031 * meaning that P will be scaled by 2^(-10) compared to floating point version
1032 */
1033 in_re = ((abs(RE(inputLeft)) + (1 << (REAL_BITS - 1))) >> REAL_BITS);
1034 in_im = ((abs(IM(inputLeft)) + (1 << (REAL_BITS - 1))) >> REAL_BITS);
1035 P[n][bk] += in_re * in_re + in_im * in_im;
1036#else
1037 P[n][bk] += MUL_R(RE(inputLeft), RE(inputLeft)) + MUL_R(IM(inputLeft), IM(inputLeft));
1038#endif
1039 }
1040 }
1041 }
1042
1043#if 0
1044 for (n = 0; n < 32; n++) {
1045 for (bk = 0; bk < 34; bk++) {
1046#ifdef FIXED_POINT
1047 printf("%d %d: %d\n", n, bk, P[n][bk] /*/(float)REAL_PRECISION*/);
1048#else
1049 printf("%d %d: %f\n", n, bk, P[n][bk] / 1024.0);
1050#endif
1051 }
1052 }
1053#endif
1054
1055 /* calculate transient reduction ratio for each parameter band b(k) */
1056 for (bk = 0; bk < ps->nr_par_bands; bk++) {
1057 for (n = ps->border_position[0]; n < ps->border_position[ps->num_env]; n++) {
1058 const real_t gamma = COEF_CONST(1.5);
1059
1060 ps->P_PeakDecayNrg[bk] = MUL_F(ps->P_PeakDecayNrg[bk], ps->alpha_decay);
1061 if (ps->P_PeakDecayNrg[bk] < P[n][bk]) {
1062 ps->P_PeakDecayNrg[bk] = P[n][bk];
1063 }
1064
1065 /* apply smoothing filter to peak decay energy */
1066 P_SmoothPeakDecayDiffNrg = ps->P_SmoothPeakDecayDiffNrg_prev[bk];
1067 P_SmoothPeakDecayDiffNrg += MUL_F((ps->P_PeakDecayNrg[bk] - P[n][bk] - ps->P_SmoothPeakDecayDiffNrg_prev[bk]), ps->alpha_smooth);
1068 ps->P_SmoothPeakDecayDiffNrg_prev[bk] = P_SmoothPeakDecayDiffNrg;
1069
1070 /* apply smoothing filter to energy */
1071 nrg = ps->P_prev[bk];
1072 nrg += MUL_F((P[n][bk] - ps->P_prev[bk]), ps->alpha_smooth);
1073 ps->P_prev[bk] = nrg;
1074
1075 /* calculate transient ratio */
1076 if (MUL_C(P_SmoothPeakDecayDiffNrg, gamma) <= nrg) {
1077 G_TransientRatio[n][bk] = REAL_CONST(1.0);
1078 } else {
1079 G_TransientRatio[n][bk] = DIV_R(nrg, (MUL_C(P_SmoothPeakDecayDiffNrg, gamma)));
1080 }
1081 }
1082 }
1083
1084#if 0
1085 for (n = 0; n < 32; n++) {
1086 for (bk = 0; bk < 34; bk++) {
1087#ifdef FIXED_POINT
1088 printf("%d %d: %f\n", n, bk, G_TransientRatio[n][bk] / (float)REAL_PRECISION);
1089#else
1090 printf("%d %d: %f\n", n, bk, G_TransientRatio[n][bk]);
1091#endif
1092 }
1093 }
1094#endif
1095
1096 /* apply stereo decorrelation filter to the signal */
1097 for (gr = 0; gr < ps->num_groups; gr++) {
1098 if (gr < ps->num_hybrid_groups) {
1099 maxsb = ps->group_border[gr] + 1;
1100 } else {
1101 maxsb = ps->group_border[gr + 1];
1102 }
1103
1104 /* QMF channel */
1105 for (sb = ps->group_border[gr]; sb < maxsb; sb++) {
1106 real_t g_DecaySlope;
1107 real_t g_DecaySlope_filt[NO_ALLPASS_LINKS];
1108
1109 /* g_DecaySlope: [0..1] */
1110 if (gr < ps->num_hybrid_groups || sb <= ps->decay_cutoff) {
1111 g_DecaySlope = FRAC_CONST(1.0);
1112 } else {
1113 int8_t decay = ps->decay_cutoff - sb;
1114 if (decay <= -20 /* -1/DECAY_SLOPE */) {
1115 g_DecaySlope = 0;
1116 } else {
1117 /* decay(int)*decay_slope(frac) = g_DecaySlope(frac) */
1118 g_DecaySlope = FRAC_CONST(1.0) + DECAY_SLOPE * decay;
1119 }
1120 }
1121
1122 /* calculate g_DecaySlope_filt for every m multiplied by filter_a[m] */
1123 for (m = 0; m < NO_ALLPASS_LINKS; m++) {
1124 g_DecaySlope_filt[m] = MUL_F(g_DecaySlope, filter_a[m]);
1125 }
1126
1127
1128 /* set delay indices */
1129 temp_delay = ps->saved_delay;
1130 for (n = 0; n < NO_ALLPASS_LINKS; n++) {
1131 temp_delay_ser[n] = ps->delay_buf_index_ser[n];
1132 }
1133
1134 for (n = ps->border_position[0]; n < ps->border_position[ps->num_env]; n++) {
1135 complex_t tmp, tmp0, R0;
1136
1137 if (gr < ps->num_hybrid_groups) {
1138 /* hybrid filterbank input */
1139 RE(inputLeft) = QMF_RE(X_hybrid_left[n][sb]);
1140 IM(inputLeft) = QMF_IM(X_hybrid_left[n][sb]);
1141 } else {
1142 /* QMF filterbank input */
1143 RE(inputLeft) = QMF_RE(X_left[n][sb]);
1144 IM(inputLeft) = QMF_IM(X_left[n][sb]);
1145 }
1146
1147 if (sb > ps->nr_allpass_bands && gr >= ps->num_hybrid_groups) {
1148 /* delay */
1149
1150 /* never hybrid subbands here, always QMF subbands */
1151 RE(tmp) = RE(ps->delay_Qmf[ps->delay_buf_index_delay[sb]][sb]);
1152 IM(tmp) = IM(ps->delay_Qmf[ps->delay_buf_index_delay[sb]][sb]);
1153 RE(R0) = RE(tmp);
1154 IM(R0) = IM(tmp);
1155 RE(ps->delay_Qmf[ps->delay_buf_index_delay[sb]][sb]) = RE(inputLeft);
1156 IM(ps->delay_Qmf[ps->delay_buf_index_delay[sb]][sb]) = IM(inputLeft);
1157 } else {
1158 /* allpass filter */
1159 uint8_t m;
1160 complex_t Phi_Fract;
1161
1162 /* fetch parameters */
1163 if (gr < ps->num_hybrid_groups) {
1164 /* select data from the hybrid subbands */
1165 RE(tmp0) = RE(ps->delay_SubQmf[temp_delay][sb]);
1166 IM(tmp0) = IM(ps->delay_SubQmf[temp_delay][sb]);
1167
1168 RE(ps->delay_SubQmf[temp_delay][sb]) = RE(inputLeft);
1169 IM(ps->delay_SubQmf[temp_delay][sb]) = IM(inputLeft);
1170
1171 RE(Phi_Fract) = RE(Phi_Fract_SubQmf[sb]);
1172 IM(Phi_Fract) = IM(Phi_Fract_SubQmf[sb]);
1173 } else {
1174 /* select data from the QMF subbands */
1175 RE(tmp0) = RE(ps->delay_Qmf[temp_delay][sb]);
1176 IM(tmp0) = IM(ps->delay_Qmf[temp_delay][sb]);
1177
1178 RE(ps->delay_Qmf[temp_delay][sb]) = RE(inputLeft);
1179 IM(ps->delay_Qmf[temp_delay][sb]) = IM(inputLeft);
1180
1181 RE(Phi_Fract) = RE(Phi_Fract_Qmf[sb]);
1182 IM(Phi_Fract) = IM(Phi_Fract_Qmf[sb]);
1183 }
1184
1185 /* z^(-2) * Phi_Fract[k] */
1186 ComplexMult(&RE(tmp), &IM(tmp), RE(tmp0), IM(tmp0), RE(Phi_Fract), IM(Phi_Fract));
1187
1188 RE(R0) = RE(tmp);
1189 IM(R0) = IM(tmp);
1190 for (m = 0; m < NO_ALLPASS_LINKS; m++) {
1191 complex_t Q_Fract_allpass, tmp2;
1192
1193 /* fetch parameters */
1194 if (gr < ps->num_hybrid_groups) {
1195 /* select data from the hybrid subbands */
1196 RE(tmp0) = RE(ps->delay_SubQmf_ser[m][temp_delay_ser[m]][sb]);
1197 IM(tmp0) = IM(ps->delay_SubQmf_ser[m][temp_delay_ser[m]][sb]);
1198
1199 if (ps->use34hybrid_bands) {
1200 RE(Q_Fract_allpass) = RE(Q_Fract_allpass_SubQmf34[sb][m]);
1201 IM(Q_Fract_allpass) = IM(Q_Fract_allpass_SubQmf34[sb][m]);
1202 } else {
1203 RE(Q_Fract_allpass) = RE(Q_Fract_allpass_SubQmf20[sb][m]);
1204 IM(Q_Fract_allpass) = IM(Q_Fract_allpass_SubQmf20[sb][m]);
1205 }
1206 } else {
1207 /* select data from the QMF subbands */
1208 RE(tmp0) = RE(ps->delay_Qmf_ser[m][temp_delay_ser[m]][sb]);
1209 IM(tmp0) = IM(ps->delay_Qmf_ser[m][temp_delay_ser[m]][sb]);
1210
1211 RE(Q_Fract_allpass) = RE(Q_Fract_allpass_Qmf[sb][m]);
1212 IM(Q_Fract_allpass) = IM(Q_Fract_allpass_Qmf[sb][m]);
1213 }
1214
1215 /* delay by a fraction */
1216 /* z^(-d(m)) * Q_Fract_allpass[k,m] */
1217 ComplexMult(&RE(tmp), &IM(tmp), RE(tmp0), IM(tmp0), RE(Q_Fract_allpass), IM(Q_Fract_allpass));
1218
1219 /* -a(m) * g_DecaySlope[k] */
1220 RE(tmp) += -MUL_F(g_DecaySlope_filt[m], RE(R0));
1221 IM(tmp) += -MUL_F(g_DecaySlope_filt[m], IM(R0));
1222
1223 /* -a(m) * g_DecaySlope[k] * Q_Fract_allpass[k,m] * z^(-d(m)) */
1224 RE(tmp2) = RE(R0) + MUL_F(g_DecaySlope_filt[m], RE(tmp));
1225 IM(tmp2) = IM(R0) + MUL_F(g_DecaySlope_filt[m], IM(tmp));
1226
1227 /* store sample */
1228 if (gr < ps->num_hybrid_groups) {
1229 RE(ps->delay_SubQmf_ser[m][temp_delay_ser[m]][sb]) = RE(tmp2);
1230 IM(ps->delay_SubQmf_ser[m][temp_delay_ser[m]][sb]) = IM(tmp2);
1231 } else {
1232 RE(ps->delay_Qmf_ser[m][temp_delay_ser[m]][sb]) = RE(tmp2);
1233 IM(ps->delay_Qmf_ser[m][temp_delay_ser[m]][sb]) = IM(tmp2);
1234 }
1235
1236 /* store for next iteration (or as output value if last iteration) */
1237 RE(R0) = RE(tmp);
1238 IM(R0) = IM(tmp);
1239 }
1240 }
1241
1242 /* select b(k) for reading the transient ratio */
1243 bk = (~NEGATE_IPD_MASK) & ps->map_group2bk[gr];
1244
1245 /* duck if a past transient is found */
1246 RE(R0) = MUL_R(G_TransientRatio[n][bk], RE(R0));
1247 IM(R0) = MUL_R(G_TransientRatio[n][bk], IM(R0));
1248
1249 if (gr < ps->num_hybrid_groups) {
1250 /* hybrid */
1251 QMF_RE(X_hybrid_right[n][sb]) = RE(R0);
1252 QMF_IM(X_hybrid_right[n][sb]) = IM(R0);
1253 } else {
1254 /* QMF */
1255 QMF_RE(X_right[n][sb]) = RE(R0);
1256 QMF_IM(X_right[n][sb]) = IM(R0);
1257 }
1258
1259 /* Update delay buffer index */
1260 if (++temp_delay >= 2) {
1261 temp_delay = 0;
1262 }
1263
1264 /* update delay indices */
1265 if (sb > ps->nr_allpass_bands && gr >= ps->num_hybrid_groups) {
1266 /* delay_D depends on the samplerate, it can hold the values 14 and 1 */
1267 if (++ps->delay_buf_index_delay[sb] >= ps->delay_D[sb]) {
1268 ps->delay_buf_index_delay[sb] = 0;
1269 }
1270 }
1271
1272 for (m = 0; m < NO_ALLPASS_LINKS; m++) {
1273 if (++temp_delay_ser[m] >= ps->num_sample_delay_ser[m]) {
1274 temp_delay_ser[m] = 0;
1275 }
1276 }
1277 }
1278 }
1279 }
1280
1281 /* update delay indices */
1282 ps->saved_delay = temp_delay;
1283 for (m = 0; m < NO_ALLPASS_LINKS; m++) {
1284 ps->delay_buf_index_ser[m] = temp_delay_ser[m];
1285 }
1286}
1287
1288#ifdef FIXED_POINT
1289#define step(shift) \
1290 if ((0x40000000l >> shift) + root <= value) \
1291 { \
1292 value -= (0x40000000l >> shift) + root; \
1293 root = (root >> 1) | (0x40000000l >> shift); \
1294 } else { \
1295 root = root >> 1; \
1296 }
1297
1298/* fixed point square root approximation */
1299static real_t ps_sqrt(real_t value)
1300{
1301 real_t root = 0;
1302
1303 step(0);
1304 step(2);
1305 step(4);
1306 step(6);
1307 step(8);
1308 step(10);
1309 step(12);
1310 step(14);
1311 step(16);
1312 step(18);
1313 step(20);
1314 step(22);
1315 step(24);
1316 step(26);
1317 step(28);
1318 step(30);
1319
1320 if (root < value) {
1321 ++root;
1322 }
1323
1324 root <<= (REAL_BITS / 2);
1325
1326 return root;
1327}
1328#else
1329#define ps_sqrt(A) sqrt(A)
1330#endif
1331
1332static const real_t ipdopd_cos_tab[] = {
1333 FRAC_CONST(1.000000000000000),
1334 FRAC_CONST(0.707106781186548),
1335 FRAC_CONST(0.000000000000000),
1336 FRAC_CONST(-0.707106781186547),
1337 FRAC_CONST(-1.000000000000000),
1338 FRAC_CONST(-0.707106781186548),
1339 FRAC_CONST(-0.000000000000000),
1340 FRAC_CONST(0.707106781186547),
1341 FRAC_CONST(1.000000000000000)
1342};
1343
1344static const real_t ipdopd_sin_tab[] = {
1345 FRAC_CONST(0.000000000000000),
1346 FRAC_CONST(0.707106781186547),
1347 FRAC_CONST(1.000000000000000),
1348 FRAC_CONST(0.707106781186548),
1349 FRAC_CONST(0.000000000000000),
1350 FRAC_CONST(-0.707106781186547),
1351 FRAC_CONST(-1.000000000000000),
1352 FRAC_CONST(-0.707106781186548),
1353 FRAC_CONST(-0.000000000000000)
1354};
1355
1356static real_t magnitude_c(complex_t c)
1357{
1358#ifdef FIXED_POINT
1359#define ps_abs(A) (((A) > 0) ? (A) : (-(A)))
1360#define ALPHA FRAC_CONST(0.948059448969)
1361#define BETA FRAC_CONST(0.392699081699)
1362
1363 real_t abs_inphase = ps_abs(RE(c));
1364 real_t abs_quadrature = ps_abs(IM(c));
1365
1366 if (abs_inphase > abs_quadrature) {
1367 return MUL_F(abs_inphase, ALPHA) + MUL_F(abs_quadrature, BETA);
1368 } else {
1369 return MUL_F(abs_quadrature, ALPHA) + MUL_F(abs_inphase, BETA);
1370 }
1371#else
1372 return sqrt(RE(c) * RE(c) + IM(c) * IM(c));
1373#endif
1374}
1375
1376static void ps_mix_phase(ps_info *ps, qmf_t X_left[38][64], qmf_t X_right[38][64],
1377 qmf_t X_hybrid_left[32][32], qmf_t X_hybrid_right[32][32])
1378{
1379 uint8_t n;
1380 uint8_t gr;
1381 uint8_t bk = 0;
1382 uint8_t sb, maxsb;
1383 uint8_t env;
1384 uint8_t nr_ipdopd_par;
1385 complex_t h11, h12, h21, h22;
1386 complex_t H11, H12, H21, H22;
1387 complex_t deltaH11, deltaH12, deltaH21, deltaH22;
1388 complex_t tempLeft;
1389 complex_t tempRight;
1390 complex_t phaseLeft;
1391 complex_t phaseRight;
1392 real_t L;
1393 const real_t *sf_iid;
1394 uint8_t no_iid_steps;
1395
1396 if (ps->iid_mode >= 3) {
1397 no_iid_steps = 15;
1398 sf_iid = sf_iid_fine;
1399 } else {
1400 no_iid_steps = 7;
1401 sf_iid = sf_iid_normal;
1402 }
1403
1404 if (ps->ipd_mode == 0 || ps->ipd_mode == 3) {
1405 nr_ipdopd_par = 11; /* resolution */
1406 } else {
1407 nr_ipdopd_par = ps->nr_ipdopd_par;
1408 }
1409
1410 for (gr = 0; gr < ps->num_groups; gr++) {
1411 bk = (~NEGATE_IPD_MASK) & ps->map_group2bk[gr];
1412
1413 /* use one channel per group in the subqmf domain */
1414 maxsb = (gr < ps->num_hybrid_groups) ? ps->group_border[gr] + 1 : ps->group_border[gr + 1];
1415
1416 for (env = 0; env < ps->num_env; env++) {
1417 if (ps->icc_mode < 3) {
1418 /* type 'A' mixing as described in 8.6.4.6.2.1 */
1419 real_t c_1, c_2;
1420 real_t cosa, sina;
1421 real_t cosb, sinb;
1422 real_t ab1, ab2;
1423 real_t ab3, ab4;
1424
1425 /*
1426 c_1 = sqrt(2.0 / (1.0 + pow(10.0, quant_iid[no_iid_steps + iid_index] / 10.0)));
1427 c_2 = sqrt(2.0 / (1.0 + pow(10.0, quant_iid[no_iid_steps - iid_index] / 10.0)));
1428 alpha = 0.5 * acos(quant_rho[icc_index]);
1429 beta = alpha * ( c_1 - c_2 ) / sqrt(2.0);
1430 */
1431
1432 //printf("%d\n", ps->iid_index[env][bk]);
1433
1434 /* calculate the scalefactors c_1 and c_2 from the intensity differences */
1435 c_1 = sf_iid[no_iid_steps + ps->iid_index[env][bk]];
1436 c_2 = sf_iid[no_iid_steps - ps->iid_index[env][bk]];
1437
1438 /* calculate alpha and beta using the ICC parameters */
1439 cosa = cos_alphas[ps->icc_index[env][bk]];
1440 sina = sin_alphas[ps->icc_index[env][bk]];
1441
1442 if (ps->iid_mode >= 3) {
1443 if (ps->iid_index[env][bk] < 0) {
1444 cosb = cos_betas_fine[-ps->iid_index[env][bk]][ps->icc_index[env][bk]];
1445 sinb = -sin_betas_fine[-ps->iid_index[env][bk]][ps->icc_index[env][bk]];
1446 } else {
1447 cosb = cos_betas_fine[ps->iid_index[env][bk]][ps->icc_index[env][bk]];
1448 sinb = sin_betas_fine[ps->iid_index[env][bk]][ps->icc_index[env][bk]];
1449 }
1450 } else {
1451 if (ps->iid_index[env][bk] < 0) {
1452 cosb = cos_betas_normal[-ps->iid_index[env][bk]][ps->icc_index[env][bk]];
1453 sinb = -sin_betas_normal[-ps->iid_index[env][bk]][ps->icc_index[env][bk]];
1454 } else {
1455 cosb = cos_betas_normal[ps->iid_index[env][bk]][ps->icc_index[env][bk]];
1456 sinb = sin_betas_normal[ps->iid_index[env][bk]][ps->icc_index[env][bk]];
1457 }
1458 }
1459
1460 ab1 = MUL_C(cosb, cosa);
1461 ab2 = MUL_C(sinb, sina);
1462 ab3 = MUL_C(sinb, cosa);
1463 ab4 = MUL_C(cosb, sina);
1464
1465 /* h_xy: COEF */
1466 RE(h11) = MUL_C(c_2, (ab1 - ab2));
1467 RE(h12) = MUL_C(c_1, (ab1 + ab2));
1468 RE(h21) = MUL_C(c_2, (ab3 + ab4));
1469 RE(h22) = MUL_C(c_1, (ab3 - ab4));
1470 } else {
1471 /* type 'B' mixing as described in 8.6.4.6.2.2 */
1472 real_t sina, cosa;
1473 real_t cosg, sing;
1474
1475 /*
1476 real_t c, rho, mu, alpha, gamma;
1477 uint8_t i;
1478
1479 i = ps->iid_index[env][bk];
1480 c = (real_t)pow(10.0, ((i)?(((i>0)?1:-1)*quant_iid[((i>0)?i:-i)-1]):0.)/20.0);
1481 rho = quant_rho[ps->icc_index[env][bk]];
1482
1483 if (rho == 0.0f && c == 1.)
1484 {
1485 alpha = (real_t)M_PI/4.0f;
1486 rho = 0.05f;
1487 } else {
1488 if (rho <= 0.05f)
1489 {
1490 rho = 0.05f;
1491 }
1492 alpha = 0.5f*(real_t)atan( (2.0f*c*rho) / (c*c-1.0f) );
1493
1494 if (alpha < 0.)
1495 {
1496 alpha += (real_t)M_PI/2.0f;
1497 }
1498 if (rho < 0.)
1499 {
1500 alpha += (real_t)M_PI;
1501 }
1502 }
1503 mu = c+1.0f/c;
1504 mu = 1+(4.0f*rho*rho-4.0f)/(mu*mu);
1505 gamma = (real_t)atan(sqrt((1.0f-sqrt(mu))/(1.0f+sqrt(mu))));
1506 */
1507
1508 if (ps->iid_mode >= 3) {
1509 uint8_t abs_iid = abs(ps->iid_index[env][bk]);
1510
1511 cosa = sincos_alphas_B_fine[no_iid_steps + ps->iid_index[env][bk]][ps->icc_index[env][bk]];
1512 sina = sincos_alphas_B_fine[30 - (no_iid_steps + ps->iid_index[env][bk])][ps->icc_index[env][bk]];
1513 cosg = cos_gammas_fine[abs_iid][ps->icc_index[env][bk]];
1514 sing = sin_gammas_fine[abs_iid][ps->icc_index[env][bk]];
1515 } else {
1516 uint8_t abs_iid = abs(ps->iid_index[env][bk]);
1517
1518 cosa = sincos_alphas_B_normal[no_iid_steps + ps->iid_index[env][bk]][ps->icc_index[env][bk]];
1519 sina = sincos_alphas_B_normal[14 - (no_iid_steps + ps->iid_index[env][bk])][ps->icc_index[env][bk]];
1520 cosg = cos_gammas_normal[abs_iid][ps->icc_index[env][bk]];
1521 sing = sin_gammas_normal[abs_iid][ps->icc_index[env][bk]];
1522 }
1523
1524 RE(h11) = MUL_C(COEF_SQRT2, MUL_C(cosa, cosg));
1525 RE(h12) = MUL_C(COEF_SQRT2, MUL_C(sina, cosg));
1526 RE(h21) = MUL_C(COEF_SQRT2, MUL_C(-cosa, sing));
1527 RE(h22) = MUL_C(COEF_SQRT2, MUL_C(sina, sing));
1528 }
1529
1530 /* calculate phase rotation parameters H_xy */
1531 /* note that the imaginary part of these parameters are only calculated when
1532 IPD and OPD are enabled
1533 */
1534 if ((ps->enable_ipdopd) && (bk < nr_ipdopd_par)) {
1535 int8_t i;
1536 real_t xy, pq, xypq;
1537
1538 /* ringbuffer index */
1539 i = ps->phase_hist;
1540
1541 /* previous value */
1542#ifdef FIXED_POINT
1543 /* divide by 4, shift right 2 bits */
1544 RE(tempLeft) = RE(ps->ipd_prev[bk][i]) >> 2;
1545 IM(tempLeft) = IM(ps->ipd_prev[bk][i]) >> 2;
1546 RE(tempRight) = RE(ps->opd_prev[bk][i]) >> 2;
1547 IM(tempRight) = IM(ps->opd_prev[bk][i]) >> 2;
1548#else
1549 RE(tempLeft) = MUL_F(RE(ps->ipd_prev[bk][i]), FRAC_CONST(0.25));
1550 IM(tempLeft) = MUL_F(IM(ps->ipd_prev[bk][i]), FRAC_CONST(0.25));
1551 RE(tempRight) = MUL_F(RE(ps->opd_prev[bk][i]), FRAC_CONST(0.25));
1552 IM(tempRight) = MUL_F(IM(ps->opd_prev[bk][i]), FRAC_CONST(0.25));
1553#endif
1554
1555 /* save current value */
1556 RE(ps->ipd_prev[bk][i]) = ipdopd_cos_tab[abs(ps->ipd_index[env][bk])];
1557 IM(ps->ipd_prev[bk][i]) = ipdopd_sin_tab[abs(ps->ipd_index[env][bk])];
1558 RE(ps->opd_prev[bk][i]) = ipdopd_cos_tab[abs(ps->opd_index[env][bk])];
1559 IM(ps->opd_prev[bk][i]) = ipdopd_sin_tab[abs(ps->opd_index[env][bk])];
1560
1561 /* add current value */
1562 RE(tempLeft) += RE(ps->ipd_prev[bk][i]);
1563 IM(tempLeft) += IM(ps->ipd_prev[bk][i]);
1564 RE(tempRight) += RE(ps->opd_prev[bk][i]);
1565 IM(tempRight) += IM(ps->opd_prev[bk][i]);
1566
1567 /* ringbuffer index */
1568 if (i == 0) {
1569 i = 2;
1570 }
1571 i--;
1572
1573 /* get value before previous */
1574#ifdef FIXED_POINT
1575 /* dividing by 2, shift right 1 bit */
1576 RE(tempLeft) += (RE(ps->ipd_prev[bk][i]) >> 1);
1577 IM(tempLeft) += (IM(ps->ipd_prev[bk][i]) >> 1);
1578 RE(tempRight) += (RE(ps->opd_prev[bk][i]) >> 1);
1579 IM(tempRight) += (IM(ps->opd_prev[bk][i]) >> 1);
1580#else
1581 RE(tempLeft) += MUL_F(RE(ps->ipd_prev[bk][i]), FRAC_CONST(0.5));
1582 IM(tempLeft) += MUL_F(IM(ps->ipd_prev[bk][i]), FRAC_CONST(0.5));
1583 RE(tempRight) += MUL_F(RE(ps->opd_prev[bk][i]), FRAC_CONST(0.5));
1584 IM(tempRight) += MUL_F(IM(ps->opd_prev[bk][i]), FRAC_CONST(0.5));
1585#endif
1586
1587#if 0 /* original code */
1588 ipd = (float)atan2(IM(tempLeft), RE(tempLeft));
1589 opd = (float)atan2(IM(tempRight), RE(tempRight));
1590
1591 /* phase rotation */
1592 RE(phaseLeft) = (float)cos(opd);
1593 IM(phaseLeft) = (float)sin(opd);
1594 opd -= ipd;
1595 RE(phaseRight) = (float)cos(opd);
1596 IM(phaseRight) = (float)sin(opd);
1597#else
1598
1599 // x = IM(tempLeft)
1600 // y = RE(tempLeft)
1601 // p = IM(tempRight)
1602 // q = RE(tempRight)
1603 // cos(atan2(x,y)) = y/sqrt((x*x) + (y*y))
1604 // sin(atan2(x,y)) = x/sqrt((x*x) + (y*y))
1605 // cos(atan2(x,y)-atan2(p,q)) = (y*q + x*p) / ( sqrt((x*x) + (y*y)) * sqrt((p*p) + (q*q)) );
1606 // sin(atan2(x,y)-atan2(p,q)) = (x*q - y*p) / ( sqrt((x*x) + (y*y)) * sqrt((p*p) + (q*q)) );
1607
1608 xy = magnitude_c(tempRight);
1609 pq = magnitude_c(tempLeft);
1610
1611 if (xy != 0) {
1612 RE(phaseLeft) = DIV_R(RE(tempRight), xy);
1613 IM(phaseLeft) = DIV_R(IM(tempRight), xy);
1614 } else {
1615 RE(phaseLeft) = 0;
1616 IM(phaseLeft) = 0;
1617 }
1618
1619 xypq = MUL_R(xy, pq);
1620
1621 if (xypq != 0) {
1622 real_t tmp1 = MUL_R(RE(tempRight), RE(tempLeft)) + MUL_R(IM(tempRight), IM(tempLeft));
1623 real_t tmp2 = MUL_R(IM(tempRight), RE(tempLeft)) - MUL_R(RE(tempRight), IM(tempLeft));
1624
1625 RE(phaseRight) = DIV_R(tmp1, xypq);
1626 IM(phaseRight) = DIV_R(tmp2, xypq);
1627 } else {
1628 RE(phaseRight) = 0;
1629 IM(phaseRight) = 0;
1630 }
1631
1632#endif
1633
1634 /* MUL_F(COEF, REAL) = COEF */
1635 IM(h11) = MUL_R(RE(h11), IM(phaseLeft));
1636 IM(h12) = MUL_R(RE(h12), IM(phaseRight));
1637 IM(h21) = MUL_R(RE(h21), IM(phaseLeft));
1638 IM(h22) = MUL_R(RE(h22), IM(phaseRight));
1639
1640 RE(h11) = MUL_R(RE(h11), RE(phaseLeft));
1641 RE(h12) = MUL_R(RE(h12), RE(phaseRight));
1642 RE(h21) = MUL_R(RE(h21), RE(phaseLeft));
1643 RE(h22) = MUL_R(RE(h22), RE(phaseRight));
1644 }
1645
1646 /* length of the envelope n_e+1 - n_e (in time samples) */
1647 /* 0 < L <= 32: integer */
1648 L = (real_t)(ps->border_position[env + 1] - ps->border_position[env]);
1649
1650 /* obtain final H_xy by means of linear interpolation */
1651 RE(deltaH11) = (RE(h11) - RE(ps->h11_prev[gr])) / L;
1652 RE(deltaH12) = (RE(h12) - RE(ps->h12_prev[gr])) / L;
1653 RE(deltaH21) = (RE(h21) - RE(ps->h21_prev[gr])) / L;
1654 RE(deltaH22) = (RE(h22) - RE(ps->h22_prev[gr])) / L;
1655
1656 RE(H11) = RE(ps->h11_prev[gr]);
1657 RE(H12) = RE(ps->h12_prev[gr]);
1658 RE(H21) = RE(ps->h21_prev[gr]);
1659 RE(H22) = RE(ps->h22_prev[gr]);
1660
1661 RE(ps->h11_prev[gr]) = RE(h11);
1662 RE(ps->h12_prev[gr]) = RE(h12);
1663 RE(ps->h21_prev[gr]) = RE(h21);
1664 RE(ps->h22_prev[gr]) = RE(h22);
1665
1666 /* only calculate imaginary part when needed */
1667 if ((ps->enable_ipdopd) && (bk < nr_ipdopd_par)) {
1668 /* obtain final H_xy by means of linear interpolation */
1669 IM(deltaH11) = (IM(h11) - IM(ps->h11_prev[gr])) / L;
1670 IM(deltaH12) = (IM(h12) - IM(ps->h12_prev[gr])) / L;
1671 IM(deltaH21) = (IM(h21) - IM(ps->h21_prev[gr])) / L;
1672 IM(deltaH22) = (IM(h22) - IM(ps->h22_prev[gr])) / L;
1673
1674 IM(H11) = IM(ps->h11_prev[gr]);
1675 IM(H12) = IM(ps->h12_prev[gr]);
1676 IM(H21) = IM(ps->h21_prev[gr]);
1677 IM(H22) = IM(ps->h22_prev[gr]);
1678
1679 if ((NEGATE_IPD_MASK & ps->map_group2bk[gr]) != 0) {
1680 IM(deltaH11) = -IM(deltaH11);
1681 IM(deltaH12) = -IM(deltaH12);
1682 IM(deltaH21) = -IM(deltaH21);
1683 IM(deltaH22) = -IM(deltaH22);
1684
1685 IM(H11) = -IM(H11);
1686 IM(H12) = -IM(H12);
1687 IM(H21) = -IM(H21);
1688 IM(H22) = -IM(H22);
1689 }
1690
1691 IM(ps->h11_prev[gr]) = IM(h11);
1692 IM(ps->h12_prev[gr]) = IM(h12);
1693 IM(ps->h21_prev[gr]) = IM(h21);
1694 IM(ps->h22_prev[gr]) = IM(h22);
1695 }
1696
1697 /* apply H_xy to the current envelope band of the decorrelated subband */
1698 for (n = ps->border_position[env]; n < ps->border_position[env + 1]; n++) {
1699 /* addition finalises the interpolation over every n */
1700 RE(H11) += RE(deltaH11);
1701 RE(H12) += RE(deltaH12);
1702 RE(H21) += RE(deltaH21);
1703 RE(H22) += RE(deltaH22);
1704 if ((ps->enable_ipdopd) && (bk < nr_ipdopd_par)) {
1705 IM(H11) += IM(deltaH11);
1706 IM(H12) += IM(deltaH12);
1707 IM(H21) += IM(deltaH21);
1708 IM(H22) += IM(deltaH22);
1709 }
1710
1711 /* channel is an alias to the subband */
1712 for (sb = ps->group_border[gr]; sb < maxsb; sb++) {
1713 complex_t inLeft, inRight;
1714
1715 /* load decorrelated samples */
1716 if (gr < ps->num_hybrid_groups) {
1717 RE(inLeft) = RE(X_hybrid_left[n][sb]);
1718 IM(inLeft) = IM(X_hybrid_left[n][sb]);
1719 RE(inRight) = RE(X_hybrid_right[n][sb]);
1720 IM(inRight) = IM(X_hybrid_right[n][sb]);
1721 } else {
1722 RE(inLeft) = RE(X_left[n][sb]);
1723 IM(inLeft) = IM(X_left[n][sb]);
1724 RE(inRight) = RE(X_right[n][sb]);
1725 IM(inRight) = IM(X_right[n][sb]);
1726 }
1727
1728 /* apply mixing */
1729 RE(tempLeft) = MUL_C(RE(H11), RE(inLeft)) + MUL_C(RE(H21), RE(inRight));
1730 IM(tempLeft) = MUL_C(RE(H11), IM(inLeft)) + MUL_C(RE(H21), IM(inRight));
1731 RE(tempRight) = MUL_C(RE(H12), RE(inLeft)) + MUL_C(RE(H22), RE(inRight));
1732 IM(tempRight) = MUL_C(RE(H12), IM(inLeft)) + MUL_C(RE(H22), IM(inRight));
1733
1734 /* only perform imaginary operations when needed */
1735 if ((ps->enable_ipdopd) && (bk < nr_ipdopd_par)) {
1736 /* apply rotation */
1737 RE(tempLeft) -= MUL_C(IM(H11), IM(inLeft)) + MUL_C(IM(H21), IM(inRight));
1738 IM(tempLeft) += MUL_C(IM(H11), RE(inLeft)) + MUL_C(IM(H21), RE(inRight));
1739 RE(tempRight) -= MUL_C(IM(H12), IM(inLeft)) + MUL_C(IM(H22), IM(inRight));
1740 IM(tempRight) += MUL_C(IM(H12), RE(inLeft)) + MUL_C(IM(H22), RE(inRight));
1741 }
1742
1743 /* store final samples */
1744 if (gr < ps->num_hybrid_groups) {
1745 RE(X_hybrid_left[n][sb]) = RE(tempLeft);
1746 IM(X_hybrid_left[n][sb]) = IM(tempLeft);
1747 RE(X_hybrid_right[n][sb]) = RE(tempRight);
1748 IM(X_hybrid_right[n][sb]) = IM(tempRight);
1749 } else {
1750 RE(X_left[n][sb]) = RE(tempLeft);
1751 IM(X_left[n][sb]) = IM(tempLeft);
1752 RE(X_right[n][sb]) = RE(tempRight);
1753 IM(X_right[n][sb]) = IM(tempRight);
1754 }
1755 }
1756 }
1757
1758 /* shift phase smoother's circular buffer index */
1759 ps->phase_hist++;
1760 if (ps->phase_hist == 2) {
1761 ps->phase_hist = 0;
1762 }
1763 }
1764 }
1765}
1766
1767void ps_free(ps_info *ps)
1768{
1769 /* free hybrid filterbank structures */
1770 hybrid_free(ps->hyb);
yuliang.hu585b7202024-09-06 17:28:13 +0800[diff] [blame^]1771 if (ps->process_buf) {
1772 free(ps->process_buf);
1773 ps->process_buf = NULL;
1774 }
1775 if (ps->process_buf2) {
1776 free(ps->process_buf2);
1777 ps->process_buf2 = NULL;
1778 }
xingri.gaoc18d4472023-02-28 02:51:02 +0000[diff] [blame]1779 faad_free(ps);
1780}
1781
1782ps_info *ps_init(uint8_t sr_index __unused, uint8_t numTimeSlotsRate)
1783{
1784 uint8_t i;
1785 uint8_t short_delay_band;
1786
1787 ps_info *ps = (ps_info*)faad_malloc(sizeof(ps_info));
1788 memset(ps, 0, sizeof(ps_info));
1789
1790 ps->hyb = hybrid_init(numTimeSlotsRate);
1791 ps->numTimeSlotsRate = numTimeSlotsRate;
1792
1793 ps->ps_data_available = 0;
1794
1795 /* delay stuff*/
1796 ps->saved_delay = 0;
1797
1798 for (i = 0; i < 64; i++) {
1799 ps->delay_buf_index_delay[i] = 0;
1800 }
1801
1802 for (i = 0; i < NO_ALLPASS_LINKS; i++) {
1803 ps->delay_buf_index_ser[i] = 0;
1804#ifdef PARAM_32KHZ
1805 if (sr_index <= 5) { /* >= 32 kHz*/
1806 ps->num_sample_delay_ser[i] = delay_length_d[1][i];
1807 } else {
1808 ps->num_sample_delay_ser[i] = delay_length_d[0][i];
1809 }
1810#else
1811 /* THESE ARE CONSTANTS NOW */
1812 ps->num_sample_delay_ser[i] = delay_length_d[i];
1813#endif
1814 }
1815
1816#ifdef PARAM_32KHZ
1817 if (sr_index <= 5) { /* >= 32 kHz*/
1818 short_delay_band = 35;
1819 ps->nr_allpass_bands = 22;
1820 ps->alpha_decay = FRAC_CONST(0.76592833836465);
1821 ps->alpha_smooth = FRAC_CONST(0.25);
1822 } else {
1823 short_delay_band = 64;
1824 ps->nr_allpass_bands = 45;
1825 ps->alpha_decay = FRAC_CONST(0.58664621951003);
1826 ps->alpha_smooth = FRAC_CONST(0.6);
1827 }
1828#else
1829 /* THESE ARE CONSTANTS NOW */
1830 short_delay_band = 35;
1831 ps->nr_allpass_bands = 22;
1832 ps->alpha_decay = FRAC_CONST(0.76592833836465);
1833 ps->alpha_smooth = FRAC_CONST(0.25);
1834#endif
1835
1836 /* THESE ARE CONSTANT NOW IF PS IS INDEPENDANT OF SAMPLERATE */
1837 for (i = 0; i < short_delay_band; i++) {
1838 ps->delay_D[i] = 14;
1839 }
1840 for (i = short_delay_band; i < 64; i++) {
1841 ps->delay_D[i] = 1;
1842 }
1843
1844 /* mixing and phase */
1845 for (i = 0; i < 50; i++) {
1846 RE(ps->h11_prev[i]) = 1;
1847 IM(ps->h12_prev[i]) = 1;
1848 RE(ps->h11_prev[i]) = 1;
1849 IM(ps->h12_prev[i]) = 1;
1850 }
1851
1852 ps->phase_hist = 0;
1853
1854 for (i = 0; i < 20; i++) {
1855 RE(ps->ipd_prev[i][0]) = 0;
1856 IM(ps->ipd_prev[i][0]) = 0;
1857 RE(ps->ipd_prev[i][1]) = 0;
1858 IM(ps->ipd_prev[i][1]) = 0;
1859 RE(ps->opd_prev[i][0]) = 0;
1860 IM(ps->opd_prev[i][0]) = 0;
1861 RE(ps->opd_prev[i][1]) = 0;
1862 IM(ps->opd_prev[i][1]) = 0;
1863 }
1864
yuliang.hu585b7202024-09-06 17:28:13 +0800[diff] [blame^]1865 ps->process_buf = calloc(32 * 32, sizeof(qmf_t));
1866 if (ps->process_buf) {
1867 ps->process_buf_size = 32 * 32 * sizeof(qmf_t);
1868 } else {
1869 ps->process_buf_size = 0;
1870 }
1871 ps->process_buf2 = calloc(32 * 32, sizeof(qmf_t));
1872 if (ps->process_buf2) {
1873 ps->process_buf_size2 = 32 * 32 * sizeof(qmf_t);
1874 } else {
1875 ps->process_buf_size2 = 0;
1876 }
xingri.gaoc18d4472023-02-28 02:51:02 +0000[diff] [blame]1877 return ps;
1878}
1879
1880/* main Parametric Stereo decoding function */
1881uint8_t ps_decode(ps_info *ps, qmf_t X_left[38][64], qmf_t X_right[38][64])
1882{
yuliang.hu585b7202024-09-06 17:28:13 +0800[diff] [blame^]1883 qmf_t (*X_hybrid_left)[32] = NULL;
1884 qmf_t (*X_hybrid_right)[32] = NULL;
1885 memset(ps->process_buf, 0, sizeof(ps->process_buf_size));
1886 memset(ps->process_buf2, 0, sizeof(ps->process_buf_size2));
1887 X_hybrid_left = (qmf_t (*)[32])ps->process_buf;
1888 X_hybrid_right = (qmf_t (*)[32])ps->process_buf2;
xingri.gaoc18d4472023-02-28 02:51:02 +0000[diff] [blame]1889
1890 /* delta decoding of the bitstream data */
1891 ps_data_decode(ps);
1892
1893 /* set up some parameters depending on filterbank type */
1894 if (ps->use34hybrid_bands) {
1895 ps->group_border = (uint8_t*)group_border34;
1896 ps->map_group2bk = (uint16_t*)map_group2bk34;
1897 ps->num_groups = 32 + 18;
1898 ps->num_hybrid_groups = 32;
1899 ps->nr_par_bands = 34;
1900 ps->decay_cutoff = 5;
1901 } else {
1902 ps->group_border = (uint8_t*)group_border20;
1903 ps->map_group2bk = (uint16_t*)map_group2bk20;
1904 ps->num_groups = 10 + 12;
1905 ps->num_hybrid_groups = 10;
1906 ps->nr_par_bands = 20;
1907 ps->decay_cutoff = 3;
1908 }
1909
1910 /* Perform further analysis on the lowest subbands to get a higher
1911 * frequency resolution
1912 */
1913 hybrid_analysis((hyb_info*)ps->hyb, X_left, X_hybrid_left,
1914 ps->use34hybrid_bands, ps->numTimeSlotsRate);
1915
1916 /* decorrelate mono signal */
1917 ps_decorrelate(ps, X_left, X_right, X_hybrid_left, X_hybrid_right);
1918
1919 /* apply mixing and phase parameters */
1920 ps_mix_phase(ps, X_left, X_right, X_hybrid_left, X_hybrid_right);
1921
1922 /* hybrid synthesis, to rebuild the SBR QMF matrices */
1923 hybrid_synthesis((hyb_info*)ps->hyb, X_left, X_hybrid_left,
1924 ps->use34hybrid_bands, ps->numTimeSlotsRate);
1925
1926 hybrid_synthesis((hyb_info*)ps->hyb, X_right, X_hybrid_right,
1927 ps->use34hybrid_bands, ps->numTimeSlotsRate);
1928
1929 return 0;
1930}
1931
1932#endif
1933