Line data Source code
1 : /*====================================================================================
2 : EVS Codec 3GPP TS26.452 Aug 12, 2021. Version 16.3.0
3 : ====================================================================================*/
4 :
5 : #include <assert.h>
6 : #include <stdint.h>
7 : #include "options.h"
8 : #include "basop_util.h"
9 : #include "options.h"
10 : #include "rom_basop_util.h"
11 : #include "rom_com.h"
12 : #include "prot_fx.h"
13 : #include "prot_fx_enc.h"
14 : #include "ivas_prot_fx.h"
15 :
16 : #define FFT_SCALING_512 1073741824 // Q22
17 : #define FFT_SCALING_640 1342177280 // Q22
18 :
19 : #define DELTA_SHIFT 2
20 : #define DELTA_SHIFT_LD64 67108864l /*DELTA_SHIFT/64.0 Q31*/
21 : #define CNG_HS 4 /* 4 bit headroom for dot product */
22 : #define CNG_S 6 /* 1 sign bit, 6 bit integer part, 9 bit frational part for input and output data */
23 :
24 : /*-------------------------------------------------------------------
25 : * Local function prototypes
26 : *-------------------------------------------------------------------*/
27 : static void getmidbands( const Word16 *part, const Word16 npart, Word16 *midband, Word16 *psize, Word16 *psize_norm, Word16 *psize_norm_exp, Word16 *psize_inv );
28 :
29 :
30 : /*-------------------------------------------------------------------
31 : * createFdCngCom()
32 : *
33 : * Create an instance of type FD_CNG_COM
34 : *-------------------------------------------------------------------*/
35 6430 : ivas_error createFdCngCom_fx(
36 : HANDLE_FD_CNG_COM *hFdCngCom /* i/o: FD_CNG structure containing all buffers and variables */
37 : )
38 : {
39 : HANDLE_FD_CNG_COM hs;
40 :
41 : /* Allocate memory */
42 6430 : hs = (HANDLE_FD_CNG_COM) malloc( sizeof( FD_CNG_COM ) );
43 :
44 6430 : if ( hs == NULL )
45 : {
46 0 : return IVAS_ERROR( IVAS_ERR_FAILED_ALLOC, "Can not allocate memory for FD CNG COM" );
47 : }
48 :
49 6430 : *hFdCngCom = hs;
50 :
51 6430 : return IVAS_ERR_OK;
52 : }
53 : /*-------------------------------------------------------------------
54 : * initFdCngCom()
55 : *
56 : *
57 : *-------------------------------------------------------------------*/
58 :
59 6427 : void ivas_initFdCngCom_fx( HANDLE_FD_CNG_COM hFdCngCom, Word16 scale )
60 : {
61 : /* Calculate CLDFB scaling factor */
62 : /* shl(i_mult2(scale, scale), 3) does not fit in 16 bit */
63 : /*hFdCngCom->scalingFactor = div_s(1, shl(i_mult2(scale, scale), 3));*/
64 6427 : assert( 2048 /*1.0/(1<<4) Q15*/ < mult( scale, scale ) );
65 : /* Exponent invScalingFactor: -16 = -(2*7 (scale) + 2 (8.0) */
66 6427 : hFdCngCom->invScalingFactor = shl( mult( scale, scale ), 1 );
67 : /* Exponent scalingFactor: -15 = -(2*7 (scale) + 2 (8.0) - 1 (1.0)) */
68 6427 : hFdCngCom->scalingFactor = div_s( 0x4000, hFdCngCom->invScalingFactor );
69 :
70 : /* Initialize the overlap-add */
71 6427 : set16_fx( hFdCngCom->timeDomainBuffer, 0, L_FRAME16k );
72 6427 : hFdCngCom->olapBufferAna = NULL;
73 6427 : set16_fx( hFdCngCom->olapBufferAna_fx, 0, FFTLEN );
74 :
75 6427 : move16();
76 6427 : set16_fx( hFdCngCom->olapBufferSynth, 0, FFTLEN );
77 6427 : hFdCngCom->olapBufferSynth2 = NULL;
78 6427 : move16();
79 :
80 : /* Initialize the comfort noise generation */
81 6427 : set32_fx( hFdCngCom->fftBuffer, 0, FFTLEN );
82 6427 : set32_fx( hFdCngCom->cngNoiseLevel, 0, FFTCLDFBLEN );
83 6427 : set16_fx( &hFdCngCom->cngNoiseLevelExp, 0, 1 );
84 :
85 : /* Initialize quantizer */
86 6427 : set32_fx( hFdCngCom->sidNoiseEst, 0, NPART );
87 6427 : set16_fx( &hFdCngCom->sidNoiseEstExp, 0, 1 );
88 6427 : set16_fx( hFdCngCom->A_cng, 0, M + 1 );
89 6427 : hFdCngCom->A_cng[0] = 4096 /*1.f Q12*/; /* 3Q12 */
90 6427 : move16();
91 :
92 : /* Set some counters and flags */
93 6427 : hFdCngCom->inactive_frame_counter = 0; /* Either SID or zero frames */
94 6427 : move16();
95 6427 : hFdCngCom->active_frame_counter = 0;
96 6427 : move16();
97 6427 : hFdCngCom->frame_type_previous = ACTIVE_FRAME;
98 6427 : move16();
99 6427 : hFdCngCom->flag_noisy_speech = 0;
100 6427 : move16();
101 6427 : hFdCngCom->likelihood_noisy_speech = 0;
102 6427 : move16();
103 6427 : hFdCngCom->numCoreBands = 0;
104 6427 : move16();
105 6427 : hFdCngCom->stopBand = 0;
106 6427 : move16();
107 6427 : hFdCngCom->startBand = 0;
108 6427 : move16();
109 6427 : hFdCngCom->stopFFTbin = 0;
110 6427 : move16();
111 6427 : hFdCngCom->frameSize = 0;
112 6427 : move16();
113 6427 : hFdCngCom->fftlen = 0;
114 6427 : move16();
115 6427 : hFdCngCom->seed = 0;
116 6427 : move16();
117 6427 : hFdCngCom->seed2 = 1;
118 6427 : move16();
119 6427 : hFdCngCom->seed3 = 2;
120 6427 : move16();
121 6427 : hFdCngCom->CngBitrate = -1;
122 6427 : move16();
123 :
124 : /* Initialize noise estimation algorithm */
125 6427 : set32_fx( hFdCngCom->periodog, 0, PERIODOGLEN );
126 6427 : mhvals( MSNUMSUBFR * MSSUBFRLEN, &( hFdCngCom->msM_win ) );
127 6427 : mhvals( MSSUBFRLEN, &( hFdCngCom->msM_subwin ) );
128 6427 : set32_fx( hFdCngCom->msPeriodogSum, 0, 2 );
129 6427 : hFdCngCom->msPeriodogSum_exp[0] = 0;
130 6427 : move16();
131 6427 : hFdCngCom->msPeriodogSum_exp[1] = 0;
132 6427 : move16();
133 6427 : set32_fx( hFdCngCom->msPsdSum, 0, 2 );
134 6427 : set16_fx( hFdCngCom->msSlope, 0, 2 );
135 6427 : set32_fx( hFdCngCom->msQeqInvAv, 0, 2 );
136 6427 : hFdCngCom->msQeqInvAv_exp[0] = 0;
137 6427 : move16();
138 6427 : hFdCngCom->msQeqInvAv_exp[1] = 0;
139 6427 : move16();
140 6427 : hFdCngCom->msFrCnt_init_counter = 0;
141 6427 : move16();
142 6427 : hFdCngCom->msFrCnt_init_thresh = 1;
143 6427 : move16();
144 6427 : hFdCngCom->init_old = 0;
145 6427 : move16();
146 6427 : hFdCngCom->offsetflag = 0;
147 6427 : move16();
148 6427 : hFdCngCom->msFrCnt = MSSUBFRLEN;
149 6427 : move16();
150 6427 : hFdCngCom->msMinBufferPtr = 0;
151 6427 : move16();
152 6427 : hFdCngCom->msAlphaCor[0] = 644245120l /*0.3f Q31*/;
153 6427 : move16();
154 6427 : hFdCngCom->msAlphaCor[1] = 644245120l /*0.3f Q31*/;
155 6427 : move16();
156 6427 : set16_fx( hFdCngCom->psize, 0, NPART );
157 : /* Initialize exponents */
158 6427 : hFdCngCom->exp_cldfb_periodog = 0;
159 6427 : move16();
160 :
161 6427 : hFdCngCom->coherence_fx = 16384; /* 0.5 in Q15 */
162 6427 : move16();
163 :
164 6427 : set32_fx( hFdCngCom->olapBufferSynth_fx, 0, FFTLEN );
165 6427 : set32_fx( hFdCngCom->olapBufferSynth2_fx, 0, FFTLEN );
166 6427 : set32_fx( hFdCngCom->exc_cng_32fx, 0, L_FRAME16k );
167 6427 : set16_fx( hFdCngCom->exc_cng, 0, L_FRAME16k );
168 :
169 6427 : return;
170 : }
171 3 : void initFdCngCom( HANDLE_FD_CNG_COM hFdCngCom, Word16 scale )
172 : {
173 : /* Calculate CLDFB scaling factor */
174 : /* shl(i_mult2(scale, scale), 3) does not fit in 16 bit */
175 : /*hFdCngCom->scalingFactor = div_s(1, shl(i_mult2(scale, scale), 3));*/
176 3 : assert( 2048 /*1.0/(1<<4) Q15*/ < mult( scale, scale ) );
177 : /* Exponent invScalingFactor: -16 = -(2*7 (scale) + 2 (8.0) */
178 3 : hFdCngCom->invScalingFactor = shl( mult( scale, scale ), 1 );
179 3 : move16();
180 : /* Exponent scalingFactor: -15 = -(2*7 (scale) + 2 (8.0) - 1 (1.0)) */
181 3 : hFdCngCom->scalingFactor = div_s( 0x4000, hFdCngCom->invScalingFactor );
182 3 : move16();
183 :
184 : /* Initialize the overlap-add */
185 3 : set16_fx( hFdCngCom->timeDomainBuffer, 0, L_FRAME16k );
186 3 : hFdCngCom->olapBufferAna = NULL;
187 3 : set16_fx( hFdCngCom->olapBufferAna_fx, 0, FFTLEN );
188 :
189 3 : set16_fx( hFdCngCom->olapBufferSynth, 0, FFTLEN );
190 3 : hFdCngCom->olapBufferSynth2 = NULL;
191 :
192 : /* Initialize the comfort noise generation */
193 3 : set32_fx( hFdCngCom->fftBuffer, 0, FFTLEN );
194 3 : set32_fx( hFdCngCom->cngNoiseLevel, 0, FFTCLDFBLEN );
195 3 : set16_fx( &hFdCngCom->cngNoiseLevelExp, 0, 1 );
196 :
197 : /* Initialize quantizer */
198 3 : set32_fx( hFdCngCom->sidNoiseEst, 0, NPART );
199 3 : set16_fx( &hFdCngCom->sidNoiseEstExp, 0, 1 );
200 3 : set16_fx( hFdCngCom->A_cng, 0, M + 1 );
201 3 : hFdCngCom->A_cng[0] = 4096 /*1.f Q12*/; /* 3Q12 */
202 3 : move16();
203 :
204 : /* Set some counters and flags */
205 3 : hFdCngCom->inactive_frame_counter = 0; /* Either SID or zero frames */
206 3 : move16();
207 3 : hFdCngCom->active_frame_counter = 0;
208 3 : move16();
209 3 : hFdCngCom->frame_type_previous = ACTIVE_FRAME;
210 3 : move16();
211 3 : hFdCngCom->flag_noisy_speech = 0;
212 3 : move16();
213 3 : hFdCngCom->likelihood_noisy_speech = 0;
214 3 : move16();
215 3 : hFdCngCom->numCoreBands = 0;
216 3 : move16();
217 3 : hFdCngCom->stopBand = 0;
218 3 : move16();
219 3 : hFdCngCom->startBand = 0;
220 3 : move16();
221 3 : hFdCngCom->stopFFTbin = 0;
222 3 : move16();
223 3 : hFdCngCom->frameSize = 0;
224 3 : move16();
225 3 : hFdCngCom->fftlen = 0;
226 3 : move16();
227 3 : hFdCngCom->seed = 0;
228 3 : move16();
229 : // hFdCngCom->seed2 = 1;
230 : // move16();
231 : // hFdCngCom->seed3 = 2;
232 : // move16();
233 3 : hFdCngCom->CngBitrate = -1;
234 3 : move16();
235 :
236 : /* Initialize noise estimation algorithm */
237 3 : set32_fx( hFdCngCom->periodog, 0, PERIODOGLEN );
238 3 : mhvals( MSNUMSUBFR * MSSUBFRLEN, &( hFdCngCom->msM_win ) );
239 3 : mhvals( MSSUBFRLEN, &( hFdCngCom->msM_subwin ) );
240 3 : set32_fx( hFdCngCom->msPeriodogSum, 0, 2 );
241 3 : hFdCngCom->msPeriodogSum_exp[0] = 0;
242 3 : move16();
243 3 : hFdCngCom->msPeriodogSum_exp[1] = 0;
244 3 : move16();
245 3 : set32_fx( hFdCngCom->msPsdSum, 0, 2 );
246 3 : set16_fx( hFdCngCom->msSlope, 0, 2 );
247 3 : set32_fx( hFdCngCom->msQeqInvAv, 0, 2 );
248 3 : hFdCngCom->msQeqInvAv_exp[0] = 0;
249 3 : move16();
250 3 : hFdCngCom->msQeqInvAv_exp[1] = 0;
251 3 : move16();
252 3 : hFdCngCom->msFrCnt_init_counter = 0;
253 3 : move16();
254 3 : hFdCngCom->msFrCnt_init_thresh = 1;
255 3 : move16();
256 3 : hFdCngCom->init_old = 0;
257 3 : move16();
258 3 : hFdCngCom->offsetflag = 0;
259 3 : move16();
260 3 : hFdCngCom->msFrCnt = MSSUBFRLEN;
261 3 : move16();
262 3 : hFdCngCom->msMinBufferPtr = 0;
263 3 : move16();
264 3 : hFdCngCom->msAlphaCor[0] = 644245120l /*0.3f Q31*/;
265 3 : move16();
266 3 : hFdCngCom->msAlphaCor[1] = 644245120l /*0.3f Q31*/;
267 3 : move16();
268 3 : set16_fx( hFdCngCom->psize, 0, NPART );
269 : /* Initialize exponents */
270 3 : hFdCngCom->exp_cldfb_periodog = 0;
271 3 : move16();
272 :
273 3 : return;
274 : }
275 :
276 : /*-------------------------------------------------------------------
277 : * deleteFdCngCom()
278 : *
279 : * Delete an instance of type FD_CNG_COM
280 : *-------------------------------------------------------------------*/
281 6430 : void deleteFdCngCom_fx( HANDLE_FD_CNG_COM *hFdCngCom ) /* i/o: Contains the variables related to the CLDFB-based CNG process */
282 : {
283 : HANDLE_FD_CNG_COM hsCom;
284 6430 : hsCom = *hFdCngCom;
285 6430 : move16();
286 6430 : IF( hsCom != NULL )
287 : {
288 6430 : free( hsCom );
289 6430 : *hFdCngCom = NULL;
290 6430 : move16();
291 : }
292 6430 : }
293 :
294 : /*-------------------------------------------------------------------
295 : * initPartitions()
296 : *
297 : * Initialize the spectral partitioning
298 : *-------------------------------------------------------------------*/
299 1060682 : void initPartitions(
300 : const Word16 *part_in,
301 : const Word16 npart_in,
302 : const Word16 startBand,
303 : const Word16 stopBand,
304 : Word16 *part_out,
305 : Word16 *npart_out,
306 : Word16 *midband,
307 : Word16 *psize,
308 : Word16 *psize_norm,
309 : Word16 *psize_norm_exp,
310 : Word16 *psize_inv,
311 : const Word16 stopBandFR )
312 : {
313 : Word16 i, j, len_out, tmp16;
314 :
315 :
316 1060682 : IF( part_in != NULL )
317 : {
318 1060682 : len_out = 0;
319 1060682 : move16();
320 1060682 : IF( GT_16( stopBandFR, startBand ) )
321 : {
322 530341 : len_out = sub( stopBandFR, startBand );
323 20683299 : FOR( i = 0; i < len_out; i++ )
324 : {
325 20152958 : part_out[i] = i;
326 20152958 : move16();
327 : }
328 : }
329 34469360 : FOR( j = 0; j < npart_in; j++ )
330 : {
331 33408678 : IF( GE_16( part_in[j], stopBand ) )
332 : {
333 0 : BREAK;
334 : }
335 33408678 : tmp16 = sub( part_in[j], startBand );
336 33408678 : test();
337 33408678 : if ( GE_16( part_in[j], stopBandFR ) && tmp16 >= 0 )
338 : {
339 24923222 : part_out[len_out++] = tmp16;
340 24923222 : move16();
341 : }
342 : }
343 : }
344 : ELSE
345 : {
346 0 : len_out = sub( stopBand, startBand );
347 0 : FOR( i = 0; i < len_out; i++ )
348 : {
349 0 : part_out[i] = i;
350 0 : move16();
351 : }
352 : }
353 :
354 1060682 : *npart_out = len_out;
355 1060682 : move16();
356 1060682 : getmidbands( part_out, len_out, midband, psize, psize_norm, psize_norm_exp, psize_inv );
357 1060682 : }
358 :
359 :
360 : /*-------------------------------------------------------------------
361 : * compress_range()
362 : *
363 : * Apply some dynamic range compression based on the log
364 : *-------------------------------------------------------------------*/
365 42153 : void compress_range(
366 : Word32 *in, // Q(31 - in_exp)
367 : Word16 in_exp,
368 : Word16 *out,
369 : const Word16 len )
370 : {
371 : Word16 i;
372 : Word32 in_s;
373 : Word32 one_s;
374 : Word32 in_exp32;
375 : Word32 L_tmp;
376 :
377 :
378 : /* out = log2( 1 + in ) */
379 42153 : IF( in_exp >= 0 )
380 : {
381 32989 : one_s = L_shr( 1073741824l /*0.5 Q31*/, in_exp ); // Q(31 - in_exp)
382 32989 : in_exp32 = L_shl( L_deposit_l( add( in_exp, 1 ) ), WORD32_BITS - 1 - LD_DATA_SCALE );
383 2058735 : FOR( i = 0; i < len; i++ )
384 : {
385 2025746 : in_s = L_add( L_shr( in[i], 1 ), one_s ); // Q(31 - in_exp)
386 2025746 : L_tmp = L_add( BASOP_Util_Log2( in_s ), in_exp32 );
387 2025746 : if ( in_s == 0 )
388 : {
389 42580 : out[i] = 0;
390 42580 : move16();
391 : }
392 2025746 : if ( in_s != 0 )
393 : {
394 1983166 : out[i] = extract_h( L_tmp );
395 1983166 : move16();
396 : }
397 2025746 : if ( out[i] == 0 )
398 : {
399 86213 : out[i] = 1;
400 86213 : move16();
401 : }
402 : }
403 : }
404 : ELSE
405 : {
406 9164 : in_exp = sub( in_exp, 1 );
407 9164 : in_exp32 = L_shl( L_deposit_l( 1 ), WORD32_BITS - 1 - LD_DATA_SCALE );
408 574017 : FOR( i = 0; i < len; i++ )
409 : {
410 564853 : L_tmp = L_add( BASOP_Util_Log2( L_add( L_shl( in[i], in_exp ), 1073741824l /*0.5 Q31*/ ) ), in_exp32 );
411 564853 : if ( in[i] == 0 )
412 : {
413 7027 : out[i] = 0;
414 7027 : move16();
415 : }
416 564853 : if ( in[i] != 0 )
417 : {
418 557826 : out[i] = extract_h( L_tmp );
419 557826 : move16();
420 : }
421 564853 : if ( out[i] == 0 )
422 : {
423 113771 : out[i] = 1;
424 113771 : move16();
425 : }
426 : }
427 : }
428 42153 : }
429 :
430 : /*-------------------------------------------------------------------
431 : * expand_range()
432 : *
433 : * Apply some dynamic range expansion to undo the compression
434 : *-------------------------------------------------------------------*/
435 42153 : void expand_range(
436 : Word16 *in, // Q15
437 : Word32 *out, // Q(31 - out_exp)
438 : Word16 *out_exp,
439 : const Word16 len )
440 : {
441 : Word16 i;
442 : Word16 s;
443 : Word16 tmp;
444 : Word32 one_s, tmp32;
445 : Word16 maxVal;
446 : Word16 maxVal2;
447 :
448 :
449 42153 : maxVal = 0;
450 42153 : move16();
451 2632752 : FOR( i = 0; i < len; i++ )
452 : {
453 2590599 : maxVal = s_max( maxVal, in[i] ); // Q15
454 : }
455 :
456 42153 : maxVal2 = maxVal;
457 42153 : move16();
458 42153 : s = 0;
459 42153 : move16();
460 438726 : WHILE( maxVal >= 0 )
461 : {
462 396573 : maxVal = sub( maxVal, 512 /*0.015625 Q15*/ ); // Q15
463 396573 : s = add( s, 1 );
464 : }
465 42153 : tmp = sub( maxVal2, maxVal );
466 :
467 42153 : one_s = L_shr( 1073741824l /*0.5 Q31*/, sub( s, 1 ) ); // Q(31 - sub( s, 1 ))
468 42153 : tmp32 = L_shr( 726941l /*0.0003385080526823181 Q31*/, s );
469 : /* out = (2^(in) - 1) */
470 2632752 : FOR( i = 0; i < len; i++ )
471 : {
472 2590599 : out[i] = L_sub( BASOP_Util_InvLog2( L_deposit_h( sub( in[i], tmp ) ) ), one_s );
473 2590599 : move32();
474 2590599 : if ( out[i] == 0 )
475 : {
476 34981 : out[i] = tmp32;
477 34981 : move32();
478 : }
479 : }
480 42153 : *out_exp = s;
481 42153 : move16();
482 42153 : }
483 :
484 : /*-------------------------------------------------------------------
485 : * expand_range_var_e()
486 : *
487 : * Apply some dynamic range expansion to undo the compression. Similar
488 : * to expand_range but has variable input exponent.
489 : *
490 : *-------------------------------------------------------------------*/
491 0 : void expand_range_var_exp(
492 : Word16 *in,
493 : Word16 in_exp,
494 : Word32 *out,
495 : Word16 *out_exp,
496 : const Word16 len )
497 : {
498 : Word16 i, tmp_e;
499 : Word32 tmp32;
500 : Word16 maxVal, maxOutExp;
501 :
502 0 : const Word32 low_lim = 726940; /* 0.0003385080526823181f in Q31 */
503 0 : move32();
504 :
505 : /* Find max possible output exponent. */
506 0 : maxVal = 0;
507 0 : move16();
508 0 : FOR( i = 0; i < len; i++ )
509 : {
510 0 : maxVal = s_max( maxVal, in[i] );
511 : }
512 0 : tmp32 = BASOP_util_Pow2( L_deposit_h( maxVal ), in_exp, &tmp_e );
513 0 : maxOutExp = tmp_e;
514 0 : move16();
515 :
516 : /* out = (2^(in) - 1) */
517 0 : FOR( i = 0; i < len; i++ )
518 : {
519 0 : tmp32 = BASOP_util_Pow2( L_deposit_h( in[i] ), in_exp, &tmp_e ); // 2^x
520 0 : tmp32 = L_sub( tmp32, L_shl( 1, sub( 31, tmp_e ) ) ); // 2^x - 1
521 0 : tmp32 = L_shr( tmp32, sub( maxOutExp, tmp_e ) ); // make exp same as maxExpOut
522 :
523 0 : out[i] = tmp32;
524 0 : move32();
525 :
526 0 : Word32 tmp_low_lim = L_shr( low_lim, maxOutExp );
527 0 : if ( LT_32( out[i], tmp_low_lim ) )
528 : {
529 0 : out[i] = tmp_low_lim;
530 0 : move32();
531 : }
532 : }
533 :
534 0 : *out_exp = maxOutExp;
535 0 : move16();
536 :
537 0 : return;
538 : }
539 :
540 : /*-------------------------------------------------------------------
541 : * minimum_statistics()
542 : *
543 : * Noise estimation using Minimum Statistics (MS)
544 : *-------------------------------------------------------------------*/
545 :
546 15 : void minimum_statistics(
547 : Word16 len, /* i : Total number of partitions (CLDFB or FFT) */
548 : Word16 lenFFT, /* i : Number of FFT partitions */
549 : Word16 *psize, /* i : Partition sizes, fractional */
550 : Word16 *msPeriodog, /* i : Periodogram (energies) */
551 : Word16 *msNoiseFloor, /* i/o: Noise floors (energies) */
552 : Word16 *msNoiseEst, /* i/o: Noise estimates (energies) */
553 : Word32 *msAlpha, /* i/o: Forgetting factors */
554 : Word16 *msPsd, /* i/o: Power Spectral Density (smoothed periodogram => energies) */
555 : Word16 *msPsdFirstMoment, /* i/o: PSD statistics of 1st order (energy means) */
556 : Word32 *msPsdSecondMoment, /* i/o: PSD statistics of 2nd order (energy variances) */
557 : Word32 *msMinBuf, /* i/o: Buffer of minima (energies) */
558 : Word32 *msBminWin, /* o : Bias correction factors */
559 : Word32 *msBminSubWin, /* o : Bias correction factors */
560 : Word32 *msCurrentMin, /* i/o: Local minima (energies) */
561 : Word32 *msCurrentMinOut, /* i/o: Local minima (energies) */
562 : Word32 *msCurrentMinSubWindow, /* i/o: Local minima (energies) */
563 : Word16 *msLocalMinFlag, /* i : Binary flag */
564 : Word16 *msNewMinFlag, /* i : Binary flag */
565 : Word16 *msPeriodogBuf, /* i/o: Buffer of periodograms (energies) */
566 : Word16 *msPeriodogBufPtr, /* i/o: Counter */
567 : HANDLE_FD_CNG_COM hFdCngCom /* i/o: FD_CNG structure containing all buffers and variables */
568 : )
569 : {
570 : Word16 i, j, k, s, s1, s2, s3;
571 : Word16 len2;
572 : Word16 current_len;
573 : Word16 start, stop, cnt;
574 : Word16 totsize;
575 : Word16 inv_totsize;
576 :
577 : Word32 tmp, tmp0, tmp1;
578 : Word32 scalar, scalar2, scalar3;
579 : Word32 snr;
580 : Word32 msAlphaHatMin2;
581 : Word32 BminCorr;
582 : Word32 QeqInvAv;
583 : Word32 *ptr;
584 : Word32 *msPsdSum;
585 : Word32 *msPeriodogSum;
586 :
587 : Word16 tmp16, tmp16_1;
588 : Word16 beta;
589 : Word16 slope;
590 : Word16 QeqInv;
591 : Word16 scalar16;
592 : Word16 scalar216;
593 : Word16 scalar316;
594 : Word16 msM_win;
595 : Word16 msM_subwin;
596 : Word16 msAlphaCorAlpha;
597 : Word16 msAlphaCorAlpha2;
598 : Word16 msPeriodogSum16;
599 : Word16 msNoiseFloor16;
600 : #ifdef BASOP_NOGLOB_DECLARE_LOCAL
601 15 : Flag Overflow = 0;
602 15 : move32();
603 : #endif
604 :
605 :
606 15 : len2 = i_mult( MSNUMSUBFR, len );
607 :
608 15 : msM_win = hFdCngCom->msM_win;
609 15 : move16();
610 15 : msM_subwin = hFdCngCom->msM_subwin;
611 15 : move16();
612 :
613 15 : msAlphaCorAlpha = MSALPHACORALPHA;
614 15 : move16();
615 15 : msAlphaCorAlpha2 = MSALPHACORALPHA2;
616 15 : move16();
617 :
618 15 : msPsdSum = hFdCngCom->msPsdSum;
619 15 : msPeriodogSum = hFdCngCom->msPeriodogSum;
620 :
621 : /* No minimum statistics at initialization */
622 15 : IF( LT_16( hFdCngCom->msFrCnt_init_counter, hFdCngCom->msFrCnt_init_thresh ) )
623 : {
624 6 : Copy( msPeriodog, msPsd, len ); /* 6Q9 */
625 6 : Copy( msPeriodog, msNoiseFloor, len ); /* 6Q9 */
626 6 : Copy( msPeriodog, msNoiseEst, len ); /* 6Q9 */
627 6 : Copy( msPeriodog, msPsdFirstMoment, len ); /* 6Q9 */
628 :
629 6 : set32_fx( msPsdSecondMoment, 0l /*0.0 Q31*/, len );
630 6 : msPeriodogSum[0] = dotp_s_fx( msPeriodog, psize, lenFFT, CNG_HS );
631 6 : move32();
632 6 : msPsdSum[0] = msPeriodogSum[0]; /* 16Q15 */
633 6 : move32();
634 :
635 6 : IF( LT_16( lenFFT, len ) )
636 : {
637 0 : msPeriodogSum[1] = dotp_s_fx( msPeriodog + lenFFT, psize + lenFFT, sub( len, lenFFT ), CNG_HS );
638 0 : move32();
639 0 : msPsdSum[1] = msPeriodogSum[1]; /* 16Q15 */
640 0 : move32();
641 : }
642 :
643 : /* Increment frame counter at initialization */
644 : /* Some frames are sometimes zero at initialization => ignore them */
645 6 : IF( LT_16( msPeriodog[0], hFdCngCom->init_old ) )
646 : {
647 2 : set32_fx( msCurrentMinOut, 2147483647l /*1.0 Q31*/, len ); /* 16Q15 */
648 2 : set32_fx( msCurrentMin, 2147483647l /*1.0 Q31*/, len ); /* 16Q15 */
649 2 : set32_fx( msMinBuf, 2147483647l /*1.0 Q31*/, len2 ); /* 16Q15 */
650 2 : set32_fx( msCurrentMinSubWindow, 2147483647l /*1.0 Q31*/, len ); /* 16Q15 */
651 :
652 2 : hFdCngCom->msFrCnt_init_counter = add( hFdCngCom->msFrCnt_init_counter, 1 );
653 2 : move16();
654 : }
655 6 : hFdCngCom->init_old = msPeriodog[0]; /* 6Q9 */
656 6 : move16();
657 : }
658 : ELSE
659 : {
660 : /* Consider the FFT and CLDFB bands separately
661 : - first iteration for FFT bins,
662 : - second one for CLDFB bands in SWB mode */
663 9 : cnt = 0;
664 9 : move16();
665 9 : start = 0;
666 9 : move16();
667 9 : stop = lenFFT;
668 9 : move16();
669 9 : totsize = sub( hFdCngCom->stopFFTbin, hFdCngCom->startBand );
670 18 : WHILE( GT_16( stop, start ) )
671 : {
672 9 : current_len = sub( stop, start );
673 :
674 : /* Compute smoothed correction factor for PSD smoothing */
675 :
676 : /* msPeriodogSum[cnt] with format 16Q15 */
677 9 : msPeriodogSum[cnt] = dotp_s_fx( msPeriodog + start, psize + start, current_len, CNG_HS );
678 9 : move32();
679 :
680 9 : IF( msPeriodogSum[cnt] == 0 )
681 : {
682 0 : hFdCngCom->msAlphaCor[cnt] = Mpy_32_16_1( hFdCngCom->msAlphaCor[cnt], msAlphaCorAlpha );
683 0 : move32();
684 : }
685 : ELSE
686 : {
687 : /* calculate scalar with normalized msPeriodogSum[cnt], exponent -2*s1 */
688 9 : s1 = norm_l( msPeriodogSum[cnt] );
689 9 : msPeriodogSum16 = round_fx_sat( L_shl_sat( msPeriodogSum[cnt], s1 ) );
690 9 : scalar = L_mult( msPeriodogSum16, msPeriodogSum16 );
691 :
692 : /* calculate difference, both elements in 16Q15 format, use absolute value
693 : to avoid -1.0 x -1.0 multiplications later */
694 9 : scalar2 = L_abs( L_sub( msPsdSum[cnt], msPeriodogSum[cnt] ) );
695 :
696 9 : s2 = WORD32_BITS - 1;
697 9 : move16();
698 9 : if ( scalar2 != 0 )
699 : {
700 : /* use absolute value to avoid -1.0 x -1.0 multiplications */
701 9 : s2 = norm_l( scalar2 );
702 : }
703 9 : scalar216 = round_fx_o( L_shl_o( scalar2, s2, &Overflow ), &Overflow );
704 9 : scalar2 = L_mult( scalar216, scalar216 );
705 :
706 : /* find common exponent */
707 9 : tmp16_1 = sub( s1, s2 );
708 9 : tmp16 = s_min( shl( abs_s( tmp16_1 ), 1 ), WORD32_BITS - 1 );
709 9 : if ( tmp16_1 < 0 )
710 : {
711 8 : scalar2 = L_shr( scalar2, tmp16 );
712 : }
713 9 : if ( tmp16_1 > 0 )
714 : {
715 0 : scalar = L_shr( scalar, tmp16 );
716 : }
717 :
718 :
719 : /* add scalar and scalar2, avoid overflows */
720 9 : scalar = L_shr( scalar, 1 );
721 9 : scalar2 = L_shr( scalar2, 1 );
722 9 : scalar3 = L_add( scalar, scalar2 );
723 :
724 : /* calculate division */
725 9 : scalar16 = BASOP_Util_Divide3232_uu_1616_Scale( scalar, scalar3, &s3 );
726 9 : s3 = s_max( s3, -( WORD16_BITS - 1 ) );
727 9 : scalar16 = shl( scalar16, s3 );
728 9 : scalar16 = s_max( scalar16, MSALPHACORMAX );
729 :
730 9 : hFdCngCom->msAlphaCor[cnt] = L_add( Mpy_32_16_1( hFdCngCom->msAlphaCor[cnt], msAlphaCorAlpha ),
731 : L_mult( scalar16, msAlphaCorAlpha2 ) );
732 9 : move32();
733 : }
734 :
735 : /* Compute SNR */
736 :
737 : /* msPeriodogSum[cnt] with format 16Q15 */
738 9 : snr = dotp_s_fx( msNoiseFloor + start, psize + start, current_len, CNG_HS );
739 :
740 9 : IF( GT_32( L_shr( Mpy_32_16_1( msPsdSum[cnt], 18431 /*0.56246299817 Q15*/ ), 13 ), snr ) )
741 : {
742 0 : tmp0 = BASOP_Util_Log2( msPsdSum[cnt] );
743 0 : tmp1 = BASOP_Util_Log2( snr );
744 0 : tmp1 = L_sub( tmp0, tmp1 );
745 0 : tmp1 = Mpy_32_16_1( tmp1, MSSNREXP );
746 0 : msAlphaHatMin2 = BASOP_Util_InvLog2( tmp1 );
747 : }
748 : ELSE
749 : {
750 9 : msAlphaHatMin2 = MSALPHAHATMIN;
751 9 : move32();
752 : }
753 9 : scalar = Mpy_32_16_1( hFdCngCom->msAlphaCor[cnt], MSALPHAMAX );
754 :
755 560 : FOR( j = start; j < stop; j++ )
756 : {
757 : /* Compute optimal smoothing parameter for PSD estimation */
758 551 : test();
759 551 : IF( ( scalar == 0 ) || ( msNoiseFloor[j] == 0 ) )
760 : {
761 0 : msAlpha[j] = msAlphaHatMin2;
762 0 : move32();
763 : }
764 : ELSE
765 : {
766 : /* calculate scalar2 with normalized msNoiseFloor[j], exponent -2*s1 */
767 551 : s1 = WORD16_BITS - 1;
768 551 : move16();
769 551 : if ( msNoiseFloor[j] != 0 )
770 : {
771 551 : s1 = norm_s( msNoiseFloor[j] );
772 : }
773 551 : msNoiseFloor16 = shl( msNoiseFloor[j], s1 );
774 551 : scalar2 = L_mult( msNoiseFloor16, msNoiseFloor16 );
775 :
776 : /* calculate difference, both elements in 6Q9 format, use absolute value
777 : to avoid -1.0 x -1.0 multiplications later */
778 551 : scalar316 = abs_s( sub( msPsd[j], msNoiseFloor[j] ) );
779 :
780 551 : s2 = WORD16_BITS - 1;
781 551 : move16();
782 551 : if ( scalar316 != 0 )
783 : {
784 : /* use absolute value to avoid -1.0 x -1.0 multiplications */
785 407 : s2 = norm_s( scalar316 );
786 : }
787 551 : scalar316 = shl( scalar316, s2 );
788 551 : scalar3 = L_mult( scalar316, scalar316 );
789 :
790 : /* find common exponent */
791 551 : tmp16_1 = sub( s1, s2 );
792 551 : tmp16 = s_min( shl( abs_s( tmp16_1 ), 1 ), WORD32_BITS - 1 );
793 551 : if ( tmp16_1 < 0 )
794 : {
795 443 : scalar3 = L_shr( scalar3, tmp16 );
796 : }
797 551 : if ( tmp16_1 > 0 )
798 : {
799 5 : scalar2 = L_shr( scalar2, tmp16 );
800 : }
801 :
802 : /* add scalar2 and scalar3, avoid overflows */
803 551 : scalar2 = L_shr( scalar2, 1 );
804 551 : scalar3 = L_shr( scalar3, 1 );
805 551 : scalar3 = L_add( scalar2, scalar3 );
806 :
807 : /* calculate division */
808 551 : tmp16 = BASOP_Util_Divide3232_uu_1616_Scale( scalar2, scalar3, &s3 );
809 551 : scalar2 = Mpy_32_16_1( scalar, tmp16 );
810 551 : s3 = s_max( s3, -( WORD32_BITS - 1 ) );
811 551 : scalar2 = L_shl( scalar2, s3 );
812 551 : msAlpha[j] = L_max( scalar2, msAlphaHatMin2 );
813 551 : move32();
814 : }
815 :
816 : /* Compute the PSD (smoothed periodogram) in each band */
817 551 : msPsd[j] = round_fx( L_add( Mpy_32_16_1( msAlpha[j], msPsd[j] ),
818 551 : Mpy_32_16_1( L_sub( 2147483647l /*1.0 Q31*/, msAlpha[j] ), msPeriodog[j] ) ) );
819 551 : move16();
820 : }
821 9 : msPsdSum[cnt] = dotp_s_fx( msPsd + start, psize + start, current_len, CNG_HS );
822 9 : move32();
823 :
824 9 : QeqInvAv = 0l /*0.0 Q31*/;
825 9 : move32();
826 :
827 : /* scalar: 7Q24 format */
828 9 : tmp = 1191182336l /*(float)(MSNUMSUBFR*MSSUBFRLEN)-1.0 Q24*/;
829 9 : move32();
830 9 : scalar = L_sub( tmp, L_mult( round_fx( tmp ), msM_win ) );
831 :
832 : /* scalar2: 4Q27 format */
833 9 : tmp = 1476395008l /*(float)MSSUBFRLEN-1.0 Q27*/;
834 9 : move32();
835 9 : scalar2 = L_sub( tmp, L_mult( round_fx( tmp ), msM_subwin ) );
836 :
837 560 : FOR( j = start; j < stop; j++ )
838 : {
839 : /* Compute variance of PSD */
840 551 : tmp = L_min( msAlpha[j], MSBETAMAX_SQRT );
841 :
842 551 : s1 = WORD32_BITS - 1;
843 551 : move16();
844 551 : if ( tmp != 0 )
845 : {
846 551 : s1 = norm_l( tmp );
847 : }
848 551 : s2 = shl( s1, 1 );
849 :
850 551 : s2 = s_min( s2, WORD32_BITS - 1 );
851 :
852 : /* beta: scaled by s2 */
853 551 : tmp16 = round_fx_o( L_shl_o( tmp, s1, &Overflow ), &Overflow );
854 551 : beta = mult_r( tmp16, tmp16 );
855 :
856 : /* scalar3: scaled by s3 */
857 551 : scalar3 = L_sub( L_deposit_l( msPsd[j] ), L_deposit_l( msPsdFirstMoment[j] ) );
858 551 : s3 = norm_l( scalar3 );
859 551 : scalar3 = L_shl( scalar3, s3 );
860 :
861 : /* msPsdFirstMoment: 6Q9 */
862 551 : tmp = L_msu( L_mult( beta, msPsdFirstMoment[j] ), beta, msPsd[j] );
863 551 : msPsdFirstMoment[j] = add( round_fx( L_shr( tmp, s2 ) ), msPsd[j] );
864 551 : move16();
865 :
866 : /* msPsdSecondMoment: 12Q19 */
867 551 : tmp0 = L_shr( Mpy_32_16_r( msPsdSecondMoment[j], beta ), s2 );
868 551 : tmp1 = Mpy_32_32( scalar3, scalar3 );
869 551 : tmp1 = L_shr( L_sub( tmp1, L_shr( Mpy_32_16_r( tmp1, beta ), s2 ) ), sub( shl( s3, 1 ), 32 ) );
870 551 : msPsdSecondMoment[j] = L_add( tmp0, tmp1 );
871 551 : move32();
872 :
873 : /* Compute inverse of amount of degrees of freedom */
874 551 : QeqInv = MSQEQINVMAX;
875 551 : move16();
876 :
877 551 : IF( msNoiseFloor[j] != 0 /*0.0 Q15*/ )
878 : {
879 551 : tmp = L_mult( msNoiseFloor[j], msNoiseFloor[j] );
880 551 : tmp16 = BASOP_Util_Divide3232_uu_1616_Scale( msPsdSecondMoment[j], tmp, &s );
881 : /* consider factor of 2 */
882 551 : s = s_min( s_max( sub( s, 1 ), -( WORD16_BITS - 1 ) ), ( WORD16_BITS - 1 ) );
883 551 : if ( s < 0 )
884 : {
885 496 : QeqInv = shl( tmp16, s );
886 : }
887 551 : QeqInv = s_min( QeqInv, MSQEQINVMAX );
888 : }
889 551 : QeqInvAv = L_add( QeqInvAv, L_mult( QeqInv, psize[j] ) );
890 :
891 : /* Compute bias correction Bmin */
892 551 : tmp0 = Mpy_32_16_1( scalar, QeqInv );
893 551 : tmp1 = L_sub( 1073741824l /*0.5 Q31*/, L_mult( msM_win, QeqInv ) );
894 551 : tmp16 = BASOP_Util_Divide3232_uu_1616_Scale( tmp0, tmp1, &s );
895 551 : msBminWin[j] = L_add( 134217728l /*1.0 Q27*/, L_shl( L_deposit_h( tmp16 ), add( s, 7 - 4 ) ) );
896 551 : move32();
897 :
898 551 : tmp0 = Mpy_32_16_1( scalar2, QeqInv );
899 551 : tmp1 = L_sub( 1073741824l /*0.5 Q31*/, L_mult( msM_subwin, QeqInv ) );
900 551 : tmp16 = BASOP_Util_Divide3232_uu_1616_Scale( tmp0, tmp1, &s );
901 551 : msBminSubWin[j] = L_add( 134217728l /*1.0 Q27*/, L_shl( L_deposit_h( tmp16 ), s ) );
902 551 : move32();
903 : }
904 :
905 9 : inv_totsize = BASOP_Util_Divide3232_uu_1616_Scale( 1, totsize, &s );
906 9 : QeqInvAv = Mpy_32_16_1( QeqInvAv, inv_totsize );
907 9 : QeqInvAv = L_shl( QeqInvAv, s );
908 9 : hFdCngCom->msQeqInvAv[cnt] = QeqInvAv;
909 9 : move32();
910 :
911 : /* New minimum? */
912 :
913 : /* exponent QeqInvAv: CNG_S, exponent MSAV: (4>>1) */
914 9 : s = CNG_S + 2 * MSAV_EXP;
915 9 : move16();
916 9 : BminCorr = Mpy_32_16_1( Sqrt32( QeqInvAv, &s ), MSAV );
917 9 : BminCorr = L_shl( BminCorr, sub( s, 1 ) );
918 :
919 : /* exponent BminCorr: 1 */
920 9 : BminCorr = L_add( BminCorr, 1073741824l /*0.5 Q31*/ );
921 :
922 560 : FOR( j = start; j < stop; j++ )
923 : {
924 : /* exponent scalar: CNG_S+1 */
925 551 : scalar = Mpy_32_16_1( BminCorr, msPsd[j] );
926 :
927 : /* exponent scalar2: CNG_S+1+4 */
928 551 : scalar2 = Mpy_32_32( scalar, msBminWin[j] );
929 :
930 551 : msNewMinFlag[j] = 0;
931 551 : move16();
932 551 : IF( LT_32( scalar2, msCurrentMin[j] ) /*0.0 Q31*/ )
933 : {
934 317 : msNewMinFlag[j] = 1;
935 317 : move16();
936 : /* exponent msCurrentMin[j]: CNG_S+1+4 */
937 317 : msCurrentMin[j] = scalar2;
938 317 : move32();
939 : /* exponent msCurrentMinSubWindow[j]: CNG_S */
940 : BASOP_SATURATE_WARNING_OFF_EVS;
941 317 : msCurrentMinSubWindow[j] = L_shl_o( Mpy_32_32( scalar, msBminSubWin[j] ), 5, &Overflow );
942 317 : move32();
943 : BASOP_SATURATE_WARNING_ON_EVS;
944 : }
945 : }
946 :
947 : /* This is used later to identify local minima */
948 9 : IF( GE_16( hFdCngCom->msFrCnt, MSSUBFRLEN ) )
949 : {
950 8 : FOR( i = 0; i < 3; i++ )
951 : {
952 6 : IF( LT_32( hFdCngCom->msQeqInvAv[cnt], L_shr( L_deposit_h( msQeqInvAv_thresh[i] ), CNG_S ) ) /*0.0 Q31*/ )
953 : {
954 0 : BREAK;
955 : }
956 : }
957 : /* format 1Q30 */
958 2 : hFdCngCom->msSlope[cnt] = msNoiseSlopeMax[i];
959 2 : move32();
960 : }
961 :
962 : /* Consider the FFT and CLDFB bands separately */
963 9 : start = stop;
964 9 : move16();
965 9 : stop = len;
966 9 : move16();
967 9 : totsize = sub( hFdCngCom->stopBand, hFdCngCom->stopFFTbin );
968 9 : cnt = add( cnt, 1 );
969 : } /*while (stop > start)*/
970 :
971 : /* Update minimum between sub windows */
972 9 : test();
973 9 : IF( GT_16( hFdCngCom->msFrCnt, 1 ) && LT_16( hFdCngCom->msFrCnt, MSSUBFRLEN ) )
974 : {
975 310 : FOR( j = 0; j < len; j++ )
976 : {
977 305 : if ( msNewMinFlag[j] > 0 )
978 : {
979 71 : msLocalMinFlag[j] = 1;
980 71 : move16();
981 : }
982 305 : if ( LT_32( msCurrentMinSubWindow[j], msCurrentMinOut[j] ) /*0.0 Q31*/ )
983 : {
984 : /* msCurrentMinOut[j] scaled with CNG_S */
985 85 : msCurrentMinOut[j] = msCurrentMinSubWindow[j];
986 85 : move32();
987 : }
988 : }
989 : /* Get the current noise floor */
990 5 : Copy_Scale_sig_32_16( msCurrentMinOut, msNoiseFloor, len, -16 );
991 : }
992 : ELSE /* sub window complete */
993 : {
994 4 : IF( GE_16( hFdCngCom->msFrCnt, MSSUBFRLEN ) )
995 : {
996 : /* Collect buffers */
997 2 : Copy32( msCurrentMinSubWindow, msMinBuf + len * hFdCngCom->msMinBufferPtr, len );
998 :
999 : /* Compute minimum among all buffers */
1000 2 : Copy32( msMinBuf, msCurrentMinOut, len );
1001 2 : ptr = msMinBuf + len;
1002 12 : FOR( i = 1; i < MSNUMSUBFR; i++ )
1003 : {
1004 625 : FOR( j = 0; j < len; j++ )
1005 : {
1006 615 : if ( LT_32( *ptr, msCurrentMinOut[j] ) /*0.0 Q31*/ )
1007 : {
1008 0 : msCurrentMinOut[j] = *ptr;
1009 0 : move32();
1010 : }
1011 615 : ptr++;
1012 : }
1013 : }
1014 :
1015 : /* Take over local minima */
1016 2 : slope = hFdCngCom->msSlope[0];
1017 2 : move16();
1018 125 : FOR( j = 0; j < len; j++ )
1019 : {
1020 123 : if ( EQ_16( j, lenFFT ) )
1021 : {
1022 0 : slope = hFdCngCom->msSlope[1];
1023 0 : move16();
1024 : }
1025 123 : test();
1026 123 : test();
1027 123 : test();
1028 123 : IF( ( msLocalMinFlag[j] != 0 ) && ( msNewMinFlag[j] == 0 ) && ( LT_32( L_shr( msCurrentMinSubWindow[j], 1 ), Mpy_32_16_1( msCurrentMinOut[j], slope ) ) /*0.0 Q31*/ ) && ( GT_32( msCurrentMinSubWindow[j], msCurrentMinOut[j] ) /*0.0 Q31*/ ) )
1029 : {
1030 0 : msCurrentMinOut[j] = msCurrentMinSubWindow[j];
1031 0 : move32();
1032 0 : i = j;
1033 0 : move16();
1034 0 : FOR( k = 0; k < MSNUMSUBFR; k++ )
1035 : {
1036 0 : msMinBuf[i] = msCurrentMinOut[j];
1037 0 : move32();
1038 0 : i = add( i, len );
1039 : }
1040 : }
1041 : }
1042 :
1043 : /* Reset */
1044 2 : set16_fx( msLocalMinFlag, 0, len );
1045 2 : set32_fx( msCurrentMin, 2147483647l /*1.0 Q31*/, len );
1046 :
1047 : /* Get the current noise floor */
1048 2 : Copy_Scale_sig_32_16( msCurrentMinOut, msNoiseFloor, len, -16 );
1049 : }
1050 : }
1051 :
1052 :
1053 : /* Detect sudden offsets based on the FFT bins (core bandwidth) */
1054 9 : IF( GT_32( Mpy_32_16_1( msPsdSum[0], 655 /*0.02 Q15*/ ), msPeriodogSum[0] ) /*0.0 Q31*/ )
1055 : {
1056 0 : IF( hFdCngCom->offsetflag > 0 )
1057 : {
1058 0 : Copy( msPeriodog, msPsd, len );
1059 0 : FOR( j = 0; j < len; j++ )
1060 : {
1061 0 : msCurrentMinOut[j] = L_deposit_h( msPeriodog[j] );
1062 0 : move32();
1063 : }
1064 0 : set32_fx( hFdCngCom->msAlphaCor, 2147483647l /*1.0 Q31*/, cnt );
1065 0 : set32_fx( msAlpha, 0l /*0.0 Q31*/, len );
1066 0 : Copy( msPeriodog, msPsdFirstMoment, len );
1067 0 : set32_fx( msPsdSecondMoment, 0l /*0.0 Q31*/, len );
1068 :
1069 0 : msPsdSum[0] = dotp_s_fx( msPeriodog, psize, lenFFT, CNG_HS );
1070 0 : move32();
1071 0 : IF( LT_16( lenFFT, len ) )
1072 : {
1073 0 : msPsdSum[1] = dotp_s_fx( msPeriodog + lenFFT, psize + lenFFT, sub( len, lenFFT ), CNG_HS );
1074 0 : move32();
1075 : }
1076 : }
1077 0 : hFdCngCom->offsetflag = 1;
1078 0 : move16();
1079 : }
1080 : ELSE
1081 : {
1082 9 : hFdCngCom->offsetflag = 0;
1083 9 : move16();
1084 : }
1085 :
1086 :
1087 : /* Increment frame counter */
1088 9 : IF( EQ_16( hFdCngCom->msFrCnt, MSSUBFRLEN ) )
1089 : {
1090 2 : hFdCngCom->msFrCnt = 1;
1091 2 : move16();
1092 2 : hFdCngCom->msMinBufferPtr = add( hFdCngCom->msMinBufferPtr, 1 );
1093 2 : move16();
1094 2 : if ( EQ_16( hFdCngCom->msMinBufferPtr, MSNUMSUBFR ) )
1095 : {
1096 0 : hFdCngCom->msMinBufferPtr = 0;
1097 0 : move16();
1098 : }
1099 : }
1100 : ELSE
1101 : {
1102 7 : hFdCngCom->msFrCnt = add( hFdCngCom->msFrCnt, 1 );
1103 : }
1104 :
1105 : /* Smooth noise estimate during CNG phases */
1106 560 : FOR( j = 0; j < len; j++ )
1107 : {
1108 551 : msNoiseEst[j] = round_fx( L_mac( L_mult( 31130 /*0.95 Q15*/, msNoiseEst[j] ), 1638 /*0.05 Q15*/, msNoiseFloor[j] ) );
1109 551 : move16();
1110 : }
1111 : }
1112 : /* Collect buffers */
1113 15 : Copy( msPeriodog, msPeriodogBuf + len * ( *msPeriodogBufPtr ), len );
1114 :
1115 15 : *msPeriodogBufPtr = add( *msPeriodogBufPtr, 1 );
1116 15 : move16();
1117 15 : if ( EQ_16( *msPeriodogBufPtr, MSBUFLEN ) )
1118 : {
1119 2 : ( *msPeriodogBufPtr ) = 0;
1120 2 : move16();
1121 : }
1122 :
1123 : /* Upper limit the noise floors with the averaged input energy */
1124 937 : FOR( j = 0; j < len; j++ )
1125 : {
1126 922 : scalar = L_mult( msPeriodogBuf[j], 6554 /*1.0/MSBUFLEN Q15*/ );
1127 :
1128 4610 : FOR( i = j + len; i < MSBUFLEN * len; i += len )
1129 : {
1130 3688 : scalar = L_mac( scalar, msPeriodogBuf[i], 6554 /*1.0/MSBUFLEN Q15*/ );
1131 : }
1132 922 : scalar16 = round_fx( scalar );
1133 922 : if ( GT_16( msNoiseEst[j], scalar16 ) /*0.0 Q15*/ )
1134 : {
1135 378 : msNoiseEst[j] = scalar16;
1136 378 : move16();
1137 : }
1138 :
1139 922 : assert( msNoiseEst[j] >= 0 /*0.0 Q15*/ );
1140 : }
1141 15 : }
1142 :
1143 42138 : void minimum_statistics_fx(
1144 : Word16 len, /* i : Total number of partitions (CLDFB or FFT) */
1145 : Word16 lenFFT, /* i : Number of FFT partitions */
1146 : Word16 *psize, /* i : Partition sizes, fractional */
1147 : Word16 *msPeriodog, /* i : Periodogram (energies) */
1148 : Word16 *msNoiseFloor, /* i/o: Noise floors (energies) */
1149 : Word16 *msNoiseEst, /* i/o: Noise estimates (energies) */
1150 : Word32 *msAlpha, /* i/o: Forgetting factors */
1151 : Word16 *msPsd, /* i/o: Power Spectral Density (smoothed periodogram => energies) */
1152 : Word16 *msPsdFirstMoment, /* i/o: PSD statistics of 1st order (energy means) */
1153 : Word32 *msPsdSecondMoment, /* i/o: PSD statistics of 2nd order (energy variances) */
1154 : Word32 *msMinBuf, /* i/o: Buffer of minima (energies) */
1155 : Word32 *msBminWin, /* o : Bias correction factors */
1156 : Word32 *msBminSubWin, /* o : Bias correction factors */
1157 : Word32 *msCurrentMin, /* i/o: Local minima (energies) */
1158 : Word32 *msCurrentMinOut, /* i/o: Local minima (energies) */
1159 : Word32 *msCurrentMinSubWindow, /* i/o: Local minima (energies) */
1160 : Word16 *msLocalMinFlag, /* i : Binary flag */
1161 : Word16 *msNewMinFlag, /* i : Binary flag */
1162 : Word16 *msPeriodogBuf, /* i/o: Buffer of periodograms (energies) */
1163 : Word16 *msPeriodogBufPtr, /* i/o: Counter */
1164 : HANDLE_FD_CNG_COM hFdCngCom, /* i/o: FD_CNG structure containing all buffers and variables */
1165 : const Word16 enc_dec, /* i : encoder/decoder indicator */
1166 : const Word16 element_mode /* i : IVAS element mode type */
1167 : )
1168 : {
1169 : Word16 i, j, k, s, s1, s2, s3;
1170 : Word16 len2;
1171 : Word16 current_len;
1172 : Word16 start, stop, cnt;
1173 : Word16 totsize;
1174 : Word16 inv_totsize;
1175 :
1176 : Word32 tmp, tmp0, tmp1;
1177 : Word32 scalar, scalar2, scalar3;
1178 : Word32 snr;
1179 : Word32 msAlphaHatMin2;
1180 : Word32 BminCorr;
1181 : Word32 QeqInvAv;
1182 : Word32 *ptr;
1183 : Word32 *msPsdSum;
1184 : Word32 *msPeriodogSum;
1185 :
1186 : Word16 tmp16, tmp16_1;
1187 : Word16 beta;
1188 : Word16 slope;
1189 : Word16 QeqInv;
1190 : Word16 scalar16;
1191 : Word16 scalar216;
1192 : Word16 scalar316;
1193 : Word16 msM_win;
1194 : Word16 msM_subwin;
1195 : Word16 msAlphaCorAlpha;
1196 : Word16 msAlphaCorAlpha2;
1197 : Word16 msPeriodogSum16;
1198 : Word16 msNoiseFloor16;
1199 : #ifdef BASOP_NOGLOB_DECLARE_LOCAL
1200 42138 : Flag Overflow = 0;
1201 42138 : move32();
1202 : #endif
1203 :
1204 :
1205 42138 : len2 = i_mult( MSNUMSUBFR, len );
1206 :
1207 42138 : msM_win = hFdCngCom->msM_win;
1208 42138 : move16();
1209 42138 : msM_subwin = hFdCngCom->msM_subwin;
1210 42138 : move16();
1211 :
1212 42138 : msAlphaCorAlpha = MSALPHACORALPHA;
1213 42138 : move16();
1214 42138 : msAlphaCorAlpha2 = MSALPHACORALPHA2;
1215 42138 : move16();
1216 :
1217 42138 : msPsdSum = hFdCngCom->msPsdSum;
1218 42138 : msPeriodogSum = hFdCngCom->msPeriodogSum;
1219 :
1220 : /* No minimum statistics at initialization */
1221 42138 : IF( LT_16( hFdCngCom->msFrCnt_init_counter, hFdCngCom->msFrCnt_init_thresh ) )
1222 : {
1223 3786 : Copy( msPeriodog, msPsd, len ); /* 6Q9 */
1224 3786 : Copy( msPeriodog, msNoiseFloor, len ); /* 6Q9 */
1225 3786 : Copy( msPeriodog, msNoiseEst, len ); /* 6Q9 */
1226 3786 : Copy( msPeriodog, msPsdFirstMoment, len ); /* 6Q9 */
1227 :
1228 3786 : set32_fx( msPsdSecondMoment, 0l /*0.0 Q31*/, len );
1229 3786 : msPeriodogSum[0] = dotp_s_fx( msPeriodog, psize, lenFFT, CNG_HS );
1230 3786 : move32();
1231 3786 : msPsdSum[0] = msPeriodogSum[0]; /* 16Q15 */
1232 3786 : move32();
1233 :
1234 3786 : IF( LT_16( lenFFT, len ) )
1235 : {
1236 0 : msPeriodogSum[1] = dotp_s_fx( msPeriodog + lenFFT, psize + lenFFT, sub( len, lenFFT ), CNG_HS );
1237 0 : move32();
1238 0 : msPsdSum[1] = msPeriodogSum[1]; /* 16Q15 */
1239 0 : move32();
1240 : }
1241 :
1242 : /* Increment frame counter at initialization */
1243 : /* Some frames are sometimes zero at initialization => ignore them */
1244 3786 : IF( LT_16( msPeriodog[0], hFdCngCom->init_old ) )
1245 : {
1246 697 : set32_fx( msCurrentMinOut, 2147483647l /*1.0 Q31*/, len ); /* 16Q15 */
1247 697 : set32_fx( msCurrentMin, 2147483647l /*1.0 Q31*/, len ); /* 16Q15 */
1248 697 : set32_fx( msMinBuf, 2147483647l /*1.0 Q31*/, len2 ); /* 16Q15 */
1249 697 : set32_fx( msCurrentMinSubWindow, 2147483647l /*1.0 Q31*/, len ); /* 16Q15 */
1250 :
1251 697 : hFdCngCom->msFrCnt_init_counter = add( hFdCngCom->msFrCnt_init_counter, 1 );
1252 697 : move16();
1253 : }
1254 3786 : hFdCngCom->init_old = msPeriodog[0]; /* 6Q9 */
1255 3786 : move16();
1256 : }
1257 : ELSE
1258 : {
1259 : /* Consider the FFT and CLDFB bands separately
1260 : - first iteration for FFT bins,
1261 : - second one for CLDFB bands in SWB mode */
1262 38352 : cnt = 0;
1263 38352 : move16();
1264 38352 : start = 0;
1265 38352 : move16();
1266 38352 : stop = lenFFT;
1267 38352 : move16();
1268 38352 : totsize = sub( hFdCngCom->stopFFTbin, hFdCngCom->startBand );
1269 76704 : WHILE( GT_16( stop, start ) )
1270 : {
1271 38352 : current_len = sub( stop, start );
1272 :
1273 : /* Compute smoothed correction factor for PSD smoothing */
1274 :
1275 : /* msPeriodogSum[cnt] with format 16Q15 */
1276 38352 : msPeriodogSum[cnt] = dotp_s_fx( msPeriodog + start, psize + start, current_len, CNG_HS );
1277 38352 : move32();
1278 :
1279 38352 : IF( msPeriodogSum[cnt] == 0 )
1280 : {
1281 0 : hFdCngCom->msAlphaCor[cnt] = Mpy_32_16_1( hFdCngCom->msAlphaCor[cnt], msAlphaCorAlpha );
1282 0 : move32();
1283 : }
1284 : ELSE
1285 : {
1286 : /* calculate scalar with normalized msPeriodogSum[cnt], exponent -2*s1 */
1287 38352 : s1 = norm_l( msPeriodogSum[cnt] );
1288 38352 : msPeriodogSum16 = round_fx_sat( L_shl( msPeriodogSum[cnt], s1 ) );
1289 38352 : scalar = L_mult( msPeriodogSum16, msPeriodogSum16 );
1290 :
1291 : /* calculate difference, both elements in 16Q15 format, use absolute value
1292 : to avoid -1.0 x -1.0 multiplications later */
1293 38352 : scalar2 = L_abs( L_sub( msPsdSum[cnt], msPeriodogSum[cnt] ) );
1294 :
1295 38352 : s2 = WORD32_BITS - 1;
1296 38352 : move16();
1297 38352 : if ( scalar2 != 0 )
1298 : {
1299 : /* use absolute value to avoid -1.0 x -1.0 multiplications */
1300 38287 : s2 = norm_l( scalar2 );
1301 : }
1302 38352 : scalar216 = round_fx_o( L_shl_o( scalar2, s2, &Overflow ), &Overflow );
1303 38352 : scalar2 = L_mult( scalar216, scalar216 );
1304 :
1305 : /* find common exponent */
1306 38352 : tmp16_1 = sub( s1, s2 );
1307 38352 : tmp16 = s_min( shl( abs_s( tmp16_1 ), 1 ), WORD32_BITS - 1 );
1308 38352 : if ( tmp16_1 < 0 )
1309 : {
1310 37001 : scalar2 = L_shr( scalar2, tmp16 );
1311 : }
1312 38352 : if ( tmp16_1 > 0 )
1313 : {
1314 560 : scalar = L_shr( scalar, tmp16 );
1315 : }
1316 :
1317 :
1318 : /* add scalar and scalar2, avoid overflows */
1319 38352 : scalar = L_shr( scalar, 1 );
1320 38352 : scalar2 = L_shr( scalar2, 1 );
1321 38352 : scalar3 = L_add( scalar, scalar2 );
1322 :
1323 : /* calculate division */
1324 38352 : scalar16 = BASOP_Util_Divide3232_uu_1616_Scale( scalar, scalar3, &s3 );
1325 38352 : s3 = s_max( s3, -( WORD16_BITS - 1 ) );
1326 38352 : scalar16 = shl( scalar16, s3 );
1327 38352 : scalar16 = s_max( scalar16, MSALPHACORMAX );
1328 :
1329 38352 : hFdCngCom->msAlphaCor[cnt] = Madd_32_16( L_mult( scalar16, msAlphaCorAlpha2 ), hFdCngCom->msAlphaCor[cnt], msAlphaCorAlpha );
1330 38352 : move32();
1331 : }
1332 :
1333 : /* Compute SNR */
1334 :
1335 : /* msPeriodogSum[cnt] with format 16Q15 */
1336 38352 : snr = dotp_s_fx( msNoiseFloor + start, psize + start, current_len, CNG_HS );
1337 :
1338 38352 : IF( GT_32( L_shr( Mpy_32_16_1( msPsdSum[cnt], 18431 /*0.56246299817 Q15*/ ), 13 ), snr ) )
1339 : {
1340 0 : tmp0 = BASOP_Util_Log2( msPsdSum[cnt] );
1341 0 : tmp1 = BASOP_Util_Log2( snr );
1342 0 : tmp1 = L_sub( tmp0, tmp1 );
1343 0 : tmp1 = Mpy_32_16_1( tmp1, MSSNREXP );
1344 0 : msAlphaHatMin2 = BASOP_Util_InvLog2( tmp1 );
1345 : }
1346 : ELSE
1347 : {
1348 38352 : msAlphaHatMin2 = MSALPHAHATMIN;
1349 38352 : move32();
1350 : }
1351 38352 : scalar = Mpy_32_16_1( hFdCngCom->msAlphaCor[cnt], MSALPHAMAX );
1352 :
1353 2395758 : FOR( j = start; j < stop; j++ )
1354 : {
1355 : /* Compute optimal smoothing parameter for PSD estimation */
1356 2357406 : test();
1357 2357406 : IF( ( scalar == 0 ) || ( msNoiseFloor[j] == 0 ) )
1358 : {
1359 84 : msAlpha[j] = msAlphaHatMin2;
1360 84 : move32();
1361 : }
1362 : ELSE
1363 : {
1364 : /* calculate scalar2 with normalized msNoiseFloor[j], exponent -2*s1 */
1365 2357322 : s1 = WORD16_BITS - 1;
1366 2357322 : move16();
1367 2357322 : if ( msNoiseFloor[j] != 0 )
1368 : {
1369 2357322 : s1 = norm_s( msNoiseFloor[j] );
1370 : }
1371 2357322 : msNoiseFloor16 = shl( msNoiseFloor[j], s1 );
1372 2357322 : scalar2 = L_mult( msNoiseFloor16, msNoiseFloor16 );
1373 :
1374 : /* calculate difference, both elements in 6Q9 format, use absolute value
1375 : to avoid -1.0 x -1.0 multiplications later */
1376 2357322 : scalar316 = abs_s( sub( msPsd[j], msNoiseFloor[j] ) );
1377 :
1378 2357322 : s2 = WORD16_BITS - 1;
1379 2357322 : move16();
1380 2357322 : if ( scalar316 != 0 )
1381 : {
1382 : /* use absolute value to avoid -1.0 x -1.0 multiplications */
1383 1999023 : s2 = norm_s( scalar316 );
1384 : }
1385 2357322 : scalar316 = shl( scalar316, s2 );
1386 2357322 : scalar3 = L_mult( scalar316, scalar316 );
1387 :
1388 : /* find common exponent */
1389 2357322 : tmp16_1 = sub( s1, s2 );
1390 2357322 : tmp16 = s_min( shl( abs_s( tmp16_1 ), 1 ), WORD32_BITS - 1 );
1391 2357322 : if ( tmp16_1 < 0 )
1392 : {
1393 2066520 : scalar3 = L_shr( scalar3, tmp16 );
1394 : }
1395 2357322 : if ( tmp16_1 > 0 )
1396 : {
1397 174822 : scalar2 = L_shr( scalar2, tmp16 );
1398 : }
1399 :
1400 : /* add scalar2 and scalar3, avoid overflows */
1401 2357322 : scalar2 = L_shr( scalar2, 1 );
1402 2357322 : scalar3 = L_shr( scalar3, 1 );
1403 2357322 : scalar3 = L_add( scalar2, scalar3 );
1404 :
1405 : /* calculate division */
1406 2357322 : tmp16 = BASOP_Util_Divide3232_uu_1616_Scale( scalar2, scalar3, &s3 );
1407 2357322 : scalar2 = Mpy_32_16_1( scalar, tmp16 );
1408 2357322 : s3 = s_max( s3, -( WORD32_BITS - 1 ) );
1409 2357322 : scalar2 = L_shl( scalar2, s3 );
1410 2357322 : msAlpha[j] = L_max( scalar2, msAlphaHatMin2 );
1411 2357322 : move32();
1412 : }
1413 :
1414 : /* Compute the PSD (smoothed periodogram) in each band */
1415 2357406 : msPsd[j] = round_fx( Madd_32_16( Mpy_32_16_1( msAlpha[j], msPsd[j] ), L_sub( 2147483647l /*1.0 Q31*/, msAlpha[j] ), msPeriodog[j] ) );
1416 2357406 : move16();
1417 : }
1418 38352 : msPsdSum[cnt] = dotp_s_fx( msPsd + start, psize + start, current_len, CNG_HS );
1419 38352 : move32();
1420 :
1421 38352 : QeqInvAv = 0l /*0.0 Q31*/;
1422 38352 : move32();
1423 :
1424 : /* scalar: 7Q24 format */
1425 38352 : tmp = 1191182336l /*(float)(MSNUMSUBFR*MSSUBFRLEN)-1.0 Q24*/;
1426 38352 : move32();
1427 38352 : scalar = L_sub( tmp, L_mult( round_fx( tmp ), msM_win ) );
1428 :
1429 : /* scalar2: 4Q27 format */
1430 38352 : tmp = 1476395008l /*(float)MSSUBFRLEN-1.0 Q27*/;
1431 38352 : move32();
1432 38352 : scalar2 = L_sub( tmp, L_mult( round_fx( tmp ), msM_subwin ) );
1433 :
1434 2395758 : FOR( j = start; j < stop; j++ )
1435 : {
1436 : /* Compute variance of PSD */
1437 2357406 : tmp = L_min( msAlpha[j], MSBETAMAX_SQRT );
1438 :
1439 2357406 : s1 = WORD32_BITS - 1;
1440 2357406 : move16();
1441 2357406 : if ( tmp != 0 )
1442 : {
1443 2357406 : s1 = norm_l( tmp );
1444 : }
1445 2357406 : s2 = shl( s1, 1 );
1446 :
1447 2357406 : s2 = s_min( s2, WORD32_BITS - 1 );
1448 :
1449 : /* beta: scaled by s2 */
1450 2357406 : tmp16 = round_fx_o( L_shl_o( tmp, s1, &Overflow ), &Overflow );
1451 2357406 : beta = mult_r( tmp16, tmp16 );
1452 :
1453 : /* scalar3: scaled by s3 */
1454 2357406 : scalar3 = L_sub( L_deposit_l( msPsd[j] ), L_deposit_l( msPsdFirstMoment[j] ) );
1455 2357406 : s3 = norm_l( scalar3 );
1456 2357406 : scalar3 = L_shl( scalar3, s3 );
1457 :
1458 : /* msPsdFirstMoment: 6Q9 */
1459 2357406 : tmp = L_msu( L_mult( beta, msPsdFirstMoment[j] ), beta, msPsd[j] );
1460 2357406 : msPsdFirstMoment[j] = add( round_fx( L_shr( tmp, s2 ) ), msPsd[j] );
1461 2357406 : move16();
1462 :
1463 : /* msPsdSecondMoment: 12Q19 */
1464 2357406 : tmp0 = L_shr( Mpy_32_16_r( msPsdSecondMoment[j], beta ), s2 );
1465 2357406 : tmp1 = Mpy_32_32( scalar3, scalar3 );
1466 2357406 : tmp1 = L_shr( L_sub( tmp1, L_shr( Mpy_32_16_r( tmp1, beta ), s2 ) ), sub( shl( s3, 1 ), 32 ) );
1467 2357406 : msPsdSecondMoment[j] = L_add( tmp0, tmp1 );
1468 2357406 : move32();
1469 :
1470 : /* Compute inverse of amount of degrees of freedom */
1471 2357406 : QeqInv = MSQEQINVMAX;
1472 2357406 : move16();
1473 :
1474 2357406 : IF( msNoiseFloor[j] != 0 /*0.0 Q15*/ )
1475 : {
1476 2357322 : tmp = L_mult( msNoiseFloor[j], msNoiseFloor[j] );
1477 2357322 : tmp16 = BASOP_Util_Divide3232_uu_1616_Scale( msPsdSecondMoment[j], tmp, &s );
1478 : /* consider factor of 2 */
1479 2357322 : s = s_min( s_max( sub( s, 1 ), -( WORD16_BITS - 1 ) ), ( WORD16_BITS - 1 ) );
1480 2357322 : if ( s < 0 )
1481 : {
1482 2199215 : QeqInv = shl( tmp16, s );
1483 : }
1484 2357322 : QeqInv = s_min( QeqInv, MSQEQINVMAX );
1485 : }
1486 2357406 : QeqInvAv = L_add( QeqInvAv, L_mult( QeqInv, psize[j] ) );
1487 :
1488 : /* Compute bias correction Bmin */
1489 2357406 : tmp0 = Mpy_32_16_1( scalar, QeqInv );
1490 2357406 : tmp1 = L_msu( 1073741824l /*0.5 Q31*/, msM_win, QeqInv );
1491 2357406 : tmp16 = BASOP_Util_Divide3232_uu_1616_Scale( tmp0, tmp1, &s );
1492 2357406 : msBminWin[j] = L_add( 134217728l /*1.0 Q27*/, L_shl( L_deposit_h( tmp16 ), add( s, 7 - 4 ) ) );
1493 2357406 : move32();
1494 :
1495 2357406 : tmp0 = Mpy_32_16_1( scalar2, QeqInv );
1496 2357406 : tmp1 = L_msu( 1073741824l /*0.5 Q31*/, msM_subwin, QeqInv );
1497 2357406 : tmp16 = BASOP_Util_Divide3232_uu_1616_Scale( tmp0, tmp1, &s );
1498 2357406 : msBminSubWin[j] = L_add( 134217728l /*1.0 Q27*/, L_shl( L_deposit_h( tmp16 ), s ) );
1499 2357406 : move32();
1500 : }
1501 :
1502 38352 : inv_totsize = BASOP_Util_Divide3232_uu_1616_Scale( 1, totsize, &s );
1503 38352 : QeqInvAv = Mpy_32_16_1( QeqInvAv, inv_totsize );
1504 38352 : QeqInvAv = L_shl( QeqInvAv, s );
1505 38352 : hFdCngCom->msQeqInvAv[cnt] = QeqInvAv;
1506 38352 : move32();
1507 :
1508 : /* New minimum? */
1509 :
1510 : /* exponent QeqInvAv: CNG_S, exponent MSAV: (4>>1) */
1511 38352 : s = CNG_S + 2 * MSAV_EXP;
1512 38352 : move16();
1513 38352 : BminCorr = Mpy_32_16_1( Sqrt32( QeqInvAv, &s ), MSAV );
1514 38352 : BminCorr = L_shl( BminCorr, sub( s, 1 ) );
1515 :
1516 : /* exponent BminCorr: 1 */
1517 38352 : BminCorr = L_add( BminCorr, 1073741824l /*0.5 Q31*/ );
1518 :
1519 2395758 : FOR( j = start; j < stop; j++ )
1520 : {
1521 : /* exponent scalar: CNG_S+1 */
1522 2357406 : scalar = Mpy_32_16_1( BminCorr, msPsd[j] );
1523 :
1524 : /* exponent scalar2: CNG_S+1+4 */
1525 2357406 : scalar2 = Mpy_32_32( scalar, msBminWin[j] );
1526 :
1527 2357406 : msNewMinFlag[j] = 0;
1528 2357406 : move16();
1529 2357406 : IF( LT_32( scalar2, msCurrentMin[j] ) /*0.0 Q31*/ )
1530 : {
1531 860523 : msNewMinFlag[j] = 1;
1532 860523 : move16();
1533 : /* exponent msCurrentMin[j]: CNG_S+1+4 */
1534 860523 : msCurrentMin[j] = scalar2;
1535 860523 : move32();
1536 : /* exponent msCurrentMinSubWindow[j]: CNG_S */
1537 : BASOP_SATURATE_WARNING_OFF_EVS;
1538 860523 : msCurrentMinSubWindow[j] = L_shl_o( Mpy_32_32( scalar, msBminSubWin[j] ), 5, &Overflow );
1539 860523 : move32();
1540 : BASOP_SATURATE_WARNING_ON_EVS;
1541 : }
1542 : }
1543 :
1544 : /* This is used later to identify local minima */
1545 38352 : IF( GE_16( hFdCngCom->msFrCnt, MSSUBFRLEN ) )
1546 : {
1547 8946 : FOR( i = 0; i < 3; i++ )
1548 : {
1549 7713 : IF( LT_32( hFdCngCom->msQeqInvAv[cnt], L_shr( L_deposit_h( msQeqInvAv_thresh[i] ), CNG_S ) ) /*0.0 Q31*/ )
1550 : {
1551 2309 : BREAK;
1552 : }
1553 : }
1554 : /* format 1Q30 */
1555 3542 : hFdCngCom->msSlope[cnt] = msNoiseSlopeMax[i];
1556 3542 : move32();
1557 : }
1558 :
1559 : /* Consider the FFT and CLDFB bands separately */
1560 38352 : start = stop;
1561 38352 : move16();
1562 38352 : stop = len;
1563 38352 : move16();
1564 38352 : totsize = sub( hFdCngCom->stopBand, hFdCngCom->stopFFTbin );
1565 38352 : cnt = add( cnt, 1 );
1566 : } /*while (stop > start)*/
1567 :
1568 : /* Update minimum between sub windows */
1569 38352 : test();
1570 38352 : IF( GT_16( hFdCngCom->msFrCnt, 1 ) && LT_16( hFdCngCom->msFrCnt, MSSUBFRLEN ) )
1571 : {
1572 1957882 : FOR( j = 0; j < len; j++ )
1573 : {
1574 1926540 : if ( msNewMinFlag[j] > 0 )
1575 : {
1576 575480 : msLocalMinFlag[j] = 1;
1577 575480 : move16();
1578 : }
1579 1926540 : if ( LT_32( msCurrentMinSubWindow[j], msCurrentMinOut[j] ) /*0.0 Q31*/ )
1580 : {
1581 : /* msCurrentMinOut[j] scaled with CNG_S */
1582 196103 : msCurrentMinOut[j] = msCurrentMinSubWindow[j];
1583 196103 : move32();
1584 : }
1585 : }
1586 : /* Get the current noise floor */
1587 31342 : Copy_Scale_sig_32_16( msCurrentMinOut, msNoiseFloor, len, -16 );
1588 : }
1589 : ELSE /* sub window complete */
1590 : {
1591 7010 : IF( GE_16( hFdCngCom->msFrCnt, MSSUBFRLEN ) )
1592 : {
1593 : /* Collect buffers */
1594 3542 : Copy32( msCurrentMinSubWindow, msMinBuf + len * hFdCngCom->msMinBufferPtr, len );
1595 :
1596 : /* Compute minimum among all buffers */
1597 3542 : Copy32( msMinBuf, msCurrentMinOut, len );
1598 3542 : ptr = msMinBuf + len;
1599 21252 : FOR( i = 1; i < MSNUMSUBFR; i++ )
1600 : {
1601 1106185 : FOR( j = 0; j < len; j++ )
1602 : {
1603 1088475 : if ( LT_32( *ptr, msCurrentMinOut[j] ) /*0.0 Q31*/ )
1604 : {
1605 204920 : msCurrentMinOut[j] = *ptr;
1606 204920 : move32();
1607 : }
1608 1088475 : ptr++;
1609 : }
1610 : }
1611 :
1612 : /* Take over local minima */
1613 3542 : slope = hFdCngCom->msSlope[0];
1614 3542 : move16();
1615 221237 : FOR( j = 0; j < len; j++ )
1616 : {
1617 217695 : if ( EQ_16( j, lenFFT ) )
1618 : {
1619 0 : slope = hFdCngCom->msSlope[1];
1620 0 : move16();
1621 : }
1622 217695 : test();
1623 217695 : test();
1624 217695 : test();
1625 217695 : IF( ( msLocalMinFlag[j] != 0 ) && ( msNewMinFlag[j] == 0 ) && ( LT_32( L_shr( msCurrentMinSubWindow[j], 1 ), Mpy_32_16_1( msCurrentMinOut[j], slope ) ) /*0.0 Q31*/ ) && ( GT_32( msCurrentMinSubWindow[j], msCurrentMinOut[j] ) /*0.0 Q31*/ ) )
1626 : {
1627 30345 : msCurrentMinOut[j] = msCurrentMinSubWindow[j];
1628 30345 : move32();
1629 30345 : i = j;
1630 30345 : move16();
1631 212415 : FOR( k = 0; k < MSNUMSUBFR; k++ )
1632 : {
1633 182070 : msMinBuf[i] = msCurrentMinOut[j];
1634 182070 : move32();
1635 182070 : i = add( i, len );
1636 : }
1637 : }
1638 : }
1639 :
1640 : /* Reset */
1641 3542 : set16_fx( msLocalMinFlag, 0, len );
1642 3542 : set32_fx( msCurrentMin, 2147483647l /*1.0 Q31*/, len );
1643 :
1644 : /* Get the current noise floor */
1645 3542 : Copy_Scale_sig_32_16( msCurrentMinOut, msNoiseFloor, len, -16 );
1646 : }
1647 : }
1648 :
1649 :
1650 : /* Detect sudden offsets based on the FFT bins (core bandwidth) */
1651 38352 : IF( GT_32( Mpy_32_16_1( msPsdSum[0], 655 /*0.02 Q15*/ ), msPeriodogSum[0] ) /*0.0 Q31*/ )
1652 : {
1653 113 : IF( hFdCngCom->offsetflag > 0 )
1654 : {
1655 47 : Copy( msPeriodog, msPsd, len );
1656 2927 : FOR( j = 0; j < len; j++ )
1657 : {
1658 2880 : msCurrentMinOut[j] = L_deposit_h( msPeriodog[j] );
1659 2880 : move32();
1660 : }
1661 47 : set32_fx( hFdCngCom->msAlphaCor, 2147483647l /*1.0 Q31*/, cnt );
1662 47 : set32_fx( msAlpha, 0l /*0.0 Q31*/, len );
1663 47 : Copy( msPeriodog, msPsdFirstMoment, len );
1664 47 : set32_fx( msPsdSecondMoment, 0l /*0.0 Q31*/, len );
1665 :
1666 47 : msPsdSum[0] = dotp_s_fx( msPeriodog, psize, lenFFT, CNG_HS );
1667 47 : move32();
1668 47 : IF( LT_16( lenFFT, len ) )
1669 : {
1670 0 : msPsdSum[1] = dotp_s_fx( msPeriodog + lenFFT, psize + lenFFT, sub( len, lenFFT ), CNG_HS );
1671 0 : move32();
1672 : }
1673 : }
1674 113 : hFdCngCom->offsetflag = 1;
1675 113 : move16();
1676 : }
1677 : ELSE
1678 : {
1679 38239 : hFdCngCom->offsetflag = 0;
1680 38239 : move16();
1681 : }
1682 :
1683 :
1684 : /* Increment frame counter */
1685 38352 : IF( EQ_16( hFdCngCom->msFrCnt, MSSUBFRLEN ) )
1686 : {
1687 3542 : hFdCngCom->msFrCnt = 1;
1688 3542 : move16();
1689 3542 : hFdCngCom->msMinBufferPtr = add( hFdCngCom->msMinBufferPtr, 1 );
1690 3542 : move16();
1691 3542 : if ( EQ_16( hFdCngCom->msMinBufferPtr, MSNUMSUBFR ) )
1692 : {
1693 392 : hFdCngCom->msMinBufferPtr = 0;
1694 392 : move16();
1695 : }
1696 : }
1697 : ELSE
1698 : {
1699 34810 : hFdCngCom->msFrCnt = add( hFdCngCom->msFrCnt, 1 );
1700 34810 : move16();
1701 : }
1702 :
1703 : /* Smooth noise estimate during CNG phases */
1704 2395758 : FOR( j = 0; j < len; j++ )
1705 : {
1706 2357406 : msNoiseEst[j] = round_fx( L_mac( L_mult( 31130 /*0.95 Q15*/, msNoiseEst[j] ), 1638 /*0.05 Q15*/, msNoiseFloor[j] ) );
1707 2357406 : move16();
1708 : }
1709 : }
1710 42138 : IF( EQ_16( enc_dec, DEC ) && EQ_16( element_mode, IVAS_CPE_TD ) )
1711 : {
1712 : // v_multc(msNoiseEst, 1.4125f, msNoiseEst, NPART_SHAPING);
1713 0 : v_multc_att( msNoiseEst, 23142, msNoiseEst, NPART_SHAPING );
1714 : }
1715 : /* Collect buffers */
1716 42138 : Copy( msPeriodog, msPeriodogBuf + len * ( *msPeriodogBufPtr ), len );
1717 :
1718 42138 : *msPeriodogBufPtr = add( *msPeriodogBufPtr, 1 );
1719 42138 : move16();
1720 42138 : if ( EQ_16( *msPeriodogBufPtr, MSBUFLEN ) )
1721 : {
1722 8070 : ( *msPeriodogBufPtr ) = 0;
1723 8070 : move16();
1724 : }
1725 :
1726 : /* Upper limit the noise floors with the averaged input energy */
1727 2631815 : FOR( j = 0; j < len; j++ )
1728 : {
1729 2589677 : scalar = L_mult( msPeriodogBuf[j], 6554 /*1.0/MSBUFLEN Q15*/ );
1730 :
1731 12948385 : FOR( i = j + len; i < MSBUFLEN * len; i += len )
1732 : {
1733 10358708 : scalar = L_mac( scalar, msPeriodogBuf[i], 6554 /*1.0/MSBUFLEN Q15*/ );
1734 : }
1735 2589677 : scalar16 = round_fx( scalar );
1736 2589677 : if ( GT_16( msNoiseEst[j], scalar16 ) /*0.0 Q15*/ )
1737 : {
1738 478202 : msNoiseEst[j] = scalar16;
1739 478202 : move16();
1740 : }
1741 :
1742 2589677 : assert( msNoiseEst[j] >= 0 /*0.0 Q15*/ );
1743 : }
1744 42138 : }
1745 :
1746 : /*-------------------------------------------------------------------
1747 : * apply_scale()
1748 : *
1749 : * Apply bitrate-dependent scale
1750 : *-------------------------------------------------------------------*/
1751 0 : void apply_scale(
1752 : Word32 *scale, /* o : scalefactor */
1753 : const Word16 bwmode, /* i : audio bandwidth */
1754 : const Word32 bitrate, /* i : Bit rate */
1755 : const SCALE_SETUP *scaleTable, /* i : Scale table Q7 */
1756 : const Word16 scaleTableSize /* i : Size of scale table */
1757 : )
1758 : {
1759 : Word16 i;
1760 : // PMT("Verifiy if the basop ued are ok for stereo too")
1761 0 : FOR( i = 0; i < scaleTableSize; i++ )
1762 : {
1763 0 : cast16();
1764 0 : IF( s_and( (Word16) EQ_16( bwmode, (Word16) scaleTable[i].bwmode ),
1765 : s_and( L_sub( bitrate, scaleTable[i].bitrateFrom ) >= 0,
1766 : L_sub( bitrate, scaleTable[i].bitrateTo ) < 0 ) ) )
1767 : {
1768 0 : BREAK;
1769 : }
1770 : }
1771 :
1772 : {
1773 0 : *scale = L_add( *scale, L_deposit_h( scaleTable[i].scale ) );
1774 0 : move32();
1775 : }
1776 0 : }
1777 :
1778 : /*-------------------------------------------------------------------
1779 : * apply_scale_ind()
1780 : *
1781 : * Apply bitrate-dependent scale
1782 : * returns index of scaleTable
1783 : *-------------------------------------------------------------------*/
1784 0 : Word16 apply_scale_ind(
1785 : Word32 *scale, /* o : scalefactor */
1786 : const Word16 bwmode, /* i : audio bandwidth */
1787 : const Word32 bitrate, /* i : Bit rate */
1788 : const SCALE_SETUP *scaleTable, /* i : Scale table */
1789 : const Word16 scaleTableSize /* i : Size of scale table */
1790 : )
1791 : {
1792 : Word16 i;
1793 : // PMT("Verifiy if the basop ued are ok for stereo too")
1794 0 : FOR( i = 0; i < scaleTableSize; i++ )
1795 : {
1796 0 : cast16();
1797 0 : IF( s_and( (Word16) EQ_16( bwmode, (Word16) scaleTable[i].bwmode ),
1798 : s_and( L_sub( bitrate, scaleTable[i].bitrateFrom ) >= 0,
1799 : L_sub( bitrate, scaleTable[i].bitrateTo ) < 0 ) ) )
1800 : {
1801 0 : BREAK;
1802 : }
1803 : }
1804 :
1805 : {
1806 0 : *scale = L_add( *scale, L_deposit_h( scaleTable[i].scale ) );
1807 0 : move32();
1808 : }
1809 0 : return i;
1810 : }
1811 :
1812 166000 : void apply_scale_ivas_fx(
1813 : Word32 *scale, /* o : scalefactor */
1814 : const Word16 bwmode, /* i : audio bandwidth */
1815 : const Word32 bitrate, /* i : Bit rate */
1816 : const SCALE_SETUP *scaleTable, /* i : Scale table Q7 */
1817 : const Word16 scaleTableSize, /* i : Size of scale table */
1818 : Word16 *index )
1819 : {
1820 : Word16 i;
1821 : // PMT("Verifiy if the basop ued are ok for stereo too")
1822 1749124 : FOR( i = 0; i < scaleTableSize; i++ )
1823 : {
1824 1749124 : cast16();
1825 1749124 : IF( s_and( (Word16) EQ_16( bwmode, (Word16) scaleTable[i].bwmode ),
1826 : s_and( L_sub( bitrate, scaleTable[i].bitrateFrom ) >= 0,
1827 : L_sub( bitrate, scaleTable[i].bitrateTo ) < 0 ) ) )
1828 : {
1829 166000 : BREAK;
1830 : }
1831 : }
1832 166000 : assert( i < scaleTableSize );
1833 :
1834 166000 : *scale = L_add( *scale, L_deposit_h( scaleTable[i].scale ) );
1835 166000 : move32();
1836 166000 : *index = i;
1837 166000 : move16();
1838 166000 : }
1839 : /*-------------------------------------------------------------------
1840 : * bandcombinepow()
1841 : *
1842 : * Compute the power for each partition
1843 : *-------------------------------------------------------------------*/
1844 51429 : void bandcombinepow(
1845 : const Word32 *bandpow, /* i : Power for each band */
1846 : const Word16 exp_bandpow, /* i : exponent of bandpow */
1847 : const Word16 nband, /* i : Number of bands */
1848 : Word16 *part, /* i : Partition upper boundaries (band indices starting from 0) */
1849 : const Word16 npart, /* i : Number of partitions */
1850 : const Word16 *psize_inv, /* i : Inverse partition sizes */
1851 : Word32 *partpow, /* o : Power for each partition */
1852 : Word16 *exp_partpow )
1853 : {
1854 :
1855 : Word16 i, p;
1856 : Word32 temp;
1857 : Word16 smin, len, prev_part;
1858 : Word16 facTabExp[NPART_SHAPING];
1859 :
1860 :
1861 51429 : IF( EQ_16( nband, npart ) )
1862 : {
1863 0 : Copy32( bandpow, partpow, nband ); // Q(31 - exp_bandpow)
1864 0 : smin = 0;
1865 0 : move16();
1866 : }
1867 : ELSE
1868 : {
1869 : /* Compute the power in each partition */
1870 51429 : prev_part = -1;
1871 51429 : move16();
1872 3056131 : FOR( p = 0; p < npart; p++ )
1873 : {
1874 3004702 : len = sub( part[p], prev_part );
1875 3004702 : facTabExp[p] = getScaleFactor32( &bandpow[prev_part + 1], len );
1876 3004702 : move16();
1877 3004702 : prev_part = part[p];
1878 3004702 : move16();
1879 : }
1880 :
1881 51429 : smin = WORD32_BITS - 1;
1882 51429 : move16();
1883 3056131 : FOR( p = 0; p < npart; p++ )
1884 : {
1885 3004702 : smin = s_min( smin, facTabExp[p] );
1886 : }
1887 :
1888 51429 : i = 0;
1889 51429 : move16();
1890 3056131 : FOR( p = 0; p < npart; p++ )
1891 : {
1892 : /* Arithmetic averaging of power for all bins in partition */
1893 3004702 : temp = 0;
1894 3004702 : move32();
1895 17840724 : FOR( ; i <= part[p]; i++ )
1896 : {
1897 14836022 : temp = L_add( temp, Mpy_32_16_1( L_shl( bandpow[i], facTabExp[p] ), psize_inv[p] ) );
1898 : }
1899 3004702 : partpow[p] = L_shr( temp, sub( facTabExp[p], smin ) );
1900 3004702 : move32();
1901 : }
1902 : }
1903 :
1904 51429 : *exp_partpow = sub( exp_bandpow, smin );
1905 51429 : move16();
1906 51429 : }
1907 :
1908 : /*-------------------------------------------------------------------
1909 : * scalebands()
1910 : *
1911 : * Scale partitions (with smoothing)
1912 : *-------------------------------------------------------------------*/
1913 3082 : void scalebands(
1914 : const Word32 *partpow, /* i : Power for each partition Qx*/
1915 : Word16 *part, /* i : Partition upper boundaries (band indices starting from 0) */
1916 : const Word16 npart, /* i : Number of partitions */
1917 : Word16 *midband, /* i : Central band of each partition */
1918 : const Word16 nFFTpart, /* i : Number of FFT partitions */
1919 : const Word16 nband, /* i : Number of bands */
1920 : Word32 *bandpow, /* o : Power for each band Qx*/
1921 : const Word16 flag_fft_en )
1922 : {
1923 : Word16 i, j, s, s1, nint, delta, delta_cmp, delta_s;
1924 : Word16 startBand, startPart, stopPart, stopPartM1;
1925 : Word32 tmp, val, partpowLD64, partpowLD64M1;
1926 :
1927 :
1928 3082 : j = 0;
1929 3082 : move16();
1930 3082 : delta = 0;
1931 3082 : move16();
1932 3082 : partpowLD64M1 = 0L; /* to avoid compilation warnings */
1933 3082 : move32();
1934 :
1935 : /* Interpolate the bin/band-wise levels from the partition levels */
1936 3082 : IF( EQ_16( nband, npart ) )
1937 : {
1938 0 : Copy32( partpow, bandpow, npart );
1939 : }
1940 : ELSE
1941 : {
1942 3082 : startBand = 0;
1943 3082 : move16();
1944 3082 : startPart = 0;
1945 3082 : move16();
1946 3082 : stopPart = nFFTpart;
1947 3082 : move16();
1948 :
1949 8624 : WHILE( startBand < nband )
1950 : {
1951 5542 : stopPartM1 = sub( stopPart, 1 );
1952 5542 : test();
1953 5542 : IF( ( flag_fft_en != 0 ) || ( GE_16( startPart, nFFTpart ) ) )
1954 : {
1955 : /* first half partition */
1956 5542 : j = startPart;
1957 5542 : move16();
1958 :
1959 21063 : FOR( i = startBand; i <= midband[j]; i++ )
1960 : {
1961 15521 : bandpow[i] = partpow[j];
1962 15521 : move32();
1963 : }
1964 5542 : j = add( j, 1 );
1965 :
1966 : /* inner partitions */
1967 5542 : IF( LT_16( j, stopPart ) )
1968 : {
1969 5394 : partpowLD64M1 = BASOP_Util_Log2( partpow[j - 1] );
1970 : }
1971 :
1972 : /* Debug values to check this variable is set. */
1973 5542 : delta = 0x4000; // 1.Q14
1974 5542 : move16();
1975 5542 : delta_cmp = 0x4000; // 1.Q14
1976 5542 : move16();
1977 5542 : s1 = 1;
1978 5542 : move16();
1979 5542 : s = 1;
1980 5542 : move16();
1981 :
1982 72263 : FOR( ; j < stopPart; j++ )
1983 : {
1984 66721 : nint = sub( midband[j], midband[j - 1] );
1985 :
1986 : /* log-linear interpolation */
1987 66721 : partpowLD64 = BASOP_Util_Log2( partpow[j] );
1988 66721 : tmp = L_sub( partpowLD64, partpowLD64M1 );
1989 66721 : tmp = Mpy_32_16_1( tmp, getNormReciprocalWord16( nint ) );
1990 :
1991 : /* scale logarithmic value */
1992 66721 : tmp = L_sub( tmp, DELTA_SHIFT_LD64 );
1993 66721 : delta_s = DELTA_SHIFT;
1994 66721 : move16();
1995 :
1996 66769 : WHILE( tmp > 0 )
1997 : {
1998 48 : tmp = L_sub( tmp, 33554432l /*0.015625 Q31*/ ); // Q31
1999 48 : delta_s = add( delta_s, 1 );
2000 : }
2001 66721 : delta_cmp = shl( 1, s_max( -15, sub( WORD16_BITS - 1, delta_s ) ) );
2002 :
2003 66721 : tmp = BASOP_Util_InvLog2( tmp );
2004 66721 : s = norm_l( tmp );
2005 66721 : s1 = sub( delta_s, s );
2006 :
2007 66721 : delta = round_fx_sat( L_shl_sat( tmp, s ) ); // Q(14+s)
2008 : /* Choose scale such that the interpolation start and end point both are representable and add 1 additional bit hr. */
2009 66721 : delta_s = sub( s_min( norm_l( partpow[j - 1] ), norm_l( partpow[j] ) ), 1 );
2010 66721 : val = L_shl( partpow[j - 1], delta_s );
2011 845595 : FOR( ; i < midband[j]; i++ )
2012 : {
2013 778874 : val = L_shl( Mpy_32_16_1( val, delta ), s1 );
2014 778874 : bandpow[i] = L_shr( val, delta_s );
2015 778874 : move32();
2016 : }
2017 66721 : bandpow[i++] = partpow[j];
2018 66721 : move32();
2019 66721 : partpowLD64M1 = partpowLD64;
2020 66721 : move32();
2021 : }
2022 :
2023 5542 : IF( GT_16( shr( delta, s ), delta_cmp ) )
2024 : {
2025 1382 : delta = 0x4000; // 1.Q14
2026 1382 : move16();
2027 1382 : s1 = 1;
2028 1382 : move16();
2029 : }
2030 :
2031 : /* last half partition */
2032 5542 : val = partpow[stopPartM1]; // Qx
2033 5542 : move32();
2034 101031 : FOR( ; i <= part[stopPartM1]; i++ )
2035 : {
2036 95489 : val = L_shl( Mpy_32_16_1( val, delta ), s1 );
2037 95489 : bandpow[i] = val; // Qx
2038 95489 : move32();
2039 : }
2040 : }
2041 5542 : startBand = add( part[stopPartM1], 1 );
2042 5542 : startPart = stopPart;
2043 5542 : move16();
2044 5542 : stopPart = npart;
2045 5542 : move16();
2046 : }
2047 : }
2048 3082 : }
2049 :
2050 96040 : void scalebands_fx(
2051 : const Word32 *partpow, /* i : Power for each partition Qx*/
2052 : Word16 *part, /* i : Partition upper boundaries (band indices starting from 0) */
2053 : const Word16 npart, /* i : Number of partitions */
2054 : Word16 *midband, /* i : Central band of each partition */
2055 : const Word16 nFFTpart, /* i : Number of FFT partitions */
2056 : const Word16 nband, /* i : Number of bands */
2057 : Word32 *bandpow, /* o : Power for each band Qx*/
2058 : const Word16 flag_fft_en )
2059 : {
2060 : Word16 i, j, s, s1, nint, delta, delta_cmp, delta_s;
2061 : Word16 startBand, startPart, stopPart, stopPartM1;
2062 : Word32 tmp, val, partpowLD64, partpowLD64M1;
2063 :
2064 :
2065 96040 : j = 0;
2066 96040 : move16();
2067 96040 : delta = 0;
2068 96040 : move16();
2069 96040 : partpowLD64M1 = 0L; /* to avoid compilation warnings */
2070 96040 : move32();
2071 :
2072 : /* Interpolate the bin/band-wise levels from the partition levels */
2073 96040 : IF( EQ_16( nband, npart ) )
2074 : {
2075 0 : Copy32( partpow, bandpow, npart );
2076 : }
2077 : ELSE
2078 : {
2079 96040 : startBand = 0;
2080 96040 : move16();
2081 96040 : startPart = 0;
2082 96040 : move16();
2083 96040 : stopPart = nFFTpart;
2084 96040 : move16();
2085 :
2086 194097 : WHILE( LT_16( startBand, nband ) )
2087 : {
2088 98057 : stopPartM1 = sub( stopPart, 1 );
2089 98057 : test();
2090 98057 : IF( ( flag_fft_en != 0 ) || ( GE_16( startPart, nFFTpart ) ) )
2091 : {
2092 : /* first half partition */
2093 96062 : j = startPart;
2094 96062 : move16();
2095 :
2096 195374 : FOR( i = startBand; i <= midband[j]; i++ )
2097 : {
2098 99312 : bandpow[i] = partpow[j]; // Qx
2099 99312 : move32();
2100 : }
2101 96062 : j = add( j, 1 );
2102 :
2103 : /* inner partitions */
2104 96062 : IF( j < stopPart )
2105 : {
2106 95837 : partpowLD64M1 = BASOP_Util_Log2( L_add( partpow[j - 1], DELTA_FX ) );
2107 : }
2108 :
2109 : /* Debug values to check this variable is set. */
2110 96062 : delta = 0x4000; // 1.Q14
2111 96062 : move16();
2112 96062 : delta_cmp = 0x4000; // 1.Q14
2113 96062 : move16();
2114 96062 : s1 = 1;
2115 96062 : move16();
2116 96062 : s = 1;
2117 96062 : move16();
2118 :
2119 5783865 : FOR( ; j < stopPart; j++ )
2120 : {
2121 5687803 : nint = sub( midband[j], midband[j - 1] );
2122 :
2123 : /* log-linear interpolation */
2124 5687803 : IF( NE_32( partpow[j], MAX_32 ) )
2125 : {
2126 5687803 : partpowLD64 = BASOP_Util_Log2( L_add( partpow[j], DELTA_FX ) );
2127 : }
2128 : ELSE
2129 : {
2130 0 : partpowLD64 = BASOP_Util_Log2( partpow[j] );
2131 : }
2132 5687803 : tmp = L_sub( partpowLD64, partpowLD64M1 );
2133 5687803 : tmp = Mpy_32_16_1( tmp, getNormReciprocalWord16( nint ) );
2134 :
2135 : /* scale logarithmic value */
2136 5687803 : tmp = L_sub( tmp, DELTA_SHIFT_LD64 );
2137 5687803 : delta_s = DELTA_SHIFT;
2138 5687803 : move16();
2139 :
2140 5844213 : WHILE( tmp > 0 )
2141 : {
2142 156410 : tmp = L_sub( tmp, 33554432l /*0.015625 Q31*/ ); // Q31
2143 156410 : delta_s = add( delta_s, 1 );
2144 : }
2145 5687803 : delta_cmp = shl( 1, s_max( -15, sub( WORD16_BITS - 1, delta_s ) ) );
2146 :
2147 5687803 : tmp = BASOP_Util_InvLog2( tmp );
2148 5687803 : s = norm_l( tmp );
2149 5687803 : s1 = sub( delta_s, s );
2150 :
2151 5687803 : delta = round_fx_sat( L_shl_sat( tmp, s ) ); // Q(14+s)
2152 : /* Choose scale such that the interpolation start and end point both are representable and add 1 additional bit hr. */
2153 5687803 : delta_s = sub( s_min( norm_l( partpow[j - 1] ), norm_l( partpow[j] ) ), 1 );
2154 5687803 : val = L_shl( partpow[j - 1], delta_s );
2155 26345313 : FOR( ; i < midband[j]; i++ )
2156 : {
2157 20657510 : val = L_shl( Mpy_32_16_1( val, delta ), s1 );
2158 20657510 : bandpow[i] = L_shr( val, delta_s );
2159 20657510 : move32();
2160 : }
2161 5687803 : bandpow[i++] = partpow[j]; // Qx
2162 5687803 : move32();
2163 5687803 : partpowLD64M1 = partpowLD64;
2164 5687803 : move32();
2165 : }
2166 :
2167 96062 : IF( GT_16( shr( delta, s ), delta_cmp ) )
2168 : {
2169 12184 : delta = 0x4000; // 1.Q14
2170 12184 : move16();
2171 12184 : s1 = 1;
2172 12184 : move16();
2173 : }
2174 :
2175 : /* last half partition */
2176 96062 : val = partpow[stopPartM1]; // Qx
2177 96062 : move32();
2178 833624 : FOR( ; i <= part[stopPartM1]; i++ )
2179 : {
2180 737562 : test();
2181 737562 : IF( val != 0 && delta != 0 )
2182 : {
2183 478303 : val = L_shl( Mpy_32_16_1( val, delta ), s1 );
2184 478303 : IF( val == 0 )
2185 : {
2186 10759 : val = 1;
2187 10759 : move32();
2188 : }
2189 : }
2190 : ELSE
2191 : {
2192 259259 : val = 0;
2193 259259 : move32();
2194 : }
2195 737562 : bandpow[i] = val; // Qx
2196 737562 : move32();
2197 : }
2198 : }
2199 98057 : startBand = add( part[stopPartM1], 1 );
2200 98057 : startPart = stopPart;
2201 98057 : move16();
2202 98057 : stopPart = npart;
2203 98057 : move16();
2204 : }
2205 : }
2206 96040 : }
2207 :
2208 : /*-------------------------------------------------------------------
2209 : * getmidbands()
2210 : *
2211 : * Get central band for each partition
2212 : *-------------------------------------------------------------------*/
2213 1060682 : static void getmidbands(
2214 : const Word16 *part, /* i : Partition upper boundaries (band indices starting from 0) */
2215 : const Word16 npart, /* i : Number of partitions */
2216 : Word16 *midband, /* o : Central band of each partition */
2217 : Word16 *psize, /* o : Partition sizes */
2218 : Word16 *psize_norm, /* o : Partition sizes, fractional values */
2219 : Word16 *psize_norm_exp, /* o : Exponent for fractional partition sizes */
2220 : Word16 *psize_inv /* o : Inverse of partition sizes */
2221 : )
2222 : {
2223 : Word16 j, max_psize, shift;
2224 : #ifdef BASOP_NOGLOB_DECLARE_LOCAL
2225 1060682 : Flag Overflow = 0;
2226 1060682 : move32();
2227 : #endif
2228 :
2229 :
2230 1060682 : max_psize = psize[0];
2231 1060682 : move16();
2232 1060682 : /* first half partition */ move16();
2233 1060682 : midband[0] = part[0];
2234 1060682 : move16();
2235 1060682 : psize[0] = add( part[0], 1 );
2236 1060682 : move16();
2237 1060682 : psize_inv[0] = getNormReciprocalWord16( psize[0] );
2238 1060682 : move16();
2239 : /* inner partitions */
2240 45076180 : FOR( j = 1; j < npart; j++ )
2241 : {
2242 44015498 : midband[j] = shr( add( add( part[j - 1], 1 ), part[j] ), 1 );
2243 44015498 : move16();
2244 44015498 : psize[j] = sub( part[j], part[j - 1] );
2245 44015498 : move16();
2246 44015498 : psize_inv[j] = getNormReciprocalWord16( psize[j] );
2247 44015498 : move16();
2248 44015498 : if ( GT_16( psize[j], max_psize ) )
2249 : {
2250 13253699 : max_psize = psize[j];
2251 13253699 : move16();
2252 : }
2253 : }
2254 :
2255 1060682 : shift = 9;
2256 1060682 : move16();
2257 1060682 : *psize_norm_exp = sub( 15, shift );
2258 1060682 : move16();
2259 46136862 : FOR( j = 0; j < npart; j++ )
2260 : {
2261 45076180 : psize_norm[j] = shl_o( psize[j], shift, &Overflow ); // Q(15 - psize_norm_exp)
2262 45076180 : move16();
2263 : }
2264 : /* minimum_statistics needs fixed exponent of 6 */
2265 1060682 : assert( norm_s( -max_psize ) >= 9 );
2266 1060682 : }
2267 :
2268 : /*-------------------------------------------------------------------
2269 : * AnalysisSTFT()
2270 : *
2271 : * STFT analysis filterbank
2272 : *-------------------------------------------------------------------*/
2273 :
2274 15 : void AnalysisSTFT(
2275 : const Word16 *timeDomainInput, /* i : pointer to time signal Q(Q)*/
2276 : Word16 Q,
2277 : Word32 *fftBuffer, /* o : FFT bins */
2278 : Word16 *fftBuffer_exp, /* i : exponent of FFT bins */
2279 : HANDLE_FD_CNG_COM hFdCngCom /* i/o: FD_CNG structure containing all buffers and variables */
2280 : )
2281 : {
2282 : Word16 i, len;
2283 : Word16 len2;
2284 : const PWord16 *olapWin;
2285 : Word16 *olapBuffer;
2286 :
2287 :
2288 15 : assert( ( hFdCngCom->fftlen >> 1 ) == hFdCngCom->frameSize );
2289 :
2290 : /* pointer inititialization */
2291 15 : assert( hFdCngCom->olapBufferAna != NULL );
2292 15 : olapBuffer = hFdCngCom->olapBufferAna;
2293 15 : olapWin = hFdCngCom->olapWinAna;
2294 :
2295 : /* olapWin factor is scaled with one bit */
2296 15 : *fftBuffer_exp = 1;
2297 15 : move16();
2298 15 : len = sub( hFdCngCom->fftlen, hFdCngCom->frameSize );
2299 15 : assert( len <= 320 ); /* see size of olapBuffer */
2300 :
2301 : /* Window the signal */
2302 15 : len2 = shr( len, 1 );
2303 2191 : FOR( i = 0; i < len2; i++ )
2304 : {
2305 2176 : move32();
2306 2176 : move32();
2307 2176 : fftBuffer[i] = L_mult( olapBuffer[i], mult_r( olapWin[i].v.im, 23170 /*1.4142135623730950488016887242097 Q14*/ ) );
2308 2176 : fftBuffer[i + len2] = L_mult( olapBuffer[i + len2], mult_r( olapWin[len2 - 1 - i].v.re, 23170 /*1.4142135623730950488016887242097 Q14*/ ) );
2309 : }
2310 15 : len2 = shr( hFdCngCom->frameSize, 1 );
2311 2191 : FOR( i = 0; i < len2; i++ )
2312 : {
2313 2176 : move32();
2314 2176 : move32();
2315 2176 : fftBuffer[i + len] = L_mult( shr_sat( timeDomainInput[i], Q ), mult_r( olapWin[i].v.re, 23170 /*1.4142135623730950488016887242097 Q14*/ ) );
2316 2176 : fftBuffer[i + len + len2] = L_mult( shr_sat( timeDomainInput[i + len2], Q ), mult_r( olapWin[len2 - 1 - i].v.im, 23170 /*1.4142135623730950488016887242097 Q14*/ ) );
2317 : }
2318 :
2319 : /* Perform FFT */
2320 15 : BASOP_rfft( fftBuffer, hFdCngCom->fftlen, fftBuffer_exp, -1 );
2321 :
2322 4367 : FOR( i = 0; i < len; i++ )
2323 : {
2324 4352 : olapBuffer[i] = shr_sat( timeDomainInput[sub( hFdCngCom->frameSize, len ) + i], Q );
2325 4352 : move16();
2326 : }
2327 15 : }
2328 :
2329 : /*-------------------------------------------------------------------
2330 : * AnalysisSTFT_fx()
2331 : *
2332 : * STFT analysis filterbank
2333 : *-------------------------------------------------------------------*/
2334 87961 : void AnalysisSTFT_fx(
2335 : const Word16 *timeDomainInput, // Q(Q)
2336 : Word16 Q,
2337 : Word32 *fftBuffer, /* o : FFT bins */
2338 : Word16 *fftBuffer_exp, /* i : exponent of FFT bins */
2339 : HANDLE_FD_CNG_COM hFdCngCom /* i/o: FD_CNG structure containing all buffers and variables */
2340 : )
2341 : {
2342 : Word16 i, len;
2343 : const Word32 *olapWin;
2344 : Word16 *olapBuffer;
2345 :
2346 87961 : assert( ( hFdCngCom->fftlen >> 1 ) == hFdCngCom->frameSize );
2347 :
2348 : /* pointer inititialization */
2349 87961 : assert( hFdCngCom->olapBufferAna != NULL );
2350 87961 : olapBuffer = hFdCngCom->olapBufferAna_fx;
2351 87961 : olapWin = hFdCngCom->olapWinAna_fx;
2352 :
2353 : /* olapWin factor is scaled with one bit */
2354 87961 : *fftBuffer_exp = 1;
2355 87961 : move16();
2356 87961 : len = sub( hFdCngCom->fftlen, hFdCngCom->frameSize );
2357 87961 : assert( len <= 320 ); /* see size of olapBuffer */
2358 :
2359 : /* Window the signal */
2360 87961 : Copy( olapBuffer + hFdCngCom->frameSize, olapBuffer, len );
2361 :
2362 25480345 : FOR( i = 0; i < hFdCngCom->frameSize; i++ )
2363 : {
2364 25392384 : olapBuffer[i + len] = shr_sat( timeDomainInput[i], Q ); // Values above MAX_16 seen with 10dB tests in float code
2365 25392384 : move16();
2366 : }
2367 :
2368 : /* Window the signal */
2369 87961 : v_L_mult_3216( olapWin, olapBuffer, fftBuffer, hFdCngCom->fftlen );
2370 :
2371 50872729 : FOR( i = 0; i < hFdCngCom->fftlen; i++ )
2372 : {
2373 50784768 : fftBuffer[i] = L_shr( fftBuffer[i], 11 );
2374 50784768 : move32();
2375 : }
2376 87961 : *fftBuffer_exp = WORD16_BITS + 11;
2377 87961 : move16();
2378 :
2379 : /* Perform FFT */
2380 87961 : RFFTN_fx( fftBuffer, hFdCngCom->fftSineTab_fx, hFdCngCom->fftlen, -1 );
2381 :
2382 :
2383 87961 : return;
2384 : }
2385 :
2386 : /*-------------------------------------------------------------------
2387 : * SynthesisSTFT_enc_ivas_fx()
2388 : *
2389 : * STFT synthesis filterbank
2390 : *-------------------------------------------------------------------*/
2391 :
2392 0 : void SynthesisSTFT_enc_ivas_fx(
2393 : Word32 *fftBuffer, /* i : pointer to FFT bins */
2394 : Word16 fftBufferExp, /* i : exponent of FFT bins */
2395 : Word16 *timeDomainOutput, /* o : pointer to time domain signal */
2396 : Word16 *olapBuffer, /* i/o : pointer to overlap buffer */
2397 : const PWord16 *olapWin, /* i : pointer to overlap window */
2398 : Word16 tcx_transition,
2399 : HANDLE_FD_CNG_COM hFdCngCom, /* i/o : pointer to FD_CNG structure containing all buffers and variables */
2400 : Word16 gen_exc,
2401 : Word16 *Q_new, /* i : Q of generated exc_cng */
2402 : const Word16 element_mode, /* i : element mode */
2403 : const Word16 nchan_out /* i : number of output channels */
2404 : )
2405 : {
2406 : Word16 i, len, scale, tmp;
2407 : Word16 len2, len3, len4;
2408 : Word16 buf[M + 1 + L_FRAME16k];
2409 :
2410 :
2411 : /* Perform IFFT */
2412 0 : scale = 0;
2413 0 : move16();
2414 0 : BASOP_rfft( fftBuffer, hFdCngCom->fftlen, &scale, 1 );
2415 0 : fftBufferExp = add( fftBufferExp, scale );
2416 0 : hFdCngCom->fftBuffer_exp = fftBufferExp;
2417 0 : move16();
2418 :
2419 0 : fftBufferExp = add( fftBufferExp, hFdCngCom->fftlenShift );
2420 :
2421 : /* Perform overlap-add */
2422 : /* Handle overlap in P/S domain for stereo */
2423 0 : test();
2424 0 : test();
2425 0 : IF( ( EQ_16( element_mode, IVAS_CPE_TD ) || EQ_16( element_mode, IVAS_CPE_DFT ) ) && EQ_16( nchan_out, 2 ) )
2426 : {
2427 0 : Copy( olapBuffer + 3 * hFdCngCom->frameSize / 4 - ( M + 1 ), buf, hFdCngCom->frameSize + M + 1 );
2428 0 : set16_fx( olapBuffer, 0, hFdCngCom->fftlen );
2429 : }
2430 : ELSE
2431 : {
2432 0 : Copy( olapBuffer + hFdCngCom->frameSize, olapBuffer, hFdCngCom->frameSize );
2433 0 : set16_fx( olapBuffer + hFdCngCom->frameSize, 0, hFdCngCom->frameSize );
2434 : }
2435 0 : len2 = shr( hFdCngCom->fftlen, 2 );
2436 0 : len4 = shr( hFdCngCom->fftlen, 3 );
2437 0 : len3 = add( len2, len4 );
2438 0 : len = add( hFdCngCom->frameSize, len4 );
2439 0 : IF( tcx_transition )
2440 : {
2441 0 : FOR( i = 0; i < len; i++ )
2442 : {
2443 0 : olapBuffer[i] = round_fx_sat( L_shl_sat( fftBuffer[i], sub( fftBufferExp, 15 ) ) ); // Q(15 - fftBufferExp)
2444 0 : move16();
2445 : }
2446 : }
2447 : ELSE
2448 : {
2449 0 : FOR( i = 0; i < len4; i++ )
2450 : {
2451 0 : olapBuffer[i + 1 * len4] = add_sat( olapBuffer[i + 1 * len4], mult_r( round_fx_sat( L_shl_sat( fftBuffer[i + 1 * len4], sub( fftBufferExp, 15 ) ) ), olapWin[i].v.im ) );
2452 0 : move16();
2453 0 : olapBuffer[i + 2 * len4] = add_sat( olapBuffer[i + 2 * len4], mult_r( round_fx_sat( L_shl_sat( fftBuffer[i + 2 * len4], sub( fftBufferExp, 15 ) ) ), olapWin[len4 - 1 - i].v.re ) );
2454 0 : move16();
2455 : }
2456 0 : FOR( i = len3; i < len; i++ )
2457 : {
2458 0 : olapBuffer[i] = round_fx_sat( L_shl_sat( fftBuffer[i], sub( fftBufferExp, 15 ) ) ); // Q(15 - fftBufferExp)
2459 0 : move16();
2460 : }
2461 : }
2462 :
2463 0 : FOR( i = 0; i < len4; i++ )
2464 : {
2465 0 : olapBuffer[i + 5 * len4] = mult_r( round_fx_sat( L_shl_sat( fftBuffer[i + 5 * len4], sub( fftBufferExp, 15 ) ) ), olapWin[i].v.re );
2466 0 : move16();
2467 0 : olapBuffer[i + 6 * len4] = mult_r( round_fx_sat( L_shl_sat( fftBuffer[i + 6 * len4], sub( fftBufferExp, 15 ) ) ), olapWin[len4 - 1 - i].v.im );
2468 0 : move16();
2469 : }
2470 :
2471 0 : len = add( len, len2 );
2472 0 : FOR( i = len; i < hFdCngCom->fftlen; i++ )
2473 : {
2474 0 : olapBuffer[i] = 0;
2475 0 : move16();
2476 : }
2477 :
2478 : /* Get time-domain signal */
2479 0 : FOR( i = 0; i < hFdCngCom->frameSize; i++ )
2480 : {
2481 0 : timeDomainOutput[i] = mult_r( olapBuffer[i + len4], hFdCngCom->fftlenFac );
2482 0 : move16();
2483 : }
2484 : /* Generate excitation */
2485 0 : test();
2486 0 : test();
2487 0 : IF( ( EQ_16( element_mode, IVAS_CPE_TD ) || EQ_16( element_mode, IVAS_CPE_DFT ) ) && EQ_16( nchan_out, 2 ) )
2488 : {
2489 0 : FOR( i = 0; i < hFdCngCom->frameSize / 2; i++ )
2490 : {
2491 0 : buf[i + ( M + 1 )] = add( buf[i + ( M + 1 )], olapBuffer[i + hFdCngCom->frameSize / 4] );
2492 0 : move16();
2493 : }
2494 :
2495 0 : FOR( i = 0; i < M + 1 + hFdCngCom->frameSize; i++ )
2496 : {
2497 0 : buf[i] = mult_r( buf[i], hFdCngCom->fftlenFac );
2498 0 : move16();
2499 : }
2500 : }
2501 : ELSE
2502 : {
2503 0 : FOR( i = 0; i < M + 1 + hFdCngCom->frameSize; i++ )
2504 : {
2505 0 : buf[i] = mult_r( olapBuffer[i + len4 - M - 1], hFdCngCom->fftlenFac );
2506 0 : move16();
2507 : }
2508 0 : tmp = buf[0];
2509 0 : move16();
2510 : }
2511 0 : IF( EQ_16( gen_exc, 1 ) )
2512 : {
2513 :
2514 0 : E_UTIL_f_preemph2( sub( *Q_new, 1 ), buf + 1, PREEMPH_FAC, M + hFdCngCom->frameSize, &tmp );
2515 0 : Residu3_fx( hFdCngCom->A_cng, buf + 1 + M, hFdCngCom->exc_cng, hFdCngCom->frameSize, 1 );
2516 : }
2517 0 : IF( EQ_16( gen_exc, 2 ) )
2518 : {
2519 0 : *Q_new = E_UTIL_f_preemph3( buf + 1, PREEMPH_FAC, M + hFdCngCom->frameSize, &tmp, 1 );
2520 0 : move16();
2521 0 : Residu3_fx( hFdCngCom->A_cng, buf + 1 + M, hFdCngCom->exc_cng, hFdCngCom->frameSize, 1 );
2522 : }
2523 0 : }
2524 :
2525 : /*-------------------------------------------------------------------
2526 : * SynthesisSTFT()
2527 : *
2528 : * STFT synthesis filterbank
2529 : *-------------------------------------------------------------------*/
2530 :
2531 1023 : void SynthesisSTFT(
2532 : Word32 *fftBuffer, /* i : pointer to FFT bins */
2533 : Word16 fftBufferExp, /* i : exponent of FFT bins */
2534 : Word16 *timeDomainOutput, /* o : pointer to time domain signal Qx*/
2535 : Word16 *olapBuffer, /* i/o : pointer to overlap buffer */
2536 : const PWord16 *olapWin, /* i : pointer to overlap window */
2537 : Word16 tcx_transition,
2538 : HANDLE_FD_CNG_COM hFdCngCom, /* i/o : pointer to FD_CNG structure containing all buffers and variables */
2539 : Word16 gen_exc,
2540 : Word16 *Q_new,
2541 : const Word16 element_mode, /* i : element mode */
2542 : const Word16 nchan_out /* i : number of output channels */
2543 : )
2544 : {
2545 : Word16 i, len, scale, tmp;
2546 : Word16 len2, len3, len4;
2547 : Word16 buf[M + 1 + L_FRAME16k];
2548 :
2549 :
2550 : /* Perform IFFT */
2551 1023 : scale = 0;
2552 1023 : move16();
2553 1023 : BASOP_rfft( fftBuffer, hFdCngCom->fftlen, &scale, 1 );
2554 1023 : fftBufferExp = add( fftBufferExp, scale );
2555 1023 : hFdCngCom->fftBuffer_exp = fftBufferExp;
2556 1023 : move16();
2557 :
2558 1023 : fftBufferExp = add( fftBufferExp, hFdCngCom->fftlenShift );
2559 :
2560 : /* Perform overlap-add */
2561 : /* Handle overlap in P/S domain for stereo */
2562 1023 : IF( ( EQ_16( element_mode, IVAS_CPE_TD ) || EQ_16( element_mode, IVAS_CPE_DFT ) ) && EQ_16( nchan_out, 2 ) )
2563 : {
2564 0 : Copy( olapBuffer + 3 * hFdCngCom->frameSize / 4 - ( M + 1 ), buf, hFdCngCom->frameSize + M + 1 );
2565 0 : set16_fx( olapBuffer, 0, hFdCngCom->fftlen );
2566 : }
2567 : ELSE
2568 : {
2569 1023 : Copy( olapBuffer + hFdCngCom->frameSize, olapBuffer, hFdCngCom->frameSize );
2570 1023 : set16_fx( olapBuffer + hFdCngCom->frameSize, 0, hFdCngCom->frameSize );
2571 : }
2572 1023 : len2 = shr( hFdCngCom->fftlen, 2 );
2573 1023 : len4 = shr( hFdCngCom->fftlen, 3 );
2574 1023 : len3 = add( len2, len4 );
2575 1023 : len = add( hFdCngCom->frameSize, len4 );
2576 1023 : IF( tcx_transition )
2577 : {
2578 0 : FOR( i = 0; i < len; i++ )
2579 : {
2580 0 : olapBuffer[i] = round_fx_sat( L_shl_sat( fftBuffer[i], fftBufferExp - 15 ) ); // Q(15 - fftBufferExp)
2581 0 : move16();
2582 : }
2583 : }
2584 : ELSE
2585 : {
2586 66495 : FOR( i = 0; i < len4; i++ )
2587 : {
2588 65472 : olapBuffer[i + 1 * len4] = add_sat( olapBuffer[i + 1 * len4], mult_r( round_fx_sat( L_shl_sat( fftBuffer[i + 1 * len4], fftBufferExp - 15 ) ), olapWin[i].v.im ) );
2589 65472 : move16();
2590 65472 : olapBuffer[i + 2 * len4] = add_sat( olapBuffer[i + 2 * len4], mult_r( round_fx_sat( L_shl_sat( fftBuffer[i + 2 * len4], fftBufferExp - 15 ) ), olapWin[len4 - 1 - i].v.re ) );
2591 65472 : move16();
2592 : }
2593 131967 : FOR( i = len3; i < len; i++ )
2594 : {
2595 130944 : olapBuffer[i] = round_fx_sat( L_shl_sat( fftBuffer[i], fftBufferExp - 15 ) );
2596 : }
2597 : }
2598 :
2599 66495 : FOR( i = 0; i < len4; i++ )
2600 : {
2601 65472 : olapBuffer[i + 5 * len4] = mult_r( round_fx_sat( L_shl_sat( fftBuffer[i + 5 * len4], fftBufferExp - 15 ) ), olapWin[i].v.re );
2602 65472 : move16();
2603 65472 : olapBuffer[i + 6 * len4] = mult_r( round_fx_sat( L_shl_sat( fftBuffer[i + 6 * len4], fftBufferExp - 15 ) ), olapWin[len4 - 1 - i].v.im );
2604 65472 : move16();
2605 : }
2606 :
2607 1023 : len = add( len, len2 );
2608 66495 : FOR( i = len; i < hFdCngCom->fftlen; i++ )
2609 : {
2610 65472 : olapBuffer[i] = 0;
2611 65472 : move16();
2612 : }
2613 :
2614 : /* Get time-domain signal */
2615 262911 : FOR( i = 0; i < hFdCngCom->frameSize; i++ )
2616 : {
2617 261888 : timeDomainOutput[i] = mult_r( olapBuffer[i + len4], hFdCngCom->fftlenFac );
2618 261888 : move16();
2619 : }
2620 : /* Generate excitation */
2621 : {
2622 280302 : FOR( i = 0; i < M + 1 + hFdCngCom->frameSize; i++ )
2623 : {
2624 279279 : buf[i] = mult_r( olapBuffer[i + len4 - M - 1], hFdCngCom->fftlenFac );
2625 279279 : move16();
2626 : }
2627 1023 : tmp = buf[0];
2628 : }
2629 1023 : IF( EQ_16( gen_exc, 1 ) )
2630 : {
2631 :
2632 0 : E_UTIL_f_preemph2( *Q_new - 1, buf + 1, PREEMPH_FAC, M + hFdCngCom->frameSize, &tmp );
2633 0 : Residu3_fx( hFdCngCom->A_cng, buf + 1 + M, hFdCngCom->exc_cng, hFdCngCom->frameSize, 1 );
2634 : }
2635 1023 : IF( EQ_16( gen_exc, 2 ) )
2636 : {
2637 0 : *Q_new = E_UTIL_f_preemph3( buf + 1, PREEMPH_FAC, M + hFdCngCom->frameSize, &tmp, 1 );
2638 0 : Residu3_fx( hFdCngCom->A_cng, buf + 1 + M, hFdCngCom->exc_cng, hFdCngCom->frameSize, 1 );
2639 : }
2640 1023 : }
2641 :
2642 18381 : void SynthesisSTFT_ivas_fx(
2643 : Word32 *fftBuffer, /* i : pointer to FFT bins */
2644 : Word16 fftBufferExp, /* i : exponent of FFT bins */
2645 : Word16 *timeDomainOutput, /* o : pointer to time domain signal Qx*/
2646 : Word16 *olapBuffer, /* i/o : pointer to overlap buffer */
2647 : const PWord16 *olapWin, /* i : pointer to overlap window */
2648 : Word16 tcx_transition,
2649 : HANDLE_FD_CNG_COM hFdCngCom, /* i/o : pointer to FD_CNG structure containing all buffers and variables */
2650 : Word16 gen_exc,
2651 : Word16 *Q_new,
2652 : const Word16 element_mode, /* i : element mode */
2653 : const Word16 nchan_out /* i : number of output channels */
2654 : )
2655 : {
2656 : Word16 i, len, scale, tmp;
2657 : Word16 len2, len3, len4;
2658 : Word16 buf[M + 1 + L_FRAME16k];
2659 :
2660 :
2661 : /* Perform IFFT */
2662 18381 : scale = 0;
2663 18381 : move16();
2664 18381 : BASOP_rfft( fftBuffer, hFdCngCom->fftlen, &scale, 1 );
2665 18381 : fftBufferExp = add( fftBufferExp, scale );
2666 18381 : hFdCngCom->fftBuffer_exp = fftBufferExp;
2667 18381 : move16();
2668 :
2669 18381 : fftBufferExp = add( fftBufferExp, hFdCngCom->fftlenShift );
2670 :
2671 : /* Perform overlap-add */
2672 : /* Handle overlap in P/S domain for stereo */
2673 18381 : test();
2674 18381 : test();
2675 18381 : IF( ( EQ_16( element_mode, IVAS_CPE_TD ) || EQ_16( element_mode, IVAS_CPE_DFT ) ) && EQ_16( nchan_out, 2 ) )
2676 : {
2677 0 : Copy( olapBuffer + 3 * hFdCngCom->frameSize / 4 - ( M + 1 ), buf, hFdCngCom->frameSize + M + 1 );
2678 0 : set16_fx( olapBuffer, 0, hFdCngCom->fftlen );
2679 : }
2680 : ELSE
2681 : {
2682 18381 : Copy( olapBuffer + hFdCngCom->frameSize, olapBuffer, hFdCngCom->frameSize );
2683 18381 : set16_fx( olapBuffer + hFdCngCom->frameSize, 0, hFdCngCom->frameSize );
2684 : }
2685 18381 : len2 = shr( hFdCngCom->fftlen, 2 );
2686 18381 : len4 = shr( hFdCngCom->fftlen, 3 );
2687 18381 : len3 = add( len2, len4 );
2688 18381 : len = add( hFdCngCom->frameSize, len4 );
2689 18381 : IF( tcx_transition )
2690 : {
2691 606892 : FOR( i = 0; i < len; i++ )
2692 : {
2693 605120 : olapBuffer[i] = round_fx_sat( L_shl_sat( fftBuffer[i], fftBufferExp - 15 ) ); // Q(15 - fftBufferExp)
2694 605120 : move16();
2695 : }
2696 : }
2697 : ELSE
2698 : {
2699 1256353 : FOR( i = 0; i < len4; i++ )
2700 : {
2701 1239744 : olapBuffer[i + 1 * len4] = add_sat( olapBuffer[i + 1 * len4], mult_r( round_fx_sat( L_shl_sat( fftBuffer[i + 1 * len4], fftBufferExp - 15 ) ), olapWin[i].v.im ) );
2702 1239744 : move16();
2703 1239744 : olapBuffer[i + 2 * len4] = add_sat( olapBuffer[i + 2 * len4], mult_r( round_fx_sat( L_shl_sat( fftBuffer[i + 2 * len4], fftBufferExp - 15 ) ), olapWin[len4 - 1 - i].v.re ) );
2704 1239744 : move16();
2705 : }
2706 2496097 : FOR( i = len3; i < len; i++ )
2707 : {
2708 2479488 : olapBuffer[i] = round_fx_sat( L_shl_sat( fftBuffer[i], fftBufferExp - 15 ) );
2709 : }
2710 : }
2711 :
2712 1379149 : FOR( i = 0; i < len4; i++ )
2713 : {
2714 1360768 : olapBuffer[i + 5 * len4] = mult_r( round_fx_sat( L_shl_sat( fftBuffer[i + 5 * len4], fftBufferExp - 15 ) ), olapWin[i].v.re );
2715 1360768 : move16();
2716 1360768 : olapBuffer[i + 6 * len4] = mult_r( round_fx_sat( L_shl_sat( fftBuffer[i + 6 * len4], fftBufferExp - 15 ) ), olapWin[len4 - 1 - i].v.im );
2717 1360768 : move16();
2718 : }
2719 :
2720 18381 : len = add( len, len2 );
2721 1379149 : FOR( i = len; i < hFdCngCom->fftlen; i++ )
2722 : {
2723 1360768 : olapBuffer[i] = 0;
2724 1360768 : move16();
2725 : }
2726 :
2727 : /* Get time-domain signal */
2728 5461453 : FOR( i = 0; i < hFdCngCom->frameSize; i++ )
2729 : {
2730 5443072 : timeDomainOutput[i] = mult_r( olapBuffer[i + len4], hFdCngCom->fftlenFac );
2731 5443072 : move16();
2732 : }
2733 : /* Generate excitation */
2734 : {
2735 5773930 : FOR( i = 0; i < M + 1 + hFdCngCom->frameSize; i++ )
2736 : {
2737 5755549 : buf[i] = mult_r( olapBuffer[i + len4 - M - 1], hFdCngCom->fftlenFac );
2738 5755549 : move16();
2739 : }
2740 18381 : tmp = buf[0];
2741 : }
2742 18381 : IF( EQ_16( gen_exc, 1 ) )
2743 : {
2744 18381 : Word16 s = getScaleFactor16( buf + 1, M + hFdCngCom->frameSize );
2745 18381 : if ( GT_16( *Q_new, s ) )
2746 : {
2747 6413 : *Q_new = s;
2748 6413 : move16();
2749 : }
2750 :
2751 18381 : E_UTIL_f_preemph2( *Q_new - 1, buf + 1, PREEMPH_FAC, M + hFdCngCom->frameSize, &tmp );
2752 18381 : Residu3_fx( hFdCngCom->A_cng, buf + 1 + M, hFdCngCom->exc_cng, hFdCngCom->frameSize, 1 );
2753 : }
2754 18381 : IF( EQ_16( gen_exc, 2 ) )
2755 : {
2756 0 : *Q_new = E_UTIL_f_preemph3( buf + 1, PREEMPH_FAC, M + hFdCngCom->frameSize, &tmp, 1 );
2757 0 : Residu3_fx( hFdCngCom->A_cng, buf + 1 + M, hFdCngCom->exc_cng, hFdCngCom->frameSize, 1 );
2758 : }
2759 18381 : }
2760 :
2761 : /**************************************************************************************
2762 : * Compute some values used in the bias correction of the minimum statistics algorithm *
2763 : **************************************************************************************/
2764 12860 : void mhvals(
2765 : const Word16 d,
2766 : Word16 *m /*, Word16 * h*/
2767 : )
2768 : {
2769 : Word16 i, j;
2770 12860 : Word16 len = SIZE_SCALE_TABLE_CN;
2771 12860 : move16();
2772 :
2773 :
2774 12860 : assert( d == 72 || d == 12 ); /* function only tested for d==72 and d==12) */
2775 12860 : i = 0;
2776 12860 : move16();
2777 109310 : FOR( i = 0; i < len; i++ )
2778 : {
2779 109310 : IF( LE_16( d, d_array[i] ) )
2780 : {
2781 12860 : BREAK;
2782 : }
2783 : }
2784 12860 : IF( EQ_16( i, len ) )
2785 : {
2786 0 : i = sub( len, 1 );
2787 0 : j = i;
2788 0 : move16();
2789 : }
2790 : ELSE
2791 : {
2792 12860 : j = sub( i, 1 );
2793 : }
2794 12860 : IF( EQ_16( d, d_array[i] ) )
2795 : {
2796 0 : *m = m_array[i];
2797 0 : move16();
2798 : }
2799 : ELSE
2800 : {
2801 : Word32 qi_m, qj_m, q_m, tmp1_m, tmp2_m;
2802 : Word16 qi_e, qj_e, q_e, tmp1_e, tmp2_e, tmp1_w16_m, tmp1_w16_e, shift;
2803 :
2804 :
2805 : /* d_array has exponent 15 */
2806 12860 : qj_e = 15;
2807 12860 : move16();
2808 12860 : qj_m = L_deposit_h( d_array[i - 1] );
2809 :
2810 12860 : qi_e = 15;
2811 12860 : move16();
2812 12860 : qi_m = L_deposit_h( d_array[i] );
2813 :
2814 12860 : q_e = 15;
2815 12860 : move16();
2816 12860 : q_m = L_deposit_h( d );
2817 :
2818 12860 : qj_m = Sqrt32( qj_m, &qj_e );
2819 12860 : qi_m = Sqrt32( qi_m, &qi_e );
2820 12860 : q_m = Sqrt32( q_m, &q_e );
2821 :
2822 12860 : tmp1_m = Mpy_32_32( qi_m, qj_m );
2823 12860 : tmp1_e = add( qi_e, qj_e );
2824 12860 : tmp1_m = L_deposit_h( BASOP_Util_Divide3232_Scale( tmp1_m, q_m, &shift ) );
2825 12860 : tmp1_e = sub( tmp1_e, q_e );
2826 12860 : tmp1_e = add( tmp1_e, shift );
2827 12860 : tmp1_m = BASOP_Util_Add_Mant32Exp( tmp1_m, tmp1_e, L_negate( qj_m ), qj_e, &tmp1_e );
2828 :
2829 12860 : tmp2_m = BASOP_Util_Add_Mant32Exp( qi_m, qi_e, L_negate( qj_m ), qj_e, &tmp2_e );
2830 12860 : tmp1_w16_m = round_fx( tmp2_m );
2831 12860 : tmp1_w16_e = tmp2_e;
2832 12860 : move16();
2833 12860 : BASOP_Util_Divide_MantExp( sub( m_array[j], m_array[i] ), 0, tmp1_w16_m, tmp1_w16_e, &tmp1_w16_m, &tmp1_w16_e );
2834 :
2835 12860 : tmp2_m = Mpy_32_16_1( tmp1_m, tmp1_w16_m );
2836 12860 : tmp2_e = add( tmp1_e, tmp1_w16_e );
2837 :
2838 12860 : tmp2_m = BASOP_Util_Add_Mant32Exp( tmp2_m, tmp2_e, L_deposit_h( m_array[i] ), 0, &tmp2_e );
2839 12860 : assert( tmp2_e == 0 );
2840 12860 : *m = extract_h( tmp2_m );
2841 12860 : move32();
2842 : }
2843 12860 : }
2844 :
2845 : /*-------------------------------------------------------------------
2846 : * rand_gauss()
2847 : *
2848 : * Random generator with Gaussian distribution with mean 0 and std 1
2849 : * Returns:
2850 : * random signal format Q3.29
2851 : *-------------------------------------------------------------------*/
2852 40269655 : Word32 rand_gauss( Word16 *seed )
2853 : {
2854 : Word32 temp;
2855 : Word16 loc_seed;
2856 :
2857 : /* This unrolled version reduces the cycles from 17 to 10 */
2858 40269655 : loc_seed = extract_l( L_mac0( 13849, *seed, 31821 ) );
2859 40269655 : temp = L_deposit_l( loc_seed );
2860 :
2861 40269655 : loc_seed = extract_l( L_mac0( 13849, loc_seed, 31821 ) );
2862 40269655 : temp = L_msu0( temp, loc_seed, -1 );
2863 :
2864 40269655 : loc_seed = extract_l( L_mac0( 13849, loc_seed, 31821 ) );
2865 40269655 : temp = L_msu0( temp, loc_seed, -1 );
2866 :
2867 40269655 : *seed = loc_seed;
2868 40269655 : move16();
2869 40269655 : return L_shl( temp, WORD16_BITS - CNG_RAND_GAUSS_SHIFT );
2870 : }
2871 :
2872 : /*-------------------------------------------------------------------
2873 : * lpc_from_spectrum()
2874 : *
2875 : *
2876 : *-------------------------------------------------------------------*/
2877 :
2878 3802 : void lpc_from_spectrum(
2879 : HANDLE_FD_CNG_COM hFdCngCom,
2880 : const Word16 start, /*i : start band*/
2881 : const Word16 stop, /*i : stop band*/
2882 : const Word16 preemph_fac /*i : preemphase factor format Q1.15*/
2883 : )
2884 : {
2885 : Word16 i, s1, s2, s3, fftlen2, scale, fftlen4, fftlen8, len, step, preemph_fac2;
2886 : Word32 maxVal, r[32], fftBuffer[FFTLEN], *ptr, *pti, nf;
2887 : Word16 tmp, r_h[32], r_l[32];
2888 : const PWord16 *table;
2889 :
2890 3802 : Word32 *powspec = hFdCngCom->cngNoiseLevel; /*i : pointer to noise levels format Q5.27*/
2891 3802 : Word16 powspec_exp = hFdCngCom->cngNoiseLevelExp;
2892 3802 : move16();
2893 3802 : Word16 fftlen = hFdCngCom->fftlen; /*i : size of fft*/
2894 3802 : Word16 *A = hFdCngCom->A_cng; /*o : lpc coefficients format Q3.12*/
2895 3802 : move16();
2896 3802 : Word16 lpcorder = M;
2897 3802 : move16();
2898 :
2899 3802 : scale = 0;
2900 3802 : move16();
2901 3802 : fftlen2 = shr( fftlen, 1 );
2902 3802 : fftlen4 = shr( fftlen, 2 );
2903 3802 : fftlen8 = shr( fftlen, 3 );
2904 :
2905 : /* Power Spectrum */
2906 3802 : maxVal = 0;
2907 3802 : move32();
2908 3802 : len = sub( stop, start );
2909 1121574 : FOR( i = 0; i < len; i++ )
2910 : {
2911 1117772 : maxVal = L_max( maxVal, L_abs( powspec[i] ) );
2912 : }
2913 3802 : s1 = norm_l( maxVal );
2914 3802 : nf = L_shr_r_sat( 1099511680l /*1e-3f Q40*/, add( sub( powspec_exp, s1 ), 9 ) );
2915 3802 : ptr = fftBuffer;
2916 3802 : pti = fftBuffer + 1;
2917 :
2918 11406 : FOR( i = 0; i < start; i++ )
2919 : {
2920 7604 : *ptr = nf;
2921 7604 : move32();
2922 7604 : *pti = L_deposit_l( 0 );
2923 7604 : move32();
2924 7604 : ptr += 2;
2925 7604 : pti += 2;
2926 : }
2927 :
2928 1121574 : FOR( ; i < stop; i++ )
2929 : {
2930 1117772 : *ptr = L_max( nf, L_shl( powspec[i - start], s1 ) );
2931 1117772 : move32();
2932 1117772 : *pti = L_deposit_l( 0 );
2933 1117772 : move32();
2934 1117772 : ptr += 2;
2935 1117772 : pti += 2;
2936 : }
2937 :
2938 3994 : FOR( ; i < fftlen2; i++ )
2939 : {
2940 192 : *ptr = nf;
2941 192 : move32();
2942 192 : *pti = L_deposit_l( 0 );
2943 192 : move32();
2944 192 : ptr += 2;
2945 192 : pti += 2;
2946 : }
2947 :
2948 3802 : fftBuffer[1] = nf;
2949 3802 : move32();
2950 :
2951 : /* Pre-emphasis */
2952 :
2953 3802 : BASOP_getTables( &table, NULL, &step, fftlen4 );
2954 3802 : tmp = round_fx( L_shr( L_add( 0x40000000, L_mult0( preemph_fac, preemph_fac ) ), 1 ) );
2955 3802 : preemph_fac2 = shr( preemph_fac, 1 );
2956 3802 : ptr = fftBuffer;
2957 3802 : *ptr = Mpy_32_16_1( *ptr, sub( tmp, preemph_fac2 ) );
2958 3802 : move32();
2959 3802 : ptr += 2;
2960 281392 : FOR( i = 1; i < fftlen8; i++ )
2961 : {
2962 277590 : move32();
2963 277590 : *ptr = Mpy_32_16_1( *ptr, sub( tmp, mult_r( preemph_fac2, add( shr( table[i - 1].v.re, 1 ), shr( table[i].v.re, 1 ) ) ) ) );
2964 277590 : ptr += 2;
2965 : }
2966 3802 : move32();
2967 3802 : *ptr = Mpy_32_16_1( *ptr, sub( tmp, mult_r( preemph_fac2, add( shr( table[fftlen8 - 1].v.re, 1 ), shr( table[fftlen8 - 1].v.im, 1 ) ) ) ) );
2968 3802 : ptr += 2;
2969 281392 : FOR( i = 1; i < fftlen8; i++ )
2970 : {
2971 277590 : move32();
2972 277590 : *ptr = Mpy_32_16_1( *ptr, sub( tmp, mult_r( preemph_fac2, add( shr( table[fftlen8 - i - 1].v.im, 1 ), shr( table[fftlen8 - i].v.im, 1 ) ) ) ) );
2973 277590 : ptr += 2;
2974 : }
2975 3802 : move32();
2976 3802 : *ptr = Mpy_32_16_1( *ptr, tmp );
2977 3802 : ptr += 2;
2978 281392 : FOR( i = 1; i < fftlen8; i++ )
2979 : {
2980 277590 : move32();
2981 277590 : *ptr = Mpy_32_16_1( *ptr, add( tmp, mult_r( preemph_fac2, add( shr( table[i - 1].v.im, 1 ), shr( table[i].v.im, 1 ) ) ) ) );
2982 277590 : ptr += 2;
2983 : }
2984 3802 : move32();
2985 3802 : *ptr = Mpy_32_16_1( *ptr, add( tmp, mult_r( preemph_fac2, add( shr( table[fftlen8 - 1].v.re, 1 ), shr( table[fftlen8 - 1].v.im, 1 ) ) ) ) );
2986 3802 : ptr += 2;
2987 281392 : FOR( i = 1; i < fftlen8; i++ )
2988 : {
2989 277590 : move32();
2990 277590 : *ptr = Mpy_32_16_1( *ptr, add( tmp, mult_r( preemph_fac2, add( shr( table[fftlen8 - i - 1].v.re, 1 ), shr( table[fftlen8 - i].v.re, 1 ) ) ) ) );
2991 277590 : ptr += 2;
2992 : }
2993 3802 : move32();
2994 3802 : fftBuffer[1] = Mpy_32_16_1( fftBuffer[1], add( tmp, preemph_fac2 ) );
2995 3802 : maxVal = 0;
2996 3802 : move32();
2997 2254938 : FOR( i = 0; i < fftlen; i++ )
2998 : {
2999 2251136 : maxVal = L_max( maxVal, L_abs( fftBuffer[i] ) );
3000 : }
3001 3802 : s2 = norm_l( maxVal );
3002 2254938 : FOR( i = 0; i < fftlen; i++ )
3003 : {
3004 2251136 : fftBuffer[i] = L_shl( fftBuffer[i], s2 );
3005 2251136 : move32();
3006 : }
3007 :
3008 : /* Autocorrelation */
3009 :
3010 3802 : BASOP_rfft( fftBuffer, fftlen, &scale, 1 );
3011 :
3012 3802 : s3 = getScaleFactor32( fftBuffer, add( lpcorder, 1 ) );
3013 :
3014 68436 : FOR( i = 0; i <= lpcorder; i++ )
3015 : {
3016 64634 : r[i] = L_shl( fftBuffer[i], s3 );
3017 64634 : move32();
3018 : }
3019 :
3020 3802 : r[0] = Mpy_32_32( r[0], 1074278656l /*1.0005f Q30*/ );
3021 3802 : move32();
3022 64634 : FOR( i = 1; i <= lpcorder; i++ )
3023 : {
3024 60832 : r[i] = Mpy_32_32( r[i], 1073741824l /*1.f Q30*/ );
3025 60832 : move32();
3026 : }
3027 3802 : s3 = getScaleFactor32( r, add( lpcorder, 1 ) );
3028 :
3029 68436 : FOR( i = 0; i <= lpcorder; i++ )
3030 : {
3031 64634 : r[i] = L_shl( r[i], s3 );
3032 64634 : move32();
3033 : }
3034 :
3035 68436 : FOR( i = 0; i <= lpcorder; i++ )
3036 : {
3037 64634 : L_Extract( r[i], &r_h[i], &r_l[i] );
3038 : }
3039 :
3040 : /* LPC */
3041 :
3042 3802 : E_LPC_lev_dur( r_h, r_l, A, NULL, lpcorder, NULL );
3043 3802 : }
3044 :
3045 : /*
3046 : msvq_decoder
3047 :
3048 : Parameters:
3049 :
3050 : cb i : Codebook (indexed cb[stages][levels][p]) format Q9.7
3051 : stages i : Number of stages
3052 : N i : Vector dimension
3053 : maxN i : Codebook vector dimension
3054 : Idx o : Indices
3055 : uq[] i : Quantized vector format Q9.7
3056 :
3057 :
3058 : Function:
3059 : multi stage vector dequantisation
3060 :
3061 : Returns:
3062 : void
3063 : */
3064 0 : void msvq_decoder(
3065 : const Word16 *const cb[], /* i : Codebook (indexed cb[*stages][levels][p]) */
3066 : const Word16 stages, /* i : Number of stages */
3067 : const Word16 N, /* i : Vector dimension */
3068 : const Word16 maxN, /* i : Codebook vector dimension */
3069 : const Word16 Idx[], /* i : Indices */
3070 : Word16 *uq /* o : quantized vector */
3071 : )
3072 : {
3073 : Word16 s, i, offset;
3074 :
3075 : // PMT("msvq_decoder Not verified")
3076 :
3077 0 : offset = i_mult( Idx[0], maxN );
3078 0 : FOR( i = 0; i < N; i++ )
3079 : {
3080 0 : uq[i] = cb[0][offset + i];
3081 0 : move16();
3082 : }
3083 :
3084 0 : FOR( s = 1; s < stages; s++ )
3085 : {
3086 0 : offset = i_mult( Idx[s], maxN );
3087 :
3088 0 : FOR( i = 0; i < N; i++ )
3089 : {
3090 0 : uq[i] = add( uq[i], cb[s][offset + i] );
3091 0 : move16();
3092 : }
3093 : }
3094 0 : }
3095 : /*-------------------------------------------------------------------
3096 : * FdCng_exc()
3097 : *
3098 : * Generate FD-CNG as LP excitation
3099 : *-------------------------------------------------------------------*/
3100 :
3101 18381 : void FdCng_exc(
3102 : HANDLE_FD_CNG_COM hFdCngCom,
3103 : Word16 *CNG_mode,
3104 : Word16 L_frame,
3105 : Word16 *lsp_old,
3106 : Word16 first_CNG,
3107 : Word16 *lspCNG,
3108 : Word16 *Aq, /* o: LPC coeffs Q12*/
3109 : Word16 *lsp_new, /* o: lsp Q15 */
3110 : Word16 *lsf_new, /* o: lsf Qlog2(2.56) */
3111 : Word16 *exc, /* o: LP excitation Q12 */
3112 : Word16 *exc2, /* o: LP excitation Q12 */
3113 : Word16 *bwe_exc /* o: LP excitation for BWE Q12*/
3114 : )
3115 : {
3116 : Word16 i;
3117 :
3118 18381 : *CNG_mode = -1;
3119 18381 : move16();
3120 :
3121 103429 : FOR( i = 0; i < L_frame / L_SUBFR; i++ )
3122 : {
3123 85048 : Copy( hFdCngCom->A_cng, Aq + i * ( M + 1 ), M + 1 );
3124 : }
3125 :
3126 18381 : E_LPC_a_lsp_conversion( Aq, lsp_new, lsp_old, M );
3127 :
3128 18381 : IF( first_CNG == 0 )
3129 : {
3130 7447 : Copy( lsp_old, lspCNG, M );
3131 : }
3132 312477 : FOR( i = 0; i < M; i++ )
3133 : {
3134 : /* AR low-pass filter */
3135 294096 : lspCNG[i] = mac_r( L_mult( CNG_ISF_FACT_FX, lspCNG[i] ), 32768 - CNG_ISF_FACT_FX, lsp_new[i] );
3136 294096 : move16(); /* Q15 (15+15+1-16) */
3137 : }
3138 :
3139 18381 : IF( EQ_16( L_frame, L_FRAME16k ) )
3140 : {
3141 11524 : lsp2lsf_fx( lsp_new, lsf_new, M, INT_FS_16k_FX );
3142 : }
3143 : ELSE
3144 : {
3145 6857 : E_LPC_lsp_lsf_conversion( lsp_new, lsf_new, M );
3146 : }
3147 18381 : Copy( hFdCngCom->exc_cng, exc, L_frame );
3148 18381 : Copy( hFdCngCom->exc_cng, exc2, L_frame );
3149 18381 : IF( bwe_exc != NULL )
3150 : {
3151 11158 : IF( EQ_16( L_frame, L_FRAME ) )
3152 : {
3153 6857 : interp_code_5over2_fx( exc2, bwe_exc, L_frame );
3154 : }
3155 : ELSE
3156 : {
3157 4301 : interp_code_4over2_fx( exc2, bwe_exc, L_frame );
3158 : }
3159 : }
3160 18381 : }
3161 :
3162 : /*-------------------------------------------------------------------
3163 : * SynthesisSTFT_fx()
3164 : *
3165 : * STFT synthesis filterbank
3166 : *-------------------------------------------------------------------*/
3167 :
3168 59366 : void SynthesisSTFT_fx(
3169 : Word32 *fftBuffer,
3170 : /* i : FFT bins */ // Q15
3171 : Word32 *timeDomainOutput,
3172 : Word32 *olapBuffer, // Qin
3173 : const Word16 *olapWin,
3174 : const Word16 tcx_transition,
3175 : HANDLE_FD_CNG_COM hFdCngCom, /* i/o: FD_CNG structure containing all buffers and variables */
3176 : const Word16 element_mode, /* i : element mode */
3177 : const Word16 nchan_out /* i : number of output channels */
3178 : )
3179 : {
3180 : Word16 i;
3181 : Word32 buf_fx[M + 1 + 320], tmp_fx;
3182 :
3183 : /* Perform IFFT */
3184 59366 : RFFTN_fx( fftBuffer, hFdCngCom->fftSineTab_fx, hFdCngCom->fftlen, 1 );
3185 :
3186 : /* Handle overlap in P/S domain for stereo */
3187 59366 : test();
3188 59366 : test();
3189 59366 : IF( ( EQ_16( element_mode, IVAS_CPE_TD ) || EQ_16( element_mode, IVAS_CPE_DFT ) ) && EQ_16( nchan_out, 2 ) )
3190 : {
3191 876 : Copy32( olapBuffer + sub( i_mult( 3, shr( hFdCngCom->frameSize, 2 ) ), ( M + 1 ) ), buf_fx, add( hFdCngCom->frameSize, M + 1 ) ); // Qin
3192 876 : set32_fx( olapBuffer, 0, hFdCngCom->fftlen );
3193 : }
3194 : ELSE
3195 : {
3196 58490 : Copy32( olapBuffer + hFdCngCom->frameSize, olapBuffer, hFdCngCom->frameSize );
3197 58490 : set32_fx( olapBuffer + hFdCngCom->frameSize, 0, hFdCngCom->frameSize ); /*olapBuffer, fftBuffer, olapWin*/
3198 : }
3199 :
3200 59366 : IF( tcx_transition )
3201 : {
3202 0 : FOR( i = 0; i < ( 5 * ( hFdCngCom->frameSize >> 2 ) ); i++ )
3203 : {
3204 0 : olapBuffer[i] = fftBuffer[i]; // Q15
3205 0 : move32();
3206 : }
3207 : }
3208 : ELSE
3209 : {
3210 8535718 : FOR( i = hFdCngCom->frameSize / 4; i < ( 3 * ( hFdCngCom->frameSize >> 2 ) ); i++ )
3211 : {
3212 8476352 : olapBuffer[i] = L_add( olapBuffer[i], Mpy_32_16_1( fftBuffer[i], olapWin[( i - ( hFdCngCom->frameSize >> 2 ) )] ) );
3213 8476352 : move32();
3214 : }
3215 8535718 : FOR( ; i < ( 5 * ( hFdCngCom->frameSize >> 2 ) ); i++ )
3216 : {
3217 8476352 : olapBuffer[i] = fftBuffer[i];
3218 8476352 : move32();
3219 : }
3220 : }
3221 8535718 : FOR( ; i < ( 7 * ( hFdCngCom->frameSize >> 2 ) ); i++ )
3222 : {
3223 8476352 : olapBuffer[i] = Mpy_32_16_1( fftBuffer[i], olapWin[( i - ( 3 * ( hFdCngCom->frameSize >> 2 ) ) )] );
3224 8476352 : move32();
3225 : }
3226 :
3227 4297542 : FOR( ; i < hFdCngCom->fftlen; i++ )
3228 : {
3229 4238176 : olapBuffer[i] = 0;
3230 4238176 : move32();
3231 : }
3232 :
3233 59366 : Word32 fftScale = 0;
3234 59366 : SWITCH( hFdCngCom->fftlen )
3235 : {
3236 27422 : case 640:
3237 27422 : fftScale = FFT_SCALING_640;
3238 27422 : move32();
3239 27422 : BREAK;
3240 31944 : case 512:
3241 31944 : fftScale = FFT_SCALING_512;
3242 31944 : move32();
3243 31944 : BREAK;
3244 0 : default:
3245 0 : assert( !"Not supported FFT length!" );
3246 : }
3247 : /* Get time-domain signal */
3248 : // v_multc(olapBuffer + hFdCngCom->frameSize / 4, (float)(hFdCngCom->fftlen / 2), timeDomainOutput, hFdCngCom->frameSize);
3249 59366 : v_multc_fixed( olapBuffer + hFdCngCom->frameSize / 4, fftScale, timeDomainOutput, hFdCngCom->frameSize ); // Q_in - 9
3250 : /* Get excitation */
3251 59366 : test();
3252 59366 : test();
3253 59366 : IF( ( EQ_16( element_mode, IVAS_CPE_TD ) || EQ_16( element_mode, IVAS_CPE_DFT ) ) && EQ_16( nchan_out, 2 ) )
3254 : {
3255 113004 : FOR( i = 0; i < hFdCngCom->frameSize / 2; i++ )
3256 : {
3257 112128 : buf_fx[i + ( M + 1 )] = L_add( buf_fx[i + ( M + 1 )], olapBuffer[( i + ( hFdCngCom->frameSize >> 2 ) )] );
3258 112128 : move32();
3259 : }
3260 : // v_multc(buf, (float)(hFdCngCom->fftlen / 2), buf, M + 1 + hFdCngCom->frameSize);
3261 876 : v_multc_fixed( buf_fx, fftScale, buf_fx, add( M + 1, hFdCngCom->frameSize ) );
3262 : }
3263 : ELSE
3264 : {
3265 : // v_multc(olapBuffer + hFdCngCom->frameSize / 4 - (M + 1), (float)(hFdCngCom->fftlen / 2), buf, M + 1 + hFdCngCom->frameSize);
3266 58490 : v_multc_fixed( olapBuffer + sub( shr( hFdCngCom->frameSize, 2 ), ( M + 1 ) ), fftScale, buf_fx, add( M + 1, hFdCngCom->frameSize ) );
3267 : }
3268 :
3269 59366 : tmp_fx = buf_fx[0];
3270 59366 : move32();
3271 : // preemph(buf + 1, PREEMPH_FAC_FLT, M + hFdCngCom->frameSize, &tmp);
3272 59366 : preemph_ivas_fx( buf_fx + 1, PREEMPH_FAC, add( M, hFdCngCom->frameSize ), &tmp_fx );
3273 : // residu(hFdCngCom->A_cng_flt, M, buf + 1 + M, hFdCngCom->exc_cng_flt, hFdCngCom->frameSize);
3274 :
3275 : // residu_ivas_fx( hFdCngCom->A_cng, Q13, M, buf_fx + 1 + M, hFdCngCom->exc_cng_32fx, hFdCngCom->frameSize );
3276 59366 : residu_ivas_fx( hFdCngCom->A_cng, ( 15 - norm_s( hFdCngCom->A_cng[0] - 1 ) ), M, buf_fx + 1 + M, hFdCngCom->exc_cng_32fx, hFdCngCom->frameSize );
3277 :
3278 59366 : return;
3279 : }
3280 :
3281 4095 : void SynthesisSTFT_dirac_fx(
3282 : Word32 *fftBuffer,
3283 : /* i : FFT bins */ // hFdCngCom->fftBuffer_exp
3284 : Word32 *timeDomainOutput,
3285 : Word32 *olapBuffer, // Q_in
3286 : const Word16 *olapWin,
3287 : const Word16 samples_out,
3288 : HANDLE_FD_CNG_COM hFdCngCom /* i/o: FD_CNG structure containing all buffers and variables */
3289 : )
3290 : {
3291 : Word16 i;
3292 : Word32 buf[M + 1 + 320], tmp;
3293 :
3294 : /* Perform IFFT */
3295 4095 : RFFTN_fx( fftBuffer, hFdCngCom->fftSineTab_fx, hFdCngCom->fftlen, 1 );
3296 :
3297 : /* Handle overlap in P/S domain for stereo */
3298 4095 : Copy32( olapBuffer + hFdCngCom->frameSize, olapBuffer, hFdCngCom->frameSize );
3299 4095 : set32_fx( olapBuffer + hFdCngCom->frameSize, 0, hFdCngCom->frameSize ); /*olapBuffer, fftBuffer, olapWin*/
3300 :
3301 611711 : FOR( i = ( hFdCngCom->frameSize >> 2 ); i < ( 3 * ( hFdCngCom->frameSize >> 2 ) ); i++ )
3302 : {
3303 607616 : olapBuffer[i] = L_add( olapBuffer[i], Mpy_32_16_1( fftBuffer[i], olapWin[( i - ( hFdCngCom->frameSize >> 2 ) )] ) );
3304 607616 : move32();
3305 : }
3306 611711 : FOR( ; i < ( 5 * ( hFdCngCom->frameSize >> 2 ) ); i++ )
3307 : {
3308 607616 : olapBuffer[i] = fftBuffer[i];
3309 607616 : move32();
3310 : }
3311 :
3312 611711 : FOR( ; i < ( 7 * ( hFdCngCom->frameSize >> 2 ) ); i++ )
3313 : {
3314 607616 : olapBuffer[i] = fftBuffer[i];
3315 607616 : move32();
3316 : }
3317 :
3318 307903 : FOR( ; i < hFdCngCom->fftlen; i++ )
3319 : {
3320 303808 : olapBuffer[i] = 0;
3321 303808 : move32();
3322 : }
3323 :
3324 4095 : Word32 fftScale = 0;
3325 4095 : move32();
3326 4095 : SWITCH( hFdCngCom->fftlen )
3327 : {
3328 2608 : case 640:
3329 2608 : fftScale = FFT_SCALING_640;
3330 2608 : move32();
3331 2608 : break;
3332 1487 : case 512:
3333 1487 : fftScale = FFT_SCALING_512;
3334 1487 : move32();
3335 1487 : break;
3336 0 : default:
3337 0 : assert( !"Not supported FFT length!" );
3338 : }
3339 :
3340 : /* Get time-domain signal */
3341 4095 : v_multc_fixed( olapBuffer + shr( hFdCngCom->frameSize, 2 ), fftScale, timeDomainOutput, samples_out ); // Q_in - 9
3342 :
3343 : /* Get excitation */
3344 4095 : v_multc_fixed( olapBuffer + sub( shr( hFdCngCom->frameSize, 2 ), ( M + 1 ) ), fftScale, buf, add( M + 1, hFdCngCom->frameSize ) );
3345 4095 : tmp = buf[0];
3346 4095 : move32();
3347 4095 : preemph_ivas_fx( buf + 1, PREEMPH_FAC, M + hFdCngCom->frameSize, &tmp );
3348 : // residu_ivas_fx( hFdCngCom->A_cng, Q13, M, buf + 1 + M, hFdCngCom->exc_cng_32fx, hFdCngCom->frameSize );
3349 4095 : residu_ivas_fx( hFdCngCom->A_cng, sub( 15, norm_s( hFdCngCom->A_cng[0] - 1 ) ), M, buf + 1 + M, hFdCngCom->exc_cng_32fx, hFdCngCom->frameSize );
3350 :
3351 : /* update and window olapBuf if we have a output frame that is shorter than the default frame size...*/
3352 4095 : IF( LT_16( samples_out, hFdCngCom->frameSize ) )
3353 : {
3354 0 : Copy32( olapBuffer + samples_out, olapBuffer + hFdCngCom->frameSize, i_mult( 3, shr( hFdCngCom->frameSize, 2 ) ) );
3355 : }
3356 611711 : FOR( i = ( 5 * ( hFdCngCom->frameSize >> 2 ) ); i < ( 7 * ( hFdCngCom->frameSize >> 2 ) ); i++ )
3357 : {
3358 607616 : olapBuffer[i] = Mpy_32_16_1( olapBuffer[i], olapWin[( i - ( 3 * ( hFdCngCom->frameSize >> 2 ) ) )] );
3359 607616 : move32();
3360 : }
3361 :
3362 4095 : return;
3363 : }
3364 :
3365 46741998 : Word32 rand_gauss_fx(
3366 : Word32 *x,
3367 : Word16 *seed,
3368 : Word16 q )
3369 : {
3370 : Word32 temp;
3371 :
3372 46741998 : temp = own_random( seed );
3373 46741998 : temp = L_add( temp, own_random( seed ) );
3374 46741998 : temp = L_add( temp, own_random( seed ) );
3375 46741998 : temp = L_shr( temp, sub( 15, q ) );
3376 :
3377 46741998 : *x = temp;
3378 46741998 : move32();
3379 :
3380 46741998 : return temp;
3381 : }
3382 :
3383 : /*-------------------------------------------------------------------
3384 : * rand_gauss_fix()
3385 : *
3386 : * Random generator with Gaussian distribution with mean 0 and std 1
3387 : *-------------------------------------------------------------------*/
3388 :
3389 50111712 : Word16 rand_gauss_fix(
3390 : Word16 *x,
3391 : Word16 *seed )
3392 : {
3393 : Word32 temp;
3394 :
3395 50111712 : temp = shr( own_random( seed ), Q2 );
3396 50111712 : temp = L_add( temp, shr( own_random( seed ), Q2 ) );
3397 50111712 : temp = L_add( temp, shr( own_random( seed ), Q2 ) );
3398 :
3399 50111712 : *x = (Word16) temp;
3400 50111712 : move32();
3401 :
3402 50111712 : return (Word16) temp;
3403 : }
|
