1 /* ----------------------------------------------------------------------------- 2 Software License for The Fraunhofer FDK AAC Codec Library for Android 3 4 Copyright 1995 - 2018 Fraunhofer-Gesellschaft zur Frderung 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 (at) 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 void LdDataVector(FIXP_DBL *srcVector, FIXP_DBL *destVector, INT n) { 118 INT i; 119 for (i = 0; i < n; i++) { 120 destVector[i] = fLog2(srcVector[i], 0); 121 } 122 } 123 124 #define MAX_POW2_PRECISION 8 125 #ifndef SINETABLE_16BIT 126 #define POW2_PRECISION MAX_POW2_PRECISION 127 #else 128 #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 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 void InitLdInt() { /* nothing to do! Use preinitialized logarithm table */ 373 } 374 375 #if (LD_INT_TAB_LEN != 0) 376 377 LNK_SECTION_CODE_L1 378 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 FDK_ASSERT((LD_INT_TAB_LEN > 0) && 383 ((FIXP_DBL)ldIntCoeff[0] == 384 (FIXP_DBL)0x80000001)); /* tab has to be initialized */ 385 386 if ((i > 0) && (i < LD_INT_TAB_LEN)) 387 return ldIntCoeff[i]; 388 else { 389 return (0); 390 } 391 } 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 FIXP_DBL fMultNorm(FIXP_DBL f1, FIXP_DBL f2, INT *result_e) { 428 INT product = 0; 429 INT norm_f1, norm_f2; 430 431 if ((f1 == (FIXP_DBL)0) || (f2 == (FIXP_DBL)0)) { 432 *result_e = 0; 433 return (FIXP_DBL)0; 434 } 435 norm_f1 = CountLeadingBits(f1); 436 f1 = f1 << norm_f1; 437 norm_f2 = CountLeadingBits(f2); 438 f2 = f2 << norm_f2; 439 440 if ((f1 == (FIXP_DBL)MINVAL_DBL) && (f2 == (FIXP_DBL)MINVAL_DBL)) { 441 product = -((FIXP_DBL)MINVAL_DBL >> 1); 442 *result_e = -(norm_f1 + norm_f2 - 1); 443 } else { 444 product = fMult(f1, f2); 445 *result_e = -(norm_f1 + norm_f2); 446 } 447 448 return (FIXP_DBL)product; 449 } 450 #endif 451 452 #ifndef FUNCTION_fDivNorm 453 FIXP_DBL fDivNorm(FIXP_DBL L_num, FIXP_DBL L_denum, INT *result_e) { 454 FIXP_DBL div; 455 INT norm_num, norm_den; 456 457 FDK_ASSERT(L_num >= (FIXP_DBL)0); 458 FDK_ASSERT(L_denum > (FIXP_DBL)0); 459 460 if (L_num == (FIXP_DBL)0) { 461 *result_e = 0; 462 return ((FIXP_DBL)0); 463 } 464 465 norm_num = CountLeadingBits(L_num); 466 L_num = L_num << norm_num; 467 L_num = L_num >> 1; 468 *result_e = -norm_num + 1; 469 470 norm_den = CountLeadingBits(L_denum); 471 L_denum = L_denum << norm_den; 472 *result_e -= -norm_den; 473 474 div = schur_div(L_num, L_denum, FRACT_BITS); 475 476 return div; 477 } 478 #endif /* !FUNCTION_fDivNorm */ 479 480 #ifndef FUNCTION_fDivNorm 481 FIXP_DBL fDivNorm(FIXP_DBL num, FIXP_DBL denom) { 482 INT e; 483 FIXP_DBL res; 484 485 FDK_ASSERT(denom >= num); 486 487 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 if (res == (FIXP_DBL)(1 << (DFRACT_BITS - 2)) && e == 1) { 492 res = (FIXP_DBL)MAXVAL_DBL; 493 } else { 494 res = scaleValue(res, e); 495 } 496 497 return res; 498 } 499 #endif /* !FUNCTION_fDivNorm */ 500 501 #ifndef FUNCTION_fDivNormSigned 502 FIXP_DBL fDivNormSigned(FIXP_DBL num, FIXP_DBL denom) { 503 INT e; 504 FIXP_DBL res; 505 int sign; 506 507 if (denom == (FIXP_DBL)0) { 508 return (FIXP_DBL)MAXVAL_DBL; 509 } 510 511 sign = ((num >= (FIXP_DBL)0) != (denom >= (FIXP_DBL)0)); 512 res = fDivNormSigned(num, denom, &e); 513 514 /* Saturate since we must output a value with exponent 0 */ 515 if ((e > 0) && (fAbs(res) >= FL2FXCONST_DBL(0.5))) { 516 if (sign) { 517 res = (FIXP_DBL)MINVAL_DBL; 518 } else { 519 res = (FIXP_DBL)MAXVAL_DBL; 520 } 521 } else { 522 res = scaleValue(res, e); 523 } 524 525 return res; 526 } 527 FIXP_DBL fDivNormSigned(FIXP_DBL L_num, FIXP_DBL L_denum, INT *result_e) { 528 FIXP_DBL div; 529 INT norm_num, norm_den; 530 int sign; 531 532 sign = ((L_num >= (FIXP_DBL)0) != (L_denum >= (FIXP_DBL)0)); 533 534 if (L_num == (FIXP_DBL)0) { 535 *result_e = 0; 536 return ((FIXP_DBL)0); 537 } 538 if (L_denum == (FIXP_DBL)0) { 539 *result_e = 14; 540 return ((FIXP_DBL)MAXVAL_DBL); 541 } 542 543 norm_num = CountLeadingBits(L_num); 544 L_num = L_num << norm_num; 545 L_num = L_num >> 2; 546 L_num = fAbs(L_num); 547 *result_e = -norm_num + 1; 548 549 norm_den = CountLeadingBits(L_denum); 550 L_denum = L_denum << norm_den; 551 L_denum = L_denum >> 1; 552 L_denum = fAbs(L_denum); 553 *result_e -= -norm_den; 554 555 div = schur_div(L_num, L_denum, FRACT_BITS); 556 557 if (sign) { 558 div = -div; 559 } 560 561 return div; 562 } 563 #endif /* FUNCTION_fDivNormSigned */ 564 565 #ifndef FUNCTION_fDivNormHighPrec 566 FIXP_DBL fDivNormHighPrec(FIXP_DBL num, FIXP_DBL denom, INT *result_e) { 567 FIXP_DBL div; 568 INT norm_num, norm_den; 569 570 FDK_ASSERT(num >= (FIXP_DBL)0); 571 FDK_ASSERT(denom > (FIXP_DBL)0); 572 573 if (num == (FIXP_DBL)0) { 574 *result_e = 0; 575 return ((FIXP_DBL)0); 576 } 577 578 norm_num = CountLeadingBits(num); 579 num = num << norm_num; 580 num = num >> 1; 581 *result_e = -norm_num + 1; 582 583 norm_den = CountLeadingBits(denom); 584 denom = denom << norm_den; 585 *result_e -= -norm_den; 586 587 div = schur_div(num, denom, 31); 588 return div; 589 } 590 #endif /* !FUNCTION_fDivNormHighPrec */ 591 592 #ifndef FUNCTION_fPow 593 FIXP_DBL f2Pow(const FIXP_DBL exp_m, const INT exp_e, INT *result_e) { 594 FIXP_DBL frac_part, result_m; 595 INT int_part; 596 597 if (exp_e > 0) { 598 INT exp_bits = DFRACT_BITS - 1 - exp_e; 599 int_part = exp_m >> exp_bits; 600 frac_part = exp_m - (FIXP_DBL)(int_part << exp_bits); 601 frac_part = frac_part << exp_e; 602 } else { 603 int_part = 0; 604 frac_part = exp_m >> -exp_e; 605 } 606 607 /* Best accuracy is around 0, so try to get there with the fractional part. */ 608 if (frac_part > FL2FXCONST_DBL(0.5f)) { 609 int_part = int_part + 1; 610 frac_part = frac_part + FL2FXCONST_DBL(-1.0f); 611 } 612 if (frac_part < FL2FXCONST_DBL(-0.5f)) { 613 int_part = int_part - 1; 614 frac_part = -(FL2FXCONST_DBL(-1.0f) - frac_part); 615 } 616 617 /* "+ 1" compensates fMultAddDiv2() of the polynomial evaluation below. */ 618 *result_e = int_part + 1; 619 620 /* Evaluate taylor polynomial which approximates 2^x */ 621 { 622 FIXP_DBL p; 623 624 /* result_m ~= 2^frac_part */ 625 p = frac_part; 626 /* First taylor series coefficient a_0 = 1.0, scaled by 0.5 due to 627 * fMultDiv2(). */ 628 result_m = FL2FXCONST_DBL(1.0f / 2.0f); 629 for (INT i = 0; i < POW2_PRECISION; i++) { 630 /* next taylor series term: a_i * x^i, x=0 */ 631 result_m = fMultAddDiv2(result_m, pow2Coeff[i], p); 632 p = fMult(p, frac_part); 633 } 634 } 635 return result_m; 636 } 637 638 FIXP_DBL f2Pow(const FIXP_DBL exp_m, const INT exp_e) { 639 FIXP_DBL result_m; 640 INT result_e; 641 642 result_m = f2Pow(exp_m, exp_e, &result_e); 643 result_e = fixMin(DFRACT_BITS - 1, fixMax(-(DFRACT_BITS - 1), result_e)); 644 645 return scaleValue(result_m, result_e); 646 } 647 648 FIXP_DBL fPow(FIXP_DBL base_m, INT base_e, FIXP_DBL exp_m, INT exp_e, 649 INT *result_e) { 650 INT ans_lg2_e, baselg2_e; 651 FIXP_DBL base_lg2, ans_lg2, result; 652 653 /* Calc log2 of base */ 654 base_lg2 = fLog2(base_m, base_e, &baselg2_e); 655 656 /* Prepare exp */ 657 { 658 INT leadingBits; 659 660 leadingBits = CountLeadingBits(fAbs(exp_m)); 661 exp_m = exp_m << leadingBits; 662 exp_e -= leadingBits; 663 } 664 665 /* Calc base pow exp */ 666 ans_lg2 = fMult(base_lg2, exp_m); 667 ans_lg2_e = exp_e + baselg2_e; 668 669 /* Calc antilog */ 670 result = f2Pow(ans_lg2, ans_lg2_e, result_e); 671 672 return result; 673 } 674 675 FIXP_DBL fLdPow(FIXP_DBL baseLd_m, INT baseLd_e, FIXP_DBL exp_m, INT exp_e, 676 INT *result_e) { 677 INT ans_lg2_e; 678 FIXP_DBL ans_lg2, result; 679 680 /* Prepare exp */ 681 { 682 INT leadingBits; 683 684 leadingBits = CountLeadingBits(fAbs(exp_m)); 685 exp_m = exp_m << leadingBits; 686 exp_e -= leadingBits; 687 } 688 689 /* Calc base pow exp */ 690 ans_lg2 = fMult(baseLd_m, exp_m); 691 ans_lg2_e = exp_e + baseLd_e; 692 693 /* Calc antilog */ 694 result = f2Pow(ans_lg2, ans_lg2_e, result_e); 695 696 return result; 697 } 698 699 FIXP_DBL fLdPow(FIXP_DBL baseLd_m, INT baseLd_e, FIXP_DBL exp_m, INT exp_e) { 700 FIXP_DBL result_m; 701 int result_e; 702 703 result_m = fLdPow(baseLd_m, baseLd_e, exp_m, exp_e, &result_e); 704 705 return SATURATE_SHIFT(result_m, -result_e, DFRACT_BITS); 706 } 707 708 FIXP_DBL fPowInt(FIXP_DBL base_m, INT base_e, INT exp, INT *pResult_e) { 709 FIXP_DBL result; 710 711 if (exp != 0) { 712 INT result_e = 0; 713 714 if (base_m != (FIXP_DBL)0) { 715 { 716 INT leadingBits; 717 leadingBits = CountLeadingBits(base_m); 718 base_m <<= leadingBits; 719 base_e -= leadingBits; 720 } 721 722 result = base_m; 723 724 { 725 int i; 726 for (i = 1; i < fAbs(exp); i++) { 727 result = fMult(result, base_m); 728 } 729 } 730 731 if (exp < 0) { 732 /* 1.0 / ans */ 733 result = fDivNorm(FL2FXCONST_DBL(0.5f), result, &result_e); 734 result_e++; 735 } else { 736 int ansScale = CountLeadingBits(result); 737 result <<= ansScale; 738 result_e -= ansScale; 739 } 740 741 result_e += exp * base_e; 742 743 } else { 744 result = (FIXP_DBL)0; 745 } 746 *pResult_e = result_e; 747 } else { 748 result = FL2FXCONST_DBL(0.5f); 749 *pResult_e = 1; 750 } 751 752 return result; 753 } 754 #endif /* FUNCTION_fPow */ 755 756 #ifndef FUNCTION_fLog2 757 FIXP_DBL CalcLog2(FIXP_DBL base_m, INT base_e, INT *result_e) { 758 return fLog2(base_m, base_e, result_e); 759 } 760 #endif /* FUNCTION_fLog2 */ 761 762 INT fixp_floorToInt(FIXP_DBL f_inp, INT sf) { 763 FDK_ASSERT(sf >= 0); 764 INT floorInt = (INT)(f_inp >> ((DFRACT_BITS - 1) - sf)); 765 return floorInt; 766 } 767 768 FIXP_DBL fixp_floor(FIXP_DBL f_inp, INT sf) { 769 FDK_ASSERT(sf >= 0); 770 INT floorInt = fixp_floorToInt(f_inp, sf); 771 FIXP_DBL f_floor = (FIXP_DBL)(floorInt << ((DFRACT_BITS - 1) - sf)); 772 return f_floor; 773 } 774 775 INT fixp_ceilToInt(FIXP_DBL f_inp, INT sf) // sf mantissaBits left of dot 776 { 777 FDK_ASSERT(sf >= 0); 778 INT sx = (DFRACT_BITS - 1) - sf; // sx mantissaBits right of dot 779 INT inpINT = (INT)f_inp; 780 781 INT mask = (0x1 << sx) - 1; 782 INT ceilInt = (INT)(f_inp >> sx); 783 784 if (inpINT & mask) { 785 ceilInt++; // increment only, if there is at least one set mantissaBit 786 // right of dot [in inpINT] 787 } 788 789 return ceilInt; 790 } 791 792 FIXP_DBL fixp_ceil(FIXP_DBL f_inp, INT sf) { 793 FDK_ASSERT(sf >= 0); 794 INT sx = (DFRACT_BITS - 1) - sf; 795 INT ceilInt = fixp_ceilToInt(f_inp, sf); 796 ULONG mask = (ULONG)0x1 << (DFRACT_BITS - 1); // 0x80000000 797 ceilInt = ceilInt 798 << sx; // no fract warn bec. shift into saturation done on int side 799 800 if ((f_inp > FL2FXCONST_DBL(0.0f)) && (ceilInt & mask)) { 801 --ceilInt; 802 } 803 FIXP_DBL f_ceil = (FIXP_DBL)ceilInt; 804 805 return f_ceil; 806 } 807 808 /***************************************************************************** 809 fixp_truncateToInt() 810 Just remove the fractional part which is located right of decimal point 811 Same as which is done when a float is casted to (INT) : 812 result_INTtype = (INT)b_floatTypeInput; 813 814 returns INT 815 *****************************************************************************/ 816 INT fixp_truncateToInt(FIXP_DBL f_inp, INT sf) // sf mantissaBits left of dot 817 // (without sign) e.g. at width 818 // 32 this would be [sign]7. 819 // supposed sf equals 8. 820 { 821 FDK_ASSERT(sf >= 0); 822 INT sx = (DFRACT_BITS - 1) - sf; // sx mantissaBits right of dot 823 // at width 32 this would be .24 824 // supposed sf equals 8. 825 INT fbaccu = (INT)f_inp; 826 INT mask = (0x1 << sx); 827 828 if ((fbaccu < 0) && (fbaccu & (mask - 1))) { 829 fbaccu = fbaccu + mask; 830 } 831 832 fbaccu = fbaccu >> sx; 833 return fbaccu; 834 } 835 836 /***************************************************************************** 837 fixp_truncate() 838 Just remove the fractional part which is located right of decimal point 839 840 returns FIXP_DBL 841 *****************************************************************************/ 842 FIXP_DBL fixp_truncate(FIXP_DBL f_inp, INT sf) { 843 FDK_ASSERT(sf >= 0); 844 INT truncateInt = fixp_truncateToInt(f_inp, sf); 845 FIXP_DBL f_truncate = (FIXP_DBL)(truncateInt << ((DFRACT_BITS - 1) - sf)); 846 return f_truncate; 847 } 848 849 /***************************************************************************** 850 fixp_roundToInt() 851 round [typical rounding] 852 853 See fct roundRef() [which is the reference] 854 returns INT 855 *****************************************************************************/ 856 INT fixp_roundToInt(FIXP_DBL f_inp, INT sf) { 857 FDK_ASSERT(sf >= 0); 858 INT sx = DFRACT_BITS - 1 - sf; 859 INT inp = (INT)f_inp; 860 INT mask1 = (0x1 << (sx - 1)); 861 INT mask2 = (0x1 << (sx)) - 1; 862 INT mask3 = 0x7FFFFFFF; 863 INT iam = inp & mask2; 864 INT rnd; 865 866 if ((inp < 0) && !(iam == mask1)) 867 rnd = inp + mask1; 868 else if ((inp > 0) && !(inp == mask3)) 869 rnd = inp + mask1; 870 else 871 rnd = inp; 872 873 rnd = rnd >> sx; 874 875 if (inp == mask3) rnd++; 876 877 return rnd; 878 } 879 880 /***************************************************************************** 881 fixp_round() 882 round [typical rounding] 883 884 See fct roundRef() [which is the reference] 885 returns FIXP_DBL 886 *****************************************************************************/ 887 FIXP_DBL fixp_round(FIXP_DBL f_inp, INT sf) { 888 FDK_ASSERT(sf >= 0); 889 INT sx = DFRACT_BITS - 1 - sf; 890 INT r = fixp_roundToInt(f_inp, sf); 891 ULONG mask = (ULONG)0x1 << (DFRACT_BITS - 1); // 0x80000000 892 r = r << sx; 893 894 if ((f_inp > FL2FXCONST_DBL(0.0f)) && (r & mask)) { 895 --r; 896 } 897 898 FIXP_DBL f_round = (FIXP_DBL)r; 899 return f_round; 900 } 901