[vtkusers] vtkDelaunay3D: wrong convex hull

[Miguel Sotaquirá]

Hi David!

Here's the hardcoded version. Just to clarify: I just want to know if the
current 3-D vertex ("c" in the code) is contained in the convex hull of its
3-D neighbors ("pts" in the code). I've used vtkDelaunay3D and vtkDelaunay2D
currently without any success.
Miguel

#include <vtkSmartPointer.h>
#include <vtkPolyData.h>
#include <vtkPoints.h>
#include <vtkUnstructuredGrid.h>
#include <vtkDelaunay3D.h>
#include <vtkPlane.h>
#include <vtkDelaunay2D.h>
#include <vtkXMLPolyDataWriter.h>
#include <vtkXMLUnstructuredGridWriter.h>

int main (int argc, char* argv[])
{
// Current vertex coordinates (this data is not used. serves only as a
reference)
double c[3] = {-0.59842, 2.12762, 5.17183};
 // Current vertex normal
double N[3] = {-0.604189, 0.796814, -0.00654253};
 // Original 3-D neighbors coordinates
double p0[3] = {-0.53592, 2.21902, 5.10933};
double p1[3] = {-0.59842, 2.17465, 5.10933};
double p2[3] = {-0.53592, 2.17465, 5.17183};
double p3[3] = {-0.66092, 2.07811, 5.17183};
double p4[3] = {-0.59842, 2.17465, 5.23433};
double p5[3] = {-0.66092, 2.12762, 5.23433};
 // Insert 3-D neighbors in a polydata
vtkSmartPointer<vtkPoints> pts = vtkSmartPointer<vtkPoints>::New();
pts->InsertNextPoint( p0 );
pts->InsertNextPoint( p1 );
pts->InsertNextPoint( p2 );
pts->InsertNextPoint( p3 );
pts->InsertNextPoint( p4 );
pts->InsertNextPoint( p5 );
/*
// 3-D convex hull: NOT WORKING
// point cloud of neighbors
vtkSmartPointer<vtkUnstructuredGrid> points =
vtkSmartPointer<vtkUnstructuredGrid>::New();
points->SetPoints( pts );
points->Update();
 // convex hull
vtkSmartPointer<vtkDelaunay3D> delaunay3D =
vtkSmartPointer<vtkDelaunay3D>::New();
delaunay3D->SetInput( points );
delaunay3D->Update();
*/
 //
// A 2D triangulation (using vtkDelaunay2D) will be obtained by projecting
// 3-D points onto a plane
//
 // create a plane, seting its origin in (0,0,0) and its normal N
vtkSmartPointer<vtkPlane> plane = vtkSmartPointer<vtkPlane>::New();
plane->SetOrigin(0.0, 0.0, 0.0);
plane->SetNormal(N[0], N[1], N[2]);
 // project each point in the point cloud onto the plane
vtkSmartPointer<vtkPoints> ptsProj = vtkSmartPointer<vtkPoints>::New();
double origin[3] = {0.0, 0.0, 0.0};
double normal[3] = {N[0], N[1], N[2]};
double projected[3];
double p[3];
for ( int i = 0; i < pts->GetNumberOfPoints(); i++ ) {
pts->GetPoint(i,p);
plane->ProjectPoint( p, origin, normal, projected );
ptsProj->InsertNextPoint(projected);
}
 // apply delaunay2D on the projected points
vtkSmartPointer<vtkPolyData> ptsForTriang =
vtkSmartPointer<vtkPolyData>::New();
ptsForTriang->SetPoints(ptsProj);
 vtkSmartPointer<vtkDelaunay2D> delaunay2D =
vtkSmartPointer<vtkDelaunay2D>::New();
delaunay2D->SetInput(ptsForTriang);
delaunay2D->Update();

// export points (original and projected) and triangulation for viewing in
Paraview
vtkSmartPointer<vtkUnstructuredGrid> points2D =
vtkSmartPointer<vtkUnstructuredGrid>::New();
vtkSmartPointer<vtkUnstructuredGrid> points3D =
vtkSmartPointer<vtkUnstructuredGrid>::New();
points2D->SetPoints(ptsProj);
points3D->SetPoints(pts);
 vtkSmartPointer<vtkXMLUnstructuredGridWriter> writerpts2D =
vtkSmartPointer<vtkXMLUnstructuredGridWriter>::New();
writerpts2D->SetFileName("points2D.vtu");
writerpts2D->SetInput(points2D);
writerpts2D->Write();
 vtkSmartPointer<vtkXMLUnstructuredGridWriter> writerpts3D =
vtkSmartPointer<vtkXMLUnstructuredGridWriter>::New();
writerpts3D->SetFileName("points3D.vtu");
writerpts3D->SetInput(points3D);
writerpts3D->Write();
 vtkSmartPointer<vtkXMLPolyDataWriter> writer =
vtkSmartPointer<vtkXMLPolyDataWriter>::New();
writer->SetFileName("delaunay2D.vtp");
writer->SetInput(delaunay2D->GetOutput());
writer->Update();
 return 0;
}

2011/7/26 David Doria <daviddoria at gmail.com>

> 2011/7/26 Miguel Sotaquirá <msotaquira at gmail.com>:
> > Hi David,
> > You're right, vtkDelaunay3D is not the correct approach in this case. Now
> > I'm trying to obtain a 2D convex hull (vtkDelaunay2D) by projecting the
> 3-D
> > points onto a plane. Here's the code:
> > // create a plane, seting its origin in (0,0,0) and its normal N (N is
> the
> > normal to  the current vertex, the center of the neighborhood)
> > vtkSmartPointer<vtkPlane> plane = vtkSmartPointer<vtkPlane>::New();
> > plane->SetOrigin(0.0, 0.0, 0.0);
> > plane->SetNormal(N[0], N[1], N[2]);
> >  // project each point in the neighborhood onto the plane
> > vtkSmartPointer<vtkPoints> ptsProj = vtkSmartPointer<vtkPoints>::New();
> > double origin[3] = {0.0, 0.0, 0.0};
> > double normal[3] = {N[0], N[1], N[2]};
> > double projected[3];
> > for ( int i = 0; i < pts->GetNumberOfPoints(); i++ ) {
> >   pts->GetPoint(i,p);
> >   plane->ProjectPoint( p, origin, normal, projected );
> >   ptsProj->InsertNextPoint(projected);}
> > // store projected points in a polydata
> > vtkSmartPointer<vtkPolyData> ptsForTriang =
> > vtkSmartPointer<vtkPolyData>::New();
> > ptsForTriang->SetPoints(ptsProj);
> > // apply delaunay2D
> > vtkSmartPointer<vtkDelaunay2D> delaunay2D =
> > vtkSmartPointer<vtkDelaunay2D>::New();
> > delaunay2D->SetInput(ptsForTriang);
> > delaunay2D->Update();
> > When using this approach to my initial problem I get this:
> > http://tinypic.com/r/2mre8ma/7
> > again, some vertices were not included in the convex hull . When I check
> the
> > 2D convex hull info it says that has 6 points (which is correct) but only
> 2
> > cells (the ones that can be seen in the image)...
> > Is this the correct approach for computing the 2D convex hull (using
> > projections)? Am I missing something?
> > Thanks again,
> > Miguel
>
> Can you hardcode those points into a compilable demo of this?
>
> David
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.vtk.org/pipermail/vtkusers/attachments/20110727/39f832ba/attachment.htm>