VTK
dox/Common/vtkMath.h
Go to the documentation of this file.
00001 /*=========================================================================
00002 
00003   Program:   Visualization Toolkit
00004   Module:    vtkMath.h
00005 
00006   Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
00007   All rights reserved.
00008   See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
00009 
00010      This software is distributed WITHOUT ANY WARRANTY; without even
00011      the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
00012      PURPOSE.  See the above copyright notice for more information.
00013 
00014 =========================================================================
00015   Copyright 2008 Sandia Corporation.
00016   Under the terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00017   license for use of this work by or on behalf of the
00018   U.S. Government. Redistribution and use in source and binary forms, with
00019   or without modification, are permitted provided that this Notice and any
00020   statement of authorship are reproduced on all copies.
00021 
00022   Contact: pppebay@sandia.gov,dcthomp@sandia.gov
00023 
00024 =========================================================================*/
00043 #ifndef __vtkMath_h
00044 #define __vtkMath_h
00045 
00046 #include "vtkObject.h"
00047 #include "vtkPolynomialSolversUnivariate.h" // For backwards compatibility of old solvers
00048 
00049 #include <assert.h> // assert() in inline implementations.
00050 
00051 #ifndef DBL_EPSILON
00052 #  define VTK_DBL_EPSILON    2.2204460492503131e-16
00053 #else  // DBL_EPSILON
00054 #  define VTK_DBL_EPSILON    DBL_EPSILON
00055 #endif  // DBL_EPSILON
00056 
00057 class vtkDataArray;
00058 class vtkPoints;
00059 class vtkMathInternal;
00060 class vtkMinimalStandardRandomSequence;
00061 class vtkBoxMuellerRandomSequence;
00062 
00063 class VTK_COMMON_EXPORT vtkMath : public vtkObject
00064 {
00065 public:
00066   static vtkMath *New();
00067   vtkTypeMacro(vtkMath,vtkObject);
00068   void PrintSelf(ostream& os, vtkIndent indent);
00069 
00071   static float Pi() { return 3.14159265358979f; };
00072 
00075   static double DoubleTwoPi() { return  6.283185307179586; };
00076 
00079   static double DoublePi() { return 3.1415926535897932384626; };
00080 
00082 
00083   static float RadiansFromDegrees( float degrees);
00084   static double RadiansFromDegrees( double degrees);
00086 
00088 
00089   static float DegreesFromRadians( float radians);
00090   static double DegreesFromRadians( double radians);
00092 
00094 
00095   VTK_LEGACY(static float DegreesToRadians());
00096   VTK_LEGACY(static double DoubleDegreesToRadians());
00098 
00100 
00101   VTK_LEGACY(static float RadiansToDegrees());
00102   VTK_LEGACY(static double DoubleRadiansToDegrees());
00104 
00106 
00107   static int Round(float f) {
00108     return static_cast<int>( f + ( f >= 0 ? 0.5 : -0.5 ) ); }
00109   static int Round(double f) {
00110     return static_cast<int>( f + ( f >= 0 ? 0.5 : -0.5 ) ); }
00112 
00114   static int Floor(double x);
00115   
00118   static vtkTypeInt64 Factorial( int N );
00119 
00123   static vtkTypeInt64 Binomial( int m, int n );
00124 
00131   static int* BeginCombination( int m, int n );
00132 
00139   static int NextCombination( int m, int n, int* combination );
00140 
00142   static void FreeCombination( int* combination);
00143 
00155   static void RandomSeed(int s);  
00156 
00165   static int GetSeed();
00166   
00176   static double Random();  
00177 
00186   static double Random( double min, double max );
00187 
00196   static double Gaussian();  
00197 
00207   static double Gaussian( double mean, double std );
00208 
00210 
00211   static void Add(const float a[3], const float b[3], float c[3]) {
00212     for (int i = 0; i < 3; ++i)
00213       c[i] = a[i] + b[i];
00214   }
00216 
00218 
00219   static void Add(const double a[3], const double b[3], double c[3]) {
00220     for (int i = 0; i < 3; ++i)
00221       c[i] = a[i] + b[i];
00222   }
00224 
00226 
00228   static void Subtract(const float a[3], const float b[3], float c[3]) {
00229     for (int i = 0; i < 3; ++i)
00230       c[i] = a[i] - b[i];
00231   }
00233 
00235 
00237   static void Subtract(const double a[3], const double b[3], double c[3]) {
00238     for (int i = 0; i < 3; ++i)
00239       c[i] = a[i] - b[i];
00240   }
00242 
00244 
00246   static void MultiplyScalar(float a[3], float s) {
00247     for (int i = 0; i < 3; ++i)
00248       a[i] *= s;
00249   }
00251 
00253 
00255   static void MultiplyScalar(double a[3], double s) {
00256     for (int i = 0; i < 3; ++i)
00257       a[i] *= s;
00258   }
00260 
00262 
00263   static float Dot(const float x[3], const float y[3]) {
00264     return ( x[0] * y[0] + x[1] * y[1] + x[2] * y[2] );};
00266 
00268 
00269   static double Dot(const double x[3], const double y[3]) {
00270     return ( x[0] * y[0] + x[1] * y[1] + x[2] * y[2] );};
00272   
00274 
00275   static void Outer(const float x[3], const float y[3], float A[3][3]) {
00276     for (int i=0; i < 3; i++)
00277       for (int j=0; j < 3; j++)
00278         A[i][j] = x[i] * y[j];
00279   }
00280   // Description:
00281   // Outer product of two 3-vectors (double-precision version).
00282   static void Outer(const double x[3], const double y[3], double A[3][3]) {
00283     for (int i=0; i < 3; i++)
00284       for (int j=0; j < 3; j++)
00285         A[i][j] = x[i] * y[j];
00286   }
00288 
00290   static void Cross(const float x[3], const float y[3], float z[3]);
00291 
00294   static void Cross(const double x[3], const double y[3], double z[3]);
00295 
00297 
00298   static float Norm(const float* x, int n); 
00299   static double Norm(const double* x, int n); 
00301 
00303 
00304   static float Norm(const float x[3]) {
00305     return static_cast<float> (sqrt( x[0] * x[0] + x[1] * x[1] + x[2] * x[2] ) );};
00307   
00309 
00310   static double Norm(const double x[3]) {
00311     return sqrt( x[0] * x[0] + x[1] * x[1] + x[2] * x[2] );};
00313   
00315   static float Normalize(float x[3]);
00316 
00319   static double Normalize(double x[3]);
00320 
00322 
00327   static void Perpendiculars(const double x[3], double y[3], double z[3], 
00328                              double theta);
00329   static void Perpendiculars(const float x[3], float y[3], float z[3],
00330                              double theta);
00332 
00334   static float Distance2BetweenPoints(const float x[3], const float y[3]);
00335 
00338   static double Distance2BetweenPoints(const double x[3], const double y[3]);
00339 
00341 
00342   static float Dot2D(const float x[2], const float y[2]) {
00343     return ( x[0] * y[0] + x[1] * y[1] );};
00345   
00347 
00348   static double Dot2D(const double x[2], const double y[2]) {
00349     return ( x[0] * y[0] + x[1] * y[1] );};
00351 
00353 
00354   static void Outer2D(const float x[2], const float y[2], float A[2][2]) 
00355     {
00356     for (int i=0; i < 2; i++)
00357       {
00358       for (int j=0; j < 2; j++)
00359         {
00360         A[i][j] = x[i] * y[j];
00361         }
00362       }
00363     }
00364   // Description:
00365   // Outer product of two 2-vectors (float version).
00366   static void Outer2D(const double x[2], const double y[2], double A[2][2]) 
00367     {
00368     for (int i=0; i < 2; i++)
00369       {
00370       for (int j=0; j < 2; j++)
00371         {
00372         A[i][j] = x[i] * y[j];
00373         }
00374       }
00375     }
00377 
00379 
00380   static float Norm2D(const float x[2]) {
00381     return static_cast<float> (sqrt( x[0] * x[0] + x[1] * x[1] ) );};
00383 
00385 
00386   static double Norm2D(const double x[2]) {
00387     return sqrt( x[0] * x[0] + x[1] * x[1] );};
00389 
00391   static float Normalize2D(float x[2]);
00392 
00395   static double Normalize2D(double x[2]);
00396 
00398 
00399   static float Determinant2x2(const float c1[2], const float c2[2]) {
00400     return (c1[0] * c2[1] - c2[0] * c1[1] );};
00402 
00404 
00405   static double Determinant2x2(double a, double b, double c, double d) {
00406     return (a * d - b * c);};
00407   static double Determinant2x2(const double c1[2], const double c2[2]) {
00408     return (c1[0] * c2[1] - c2[0] * c1[1] );};
00410 
00412 
00413   static void LUFactor3x3(float A[3][3], int index[3]);
00414   static void LUFactor3x3(double A[3][3], int index[3]);
00416 
00418 
00419   static void LUSolve3x3(const float A[3][3], const int index[3], 
00420                          float x[3]);
00421   static void LUSolve3x3(const double A[3][3], const int index[3], 
00422                          double x[3]);
00424 
00426 
00428   static void LinearSolve3x3(const float A[3][3], const float x[3], 
00429                              float y[3]);
00430   static void LinearSolve3x3(const double A[3][3], const double x[3], 
00431                              double y[3]);
00433 
00435 
00436   static void Multiply3x3(const float A[3][3], const float in[3], 
00437                           float out[3]);
00438   static void Multiply3x3(const double A[3][3], const double in[3], 
00439                           double out[3]);
00441   
00443 
00444   static void Multiply3x3(const float A[3][3], const float B[3][3], 
00445                           float C[3][3]);
00446   static void Multiply3x3(const double A[3][3], const double B[3][3], 
00447                           double C[3][3]);
00449 
00451 
00453   static void MultiplyMatrix(const double **A, const double **B,
00454                              unsigned int rowA, unsigned int colA, 
00455                              unsigned int rowB, unsigned int colB,
00456                              double **C);
00458 
00460 
00462   static void Transpose3x3(const float A[3][3], float AT[3][3]);
00463   static void Transpose3x3(const double A[3][3], double AT[3][3]);
00465 
00467 
00469   static void Invert3x3(const float A[3][3], float AI[3][3]);
00470   static void Invert3x3(const double A[3][3], double AI[3][3]);
00472 
00474 
00475   static void Identity3x3(float A[3][3]);
00476   static void Identity3x3(double A[3][3]);
00478 
00480 
00481   static double Determinant3x3(float A[3][3]);
00482   static double Determinant3x3(double A[3][3]);
00484 
00486 
00487   static float Determinant3x3(const float c1[3], 
00488                               const float c2[3], 
00489                               const float c3[3]);
00491 
00493 
00494   static double Determinant3x3(const double c1[3], 
00495                                const double c2[3], 
00496                                const double c3[3]);
00498 
00500 
00502   static double Determinant3x3(double a1, double a2, double a3, 
00503                                double b1, double b2, double b3, 
00504                                double c1, double c2, double c3);
00506 
00508 
00510   static void QuaternionToMatrix3x3(const float quat[4], float A[3][3]); 
00511   static void QuaternionToMatrix3x3(const double quat[4], double A[3][3]); 
00513 
00515 
00518   static void Matrix3x3ToQuaternion(const float A[3][3], float quat[4]);
00519   static void Matrix3x3ToQuaternion(const double A[3][3], double quat[4]);
00521   
00523 
00526   static void Orthogonalize3x3(const float A[3][3], float B[3][3]);
00527   static void Orthogonalize3x3(const double A[3][3], double B[3][3]);
00529 
00531 
00535   static void Diagonalize3x3(const float A[3][3], float w[3], float V[3][3]);
00536   static void Diagonalize3x3(const double A[3][3],double w[3],double V[3][3]);
00538 
00540 
00547   static void SingularValueDecomposition3x3(const float A[3][3],
00548                                             float U[3][3], float w[3],
00549                                             float VT[3][3]);
00550   static void SingularValueDecomposition3x3(const double A[3][3],
00551                                             double U[3][3], double w[3],
00552                                             double VT[3][3]);
00554 
00559   static int SolveLinearSystem(double **A, double *x, int size);
00560 
00564   static int InvertMatrix(double **A, double **AI, int size);
00565 
00567 
00569   static int InvertMatrix(double **A, double **AI, int size,
00570                           int *tmp1Size, double *tmp2Size);
00572 
00578   static int LUFactorLinearSystem(double **A, int *index, int size);
00579 
00581 
00583   static int LUFactorLinearSystem(double **A, int *index, int size,
00584                                   double *tmpSize);
00586 
00588 
00594   static void LUSolveLinearSystem(double **A, int *index, 
00595                                   double *x, int size);
00597 
00605   static double EstimateMatrixCondition(double **A, int size);
00606 
00608 
00612   static int Jacobi(float **a, float *w, float **v);
00613   static int Jacobi(double **a, double *w, double **v);
00615 
00617 
00622   static int JacobiN(float **a, int n, float *w, float **v);
00623   static int JacobiN(double **a, int n, double *w, double **v);
00625 
00632   static double* SolveCubic(double c0, double c1, double c2, double c3);
00633 
00640   static double* SolveQuadratic(double c0, double c1, double c2);
00641 
00645   static double* SolveLinear(double c0, double c1);
00646 
00648 
00659   static int SolveCubic(double c0, double c1, double c2, double c3, 
00660                         double *r1, double *r2, double *r3, int *num_roots);
00662 
00664 
00668   static int SolveQuadratic(double c0, double c1, double c2, 
00669                             double *r1, double *r2, int *num_roots);
00671   
00677   static int SolveQuadratic( double* c, double* r, int* m );
00678 
00683   static int SolveLinear(double c0, double c1, double *r1, int *num_roots);
00684 
00686 
00696   static int SolveHomogeneousLeastSquares(int numberOfSamples, double **xt, int xOrder,
00697                                 double **mt);
00699 
00700 
00702 
00713   static int SolveLeastSquares(int numberOfSamples, double **xt, int xOrder,
00714                                double **yt, int yOrder, double **mt, int checkHomogeneous=1);
00716 
00718 
00720   static void RGBToHSV(const float rgb[3], float hsv[3])
00721     { RGBToHSV(rgb[0], rgb[1], rgb[2], hsv, hsv+1, hsv+2); }
00722   static void RGBToHSV(float r, float g, float b, float *h, float *s, float *v);
00723   static double* RGBToHSV(const double rgb[3]);
00724   static double* RGBToHSV(double r, double g, double b);
00725   static void RGBToHSV(const double rgb[3], double hsv[3])
00726     { RGBToHSV(rgb[0], rgb[1], rgb[2], hsv, hsv+1, hsv+2); }
00727   static void RGBToHSV(double r, double g, double b, double *h, double *s, double *v);
00729 
00731 
00733   static void HSVToRGB(const float hsv[3], float rgb[3])
00734     { HSVToRGB(hsv[0], hsv[1], hsv[2], rgb, rgb+1, rgb+2); }
00735   static void HSVToRGB(float h, float s, float v, float *r, float *g, float *b);
00736   static double* HSVToRGB(const double hsv[3]);
00737   static double* HSVToRGB(double h, double s, double v);
00738   static void HSVToRGB(const double hsv[3], double rgb[3])
00739     { HSVToRGB(hsv[0], hsv[1], hsv[2], rgb, rgb+1, rgb+2); }
00740   static void HSVToRGB(double h, double s, double v, double *r, double *g, double *b);
00742 
00744 
00745   static void LabToXYZ(const double lab[3], double xyz[3]) {
00746     LabToXYZ(lab[0], lab[1], lab[2], xyz+0, xyz+1, xyz+2);
00747   }
00748   static void LabToXYZ(double L, double a, double b,
00749                        double *x, double *y, double *z);
00750   static double *LabToXYZ(const double lab[3]);
00752 
00754 
00755   static void XYZToLab(const double xyz[3], double lab[3]) {
00756     XYZToLab(xyz[0], xyz[1], xyz[2], lab+0, lab+1, lab+2);
00757   }
00758   static void XYZToLab(double x, double y, double z,
00759                        double *L, double *a, double *b);
00760   static double *XYZToLab(const double xyz[3]);
00762 
00764 
00765   static void XYZToRGB(const double xyz[3], double rgb[3]) {
00766     XYZToRGB(xyz[0], xyz[1], xyz[2], rgb+0, rgb+1, rgb+2);
00767   }
00768   static void XYZToRGB(double x, double y, double z,
00769                        double *r, double *g, double *b);
00770   static double *XYZToRGB(const double xyz[3]);
00772 
00774 
00775   static void RGBToXYZ(const double rgb[3], double xyz[3]) {
00776     RGBToXYZ(rgb[0], rgb[1], rgb[2], xyz+0, xyz+1, xyz+2);
00777   }
00778   static void RGBToXYZ(double r, double g, double b,
00779                        double *x, double *y, double *z);
00780   static double *RGBToXYZ(const double rgb[3]);
00782 
00784 
00785   static void RGBToLab(const double rgb[3], double lab[3]) {
00786     RGBToLab(rgb[0], rgb[1], rgb[2], lab+0, lab+1, lab+2);
00787   }
00788   static void RGBToLab(double red, double green, double blue,
00789                        double *L, double *a, double *b);
00790   static double *RGBToLab(const double rgb[3]);
00792 
00794 
00795   static void LabToRGB(const double lab[3], double rgb[3]) {
00796     LabToRGB(lab[0], lab[1], lab[2], rgb+0, rgb+1, rgb+2);
00797   }
00798   static void LabToRGB(double L, double a, double b,
00799                        double *red, double *green, double *blue);
00800   static double *LabToRGB(const double lab[3]);
00802 
00804 
00805   static void UninitializeBounds(double bounds[6]){
00806     bounds[0] = 1.0;
00807     bounds[1] = -1.0;
00808     bounds[2] = 1.0;
00809     bounds[3] = -1.0;
00810     bounds[4] = 1.0;
00811     bounds[5] = -1.0;
00812   }
00814   
00816 
00817   static int AreBoundsInitialized(double bounds[6]){
00818     if ( bounds[1]-bounds[0]<0.0 )
00819       {
00820       return 0;
00821       }
00822     return 1;
00823   }
00825 
00827 
00829   static void ClampValue(double *value, const double range[2]);
00830   static void ClampValue(double value, const double range[2], double *clamped_value);
00831   static void ClampValues(
00832     double *values, int nb_values, const double range[2]);
00833   static void ClampValues(
00834     const double *values, int nb_values, const double range[2], double *clamped_values);
00836 
00838 
00841   static double ClampAndNormalizeValue(double value,
00842                                        const double range[2]);
00844 
00846 
00852   static int GetScalarTypeFittingRange(
00853     double range_min, double range_max, 
00854     double scale = 1.0, double shift = 0.0);
00856 
00858 
00864   static int GetAdjustedScalarRange(
00865     vtkDataArray *array, int comp, double range[2]);
00867  
00870   static int ExtentIsWithinOtherExtent(int extent1[6], int extent2[6]);
00871   
00875   static int BoundsIsWithinOtherBounds(double bounds1[6], double bounds2[6], double delta[3]);
00876   
00880   static int PointIsWithinBounds(double point[3], double bounds[6], double delta[3]);
00881 
00885   static void SpiralPoints(vtkIdType num, vtkPoints * offsets);
00886   
00894   static double Solve3PointCircle(const double p1[3], const double p2[3], const double p3[3], double center[3]);
00895 
00897   static double Inf();
00898 
00900   static double NegInf();
00901 
00903   static double Nan();
00904 
00907   static int IsInf(double x);
00908 
00909   // Test if a number is equal to the special floating point value Not-A-Number (Nan).
00910   static int IsNan(double x);
00911 
00912 protected:
00913   vtkMath() {};
00914   ~vtkMath() {};
00915   
00916   static vtkMathInternal Internal;
00917 private:
00918   vtkMath(const vtkMath&);  // Not implemented.
00919   void operator=(const vtkMath&);  // Not implemented.
00920 };
00921 
00922 //BTX
00923 class vtkMathInternal
00924 {
00925 public:
00926   vtkMathInternal();
00927   ~vtkMathInternal();
00928   vtkMinimalStandardRandomSequence *Uniform;
00929   vtkBoxMuellerRandomSequence *Gaussian;
00930 };
00931 //ETX
00932 
00933 //----------------------------------------------------------------------------
00934 inline float vtkMath::RadiansFromDegrees( float x )
00935 {
00936   return x * 0.017453292f;
00937 }
00938 
00939 //----------------------------------------------------------------------------
00940 inline double vtkMath::RadiansFromDegrees( double x )
00941 {
00942   return x * 0.017453292519943295;
00943 }
00944 
00945 //----------------------------------------------------------------------------
00946 inline float vtkMath::DegreesFromRadians( float x )
00947 {
00948   return x * 57.2957795131f;
00949 }
00950 
00951 //----------------------------------------------------------------------------
00952 inline double vtkMath::DegreesFromRadians( double x )
00953 {
00954   return x * 57.29577951308232;
00955 }
00956 
00957 #ifndef VTK_LEGACY_REMOVE
00958 //----------------------------------------------------------------------------
00959 inline float vtkMath::DegreesToRadians()
00960 {
00961   VTK_LEGACY_REPLACED_BODY(vtkMath::DegreesToRadians, "VTK 5.4",
00962                            vtkMath::RadiansFromDegrees);
00963 
00964   return vtkMath::RadiansFromDegrees( 1.f );
00965 }
00966 
00967 //----------------------------------------------------------------------------
00968 inline double vtkMath::DoubleDegreesToRadians()
00969 {
00970   VTK_LEGACY_REPLACED_BODY(vtkMath::DoubleDegreesToRadians, "VTK 5.4",
00971                            vtkMath::RadiansFromDegrees);
00972 
00973 
00974   return vtkMath::RadiansFromDegrees( 1. );
00975 }
00976 
00977 //----------------------------------------------------------------------------
00978 inline float vtkMath::RadiansToDegrees()
00979 {
00980   VTK_LEGACY_REPLACED_BODY(vtkMath::RadiansToDegrees, "VTK 5.4",
00981                            vtkMath::DegreesFromRadians);
00982 
00983   return vtkMath::DegreesFromRadians( 1.f );
00984 }
00985 
00986 //----------------------------------------------------------------------------
00987 inline double vtkMath::DoubleRadiansToDegrees()
00988 {
00989   VTK_LEGACY_REPLACED_BODY(vtkMath::DoubleRadiansToDegrees, "VTK 5.4",
00990                            vtkMath::DegreesFromRadians);
00991 
00992   return vtkMath::DegreesFromRadians( 1. );
00993 }
00994 #endif // #ifndef VTK_LEGACY_REMOVE
00995 
00996 //----------------------------------------------------------------------------
00997 inline vtkTypeInt64 vtkMath::Factorial( int N )
00998 {
00999   vtkTypeInt64 r = 1;
01000   while ( N > 1 )
01001     {
01002     r *= N--;
01003     }
01004   return r;
01005 }
01006 
01007 //----------------------------------------------------------------------------
01008 inline int vtkMath::Floor(double x)
01009 {
01010 #if defined i386 || defined _M_IX86
01011   union { int i[2]; double d; } u;
01012   // use 52-bit precision of IEEE double to round (x - 0.25) to 
01013   // the nearest multiple of 0.5, according to prevailing rounding
01014   // mode which is IEEE round-to-nearest,even
01015   u.d = (x - 0.25) + 3377699720527872.0; // (2**51)*1.5
01016   // extract mantissa, use shift to divide by 2 and hence get rid
01017   // of the bit that gets messed up because the FPU uses
01018   // round-to-nearest,even mode instead of round-to-nearest,+infinity
01019   return u.i[0] >> 1;
01020 #else
01021   return static_cast<int>(floor( x ) );
01022 #endif
01023 }
01024 
01025 //----------------------------------------------------------------------------
01026 inline float vtkMath::Normalize(float x[3])
01027 {
01028   float den; 
01029   if ( ( den = vtkMath::Norm( x ) ) != 0.0 )
01030     {
01031     for (int i=0; i < 3; i++)
01032       {
01033       x[i] /= den;
01034       }
01035     }
01036   return den;
01037 }
01038 
01039 //----------------------------------------------------------------------------
01040 inline double vtkMath::Normalize(double x[3])
01041 {
01042   double den; 
01043   if ( ( den = vtkMath::Norm( x ) ) != 0.0 )
01044     {
01045     for (int i=0; i < 3; i++)
01046       {
01047       x[i] /= den;
01048       }
01049     }
01050   return den;
01051 }
01052 
01053 //----------------------------------------------------------------------------
01054 inline float vtkMath::Normalize2D(float x[3])
01055 {
01056   float den; 
01057   if ( ( den = vtkMath::Norm2D( x ) ) != 0.0 )
01058     {
01059     for (int i=0; i < 2; i++)
01060       {
01061       x[i] /= den;
01062       }
01063     }
01064   return den;
01065 }
01066 
01067 //----------------------------------------------------------------------------
01068 inline double vtkMath::Normalize2D(double x[3])
01069 {
01070   double den; 
01071   if ( ( den = vtkMath::Norm2D( x ) ) != 0.0 )
01072     {
01073     for (int i=0; i < 2; i++)
01074       {
01075       x[i] /= den;
01076       }
01077     }
01078   return den;
01079 }
01080 
01081 //----------------------------------------------------------------------------
01082 inline float vtkMath::Determinant3x3(const float c1[3], 
01083                                      const float c2[3], 
01084                                      const float c3[3])
01085 {
01086   return c1[0] * c2[1] * c3[2] + c2[0] * c3[1] * c1[2] + c3[0] * c1[1] * c2[2] -
01087          c1[0] * c3[1] * c2[2] - c2[0] * c1[1] * c3[2] - c3[0] * c2[1] * c1[2];
01088 }
01089 
01090 //----------------------------------------------------------------------------
01091 inline double vtkMath::Determinant3x3(const double c1[3], 
01092                                       const double c2[3], 
01093                                       const double c3[3])
01094 {
01095   return c1[0] * c2[1] * c3[2] + c2[0] * c3[1] * c1[2] + c3[0] * c1[1] * c2[2] -
01096          c1[0] * c3[1] * c2[2] - c2[0] * c1[1] * c3[2] - c3[0] * c2[1] * c1[2];
01097 }
01098 
01099 //----------------------------------------------------------------------------
01100 inline double vtkMath::Determinant3x3(double a1, double a2, double a3, 
01101                                       double b1, double b2, double b3, 
01102                                       double c1, double c2, double c3)
01103 {
01104     return ( a1 * vtkMath::Determinant2x2( b2, b3, c2, c3 )
01105            - b1 * vtkMath::Determinant2x2( a2, a3, c2, c3 )
01106            + c1 * vtkMath::Determinant2x2( a2, a3, b2, b3 ) );
01107 }
01108 
01109 //----------------------------------------------------------------------------
01110 inline float vtkMath::Distance2BetweenPoints(const float x[3], 
01111                                              const float y[3])
01112 {
01113   return ( ( x[0] - y[0] ) * ( x[0] - y[0] ) 
01114            + ( x[1] - y[1] ) * ( x[1] - y[1] ) 
01115            + ( x[2] - y[2] ) * ( x[2] - y[2] ) );
01116 }
01117 
01118 //----------------------------------------------------------------------------
01119 inline double vtkMath::Distance2BetweenPoints(const double x[3], 
01120                                               const double y[3])
01121 {
01122   return ( ( x[0] - y[0] ) * ( x[0] - y[0] ) 
01123            + ( x[1] - y[1] ) * ( x[1] - y[1] ) 
01124            + ( x[2] - y[2] ) * ( x[2] - y[2] ) );
01125 }
01126 
01127 //----------------------------------------------------------------------------
01128 // Cross product of two 3-vectors. Result (a x b) is stored in z[3].
01129 inline void vtkMath::Cross(const float x[3], const float y[3], float z[3])
01130 {
01131   float Zx = x[1] * y[2] - x[2] * y[1]; 
01132   float Zy = x[2] * y[0] - x[0] * y[2];
01133   float Zz = x[0] * y[1] - x[1] * y[0];
01134   z[0] = Zx; z[1] = Zy; z[2] = Zz; 
01135 }
01136 
01137 //----------------------------------------------------------------------------
01138 // Cross product of two 3-vectors. Result (a x b) is stored in z[3].
01139 inline void vtkMath::Cross(const double x[3], const double y[3], double z[3])
01140 {
01141   double Zx = x[1] * y[2] - x[2] * y[1]; 
01142   double Zy = x[2] * y[0] - x[0] * y[2];
01143   double Zz = x[0] * y[1] - x[1] * y[0];
01144   z[0] = Zx; z[1] = Zy; z[2] = Zz; 
01145 }
01146 
01147 //BTX
01148 //----------------------------------------------------------------------------
01149 template<class T>
01150 inline double vtkDeterminant3x3(T A[3][3])
01151 {
01152   return A[0][0] * A[1][1] * A[2][2] + A[1][0] * A[2][1] * A[0][2] + 
01153          A[2][0] * A[0][1] * A[1][2] - A[0][0] * A[2][1] * A[1][2] - 
01154          A[1][0] * A[0][1] * A[2][2] - A[2][0] * A[1][1] * A[0][2];
01155 }
01156 //ETX
01157 
01158 //----------------------------------------------------------------------------
01159 inline double vtkMath::Determinant3x3(float A[3][3])
01160 {
01161   return vtkDeterminant3x3( A );
01162 }
01163 
01164 //----------------------------------------------------------------------------
01165 inline double vtkMath::Determinant3x3(double A[3][3])
01166 {
01167   return vtkDeterminant3x3( A );
01168 }
01169 
01170 //----------------------------------------------------------------------------
01171 inline double* vtkMath::SolveCubic(double c0, double c1, double c2, double c3)
01172 {
01173   return vtkPolynomialSolversUnivariate::SolveCubic( c0, c1, c2, c3 );
01174 }
01175 
01176 //----------------------------------------------------------------------------
01177 inline double* vtkMath::SolveQuadratic(double c0, double c1, double c2)
01178 {
01179   return vtkPolynomialSolversUnivariate::SolveQuadratic( c0, c1, c2 );
01180 }
01181 
01182 //----------------------------------------------------------------------------
01183 inline double* vtkMath::SolveLinear(double c0, double c1)
01184 {
01185   return vtkPolynomialSolversUnivariate::SolveLinear( c0, c1 );
01186 }
01187 
01188 //----------------------------------------------------------------------------
01189 inline int vtkMath::SolveCubic(double c0, double c1, double c2, double c3, 
01190                                double *r1, double *r2, double *r3, int *num_roots)
01191 {
01192   return vtkPolynomialSolversUnivariate::SolveCubic( c0, c1, c2, c3, r1, r2, r3, num_roots );
01193 }
01194 
01195 //----------------------------------------------------------------------------
01196 inline int vtkMath::SolveQuadratic(double c0, double c1, double c2, 
01197                                    double *r1, double *r2, int *num_roots)
01198 {
01199   return vtkPolynomialSolversUnivariate::SolveQuadratic( c0, c1, c2, r1, r2, num_roots );
01200 }
01201   
01202 //----------------------------------------------------------------------------
01203 inline int vtkMath::SolveQuadratic( double* c, double* r, int* m )
01204 {
01205   return vtkPolynomialSolversUnivariate::SolveQuadratic( c, r, m );
01206 }
01207 
01208 //----------------------------------------------------------------------------
01209 inline int vtkMath::SolveLinear(double c0, double c1, double *r1, int *num_roots)
01210 {
01211   return vtkPolynomialSolversUnivariate::SolveLinear( c0, c1, r1, num_roots );
01212 }
01213 
01214 //----------------------------------------------------------------------------
01215 inline void vtkMath::ClampValue(double *value, const double range[2])
01216 {
01217   if (value && range)
01218     {
01219     if (*value < range[0])
01220       {
01221       *value = range[0];
01222       }
01223     else if (*value > range[1])
01224       {
01225       *value = range[1];
01226       }
01227     }
01228 }
01229 
01230 //----------------------------------------------------------------------------
01231 inline void vtkMath::ClampValue(
01232   double value, const double range[2], double *clamped_value)
01233 {
01234   if (range && clamped_value)
01235     {
01236     if (value < range[0])
01237       {
01238       *clamped_value = range[0];
01239       }
01240     else if (value > range[1])
01241       {
01242       *clamped_value = range[1];
01243       }
01244     else
01245       {
01246       *clamped_value = value;
01247       }
01248     }
01249 }
01250 
01251 // ---------------------------------------------------------------------------
01252 inline double vtkMath::ClampAndNormalizeValue(double value,
01253                                               const double range[2])
01254 {
01255   assert("pre: valid_range" && range[0]<=range[1]);
01256   
01257   double result;
01258   if(range[0]==range[1])
01259     {
01260       result=0.0;
01261     }
01262   else
01263     {
01264       // clamp
01265       if(value<range[0])
01266         {
01267           result=range[0];
01268         }
01269       else
01270         {
01271           if(value>range[1])
01272             {
01273               result=range[1];
01274             }
01275           else
01276             {
01277               result=value;
01278             }
01279         }
01280 
01281       // normalize
01282       result=( result - range[0] ) / ( range[1] - range[0] );
01283     }
01284 
01285   assert("post: valid_result" && result>=0.0 && result<=1.0);
01286 
01287   return result;
01288 }
01289 
01290 #endif