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