common/vtkMath.h Source File

00001 /*=========================================================================
00002 
00003   Program:   Visualization Toolkit
00004   Module:    $RCSfile: vtkMath.h,v $
00005   Language:  C++
00006 
00007 
00008 Copyright (c) 1993-2001 Ken Martin, Will Schroeder, Bill Lorensen 
00009 All rights reserved.
00010 
00011 Redistribution and use in source and binary forms, with or without
00012 modification, are permitted provided that the following conditions are met:
00013 
00014  * Redistributions of source code must retain the above copyright notice,
00015    this list of conditions and the following disclaimer.
00016 
00017  * Redistributions in binary form must reproduce the above copyright notice,
00018    this list of conditions and the following disclaimer in the documentation
00019    and/or other materials provided with the distribution.
00020 
00021  * Neither name of Ken Martin, Will Schroeder, or Bill Lorensen nor the names
00022    of any contributors may be used to endorse or promote products derived
00023    from this software without specific prior written permission.
00024 
00025  * Modified source versions must be plainly marked as such, and must not be
00026    misrepresented as being the original software.
00027 
00028 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS ``AS IS''
00029 AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00030 IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
00031 ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHORS OR CONTRIBUTORS BE LIABLE FOR
00032 ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
00033 DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
00034 SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
00035 CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
00036 OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
00037 OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00038 
00039 =========================================================================*/
00056 #ifndef __vtkMath_h
00057 #define __vtkMath_h
00058 
00059 #include "vtkObject.h"
00060 
00061 class VTK_EXPORT vtkMath : public vtkObject
00062 {
00063 public:
00064   static vtkMath *New();
00065   vtkTypeMacro(vtkMath,vtkObject);
00066   
00068   static float Pi() {return 3.14159265358979;};
00069   static float DegreesToRadians() {return 0.017453292;};
00070 
00072   static double DoubleDegreesToRadians() {return 0.017453292519943295;};
00073 
00075   static float Dot(const float x[3], const float y[3]) {
00076     return (x[0]*y[0] + x[1]*y[1] + x[2]*y[2]);};
00077 
00079   static double Dot(const double x[3], const double y[3]) {
00080     return (x[0]*y[0] + x[1]*y[1] + x[2]*y[2]);};
00081   
00083   static void Cross(const float x[3], const float y[3], float z[3]);
00084 
00087   static void Cross(const double x[3], const double y[3], double z[3]);
00088 
00090   static float Norm(const float x[3]) {
00091     return sqrt(x[0]*x[0] + x[1]*x[1] + x[2]*x[2]);};
00092   
00094   static double Norm(const double x[3]) {
00095     return sqrt(x[0]*x[0] + x[1]*x[1] + x[2]*x[2]);};
00096   
00098   static float Normalize(float x[3]);
00099 
00102   static double Normalize(double x[3]);
00103 
00109   static void Perpendiculars(const double x[3], double y[3], double z[3], 
00110               double theta);
00111   static void Perpendiculars(const float x[3], float y[3], float z[3],
00112               double theta);
00113 
00115   static float Distance2BetweenPoints(const float x[3], const float y[3]);
00116 
00119   static double Distance2BetweenPoints(const double x[3], const double y[3]);
00120 
00122   static float Dot2D(const float x[3], const float y[3]) {
00123     return (x[0]*y[0] + x[1]*y[1]);};
00124   
00127   static double Dot2D(const double x[3], const double y[3]) {
00128     return (x[0]*y[0] + x[1]*y[1]);};
00129 
00131   static float Norm2D(const float x[3]) {
00132     return sqrt(x[0]*x[0] + x[1]*x[1]);};
00133 
00136   static double Norm2D(const double x[3]) {
00137     return sqrt(x[0]*x[0] + x[1]*x[1]);};
00138 
00141   static float Normalize2D(float x[3]);
00142 
00145   static double Normalize2D(double x[3]);
00146 
00148   static float Determinant2x2(const float c1[2], const float c2[2]) {
00149     return (c1[0]*c2[1] - c2[0]*c1[1]);};
00150 
00152   static double Determinant2x2(double a, double b, double c, double d) {
00153     return (a * d - b * c);};
00154 
00157   static void LUFactor3x3(float A[3][3], int index[3]);
00158   static void LUFactor3x3(double A[3][3], int index[3]);
00159 
00162   static void LUSolve3x3(const float A[3][3], const int index[3], 
00163           float x[3]);
00164   static void LUSolve3x3(const double A[3][3], const int index[3], 
00165           double x[3]);
00166 
00169   static void LinearSolve3x3(const float A[3][3], const float x[3], 
00170               float y[3]);
00171   static void LinearSolve3x3(const double A[3][3], const double x[3], 
00172               double y[3]);
00173 
00175   static void Multiply3x3(const float A[3][3], const float x[3], float y[3]);
00176   static void Multiply3x3(const double A[3][3], const double x[3], 
00177            double y[3]);
00178   
00181   static void Multiply3x3(const float A[3][3], const float B[3][3], 
00182            float C[3][3]);
00183   static void Multiply3x3(const double A[3][3], const double B[3][3], 
00184            double C[3][3]);
00185 
00187   static void Transpose3x3(const float A[3][3], float AT[3][3]);
00188   static void Transpose3x3(const double A[3][3], double AT[3][3]);
00189 
00191   static void Invert3x3(const float A[3][3], float AI[3][3]);
00192   static void Invert3x3(const double A[3][3], double AI[3][3]);
00193 
00195   static void Identity3x3(float A[3][3]);
00196   static void Identity3x3(double A[3][3]);
00197 
00199   static double Determinant3x3(float A[3][3]);
00200   static double Determinant3x3(double A[3][3]);
00201 
00203   static float Determinant3x3(const float c1[3], 
00204                const float c2[3], 
00205                const float c3[3]);
00206 
00209   static double Determinant3x3(double a1, double a2, double a3, 
00210                 double b1, double b2, double b3, 
00211                 double c1, double c2, double c3);
00212 
00215   static void QuaternionToMatrix3x3(const float quat[4], float A[3][3]); 
00216   static void QuaternionToMatrix3x3(const double quat[4], double A[3][3]); 
00217 
00221   static void Matrix3x3ToQuaternion(const float A[3][3], float quat[4]);
00222   static void Matrix3x3ToQuaternion(const double A[3][3], double quat[4]);
00223   
00228   static void Orthogonalize3x3(const float A[3][3], float B[3][3]);
00229   static void Orthogonalize3x3(const double A[3][3], double B[3][3]);
00230 
00235   static void Diagonalize3x3(const float A[3][3], float w[3], float V[3][3]);
00236   static void Diagonalize3x3(const double A[3][3],double w[3],double V[3][3]);
00237 
00242   static void SingularValueDecomposition3x3(const float A[3][3],
00243                    float U[3][3], float w[3],
00244                    float VT[3][3]);
00245   static void SingularValueDecomposition3x3(const double A[3][3],
00246                    double U[3][3], double w[3],
00247                    double VT[3][3]);
00248   
00253   static int SolveLinearSystem(double **A, double *x, int size);
00254 
00258   static int InvertMatrix(double **A, double **AI, int size);
00259 
00262   static int InvertMatrix(double **A, double **AI, int size,
00263            int *tmp1Size, double *tmp2Size);
00264 
00270   static int LUFactorLinearSystem(double **A, int *index, int size);
00271 
00274   static int LUFactorLinearSystem(double **A, int *index, int size,
00275               double *tmpSize);
00276 
00283   static void LUSolveLinearSystem(double **A, int *index, 
00284               double *x, int size);
00285 
00293   static double EstimateMatrixCondition(double **A, int size);
00294 
00300   static void RandomSeed(long s);  
00301 
00304   static float Random();  
00305 
00307   static float Random(float min, float max);
00308 
00313   static int Jacobi(float **a, float *w, float **v);
00314   static int Jacobi(double **a, double *w, double **v);
00315 
00321   static int JacobiN(float **a, int n, float *w, float **v);
00322   static int JacobiN(double **a, int n, double *w, double **v);
00323 
00330   static double* SolveCubic(double c0, double c1, double c2, double c3);
00331 
00338   static double* SolveQuadratic(double c0, double c1, double c2);
00339 
00343   static double* SolveLinear(double c0, double c1);
00344 
00356   static int SolveCubic(double c0, double c1, double c2, double c3, 
00357          double *r1, double *r2, double *r3, int *num_roots);
00358 
00363   static int SolveQuadratic(double c0, double c1, double c2, 
00364              double *r1, double *r2, int *num_roots);
00365   
00370   static int SolveLinear(double c0, double c1, double *r1, int *num_roots);
00371 
00372 
00380   static int SolveLeastSquares(int numberOfSamples, double **xt, int xOrder,
00381                                double **yt, int yOrder, double **mt);
00382   
00383 protected:
00384   vtkMath() {};
00385   ~vtkMath() {};
00386   vtkMath(const vtkMath&) {};
00387   void operator=(const vtkMath&) {};
00388   
00389   static long Seed;
00390 };
00391 
00392 inline float vtkMath::Normalize(float x[3])
00393 {
00394   float den; 
00395   if ( (den = vtkMath::Norm(x)) != 0.0 )
00396     {
00397     for (int i=0; i < 3; i++)
00398       {
00399       x[i] /= den;
00400       }
00401     }
00402   return den;
00403 }
00404 inline double vtkMath::Normalize(double x[3])
00405 {
00406   double den; 
00407   if ( (den = vtkMath::Norm(x)) != 0.0 )
00408     {
00409     for (int i=0; i < 3; i++)
00410       {
00411       x[i] /= den;
00412       }
00413     }
00414   return den;
00415 }
00416 
00417 inline float vtkMath::Normalize2D(float x[3])
00418 {
00419   float den; 
00420   if ( (den = vtkMath::Norm2D(x)) != 0.0 )
00421     {
00422     for (int i=0; i < 2; i++)
00423       {
00424       x[i] /= den;
00425       }
00426     }
00427   return den;
00428 }
00429 
00430 inline double vtkMath::Normalize2D(double x[3])
00431 {
00432   double den; 
00433   if ( (den = vtkMath::Norm2D(x)) != 0.0 )
00434     {
00435     for (int i=0; i < 2; i++)
00436       {
00437       x[i] /= den;
00438       }
00439     }
00440   return den;
00441 }
00442 
00443 inline float vtkMath::Determinant3x3(const float c1[3], 
00444                  const float c2[3], 
00445                  const float c3[3])
00446 {
00447   return c1[0]*c2[1]*c3[2] + c2[0]*c3[1]*c1[2] + c3[0]*c1[1]*c2[2] -
00448          c1[0]*c3[1]*c2[2] - c2[0]*c1[1]*c3[2] - c3[0]*c2[1]*c1[2];
00449 }
00450 
00451 inline double vtkMath::Determinant3x3(double a1, double a2, double a3, 
00452                   double b1, double b2, double b3, 
00453                   double c1, double c2, double c3)
00454 {
00455     return ( a1 * vtkMath::Determinant2x2( b2, b3, c2, c3 )
00456       - b1 * vtkMath::Determinant2x2( a2, a3, c2, c3 )
00457            + c1 * vtkMath::Determinant2x2( a2, a3, b2, b3 ) );
00458 }
00459 
00460 inline float vtkMath::Distance2BetweenPoints(const float x[3], 
00461                     const float y[3])
00462 {
00463   return ((x[0]-y[0])*(x[0]-y[0]) + (x[1]-y[1])*(x[1]-y[1]) +
00464           (x[2]-y[2])*(x[2]-y[2]));
00465 }
00466 inline double vtkMath::Distance2BetweenPoints(const double x[3], 
00467                      const double y[3])
00468 {
00469   return ((x[0]-y[0])*(x[0]-y[0]) + (x[1]-y[1])*(x[1]-y[1]) +
00470           (x[2]-y[2])*(x[2]-y[2]));
00471 }
00472 
00473 inline float vtkMath::Random(float min, float max)
00474 {
00475   return (min + vtkMath::Random()*(max-min));
00476 }
00477 
00478 // Cross product of two 3-vectors. Result vector in z[3].
00479 inline void vtkMath::Cross(const float x[3], const float y[3], float z[3])
00480 {
00481   float Zx = x[1]*y[2] - x[2]*y[1]; 
00482   float Zy = x[2]*y[0] - x[0]*y[2];
00483   float Zz = x[0]*y[1] - x[1]*y[0];
00484   z[0] = Zx; z[1] = Zy; z[2] = Zz; 
00485 }
00486 
00487 // Cross product of two 3-vectors. Result vector in z[3].
00488 inline void vtkMath::Cross(const double x[3], const double y[3], double z[3])
00489 {
00490   double Zx = x[1]*y[2] - x[2]*y[1]; 
00491   double Zy = x[2]*y[0] - x[0]*y[2];
00492   double Zz = x[0]*y[1] - x[1]*y[0];
00493   z[0] = Zx; z[1] = Zy; z[2] = Zz; 
00494 }
00495 
00496 //BTX
00497 //----------------------------------------------------------------------------
00498 template<class T>
00499 static inline double vtkDeterminant3x3(T A[3][3])
00500 {
00501   return A[0][0]*A[1][1]*A[2][2] + A[1][0]*A[2][1]*A[0][2] + 
00502          A[2][0]*A[0][1]*A[1][2] - A[0][0]*A[2][1]*A[1][2] - 
00503          A[1][0]*A[0][1]*A[2][2] - A[2][0]*A[1][1]*A[0][2];
00504 }
00505 //ETX
00506 
00507 inline double vtkMath::Determinant3x3(float A[3][3])
00508 {
00509   return vtkDeterminant3x3(A);
00510 }
00511 
00512 inline double vtkMath::Determinant3x3(double A[3][3])
00513 {
00514   return vtkDeterminant3x3(A);
00515 }
00516 
00517 
00518 #endif