common/alnumeric.h

   1 #ifndef AL_NUMERIC_H
   2 #define AL_NUMERIC_H
   3
   4 #include "config_simd.h"
   5
   6 #include <algorithm>
   7 #include <array>
   8 #include <cmath>
   9 #include <cstddef>
  10 #include <cstdint>
  11 #include <iterator>
  12 #include <type_traits>
  13 #ifdef HAVE_INTRIN_H
  14 #include <intrin.h>
  15 #endif
  16 #if HAVE_SSE_INTRINSICS
  17 #include <xmmintrin.h>
  18 #endif
  19
  20 #include "albit.h"
  21 #include "altraits.h"
  22 #include "opthelpers.h"
  23
  24
  25 constexpr auto operator "" _i64(unsigned long long n) noexcept { return static_cast<std::int64_t>(n); }
  26 constexpr auto operator "" _u64(unsigned long long n) noexcept { return static_cast<std::uint64_t>(n); }
  27
  28 constexpr auto operator "" _z(unsigned long long n) noexcept
  29 { return static_cast<std::make_signed_t<std::size_t>>(n); }
  30 constexpr auto operator "" _uz(unsigned long long n) noexcept { return static_cast<std::size_t>(n); }
  31 constexpr auto operator "" _zu(unsigned long long n) noexcept { return static_cast<std::size_t>(n); }
  32
  33
  34 constexpr auto GetCounterSuffix(size_t count) noexcept -> const char*
  35 {
  36     auto &suffix = (((count%100)/10) == 1) ? "th" :
  37         ((count%10) == 1) ? "st" :
  38         ((count%10) == 2) ? "nd" :
  39         ((count%10) == 3) ? "rd" : "th";
  40     return std::data(suffix);
  41 }
  42
  43
  44 constexpr inline float lerpf(float val1, float val2, float mu) noexcept
  45 { return val1 + (val2-val1)*mu; }
  46 constexpr inline double lerpd(double val1, double val2, double mu) noexcept
  47 { return val1 + (val2-val1)*mu; }
  48
  49
  50 /** Find the next power-of-2 for non-power-of-2 numbers. */
  51 inline uint32_t NextPowerOf2(uint32_t value) noexcept
  52 {
  53     if(value > 0)
  54     {
  55         value--;
  56         value |= value>>1;
  57         value |= value>>2;
  58         value |= value>>4;
  59         value |= value>>8;
  60         value |= value>>16;
  61     }
  62     return value+1;
  63 }
  64
  65 /**
  66  * If the value is not already a multiple of r, round down to the next
  67  * multiple.
  68  */
  69 template<typename T>
  70 constexpr T RoundDown(T value, al::type_identity_t<T> r) noexcept
  71 { return value - (value%r); }
  72
  73 /**
  74  * If the value is not already a multiple of r, round up to the next multiple.
  75  */
  76 template<typename T>
  77 constexpr T RoundUp(T value, al::type_identity_t<T> r) noexcept
  78 { return RoundDown(value + r-1, r); }
  79
  80
  81 /**
  82  * Fast float-to-int conversion. No particular rounding mode is assumed; the
  83  * IEEE-754 default is round-to-nearest with ties-to-even, though an app could
  84  * change it on its own threads. On some systems, a truncating conversion may
  85  * always be the fastest method.
  86  */
  87 inline int fastf2i(float f) noexcept
  88 {
  89 #if HAVE_SSE_INTRINSICS
  90     return _mm_cvt_ss2si(_mm_set_ss(f));
  91
  92 #elif defined(_MSC_VER) && defined(_M_IX86_FP) && _M_IX86_FP == 0
  93
  94     int i;
  95     __asm fld f
  96     __asm fistp i
  97     return i;
  98
  99 #elif (defined(__GNUC__) || defined(__clang__)) && (defined(__i386__) || defined(__x86_64__)) \
 100     && !defined(__SSE_MATH__)
 101
 102     int i;
 103     __asm__ __volatile__("fistpl %0" : "=m"(i) : "t"(f) : "st");
 104     return i;
 105
 106 #else
 107
 108     return static_cast<int>(f);
 109 #endif
 110 }
 111 inline unsigned int fastf2u(float f) noexcept
 112 { return static_cast<unsigned int>(fastf2i(f)); }
 113
 114 /** Converts float-to-int using standard behavior (truncation). */
 115 inline int float2int(float f) noexcept
 116 {
 117 #if HAVE_SSE_INTRINSICS
 118     return _mm_cvtt_ss2si(_mm_set_ss(f));
 119
 120 #elif (defined(_MSC_VER) && defined(_M_IX86_FP) && _M_IX86_FP == 0) \
 121     || ((defined(__GNUC__) || defined(__clang__)) && (defined(__i386__) || defined(__x86_64__)) \
 122         && !defined(__SSE_MATH__))
 123     const int conv_i{al::bit_cast<int>(f)};
 124
 125     const int sign{(conv_i>>31) | 1};
 126     const int shift{((conv_i>>23)&0xff) - (127+23)};
 127
 128     /* Over/underflow */
 129     if(shift >= 31 || shift < -23) UNLIKELY
 130         return 0;
 131
 132     const int mant{(conv_i&0x7fffff) | 0x800000};
 133     if(shift < 0) LIKELY
 134         return (mant >> -shift) * sign;
 135     return (mant << shift) * sign;
 136
 137 #else
 138
 139     return static_cast<int>(f);
 140 #endif
 141 }
 142 inline unsigned int float2uint(float f) noexcept
 143 { return static_cast<unsigned int>(float2int(f)); }
 144
 145 /** Converts double-to-int using standard behavior (truncation). */
 146 inline int double2int(double d) noexcept
 147 {
 148 #if HAVE_SSE_INTRINSICS
 149     return _mm_cvttsd_si32(_mm_set_sd(d));
 150
 151 #elif (defined(_MSC_VER) && defined(_M_IX86_FP) && _M_IX86_FP < 2) \
 152     || ((defined(__GNUC__) || defined(__clang__)) && (defined(__i386__) || defined(__x86_64__)) \
 153         && !defined(__SSE2_MATH__))
 154     const int64_t conv_i64{al::bit_cast<int64_t>(d)};
 155
 156     const int sign{static_cast<int>(conv_i64 >> 63) | 1};
 157     const int shift{(static_cast<int>(conv_i64 >> 52) & 0x7ff) - (1023 + 52)};
 158
 159     /* Over/underflow */
 160     if(shift >= 63 || shift < -52) UNLIKELY
 161         return 0;
 162
 163     const int64_t mant{(conv_i64 & 0xfffffffffffff_i64) | 0x10000000000000_i64};
 164     if(shift < 0) LIKELY
 165         return static_cast<int>(mant >> -shift) * sign;
 166     return static_cast<int>(mant << shift) * sign;
 167
 168 #else
 169
 170     return static_cast<int>(d);
 171 #endif
 172 }
 173
 174 /**
 175  * Rounds a float to the nearest integral value, according to the current
 176  * rounding mode. This is essentially an inlined version of rintf, although
 177  * makes fewer promises (e.g. -0 or -0.25 rounded to 0 may result in +0).
 178  */
 179 inline float fast_roundf(float f) noexcept
 180 {
 181 #if (defined(__GNUC__) || defined(__clang__)) && (defined(__i386__) || defined(__x86_64__)) \
 182     && !defined(__SSE_MATH__)
 183
 184     float out;
 185     __asm__ __volatile__("frndint" : "=t"(out) : "0"(f));
 186     return out;
 187
 188 #elif (defined(__GNUC__) || defined(__clang__)) && defined(__aarch64__)
 189
 190     float out;
 191     __asm__ volatile("frintx %s0, %s1" : "=w"(out) : "w"(f));
 192     return out;
 193
 194 #else
 195
 196     /* Integral limit, where sub-integral precision is not available for
 197      * floats.
 198      */
 199     static constexpr std::array ilim{
 200          8388608.0f /*  0x1.0p+23 */,
 201         -8388608.0f /* -0x1.0p+23 */
 202     };
 203     const unsigned int conv_i{al::bit_cast<unsigned int>(f)};
 204
 205     const unsigned int sign{(conv_i>>31)&0x01};
 206     const unsigned int expo{(conv_i>>23)&0xff};
 207
 208     if(expo >= 150/*+23*/) UNLIKELY
 209     {
 210         /* An exponent (base-2) of 23 or higher is incapable of sub-integral
 211          * precision, so it's already an integral value. We don't need to worry
 212          * about infinity or NaN here.
 213          */
 214         return f;
 215     }
 216     /* Adding the integral limit to the value (with a matching sign) forces a
 217      * result that has no sub-integral precision, and is consequently forced to
 218      * round to an integral value. Removing the integral limit then restores
 219      * the initial value rounded to the integral. The compiler should not
 220      * optimize this out because of non-associative rules on floating-point
 221      * math (as long as you don't use -fassociative-math,
 222      * -funsafe-math-optimizations, -ffast-math, or -Ofast, in which case this
 223      * may break without __builtin_assoc_barrier support).
 224      */
 225 #if HAS_BUILTIN(__builtin_assoc_barrier)
 226     return __builtin_assoc_barrier(f + ilim[sign]) - ilim[sign];
 227 #else
 228     f += ilim[sign];
 229     return f - ilim[sign];
 230 #endif
 231 #endif
 232 }
 233
 234
 235 // Converts level (mB) to gain.
 236 inline float level_mb_to_gain(float x)
 237 {
 238     if(x <= -10'000.0f)
 239         return 0.0f;
 240     return std::pow(10.0f, x / 2'000.0f);
 241 }
 242
 243 // Converts gain to level (mB).
 244 inline float gain_to_level_mb(float x)
 245 {
 246     if(x <= 1e-05f)
 247         return -10'000.0f;
 248     return std::max(std::log10(x) * 2'000.0f, -10'000.0f);
 249 }
 250
 251 #endif /* AL_NUMERIC_H */