The Tor Browser / file comparison
comparison: media/libopus/celt/mathops.h
media/libopus/celt/mathops.h
- branch
- TOR_BUG_9701
- changeset 10
- ac0c01689b40
equal
deleted
inserted
replaced
| -1:000000000000 | 0:ed44997f4202 |
|---|---|
| 1 /* Copyright (c) 2002-2008 Jean-Marc Valin | |
| 2 Copyright (c) 2007-2008 CSIRO | |
| 3 Copyright (c) 2007-2009 Xiph.Org Foundation | |
| 4 Written by Jean-Marc Valin */ | |
| 5 /** | |
| 6 @file mathops.h | |
| 7 @brief Various math functions | |
| 8 */ | |
| 9 /* | |
| 10 Redistribution and use in source and binary forms, with or without | |
| 11 modification, are permitted provided that the following conditions | |
| 12 are met: | |
| 13 | |
| 14 - Redistributions of source code must retain the above copyright | |
| 15 notice, this list of conditions and the following disclaimer. | |
| 16 | |
| 17 - Redistributions in binary form must reproduce the above copyright | |
| 18 notice, this list of conditions and the following disclaimer in the | |
| 19 documentation and/or other materials provided with the distribution. | |
| 20 | |
| 21 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS | |
| 22 ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT | |
| 23 LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR | |
| 24 A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER | |
| 25 OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, | |
| 26 EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, | |
| 27 PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR | |
| 28 PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF | |
| 29 LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING | |
| 30 NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS | |
| 31 SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. | |
| 32 */ | |
| 33 | |
| 34 #ifndef MATHOPS_H | |
| 35 #define MATHOPS_H | |
| 36 | |
| 37 #include "arch.h" | |
| 38 #include "entcode.h" | |
| 39 #include "os_support.h" | |
| 40 | |
| 41 /* Multiplies two 16-bit fractional values. Bit-exactness of this macro is important */ | |
| 42 #define FRAC_MUL16(a,b) ((16384+((opus_int32)(opus_int16)(a)*(opus_int16)(b)))>>15) | |
| 43 | |
| 44 unsigned isqrt32(opus_uint32 _val); | |
| 45 | |
| 46 #ifndef OVERRIDE_CELT_MAXABS16 | |
| 47 static OPUS_INLINE opus_val32 celt_maxabs16(const opus_val16 *x, int len) | |
| 48 { | |
| 49 int i; | |
| 50 opus_val16 maxval = 0; | |
| 51 opus_val16 minval = 0; | |
| 52 for (i=0;i<len;i++) | |
| 53 { | |
| 54 maxval = MAX16(maxval, x[i]); | |
| 55 minval = MIN16(minval, x[i]); | |
| 56 } | |
| 57 return MAX32(EXTEND32(maxval),-EXTEND32(minval)); | |
| 58 } | |
| 59 #endif | |
| 60 | |
| 61 #ifndef OVERRIDE_CELT_MAXABS32 | |
| 62 #ifdef FIXED_POINT | |
| 63 static OPUS_INLINE opus_val32 celt_maxabs32(const opus_val32 *x, int len) | |
| 64 { | |
| 65 int i; | |
| 66 opus_val32 maxval = 0; | |
| 67 opus_val32 minval = 0; | |
| 68 for (i=0;i<len;i++) | |
| 69 { | |
| 70 maxval = MAX32(maxval, x[i]); | |
| 71 minval = MIN32(minval, x[i]); | |
| 72 } | |
| 73 return MAX32(maxval, -minval); | |
| 74 } | |
| 75 #else | |
| 76 #define celt_maxabs32(x,len) celt_maxabs16(x,len) | |
| 77 #endif | |
| 78 #endif | |
| 79 | |
| 80 | |
| 81 #ifndef FIXED_POINT | |
| 82 | |
| 83 #define PI 3.141592653f | |
| 84 #define celt_sqrt(x) ((float)sqrt(x)) | |
| 85 #define celt_rsqrt(x) (1.f/celt_sqrt(x)) | |
| 86 #define celt_rsqrt_norm(x) (celt_rsqrt(x)) | |
| 87 #define celt_cos_norm(x) ((float)cos((.5f*PI)*(x))) | |
| 88 #define celt_rcp(x) (1.f/(x)) | |
| 89 #define celt_div(a,b) ((a)/(b)) | |
| 90 #define frac_div32(a,b) ((float)(a)/(b)) | |
| 91 | |
| 92 #ifdef FLOAT_APPROX | |
| 93 | |
| 94 /* Note: This assumes radix-2 floating point with the exponent at bits 23..30 and an offset of 127 | |
| 95 denorm, +/- inf and NaN are *not* handled */ | |
| 96 | |
| 97 /** Base-2 log approximation (log2(x)). */ | |
| 98 static OPUS_INLINE float celt_log2(float x) | |
| 99 { | |
| 100 int integer; | |
| 101 float frac; | |
| 102 union { | |
| 103 float f; | |
| 104 opus_uint32 i; | |
| 105 } in; | |
| 106 in.f = x; | |
| 107 integer = (in.i>>23)-127; | |
| 108 in.i -= integer<<23; | |
| 109 frac = in.f - 1.5f; | |
| 110 frac = -0.41445418f + frac*(0.95909232f | |
| 111 + frac*(-0.33951290f + frac*0.16541097f)); | |
| 112 return 1+integer+frac; | |
| 113 } | |
| 114 | |
| 115 /** Base-2 exponential approximation (2^x). */ | |
| 116 static OPUS_INLINE float celt_exp2(float x) | |
| 117 { | |
| 118 int integer; | |
| 119 float frac; | |
| 120 union { | |
| 121 float f; | |
| 122 opus_uint32 i; | |
| 123 } res; | |
| 124 integer = floor(x); | |
| 125 if (integer < -50) | |
| 126 return 0; | |
| 127 frac = x-integer; | |
| 128 /* K0 = 1, K1 = log(2), K2 = 3-4*log(2), K3 = 3*log(2) - 2 */ | |
| 129 res.f = 0.99992522f + frac * (0.69583354f | |
| 130 + frac * (0.22606716f + 0.078024523f*frac)); | |
| 131 res.i = (res.i + (integer<<23)) & 0x7fffffff; | |
| 132 return res.f; | |
| 133 } | |
| 134 | |
| 135 #else | |
| 136 #define celt_log2(x) ((float)(1.442695040888963387*log(x))) | |
| 137 #define celt_exp2(x) ((float)exp(0.6931471805599453094*(x))) | |
| 138 #endif | |
| 139 | |
| 140 #endif | |
| 141 | |
| 142 #ifdef FIXED_POINT | |
| 143 | |
| 144 #include "os_support.h" | |
| 145 | |
| 146 #ifndef OVERRIDE_CELT_ILOG2 | |
| 147 /** Integer log in base2. Undefined for zero and negative numbers */ | |
| 148 static OPUS_INLINE opus_int16 celt_ilog2(opus_int32 x) | |
| 149 { | |
| 150 celt_assert2(x>0, "celt_ilog2() only defined for strictly positive numbers"); | |
| 151 return EC_ILOG(x)-1; | |
| 152 } | |
| 153 #endif | |
| 154 | |
| 155 | |
| 156 /** Integer log in base2. Defined for zero, but not for negative numbers */ | |
| 157 static OPUS_INLINE opus_int16 celt_zlog2(opus_val32 x) | |
| 158 { | |
| 159 return x <= 0 ? 0 : celt_ilog2(x); | |
| 160 } | |
| 161 | |
| 162 opus_val16 celt_rsqrt_norm(opus_val32 x); | |
| 163 | |
| 164 opus_val32 celt_sqrt(opus_val32 x); | |
| 165 | |
| 166 opus_val16 celt_cos_norm(opus_val32 x); | |
| 167 | |
| 168 /** Base-2 logarithm approximation (log2(x)). (Q14 input, Q10 output) */ | |
| 169 static OPUS_INLINE opus_val16 celt_log2(opus_val32 x) | |
| 170 { | |
| 171 int i; | |
| 172 opus_val16 n, frac; | |
| 173 /* -0.41509302963303146, 0.9609890551383969, -0.31836011537636605, | |
| 174 0.15530808010959576, -0.08556153059057618 */ | |
| 175 static const opus_val16 C[5] = {-6801+(1<<(13-DB_SHIFT)), 15746, -5217, 2545, -1401}; | |
| 176 if (x==0) | |
| 177 return -32767; | |
| 178 i = celt_ilog2(x); | |
| 179 n = VSHR32(x,i-15)-32768-16384; | |
| 180 frac = ADD16(C[0], MULT16_16_Q15(n, ADD16(C[1], MULT16_16_Q15(n, ADD16(C[2], MULT16_16_Q15(n, ADD16(C[3], MULT16_16_Q15(n, C[4])))))))); | |
| 181 return SHL16(i-13,DB_SHIFT)+SHR16(frac,14-DB_SHIFT); | |
| 182 } | |
| 183 | |
| 184 /* | |
| 185 K0 = 1 | |
| 186 K1 = log(2) | |
| 187 K2 = 3-4*log(2) | |
| 188 K3 = 3*log(2) - 2 | |
| 189 */ | |
| 190 #define D0 16383 | |
| 191 #define D1 22804 | |
| 192 #define D2 14819 | |
| 193 #define D3 10204 | |
| 194 | |
| 195 static OPUS_INLINE opus_val32 celt_exp2_frac(opus_val16 x) | |
| 196 { | |
| 197 opus_val16 frac; | |
| 198 frac = SHL16(x, 4); | |
| 199 return ADD16(D0, MULT16_16_Q15(frac, ADD16(D1, MULT16_16_Q15(frac, ADD16(D2 , MULT16_16_Q15(D3,frac)))))); | |
| 200 } | |
| 201 /** Base-2 exponential approximation (2^x). (Q10 input, Q16 output) */ | |
| 202 static OPUS_INLINE opus_val32 celt_exp2(opus_val16 x) | |
| 203 { | |
| 204 int integer; | |
| 205 opus_val16 frac; | |
| 206 integer = SHR16(x,10); | |
| 207 if (integer>14) | |
| 208 return 0x7f000000; | |
| 209 else if (integer < -15) | |
| 210 return 0; | |
| 211 frac = celt_exp2_frac(x-SHL16(integer,10)); | |
| 212 return VSHR32(EXTEND32(frac), -integer-2); | |
| 213 } | |
| 214 | |
| 215 opus_val32 celt_rcp(opus_val32 x); | |
| 216 | |
| 217 #define celt_div(a,b) MULT32_32_Q31((opus_val32)(a),celt_rcp(b)) | |
| 218 | |
| 219 opus_val32 frac_div32(opus_val32 a, opus_val32 b); | |
| 220 | |
| 221 #define M1 32767 | |
| 222 #define M2 -21 | |
| 223 #define M3 -11943 | |
| 224 #define M4 4936 | |
| 225 | |
| 226 /* Atan approximation using a 4th order polynomial. Input is in Q15 format | |
| 227 and normalized by pi/4. Output is in Q15 format */ | |
| 228 static OPUS_INLINE opus_val16 celt_atan01(opus_val16 x) | |
| 229 { | |
| 230 return MULT16_16_P15(x, ADD32(M1, MULT16_16_P15(x, ADD32(M2, MULT16_16_P15(x, ADD32(M3, MULT16_16_P15(M4, x))))))); | |
| 231 } | |
| 232 | |
| 233 #undef M1 | |
| 234 #undef M2 | |
| 235 #undef M3 | |
| 236 #undef M4 | |
| 237 | |
| 238 /* atan2() approximation valid for positive input values */ | |
| 239 static OPUS_INLINE opus_val16 celt_atan2p(opus_val16 y, opus_val16 x) | |
| 240 { | |
| 241 if (y < x) | |
| 242 { | |
| 243 opus_val32 arg; | |
| 244 arg = celt_div(SHL32(EXTEND32(y),15),x); | |
| 245 if (arg >= 32767) | |
| 246 arg = 32767; | |
| 247 return SHR16(celt_atan01(EXTRACT16(arg)),1); | |
| 248 } else { | |
| 249 opus_val32 arg; | |
| 250 arg = celt_div(SHL32(EXTEND32(x),15),y); | |
| 251 if (arg >= 32767) | |
| 252 arg = 32767; | |
| 253 return 25736-SHR16(celt_atan01(EXTRACT16(arg)),1); | |
| 254 } | |
| 255 } | |
| 256 | |
| 257 #endif /* FIXED_POINT */ | |
| 258 #endif /* MATHOPS_H */ |
