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