[vtkusers] computing the volume of a simple multi triangles polydata

Michele michele.conconi at unibo.it
Thu Jun 6 10:33:22 EDT 2013 

Hi everybody,

I am tring to build a polydata starting from an ordered set of points in
space in order to evaluate its volume.

I built an exapmple (see below) which generate a correct visual
representation of the polydata, but the computed volume is wrong (6.8 rather
then 48 mm^3).

I tried aldo the vtkDelaunay2D, with no better results.

What am I doing wrong?
Thanks in advance for the help

Michele

------------------------------


#include <vtkVersion.h>
#include <vtkSmartPointer.h>
#include <vtkCellArray.h>
#include <vtkCellData.h>
#include <vtkDoubleArray.h>
#include <vtkPoints.h>
#include <vtkPolyLine.h>
#include <vtkPolyData.h>
#include <vtkPolyDataMapper.h>
#include <vtkActor.h>
#include <vtkRenderWindow.h>
#include <vtkRenderer.h>
#include <vtkRenderWindowInteractor.h>
#include <vtkMassProperties.h>
#include <vtkTriangle.h>



int main( int argc, char *argv[] )
{
	// build the points set to be triangulated
	int Np = 5;
	double **P1, **P2, **P3;

	P1 = new double*[Np];
	P2 = new double*[Np];
	P3 = new double*[Np];

	for(int i = 0; i < Np; i++)
	{
		P1[i] = new double[4];
		P2[i] = new double[4];
		P3[i] = new double[4];
	}

	for(int i=0;i<Np;i++)
	{
		P1[i][0]	= 0+i;
		P1[i][1]	= 1;
		P1[i][2]	= 0;

		P2[i][0]	= 0+i;
		P2[i][1]	= 3;
		P2[i][2]	= 2;

		P3[i][0]	= 0+i;
		P3[i][1]	= 5;
		P3[i][2]	= 4;
	}


	vtkSmartPointer&lt;vtkPoints> points = vtkSmartPointer<vtkPoints>::New();

	for(int i=0;i<Np;i++)
		points->InsertNextPoint(P1[i]);

	for(int i=0;i<Np;i++)
		points->InsertNextPoint(P2[i]);


	for(int i=0;i<Np;i++)
		points->InsertNextPoint(P3[i]);

	for(int i=0;i<Np;i++)
		points->InsertNextPoint(P1[i][0],0,P1[i][2]);


	for(int i=0;i<Np;i++)
		points->InsertNextPoint(P2[i][0],0,P2[i][2]);

	for(int i=0;i<Np;i++)
		points->InsertNextPoint(P3[i][0],0,P3[i][2]);


	//build manually the triangle

	vtkSmartPointer<vtkCellArray> triangles =
vtkSmartPointer<vtkCellArray>::New();
	vtkSmartPointer<vtkTriangle> triangle =
vtkSmartPointer<vtkTriangle>::New();
    
	// ------------------------------------- base
	for(int j=3;j<5;j++)
		for(int i = (j*Np)+1; i < (j+1)*Np; i++)
		{
			triangle->Initialize();
			triangle->GetPointIds()->SetId(0,i-1);
			triangle->GetPointIds()->SetId(1,i);
			triangle->GetPointIds()->SetId(2,Np+i-1);
			triangles ->InsertNextCell(triangle);
		}

	for(int j=4;j<6;j++)
		for(int i = j*Np; i < (j+1)*Np-1; i++)
		{
			triangle->Initialize();
			triangle->GetPointIds()->SetId(0,i);
			triangle->GetPointIds()->SetId(1,i+1);
			triangle->GetPointIds()->SetId(2,i-Np+1);
			triangles ->InsertNextCell(triangle);
		}


	// ---------------------------------------- side 1 

	for(int i = 1; i < Np; i++)
		{
			triangle->Initialize();
			triangle->GetPointIds()->SetId(0,i-1);
			triangle->GetPointIds()->SetId(1,i);
			triangle->GetPointIds()->SetId(2,3*Np+i-1);
			triangles ->InsertNextCell(triangle);
		}

	for(int i = 3*Np; i < 4*Np-1; i++)
		{
			triangle->Initialize();
			triangle->GetPointIds()->SetId(0,i);
			triangle->GetPointIds()->SetId(1,i+1);
			triangle->GetPointIds()->SetId(2,i-3*Np+1);
			triangles ->InsertNextCell(triangle);
		}

	// ---------------------------------------- side 2 
	for(int i = 2*Np+1; i < 3*Np; i++)
		{
			triangle->Initialize();
			triangle->GetPointIds()->SetId(0,i-1);
			triangle->GetPointIds()->SetId(1,i);
			triangle->GetPointIds()->SetId(2,3*Np+i-1);
			triangles ->InsertNextCell(triangle);
		}

	for(int i = 5*Np; i < 6*Np-1; i++)
		{
			triangle->Initialize();
			triangle->GetPointIds()->SetId(0,i);
			triangle->GetPointIds()->SetId(1,i+1);
			triangle->GetPointIds()->SetId(2,i-3*Np+1);
			triangles ->InsertNextCell(triangle);
		}

	// ---------------------------------------- side 3 

			triangle->Initialize();
			triangle->GetPointIds()->SetId(0,0);
			triangle->GetPointIds()->SetId(1,3*Np);
			triangle->GetPointIds()->SetId(2,4*Np);
			triangles ->InsertNextCell(triangle);

			triangle->Initialize();
			triangle->GetPointIds()->SetId(0,Np);
			triangle->GetPointIds()->SetId(1,4*Np);
			triangle->GetPointIds()->SetId(2,5*Np);
			triangles ->InsertNextCell(triangle);

			triangle->Initialize();
			triangle->GetPointIds()->SetId(0,4*Np);
			triangle->GetPointIds()->SetId(1,Np);
			triangle->GetPointIds()->SetId(2,0);
			triangles ->InsertNextCell(triangle);

			triangle->Initialize();
			triangle->GetPointIds()->SetId(0,5*Np);
			triangle->GetPointIds()->SetId(1,2*Np);
			triangle->GetPointIds()->SetId(2,Np);
			triangles ->InsertNextCell(triangle);

	// ---------------------------------------- side 4 

			triangle->Initialize();
			triangle->GetPointIds()->SetId(0,Np-1);
			triangle->GetPointIds()->SetId(1,4*Np-1);
			triangle->GetPointIds()->SetId(2,5*Np-1);
			triangles ->InsertNextCell(triangle);

			triangle->Initialize();
			triangle->GetPointIds()->SetId(0,2*Np-1);
			triangle->GetPointIds()->SetId(1,5*Np-1);
			triangle->GetPointIds()->SetId(2,6*Np-1);
			triangles ->InsertNextCell(triangle);

			triangle->Initialize();
			triangle->GetPointIds()->SetId(0,5*Np-1);
			triangle->GetPointIds()->SetId(1,2*Np-1);
			triangle->GetPointIds()->SetId(2,Np-1);
			triangles ->InsertNextCell(triangle);

			triangle->Initialize();
			triangle->GetPointIds()->SetId(0,6*Np-1);
			triangle->GetPointIds()->SetId(1,3*Np-1);
			triangle->GetPointIds()->SetId(2,2*Np-1);
			triangles ->InsertNextCell(triangle);
	
		
	// -------------------- top

	for(int j=0;j<2;j++)
		for(int i = (j*Np)+1; i < (j+1)*Np; i++)
		{
			triangle->Initialize();
			triangle->GetPointIds()->SetId(0,i-1);
			triangle->GetPointIds()->SetId(1,i);
			triangle->GetPointIds()->SetId(2,Np+i-1);
			triangles ->InsertNextCell(triangle);
		}

	for(int j=1;j<3;j++)
		for(int i = j*Np; i < (j+1)*Np-1; i++)
		{
			triangle->Initialize();
			triangle->GetPointIds()->SetId(0,i);
			triangle->GetPointIds()->SetId(1,i+1);
			triangle->GetPointIds()->SetId(2,i-Np+1);
			triangles ->InsertNextCell(triangle);
		}	
	

	// build the polydata

	vtkPolyData *polyData = vtkPolyData::New();

	polyData->SetPoints(points);
	polyData->SetPolys(triangles );


	//evaluate the volume : expected outcome for the current points set 48

	vtkMassProperties *mass = vtkMassProperties::New();

	mass->SetInput(polyData);
	cout<<"volume "<<mass->GetVolume()<<endl;

	  // Setup actor and mapper
	  vtkSmartPointer&lt;vtkPolyDataMapper> mapper =   
vtkSmartPointer<vtkPolyDataMapper>::New();
	  mapper->SetInput(polyData);

	 
	  vtkSmartPointer<vtkActor> actor = vtkSmartPointer<vtkActor>::New();
	  actor->SetMapper(mapper);
	 
	  // Setup render window, renderer, and interactor
	  vtkSmartPointer<vtkRenderer> renderer = 
vtkSmartPointer<vtkRenderer>::New();
	  vtkSmartPointer<vtkRenderWindow> renderWindow = 
vtkSmartPointer<vtkRenderWindow>::New();
	  renderWindow->AddRenderer(renderer);
	  vtkSmartPointer<vtkRenderWindowInteractor> renderWindowInteractor = 
vtkSmartPointer<vtkRenderWindowInteractor>::New();
	  renderWindowInteractor->SetRenderWindow(renderWindow);
	  renderer->AddActor(actor);
	 
	  renderWindow->Render();
	  renderWindowInteractor->Start();

return 0;
}



--
View this message in context: http://vtk.1045678.n5.nabble.com/computing-the-volume-of-a-simple-multi-triangles-polydata-tp5721243.html
Sent from the VTK - Users mailing list archive at Nabble.com.

More information about the vtkusers mailing list