1 /**
2 * Copyright (C) 2022 The Android Open Source Project
3 *
4 * Licensed under the Apache License, Version 2.0 (the "License");
5 * you may not use this file except in compliance with the License.
6 * You may obtain a copy of the License at
7 *
8 * http://www.apache.org/licenses/LICENSE-2.0
9 *
10 * Unless required by applicable law or agreed to in writing, software
11 * distributed under the License is distributed on an "AS IS" BASIS,
12 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13 * See the License for the specific language governing permissions and
14 * limitations under the License.
15 */
16 /*------------------------------------------------------------------------------
17 *
18 * Subband processing consists of:
19 * inverse quantisation (defined in a separate file),
20 * predictor coefficient update (Pole and Zero Coeff update),
21 * predictor filtering.
22 *
23 *----------------------------------------------------------------------------*/
24
25 #ifndef SUBBANDFUNCTIONSCOMMON_H
26 #define SUBBANDFUNCTIONSCOMMON_H
27
28 enum reg64_reg { reg64_H = 1, reg64_L = 0 };
29
30 void processSubband(const int32_t qCode, const int32_t ditherVal, Subband_data* SubbandDataPt,
31 IQuantiser_data* iqDataPt);
32 void processSubbandLL(const int32_t qCode, const int32_t ditherVal, Subband_data* SubbandDataPt,
33 IQuantiser_data* iqDataPt);
34 void processSubbandHL(const int32_t qCode, const int32_t ditherVal, Subband_data* SubbandDataPt,
35 IQuantiser_data* iqDataPt);
36
37 /* Function to carry out inverse quantisation for LL, LH and HH subband types */
invertQuantisation(const int32_t qCode,const int32_t ditherVal,IQuantiser_data * iqdata_pt)38 XBT_INLINE_ void invertQuantisation(const int32_t qCode, const int32_t ditherVal,
39 IQuantiser_data* iqdata_pt) {
40 int32_t invQ;
41 int32_t index;
42 int32_t acc;
43 reg64_t tmp_r64;
44 int64_t tmp_acc;
45 int32_t tmp_accL;
46 int32_t tmp_accH;
47 uint32_t tmp_round0;
48 uint32_t tmp_round1;
49
50 unsigned u16t;
51 /* log delta leak value (Q23) */
52 const uint32_t logDeltaLeakVal = 0x7F6CL;
53
54 /* Turn the quantised code back into an index into the threshold table. This
55 * involves bitwise inversion of the code (if -ve) and adding 1 (phantom
56 * element at table base). Then set invQ to be +/- the threshold value,
57 * depending on the code sign. */
58 index = qCode;
59 if (qCode < 0) {
60 index = (~index);
61 }
62 index = index + 1;
63 invQ = iqdata_pt->thresholdTablePtr_sl1[index];
64 if (qCode < 0) {
65 invQ = -invQ;
66 }
67
68 /* Load invQ into the accumulator. Add the product of the dither value times
69 * the indexed dither table value. Then get the result back from the
70 * accumulator as an updated invQ. */
71 tmp_r64.s64 = ((int64_t)ditherVal * iqdata_pt->ditherTablePtr_sf1[index]);
72 tmp_r64.s32.h += invQ >> 1;
73
74 acc = tmp_r64.s32.h;
75
76 tmp_round1 = tmp_r64.s32.h & 0x00000001L;
77 if (tmp_r64.u32.l >= 0x80000000) {
78 acc++;
79 }
80 if (tmp_round1 == 0 && tmp_r64.s32.l == (int32_t)0x80000000L) {
81 acc--;
82 }
83 acc = ssat24(acc);
84
85 invQ = acc;
86
87 /* Scale invQ by the current delta value. Left-shift the result (in the
88 * accumulator) by 4 positions for the delta scaling. Get the updated invQ
89 * back from the accumulator. */
90
91 u16t = iqdata_pt->logDelta;
92 tmp_acc = ((int64_t)invQ * iqdata_pt->delta);
93 tmp_accL = u16t * logDeltaLeakVal;
94 tmp_accH = iqdata_pt->incrTablePtr[index];
95 acc = (int32_t)(tmp_acc >> (23 - deltaScale));
96 invQ = ssat24(acc);
97
98 /* Now update the value of logDelta. Load the accumulator with the index
99 * value of the logDelta increment table. Add the product of the current
100 * logDelta scaled by a leaky coefficient (16310 in Q14). Get the value back
101 * from the accumulator. */
102 tmp_accH += tmp_accL >> (32 - 17);
103
104 acc = tmp_accH;
105
106 tmp_r64.u32.l = ((uint32_t)tmp_accL << 17);
107 tmp_r64.s32.h = tmp_accH;
108
109 tmp_round0 = tmp_r64.u32.l;
110 tmp_round1 = (int32_t)(tmp_r64.u64 >> 1);
111 if (tmp_round0 >= 0x80000000L) {
112 acc++;
113 }
114 if (tmp_round1 == 0x40000000L) {
115 acc--;
116 }
117
118 /* Limit the updated logDelta between 0 and its subband-specific maximum. */
119 if (acc < 0) {
120 acc = 0;
121 }
122 if (acc > iqdata_pt->maxLogDelta) {
123 acc = iqdata_pt->maxLogDelta;
124 }
125
126 iqdata_pt->logDelta = (uint16_t)acc;
127
128 /* The updated value of delta is the logTable output (indexed by 5 bits from
129 * the updated logDelta) shifted by a value involving the logDelta minimum
130 * and the updated logDelta itself. */
131 iqdata_pt->delta = iqdata_pt->iquantTableLogPtr[(acc >> 3) & 0x1f] >>
132 (22 - 25 - iqdata_pt->minLogDelta - (acc >> 8));
133
134 iqdata_pt->invQ = invQ;
135 }
136
137 /* Function to carry out inverse quantisation for a HL subband type */
invertQuantisationHL(const int32_t qCode,const int32_t ditherVal,IQuantiser_data * iqdata_pt)138 XBT_INLINE_ void invertQuantisationHL(const int32_t qCode, const int32_t ditherVal,
139 IQuantiser_data* iqdata_pt) {
140 int32_t invQ;
141 int32_t index;
142 int32_t acc;
143 reg64_t tmp_r64;
144 int64_t tmp_acc;
145 int32_t tmp_accL;
146 int32_t tmp_accH;
147 uint32_t tmp_round0;
148 uint32_t tmp_round1;
149
150 unsigned u16t;
151 /* log delta leak value (Q23) */
152 const uint32_t logDeltaLeakVal = 0x7F6CL;
153
154 /* Turn the quantised code back into an index into the threshold table. This
155 * involves bitwise inversion of the code (if -ve) and adding 1 (phantom
156 * element at table base). Then set invQ to be +/- the threshold value,
157 * depending on the code sign. */
158 index = qCode;
159 if (qCode < 0) {
160 index = (~index);
161 }
162 index = index + 1;
163 invQ = iqdata_pt->thresholdTablePtr_sl1[index];
164 if (qCode < 0) {
165 invQ = -invQ;
166 }
167
168 /* Load invQ into the accumulator. Add the product of the dither value times
169 * the indexed dither table value. Then get the result back from the
170 * accumulator as an updated invQ. */
171 tmp_r64.s64 = ((int64_t)ditherVal * iqdata_pt->ditherTablePtr_sf1[index]);
172 tmp_r64.s32.h += invQ >> 1;
173
174 acc = tmp_r64.s32.h;
175
176 tmp_round1 = tmp_r64.s32.h & 0x00000001L;
177 if (tmp_r64.u32.l >= 0x80000000) {
178 acc++;
179 }
180 if (tmp_round1 == 0 && tmp_r64.u32.l == 0x80000000L) {
181 acc--;
182 }
183 acc = ssat24(acc);
184
185 invQ = acc;
186
187 /* Scale invQ by the current delta value. Left-shift the result (in the
188 * accumulator) by 4 positions for the delta scaling. Get the updated invQ
189 * back from the accumulator. */
190 u16t = iqdata_pt->logDelta;
191 tmp_acc = ((int64_t)invQ * iqdata_pt->delta);
192 tmp_accL = u16t * logDeltaLeakVal;
193 tmp_accH = iqdata_pt->incrTablePtr[index];
194 acc = (int32_t)(tmp_acc >> (23 - deltaScale));
195 invQ = acc;
196
197 /* Now update the value of logDelta. Load the accumulator with the index
198 * value of the logDelta increment table. Add the product of the current
199 * logDelta scaled by a leaky coefficient (16310 in Q14). Get the value back
200 * from the accumulator. */
201 tmp_accH += tmp_accL >> (32 - 17);
202
203 acc = tmp_accH;
204
205 tmp_r64.u32.l = ((uint32_t)tmp_accL << 17);
206 tmp_r64.s32.h = tmp_accH;
207
208 tmp_round0 = tmp_r64.u32.l;
209 tmp_round1 = (int32_t)(tmp_r64.u64 >> 1);
210 if (tmp_round0 >= 0x80000000L) {
211 acc++;
212 }
213 if (tmp_round1 == 0x40000000L) {
214 acc--;
215 }
216
217 /* Limit the updated logDelta between 0 and its subband-specific maximum. */
218 if (acc < 0) {
219 acc = 0;
220 }
221 if (acc > iqdata_pt->maxLogDelta) {
222 acc = iqdata_pt->maxLogDelta;
223 }
224
225 iqdata_pt->logDelta = (uint16_t)acc;
226
227 /* The updated value of delta is the logTable output (indexed by 5 bits from
228 * the updated logDelta) shifted by a value involving the logDelta minimum
229 * and the updated logDelta itself. */
230 iqdata_pt->delta = iqdata_pt->iquantTableLogPtr[(acc >> 3) & 0x1f] >>
231 (22 - 25 - iqdata_pt->minLogDelta - (acc >> 8));
232
233 iqdata_pt->invQ = invQ;
234 }
235
236 /* Function to carry out prediction ARMA filtering for the current subband
237 * performPredictionFiltering should only be used for HH and LH subband! */
performPredictionFiltering(const int32_t invQ,Subband_data * SubbandDataPt)238 XBT_INLINE_ void performPredictionFiltering(const int32_t invQ, Subband_data* SubbandDataPt) {
239 int32_t poleVal;
240 int32_t acc;
241 int64_t accL;
242 uint32_t pointer;
243 int32_t poleDelayLine;
244 int32_t predVal;
245 int32_t* zeroCoeffPt = SubbandDataPt->m_ZeroCoeffData.m_zeroCoeff;
246 int32_t* poleCoeff = SubbandDataPt->m_PoleCoeffData.m_poleCoeff;
247 int32_t zData0;
248 int32_t* cbuf_pt;
249 int32_t invQincr_pos;
250 int32_t invQincr_neg;
251 int32_t k;
252 int32_t oldZData;
253 /* Pole coefficient and data indices */
254 enum { a1 = 0, a2 = 1 };
255
256 /* Write the newest pole input sample to the pole delay line.
257 * Ensure the sum of the current dequantised error and the previous
258 * predictor output is saturated if necessary. */
259 poleDelayLine = invQ + SubbandDataPt->m_predData.m_predVal;
260
261 poleDelayLine = ssat24(poleDelayLine);
262
263 /* Pole filter convolution. Shift convolution result 1 place to the left
264 * before retrieving it, since the pole coefficients are Q22 (data is Q23)
265 * and we want a Q23 result */
266 accL = ((int64_t)poleCoeff[a2] * (int64_t)SubbandDataPt->m_predData.m_poleDelayLine[a2]);
267 /* Update the pole delay line for the next pass by writing the new input
268 * sample into the 2nd element */
269 SubbandDataPt->m_predData.m_poleDelayLine[a2] = poleDelayLine;
270 accL += ((int64_t)poleCoeff[a1] * (int64_t)poleDelayLine);
271 poleVal = (int32_t)(accL >> 22);
272 poleVal = ssat24(poleVal);
273
274 /* Create (2^(-7)) * sgn(invQ) in Q22 format. */
275 if (invQ == 0) {
276 invQincr_pos = 0L;
277 } else {
278 invQincr_pos = 0x800000;
279 }
280 if (invQ < 0L) {
281 invQincr_pos = -invQincr_pos;
282 }
283
284 invQincr_neg = 0x0080 - invQincr_pos;
285 invQincr_pos += 0x0080;
286
287 pointer = (SubbandDataPt->m_predData.m_zeroDelayLine.pointer++) + 12;
288 cbuf_pt = &SubbandDataPt->m_predData.m_zeroDelayLine.buffer[pointer];
289 /* partial manual unrolling to improve performance */
290 if (SubbandDataPt->m_predData.m_zeroDelayLine.pointer >= 12) {
291 SubbandDataPt->m_predData.m_zeroDelayLine.pointer = 0;
292 }
293
294 SubbandDataPt->m_predData.m_zeroDelayLine.modulo = invQ;
295
296 /* Iterate over the number of coefficients for this subband */
297 oldZData = invQ;
298 accL = 0;
299 for (k = 0; k < 12; k++) {
300 uint32_t tmp_round0;
301 int32_t coeffValue;
302
303 zData0 = (*(cbuf_pt--));
304 coeffValue = *(zeroCoeffPt + k);
305 if (zData0 < 0L) {
306 acc = invQincr_neg - coeffValue;
307 } else {
308 acc = invQincr_pos - coeffValue;
309 }
310 tmp_round0 = acc;
311 acc = (acc >> 8) + coeffValue;
312 if (((tmp_round0 << 23) ^ 0x80000000) == 0) {
313 acc--;
314 }
315 accL += (int64_t)acc * (int64_t)(oldZData);
316 oldZData = zData0;
317 *(zeroCoeffPt + k) = acc;
318 }
319
320 acc = (int32_t)(accL >> 22);
321 acc = ssat24(acc);
322 /* Predictor output is the sum of the pole and zero filter outputs. Ensure
323 * this is saturated, if necessary. */
324 predVal = acc + poleVal;
325 predVal = ssat24(predVal);
326 SubbandDataPt->m_predData.m_zeroVal = acc;
327 SubbandDataPt->m_predData.m_predVal = predVal;
328
329 /* Update the zero filter delay line by writing the new input sample to the
330 * circular buffer. */
331 SubbandDataPt->m_predData.m_zeroDelayLine
332 .buffer[SubbandDataPt->m_predData.m_zeroDelayLine.pointer] =
333 SubbandDataPt->m_predData.m_zeroDelayLine.modulo;
334 SubbandDataPt->m_predData.m_zeroDelayLine
335 .buffer[SubbandDataPt->m_predData.m_zeroDelayLine.pointer + 12] =
336 SubbandDataPt->m_predData.m_zeroDelayLine.modulo;
337 }
338
performPredictionFilteringLL(const int32_t invQ,Subband_data * SubbandDataPt)339 XBT_INLINE_ void performPredictionFilteringLL(const int32_t invQ, Subband_data* SubbandDataPt) {
340 int32_t poleVal;
341 int32_t acc;
342 int64_t accL;
343 uint32_t pointer;
344 int32_t poleDelayLine;
345 int32_t predVal;
346 int32_t* zeroCoeffPt = SubbandDataPt->m_ZeroCoeffData.m_zeroCoeff;
347 int32_t* poleCoeff = SubbandDataPt->m_PoleCoeffData.m_poleCoeff;
348 int32_t* cbuf_pt;
349 int32_t invQincr_pos;
350 int32_t invQincr_neg;
351 int32_t k;
352 int32_t oldZData;
353 /* Pole coefficient and data indices */
354 enum { a1 = 0, a2 = 1 };
355
356 /* Write the newest pole input sample to the pole delay line.
357 * Ensure the sum of the current dequantised error and the previous
358 * predictor output is saturated if necessary. */
359 poleDelayLine = invQ + SubbandDataPt->m_predData.m_predVal;
360
361 poleDelayLine = ssat24(poleDelayLine);
362
363 /* Pole filter convolution. Shift convolution result 1 place to the left
364 * before retrieving it, since the pole coefficients are Q22 (data is Q23)
365 * and we want a Q23 result */
366 accL = ((int64_t)poleCoeff[a2] * (int64_t)SubbandDataPt->m_predData.m_poleDelayLine[a2]);
367 /* Update the pole delay line for the next pass by writing the new input
368 * sample into the 2nd element */
369 SubbandDataPt->m_predData.m_poleDelayLine[a2] = poleDelayLine;
370 accL += ((int64_t)poleCoeff[a1] * (int64_t)poleDelayLine);
371 poleVal = (int32_t)(accL >> 22);
372 poleVal = ssat24(poleVal);
373 // store poleVal to free one register.
374 SubbandDataPt->m_predData.m_predVal = poleVal;
375
376 /* Create (2^(-7)) * sgn(invQ) in Q22 format. */
377 if (invQ == 0) {
378 invQincr_pos = 0L;
379 } else {
380 invQincr_pos = 0x800000;
381 }
382 if (invQ < 0L) {
383 invQincr_pos = -invQincr_pos;
384 }
385
386 invQincr_neg = 0x0080 - invQincr_pos;
387 invQincr_pos += 0x0080;
388
389 pointer = (SubbandDataPt->m_predData.m_zeroDelayLine.pointer++) + 24;
390 cbuf_pt = &SubbandDataPt->m_predData.m_zeroDelayLine.buffer[pointer];
391 /* partial manual unrolling to improve performance */
392 if (SubbandDataPt->m_predData.m_zeroDelayLine.pointer >= 24) {
393 SubbandDataPt->m_predData.m_zeroDelayLine.pointer = 0;
394 }
395
396 SubbandDataPt->m_predData.m_zeroDelayLine.modulo = invQ;
397
398 /* Iterate over the number of coefficients for this subband */
399
400 oldZData = invQ;
401 accL = 0;
402 for (k = 0; k < 24; k++) {
403 int32_t zData0;
404 int32_t coeffValue;
405
406 zData0 = (*(cbuf_pt--));
407 coeffValue = *(zeroCoeffPt + k);
408 if (zData0 < 0L) {
409 acc = invQincr_neg - coeffValue;
410 } else {
411 acc = invQincr_pos - coeffValue;
412 }
413 if (((acc << 23) ^ 0x80000000) == 0) {
414 coeffValue--;
415 }
416 acc = (acc >> 8) + coeffValue;
417 accL += (int64_t)acc * (int64_t)(oldZData);
418 oldZData = zData0;
419 *(zeroCoeffPt + k) = acc;
420 }
421
422 acc = (int32_t)(accL >> 22);
423 acc = ssat24(acc);
424 /* Predictor output is the sum of the pole and zero filter outputs. Ensure
425 * this is saturated, if necessary. */
426 // recover value of PoleVal stored at beginning of routine...
427 predVal = acc + SubbandDataPt->m_predData.m_predVal;
428 predVal = ssat24(predVal);
429 SubbandDataPt->m_predData.m_zeroVal = acc;
430 SubbandDataPt->m_predData.m_predVal = predVal;
431
432 /* Update the zero filter delay line by writing the new input sample to the
433 * circular buffer. */
434 SubbandDataPt->m_predData.m_zeroDelayLine
435 .buffer[SubbandDataPt->m_predData.m_zeroDelayLine.pointer] =
436 SubbandDataPt->m_predData.m_zeroDelayLine.modulo;
437 SubbandDataPt->m_predData.m_zeroDelayLine
438 .buffer[SubbandDataPt->m_predData.m_zeroDelayLine.pointer + 24] =
439 SubbandDataPt->m_predData.m_zeroDelayLine.modulo;
440 }
441
performPredictionFilteringHL(const int32_t invQ,Subband_data * SubbandDataPt)442 XBT_INLINE_ void performPredictionFilteringHL(const int32_t invQ, Subband_data* SubbandDataPt) {
443 int32_t poleVal;
444 int32_t acc;
445 int64_t accL;
446 uint32_t pointer;
447 int32_t poleDelayLine;
448 int32_t predVal;
449 int32_t* zeroCoeffPt = SubbandDataPt->m_ZeroCoeffData.m_zeroCoeff;
450 int32_t* poleCoeff = SubbandDataPt->m_PoleCoeffData.m_poleCoeff;
451 int32_t zData0;
452 int32_t* cbuf_pt;
453 int32_t invQincr_pos;
454 int32_t invQincr_neg;
455 int32_t k;
456 int32_t oldZData;
457 const int32_t roundCte = 0x80000000;
458 /* Pole coefficient and data indices */
459 enum { a1 = 0, a2 = 1 };
460
461 /* Write the newest pole input sample to the pole delay line.
462 * Ensure the sum of the current dequantised error and the previous
463 * predictor output is saturated if necessary. */
464 poleDelayLine = invQ + SubbandDataPt->m_predData.m_predVal;
465
466 poleDelayLine = ssat24(poleDelayLine);
467
468 /* Pole filter convolution. Shift convolution result 1 place to the left
469 * before retrieving it, since the pole coefficients are Q22 (data is Q23)
470 * and we want a Q23 result */
471 accL = ((int64_t)poleCoeff[a2] * (int64_t)SubbandDataPt->m_predData.m_poleDelayLine[a2]);
472 /* Update the pole delay line for the next pass by writing the new input
473 * sample into the 2nd element */
474 SubbandDataPt->m_predData.m_poleDelayLine[a2] = poleDelayLine;
475 accL += ((int64_t)poleCoeff[a1] * (int64_t)poleDelayLine);
476 poleVal = (int32_t)(accL >> 22);
477 poleVal = ssat24(poleVal);
478
479 /* Create (2^(-7)) * sgn(invQ) in Q22 format. */
480 invQincr_pos = 0L;
481 if (invQ != 0) {
482 invQincr_pos = 0x800000;
483 }
484 if (invQ < 0L) {
485 invQincr_pos = -invQincr_pos;
486 }
487
488 invQincr_neg = 0x0080 - invQincr_pos;
489 invQincr_pos += 0x0080;
490
491 pointer = (SubbandDataPt->m_predData.m_zeroDelayLine.pointer++) + 6;
492 cbuf_pt = &SubbandDataPt->m_predData.m_zeroDelayLine.buffer[pointer];
493 /* partial manual unrolling to improve performance */
494 if (SubbandDataPt->m_predData.m_zeroDelayLine.pointer >= 6) {
495 SubbandDataPt->m_predData.m_zeroDelayLine.pointer = 0;
496 }
497
498 SubbandDataPt->m_predData.m_zeroDelayLine.modulo = invQ;
499
500 /* Iterate over the number of coefficients for this subband */
501 oldZData = invQ;
502 accL = 0;
503
504 for (k = 0; k < 6; k++) {
505 uint32_t tmp_round0;
506 int32_t coeffValue;
507
508 zData0 = (*(cbuf_pt--));
509 coeffValue = *(zeroCoeffPt + k);
510 if (zData0 < 0L) {
511 acc = invQincr_neg - coeffValue;
512 } else {
513 acc = invQincr_pos - coeffValue;
514 }
515 tmp_round0 = acc;
516 acc = (acc >> 8) + coeffValue;
517 if (((tmp_round0 << 23) ^ roundCte) == 0) {
518 acc--;
519 }
520 accL += (int64_t)acc * (int64_t)(oldZData);
521 oldZData = zData0;
522 *(zeroCoeffPt + k) = acc;
523 }
524
525 acc = (int32_t)(accL >> 22);
526 acc = ssat24(acc);
527 /* Predictor output is the sum of the pole and zero filter outputs. Ensure
528 * this is saturated, if necessary. */
529 predVal = acc + poleVal;
530 predVal = ssat24(predVal);
531 SubbandDataPt->m_predData.m_zeroVal = acc;
532 SubbandDataPt->m_predData.m_predVal = predVal;
533
534 /* Update the zero filter delay line by writing the new input sample to the
535 * circular buffer. */
536 SubbandDataPt->m_predData.m_zeroDelayLine
537 .buffer[SubbandDataPt->m_predData.m_zeroDelayLine.pointer] =
538 SubbandDataPt->m_predData.m_zeroDelayLine.modulo;
539 SubbandDataPt->m_predData.m_zeroDelayLine
540 .buffer[SubbandDataPt->m_predData.m_zeroDelayLine.pointer + 6] =
541 SubbandDataPt->m_predData.m_zeroDelayLine.modulo;
542 }
543
544 #endif // SUBBANDFUNCTIONSCOMMON_H
545