common/alnumeric.h

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