[vtkusers] question about vtkMeshQuality.cxx

Mathieu Malaterre mathieu.malaterre at kitware.com
Tue Jul 20 09:59:14 EDT 2004


Alexandre,

	What is the license attach to this code ?

Thanks
Mathieu

Alexandre gouaillard wrote:
> Indeed , vtkMeshQuality was not enough for us too.
> 
> The problem with the "quality" of a cell is that many measures exist in 
> the litterature depending on the application.
> As philippe seems interested in making something additions, I 
> implemented one of the triangle shape measures existing :
> 1/2 perimeter * longest edge / area
> normalized by the value for an equilateral triangle.
> It measure the "stretching factor" of a triangle.
> 
> please find the code attached.
> it's not perfect, some checking steps were never implemented, but in 
> those cases, a comment line exists to remind that a check should be made.
> If it can be of any help, use it.
> 
> alex.
> 
> 
> A 11:05 19/07/2004 -0700, Philippe P. Pebay a ecrit :
> 
>> Hi
>>
>> I was looking at using vtkMeshQuality for some paper results, but it 
>> seems that the name of this class is too generic, because:
>> 1. it only addresses tetrahedral meshes
>> 2. in the case of tetrahedral meshes, there are many quality measures 
>> that can be used, but this class only provides the "radius-ratio" 
>> (normalized circumradius / inradius).
>>
>> I haven't seen any other quality metrics in VTK; would anyone mind if 
>> we revise vtkMeshQuality to compute various quality measures on 
>> several different types of cells ?  We would like to compute the 
>> following quality measures:
>> 1. tetrahedral meshes: radius-ratio, aspect-ratio, shape measures 
>> (Frobenius-like), etc.
>> 2. hexahedral meshes
>> 3. triangular and quadrilateral elements: various measures in the 
>> literature.
>>
>> Also, I noticed that the Execute method in vtkMeshQuality.cxx declares 
>> p# and dp# (#= 1,2,3,4) as double arrays, and then copies p's into 
>> dp's. But it appears that dp's are never altered, so why do 2 copies 
>> of the points exist ?
>>
>> Cheers
>> -- 
>> ===============================================================
>> Philippe P. Pebay             |  pppebay at ca.sandia.gov
>> Sandia National Laboratories  |  Phone: +1 925 294 2024
>> PO Box 969, MS 9051           |  Cell : +1 925 784 3255
>> Livermore, CA 94550 U.S.A.    |  Fax  : +1 925 294 2595
>>            www.ca.sandia.gov/CRF/staff/pebay.html
>> ===============================================================
>> "Success is the ability to go from one failure
>> to another with no loss of enthusiasm."
>> -- Sir Winston Churchill
>>
>>
>> _______________________________________________
>> This is the private VTK discussion list. Please keep messages 
>> on-topic. Check the FAQ at: <http://public.kitware.com/cgi-bin/vtkfaq>
>> Follow this link to subscribe/unsubscribe:
>> http://www.vtk.org/mailman/listinfo/vtkusers
> 
> 
> ------------------------------------------------------------------------
> 
> /*=========================================================================
> 
>   Program:   Mailleur 3D multi-résolution (Creatis PFE  2001~)
>   Module:    vtkTriangleShapeFactor.h
>   Language:  C++
>   Date:      2003/03
>   Auteurs:   ALex. GOUAILLARD
> 
> =========================================================================*/
> // .NAME vtkTriangleShapeFactor 
> // .SECTION Description
> // This class compute the mean shape factor 
> // of all the triangles of a given PolyData
> 
> #ifndef __vtkTriangleShapeFactor_h
> #define __vtkTriangleShapeFactor_h
> 
> #include <vtkObject.h>
> #include <vtkPolyData.h>
> #include <vtkMath.h>
> 
> class VTK_EXPORT vtkTriangleShapeFactor : public vtkObject 
> {
> 
> public:
> 
>   // VTK Constructor 
>   static vtkTriangleShapeFactor *New();
>   
>   // VTK inheritance scheme
>   vtkTypeMacro(vtkTriangleShapeFactor,vtkObject);
> 
>   // To get the result as a float
>   float GetMeanShapeFactor() {return this->MeanShapeFactor;};
> 
>   //To get thie min angle of 
>   static double GetMinAngle(vtkPolyData *input,int IdCell);
> 
>   // Set a PolyData as input
>   void SetInput(vtkPolyData *UserSetInput) {this->Input = UserSetInput;};  
> 
>   // compute the value and put it in this->MeanShapeFactor
>   void Update();
> 
> protected:
>   vtkTriangleShapeFactor(); //constructeur
>   ~vtkTriangleShapeFactor() {}; // destructeur
>   vtkTriangleShapeFactor(const vtkTriangleShapeFactor&) {};
>   void operator=(const vtkTriangleShapeFactor&) {};
> 
>   float MeanShapeFactor;
> 
>   vtkPolyData *Input;
> 
> };
> 
> #endif
> 
> 
> ------------------------------------------------------------------------
> 
> /*=========================================================================
> 
>   Program:   Mailleur 3D multi-résolution (Creatis PFE  2001~)
>   Module:    vtkTriangleShapeFactor.cxx
>   Language:  C++
>   Date:      2003/03
>   Auteurs:   Alex. GOUAILLARD
> 
> =========================================================================*/
> #include <vtkObjectFactory.h>
> #include <vtkMath.h>
> #include <vtkTriangle.h>
> 
> #include "vtkTriangleShapeFactor.h"
> #include "vtkParameterization.h"
> //--------------------------------------------------------------------------
> 
> vtkTriangleShapeFactor* vtkTriangleShapeFactor::New()
> {
>   // First try to create the object from the vtkObjectFactory
>   vtkObject* ret = vtkObjectFactory::CreateInstance("vtkTriangleShapeFactor");
>   if(ret)
>     {
>     return (vtkTriangleShapeFactor*)ret;
>     }
>   // If the factory was unable to create the object, then create it here.
>   return new vtkTriangleShapeFactor;
> }
> 
> 
> 
> 
> 
> vtkTriangleShapeFactor::vtkTriangleShapeFactor()
> {
> 	this->Input = NULL;
> 	this->MeanShapeFactor = 0.0;
> }
> 
> 
> 
> 
> 
> void vtkTriangleShapeFactor::Update() 
> {
> 	// check inputs
> 	if (!this->Input) return;
> 	
> 	// check if it's a triangular mesh
> 
> 	// name variables
> 	double QualitySum;
> 	double RegularTriangleShapeFactor;
> 	double LongestEdgeLength;
> 	double TempLength;
> 	double Perimeter;
> 
> 	int NbCells;
> 	int IdVertex1;
> 	int IdVertex2;
> 	int IdVertex3;
> 	
> 	register int i;
> 
> 	double Vertex1[3];
> 	double Vertex2[3];
> 	double Vertex3[3];
> 
> 	// initialization
> 	NbCells = this->Input->GetNumberOfCells();
> 	RegularTriangleShapeFactor = 2.0 * sqrt(3.0);
> 	QualitySum = 0.0;
> 
> 	// loop on all the cells and compute the unnormalyzed shape factor
> 	for (i=0;i<NbCells;i++)
> 	{
> 		// Initialyze
> 		LongestEdgeLength = 0.0;
> 		Perimeter = 0.0;
> 
> 		// get the points coordinates
> 		vtkParameterization::GetFaceVertices(this->Input,i,IdVertex1,IdVertex2,IdVertex3);
> 		this->Input->GetPoint(IdVertex1,Vertex1);
> 		this->Input->GetPoint(IdVertex2,Vertex2);
> 		this->Input->GetPoint(IdVertex3,Vertex3);
> 
> 		// compute the length of each edges
> 		TempLength = sqrt(vtkMath::Distance2BetweenPoints(Vertex1,Vertex2));
> 		if(TempLength > LongestEdgeLength) LongestEdgeLength = TempLength;
> 		Perimeter += TempLength;
> 		
> 		TempLength = sqrt(vtkMath::Distance2BetweenPoints(Vertex1,Vertex3));
> 		if(TempLength > LongestEdgeLength) LongestEdgeLength = TempLength;
> 		Perimeter += TempLength;
> 
> 		TempLength = sqrt(vtkMath::Distance2BetweenPoints(Vertex3,Vertex2));
> 		if(TempLength > LongestEdgeLength) LongestEdgeLength = TempLength;
> 		Perimeter += TempLength;
> 
> 		// shape factor = half perimeter * longest edge / area.
> 		QualitySum += 0.5 * Perimeter * LongestEdgeLength / vtkTriangle::TriangleArea(Vertex1,Vertex2,Vertex3);
> 
> 	}
> 
> 	// normalize
> 	this->MeanShapeFactor = NbCells * RegularTriangleShapeFactor / QualitySum;
> 
> }
> 
> double vtkTriangleShapeFactor::GetMinAngle(vtkPolyData *input,int IdCell)
> {
> 	
> 	register int i;
> 
> 	int Id[3];
> 	
> 	double A[3];
> 	double B[3];
> 	double C[3];
> 
> 	double BA[3];
> 	double BC[3];
> 	double AC[3];
> 
> 	double dot[3];
> 	double cross[3][3];
> 	double tan_teta[3];
> 
> 	double min;
> 
> 	for(i=0;i<3;i++)
> 		tan_teta[i] = -1;
> 	
> 	vtkParameterization::GetFaceVertices(input,IdCell,Id[0],Id[1],Id[2]);
> 	input->GetPoint(Id[0],A);
> 	input->GetPoint(Id[1],B);
> 	input->GetPoint(Id[2],C);
> 
> 	for(i=0;i<3;i++)
> 	{
> 		BA[i] = A[i] - B[i];
> 		BC[i] = C[i] - B[i];
> 		AC[i] = C[i] - A[i];
> 	}
> 
> 	dot[0] = -1.0 * vtkMath::Dot(BA,AC);
> 	dot[1] = vtkMath::Dot(BA,BC);
> 	dot[2] = vtkMath::Dot(AC,BC);
> 
> 	vtkMath::Cross(BA,AC,cross[0]);
> 	vtkMath::Cross(BA,BC,cross[1]);
> 	vtkMath::Cross(AC,BC,cross[2]);
> 
> 	for(i=0;i<3;i++)
> 	{
> 		if(dot[i] < 0)
> 			tan_teta[i] = -1.0 * vtkMath::Norm(cross[i]) / dot[i];
> 		else if(dot[i] > 0)
> 			tan_teta[i] = vtkMath::Norm(cross[i]) / dot[i];
> 	}
> 
> 	if(tan_teta[0] != -1) //if(teta[0] != pi/2)
> 		min = tan_teta[0];
> 		
> 	else 
> 		min = tan_teta[1];
> 
> 	for(i=0;i<3;i++)
> 		if((tan_teta[i] < min) && (tan_teta[i] > 0))
> 			min = tan_teta[i];
> 
> 	return(min);
> }
> 
> 		
> 	
> 
> 	
> 
> 
> ------------------------------------------------------------------------
> 
> _______________________________________________
> This is the private VTK discussion list. 
> Please keep messages on-topic. Check the FAQ at: <http://public.kitware.com/cgi-bin/vtkfaq>
> Follow this link to subscribe/unsubscribe:
> http://www.vtk.org/mailman/listinfo/vtkusers






More information about the vtkusers mailing list