VTK  9.2.6
vtkMath.h
Go to the documentation of this file.
1 /*=========================================================================
2 
3  Program: Visualization Toolkit
4  Module: vtkMath.h
5 
6  Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
7  All rights reserved.
8  See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
9 
10  This software is distributed WITHOUT ANY WARRANTY; without even
11  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
12  PURPOSE. See the above copyright notice for more information.
13 
14 =========================================================================
15  Copyright 2011 Sandia Corporation.
16  Under the terms of Contract DE-AC04-94AL85000, there is a non-exclusive
17  license for use of this work by or on behalf of the
18  U.S. Government. Redistribution and use in source and binary forms, with
19  or without modification, are permitted provided that this Notice and any
20  statement of authorship are reproduced on all copies.
21 
22  Contact: pppebay@sandia.gov,dcthomp@sandia.gov
23 
24 =========================================================================*/
154 #ifndef vtkMath_h
155 #define vtkMath_h
156 
157 #include "vtkCommonCoreModule.h" // For export macro
158 #include "vtkMathPrivate.hxx" // For Matrix meta-class helpers
159 #include "vtkMatrixUtilities.h" // For Matrix wrapping / mapping
160 #include "vtkObject.h"
161 #include "vtkSmartPointer.h" // For vtkSmartPointer.
162 #include "vtkTypeTraits.h" // For type traits
163 
164 #include "vtkMathConfigure.h" // For <cmath> and VTK_HAS_ISNAN etc.
165 
166 #include <algorithm> // for std::clamp
167 #include <cassert> // assert() in inline implementations.
168 
169 #ifndef DBL_MIN
170 #define VTK_DBL_MIN 2.2250738585072014e-308
171 #else // DBL_MIN
172 #define VTK_DBL_MIN DBL_MIN
173 #endif // DBL_MIN
174 
175 #ifndef DBL_EPSILON
176 #define VTK_DBL_EPSILON 2.2204460492503131e-16
177 #else // DBL_EPSILON
178 #define VTK_DBL_EPSILON DBL_EPSILON
179 #endif // DBL_EPSILON
180 
181 #ifndef VTK_DBL_EPSILON
182 #ifndef DBL_EPSILON
183 #define VTK_DBL_EPSILON 2.2204460492503131e-16
184 #else // DBL_EPSILON
185 #define VTK_DBL_EPSILON DBL_EPSILON
186 #endif // DBL_EPSILON
187 #endif // VTK_DBL_EPSILON
188 
189 class vtkDataArray;
190 class vtkPoints;
191 class vtkMathInternal;
194 
195 namespace vtk_detail
196 {
197 // forward declaration
198 template <typename OutT>
199 void RoundDoubleToIntegralIfNecessary(double val, OutT* ret);
200 } // end namespace vtk_detail
201 
202 class VTKCOMMONCORE_EXPORT vtkMath : public vtkObject
203 {
204 public:
205  static vtkMath* New();
206  vtkTypeMacro(vtkMath, vtkObject);
207  void PrintSelf(ostream& os, vtkIndent indent) override;
208 
212  static constexpr double Pi() { return 3.141592653589793; }
213 
215 
218  static float RadiansFromDegrees(float degrees);
219  static double RadiansFromDegrees(double degrees);
221 
223 
226  static float DegreesFromRadians(float radians);
227  static double DegreesFromRadians(double radians);
229 
233 #if 1
234  static int Round(float f) { return static_cast<int>(f + (f >= 0.0 ? 0.5 : -0.5)); }
235  static int Round(double f) { return static_cast<int>(f + (f >= 0.0 ? 0.5 : -0.5)); }
236 #endif
237 
242  template <typename OutT>
243  static void RoundDoubleToIntegralIfNecessary(double val, OutT* ret)
244  {
245  // Can't specialize template methods in a template class, so we move the
246  // implementations to a external namespace.
248  }
249 
255  static int Floor(double x);
256 
262  static int Ceil(double x);
263 
269  static int CeilLog2(vtkTypeUInt64 x);
270 
275  template <class T>
276  static T Min(const T& a, const T& b);
277 
282  template <class T>
283  static T Max(const T& a, const T& b);
284 
288  static bool IsPowerOfTwo(vtkTypeUInt64 x);
289 
295  static int NearestPowerOfTwo(int x);
296 
301  static vtkTypeInt64 Factorial(int N);
302 
308  static vtkTypeInt64 Binomial(int m, int n);
309 
321  static int* BeginCombination(int m, int n);
322 
333  static int NextCombination(int m, int n, int* combination);
334 
338  static void FreeCombination(int* combination);
339 
355  static void RandomSeed(int s);
356 
368  static int GetSeed();
369 
383  static double Random();
384 
397  static double Random(double min, double max);
398 
411  static double Gaussian();
412 
425  static double Gaussian(double mean, double std);
426 
431  template <class VectorT1, class VectorT2>
432  static void Assign(const VectorT1& a, VectorT2&& b)
433  {
434  b[0] = a[0];
435  b[1] = a[1];
436  b[2] = a[2];
437  }
438 
442  static void Assign(const double a[3], double b[3]) { vtkMath::Assign<>(a, b); }
443 
447  static void Add(const float a[3], const float b[3], float c[3])
448  {
449  for (int i = 0; i < 3; ++i)
450  {
451  c[i] = a[i] + b[i];
452  }
453  }
454 
458  static void Add(const double a[3], const double b[3], double c[3])
459  {
460  for (int i = 0; i < 3; ++i)
461  {
462  c[i] = a[i] + b[i];
463  }
464  }
465 
469  static void Subtract(const float a[3], const float b[3], float c[3])
470  {
471  for (int i = 0; i < 3; ++i)
472  {
473  c[i] = a[i] - b[i];
474  }
475  }
476 
480  static void Subtract(const double a[3], const double b[3], double c[3])
481  {
482  for (int i = 0; i < 3; ++i)
483  {
484  c[i] = a[i] - b[i];
485  }
486  }
487 
493  template <class VectorT1, class VectorT2, class VectorT3>
494  static void Subtract(const VectorT1& a, const VectorT2& b, VectorT3&& c)
495  {
496  c[0] = a[0] - b[0];
497  c[1] = a[1] - b[1];
498  c[2] = a[2] - b[2];
499  }
500 
505  static void MultiplyScalar(float a[3], float s)
506  {
507  for (int i = 0; i < 3; ++i)
508  {
509  a[i] *= s;
510  }
511  }
512 
517  static void MultiplyScalar2D(float a[2], float s)
518  {
519  for (int i = 0; i < 2; ++i)
520  {
521  a[i] *= s;
522  }
523  }
524 
529  static void MultiplyScalar(double a[3], double s)
530  {
531  for (int i = 0; i < 3; ++i)
532  {
533  a[i] *= s;
534  }
535  }
536 
541  static void MultiplyScalar2D(double a[2], double s)
542  {
543  for (int i = 0; i < 2; ++i)
544  {
545  a[i] *= s;
546  }
547  }
548 
552  static float Dot(const float a[3], const float b[3])
553  {
554  return a[0] * b[0] + a[1] * b[1] + a[2] * b[2];
555  }
556 
560  static double Dot(const double a[3], const double b[3])
561  {
562  return a[0] * b[0] + a[1] * b[1] + a[2] * b[2];
563  }
564 
580  template <typename ReturnTypeT = double, typename TupleRangeT1, typename TupleRangeT2,
581  typename EnableT = typename std::conditional<!std::is_pointer<TupleRangeT1>::value &&
583  TupleRangeT1, TupleRangeT2>::type::value_type>
584  static ReturnTypeT Dot(const TupleRangeT1& a, const TupleRangeT2& b)
585  {
586  return a[0] * b[0] + a[1] * b[1] + a[2] * b[2];
587  }
588 
592  static void Outer(const float a[3], const float b[3], float c[3][3])
593  {
594  for (int i = 0; i < 3; ++i)
595  {
596  for (int j = 0; j < 3; ++j)
597  {
598  c[i][j] = a[i] * b[j];
599  }
600  }
601  }
602 
606  static void Outer(const double a[3], const double b[3], double c[3][3])
607  {
608  for (int i = 0; i < 3; ++i)
609  {
610  for (int j = 0; j < 3; ++j)
611  {
612  c[i][j] = a[i] * b[j];
613  }
614  }
615  }
616 
621  static void Cross(const float a[3], const float b[3], float c[3]);
622 
627  static void Cross(const double a[3], const double b[3], double c[3]);
628 
630 
633  static float Norm(const float* x, int n);
634  static double Norm(const double* x, int n);
636 
640  static float Norm(const float v[3]) { return std::sqrt(v[0] * v[0] + v[1] * v[1] + v[2] * v[2]); }
641 
645  static double Norm(const double v[3])
646  {
647  return std::sqrt(v[0] * v[0] + v[1] * v[1] + v[2] * v[2]);
648  }
649 
659  template <typename ReturnTypeT = double, typename TupleRangeT>
660  static ReturnTypeT SquaredNorm(const TupleRangeT& v)
661  {
662  return v[0] * v[0] + v[1] * v[1] + v[2] * v[2];
663  }
664 
669  static float Normalize(float v[3]);
670 
675  static double Normalize(double v[3]);
676 
678 
685  static void Perpendiculars(const double v1[3], double v2[3], double v3[3], double theta);
686  static void Perpendiculars(const float v1[3], float v2[3], float v3[3], double theta);
688 
690 
695  static bool ProjectVector(const float a[3], const float b[3], float projection[3]);
696  static bool ProjectVector(const double a[3], const double b[3], double projection[3]);
698 
700 
706  static bool ProjectVector2D(const float a[2], const float b[2], float projection[2]);
707  static bool ProjectVector2D(const double a[2], const double b[2], double projection[2]);
709 
725  template <typename ReturnTypeT = double, typename TupleRangeT1, typename TupleRangeT2,
726  typename EnableT = typename std::conditional<!std::is_pointer<TupleRangeT1>::value &&
728  TupleRangeT1, TupleRangeT2>::type::value_type>
729  static ReturnTypeT Distance2BetweenPoints(const TupleRangeT1& p1, const TupleRangeT2& p2);
730 
735  static float Distance2BetweenPoints(const float p1[3], const float p2[3]);
736 
741  static double Distance2BetweenPoints(const double p1[3], const double p2[3]);
742 
746  static double AngleBetweenVectors(const double v1[3], const double v2[3]);
747 
752  const double v1[3], const double v2[3], const double vn[3]);
753 
758  static double GaussianAmplitude(const double variance, const double distanceFromMean);
759 
764  static double GaussianAmplitude(const double mean, const double variance, const double position);
765 
771  static double GaussianWeight(const double variance, const double distanceFromMean);
772 
778  static double GaussianWeight(const double mean, const double variance, const double position);
779 
783  static float Dot2D(const float x[2], const float y[2]) { return x[0] * y[0] + x[1] * y[1]; }
784 
788  static double Dot2D(const double x[2], const double y[2]) { return x[0] * y[0] + x[1] * y[1]; }
789 
793  static void Outer2D(const float x[2], const float y[2], float A[2][2])
794  {
795  for (int i = 0; i < 2; ++i)
796  {
797  for (int j = 0; j < 2; ++j)
798  {
799  A[i][j] = x[i] * y[j];
800  }
801  }
802  }
803 
807  static void Outer2D(const double x[2], const double y[2], double A[2][2])
808  {
809  for (int i = 0; i < 2; ++i)
810  {
811  for (int j = 0; j < 2; ++j)
812  {
813  A[i][j] = x[i] * y[j];
814  }
815  }
816  }
817 
822  static float Norm2D(const float x[2]) { return std::sqrt(x[0] * x[0] + x[1] * x[1]); }
823 
828  static double Norm2D(const double x[2]) { return std::sqrt(x[0] * x[0] + x[1] * x[1]); }
829 
834  static float Normalize2D(float v[2]);
835 
840  static double Normalize2D(double v[2]);
841 
845  static float Determinant2x2(const float c1[2], const float c2[2])
846  {
847  return c1[0] * c2[1] - c2[0] * c1[1];
848  }
849 
851 
854  static double Determinant2x2(double a, double b, double c, double d) { return a * d - b * c; }
855  static double Determinant2x2(const double c1[2], const double c2[2])
856  {
857  return c1[0] * c2[1] - c2[0] * c1[1];
858  }
860 
862 
865  static void LUFactor3x3(float A[3][3], int index[3]);
866  static void LUFactor3x3(double A[3][3], int index[3]);
868 
870 
873  static void LUSolve3x3(const float A[3][3], const int index[3], float x[3]);
874  static void LUSolve3x3(const double A[3][3], const int index[3], double x[3]);
876 
878 
882  static void LinearSolve3x3(const float A[3][3], const float x[3], float y[3]);
883  static void LinearSolve3x3(const double A[3][3], const double x[3], double y[3]);
885 
887 
890  static void Multiply3x3(const float A[3][3], const float v[3], float u[3]);
891  static void Multiply3x3(const double A[3][3], const double v[3], double u[3]);
893 
895 
898  static void Multiply3x3(const float A[3][3], const float B[3][3], float C[3][3]);
899  static void Multiply3x3(const double A[3][3], const double B[3][3], double C[3][3]);
901 
925  template <int RowsT, int MidDimT, int ColsT,
926  class LayoutT1 = vtkMatrixUtilities::Layout::Identity,
927  class LayoutT2 = vtkMatrixUtilities::Layout::Identity, class MatrixT1, class MatrixT2,
928  class MatrixT3>
929  static void MultiplyMatrix(const MatrixT1& M1, const MatrixT2& M2, MatrixT3&& M3)
930  {
931  vtkMathPrivate::MultiplyMatrix<RowsT, MidDimT, ColsT, LayoutT1, LayoutT2>::Compute(M1, M2, M3);
932  }
933 
954  template <int RowsT, int ColsT, class LayoutT = vtkMatrixUtilities::Layout::Identity,
955  class MatrixT, class VectorT1, class VectorT2>
956  static void MultiplyMatrixWithVector(const MatrixT& M, const VectorT1& X, VectorT2&& Y)
957  {
958  vtkMathPrivate::MultiplyMatrix<RowsT, ColsT, 1, LayoutT>::Compute(M, X, Y);
959  }
960 
966  template <class ScalarT, int SizeT, class VectorT1, class VectorT2>
967  static ScalarT Dot(const VectorT1& x, const VectorT2& y)
968  {
969  return vtkMathPrivate::ContractRowWithCol<ScalarT, 1, SizeT, 1, 0, 0,
970  vtkMatrixUtilities::Layout::Identity, vtkMatrixUtilities::Layout::Transpose>::Compute(x, y);
971  }
972 
989  template <int SizeT, class LayoutT = vtkMatrixUtilities::Layout::Identity, class MatrixT>
991  const MatrixT& M)
992  {
993  return vtkMathPrivate::Determinant<SizeT, LayoutT>::Compute(M);
994  }
995 
1011  template <int SizeT, class LayoutT = vtkMatrixUtilities::Layout::Identity, class MatrixT1,
1012  class MatrixT2>
1013  static void InvertMatrix(const MatrixT1& M1, MatrixT2&& M2)
1014  {
1015  vtkMathPrivate::InvertMatrix<SizeT, LayoutT>::Compute(M1, M2);
1016  }
1017 
1031  template <int RowsT, int ColsT, class LayoutT = vtkMatrixUtilities::Layout::Identity,
1032  class MatrixT, class VectorT1, class VectorT2>
1033  static void LinearSolve(const MatrixT& M, const VectorT1& x, VectorT2& y)
1034  {
1035  vtkMathPrivate::LinearSolve<RowsT, ColsT, LayoutT>::Compute(M, x, y);
1036  }
1037 
1052  template <class ScalarT, int SizeT, class LayoutT = vtkMatrixUtilities::Layout::Identity,
1053  class VectorT1, class MatrixT, class VectorT2>
1054  static ScalarT Dot(const VectorT1& x, const MatrixT& M, const VectorT2& y)
1055  {
1056  ScalarT tmp[SizeT];
1057  vtkMathPrivate::MultiplyMatrix<SizeT, SizeT, 1, LayoutT>::Compute(M, y, tmp);
1058  return vtkMathPrivate::ContractRowWithCol<ScalarT, 1, SizeT, 1, 0, 0,
1059  vtkMatrixUtilities::Layout::Identity, vtkMatrixUtilities::Layout::Transpose>::Compute(x, tmp);
1060  }
1061 
1067  static void MultiplyMatrix(const double* const* A, const double* const* B, unsigned int rowA,
1068  unsigned int colA, unsigned int rowB, unsigned int colB, double** C);
1069 
1071 
1075  static void Transpose3x3(const float A[3][3], float AT[3][3]);
1076  static void Transpose3x3(const double A[3][3], double AT[3][3]);
1078 
1080 
1084  static void Invert3x3(const float A[3][3], float AI[3][3]);
1085  static void Invert3x3(const double A[3][3], double AI[3][3]);
1087 
1089 
1092  static void Identity3x3(float A[3][3]);
1093  static void Identity3x3(double A[3][3]);
1095 
1097 
1100  static double Determinant3x3(const float A[3][3]);
1101  static double Determinant3x3(const double A[3][3]);
1103 
1107  static float Determinant3x3(const float c1[3], const float c2[3], const float c3[3]);
1108 
1112  static double Determinant3x3(const double c1[3], const double c2[3], const double c3[3]);
1113 
1120  static double Determinant3x3(double a1, double a2, double a3, double b1, double b2, double b3,
1121  double c1, double c2, double c3);
1122 
1124 
1131  static void QuaternionToMatrix3x3(const float quat[4], float A[3][3]);
1132  static void QuaternionToMatrix3x3(const double quat[4], double A[3][3]);
1133  template <class QuaternionT, class MatrixT,
1134  class EnableT = typename std::enable_if<!vtkMatrixUtilities::MatrixIs2DArray<MatrixT>()>::type>
1135  static void QuaternionToMatrix3x3(const QuaternionT& q, MatrixT&& A);
1137 
1139 
1148  static void Matrix3x3ToQuaternion(const float A[3][3], float quat[4]);
1149  static void Matrix3x3ToQuaternion(const double A[3][3], double quat[4]);
1150  template <class MatrixT, class QuaternionT,
1151  class EnableT = typename std::enable_if<!vtkMatrixUtilities::MatrixIs2DArray<MatrixT>()>::type>
1152  static void Matrix3x3ToQuaternion(const MatrixT& A, QuaternionT&& q);
1154 
1156 
1162  static void MultiplyQuaternion(const float q1[4], const float q2[4], float q[4]);
1163  static void MultiplyQuaternion(const double q1[4], const double q2[4], double q[4]);
1165 
1167 
1171  static void RotateVectorByNormalizedQuaternion(const float v[3], const float q[4], float r[3]);
1172  static void RotateVectorByNormalizedQuaternion(const double v[3], const double q[4], double r[3]);
1174 
1176 
1180  static void RotateVectorByWXYZ(const float v[3], const float q[4], float r[3]);
1181  static void RotateVectorByWXYZ(const double v[3], const double q[4], double r[3]);
1183 
1185 
1190  static void Orthogonalize3x3(const float A[3][3], float B[3][3]);
1191  static void Orthogonalize3x3(const double A[3][3], double B[3][3]);
1193 
1195 
1201  static void Diagonalize3x3(const float A[3][3], float w[3], float V[3][3]);
1202  static void Diagonalize3x3(const double A[3][3], double w[3], double V[3][3]);
1204 
1206 
1216  const float A[3][3], float U[3][3], float w[3], float VT[3][3]);
1218  const double A[3][3], double U[3][3], double w[3], double VT[3][3]);
1220 
1229  double a00, double a01, double a10, double a11, double b0, double b1, double& x0, double& x1);
1230 
1239  static vtkTypeBool SolveLinearSystem(double** A, double* x, int size);
1240 
1247  static vtkTypeBool InvertMatrix(double** A, double** AI, int size);
1248 
1255  double** A, double** AI, int size, int* tmp1Size, double* tmp2Size);
1256 
1279  static vtkTypeBool LUFactorLinearSystem(double** A, int* index, int size);
1280 
1286  static vtkTypeBool LUFactorLinearSystem(double** A, int* index, int size, double* tmpSize);
1287 
1296  static void LUSolveLinearSystem(double** A, int* index, double* x, int size);
1297 
1306  static double EstimateMatrixCondition(const double* const* A, int size);
1307 
1309 
1317  static vtkTypeBool Jacobi(float** a, float* w, float** v);
1318  static vtkTypeBool Jacobi(double** a, double* w, double** v);
1320 
1322 
1331  static vtkTypeBool JacobiN(float** a, int n, float* w, float** v);
1332  static vtkTypeBool JacobiN(double** a, int n, double* w, double** v);
1334 
1349  int numberOfSamples, double** xt, int xOrder, double** mt);
1350 
1365  static vtkTypeBool SolveLeastSquares(int numberOfSamples, double** xt, int xOrder, double** yt,
1366  int yOrder, double** mt, int checkHomogeneous = 1);
1367 
1369 
1376  static void RGBToHSV(const float rgb[3], float hsv[3])
1377  {
1378  RGBToHSV(rgb[0], rgb[1], rgb[2], hsv, hsv + 1, hsv + 2);
1379  }
1380  static void RGBToHSV(float r, float g, float b, float* h, float* s, float* v);
1381  static void RGBToHSV(const double rgb[3], double hsv[3])
1382  {
1383  RGBToHSV(rgb[0], rgb[1], rgb[2], hsv, hsv + 1, hsv + 2);
1384  }
1385  static void RGBToHSV(double r, double g, double b, double* h, double* s, double* v);
1387 
1389 
1396  static void HSVToRGB(const float hsv[3], float rgb[3])
1397  {
1398  HSVToRGB(hsv[0], hsv[1], hsv[2], rgb, rgb + 1, rgb + 2);
1399  }
1400  static void HSVToRGB(float h, float s, float v, float* r, float* g, float* b);
1401  static void HSVToRGB(const double hsv[3], double rgb[3])
1402  {
1403  HSVToRGB(hsv[0], hsv[1], hsv[2], rgb, rgb + 1, rgb + 2);
1404  }
1405  static void HSVToRGB(double h, double s, double v, double* r, double* g, double* b);
1407 
1409 
1412  static void LabToXYZ(const double lab[3], double xyz[3])
1413  {
1414  LabToXYZ(lab[0], lab[1], lab[2], xyz + 0, xyz + 1, xyz + 2);
1415  }
1416  static void LabToXYZ(double L, double a, double b, double* x, double* y, double* z);
1418 
1420 
1423  static void XYZToLab(const double xyz[3], double lab[3])
1424  {
1425  XYZToLab(xyz[0], xyz[1], xyz[2], lab + 0, lab + 1, lab + 2);
1426  }
1427  static void XYZToLab(double x, double y, double z, double* L, double* a, double* b);
1429 
1431 
1434  static void XYZToRGB(const double xyz[3], double rgb[3])
1435  {
1436  XYZToRGB(xyz[0], xyz[1], xyz[2], rgb + 0, rgb + 1, rgb + 2);
1437  }
1438  static void XYZToRGB(double x, double y, double z, double* r, double* g, double* b);
1440 
1442 
1445  static void RGBToXYZ(const double rgb[3], double xyz[3])
1446  {
1447  RGBToXYZ(rgb[0], rgb[1], rgb[2], xyz + 0, xyz + 1, xyz + 2);
1448  }
1449  static void RGBToXYZ(double r, double g, double b, double* x, double* y, double* z);
1451 
1453 
1459  static void RGBToLab(const double rgb[3], double lab[3])
1460  {
1461  RGBToLab(rgb[0], rgb[1], rgb[2], lab + 0, lab + 1, lab + 2);
1462  }
1463  static void RGBToLab(double red, double green, double blue, double* L, double* a, double* b);
1465 
1467 
1470  static void LabToRGB(const double lab[3], double rgb[3])
1471  {
1472  LabToRGB(lab[0], lab[1], lab[2], rgb + 0, rgb + 1, rgb + 2);
1473  }
1474  static void LabToRGB(double L, double a, double b, double* red, double* green, double* blue);
1476 
1478 
1481  static void UninitializeBounds(double bounds[6])
1482  {
1483  bounds[0] = 1.0;
1484  bounds[1] = -1.0;
1485  bounds[2] = 1.0;
1486  bounds[3] = -1.0;
1487  bounds[4] = 1.0;
1488  bounds[5] = -1.0;
1489  }
1491 
1493 
1496  static vtkTypeBool AreBoundsInitialized(const double bounds[6])
1497  {
1498  if (bounds[1] - bounds[0] < 0.0)
1499  {
1500  return 0;
1501  }
1502  return 1;
1503  }
1505 
1510  template <class T>
1511  static T ClampValue(const T& value, const T& min, const T& max);
1512 
1514 
1518  static void ClampValue(double* value, const double range[2]);
1519  static void ClampValue(double value, const double range[2], double* clamped_value);
1520  static void ClampValues(double* values, int nb_values, const double range[2]);
1521  static void ClampValues(
1522  const double* values, int nb_values, const double range[2], double* clamped_values);
1524 
1531  static double ClampAndNormalizeValue(double value, const double range[2]);
1532 
1537  template <class T1, class T2>
1538  static void TensorFromSymmetricTensor(const T1 symmTensor[6], T2 tensor[9]);
1539 
1545  template <class T>
1546  static void TensorFromSymmetricTensor(T tensor[9]);
1547 
1557  double range_min, double range_max, double scale = 1.0, double shift = 0.0);
1558 
1567  static vtkTypeBool GetAdjustedScalarRange(vtkDataArray* array, int comp, double range[2]);
1568 
1573  static vtkTypeBool ExtentIsWithinOtherExtent(const int extent1[6], const int extent2[6]);
1574 
1581  const double bounds1[6], const double bounds2[6], const double delta[3]);
1582 
1589  const double point[3], const double bounds[6], const double delta[3]);
1590 
1601  const double bounds[6], const double normal[3], const double point[3]);
1602 
1612  static double Solve3PointCircle(
1613  const double p1[3], const double p2[3], const double p3[3], double center[3]);
1614 
1618  static double Inf();
1619 
1623  static double NegInf();
1624 
1628  static double Nan();
1629 
1633  static vtkTypeBool IsInf(double x);
1634 
1638  static vtkTypeBool IsNan(double x);
1639 
1644  static bool IsFinite(double x);
1645 
1650  static int QuadraticRoot(double a, double b, double c, double min, double max, double* u);
1651 
1652  enum class ConvolutionMode
1653  {
1654  FULL,
1655  SAME,
1656  VALID
1657  };
1658 
1681  template <class Iter1, class Iter2, class Iter3>
1682  static void Convolve1D(Iter1 beginSample, Iter1 endSample, Iter2 beginKernel, Iter2 endKernel,
1683  Iter3 beginOut, Iter3 endOut, ConvolutionMode mode = ConvolutionMode::FULL)
1684  {
1685  int sampleSize = std::distance(beginSample, endSample);
1686  int kernelSize = std::distance(beginKernel, endKernel);
1687  int outSize = std::distance(beginOut, endOut);
1688 
1689  if (sampleSize <= 0 || kernelSize <= 0 || outSize <= 0)
1690  {
1691  return;
1692  }
1693 
1694  int begin = 0;
1695  int end = outSize;
1696 
1697  switch (mode)
1698  {
1699  case ConvolutionMode::SAME:
1700  begin = static_cast<int>(std::ceil(std::min(sampleSize, kernelSize) / 2.0)) - 1;
1701  end = begin + std::max(sampleSize, kernelSize);
1702  break;
1703  case ConvolutionMode::VALID:
1704  begin = std::min(sampleSize, kernelSize) - 1;
1705  end = begin + std::abs(sampleSize - kernelSize) + 1;
1706  break;
1707  case ConvolutionMode::FULL:
1708  default:
1709  break;
1710  }
1711 
1712  for (int i = begin; i < end; i++)
1713  {
1714  Iter3 out = beginOut + i - begin;
1715  *out = 0;
1716  for (int j = std::max(i - sampleSize + 1, 0); j <= std::min(i, kernelSize - 1); j++)
1717  {
1718  *out += *(beginSample + (i - j)) * *(beginKernel + j);
1719  }
1720  }
1721  }
1722 
1723 protected:
1724  vtkMath() = default;
1725  ~vtkMath() override = default;
1726 
1728 
1729 private:
1730  vtkMath(const vtkMath&) = delete;
1731  void operator=(const vtkMath&) = delete;
1732 };
1733 
1734 //----------------------------------------------------------------------------
1735 inline float vtkMath::RadiansFromDegrees(float x)
1736 {
1737  return x * 0.017453292f;
1738 }
1739 
1740 //----------------------------------------------------------------------------
1741 inline double vtkMath::RadiansFromDegrees(double x)
1742 {
1743  return x * 0.017453292519943295;
1744 }
1745 
1746 //----------------------------------------------------------------------------
1747 inline float vtkMath::DegreesFromRadians(float x)
1748 {
1749  return x * 57.2957795131f;
1750 }
1751 
1752 //----------------------------------------------------------------------------
1753 inline double vtkMath::DegreesFromRadians(double x)
1754 {
1755  return x * 57.29577951308232;
1756 }
1757 
1758 //----------------------------------------------------------------------------
1759 inline bool vtkMath::IsPowerOfTwo(vtkTypeUInt64 x)
1760 {
1761  return ((x != 0) & ((x & (x - 1)) == 0));
1762 }
1763 
1764 //----------------------------------------------------------------------------
1765 // Credit goes to Peter Hart and William Lewis on comp.lang.python 1997
1767 {
1768  unsigned int z = static_cast<unsigned int>(((x > 0) ? x - 1 : 0));
1769  z |= z >> 1;
1770  z |= z >> 2;
1771  z |= z >> 4;
1772  z |= z >> 8;
1773  z |= z >> 16;
1774  return static_cast<int>(z + 1);
1775 }
1776 
1777 //----------------------------------------------------------------------------
1778 // Modify the trunc() operation provided by static_cast<int>() to get floor(),
1779 // Note that in C++ conditions evaluate to values of 1 or 0 (true or false).
1780 inline int vtkMath::Floor(double x)
1781 {
1782  int i = static_cast<int>(x);
1783  return i - (i > x);
1784 }
1785 
1786 //----------------------------------------------------------------------------
1787 // Modify the trunc() operation provided by static_cast<int>() to get ceil(),
1788 // Note that in C++ conditions evaluate to values of 1 or 0 (true or false).
1789 inline int vtkMath::Ceil(double x)
1790 {
1791  int i = static_cast<int>(x);
1792  return i + (i < x);
1793 }
1794 
1795 //----------------------------------------------------------------------------
1796 template <class T>
1797 inline T vtkMath::Min(const T& a, const T& b)
1798 {
1799  return (b <= a ? b : a);
1800 }
1801 
1802 //----------------------------------------------------------------------------
1803 template <class T>
1804 inline T vtkMath::Max(const T& a, const T& b)
1805 {
1806  return (b > a ? b : a);
1807 }
1808 
1809 //----------------------------------------------------------------------------
1810 inline float vtkMath::Normalize(float v[3])
1811 {
1812  float den = vtkMath::Norm(v);
1813  if (den != 0.0)
1814  {
1815  for (int i = 0; i < 3; ++i)
1816  {
1817  v[i] /= den;
1818  }
1819  }
1820  return den;
1821 }
1822 
1823 //----------------------------------------------------------------------------
1824 inline double vtkMath::Normalize(double v[3])
1825 {
1826  double den = vtkMath::Norm(v);
1827  if (den != 0.0)
1828  {
1829  for (int i = 0; i < 3; ++i)
1830  {
1831  v[i] /= den;
1832  }
1833  }
1834  return den;
1835 }
1836 
1837 //----------------------------------------------------------------------------
1838 inline float vtkMath::Normalize2D(float v[2])
1839 {
1840  float den = vtkMath::Norm2D(v);
1841  if (den != 0.0)
1842  {
1843  for (int i = 0; i < 2; ++i)
1844  {
1845  v[i] /= den;
1846  }
1847  }
1848  return den;
1849 }
1850 
1851 //----------------------------------------------------------------------------
1852 inline double vtkMath::Normalize2D(double v[2])
1853 {
1854  double den = vtkMath::Norm2D(v);
1855  if (den != 0.0)
1856  {
1857  for (int i = 0; i < 2; ++i)
1858  {
1859  v[i] /= den;
1860  }
1861  }
1862  return den;
1863 }
1864 
1865 //----------------------------------------------------------------------------
1866 inline float vtkMath::Determinant3x3(const float c1[3], const float c2[3], const float c3[3])
1867 {
1868  return c1[0] * c2[1] * c3[2] + c2[0] * c3[1] * c1[2] + c3[0] * c1[1] * c2[2] -
1869  c1[0] * c3[1] * c2[2] - c2[0] * c1[1] * c3[2] - c3[0] * c2[1] * c1[2];
1870 }
1871 
1872 //----------------------------------------------------------------------------
1873 inline double vtkMath::Determinant3x3(const double c1[3], const double c2[3], const double c3[3])
1874 {
1875  return c1[0] * c2[1] * c3[2] + c2[0] * c3[1] * c1[2] + c3[0] * c1[1] * c2[2] -
1876  c1[0] * c3[1] * c2[2] - c2[0] * c1[1] * c3[2] - c3[0] * c2[1] * c1[2];
1877 }
1878 
1879 //----------------------------------------------------------------------------
1881  double a1, double a2, double a3, double b1, double b2, double b3, double c1, double c2, double c3)
1882 {
1883  return (a1 * vtkMath::Determinant2x2(b2, b3, c2, c3) -
1884  b1 * vtkMath::Determinant2x2(a2, a3, c2, c3) + c1 * vtkMath::Determinant2x2(a2, a3, b2, b3));
1885 }
1886 
1887 //----------------------------------------------------------------------------
1888 inline float vtkMath::Distance2BetweenPoints(const float p1[3], const float p2[3])
1889 {
1890  return ((p1[0] - p2[0]) * (p1[0] - p2[0]) + (p1[1] - p2[1]) * (p1[1] - p2[1]) +
1891  (p1[2] - p2[2]) * (p1[2] - p2[2]));
1892 }
1893 
1894 //----------------------------------------------------------------------------
1895 inline double vtkMath::Distance2BetweenPoints(const double p1[3], const double p2[3])
1896 {
1897  return ((p1[0] - p2[0]) * (p1[0] - p2[0]) + (p1[1] - p2[1]) * (p1[1] - p2[1]) +
1898  (p1[2] - p2[2]) * (p1[2] - p2[2]));
1899 }
1900 
1901 //----------------------------------------------------------------------------
1902 template <typename ReturnTypeT, typename TupleRangeT1, typename TupleRangeT2, typename EnableT>
1903 inline ReturnTypeT vtkMath::Distance2BetweenPoints(const TupleRangeT1& p1, const TupleRangeT2& p2)
1904 {
1905  return ((p1[0] - p2[0]) * (p1[0] - p2[0]) + (p1[1] - p2[1]) * (p1[1] - p2[1]) +
1906  (p1[2] - p2[2]) * (p1[2] - p2[2]));
1907 }
1908 
1909 //----------------------------------------------------------------------------
1910 // Cross product of two 3-vectors. Result (a x b) is stored in c[3].
1911 inline void vtkMath::Cross(const float a[3], const float b[3], float c[3])
1912 {
1913  float Cx = a[1] * b[2] - a[2] * b[1];
1914  float Cy = a[2] * b[0] - a[0] * b[2];
1915  float Cz = a[0] * b[1] - a[1] * b[0];
1916  c[0] = Cx;
1917  c[1] = Cy;
1918  c[2] = Cz;
1919 }
1920 
1921 //----------------------------------------------------------------------------
1922 // Cross product of two 3-vectors. Result (a x b) is stored in c[3].
1923 inline void vtkMath::Cross(const double a[3], const double b[3], double c[3])
1924 {
1925  double Cx = a[1] * b[2] - a[2] * b[1];
1926  double Cy = a[2] * b[0] - a[0] * b[2];
1927  double Cz = a[0] * b[1] - a[1] * b[0];
1928  c[0] = Cx;
1929  c[1] = Cy;
1930  c[2] = Cz;
1931 }
1932 
1933 //----------------------------------------------------------------------------
1934 template <class T>
1935 inline double vtkDeterminant3x3(const T A[3][3])
1936 {
1937  return A[0][0] * A[1][1] * A[2][2] + A[1][0] * A[2][1] * A[0][2] + A[2][0] * A[0][1] * A[1][2] -
1938  A[0][0] * A[2][1] * A[1][2] - A[1][0] * A[0][1] * A[2][2] - A[2][0] * A[1][1] * A[0][2];
1939 }
1940 
1941 //----------------------------------------------------------------------------
1942 inline double vtkMath::Determinant3x3(const float A[3][3])
1943 {
1944  return vtkDeterminant3x3(A);
1945 }
1946 
1947 //----------------------------------------------------------------------------
1948 inline double vtkMath::Determinant3x3(const double A[3][3])
1949 {
1950  return vtkDeterminant3x3(A);
1951 }
1952 
1953 //----------------------------------------------------------------------------
1954 template <class T>
1955 inline T vtkMath::ClampValue(const T& value, const T& min, const T& max)
1956 {
1957  assert("pre: valid_range" && min <= max);
1958 
1959 #if __cplusplus >= 201703L
1960  return std::clamp(value, min, max);
1961 #else
1962  // compilers are good at optimizing the ternary operator,
1963  // use '<' since it is preferred by STL for custom types
1964  T v = (min < value ? value : min);
1965  return (v < max ? v : max);
1966 #endif
1967 }
1968 
1969 //----------------------------------------------------------------------------
1970 inline void vtkMath::ClampValue(double* value, const double range[2])
1971 {
1972  if (value && range)
1973  {
1974  assert("pre: valid_range" && range[0] <= range[1]);
1975 
1976  *value = vtkMath::ClampValue(*value, range[0], range[1]);
1977  }
1978 }
1979 
1980 //----------------------------------------------------------------------------
1981 inline void vtkMath::ClampValue(double value, const double range[2], double* clamped_value)
1982 {
1983  if (range && clamped_value)
1984  {
1985  assert("pre: valid_range" && range[0] <= range[1]);
1986 
1987  *clamped_value = vtkMath::ClampValue(value, range[0], range[1]);
1988  }
1989 }
1990 
1991 // ---------------------------------------------------------------------------
1992 inline double vtkMath::ClampAndNormalizeValue(double value, const double range[2])
1993 {
1994  assert("pre: valid_range" && range[0] <= range[1]);
1995 
1996  double result;
1997  if (range[0] == range[1])
1998  {
1999  result = 0.0;
2000  }
2001  else
2002  {
2003  // clamp
2004  result = vtkMath::ClampValue(value, range[0], range[1]);
2005 
2006  // normalize
2007  result = (result - range[0]) / (range[1] - range[0]);
2008  }
2009 
2010  assert("post: valid_result" && result >= 0.0 && result <= 1.0);
2011 
2012  return result;
2013 }
2014 
2015 //-----------------------------------------------------------------------------
2016 template <class T1, class T2>
2017 inline void vtkMath::TensorFromSymmetricTensor(const T1 symmTensor[9], T2 tensor[9])
2018 {
2019  for (int i = 0; i < 3; ++i)
2020  {
2021  tensor[4 * i] = symmTensor[i];
2022  }
2023  tensor[1] = tensor[3] = symmTensor[3];
2024  tensor[2] = tensor[6] = symmTensor[5];
2025  tensor[5] = tensor[7] = symmTensor[4];
2026 }
2027 
2028 //-----------------------------------------------------------------------------
2029 template <class T>
2030 inline void vtkMath::TensorFromSymmetricTensor(T tensor[9])
2031 {
2032  tensor[6] = tensor[5]; // XZ
2033  tensor[7] = tensor[4]; // YZ
2034  tensor[8] = tensor[2]; // ZZ
2035  tensor[4] = tensor[1]; // YY
2036  tensor[5] = tensor[7]; // YZ
2037  tensor[2] = tensor[6]; // XZ
2038  tensor[1] = tensor[3]; // XY
2039 }
2040 
2041 namespace
2042 {
2043 template <class QuaternionT, class MatrixT>
2044 inline void vtkQuaternionToMatrix3x3(const QuaternionT& quat, MatrixT& A)
2045 {
2047 
2048  Scalar ww = quat[0] * quat[0];
2049  Scalar wx = quat[0] * quat[1];
2050  Scalar wy = quat[0] * quat[2];
2051  Scalar wz = quat[0] * quat[3];
2052 
2053  Scalar xx = quat[1] * quat[1];
2054  Scalar yy = quat[2] * quat[2];
2055  Scalar zz = quat[3] * quat[3];
2056 
2057  Scalar xy = quat[1] * quat[2];
2058  Scalar xz = quat[1] * quat[3];
2059  Scalar yz = quat[2] * quat[3];
2060 
2061  Scalar rr = xx + yy + zz;
2062  // normalization factor, just in case quaternion was not normalized
2063  Scalar f = 1 / (ww + rr);
2064  Scalar s = (ww - rr) * f;
2065  f *= 2;
2066 
2068 
2069  Wrapper::template Get<0, 0>(A) = xx * f + s;
2070  Wrapper::template Get<1, 0>(A) = (xy + wz) * f;
2071  Wrapper::template Get<2, 0>(A) = (xz - wy) * f;
2072 
2073  Wrapper::template Get<0, 1>(A) = (xy - wz) * f;
2074  Wrapper::template Get<1, 1>(A) = yy * f + s;
2075  Wrapper::template Get<2, 1>(A) = (yz + wx) * f;
2076 
2077  Wrapper::template Get<0, 2>(A) = (xz + wy) * f;
2078  Wrapper::template Get<1, 2>(A) = (yz - wx) * f;
2079  Wrapper::template Get<2, 2>(A) = zz * f + s;
2080 }
2081 } // anonymous namespace
2082 
2083 //------------------------------------------------------------------------------
2084 inline void vtkMath::QuaternionToMatrix3x3(const float quat[4], float A[3][3])
2085 {
2086  vtkQuaternionToMatrix3x3(quat, A);
2087 }
2088 
2089 //------------------------------------------------------------------------------
2090 inline void vtkMath::QuaternionToMatrix3x3(const double quat[4], double A[3][3])
2091 {
2092  vtkQuaternionToMatrix3x3(quat, A);
2093 }
2094 
2095 //-----------------------------------------------------------------------------
2096 template <class QuaternionT, class MatrixT, class EnableT>
2097 inline void vtkMath::QuaternionToMatrix3x3(const QuaternionT& q, MatrixT&& A)
2098 {
2099  vtkQuaternionToMatrix3x3(q, A);
2100 }
2101 
2102 namespace
2103 {
2104 //------------------------------------------------------------------------------
2105 // The solution is based on
2106 // Berthold K. P. Horn (1987),
2107 // "Closed-form solution of absolute orientation using unit quaternions,"
2108 // Journal of the Optical Society of America A, 4:629-642
2109 template <class MatrixT, class QuaternionT>
2110 inline void vtkMatrix3x3ToQuaternion(const MatrixT& A, QuaternionT& quat)
2111 {
2113 
2114  Scalar N[4][4];
2115 
2117 
2118  // on-diagonal elements
2119  N[0][0] = Wrapper::template Get<0, 0>(A) + Wrapper::template Get<1, 1>(A) +
2120  Wrapper::template Get<2, 2>(A);
2121  N[1][1] = Wrapper::template Get<0, 0>(A) - Wrapper::template Get<1, 1>(A) -
2122  Wrapper::template Get<2, 2>(A);
2123  N[2][2] = -Wrapper::template Get<0, 0>(A) + Wrapper::template Get<1, 1>(A) -
2124  Wrapper::template Get<2, 2>(A);
2125  N[3][3] = -Wrapper::template Get<0, 0>(A) - Wrapper::template Get<1, 1>(A) +
2126  Wrapper::template Get<2, 2>(A);
2127 
2128  // off-diagonal elements
2129  N[0][1] = N[1][0] = Wrapper::template Get<2, 1>(A) - Wrapper::template Get<1, 2>(A);
2130  N[0][2] = N[2][0] = Wrapper::template Get<0, 2>(A) - Wrapper::template Get<2, 0>(A);
2131  N[0][3] = N[3][0] = Wrapper::template Get<1, 0>(A) - Wrapper::template Get<0, 1>(A);
2132 
2133  N[1][2] = N[2][1] = Wrapper::template Get<1, 0>(A) + Wrapper::template Get<0, 1>(A);
2134  N[1][3] = N[3][1] = Wrapper::template Get<0, 2>(A) + Wrapper::template Get<2, 0>(A);
2135  N[2][3] = N[3][2] = Wrapper::template Get<2, 1>(A) + Wrapper::template Get<1, 2>(A);
2136 
2137  Scalar eigenvectors[4][4], eigenvalues[4];
2138 
2139  // convert into format that JacobiN can use,
2140  // then use Jacobi to find eigenvalues and eigenvectors
2141  Scalar *NTemp[4], *eigenvectorsTemp[4];
2142  for (int i = 0; i < 4; ++i)
2143  {
2144  NTemp[i] = N[i];
2145  eigenvectorsTemp[i] = eigenvectors[i];
2146  }
2147  vtkMath::JacobiN(NTemp, 4, eigenvalues, eigenvectorsTemp);
2148 
2149  // the first eigenvector is the one we want
2150  quat[0] = eigenvectors[0][0];
2151  quat[1] = eigenvectors[1][0];
2152  quat[2] = eigenvectors[2][0];
2153  quat[3] = eigenvectors[3][0];
2154 }
2155 } // anonymous namespace
2156 
2157 //------------------------------------------------------------------------------
2158 inline void vtkMath::Matrix3x3ToQuaternion(const float A[3][3], float quat[4])
2159 {
2160  vtkMatrix3x3ToQuaternion(A, quat);
2161 }
2162 
2163 //------------------------------------------------------------------------------
2164 inline void vtkMath::Matrix3x3ToQuaternion(const double A[3][3], double quat[4])
2165 {
2166  vtkMatrix3x3ToQuaternion(A, quat);
2167 }
2168 
2169 //-----------------------------------------------------------------------------
2170 template <class MatrixT, class QuaternionT, class EnableT>
2171 inline void vtkMath::Matrix3x3ToQuaternion(const MatrixT& A, QuaternionT&& q)
2172 {
2173  vtkMatrix3x3ToQuaternion(A, q);
2174 }
2175 
2176 namespace vtk_detail
2177 {
2178 // Can't specialize templates inside a template class, so we move the impl here.
2179 template <typename OutT>
2180 void RoundDoubleToIntegralIfNecessary(double val, OutT* ret)
2181 { // OutT is integral -- clamp and round
2182  if (!vtkMath::IsNan(val))
2183  {
2184  double min = static_cast<double>(vtkTypeTraits<OutT>::Min());
2185  double max = static_cast<double>(vtkTypeTraits<OutT>::Max());
2186  val = vtkMath::ClampValue(val, min, max);
2187  *ret = static_cast<OutT>((val >= 0.0) ? (val + 0.5) : (val - 0.5));
2188  }
2189  else
2190  *ret = 0;
2191 }
2192 template <>
2193 inline void RoundDoubleToIntegralIfNecessary(double val, double* retVal)
2194 { // OutT is double: passthrough
2195  *retVal = val;
2196 }
2197 template <>
2198 inline void RoundDoubleToIntegralIfNecessary(double val, float* retVal)
2199 { // OutT is float -- just clamp (as doubles, then the cast to float is well-defined.)
2200  if (!vtkMath::IsNan(val))
2201  {
2202  double min = static_cast<double>(vtkTypeTraits<float>::Min());
2203  double max = static_cast<double>(vtkTypeTraits<float>::Max());
2204  val = vtkMath::ClampValue(val, min, max);
2205  }
2206 
2207  *retVal = static_cast<float>(val);
2208 }
2209 } // end namespace vtk_detail
2210 
2211 //-----------------------------------------------------------------------------
2212 #if defined(VTK_HAS_ISINF) || defined(VTK_HAS_STD_ISINF)
2213 #define VTK_MATH_ISINF_IS_INLINE
2214 inline vtkTypeBool vtkMath::IsInf(double x)
2215 {
2216 #if defined(VTK_HAS_STD_ISINF)
2217  return std::isinf(x);
2218 #else
2219  return (isinf(x) != 0); // Force conversion to bool
2220 #endif
2221 }
2222 #endif
2223 
2224 //-----------------------------------------------------------------------------
2225 #if defined(VTK_HAS_ISNAN) || defined(VTK_HAS_STD_ISNAN)
2226 #define VTK_MATH_ISNAN_IS_INLINE
2227 inline vtkTypeBool vtkMath::IsNan(double x)
2228 {
2229 #if defined(VTK_HAS_STD_ISNAN)
2230  return std::isnan(x);
2231 #else
2232  return (isnan(x) != 0); // Force conversion to bool
2233 #endif
2234 }
2235 #endif
2236 
2237 //-----------------------------------------------------------------------------
2238 #if defined(VTK_HAS_ISFINITE) || defined(VTK_HAS_STD_ISFINITE) || defined(VTK_HAS_FINITE)
2239 #define VTK_MATH_ISFINITE_IS_INLINE
2240 inline bool vtkMath::IsFinite(double x)
2241 {
2242 #if defined(VTK_HAS_STD_ISFINITE)
2243  return std::isfinite(x);
2244 #elif defined(VTK_HAS_ISFINITE)
2245  return (isfinite(x) != 0); // Force conversion to bool
2246 #else
2247  return (finite(x) != 0); // Force conversion to bool
2248 #endif
2249 }
2250 #endif
2251 
2252 #endif
Gaussian sequence of pseudo random numbers implemented with the Box-Mueller transform.
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:165
a simple class to control print indentation
Definition: vtkIndent.h:119
performs common math operations
Definition: vtkMath.h:203
static ReturnTypeT Distance2BetweenPoints(const TupleRangeT1 &p1, const TupleRangeT2 &p2)
Compute distance squared between two points p1 and p2.
Definition: vtkMath.h:1903
static void Multiply3x3(const float A[3][3], const float B[3][3], float C[3][3])
Multiply one 3x3 matrix by another according to C = AB.
static double Dot(const double a[3], const double b[3])
Dot product of two 3-vectors (double version).
Definition: vtkMath.h:560
static int GetScalarTypeFittingRange(double range_min, double range_max, double scale=1.0, double shift=0.0)
Return the scalar type that is most likely to have enough precision to store a given range of data on...
static void RGBToXYZ(double r, double g, double b, double *x, double *y, double *z)
Convert color from the RGB system to CIE XYZ.
static void Multiply3x3(const double A[3][3], const double B[3][3], double C[3][3])
Multiply one 3x3 matrix by another according to C = AB.
static double Norm(const double *x, int n)
Compute the norm of n-vector.
static int Round(float f)
Rounds a float to the nearest integer.
Definition: vtkMath.h:234
static int * BeginCombination(int m, int n)
Start iterating over "m choose n" objects.
static void MultiplyMatrixWithVector(const MatrixT &M, const VectorT1 &X, VectorT2 &&Y)
Multiply matrix M with vector Y such that Y = M x X.
Definition: vtkMath.h:956
static double Norm2D(const double x[2])
Compute the norm of a 2-vector.
Definition: vtkMath.h:828
static void XYZToRGB(double x, double y, double z, double *r, double *g, double *b)
Convert color from the CIE XYZ system to RGB.
static void Subtract(const float a[3], const float b[3], float c[3])
Subtraction of two 3-vectors (float version).
Definition: vtkMath.h:469
static void LUSolve3x3(const double A[3][3], const int index[3], double x[3])
LU back substitution for a 3x3 matrix.
static vtkTypeBool SolveHomogeneousLeastSquares(int numberOfSamples, double **xt, int xOrder, double **mt)
Solves for the least squares best fit matrix for the homogeneous equation X'M' = 0'.
static void Outer2D(const float x[2], const float y[2], float A[2][2])
Outer product of two 2-vectors (float version).
Definition: vtkMath.h:793
static bool ProjectVector(const double a[3], const double b[3], double projection[3])
Compute the projection of vector a on vector b and return it in projection[3].
static vtkSmartPointer< vtkMathInternal > Internal
Definition: vtkMath.h:1727
static float Norm(const float *x, int n)
Compute the norm of n-vector.
static vtkTypeBool ExtentIsWithinOtherExtent(const int extent1[6], const int extent2[6])
Return true if first 3D extent is within second 3D extent Extent is x-min, x-max, y-min,...
static void Add(const double a[3], const double b[3], double c[3])
Addition of two 3-vectors (double version).
Definition: vtkMath.h:458
static void RGBToHSV(float r, float g, float b, float *h, float *s, float *v)
Convert color in RGB format (Red, Green, Blue) to HSV format (Hue, Saturation, Value).
static float Norm(const float v[3])
Compute the norm of 3-vector (float version).
Definition: vtkMath.h:640
static ReturnTypeT Dot(const TupleRangeT1 &a, const TupleRangeT2 &b)
Compute dot product between two points p1 and p2.
Definition: vtkMath.h:584
static vtkTypeBool Jacobi(double **a, double *w, double **v)
Jacobi iteration for the solution of eigenvectors/eigenvalues of a 3x3 real symmetric matrix.
static void XYZToLab(const double xyz[3], double lab[3])
Convert Color from the CIE XYZ system to CIE-L*ab.
Definition: vtkMath.h:1423
void PrintSelf(ostream &os, vtkIndent indent) override
Methods invoked by print to print information about the object including superclasses.
static vtkTypeInt64 Factorial(int N)
Compute N factorial, N! = N*(N-1) * (N-2)...*3*2*1.
static vtkTypeInt64 Binomial(int m, int n)
The number of combinations of n objects from a pool of m objects (m>n).
static double Random()
Generate pseudo-random numbers distributed according to the uniform distribution between 0....
static void LinearSolve(const MatrixT &M, const VectorT1 &x, VectorT2 &y)
This method solves linear systems M * x = y.
Definition: vtkMath.h:1033
static void Identity3x3(float A[3][3])
Set A to the identity matrix.
static double GaussianAmplitude(const double variance, const double distanceFromMean)
Compute the amplitude of a Gaussian function with mean=0 and specified variance.
static void SingularValueDecomposition3x3(const float A[3][3], float U[3][3], float w[3], float VT[3][3])
Perform singular value decomposition on a 3x3 matrix.
static double Nan()
Special IEEE-754 number used to represent Not-A-Number (Nan).
static void Perpendiculars(const float v1[3], float v2[3], float v3[3], double theta)
Given a unit vector v1, find two unit vectors v2 and v3 such that v1 cross v2 = v3 (i....
static double Gaussian(double mean, double std)
Generate pseudo-random numbers distributed according to the Gaussian distribution with mean mean and ...
static bool IsFinite(double x)
Test if a number has finite value i.e.
static void LUSolveLinearSystem(double **A, int *index, double *x, int size)
Solve linear equations Ax = b using LU decomposition A = LU where L is lower triangular matrix and U ...
static double EstimateMatrixCondition(const double *const *A, int size)
Estimate the condition number of a LU factored matrix.
static void LUFactor3x3(float A[3][3], int index[3])
LU Factorization of a 3x3 matrix.
static void FreeCombination(int *combination)
Free the "iterator" array created by vtkMath::BeginCombination.
static double Random(double min, double max)
Generate pseudo-random numbers distributed according to the uniform distribution between min and max.
static void TensorFromSymmetricTensor(const T1 symmTensor[6], T2 tensor[9])
Convert a 6-Component symmetric tensor into a 9-Component tensor, no allocation performed.
static void LabToXYZ(const double lab[3], double xyz[3])
Convert color from the CIE-L*ab system to CIE XYZ.
Definition: vtkMath.h:1412
static double Solve3PointCircle(const double p1[3], const double p2[3], const double p3[3], double center[3])
In Euclidean space, there is a unique circle passing through any given three non-collinear points P1,...
static vtkTypeBool PointIsWithinBounds(const double point[3], const double bounds[6], const double delta[3])
Return true if point is within the given 3D bounds Bounds is x-min, x-max, y-min, y-max,...
static float Dot(const float a[3], const float b[3])
Dot product of two 3-vectors (float version).
Definition: vtkMath.h:552
static void Diagonalize3x3(const float A[3][3], float w[3], float V[3][3])
Diagonalize a symmetric 3x3 matrix and return the eigenvalues in w and the eigenvectors in the column...
static void LabToXYZ(double L, double a, double b, double *x, double *y, double *z)
Convert color from the CIE-L*ab system to CIE XYZ.
static vtkMatrixUtilities::ScalarTypeExtractor< MatrixT >::value_type Determinant(const MatrixT &M)
Computes the determinant of input square SizeT x SizeT matrix M.
Definition: vtkMath.h:990
static vtkTypeBool GetAdjustedScalarRange(vtkDataArray *array, int comp, double range[2])
Get a vtkDataArray's scalar range for a given component.
static bool ProjectVector(const float a[3], const float b[3], float projection[3])
Compute the projection of vector a on vector b and return it in projection[3].
static void Cross(const float a[3], const float b[3], float c[3])
Cross product of two 3-vectors.
Definition: vtkMath.h:1911
static void MultiplyScalar2D(float a[2], float s)
Multiplies a 2-vector by a scalar (float version).
Definition: vtkMath.h:517
static void HSVToRGB(const float hsv[3], float rgb[3])
Convert color in HSV format (Hue, Saturation, Value) to RGB format (Red, Green, Blue).
Definition: vtkMath.h:1396
static void Assign(const double a[3], double b[3])
Assign values to a 3-vector (double version).
Definition: vtkMath.h:442
static double Determinant2x2(const double c1[2], const double c2[2])
Calculate the determinant of a 2x2 matrix: | a b | | c d |.
Definition: vtkMath.h:855
static T Max(const T &a, const T &b)
Returns the maximum of the two arguments provided.
Definition: vtkMath.h:1804
static void Outer2D(const double x[2], const double y[2], double A[2][2])
Outer product of two 2-vectors (double version).
Definition: vtkMath.h:807
static void RandomSeed(int s)
Initialize seed value.
static double NegInf()
Special IEEE-754 number used to represent negative infinity.
static void MultiplyScalar2D(double a[2], double s)
Multiplies a 2-vector by a scalar (double version).
Definition: vtkMath.h:541
static void LabToRGB(double L, double a, double b, double *red, double *green, double *blue)
Convert color from the CIE-L*ab system to RGB.
static double Gaussian()
Generate pseudo-random numbers distributed according to the standard normal distribution.
static int Ceil(double x)
Rounds a double to the nearest integer not less than itself.
Definition: vtkMath.h:1789
static void HSVToRGB(const double hsv[3], double rgb[3])
Convert color in HSV format (Hue, Saturation, Value) to RGB format (Red, Green, Blue).
Definition: vtkMath.h:1401
~vtkMath() override=default
static ScalarT Dot(const VectorT1 &x, const VectorT2 &y)
Computes the dot product between 2 vectors x and y.
Definition: vtkMath.h:967
static double Inf()
Special IEEE-754 number used to represent positive infinity.
static vtkMath * New()
static double GaussianAmplitude(const double mean, const double variance, const double position)
Compute the amplitude of a Gaussian function with specified mean and variance.
static vtkTypeBool Jacobi(float **a, float *w, float **v)
Jacobi iteration for the solution of eigenvectors/eigenvalues of a 3x3 real symmetric matrix.
static int PlaneIntersectsAABB(const double bounds[6], const double normal[3], const double point[3])
Implements Plane / Axis-Aligned Bounding-Box intersection as described in Graphics Gems IV,...
static void RGBToXYZ(const double rgb[3], double xyz[3])
Convert color from the RGB system to CIE XYZ.
Definition: vtkMath.h:1445
static void QuaternionToMatrix3x3(const float quat[4], float A[3][3])
Convert a quaternion to a 3x3 rotation matrix.
Definition: vtkMath.h:2084
static int NearestPowerOfTwo(int x)
Compute the nearest power of two that is not less than x.
Definition: vtkMath.h:1766
static void HSVToRGB(double h, double s, double v, double *r, double *g, double *b)
Convert color in HSV format (Hue, Saturation, Value) to RGB format (Red, Green, Blue).
static void SingularValueDecomposition3x3(const double A[3][3], double U[3][3], double w[3], double VT[3][3])
Perform singular value decomposition on a 3x3 matrix.
static double SignedAngleBetweenVectors(const double v1[3], const double v2[3], const double vn[3])
Compute signed angle in radians between two vectors with regard to a third orthogonal vector.
static float Normalize2D(float v[2])
Normalize (in place) a 2-vector.
Definition: vtkMath.h:1838
static void Invert3x3(const double A[3][3], double AI[3][3])
Invert a 3x3 matrix.
static void HSVToRGB(float h, float s, float v, float *r, float *g, float *b)
Convert color in HSV format (Hue, Saturation, Value) to RGB format (Red, Green, Blue).
static void MultiplyQuaternion(const double q1[4], const double q2[4], double q[4])
Multiply two quaternions.
static void Multiply3x3(const double A[3][3], const double v[3], double u[3])
Multiply a vector by a 3x3 matrix.
static void Outer(const double a[3], const double b[3], double c[3][3])
Outer product of two 3-vectors (double version).
Definition: vtkMath.h:606
static vtkTypeBool InvertMatrix(double **A, double **AI, int size, int *tmp1Size, double *tmp2Size)
Thread safe version of InvertMatrix method.
static vtkTypeBool InvertMatrix(double **A, double **AI, int size)
Invert input square matrix A into matrix AI.
static void LUSolve3x3(const float A[3][3], const int index[3], float x[3])
LU back substitution for a 3x3 matrix.
static int GetSeed()
Return the current seed used by the random number generator.
static void Assign(const VectorT1 &a, VectorT2 &&b)
Assign values to a 3-vector (templated version).
Definition: vtkMath.h:432
static float RadiansFromDegrees(float degrees)
Convert degrees into radians.
Definition: vtkMath.h:1735
static void Convolve1D(Iter1 beginSample, Iter1 endSample, Iter2 beginKernel, Iter2 endKernel, Iter3 beginOut, Iter3 endOut, ConvolutionMode mode=ConvolutionMode::FULL)
Compute the convolution of a sampled 1D signal by a given kernel.
Definition: vtkMath.h:1682
static void RotateVectorByWXYZ(const double v[3], const double q[4], double r[3])
rotate a vector by WXYZ using // https://en.wikipedia.org/wiki/Rodrigues%27_rotation_formula
static void Add(const float a[3], const float b[3], float c[3])
Addition of two 3-vectors (float version).
Definition: vtkMath.h:447
static int CeilLog2(vtkTypeUInt64 x)
Gives the exponent of the lowest power of two not less than x.
static vtkTypeBool AreBoundsInitialized(const double bounds[6])
Are the bounds initialized?
Definition: vtkMath.h:1496
static bool ProjectVector2D(const double a[2], const double b[2], double projection[2])
Compute the projection of 2D vector a on 2D vector b and returns the result in projection[2].
static vtkTypeBool JacobiN(float **a, int n, float *w, float **v)
JacobiN iteration for the solution of eigenvectors/eigenvalues of a nxn real symmetric matrix.
static int NextCombination(int m, int n, int *combination)
Given m, n, and a valid combination of n integers in the range [0,m[, this function alters the intege...
static constexpr double Pi()
A mathematical constant.
Definition: vtkMath.h:212
static double GaussianWeight(const double variance, const double distanceFromMean)
Compute the amplitude of an unnormalized Gaussian function with mean=0 and specified variance.
static void Multiply3x3(const float A[3][3], const float v[3], float u[3])
Multiply a vector by a 3x3 matrix.
static void Subtract(const double a[3], const double b[3], double c[3])
Subtraction of two 3-vectors (double version).
Definition: vtkMath.h:480
static void Matrix3x3ToQuaternion(const float A[3][3], float quat[4])
Convert a 3x3 matrix into a quaternion.
Definition: vtkMath.h:2158
static void Orthogonalize3x3(const double A[3][3], double B[3][3])
Orthogonalize a 3x3 matrix and put the result in B.
static void XYZToRGB(const double xyz[3], double rgb[3])
Convert color from the CIE XYZ system to RGB.
Definition: vtkMath.h:1434
static double ClampAndNormalizeValue(double value, const double range[2])
Clamp a value against a range and then normalize it between 0 and 1.
Definition: vtkMath.h:1992
static void MultiplyScalar(double a[3], double s)
Multiplies a 3-vector by a scalar (double version).
Definition: vtkMath.h:529
static double Dot2D(const double x[2], const double y[2])
Dot product of two 2-vectors.
Definition: vtkMath.h:788
static void LinearSolve3x3(const float A[3][3], const float x[3], float y[3])
Solve Ay = x for y and place the result in y.
static void MultiplyMatrix(const MatrixT1 &M1, const MatrixT2 &M2, MatrixT3 &&M3)
Multiply matrices such that M3 = M1 x M2.
Definition: vtkMath.h:929
static vtkTypeBool IsNan(double x)
Test if a number is equal to the special floating point value Not-A-Number (Nan).
static void Diagonalize3x3(const double A[3][3], double w[3], double V[3][3])
Diagonalize a symmetric 3x3 matrix and return the eigenvalues in w and the eigenvectors in the column...
static void RGBToLab(const double rgb[3], double lab[3])
Convert color from the RGB system to CIE-L*ab.
Definition: vtkMath.h:1459
static int Floor(double x)
Rounds a double to the nearest integer not greater than itself.
Definition: vtkMath.h:1780
static void RotateVectorByNormalizedQuaternion(const double v[3], const double q[4], double r[3])
rotate a vector by a normalized quaternion using // https://en.wikipedia.org/wiki/Rodrigues%27_rotati...
static void Subtract(const VectorT1 &a, const VectorT2 &b, VectorT3 &&c)
Subtraction of two 3-vectors (templated version).
Definition: vtkMath.h:494
static vtkTypeBool BoundsIsWithinOtherBounds(const double bounds1[6], const double bounds2[6], const double delta[3])
Return true if first 3D bounds is within the second 3D bounds Bounds is x-min, x-max,...
static double Determinant2x2(double a, double b, double c, double d)
Calculate the determinant of a 2x2 matrix: | a b | | c d |.
Definition: vtkMath.h:854
static void RGBToHSV(const double rgb[3], double hsv[3])
Convert color in RGB format (Red, Green, Blue) to HSV format (Hue, Saturation, Value).
Definition: vtkMath.h:1381
static vtkTypeBool JacobiN(double **a, int n, double *w, double **v)
JacobiN iteration for the solution of eigenvectors/eigenvalues of a nxn real symmetric matrix.
static double AngleBetweenVectors(const double v1[3], const double v2[3])
Compute angle in radians between two vectors.
static void MultiplyMatrix(const double *const *A, const double *const *B, unsigned int rowA, unsigned int colA, unsigned int rowB, unsigned int colB, double **C)
General matrix multiplication.
static float DegreesFromRadians(float radians)
Convert radians into degrees.
Definition: vtkMath.h:1747
static float Determinant2x2(const float c1[2], const float c2[2])
Compute determinant of 2x2 matrix.
Definition: vtkMath.h:845
static int Round(double f)
Definition: vtkMath.h:235
static vtkTypeBool IsInf(double x)
Test if a number is equal to the special floating point value infinity.
static void UninitializeBounds(double bounds[6])
Set the bounds to an uninitialized state.
Definition: vtkMath.h:1481
vtkMath()=default
static void RGBToHSV(double r, double g, double b, double *h, double *s, double *v)
Convert color in RGB format (Red, Green, Blue) to HSV format (Hue, Saturation, Value).
static void Outer(const float a[3], const float b[3], float c[3][3])
Outer product of two 3-vectors (float version).
Definition: vtkMath.h:592
static double Norm(const double v[3])
Compute the norm of 3-vector (double version).
Definition: vtkMath.h:645
static void RoundDoubleToIntegralIfNecessary(double val, OutT *ret)
Round a double to type OutT if OutT is integral, otherwise simply clamp the value to the output range...
Definition: vtkMath.h:243
static void RotateVectorByWXYZ(const float v[3], const float q[4], float r[3])
rotate a vector by WXYZ using // https://en.wikipedia.org/wiki/Rodrigues%27_rotation_formula
static bool IsPowerOfTwo(vtkTypeUInt64 x)
Returns true if integer is a power of two.
Definition: vtkMath.h:1759
static void Invert3x3(const float A[3][3], float AI[3][3])
Invert a 3x3 matrix.
static float Normalize(float v[3])
Normalize (in place) a 3-vector.
Definition: vtkMath.h:1810
static void Transpose3x3(const double A[3][3], double AT[3][3])
Transpose a 3x3 matrix.
static ReturnTypeT SquaredNorm(const TupleRangeT &v)
Compute the squared norm of a 3-vector.
Definition: vtkMath.h:660
static double Determinant3x3(const float A[3][3])
Return the determinant of a 3x3 matrix.
Definition: vtkMath.h:1942
static float Dot2D(const float x[2], const float y[2])
Dot product of two 2-vectors.
Definition: vtkMath.h:783
ConvolutionMode
Definition: vtkMath.h:1653
static void RotateVectorByNormalizedQuaternion(const float v[3], const float q[4], float r[3])
rotate a vector by a normalized quaternion using // https://en.wikipedia.org/wiki/Rodrigues%27_rotati...
static ScalarT Dot(const VectorT1 &x, const MatrixT &M, const VectorT2 &y)
Computes the dot product x^T M y, where x and y are vectors and M is a metric matrix.
Definition: vtkMath.h:1054
static void RGBToHSV(const float rgb[3], float hsv[3])
Convert color in RGB format (Red, Green, Blue) to HSV format (Hue, Saturation, Value).
Definition: vtkMath.h:1376
static void Orthogonalize3x3(const float A[3][3], float B[3][3])
Orthogonalize a 3x3 matrix and put the result in B.
static bool ProjectVector2D(const float a[2], const float b[2], float projection[2])
Compute the projection of 2D vector a on 2D vector b and returns the result in projection[2].
static vtkTypeBool SolveLinearSystemGEPP2x2(double a00, double a01, double a10, double a11, double b0, double b1, double &x0, double &x1)
Solve linear equation Ax = b using Gaussian Elimination with Partial Pivoting for a 2x2 system.
static vtkTypeBool SolveLinearSystem(double **A, double *x, int size)
Solve linear equations Ax = b using Crout's method.
static void LabToRGB(const double lab[3], double rgb[3])
Convert color from the CIE-L*ab system to RGB.
Definition: vtkMath.h:1470
static float Norm2D(const float x[2])
Compute the norm of a 2-vector.
Definition: vtkMath.h:822
static double GaussianWeight(const double mean, const double variance, const double position)
Compute the amplitude of an unnormalized Gaussian function with specified mean and variance.
static vtkTypeBool LUFactorLinearSystem(double **A, int *index, int size, double *tmpSize)
Thread safe version of LUFactorLinearSystem method.
static void LinearSolve3x3(const double A[3][3], const double x[3], double y[3])
Solve Ay = x for y and place the result in y.
static void XYZToLab(double x, double y, double z, double *L, double *a, double *b)
Convert Color from the CIE XYZ system to CIE-L*ab.
static void InvertMatrix(const MatrixT1 &M1, MatrixT2 &&M2)
Computes the inverse of input matrix M1 into M2.
Definition: vtkMath.h:1013
static void MultiplyScalar(float a[3], float s)
Multiplies a 3-vector by a scalar (float version).
Definition: vtkMath.h:505
static T Min(const T &a, const T &b)
Returns the minimum of the two arguments provided.
Definition: vtkMath.h:1797
static void Perpendiculars(const double v1[3], double v2[3], double v3[3], double theta)
Given a unit vector v1, find two unit vectors v2 and v3 such that v1 cross v2 = v3 (i....
static T ClampValue(const T &value, const T &min, const T &max)
Clamp some value against a range, return the result.
Definition: vtkMath.h:1955
static vtkTypeBool SolveLeastSquares(int numberOfSamples, double **xt, int xOrder, double **yt, int yOrder, double **mt, int checkHomogeneous=1)
Solves for the least squares best fit matrix for the equation X'M' = Y'.
static void Identity3x3(double A[3][3])
Set A to the identity matrix.
static void LUFactor3x3(double A[3][3], int index[3])
LU Factorization of a 3x3 matrix.
static vtkTypeBool LUFactorLinearSystem(double **A, int *index, int size)
Factor linear equations Ax = b using LU decomposition into the form A = LU where L is a unit lower tr...
static void RGBToLab(double red, double green, double blue, double *L, double *a, double *b)
Convert color from the RGB system to CIE-L*ab.
static void MultiplyQuaternion(const float q1[4], const float q2[4], float q[4])
Multiply two quaternions.
static void ClampValues(const double *values, int nb_values, const double range[2], double *clamped_values)
Clamp some values against a range The method without 'clamped_values' will perform in-place clamping.
static void Transpose3x3(const float A[3][3], float AT[3][3])
Transpose a 3x3 matrix.
static void ClampValues(double *values, int nb_values, const double range[2])
Clamp some values against a range The method without 'clamped_values' will perform in-place clamping.
static int QuadraticRoot(double a, double b, double c, double min, double max, double *u)
find roots of ax^2+bx+c=0 in the interval min,max.
Matrix wrapping class.
Park and Miller Sequence of pseudo random numbers.
abstract base class for most VTK objects
Definition: vtkObject.h:82
represent and manipulate 3D points
Definition: vtkPoints.h:149
Computes the portion of a dataset which is inside a selection.
@ point
Definition: vtkX3D.h:242
@ mode
Definition: vtkX3D.h:253
@ value
Definition: vtkX3D.h:226
@ scale
Definition: vtkX3D.h:235
@ range
Definition: vtkX3D.h:244
@ center
Definition: vtkX3D.h:236
@ type
Definition: vtkX3D.h:522
@ position
Definition: vtkX3D.h:267
@ size
Definition: vtkX3D.h:259
@ index
Definition: vtkX3D.h:252
void RoundDoubleToIntegralIfNecessary(double val, OutT *ret)
Definition: vtkMath.h:2180
detail::ScalarTypeExtractor< std::is_array< DerefContainer >::value||std::is_pointer< DerefContainer >::value, ContainerT >::value_type value_type
Template defining traits of native types used by VTK.
Definition: vtkTypeTraits.h:33
int vtkTypeBool
Definition: vtkABI.h:69
double vtkDeterminant3x3(const T A[3][3])
Definition: vtkMath.h:1935
#define max(a, b)