Line data Source code
1 : /******************************************************************************************************
2 :
3 : (C) 2022-2025 IVAS codec Public Collaboration with portions copyright Dolby International AB, Ericsson AB,
4 : Fraunhofer-Gesellschaft zur Foerderung der angewandten Forschung e.V., Huawei Technologies Co. LTD.,
5 : Koninklijke Philips N.V., Nippon Telegraph and Telephone Corporation, Nokia Technologies Oy, Orange,
6 : Panasonic Holdings Corporation, Qualcomm Technologies, Inc., VoiceAge Corporation, and other
7 : contributors to this repository. All Rights Reserved.
8 :
9 : This software is protected by copyright law and by international treaties.
10 : The IVAS codec Public Collaboration consisting of Dolby International AB, Ericsson AB,
11 : Fraunhofer-Gesellschaft zur Foerderung der angewandten Forschung e.V., Huawei Technologies Co. LTD.,
12 : Koninklijke Philips N.V., Nippon Telegraph and Telephone Corporation, Nokia Technologies Oy, Orange,
13 : Panasonic Holdings Corporation, Qualcomm Technologies, Inc., VoiceAge Corporation, and other
14 : contributors to this repository retain full ownership rights in their respective contributions in
15 : the software. This notice grants no license of any kind, including but not limited to patent
16 : license, nor is any license granted by implication, estoppel or otherwise.
17 :
18 : Contributors are required to enter into the IVAS codec Public Collaboration agreement before making
19 : contributions.
20 :
21 : This software is provided "AS IS", without any express or implied warranties. The software is in the
22 : development stage. It is intended exclusively for experts who have experience with such software and
23 : solely for the purpose of inspection. All implied warranties of non-infringement, merchantability
24 : and fitness for a particular purpose are hereby disclaimed and excluded.
25 :
26 : Any dispute, controversy or claim arising under or in relation to providing this software shall be
27 : submitted to and settled by the final, binding jurisdiction of the courts of Munich, Germany in
28 : accordance with the laws of the Federal Republic of Germany excluding its conflict of law rules and
29 : the United Nations Convention on Contracts on the International Sales of Goods.
30 :
31 : *******************************************************************************************************/
32 :
33 : /*====================================================================================
34 : EVS Codec 3GPP TS26.452 Aug 12, 2021. Version 16.3.0
35 : ====================================================================================*/
36 :
37 : #include <assert.h>
38 : #include <stdint.h>
39 : #include "options.h"
40 : #include "complex_basop.h"
41 : #include "basop_util.h"
42 : #include "basop32.h"
43 : #include "rom_com.h"
44 : #include "rom_basic_math.h"
45 : #include "basop_settings.h"
46 : #include "cnst.h"
47 :
48 : extern const Word32 SqrtTable[32]; // Q31
49 : extern const Word16 SqrtDiffTable[32]; /* Q15 */
50 :
51 : extern const Word32 ISqrtTable[32];
52 : extern const Word16 ISqrtDiffTable[32];
53 :
54 : extern const Word32 InvTable[32];
55 : extern const Word16 InvDiffTable[32];
56 :
57 854420549 : Word32 BASOP_Util_Log2( Word32 x )
58 : {
59 : Word32 exp;
60 : Word16 exp_e;
61 : Word16 nIn;
62 : Word16 accuSqr;
63 : Word32 accuRes;
64 :
65 :
66 854420549 : assert( x >= 0 );
67 :
68 854420549 : if ( x == 0 )
69 : {
70 :
71 72200 : return ( (Word32) MIN_32 );
72 : }
73 :
74 : /* normalize input, calculate integer part */
75 854348349 : exp_e = norm_l( x );
76 854348349 : x = L_shl( x, exp_e );
77 854348349 : exp = L_deposit_l( exp_e );
78 :
79 : /* calculate (1-normalized_input) */
80 854348349 : nIn = extract_h( L_sub( MAX_32, x ) );
81 :
82 : /* approximate ln() for fractional part (nIn *c0 + nIn^2*c1 + nIn^3*c2 + ... + nIn^8 *c7) */
83 :
84 : /* iteration 1, no need for accumulation */
85 854348349 : accuRes = L_mult( nIn, ldCoeff[0] ); /* nIn^i * coeff[0] */
86 854348349 : accuSqr = mult( nIn, nIn ); /* nIn^2, nIn^3 .... */
87 :
88 : /* iteration 2 */
89 854348349 : accuRes = L_mac( accuRes, accuSqr, ldCoeff[1] ); /* nIn^i * coeff[1] */
90 854348349 : accuSqr = mult( accuSqr, nIn ); /* nIn^2, nIn^3 .... */
91 :
92 : /* iteration 3 */
93 854348349 : accuRes = L_mac( accuRes, accuSqr, ldCoeff[2] ); /* nIn^i * coeff[2] */
94 854348349 : accuSqr = mult( accuSqr, nIn ); /* nIn^2, nIn^3 .... */
95 :
96 : /* iteration 4 */
97 854348349 : accuRes = L_mac( accuRes, accuSqr, ldCoeff[3] ); /* nIn^i * coeff[3] */
98 854348349 : accuSqr = mult( accuSqr, nIn ); /* nIn^2, nIn^3 .... */
99 :
100 : /* iteration 5 */
101 854348349 : accuRes = L_mac( accuRes, accuSqr, ldCoeff[4] ); /* nIn^i * coeff[4] */
102 854348349 : accuSqr = mult( accuSqr, nIn ); /* nIn^2, nIn^3 .... */
103 :
104 : /* iteration 6 */
105 854348349 : accuRes = L_mac( accuRes, accuSqr, ldCoeff[5] ); /* nIn^i * coeff[5] */
106 854348349 : accuSqr = mult( accuSqr, nIn ); /* nIn^2, nIn^3 .... */
107 :
108 : /* iteration 7, no need to calculate accuSqr any more */
109 854348349 : accuRes = L_mac( accuRes, accuSqr, ldCoeff[6] ); /* nIn^i * coeff[6] */
110 :
111 : /* ld(fractional part) = ln(fractional part)/ln(2), 1/ln(2) = (1 + 0.44269504) */
112 854348349 : accuRes = L_mac0( L_shr( accuRes, 1 ), extract_h( accuRes ), 14506 );
113 :
114 854348349 : accuRes = L_shr( accuRes, LD_DATA_SCALE - 1 ); /* fractional part/LD_DATA_SCALE */
115 854348349 : exp = L_shl( exp, ( 31 - LD_DATA_SCALE ) ); /* integer part/LD_DATA_SCALE */
116 854348349 : accuRes = L_sub( accuRes, exp ); /* result = integer part + fractional part */
117 :
118 :
119 854348349 : return ( accuRes );
120 : }
121 :
122 11230835 : Word32 BASOP_Util_InvLog2( Word32 x )
123 : {
124 : Word16 frac;
125 : Word16 exp;
126 : Word32 retVal;
127 : UWord32 index3;
128 : UWord32 index2;
129 : UWord32 index1;
130 : UWord32 lookup3f;
131 : UWord32 lookup12;
132 : UWord32 lookup;
133 :
134 :
135 11230835 : if ( x < -1040187392l /*-31.0/64.0 Q31*/ )
136 : {
137 :
138 71 : return 0;
139 : }
140 11230764 : test();
141 11230764 : if ( ( GE_32( x, 1040187392l /*31.0/64.0 Q31*/ ) ) || ( x == 0 ) )
142 : {
143 :
144 0 : return 0x7FFFFFFF;
145 : }
146 :
147 11230764 : frac = extract_l( L_and( x, 0x3FF ) );
148 :
149 11230764 : index3 = L_and( L_shr( x, 10 ), 0x1F );
150 11230764 : index2 = L_and( L_shr( x, 15 ), 0x1F );
151 11230764 : index1 = L_and( L_shr( x, 20 ), 0x1F );
152 :
153 11230764 : exp = extract_l( L_shr( x, 25 ) );
154 11230764 : if ( x > 0 )
155 : {
156 112 : exp = sub( 31, exp );
157 : }
158 11230764 : if ( x < 0 )
159 : {
160 11230652 : exp = negate( exp );
161 : }
162 :
163 11230764 : lookup3f = L_add( exp2x_tab_long[index3], L_shr( Mpy_32_16_1( 0x0016302F, frac ), 1 ) );
164 11230764 : lookup12 = Mpy_32_32( exp2_tab_long[index1], exp2w_tab_long[index2] );
165 11230764 : lookup = Mpy_32_32( lookup12, lookup3f );
166 :
167 11230764 : retVal = L_shr( lookup, sub( exp, 3 ) );
168 :
169 :
170 11230764 : return retVal;
171 : }
172 :
173 4608922 : Word32 BASOP_Util_Log10( Word32 x, Word16 e )
174 : {
175 4608922 : test();
176 4608922 : IF( e >= 0 && LE_16( e, 31 ) )
177 : {
178 4200834 : IF( EQ_32( x, L_shl_sat( 1, sub( 31, e ) ) ) )
179 : {
180 119438 : return 0;
181 : }
182 : }
183 4489484 : Word32 res = BASOP_Util_Log2( x );
184 4489484 : res = L_add( Mpy_32_32( res, 646456993 /* log10(2) in Q31 */ ), Mpy_32_32( L_shl( e, 24 ), 1292913986 /* log10(2) in Q32 */ ) ); // Adjusting for the exponent mismatch: multiplying first so as to avoid saturation
185 : /* log10(2) is used in Q32 to support exponent till 127 in Mpy_32_32( L_shl( e, 24 ), 1292913986 )*/
186 4489484 : return res;
187 : }
188 :
189 29541431 : Word32 BASOP_Util_Loge( Word32 x, Word16 e )
190 : {
191 29541431 : Word32 res = BASOP_Util_Log2( x );
192 29541431 : res = L_add( Mpy_32_32( res, 1488522235 /* loge(2) in Q31 */ ), Mpy_32_32( L_shl( e, 25 ), 1488522235 /* loge(2) in Q31 */ ) ); // Adjusting for the exponent mismatch: multiplying first so as to avoid saturation
193 29541431 : return res;
194 : }
195 :
196 64209577 : Word16 BASOP_Util_Add_MantExp /*!< Exponent of result */
197 : ( Word16 a_m, /*!< Mantissa of 1st operand a */
198 : Word16 a_e, /*!< Exponent of 1st operand a */
199 : Word16 b_m, /*!< Mantissa of 2nd operand b */
200 : Word16 b_e, /*!< Exponent of 2nd operand b */
201 : Word16 *ptrSum_m ) /*!< Mantissa of result */
202 : {
203 : Word32 L_lm, L_hm;
204 : Word16 shift;
205 :
206 :
207 : /* Compare exponents: the difference is limited to +/- 15
208 : The Word16 mantissa of the operand with higher exponent is moved into the low
209 : part of a Word32 and shifted left by the exponent difference. Then, the
210 : unshifted mantissa of the operand with the lower exponent is added to the lower
211 : 16 bits. The addition result is normalized and the upper Word16 of the result represents
212 : the mantissa to return. The returned exponent takes into account all shift operations
213 : including the final 16-bit extraction.
214 : Note: The resulting mantissa may be inaccurate in the case, where the mantissa of the operand
215 : with higher exponent is not really left-aligned, while the mantissa of the operand with
216 : lower exponent is so. If in such a case, the difference in exponents is more than 15,
217 : an inaccuracy is introduced.
218 : Example:
219 : A: a_e = 20, a_m = 0x0001
220 : B: b_e = 0, b_m = 0x4000
221 : correct: A+B=1*2^20+1*2^14=0x0010.0000+0x0000.4000=0x0010.4000=0x4100*2^6
222 : previously: A+B=1*2^20+1*2^14=0x0001+0x0000=0x0001*2^20
223 : this version: A+B=1*2^20+1*2^14=0x0000.8000+0x0000.4000=0x6000*2^6
224 : */
225 :
226 64209577 : shift = sub( a_e, b_e );
227 64209577 : if ( shift >= 0 )
228 43843243 : shift = s_min( 15, shift );
229 :
230 64209577 : if ( shift < 0 )
231 20366334 : shift = s_max( -15, shift );
232 64209577 : a_e = s_max( a_e, b_e );
233 64209577 : L_hm = L_deposit_l( a_m ); /* mantissa belonging to higher exponent */
234 64209577 : L_lm = L_deposit_l( a_m ); /* mantissa belonging to lower exponent */
235 64209577 : if ( shift >= 0 )
236 43843243 : L_lm = L_deposit_l( b_m );
237 64209577 : if ( shift < 0 )
238 20366334 : L_hm = L_deposit_l( b_m );
239 :
240 64209577 : if ( shift > 0 )
241 33801645 : shift = negate( shift );
242 :
243 64209577 : L_hm = L_shr( L_hm, shift ); /* shift left due to negative shift parameter */
244 64209577 : a_e = add( a_e, shift );
245 64209577 : L_hm = L_add( L_hm, L_lm );
246 64209577 : shift = norm_l( L_hm );
247 64209577 : L_hm = L_shl( L_hm, shift );
248 64209577 : *ptrSum_m = extract_h( L_hm );
249 64209577 : move16();
250 :
251 64209577 : a_e = sub( a_e, shift );
252 64209577 : if ( L_hm )
253 63467413 : a_e = add( a_e, 16 );
254 64209577 : return ( a_e );
255 : }
256 :
257 :
258 2126845 : void BASOP_Util_Divide_MantExp( Word16 a_m, /*!< Mantissa of dividend a */
259 : Word16 a_e, /*!< Exponent of dividend a */
260 : Word16 b_m, /*!< Mantissa of divisor b */
261 : Word16 b_e, /*!< Exponent of divisor b */
262 : Word16 *ptrResult_m, /*!< Mantissa of quotient a/b */
263 : Word16 *ptrResult_e /*!< Exponent of quotient a/b */
264 : )
265 : {
266 : Word16 index, frac;
267 : Word16 preShift, postShift;
268 : Word16 m;
269 : Word32 m32;
270 : #ifndef ISSUE_1836_replace_overflow_libcom
271 : #ifdef BASOP_NOGLOB_DECLARE_LOCAL
272 : Flag Overflow = 0;
273 : #endif
274 : #endif
275 :
276 :
277 2126845 : assert( b_m != 0 );
278 :
279 : /* normalize b */
280 2126845 : preShift = norm_s( b_m );
281 2126845 : m = shl( b_m, preShift );
282 :
283 : /* make b positive */
284 : BASOP_SATURATE_WARNING_OFF_EVS;
285 2126845 : if ( m < 0 )
286 0 : m = negate( m );
287 : BASOP_SATURATE_WARNING_ON_EVS;
288 :
289 : /* get table index (upper 6 bits minus 16) */
290 : /* index = (m >> 9) - 32; */
291 2126845 : index = mac_r( -32768 - ( 32 << 16 ), m, 1 << 6 );
292 :
293 : /* get fractional part for interpolation (lower 9 bits) */
294 2126845 : frac = shl( s_and( m, 0x1FF ), 1 ); /* Q10 */
295 :
296 : /* interpolate 1/b */
297 2126845 : m = msu_r( InvTable[index], InvDiffTable[index], frac );
298 :
299 : /* restore sign */
300 2126845 : if ( b_m < 0 )
301 0 : m = negate( m );
302 :
303 : /* multiply with a */
304 2126845 : m32 = L_mult( a_m, m );
305 :
306 : /* normalize result */
307 2126845 : postShift = norm_l( m32 );
308 : #ifdef ISSUE_1836_replace_overflow_libcom
309 2126845 : m = round_fx_sat( L_shl( m32, postShift ) );
310 : #else
311 : m = round_fx_o( L_shl( m32, postShift ), &Overflow );
312 : #endif
313 :
314 : /* exponent */
315 2126845 : *ptrResult_e = sub( add( add( a_e, sub( 1, b_e ) ), preShift ), postShift );
316 2126845 : move16();
317 :
318 2126845 : *ptrResult_m = m;
319 2126845 : move16();
320 2126845 : }
321 :
322 :
323 : /* local function for Sqrt16 */
324 97421326 : static Word16 Sqrt16_common( Word16 m,
325 : Word16 e )
326 : {
327 : Word16 index, frac;
328 : #ifndef ISSUE_1836_replace_overflow_libcom
329 : Flag Overflow;
330 : #endif
331 :
332 :
333 97421326 : assert( ( m >= 0x4000 ) || ( m == 0 ) );
334 :
335 : /* get table index (upper 6 bits minus 32) */
336 : /* index = (m >> 9) - 32; */
337 97421326 : index = mac_r( -32768 - ( 32 << 16 ), m, 1 << 6 );
338 :
339 : /* get fractional part for interpolation (lower 9 bits) */
340 97421326 : frac = s_and( m, 0x1FF ); /* Q9 */
341 :
342 : /* interpolate */
343 97421326 : if ( m != 0 )
344 : {
345 : #ifdef ISSUE_1836_replace_overflow_libcom
346 95590131 : m = mac_r_sat( SqrtTable[index], SqrtDiffTable[index], frac );
347 : #else
348 : BASOP_SATURATE_WARNING_OFF_EVS;
349 : m = mac_ro( SqrtTable[index], SqrtDiffTable[index], frac, &Overflow );
350 : BASOP_SATURATE_WARNING_ON_EVS;
351 : #endif
352 : }
353 :
354 : /* handle odd exponents */
355 97421326 : if ( s_and( e, 1 ) != 0 )
356 42331616 : m = mult_r( m, 0x5a82 );
357 :
358 97421326 : return m;
359 : }
360 :
361 :
362 : /* local function for Sqrt32 and Sqrt32norm */
363 945157343 : static Word32 Sqrt32_common( Word32 m,
364 : Word16 e )
365 : {
366 : Word16 m16, index, frac;
367 : #ifndef ISSUE_1836_replace_overflow_libcom
368 : #ifdef BASOP_NOGLOB_DECLARE_LOCAL
369 : Flag Overflow = 0;
370 : #endif
371 : #endif
372 :
373 945157343 : assert( ( m >= 0x40000000 ) || ( m == 0 ) );
374 :
375 : #ifdef ISSUE_1836_replace_overflow_libcom
376 945157343 : m16 = round_fx_sat( m );
377 : #else
378 : m16 = round_fx_o( m, &Overflow );
379 : #endif
380 :
381 : /* get table index (upper 6 bits minus 32) */
382 : /* index = (m16 >> 9) - 32; */
383 945157343 : index = mac_r( -32768 - ( 32 << 16 ), m16, 1 << 6 );
384 :
385 : /* get fractional part for interpolation (lower 9 bits) */
386 945157343 : frac = s_and( m16, 0x1FF ); /* Q9 */
387 :
388 : /* interpolate */
389 945157343 : if ( m != 0 )
390 : {
391 : BASOP_SATURATE_WARNING_OFF_EVS;
392 895919453 : m = L_mac_sat( SqrtTable[index], SqrtDiffTable[index], frac );
393 : BASOP_SATURATE_WARNING_ON_EVS;
394 : }
395 :
396 : /* handle odd exponents */
397 945157343 : if ( s_and( e, 1 ) != 0 )
398 446272782 : m = Mpy_32_16_1( m, 0x5a82 );
399 :
400 945157343 : return m;
401 : }
402 :
403 : /* local function for ISqrt16 and ISqrt16norm */
404 51294580 : static Word16 ISqrt16_common( Word16 m,
405 : Word16 e )
406 : {
407 : Word16 index, frac;
408 :
409 51294580 : assert( m >= 0x4000 );
410 :
411 : /* get table index (upper 6 bits minus 32) */
412 : /* index = (m >> 9) - 32; */
413 51294580 : index = mac_r( -32768 - ( 32 << 16 ), m, 1 << 6 );
414 :
415 : /* get fractional part for interpolation (lower 9 bits) */
416 51294580 : frac = s_and( m, 0x1FF ); /* Q9 */
417 :
418 : /* interpolate */
419 51294580 : m = msu_r( ISqrtTable[index], ISqrtDiffTable[index], frac );
420 :
421 : /* handle even exponents */
422 51294580 : if ( s_and( e, 1 ) == 0 )
423 : {
424 22421830 : m = mult_r( m, 0x5a82 );
425 : }
426 :
427 51294580 : return m;
428 : }
429 :
430 : /* local function for ISqrt32 and ISqrt32norm */
431 99420423 : static Word32 ISqrt32_common( Word32 m,
432 : Word16 e )
433 : {
434 : Word16 m16, index, frac;
435 : #ifndef ISSUE_1836_replace_overflow_libcom
436 : #ifdef BASOP_NOGLOB_DECLARE_LOCAL
437 : Flag Overflow = 0;
438 : #endif
439 : #endif
440 :
441 99420423 : assert( m >= 0x40000000 );
442 : #ifdef ISSUE_1836_replace_overflow_libcom
443 99420423 : m16 = round_fx_sat( m );
444 : #else
445 : #ifdef BASOP_NOGLOB_DECLARE_LOCAL
446 : m16 = round_fx_o( m, &Overflow );
447 : #else
448 : m16 = round_fx( m );
449 : #endif
450 : #endif
451 :
452 :
453 : /* get table index (upper 6 bits minus 32) */
454 : /* index = (m16 >> 25) - 32; */
455 99420423 : index = mac_r( -32768 - ( 32 << 16 ), m16, 1 << 6 );
456 :
457 : /* get fractional part for interpolation (lower 9 bits) */
458 99420423 : frac = s_and( m16, 0x1FF ); /* Q9 */
459 :
460 : /* interpolate */
461 99420423 : m = L_msu( ISqrtTable[index], ISqrtDiffTable[index], frac );
462 :
463 : /* handle even exponents */
464 99420423 : if ( s_and( e, 1 ) == 0 )
465 : {
466 53741875 : m = Mpy_32_16_1( m, 0x5a82 );
467 : }
468 :
469 99420423 : return m;
470 : }
471 :
472 96843690 : Word16 Sqrt16( /*!< output mantissa */
473 : Word16 mantissa, /*!< input mantissa */
474 : Word16 *exponent /*!< pointer to exponent */
475 : )
476 : {
477 : Word16 preShift, e;
478 :
479 96843690 : assert( mantissa >= 0 );
480 :
481 : /* normalize */
482 96843690 : preShift = norm_s( mantissa );
483 :
484 96843690 : e = sub( *exponent, preShift );
485 96843690 : mantissa = shl( mantissa, preShift );
486 :
487 : /* calc mantissa */
488 96843690 : mantissa = Sqrt16_common( mantissa, e );
489 :
490 : /* e = (e + 1) >> 1 */
491 96843690 : *exponent = mult_r( e, 1 << 14 );
492 96843690 : move16();
493 :
494 96843690 : return mantissa;
495 : }
496 :
497 :
498 577636 : Word16 Sqrt16norm( /*!< output mantissa */
499 : Word16 mantissa, /*!< normalized input mantissa */
500 : Word16 *exponent /*!< pointer to exponent */
501 : )
502 : {
503 :
504 577636 : assert( ( mantissa >= 0x4000 ) || ( mantissa == 0 ) );
505 :
506 : /* calc mantissa */
507 577636 : mantissa = Sqrt16_common( mantissa, *exponent );
508 :
509 : /* e = (e + 1) >> 1 */
510 577636 : *exponent = mult_r( *exponent, 1 << 14 );
511 577636 : move16();
512 :
513 577636 : return mantissa;
514 : }
515 :
516 :
517 51294580 : Word16 ISqrt16( /*!< output mantissa */
518 : Word16 mantissa, /*!< input mantissa */
519 : Word16 *exponent /*!< pointer to exponent */
520 : )
521 : {
522 : Word16 preShift, e;
523 :
524 51294580 : assert( mantissa > 0 );
525 :
526 : /* normalize */
527 51294580 : preShift = norm_s( mantissa );
528 :
529 51294580 : e = sub( *exponent, preShift );
530 51294580 : mantissa = shl( mantissa, preShift );
531 :
532 : /* calc mantissa */
533 51294580 : mantissa = ISqrt16_common( mantissa, e );
534 :
535 : /* e = (2 - e) >> 1 */
536 51294580 : *exponent = msu_r( 1L << 15, e, 1 << 14 );
537 51294580 : move16();
538 :
539 51294580 : return mantissa;
540 : }
541 :
542 945157343 : Word32 Sqrt32( /*!< output mantissa */
543 : Word32 mantissa, /*!< input mantissa */
544 : Word16 *exponent /*!< pointer to exponent */
545 : )
546 : {
547 : Word16 preShift, e;
548 :
549 945157343 : assert( mantissa >= 0 );
550 :
551 : /* normalize */
552 945157343 : preShift = norm_l( mantissa );
553 :
554 945157343 : e = sub( *exponent, preShift );
555 945157343 : mantissa = L_shl( mantissa, preShift );
556 :
557 : /* calc mantissa */
558 945157343 : mantissa = Sqrt32_common( mantissa, e );
559 :
560 : /* e = (e + 1) >> 1 */
561 945157343 : *exponent = mult_r( e, 1 << 14 );
562 945157343 : move16();
563 :
564 945157343 : return mantissa;
565 : }
566 :
567 :
568 99387343 : Word32 ISqrt32( /*!< output mantissa */
569 : Word32 mantissa, /*!< input mantissa */
570 : Word16 *exponent /*!< pointer to exponent */
571 : )
572 : {
573 : Word16 preShift, e;
574 :
575 99387343 : assert( mantissa > 0 );
576 :
577 : /* normalize */
578 99387343 : preShift = norm_l( mantissa );
579 :
580 99387343 : e = sub( *exponent, preShift );
581 99387343 : mantissa = L_shl( mantissa, preShift );
582 :
583 : /* calc mantissa */
584 99387343 : mantissa = ISqrt32_common( mantissa, e );
585 :
586 : /* e = (2 - e) >> 1 */
587 99387343 : *exponent = msu_r( 1L << 15, e, 1 << 14 );
588 99387343 : move16();
589 :
590 99387343 : return mantissa;
591 : }
592 :
593 33080 : Word32 ISqrt32norm( /*!< output mantissa */
594 : Word32 mantissa, /*!< normalized input mantissa */
595 : Word16 *exponent /*!< pointer to exponent */
596 : )
597 : {
598 :
599 33080 : assert( mantissa >= 0x40000000 );
600 :
601 : /* calc mantissa */
602 33080 : mantissa = ISqrt32_common( mantissa, *exponent );
603 :
604 : /* e = (2 - e) >> 1 */
605 33080 : *exponent = msu_r( 1L << 15, *exponent, 1 << 14 );
606 33080 : move16();
607 :
608 33080 : return mantissa;
609 : }
610 :
611 29461420 : Word16 Inv16( /*!< output mantissa */
612 : Word16 mantissa, /*!< input mantissa */
613 : Word16 *exponent /*!< pointer to exponent */
614 : )
615 : {
616 : Word16 index, frac;
617 : Word16 preShift;
618 : Word16 m, e;
619 :
620 :
621 29461420 : assert( mantissa != 0 );
622 :
623 : /* absolute */
624 : BASOP_SATURATE_WARNING_OFF_EVS;
625 29461420 : m = abs_s( mantissa );
626 : BASOP_SATURATE_WARNING_ON_EVS;
627 :
628 : /* normalize */
629 29461420 : preShift = norm_s( m );
630 :
631 29461420 : e = sub( *exponent, preShift );
632 29461420 : m = shl( m, preShift );
633 :
634 : /* get table index (upper 6 bits minus 32) */
635 : /* index = (m >> 9) - 32; */
636 29461420 : index = mac_r( -32768 - ( 32 << 16 ), m, 1 << 6 );
637 :
638 : /* get fractional part for interpolation (lower 9 bits) */
639 29461420 : frac = shl( s_and( m, 0x1FF ), 1 ); /* Q10 */
640 :
641 : /* interpolate */
642 29461420 : m = msu_r( InvTable[index], InvDiffTable[index], frac );
643 :
644 : /* restore sign */
645 29461420 : if ( mantissa < 0 )
646 2511801 : m = negate( m );
647 :
648 : /* e = 1 - e */
649 29461420 : *exponent = sub( 1, e );
650 29461420 : move16();
651 :
652 29461420 : return m;
653 : }
654 :
655 :
656 56336247 : void BASOP_Util_Sqrt_InvSqrt_MantExp( Word16 mantissa, /*!< mantissa */
657 : Word16 exponent, /*!< expoinent */
658 : Word16 *sqrt_mant, /*!< Pointer to sqrt mantissa */
659 : Word16 *sqrt_exp, /*!< Pointer to sqrt exponent */
660 : Word16 *isqrt_mant, /*!< Pointer to 1/sqrt mantissa */
661 : Word16 *isqrt_exp /*!< Pointer to 1/sqrt exponent */
662 : )
663 : {
664 : Word16 index, frac;
665 : Word16 preShift;
666 : Word16 m, mi, e_odd;
667 :
668 :
669 56336247 : assert( mantissa > 0 );
670 :
671 : /* normalize */
672 56336247 : preShift = norm_s( mantissa );
673 :
674 56336247 : exponent = sub( exponent, preShift );
675 56336247 : mantissa = shl( mantissa, preShift );
676 :
677 : /* get table index (upper 6 bits minus 32) */
678 : /* index = (m >> 9) - 32; */
679 56336247 : index = mac_r( -32768 - ( 32 << 16 ), mantissa, 1 << 6 );
680 :
681 : /* get fractional part for interpolation (lower 9 bits) */
682 56336247 : frac = s_and( mantissa, 0x1FF ); /* Q9 */
683 :
684 : /* interpolate */
685 : BASOP_SATURATE_WARNING_OFF_EVS;
686 56336247 : m = mac_r_sat( SqrtTable[index], SqrtDiffTable[index], frac );
687 56336247 : mi = msu_r_sat( ISqrtTable[index], ISqrtDiffTable[index], frac );
688 : BASOP_SATURATE_WARNING_ON_EVS;
689 :
690 : /* handle even/odd exponents */
691 56336247 : e_odd = s_and( exponent, 1 );
692 56336247 : if ( e_odd != 0 )
693 28908600 : m = mult_r( m, 0x5a82 );
694 56336247 : if ( e_odd == 0 )
695 27427647 : mi = mult_r( mi, 0x5a82 );
696 :
697 : /* e = (e + 1) >> 1 */
698 56336247 : *sqrt_exp = mult_r( exponent, 1 << 14 );
699 56336247 : move16();
700 :
701 : /* e = (2 - e) >> 1 */
702 56336247 : *isqrt_exp = msu_r( 1L << 15, exponent, 1 << 14 );
703 56336247 : move16();
704 :
705 : /* Write result */
706 56336247 : *sqrt_mant = m;
707 56336247 : move16();
708 56336247 : *isqrt_mant = mi;
709 56336247 : move16();
710 56336247 : }
711 :
712 :
713 : /********************************************************************/
714 : /*!
715 : \brief Calculates the scalefactor needed to normalize input array
716 :
717 : The scalefactor needed to normalize the Word16 input array is returned <br>
718 : If the input array contains only '0', a scalefactor 0 is returned <br>
719 : Scaling factor is determined wrt a normalized target x: 16384 <= x <= 32767 for positive x <br>
720 : and -32768 <= x <= -16384 for negative x
721 : */
722 :
723 13448387 : Word16 getScaleFactor16( /* o: measured headroom in range [0..15], 0 if all x[i] == 0 */
724 : const Word16 *x, /* i: array containing 16-bit data */
725 : const Word16 len_x ) /* i: length of the array to scan */
726 : {
727 : Word16 i, i_min, i_max;
728 : Word16 x_min, x_max;
729 :
730 :
731 13448387 : x_max = 0;
732 13448387 : move16();
733 13448387 : x_min = 0;
734 13448387 : move16();
735 10085236969 : FOR( i = 0; i < len_x; i++ )
736 : {
737 10071788582 : if ( x[i] >= 0 )
738 6590216176 : x_max = s_max( x_max, x[i] );
739 10071788582 : if ( x[i] < 0 )
740 3481572406 : x_min = s_min( x_min, x[i] );
741 : }
742 :
743 13448387 : i_max = 0x10;
744 13448387 : move16();
745 13448387 : i_min = 0x10;
746 13448387 : move16();
747 :
748 13448387 : if ( x_max != 0 )
749 11110182 : i_max = norm_s( x_max );
750 :
751 13448387 : if ( x_min != 0 )
752 10867214 : i_min = norm_s( x_min );
753 :
754 13448387 : i = s_and( s_min( i_max, i_min ), 0xF );
755 :
756 :
757 13448387 : return i;
758 : }
759 :
760 :
761 : /********************************************************************/
762 : /*!
763 : \brief Calculates the scalefactor needed to normalize input array
764 :
765 : The scalefactor needed to normalize the Word32 input array is returned <br>
766 : If the input array contains only '0', a scalefactor 0 is returned <br>
767 : Scaling factor is determined wrt a normalized target x: 1073741824 <= x <= 2147483647 for positive x <br>
768 : and -2147483648 <= x <= -1073741824 for negative x
769 : */
770 :
771 123944297 : Word16 getScaleFactor32( /* o: measured headroom in range [0..31], 0 if all x[i] == 0 */
772 : const Word32 *x, /* i: array containing 32-bit data */
773 : const Word16 len_x ) /* i: length of the array to scan */
774 : {
775 : Word16 i, i_min, i_max;
776 : Word32 x_min, x_max;
777 :
778 :
779 123944297 : x_max = 0;
780 123944297 : move32();
781 123944297 : x_min = 0;
782 123944297 : move32();
783 15197422242 : FOR( i = 0; i < len_x; i++ )
784 : {
785 15073477945 : if ( x[i] >= 0 )
786 9910614675 : x_max = L_max( x_max, x[i] );
787 15073477945 : if ( x[i] < 0 )
788 5162863270 : x_min = L_min( x_min, x[i] );
789 : }
790 :
791 123944297 : i_max = 0x20;
792 123944297 : move16();
793 123944297 : i_min = 0x20;
794 123944297 : move16();
795 :
796 123944297 : if ( x_max != 0 )
797 114190173 : i_max = norm_l( x_max );
798 :
799 123944297 : if ( x_min != 0 )
800 101240564 : i_min = norm_l( x_min );
801 :
802 123944297 : i = s_and( s_min( i_max, i_min ), 0x1F );
803 :
804 :
805 123944297 : return i;
806 : }
807 :
808 19551 : Word16 getScaleFactor32_copy( /* o: measured headroom in range [0..31], 0 if all x[i] == 0 */
809 : const Word32 *x, /* i: array containing 32-bit data */
810 : const Word32 len_x ) /* i: length of the array to scan */
811 : {
812 : Word32 i;
813 : Word16 i_min, i_max;
814 : Word32 x_min, x_max;
815 :
816 :
817 19551 : x_max = 0;
818 19551 : move32();
819 19551 : x_min = 0;
820 19551 : move32();
821 129305631 : FOR( i = 0; i < len_x; i++ )
822 : {
823 129286080 : if ( x[i] >= 0 )
824 107800501 : x_max = L_max( x_max, x[i] );
825 129286080 : if ( x[i] < 0 )
826 21485579 : x_min = L_min( x_min, x[i] );
827 : }
828 :
829 19551 : i_max = 0x20;
830 19551 : move16();
831 19551 : i_min = 0x20;
832 19551 : move16();
833 :
834 19551 : if ( x_max != 0 )
835 19542 : i_max = norm_l( x_max );
836 :
837 19551 : if ( x_min != 0 )
838 19551 : i_min = norm_l( x_min );
839 :
840 19551 : i_max = s_and( s_min( i_max, i_min ), 0x1F );
841 :
842 :
843 19551 : return i_max;
844 : }
845 :
846 1124378 : Word16 normalize16( Word16 mantissa, Word16 *pexponent )
847 : {
848 : Word16 tmp;
849 :
850 1124378 : tmp = norm_s( mantissa );
851 1124378 : mantissa = shl( mantissa, tmp );
852 1124378 : move16();
853 1124378 : *pexponent = sub( *pexponent, tmp );
854 :
855 :
856 1124378 : return mantissa;
857 : }
858 556152 : Word16 divide3216( Word32 x, Word16 y )
859 : {
860 : Word16 z;
861 :
862 :
863 556152 : z = 0;
864 556152 : move16();
865 556152 : if ( 0 == y )
866 : {
867 0 : return 0x7fff;
868 : }
869 :
870 556152 : IF( x != 0 )
871 : {
872 : Word16 den, sign;
873 : Word32 num;
874 556152 : num = L_abs( x );
875 556152 : den = abs_s( y );
876 :
877 556152 : sign = s_and( s_xor( extract_h( x ), y ), -32768 /* 0x8000 */ );
878 :
879 556152 : z = div_l( num, den );
880 556152 : if ( 0 != sign )
881 : {
882 0 : z = negate( z );
883 : }
884 : }
885 :
886 :
887 556152 : return z;
888 : }
889 :
890 6042211 : Word16 divide1616( Word16 x, Word16 y )
891 : {
892 : Word16 z, num, den, sign;
893 :
894 :
895 6042211 : num = abs_s( x );
896 6042211 : den = abs_s( y );
897 :
898 6042211 : sign = s_and( s_xor( x, y ), -32768 /* 0x8000 */ );
899 :
900 6042211 : move16();
901 6042211 : z = 0x7fff;
902 6042211 : if ( LT_16( num, den ) )
903 6041753 : z = div_s( num, den );
904 :
905 6042211 : if ( 0 != sign )
906 : {
907 581643 : z = negate( z );
908 : }
909 :
910 :
911 6042211 : return z;
912 : }
913 :
914 34263060 : Word16 divide3232( Word32 L_num, Word32 L_denom )
915 : {
916 : Word16 z;
917 : Word32 sign;
918 : #ifndef ISSUE_1836_replace_overflow_libcom
919 : #ifdef BASOP_NOGLOB_DECLARE_LOCAL
920 : Flag Overflow = 0;
921 : #endif
922 : #endif
923 :
924 :
925 34263060 : sign = L_and( L_xor( L_num, L_denom ), (Word32) 0x80000000 );
926 :
927 34263060 : L_num = L_abs( L_num );
928 34263060 : L_denom = L_abs( L_denom );
929 :
930 : /* limit the range of denominator to Word16 */
931 34263060 : z = s_min( norm_l( L_num ), norm_l( L_denom ) );
932 34263060 : L_num = L_shl( L_num, z );
933 34263060 : L_denom = L_shl( L_denom, z );
934 :
935 : /* round_fx instead of extract_h improves spectral distortion in E_UTIL_lev_dur (schur version). */
936 : #ifdef ISSUE_1836_replace_overflow_libcom
937 34263060 : z = div_l( L_num, round_fx_sat( L_denom ) );
938 : #else
939 : z = div_l( L_num, round_fx_o( L_denom, &Overflow ) );
940 : #endif
941 34263060 : if ( 0 != sign )
942 : {
943 14993169 : z = negate( z );
944 : }
945 :
946 :
947 34263060 : return z;
948 : }
949 :
950 77821310 : Word16 BASOP_Util_Divide3232_uu_1616_Scale( Word32 x, Word32 y, Word16 *s )
951 : {
952 : Word16 z;
953 : Word16 sx;
954 : Word16 sy;
955 : Word16 x16;
956 : Word16 y16;
957 :
958 :
959 77821310 : assert( x >= 0 );
960 77821310 : assert( y > 0 );
961 :
962 77821310 : if ( x == 0 )
963 : {
964 6669818 : *s = 0;
965 6669818 : move16();
966 :
967 :
968 6669818 : return ( 0 );
969 : }
970 :
971 71151492 : sx = norm_l( x );
972 71151492 : sy = norm_l( y );
973 :
974 71151492 : x16 = extract_h( L_shl( x, sx ) );
975 71151492 : y16 = extract_h( L_shl( y, sy ) );
976 :
977 71151492 : if ( GT_16( x16, y16 ) )
978 : {
979 25580386 : sx = sub( sx, 1 );
980 : }
981 :
982 71151492 : if ( LT_16( y16, x16 ) )
983 : {
984 25580386 : x16 = mult_r( x16, 0x4000 );
985 : }
986 :
987 :
988 71151492 : z = div_s( x16, y16 );
989 71151492 : move16();
990 71151492 : *s = sub( sy, sx );
991 :
992 :
993 71151492 : return ( z );
994 : }
995 :
996 27563476 : Word32 div_w( Word32 L_num, Word32 L_den )
997 : {
998 27563476 : Word32 L_var_out = 0;
999 : Word16 iteration;
1000 27563476 : move32();
1001 :
1002 :
1003 27563476 : IF( L_den == 0 )
1004 : {
1005 : /* printf("Division by 0 in div_l, Fatal error in "); printStack(); */
1006 0 : return ( 0 );
1007 : }
1008 :
1009 27563476 : test();
1010 27563476 : IF( ( L_num < 0 ) || ( L_den < 0 ) )
1011 : {
1012 : /* printf("Division Error in div_l, Fatal error in "); printStack(); */
1013 0 : return ( 0 );
1014 : }
1015 : Word64 W_num, W_den;
1016 27563476 : W_num = W_deposit32_h( L_num );
1017 27563476 : W_den = W_deposit32_h( L_den );
1018 :
1019 27563476 : IF( W_sub( W_num, W_den ) >= 0 )
1020 : {
1021 0 : return MAX_32;
1022 : }
1023 : ELSE
1024 : {
1025 27563476 : W_num = W_shr( W_num, 1 );
1026 27563476 : W_den = W_shr( W_den, 1 );
1027 :
1028 882031232 : FOR( iteration = 0; iteration < 31; iteration++ )
1029 : {
1030 854467756 : L_var_out = L_shl( L_var_out, 1 );
1031 854467756 : W_num = W_shl( W_num, 1 );
1032 :
1033 854467756 : IF( W_sub( W_num, W_den ) >= 0 )
1034 : {
1035 506714801 : W_num = W_sub( W_num, W_den );
1036 506714801 : L_var_out = L_add( L_var_out, 1 );
1037 : }
1038 : }
1039 :
1040 27563476 : return L_var_out;
1041 : }
1042 : }
1043 :
1044 28038160 : Word32 BASOP_Util_Divide3232_Scale_cadence( Word32 x, Word32 y, Word16 *s )
1045 : {
1046 : Word32 z;
1047 : Word16 sx;
1048 : Word16 sy;
1049 : Word32 sign;
1050 :
1051 : /* assert (x >= (Word32)0); */
1052 28038160 : assert( y != (Word32) 0 );
1053 :
1054 28038160 : sign = 0;
1055 28038160 : move16();
1056 :
1057 28038160 : IF( x < 0 )
1058 : {
1059 0 : x = L_negate( x );
1060 0 : sign = L_xor( sign, 1 );
1061 : }
1062 :
1063 28038160 : IF( y < 0 )
1064 : {
1065 0 : y = L_negate( y );
1066 0 : sign = L_xor( sign, 1 );
1067 : }
1068 :
1069 28038160 : IF( x == (Word32) 0 )
1070 : {
1071 475309 : *s = 0;
1072 475309 : return ( (Word32) 0 );
1073 : }
1074 :
1075 27562851 : sx = norm_l( x );
1076 27562851 : x = L_shl( x, sx );
1077 27562851 : x = L_shr( x, 1 );
1078 27562851 : move16();
1079 27562851 : *s = sub( 1, sx );
1080 :
1081 27562851 : sy = norm_l( y );
1082 27562851 : y = L_shl( y, sy );
1083 27562851 : move16();
1084 27562851 : *s = add( *s, sy );
1085 :
1086 27562851 : z = div_w( x, y );
1087 :
1088 27562851 : if ( sign != 0 )
1089 : {
1090 0 : z = L_negate( z );
1091 : }
1092 :
1093 27562851 : return z;
1094 : }
1095 :
1096 : Word32 div_w_newton( Word32 num, Word32 den );
1097 : /*
1098 : Table of 256 precalculated estimates to be used by the "div_w_newton"
1099 : function using the Newton/Raphson method.
1100 : Note: The first table value (for denominator near 0x40000000) is not fully
1101 : accurate and should not be used.
1102 : */
1103 : Word32 div_w_newton_lookup[256] = {
1104 : /* Precalculated rounded results for 0x40000000 / b with b in [0x40000000 ... 0x7FFFFFFF] */
1105 : 0x7FFFFFFF, // 1.000000000000000 i=0 0.5 / 0.5+0/512 (b=0x40000000)
1106 : 0x7F807F80, // 0.996108949416342 i=1 0.5 / 0.5+1/512 (b=0x40400000)
1107 : 0x7F01FC07, // 0.992248062015504 i=2 0.5 / 0.5+2/512 (b=0x40800000)
1108 : 0x7E8472A8, // 0.988416988416988 i=3 0.5 / 0.5+3/512 (b=0x40C00000)
1109 : 0x7E07E07E, // 0.984615384615385 i=4 0.5 / 0.5+4/512 (b=0x41000000)
1110 : 0x7D8C42B2, // 0.980842911877395 i=5 0.5 / 0.5+5/512 (b=0x41400000)
1111 : 0x7D119679, // 0.977099236641221 i=6 0.5 / 0.5+6/512 (b=0x41800000)
1112 : 0x7C97D910, // 0.973384030418251 i=7 0.5 / 0.5+7/512 (b=0x41C00000)
1113 : 0x7C1F07C1, // 0.969696969696970 i=8 0.5 / 0.5+8/512 (b=0x42000000)
1114 : 0x7BA71FE1, // 0.966037735849057 i=9 0.5 / 0.5+9/512 (b=0x42400000)
1115 : 0x7B301ECC, // 0.962406015037594 i=10 0.5 / 0.5+10/512 (b=0x42800000)
1116 : 0x7ABA01EA, // 0.958801498127341 i=11 0.5 / 0.5+11/512 (b=0x42C00000)
1117 : 0x7A44C6AF, // 0.955223880597015 i=12 0.5 / 0.5+12/512 (b=0x43000000)
1118 : 0x79D06A96, // 0.951672862453532 i=13 0.5 / 0.5+13/512 (b=0x43400000)
1119 : 0x795CEB24, // 0.948148148148148 i=14 0.5 / 0.5+14/512 (b=0x43800000)
1120 : 0x78EA45E7, // 0.944649446494465 i=15 0.5 / 0.5+15/512 (b=0x43C00000)
1121 : 0x78787878, // 0.941176470588235 i=16 0.5 / 0.5+16/512 (b=0x44000000)
1122 : 0x78078078, // 0.937728937728938 i=17 0.5 / 0.5+17/512 (b=0x44400000)
1123 : 0x77975B8F, // 0.934306569343066 i=18 0.5 / 0.5+18/512 (b=0x44800000)
1124 : 0x77280772, // 0.930909090909091 i=19 0.5 / 0.5+19/512 (b=0x44C00000)
1125 : 0x76B981DA, // 0.927536231884058 i=20 0.5 / 0.5+20/512 (b=0x45000000)
1126 : 0x764BC88C, // 0.924187725631769 i=21 0.5 / 0.5+21/512 (b=0x45400000)
1127 : 0x75DED952, // 0.920863309352518 i=22 0.5 / 0.5+22/512 (b=0x45800000)
1128 : 0x7572B201, // 0.917562724014337 i=23 0.5 / 0.5+23/512 (b=0x45C00000)
1129 : 0x75075075, // 0.914285714285714 i=24 0.5 / 0.5+24/512 (b=0x46000000)
1130 : 0x749CB28F, // 0.911032028469751 i=25 0.5 / 0.5+25/512 (b=0x46400000)
1131 : 0x7432D63D, // 0.907801418439716 i=26 0.5 / 0.5+26/512 (b=0x46800000)
1132 : 0x73C9B971, // 0.904593639575972 i=27 0.5 / 0.5+27/512 (b=0x46C00000)
1133 : 0x73615A24, // 0.901408450704225 i=28 0.5 / 0.5+28/512 (b=0x47000000)
1134 : 0x72F9B658, // 0.898245614035088 i=29 0.5 / 0.5+29/512 (b=0x47400000)
1135 : 0x7292CC15, // 0.895104895104895 i=30 0.5 / 0.5+30/512 (b=0x47800000)
1136 : 0x722C996B, // 0.891986062717770 i=31 0.5 / 0.5+31/512 (b=0x47C00000)
1137 : 0x71C71C71, // 0.888888888888889 i=32 0.5 / 0.5+32/512 (b=0x48000000)
1138 : 0x71625344, // 0.885813148788927 i=33 0.5 / 0.5+33/512 (b=0x48400000)
1139 : 0x70FE3C07, // 0.882758620689655 i=34 0.5 / 0.5+34/512 (b=0x48800000)
1140 : 0x709AD4E4, // 0.879725085910653 i=35 0.5 / 0.5+35/512 (b=0x48C00000)
1141 : 0x70381C0E, // 0.876712328767123 i=36 0.5 / 0.5+36/512 (b=0x49000000)
1142 : 0x6FD60FBA, // 0.873720136518771 i=37 0.5 / 0.5+37/512 (b=0x49400000)
1143 : 0x6F74AE26, // 0.870748299319728 i=38 0.5 / 0.5+38/512 (b=0x49800000)
1144 : 0x6F13F596, // 0.867796610169492 i=39 0.5 / 0.5+39/512 (b=0x49C00000)
1145 : 0x6EB3E453, // 0.864864864864865 i=40 0.5 / 0.5+40/512 (b=0x4A000000)
1146 : 0x6E5478AC, // 0.861952861952862 i=41 0.5 / 0.5+41/512 (b=0x4A400000)
1147 : 0x6DF5B0F7, // 0.859060402684564 i=42 0.5 / 0.5+42/512 (b=0x4A800000)
1148 : 0x6D978B8E, // 0.856187290969900 i=43 0.5 / 0.5+43/512 (b=0x4AC00000)
1149 : 0x6D3A06D3, // 0.853333333333333 i=44 0.5 / 0.5+44/512 (b=0x4B000000)
1150 : 0x6CDD212B, // 0.850498338870432 i=45 0.5 / 0.5+45/512 (b=0x4B400000)
1151 : 0x6C80D901, // 0.847682119205298 i=46 0.5 / 0.5+46/512 (b=0x4B800000)
1152 : 0x6C252CC7, // 0.844884488448845 i=47 0.5 / 0.5+47/512 (b=0x4BC00000)
1153 : 0x6BCA1AF2, // 0.842105263157895 i=48 0.5 / 0.5+48/512 (b=0x4C000000)
1154 : 0x6B6FA1FE, // 0.839344262295082 i=49 0.5 / 0.5+49/512 (b=0x4C400000)
1155 : 0x6B15C06B, // 0.836601307189543 i=50 0.5 / 0.5+50/512 (b=0x4C800000)
1156 : 0x6ABC74BE, // 0.833876221498371 i=51 0.5 / 0.5+51/512 (b=0x4CC00000)
1157 : 0x6A63BD81, // 0.831168831168831 i=52 0.5 / 0.5+52/512 (b=0x4D000000)
1158 : 0x6A0B9944, // 0.828478964401295 i=53 0.5 / 0.5+53/512 (b=0x4D400000)
1159 : 0x69B4069B, // 0.825806451612903 i=54 0.5 / 0.5+54/512 (b=0x4D800000)
1160 : 0x695D041D, // 0.823151125401929 i=55 0.5 / 0.5+55/512 (b=0x4DC00000)
1161 : 0x69069069, // 0.820512820512820 i=56 0.5 / 0.5+56/512 (b=0x4E000000)
1162 : 0x68B0AA1F, // 0.817891373801917 i=57 0.5 / 0.5+57/512 (b=0x4E400000)
1163 : 0x685B4FE5, // 0.815286624203822 i=58 0.5 / 0.5+58/512 (b=0x4E800000)
1164 : 0x68068068, // 0.812698412698413 i=59 0.5 / 0.5+59/512 (b=0x4EC00000)
1165 : 0x67B23A54, // 0.810126582278481 i=60 0.5 / 0.5+60/512 (b=0x4F000000)
1166 : 0x675E7C5D, // 0.807570977917981 i=61 0.5 / 0.5+61/512 (b=0x4F400000)
1167 : 0x670B453B, // 0.805031446540881 i=62 0.5 / 0.5+62/512 (b=0x4F800000)
1168 : 0x66B893A9, // 0.802507836990596 i=63 0.5 / 0.5+63/512 (b=0x4FC00000)
1169 : 0x66666666, // 0.800000000000000 i=64 0.5 / 0.5+64/512 (b=0x50000000)
1170 : 0x6614BC36, // 0.797507788161994 i=65 0.5 / 0.5+65/512 (b=0x50400000)
1171 : 0x65C393E0, // 0.795031055900621 i=66 0.5 / 0.5+66/512 (b=0x50800000)
1172 : 0x6572EC2F, // 0.792569659442725 i=67 0.5 / 0.5+67/512 (b=0x50C00000)
1173 : 0x6522C3F3, // 0.790123456790123 i=68 0.5 / 0.5+68/512 (b=0x51000000)
1174 : 0x64D319FE, // 0.787692307692308 i=69 0.5 / 0.5+69/512 (b=0x51400000)
1175 : 0x6483ED27, // 0.785276073619632 i=70 0.5 / 0.5+70/512 (b=0x51800000)
1176 : 0x64353C48, // 0.782874617737003 i=71 0.5 / 0.5+71/512 (b=0x51C00000)
1177 : 0x63E7063E, // 0.780487804878049 i=72 0.5 / 0.5+72/512 (b=0x52000000)
1178 : 0x639949EB, // 0.778115501519757 i=73 0.5 / 0.5+73/512 (b=0x52400000)
1179 : 0x634C0634, // 0.775757575757576 i=74 0.5 / 0.5+74/512 (b=0x52800000)
1180 : 0x62FF3A01, // 0.773413897280967 i=75 0.5 / 0.5+75/512 (b=0x52C00000)
1181 : 0x62B2E43D, // 0.771084337349398 i=76 0.5 / 0.5+76/512 (b=0x53000000)
1182 : 0x626703D8, // 0.768768768768769 i=77 0.5 / 0.5+77/512 (b=0x53400000)
1183 : 0x621B97C2, // 0.766467065868264 i=78 0.5 / 0.5+78/512 (b=0x53800000)
1184 : 0x61D09EF3, // 0.764179104477612 i=79 0.5 / 0.5+79/512 (b=0x53C00000)
1185 : 0x61861861, // 0.761904761904762 i=80 0.5 / 0.5+80/512 (b=0x54000000)
1186 : 0x613C0309, // 0.759643916913947 i=81 0.5 / 0.5+81/512 (b=0x54400000)
1187 : 0x60F25DEA, // 0.757396449704142 i=82 0.5 / 0.5+82/512 (b=0x54800000)
1188 : 0x60A92806, // 0.755162241887906 i=83 0.5 / 0.5+83/512 (b=0x54C00000)
1189 : 0x60606060, // 0.752941176470588 i=84 0.5 / 0.5+84/512 (b=0x55000000)
1190 : 0x60180601, // 0.750733137829912 i=85 0.5 / 0.5+85/512 (b=0x55400000)
1191 : 0x5FD017F4, // 0.748538011695906 i=86 0.5 / 0.5+86/512 (b=0x55800000)
1192 : 0x5F889545, // 0.746355685131195 i=87 0.5 / 0.5+87/512 (b=0x55C00000)
1193 : 0x5F417D05, // 0.744186046511628 i=88 0.5 / 0.5+88/512 (b=0x56000000)
1194 : 0x5EFACE48, // 0.742028985507246 i=89 0.5 / 0.5+89/512 (b=0x56400000)
1195 : 0x5EB48823, // 0.739884393063584 i=90 0.5 / 0.5+90/512 (b=0x56800000)
1196 : 0x5E6EA9AE, // 0.737752161383285 i=91 0.5 / 0.5+91/512 (b=0x56C00000)
1197 : 0x5E293205, // 0.735632183908046 i=92 0.5 / 0.5+92/512 (b=0x57000000)
1198 : 0x5DE42046, // 0.733524355300860 i=93 0.5 / 0.5+93/512 (b=0x57400000)
1199 : 0x5D9F7390, // 0.731428571428571 i=94 0.5 / 0.5+94/512 (b=0x57800000)
1200 : 0x5D5B2B08, // 0.729344729344729 i=95 0.5 / 0.5+95/512 (b=0x57C00000)
1201 : 0x5D1745D1, // 0.727272727272727 i=96 0.5 / 0.5+96/512 (b=0x58000000)
1202 : 0x5CD3C315, // 0.725212464589235 i=97 0.5 / 0.5+97/512 (b=0x58400000)
1203 : 0x5C90A1FD, // 0.723163841807910 i=98 0.5 / 0.5+98/512 (b=0x58800000)
1204 : 0x5C4DE1B6, // 0.721126760563380 i=99 0.5 / 0.5+99/512 (b=0x58C00000)
1205 : 0x5C0B8170, // 0.719101123595506 i=100 0.5 / 0.5+100/512 (b=0x59000000)
1206 : 0x5BC9805B, // 0.717086834733894 i=101 0.5 / 0.5+101/512 (b=0x59400000)
1207 : 0x5B87DDAD, // 0.715083798882682 i=102 0.5 / 0.5+102/512 (b=0x59800000)
1208 : 0x5B46989A, // 0.713091922005571 i=103 0.5 / 0.5+103/512 (b=0x59C00000)
1209 : 0x5B05B05B, // 0.711111111111111 i=104 0.5 / 0.5+104/512 (b=0x5A000000)
1210 : 0x5AC5242A, // 0.709141274238227 i=105 0.5 / 0.5+105/512 (b=0x5A400000)
1211 : 0x5A84F345, // 0.707182320441989 i=106 0.5 / 0.5+106/512 (b=0x5A800000)
1212 : 0x5A451CEA, // 0.705234159779614 i=107 0.5 / 0.5+107/512 (b=0x5AC00000)
1213 : 0x5A05A05A, // 0.703296703296703 i=108 0.5 / 0.5+108/512 (b=0x5B000000)
1214 : 0x59C67CD8, // 0.701369863013699 i=109 0.5 / 0.5+109/512 (b=0x5B400000)
1215 : 0x5987B1A9, // 0.699453551912568 i=110 0.5 / 0.5+110/512 (b=0x5B800000)
1216 : 0x59493E14, // 0.697547683923706 i=111 0.5 / 0.5+111/512 (b=0x5BC00000)
1217 : 0x590B2164, // 0.695652173913043 i=112 0.5 / 0.5+112/512 (b=0x5C000000)
1218 : 0x58CD5AE2, // 0.693766937669377 i=113 0.5 / 0.5+113/512 (b=0x5C400000)
1219 : 0x588FE9DC, // 0.691891891891892 i=114 0.5 / 0.5+114/512 (b=0x5C800000)
1220 : 0x5852CDA0, // 0.690026954177898 i=115 0.5 / 0.5+115/512 (b=0x5CC00000)
1221 : 0x58160581, // 0.688172043010753 i=116 0.5 / 0.5+116/512 (b=0x5D000000)
1222 : 0x57D990D0, // 0.686327077747989 i=117 0.5 / 0.5+117/512 (b=0x5D400000)
1223 : 0x579D6EE3, // 0.684491978609626 i=118 0.5 / 0.5+118/512 (b=0x5D800000)
1224 : 0x57619F0F, // 0.682666666666667 i=119 0.5 / 0.5+119/512 (b=0x5DC00000)
1225 : 0x572620AE, // 0.680851063829787 i=120 0.5 / 0.5+120/512 (b=0x5E000000)
1226 : 0x56EAF319, // 0.679045092838196 i=121 0.5 / 0.5+121/512 (b=0x5E400000)
1227 : 0x56B015AC, // 0.677248677248677 i=122 0.5 / 0.5+122/512 (b=0x5E800000)
1228 : 0x567587C4, // 0.675461741424802 i=123 0.5 / 0.5+123/512 (b=0x5EC00000)
1229 : 0x563B48C2, // 0.673684210526316 i=124 0.5 / 0.5+124/512 (b=0x5F000000)
1230 : 0x56015805, // 0.671916010498688 i=125 0.5 / 0.5+125/512 (b=0x5F400000)
1231 : 0x55C7B4F1, // 0.670157068062827 i=126 0.5 / 0.5+126/512 (b=0x5F800000)
1232 : 0x558E5EE9, // 0.668407310704961 i=127 0.5 / 0.5+127/512 (b=0x5FC00000)
1233 : 0x55555555, // 0.666666666666667 i=128 0.5 / 0.5+128/512 (b=0x60000000)
1234 : 0x551C979A, // 0.664935064935065 i=129 0.5 / 0.5+129/512 (b=0x60400000)
1235 : 0x54E42523, // 0.663212435233161 i=130 0.5 / 0.5+130/512 (b=0x60800000)
1236 : 0x54ABFD5A, // 0.661498708010336 i=131 0.5 / 0.5+131/512 (b=0x60C00000)
1237 : 0x54741FAB, // 0.659793814432990 i=132 0.5 / 0.5+132/512 (b=0x61000000)
1238 : 0x543C8B84, // 0.658097686375321 i=133 0.5 / 0.5+133/512 (b=0x61400000)
1239 : 0x54054054, // 0.656410256410256 i=134 0.5 / 0.5+134/512 (b=0x61800000)
1240 : 0x53CE3D8B, // 0.654731457800512 i=135 0.5 / 0.5+135/512 (b=0x61C00000)
1241 : 0x5397829C, // 0.653061224489796 i=136 0.5 / 0.5+136/512 (b=0x62000000)
1242 : 0x53610EFB, // 0.651399491094148 i=137 0.5 / 0.5+137/512 (b=0x62400000)
1243 : 0x532AE21C, // 0.649746192893401 i=138 0.5 / 0.5+138/512 (b=0x62800000)
1244 : 0x52F4FB76, // 0.648101265822785 i=139 0.5 / 0.5+139/512 (b=0x62C00000)
1245 : 0x52BF5A81, // 0.646464646464647 i=140 0.5 / 0.5+140/512 (b=0x63000000)
1246 : 0x5289FEB5, // 0.644836272040302 i=141 0.5 / 0.5+141/512 (b=0x63400000)
1247 : 0x5254E78E, // 0.643216080402010 i=142 0.5 / 0.5+142/512 (b=0x63800000)
1248 : 0x52201488, // 0.641604010025063 i=143 0.5 / 0.5+143/512 (b=0x63C00000)
1249 : 0x51EB851E, // 0.640000000000000 i=144 0.5 / 0.5+144/512 (b=0x64000000)
1250 : 0x51B738D1, // 0.638403990024938 i=145 0.5 / 0.5+145/512 (b=0x64400000)
1251 : 0x51832F1F, // 0.636815920398010 i=146 0.5 / 0.5+146/512 (b=0x64800000)
1252 : 0x514F678B, // 0.635235732009926 i=147 0.5 / 0.5+147/512 (b=0x64C00000)
1253 : 0x511BE195, // 0.633663366336634 i=148 0.5 / 0.5+148/512 (b=0x65000000)
1254 : 0x50E89CC2, // 0.632098765432099 i=149 0.5 / 0.5+149/512 (b=0x65400000)
1255 : 0x50B59897, // 0.630541871921182 i=150 0.5 / 0.5+150/512 (b=0x65800000)
1256 : 0x5082D499, // 0.628992628992629 i=151 0.5 / 0.5+151/512 (b=0x65C00000)
1257 : 0x50505050, // 0.627450980392157 i=152 0.5 / 0.5+152/512 (b=0x66000000)
1258 : 0x501E0B44, // 0.625916870415648 i=153 0.5 / 0.5+153/512 (b=0x66400000)
1259 : 0x4FEC04FE, // 0.624390243902439 i=154 0.5 / 0.5+154/512 (b=0x66800000)
1260 : 0x4FBA3D0A, // 0.622871046228710 i=155 0.5 / 0.5+155/512 (b=0x66C00000)
1261 : 0x4F88B2F3, // 0.621359223300971 i=156 0.5 / 0.5+156/512 (b=0x67000000)
1262 : 0x4F576646, // 0.619854721549637 i=157 0.5 / 0.5+157/512 (b=0x67400000)
1263 : 0x4F265691, // 0.618357487922705 i=158 0.5 / 0.5+158/512 (b=0x67800000)
1264 : 0x4EF58364, // 0.616867469879518 i=159 0.5 / 0.5+159/512 (b=0x67C00000)
1265 : 0x4EC4EC4E, // 0.615384615384615 i=160 0.5 / 0.5+160/512 (b=0x68000000)
1266 : 0x4E9490E1, // 0.613908872901679 i=161 0.5 / 0.5+161/512 (b=0x68400000)
1267 : 0x4E6470B0, // 0.612440191387560 i=162 0.5 / 0.5+162/512 (b=0x68800000)
1268 : 0x4E348B4D, // 0.610978520286396 i=163 0.5 / 0.5+163/512 (b=0x68C00000)
1269 : 0x4E04E04E, // 0.609523809523810 i=164 0.5 / 0.5+164/512 (b=0x69000000)
1270 : 0x4DD56F47, // 0.608076009501188 i=165 0.5 / 0.5+165/512 (b=0x69400000)
1271 : 0x4DA637CF, // 0.606635071090047 i=166 0.5 / 0.5+166/512 (b=0x69800000)
1272 : 0x4D77397E, // 0.605200945626478 i=167 0.5 / 0.5+167/512 (b=0x69C00000)
1273 : 0x4D4873EC, // 0.603773584905660 i=168 0.5 / 0.5+168/512 (b=0x6A000000)
1274 : 0x4D19E6B3, // 0.602352941176471 i=169 0.5 / 0.5+169/512 (b=0x6A400000)
1275 : 0x4CEB916D, // 0.600938967136150 i=170 0.5 / 0.5+170/512 (b=0x6A800000)
1276 : 0x4CBD73B5, // 0.599531615925059 i=171 0.5 / 0.5+171/512 (b=0x6AC00000)
1277 : 0x4C8F8D28, // 0.598130841121495 i=172 0.5 / 0.5+172/512 (b=0x6B000000)
1278 : 0x4C61DD63, // 0.596736596736597 i=173 0.5 / 0.5+173/512 (b=0x6B400000)
1279 : 0x4C346404, // 0.595348837209302 i=174 0.5 / 0.5+174/512 (b=0x6B800000)
1280 : 0x4C0720AB, // 0.593967517401392 i=175 0.5 / 0.5+175/512 (b=0x6BC00000)
1281 : 0x4BDA12F6, // 0.592592592592593 i=176 0.5 / 0.5+176/512 (b=0x6C000000)
1282 : 0x4BAD3A87, // 0.591224018475751 i=177 0.5 / 0.5+177/512 (b=0x6C400000)
1283 : 0x4B809701, // 0.589861751152074 i=178 0.5 / 0.5+178/512 (b=0x6C800000)
1284 : 0x4B542804, // 0.588505747126437 i=179 0.5 / 0.5+179/512 (b=0x6CC00000)
1285 : 0x4B27ED36, // 0.587155963302752 i=180 0.5 / 0.5+180/512 (b=0x6D000000)
1286 : 0x4AFBE639, // 0.585812356979405 i=181 0.5 / 0.5+181/512 (b=0x6D400000)
1287 : 0x4AD012B4, // 0.584474885844749 i=182 0.5 / 0.5+182/512 (b=0x6D800000)
1288 : 0x4AA4724B, // 0.583143507972665 i=183 0.5 / 0.5+183/512 (b=0x6DC00000)
1289 : 0x4A7904A7, // 0.581818181818182 i=184 0.5 / 0.5+184/512 (b=0x6E000000)
1290 : 0x4A4DC96E, // 0.580498866213152 i=185 0.5 / 0.5+185/512 (b=0x6E400000)
1291 : 0x4A22C04A, // 0.579185520361991 i=186 0.5 / 0.5+186/512 (b=0x6E800000)
1292 : 0x49F7E8E2, // 0.577878103837472 i=187 0.5 / 0.5+187/512 (b=0x6EC00000)
1293 : 0x49CD42E2, // 0.576576576576577 i=188 0.5 / 0.5+188/512 (b=0x6F000000)
1294 : 0x49A2CDF3, // 0.575280898876405 i=189 0.5 / 0.5+189/512 (b=0x6F400000)
1295 : 0x497889C2, // 0.573991031390135 i=190 0.5 / 0.5+190/512 (b=0x6F800000)
1296 : 0x494E75FA, // 0.572706935123042 i=191 0.5 / 0.5+191/512 (b=0x6FC00000)
1297 : 0x49249249, // 0.571428571428571 i=192 0.5 / 0.5+192/512 (b=0x70000000)
1298 : 0x48FADE5C, // 0.570155902004454 i=193 0.5 / 0.5+193/512 (b=0x70400000)
1299 : 0x48D159E2, // 0.568888888888889 i=194 0.5 / 0.5+194/512 (b=0x70800000)
1300 : 0x48A8048A, // 0.567627494456763 i=195 0.5 / 0.5+195/512 (b=0x70C00000)
1301 : 0x487EDE04, // 0.566371681415929 i=196 0.5 / 0.5+196/512 (b=0x71000000)
1302 : 0x4855E601, // 0.565121412803532 i=197 0.5 / 0.5+197/512 (b=0x71400000)
1303 : 0x482D1C31, // 0.563876651982379 i=198 0.5 / 0.5+198/512 (b=0x71800000)
1304 : 0x48048048, // 0.562637362637363 i=199 0.5 / 0.5+199/512 (b=0x71C00000)
1305 : 0x47DC11F7, // 0.561403508771930 i=200 0.5 / 0.5+200/512 (b=0x72000000)
1306 : 0x47B3D0F1, // 0.560175054704595 i=201 0.5 / 0.5+201/512 (b=0x72400000)
1307 : 0x478BBCEC, // 0.558951965065502 i=202 0.5 / 0.5+202/512 (b=0x72800000)
1308 : 0x4763D59C, // 0.557734204793028 i=203 0.5 / 0.5+203/512 (b=0x72C00000)
1309 : 0x473C1AB6, // 0.556521739130435 i=204 0.5 / 0.5+204/512 (b=0x73000000)
1310 : 0x47148BF0, // 0.555314533622560 i=205 0.5 / 0.5+205/512 (b=0x73400000)
1311 : 0x46ED2901, // 0.554112554112554 i=206 0.5 / 0.5+206/512 (b=0x73800000)
1312 : 0x46C5F19F, // 0.552915766738661 i=207 0.5 / 0.5+207/512 (b=0x73C00000)
1313 : 0x469EE584, // 0.551724137931034 i=208 0.5 / 0.5+208/512 (b=0x74000000)
1314 : 0x46780467, // 0.550537634408602 i=209 0.5 / 0.5+209/512 (b=0x74400000)
1315 : 0x46514E02, // 0.549356223175966 i=210 0.5 / 0.5+210/512 (b=0x74800000)
1316 : 0x462AC20E, // 0.548179871520343 i=211 0.5 / 0.5+211/512 (b=0x74C00000)
1317 : 0x46046046, // 0.547008547008547 i=212 0.5 / 0.5+212/512 (b=0x75000000)
1318 : 0x45DE2864, // 0.545842217484009 i=213 0.5 / 0.5+213/512 (b=0x75400000)
1319 : 0x45B81A25, // 0.544680851063830 i=214 0.5 / 0.5+214/512 (b=0x75800000)
1320 : 0x45923543, // 0.543524416135881 i=215 0.5 / 0.5+215/512 (b=0x75C00000)
1321 : 0x456C797D, // 0.542372881355932 i=216 0.5 / 0.5+216/512 (b=0x76000000)
1322 : 0x4546E68F, // 0.541226215644820 i=217 0.5 / 0.5+217/512 (b=0x76400000)
1323 : 0x45217C38, // 0.540084388185654 i=218 0.5 / 0.5+218/512 (b=0x76800000)
1324 : 0x44FC3A34, // 0.538947368421053 i=219 0.5 / 0.5+219/512 (b=0x76C00000)
1325 : 0x44D72044, // 0.537815126050420 i=220 0.5 / 0.5+220/512 (b=0x77000000)
1326 : 0x44B22E27, // 0.536687631027254 i=221 0.5 / 0.5+221/512 (b=0x77400000)
1327 : 0x448D639D, // 0.535564853556485 i=222 0.5 / 0.5+222/512 (b=0x77800000)
1328 : 0x4468C066, // 0.534446764091858 i=223 0.5 / 0.5+223/512 (b=0x77C00000)
1329 : 0x44444444, // 0.533333333333333 i=224 0.5 / 0.5+224/512 (b=0x78000000)
1330 : 0x441FEEF8, // 0.532224532224532 i=225 0.5 / 0.5+225/512 (b=0x78400000)
1331 : 0x43FBC043, // 0.531120331950207 i=226 0.5 / 0.5+226/512 (b=0x78800000)
1332 : 0x43D7B7EA, // 0.530020703933747 i=227 0.5 / 0.5+227/512 (b=0x78C00000)
1333 : 0x43B3D5AF, // 0.528925619834711 i=228 0.5 / 0.5+228/512 (b=0x79000000)
1334 : 0x43901956, // 0.527835051546392 i=229 0.5 / 0.5+229/512 (b=0x79400000)
1335 : 0x436C82A2, // 0.526748971193416 i=230 0.5 / 0.5+230/512 (b=0x79800000)
1336 : 0x43491158, // 0.525667351129363 i=231 0.5 / 0.5+231/512 (b=0x79C00000)
1337 : 0x4325C53E, // 0.524590163934426 i=232 0.5 / 0.5+232/512 (b=0x7A000000)
1338 : 0x43029E1A, // 0.523517382413088 i=233 0.5 / 0.5+233/512 (b=0x7A400000)
1339 : 0x42DF9BB0, // 0.522448979591837 i=234 0.5 / 0.5+234/512 (b=0x7A800000)
1340 : 0x42BCBDC8, // 0.521384928716904 i=235 0.5 / 0.5+235/512 (b=0x7AC00000)
1341 : 0x429A0429, // 0.520325203252033 i=236 0.5 / 0.5+236/512 (b=0x7B000000)
1342 : 0x42776E9A, // 0.519269776876268 i=237 0.5 / 0.5+237/512 (b=0x7B400000)
1343 : 0x4254FCE4, // 0.518218623481781 i=238 0.5 / 0.5+238/512 (b=0x7B800000)
1344 : 0x4232AECD, // 0.517171717171717 i=239 0.5 / 0.5+239/512 (b=0x7BC00000)
1345 : 0x42108421, // 0.516129032258065 i=240 0.5 / 0.5+240/512 (b=0x7C000000)
1346 : 0x41EE7CA6, // 0.515090543259557 i=241 0.5 / 0.5+241/512 (b=0x7C400000)
1347 : 0x41CC9829, // 0.514056224899598 i=242 0.5 / 0.5+242/512 (b=0x7C800000)
1348 : 0x41AAD671, // 0.513026052104208 i=243 0.5 / 0.5+243/512 (b=0x7CC00000)
1349 : 0x4189374B, // 0.512000000000000 i=244 0.5 / 0.5+244/512 (b=0x7D000000)
1350 : 0x4167BA81, // 0.510978043912176 i=245 0.5 / 0.5+245/512 (b=0x7D400000)
1351 : 0x41465FDF, // 0.509960159362550 i=246 0.5 / 0.5+246/512 (b=0x7D800000)
1352 : 0x41252730, // 0.508946322067594 i=247 0.5 / 0.5+247/512 (b=0x7DC00000)
1353 : 0x41041041, // 0.507936507936508 i=248 0.5 / 0.5+248/512 (b=0x7E000000)
1354 : 0x40E31ADE, // 0.506930693069307 i=249 0.5 / 0.5+249/512 (b=0x7E400000)
1355 : 0x40C246D4, // 0.505928853754941 i=250 0.5 / 0.5+250/512 (b=0x7E800000)
1356 : 0x40A193F1, // 0.504930966469428 i=251 0.5 / 0.5+251/512 (b=0x7EC00000)
1357 : 0x40810204, // 0.503937007874016 i=252 0.5 / 0.5+252/512 (b=0x7F000000)
1358 : 0x406090D9, // 0.502946954813359 i=253 0.5 / 0.5+253/512 (b=0x7F400000)
1359 : 0x40404040, // 0.501960784313725 i=254 0.5 / 0.5+254/512 (b=0x7F800000)
1360 : 0x40201008 // 0.500978473581213 i=255 0.5 / 0.5+255/512 (b=0x7FC00000)
1361 : };
1362 :
1363 :
1364 : /*
1365 : * Fractional multiplication of signed a and b, both in Q31. The result is doubled.
1366 : * Note: in this test, saturation is not needed.
1367 : * BASOP weights: 3
1368 : */
1369 :
1370 1825107222 : static Word32 L_dmult( Word32 L_var1, Word32 L_var2 )
1371 : {
1372 1825107222 : Word64 L64_var1 = W_mult0_32_32( L_var1, L_var2 );
1373 1825107222 : L64_var1 = W_shr( L64_var1, 30 );
1374 1825107222 : return W_extract_l( L64_var1 );
1375 : }
1376 :
1377 : /*
1378 : * 32 by 32 bit division, following the Newton / Raphson method.
1379 : * Usage of this low-level procedure by the caller:
1380 : * 1. Numerator <num> can use the full range of signed 32-bit datatypes: 0x80000000...0x7FFFFFFF
1381 : * Note: Since <num> is not normalized here, but it is multplied to the reciprocal of the
1382 : * denominator, the caller should use a normalized <num> and handle any exponent outside.
1383 : * 2. Denominator <den> must be normalized into range 0x40000001 to 0x7FFFFFFF (all positive)
1384 : * Note: In case of den=0x40000000, the caller is not allowed to call the division routine,
1385 : since the result is known.
1386 : * Note: num / 0x40000000 equals to num with exp += 1.
1387 : * 3. The result is in range 0x40000000 to 0x7FFFFFFF, finally multiplied by <num>.
1388 : * BASOP weights: 24 (incl. L_dmult)
1389 : */
1390 :
1391 608369074 : Word32 div_w_newton( Word32 num, Word32 den )
1392 : {
1393 : Word32 x0, x1, x2, x3, diff, result;
1394 :
1395 608369074 : x0 = div_w_newton_lookup[sub( extract_l( L_shr( den, 22 ) ), 256 )];
1396 608369074 : move32();
1397 :
1398 608369074 : diff = L_sub( 0x40000000, Mpy_32_32( den, x0 ) );
1399 :
1400 608369074 : x1 = L_add( x0, L_dmult( x0, diff ) );
1401 608369074 : diff = L_sub( 0x40000000, Mpy_32_32( den, x1 ) );
1402 :
1403 608369074 : x2 = L_add( x1, L_dmult( x1, diff ) );
1404 608369074 : diff = L_sub( 0x40000000, Mpy_32_32( den, x2 ) );
1405 :
1406 608369074 : x3 = L_add( x2, L_dmult( x2, diff ) );
1407 :
1408 608369074 : result = Mpy_32_32( num, x3 );
1409 :
1410 608369074 : return result;
1411 : }
1412 :
1413 : /*
1414 : * 32 / 32 division
1415 : * Usage of this global procedure by the caller:
1416 : * 1. Numerator <num> can use the full range of signed 32-bit datatypes: 0x80000000...0x7FFFFFFF
1417 : * Note: Since <num> is not normalized here, but it is multplied to the reciprocal of the
1418 : * denominator, the caller should use a normalized <num> and handle any exponent outside.
1419 : * 2. Denominator can use the full range of signed 32-bit datatypes: 0x80000000...0x7FFFFFFF
1420 : * except den=0x00000000. In case of 0x80000000, it becomes internally negated to 0x7FFFFFFF.
1421 : * 3. The result is 0x00000000 (*s=0)for <num> equals 0x00000000.
1422 : * The result is rather left aligned, with up to 1 bit headroom.
1423 : * 4. The result exponent is stored in s[0]
1424 : * BASOP weights: 41 (incl. div_w_newton)
1425 : */
1426 :
1427 633426114 : Word32 BASOP_Util_Divide3232_Scale_newton( Word32 x, Word32 y, Word16 *s )
1428 : {
1429 : Word32 z;
1430 : Word16 sx;
1431 : Word16 sy;
1432 : Word32 sign;
1433 :
1434 633426114 : assert( y != (Word32) 0 );
1435 :
1436 : /* Early exit, if numerator is zero */
1437 633426114 : IF( x == (Word32) 0 )
1438 : {
1439 2580471 : *s = 0;
1440 2580471 : return ( (Word32) 0 );
1441 : }
1442 630845643 : IF( EQ_32( y, (Word32) 0x80000000 ) )
1443 : {
1444 : /* Division by -1.0: same as negation of numerator */
1445 : /* Return normalized negated numerator */
1446 98 : sx = norm_l( x );
1447 98 : x = L_shl( x, sx );
1448 98 : *s = negate( sx );
1449 98 : return L_negate( x );
1450 : }
1451 630845545 : sign = y;
1452 630845545 : move32();
1453 630845545 : if ( y < 0 )
1454 : {
1455 12911554 : y = L_negate( y );
1456 : }
1457 :
1458 : /* Normalize numerator */
1459 630845545 : sx = norm_l( x );
1460 630845545 : x = L_shl( x, sx );
1461 :
1462 : /* Normalize denominator */
1463 630845545 : sy = norm_l( y );
1464 630845545 : y = L_shl( y, sy );
1465 :
1466 : /* Store exponent: + 1 for div_w_newton computing 0.5*num/den */
1467 630845545 : *s = sub( add( sy, 1 ), sx );
1468 630845545 : move16();
1469 :
1470 : /* Special treatment for den=0x40000000 */
1471 : /* Result is known: z=2*num */
1472 630845545 : IF( EQ_32( y, 0x40000000 ) )
1473 : {
1474 22476471 : if ( sign < 0 )
1475 : {
1476 0 : x = L_negate( x );
1477 : }
1478 22476471 : return x;
1479 : }
1480 :
1481 : /* Invoke division applying Newton/Raphson-Algorithm */
1482 608369074 : z = div_w_newton( x, y );
1483 :
1484 608369074 : if ( sign < 0 )
1485 : {
1486 12911554 : z = L_negate( z );
1487 : }
1488 :
1489 608369074 : return z;
1490 : }
1491 :
1492 382243950 : Word16 BASOP_Util_Divide3232_Scale( Word32 x, Word32 y, Word16 *s )
1493 : {
1494 : Word16 z;
1495 : Word16 sy;
1496 :
1497 :
1498 382243950 : sy = norm_l( y );
1499 382243950 : if ( sy > 0 )
1500 : {
1501 296346493 : sy = sub( sy, 1 );
1502 : }
1503 382243950 : y = L_shl( y, sy );
1504 :
1505 382243950 : z = BASOP_Util_Divide3216_Scale( x, extract_h( y ), s );
1506 382243950 : move16();
1507 382243950 : *s = add( *s, sy );
1508 :
1509 :
1510 382243950 : return ( z );
1511 : }
1512 :
1513 :
1514 39445165 : Word16 BASOP_Util_Divide1616_Scale( Word16 x, Word16 y, Word16 *s )
1515 : {
1516 : Word16 z;
1517 : Word16 sx;
1518 : Word16 sy;
1519 : Word16 sign;
1520 :
1521 :
1522 : /* assert (x >= (Word16)0); */
1523 39445165 : assert( y != (Word16) 0 );
1524 :
1525 39445165 : sign = 0;
1526 39445165 : move16();
1527 :
1528 39445165 : IF( x < 0 )
1529 : {
1530 1081665 : x = negate( x );
1531 1081665 : sign = s_xor( sign, 1 );
1532 : }
1533 :
1534 39445165 : IF( y < 0 )
1535 : {
1536 6423 : y = negate( y );
1537 6423 : sign = s_xor( sign, 1 );
1538 : }
1539 :
1540 39445165 : IF( x == (Word16) 0 )
1541 : {
1542 1194881 : move16();
1543 1194881 : *s = 0;
1544 :
1545 :
1546 1194881 : return ( (Word16) 0 );
1547 : }
1548 :
1549 38250284 : sx = norm_s( x );
1550 38250284 : x = shl( x, sx );
1551 38250284 : x = shr( x, 1 );
1552 38250284 : move16();
1553 38250284 : *s = sub( 1, sx );
1554 :
1555 38250284 : sy = norm_s( y );
1556 38250284 : y = shl( y, sy );
1557 38250284 : move16();
1558 38250284 : *s = add( *s, sy );
1559 :
1560 38250284 : z = div_s( x, y );
1561 :
1562 38250284 : if ( sign != 0 )
1563 : {
1564 1080826 : z = negate( z );
1565 : }
1566 :
1567 :
1568 38250284 : return z;
1569 : }
1570 :
1571 :
1572 51029 : void set_val_Word16(
1573 : Word16 X[],
1574 : const Word16 val,
1575 : Word16 n )
1576 : {
1577 : Word16 i;
1578 :
1579 :
1580 46987020 : FOR( i = 0; i < n; i++ )
1581 : {
1582 46935991 : X[i] = val;
1583 46935991 : move16();
1584 : }
1585 :
1586 :
1587 51029 : return;
1588 : }
1589 :
1590 969942 : void set_val_Word32(
1591 : Word32 X[],
1592 : const Word32 val,
1593 : Word16 n )
1594 : {
1595 : Word16 i;
1596 :
1597 :
1598 288831951 : FOR( i = 0; i < n; i++ )
1599 : {
1600 287862009 : X[i] = val;
1601 287862009 : move32();
1602 : }
1603 :
1604 :
1605 969942 : return;
1606 : }
1607 :
1608 20475180 : Word16 mult0(
1609 : Word16 x,
1610 : Word16 y )
1611 : {
1612 20475180 : return extract_l( L_mult0( x, y ) );
1613 : }
1614 :
1615 2172585 : void copyWord8( const Word8 *src, Word8 *dst, const Word32 n )
1616 : {
1617 : Word32 i;
1618 :
1619 :
1620 23179822 : FOR( i = 0; i < n; i++ )
1621 : {
1622 21007237 : dst[i] = src[i];
1623 21007237 : move16();
1624 : }
1625 2172585 : }
1626 :
1627 :
1628 0 : void set_zero_Word8( Word8 X[], Word32 n )
1629 : {
1630 : Word32 i;
1631 :
1632 :
1633 0 : FOR( i = 0; i < n; i++ )
1634 : {
1635 0 : X[i] = 0;
1636 0 : move16();
1637 : }
1638 0 : }
1639 :
1640 :
1641 0 : Word32 L_mult0_3216( Word32 x, Word16 y )
1642 : {
1643 : UWord16 mpy_low16;
1644 : Word32 mpy_high32;
1645 :
1646 :
1647 0 : Mpy_32_16_ss( x, y, &mpy_high32, &mpy_low16 );
1648 :
1649 0 : mpy_high32 = L_add( L_shl( mpy_high32, 15 ), L_lshr( L_deposit_h( mpy_low16 ), 17 ) );
1650 :
1651 :
1652 0 : return mpy_high32;
1653 : }
1654 :
1655 397970 : Word16 BASOP_util_norm_l_dim2_cplx( const Word32 *const *re, /*!< Real part of 32 Bit input */
1656 : const Word32 *const *im, /*!< Imag part if 32 Bit input */
1657 : Word16 startBand, /*!< start band of cplx data */
1658 : Word16 stopBand, /*!< stop band of cplx data */
1659 : Word16 startSlot, /*!< start slot of cplx data */
1660 : Word16 stopSlot /*!< stop slot of cplx data */
1661 : )
1662 : {
1663 : Word16 col;
1664 : Word16 band;
1665 : Word16 maxShift;
1666 : Word32 maxVal;
1667 :
1668 :
1669 397970 : maxVal = L_deposit_l( 1 );
1670 :
1671 6765490 : FOR( col = startSlot; col < stopSlot; col++ )
1672 : {
1673 364890720 : FOR( band = startBand; band < stopBand; band++ )
1674 : {
1675 358523200 : maxVal = L_max( maxVal, L_abs( re[col][band] ) );
1676 358523200 : maxVal = L_max( maxVal, L_abs( im[col][band] ) );
1677 : }
1678 : }
1679 397970 : maxShift = norm_l( maxVal );
1680 :
1681 :
1682 397970 : return ( maxShift );
1683 : }
1684 :
1685 763948 : Word16 BASOP_util_norm_s_bands2shift( Word16 x )
1686 : {
1687 : Word16 shift;
1688 :
1689 763948 : shift = sub( WORD16_BITS - 1, norm_s( negate( x ) ) );
1690 :
1691 763948 : return ( shift );
1692 : }
1693 :
1694 : #define SINETAB SineTable512_fx
1695 : #define LD 9
1696 :
1697 : /*
1698 : * Calculates coarse lookup values for sine/cosine and residual angle.
1699 : * \param x angle in radians with exponent = 2 or as radix 2 with exponent = 0.
1700 : * \param scale shall always be 2
1701 : * \param sine pointer to where the sine lookup value is stored into
1702 : * \param cosine pointer to where the cosine lookup value is stored into
1703 : * \param flag_radix2 flag indicating radix 2 angle if non-zero.
1704 : */
1705 135696170 : static Word16 fixp_sin_cos_residual_16(
1706 : Word16 x,
1707 : const Word16 scale,
1708 : Word16 *sine,
1709 : Word16 *cosine,
1710 : Word8 flag_radix2 )
1711 : {
1712 : Word16 residual;
1713 : Word16 s;
1714 : Word16 ssign;
1715 : Word16 csign;
1716 135696170 : Word16 tmp, cl = 0, sl = 0;
1717 135696170 : const Word16 shift = 15 - LD - 1 - scale;
1718 :
1719 135696170 : if ( flag_radix2 == 0 )
1720 : {
1721 43124984 : x = mult_r( x, FL2WORD16( 1.0 / EVS_PI ) );
1722 : }
1723 135696170 : s = shr( x, shift );
1724 :
1725 135696170 : residual = s_and( x, ( 1 << shift ) - 1 );
1726 : /* We assume "2+scale" is a constant */
1727 135696170 : residual = shl( residual, 2 + scale );
1728 135696170 : residual = mult_r( residual, FL2WORD16( EVS_PI / 4.0 ) );
1729 :
1730 : /* Sine sign symmetry */
1731 135696170 : ssign = s_and( s, ( 1 << LD ) << 1 );
1732 :
1733 : /* Cosine sign symmetry */
1734 135696170 : csign = s_and( add( s, ( 1 << LD ) ), ( 1 << LD ) << 1 );
1735 :
1736 : /* Modulo EVS_PI */
1737 135696170 : s = s_and( s, ( 2 << LD ) - 1 );
1738 :
1739 : /* EVS_PI/2 symmetry */
1740 135696170 : s = s_min( s, sub( 2 << LD, s ) );
1741 :
1742 : {
1743 135696170 : tmp = s_min( sub( 1 << LD, s ), s );
1744 135696170 : s = sub( tmp, s );
1745 :
1746 135696170 : if ( !s )
1747 : {
1748 81384245 : move16();
1749 81384245 : sl = SINETAB[tmp].v.im;
1750 : }
1751 135696170 : if ( !s )
1752 : {
1753 81384245 : move16();
1754 81384245 : cl = SINETAB[tmp].v.re;
1755 : }
1756 135696170 : if ( s )
1757 : {
1758 54311925 : move16();
1759 54311925 : sl = SINETAB[tmp].v.re;
1760 : }
1761 135696170 : if ( s )
1762 : {
1763 54311925 : move16();
1764 54311925 : cl = SINETAB[tmp].v.im;
1765 : }
1766 :
1767 135696170 : if ( ssign )
1768 : {
1769 41167685 : sl = negate( sl );
1770 : }
1771 135696170 : if ( csign )
1772 : {
1773 38620821 : cl = negate( cl );
1774 : }
1775 :
1776 135696170 : move16();
1777 135696170 : move16();
1778 135696170 : *sine = sl;
1779 135696170 : *cosine = cl;
1780 : }
1781 :
1782 135696170 : return residual;
1783 : }
1784 :
1785 43124984 : Word16 getCosWord16( Word16 theta )
1786 : {
1787 : Word16 result, residual, sine, cosine;
1788 :
1789 43124984 : residual = fixp_sin_cos_residual_16( theta, 2, &sine, &cosine, 0 );
1790 : /* This negation prevents the subsequent addition from overflow */
1791 : /* The negation cannot overflow, sine is in range [0x0..0x7FFF] */
1792 43124984 : sine = negate( sine );
1793 43124984 : result = mac_r( L_mult0( sine, residual ), cosine, 16384 );
1794 :
1795 :
1796 43124984 : return result;
1797 : }
1798 :
1799 18228054 : Word16 getSinWord16( Word16 theta )
1800 : {
1801 : Word16 sine;
1802 18228054 : Word32 theta_new = L_sub( EVS_PI_BY_2_FX, theta );
1803 : Word16 l_theta;
1804 18228054 : IF( GT_32( theta_new, EVS_PI_FX ) )
1805 : {
1806 982003 : l_theta = extract_l( L_sub( L_sub( theta_new, EVS_PI_FX ), EVS_PI_FX ) );
1807 : }
1808 17246051 : ELSE IF( LT_32( theta_new, -EVS_PI_FX ) )
1809 : {
1810 0 : l_theta = extract_l( L_add( L_add( theta_new, EVS_PI_FX ), EVS_PI_FX ) );
1811 : }
1812 : ELSE
1813 : {
1814 17246051 : l_theta = extract_l( theta_new );
1815 : }
1816 18228054 : sine = getCosWord16( l_theta );
1817 18228054 : IF( EQ_16( sine, ONE_IN_Q14 ) )
1818 : {
1819 21393 : sine = MAX_16;
1820 : }
1821 : ELSE
1822 : {
1823 18206661 : sine = shl( sine, 1 );
1824 : }
1825 18228054 : return sine;
1826 : }
1827 :
1828 92571186 : Word16 getCosWord16R2(
1829 : Word16 theta )
1830 : {
1831 : Word16 result, residual, sine, cosine;
1832 :
1833 92571186 : residual = fixp_sin_cos_residual_16( theta, 1, &sine, &cosine, 1 );
1834 : /* This negation prevents the subsequent addition from overflow */
1835 : /* The negation cannot overflow, sine is in range [0x0..0x7FFF] */
1836 : BASOP_SATURATE_WARNING_OFF
1837 92571186 : sine = negate( sine );
1838 : /* Saturation has been included based on the recommendation from the PC group */
1839 92571186 : result = msu_r_sat( L_mult( sine, residual ), cosine, -32768 );
1840 : BASOP_SATURATE_WARNING_ON
1841 :
1842 92571186 : return result;
1843 : }
1844 :
1845 :
1846 39086052 : Word16 getSineWord16R2( Word16 theta )
1847 : {
1848 39086052 : IF( EQ_16( theta, (Word16) 0 ) )
1849 : {
1850 10517849 : return 0;
1851 : }
1852 28568203 : ELSE IF( LT_16( theta, (Word16) -24576 ) )
1853 : {
1854 120594 : theta = add( add( theta, (Word16) 32767 ), (Word16) 1 );
1855 : }
1856 28568203 : return getCosWord16R2( sub( 8192, theta ) );
1857 : }
1858 :
1859 : /*
1860 : * Calculate Integer Square Root of 'val'. This is the equivalent of (int)sqrt(val);
1861 : * The return value will be truncated to the lowest integer (throwing away the fractionnal part.
1862 : *
1863 : * There are many ways to do this. The approach here is to use a simple function to get a
1864 : * 1st estimate of (int)sqrt(val) and then correct this estimate if it is too low or too high.
1865 : *
1866 : * Using Word16, the range of 'val' is limited to roughly 2^30.
1867 : *
1868 : * Complexity: Worst=31Clks, Best=27Clks
1869 : */
1870 5250704 : Word16 getSqrtWord32( Word32 val )
1871 : {
1872 : Word32 L_temp, L_temp2;
1873 : Word16 temp, temp2;
1874 : Word16 exp, exp2;
1875 :
1876 : /* Calc Approximation */
1877 5250704 : exp2 = norm_l( val );
1878 5250704 : L_temp2 = L_shl( val, exp2 );
1879 5250704 : exp = sub( 31 - 32, exp2 );
1880 5250704 : L_temp = Isqrt_lc( L_temp2, &exp ); /* 12 clks */
1881 :
1882 5250704 : temp = round_fx_sat( L_temp );
1883 5250704 : L_temp = Mpy_32_16_1( L_temp2, temp ); /* 2 clks */
1884 :
1885 5250704 : L_temp = L_shl( L_temp, sub( exp, exp2 ) );
1886 :
1887 : /* The Approximation Error Range is -1..+7, so Too Low by 1 or Up to Too High by 7 */
1888 5250704 : temp = round_fx( L_temp );
1889 :
1890 : /* Too High? */
1891 5250704 : if ( L_msu0( val, temp, temp ) < 0 )
1892 : {
1893 : /* Reduce by 2 */
1894 4215893 : temp = sub( temp, 2 );
1895 : }
1896 : /* Too High? */
1897 5250704 : if ( L_msu0( val, temp, temp ) < 0 )
1898 : {
1899 : /* Reduce by 2 */
1900 478609 : temp = sub( temp, 2 );
1901 : }
1902 : /* Too High? */
1903 5250704 : if ( L_msu0( val, temp, temp ) < 0 )
1904 : {
1905 : /* Reduce by 2 */
1906 586 : temp = sub( temp, 2 );
1907 : }
1908 : /* Too High? */
1909 5250704 : if ( L_msu0( val, temp, temp ) < 0 )
1910 : {
1911 : /* Reduce by 1 */
1912 6 : temp = sub( temp, 1 );
1913 : }
1914 :
1915 : /* Try +1 */
1916 5250704 : temp2 = add( temp, 1 );
1917 : /* It fits? */
1918 5250704 : if ( L_msu0( val, temp2, temp2 ) >= 0 )
1919 : {
1920 : /* Yes */
1921 3452881 : temp = temp2;
1922 3452881 : move16();
1923 : }
1924 5250704 : return temp;
1925 : }
1926 0 : Word16 findIndexOfMinWord32( Word32 *x, const Word16 len )
1927 : {
1928 : Word16 i, indx;
1929 :
1930 :
1931 0 : indx = 0;
1932 0 : move16();
1933 0 : FOR( i = 1; i < len; i++ )
1934 : {
1935 0 : if ( LT_32( x[i], x[indx] ) )
1936 : {
1937 0 : indx = i;
1938 0 : move16();
1939 : }
1940 : }
1941 :
1942 :
1943 0 : return indx;
1944 : }
1945 :
1946 34888 : Word16 findIndexOfMinWord64( Word64 *x, const Word16 len )
1947 : {
1948 : Word16 i, indx;
1949 :
1950 :
1951 34888 : indx = 0;
1952 34888 : move16();
1953 279104 : FOR( i = 1; i < len; i++ )
1954 : {
1955 244216 : if ( LT_64( x[i], x[indx] ) )
1956 : {
1957 43905 : indx = i;
1958 43905 : move16();
1959 : }
1960 : }
1961 :
1962 :
1963 34888 : return indx;
1964 : }
1965 :
1966 :
1967 2575169806 : Word16 imult1616( Word16 x, Word16 y )
1968 : {
1969 2575169806 : assert( (int) x * (int) y < 32768 && (int) x * (int) y >= -32768 );
1970 2575169806 : return i_mult( x, y );
1971 : }
1972 :
1973 1860921482 : Word32 imult3216( Word32 x, Word16 y )
1974 : {
1975 : Word32 mh;
1976 : UWord16 ml;
1977 :
1978 1860921482 : Mpy_32_16_ss( x, y, &mh, &ml );
1979 :
1980 1860921482 : mh = L_shl( mh, 15 );
1981 1860921482 : ml = (UWord16) lshr( (Word16) ml, 1 );
1982 :
1983 1860921482 : return L_or( mh, L_deposit_l( ml ) );
1984 : }
1985 :
1986 :
1987 0 : Word16 idiv1616U_IVAS(
1988 : Word16 x,
1989 : Word16 y )
1990 : {
1991 : Word16 s;
1992 : Word16 tmp;
1993 :
1994 : /* make y > x */
1995 0 : s = add( sub( norm_s( y ), norm_s( x ) ), 1 );
1996 0 : s = s_max( s, 0 );
1997 :
1998 : BASOP_SATURATE_WARNING_OFF
1999 0 : y = shl( y, s );
2000 : BASOP_SATURATE_WARNING_ON
2001 :
2002 : /* divide and shift */
2003 0 : tmp = div_s( x, y );
2004 0 : y = shr( tmp, sub( 15, s ) );
2005 :
2006 0 : return y;
2007 : }
2008 32858362 : Word16 idiv1616U( Word16 x, Word16 y )
2009 : {
2010 : Word16 sx, sy;
2011 :
2012 : /* make y > x to meet the requirements for div_s parameters */
2013 32858362 : sx = norm_s( x );
2014 32858362 : sy = norm_s( y );
2015 32858362 : x = shl( x, sx );
2016 32858362 : y = shl( y, sy );
2017 :
2018 32858362 : if ( x >= y )
2019 : {
2020 9672739 : x = shr( x, 1 );
2021 9672739 : sx = sub( sx, 1 );
2022 : }
2023 :
2024 : /* divide and shift */
2025 32858362 : y = shr( div_s( x, y ), sub( 15, sub( sy, sx ) ) );
2026 :
2027 32858362 : return y;
2028 : }
2029 :
2030 :
2031 98974881 : Word16 idiv1616( Word16 x, Word16 y )
2032 : {
2033 : Word16 s, num, den, sign;
2034 :
2035 :
2036 98974881 : num = abs_s( x );
2037 98974881 : den = abs_s( y );
2038 :
2039 98974881 : sign = s_and( s_xor( x, y ), -32768 /* 0x8000 */ );
2040 :
2041 : /* make num > den */
2042 98974881 : s = add( sub( norm_s( den ), norm_s( num ) ), 1 );
2043 98974881 : s = s_max( s, 0 );
2044 :
2045 98974881 : den = shl( den, s );
2046 :
2047 : /* divide and shift */
2048 98974881 : y = shr( div_s( num, den ), sub( 15, s ) );
2049 :
2050 98974881 : if ( 0 != sign )
2051 : {
2052 150342 : y = negate( y );
2053 : }
2054 :
2055 :
2056 98974881 : return y;
2057 : }
2058 :
2059 3 : Word16 idiv1616_1( Word16 x, Word16 y )
2060 : {
2061 3 : IF( L_mult0( x, y ) < 0 )
2062 : {
2063 0 : return negate( idiv1616( abs_s( x ), abs_s( y ) ) );
2064 : }
2065 3 : ELSE IF( L_mult0( x, y ) > 0 )
2066 : {
2067 3 : return idiv1616( x, y );
2068 : }
2069 : ELSE
2070 : {
2071 0 : return 0;
2072 : }
2073 : }
2074 :
2075 : #ifndef FIX_ISSUE_1965_REPLACE_CARRY_OVERFLOW
2076 : Word32 norm_llQ31( /* o : normalized result Q31 */
2077 : Word32 L_c, /* i : upper bits of accu Q-1 */
2078 : Word32 L_sum, /* i : lower bits of accu, unsigned Q31 */
2079 : Word16 *exp /* o : exponent of result in [-32,31] Q0 */
2080 : )
2081 : {
2082 : #ifdef FIX_ISSUE_1965_REPLACE_CARRY_OVERFLOW
2083 : move32(); /* compensate missing instrumentation */
2084 : move32(); /* compensate missing instrumentation */
2085 : return w_norm_llQ31( ( (Word64) L_c << 32 ) | (UWord32) L_sum, exp );
2086 : #else
2087 : Word16 i;
2088 : Word32 L_tmp;
2089 : #ifdef BASOP_NOGLOB_DECLARE_LOCAL
2090 : #ifndef ISSUE_1836_replace_overflow_libcom
2091 : Flag Overflow = 0;
2092 : #endif
2093 : Flag Carry = 0;
2094 : #endif /* BASOP_NOGLOB */
2095 :
2096 : /* Move MSBit of L_sum into L_c */
2097 : #ifdef ISSUE_1836_replace_overflow_libcom
2098 : Carry = 0;
2099 : L_tmp = L_add_c( L_sum, L_sum, &Carry ); /* L_tmp = L_sum << 1 */
2100 : L_c = L_add_c( L_c, L_c, &Carry );
2101 : #else
2102 : L_tmp = L_add_co( L_sum, L_sum, &Carry, &Overflow ); /* L_tmp = L_sum << 1 */
2103 : L_c = L_add_co( L_c, L_c, &Carry, &Overflow );
2104 : #endif
2105 : L_add( 0, 0 );
2106 : test();
2107 : IF( ( L_c != (Word32) 0L ) && ( L_c != (Word32) 0xFFFFFFFFL ) )
2108 : {
2109 : i = norm_l( L_c );
2110 : L_c = L_shl( L_c, i );
2111 : i = sub( 31, i ); /* positive exponent */
2112 : L_sum = L_lshr( L_tmp, 1 ); /* L_sum with MSBit=0 */
2113 : L_sum = L_lshr( L_sum, i );
2114 : L_sum = L_add( L_c, L_sum );
2115 : }
2116 : ELSE
2117 : {
2118 : i = -32;
2119 : move16(); /* default exponent, if total sum=0 */
2120 : IF( L_sum )
2121 : {
2122 : i = norm_l( L_sum );
2123 : L_sum = L_shl( L_sum, i );
2124 : i = negate( i ); /* negative or zero exponent */
2125 : }
2126 : }
2127 : *exp = i;
2128 : move16();
2129 : return L_sum;
2130 : #endif /* #ifdef FIX_ISSUE_1965_REPLACE_CARRY_OVERFLOW */
2131 : }
2132 : #endif /* #ifndef FIX_ISSUE_1965_REPLACE_CARRY_OVERFLOW */
2133 :
2134 :
2135 1489740 : Word32 Dot_product16HQ( /* o : normalized result Q31 */
2136 : const Word32 L_off, /* i : initial sum value Qn */
2137 : const Word16 x[], /* i : x vector Qn */
2138 : const Word16 y[], /* i : y vector Qn */
2139 : const Word16 lg, /* i : vector length, range [0..7FFF] Q0 */
2140 : Word16 *exp /* o : exponent of result in [-32,31] Q0 */
2141 : )
2142 : {
2143 : Word16 i;
2144 : Word32 L_sum;
2145 : Word64 L_sum64;
2146 :
2147 1489740 : L_sum64 = W_deposit32_l( L_off );
2148 :
2149 93024993 : FOR( i = 0; i < lg; i++ )
2150 : {
2151 91535253 : L_sum64 = W_mac_16_16( L_sum64, x[i], y[i] );
2152 : }
2153 :
2154 1489740 : L_sum = w_norm_llQ31( L_sum64, exp );
2155 1489740 : return L_sum;
2156 : }
2157 :
2158 11853836 : Word32 sum_array_norm( /* o : normalized result Q31 */
2159 : const Word32 x[], /* i : x vector Qn */
2160 : const Word16 lg, /* i : vector length, range [0..7FFF] Q0 */
2161 : Word16 *exp /* o : exponent of result in [-32,31] Q0 */
2162 : )
2163 : {
2164 : Word16 i;
2165 11853836 : Word64 W_tmp = 0;
2166 : Word32 L_tmp;
2167 391641658 : FOR( i = 0; i < lg; i++ )
2168 : {
2169 379787822 : W_tmp = W_add( W_tmp, x[i] ); /*Q31*/
2170 : }
2171 11853836 : L_tmp = w_norm_llQ31( W_tmp, exp ); /*Q31 - *exp*/
2172 11853836 : return L_tmp;
2173 : }
2174 :
2175 :
2176 613675 : Word32 Norm32Norm( const Word32 *x, const Word16 headroom, const Word16 length, Word16 *result_e )
2177 : {
2178 : Word32 L_tmp, L_tmp2;
2179 : Word16 i, shift, tmp;
2180 :
2181 613675 : move16();
2182 613675 : shift = headroom;
2183 :
2184 613675 : L_tmp = L_deposit_l( 0 );
2185 :
2186 122466899 : FOR( i = 0; i < length; i++ )
2187 : {
2188 121853224 : L_tmp2 = L_sub( L_tmp, 0x40000000 );
2189 121853224 : if ( L_tmp2 >= 0 )
2190 1027081 : shift = sub( shift, 1 );
2191 121853224 : if ( L_tmp2 >= 0 )
2192 1027081 : L_tmp = L_shr( L_tmp, 2 );
2193 :
2194 121853224 : tmp = round_fx_sat( L_shl_sat( x[i], shift ) );
2195 121853224 : L_tmp = L_mac0_sat( L_tmp, tmp, tmp ); /* exponent = (1-shift*2) , Q(30+shift*2) */
2196 : }
2197 :
2198 613675 : move16();
2199 613675 : *result_e = sub( 1, shl( shift, 1 ) );
2200 :
2201 613675 : return L_tmp;
2202 : }
2203 :
2204 1868893 : Word32 Dot_productSq16HQ( /* o : normalized result Q31 */
2205 : const Word32 L_off, /* i : initial sum value Qn */
2206 : const Word16 x[], /* i : x vector Qn */
2207 : const Word16 lg, /* i : vector length, range [0..7FFF] Q0 */
2208 : Word16 *exp /* o : exponent of result in [-32,31] Q0 */
2209 : )
2210 : {
2211 : Word16 i;
2212 : Word32 L_sum;
2213 : Word64 L_sum64;
2214 :
2215 1868893 : L_sum64 = W_deposit32_l( L_off );
2216 :
2217 138906302 : FOR( i = 0; i < lg; i++ )
2218 : {
2219 137037409 : L_sum64 = W_mac_16_16( L_sum64, x[i], x[i] );
2220 : }
2221 1868893 : L_sum = w_norm_llQ31( L_sum64, exp );
2222 :
2223 1868893 : return L_sum;
2224 : }
2225 :
2226 863600 : Word32 dotp_s_fx( const Word16 *x, const Word16 *y, const Word16 n, Word16 s )
2227 : {
2228 : Word16 i;
2229 : Word16 n2;
2230 : Word32 L_tmp;
2231 : Word32 L_sum;
2232 :
2233 :
2234 863600 : L_sum = 0;
2235 863600 : move32();
2236 :
2237 863600 : n2 = shr( n, 1 );
2238 :
2239 863600 : s = sub( s, 1 );
2240 :
2241 9082489 : FOR( i = 0; i < n2; i++ )
2242 : {
2243 8218889 : L_tmp = L_mult0( x[2 * i], y[2 * i] );
2244 8218889 : L_tmp = L_mac0( L_tmp, x[2 * i + 1], y[2 * i + 1] );
2245 8218889 : L_sum = L_add( L_sum, L_shr( L_tmp, s ) );
2246 : }
2247 :
2248 863600 : IF( s_and( n, 1 ) )
2249 : {
2250 124312 : L_tmp = L_mult0( x[n - 1], y[n - 1] );
2251 124312 : L_sum = L_add( L_sum, L_shr( L_tmp, s ) );
2252 : }
2253 :
2254 :
2255 863600 : return L_sum;
2256 : }
2257 :
2258 :
2259 214967641 : Word32 BASOP_util_Pow2(
2260 : const Word32 exp_m,
2261 : const Word16 exp_e,
2262 : Word16 *result_e )
2263 : {
2264 : static const Word16 pow2Coeff[8] = {
2265 : 22713 /*0.693147180559945309417232121458177 Q15*/, /* ln(2)^1 /1! */
2266 : 7872 /*0.240226506959100712333551263163332 Q15*/, /* ln(2)^2 /2! */
2267 : 1819 /*0.0555041086648215799531422637686218 Q15*/, /* ln(2)^3 /3! */
2268 : 315 /*0.00961812910762847716197907157365887 Q15*/, /* ln(2)^4 /4! */
2269 : 44 /*0.00133335581464284434234122219879962 Q15*/, /* ln(2)^5 /5! */
2270 : 5 /*1.54035303933816099544370973327423e-4 Q15*/, /* ln(2)^6 /6! */
2271 : 0 /*1.52527338040598402800254390120096e-5 Q15*/, /* ln(2)^7 /7! */
2272 : 0 /*1.32154867901443094884037582282884e-6 Q15*/ /* ln(2)^8 /8! */
2273 : };
2274 :
2275 : Word32 frac_part, tmp_frac, result_m;
2276 : Word16 int_part;
2277 :
2278 214967641 : int_part = 0; /* to avoid compilation warnings */
2279 214967641 : frac_part = 0; /* to avoid compilation warnings */
2280 :
2281 214967641 : IF( exp_e > 0 )
2282 : {
2283 : /* "+ 1" compensates L_shr(,1) of the polynomial evaluation at the loop end. */
2284 :
2285 214835721 : int_part = add( 1, extract_l( L_shr( exp_m, sub( 31, exp_e ) ) ) );
2286 214835721 : frac_part = L_lshl( exp_m, exp_e );
2287 214835721 : frac_part = L_and( 0x7FFFFFFF, frac_part );
2288 : }
2289 214967641 : if ( exp_e <= 0 )
2290 131920 : frac_part = L_shl( exp_m, exp_e );
2291 214967641 : if ( exp_e <= 0 )
2292 : {
2293 131920 : int_part = 1;
2294 131920 : move16();
2295 : }
2296 :
2297 : /* Best accuracy is around 0, so try to get there with the fractional part. */
2298 214967641 : IF( ( tmp_frac = L_sub( frac_part, 1073741824l /*0.5 Q31*/ ) ) >= 0 )
2299 : {
2300 105838564 : int_part = add( int_part, 1 );
2301 105838564 : frac_part = L_sub( tmp_frac, 1073741824l /*0.5 Q31*/ );
2302 : }
2303 109129077 : ELSE IF( ( tmp_frac = L_add( frac_part, 1073741824l /*0.5 Q31*/ ) ) < 0 )
2304 : {
2305 0 : int_part = sub( int_part, 1 );
2306 0 : frac_part = L_add( tmp_frac, 1073741824l /*0.5 Q31*/ );
2307 : }
2308 :
2309 : /* Evaluate taylor polynomial which approximates 2^x */
2310 : {
2311 : Word32 p;
2312 : Word16 i;
2313 :
2314 :
2315 : /* First taylor series coefficient a_0 = 1.0, scaled by 0.5 due to L_shr(,1). */
2316 214967641 : result_m = L_add( 1073741824l /*1.0/2.0 Q31*/, L_shr( Mpy_32_16_1( frac_part, pow2Coeff[0] ), 1 ) );
2317 214967641 : p = Mpy_32_32( frac_part, frac_part );
2318 1504773487 : FOR( i = 1; i < 7; i++ )
2319 : {
2320 : /* next taylor series term: a_i * x^i, x=0 */
2321 1289805846 : result_m = L_add( result_m, L_shr( Mpy_32_16_1( p, pow2Coeff[i] ), 1 ) );
2322 1289805846 : p = Mpy_32_32( p, frac_part );
2323 : }
2324 214967641 : result_m = L_add( result_m, L_shr( Mpy_32_16_1( p, pow2Coeff[i] ), 1 ) );
2325 : }
2326 214967641 : *result_e = int_part;
2327 214967641 : move16();
2328 214967641 : return result_m;
2329 : }
2330 :
2331 0 : Word16 findIndexOfMaxWord32( Word32 *x, const Word16 len )
2332 : {
2333 : Word16 i, indx;
2334 :
2335 :
2336 0 : indx = 0;
2337 0 : move16();
2338 0 : FOR( i = 1; i < len; i++ )
2339 : {
2340 0 : if ( GT_32( x[i], x[indx] ) )
2341 : {
2342 0 : indx = i;
2343 0 : move16();
2344 : }
2345 : }
2346 :
2347 :
2348 0 : return indx;
2349 : }
2350 :
2351 55252939 : Word16 getNormReciprocalWord16( Word16 x )
2352 : {
2353 :
2354 55252939 : assert( x < (Word16) ( sizeof( BASOP_util_normReciprocal ) / sizeof( BASOP_util_normReciprocal[0] ) ) );
2355 :
2356 55252939 : return extract_h( BASOP_util_normReciprocal[x] );
2357 : }
2358 64 : Word16 getNormReciprocalWord16Scale( Word16 x, Word16 s )
2359 : {
2360 :
2361 64 : assert( x < (Word16) ( sizeof( BASOP_util_normReciprocal ) / sizeof( BASOP_util_normReciprocal[0] ) ) );
2362 :
2363 64 : return round_fx( L_shl( BASOP_util_normReciprocal[x], s ) );
2364 : }
2365 :
2366 :
2367 : /*! r: result of division x/y, not normalized */
2368 468656153 : Word16 BASOP_Util_Divide3216_Scale(
2369 : Word32 x, /* i : numerator, signed */
2370 : Word16 y, /* i : denominator, signed */
2371 : Word16 *s ) /* o : scaling, 0, if x==0 */
2372 : {
2373 : Word16 z;
2374 : Word16 sx;
2375 : Word16 sy;
2376 : Word16 sign;
2377 :
2378 : /*assert (x > (Word32)0);
2379 : assert (y >= (Word16)0);*/
2380 :
2381 : /* check, if numerator equals zero, return zero then */
2382 468656153 : IF( x == (Word32) 0 )
2383 : {
2384 31867198 : *s = 0;
2385 31867198 : move16();
2386 :
2387 31867198 : return ( (Word16) 0 );
2388 : }
2389 :
2390 436788955 : sign = s_xor( extract_h( x ), y ); /* just to exor the sign bits */
2391 : BASOP_SATURATE_WARNING_OFF
2392 436788955 : x = L_abs( x );
2393 436788955 : y = abs_s( y );
2394 : BASOP_SATURATE_WARNING_ON
2395 436788955 : sx = sub( norm_l( x ), 1 );
2396 436788955 : x = L_shl( x, sx );
2397 436788955 : sy = norm_s( y );
2398 436788955 : y = shl( y, sy );
2399 436788955 : *s = sub( sy, sx );
2400 436788955 : move16();
2401 :
2402 436788955 : z = div_s( round_fx( x ), y );
2403 :
2404 436788955 : if ( sign < 0 ) /* if sign bits differ, negate the result */
2405 : {
2406 21379912 : z = negate( z );
2407 : }
2408 :
2409 436788955 : return z;
2410 : }
2411 :
2412 :
2413 : /*************************************************************************
2414 : *
2415 : * FUNCTION: Log2_norm()
2416 : *
2417 : * PURPOSE: Computes log2(L_x, exp), where L_x is positive and
2418 : * normalized, and exp is the normalisation exponent
2419 : * If L_x is negative or zero, the result is 0.
2420 : *
2421 : * DESCRIPTION:
2422 : * The function Log2(L_x) is approximated by a table and linear
2423 : * interpolation. The following steps are used to compute Log2(L_x)
2424 : *
2425 : * 1- exponent = 30-norm_exponent
2426 : * 2- i = bit25-b31 of L_x; 32<=i<=63 (because of normalization).
2427 : * 3- a = bit10-b24
2428 : * 4- i -=32
2429 : * 5- fraction = table[i]<<16 - (table[i] - table[i+1]) * a * 2
2430 : *
2431 : *************************************************************************/
2432 :
2433 43718631 : Word32 BASOP_Util_fPow(
2434 : Word32 base_m,
2435 : Word16 base_e,
2436 : Word32 exp_m,
2437 : Word16 exp_e,
2438 : Word16 *result_e )
2439 : {
2440 :
2441 : Word16 ans_lg2_e, base_lg2_e;
2442 : Word32 base_lg2_m, ans_lg2_m, result_m;
2443 : Word16 shift;
2444 :
2445 43718631 : test();
2446 43718631 : IF( ( base_m == 0 ) && ( exp_m != 0 ) )
2447 : {
2448 0 : *result_e = 0;
2449 0 : move16();
2450 0 : return 0;
2451 : }
2452 : /* Calc log2 of base */
2453 43718631 : shift = norm_l( base_m );
2454 43718631 : base_m = L_shl( base_m, shift );
2455 43718631 : base_e = sub( base_e, shift );
2456 43718631 : base_lg2_m = BASOP_Util_Log2( base_m );
2457 :
2458 : /* shift: max left shift such that neither base_e or base_lg2_m saturate. */
2459 43718631 : shift = sub( s_min( norm_s( base_e ), WORD16_BITS - 1 - LD_DATA_SCALE ), 1 );
2460 : /* Compensate shift into exponent of result. */
2461 43718631 : base_lg2_e = sub( WORD16_BITS - 1, shift );
2462 43718631 : base_lg2_m = L_add( L_shr( base_lg2_m, sub( WORD16_BITS - 1 - LD_DATA_SCALE, shift ) ), L_deposit_h( shl( base_e, shift ) ) );
2463 :
2464 : /* Prepare exp */
2465 43718631 : shift = norm_l( exp_m );
2466 43718631 : exp_m = L_shl( exp_m, shift );
2467 43718631 : exp_e = sub( exp_e, shift );
2468 :
2469 : /* Calc base pow exp */
2470 43718631 : ans_lg2_m = Mpy_32_32( base_lg2_m, exp_m );
2471 43718631 : ans_lg2_e = add( exp_e, base_lg2_e );
2472 :
2473 : /* Calc antilog */
2474 43718631 : result_m = BASOP_util_Pow2( ans_lg2_m, ans_lg2_e, result_e );
2475 :
2476 43718631 : return result_m;
2477 : }
2478 :
2479 : /*___________________________________________________________________________
2480 : | |
2481 : | Function Name : Dot_product12_offs() |
2482 : | |
2483 : | Compute scalar product of <x[],y[]> using accumulator. |
2484 : | The parameter 'L_off' is added to the accumulation result. |
2485 : | The result is normalized (in Q31) with exponent (0..30). |
2486 : | Notes: |
2487 : | o data in x[],y[] must provide enough headroom for accumulation |
2488 : | o L_off must correspond in format with product of x,y |
2489 : | Example: 0.01f for Q9 x Q9: 0x0000147B in Q19 |
2490 : | means: L_off = FL2WORD32_SCALE(0.01,31-19) |
2491 : |---------------------------------------------------------------------------|
2492 : | Algorithm: |
2493 : | |
2494 : | dot_product = L_off + sum(x[i]*y[i]) i=0..N-1 |
2495 : |___________________________________________________________________________|
2496 : */
2497 336358 : Word32 Dot_product12_offs( /* (o) Q31: normalized result (1 < val <= -1) */
2498 : const Word16 x[], /* (i) 12bits: x vector */
2499 : const Word16 y[], /* (i) 12bits: y vector */
2500 : const Word16 lg, /* (i) : vector length in range [1..256] */
2501 : Word16 *exp, /* (o) : exponent of result (0..+30) */
2502 : Word32 L_off /* (i) initial summation offset / 2 */
2503 : )
2504 : {
2505 : Word16 i, sft;
2506 : Word32 L_sum;
2507 :
2508 : Word64 L_sum64;
2509 :
2510 336358 : L_sum64 = W_deposit32_l( L_off );
2511 17074438 : FOR( i = 0; i < lg; i++ )
2512 : {
2513 16738080 : L_sum64 = W_mac0_16_16( L_sum64, x[i], y[i] );
2514 : }
2515 336358 : L_sum = W_sat_l( L_sum64 );
2516 :
2517 : /* Normalize acc in Q31 */
2518 :
2519 336358 : sft = norm_l( L_sum );
2520 336358 : if ( exp != NULL )
2521 : {
2522 170444 : L_sum = L_shl( L_sum, sft );
2523 : }
2524 :
2525 : /* exponent = 0..30, when L_sum != 0 */
2526 336358 : if ( L_sum != 0 )
2527 : {
2528 336357 : sft = sub( 31, sft );
2529 : }
2530 :
2531 336358 : if ( exp != NULL )
2532 : {
2533 170444 : *exp = sft;
2534 170444 : move16();
2535 : }
2536 :
2537 336358 : return L_sum;
2538 : }
2539 :
2540 15530 : Word32 Dot_product15_offs( /* (o) Q31: normalized result (1 < val <= -1) */
2541 : const Word16 x[], /* (i) 15bits: x vector */
2542 : const Word16 y[], /* (i) 15bits: y vector */
2543 : const Word16 lg, /* (i) : vector length in range [1..256] */
2544 : Word16 *exp, /* (o) : exponent of result (0..+30) */
2545 : Word32 L_off /* (i) initial summation offset */
2546 : )
2547 : {
2548 : Word16 i, sft, fac, ld;
2549 : Word32 L_sum;
2550 :
2551 15530 : ld = sub( 14, norm_s( lg ) );
2552 15530 : fac = shr( -32768, ld );
2553 15530 : L_sum = L_shr( L_off, ld );
2554 :
2555 1009450 : FOR( i = 0; i < lg; i++ )
2556 : {
2557 993920 : L_sum = L_add( L_sum, Mpy_32_16_1( L_msu( 0, y[i], fac ), x[i] ) );
2558 : }
2559 :
2560 : /* Avoid returning 0 */
2561 15530 : if ( L_sum == 0 )
2562 : {
2563 18 : L_sum = L_add( L_sum, 1 );
2564 : }
2565 :
2566 : /* Normalize acc in Q31 */
2567 15530 : sft = norm_l( L_sum );
2568 15530 : L_sum = L_shl( L_sum, sft );
2569 :
2570 : /* exponent = 0..30, when L_sum != 0 */
2571 15530 : if ( L_sum != 0 )
2572 : {
2573 15530 : sft = add( ld, sub( 30, sft ) );
2574 : }
2575 :
2576 15530 : *exp = sft;
2577 15530 : move16();
2578 :
2579 15530 : return L_sum;
2580 : }
2581 :
2582 158758419 : Word16 BASOP_Util_Cmp_Mant32Exp /*!< o: flag: result of comparison */
2583 : /* 0, if a == b */
2584 : /* 1, if a > b */
2585 : /* -1, if a < b */
2586 : ( Word32 a_m, /*!< i: Mantissa of 1st operand a */
2587 : Word16 a_e, /*!< i: Exponent of 1st operand a */
2588 : Word32 b_m, /*!< i: Mantissa of 2nd operand b */
2589 : Word16 b_e ) /*!< i: Exponent of 2nd operand b */
2590 :
2591 : {
2592 : Word32 diff_m;
2593 : Word16 diff_e, shift, result;
2594 :
2595 :
2596 : /*
2597 : This function compares two input parameters, both represented by a 32-bit mantissa and a 16-bit exponent.
2598 : If both values are identical, 0 is returned.
2599 : If a is greater b, 1 is returned.
2600 : If a is less than b, -1 is returned.
2601 : */
2602 :
2603 : /* Check, if both mantissa and exponents are identical, when normalized: return 0 */
2604 158758419 : shift = norm_l( a_m );
2605 158758419 : if ( shift )
2606 79428154 : a_m = L_shl( a_m, shift );
2607 158758419 : if ( shift )
2608 79428154 : a_e = sub( a_e, shift );
2609 :
2610 158758419 : shift = norm_l( b_m );
2611 158758419 : if ( shift )
2612 101776943 : b_m = L_shl( b_m, shift );
2613 158758419 : if ( shift )
2614 101776943 : b_e = sub( b_e, shift );
2615 :
2616 : /* align exponent, if any mantissa is zero */
2617 158758419 : if ( !a_m )
2618 : {
2619 2472916 : a_e = b_e;
2620 2472916 : move16();
2621 : }
2622 158758419 : if ( !b_m )
2623 : {
2624 3347806 : b_e = a_e;
2625 3347806 : move16();
2626 : }
2627 :
2628 : BASOP_SATURATE_WARNING_OFF_EVS
2629 158758419 : diff_m = L_sub_sat( a_m, b_m );
2630 : BASOP_SATURATE_WARNING_ON_EVS
2631 158758419 : diff_e = sub( a_e, b_e );
2632 :
2633 158758419 : test();
2634 158758419 : IF( diff_m == 0 && diff_e == 0 )
2635 : {
2636 891765 : return 0;
2637 : }
2638 :
2639 : /* Check sign, exponent and mantissa to identify, whether a is greater b or not */
2640 157866654 : result = sub( 0, 1 );
2641 :
2642 157866654 : IF( a_m >= 0 )
2643 : {
2644 : /* a is positive */
2645 157233598 : if ( b_m < 0 )
2646 : {
2647 115104 : result = 1;
2648 115104 : move16();
2649 : }
2650 :
2651 157233598 : test();
2652 157233598 : test();
2653 157233598 : test();
2654 157233598 : if ( ( b_m >= 0 ) && ( ( diff_e > 0 ) || ( diff_e == 0 && diff_m > 0 ) ) )
2655 : {
2656 100140295 : result = 1;
2657 100140295 : move16();
2658 : }
2659 : }
2660 : ELSE
2661 : {
2662 : /* a is negative */
2663 633056 : test();
2664 633056 : test();
2665 633056 : test();
2666 633056 : if ( ( b_m < 0 ) && ( ( diff_e < 0 ) || ( diff_e == 0 && diff_m > 0 ) ) )
2667 : {
2668 142823 : result = 1;
2669 142823 : move16();
2670 : }
2671 : }
2672 157866654 : return result;
2673 : }
2674 :
2675 : /*
2676 :
2677 : headroom is introduced into acc
2678 : */
2679 :
2680 :
2681 9545643588 : Word32 BASOP_Util_Add_Mant32Exp /* o : normalized result mantissa */
2682 : ( Word32 a_m, /* i : Mantissa of 1st operand a */
2683 : Word16 a_e, /* i : Exponent of 1st operand a */
2684 : Word32 b_m, /* i : Mantissa of 2nd operand b */
2685 : Word16 b_e, /* i : Exponent of 2nd operand b */
2686 : Word16 *ptr_e ) /* o : exponent of result */
2687 : {
2688 : Word32 L_tmp;
2689 : Word16 shift;
2690 :
2691 : /* Compare exponents: the difference is limited to +/- 30
2692 : The Word32 mantissa of the operand with lower exponent is shifted right by the exponent difference.
2693 : Then, the unshifted mantissa of the operand with the higher exponent is added. The addition result
2694 : is normalized and the result represents the mantissa to return. The returned exponent takes into
2695 : account all shift operations.
2696 : */
2697 :
2698 9545643588 : if ( !a_m )
2699 : {
2700 1788897457 : a_e = b_e;
2701 1788897457 : move16();
2702 : }
2703 :
2704 9545643588 : if ( !b_m )
2705 : {
2706 869894809 : b_e = a_e;
2707 869894809 : move16();
2708 : }
2709 :
2710 9545643588 : shift = sub( a_e, b_e );
2711 9545643588 : shift = s_max( -31, shift );
2712 9545643588 : shift = s_min( 31, shift );
2713 9545643588 : if ( shift < 0 )
2714 : {
2715 : /* exponent of b is greater than exponent of a, shr a_m */
2716 2958037069 : a_m = L_shl( a_m, shift );
2717 : }
2718 9545643588 : if ( shift > 0 )
2719 : {
2720 : /* exponent of a is greater than exponent of b */
2721 3959151109 : b_m = L_shr( b_m, shift );
2722 : }
2723 9545643588 : a_e = add( s_max( a_e, b_e ), 1 );
2724 9545643588 : L_tmp = L_add( L_shr( a_m, 1 ), L_shr( b_m, 1 ) );
2725 9545643588 : shift = norm_l( L_tmp );
2726 9545643588 : if ( shift )
2727 8409494572 : L_tmp = L_shl( L_tmp, shift );
2728 9545643588 : if ( L_tmp == 0 )
2729 : {
2730 627650632 : a_e = 0;
2731 627650632 : move16();
2732 : }
2733 9545643588 : if ( L_tmp != 0 )
2734 8917992956 : a_e = sub( a_e, shift );
2735 9545643588 : *ptr_e = a_e;
2736 :
2737 9545643588 : return ( L_tmp );
2738 : }
2739 :
2740 :
2741 : static const Word16 shift_lc[] = { 9, 10 };
2742 :
2743 0 : Word32 Isqrt_lc1(
2744 : Word32 frac, /* i : Q31: normalized value (1.0 < frac <= 0.5) */
2745 : Word16 *exp /* i/o: exponent (value = frac x 2^exponent) */
2746 : )
2747 : {
2748 : Word16 i, a;
2749 : Word32 L_tmp;
2750 :
2751 0 : IF( frac <= (Word32) 0 )
2752 : {
2753 0 : *exp = 0;
2754 0 : move16();
2755 0 : return 0x7fffffff; /*0x7fffffff*/
2756 : }
2757 :
2758 : /* If exponant odd -> shift right by 10 (otherwise 9) */
2759 0 : L_tmp = L_shr( frac, shift_lc[s_and( *exp, 1 )] );
2760 :
2761 : /* 1) -16384 to shift left and change sign */
2762 : /* 2) 32768 to Add 1 to Exponent like it was divided by 2 */
2763 : /* 3) We let the mac_r add another 0.5 because it imitates */
2764 : /* the behavior of shr on negative number that should */
2765 : /* not be rounded towards negative infinity. */
2766 : /* It replaces: */
2767 : /* *exp = negate(shr(sub(*exp, 1), 1)); move16(); */
2768 0 : *exp = mac_r( 32768, *exp, -16384 );
2769 0 : move16();
2770 :
2771 0 : a = extract_l( L_tmp ); /* Extract b10-b24 */
2772 0 : a = lshr( a, 1 );
2773 :
2774 0 : i = mac_r( L_tmp, -16 * 2 - 1, 16384 ); /* Extract b25-b31 minus 16 */
2775 :
2776 0 : L_tmp = L_msu( L_table_isqrt[i], table_isqrt_diff[i], a ); /* table[i] << 16 - diff*a*2 */
2777 :
2778 0 : return L_tmp;
2779 : }
2780 :
2781 :
2782 12087 : void bufferCopyFx(
2783 : Word16 *src, /*<! Qx pointer to input buffer */
2784 : Word16 *dest, /*<! Qx pointer to output buffer */
2785 : Word16 length, /*<! Q0 length of buffer to copy */
2786 : Word16 Qf_src, /*<! Q0 Q format (frac-bits) of source buffer */
2787 : Word16 Qf_dest, /*<! Q0 Q format (frac-bits) of dest buffer */
2788 : Word16 Q_src, /*<! Q0 exponent of source buffer */
2789 : Word16 Q_dest /*<! Q0 exponent of destination buffer */
2790 : )
2791 : {
2792 : Word16 tmp_16, i;
2793 :
2794 : /*Copy( st->old_exc, exc_buf, st->old_exc_len);*/
2795 12087 : tmp_16 = sub( sub( Qf_src, Qf_dest ), sub( Q_src, Q_dest ) );
2796 12087 : IF( tmp_16 > 0 ) /*if value will be shifted right, do a multiplication with rounding ->preserves more accuracy*/
2797 : {
2798 3352 : tmp_16 = shl( 1, sub( 15, tmp_16 ) );
2799 957456 : FOR( i = 0; i < length; i++ )
2800 : {
2801 954104 : *( dest + i ) = mult_r( *( src + i ), tmp_16 );
2802 954104 : move16();
2803 : }
2804 : }
2805 8735 : ELSE IF( tmp_16 < 0 ) /*leftshift - no accuracy preservation needed*/
2806 : {
2807 1168597 : FOR( i = 0; i < length; i++ )
2808 : {
2809 1165031 : *( dest + i ) = shr_sat( *( src + i ), tmp_16 );
2810 1165031 : move16();
2811 : }
2812 : }
2813 : ELSE /*no shift, simply copy*/
2814 : {
2815 1932536 : FOR( i = 0; i < length; i++ )
2816 : {
2817 1927367 : *( dest + i ) = *( src + i );
2818 1927367 : move16();
2819 : }
2820 : }
2821 12087 : }
2822 :
2823 39654 : Word32 dotWord32_16_Mant32Exp( const Word32 *bufX32, /* i: 32-bit buffer with unknown headroom */
2824 : Word16 bufX32_exp, /* i: exponent of buffer bufX32 */
2825 : const Word16 *bufY16, /* i: 16-bit buffer quite right-aligned */
2826 : Word16 bufY16_exp, /* i: exponent of buffer bufY16 */
2827 : Word16 len, /* i: buffer len to process */
2828 : Word16 *exp ) /* o: result exponent */
2829 : {
2830 : Word32 L_sum;
2831 : Word16 shift, shift1, i;
2832 :
2833 :
2834 39654 : shift = getScaleFactor32( bufX32, len ); /* current available headroom */
2835 39654 : shift = sub( shift, sub( 14, norm_s( len ) ) ); /* reduced required headroom */
2836 39654 : L_sum = 0; /* Clear accu */
2837 39654 : move32();
2838 2478006 : FOR( i = 0; i < len; i++ )
2839 : {
2840 2438352 : L_sum = L_mac0( L_sum, round_fx( L_shl( bufX32[i], shift ) ), bufY16[i] );
2841 : }
2842 39654 : shift1 = norm_l( L_sum );
2843 39654 : L_sum = L_shl( L_sum, shift1 ); /* return value */
2844 :
2845 39654 : shift = sub( add( bufX32_exp, bufY16_exp ), add( shift, shift1 ) );
2846 39654 : shift = add( shift, 1 ); /* compensate for factor of 2 introduced by L_mac0 */
2847 : /* In case of NULL result, we want to have a 0 exponent too */
2848 39654 : if ( L_sum == 0 )
2849 494 : shift = 0;
2850 39654 : *exp = shift;
2851 39654 : move16();
2852 :
2853 :
2854 39654 : return L_sum;
2855 : }
2856 :
2857 11592 : Word16 BASOP_Util_lin2dB( Word32 x, Word16 x_e, Word16 fEnergy )
2858 : {
2859 11592 : assert( x >= 0 );
2860 :
2861 : /* log2 */
2862 11592 : x = L_shr( BASOP_Util_Log2( x ), 25 - 16 ); /* Q16 */
2863 :
2864 : /* add exponent */
2865 11592 : x = L_msu( x, x_e, -32768 /* 0x8000 */ );
2866 :
2867 : /* convert log2 to 20*log10 */
2868 11592 : x = Mpy_32_16_1( x, 24660 /*6.0206f Q12*/ ); /* Q13 */
2869 :
2870 : /* if energy divide by 2 (->10*log10) */
2871 11592 : if ( fEnergy != 0 )
2872 3532 : x = L_shr( x, 1 );
2873 :
2874 : /* return dB as 7Q8 */
2875 11592 : return round_fx( L_shl( x, 8 - 13 + 16 ) ); /* Q8 */
2876 : }
2877 :
2878 : /* --- fixp_atan() ---- */
2879 : #define Q_ATANINP ( 25 ) /* Input in q25, Output in q14 */
2880 : #define Q_ATANOUT ( 14 )
2881 : #define ATI_SF ( ( 32 - 1 ) - Q_ATANINP ) /* 6 */
2882 : #define ATO_SF ( ( 16 - 1 ) - Q_ATANOUT ) /* 1 ] -pi/2 .. pi/2 [ */
2883 : /* --- fixp_atan2() --- */
2884 : #define Q_ATAN2OUT ( 13 )
2885 : #define AT2O_SF ( ( 16 - 1 ) - Q_ATAN2OUT ) /* 2 ] -pi .. pi ] */
2886 :
2887 :
2888 55633440 : Word16 BASOP_util_atan2( /* o: atan2(y,x) [-pi,pi] Q13 */
2889 : Word32 y, /* i: */
2890 : Word32 x, /* i: */
2891 : Word16 e /* i: exponent difference (exp_y - exp_x) */
2892 : )
2893 : {
2894 : Word16 q;
2895 : Word32 at;
2896 55633440 : Word16 ret = -32768 /*-1.0f Q15*/;
2897 : Word16 sf, sfo, stf;
2898 : Word32 L_sign;
2899 :
2900 55633440 : if ( L_or( y, x ) == 0 )
2901 : {
2902 395953 : return 0;
2903 : }
2904 :
2905 55237487 : IF( x == 0 )
2906 : {
2907 472235 : ret = 12868 /*+EVS_PI/2 Q13*/;
2908 472235 : move16();
2909 472235 : if ( y < 0 )
2910 : {
2911 272265 : ret = negate( ret );
2912 : }
2913 :
2914 472235 : return ret;
2915 : }
2916 :
2917 : /* --- division */
2918 54765252 : L_sign = L_and( L_xor( x, y ), (Word32) 0x80000000 );
2919 :
2920 54765252 : q = 32767 /*1.0f Q15*/; /* y/x = neg/zero = -Inf */
2921 54765252 : sf = 0;
2922 : BASOP_SATURATE_WARNING_OFF_EVS
2923 54765252 : q = BASOP_Util_Divide3232_uu_1616_Scale( L_abs( y ), L_abs( x ), &sf );
2924 : BASOP_SATURATE_WARNING_ON_EVS
2925 :
2926 : BASOP_SATURATE_WARNING_OFF_EVS
2927 54765252 : if ( L_sign < 0 )
2928 24011000 : q = negate( q );
2929 : BASOP_SATURATE_WARNING_ON_EVS
2930 :
2931 54765252 : sfo = add( sf, e );
2932 :
2933 : /* --- atan() */
2934 54765252 : IF( GT_16( sfo, ATI_SF ) )
2935 : {
2936 : /* --- could not calc fixp_atan() here bec of input data out of range
2937 : ==> therefore give back boundary values */
2938 :
2939 2089337 : sfo = s_min( sfo, MAXSFTAB );
2940 :
2941 : /*q = FL2WORD16( 0.0f ); move16();*/
2942 :
2943 2089337 : if ( q > 0 )
2944 : {
2945 215222 : move16();
2946 215222 : q = +f_atan_expand_range[sfo - ATI_SF - 1];
2947 : }
2948 2089337 : if ( q < 0 )
2949 : {
2950 309013 : move16();
2951 309013 : q = -f_atan_expand_range[sfo - ATI_SF - 1];
2952 : }
2953 : }
2954 : ELSE
2955 : {
2956 : /* --- calc of fixp_atan() is possible; input data within range
2957 : ==> set q on fixed scale level as desired from fixp_atan() */
2958 52675915 : stf = sub( sfo, ATI_SF );
2959 :
2960 52675915 : at = L_deposit_h( q );
2961 52675915 : if ( stf < 0 )
2962 52270423 : at = L_shl( at, stf );
2963 :
2964 52675915 : q = BASOP_util_atan( at ); /* ATO_SF*/
2965 : }
2966 :
2967 :
2968 : /* --- atan2() */
2969 :
2970 54765252 : ret = shr( q, ( AT2O_SF - ATO_SF ) ); /* now AT2O_SF for atan2 */
2971 54765252 : IF( x < 0 )
2972 : {
2973 9287778 : if ( y >= 0 )
2974 : {
2975 5396881 : ret = add( ret, 25736 /*EVS_PI Q13*/ );
2976 : }
2977 9287778 : if ( y < 0 )
2978 : {
2979 3890897 : ret = sub( ret, 25736 /* EVS_PI Q13*/ );
2980 : }
2981 : }
2982 :
2983 54765252 : return ret;
2984 : }
2985 :
2986 : /* SNR of fixp_atan() = 56 dB*/
2987 : #define ONEBY3P56 0x26800000 /* 1.0/3.56 in q31*/
2988 : #define P281 0x00026000 /* 0.281 in q19*/
2989 : #define ONEP571 0x6487 /* 1.571 in q14*/
2990 :
2991 52750264 : Word16 BASOP_util_atan( /* o: atan(x) [-pi/2;pi/2] 1Q14 */
2992 : Word32 x /* i: input data (-64;64) 6Q25 */
2993 : )
2994 : {
2995 : Word16 sign, result, exp;
2996 : Word16 res_e;
2997 : Word16 tmp, xx;
2998 :
2999 :
3000 52750264 : sign = 0;
3001 52750264 : move16();
3002 52750264 : if ( x < 0 )
3003 : {
3004 23684843 : sign = 1;
3005 23684843 : move16();
3006 : }
3007 52750264 : x = L_abs( x );
3008 :
3009 : /* calc of arctan */
3010 52750264 : IF( LT_32( x, 1509950l /*0.045f/64.0f Q31*/ ) )
3011 : {
3012 6137474 : result = round_fx( L_shl( x, 5 ) ); /*Q14*/
3013 : /*BASOP_util_atan_16(0.0444059968): max error 0.0000567511, mean 0.000017, abs mean 0.000017*/
3014 : }
3015 : ELSE
3016 46612790 : IF( LT_32( x, ( L_shl( 1, Q_ATANINP ) - 8482560l /*0.00395 Q31*/ ) ) )
3017 : {
3018 22783652 : xx = round_fx( L_shl( x, 6 ) );
3019 22783652 : tmp = mult_r( xx, xx ); /* q15 * q15 - (16-1) = q15*/
3020 22783652 : tmp = mult_r( tmp, 0x1340 ); /* 15 * (ONEBY3P56) q14 - (16-1) = q14*/
3021 22783652 : tmp = add( tmp, 0x4000 ); /*L_shl(1,14) = 524288*/ /* q14 + q14 = q14 */
3022 22783652 : res_e = Q_ATANOUT - 15 + 14 - 16 + 1;
3023 22783652 : move16();
3024 22783652 : if ( GT_16( xx, tmp ) )
3025 : {
3026 4944264 : res_e = add( res_e, 1 );
3027 : }
3028 22783652 : if ( GT_16( xx, tmp ) )
3029 : {
3030 4944264 : xx = shr( xx, 1 );
3031 : }
3032 22783652 : result = div_s( xx, tmp );
3033 22783652 : result = msu_r( 0, result, shl( -32768, res_e ) );
3034 : /*BASOP_util_atan_16(0.7471138239): max error 0.0020029545, mean 0.000715, abs mean 0.000715*/
3035 : }
3036 23829138 : ELSE IF( LT_32( x, 42949673l /*1.28/64.0 Q31*/ ) )
3037 : {
3038 : Word16 delta_fix;
3039 8466216 : Word32 PI_BY_4 = 1686629684l /*3.1415926/4.0 Q31*/ / 2; /* pi/4 in q30 */
3040 :
3041 8466216 : delta_fix = round_fx( L_shl( L_sub( x, 33554432l /*1.0/64.0 Q31*/ ), 5 ) ); /* q30 */
3042 8466216 : result = round_fx( L_sub( L_add( PI_BY_4, L_msu( 0, delta_fix, -16384 ) ), ( L_mult0( delta_fix, delta_fix ) ) ) );
3043 : /* BASOP_Util_fPow(0.7472000122): max error 0.0020237688, mean 0.000026, abs mean 0.000520 */
3044 : }
3045 : ELSE
3046 : {
3047 15362922 : exp = sub( norm_l( x ), 1 );
3048 15362922 : xx = round_fx( L_shl( x, exp ) );
3049 : /* q25+exp * q25+exp - (16-1) = q19+2*exp*/
3050 15362922 : tmp = mac_r( L_shl( P281, shl( exp, 1 ) ), xx, xx ); /* q19+2*exp + q19+2*exp = q19+2*exp*/
3051 15362922 : res_e = norm_s( tmp );
3052 15362922 : result = div_s( xx, shl( tmp, res_e ) );
3053 15362922 : result = shl( result, add( add( Q_ATANOUT - Q_ATANINP /*-exp*/ + 19 /*+2*exp*/ - 16 + 1, res_e ), exp ) );
3054 15362922 : result = sub( ONEP571, result ); /* q14 + q14 = q14*/
3055 : /*BASOP_Util_fPow(1.2799999714): max error 0.0020168927, mean 0.000066, abs mean 0.000072*/
3056 : }
3057 :
3058 52750264 : if ( sign )
3059 : {
3060 23684843 : result = negate( result );
3061 : }
3062 :
3063 52750264 : return ( result );
3064 : }
3065 :
3066 : /* compare two positive normalized 16 bit mantissa/exponent values */
3067 : /* return value: positive if first value greater, negative if second value greater, zero if equal */
3068 1730773 : Word16 compMantExp16Unorm( Word16 m1, Word16 e1, Word16 m2, Word16 e2 )
3069 : {
3070 : Word16 tmp;
3071 :
3072 1730773 : assert( ( m1 >= 0x4000 ) && ( m2 >= 0x4000 ) ); /* comparisons below work only for normalized mantissas */
3073 :
3074 1730773 : tmp = sub( e1, e2 );
3075 1730773 : if ( tmp == 0 )
3076 958527 : tmp = sub( m1, m2 );
3077 :
3078 1730773 : return tmp;
3079 : }
3080 :
3081 3619376666 : cmplx CL_scale_t( cmplx x, Word16 y )
3082 : {
3083 : cmplx result;
3084 3619376666 : result.re = Mpy_32_16_1( x.re, y );
3085 3619376666 : result.im = Mpy_32_16_1( x.im, y );
3086 : #ifdef WMOPS
3087 : multiCounter[currCounter].Mpy_32_16_1--;
3088 : multiCounter[currCounter].Mpy_32_16_1--;
3089 : multiCounter[currCounter].CL_scale++;
3090 : #endif
3091 3619376666 : return ( result );
3092 : }
3093 :
3094 21361658 : cmplx CL_dscale_t( cmplx x, Word16 y1, Word16 y2 )
3095 : {
3096 : cmplx result;
3097 21361658 : result.re = Mpy_32_16_1( x.re, y1 );
3098 21361658 : result.im = Mpy_32_16_1( x.im, y2 );
3099 : #ifdef WMOPS
3100 : multiCounter[currCounter].Mpy_32_16_1--;
3101 : multiCounter[currCounter].Mpy_32_16_1--;
3102 : multiCounter[currCounter].CL_dscale++;
3103 : #endif /* #ifdef WMOPS */
3104 21361658 : return ( result );
3105 : }
3106 :
3107 485496355 : cmplx CL_mult_32x16( cmplx input, cmplx_s coeff )
3108 : {
3109 : cmplx result;
3110 485496355 : result.re = L_sub( Mpy_32_16_1( input.re, coeff.re ), Mpy_32_16_1( input.im, coeff.im ) );
3111 485496355 : result.im = L_add( Mpy_32_16_1( input.re, coeff.im ), Mpy_32_16_1( input.im, coeff.re ) );
3112 : #ifdef WMOPS
3113 : multiCounter[currCounter].CL_multr_32x16++;
3114 : multiCounter[currCounter].Mpy_32_16_1--;
3115 : multiCounter[currCounter].Mpy_32_16_1--;
3116 : multiCounter[currCounter].Mpy_32_16_1--;
3117 : multiCounter[currCounter].Mpy_32_16_1--;
3118 : multiCounter[currCounter].L_sub--;
3119 : multiCounter[currCounter].L_add--;
3120 : #endif
3121 485496355 : return result;
3122 : }
|
