xbmc/utils/MathUtils.h

   1 /*
   2  *  Copyright (C) 2005-2018 Team Kodi
   3  *  This file is part of Kodi - https://kodi.tv
   4  *
   5  *  SPDX-License-Identifier: GPL-2.0-or-later
   6  *  See LICENSES/README.md for more information.
   7  */
   8
   9 #pragma once
  10
  11 #include <assert.h>
  12 #include <climits>
  13 #include <cmath>
  14 #include <stdint.h>
  15 #include <type_traits>
  16
  17 #if defined(HAVE_SSE2) && defined(__SSE2__)
  18 #include <emmintrin.h>
  19 #endif
  20
  21 // use real compiler defines in here as we want to
  22 // avoid including system.h or other magic includes.
  23 // use 'gcc -dM -E - < /dev/null' or similar to find them.
  24
  25 // clang-format off
  26 #if defined(__aarch64__) || \
  27     defined(__alpha__) || \
  28     defined(__arc__) || \
  29     defined(__arm__) || \
  30     defined(__loongarch__) || \
  31     defined(_M_ARM) || \
  32     defined(__mips__) || \
  33     defined(__or1k__) || \
  34     defined(__powerpc__) || \
  35     defined(__ppc__) || \
  36     defined(__riscv) || \
  37     defined(__SH4__) || \
  38     defined(__s390x__) || \
  39     defined(__sparc__) || \
  40     defined(__xtensa__)
  41 #define DISABLE_MATHUTILS_ASM_ROUND_INT
  42 #endif
  43 // clang-format on
  44
  45 /*! \brief Math utility class.
  46  Note that the test() routine should return true for all implementations
  47
  48  See http://ldesoras.free.fr/doc/articles/rounding_en.pdf for an explanation
  49  of the technique used on x86.
  50  */
  51 namespace MathUtils
  52 {
  53   // GCC does something stupid with optimization on release builds if we try
  54   // to assert in these functions
  55
  56   /*! \brief Round to nearest integer.
  57    This routine does fast rounding to the nearest integer.
  58    In the case (k + 0.5 for any integer k) we round up to k+1, and in all other
  59    instances we should return the nearest integer.
  60    Thus, { -1.5, -0.5, 0.5, 1.5 } is rounded to { -1, 0, 1, 2 }.
  61    It preserves the property that round(k) - round(k-1) = 1 for all doubles k.
  62
  63    Make sure MathUtils::test() returns true for each implementation.
  64    \sa truncate_int, test
  65   */
  66   inline int round_int(double x)
  67   {
  68     assert(x > static_cast<double>((int) (INT_MIN / 2)) - 1.0);
  69     assert(x < static_cast<double>((int) (INT_MAX / 2)) + 1.0);
  70
  71 #if defined(DISABLE_MATHUTILS_ASM_ROUND_INT)
  72     /* This implementation warrants some further explanation.
  73      *
  74      * First, a couple of notes on rounding:
  75      * 1) C casts from float/double to integer round towards zero.
  76      * 2) Float/double additions are rounded according to the normal rules,
  77      *    in other words: on some architectures, it's fixed at compile-time,
  78      *    and on others it can be set using fesetround()). The following
  79      *    analysis assumes round-to-nearest with ties rounding to even. This
  80      *    is a fairly sensible choice, and is the default with ARM VFP.
  81      *
  82      * What this function wants is round-to-nearest with ties rounding to
  83      * +infinity. This isn't an IEEE rounding mode, even if we could guarantee
  84      * that all architectures supported fesetround(), which they don't. Instead,
  85      * this adds an offset of 2147483648.5 (= 0x80000000.8p0), then casts to
  86      * an unsigned int (crucially, all possible inputs are now in a range where
  87      * round to zero acts the same as round to -infinity) and then subtracts
  88      * 0x80000000 in the integer domain. The 0.5 component of the offset
  89      * converts what is effectively a round down into a round to nearest, with
  90      * ties rounding up, as desired.
  91      *
  92      * There is a catch, that because there is a double rounding, there is a
  93      * small region where the input falls just *below* a tie, where the addition
  94      * of the offset causes a round *up* to an exact integer, due to the finite
  95      * level of precision available in floating point. You need to be aware of
  96      * this when calling this function, although at present it is not believed
  97      * that XBMC ever attempts to round numbers in this window.
  98      *
  99      * It is worth proving the size of the affected window. Recall that double
 100      * precision employs a mantissa of 52 bits.
 101      * 1) For all inputs -0.5 <= x <= INT_MAX
 102      *    Once the offset is applied, the most significant binary digit in the
 103      *    floating-point representation is +2^31.
 104      *    At this magnitude, the smallest step representable in double precision
 105      *    is 2^31 / 2^52 = 0.000000476837158203125
 106      *    So the size of the range which is rounded up due to the addition is
 107      *    half the size of this step, or 0.0000002384185791015625
 108      *
 109      * 2) For all inputs INT_MIN/2 < x < -0.5
 110      *    Once the offset is applied, the most significant binary digit in the
 111      *    floating-point representation is +2^30.
 112      *    At this magnitude, the smallest step representable in double precision
 113      *    is 2^30 / 2^52 = 0.0000002384185791015625
 114      *    So the size of the range which is rounded up due to the addition is
 115      *    half the size of this step, or 0.00000011920928955078125
 116      *
 117      * 3) For all inputs INT_MIN <= x <= INT_MIN/2
 118      *    The representation once the offset is applied has equal or greater
 119      *    precision than the input, so the addition does not cause rounding.
 120      */
 121     return ((unsigned int) (x + 2147483648.5)) - 0x80000000;
 122
 123 #else
 124     const float round_to_nearest = 0.5f;
 125     int i;
 126 #if defined(HAVE_SSE2) && defined(__SSE2__)
 127     const float round_dn_to_nearest = 0.4999999f;
 128     i = (x > 0) ? _mm_cvttsd_si32(_mm_set_sd(x + static_cast<double>(round_to_nearest)))
 129                 : _mm_cvttsd_si32(_mm_set_sd(x - static_cast<double>(round_dn_to_nearest)));
 130
 131 #elif defined(TARGET_WINDOWS)
 132     __asm
 133     {
 134       fld x
 135       fadd st, st (0)
 136       fadd round_to_nearest
 137       fistp i
 138       sar i, 1
 139     }
 140
 141 #else
 142     __asm__ __volatile__ (
 143       "fadd %%st\n\t"
 144       "fadd %%st(1)\n\t"
 145       "fistpl %0\n\t"
 146       "sarl $1, %0\n"
 147       : "=m"(i) : "u"(round_to_nearest), "t"(x) : "st"
 148     );
 149
 150 #endif
 151     return i;
 152 #endif
 153   }
 154
 155   /*! \brief Round to nearest integer.
 156    \sa truncate_int, test
 157   */
 158   inline int round_int(float x)
 159   {
 160     return round_int(static_cast<double>(x));
 161   }
 162
 163   /*! \brief Truncate to nearest integer.
 164    This routine does fast truncation to an integer.
 165    It should simply drop the fractional portion of the floating point number.
 166
 167    Make sure MathUtils::test() returns true for each implementation.
 168    \sa round_int, test
 169   */
 170   inline int truncate_int(double x)
 171   {
 172     assert(x > static_cast<double>(INT_MIN / 2) - 1.0);
 173     assert(x < static_cast<double>(INT_MAX / 2) + 1.0);
 174     return static_cast<int>(x);
 175   }
 176
 177   inline int64_t abs(int64_t a)
 178   {
 179     return (a < 0) ? -a : a;
 180   }
 181
 182   inline unsigned bitcount(unsigned v)
 183   {
 184     unsigned c = 0;
 185     for (c = 0; v; c++)
 186       v &= v - 1; // clear the least significant bit set
 187     return c;
 188   }
 189
 190   inline void hack()
 191   {
 192     // stupid hack to keep compiler from dropping these
 193     // functions as unused
 194     MathUtils::round_int(0.0);
 195     MathUtils::truncate_int(0.0);
 196     MathUtils::abs(0);
 197   }
 198
 199   /**
 200    * Compare two floating-point numbers for equality and regard them
 201    * as equal if their difference is below a given threshold.
 202    *
 203    * It is usually not useful to compare float numbers for equality with
 204    * the standard operator== since very close numbers might have different
 205    * representations.
 206    */
 207   template<typename FloatT>
 208   inline bool FloatEquals(FloatT f1, FloatT f2, FloatT maxDelta)
 209   {
 210     return (std::abs(f2 - f1) < maxDelta);
 211   }
 212
 213   /*!
 214    * \brief Round a floating point number to nearest multiple
 215    * \param value The value to round
 216    * \param multiple The multiple
 217    * \return The rounded value
 218    */
 219   template<typename T, std::enable_if_t<std::is_floating_point<T>::value, bool> = true>
 220   inline T RoundF(const T value, const T multiple)
 221   {
 222     if (multiple == 0)
 223       return value;
 224
 225     return static_cast<T>(std::round(static_cast<double>(value) / static_cast<double>(multiple)) *
 226                           static_cast<double>(multiple));
 227   }
 228
 229 #if 0
 230   /*! \brief test routine for round_int and truncate_int
 231    Must return true on all platforms.
 232    */
 233   inline bool test()
 234   {
 235     for (int i = -8; i < 8; ++i)
 236     {
 237       double d = 0.25*i;
 238       int r = (i < 0) ? (i - 1) / 4 : (i + 2) / 4;
 239       int t = i / 4;
 240       if (round_int(d) != r || truncate_int(d) != t)
 241         return false;
 242     }
 243     return true;
 244   }
 245 #endif
 246 } // namespace MathUtils
 247