/src/aac/libFDK/src/fixpoint_math.cpp
Line | Count | Source (jump to first uncovered line) |
1 | | /* ----------------------------------------------------------------------------- |
2 | | Software License for The Fraunhofer FDK AAC Codec Library for Android |
3 | | |
4 | | © Copyright 1995 - 2020 Fraunhofer-Gesellschaft zur Förderung der angewandten |
5 | | Forschung e.V. All rights reserved. |
6 | | |
7 | | 1. INTRODUCTION |
8 | | The Fraunhofer FDK AAC Codec Library for Android ("FDK AAC Codec") is software |
9 | | that implements the MPEG Advanced Audio Coding ("AAC") encoding and decoding |
10 | | scheme for digital audio. This FDK AAC Codec software is intended to be used on |
11 | | a wide variety of Android devices. |
12 | | |
13 | | AAC's HE-AAC and HE-AAC v2 versions are regarded as today's most efficient |
14 | | general perceptual audio codecs. AAC-ELD is considered the best-performing |
15 | | full-bandwidth communications codec by independent studies and is widely |
16 | | deployed. AAC has been standardized by ISO and IEC as part of the MPEG |
17 | | specifications. |
18 | | |
19 | | Patent licenses for necessary patent claims for the FDK AAC Codec (including |
20 | | those of Fraunhofer) may be obtained through Via Licensing |
21 | | (www.vialicensing.com) or through the respective patent owners individually for |
22 | | the purpose of encoding or decoding bit streams in products that are compliant |
23 | | with the ISO/IEC MPEG audio standards. Please note that most manufacturers of |
24 | | Android devices already license these patent claims through Via Licensing or |
25 | | directly from the patent owners, and therefore FDK AAC Codec software may |
26 | | already be covered under those patent licenses when it is used for those |
27 | | licensed purposes only. |
28 | | |
29 | | Commercially-licensed AAC software libraries, including floating-point versions |
30 | | with enhanced sound quality, are also available from Fraunhofer. Users are |
31 | | encouraged to check the Fraunhofer website for additional applications |
32 | | information and documentation. |
33 | | |
34 | | 2. COPYRIGHT LICENSE |
35 | | |
36 | | Redistribution and use in source and binary forms, with or without modification, |
37 | | are permitted without payment of copyright license fees provided that you |
38 | | satisfy the following conditions: |
39 | | |
40 | | You must retain the complete text of this software license in redistributions of |
41 | | the FDK AAC Codec or your modifications thereto in source code form. |
42 | | |
43 | | You must retain the complete text of this software license in the documentation |
44 | | and/or other materials provided with redistributions of the FDK AAC Codec or |
45 | | your modifications thereto in binary form. You must make available free of |
46 | | charge copies of the complete source code of the FDK AAC Codec and your |
47 | | modifications thereto to recipients of copies in binary form. |
48 | | |
49 | | The name of Fraunhofer may not be used to endorse or promote products derived |
50 | | from this library without prior written permission. |
51 | | |
52 | | You may not charge copyright license fees for anyone to use, copy or distribute |
53 | | the FDK AAC Codec software or your modifications thereto. |
54 | | |
55 | | Your modified versions of the FDK AAC Codec must carry prominent notices stating |
56 | | that you changed the software and the date of any change. For modified versions |
57 | | of the FDK AAC Codec, the term "Fraunhofer FDK AAC Codec Library for Android" |
58 | | must be replaced by the term "Third-Party Modified Version of the Fraunhofer FDK |
59 | | AAC Codec Library for Android." |
60 | | |
61 | | 3. NO PATENT LICENSE |
62 | | |
63 | | NO EXPRESS OR IMPLIED LICENSES TO ANY PATENT CLAIMS, including without |
64 | | limitation the patents of Fraunhofer, ARE GRANTED BY THIS SOFTWARE LICENSE. |
65 | | Fraunhofer provides no warranty of patent non-infringement with respect to this |
66 | | software. |
67 | | |
68 | | You may use this FDK AAC Codec software or modifications thereto only for |
69 | | purposes that are authorized by appropriate patent licenses. |
70 | | |
71 | | 4. DISCLAIMER |
72 | | |
73 | | This FDK AAC Codec software is provided by Fraunhofer on behalf of the copyright |
74 | | holders and contributors "AS IS" and WITHOUT ANY EXPRESS OR IMPLIED WARRANTIES, |
75 | | including but not limited to the implied warranties of merchantability and |
76 | | fitness for a particular purpose. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR |
77 | | CONTRIBUTORS BE LIABLE for any direct, indirect, incidental, special, exemplary, |
78 | | or consequential damages, including but not limited to procurement of substitute |
79 | | goods or services; loss of use, data, or profits, or business interruption, |
80 | | however caused and on any theory of liability, whether in contract, strict |
81 | | liability, or tort (including negligence), arising in any way out of the use of |
82 | | this software, even if advised of the possibility of such damage. |
83 | | |
84 | | 5. CONTACT INFORMATION |
85 | | |
86 | | Fraunhofer Institute for Integrated Circuits IIS |
87 | | Attention: Audio and Multimedia Departments - FDK AAC LL |
88 | | Am Wolfsmantel 33 |
89 | | 91058 Erlangen, Germany |
90 | | |
91 | | www.iis.fraunhofer.de/amm |
92 | | amm-info@iis.fraunhofer.de |
93 | | ----------------------------------------------------------------------------- */ |
94 | | |
95 | | /******************* Library for basic calculation routines ******************** |
96 | | |
97 | | Author(s): M. Gayer |
98 | | |
99 | | Description: Fixed point specific mathematical functions |
100 | | |
101 | | *******************************************************************************/ |
102 | | |
103 | | #include "fixpoint_math.h" |
104 | | |
105 | | /* |
106 | | * Hardware specific implementations |
107 | | */ |
108 | | |
109 | | /* |
110 | | * Fallback implementations |
111 | | */ |
112 | | |
113 | | /***************************************************************************** |
114 | | functionname: LdDataVector |
115 | | *****************************************************************************/ |
116 | | LNK_SECTION_CODE_L1 |
117 | 0 | void LdDataVector(FIXP_DBL *srcVector, FIXP_DBL *destVector, INT n) { |
118 | 0 | INT i; |
119 | 0 | for (i = 0; i < n; i++) { |
120 | 0 | destVector[i] = fLog2(srcVector[i], 0); |
121 | 0 | } |
122 | 0 | } |
123 | | |
124 | | #define MAX_POW2_PRECISION 8 |
125 | | #ifndef SINETABLE_16BIT |
126 | | #define POW2_PRECISION MAX_POW2_PRECISION |
127 | | #else |
128 | 90.7M | #define POW2_PRECISION 5 |
129 | | #endif |
130 | | |
131 | | /* |
132 | | Taylor series coefficients of the function x^2. The first coefficient is |
133 | | ommited (equal to 1.0). |
134 | | |
135 | | pow2Coeff[i-1] = (1/i!) d^i(2^x)/dx^i, i=1..MAX_POW2_PRECISION |
136 | | To evaluate the taylor series around x = 0, the coefficients are: 1/!i * |
137 | | ln(2)^i |
138 | | */ |
139 | | #ifndef POW2COEFF_16BIT |
140 | | RAM_ALIGN |
141 | | LNK_SECTION_CONSTDATA_L1 |
142 | | static const FIXP_DBL pow2Coeff[MAX_POW2_PRECISION] = { |
143 | | FL2FXCONST_DBL(0.693147180559945309417232121458177), /* ln(2)^1 /1! */ |
144 | | FL2FXCONST_DBL(0.240226506959100712333551263163332), /* ln(2)^2 /2! */ |
145 | | FL2FXCONST_DBL(0.0555041086648215799531422637686218), /* ln(2)^3 /3! */ |
146 | | FL2FXCONST_DBL(0.00961812910762847716197907157365887), /* ln(2)^4 /4! */ |
147 | | FL2FXCONST_DBL(0.00133335581464284434234122219879962), /* ln(2)^5 /5! */ |
148 | | FL2FXCONST_DBL(1.54035303933816099544370973327423e-4), /* ln(2)^6 /6! */ |
149 | | FL2FXCONST_DBL(1.52527338040598402800254390120096e-5), /* ln(2)^7 /7! */ |
150 | | FL2FXCONST_DBL(1.32154867901443094884037582282884e-6) /* ln(2)^8 /8! */ |
151 | | }; |
152 | | #else |
153 | | RAM_ALIGN |
154 | | LNK_SECTION_CONSTDATA_L1 |
155 | | static const FIXP_SGL pow2Coeff[MAX_POW2_PRECISION] = { |
156 | | FL2FXCONST_SGL(0.693147180559945309417232121458177), /* ln(2)^1 /1! */ |
157 | | FL2FXCONST_SGL(0.240226506959100712333551263163332), /* ln(2)^2 /2! */ |
158 | | FL2FXCONST_SGL(0.0555041086648215799531422637686218), /* ln(2)^3 /3! */ |
159 | | FL2FXCONST_SGL(0.00961812910762847716197907157365887), /* ln(2)^4 /4! */ |
160 | | FL2FXCONST_SGL(0.00133335581464284434234122219879962), /* ln(2)^5 /5! */ |
161 | | FL2FXCONST_SGL(1.54035303933816099544370973327423e-4), /* ln(2)^6 /6! */ |
162 | | FL2FXCONST_SGL(1.52527338040598402800254390120096e-5), /* ln(2)^7 /7! */ |
163 | | FL2FXCONST_SGL(1.32154867901443094884037582282884e-6) /* ln(2)^8 /8! */ |
164 | | }; |
165 | | #endif |
166 | | |
167 | | /***************************************************************************** |
168 | | |
169 | | functionname: CalcInvLdData |
170 | | description: Delivers the inverse of function CalcLdData(). |
171 | | Delivers 2^(op*LD_DATA_SCALING) |
172 | | input: Input op is assumed to be fractional -1.0 < op < 1.0 |
173 | | output: For op == 0, the result is MAXVAL_DBL (almost 1.0). |
174 | | For negative input values the output should be treated as a |
175 | | positive fractional value. For positive input values the output should be |
176 | | treated as a positive integer value. This function does not output negative |
177 | | values. |
178 | | |
179 | | *****************************************************************************/ |
180 | | /* Date: 06-JULY-2012 Arthur Tritthart, IIS Fraunhofer Erlangen */ |
181 | | /* Version with 3 table lookup and 1 linear interpolations */ |
182 | | /* Algorithm: compute power of 2, argument x is in Q7.25 format */ |
183 | | /* result = 2^(x/64) */ |
184 | | /* We split exponent (x/64) into 5 components: */ |
185 | | /* integer part: represented by b31..b25 (exp) */ |
186 | | /* fractional part 1: represented by b24..b20 (lookup1) */ |
187 | | /* fractional part 2: represented by b19..b15 (lookup2) */ |
188 | | /* fractional part 3: represented by b14..b10 (lookup3) */ |
189 | | /* fractional part 4: represented by b09..b00 (frac) */ |
190 | | /* => result = (lookup1*lookup2*(lookup3+C1*frac)<<3)>>exp */ |
191 | | /* Due to the fact, that all lookup values contain a factor 0.5 */ |
192 | | /* the result has to be shifted by 3 to the right also. */ |
193 | | /* Table exp2_tab_long contains the log2 for 0 to 1.0 in steps */ |
194 | | /* of 1/32, table exp2w_tab_long the log2 for 0 to 1/32 in steps*/ |
195 | | /* of 1/1024, table exp2x_tab_long the log2 for 0 to 1/1024 in */ |
196 | | /* steps of 1/32768. Since the 2-logarithm of very very small */ |
197 | | /* negative value is rather linear, we can use interpolation. */ |
198 | | /* Limitations: */ |
199 | | /* For x <= 0, the result is fractional positive */ |
200 | | /* For x > 0, the result is integer in range 1...7FFF.FFFF */ |
201 | | /* For x < -31/64, we have to clear the result */ |
202 | | /* For x = 0, the result is ~1.0 (0x7FFF.FFFF) */ |
203 | | /* For x >= 31/64, the result is 0x7FFF.FFFF */ |
204 | | |
205 | | /* This table is used for lookup 2^x with */ |
206 | | /* x in range [0...1.0[ in steps of 1/32 */ |
207 | | LNK_SECTION_DATA_L1 |
208 | | const UINT exp2_tab_long[32] = { |
209 | | 0x40000000, 0x4166C34C, 0x42D561B4, 0x444C0740, 0x45CAE0F2, 0x47521CC6, |
210 | | 0x48E1E9BA, 0x4A7A77D4, 0x4C1BF829, 0x4DC69CDD, 0x4F7A9930, 0x51382182, |
211 | | 0x52FF6B55, 0x54D0AD5A, 0x56AC1F75, 0x5891FAC1, 0x5A82799A, 0x5C7DD7A4, |
212 | | 0x5E8451D0, 0x60962665, 0x62B39509, 0x64DCDEC3, 0x6712460B, 0x69540EC9, |
213 | | 0x6BA27E65, 0x6DFDDBCC, 0x70666F76, 0x72DC8374, 0x75606374, 0x77F25CCE, |
214 | | 0x7A92BE8B, 0x7D41D96E |
215 | | // 0x80000000 |
216 | | }; |
217 | | |
218 | | /* This table is used for lookup 2^x with */ |
219 | | /* x in range [0...1/32[ in steps of 1/1024 */ |
220 | | LNK_SECTION_DATA_L1 |
221 | | const UINT exp2w_tab_long[32] = { |
222 | | 0x40000000, 0x400B1818, 0x4016321B, 0x40214E0C, 0x402C6BE9, 0x40378BB4, |
223 | | 0x4042AD6D, 0x404DD113, 0x4058F6A8, 0x40641E2B, 0x406F479E, 0x407A7300, |
224 | | 0x4085A051, 0x4090CF92, 0x409C00C4, 0x40A733E6, 0x40B268FA, 0x40BD9FFF, |
225 | | 0x40C8D8F5, 0x40D413DD, 0x40DF50B8, 0x40EA8F86, 0x40F5D046, 0x410112FA, |
226 | | 0x410C57A2, 0x41179E3D, 0x4122E6CD, 0x412E3152, 0x41397DCC, 0x4144CC3B, |
227 | | 0x41501CA0, 0x415B6EFB, |
228 | | // 0x4166C34C, |
229 | | }; |
230 | | /* This table is used for lookup 2^x with */ |
231 | | /* x in range [0...1/1024[ in steps of 1/32768 */ |
232 | | LNK_SECTION_DATA_L1 |
233 | | const UINT exp2x_tab_long[32] = { |
234 | | 0x40000000, 0x400058B9, 0x4000B173, 0x40010A2D, 0x400162E8, 0x4001BBA3, |
235 | | 0x4002145F, 0x40026D1B, 0x4002C5D8, 0x40031E95, 0x40037752, 0x4003D011, |
236 | | 0x400428CF, 0x4004818E, 0x4004DA4E, 0x4005330E, 0x40058BCE, 0x4005E48F, |
237 | | 0x40063D51, 0x40069613, 0x4006EED5, 0x40074798, 0x4007A05B, 0x4007F91F, |
238 | | 0x400851E4, 0x4008AAA8, 0x4009036E, 0x40095C33, 0x4009B4FA, 0x400A0DC0, |
239 | | 0x400A6688, 0x400ABF4F, |
240 | | // 0x400B1818 |
241 | | }; |
242 | | |
243 | | /***************************************************************************** |
244 | | functionname: InitLdInt and CalcLdInt |
245 | | description: Create and access table with integer LdData (0 to |
246 | | LD_INT_TAB_LEN) |
247 | | *****************************************************************************/ |
248 | | #ifndef LD_INT_TAB_LEN |
249 | | #define LD_INT_TAB_LEN \ |
250 | 1.38M | 193 /* Default tab length. Lower value should be set in fix.h */ |
251 | | #endif |
252 | | |
253 | | #if (LD_INT_TAB_LEN <= 120) |
254 | | LNK_SECTION_CONSTDATA_L1 |
255 | | static const FIXP_DBL ldIntCoeff[] = { |
256 | | (FIXP_DBL)0x80000001, (FIXP_DBL)0x00000000, (FIXP_DBL)0x02000000, |
257 | | (FIXP_DBL)0x032b8034, (FIXP_DBL)0x04000000, (FIXP_DBL)0x04a4d3c2, |
258 | | (FIXP_DBL)0x052b8034, (FIXP_DBL)0x059d5da0, (FIXP_DBL)0x06000000, |
259 | | (FIXP_DBL)0x06570069, (FIXP_DBL)0x06a4d3c2, (FIXP_DBL)0x06eb3a9f, |
260 | | (FIXP_DBL)0x072b8034, (FIXP_DBL)0x0766a009, (FIXP_DBL)0x079d5da0, |
261 | | (FIXP_DBL)0x07d053f7, (FIXP_DBL)0x08000000, (FIXP_DBL)0x082cc7ee, |
262 | | (FIXP_DBL)0x08570069, (FIXP_DBL)0x087ef05b, (FIXP_DBL)0x08a4d3c2, |
263 | | (FIXP_DBL)0x08c8ddd4, (FIXP_DBL)0x08eb3a9f, (FIXP_DBL)0x090c1050, |
264 | | (FIXP_DBL)0x092b8034, (FIXP_DBL)0x0949a785, (FIXP_DBL)0x0966a009, |
265 | | (FIXP_DBL)0x0982809d, (FIXP_DBL)0x099d5da0, (FIXP_DBL)0x09b74949, |
266 | | (FIXP_DBL)0x09d053f7, (FIXP_DBL)0x09e88c6b, (FIXP_DBL)0x0a000000, |
267 | | (FIXP_DBL)0x0a16bad3, (FIXP_DBL)0x0a2cc7ee, (FIXP_DBL)0x0a423162, |
268 | | (FIXP_DBL)0x0a570069, (FIXP_DBL)0x0a6b3d79, (FIXP_DBL)0x0a7ef05b, |
269 | | (FIXP_DBL)0x0a92203d, (FIXP_DBL)0x0aa4d3c2, (FIXP_DBL)0x0ab7110e, |
270 | | (FIXP_DBL)0x0ac8ddd4, (FIXP_DBL)0x0ada3f60, (FIXP_DBL)0x0aeb3a9f, |
271 | | (FIXP_DBL)0x0afbd42b, (FIXP_DBL)0x0b0c1050, (FIXP_DBL)0x0b1bf312, |
272 | | (FIXP_DBL)0x0b2b8034, (FIXP_DBL)0x0b3abb40, (FIXP_DBL)0x0b49a785, |
273 | | (FIXP_DBL)0x0b584822, (FIXP_DBL)0x0b66a009, (FIXP_DBL)0x0b74b1fd, |
274 | | (FIXP_DBL)0x0b82809d, (FIXP_DBL)0x0b900e61, (FIXP_DBL)0x0b9d5da0, |
275 | | (FIXP_DBL)0x0baa708f, (FIXP_DBL)0x0bb74949, (FIXP_DBL)0x0bc3e9ca, |
276 | | (FIXP_DBL)0x0bd053f7, (FIXP_DBL)0x0bdc899b, (FIXP_DBL)0x0be88c6b, |
277 | | (FIXP_DBL)0x0bf45e09, (FIXP_DBL)0x0c000000, (FIXP_DBL)0x0c0b73cb, |
278 | | (FIXP_DBL)0x0c16bad3, (FIXP_DBL)0x0c21d671, (FIXP_DBL)0x0c2cc7ee, |
279 | | (FIXP_DBL)0x0c379085, (FIXP_DBL)0x0c423162, (FIXP_DBL)0x0c4caba8, |
280 | | (FIXP_DBL)0x0c570069, (FIXP_DBL)0x0c6130af, (FIXP_DBL)0x0c6b3d79, |
281 | | (FIXP_DBL)0x0c7527b9, (FIXP_DBL)0x0c7ef05b, (FIXP_DBL)0x0c88983f, |
282 | | (FIXP_DBL)0x0c92203d, (FIXP_DBL)0x0c9b8926, (FIXP_DBL)0x0ca4d3c2, |
283 | | (FIXP_DBL)0x0cae00d2, (FIXP_DBL)0x0cb7110e, (FIXP_DBL)0x0cc0052b, |
284 | | (FIXP_DBL)0x0cc8ddd4, (FIXP_DBL)0x0cd19bb0, (FIXP_DBL)0x0cda3f60, |
285 | | (FIXP_DBL)0x0ce2c97d, (FIXP_DBL)0x0ceb3a9f, (FIXP_DBL)0x0cf39355, |
286 | | (FIXP_DBL)0x0cfbd42b, (FIXP_DBL)0x0d03fda9, (FIXP_DBL)0x0d0c1050, |
287 | | (FIXP_DBL)0x0d140ca0, (FIXP_DBL)0x0d1bf312, (FIXP_DBL)0x0d23c41d, |
288 | | (FIXP_DBL)0x0d2b8034, (FIXP_DBL)0x0d3327c7, (FIXP_DBL)0x0d3abb40, |
289 | | (FIXP_DBL)0x0d423b08, (FIXP_DBL)0x0d49a785, (FIXP_DBL)0x0d510118, |
290 | | (FIXP_DBL)0x0d584822, (FIXP_DBL)0x0d5f7cff, (FIXP_DBL)0x0d66a009, |
291 | | (FIXP_DBL)0x0d6db197, (FIXP_DBL)0x0d74b1fd, (FIXP_DBL)0x0d7ba190, |
292 | | (FIXP_DBL)0x0d82809d, (FIXP_DBL)0x0d894f75, (FIXP_DBL)0x0d900e61, |
293 | | (FIXP_DBL)0x0d96bdad, (FIXP_DBL)0x0d9d5da0, (FIXP_DBL)0x0da3ee7f, |
294 | | (FIXP_DBL)0x0daa708f, (FIXP_DBL)0x0db0e412, (FIXP_DBL)0x0db74949, |
295 | | (FIXP_DBL)0x0dbda072, (FIXP_DBL)0x0dc3e9ca, (FIXP_DBL)0x0dca258e}; |
296 | | |
297 | | #elif (LD_INT_TAB_LEN <= 193) |
298 | | LNK_SECTION_CONSTDATA_L1 |
299 | | static const FIXP_DBL ldIntCoeff[] = { |
300 | | (FIXP_DBL)0x80000001, (FIXP_DBL)0x00000000, (FIXP_DBL)0x02000000, |
301 | | (FIXP_DBL)0x032b8034, (FIXP_DBL)0x04000000, (FIXP_DBL)0x04a4d3c2, |
302 | | (FIXP_DBL)0x052b8034, (FIXP_DBL)0x059d5da0, (FIXP_DBL)0x06000000, |
303 | | (FIXP_DBL)0x06570069, (FIXP_DBL)0x06a4d3c2, (FIXP_DBL)0x06eb3a9f, |
304 | | (FIXP_DBL)0x072b8034, (FIXP_DBL)0x0766a009, (FIXP_DBL)0x079d5da0, |
305 | | (FIXP_DBL)0x07d053f7, (FIXP_DBL)0x08000000, (FIXP_DBL)0x082cc7ee, |
306 | | (FIXP_DBL)0x08570069, (FIXP_DBL)0x087ef05b, (FIXP_DBL)0x08a4d3c2, |
307 | | (FIXP_DBL)0x08c8ddd4, (FIXP_DBL)0x08eb3a9f, (FIXP_DBL)0x090c1050, |
308 | | (FIXP_DBL)0x092b8034, (FIXP_DBL)0x0949a785, (FIXP_DBL)0x0966a009, |
309 | | (FIXP_DBL)0x0982809d, (FIXP_DBL)0x099d5da0, (FIXP_DBL)0x09b74949, |
310 | | (FIXP_DBL)0x09d053f7, (FIXP_DBL)0x09e88c6b, (FIXP_DBL)0x0a000000, |
311 | | (FIXP_DBL)0x0a16bad3, (FIXP_DBL)0x0a2cc7ee, (FIXP_DBL)0x0a423162, |
312 | | (FIXP_DBL)0x0a570069, (FIXP_DBL)0x0a6b3d79, (FIXP_DBL)0x0a7ef05b, |
313 | | (FIXP_DBL)0x0a92203d, (FIXP_DBL)0x0aa4d3c2, (FIXP_DBL)0x0ab7110e, |
314 | | (FIXP_DBL)0x0ac8ddd4, (FIXP_DBL)0x0ada3f60, (FIXP_DBL)0x0aeb3a9f, |
315 | | (FIXP_DBL)0x0afbd42b, (FIXP_DBL)0x0b0c1050, (FIXP_DBL)0x0b1bf312, |
316 | | (FIXP_DBL)0x0b2b8034, (FIXP_DBL)0x0b3abb40, (FIXP_DBL)0x0b49a785, |
317 | | (FIXP_DBL)0x0b584822, (FIXP_DBL)0x0b66a009, (FIXP_DBL)0x0b74b1fd, |
318 | | (FIXP_DBL)0x0b82809d, (FIXP_DBL)0x0b900e61, (FIXP_DBL)0x0b9d5da0, |
319 | | (FIXP_DBL)0x0baa708f, (FIXP_DBL)0x0bb74949, (FIXP_DBL)0x0bc3e9ca, |
320 | | (FIXP_DBL)0x0bd053f7, (FIXP_DBL)0x0bdc899b, (FIXP_DBL)0x0be88c6b, |
321 | | (FIXP_DBL)0x0bf45e09, (FIXP_DBL)0x0c000000, (FIXP_DBL)0x0c0b73cb, |
322 | | (FIXP_DBL)0x0c16bad3, (FIXP_DBL)0x0c21d671, (FIXP_DBL)0x0c2cc7ee, |
323 | | (FIXP_DBL)0x0c379085, (FIXP_DBL)0x0c423162, (FIXP_DBL)0x0c4caba8, |
324 | | (FIXP_DBL)0x0c570069, (FIXP_DBL)0x0c6130af, (FIXP_DBL)0x0c6b3d79, |
325 | | (FIXP_DBL)0x0c7527b9, (FIXP_DBL)0x0c7ef05b, (FIXP_DBL)0x0c88983f, |
326 | | (FIXP_DBL)0x0c92203d, (FIXP_DBL)0x0c9b8926, (FIXP_DBL)0x0ca4d3c2, |
327 | | (FIXP_DBL)0x0cae00d2, (FIXP_DBL)0x0cb7110e, (FIXP_DBL)0x0cc0052b, |
328 | | (FIXP_DBL)0x0cc8ddd4, (FIXP_DBL)0x0cd19bb0, (FIXP_DBL)0x0cda3f60, |
329 | | (FIXP_DBL)0x0ce2c97d, (FIXP_DBL)0x0ceb3a9f, (FIXP_DBL)0x0cf39355, |
330 | | (FIXP_DBL)0x0cfbd42b, (FIXP_DBL)0x0d03fda9, (FIXP_DBL)0x0d0c1050, |
331 | | (FIXP_DBL)0x0d140ca0, (FIXP_DBL)0x0d1bf312, (FIXP_DBL)0x0d23c41d, |
332 | | (FIXP_DBL)0x0d2b8034, (FIXP_DBL)0x0d3327c7, (FIXP_DBL)0x0d3abb40, |
333 | | (FIXP_DBL)0x0d423b08, (FIXP_DBL)0x0d49a785, (FIXP_DBL)0x0d510118, |
334 | | (FIXP_DBL)0x0d584822, (FIXP_DBL)0x0d5f7cff, (FIXP_DBL)0x0d66a009, |
335 | | (FIXP_DBL)0x0d6db197, (FIXP_DBL)0x0d74b1fd, (FIXP_DBL)0x0d7ba190, |
336 | | (FIXP_DBL)0x0d82809d, (FIXP_DBL)0x0d894f75, (FIXP_DBL)0x0d900e61, |
337 | | (FIXP_DBL)0x0d96bdad, (FIXP_DBL)0x0d9d5da0, (FIXP_DBL)0x0da3ee7f, |
338 | | (FIXP_DBL)0x0daa708f, (FIXP_DBL)0x0db0e412, (FIXP_DBL)0x0db74949, |
339 | | (FIXP_DBL)0x0dbda072, (FIXP_DBL)0x0dc3e9ca, (FIXP_DBL)0x0dca258e, |
340 | | (FIXP_DBL)0x0dd053f7, (FIXP_DBL)0x0dd6753e, (FIXP_DBL)0x0ddc899b, |
341 | | (FIXP_DBL)0x0de29143, (FIXP_DBL)0x0de88c6b, (FIXP_DBL)0x0dee7b47, |
342 | | (FIXP_DBL)0x0df45e09, (FIXP_DBL)0x0dfa34e1, (FIXP_DBL)0x0e000000, |
343 | | (FIXP_DBL)0x0e05bf94, (FIXP_DBL)0x0e0b73cb, (FIXP_DBL)0x0e111cd2, |
344 | | (FIXP_DBL)0x0e16bad3, (FIXP_DBL)0x0e1c4dfb, (FIXP_DBL)0x0e21d671, |
345 | | (FIXP_DBL)0x0e275460, (FIXP_DBL)0x0e2cc7ee, (FIXP_DBL)0x0e323143, |
346 | | (FIXP_DBL)0x0e379085, (FIXP_DBL)0x0e3ce5d8, (FIXP_DBL)0x0e423162, |
347 | | (FIXP_DBL)0x0e477346, (FIXP_DBL)0x0e4caba8, (FIXP_DBL)0x0e51daa8, |
348 | | (FIXP_DBL)0x0e570069, (FIXP_DBL)0x0e5c1d0b, (FIXP_DBL)0x0e6130af, |
349 | | (FIXP_DBL)0x0e663b74, (FIXP_DBL)0x0e6b3d79, (FIXP_DBL)0x0e7036db, |
350 | | (FIXP_DBL)0x0e7527b9, (FIXP_DBL)0x0e7a1030, (FIXP_DBL)0x0e7ef05b, |
351 | | (FIXP_DBL)0x0e83c857, (FIXP_DBL)0x0e88983f, (FIXP_DBL)0x0e8d602e, |
352 | | (FIXP_DBL)0x0e92203d, (FIXP_DBL)0x0e96d888, (FIXP_DBL)0x0e9b8926, |
353 | | (FIXP_DBL)0x0ea03232, (FIXP_DBL)0x0ea4d3c2, (FIXP_DBL)0x0ea96df0, |
354 | | (FIXP_DBL)0x0eae00d2, (FIXP_DBL)0x0eb28c7f, (FIXP_DBL)0x0eb7110e, |
355 | | (FIXP_DBL)0x0ebb8e96, (FIXP_DBL)0x0ec0052b, (FIXP_DBL)0x0ec474e4, |
356 | | (FIXP_DBL)0x0ec8ddd4, (FIXP_DBL)0x0ecd4012, (FIXP_DBL)0x0ed19bb0, |
357 | | (FIXP_DBL)0x0ed5f0c4, (FIXP_DBL)0x0eda3f60, (FIXP_DBL)0x0ede8797, |
358 | | (FIXP_DBL)0x0ee2c97d, (FIXP_DBL)0x0ee70525, (FIXP_DBL)0x0eeb3a9f, |
359 | | (FIXP_DBL)0x0eef69ff, (FIXP_DBL)0x0ef39355, (FIXP_DBL)0x0ef7b6b4, |
360 | | (FIXP_DBL)0x0efbd42b, (FIXP_DBL)0x0effebcd, (FIXP_DBL)0x0f03fda9, |
361 | | (FIXP_DBL)0x0f0809cf, (FIXP_DBL)0x0f0c1050, (FIXP_DBL)0x0f10113b, |
362 | | (FIXP_DBL)0x0f140ca0, (FIXP_DBL)0x0f18028d, (FIXP_DBL)0x0f1bf312, |
363 | | (FIXP_DBL)0x0f1fde3d, (FIXP_DBL)0x0f23c41d, (FIXP_DBL)0x0f27a4c0, |
364 | | (FIXP_DBL)0x0f2b8034}; |
365 | | |
366 | | #else |
367 | | #error "ldInt table size too small" |
368 | | |
369 | | #endif |
370 | | |
371 | | LNK_SECTION_INITCODE |
372 | 0 | void InitLdInt() { /* nothing to do! Use preinitialized logarithm table */ |
373 | 0 | } |
374 | | |
375 | | #if (LD_INT_TAB_LEN != 0) |
376 | | |
377 | | LNK_SECTION_CODE_L1 |
378 | 1.38M | FIXP_DBL CalcLdInt(INT i) { |
379 | | /* calculates ld(op)/LD_DATA_SCALING */ |
380 | | /* op is assumed to be an integer value between 1 and LD_INT_TAB_LEN */ |
381 | | |
382 | 1.38M | FDK_ASSERT((LD_INT_TAB_LEN > 0) && |
383 | 1.38M | ((FIXP_DBL)ldIntCoeff[0] == |
384 | 1.38M | (FIXP_DBL)0x80000001)); /* tab has to be initialized */ |
385 | | |
386 | 1.38M | if ((i > 0) && (i < LD_INT_TAB_LEN)) |
387 | 1.38M | return ldIntCoeff[i]; |
388 | 0 | else { |
389 | 0 | return (0); |
390 | 0 | } |
391 | 1.38M | } |
392 | | #endif /* (LD_INT_TAB_LEN!=0) */ |
393 | | |
394 | | #if !defined(FUNCTION_schur_div) |
395 | | /***************************************************************************** |
396 | | |
397 | | functionname: schur_div |
398 | | description: delivers op1/op2 with op3-bit accuracy |
399 | | |
400 | | *****************************************************************************/ |
401 | | |
402 | | FIXP_DBL schur_div(FIXP_DBL num, FIXP_DBL denum, INT count) { |
403 | | INT L_num = (LONG)num >> 1; |
404 | | INT L_denum = (LONG)denum >> 1; |
405 | | INT div = 0; |
406 | | INT k = count; |
407 | | |
408 | | FDK_ASSERT(num >= (FIXP_DBL)0); |
409 | | FDK_ASSERT(denum > (FIXP_DBL)0); |
410 | | FDK_ASSERT(num <= denum); |
411 | | |
412 | | if (L_num != 0) |
413 | | while (--k) { |
414 | | div <<= 1; |
415 | | L_num <<= 1; |
416 | | if (L_num >= L_denum) { |
417 | | L_num -= L_denum; |
418 | | div++; |
419 | | } |
420 | | } |
421 | | return (FIXP_DBL)(div << (DFRACT_BITS - count)); |
422 | | } |
423 | | |
424 | | #endif /* !defined(FUNCTION_schur_div) */ |
425 | | |
426 | | #ifndef FUNCTION_fMultNorm |
427 | 13.1M | FIXP_DBL fMultNorm(FIXP_DBL f1, FIXP_DBL f2, INT *result_e) { |
428 | 13.1M | INT product = 0; |
429 | 13.1M | INT norm_f1, norm_f2; |
430 | | |
431 | 13.1M | if ((f1 == (FIXP_DBL)0) || (f2 == (FIXP_DBL)0)) { |
432 | 1.50M | *result_e = 0; |
433 | 1.50M | return (FIXP_DBL)0; |
434 | 1.50M | } |
435 | 11.6M | norm_f1 = CountLeadingBits(f1); |
436 | 11.6M | f1 = f1 << norm_f1; |
437 | 11.6M | norm_f2 = CountLeadingBits(f2); |
438 | 11.6M | f2 = f2 << norm_f2; |
439 | | |
440 | 11.6M | if ((f1 == (FIXP_DBL)MINVAL_DBL) && (f2 == (FIXP_DBL)MINVAL_DBL)) { |
441 | 0 | product = -((FIXP_DBL)MINVAL_DBL >> 1); |
442 | 0 | *result_e = -(norm_f1 + norm_f2 - 1); |
443 | 11.6M | } else { |
444 | 11.6M | product = fMult(f1, f2); |
445 | 11.6M | *result_e = -(norm_f1 + norm_f2); |
446 | 11.6M | } |
447 | | |
448 | 11.6M | return (FIXP_DBL)product; |
449 | 13.1M | } |
450 | | #endif |
451 | | |
452 | | #ifndef FUNCTION_fDivNorm |
453 | 350M | FIXP_DBL fDivNorm(FIXP_DBL L_num, FIXP_DBL L_denum, INT *result_e) { |
454 | 350M | FIXP_DBL div; |
455 | 350M | INT norm_num, norm_den; |
456 | | |
457 | 350M | FDK_ASSERT(L_num >= (FIXP_DBL)0); |
458 | 350M | FDK_ASSERT(L_denum > (FIXP_DBL)0); |
459 | | |
460 | 350M | if (L_num == (FIXP_DBL)0) { |
461 | 801k | *result_e = 0; |
462 | 801k | return ((FIXP_DBL)0); |
463 | 801k | } |
464 | | |
465 | 350M | norm_num = CountLeadingBits(L_num); |
466 | 350M | L_num = L_num << norm_num; |
467 | 350M | L_num = L_num >> 1; |
468 | 350M | *result_e = -norm_num + 1; |
469 | | |
470 | 350M | norm_den = CountLeadingBits(L_denum); |
471 | 350M | L_denum = L_denum << norm_den; |
472 | 350M | *result_e -= -norm_den; |
473 | | |
474 | 350M | div = schur_div(L_num, L_denum, FRACT_BITS); |
475 | | |
476 | 350M | return div; |
477 | 350M | } |
478 | | #endif /* !FUNCTION_fDivNorm */ |
479 | | |
480 | | #ifndef FUNCTION_fDivNorm |
481 | 268M | FIXP_DBL fDivNorm(FIXP_DBL num, FIXP_DBL denom) { |
482 | 268M | INT e; |
483 | 268M | FIXP_DBL res; |
484 | | |
485 | 268M | FDK_ASSERT(denom >= num); |
486 | | |
487 | 268M | res = fDivNorm(num, denom, &e); |
488 | | |
489 | | /* Avoid overflow since we must output a value with exponent 0 |
490 | | there is no other choice than saturating to almost 1.0f */ |
491 | 268M | if (res == (FIXP_DBL)(1 << (DFRACT_BITS - 2)) && e == 1) { |
492 | 262k | res = (FIXP_DBL)MAXVAL_DBL; |
493 | 268M | } else { |
494 | 268M | res = scaleValue(res, e); |
495 | 268M | } |
496 | | |
497 | 268M | return res; |
498 | 268M | } |
499 | | #endif /* !FUNCTION_fDivNorm */ |
500 | | |
501 | | #ifndef FUNCTION_fDivNormSigned |
502 | 0 | FIXP_DBL fDivNormSigned(FIXP_DBL num, FIXP_DBL denom) { |
503 | 0 | INT e; |
504 | 0 | FIXP_DBL res; |
505 | 0 | int sign; |
506 | |
|
507 | 0 | if (denom == (FIXP_DBL)0) { |
508 | 0 | return (FIXP_DBL)MAXVAL_DBL; |
509 | 0 | } |
510 | | |
511 | 0 | sign = ((num >= (FIXP_DBL)0) != (denom >= (FIXP_DBL)0)); |
512 | 0 | res = fDivNormSigned(num, denom, &e); |
513 | | |
514 | | /* Saturate since we must output a value with exponent 0 */ |
515 | 0 | if ((e > 0) && (fAbs(res) >= FL2FXCONST_DBL(0.5))) { |
516 | 0 | if (sign) { |
517 | 0 | res = (FIXP_DBL)MINVAL_DBL; |
518 | 0 | } else { |
519 | 0 | res = (FIXP_DBL)MAXVAL_DBL; |
520 | 0 | } |
521 | 0 | } else { |
522 | 0 | res = scaleValue(res, e); |
523 | 0 | } |
524 | |
|
525 | 0 | return res; |
526 | 0 | } |
527 | 313k | FIXP_DBL fDivNormSigned(FIXP_DBL L_num, FIXP_DBL L_denum, INT *result_e) { |
528 | 313k | FIXP_DBL div; |
529 | 313k | INT norm_num, norm_den; |
530 | 313k | int sign; |
531 | | |
532 | 313k | sign = ((L_num >= (FIXP_DBL)0) != (L_denum >= (FIXP_DBL)0)); |
533 | | |
534 | 313k | if (L_num == (FIXP_DBL)0) { |
535 | 5.46k | *result_e = 0; |
536 | 5.46k | return ((FIXP_DBL)0); |
537 | 5.46k | } |
538 | 308k | if (L_denum == (FIXP_DBL)0) { |
539 | 0 | *result_e = 14; |
540 | 0 | return ((FIXP_DBL)MAXVAL_DBL); |
541 | 0 | } |
542 | | |
543 | 308k | norm_num = CountLeadingBits(L_num); |
544 | 308k | L_num = L_num << norm_num; |
545 | 308k | L_num = L_num >> 2; |
546 | 308k | L_num = fAbs(L_num); |
547 | 308k | *result_e = -norm_num + 1; |
548 | | |
549 | 308k | norm_den = CountLeadingBits(L_denum); |
550 | 308k | L_denum = L_denum << norm_den; |
551 | 308k | L_denum = L_denum >> 1; |
552 | 308k | L_denum = fAbs(L_denum); |
553 | 308k | *result_e -= -norm_den; |
554 | | |
555 | 308k | div = schur_div(L_num, L_denum, FRACT_BITS); |
556 | | |
557 | 308k | if (sign) { |
558 | 77.4k | div = -div; |
559 | 77.4k | } |
560 | | |
561 | 308k | return div; |
562 | 308k | } |
563 | | #endif /* FUNCTION_fDivNormSigned */ |
564 | | |
565 | | #ifndef FUNCTION_fDivNormHighPrec |
566 | 406k | FIXP_DBL fDivNormHighPrec(FIXP_DBL num, FIXP_DBL denom, INT *result_e) { |
567 | 406k | FIXP_DBL div; |
568 | 406k | INT norm_num, norm_den; |
569 | | |
570 | 406k | FDK_ASSERT(num >= (FIXP_DBL)0); |
571 | 406k | FDK_ASSERT(denom > (FIXP_DBL)0); |
572 | | |
573 | 406k | if (num == (FIXP_DBL)0) { |
574 | 0 | *result_e = 0; |
575 | 0 | return ((FIXP_DBL)0); |
576 | 0 | } |
577 | | |
578 | 406k | norm_num = CountLeadingBits(num); |
579 | 406k | num = num << norm_num; |
580 | 406k | num = num >> 1; |
581 | 406k | *result_e = -norm_num + 1; |
582 | | |
583 | 406k | norm_den = CountLeadingBits(denom); |
584 | 406k | denom = denom << norm_den; |
585 | 406k | *result_e -= -norm_den; |
586 | | |
587 | 406k | div = schur_div(num, denom, 31); |
588 | 406k | return div; |
589 | 406k | } |
590 | | #endif /* !FUNCTION_fDivNormHighPrec */ |
591 | | |
592 | | #ifndef FUNCTION_fPow |
593 | 15.1M | FIXP_DBL f2Pow(const FIXP_DBL exp_m, const INT exp_e, INT *result_e) { |
594 | 15.1M | FIXP_DBL frac_part, result_m; |
595 | 15.1M | INT int_part; |
596 | | |
597 | 15.1M | if (exp_e > 0) { |
598 | 14.0M | INT exp_bits = DFRACT_BITS - 1 - exp_e; |
599 | 14.0M | int_part = exp_m >> exp_bits; |
600 | 14.0M | frac_part = exp_m - (FIXP_DBL)(int_part << exp_bits); |
601 | 14.0M | frac_part = frac_part << exp_e; |
602 | 14.0M | } else { |
603 | 1.09M | int_part = 0; |
604 | 1.09M | frac_part = exp_m >> -exp_e; |
605 | 1.09M | } |
606 | | |
607 | | /* Best accuracy is around 0, so try to get there with the fractional part. */ |
608 | 15.1M | if (frac_part > FL2FXCONST_DBL(0.5f)) { |
609 | 7.89M | int_part = int_part + 1; |
610 | 7.89M | frac_part = frac_part + FL2FXCONST_DBL(-1.0f); |
611 | 7.89M | } |
612 | 15.1M | if (frac_part < FL2FXCONST_DBL(-0.5f)) { |
613 | 4.05k | int_part = int_part - 1; |
614 | 4.05k | frac_part = -(FL2FXCONST_DBL(-1.0f) - frac_part); |
615 | 4.05k | } |
616 | | |
617 | | /* "+ 1" compensates fMultAddDiv2() of the polynomial evaluation below. */ |
618 | 15.1M | *result_e = int_part + 1; |
619 | | |
620 | | /* Evaluate taylor polynomial which approximates 2^x */ |
621 | 15.1M | { |
622 | 15.1M | FIXP_DBL p; |
623 | | |
624 | | /* result_m ~= 2^frac_part */ |
625 | 15.1M | p = frac_part; |
626 | | /* First taylor series coefficient a_0 = 1.0, scaled by 0.5 due to |
627 | | * fMultDiv2(). */ |
628 | 15.1M | result_m = FL2FXCONST_DBL(1.0f / 2.0f); |
629 | 90.7M | for (INT i = 0; i < POW2_PRECISION; i++) { |
630 | | /* next taylor series term: a_i * x^i, x=0 */ |
631 | 75.6M | result_m = fMultAddDiv2(result_m, pow2Coeff[i], p); |
632 | 75.6M | p = fMult(p, frac_part); |
633 | 75.6M | } |
634 | 15.1M | } |
635 | 15.1M | return result_m; |
636 | 15.1M | } |
637 | | |
638 | 0 | FIXP_DBL f2Pow(const FIXP_DBL exp_m, const INT exp_e) { |
639 | 0 | FIXP_DBL result_m; |
640 | 0 | INT result_e; |
641 | |
|
642 | 0 | result_m = f2Pow(exp_m, exp_e, &result_e); |
643 | 0 | result_e = fixMin(DFRACT_BITS - 1, fixMax(-(DFRACT_BITS - 1), result_e)); |
644 | |
|
645 | 0 | return scaleValue(result_m, result_e); |
646 | 0 | } |
647 | | |
648 | | FIXP_DBL fPow(FIXP_DBL base_m, INT base_e, FIXP_DBL exp_m, INT exp_e, |
649 | 806k | INT *result_e) { |
650 | 806k | INT ans_lg2_e, baselg2_e; |
651 | 806k | FIXP_DBL base_lg2, ans_lg2, result; |
652 | | |
653 | 806k | if (base_m <= (FIXP_DBL)0) { |
654 | 0 | result = (FIXP_DBL)0; |
655 | 0 | *result_e = 0; |
656 | 0 | return result; |
657 | 0 | } |
658 | | |
659 | | /* Calc log2 of base */ |
660 | 806k | base_lg2 = fLog2(base_m, base_e, &baselg2_e); |
661 | | |
662 | | /* Prepare exp */ |
663 | 806k | { |
664 | 806k | INT leadingBits; |
665 | | |
666 | 806k | leadingBits = CountLeadingBits(fAbs(exp_m)); |
667 | 806k | exp_m = exp_m << leadingBits; |
668 | 806k | exp_e -= leadingBits; |
669 | 806k | } |
670 | | |
671 | | /* Calc base pow exp */ |
672 | 806k | ans_lg2 = fMult(base_lg2, exp_m); |
673 | 806k | ans_lg2_e = exp_e + baselg2_e; |
674 | | |
675 | | /* Calc antilog */ |
676 | 806k | result = f2Pow(ans_lg2, ans_lg2_e, result_e); |
677 | | |
678 | 806k | return result; |
679 | 806k | } |
680 | | |
681 | | FIXP_DBL fLdPow(FIXP_DBL baseLd_m, INT baseLd_e, FIXP_DBL exp_m, INT exp_e, |
682 | 494k | INT *result_e) { |
683 | 494k | INT ans_lg2_e; |
684 | 494k | FIXP_DBL ans_lg2, result; |
685 | | |
686 | | /* Prepare exp */ |
687 | 494k | { |
688 | 494k | INT leadingBits; |
689 | | |
690 | 494k | leadingBits = CountLeadingBits(fAbs(exp_m)); |
691 | 494k | exp_m = exp_m << leadingBits; |
692 | 494k | exp_e -= leadingBits; |
693 | 494k | } |
694 | | |
695 | | /* Calc base pow exp */ |
696 | 494k | ans_lg2 = fMult(baseLd_m, exp_m); |
697 | 494k | ans_lg2_e = exp_e + baseLd_e; |
698 | | |
699 | | /* Calc antilog */ |
700 | 494k | result = f2Pow(ans_lg2, ans_lg2_e, result_e); |
701 | | |
702 | 494k | return result; |
703 | 494k | } |
704 | | |
705 | 0 | FIXP_DBL fLdPow(FIXP_DBL baseLd_m, INT baseLd_e, FIXP_DBL exp_m, INT exp_e) { |
706 | 0 | FIXP_DBL result_m; |
707 | 0 | int result_e; |
708 | |
|
709 | 0 | result_m = fLdPow(baseLd_m, baseLd_e, exp_m, exp_e, &result_e); |
710 | |
|
711 | 0 | return SATURATE_SHIFT(result_m, -result_e, DFRACT_BITS); |
712 | 0 | } |
713 | | |
714 | 4.00k | FIXP_DBL fPowInt(FIXP_DBL base_m, INT base_e, INT exp, INT *pResult_e) { |
715 | 4.00k | FIXP_DBL result; |
716 | | |
717 | 4.00k | if (exp != 0) { |
718 | 914 | INT result_e = 0; |
719 | | |
720 | 914 | if (base_m != (FIXP_DBL)0) { |
721 | 914 | { |
722 | 914 | INT leadingBits; |
723 | 914 | leadingBits = CountLeadingBits(base_m); |
724 | 914 | base_m <<= leadingBits; |
725 | 914 | base_e -= leadingBits; |
726 | 914 | } |
727 | | |
728 | 914 | result = base_m; |
729 | | |
730 | 914 | { |
731 | 914 | int i; |
732 | 5.30k | for (i = 1; i < fAbs(exp); i++) { |
733 | 4.38k | result = fMult(result, base_m); |
734 | 4.38k | } |
735 | 914 | } |
736 | | |
737 | 914 | if (exp < 0) { |
738 | | /* 1.0 / ans */ |
739 | 0 | result = fDivNorm(FL2FXCONST_DBL(0.5f), result, &result_e); |
740 | 0 | result_e++; |
741 | 914 | } else { |
742 | 914 | int ansScale = CountLeadingBits(result); |
743 | 914 | result <<= ansScale; |
744 | 914 | result_e -= ansScale; |
745 | 914 | } |
746 | | |
747 | 914 | result_e += exp * base_e; |
748 | | |
749 | 914 | } else { |
750 | 0 | result = (FIXP_DBL)0; |
751 | 0 | } |
752 | 914 | *pResult_e = result_e; |
753 | 3.09k | } else { |
754 | 3.09k | result = FL2FXCONST_DBL(0.5f); |
755 | 3.09k | *pResult_e = 1; |
756 | 3.09k | } |
757 | | |
758 | 4.00k | return result; |
759 | 4.00k | } |
760 | | #endif /* FUNCTION_fPow */ |
761 | | |
762 | | #ifndef FUNCTION_fLog2 |
763 | 5.15M | FIXP_DBL CalcLog2(FIXP_DBL base_m, INT base_e, INT *result_e) { |
764 | 5.15M | return fLog2(base_m, base_e, result_e); |
765 | 5.15M | } |
766 | | #endif /* FUNCTION_fLog2 */ |
767 | | |
768 | 0 | INT fixp_floorToInt(FIXP_DBL f_inp, INT sf) { |
769 | 0 | FDK_ASSERT(sf >= 0); |
770 | 0 | INT floorInt = (INT)(f_inp >> ((DFRACT_BITS - 1) - sf)); |
771 | 0 | return floorInt; |
772 | 0 | } |
773 | | |
774 | 0 | FIXP_DBL fixp_floor(FIXP_DBL f_inp, INT sf) { |
775 | 0 | FDK_ASSERT(sf >= 0); |
776 | 0 | INT floorInt = fixp_floorToInt(f_inp, sf); |
777 | 0 | FIXP_DBL f_floor = (FIXP_DBL)(floorInt << ((DFRACT_BITS - 1) - sf)); |
778 | 0 | return f_floor; |
779 | 0 | } |
780 | | |
781 | | INT fixp_ceilToInt(FIXP_DBL f_inp, INT sf) // sf mantissaBits left of dot |
782 | 0 | { |
783 | 0 | FDK_ASSERT(sf >= 0); |
784 | 0 | INT sx = (DFRACT_BITS - 1) - sf; // sx mantissaBits right of dot |
785 | 0 | INT inpINT = (INT)f_inp; |
786 | |
|
787 | 0 | INT mask = (0x1 << sx) - 1; |
788 | 0 | INT ceilInt = (INT)(f_inp >> sx); |
789 | |
|
790 | 0 | if (inpINT & mask) { |
791 | 0 | ceilInt++; // increment only, if there is at least one set mantissaBit |
792 | | // right of dot [in inpINT] |
793 | 0 | } |
794 | |
|
795 | 0 | return ceilInt; |
796 | 0 | } |
797 | | |
798 | 0 | FIXP_DBL fixp_ceil(FIXP_DBL f_inp, INT sf) { |
799 | 0 | FDK_ASSERT(sf >= 0); |
800 | 0 | INT sx = (DFRACT_BITS - 1) - sf; |
801 | 0 | INT ceilInt = fixp_ceilToInt(f_inp, sf); |
802 | 0 | ULONG mask = (ULONG)0x1 << (DFRACT_BITS - 1); // 0x80000000 |
803 | 0 | ceilInt = ceilInt |
804 | 0 | << sx; // no fract warn bec. shift into saturation done on int side |
805 | |
|
806 | 0 | if ((f_inp > FL2FXCONST_DBL(0.0f)) && (ceilInt & mask)) { |
807 | 0 | --ceilInt; |
808 | 0 | } |
809 | 0 | FIXP_DBL f_ceil = (FIXP_DBL)ceilInt; |
810 | |
|
811 | 0 | return f_ceil; |
812 | 0 | } |
813 | | |
814 | | /***************************************************************************** |
815 | | fixp_truncateToInt() |
816 | | Just remove the fractional part which is located right of decimal point |
817 | | Same as which is done when a float is casted to (INT) : |
818 | | result_INTtype = (INT)b_floatTypeInput; |
819 | | |
820 | | returns INT |
821 | | *****************************************************************************/ |
822 | | INT fixp_truncateToInt(FIXP_DBL f_inp, INT sf) // sf mantissaBits left of dot |
823 | | // (without sign) e.g. at width |
824 | | // 32 this would be [sign]7. |
825 | | // supposed sf equals 8. |
826 | 0 | { |
827 | 0 | FDK_ASSERT(sf >= 0); |
828 | 0 | INT sx = (DFRACT_BITS - 1) - sf; // sx mantissaBits right of dot |
829 | | // at width 32 this would be .24 |
830 | | // supposed sf equals 8. |
831 | 0 | INT fbaccu = (INT)f_inp; |
832 | 0 | INT mask = (0x1 << sx); |
833 | |
|
834 | 0 | if ((fbaccu < 0) && (fbaccu & (mask - 1))) { |
835 | 0 | fbaccu = fbaccu + mask; |
836 | 0 | } |
837 | |
|
838 | 0 | fbaccu = fbaccu >> sx; |
839 | 0 | return fbaccu; |
840 | 0 | } |
841 | | |
842 | | /***************************************************************************** |
843 | | fixp_truncate() |
844 | | Just remove the fractional part which is located right of decimal point |
845 | | |
846 | | returns FIXP_DBL |
847 | | *****************************************************************************/ |
848 | 0 | FIXP_DBL fixp_truncate(FIXP_DBL f_inp, INT sf) { |
849 | 0 | FDK_ASSERT(sf >= 0); |
850 | 0 | INT truncateInt = fixp_truncateToInt(f_inp, sf); |
851 | 0 | FIXP_DBL f_truncate = (FIXP_DBL)(truncateInt << ((DFRACT_BITS - 1) - sf)); |
852 | 0 | return f_truncate; |
853 | 0 | } |
854 | | |
855 | | /***************************************************************************** |
856 | | fixp_roundToInt() |
857 | | round [typical rounding] |
858 | | |
859 | | See fct roundRef() [which is the reference] |
860 | | returns INT |
861 | | *****************************************************************************/ |
862 | 0 | INT fixp_roundToInt(FIXP_DBL f_inp, INT sf) { |
863 | 0 | FDK_ASSERT(sf >= 0); |
864 | 0 | INT sx = DFRACT_BITS - 1 - sf; |
865 | 0 | INT inp = (INT)f_inp; |
866 | 0 | INT mask1 = (0x1 << (sx - 1)); |
867 | 0 | INT mask2 = (0x1 << (sx)) - 1; |
868 | 0 | INT mask3 = 0x7FFFFFFF; |
869 | 0 | INT iam = inp & mask2; |
870 | 0 | INT rnd; |
871 | |
|
872 | 0 | if ((inp < 0) && !(iam == mask1)) |
873 | 0 | rnd = inp + mask1; |
874 | 0 | else if ((inp > 0) && !(inp == mask3)) |
875 | 0 | rnd = inp + mask1; |
876 | 0 | else |
877 | 0 | rnd = inp; |
878 | |
|
879 | 0 | rnd = rnd >> sx; |
880 | |
|
881 | 0 | if (inp == mask3) rnd++; |
882 | |
|
883 | 0 | return rnd; |
884 | 0 | } |
885 | | |
886 | | /***************************************************************************** |
887 | | fixp_round() |
888 | | round [typical rounding] |
889 | | |
890 | | See fct roundRef() [which is the reference] |
891 | | returns FIXP_DBL |
892 | | *****************************************************************************/ |
893 | 0 | FIXP_DBL fixp_round(FIXP_DBL f_inp, INT sf) { |
894 | 0 | FDK_ASSERT(sf >= 0); |
895 | 0 | INT sx = DFRACT_BITS - 1 - sf; |
896 | 0 | INT r = fixp_roundToInt(f_inp, sf); |
897 | 0 | ULONG mask = (ULONG)0x1 << (DFRACT_BITS - 1); // 0x80000000 |
898 | 0 | r = r << sx; |
899 | |
|
900 | 0 | if ((f_inp > FL2FXCONST_DBL(0.0f)) && (r & mask)) { |
901 | 0 | --r; |
902 | 0 | } |
903 | |
|
904 | 0 | FIXP_DBL f_round = (FIXP_DBL)r; |
905 | 0 | return f_round; |
906 | 0 | } |