HDK
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Math.h
Go to the documentation of this file.
1 // Copyright Contributors to the OpenVDB Project
2 // SPDX-License-Identifier: MPL-2.0
3 //
4 /// @file Math.h
5 /// @brief General-purpose arithmetic and comparison routines, most of which
6 /// accept arbitrary value types (or at least arbitrary numeric value types)
7 
8 #ifndef OPENVDB_MATH_HAS_BEEN_INCLUDED
9 #define OPENVDB_MATH_HAS_BEEN_INCLUDED
10 
11 #include <openvdb/Platform.h>
12 #include <openvdb/version.h>
13 #include <hboost/numeric/conversion/conversion_traits.hpp>
14 #include <algorithm> // for std::max()
15 #include <cassert>
16 #include <cmath> // for std::ceil(), std::fabs(), std::pow(), std::sqrt(), etc.
17 #include <cstdlib> // for abs(int)
18 #include <random>
19 #include <string>
20 #include <type_traits> // for std::is_arithmetic
21 
22 
23 // Compile pragmas
24 
25 // Intel(r) compiler fires remark #1572: floating-point equality and inequality
26 // comparisons are unrealiable when == or != is used with floating point operands.
27 #if defined(__INTEL_COMPILER)
28  #define OPENVDB_NO_FP_EQUALITY_WARNING_BEGIN \
29  _Pragma("warning (push)") \
30  _Pragma("warning (disable:1572)")
31  #define OPENVDB_NO_FP_EQUALITY_WARNING_END \
32  _Pragma("warning (pop)")
33 #elif defined(__clang__)
34  #define OPENVDB_NO_FP_EQUALITY_WARNING_BEGIN \
35  PRAGMA(clang diagnostic push) \
36  PRAGMA(clang diagnostic ignored "-Wfloat-equal")
37  #define OPENVDB_NO_FP_EQUALITY_WARNING_END \
38  PRAGMA(clang diagnostic pop)
39 #else
40  // For GCC, #pragma GCC diagnostic ignored "-Wfloat-equal"
41  // isn't working until gcc 4.2+,
42  // Trying
43  // #pragma GCC system_header
44  // creates other problems, most notably "warning: will never be executed"
45  // in from templates, unsure of how to work around.
46  // If necessary, could use integer based comparisons for equality
47  #define OPENVDB_NO_FP_EQUALITY_WARNING_BEGIN
48  #define OPENVDB_NO_FP_EQUALITY_WARNING_END
49 #endif
50 
51 
52 #ifdef OPENVDB_IS_POD
53 #undef OPENVDB_IS_POD
54 #endif
55 #define OPENVDB_IS_POD(Type) \
56 static_assert(std::is_standard_layout<Type>::value, \
57  #Type" must be a POD type (satisfy StandardLayoutType.)"); \
58 static_assert(std::is_trivial<Type>::value, \
59  #Type" must be a POD type (satisfy TrivialType.)");
60 
61 namespace openvdb {
63 namespace OPENVDB_VERSION_NAME {
64 
65 /// @brief Return the value of type T that corresponds to zero.
66 /// @note A zeroVal<T>() specialization must be defined for each @c ValueType T
67 /// that cannot be constructed using the form @c T(0). For example, @c std::string(0)
68 /// treats 0 as @c nullptr and throws a @c std::logic_error.
69 template<typename T> inline T zeroVal() { return T(0); }
70 /// Return the @c std::string value that corresponds to zero.
71 template<> inline std::string zeroVal<std::string>() { return ""; }
72 /// Return the @c bool value that corresponds to zero.
73 template<> inline bool zeroVal<bool>() { return false; }
74 
75 namespace math {
76 
77 /// @todo These won't be needed if we eliminate StringGrids.
78 //@{
79 /// @brief Needed to support the <tt>(zeroVal<ValueType>() + val)</tt> idiom
80 /// when @c ValueType is @c std::string
81 inline std::string operator+(const std::string& s, bool) { return s; }
82 inline std::string operator+(const std::string& s, int) { return s; }
83 inline std::string operator+(const std::string& s, float) { return s; }
84 inline std::string operator+(const std::string& s, double) { return s; }
85 //@}
86 
87 /// @brief Componentwise adder for POD types.
88 template<typename Type1, typename Type2>
89 inline auto cwiseAdd(const Type1& v, const Type2 s)
90 {
92  return v + s;
94 }
95 
96 /// @brief Componentwise less than for POD types.
97 template<typename Type1, typename Type2>
98 inline bool cwiseLessThan(const Type1& a, const Type2& b)
99 {
101  return a < b;
103 }
104 
105 /// @brief Componentwise greater than for POD types.
106 template<typename Type1, typename Type2>
107 inline bool cwiseGreaterThan(const Type1& a, const Type2& b)
108 {
110  return a > b;
112 }
113 
114 
115 
116 /// @brief Pi constant taken from Boost to match old behaviour
117 /// @note Available in C++20
118 template <typename T> inline constexpr T pi() { return 3.141592653589793238462643383279502884e+00; }
119 template <> inline constexpr float pi() { return 3.141592653589793238462643383279502884e+00F; }
120 template <> inline constexpr double pi() { return 3.141592653589793238462643383279502884e+00; }
121 template <> inline constexpr long double pi() { return 3.141592653589793238462643383279502884e+00L; }
122 
123 
124 /// @brief Return the unary negation of the given value.
125 /// @note A negative<T>() specialization must be defined for each ValueType T
126 /// for which unary negation is not defined.
127 template<typename T> inline T negative(const T& val)
128 {
129 // disable unary minus on unsigned warning
130 #if defined(_MSC_VER)
131 #pragma warning(push)
132 #pragma warning(disable:4146)
133 #endif
134  return T(-val);
135 #if defined(_MSC_VER)
136 #pragma warning(pop)
137 #endif
138 }
139 /// Return the negation of the given boolean.
140 template<> inline bool negative(const bool& val) { return !val; }
141 /// Return the "negation" of the given string.
142 template<> inline std::string negative(const std::string& val) { return val; }
143 
144 
145 //@{
146 /// Tolerance for floating-point comparison
147 template<typename T> struct Tolerance { static T value() { return zeroVal<T>(); } };
148 template<> struct Tolerance<float> { static float value() { return 1e-8f; } };
149 template<> struct Tolerance<double> { static double value() { return 1e-15; } };
150 //@}
151 
152 //@{
153 /// Delta for small floating-point offsets
154 template<typename T> struct Delta { static T value() { return zeroVal<T>(); } };
155 template<> struct Delta<float> { static float value() { return 1e-5f; } };
156 template<> struct Delta<double> { static double value() { return 1e-9; } };
157 //@}
158 
159 
160 // ==========> Random Values <==================
161 
162 /// @brief Simple generator of random numbers over the range [0, 1)
163 /// @details Thread-safe as long as each thread has its own Rand01 instance
164 template<typename FloatType = double, typename EngineType = std::mt19937>
165 class Rand01
166 {
167 private:
168  EngineType mEngine;
169  std::uniform_real_distribution<FloatType> mRand;
170 
171 public:
173 
174  /// @brief Initialize the generator.
175  /// @param engine random number generator
176  Rand01(const EngineType& engine): mEngine(engine) {}
177 
178  /// @brief Initialize the generator.
179  /// @param seed seed value for the random number generator
180  Rand01(unsigned int seed): mEngine(static_cast<typename EngineType::result_type>(seed)) {}
181 
182  /// Set the seed value for the random number generator
183  void setSeed(unsigned int seed)
184  {
185  mEngine.seed(static_cast<typename EngineType::result_type>(seed));
186  }
187 
188  /// Return a const reference to the random number generator.
189  const EngineType& engine() const { return mEngine; }
190 
191  /// Return a uniformly distributed random number in the range [0, 1).
192  FloatType operator()() { return mRand(mEngine); }
193 };
194 
196 
197 
198 /// @brief Simple random integer generator
199 /// @details Thread-safe as long as each thread has its own RandInt instance
200 template<typename IntType = int, typename EngineType = std::mt19937>
201 class RandInt
202 {
203 private:
204  using Distr = std::uniform_int_distribution<IntType>;
205  EngineType mEngine;
206  Distr mRand;
207 
208 public:
209  /// @brief Initialize the generator.
210  /// @param engine random number generator
211  /// @param imin,imax generate integers that are uniformly distributed over [imin, imax]
212  RandInt(const EngineType& engine, IntType imin, IntType imax):
213  mEngine(engine),
214  mRand(std::min(imin, imax), std::max(imin, imax))
215  {}
216 
217  /// @brief Initialize the generator.
218  /// @param seed seed value for the random number generator
219  /// @param imin,imax generate integers that are uniformly distributed over [imin, imax]
220  RandInt(unsigned int seed, IntType imin, IntType imax):
221  mEngine(static_cast<typename EngineType::result_type>(seed)),
222  mRand(std::min(imin, imax), std::max(imin, imax))
223  {}
224 
225  /// Change the range over which integers are distributed to [imin, imax].
226  void setRange(IntType imin, IntType imax)
227  {
228  mRand = Distr(std::min(imin, imax), std::max(imin, imax));
229  }
230 
231  /// Set the seed value for the random number generator
232  void setSeed(unsigned int seed)
233  {
234  mEngine.seed(static_cast<typename EngineType::result_type>(seed));
235  }
236 
237  /// Return a const reference to the random number generator.
238  const EngineType& engine() const { return mEngine; }
239 
240  /// Return a randomly-generated integer in the current range.
241  IntType operator()() { return mRand(mEngine); }
242 
243  /// @brief Return a randomly-generated integer in the new range [imin, imax],
244  /// without changing the current range.
245  IntType operator()(IntType imin, IntType imax)
246  {
247  const IntType lo = std::min(imin, imax), hi = std::max(imin, imax);
248  return mRand(mEngine, typename Distr::param_type(lo, hi));
249  }
250 };
251 
253 
254 
255 // ==========> Clamp <==================
256 
257 /// Return @a x clamped to [@a min, @a max]
258 template<typename Type>
259 inline Type
261 {
262  assert( !(min>max) );
263  return x > min ? x < max ? x : max : min;
264 }
265 
266 
267 /// Return @a x clamped to [0, 1]
268 template<typename Type>
269 inline Type
270 Clamp01(Type x) { return x > Type(0) ? x < Type(1) ? x : Type(1) : Type(0); }
271 
272 
273 /// Return @c true if @a x is outside [0,1]
274 template<typename Type>
275 inline bool
277 {
278  if (x >= Type(0) && x <= Type(1)) return false;
279  x = x < Type(0) ? Type(0) : Type(1);
280  return true;
281 }
282 
283 /// @brief Return 0 if @a x < @a 0, 1 if @a x > 1 or else (3 &minus; 2 @a x) @a x&sup2;.
284 template<typename Type>
285 inline Type
287 {
288  return x > 0 ? x < 1 ? (3-2*x)*x*x : Type(1) : Type(0);
289 }
290 
291 /// @brief Return 0 if @a x < @a min, 1 if @a x > @a max or else (3 &minus; 2 @a t) @a t&sup2;,
292 /// where @a t = (@a x &minus; @a min)/(@a max &minus; @a min).
293 template<typename Type>
294 inline Type
296 {
297  assert(min < max);
298  return SmoothUnitStep((x-min)/(max-min));
299 }
300 
301 
302 // ==========> Absolute Value <==================
303 
304 
305 //@{
306 /// Return the absolute value of the given quantity.
307 inline int32_t Abs(int32_t i) { return abs(i); }
308 inline int64_t Abs(int64_t i)
309 {
310 #ifdef _MSC_VER
311  return (i < int64_t(0) ? -i : i);
312 #else
313  return labs(i);
314 #endif
315 }
316 inline float Abs(float x) { return std::fabs(x); }
317 inline double Abs(double x) { return std::fabs(x); }
318 inline long double Abs(long double x) { return std::fabs(x); }
319 inline uint32_t Abs(uint32_t i) { return i; }
320 inline uint64_t Abs(uint64_t i) { return i; }
321 inline bool Abs(bool b) { return b; }
322 // On systems like macOS and FreeBSD, size_t and uint64_t are different types
323 template <typename T>
325 Abs(T i) { return i; }
326 //@}
327 
328 
329 ////////////////////////////////////////
330 
331 
332 // ==========> Value Comparison <==================
333 
334 
335 /// Return @c true if @a x is exactly equal to zero.
336 template<typename Type>
337 inline bool
338 isZero(const Type& x)
339 {
341  return x == zeroVal<Type>();
343 }
344 
345 
346 /// @brief Return @c true if @a x is equal to zero to within
347 /// the default floating-point comparison tolerance.
348 template<typename Type>
349 inline bool
351 {
352  const Type tolerance = Type(zeroVal<Type>() + Tolerance<Type>::value());
353  return !(x > tolerance) && !(x < -tolerance);
354 }
355 
356 /// Return @c true if @a x is equal to zero to within the given tolerance.
357 template<typename Type>
358 inline bool
359 isApproxZero(const Type& x, const Type& tolerance)
360 {
361  return !(x > tolerance) && !(x < -tolerance);
362 }
363 
364 
365 /// Return @c true if @a x is less than zero.
366 template<typename Type>
367 inline bool
368 isNegative(const Type& x) { return x < zeroVal<Type>(); }
369 
370 // Return false, since bool values are never less than zero.
371 template<> inline bool isNegative<bool>(const bool&) { return false; }
372 
373 
374 /// Return @c true if @a x is finite.
375 inline bool
376 isFinite(const float x) { return std::isfinite(x); }
377 
378 /// Return @c true if @a x is finite.
380 inline bool
381 isFinite(const Type& x) { return std::isfinite(static_cast<double>(x)); }
382 
383 
384 /// Return @c true if @a x is an infinity value (either positive infinity or negative infinity).
385 inline bool
386 isInfinite(const float x) { return std::isinf(x); }
387 
388 /// Return @c true if @a x is an infinity value (either positive infinity or negative infinity).
390 inline bool
391 isInfinite(const Type& x) { return std::isinf(static_cast<double>(x)); }
392 
393 
394 /// Return @c true if @a x is a NaN (Not-A-Number) value.
395 inline bool
396 isNan(const float x) { return std::isnan(x); }
397 
398 /// Return @c true if @a x is a NaN (Not-A-Number) value.
400 inline bool
401 isNan(const Type& x) { return std::isnan(static_cast<double>(x)); }
402 
403 
404 /// Return @c true if @a a is equal to @a b to within the given tolerance.
405 template<typename Type>
406 inline bool
407 isApproxEqual(const Type& a, const Type& b, const Type& tolerance)
408 {
409  return !cwiseGreaterThan(Abs(a - b), tolerance);
410 }
411 
412 /// @brief Return @c true if @a a is equal to @a b to within
413 /// the default floating-point comparison tolerance.
414 template<typename Type>
415 inline bool
416 isApproxEqual(const Type& a, const Type& b)
417 {
418  const Type tolerance = Type(zeroVal<Type>() + Tolerance<Type>::value());
419  return isApproxEqual(a, b, tolerance);
420 }
421 
422 #define OPENVDB_EXACT_IS_APPROX_EQUAL(T) \
423  template<> inline bool isApproxEqual<T>(const T& a, const T& b) { return a == b; } \
424  template<> inline bool isApproxEqual<T>(const T& a, const T& b, const T&) { return a == b; } \
425  /**/
426 
429 
430 
431 /// @brief Return @c true if @a a is larger than @a b to within
432 /// the given tolerance, i.e., if @a b - @a a < @a tolerance.
433 template<typename Type>
434 inline bool
436 {
437  return (b - a < tolerance);
438 }
439 
440 
441 /// @brief Return @c true if @a a is exactly equal to @a b.
442 template<typename T0, typename T1>
443 inline bool
444 isExactlyEqual(const T0& a, const T1& b)
445 {
447  return a == b;
449 }
450 
451 
452 template<typename Type>
453 inline bool
454 isRelOrApproxEqual(const Type& a, const Type& b, const Type& absTol, const Type& relTol)
455 {
456  // First check to see if we are inside the absolute tolerance
457  // Necessary for numbers close to 0
458  if (!(Abs(a - b) > absTol)) return true;
459 
460  // Next check to see if we are inside the relative tolerance
461  // to handle large numbers that aren't within the abs tolerance
462  // but could be the closest floating point representation
463  double relError;
464  if (Abs(b) > Abs(a)) {
465  relError = Abs((a - b) / b);
466  } else {
467  relError = Abs((a - b) / a);
468  }
469  return (relError <= relTol);
470 }
471 
472 template<>
473 inline bool
474 isRelOrApproxEqual(const bool& a, const bool& b, const bool&, const bool&)
475 {
476  return (a == b);
477 }
478 
479 
480 // Avoid strict aliasing issues by using type punning
481 // http://cellperformance.beyond3d.com/articles/2006/06/understanding-strict-aliasing.html
482 // Using "casting through a union(2)"
483 inline int32_t
484 floatToInt32(const float aFloatValue)
485 {
486  union FloatOrInt32 { float floatValue; int32_t int32Value; };
487  const FloatOrInt32* foi = reinterpret_cast<const FloatOrInt32*>(&aFloatValue);
488  return foi->int32Value;
489 }
490 
491 
492 inline int64_t
493 doubleToInt64(const double aDoubleValue)
494 {
495  union DoubleOrInt64 { double doubleValue; int64_t int64Value; };
496  const DoubleOrInt64* dol = reinterpret_cast<const DoubleOrInt64*>(&aDoubleValue);
497  return dol->int64Value;
498 }
499 
500 
501 // aUnitsInLastPlace is the allowed difference between the least significant digits
502 // of the numbers' floating point representation
503 // Please read the reference paper before trying to use isUlpsEqual
504 // http://www.cygnus-software.com/papers/comparingfloats/comparingfloats.htm
505 inline bool
506 isUlpsEqual(const double aLeft, const double aRight, const int64_t aUnitsInLastPlace)
507 {
508  int64_t longLeft = doubleToInt64(aLeft);
509  // Because of 2's complement, must restore lexicographical order
510  if (longLeft < 0) {
511  longLeft = INT64_C(0x8000000000000000) - longLeft;
512  }
513 
514  int64_t longRight = doubleToInt64(aRight);
515  // Because of 2's complement, must restore lexicographical order
516  if (longRight < 0) {
517  longRight = INT64_C(0x8000000000000000) - longRight;
518  }
519 
520  int64_t difference = labs(longLeft - longRight);
521  return (difference <= aUnitsInLastPlace);
522 }
523 
524 inline bool
525 isUlpsEqual(const float aLeft, const float aRight, const int32_t aUnitsInLastPlace)
526 {
527  int32_t intLeft = floatToInt32(aLeft);
528  // Because of 2's complement, must restore lexicographical order
529  if (intLeft < 0) {
530  intLeft = 0x80000000 - intLeft;
531  }
532 
533  int32_t intRight = floatToInt32(aRight);
534  // Because of 2's complement, must restore lexicographical order
535  if (intRight < 0) {
536  intRight = 0x80000000 - intRight;
537  }
538 
539  int32_t difference = abs(intLeft - intRight);
540  return (difference <= aUnitsInLastPlace);
541 }
542 
543 
544 ////////////////////////////////////////
545 
546 
547 // ==========> Pow <==================
548 
549 /// Return @a x<sup>2</sup>.
550 template<typename Type>
551 inline Type Pow2(Type x) { return x*x; }
552 
553 /// Return @a x<sup>3</sup>.
554 template<typename Type>
555 inline Type Pow3(Type x) { return x*x*x; }
556 
557 /// Return @a x<sup>4</sup>.
558 template<typename Type>
559 inline Type Pow4(Type x) { return Pow2(Pow2(x)); }
560 
561 /// Return @a x<sup>n</sup>.
562 template<typename Type>
563 Type
564 Pow(Type x, int n)
565 {
566  Type ans = 1;
567  if (n < 0) {
568  n = -n;
569  x = Type(1)/x;
570  }
571  while (n--) ans *= x;
572  return ans;
573 }
574 
575 //@{
576 /// Return @a b<sup>e</sup>.
577 inline float
578 Pow(float b, float e)
579 {
580  assert( b >= 0.0f && "Pow(float,float): base is negative" );
581  return powf(b,e);
582 }
583 
584 inline double
585 Pow(double b, double e)
586 {
587  assert( b >= 0.0 && "Pow(double,double): base is negative" );
588  return std::pow(b,e);
589 }
590 //@}
591 
592 
593 // ==========> Max <==================
594 
595 /// Return the maximum of two values
596 template<typename Type>
597 inline const Type&
598 Max(const Type& a, const Type& b)
599 {
600  return std::max(a,b);
601 }
602 
603 /// Return the maximum of three values
604 template<typename Type>
605 inline const Type&
606 Max(const Type& a, const Type& b, const Type& c)
607 {
608  return std::max(std::max(a,b), c);
609 }
610 
611 /// Return the maximum of four values
612 template<typename Type>
613 inline const Type&
614 Max(const Type& a, const Type& b, const Type& c, const Type& d)
615 {
616  return std::max(std::max(a,b), std::max(c,d));
617 }
618 
619 /// Return the maximum of five values
620 template<typename Type>
621 inline const Type&
622 Max(const Type& a, const Type& b, const Type& c, const Type& d, const Type& e)
623 {
624  return std::max(std::max(a,b), Max(c,d,e));
625 }
626 
627 /// Return the maximum of six values
628 template<typename Type>
629 inline const Type&
630 Max(const Type& a, const Type& b, const Type& c, const Type& d, const Type& e, const Type& f)
631 {
632  return std::max(Max(a,b,c), Max(d,e,f));
633 }
634 
635 /// Return the maximum of seven values
636 template<typename Type>
637 inline const Type&
638 Max(const Type& a, const Type& b, const Type& c, const Type& d,
639  const Type& e, const Type& f, const Type& g)
640 {
641  return std::max(Max(a,b,c,d), Max(e,f,g));
642 }
643 
644 /// Return the maximum of eight values
645 template<typename Type>
646 inline const Type&
647 Max(const Type& a, const Type& b, const Type& c, const Type& d,
648  const Type& e, const Type& f, const Type& g, const Type& h)
649 {
650  return std::max(Max(a,b,c,d), Max(e,f,g,h));
651 }
652 
653 
654 // ==========> Min <==================
655 
656 /// Return the minimum of two values
657 template<typename Type>
658 inline const Type&
659 Min(const Type& a, const Type& b) { return std::min(a, b); }
660 
661 /// Return the minimum of three values
662 template<typename Type>
663 inline const Type&
664 Min(const Type& a, const Type& b, const Type& c) { return std::min(std::min(a, b), c); }
665 
666 /// Return the minimum of four values
667 template<typename Type>
668 inline const Type&
669 Min(const Type& a, const Type& b, const Type& c, const Type& d)
670 {
671  return std::min(std::min(a, b), std::min(c, d));
672 }
673 
674 /// Return the minimum of five values
675 template<typename Type>
676 inline const Type&
677 Min(const Type& a, const Type& b, const Type& c, const Type& d, const Type& e)
678 {
679  return std::min(std::min(a,b), Min(c,d,e));
680 }
681 
682 /// Return the minimum of six values
683 template<typename Type>
684 inline const Type&
685 Min(const Type& a, const Type& b, const Type& c, const Type& d, const Type& e, const Type& f)
686 {
687  return std::min(Min(a,b,c), Min(d,e,f));
688 }
689 
690 /// Return the minimum of seven values
691 template<typename Type>
692 inline const Type&
693 Min(const Type& a, const Type& b, const Type& c, const Type& d,
694  const Type& e, const Type& f, const Type& g)
695 {
696  return std::min(Min(a,b,c,d), Min(e,f,g));
697 }
698 
699 /// Return the minimum of eight values
700 template<typename Type>
701 inline const Type&
702 Min(const Type& a, const Type& b, const Type& c, const Type& d,
703  const Type& e, const Type& f, const Type& g, const Type& h)
704 {
705  return std::min(Min(a,b,c,d), Min(e,f,g,h));
706 }
707 
708 
709 // ============> Exp <==================
710 
711 /// Return @a e<sup>x</sup>.
712 template<typename Type>
713 inline Type Exp(const Type& x) { return std::exp(x); }
714 
715 // ============> Sin <==================
716 
717 //@{
718 /// Return sin @a x.
719 inline float Sin(const float& x) { return std::sin(x); }
720 
721 inline double Sin(const double& x) { return std::sin(x); }
722 //@}
723 
724 // ============> Cos <==================
725 
726 //@{
727 /// Return cos @a x.
728 inline float Cos(const float& x) { return std::cos(x); }
729 
730 inline double Cos(const double& x) { return std::cos(x); }
731 //@}
732 
733 
734 ////////////////////////////////////////
735 
736 
737 /// Return the sign of the given value as an integer (either -1, 0 or 1).
738 template <typename Type>
739 inline int Sign(const Type &x) { return (zeroVal<Type>() < x) - (x < zeroVal<Type>()); }
740 
741 
742 /// @brief Return @c true if @a a and @a b have different signs.
743 /// @note Zero is considered a positive number.
744 template <typename Type>
745 inline bool
746 SignChange(const Type& a, const Type& b)
747 {
748  return ( (a<zeroVal<Type>()) ^ (b<zeroVal<Type>()) );
749 }
750 
751 
752 /// @brief Return @c true if the interval [@a a, @a b] includes zero,
753 /// i.e., if either @a a or @a b is zero or if they have different signs.
754 template <typename Type>
755 inline bool
756 ZeroCrossing(const Type& a, const Type& b)
757 {
758  return a * b <= zeroVal<Type>();
759 }
760 
761 
762 //@{
763 /// Return the square root of a floating-point value.
764 inline float Sqrt(float x) { return std::sqrt(x); }
765 inline double Sqrt(double x) { return std::sqrt(x); }
766 inline long double Sqrt(long double x) { return std::sqrt(x); }
767 //@}
768 
769 
770 //@{
771 /// Return the cube root of a floating-point value.
772 inline float Cbrt(float x) { return std::cbrt(x); }
773 inline double Cbrt(double x) { return std::cbrt(x); }
774 inline long double Cbrt(long double x) { return std::cbrt(x); }
775 //@}
776 
777 
778 //@{
779 /// Return the remainder of @a x / @a y.
780 inline int Mod(int x, int y) { return (x % y); }
781 inline float Mod(float x, float y) { return std::fmod(x, y); }
782 inline double Mod(double x, double y) { return std::fmod(x, y); }
783 inline long double Mod(long double x, long double y) { return std::fmod(x, y); }
784 template<typename Type> inline Type Remainder(Type x, Type y) { return Mod(x, y); }
785 //@}
786 
787 
788 //@{
789 /// Return @a x rounded up to the nearest integer.
790 inline float RoundUp(float x) { return std::ceil(x); }
791 inline double RoundUp(double x) { return std::ceil(x); }
792 inline long double RoundUp(long double x) { return std::ceil(x); }
793 //@}
794 /// Return @a x rounded up to the nearest multiple of @a base.
795 template<typename Type>
796 inline Type
798 {
799  Type remainder = Remainder(x, base);
800  return remainder ? x-remainder+base : x;
801 }
802 
803 
804 //@{
805 /// Return @a x rounded down to the nearest integer.
806 inline float RoundDown(float x) { return std::floor(x); }
807 inline double RoundDown(double x) { return std::floor(x); }
808 inline long double RoundDown(long double x) { return std::floor(x); }
809 //@}
810 /// Return @a x rounded down to the nearest multiple of @a base.
811 template<typename Type>
812 inline Type
814 {
815  Type remainder = Remainder(x, base);
816  return remainder ? x-remainder : x;
817 }
818 
819 
820 //@{
821 /// Return @a x rounded to the nearest integer.
822 inline float Round(float x) { return RoundDown(x + 0.5f); }
823 inline double Round(double x) { return RoundDown(x + 0.5); }
824 inline long double Round(long double x) { return RoundDown(x + 0.5l); }
825 //@}
826 
827 
828 /// Return the euclidean remainder of @a x.
829 /// Note unlike % operator this will always return a positive result
830 template<typename Type>
831 inline Type
832 EuclideanRemainder(Type x) { return x - RoundDown(x); }
833 
834 
835 /// Return the integer part of @a x.
836 template<typename Type>
837 inline Type
839 {
840  return (x > 0 ? RoundDown(x) : RoundUp(x));
841 }
842 
843 /// Return the fractional part of @a x.
844 template<typename Type>
845 inline Type
846 FractionalPart(Type x) { return Mod(x,Type(1)); }
847 
848 
849 //@{
850 /// Return the floor of @a x.
851 inline int Floor(float x) { return int(RoundDown(x)); }
852 inline int Floor(double x) { return int(RoundDown(x)); }
853 inline int Floor(long double x) { return int(RoundDown(x)); }
854 //@}
855 
856 
857 //@{
858 /// Return the ceiling of @a x.
859 inline int Ceil(float x) { return int(RoundUp(x)); }
860 inline int Ceil(double x) { return int(RoundUp(x)); }
861 inline int Ceil(long double x) { return int(RoundUp(x)); }
862 //@}
863 
864 
865 /// Return @a x if it is greater or equal in magnitude than @a delta. Otherwise, return zero.
866 template<typename Type>
867 inline Type Chop(Type x, Type delta) { return (Abs(x) < delta ? zeroVal<Type>() : x); }
868 
869 
870 /// Return @a x truncated to the given number of decimal digits.
871 template<typename Type>
872 inline Type
873 Truncate(Type x, unsigned int digits)
874 {
875  Type tenth = Pow(10,digits);
876  return RoundDown(x*tenth+0.5)/tenth;
877 }
878 
879 ////////////////////////////////////////
880 
881 
882 /// @brief 8-bit integer values print to std::ostreams as characters.
883 /// Cast them so that they print as integers instead.
884 template<typename T>
886  && !std::is_same<T, uint8_t>::value, const T&>::type { return val; }
887 inline int32_t PrintCast(int8_t val) { return int32_t(val); }
888 inline uint32_t PrintCast(uint8_t val) { return uint32_t(val); }
889 
890 
891 ////////////////////////////////////////
892 
893 
894 /// Return the inverse of @a x.
895 template<typename Type>
896 inline Type
898 {
899  assert(x);
900  return Type(1)/x;
901 }
902 
903 
904 enum Axis {
905  X_AXIS = 0,
906  Y_AXIS = 1,
907  Z_AXIS = 2
908 };
909 
910 // enum values are consistent with their historical mx analogs.
920 };
921 
922 
923 template <typename S, typename T>
924 struct promote {
925  using type = typename hboost::numeric::conversion_traits<S, T>::supertype;
926 };
927 
928 
929 /// @brief Return the index [0,1,2] of the smallest value in a 3D vector.
930 /// @note This methods assumes operator[] exists and avoids branching.
931 /// @details If two components of the input vector are equal and smaller than the
932 /// third component, the largest index of the two is always returned.
933 /// If all three vector components are equal the largest index, i.e. 2, is
934 /// returned. In other words the return value corresponds to the largest index
935 /// of the of the smallest vector components.
936 template<typename Vec3T>
937 size_t
938 MinIndex(const Vec3T& v)
939 {
940  static const size_t hashTable[8] = { 2, 1, 9, 1, 2, 9, 0, 0 };//9 is a dummy value
941  const size_t hashKey =
942  ((v[0] < v[1]) << 2) + ((v[0] < v[2]) << 1) + (v[1] < v[2]);// ?*4+?*2+?*1
943  return hashTable[hashKey];
944 }
945 
946 
947 /// @brief Return the index [0,1,2] of the largest value in a 3D vector.
948 /// @note This methods assumes operator[] exists and avoids branching.
949 /// @details If two components of the input vector are equal and larger than the
950 /// third component, the largest index of the two is always returned.
951 /// If all three vector components are equal the largest index, i.e. 2, is
952 /// returned. In other words the return value corresponds to the largest index
953 /// of the largest vector components.
954 template<typename Vec3T>
955 size_t
956 MaxIndex(const Vec3T& v)
957 {
958  static const size_t hashTable[8] = { 2, 1, 9, 1, 2, 9, 0, 0 };//9 is a dummy value
959  const size_t hashKey =
960  ((v[0] > v[1]) << 2) + ((v[0] > v[2]) << 1) + (v[1] > v[2]);// ?*4+?*2+?*1
961  return hashTable[hashKey];
962 }
963 
964 } // namespace math
965 } // namespace OPENVDB_VERSION_NAME
966 } // namespace openvdb
967 
968 #endif // OPENVDB_MATH_MATH_HAS_BEEN_INCLUDED
void setSeed(unsigned int seed)
Set the seed value for the random number generator.
Definition: Math.h:183
int32_t floatToInt32(const float aFloatValue)
Definition: Math.h:484
#define OPENVDB_NO_TYPE_CONVERSION_WARNING_END
Definition: Platform.h:207
SYS_API double cos(double x)
GLboolean GLboolean GLboolean b
Definition: glcorearb.h:1221
SYS_API double fmod(double x, double y)
Type Inv(Type x)
Return the inverse of x.
Definition: Math.h:897
Type Truncate(Type x, unsigned int digits)
Return x truncated to the given number of decimal digits.
Definition: Math.h:873
int Ceil(float x)
Return the ceiling of x.
Definition: Math.h:859
Simple generator of random numbers over the range [0, 1)
Definition: Math.h:165
bool cwiseGreaterThan(const Mat< SIZE, T > &m0, const Mat< SIZE, T > &m1)
Definition: Mat.h:1051
Type Pow(Type x, int n)
Return xn.
Definition: Math.h:564
Type IntegerPart(Type x)
Return the integer part of x.
Definition: Math.h:838
bool ZeroCrossing(const Type &a, const Type &b)
Return true if the interval [a, b] includes zero, i.e., if either a or b is zero or if they have diff...
Definition: Math.h:756
bool isExactlyEqual(const T0 &a, const T1 &b)
Return true if a is exactly equal to b.
Definition: Math.h:444
Type Pow2(Type x)
Return x2.
Definition: Math.h:551
T negative(const T &val)
Return the unary negation of the given value.
Definition: Math.h:127
size_t MaxIndex(const Vec3T &v)
Return the index [0,1,2] of the largest value in a 3D vector.
Definition: Math.h:956
Type FractionalPart(Type x)
Return the fractional part of x.
Definition: Math.h:846
Mat3< Type1 > cwiseAdd(const Mat3< Type1 > &m, const Type2 s)
Definition: Mat3.h:820
const GLfloat * c
Definition: glew.h:16631
GLboolean GLboolean g
Definition: glcorearb.h:1221
vfloat4 sqrt(const vfloat4 &a)
Definition: simd.h:7458
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h:178
const EngineType & engine() const
Return a const reference to the random number generator.
Definition: Math.h:238
float Cos(const float &x)
Return cos x.
Definition: Math.h:728
SYS_API float powf(float x, float y)
Rand01(unsigned int seed)
Initialize the generator.
Definition: Math.h:180
bool SignChange(const Type &a, const Type &b)
Return true if a and b have different signs.
Definition: Math.h:746
float RoundDown(float x)
Return x rounded down to the nearest integer.
Definition: Math.h:806
RandInt(const EngineType &engine, IntType imin, IntType imax)
Initialize the generator.
Definition: Math.h:212
bool zeroVal< bool >()
Return the bool value that corresponds to zero.
Definition: Math.h:73
bool isNegative(const Type &x)
Return true if x is less than zero.
Definition: Math.h:368
int64_t doubleToInt64(const double aDoubleValue)
Definition: Math.h:493
Tolerance for floating-point comparison.
Definition: Math.h:147
float RoundUp(float x)
Return x rounded up to the nearest integer.
Definition: Math.h:790
float Cbrt(float x)
Return the cube root of a floating-point value.
Definition: Math.h:772
Simple random integer generator.
Definition: Math.h:201
typename hboost::numeric::conversion_traits< S, T >::supertype type
Definition: Math.h:925
GLdouble l
Definition: glew.h:9164
GLint GLenum GLint x
Definition: glcorearb.h:408
ImageBuf OIIO_API pow(const ImageBuf &A, cspan< float > B, ROI roi={}, int nthreads=0)
bool isNegative< bool >(const bool &)
Definition: Math.h:371
size_t MinIndex(const Vec3T &v)
Return the index [0,1,2] of the smallest value in a 3D vector.
Definition: Math.h:938
bool isApproxEqual(const Type &a, const Type &b, const Type &tolerance)
Return true if a is equal to b to within the given tolerance.
Definition: Math.h:407
GLint GLint GLsizei GLint GLenum GLenum type
Definition: glcorearb.h:107
ImageBuf OIIO_API min(Image_or_Const A, Image_or_Const B, ROI roi={}, int nthreads=0)
IMATH_INTERNAL_NAMESPACE_HEADER_ENTER T abs(T a)
Definition: ImathFun.h:55
Coord Abs(const Coord &xyz)
Definition: Coord.h:514
float Round(float x)
Return x rounded to the nearest integer.
Definition: Math.h:822
const GLdouble * v
Definition: glcorearb.h:836
bool isInfinite(const float x)
Return true if x is an infinity value (either positive infinity or negative infinity).
Definition: Math.h:386
Type Pow4(Type x)
Return x4.
Definition: Math.h:559
GLboolean GLboolean GLboolean GLboolean a
Definition: glcorearb.h:1221
GLsizei const GLchar *const * string
Definition: glcorearb.h:813
Delta for small floating-point offsets.
Definition: Math.h:154
float Sqrt(float x)
Return the square root of a floating-point value.
Definition: Math.h:764
ImageBuf OIIO_API max(Image_or_Const A, Image_or_Const B, ROI roi={}, int nthreads=0)
bool isNan(const float x)
Return true if x is a NaN (Not-A-Number) value.
Definition: Math.h:396
Type Pow3(Type x)
Return x3.
Definition: Math.h:555
bool cwiseLessThan(const Mat< SIZE, T > &m0, const Mat< SIZE, T > &m1)
Definition: Mat.h:1037
bool isRelOrApproxEqual(const Type &a, const Type &b, const Type &absTol, const Type &relTol)
Definition: Math.h:454
typedef int(WINAPI *PFNWGLRELEASEPBUFFERDCARBPROC)(HPBUFFERARB hPbuffer
#define OPENVDB_NO_FP_EQUALITY_WARNING_BEGIN
Definition: Math.h:47
bool isUlpsEqual(const double aLeft, const double aRight, const int64_t aUnitsInLastPlace)
Definition: Math.h:506
bool isApproxZero(const Type &x)
Return true if x is equal to zero to within the default floating-point comparison tolerance...
Definition: Math.h:350
int floor(T x)
Definition: ImathFun.h:150
Type Clamp(Type x, Type min, Type max)
Return x clamped to [min, max].
Definition: Math.h:260
GLdouble n
Definition: glcorearb.h:2007
#define OPENVDB_EXACT_IS_APPROX_EQUAL(T)
Definition: Math.h:422
RandInt(unsigned int seed, IntType imin, IntType imax)
Initialize the generator.
Definition: Math.h:220
GLuint GLfloat * val
Definition: glcorearb.h:1607
GLfloat GLfloat GLfloat GLfloat h
Definition: glcorearb.h:2001
int Sign(const Type &x)
Return the sign of the given value as an integer (either -1, 0 or 1).
Definition: Math.h:739
Library and file format version numbers.
void setSeed(unsigned int seed)
Set the seed value for the random number generator.
Definition: Math.h:232
#define OPENVDB_NO_TYPE_CONVERSION_WARNING_BEGIN
Bracket code with OPENVDB_NO_TYPE_CONVERSION_WARNING_BEGIN/_END, to inhibit warnings about type conve...
Definition: Platform.h:206
const Type & Min(const Type &a, const Type &b)
Return the minimum of two values.
Definition: Math.h:659
GLsizei const GLfloat * value
Definition: glcorearb.h:823
GLfloat f
Definition: glcorearb.h:1925
const EngineType & engine() const
Return a const reference to the random number generator.
Definition: Math.h:189
Type Chop(Type x, Type delta)
Return x if it is greater or equal in magnitude than delta. Otherwise, return zero.
Definition: Math.h:867
constexpr T pi()
Pi constant taken from Boost to match old behaviour.
Definition: Math.h:118
Vec3< typename promote< T, typename Coord::ValueType >::type > operator+(const Vec3< T > &v0, const Coord &v1)
Allow a Coord to be added to or subtracted from a Vec3.
Definition: Coord.h:524
INT64 INT64 INT64 remainder
Definition: wglew.h:1182
int Floor(float x)
Return the floor of x.
Definition: Math.h:851
bool ClampTest01(Type &x)
Return true if x is outside [0,1].
Definition: Math.h:276
IntType operator()()
Return a randomly-generated integer in the current range.
Definition: Math.h:241
float Sin(const float &x)
Return sin x.
Definition: Math.h:719
#define const
Definition: zconf.h:214
int ceil(T x)
Definition: ImathFun.h:158
T zeroVal()
Return the value of type T that corresponds to zero.
Definition: Math.h:69
Rand01(const EngineType &engine)
Initialize the generator.
Definition: Math.h:176
Type SmoothUnitStep(Type x)
Return 0 if x < 0, 1 if x > 1 or else (3 − 2 x) x .
Definition: Math.h:286
void setRange(IntType imin, IntType imax)
Change the range over which integers are distributed to [imin, imax].
Definition: Math.h:226
GLdouble s
Definition: glew.h:1395
#define OPENVDB_VERSION_NAME
The version namespace name for this library version.
Definition: version.h:114
Type Remainder(Type x, Type y)
Return the remainder of x / y.
Definition: Math.h:784
const Type & Max(const Type &a, const Type &b)
Return the maximum of two values.
Definition: Math.h:598
GLint y
Definition: glcorearb.h:102
SYS_API double sin(double x)
int Mod(int x, int y)
Return the remainder of x / y.
Definition: Math.h:780
auto PrintCast(const T &val) -> typename std::enable_if<!std::is_same< T, int8_t >::value &&!std::is_same< T, uint8_t >::value, const T & >::type
8-bit integer values print to std::ostreams as characters. Cast them so that they print as integers i...
Definition: Math.h:885
bool isZero(const Type &x)
Return true if x is exactly equal to zero.
Definition: Math.h:338
bool isFinite(const float x)
Return true if x is finite.
Definition: Math.h:376
Type Exp(const Type &x)
Return ex.
Definition: Math.h:713
IntType operator()(IntType imin, IntType imax)
Return a randomly-generated integer in the new range [imin, imax], without changing the current range...
Definition: Math.h:245
FloatType operator()()
Return a uniformly distributed random number in the range [0, 1).
Definition: Math.h:192
Type Clamp01(Type x)
Return x clamped to [0, 1].
Definition: Math.h:270
#define OPENVDB_NO_FP_EQUALITY_WARNING_END
Definition: Math.h:48
bool isApproxLarger(const Type &a, const Type &b, const Type &tolerance)
Return true if a is larger than b to within the given tolerance, i.e., if b - a < tolerance...
Definition: Math.h:435