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.443 Nov 04, 2021. Version 12.14.0 / 13.10.0 / 14.6.0 / 15.4.0 / 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 14843897 : 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 14843897 : assert( x >= 0 );
67 :
68 14843897 : if ( x == 0 )
69 : {
70 :
71 57047 : return ( (Word32) MIN_32 );
72 : }
73 :
74 : /* normalize input, calculate integer part */
75 14786850 : exp_e = norm_l( x );
76 14786850 : x = L_shl( x, exp_e );
77 14786850 : exp = L_deposit_l( exp_e );
78 :
79 : /* calculate (1-normalized_input) */
80 14786850 : 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 14786850 : accuRes = L_mult( nIn, ldCoeff[0] ); /* nIn^i * coeff[0] */
86 14786850 : accuSqr = mult( nIn, nIn ); /* nIn^2, nIn^3 .... */
87 :
88 : /* iteration 2 */
89 14786850 : accuRes = L_mac( accuRes, accuSqr, ldCoeff[1] ); /* nIn^i * coeff[1] */
90 14786850 : accuSqr = mult( accuSqr, nIn ); /* nIn^2, nIn^3 .... */
91 :
92 : /* iteration 3 */
93 14786850 : accuRes = L_mac( accuRes, accuSqr, ldCoeff[2] ); /* nIn^i * coeff[2] */
94 14786850 : accuSqr = mult( accuSqr, nIn ); /* nIn^2, nIn^3 .... */
95 :
96 : /* iteration 4 */
97 14786850 : accuRes = L_mac( accuRes, accuSqr, ldCoeff[3] ); /* nIn^i * coeff[3] */
98 14786850 : accuSqr = mult( accuSqr, nIn ); /* nIn^2, nIn^3 .... */
99 :
100 : /* iteration 5 */
101 14786850 : accuRes = L_mac( accuRes, accuSqr, ldCoeff[4] ); /* nIn^i * coeff[4] */
102 14786850 : accuSqr = mult( accuSqr, nIn ); /* nIn^2, nIn^3 .... */
103 :
104 : /* iteration 6 */
105 14786850 : accuRes = L_mac( accuRes, accuSqr, ldCoeff[5] ); /* nIn^i * coeff[5] */
106 14786850 : accuSqr = mult( accuSqr, nIn ); /* nIn^2, nIn^3 .... */
107 :
108 : /* iteration 7, no need to calculate accuSqr any more */
109 14786850 : 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 14786850 : accuRes = L_mac0( L_shr( accuRes, 1 ), extract_h( accuRes ), 14506 );
113 :
114 14786850 : accuRes = L_shr( accuRes, LD_DATA_SCALE - 1 ); /* fractional part/LD_DATA_SCALE */
115 14786850 : exp = L_shl( exp, ( 31 - LD_DATA_SCALE ) ); /* integer part/LD_DATA_SCALE */
116 14786850 : accuRes = L_sub( accuRes, exp ); /* result = integer part + fractional part */
117 :
118 :
119 14786850 : return ( accuRes );
120 : }
121 :
122 7638074 : 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 7638074 : if ( x < -1040187392l /*-31.0/64.0 Q31*/ )
136 : {
137 :
138 0 : return 0;
139 : }
140 7638074 : test();
141 7638074 : if ( ( GE_32( x, 1040187392l /*31.0/64.0 Q31*/ ) ) || ( x == 0 ) )
142 : {
143 :
144 0 : return 0x7FFFFFFF;
145 : }
146 :
147 7638074 : frac = extract_l( L_and( x, 0x3FF ) );
148 :
149 7638074 : index3 = L_and( L_shr( x, 10 ), 0x1F );
150 7638074 : index2 = L_and( L_shr( x, 15 ), 0x1F );
151 7638074 : index1 = L_and( L_shr( x, 20 ), 0x1F );
152 :
153 7638074 : exp = extract_l( L_shr( x, 25 ) );
154 7638074 : if ( x > 0 )
155 : {
156 144 : exp = sub( 31, exp );
157 : }
158 7638074 : if ( x < 0 )
159 : {
160 7637930 : exp = negate( exp );
161 : }
162 :
163 7638074 : lookup3f = L_add( exp2x_tab_long[index3], L_shr( Mpy_32_16_1( 0x0016302F, frac ), 1 ) );
164 7638074 : lookup12 = Mpy_32_32( exp2_tab_long[index1], exp2w_tab_long[index2] );
165 7638074 : lookup = Mpy_32_32( lookup12, lookup3f );
166 :
167 7638074 : retVal = L_shr( lookup, sub( exp, 3 ) );
168 :
169 :
170 7638074 : return retVal;
171 : }
172 :
173 0 : Word32 BASOP_Util_Log10( Word32 x, Word16 e )
174 : {
175 0 : test();
176 0 : IF( e >= 0 && LE_16( e, 31 ) )
177 : {
178 0 : IF( EQ_32( x, L_shl_sat( 1, sub( 31, e ) ) ) )
179 : {
180 0 : return 0;
181 : }
182 : }
183 0 : Word32 res = BASOP_Util_Log2( x );
184 0 : 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 0 : return res;
187 : }
188 :
189 0 : Word32 BASOP_Util_Loge( Word32 x, Word16 e )
190 : {
191 0 : Word32 res = BASOP_Util_Log2( x );
192 0 : 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 0 : return res;
194 : }
195 :
196 35709661 : 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 35709661 : shift = sub( a_e, b_e );
227 35709661 : if ( shift >= 0 )
228 19798692 : shift = s_min( 15, shift );
229 :
230 35709661 : if ( shift < 0 )
231 15910969 : shift = s_max( -15, shift );
232 35709661 : a_e = s_max( a_e, b_e );
233 35709661 : L_hm = L_deposit_l( a_m ); /* mantissa belonging to higher exponent */
234 35709661 : L_lm = L_deposit_l( a_m ); /* mantissa belonging to lower exponent */
235 35709661 : if ( shift >= 0 )
236 19798692 : L_lm = L_deposit_l( b_m );
237 35709661 : if ( shift < 0 )
238 15910969 : L_hm = L_deposit_l( b_m );
239 :
240 35709661 : if ( shift > 0 )
241 17077164 : shift = negate( shift );
242 :
243 35709661 : L_hm = L_shr( L_hm, shift ); /* shift left due to negative shift parameter */
244 35709661 : a_e = add( a_e, shift );
245 35709661 : L_hm = L_add( L_hm, L_lm );
246 35709661 : shift = norm_l( L_hm );
247 35709661 : L_hm = L_shl( L_hm, shift );
248 35709661 : *ptrSum_m = extract_h( L_hm );
249 35709661 : move16();
250 :
251 35709661 : a_e = sub( a_e, shift );
252 35709661 : if ( L_hm )
253 35473459 : a_e = add( a_e, 16 );
254 35709661 : return ( a_e );
255 : }
256 :
257 :
258 769212 : 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 769212 : assert( b_m != 0 );
278 :
279 : /* normalize b */
280 769212 : preShift = norm_s( b_m );
281 769212 : m = shl( b_m, preShift );
282 :
283 : /* make b positive */
284 : BASOP_SATURATE_WARNING_OFF_EVS;
285 769212 : 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 769212 : index = mac_r( -32768 - ( 32 << 16 ), m, 1 << 6 );
292 :
293 : /* get fractional part for interpolation (lower 9 bits) */
294 769212 : frac = shl( s_and( m, 0x1FF ), 1 ); /* Q10 */
295 :
296 : /* interpolate 1/b */
297 769212 : m = msu_r( InvTable[index], InvDiffTable[index], frac );
298 :
299 : /* restore sign */
300 769212 : if ( b_m < 0 )
301 0 : m = negate( m );
302 :
303 : /* multiply with a */
304 769212 : m32 = L_mult( a_m, m );
305 :
306 : /* normalize result */
307 769212 : postShift = norm_l( m32 );
308 : #ifdef ISSUE_1836_replace_overflow_libcom
309 769212 : 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 769212 : *ptrResult_e = sub( add( add( a_e, sub( 1, b_e ) ), preShift ), postShift );
316 769212 : move16();
317 :
318 769212 : *ptrResult_m = m;
319 769212 : move16();
320 769212 : }
321 :
322 :
323 : /* local function for Sqrt16 */
324 79395255 : 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 79395255 : assert( ( m >= 0x4000 ) || ( m == 0 ) );
334 :
335 : /* get table index (upper 6 bits minus 32) */
336 : /* index = (m >> 9) - 32; */
337 79395255 : index = mac_r( -32768 - ( 32 << 16 ), m, 1 << 6 );
338 :
339 : /* get fractional part for interpolation (lower 9 bits) */
340 79395255 : frac = s_and( m, 0x1FF ); /* Q9 */
341 :
342 : /* interpolate */
343 79395255 : if ( m != 0 )
344 : {
345 : #ifdef ISSUE_1836_replace_overflow_libcom
346 78682259 : 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 79395255 : if ( s_and( e, 1 ) != 0 )
356 35374831 : m = mult_r( m, 0x5a82 );
357 :
358 79395255 : return m;
359 : }
360 :
361 :
362 : /* local function for Sqrt32 and Sqrt32norm */
363 721830498 : 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 721830498 : assert( ( m >= 0x40000000 ) || ( m == 0 ) );
374 :
375 : #ifdef ISSUE_1836_replace_overflow_libcom
376 721830498 : 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 721830498 : index = mac_r( -32768 - ( 32 << 16 ), m16, 1 << 6 );
384 :
385 : /* get fractional part for interpolation (lower 9 bits) */
386 721830498 : frac = s_and( m16, 0x1FF ); /* Q9 */
387 :
388 : /* interpolate */
389 721830498 : if ( m != 0 )
390 : {
391 : BASOP_SATURATE_WARNING_OFF_EVS;
392 675775008 : m = L_mac_sat( SqrtTable[index], SqrtDiffTable[index], frac );
393 : BASOP_SATURATE_WARNING_ON_EVS;
394 : }
395 :
396 : /* handle odd exponents */
397 721830498 : if ( s_and( e, 1 ) != 0 )
398 341635077 : m = Mpy_32_16_1( m, 0x5a82 );
399 :
400 721830498 : return m;
401 : }
402 :
403 : /* local function for ISqrt16 and ISqrt16norm */
404 41165672 : static Word16 ISqrt16_common( Word16 m,
405 : Word16 e )
406 : {
407 : Word16 index, frac;
408 :
409 41165672 : assert( m >= 0x4000 );
410 :
411 : /* get table index (upper 6 bits minus 32) */
412 : /* index = (m >> 9) - 32; */
413 41165672 : index = mac_r( -32768 - ( 32 << 16 ), m, 1 << 6 );
414 :
415 : /* get fractional part for interpolation (lower 9 bits) */
416 41165672 : frac = s_and( m, 0x1FF ); /* Q9 */
417 :
418 : /* interpolate */
419 41165672 : m = msu_r( ISqrtTable[index], ISqrtDiffTable[index], frac );
420 :
421 : /* handle even exponents */
422 41165672 : if ( s_and( e, 1 ) == 0 )
423 : {
424 18085922 : m = mult_r( m, 0x5a82 );
425 : }
426 :
427 41165672 : return m;
428 : }
429 :
430 : /* local function for ISqrt32 and ISqrt32norm */
431 68289564 : 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 68289564 : assert( m >= 0x40000000 );
442 : #ifdef ISSUE_1836_replace_overflow_libcom
443 68289564 : 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 68289564 : index = mac_r( -32768 - ( 32 << 16 ), m16, 1 << 6 );
456 :
457 : /* get fractional part for interpolation (lower 9 bits) */
458 68289564 : frac = s_and( m16, 0x1FF ); /* Q9 */
459 :
460 : /* interpolate */
461 68289564 : m = L_msu( ISqrtTable[index], ISqrtDiffTable[index], frac );
462 :
463 : /* handle even exponents */
464 68289564 : if ( s_and( e, 1 ) == 0 )
465 : {
466 38944386 : m = Mpy_32_16_1( m, 0x5a82 );
467 : }
468 :
469 68289564 : return m;
470 : }
471 :
472 78822325 : Word16 Sqrt16( /*!< output mantissa */
473 : Word16 mantissa, /*!< input mantissa */
474 : Word16 *exponent /*!< pointer to exponent */
475 : )
476 : {
477 : Word16 preShift, e;
478 :
479 78822325 : assert( mantissa >= 0 );
480 :
481 : /* normalize */
482 78822325 : preShift = norm_s( mantissa );
483 :
484 78822325 : e = sub( *exponent, preShift );
485 78822325 : mantissa = shl( mantissa, preShift );
486 :
487 : /* calc mantissa */
488 78822325 : mantissa = Sqrt16_common( mantissa, e );
489 :
490 : /* e = (e + 1) >> 1 */
491 78822325 : *exponent = mult_r( e, 1 << 14 );
492 78822325 : move16();
493 :
494 78822325 : return mantissa;
495 : }
496 :
497 :
498 572930 : Word16 Sqrt16norm( /*!< output mantissa */
499 : Word16 mantissa, /*!< normalized input mantissa */
500 : Word16 *exponent /*!< pointer to exponent */
501 : )
502 : {
503 :
504 572930 : assert( ( mantissa >= 0x4000 ) || ( mantissa == 0 ) );
505 :
506 : /* calc mantissa */
507 572930 : mantissa = Sqrt16_common( mantissa, *exponent );
508 :
509 : /* e = (e + 1) >> 1 */
510 572930 : *exponent = mult_r( *exponent, 1 << 14 );
511 572930 : move16();
512 :
513 572930 : return mantissa;
514 : }
515 :
516 :
517 41165672 : Word16 ISqrt16( /*!< output mantissa */
518 : Word16 mantissa, /*!< input mantissa */
519 : Word16 *exponent /*!< pointer to exponent */
520 : )
521 : {
522 : Word16 preShift, e;
523 :
524 41165672 : assert( mantissa > 0 );
525 :
526 : /* normalize */
527 41165672 : preShift = norm_s( mantissa );
528 :
529 41165672 : e = sub( *exponent, preShift );
530 41165672 : mantissa = shl( mantissa, preShift );
531 :
532 : /* calc mantissa */
533 41165672 : mantissa = ISqrt16_common( mantissa, e );
534 :
535 : /* e = (2 - e) >> 1 */
536 41165672 : *exponent = msu_r( 1L << 15, e, 1 << 14 );
537 41165672 : move16();
538 :
539 41165672 : return mantissa;
540 : }
541 :
542 721830498 : Word32 Sqrt32( /*!< output mantissa */
543 : Word32 mantissa, /*!< input mantissa */
544 : Word16 *exponent /*!< pointer to exponent */
545 : )
546 : {
547 : Word16 preShift, e;
548 :
549 721830498 : assert( mantissa >= 0 );
550 :
551 : /* normalize */
552 721830498 : preShift = norm_l( mantissa );
553 :
554 721830498 : e = sub( *exponent, preShift );
555 721830498 : mantissa = L_shl( mantissa, preShift );
556 :
557 : /* calc mantissa */
558 721830498 : mantissa = Sqrt32_common( mantissa, e );
559 :
560 : /* e = (e + 1) >> 1 */
561 721830498 : *exponent = mult_r( e, 1 << 14 );
562 721830498 : move16();
563 :
564 721830498 : return mantissa;
565 : }
566 :
567 :
568 68275576 : Word32 ISqrt32( /*!< output mantissa */
569 : Word32 mantissa, /*!< input mantissa */
570 : Word16 *exponent /*!< pointer to exponent */
571 : )
572 : {
573 : Word16 preShift, e;
574 :
575 68275576 : assert( mantissa > 0 );
576 :
577 : /* normalize */
578 68275576 : preShift = norm_l( mantissa );
579 :
580 68275576 : e = sub( *exponent, preShift );
581 68275576 : mantissa = L_shl( mantissa, preShift );
582 :
583 : /* calc mantissa */
584 68275576 : mantissa = ISqrt32_common( mantissa, e );
585 :
586 : /* e = (2 - e) >> 1 */
587 68275576 : *exponent = msu_r( 1L << 15, e, 1 << 14 );
588 68275576 : move16();
589 :
590 68275576 : return mantissa;
591 : }
592 :
593 13988 : Word32 ISqrt32norm( /*!< output mantissa */
594 : Word32 mantissa, /*!< normalized input mantissa */
595 : Word16 *exponent /*!< pointer to exponent */
596 : )
597 : {
598 :
599 13988 : assert( mantissa >= 0x40000000 );
600 :
601 : /* calc mantissa */
602 13988 : mantissa = ISqrt32_common( mantissa, *exponent );
603 :
604 : /* e = (2 - e) >> 1 */
605 13988 : *exponent = msu_r( 1L << 15, *exponent, 1 << 14 );
606 13988 : move16();
607 :
608 13988 : return mantissa;
609 : }
610 :
611 5793138 : 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 5793138 : assert( mantissa != 0 );
622 :
623 : /* absolute */
624 : BASOP_SATURATE_WARNING_OFF_EVS;
625 5793138 : m = abs_s( mantissa );
626 : BASOP_SATURATE_WARNING_ON_EVS;
627 :
628 : /* normalize */
629 5793138 : preShift = norm_s( m );
630 :
631 5793138 : e = sub( *exponent, preShift );
632 5793138 : m = shl( m, preShift );
633 :
634 : /* get table index (upper 6 bits minus 32) */
635 : /* index = (m >> 9) - 32; */
636 5793138 : index = mac_r( -32768 - ( 32 << 16 ), m, 1 << 6 );
637 :
638 : /* get fractional part for interpolation (lower 9 bits) */
639 5793138 : frac = shl( s_and( m, 0x1FF ), 1 ); /* Q10 */
640 :
641 : /* interpolate */
642 5793138 : m = msu_r( InvTable[index], InvDiffTable[index], frac );
643 :
644 : /* restore sign */
645 5793138 : if ( mantissa < 0 )
646 2514064 : m = negate( m );
647 :
648 : /* e = 1 - e */
649 5793138 : *exponent = sub( 1, e );
650 5793138 : move16();
651 :
652 5793138 : return m;
653 : }
654 :
655 :
656 15943040 : 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 15943040 : assert( mantissa > 0 );
670 :
671 : /* normalize */
672 15943040 : preShift = norm_s( mantissa );
673 :
674 15943040 : exponent = sub( exponent, preShift );
675 15943040 : mantissa = shl( mantissa, preShift );
676 :
677 : /* get table index (upper 6 bits minus 32) */
678 : /* index = (m >> 9) - 32; */
679 15943040 : index = mac_r( -32768 - ( 32 << 16 ), mantissa, 1 << 6 );
680 :
681 : /* get fractional part for interpolation (lower 9 bits) */
682 15943040 : frac = s_and( mantissa, 0x1FF ); /* Q9 */
683 :
684 : /* interpolate */
685 : BASOP_SATURATE_WARNING_OFF_EVS;
686 15943040 : m = mac_r_sat( SqrtTable[index], SqrtDiffTable[index], frac );
687 15943040 : mi = msu_r_sat( ISqrtTable[index], ISqrtDiffTable[index], frac );
688 : BASOP_SATURATE_WARNING_ON_EVS;
689 :
690 : /* handle even/odd exponents */
691 15943040 : e_odd = s_and( exponent, 1 );
692 15943040 : if ( e_odd != 0 )
693 8054553 : m = mult_r( m, 0x5a82 );
694 15943040 : if ( e_odd == 0 )
695 7888487 : mi = mult_r( mi, 0x5a82 );
696 :
697 : /* e = (e + 1) >> 1 */
698 15943040 : *sqrt_exp = mult_r( exponent, 1 << 14 );
699 15943040 : move16();
700 :
701 : /* e = (2 - e) >> 1 */
702 15943040 : *isqrt_exp = msu_r( 1L << 15, exponent, 1 << 14 );
703 15943040 : move16();
704 :
705 : /* Write result */
706 15943040 : *sqrt_mant = m;
707 15943040 : move16();
708 15943040 : *isqrt_mant = mi;
709 15943040 : move16();
710 15943040 : }
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 4717598 : 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 4717598 : x_max = 0;
732 4717598 : move16();
733 4717598 : x_min = 0;
734 4717598 : move16();
735 6003255352 : FOR( i = 0; i < len_x; i++ )
736 : {
737 5998537754 : if ( x[i] >= 0 )
738 3733098652 : x_max = s_max( x_max, x[i] );
739 5998537754 : if ( x[i] < 0 )
740 2265439102 : x_min = s_min( x_min, x[i] );
741 : }
742 :
743 4717598 : i_max = 0x10;
744 4717598 : move16();
745 4717598 : i_min = 0x10;
746 4717598 : move16();
747 :
748 4717598 : if ( x_max != 0 )
749 4274692 : i_max = norm_s( x_max );
750 :
751 4717598 : if ( x_min != 0 )
752 4172082 : i_min = norm_s( x_min );
753 :
754 4717598 : i = s_and( s_min( i_max, i_min ), 0xF );
755 :
756 :
757 4717598 : 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 98260718 : 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 98260718 : x_max = 0;
780 98260718 : move32();
781 98260718 : x_min = 0;
782 98260718 : move32();
783 9417327412 : FOR( i = 0; i < len_x; i++ )
784 : {
785 9319066694 : if ( x[i] >= 0 )
786 6274105184 : x_max = L_max( x_max, x[i] );
787 9319066694 : if ( x[i] < 0 )
788 3044961510 : x_min = L_min( x_min, x[i] );
789 : }
790 :
791 98260718 : i_max = 0x20;
792 98260718 : move16();
793 98260718 : i_min = 0x20;
794 98260718 : move16();
795 :
796 98260718 : if ( x_max != 0 )
797 90744700 : i_max = norm_l( x_max );
798 :
799 98260718 : if ( x_min != 0 )
800 82877590 : i_min = norm_l( x_min );
801 :
802 98260718 : i = s_and( s_min( i_max, i_min ), 0x1F );
803 :
804 :
805 98260718 : return i;
806 : }
807 :
808 19553 : 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 19553 : x_max = 0;
818 19553 : move32();
819 19553 : x_min = 0;
820 19553 : move32();
821 129305633 : FOR( i = 0; i < len_x; i++ )
822 : {
823 129286080 : if ( x[i] >= 0 )
824 107871501 : x_max = L_max( x_max, x[i] );
825 129286080 : if ( x[i] < 0 )
826 21414579 : x_min = L_min( x_min, x[i] );
827 : }
828 :
829 19553 : i_max = 0x20;
830 19553 : move16();
831 19553 : i_min = 0x20;
832 19553 : move16();
833 :
834 19553 : if ( x_max != 0 )
835 19544 : i_max = norm_l( x_max );
836 :
837 19553 : if ( x_min != 0 )
838 19553 : i_min = norm_l( x_min );
839 :
840 19553 : i_max = s_and( s_min( i_max, i_min ), 0x1F );
841 :
842 :
843 19553 : return i_max;
844 : }
845 :
846 0 : Word16 normalize16( Word16 mantissa, Word16 *pexponent )
847 : {
848 : Word16 tmp;
849 :
850 0 : tmp = norm_s( mantissa );
851 0 : mantissa = shl( mantissa, tmp );
852 0 : move16();
853 0 : *pexponent = sub( *pexponent, tmp );
854 :
855 :
856 0 : return mantissa;
857 : }
858 162172 : Word16 divide3216( Word32 x, Word16 y )
859 : {
860 : Word16 z;
861 :
862 :
863 162172 : z = 0;
864 162172 : move16();
865 162172 : if ( 0 == y )
866 : {
867 0 : return 0x7fff;
868 : }
869 :
870 162172 : IF( x != 0 )
871 : {
872 : Word16 den, sign;
873 : Word32 num;
874 162172 : num = L_abs( x );
875 162172 : den = abs_s( y );
876 :
877 162172 : sign = s_and( s_xor( extract_h( x ), y ), -32768 /* 0x8000 */ );
878 :
879 162172 : z = div_l( num, den );
880 162172 : if ( 0 != sign )
881 : {
882 0 : z = negate( z );
883 : }
884 : }
885 :
886 :
887 162172 : return z;
888 : }
889 :
890 89327 : Word16 divide1616( Word16 x, Word16 y )
891 : {
892 : Word16 z, num, den, sign;
893 :
894 :
895 89327 : num = abs_s( x );
896 89327 : den = abs_s( y );
897 :
898 89327 : sign = s_and( s_xor( x, y ), -32768 /* 0x8000 */ );
899 :
900 89327 : move16();
901 89327 : z = 0x7fff;
902 89327 : if ( LT_16( num, den ) )
903 89327 : z = div_s( num, den );
904 :
905 89327 : if ( 0 != sign )
906 : {
907 11076 : z = negate( z );
908 : }
909 :
910 :
911 89327 : return z;
912 : }
913 :
914 4107606 : 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 4107606 : sign = L_and( L_xor( L_num, L_denom ), (Word32) 0x80000000 );
926 :
927 4107606 : L_num = L_abs( L_num );
928 4107606 : L_denom = L_abs( L_denom );
929 :
930 : /* limit the range of denominator to Word16 */
931 4107606 : z = s_min( norm_l( L_num ), norm_l( L_denom ) );
932 4107606 : L_num = L_shl( L_num, z );
933 4107606 : 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 4107606 : 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 4107606 : if ( 0 != sign )
942 : {
943 838521 : z = negate( z );
944 : }
945 :
946 :
947 4107606 : return z;
948 : }
949 :
950 44119028 : 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 44119028 : assert( x >= 0 );
960 44119028 : assert( y > 0 );
961 :
962 44119028 : if ( x == 0 )
963 : {
964 4421809 : *s = 0;
965 4421809 : move16();
966 :
967 :
968 4421809 : return ( 0 );
969 : }
970 :
971 39697219 : sx = norm_l( x );
972 39697219 : sy = norm_l( y );
973 :
974 39697219 : x16 = extract_h( L_shl( x, sx ) );
975 39697219 : y16 = extract_h( L_shl( y, sy ) );
976 :
977 39697219 : if ( GT_16( x16, y16 ) )
978 : {
979 13552949 : sx = sub( sx, 1 );
980 : }
981 :
982 39697219 : if ( LT_16( y16, x16 ) )
983 : {
984 13552949 : x16 = mult_r( x16, 0x4000 );
985 : }
986 :
987 :
988 39697219 : z = div_s( x16, y16 );
989 39697219 : move16();
990 39697219 : *s = sub( sy, sx );
991 :
992 :
993 39697219 : return ( z );
994 : }
995 :
996 27371989 : Word32 div_w( Word32 L_num, Word32 L_den )
997 : {
998 27371989 : Word32 L_var_out = 0;
999 : Word16 iteration;
1000 27371989 : move32();
1001 :
1002 :
1003 27371989 : IF( L_den == 0 )
1004 : {
1005 : /* printf("Division by 0 in div_l, Fatal error in "); printStack(); */
1006 0 : return ( 0 );
1007 : }
1008 :
1009 27371989 : test();
1010 27371989 : 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 27371989 : W_num = W_deposit32_h( L_num );
1017 27371989 : W_den = W_deposit32_h( L_den );
1018 :
1019 27371989 : IF( W_sub( W_num, W_den ) >= 0 )
1020 : {
1021 0 : return MAX_32;
1022 : }
1023 : ELSE
1024 : {
1025 27371989 : W_num = W_shr( W_num, 1 );
1026 27371989 : W_den = W_shr( W_den, 1 );
1027 :
1028 875903648 : FOR( iteration = 0; iteration < 31; iteration++ )
1029 : {
1030 848531659 : L_var_out = L_shl( L_var_out, 1 );
1031 848531659 : W_num = W_shl( W_num, 1 );
1032 :
1033 848531659 : IF( W_sub( W_num, W_den ) >= 0 )
1034 : {
1035 503175649 : W_num = W_sub( W_num, W_den );
1036 503175649 : L_var_out = L_add( L_var_out, 1 );
1037 : }
1038 : }
1039 :
1040 27371989 : return L_var_out;
1041 : }
1042 : }
1043 :
1044 27836840 : 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 27836840 : assert( y != (Word32) 0 );
1053 :
1054 27836840 : sign = 0;
1055 27836840 : move16();
1056 :
1057 27836840 : IF( x < 0 )
1058 : {
1059 0 : x = L_negate( x );
1060 0 : sign = L_xor( sign, 1 );
1061 : }
1062 :
1063 27836840 : IF( y < 0 )
1064 : {
1065 0 : y = L_negate( y );
1066 0 : sign = L_xor( sign, 1 );
1067 : }
1068 :
1069 27836840 : IF( x == (Word32) 0 )
1070 : {
1071 464851 : *s = 0;
1072 464851 : return ( (Word32) 0 );
1073 : }
1074 :
1075 27371989 : sx = norm_l( x );
1076 27371989 : x = L_shl( x, sx );
1077 27371989 : x = L_shr( x, 1 );
1078 27371989 : move16();
1079 27371989 : *s = sub( 1, sx );
1080 :
1081 27371989 : sy = norm_l( y );
1082 27371989 : y = L_shl( y, sy );
1083 27371989 : move16();
1084 27371989 : *s = add( *s, sy );
1085 :
1086 27371989 : z = div_w( x, y );
1087 :
1088 27371989 : if ( sign != 0 )
1089 : {
1090 0 : z = L_negate( z );
1091 : }
1092 :
1093 27371989 : 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 967043625 : static Word32 L_dmult( Word32 L_var1, Word32 L_var2 )
1371 : {
1372 967043625 : Word64 L64_var1 = W_mult0_32_32( L_var1, L_var2 );
1373 967043625 : L64_var1 = W_shr( L64_var1, 30 );
1374 967043625 : 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 322347875 : Word32 div_w_newton( Word32 num, Word32 den )
1392 : {
1393 : Word32 x0, x1, x2, x3, diff, result;
1394 :
1395 322347875 : x0 = div_w_newton_lookup[sub( extract_l( L_shr( den, 22 ) ), 256 )];
1396 322347875 : move32();
1397 :
1398 322347875 : diff = L_sub( 0x40000000, Mpy_32_32( den, x0 ) );
1399 :
1400 322347875 : x1 = L_add( x0, L_dmult( x0, diff ) );
1401 322347875 : diff = L_sub( 0x40000000, Mpy_32_32( den, x1 ) );
1402 :
1403 322347875 : x2 = L_add( x1, L_dmult( x1, diff ) );
1404 322347875 : diff = L_sub( 0x40000000, Mpy_32_32( den, x2 ) );
1405 :
1406 322347875 : x3 = L_add( x2, L_dmult( x2, diff ) );
1407 :
1408 322347875 : result = Mpy_32_32( num, x3 );
1409 :
1410 322347875 : 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 339967180 : 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 339967180 : assert( y != (Word32) 0 );
1435 :
1436 : /* Early exit, if numerator is zero */
1437 339967180 : IF( x == (Word32) 0 )
1438 : {
1439 1900335 : *s = 0;
1440 1900335 : return ( (Word32) 0 );
1441 : }
1442 338066845 : IF( EQ_32( y, (Word32) 0x80000000 ) )
1443 : {
1444 : /* Division by -1.0: same as negation of numerator */
1445 : /* Return normalized negated numerator */
1446 101 : sx = norm_l( x );
1447 101 : x = L_shl( x, sx );
1448 101 : *s = negate( sx );
1449 101 : return L_negate( x );
1450 : }
1451 338066744 : sign = y;
1452 338066744 : move32();
1453 338066744 : if ( y < 0 )
1454 : {
1455 12849638 : y = L_negate( y );
1456 : }
1457 :
1458 : /* Normalize numerator */
1459 338066744 : sx = norm_l( x );
1460 338066744 : x = L_shl( x, sx );
1461 :
1462 : /* Normalize denominator */
1463 338066744 : sy = norm_l( y );
1464 338066744 : y = L_shl( y, sy );
1465 :
1466 : /* Store exponent: + 1 for div_w_newton computing 0.5*num/den */
1467 338066744 : *s = sub( add( sy, 1 ), sx );
1468 338066744 : move16();
1469 :
1470 : /* Special treatment for den=0x40000000 */
1471 : /* Result is known: z=2*num */
1472 338066744 : IF( EQ_32( y, 0x40000000 ) )
1473 : {
1474 15718869 : if ( sign < 0 )
1475 : {
1476 0 : x = L_negate( x );
1477 : }
1478 15718869 : return x;
1479 : }
1480 :
1481 : /* Invoke division applying Newton/Raphson-Algorithm */
1482 322347875 : z = div_w_newton( x, y );
1483 :
1484 322347875 : if ( sign < 0 )
1485 : {
1486 12849638 : z = L_negate( z );
1487 : }
1488 :
1489 322347875 : return z;
1490 : }
1491 :
1492 266906317 : Word16 BASOP_Util_Divide3232_Scale( Word32 x, Word32 y, Word16 *s )
1493 : {
1494 : Word16 z;
1495 : Word16 sy;
1496 :
1497 :
1498 266906317 : sy = norm_l( y );
1499 266906317 : if ( sy > 0 )
1500 : {
1501 236253534 : sy = sub( sy, 1 );
1502 : }
1503 266906317 : y = L_shl( y, sy );
1504 :
1505 266906317 : z = BASOP_Util_Divide3216_Scale( x, extract_h( y ), s );
1506 266906317 : move16();
1507 266906317 : *s = add( *s, sy );
1508 :
1509 :
1510 266906317 : return ( z );
1511 : }
1512 :
1513 :
1514 14782063 : 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 14782063 : assert( y != (Word16) 0 );
1524 :
1525 14782063 : sign = 0;
1526 14782063 : move16();
1527 :
1528 14782063 : IF( x < 0 )
1529 : {
1530 665930 : x = negate( x );
1531 665930 : sign = s_xor( sign, 1 );
1532 : }
1533 :
1534 14782063 : IF( y < 0 )
1535 : {
1536 0 : y = negate( y );
1537 0 : sign = s_xor( sign, 1 );
1538 : }
1539 :
1540 14782063 : IF( x == (Word16) 0 )
1541 : {
1542 790348 : move16();
1543 790348 : *s = 0;
1544 :
1545 :
1546 790348 : return ( (Word16) 0 );
1547 : }
1548 :
1549 13991715 : sx = norm_s( x );
1550 13991715 : x = shl( x, sx );
1551 13991715 : x = shr( x, 1 );
1552 13991715 : move16();
1553 13991715 : *s = sub( 1, sx );
1554 :
1555 13991715 : sy = norm_s( y );
1556 13991715 : y = shl( y, sy );
1557 13991715 : move16();
1558 13991715 : *s = add( *s, sy );
1559 :
1560 13991715 : z = div_s( x, y );
1561 :
1562 13991715 : if ( sign != 0 )
1563 : {
1564 665930 : z = negate( z );
1565 : }
1566 :
1567 :
1568 13991715 : return z;
1569 : }
1570 :
1571 :
1572 330 : void set_val_Word16(
1573 : Word16 X[],
1574 : const Word16 val,
1575 : Word16 n )
1576 : {
1577 : Word16 i;
1578 :
1579 :
1580 317130 : FOR( i = 0; i < n; i++ )
1581 : {
1582 316800 : X[i] = val;
1583 316800 : move16();
1584 : }
1585 :
1586 :
1587 330 : return;
1588 : }
1589 :
1590 0 : void set_val_Word32(
1591 : Word32 X[],
1592 : const Word32 val,
1593 : Word16 n )
1594 : {
1595 : Word16 i;
1596 :
1597 :
1598 0 : FOR( i = 0; i < n; i++ )
1599 : {
1600 0 : X[i] = val;
1601 0 : move32();
1602 : }
1603 :
1604 :
1605 0 : return;
1606 : }
1607 :
1608 595965 : Word16 mult0(
1609 : Word16 x,
1610 : Word16 y )
1611 : {
1612 595965 : return extract_l( L_mult0( x, y ) );
1613 : }
1614 :
1615 0 : void copyWord8( const Word8 *src, Word8 *dst, const Word32 n )
1616 : {
1617 : Word32 i;
1618 :
1619 :
1620 0 : FOR( i = 0; i < n; i++ )
1621 : {
1622 0 : dst[i] = src[i];
1623 0 : move16();
1624 : }
1625 0 : }
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 0 : 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 0 : maxVal = L_deposit_l( 1 );
1670 :
1671 0 : FOR( col = startSlot; col < stopSlot; col++ )
1672 : {
1673 0 : FOR( band = startBand; band < stopBand; band++ )
1674 : {
1675 0 : maxVal = L_max( maxVal, L_abs( re[col][band] ) );
1676 0 : maxVal = L_max( maxVal, L_abs( im[col][band] ) );
1677 : }
1678 : }
1679 0 : maxShift = norm_l( maxVal );
1680 :
1681 :
1682 0 : return ( maxShift );
1683 : }
1684 :
1685 0 : Word16 BASOP_util_norm_s_bands2shift( Word16 x )
1686 : {
1687 : Word16 shift;
1688 :
1689 0 : shift = sub( WORD16_BITS - 1, norm_s( negate( x ) ) );
1690 :
1691 0 : 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 65411950 : 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 65411950 : Word16 tmp, cl = 0, sl = 0;
1717 65411950 : const Word16 shift = 15 - LD - 1 - scale;
1718 :
1719 65411950 : if ( flag_radix2 == 0 )
1720 : {
1721 3533653 : x = mult_r( x, FL2WORD16( 1.0 / EVS_PI ) );
1722 : }
1723 65411950 : s = shr( x, shift );
1724 :
1725 65411950 : residual = s_and( x, ( 1 << shift ) - 1 );
1726 : /* We assume "2+scale" is a constant */
1727 65411950 : residual = shl( residual, 2 + scale );
1728 65411950 : residual = mult_r( residual, FL2WORD16( EVS_PI / 4.0 ) );
1729 :
1730 : /* Sine sign symmetry */
1731 65411950 : ssign = s_and( s, ( 1 << LD ) << 1 );
1732 :
1733 : /* Cosine sign symmetry */
1734 65411950 : csign = s_and( add( s, ( 1 << LD ) ), ( 1 << LD ) << 1 );
1735 :
1736 : /* Modulo EVS_PI */
1737 65411950 : s = s_and( s, ( 2 << LD ) - 1 );
1738 :
1739 : /* EVS_PI/2 symmetry */
1740 65411950 : s = s_min( s, sub( 2 << LD, s ) );
1741 :
1742 : {
1743 65411950 : tmp = s_min( sub( 1 << LD, s ), s );
1744 65411950 : s = sub( tmp, s );
1745 :
1746 65411950 : if ( !s )
1747 : {
1748 41798885 : move16();
1749 41798885 : sl = SINETAB[tmp].v.im;
1750 : }
1751 65411950 : if ( !s )
1752 : {
1753 41798885 : move16();
1754 41798885 : cl = SINETAB[tmp].v.re;
1755 : }
1756 65411950 : if ( s )
1757 : {
1758 23613065 : move16();
1759 23613065 : sl = SINETAB[tmp].v.re;
1760 : }
1761 65411950 : if ( s )
1762 : {
1763 23613065 : move16();
1764 23613065 : cl = SINETAB[tmp].v.im;
1765 : }
1766 :
1767 65411950 : if ( ssign )
1768 : {
1769 19637238 : sl = negate( sl );
1770 : }
1771 65411950 : if ( csign )
1772 : {
1773 17879400 : cl = negate( cl );
1774 : }
1775 :
1776 65411950 : move16();
1777 65411950 : move16();
1778 65411950 : *sine = sl;
1779 65411950 : *cosine = cl;
1780 : }
1781 :
1782 65411950 : return residual;
1783 : }
1784 :
1785 3533653 : Word16 getCosWord16( Word16 theta )
1786 : {
1787 : Word16 result, residual, sine, cosine;
1788 :
1789 3533653 : 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 3533653 : sine = negate( sine );
1793 3533653 : result = mac_r( L_mult0( sine, residual ), cosine, 16384 );
1794 :
1795 :
1796 3533653 : return result;
1797 : }
1798 :
1799 689378 : Word16 getSinWord16( Word16 theta )
1800 : {
1801 : Word16 sine;
1802 689378 : Word32 theta_new = L_sub( EVS_PI_BY_2_FX, theta );
1803 : Word16 l_theta;
1804 689378 : IF( GT_32( theta_new, EVS_PI_FX ) )
1805 : {
1806 61221 : l_theta = extract_l( L_sub( L_sub( theta_new, EVS_PI_FX ), EVS_PI_FX ) );
1807 : }
1808 628157 : 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 628157 : l_theta = extract_l( theta_new );
1815 : }
1816 689378 : sine = getCosWord16( l_theta );
1817 689378 : IF( EQ_16( sine, ONE_IN_Q14 ) )
1818 : {
1819 6706 : sine = MAX_16;
1820 : }
1821 : ELSE
1822 : {
1823 682672 : sine = shl( sine, 1 );
1824 : }
1825 689378 : return sine;
1826 : }
1827 :
1828 61878297 : Word16 getCosWord16R2(
1829 : Word16 theta )
1830 : {
1831 : Word16 result, residual, sine, cosine;
1832 :
1833 61878297 : 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 61878297 : sine = negate( sine );
1838 : /* Saturation has been included based on the recommendation from the PC group */
1839 61878297 : result = msu_r_sat( L_mult( sine, residual ), cosine, -32768 );
1840 : BASOP_SATURATE_WARNING_ON
1841 :
1842 61878297 : return result;
1843 : }
1844 :
1845 :
1846 26189840 : Word16 getSineWord16R2( Word16 theta )
1847 : {
1848 26189840 : IF( EQ_16( theta, (Word16) 0 ) )
1849 : {
1850 7594114 : return 0;
1851 : }
1852 18595726 : ELSE IF( LT_16( theta, (Word16) -24576 ) )
1853 : {
1854 120665 : theta = add( add( theta, (Word16) 32767 ), (Word16) 1 );
1855 : }
1856 18595726 : 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 1700715 : Word16 getSqrtWord32( Word32 val )
1871 : {
1872 : Word32 L_temp, L_temp2;
1873 : Word16 temp, temp2;
1874 : Word16 exp, exp2;
1875 :
1876 : /* Calc Approximation */
1877 1700715 : exp2 = norm_l( val );
1878 1700715 : L_temp2 = L_shl( val, exp2 );
1879 1700715 : exp = sub( 31 - 32, exp2 );
1880 1700715 : L_temp = Isqrt_lc( L_temp2, &exp ); /* 12 clks */
1881 :
1882 1700715 : temp = round_fx_sat( L_temp );
1883 1700715 : L_temp = Mpy_32_16_1( L_temp2, temp ); /* 2 clks */
1884 :
1885 1700715 : 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 1700715 : temp = round_fx( L_temp );
1889 :
1890 : /* Too High? */
1891 1700715 : if ( L_msu0( val, temp, temp ) < 0 )
1892 : {
1893 : /* Reduce by 2 */
1894 1391862 : temp = sub( temp, 2 );
1895 : }
1896 : /* Too High? */
1897 1700715 : if ( L_msu0( val, temp, temp ) < 0 )
1898 : {
1899 : /* Reduce by 2 */
1900 156331 : temp = sub( temp, 2 );
1901 : }
1902 : /* Too High? */
1903 1700715 : if ( L_msu0( val, temp, temp ) < 0 )
1904 : {
1905 : /* Reduce by 2 */
1906 249 : temp = sub( temp, 2 );
1907 : }
1908 : /* Too High? */
1909 1700715 : if ( L_msu0( val, temp, temp ) < 0 )
1910 : {
1911 : /* Reduce by 1 */
1912 3 : temp = sub( temp, 1 );
1913 : }
1914 :
1915 : /* Try +1 */
1916 1700715 : temp2 = add( temp, 1 );
1917 : /* It fits? */
1918 1700715 : if ( L_msu0( val, temp2, temp2 ) >= 0 )
1919 : {
1920 : /* Yes */
1921 1171597 : temp = temp2;
1922 1171597 : move16();
1923 : }
1924 1700715 : 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 0 : Word16 findIndexOfMinWord64( Word64 *x, const Word16 len )
1947 : {
1948 : Word16 i, indx;
1949 :
1950 :
1951 0 : indx = 0;
1952 0 : move16();
1953 0 : FOR( i = 1; i < len; i++ )
1954 : {
1955 0 : if ( LT_64( x[i], x[indx] ) )
1956 : {
1957 0 : indx = i;
1958 0 : move16();
1959 : }
1960 : }
1961 :
1962 :
1963 0 : return indx;
1964 : }
1965 :
1966 :
1967 1530123834 : Word16 imult1616( Word16 x, Word16 y )
1968 : {
1969 1530123834 : assert( (int) x * (int) y < 32768 && (int) x * (int) y >= -32768 );
1970 1530123834 : return i_mult( x, y );
1971 : }
1972 :
1973 1467424165 : Word32 imult3216( Word32 x, Word16 y )
1974 : {
1975 : Word32 mh;
1976 : UWord16 ml;
1977 :
1978 1467424165 : Mpy_32_16_ss( x, y, &mh, &ml );
1979 :
1980 1467424165 : mh = L_shl( mh, 15 );
1981 1467424165 : ml = (UWord16) lshr( (Word16) ml, 1 );
1982 :
1983 1467424165 : 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 1848570 : 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 1848570 : sx = norm_s( x );
2014 1848570 : sy = norm_s( y );
2015 1848570 : x = shl( x, sx );
2016 1848570 : y = shl( y, sy );
2017 :
2018 1848570 : if ( x >= y )
2019 : {
2020 1003817 : x = shr( x, 1 );
2021 1003817 : sx = sub( sx, 1 );
2022 : }
2023 :
2024 : /* divide and shift */
2025 1848570 : y = shr( div_s( x, y ), sub( 15, sub( sy, sx ) ) );
2026 :
2027 1848570 : return y;
2028 : }
2029 :
2030 :
2031 50476773 : Word16 idiv1616( Word16 x, Word16 y )
2032 : {
2033 : Word16 s, num, den, sign;
2034 :
2035 :
2036 50476773 : num = abs_s( x );
2037 50476773 : den = abs_s( y );
2038 :
2039 50476773 : sign = s_and( s_xor( x, y ), -32768 /* 0x8000 */ );
2040 :
2041 : /* make num > den */
2042 50476773 : s = add( sub( norm_s( den ), norm_s( num ) ), 1 );
2043 50476773 : s = s_max( s, 0 );
2044 :
2045 50476773 : den = shl( den, s );
2046 :
2047 : /* divide and shift */
2048 50476773 : y = shr( div_s( num, den ), sub( 15, s ) );
2049 :
2050 50476773 : if ( 0 != sign )
2051 : {
2052 148474 : y = negate( y );
2053 : }
2054 :
2055 :
2056 50476773 : return y;
2057 : }
2058 :
2059 0 : Word16 idiv1616_1( Word16 x, Word16 y )
2060 : {
2061 0 : IF( L_mult0( x, y ) < 0 )
2062 : {
2063 0 : return negate( idiv1616( abs_s( x ), abs_s( y ) ) );
2064 : }
2065 0 : ELSE IF( L_mult0( x, y ) > 0 )
2066 : {
2067 0 : return idiv1616( x, y );
2068 : }
2069 : ELSE
2070 : {
2071 0 : return 0;
2072 : }
2073 : }
2074 :
2075 11653234 : Word32 norm_llQ31( /* o : normalized result Q31 */
2076 : Word32 L_c, /* i : upper bits of accu Q-1 */
2077 : Word32 L_sum, /* i : lower bits of accu, unsigned Q31 */
2078 : Word16 *exp /* o : exponent of result in [-32,31] Q0 */
2079 : )
2080 : {
2081 : Word16 i;
2082 : Word32 L_tmp;
2083 : #ifdef BASOP_NOGLOB_DECLARE_LOCAL
2084 : #ifndef ISSUE_1836_replace_overflow_libcom
2085 : Flag Overflow = 0;
2086 : #endif
2087 11653234 : Flag Carry = 0;
2088 : #endif /* BASOP_NOGLOB */
2089 :
2090 : /* Move MSBit of L_sum into L_c */
2091 11653234 : Carry = 0;
2092 : #ifdef ISSUE_1836_replace_overflow_libcom
2093 11653234 : L_tmp = L_add_c( L_sum, L_sum, &Carry ); /* L_tmp = L_sum << 1 */
2094 11653234 : L_c = L_add_c( L_c, L_c, &Carry );
2095 : #else
2096 : L_tmp = L_add_co( L_sum, L_sum, &Carry, &Overflow ); /* L_tmp = L_sum << 1 */
2097 : L_c = L_add_co( L_c, L_c, &Carry, &Overflow );
2098 : #endif
2099 11653234 : L_add( 0, 0 );
2100 11653234 : test();
2101 11653234 : IF( ( L_c != (Word32) 0L ) && ( L_c != (Word32) 0xFFFFFFFFL ) )
2102 : {
2103 3376657 : i = norm_l( L_c );
2104 3376657 : L_c = L_shl( L_c, i );
2105 3376657 : i = sub( 31, i ); /* positive exponent */
2106 3376657 : L_sum = L_lshr( L_tmp, 1 ); /* L_sum with MSBit=0 */
2107 3376657 : L_sum = L_lshr( L_sum, i );
2108 3376657 : L_sum = L_add( L_c, L_sum );
2109 : }
2110 : ELSE
2111 : {
2112 8276577 : i = -32;
2113 8276577 : move16(); /* default exponent, if total sum=0 */
2114 8276577 : IF( L_sum )
2115 : {
2116 2682033 : i = norm_l( L_sum );
2117 2682033 : L_sum = L_shl( L_sum, i );
2118 2682033 : i = negate( i ); /* negative or zero exponent */
2119 : }
2120 : }
2121 11653234 : *exp = i;
2122 11653234 : move16();
2123 11653234 : return L_sum;
2124 : }
2125 :
2126 : #ifndef FIX_ISSUE_1817_REPLACE_CARRY_OVERFLOW
2127 : /* note: now available in basop_util.h */
2128 : Word32 w_norm_llQ31( Word64 L_sum, Word16 *exp );
2129 : Word32 w_norm_llQ31( /* o : normalized result Q31 */
2130 : Word64 L_sum, /* i : upper and lower bits of accu, unsigned Q31 */
2131 : Word16 *exp /* o : exponent of result in [-32,31] Q0 */
2132 : )
2133 : {
2134 : Word32 L_tmp;
2135 : Word16 exp_val;
2136 : Word64 L64_inp64 = L_sum;
2137 : move64();
2138 :
2139 : L64_inp64 = W_shl( L64_inp64, 1 );
2140 : exp_val = W_norm( L64_inp64 );
2141 : L64_inp64 = W_shl( L64_inp64, exp_val );
2142 : exp_val = sub( 31, exp_val );
2143 : if ( EQ_64( L_sum, 0 ) )
2144 : {
2145 : exp_val = -32;
2146 : move16();
2147 : }
2148 : *exp = exp_val;
2149 : move16();
2150 :
2151 : L_tmp = W_extract_h( L64_inp64 );
2152 : return L_tmp;
2153 : }
2154 : #endif
2155 :
2156 29425 : Word32 Dot_product16HQ( /* o : normalized result Q31 */
2157 : const Word32 L_off, /* i : initial sum value Qn */
2158 : const Word16 x[], /* i : x vector Qn */
2159 : const Word16 y[], /* i : y 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 : Word32 L_sum;
2166 : Word64 L_sum64;
2167 :
2168 29425 : L_sum64 = W_deposit32_l( L_off );
2169 :
2170 1270304 : FOR( i = 0; i < lg; i++ )
2171 : {
2172 1240879 : L_sum64 = W_mac_16_16( L_sum64, x[i], y[i] );
2173 : }
2174 :
2175 29425 : L_sum = w_norm_llQ31( L_sum64, exp );
2176 29425 : return L_sum;
2177 : }
2178 :
2179 102675 : Word32 sum_array_norm( /* o : normalized result Q31 */
2180 : const Word32 x[], /* i : x vector Qn */
2181 : const Word16 lg, /* i : vector length, range [0..7FFF] Q0 */
2182 : Word16 *exp /* o : exponent of result in [-32,31] Q0 */
2183 : )
2184 : {
2185 : Word16 i;
2186 102675 : Word64 W_tmp = 0;
2187 : Word32 L_tmp;
2188 29879597 : FOR( i = 0; i < lg; i++ )
2189 : {
2190 29776922 : W_tmp = W_add( W_tmp, x[i] ); /*Q31*/
2191 : }
2192 102675 : L_tmp = w_norm_llQ31( W_tmp, exp ); /*Q31 - *exp*/
2193 102675 : return L_tmp;
2194 : }
2195 :
2196 :
2197 1614 : Word32 Norm32Norm( const Word32 *x, const Word16 headroom, const Word16 length, Word16 *result_e )
2198 : {
2199 : Word32 L_tmp, L_tmp2;
2200 : Word16 i, shift, tmp;
2201 :
2202 1614 : move16();
2203 1614 : shift = headroom;
2204 :
2205 1614 : L_tmp = L_deposit_l( 0 );
2206 :
2207 194148 : FOR( i = 0; i < length; i++ )
2208 : {
2209 192534 : L_tmp2 = L_sub( L_tmp, 0x40000000 );
2210 192534 : if ( L_tmp2 >= 0 )
2211 3225 : shift = sub( shift, 1 );
2212 192534 : if ( L_tmp2 >= 0 )
2213 3225 : L_tmp = L_shr( L_tmp, 2 );
2214 :
2215 192534 : tmp = round_fx_sat( L_shl_sat( x[i], shift ) );
2216 192534 : L_tmp = L_mac0_sat( L_tmp, tmp, tmp ); /* exponent = (1-shift*2) , Q(30+shift*2) */
2217 : }
2218 :
2219 1614 : move16();
2220 1614 : *result_e = sub( 1, shl( shift, 1 ) );
2221 :
2222 1614 : return L_tmp;
2223 : }
2224 :
2225 399740 : Word32 Dot_productSq16HQ( /* o : normalized result Q31 */
2226 : const Word32 L_off, /* i : initial sum value Qn */
2227 : const Word16 x[], /* i : x vector Qn */
2228 : const Word16 lg, /* i : vector length, range [0..7FFF] Q0 */
2229 : Word16 *exp /* o : exponent of result in [-32,31] Q0 */
2230 : )
2231 : {
2232 : Word16 i;
2233 : Word32 L_sum;
2234 : Word64 L_sum64;
2235 :
2236 399740 : L_sum64 = W_deposit32_l( L_off );
2237 :
2238 44792108 : FOR( i = 0; i < lg; i++ )
2239 : {
2240 44392368 : L_sum64 = W_mac_16_16( L_sum64, x[i], x[i] );
2241 : }
2242 399740 : L_sum = w_norm_llQ31( L_sum64, exp );
2243 :
2244 399740 : return L_sum;
2245 : }
2246 :
2247 125163 : Word32 dotp_s_fx( const Word16 *x, const Word16 *y, const Word16 n, Word16 s )
2248 : {
2249 : Word16 i;
2250 : Word16 n2;
2251 : Word32 L_tmp;
2252 : Word32 L_sum;
2253 :
2254 :
2255 125163 : L_sum = 0;
2256 125163 : move32();
2257 :
2258 125163 : n2 = shr( n, 1 );
2259 :
2260 125163 : s = sub( s, 1 );
2261 :
2262 3935169 : FOR( i = 0; i < n2; i++ )
2263 : {
2264 3810006 : L_tmp = L_mult0( x[2 * i], y[2 * i] );
2265 3810006 : L_tmp = L_mac0( L_tmp, x[2 * i + 1], y[2 * i + 1] );
2266 3810006 : L_sum = L_add( L_sum, L_shr( L_tmp, s ) );
2267 : }
2268 :
2269 125163 : IF( s_and( n, 1 ) )
2270 : {
2271 70047 : L_tmp = L_mult0( x[n - 1], y[n - 1] );
2272 70047 : L_sum = L_add( L_sum, L_shr( L_tmp, s ) );
2273 : }
2274 :
2275 :
2276 125163 : return L_sum;
2277 : }
2278 :
2279 :
2280 57972322 : Word32 BASOP_util_Pow2(
2281 : const Word32 exp_m,
2282 : const Word16 exp_e,
2283 : Word16 *result_e )
2284 : {
2285 : static const Word16 pow2Coeff[8] = {
2286 : 22713 /*0.693147180559945309417232121458177 Q15*/, /* ln(2)^1 /1! */
2287 : 7872 /*0.240226506959100712333551263163332 Q15*/, /* ln(2)^2 /2! */
2288 : 1819 /*0.0555041086648215799531422637686218 Q15*/, /* ln(2)^3 /3! */
2289 : 315 /*0.00961812910762847716197907157365887 Q15*/, /* ln(2)^4 /4! */
2290 : 44 /*0.00133335581464284434234122219879962 Q15*/, /* ln(2)^5 /5! */
2291 : 5 /*1.54035303933816099544370973327423e-4 Q15*/, /* ln(2)^6 /6! */
2292 : 0 /*1.52527338040598402800254390120096e-5 Q15*/, /* ln(2)^7 /7! */
2293 : 0 /*1.32154867901443094884037582282884e-6 Q15*/ /* ln(2)^8 /8! */
2294 : };
2295 :
2296 : Word32 frac_part, tmp_frac, result_m;
2297 : Word16 int_part;
2298 :
2299 57972322 : int_part = 0; /* to avoid compilation warnings */
2300 57972322 : frac_part = 0; /* to avoid compilation warnings */
2301 :
2302 57972322 : IF( exp_e > 0 )
2303 : {
2304 : /* "+ 1" compensates L_shr(,1) of the polynomial evaluation at the loop end. */
2305 :
2306 57904322 : int_part = add( 1, extract_l( L_shr( exp_m, sub( 31, exp_e ) ) ) );
2307 57904322 : frac_part = L_lshl( exp_m, exp_e );
2308 57904322 : frac_part = L_and( 0x7FFFFFFF, frac_part );
2309 : }
2310 57972322 : if ( exp_e <= 0 )
2311 68000 : frac_part = L_shl( exp_m, exp_e );
2312 57972322 : if ( exp_e <= 0 )
2313 : {
2314 68000 : int_part = 1;
2315 68000 : move16();
2316 : }
2317 :
2318 : /* Best accuracy is around 0, so try to get there with the fractional part. */
2319 57972322 : IF( ( tmp_frac = L_sub( frac_part, 1073741824l /*0.5 Q31*/ ) ) >= 0 )
2320 : {
2321 28324865 : int_part = add( int_part, 1 );
2322 28324865 : frac_part = L_sub( tmp_frac, 1073741824l /*0.5 Q31*/ );
2323 : }
2324 29647457 : ELSE IF( ( tmp_frac = L_add( frac_part, 1073741824l /*0.5 Q31*/ ) ) < 0 )
2325 : {
2326 0 : int_part = sub( int_part, 1 );
2327 0 : frac_part = L_add( tmp_frac, 1073741824l /*0.5 Q31*/ );
2328 : }
2329 :
2330 : /* Evaluate taylor polynomial which approximates 2^x */
2331 : {
2332 : Word32 p;
2333 : Word16 i;
2334 :
2335 :
2336 : /* First taylor series coefficient a_0 = 1.0, scaled by 0.5 due to L_shr(,1). */
2337 57972322 : result_m = L_add( 1073741824l /*1.0/2.0 Q31*/, L_shr( Mpy_32_16_1( frac_part, pow2Coeff[0] ), 1 ) );
2338 57972322 : p = Mpy_32_32( frac_part, frac_part );
2339 405806254 : FOR( i = 1; i < 7; i++ )
2340 : {
2341 : /* next taylor series term: a_i * x^i, x=0 */
2342 347833932 : result_m = L_add( result_m, L_shr( Mpy_32_16_1( p, pow2Coeff[i] ), 1 ) );
2343 347833932 : p = Mpy_32_32( p, frac_part );
2344 : }
2345 57972322 : result_m = L_add( result_m, L_shr( Mpy_32_16_1( p, pow2Coeff[i] ), 1 ) );
2346 : }
2347 57972322 : *result_e = int_part;
2348 57972322 : move16();
2349 57972322 : return result_m;
2350 : }
2351 :
2352 0 : Word16 findIndexOfMaxWord32( Word32 *x, const Word16 len )
2353 : {
2354 : Word16 i, indx;
2355 :
2356 :
2357 0 : indx = 0;
2358 0 : move16();
2359 0 : FOR( i = 1; i < len; i++ )
2360 : {
2361 0 : if ( GT_32( x[i], x[indx] ) )
2362 : {
2363 0 : indx = i;
2364 0 : move16();
2365 : }
2366 : }
2367 :
2368 :
2369 0 : return indx;
2370 : }
2371 :
2372 53742402 : Word16 getNormReciprocalWord16( Word16 x )
2373 : {
2374 :
2375 53742402 : assert( x < (Word16) ( sizeof( BASOP_util_normReciprocal ) / sizeof( BASOP_util_normReciprocal[0] ) ) );
2376 :
2377 53742402 : return extract_h( BASOP_util_normReciprocal[x] );
2378 : }
2379 68 : Word16 getNormReciprocalWord16Scale( Word16 x, Word16 s )
2380 : {
2381 :
2382 68 : assert( x < (Word16) ( sizeof( BASOP_util_normReciprocal ) / sizeof( BASOP_util_normReciprocal[0] ) ) );
2383 :
2384 68 : return round_fx( L_shl( BASOP_util_normReciprocal[x], s ) );
2385 : }
2386 :
2387 :
2388 : /*! r: result of division x/y, not normalized */
2389 329474102 : Word16 BASOP_Util_Divide3216_Scale(
2390 : Word32 x, /* i : numerator, signed */
2391 : Word16 y, /* i : denominator, signed */
2392 : Word16 *s ) /* o : scaling, 0, if x==0 */
2393 : {
2394 : Word16 z;
2395 : Word16 sx;
2396 : Word16 sy;
2397 : Word16 sign;
2398 :
2399 : /*assert (x > (Word32)0);
2400 : assert (y >= (Word16)0);*/
2401 :
2402 : /* check, if numerator equals zero, return zero then */
2403 329474102 : IF( x == (Word32) 0 )
2404 : {
2405 29347389 : *s = 0;
2406 29347389 : move16();
2407 :
2408 29347389 : return ( (Word16) 0 );
2409 : }
2410 :
2411 300126713 : sign = s_xor( extract_h( x ), y ); /* just to exor the sign bits */
2412 : BASOP_SATURATE_WARNING_OFF
2413 300126713 : x = L_abs( x );
2414 300126713 : y = abs_s( y );
2415 : BASOP_SATURATE_WARNING_ON
2416 300126713 : sx = sub( norm_l( x ), 1 );
2417 300126713 : x = L_shl( x, sx );
2418 300126713 : sy = norm_s( y );
2419 300126713 : y = shl( y, sy );
2420 300126713 : *s = sub( sy, sx );
2421 300126713 : move16();
2422 :
2423 300126713 : z = div_s( round_fx( x ), y );
2424 :
2425 300126713 : if ( sign < 0 ) /* if sign bits differ, negate the result */
2426 : {
2427 3840318 : z = negate( z );
2428 : }
2429 :
2430 300126713 : return z;
2431 : }
2432 :
2433 :
2434 : /*************************************************************************
2435 : *
2436 : * FUNCTION: Log2_norm()
2437 : *
2438 : * PURPOSE: Computes log2(L_x, exp), where L_x is positive and
2439 : * normalized, and exp is the normalisation exponent
2440 : * If L_x is negative or zero, the result is 0.
2441 : *
2442 : * DESCRIPTION:
2443 : * The function Log2(L_x) is approximated by a table and linear
2444 : * interpolation. The following steps are used to compute Log2(L_x)
2445 : *
2446 : * 1- exponent = 30-norm_exponent
2447 : * 2- i = bit25-b31 of L_x; 32<=i<=63 (because of normalization).
2448 : * 3- a = bit10-b24
2449 : * 4- i -=32
2450 : * 5- fraction = table[i]<<16 - (table[i] - table[i+1]) * a * 2
2451 : *
2452 : *************************************************************************/
2453 :
2454 944765 : Word32 BASOP_Util_fPow(
2455 : Word32 base_m,
2456 : Word16 base_e,
2457 : Word32 exp_m,
2458 : Word16 exp_e,
2459 : Word16 *result_e )
2460 : {
2461 :
2462 : Word16 ans_lg2_e, base_lg2_e;
2463 : Word32 base_lg2_m, ans_lg2_m, result_m;
2464 : Word16 shift;
2465 :
2466 944765 : test();
2467 944765 : IF( ( base_m == 0 ) && ( exp_m != 0 ) )
2468 : {
2469 0 : *result_e = 0;
2470 0 : move16();
2471 0 : return 0;
2472 : }
2473 : /* Calc log2 of base */
2474 944765 : shift = norm_l( base_m );
2475 944765 : base_m = L_shl( base_m, shift );
2476 944765 : base_e = sub( base_e, shift );
2477 944765 : base_lg2_m = BASOP_Util_Log2( base_m );
2478 :
2479 : /* shift: max left shift such that neither base_e or base_lg2_m saturate. */
2480 944765 : shift = sub( s_min( norm_s( base_e ), WORD16_BITS - 1 - LD_DATA_SCALE ), 1 );
2481 : /* Compensate shift into exponent of result. */
2482 944765 : base_lg2_e = sub( WORD16_BITS - 1, shift );
2483 944765 : 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 ) ) );
2484 :
2485 : /* Prepare exp */
2486 944765 : shift = norm_l( exp_m );
2487 944765 : exp_m = L_shl( exp_m, shift );
2488 944765 : exp_e = sub( exp_e, shift );
2489 :
2490 : /* Calc base pow exp */
2491 944765 : ans_lg2_m = Mpy_32_32( base_lg2_m, exp_m );
2492 944765 : ans_lg2_e = add( exp_e, base_lg2_e );
2493 :
2494 : /* Calc antilog */
2495 944765 : result_m = BASOP_util_Pow2( ans_lg2_m, ans_lg2_e, result_e );
2496 :
2497 944765 : return result_m;
2498 : }
2499 :
2500 : /*___________________________________________________________________________
2501 : | |
2502 : | Function Name : Dot_product12_offs() |
2503 : | |
2504 : | Compute scalar product of <x[],y[]> using accumulator. |
2505 : | The parameter 'L_off' is added to the accumulation result. |
2506 : | The result is normalized (in Q31) with exponent (0..30). |
2507 : | Notes: |
2508 : | o data in x[],y[] must provide enough headroom for accumulation |
2509 : | o L_off must correspond in format with product of x,y |
2510 : | Example: 0.01f for Q9 x Q9: 0x0000147B in Q19 |
2511 : | means: L_off = FL2WORD32_SCALE(0.01,31-19) |
2512 : |---------------------------------------------------------------------------|
2513 : | Algorithm: |
2514 : | |
2515 : | dot_product = L_off + sum(x[i]*y[i]) i=0..N-1 |
2516 : |___________________________________________________________________________|
2517 : */
2518 18640 : Word32 Dot_product12_offs( /* (o) Q31: normalized result (1 < val <= -1) */
2519 : const Word16 x[], /* (i) 12bits: x vector */
2520 : const Word16 y[], /* (i) 12bits: y vector */
2521 : const Word16 lg, /* (i) : vector length in range [1..256] */
2522 : Word16 *exp, /* (o) : exponent of result (0..+30) */
2523 : Word32 L_off /* (i) initial summation offset / 2 */
2524 : )
2525 : {
2526 : Word16 i, sft;
2527 : Word32 L_sum;
2528 :
2529 : Word64 L_sum64;
2530 :
2531 18640 : L_sum64 = W_deposit32_l( L_off );
2532 1717152 : FOR( i = 0; i < lg; i++ )
2533 : {
2534 1698512 : L_sum64 = W_mac0_16_16( L_sum64, x[i], y[i] );
2535 : }
2536 18640 : L_sum = W_sat_l( L_sum64 );
2537 :
2538 : /* Normalize acc in Q31 */
2539 :
2540 18640 : sft = norm_l( L_sum );
2541 18640 : if ( exp != NULL )
2542 : {
2543 18640 : L_sum = L_shl( L_sum, sft );
2544 : }
2545 :
2546 : /* exponent = 0..30, when L_sum != 0 */
2547 18640 : if ( L_sum != 0 )
2548 : {
2549 18640 : sft = sub( 31, sft );
2550 : }
2551 :
2552 18640 : if ( exp != NULL )
2553 : {
2554 18640 : *exp = sft;
2555 18640 : move16();
2556 : }
2557 :
2558 18640 : return L_sum;
2559 : }
2560 :
2561 0 : Word32 Dot_product15_offs( /* (o) Q31: normalized result (1 < val <= -1) */
2562 : const Word16 x[], /* (i) 15bits: x vector */
2563 : const Word16 y[], /* (i) 15bits: y vector */
2564 : const Word16 lg, /* (i) : vector length in range [1..256] */
2565 : Word16 *exp, /* (o) : exponent of result (0..+30) */
2566 : Word32 L_off /* (i) initial summation offset */
2567 : )
2568 : {
2569 : Word16 i, sft, fac, ld;
2570 : Word32 L_sum;
2571 :
2572 0 : ld = sub( 14, norm_s( lg ) );
2573 0 : fac = shr( -32768, ld );
2574 0 : L_sum = L_shr( L_off, ld );
2575 :
2576 0 : FOR( i = 0; i < lg; i++ )
2577 : {
2578 0 : L_sum = L_add( L_sum, Mpy_32_16_1( L_msu( 0, y[i], fac ), x[i] ) );
2579 : }
2580 :
2581 : /* Avoid returning 0 */
2582 0 : if ( L_sum == 0 )
2583 : {
2584 0 : L_sum = L_add( L_sum, 1 );
2585 : }
2586 :
2587 : /* Normalize acc in Q31 */
2588 0 : sft = norm_l( L_sum );
2589 0 : L_sum = L_shl( L_sum, sft );
2590 :
2591 : /* exponent = 0..30, when L_sum != 0 */
2592 0 : if ( L_sum != 0 )
2593 : {
2594 0 : sft = add( ld, sub( 30, sft ) );
2595 : }
2596 :
2597 0 : *exp = sft;
2598 0 : move16();
2599 :
2600 0 : return L_sum;
2601 : }
2602 :
2603 20145801 : Word16 BASOP_Util_Cmp_Mant32Exp /*!< o: flag: result of comparison */
2604 : /* 0, if a == b */
2605 : /* 1, if a > b */
2606 : /* -1, if a < b */
2607 : ( Word32 a_m, /*!< i: Mantissa of 1st operand a */
2608 : Word16 a_e, /*!< i: Exponent of 1st operand a */
2609 : Word32 b_m, /*!< i: Mantissa of 2nd operand b */
2610 : Word16 b_e ) /*!< i: Exponent of 2nd operand b */
2611 :
2612 : {
2613 : Word32 diff_m;
2614 : Word16 diff_e, shift, result;
2615 :
2616 :
2617 : /*
2618 : This function compares two input parameters, both represented by a 32-bit mantissa and a 16-bit exponent.
2619 : If both values are identical, 0 is returned.
2620 : If a is greater b, 1 is returned.
2621 : If a is less than b, -1 is returned.
2622 : */
2623 :
2624 : /* Check, if both mantissa and exponents are identical, when normalized: return 0 */
2625 20145801 : shift = norm_l( a_m );
2626 20145801 : if ( shift )
2627 6333270 : a_m = L_shl( a_m, shift );
2628 20145801 : if ( shift )
2629 6333270 : a_e = sub( a_e, shift );
2630 :
2631 20145801 : shift = norm_l( b_m );
2632 20145801 : if ( shift )
2633 9676165 : b_m = L_shl( b_m, shift );
2634 20145801 : if ( shift )
2635 9676165 : b_e = sub( b_e, shift );
2636 :
2637 : /* align exponent, if any mantissa is zero */
2638 20145801 : if ( !a_m )
2639 : {
2640 656288 : a_e = b_e;
2641 656288 : move16();
2642 : }
2643 20145801 : if ( !b_m )
2644 : {
2645 2166644 : b_e = a_e;
2646 2166644 : move16();
2647 : }
2648 :
2649 : BASOP_SATURATE_WARNING_OFF_EVS
2650 20145801 : diff_m = L_sub_sat( a_m, b_m );
2651 : BASOP_SATURATE_WARNING_ON_EVS
2652 20145801 : diff_e = sub( a_e, b_e );
2653 :
2654 20145801 : test();
2655 20145801 : IF( diff_m == 0 && diff_e == 0 )
2656 : {
2657 80780 : return 0;
2658 : }
2659 :
2660 : /* Check sign, exponent and mantissa to identify, whether a is greater b or not */
2661 20065021 : result = sub( 0, 1 );
2662 :
2663 20065021 : IF( a_m >= 0 )
2664 : {
2665 : /* a is positive */
2666 20063001 : if ( b_m < 0 )
2667 : {
2668 0 : result = 1;
2669 0 : move16();
2670 : }
2671 :
2672 20063001 : test();
2673 20063001 : test();
2674 20063001 : test();
2675 20063001 : if ( ( b_m >= 0 ) && ( ( diff_e > 0 ) || ( diff_e == 0 && diff_m > 0 ) ) )
2676 : {
2677 12931811 : result = 1;
2678 12931811 : move16();
2679 : }
2680 : }
2681 : ELSE
2682 : {
2683 : /* a is negative */
2684 2020 : test();
2685 2020 : test();
2686 2020 : test();
2687 2020 : if ( ( b_m < 0 ) && ( ( diff_e < 0 ) || ( diff_e == 0 && diff_m > 0 ) ) )
2688 : {
2689 20 : result = 1;
2690 20 : move16();
2691 : }
2692 : }
2693 20065021 : return result;
2694 : }
2695 :
2696 : /*
2697 :
2698 : headroom is introduced into acc
2699 : */
2700 :
2701 :
2702 2739786546 : Word32 BASOP_Util_Add_Mant32Exp /* o : normalized result mantissa */
2703 : ( Word32 a_m, /* i : Mantissa of 1st operand a */
2704 : Word16 a_e, /* i : Exponent of 1st operand a */
2705 : Word32 b_m, /* i : Mantissa of 2nd operand b */
2706 : Word16 b_e, /* i : Exponent of 2nd operand b */
2707 : Word16 *ptr_e ) /* o : exponent of result */
2708 : {
2709 : Word32 L_tmp;
2710 : Word16 shift;
2711 :
2712 : /* Compare exponents: the difference is limited to +/- 30
2713 : The Word32 mantissa of the operand with lower exponent is shifted right by the exponent difference.
2714 : Then, the unshifted mantissa of the operand with the higher exponent is added. The addition result
2715 : is normalized and the result represents the mantissa to return. The returned exponent takes into
2716 : account all shift operations.
2717 : */
2718 :
2719 2739786546 : if ( !a_m )
2720 : {
2721 817631582 : a_e = b_e;
2722 817631582 : move16();
2723 : }
2724 :
2725 2739786546 : if ( !b_m )
2726 : {
2727 400759512 : b_e = a_e;
2728 400759512 : move16();
2729 : }
2730 :
2731 2739786546 : shift = sub( a_e, b_e );
2732 2739786546 : shift = s_max( -31, shift );
2733 2739786546 : shift = s_min( 31, shift );
2734 2739786546 : if ( shift < 0 )
2735 : {
2736 : /* exponent of b is greater than exponent of a, shr a_m */
2737 597383682 : a_m = L_shl( a_m, shift );
2738 : }
2739 2739786546 : if ( shift > 0 )
2740 : {
2741 : /* exponent of a is greater than exponent of b */
2742 1079875459 : b_m = L_shr( b_m, shift );
2743 : }
2744 2739786546 : a_e = add( s_max( a_e, b_e ), 1 );
2745 2739786546 : L_tmp = L_add( L_shr( a_m, 1 ), L_shr( b_m, 1 ) );
2746 2739786546 : shift = norm_l( L_tmp );
2747 2739786546 : if ( shift )
2748 2213336257 : L_tmp = L_shl( L_tmp, shift );
2749 2739786546 : if ( L_tmp == 0 )
2750 : {
2751 345581927 : a_e = 0;
2752 345581927 : move16();
2753 : }
2754 2739786546 : if ( L_tmp != 0 )
2755 2394204619 : a_e = sub( a_e, shift );
2756 2739786546 : *ptr_e = a_e;
2757 :
2758 2739786546 : return ( L_tmp );
2759 : }
2760 :
2761 :
2762 : static const Word16 shift_lc[] = { 9, 10 };
2763 :
2764 0 : Word32 Isqrt_lc1(
2765 : Word32 frac, /* i : Q31: normalized value (1.0 < frac <= 0.5) */
2766 : Word16 *exp /* i/o: exponent (value = frac x 2^exponent) */
2767 : )
2768 : {
2769 : Word16 i, a;
2770 : Word32 L_tmp;
2771 :
2772 0 : IF( frac <= (Word32) 0 )
2773 : {
2774 0 : *exp = 0;
2775 0 : move16();
2776 0 : return 0x7fffffff; /*0x7fffffff*/
2777 : }
2778 :
2779 : /* If exponant odd -> shift right by 10 (otherwise 9) */
2780 0 : L_tmp = L_shr( frac, shift_lc[s_and( *exp, 1 )] );
2781 :
2782 : /* 1) -16384 to shift left and change sign */
2783 : /* 2) 32768 to Add 1 to Exponent like it was divided by 2 */
2784 : /* 3) We let the mac_r add another 0.5 because it imitates */
2785 : /* the behavior of shr on negative number that should */
2786 : /* not be rounded towards negative infinity. */
2787 : /* It replaces: */
2788 : /* *exp = negate(shr(sub(*exp, 1), 1)); move16(); */
2789 0 : *exp = mac_r( 32768, *exp, -16384 );
2790 0 : move16();
2791 :
2792 0 : a = extract_l( L_tmp ); /* Extract b10-b24 */
2793 0 : a = lshr( a, 1 );
2794 :
2795 0 : i = mac_r( L_tmp, -16 * 2 - 1, 16384 ); /* Extract b25-b31 minus 16 */
2796 :
2797 0 : L_tmp = L_msu( L_table_isqrt[i], table_isqrt_diff[i], a ); /* table[i] << 16 - diff*a*2 */
2798 :
2799 0 : return L_tmp;
2800 : }
2801 :
2802 :
2803 11005 : void bufferCopyFx(
2804 : Word16 *src, /*<! Qx pointer to input buffer */
2805 : Word16 *dest, /*<! Qx pointer to output buffer */
2806 : Word16 length, /*<! Q0 length of buffer to copy */
2807 : Word16 Qf_src, /*<! Q0 Q format (frac-bits) of source buffer */
2808 : Word16 Qf_dest, /*<! Q0 Q format (frac-bits) of dest buffer */
2809 : Word16 Q_src, /*<! Q0 exponent of source buffer */
2810 : Word16 Q_dest /*<! Q0 exponent of destination buffer */
2811 : )
2812 : {
2813 : Word16 tmp_16, i;
2814 :
2815 : /*Copy( st->old_exc, exc_buf, st->old_exc_len);*/
2816 11005 : tmp_16 = sub( sub( Qf_src, Qf_dest ), sub( Q_src, Q_dest ) );
2817 11005 : IF( tmp_16 > 0 ) /*if value will be shifted right, do a multiplication with rounding ->preserves more accuracy*/
2818 : {
2819 3045 : tmp_16 = shl( 1, sub( 15, tmp_16 ) );
2820 994458 : FOR( i = 0; i < length; i++ )
2821 : {
2822 991413 : *( dest + i ) = mult_r( *( src + i ), tmp_16 );
2823 991413 : move16();
2824 : }
2825 : }
2826 7960 : ELSE IF( tmp_16 < 0 ) /*leftshift - no accuracy preservation needed*/
2827 : {
2828 438981 : FOR( i = 0; i < length; i++ )
2829 : {
2830 437303 : *( dest + i ) = shr_sat( *( src + i ), tmp_16 );
2831 437303 : move16();
2832 : }
2833 : }
2834 : ELSE /*no shift, simply copy*/
2835 : {
2836 2373450 : FOR( i = 0; i < length; i++ )
2837 : {
2838 2367168 : *( dest + i ) = *( src + i );
2839 2367168 : move16();
2840 : }
2841 : }
2842 11005 : }
2843 :
2844 40062 : Word32 dotWord32_16_Mant32Exp( const Word32 *bufX32, /* i: 32-bit buffer with unknown headroom */
2845 : Word16 bufX32_exp, /* i: exponent of buffer bufX32 */
2846 : const Word16 *bufY16, /* i: 16-bit buffer quite right-aligned */
2847 : Word16 bufY16_exp, /* i: exponent of buffer bufY16 */
2848 : Word16 len, /* i: buffer len to process */
2849 : Word16 *exp ) /* o: result exponent */
2850 : {
2851 : Word32 L_sum;
2852 : Word16 shift, shift1, i;
2853 :
2854 :
2855 40062 : shift = getScaleFactor32( bufX32, len ); /* current available headroom */
2856 40062 : shift = sub( shift, sub( 14, norm_s( len ) ) ); /* reduced required headroom */
2857 40062 : L_sum = 0; /* Clear accu */
2858 40062 : move32();
2859 2503453 : FOR( i = 0; i < len; i++ )
2860 : {
2861 2463391 : L_sum = L_mac0( L_sum, round_fx( L_shl( bufX32[i], shift ) ), bufY16[i] );
2862 : }
2863 40062 : shift1 = norm_l( L_sum );
2864 40062 : L_sum = L_shl( L_sum, shift1 ); /* return value */
2865 :
2866 40062 : shift = sub( add( bufX32_exp, bufY16_exp ), add( shift, shift1 ) );
2867 40062 : shift = add( shift, 1 ); /* compensate for factor of 2 introduced by L_mac0 */
2868 : /* In case of NULL result, we want to have a 0 exponent too */
2869 40062 : if ( L_sum == 0 )
2870 569 : shift = 0;
2871 40062 : *exp = shift;
2872 40062 : move16();
2873 :
2874 :
2875 40062 : return L_sum;
2876 : }
2877 :
2878 454 : Word16 BASOP_Util_lin2dB( Word32 x, Word16 x_e, Word16 fEnergy )
2879 : {
2880 454 : assert( x >= 0 );
2881 :
2882 : /* log2 */
2883 454 : x = L_shr( BASOP_Util_Log2( x ), 25 - 16 ); /* Q16 */
2884 :
2885 : /* add exponent */
2886 454 : x = L_msu( x, x_e, -32768 /* 0x8000 */ );
2887 :
2888 : /* convert log2 to 20*log10 */
2889 454 : x = Mpy_32_16_1( x, 24660 /*6.0206f Q12*/ ); /* Q13 */
2890 :
2891 : /* if energy divide by 2 (->10*log10) */
2892 454 : if ( fEnergy != 0 )
2893 454 : x = L_shr( x, 1 );
2894 :
2895 : /* return dB as 7Q8 */
2896 454 : return round_fx( L_shl( x, 8 - 13 + 16 ) ); /* Q8 */
2897 : }
2898 :
2899 : /* --- fixp_atan() ---- */
2900 : #define Q_ATANINP ( 25 ) /* Input in q25, Output in q14 */
2901 : #define Q_ATANOUT ( 14 )
2902 : #define ATI_SF ( ( 32 - 1 ) - Q_ATANINP ) /* 6 */
2903 : #define ATO_SF ( ( 16 - 1 ) - Q_ATANOUT ) /* 1 ] -pi/2 .. pi/2 [ */
2904 : /* --- fixp_atan2() --- */
2905 : #define Q_ATAN2OUT ( 13 )
2906 : #define AT2O_SF ( ( 16 - 1 ) - Q_ATAN2OUT ) /* 2 ] -pi .. pi ] */
2907 :
2908 :
2909 34709744 : Word16 BASOP_util_atan2( /* o: atan2(y,x) [-pi,pi] Q13 */
2910 : Word32 y, /* i: */
2911 : Word32 x, /* i: */
2912 : Word16 e /* i: exponent difference (exp_y - exp_x) */
2913 : )
2914 : {
2915 : Word16 q;
2916 : Word32 at;
2917 34709744 : Word16 ret = -32768 /*-1.0f Q15*/;
2918 : Word16 sf, sfo, stf;
2919 : Word32 L_sign;
2920 :
2921 34709744 : if ( L_or( y, x ) == 0 )
2922 : {
2923 139157 : return 0;
2924 : }
2925 :
2926 34570587 : IF( x == 0 )
2927 : {
2928 462812 : ret = 12868 /*+EVS_PI/2 Q13*/;
2929 462812 : move16();
2930 462812 : if ( y < 0 )
2931 : {
2932 267522 : ret = negate( ret );
2933 : }
2934 :
2935 462812 : return ret;
2936 : }
2937 :
2938 : /* --- division */
2939 34107775 : L_sign = L_and( L_xor( x, y ), (Word32) 0x80000000 );
2940 :
2941 34107775 : q = 32767 /*1.0f Q15*/; /* y/x = neg/zero = -Inf */
2942 34107775 : sf = 0;
2943 : BASOP_SATURATE_WARNING_OFF_EVS
2944 34107775 : q = BASOP_Util_Divide3232_uu_1616_Scale( L_abs( y ), L_abs( x ), &sf );
2945 : BASOP_SATURATE_WARNING_ON_EVS
2946 :
2947 : BASOP_SATURATE_WARNING_OFF_EVS
2948 34107775 : if ( L_sign < 0 )
2949 16707845 : q = negate( q );
2950 : BASOP_SATURATE_WARNING_ON_EVS
2951 :
2952 34107775 : sfo = add( sf, e );
2953 :
2954 : /* --- atan() */
2955 34107775 : IF( GT_16( sfo, ATI_SF ) )
2956 : {
2957 : /* --- could not calc fixp_atan() here bec of input data out of range
2958 : ==> therefore give back boundary values */
2959 :
2960 1863856 : sfo = s_min( sfo, MAXSFTAB );
2961 :
2962 : /*q = FL2WORD16( 0.0f ); move16();*/
2963 :
2964 1863856 : if ( q > 0 )
2965 : {
2966 149913 : move16();
2967 149913 : q = +f_atan_expand_range[sfo - ATI_SF - 1];
2968 : }
2969 1863856 : if ( q < 0 )
2970 : {
2971 241046 : move16();
2972 241046 : q = -f_atan_expand_range[sfo - ATI_SF - 1];
2973 : }
2974 : }
2975 : ELSE
2976 : {
2977 : /* --- calc of fixp_atan() is possible; input data within range
2978 : ==> set q on fixed scale level as desired from fixp_atan() */
2979 32243919 : stf = sub( sfo, ATI_SF );
2980 :
2981 32243919 : at = L_deposit_h( q );
2982 32243919 : if ( stf < 0 )
2983 31983699 : at = L_shl( at, stf );
2984 :
2985 32243919 : q = BASOP_util_atan( at ); /* ATO_SF*/
2986 : }
2987 :
2988 :
2989 : /* --- atan2() */
2990 :
2991 34107775 : ret = shr( q, ( AT2O_SF - ATO_SF ) ); /* now AT2O_SF for atan2 */
2992 34107775 : IF( x < 0 )
2993 : {
2994 6251576 : if ( y >= 0 )
2995 : {
2996 3236782 : ret = add( ret, 25736 /*EVS_PI Q13*/ );
2997 : }
2998 6251576 : if ( y < 0 )
2999 : {
3000 3014794 : ret = sub( ret, 25736 /* EVS_PI Q13*/ );
3001 : }
3002 : }
3003 :
3004 34107775 : return ret;
3005 : }
3006 :
3007 : /* SNR of fixp_atan() = 56 dB*/
3008 : #define ONEBY3P56 0x26800000 /* 1.0/3.56 in q31*/
3009 : #define P281 0x00026000 /* 0.281 in q19*/
3010 : #define ONEP571 0x6487 /* 1.571 in q14*/
3011 :
3012 32243919 : Word16 BASOP_util_atan( /* o: atan(x) [-pi/2;pi/2] 1Q14 */
3013 : Word32 x /* i: input data (-64;64) 6Q25 */
3014 : )
3015 : {
3016 : Word16 sign, result, exp;
3017 : Word16 res_e;
3018 : Word16 tmp, xx;
3019 :
3020 :
3021 32243919 : sign = 0;
3022 32243919 : move16();
3023 32243919 : if ( x < 0 )
3024 : {
3025 16409501 : sign = 1;
3026 16409501 : move16();
3027 : }
3028 32243919 : x = L_abs( x );
3029 :
3030 : /* calc of arctan */
3031 32243919 : IF( LT_32( x, 1509950l /*0.045f/64.0f Q31*/ ) )
3032 : {
3033 4400724 : result = round_fx( L_shl( x, 5 ) ); /*Q14*/
3034 : /*BASOP_util_atan_16(0.0444059968): max error 0.0000567511, mean 0.000017, abs mean 0.000017*/
3035 : }
3036 : ELSE
3037 27843195 : IF( LT_32( x, ( L_shl( 1, Q_ATANINP ) - 8482560l /*0.00395 Q31*/ ) ) )
3038 : {
3039 14007356 : xx = round_fx( L_shl( x, 6 ) );
3040 14007356 : tmp = mult_r( xx, xx ); /* q15 * q15 - (16-1) = q15*/
3041 14007356 : tmp = mult_r( tmp, 0x1340 ); /* 15 * (ONEBY3P56) q14 - (16-1) = q14*/
3042 14007356 : tmp = add( tmp, 0x4000 ); /*L_shl(1,14) = 524288*/ /* q14 + q14 = q14 */
3043 14007356 : res_e = Q_ATANOUT - 15 + 14 - 16 + 1;
3044 14007356 : move16();
3045 14007356 : if ( GT_16( xx, tmp ) )
3046 : {
3047 2982497 : res_e = add( res_e, 1 );
3048 : }
3049 14007356 : if ( GT_16( xx, tmp ) )
3050 : {
3051 2982497 : xx = shr( xx, 1 );
3052 : }
3053 14007356 : result = div_s( xx, tmp );
3054 14007356 : result = msu_r( 0, result, shl( -32768, res_e ) );
3055 : /*BASOP_util_atan_16(0.7471138239): max error 0.0020029545, mean 0.000715, abs mean 0.000715*/
3056 : }
3057 13835839 : ELSE IF( LT_32( x, 42949673l /*1.28/64.0 Q31*/ ) )
3058 : {
3059 : Word16 delta_fix;
3060 5234780 : Word32 PI_BY_4 = 1686629684l /*3.1415926/4.0 Q31*/ / 2; /* pi/4 in q30 */
3061 :
3062 5234780 : delta_fix = round_fx( L_shl( L_sub( x, 33554432l /*1.0/64.0 Q31*/ ), 5 ) ); /* q30 */
3063 5234780 : result = round_fx( L_sub( L_add( PI_BY_4, L_msu( 0, delta_fix, -16384 ) ), ( L_mult0( delta_fix, delta_fix ) ) ) );
3064 : /* BASOP_Util_fPow(0.7472000122): max error 0.0020237688, mean 0.000026, abs mean 0.000520 */
3065 : }
3066 : ELSE
3067 : {
3068 8601059 : exp = sub( norm_l( x ), 1 );
3069 8601059 : xx = round_fx( L_shl( x, exp ) );
3070 : /* q25+exp * q25+exp - (16-1) = q19+2*exp*/
3071 8601059 : tmp = mac_r( L_shl( P281, shl( exp, 1 ) ), xx, xx ); /* q19+2*exp + q19+2*exp = q19+2*exp*/
3072 8601059 : res_e = norm_s( tmp );
3073 8601059 : result = div_s( xx, shl( tmp, res_e ) );
3074 8601059 : result = shl( result, add( add( Q_ATANOUT - Q_ATANINP /*-exp*/ + 19 /*+2*exp*/ - 16 + 1, res_e ), exp ) );
3075 8601059 : result = sub( ONEP571, result ); /* q14 + q14 = q14*/
3076 : /*BASOP_Util_fPow(1.2799999714): max error 0.0020168927, mean 0.000066, abs mean 0.000072*/
3077 : }
3078 :
3079 32243919 : if ( sign )
3080 : {
3081 16409501 : result = negate( result );
3082 : }
3083 :
3084 32243919 : return ( result );
3085 : }
3086 :
3087 : /* compare two positive normalized 16 bit mantissa/exponent values */
3088 : /* return value: positive if first value greater, negative if second value greater, zero if equal */
3089 539291 : Word16 compMantExp16Unorm( Word16 m1, Word16 e1, Word16 m2, Word16 e2 )
3090 : {
3091 : Word16 tmp;
3092 :
3093 539291 : assert( ( m1 >= 0x4000 ) && ( m2 >= 0x4000 ) ); /* comparisons below work only for normalized mantissas */
3094 :
3095 539291 : tmp = sub( e1, e2 );
3096 539291 : if ( tmp == 0 )
3097 307094 : tmp = sub( m1, m2 );
3098 :
3099 539291 : return tmp;
3100 : }
3101 :
3102 2113083076 : cmplx CL_scale_t( cmplx x, Word16 y )
3103 : {
3104 : cmplx result;
3105 2113083076 : result.re = Mpy_32_16_1( x.re, y );
3106 2113083076 : result.im = Mpy_32_16_1( x.im, y );
3107 : #ifdef WMOPS
3108 : multiCounter[currCounter].Mpy_32_16_1--;
3109 : multiCounter[currCounter].Mpy_32_16_1--;
3110 : multiCounter[currCounter].CL_scale++;
3111 : #endif
3112 2113083076 : return ( result );
3113 : }
3114 :
3115 7645716 : cmplx CL_dscale_t( cmplx x, Word16 y1, Word16 y2 )
3116 : {
3117 : cmplx result;
3118 7645716 : result.re = Mpy_32_16_1( x.re, y1 );
3119 7645716 : result.im = Mpy_32_16_1( x.im, y2 );
3120 : #ifdef WMOPS
3121 : multiCounter[currCounter].Mpy_32_16_1--;
3122 : multiCounter[currCounter].Mpy_32_16_1--;
3123 : multiCounter[currCounter].CL_dscale++;
3124 : #endif /* #ifdef WMOPS */
3125 7645716 : return ( result );
3126 : }
3127 :
3128 54041541 : cmplx CL_mult_32x16( cmplx input, cmplx_s coeff )
3129 : {
3130 : cmplx result;
3131 54041541 : result.re = L_sub( Mpy_32_16_1( input.re, coeff.re ), Mpy_32_16_1( input.im, coeff.im ) );
3132 54041541 : result.im = L_add( Mpy_32_16_1( input.re, coeff.im ), Mpy_32_16_1( input.im, coeff.re ) );
3133 : #ifdef WMOPS
3134 : multiCounter[currCounter].CL_multr_32x16++;
3135 : multiCounter[currCounter].Mpy_32_16_1--;
3136 : multiCounter[currCounter].Mpy_32_16_1--;
3137 : multiCounter[currCounter].Mpy_32_16_1--;
3138 : multiCounter[currCounter].Mpy_32_16_1--;
3139 : multiCounter[currCounter].L_sub--;
3140 : multiCounter[currCounter].L_add--;
3141 : #endif
3142 54041541 : return result;
3143 : }
|
