[vtkusers] Problem with vtkPolyDataToImageStencil
Roman Grothausmann
roman.grothausmann at helmholtz-berlin.de
Fri Jul 27 03:52:36 EDT 2012
Hello Matheus,
I've experienced the same problem, and thanks to David I've now a more
clear understanding when this happens. So far I've only checked if the
mesh is closed. I will add a check with tolerance >0.
With colleagues I've prepared a blender plug-in based on vtk filters for
voxelizing blender meshes. It will be published on Midas as soon as my
colleagues have given their OK. Let me know if You'd be interesting in
testing the script before.
Otherwise You could try vtkPartialVolumeModeller from Cory Quammen and
Taylor II R.M.: http://www.insight-journal.org/browse/publication/792
Corry is just upgrading the code a bit.
One other option I know of is to use ITK, which during my testing, was
not quicker than the VTK approach with vtkPolyDataToImageStencil. Below
is a little test program.
David: A more robust vtkPolyDataToImageStencil would be really great!
Let me know if there is any upgrade to test.
Happy coding
Roman
////////program to voxelize a vtk-polydata-mesh into a itk-voxel-image
////////supposingly quicker than pure vtk:
http://www.cmake.org/Wiki/VTK/Examples/PolyDataToImageData
//allocating image before voxelization
//#include <itkImage.h>
//#include <itkQuadEdgeMesh.h>
#include <itkVTKPolyDataReader.h> //can only read ASCII legacy VTK-files
containing only triangles!!!
#include <itkTriangleMeshToBinaryImageFilter.h>
#include <itkImageFileWriter.h>
#include "itkFilterWatcher2.h"
int main( int argc, char *argv[] ){
if( argc != 9 )
{
std::cerr << "Usage: " << argv[0];
std::cerr << " input_VTKpolydata-mesh";
std::cerr << " outputImage";
std::cerr << " dx dy dz ox oy oz";
std::cerr << std::endl;
return EXIT_FAILURE;
}
typedef unsigned char PixelType;
const unsigned int Dimension = 3;
//typedef itk::DefaultDynamicMeshTraits<double, 3, 3,double,double>
TriangleMeshTraits;
//typedef itk::Mesh<double,3, TriangleMeshTraits> TriangleMeshType;
//typedef itk::QuadEdgeMesh< float, Dimension > MeshType;
typedef itk::Mesh< double, 3 > MeshType;
typedef itk::Image< PixelType, Dimension > ImageType;
typedef itk::VTKPolyDataReader<MeshType> VTKmeshreaderType;
VTKmeshreaderType::Pointer meshReader = VTKmeshreaderType::New();
meshReader->SetFileName(argv[1]);
std::cout << "Reading: " << argv[1] << std::endl;
try
{
meshReader->Update();
}
catch( itk::ExceptionObject & excp )
{
std::cerr << "Exception: " << excp << std::endl;
return EXIT_FAILURE;
}
std::cout << "Read: " << argv[1] << std::endl;
// Set Size, Spacing and origin
ImageType::SizeType size;
size[ 0 ] = atoi(argv[3]);
size[ 1 ] = atoi(argv[4]);
size[ 2 ] = atoi(argv[5]);
ImageType::SpacingType spacing;
spacing[0] = 1; //100.0 / size[0];
spacing[1] = 1; //100.0 / size[1];
spacing[2] = 1; //100.0 / size[2];
ImageType::PointType origin;
origin[0]= atoi(argv[6]);
origin[1]= atoi(argv[7]);
origin[2]= atoi(argv[8]);
//allocate the output image
ImageType::Pointer output = ImageType::New();
output->SetRegions(size);
output->SetSpacing(spacing);
output->SetOrigin(origin);
//output->SetRegions(meshReader->GetOutput()->GetRequestedRegion());
output->Allocate();
std::cout << "Image allocated!" << std::endl;
//Set Inside/Outside voxel value
const PixelType empty = 0;
const PixelType fill = 255;
typedef itk::TriangleMeshToBinaryImageFilter<MeshType, ImageType>
MeshFilterType;
MeshFilterType::Pointer meshFilter = MeshFilterType::New();
meshFilter->SetInput(meshReader->GetOutput());
meshFilter->SetInfoImage(output); //
meshFilter->SetTolerance (1.0);
//meshFilter->SetSize (size);
//meshFilter->SetSpacing (spacing);
//meshFilter->SetOrigin(origin);
//meshFilter->SetIndex (index);
//meshFilter->SetUseObjectValue( true );
meshFilter->SetInsideValue(fill);
meshFilter->SetOutsideValue(empty);
FilterWatcher watcher(meshFilter, "filter");
meshFilter->Update();
// Write the image
typedef itk::ImageFileWriter< ImageType > WriterType;
WriterType::Pointer writer = WriterType::New();
writer->SetFileName(argv[2]);
writer->SetInput(meshFilter->GetOutput());
try
{
meshFilter->Update();
writer->Update();
}
catch( itk::ExceptionObject & excp )
{
std::cerr << excp << std::endl;
return EXIT_FAILURE;
}
return EXIT_SUCCESS;
}
On 27.07.2012 00:53, David Gobbi wrote:
> Hi Matheus,
>
> The vtkPolyDataToImageStencil filter is known to produce streaks like
> those under these conditions:
> 1) if the input data has free edges or is non-manifold
> 2) if the distance between points is small enough that it approaches
> the numerical tolerance of 32-bit floats
>
> You can check the first condition with vtkFeatureEdges, using the
> following code:
>
> vtkFeatureEdges *edges = vtkFeatureEdges::New();
> edges->SetInputConnection(yourdata->GetOutputPort());
> edges->FeatureEdgesOff();
> edges->NonManifoldEdgesOn();
> edges->BoundaryEdgesOn();
> edges->Update();
> cout << edges->GetOutput()->GetNumberOfCells() << endl;
>
> It should print "0" if your data is nice, well-defined closed surface.
>
> You can check the second condition with vtkCleanPolyData. If your
> polydata is changed by vtkCleanPolyData when the tolerance is set to
> 1e-7, then vtkPolyDataToImageStencil might give incorrect results due
> to roundoff error.
>
> I have plans to improve the vtkPolyDataToImageStencil code to make
> it more robust, but it will be several months (at least) before it rises
> to the top of my to-do list.
>
> - David
>
>
> On Thu, Jul 26, 2012 at 3:32 PM, matheus_viana <vianamp at gmail.com> wrote:
>> Hello guys.
>>
>> I've a polydata and I'd like to voxelize that in order to obtain an
>> ImageData. In order to do so, I'm using the following code:
>>
>> vtkImageData *blankImage = vtkImageData::New();
>> blankImage -> SetExtent(extent);
>> blankImage -> SetOrigin(0,0,0);
>> blankImage -> SetSpacing(1,1,1);
>> blankImage -> SetScalarTypeToUnsignedChar();
>> blankImage -> AllocateScalars();
>>
>> vtkPolyDataToImageStencil *pol2Stenc = vtkPolyDataToImageStencil::New();
>> pol2Stenc -> SetTolerance(0.5);
>> pol2Stenc -> SetInput(myPolyData);
>> pol2Stenc -> SetInformationInput(blankImage);
>> pol2Stenc -> Update();
>>
>> vtkImageStencil *stencil = vtkImageStencil::New();
>> stencil -> SetInput(blankImage);
>> stencil -> ReverseStencilOff();
>> stencil -> SetStencil(pol2Stenc->GetOutput());
>> stencil -> Update();
>> vtkImageData *image = stencil -> GetOutput();
>>
>> It does the job, but I've observed some artifacts in the final image, as can
>> be seen in the attached picture.
>>
>> Dos anyone have any idea about how to solve this problem? Is that a bug? Or
>> am I doing any stupid thing?
>>
>> Cheers
>> Matheus
>>
>> http://vtk.1045678.n5.nabble.com/file/n5714947/Viana-vtkPolyDataToImageStencil.jpg
>>
>>
>>
>> --
>> View this message in context: http://vtk.1045678.n5.nabble.com/Problem-with-vtkPolyDataToImageStencil-tp5714947.html
>> Sent from the VTK - Users mailing list archive at Nabble.com.
>> _______________________________________________
>> Powered by www.kitware.com
>>
>> Visit other Kitware open-source projects at http://www.kitware.com/opensource/opensource.html
>>
>> Please keep messages on-topic and check the VTK FAQ at: http://www.vtk.org/Wiki/VTK_FAQ
>>
>> Follow this link to subscribe/unsubscribe:
>> http://www.vtk.org/mailman/listinfo/vtkusers
> _______________________________________________
> Powered by www.kitware.com
>
> Visit other Kitware open-source projects at http://www.kitware.com/opensource/opensource.html
>
> Please keep messages on-topic and check the VTK FAQ at: http://www.vtk.org/Wiki/VTK_FAQ
>
> Follow this link to subscribe/unsubscribe:
> http://www.vtk.org/mailman/listinfo/vtkusers
>
--
Roman Grothausmann
Helmholtz-Zentrum Berlin für Materialien und Energie GmbH
Bereich Funktionale Materialien
Institut für angewandte Materialforschung
Hahn-Meitner-Platz 1
D-14109 Berlin
(früher Hahn-Meitner-Institut und BESSY)
Tel.: +49-(0)30-8062-42816
Fax.: +49-(0)30-8062-43059
Vorsitzender des Aufsichtsrats: Prof. Dr. Dr. h.c. mult. Joachim Treusch
Stellvertretende Vorsitzende: Dr. Beatrix Vierkorn-Rudolph
Geschäftsführer: Prof. Dr. Anke Rita Kaysser-Pyzalla, Dr. Ulrich Breuer
Sitz der Gesellschaft: Berlin
Handelsregister: AG Charlottenburg, 89 HRB 5583
More information about the vtkusers
mailing list