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