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