| 1 | /* | 
|---|
| 2 | * Copyright 2006 The Android Open Source Project | 
|---|
| 3 | * | 
|---|
| 4 | * Use of this source code is governed by a BSD-style license that can be | 
|---|
| 5 | * found in the LICENSE file. | 
|---|
| 6 | */ | 
|---|
| 7 |  | 
|---|
| 8 | #ifndef SkFloatingPoint_DEFINED | 
|---|
| 9 | #define SkFloatingPoint_DEFINED | 
|---|
| 10 |  | 
|---|
| 11 | #include "include/core/SkTypes.h" | 
|---|
| 12 | #include "include/private/SkFloatBits.h" | 
|---|
| 13 | #include "include/private/SkSafe_math.h" | 
|---|
| 14 | #include <float.h> | 
|---|
| 15 | #include <math.h> | 
|---|
| 16 | #include <cmath> | 
|---|
| 17 | #include <cstring> | 
|---|
| 18 | #include <limits> | 
|---|
| 19 |  | 
|---|
| 20 |  | 
|---|
| 21 | #if SK_CPU_SSE_LEVEL >= SK_CPU_SSE_LEVEL_SSE1 | 
|---|
| 22 | #include <xmmintrin.h> | 
|---|
| 23 | #elif defined(SK_ARM_HAS_NEON) | 
|---|
| 24 | #include <arm_neon.h> | 
|---|
| 25 | #endif | 
|---|
| 26 |  | 
|---|
| 27 | // For _POSIX_VERSION | 
|---|
| 28 | #if defined(__unix__) || (defined(__APPLE__) && defined(__MACH__)) | 
|---|
| 29 | #include <unistd.h> | 
|---|
| 30 | #endif | 
|---|
| 31 |  | 
|---|
| 32 | constexpr float SK_FloatSqrt2 = 1.41421356f; | 
|---|
| 33 | constexpr float SK_FloatPI    = 3.14159265f; | 
|---|
| 34 | constexpr double SK_DoublePI  = 3.14159265358979323846264338327950288; | 
|---|
| 35 |  | 
|---|
| 36 | // C++98 cmath std::pow seems to be the earliest portable way to get float pow. | 
|---|
| 37 | // However, on Linux including cmath undefines isfinite. | 
|---|
| 38 | // http://gcc.gnu.org/bugzilla/show_bug.cgi?id=14608 | 
|---|
| 39 | static inline float sk_float_pow(float base, float exp) { | 
|---|
| 40 | return powf(base, exp); | 
|---|
| 41 | } | 
|---|
| 42 |  | 
|---|
| 43 | #define sk_float_sqrt(x)        sqrtf(x) | 
|---|
| 44 | #define sk_float_sin(x)         sinf(x) | 
|---|
| 45 | #define sk_float_cos(x)         cosf(x) | 
|---|
| 46 | #define sk_float_tan(x)         tanf(x) | 
|---|
| 47 | #define sk_float_floor(x)       floorf(x) | 
|---|
| 48 | #define sk_float_ceil(x)        ceilf(x) | 
|---|
| 49 | #define sk_float_trunc(x)       truncf(x) | 
|---|
| 50 | #ifdef SK_BUILD_FOR_MAC | 
|---|
| 51 | #    define sk_float_acos(x)    static_cast<float>(acos(x)) | 
|---|
| 52 | #    define sk_float_asin(x)    static_cast<float>(asin(x)) | 
|---|
| 53 | #else | 
|---|
| 54 | #    define sk_float_acos(x)    acosf(x) | 
|---|
| 55 | #    define sk_float_asin(x)    asinf(x) | 
|---|
| 56 | #endif | 
|---|
| 57 | #define sk_float_atan2(y,x)     atan2f(y,x) | 
|---|
| 58 | #define sk_float_abs(x)         fabsf(x) | 
|---|
| 59 | #define sk_float_copysign(x, y) copysignf(x, y) | 
|---|
| 60 | #define sk_float_mod(x,y)       fmodf(x,y) | 
|---|
| 61 | #define sk_float_exp(x)         expf(x) | 
|---|
| 62 | #define sk_float_log(x)         logf(x) | 
|---|
| 63 |  | 
|---|
| 64 | constexpr float sk_float_degrees_to_radians(float degrees) { | 
|---|
| 65 | return degrees * (SK_FloatPI / 180); | 
|---|
| 66 | } | 
|---|
| 67 |  | 
|---|
| 68 | constexpr float sk_float_radians_to_degrees(float radians) { | 
|---|
| 69 | return radians * (180 / SK_FloatPI); | 
|---|
| 70 | } | 
|---|
| 71 |  | 
|---|
| 72 | #define sk_float_round(x) sk_float_floor((x) + 0.5f) | 
|---|
| 73 |  | 
|---|
| 74 | // can't find log2f on android, but maybe that just a tool bug? | 
|---|
| 75 | #ifdef SK_BUILD_FOR_ANDROID | 
|---|
| 76 | static inline float sk_float_log2(float x) { | 
|---|
| 77 | const double inv_ln_2 = 1.44269504088896; | 
|---|
| 78 | return (float)(log(x) * inv_ln_2); | 
|---|
| 79 | } | 
|---|
| 80 | #else | 
|---|
| 81 | #define sk_float_log2(x)        log2f(x) | 
|---|
| 82 | #endif | 
|---|
| 83 |  | 
|---|
| 84 | static inline bool sk_float_isfinite(float x) { | 
|---|
| 85 | return SkFloatBits_IsFinite(SkFloat2Bits(x)); | 
|---|
| 86 | } | 
|---|
| 87 |  | 
|---|
| 88 | static inline bool sk_floats_are_finite(float a, float b) { | 
|---|
| 89 | return sk_float_isfinite(a) && sk_float_isfinite(b); | 
|---|
| 90 | } | 
|---|
| 91 |  | 
|---|
| 92 | static inline bool sk_floats_are_finite(const float array[], int count) { | 
|---|
| 93 | float prod = 0; | 
|---|
| 94 | for (int i = 0; i < count; ++i) { | 
|---|
| 95 | prod *= array[i]; | 
|---|
| 96 | } | 
|---|
| 97 | // At this point, prod will either be NaN or 0 | 
|---|
| 98 | return prod == 0;   // if prod is NaN, this check will return false | 
|---|
| 99 | } | 
|---|
| 100 |  | 
|---|
| 101 | static inline bool sk_float_isinf(float x) { | 
|---|
| 102 | return SkFloatBits_IsInf(SkFloat2Bits(x)); | 
|---|
| 103 | } | 
|---|
| 104 |  | 
|---|
| 105 | static inline bool sk_float_isnan(float x) { | 
|---|
| 106 | return !(x == x); | 
|---|
| 107 | } | 
|---|
| 108 |  | 
|---|
| 109 | #define sk_double_isnan(a)          sk_float_isnan(a) | 
|---|
| 110 |  | 
|---|
| 111 | #define SK_MaxS32FitsInFloat    2147483520 | 
|---|
| 112 | #define SK_MinS32FitsInFloat    -SK_MaxS32FitsInFloat | 
|---|
| 113 |  | 
|---|
| 114 | #define SK_MaxS64FitsInFloat    (SK_MaxS64 >> (63-24) << (63-24))   // 0x7fffff8000000000 | 
|---|
| 115 | #define SK_MinS64FitsInFloat    -SK_MaxS64FitsInFloat | 
|---|
| 116 |  | 
|---|
| 117 | /** | 
|---|
| 118 | *  Return the closest int for the given float. Returns SK_MaxS32FitsInFloat for NaN. | 
|---|
| 119 | */ | 
|---|
| 120 | static inline int sk_float_saturate2int(float x) { | 
|---|
| 121 | x = x < SK_MaxS32FitsInFloat ? x : SK_MaxS32FitsInFloat; | 
|---|
| 122 | x = x > SK_MinS32FitsInFloat ? x : SK_MinS32FitsInFloat; | 
|---|
| 123 | return (int)x; | 
|---|
| 124 | } | 
|---|
| 125 |  | 
|---|
| 126 | /** | 
|---|
| 127 | *  Return the closest int for the given double. Returns SK_MaxS32 for NaN. | 
|---|
| 128 | */ | 
|---|
| 129 | static inline int sk_double_saturate2int(double x) { | 
|---|
| 130 | x = x < SK_MaxS32 ? x : SK_MaxS32; | 
|---|
| 131 | x = x > SK_MinS32 ? x : SK_MinS32; | 
|---|
| 132 | return (int)x; | 
|---|
| 133 | } | 
|---|
| 134 |  | 
|---|
| 135 | /** | 
|---|
| 136 | *  Return the closest int64_t for the given float. Returns SK_MaxS64FitsInFloat for NaN. | 
|---|
| 137 | */ | 
|---|
| 138 | static inline int64_t sk_float_saturate2int64(float x) { | 
|---|
| 139 | x = x < SK_MaxS64FitsInFloat ? x : SK_MaxS64FitsInFloat; | 
|---|
| 140 | x = x > SK_MinS64FitsInFloat ? x : SK_MinS64FitsInFloat; | 
|---|
| 141 | return (int64_t)x; | 
|---|
| 142 | } | 
|---|
| 143 |  | 
|---|
| 144 | #define sk_float_floor2int(x)   sk_float_saturate2int(sk_float_floor(x)) | 
|---|
| 145 | #define sk_float_round2int(x)   sk_float_saturate2int(sk_float_floor((x) + 0.5f)) | 
|---|
| 146 | #define sk_float_ceil2int(x)    sk_float_saturate2int(sk_float_ceil(x)) | 
|---|
| 147 |  | 
|---|
| 148 | #define sk_float_floor2int_no_saturate(x)   (int)sk_float_floor(x) | 
|---|
| 149 | #define sk_float_round2int_no_saturate(x)   (int)sk_float_floor((x) + 0.5f) | 
|---|
| 150 | #define sk_float_ceil2int_no_saturate(x)    (int)sk_float_ceil(x) | 
|---|
| 151 |  | 
|---|
| 152 | #define sk_double_floor(x)          floor(x) | 
|---|
| 153 | #define sk_double_round(x)          floor((x) + 0.5) | 
|---|
| 154 | #define sk_double_ceil(x)           ceil(x) | 
|---|
| 155 | #define sk_double_floor2int(x)      (int)floor(x) | 
|---|
| 156 | #define sk_double_round2int(x)      (int)floor((x) + 0.5) | 
|---|
| 157 | #define sk_double_ceil2int(x)       (int)ceil(x) | 
|---|
| 158 |  | 
|---|
| 159 | // Cast double to float, ignoring any warning about too-large finite values being cast to float. | 
|---|
| 160 | // Clang thinks this is undefined, but it's actually implementation defined to return either | 
|---|
| 161 | // the largest float or infinity (one of the two bracketing representable floats).  Good enough! | 
|---|
| 162 | SK_ATTRIBUTE(no_sanitize( "float-cast-overflow")) | 
|---|
| 163 | static inline float sk_double_to_float(double x) { | 
|---|
| 164 | return static_cast<float>(x); | 
|---|
| 165 | } | 
|---|
| 166 |  | 
|---|
| 167 | #define SK_FloatNaN                 std::numeric_limits<float>::quiet_NaN() | 
|---|
| 168 | #define SK_FloatInfinity            (+std::numeric_limits<float>::infinity()) | 
|---|
| 169 | #define SK_FloatNegativeInfinity    (-std::numeric_limits<float>::infinity()) | 
|---|
| 170 |  | 
|---|
| 171 | #define SK_DoubleNaN                std::numeric_limits<double>::quiet_NaN() | 
|---|
| 172 |  | 
|---|
| 173 | // Returns false if any of the floats are outside of [0...1] | 
|---|
| 174 | // Returns true if count is 0 | 
|---|
| 175 | bool sk_floats_are_unit(const float array[], size_t count); | 
|---|
| 176 |  | 
|---|
| 177 | static inline float sk_float_rsqrt_portable(float x) { | 
|---|
| 178 | // Get initial estimate. | 
|---|
| 179 | int i; | 
|---|
| 180 | memcpy(&i, &x, 4); | 
|---|
| 181 | i = 0x5F1FFFF9 - (i>>1); | 
|---|
| 182 | float estimate; | 
|---|
| 183 | memcpy(&estimate, &i, 4); | 
|---|
| 184 |  | 
|---|
| 185 | // One step of Newton's method to refine. | 
|---|
| 186 | const float estimate_sq = estimate*estimate; | 
|---|
| 187 | estimate *= 0.703952253f*(2.38924456f-x*estimate_sq); | 
|---|
| 188 | return estimate; | 
|---|
| 189 | } | 
|---|
| 190 |  | 
|---|
| 191 | // Fast, approximate inverse square root. | 
|---|
| 192 | // Compare to name-brand "1.0f / sk_float_sqrt(x)".  Should be around 10x faster on SSE, 2x on NEON. | 
|---|
| 193 | static inline float sk_float_rsqrt(float x) { | 
|---|
| 194 | // We want all this inlined, so we'll inline SIMD and just take the hit when we don't know we've got | 
|---|
| 195 | // it at compile time.  This is going to be too fast to productively hide behind a function pointer. | 
|---|
| 196 | // | 
|---|
| 197 | // We do one step of Newton's method to refine the estimates in the NEON and portable paths.  No | 
|---|
| 198 | // refinement is faster, but very innacurate.  Two steps is more accurate, but slower than 1/sqrt. | 
|---|
| 199 | // | 
|---|
| 200 | // Optimized constants in the portable path courtesy of http://rrrola.wz.cz/inv_sqrt.html | 
|---|
| 201 | #if SK_CPU_SSE_LEVEL >= SK_CPU_SSE_LEVEL_SSE1 | 
|---|
| 202 | return _mm_cvtss_f32(_mm_rsqrt_ss(_mm_set_ss(x))); | 
|---|
| 203 | #elif defined(SK_ARM_HAS_NEON) | 
|---|
| 204 | // Get initial estimate. | 
|---|
| 205 | const float32x2_t xx = vdup_n_f32(x);  // Clever readers will note we're doing everything 2x. | 
|---|
| 206 | float32x2_t estimate = vrsqrte_f32(xx); | 
|---|
| 207 |  | 
|---|
| 208 | // One step of Newton's method to refine. | 
|---|
| 209 | const float32x2_t estimate_sq = vmul_f32(estimate, estimate); | 
|---|
| 210 | estimate = vmul_f32(estimate, vrsqrts_f32(xx, estimate_sq)); | 
|---|
| 211 | return vget_lane_f32(estimate, 0);  // 1 will work fine too; the answer's in both places. | 
|---|
| 212 | #else | 
|---|
| 213 | return sk_float_rsqrt_portable(x); | 
|---|
| 214 | #endif | 
|---|
| 215 | } | 
|---|
| 216 |  | 
|---|
| 217 | // This is the number of significant digits we can print in a string such that when we read that | 
|---|
| 218 | // string back we get the floating point number we expect.  The minimum value C requires is 6, but | 
|---|
| 219 | // most compilers support 9 | 
|---|
| 220 | #ifdef FLT_DECIMAL_DIG | 
|---|
| 221 | #define SK_FLT_DECIMAL_DIG FLT_DECIMAL_DIG | 
|---|
| 222 | #else | 
|---|
| 223 | #define SK_FLT_DECIMAL_DIG 9 | 
|---|
| 224 | #endif | 
|---|
| 225 |  | 
|---|
| 226 | // IEEE defines how float divide behaves for non-finite values and zero-denoms, but C does not | 
|---|
| 227 | // so we have a helper that suppresses the possible undefined-behavior warnings. | 
|---|
| 228 |  | 
|---|
| 229 | SK_ATTRIBUTE(no_sanitize( "float-divide-by-zero")) | 
|---|
| 230 | static inline float sk_ieee_float_divide(float numer, float denom) { | 
|---|
| 231 | return numer / denom; | 
|---|
| 232 | } | 
|---|
| 233 |  | 
|---|
| 234 | SK_ATTRIBUTE(no_sanitize( "float-divide-by-zero")) | 
|---|
| 235 | static inline double sk_ieee_double_divide(double numer, double denom) { | 
|---|
| 236 | return numer / denom; | 
|---|
| 237 | } | 
|---|
| 238 |  | 
|---|
| 239 | // While we clean up divide by zero, we'll replace places that do divide by zero with this TODO. | 
|---|
| 240 | static inline float sk_ieee_float_divide_TODO_IS_DIVIDE_BY_ZERO_SAFE_HERE(float n, float d) { | 
|---|
| 241 | return sk_ieee_float_divide(n,d); | 
|---|
| 242 | } | 
|---|
| 243 | static inline float sk_ieee_double_divide_TODO_IS_DIVIDE_BY_ZERO_SAFE_HERE(double n, double d) { | 
|---|
| 244 | return sk_ieee_double_divide(n,d); | 
|---|
| 245 | } | 
|---|
| 246 |  | 
|---|
| 247 | static inline float sk_fmaf(float f, float m, float a) { | 
|---|
| 248 | #if defined(FP_FAST_FMA) | 
|---|
| 249 | return std::fmaf(f,m,a); | 
|---|
| 250 | #else | 
|---|
| 251 | return f*m+a; | 
|---|
| 252 | #endif | 
|---|
| 253 | } | 
|---|
| 254 |  | 
|---|
| 255 | #endif | 
|---|
| 256 |  | 
|---|