VTK
|
00001 /*========================================================================= 00002 00003 Program: Visualization Toolkit 00004 Module: $RCSfile: vtkTetra.h,v $ 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 =========================================================================*/ 00032 #ifndef __vtkTetra_h 00033 #define __vtkTetra_h 00034 00035 #include "vtkCell3D.h" 00036 00037 class vtkLine; 00038 class vtkTriangle; 00039 class vtkUnstructuredGrid; 00040 00041 class VTK_FILTERING_EXPORT vtkTetra : public vtkCell3D 00042 { 00043 public: 00044 static vtkTetra *New(); 00045 vtkTypeRevisionMacro(vtkTetra,vtkCell3D); 00046 void PrintSelf(ostream& os, vtkIndent indent); 00047 00049 00050 virtual void GetEdgePoints(int edgeId, int* &pts); 00051 virtual void GetFacePoints(int faceId, int* &pts); 00053 00055 00056 int GetCellType() {return VTK_TETRA;} 00057 int GetNumberOfEdges() {return 6;} 00058 int GetNumberOfFaces() {return 4;} 00059 vtkCell *GetEdge(int edgeId); 00060 vtkCell *GetFace(int faceId); 00061 void Contour(double value, vtkDataArray *cellScalars, 00062 vtkPointLocator *locator, vtkCellArray *verts, 00063 vtkCellArray *lines, vtkCellArray *polys, 00064 vtkPointData *inPd, vtkPointData *outPd, 00065 vtkCellData *inCd, vtkIdType cellId, vtkCellData *outCd); 00066 void Clip(double value, vtkDataArray *cellScalars, 00067 vtkPointLocator *locator, vtkCellArray *connectivity, 00068 vtkPointData *inPd, vtkPointData *outPd, 00069 vtkCellData *inCd, vtkIdType cellId, vtkCellData *outCd, 00070 int insideOut); 00071 int EvaluatePosition(double x[3], double* closestPoint, 00072 int& subId, double pcoords[3], 00073 double& dist2, double *weights); 00074 void EvaluateLocation(int& subId, double pcoords[3], double x[3], 00075 double *weights); 00076 int IntersectWithLine(double p1[3], double p2[3], double tol, double& t, 00077 double x[3], double pcoords[3], int& subId); 00078 int Triangulate(int index, vtkIdList *ptIds, vtkPoints *pts); 00079 void Derivatives(int subId, double pcoords[3], double *values, 00080 int dim, double *derivs); 00081 virtual double *GetParametricCoords(); 00083 00087 int CellBoundary(int subId, double pcoords[3], vtkIdList *pts); 00088 00090 int GetParametricCenter(double pcoords[3]); 00091 00094 double GetParametricDistance(double pcoords[3]); 00095 00097 00098 static void TetraCenter(double p1[3], double p2[3], double p3[3], double p4[3], 00099 double center[3]); 00101 00103 00106 static double Circumsphere(double p1[3], double p2[3], double p3[3], 00107 double p4[3], double center[3]); 00109 00111 00114 static double Insphere(double p1[3], double p2[3], double p3[3], 00115 double p4[3], double center[3]); 00117 00119 00130 static int BarycentricCoords(double x[3], double x1[3], double x2[3], 00131 double x3[3], double x4[3], double bcoords[4]); 00133 00135 00137 static double ComputeVolume(double p1[3], double p2[3], double p3[3], 00138 double p4[3]); 00140 00144 int JacobianInverse(double **inverse, double derivs[12]); 00145 00147 00148 static void InterpolationFunctions(double pcoords[3], double weights[4]); 00149 // Description: 00150 // @deprecated Replaced by vtkTetra::InterpolateDerivs as of VTK 5.2 00151 static void InterpolationDerivs(double pcoords[3], double derivs[12]); 00152 // Description: 00153 // Compute the interpolation functions/derivatives 00154 // (aka shape functions/derivatives) 00155 virtual void InterpolateFunctions(double pcoords[3], double weights[4]) 00156 { 00157 vtkTetra::InterpolationFunctions(pcoords,weights); 00158 } 00159 virtual void InterpolateDerivs(double pcoords[3], double derivs[12]) 00160 { 00161 vtkTetra::InterpolationDerivs(pcoords,derivs); 00162 } 00164 00166 00168 static int *GetEdgeArray(int edgeId); 00169 static int *GetFaceArray(int faceId); 00171 00172 protected: 00173 vtkTetra(); 00174 ~vtkTetra(); 00175 00176 vtkLine *Line; 00177 vtkTriangle *Triangle; 00178 00179 private: 00180 vtkTetra(const vtkTetra&); // Not implemented. 00181 void operator=(const vtkTetra&); // Not implemented. 00182 }; 00183 00184 inline int vtkTetra::GetParametricCenter(double pcoords[3]) 00185 { 00186 pcoords[0] = pcoords[1] = pcoords[2] = 0.25; 00187 return 0; 00188 } 00189 00190 #endif 00191 00192 00193