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