[vtkusers] Solid Voxelization with VTK
Berti Krüger
berti_krueger at hotmail.com
Wed Jun 13 21:21:07 EDT 2018
Hi Elvis,
thanks for jumping in. As you can see in David's response you were right or at
least close with your assumptions.
Regards
Berti
Am 11.06.2018 um 08:50 schrieb Elvis Stansvik:
Den mån 11 juni 2018 08:33Berti Krüger <berti_krueger at hotmail.com<mailto:berti_krueger at hotmail.com>> skrev:
Hi David,
thank you very much for all your help and for all the effort.
I am sorry, but i am an absolute beginner with VTK and for me it is really hard to find informations / documentation / explanations beyond the doxygen class documentation about the advanced visualization/manipulation operations in VTK (the VTK User Manual sadly only covers the very basics).
So thanks again for all your clarifications and all your help so far.
I replaced the visualization part of my pipeline with your code which makes it so much clearer, shorter and better in terms of the visualization. And it is also easy enough so that even i get it :-)
Unfortunately when i extend the simple cone source voxelization example i attached the last time by using general STL meshes via the vtkSTLReader, i get missing voxels in the visualization. To be sure that the problem is indeed in the visualization part and not in the voxelization part i used binvox, a mature voxelization program (http://www.patrickmin.com/binvox/ ) to do the voxelization and the code you suggested for the visualization of the vtkImageData via VTK.
Using this setup the Utah Teapot looks like this (the handle and other thin parts of the mesh are missing in the visualization):
[Teapot missing voxels]
For surface-only-voxelized meshes there even more voxels missing (the example here is the Stanford Bunny which has only been surface voxelized):
[bunny surface voxelization visualization]
What i don't get is that the visualization code is as simple as it gets:
vtkNew<vtkThreshold> selector;
selector->SetInputData(voxelImage);
selector->ThresholdByUpper(1.0);
selector->Update();
vtkNew<vtkDataSetMapper> mapper;
mapper->SetInputConnection(selector->GetOutputPort());
mapper->ScalarVisibilityOff();
mapper->StaticOn();
When i understand it correctly, then we simply visualize every cell (voxel) whose threshold is equal or greater than 1.0 which means we visualize (filter) all voxels which are set (in my vtkImageData voxels which are set are 1, voxels which are not set 0) and the vtkDataSetMapper maps these to graphic primitives and this does the rest of the visualization. Not much to do wrong here and that basically should do the trick but it doesn't.
Any idea why this does not work correctly all the time and has problems with thin parts?
I'm just jumping in from the side here, and I'm by no means a VTK expert. But I think what you're seeing is simply due to how VTK does volume rendering. It interpolates the given data points, which means you need at least two data points to get a voxel.
I'm guessing that at the resolution you're doing the voxelization here, your teapot handle is simply too thin to result in even a single voxel.
Think of your data points as a grid of fenceposts. You'll get a voxel in the middle of each 2x2x2 cluster of posts.
Someone correct me if I'm wrong here.
HTH,
Elvis
I therefore tried another approach by using the vtkGlyph3DMapper based on the example code from
<https://www.vtk.org/Wiki/VTK/Examples/Cxx/Visualization/Glyph3DMapper>
https://www.vtk.org/Wiki/VTK/Examples/Cxx/Visualization/Glyph3DMapper :
vtkNew<vtkPoints> points;
for (int x = 0; x < dims[0]; x++)
{
for (int y = 0; y < dims[1]; y++)
{
for (int z = 0; z < dims[2]; z++)
{
char* pixel = static_cast<char*>(voxelImage->GetScalarPointer(x, y, z));
if(pixel[0] == 1) // if the voxel is set (1) in the voxel image
points->InsertNextPoint(x, y, z);
}
}
}
vtkNew<vtkPolyData> polydata;
polydata->SetPoints(points);
vtkNew<vtkPolyData> glyph;
vtkNew<vtkCubeSource> cubeSource;
vtkNew<vtkGlyph3DMapper> glyph3Dmapper;
glyph3Dmapper->SetSourceConnection(cubeSource->GetOutputPort());
glyph3Dmapper->SetInputData(polydata);
glyph3Dmapper->Update();
actor->SetMapper(glyph3Dmapper);
The visualization with this approach is a bit slower (probably because it does no backface culling of the backfaces of the cube polygons) but it seems to work correctly even with surface only voxelized meshes:
[bunny ok][teapot]
I attached a simple working vtk demo with both visualizations and integrated voxelized meshes (no external data files or other dependencies) as a simple cmake/vtk compileable project (PolyDataToImageDataVisualizationVoxel). Maybe it might help other vtk users in the future.
Are there any other possibilities (maybe faster) to do the visualization with vtk ?
Now having a visualization that seems to work, i tried this visualization with the suggested voxelization code from the vtk examples here https://lorensen.github.io/VTKExamples/site/Cxx/PolyData/PolyDataToImageData/
and while the voxelization is lighting-fast and some example meshes did get correctly voxelized, i also had a lot of meshes whose voxelizations had misclassified voxels (e.g. voxels are set where voxels should not have been set and the opposite (holes)) like with the Utah Teapot here:
[teapot errors][teapot errors 2][teapot errors3]
The problem is that for my project i have to find paths through the voxel image. So miscIassified voxels would lead to either wrong paths found (paths which do not really exist like in the first image) or maybe to a situation where no path can be found because of holes.
I attached a compileable cmake/vtk example (PolyDataToImageDataVisualizationSTL) together with some small example stl-mesh-files (teapot etc.) using the example code from the vtk example repository for showing the problem with the given approach there.
Is it possible by tweaking the vtk PolyDataToImageData example code to always get a correct voxelization or is another approach needed?
If you or anybody else has ideas / suggestions i would really be glad.
Thank for reading if you managed to get this far and have not fallen asleep and thank you very much in advance.
Regards
Berti
Am 06.06.2018 um 01:49 schrieb David Gobbi:
Hi Berti,
The correct pad filter to use for this is "vtkImageConstantPad", not "WrapPad", but padding really isn't necessary. Neither is the use of cell scalars, it's more appropriate to use point scalars when working with image data since each voxel is represented as a point scalar.
The threshold should be set like this, for a threshold halfway between your "inval" and "outval":
selector->ThresholdByUpper(0.5*(startLabel + endLabel));
The lookup table wasn't correctly set. Calling SetRange() does not set the number of entries in the lookup table (by default, the lookup table has 256 entries).
lut->SetTableRange(0, 1);
lut->SetNumberOfColors(2);
Also, your color for "1" was transparent, whereas it has to be opaque in order to properly see a scalar-colored cone.
I stripped the visualization part of your pipeline down to the following. No padding is done, the thresholding uses point scalars, no shifting, and no lookup table (instead, I call ScalarVisibilityOff() and set the color via the actor's property object). I replaced vtkPolyDataMapper with vtkDataSetMapper, which automatically skins the data.
unsigned char startLabel = outval;
unsigned char endLabel = inval;
vtkNew<vtkThreshold> selector;
// the call to SetInputArrayToProcess() isn't necessary, the default is fine
//selector->SetInputArrayToProcess(0, 0, 0,
// vtkDataObject::FIELD_ASSOCIATION_POINTS,
// vtkDataSetAttributes::SCALARS);
selector->SetInputData(voxelImage);
selector->ThresholdByUpper(0.5*(startLabel + endLabel));
selector->Update();
vtkNew<vtkDataSetMapper> mapper;
mapper->SetInputConnection(selector->GetOutputPort());
mapper->ScalarVisibilityOff();
vtkNew<vtkActor> actor;
actor->SetMapper(mapper);
actor->GetProperty()->SetColor(0.0, 1.0, 0.0);
An image of the resulting cone is attached.
- David
On Tue, Jun 5, 2018 at 1:26 AM, Berti Krüger <berti_krueger at hotmail.com<mailto:berti_krueger at hotmail.com>> wrote:
Hi David,
thank you very much for your help.
I tried the example from your link. Unfortunately the result is same as using the vtkVoxelModeller:
[Screenshot result voxelization]
I get only a voxelization of the surface. The interior of the voxelized cone is still hollow and not filled with voxels (at least this is implied by the visualization).
For the visualization i use Bill Lorensen's example code from here:
https://lorensen.github.io/VTKExamples/site/Cxx/Medical/GenerateCubesFromLabels/
I have attached this small example which compiles with VTK 8.1 to show the issue. It basicly only consist of the example code from your link (https://lorensen.github.io/VTKExamples/site/Cxx/PolyData/PolyDataToImageData/) and the visualization code from the Bill Lorensen's example from the link above. I have not changed much.
Any idea what is going wrong or what can i do get also voxelization of the interior (solid voxelization)?
Thank you very much in advance.
Regards
Berti
________________________________
Von: David Gobbi <david.gobbi at gmail.com<mailto:david.gobbi at gmail.com>>
Gesendet: Montag, 28. Mai 2018 23:13
An: Berti Krüger
Cc: vtkusers at vtk.org<mailto:vtkusers at vtk.org>
Betreff: Re: [vtkusers] Solid Voxelization with VTK
Hi Berti,
If its a triangulated surface that you want to fill with voxels, use vtkPolyDataToImageStencil:
https://lorensen.github.io/VTKExamples/site/Cxx/PolyData/PolyDataToImageData/
If its a mesh of 3D elements and you want to sample the elements to create voxels, try vtkResampleToImage.
- David
On Mon, May 28, 2018 at 4:41 PM, Berti Krüger <berti_krueger at hotmail.com<mailto:berti_krueger at hotmail.com>> wrote:
Hello Everyone.
For my Project i have to voxelize a 3D-Triangle-Mesh. I already found the
vtkVoxelModeller which, while somewhat slow (only around 200 Triangles per
second), works, but only voxelizes the outer shell where the mesh boundary
triangles are. The inner part of the mesh stays hollow.
Is there some way to get solid voxelization of 3D-Triangle-Meshes out of the
box with VTK ?
Thank you very much in advance,
Berti
_______________________________________________
Powered by www.kitware.com<http://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
Search the list archives at: http://markmail.org/search/?q=vtkusers
Follow this link to subscribe/unsubscribe:
https://public.kitware.com/mailman/listinfo/vtkusers
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://public.kitware.com/pipermail/vtkusers/attachments/20180614/a5a3298d/attachment.html>
More information about the vtkusers
mailing list