Line data Source code
1 : /******************************************************************************************************
2 :
3 : (C) 2022-2025 IVAS codec Public Collaboration with portions copyright Dolby International AB, Ericsson AB,
4 : Fraunhofer-Gesellschaft zur Foerderung der angewandten Forschung e.V., Huawei Technologies Co. LTD.,
5 : Koninklijke Philips N.V., Nippon Telegraph and Telephone Corporation, Nokia Technologies Oy, Orange,
6 : Panasonic Holdings Corporation, Qualcomm Technologies, Inc., VoiceAge Corporation, and other
7 : contributors to this repository. All Rights Reserved.
8 :
9 : This software is protected by copyright law and by international treaties.
10 : The IVAS codec Public Collaboration consisting of Dolby International AB, Ericsson AB,
11 : Fraunhofer-Gesellschaft zur Foerderung der angewandten Forschung e.V., Huawei Technologies Co. LTD.,
12 : Koninklijke Philips N.V., Nippon Telegraph and Telephone Corporation, Nokia Technologies Oy, Orange,
13 : Panasonic Holdings Corporation, Qualcomm Technologies, Inc., VoiceAge Corporation, and other
14 : contributors to this repository retain full ownership rights in their respective contributions in
15 : the software. This notice grants no license of any kind, including but not limited to patent
16 : license, nor is any license granted by implication, estoppel or otherwise.
17 :
18 : Contributors are required to enter into the IVAS codec Public Collaboration agreement before making
19 : contributions.
20 :
21 : This software is provided "AS IS", without any express or implied warranties. The software is in the
22 : development stage. It is intended exclusively for experts who have experience with such software and
23 : solely for the purpose of inspection. All implied warranties of non-infringement, merchantability
24 : and fitness for a particular purpose are hereby disclaimed and excluded.
25 :
26 : Any dispute, controversy or claim arising under or in relation to providing this software shall be
27 : submitted to and settled by the final, binding jurisdiction of the courts of Munich, Germany in
28 : accordance with the laws of the Federal Republic of Germany excluding its conflict of law rules and
29 : the United Nations Convention on Contracts on the International Sales of Goods.
30 :
31 : *******************************************************************************************************/
32 :
33 : #include <stdint.h>
34 : #include "options.h"
35 : #include "cnst.h"
36 : #include "prot_fx.h"
37 : #include "rom_com.h"
38 : #include "ivas_rom_com.h"
39 : #include <math.h>
40 : #include <assert.h>
41 : #include "wmc_auto.h"
42 : #include "ivas_prot_fx.h"
43 :
44 : /*-------------------------------------------------------------------
45 : * sns_compute_scf_fx()
46 : *
47 : *
48 : *-------------------------------------------------------------------*/
49 :
50 637503 : void sns_compute_scf_fx(
51 : Word32 spectrum[], /* i : Spectrum Q_in */
52 : const PsychoacousticParameters *pPsychParams,
53 : const Word16 L_frame,
54 : Word32 *scf, /* o : Scalefactors (Q16)*/
55 : Word16 Q_in )
56 : {
57 : Word16 i, n, k;
58 : Word32 x[FDNS_NPTS], xs[FDNS_NPTS], mean, xl4[SNS_NPTS], nf, xl[FDNS_NPTS];
59 : Word64 x_64[FDNS_NPTS];
60 : Word64 sum;
61 : Word32 L_tmp;
62 : const Word32 *pow_tilt;
63 : Word16 q_shift, q_out, f_tmp;
64 : Word16 bw, inv_bw, exp;
65 637503 : const UWord8 nBands = pPsychParams->nBands;
66 637503 : move16();
67 637503 : const UWord8 *bandLengths = pPsychParams->bandLengths;
68 :
69 637503 : const Word16 w_0 = 2730; // (1.0f / 12.0f) in Q15
70 637503 : move16();
71 637503 : const Word16 w_1 = 5461; // (2.0f / 12.0f) in Q15
72 637503 : move16();
73 637503 : const Word16 w_2 = 8192; // 0.25f ( 3.0f / 12.0f ) in Q15
74 637503 : move16();
75 637503 : const Word16 w_3 = w_2; // q15
76 637503 : move16();
77 637503 : const Word16 w_4 = w_1; // q15
78 637503 : move16();
79 637503 : const Word16 w_5 = w_0; // q15
80 637503 : move16();
81 :
82 637503 : assert( nBands == FDNS_NPTS );
83 :
84 637503 : set64_fx( x_64, 0, FDNS_NPTS );
85 637503 : set32_fx( x, 0, FDNS_NPTS );
86 :
87 637503 : IF( bandLengths == NULL )
88 : {
89 1361 : bw = shr( L_frame, 6 );
90 1361 : move16();
91 :
92 1361 : exp = norm_s( bw );
93 1361 : inv_bw = div_s( ONE_IN_Q14, shl( bw, exp ) ); // Q:15+14-exp = 29-exp
94 1361 : inv_bw = shl( inv_bw, sub( exp, 14 ) ); // Q15
95 :
96 : /* Energy per band */
97 1361 : k = 0;
98 1361 : move16();
99 88465 : FOR( i = 0; i < nBands; ++i )
100 : {
101 87104 : sum = 0;
102 87104 : move64();
103 821184 : FOR( n = 0; n < bw; ( ++n, ++k ) )
104 : {
105 : /* x[i] += spectrum[k];
106 : inv_bw is for x[i] /= bw; */
107 734080 : sum = W_mac_32_16( sum, spectrum[k], inv_bw ); // Q_in+15+1
108 : }
109 87104 : x_64[i] = sum; // Q_in+16
110 87104 : move64();
111 : }
112 : }
113 : ELSE
114 : {
115 : /* Energy per band */
116 636142 : k = 0;
117 636142 : move16();
118 41349230 : FOR( i = 0; i < nBands; ++i )
119 : {
120 40713088 : exp = norm_s( bandLengths[i] );
121 40713088 : inv_bw = div_s( ONE_IN_Q14, shl( bandLengths[i], exp ) ); // Q:15+14-exp
122 40713088 : inv_bw = shl( inv_bw, sub( exp, 14 ) ); // Q15
123 :
124 40713088 : sum = 0;
125 40713088 : move64();
126 389654016 : FOR( n = 0; n < bandLengths[i]; ( ++n, ++k ) )
127 : {
128 : /* x[i] += spectrum[k];
129 : inv_bw is for x[i] /= bandLengths[i]; */
130 348940928 : sum = W_mac_32_16( sum, spectrum[k], inv_bw ); // Q_in+15+1
131 : }
132 40713088 : x_64[i] = sum; // Q_in+16
133 40713088 : move64();
134 : }
135 : }
136 :
137 : /* Move accumulated values to 32-bit */
138 637503 : q_shift = W_norm_arr( x_64, nBands ); // W_norm_arr return 63 when all the values of the input buffer are zeroes
139 637503 : IF( EQ_16( q_shift, 63 ) )
140 : {
141 : /* If all the values of x_64 are zeros, the scale factor (scf) values will be calculated as zeroes as per the below operations.
142 : To avoid extra computations in such a case, set scf values as zeroes and return. */
143 :
144 11514 : set_zero_fx( scf, SNS_NPTS );
145 :
146 11514 : return;
147 : }
148 :
149 40689285 : FOR( i = 0; i < nBands; ++i )
150 : {
151 40063296 : x[i] = W_extract_h( W_shl( x_64[i], q_shift ) ); // Q: Q_in+16+q_shift-32 = Q_in+q_shift-16
152 40063296 : move32();
153 : }
154 :
155 : /* Smoothing */
156 625989 : xs[0] = Madd_32_16( Mpy_32_16_1( x[0], 24576 /* 0.75 in Q15 */ ), x[1], 8192 /* 0.25 in Q15 */ ); // Q_in+q_shift-16
157 625989 : move32();
158 :
159 39437307 : FOR( i = 1; i < FDNS_NPTS - 1; i++ )
160 : {
161 38811318 : xs[i] = Madd_32_16( Madd_32_16( Mpy_32_16_1( x[i], 16384 /* 0.5 in Q15 */ ), x[i - 1], 8192 /* 0.25 in Q15 */ ), x[i + 1], 8192 /* 0.25 in Q15 */ ); // Q_in+q_shift-16
162 38811318 : move32();
163 : }
164 :
165 625989 : xs[FDNS_NPTS - 1] = Madd_32_16( Mpy_32_16_1( x[FDNS_NPTS - 1], 24576 /* 0.75 in Q15 */ ), x[FDNS_NPTS - 2], 8192 /* 0.25 in Q15 */ ); // Q_in+q_shift-16
166 625989 : move32();
167 :
168 : /* Pre-emphasis */
169 625989 : SWITCH( L_frame )
170 : {
171 56543 : case L_FRAME16k:
172 56543 : pow_tilt = pow_tilt_16k; // Q23
173 56543 : BREAK;
174 247475 : case L_FRAME25_6k:
175 247475 : pow_tilt = pow_tilt_25_6k; // Q23
176 247475 : BREAK;
177 321971 : case L_FRAME32k:
178 321971 : pow_tilt = pow_tilt_32k; // Q23
179 321971 : BREAK;
180 0 : default:
181 0 : pow_tilt = NULL;
182 0 : assert( !"illegal frame length in sns_compute_scf_fx" );
183 : }
184 :
185 40689285 : FOR( i = 0; i < FDNS_NPTS; i++ )
186 : {
187 40063296 : xs[i] = Mpy_32_32( xs[i], pow_tilt[i] ); // Q_in+q_shift-16+23-31 = Q_in+q_shift-24
188 40063296 : move32();
189 : }
190 :
191 : /* Noise floor at -40dB */
192 625989 : sum = 0;
193 625989 : move64();
194 40689285 : FOR( i = 0; i < FDNS_NPTS; i++ )
195 : {
196 40063296 : sum = W_mac_32_16( sum, xs[i], 1 ); // Q_in+q_shift-24+1
197 : }
198 :
199 625989 : q_out = sub( add( Q_in, q_shift ), 24 );
200 :
201 : /* mean = sum / FDNS_NPTS;
202 : -Q6 is for division with FDNS_NPTS and -Q1 is to reduce Q by one */
203 625989 : mean = W_shl_sat_l( sum, -Q7 ); // q_out
204 625989 : nf = Mpy_32_32( mean, 214748 /* powf( 10.0f, -4.0f ) in Q31 */ ); // q_out
205 :
206 625989 : IF( LE_32( nf, L_shl_sat( 256, sub( q_out, 40 ) ) ) ) /* powf( 2.0f, -32.0f ) in Q40 */
207 : {
208 15 : nf = 256;
209 15 : move32();
210 15 : q_out = 40;
211 15 : move16();
212 : }
213 :
214 40689285 : FOR( i = 0; i < FDNS_NPTS; i++ )
215 : {
216 40063296 : xs[i] = L_max( xs[i], nf ); // q_out
217 40063296 : move32();
218 : }
219 :
220 : /* Log-domain */
221 40689285 : FOR( i = 0; i < FDNS_NPTS; i++ )
222 : {
223 : /* xl[i] = logf( xs[i] ) * scale_log;
224 : scale_log = INV_LOG_2 * 0.5f; */
225 :
226 40063296 : exp = norm_l( xs[i] );
227 40063296 : f_tmp = Log2_norm_lc( L_shl( xs[i], exp ) ); // Q15
228 40063296 : exp = sub( sub( 30, exp ), q_out );
229 40063296 : L_tmp = L_mac( L_deposit_h( exp ), f_tmp, 1 ); // Q16
230 40063296 : xl[i] = L_shr( L_tmp, 1 ); // Q16
231 40063296 : move32();
232 : }
233 :
234 : /* Downsampling */
235 625989 : L_tmp = Madd_32_16( Mpy_32_16_1( xl[0], w_0 ), xl[0], w_1 ); // Q16
236 625989 : L_tmp = Madd_32_16( L_tmp, xl[1], w_2 ); // Q16
237 625989 : L_tmp = Madd_32_16( L_tmp, xl[2], w_3 ); // Q16
238 625989 : L_tmp = Madd_32_16( L_tmp, xl[3], w_4 ); // Q16
239 625989 : xl4[0] = Madd_32_16( L_tmp, xl[4], w_5 ); // Q16
240 :
241 :
242 9389835 : FOR( n = 1; n < SNS_NPTS - 1; n++ )
243 : {
244 8763846 : L_tmp = Mpy_32_16_1( xl[4 * n - 1], w_0 ); // Q16
245 8763846 : L_tmp = Madd_32_16( L_tmp, xl[4 * n], w_1 ); // Q16
246 8763846 : L_tmp = Madd_32_16( L_tmp, xl[4 * n + 1], w_2 ); // Q16
247 8763846 : L_tmp = Madd_32_16( L_tmp, xl[4 * n + 2], w_3 ); // Q16
248 8763846 : L_tmp = Madd_32_16( L_tmp, xl[4 * n + 3], w_4 ); // Q16
249 8763846 : xl4[n] = Madd_32_16( L_tmp, xl[4 * n + 4], w_5 ); // Q16
250 8763846 : move32();
251 : }
252 :
253 625989 : L_tmp = Mpy_32_16_1( xl[FDNS_NPTS - 5], w_0 ); // Q16
254 625989 : L_tmp = Madd_32_16( L_tmp, xl[FDNS_NPTS - 4], w_1 ); // Q16
255 625989 : L_tmp = Madd_32_16( L_tmp, xl[FDNS_NPTS - 3], w_2 ); // Q16
256 625989 : L_tmp = Madd_32_16( L_tmp, xl[FDNS_NPTS - 2], w_3 ); // Q16
257 625989 : L_tmp = Madd_32_16( L_tmp, xl[FDNS_NPTS - 1], w_4 ); // Q16
258 625989 : xl4[SNS_NPTS - 1] = Madd_32_16( L_tmp, xl[FDNS_NPTS - 1], w_5 ); // Q16
259 625989 : move32();
260 :
261 : /* Remove mean and scaling */
262 625989 : sum = 0;
263 625989 : move64();
264 10641813 : FOR( i = 0; i < SNS_NPTS; i++ )
265 : {
266 10015824 : sum = W_mac_32_16( sum, xl4[i], 1 ); // Q16+1
267 : }
268 : /* mean = sum / SNS_NPTS;
269 : -Q4 is for division with SNS_NPTS and -Q1 is to reduce Q by one */
270 625989 : mean = W_shl_sat_l( sum, -Q5 ); // Q16
271 :
272 10641813 : FOR( i = 0; i < SNS_NPTS; i++ )
273 : {
274 10015824 : scf[i] = Mpy_32_16_1( L_sub( xl4[i], mean ), 27853 /* 0.85 in in Q15 */ ); // Q16
275 10015824 : move32();
276 : }
277 :
278 625989 : return;
279 : }
280 :
281 : /*-------------------------------------------------------------------
282 : * sns_interpolate_scalefactors_fx()
283 : *
284 : *
285 : *-------------------------------------------------------------------*/
286 :
287 1904221 : void sns_interpolate_scalefactors_fx(
288 : Word32 *scf_int, /* o : interpolated scalefactors for spectrum shaping q16*/
289 : const Word32 *scf, /* i : sns scalefactors as derived from the signal or read from the bitstream Q16*/
290 : Word16 encoder_side /* i : flag, if scalefactors have to be inverted */
291 : )
292 : {
293 : Word16 n;
294 : Word32 L_tmp;
295 : Word16 exp;
296 :
297 : /* Interpolation */
298 1904221 : scf_int[0] = scf[0]; // q16
299 1904221 : scf_int[1] = scf[0]; // q16
300 1904221 : move32();
301 1904221 : move32();
302 :
303 30467536 : FOR( n = 0; n <= M - 2; n++ )
304 : {
305 28563315 : scf_int[n * 4 + 2] = Madd_32_16( scf[n], L_sub( scf[n + 1], scf[n] ), 4096 /* 4096 -> 1/8 in Q15 */ ); // q16
306 28563315 : scf_int[n * 4 + 3] = Madd_32_16( scf[n], L_sub( scf[n + 1], scf[n] ), 12288 /* 12288 -> 3/8 in Q15 */ ); // q16
307 28563315 : scf_int[n * 4 + 4] = Madd_32_16( scf[n], L_sub( scf[n + 1], scf[n] ), 20480 /* 20480 -> 5/8 in Q15 */ ); // q16
308 28563315 : scf_int[n * 4 + 5] = Madd_32_16( scf[n], L_sub( scf[n + 1], scf[n] ), 28672 /* 28672 -> 7/8 in Q15 */ ); // q16
309 28563315 : move32();
310 28563315 : move32();
311 28563315 : move32();
312 28563315 : move32();
313 : }
314 :
315 1904221 : scf_int[FDNS_NPTS - 2] = Madd_32_16( scf[M - 1], L_sub( scf[M - 1], scf[M - 2] ), 4096 /* 4096 -> 1/8 in Q15 */ ); // q16
316 1904221 : scf_int[FDNS_NPTS - 1] = Madd_32_16( scf[M - 1], L_sub( scf[M - 1], scf[M - 2] ), 12288 /* 12288 -> 3/8 in Q15 */ ); // q16
317 1904221 : move32();
318 1904221 : move32();
319 :
320 : /* Inversion at encoder-side */
321 1904221 : IF( encoder_side == ENC /*0*/ )
322 : {
323 76530220 : FOR( n = 0; n < FDNS_NPTS; n++ )
324 : {
325 75352832 : scf_int[n] = L_negate( scf_int[n] ); // q16
326 75352832 : move32();
327 : }
328 : }
329 :
330 : /* Linear domain */
331 123774365 : FOR( n = 0; n < FDNS_NPTS; n++ )
332 : {
333 121870144 : L_tmp = BASOP_util_Pow2( scf_int[n], 15 /*31-Q16*/, &exp ); // Q=31-exp
334 121870144 : exp = sub( 15, exp );
335 121870144 : scf_int[n] = L_shr( L_tmp, exp ); // q16
336 121870144 : move32();
337 : }
338 :
339 1904221 : return;
340 : }
341 :
342 :
343 : /*-------------------------------------------------------------------
344 : * sns_shape_spectrum_fx()
345 : *
346 : *
347 : *-------------------------------------------------------------------*/
348 :
349 1755443 : void sns_shape_spectrum_fx(
350 : Word32 spectrum[], /* i/o: spectrum to be shaped Input Q is q_spectrum and ouput Q is (q_spectrum-1)*/
351 : Word16 *q_spectrum, /* i/o: Q of spectrum */
352 : const PsychoacousticParameters *pPsychParams, /* i : psychoacoustic parameters used to get the frequency bands */
353 : const Word32 *scf_int, /* i : already interpolated SNS scalefactors */
354 : const Word16 q_scf_int, /* i : Q of interpolated SNS scalefactors q_scf_int*/
355 : const Word16 L_frame, /* i : frame length */
356 : Word16 *length )
357 : {
358 1755443 : Word16 i, n, k, tmp_k, bw, q_tmp = 0, shift, min_shift = 63;
359 1755443 : move16();
360 1755443 : move16();
361 : Word64 L64_tmp[L_FRAME48k];
362 1755443 : const UWord8 nBands = pPsychParams->nBands;
363 1755443 : move16();
364 1755443 : const UWord8 *bandLengths = pPsychParams->bandLengths;
365 :
366 1755443 : IF( bandLengths == NULL )
367 : {
368 3359 : bw = divide3216( shl( L_frame, 1 ), nBands ); // q0
369 :
370 : /* Shape spectrum */
371 3359 : k = 0;
372 3359 : move16();
373 218335 : FOR( i = 0; i < nBands; ++i )
374 : {
375 2446016 : FOR( n = 0; n < bw; ( ++n, ++k ) )
376 : {
377 2231040 : L64_tmp[k] = W_mult_32_32( spectrum[k], scf_int[i] ); // Q = q_spectrum + q_scf_int + 1
378 2231040 : move64();
379 2231040 : shift = W_norm( L64_tmp[k] );
380 2231040 : if ( LT_16( shift, min_shift ) )
381 : {
382 9004 : min_shift = shift;
383 9004 : move16();
384 : }
385 : }
386 : }
387 3359 : tmp_k = k;
388 3359 : move16();
389 3359 : if ( length != NULL )
390 : {
391 3359 : *length = k;
392 3359 : move16();
393 : }
394 3359 : q_tmp = sub( add( add( *q_spectrum, q_scf_int ), min_shift ), 32 );
395 3359 : if ( GT_16( q_tmp, 30 ) )
396 : {
397 34 : q_tmp = 30;
398 34 : move16();
399 : }
400 2234399 : FOR( k = 0; k < tmp_k; k++ )
401 : {
402 2231040 : L64_tmp[k] = W_shr( L64_tmp[k], sub( add( *q_spectrum, q_scf_int ), q_tmp ) ); // Q = q_tmp+1
403 2231040 : move64();
404 2231040 : spectrum[k] = W_sat_l( L64_tmp[k] ); // Q = q_tmp+1
405 2231040 : move64();
406 : }
407 3359 : *q_spectrum = q_tmp;
408 3359 : move16();
409 : }
410 : ELSE
411 : {
412 : /* Shape spectrum */
413 1752084 : k = 0;
414 1752084 : move16();
415 113885460 : FOR( i = 0; i < nBands; ++i )
416 : {
417 1074648576 : FOR( n = 0; n < bandLengths[i]; ( ++n, ++k ) )
418 : {
419 962515200 : L64_tmp[k] = W_mult_32_32( spectrum[k], scf_int[i] ); // Q = q_spectrum + q_scf_int + 1
420 962515200 : move64();
421 962515200 : shift = W_norm( L64_tmp[k] );
422 962515200 : test();
423 962515200 : if ( LT_16( shift, min_shift ) && NE_64( L64_tmp[k], 0 ) )
424 : {
425 5219594 : min_shift = shift;
426 5219594 : move16();
427 : }
428 : }
429 : }
430 1752084 : tmp_k = k;
431 1752084 : move16();
432 1752084 : if ( length != NULL )
433 : {
434 1750510 : *length = k;
435 1750510 : move16();
436 : }
437 1752084 : q_tmp = sub( add( add( *q_spectrum, q_scf_int ), min_shift ), 32 );
438 1752084 : if ( GT_16( q_tmp, 30 ) )
439 : {
440 1817 : q_tmp = 30;
441 1817 : move16();
442 : }
443 964267284 : FOR( k = 0; k < tmp_k; k++ )
444 : {
445 962515200 : L64_tmp[k] = W_shr( L64_tmp[k], sub( add( *q_spectrum, q_scf_int ), q_tmp ) ); // q_tmp+1
446 962515200 : move64();
447 962515200 : spectrum[k] = W_sat_l( L64_tmp[k] ); // Q = q_tmp+1
448 962515200 : move64();
449 : }
450 1752084 : *q_spectrum = q_tmp;
451 1752084 : move16();
452 : }
453 :
454 1755443 : return;
455 : }
|