<div dir="ltr"><div>Would it be acceptable to step out of VTK for meshing?</div><div>You could try to do something like the following:</div><div>1. Import DICOM data</div><div>2. Select region of interest</div><div>3. Export region of interest to STL</div>
<div>4. Use dedicated mesh generator like netgen</div><div> <a href="http://sourceforge.net/apps/mediawiki/netgen-mesher/index.php?title=Main_Page">http://sourceforge.net/apps/mediawiki/netgen-mesher/index.php?title=Main_Page</a></div>
<div> to mesh the geometry</div><div>5. Import mesh into VTK</div><div> </div><div>Export/Import can ofcourse also be replaced by</div><div>direct calls from C++ into the netgen library.</div><div>I am not sure if direct VTK support is available in</div>
<div>netgen.</div><div> </div><div> </div><div>Regards,</div><div> </div><div>Marco</div><div> </div></div><div class="gmail_extra"><br><br><div class="gmail_quote">On Wed, Nov 27, 2013 at 5:20 AM, Harold <span dir="ltr"><<a href="mailto:boogiedoll@yahoo.com" target="_blank">boogiedoll@yahoo.com</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div><div style="font-family:HelveticaNeue,Helvetica Neue,Helvetica,Arial,Lucida Grande,sans-serif;font-size:10pt"><div>
Hi Miro,</div><div><br></div><div>Thank you for the guide. I'm new to vtk so while messing around, I managed to apply the followings to a set of dicom images (in a very crude way) and save it as a .vtk file.</div><div>
1. Marching Cubes - get mesh</div><div>2. Decimation - reduce the number of triangles</div><div>3. Smoothing - smooth the mesh</div><div>4. Cropping - select the largest region using connectivity filter</div><div><br></div>
<div>I guess my next step is to find out how to do pre-processing and refining of the defected mesh, and also making statistical models (if there are enough training images...).</div><div><br></div><div>By the way, here's my crude code based on the the MarchingCubes Example at <a href="http://www.vtk.org/Wiki/VTK/Examples/Cxx/Modelling/MarchingCubes" target="_blank">http://www.vtk.org/Wiki/VTK/Examples/Cxx/Modelling/MarchingCubes</a>". Hopefully
it can be a small help to other novice like me out there.</div><div><br></div><div>/////////////////////////////////</div><div>int main(int argc, char *argv[])</div><div>{</div><div> vtkSmartPointer<vtkImageData> volume =</div>
<div> vtkSmartPointer<vtkImageData>::New();</div><div> double isoValue;</div><div><br></div><div> if (argc < 3) </div><div> {</div><div><span style="white-space:pre-wrap"> </span>// Create a sphere if there's no dicom input</div>
<div> vtkSmartPointer<vtkSphereSource> sphereSource = </div><div> vtkSmartPointer<vtkSphereSource>::New();</div><div> sphereSource->SetPhiResolution(20);</div><div> sphereSource->SetThetaResolution(20);</div>
<div> sphereSource->Update();</div><div> </div><div> double
bounds[6];</div><div> sphereSource->GetOutput()->GetBounds(bounds);</div><div> for (unsigned int i = 0; i < 6; i += 2)</div><div> {</div><div> double range = bounds[i+1] - bounds[i];</div><div> bounds[i] = bounds[i] - .1 * range;</div>
<div> bounds[i+1] = bounds[i+1] + .1 * range;</div><div> }</div><div> vtkSmartPointer<vtkVoxelModeller> voxelModeller = </div><div> vtkSmartPointer<vtkVoxelModeller>::New();</div><div> voxelModeller->SetSampleDimensions(50,50,50);</div>
<div> voxelModeller->SetModelBounds(bounds);</div><div> voxelModeller->SetScalarTypeToFloat();</div><div> voxelModeller->SetMaximumDistance(.1);</div><div> </div><div>
voxelModeller->SetInputConnection(sphereSource->GetOutputPort());</div><div> voxelModeller->Update();</div><div> isoValue = 0.5;</div><div> volume->DeepCopy(voxelModeller->GetOutput());</div><div>
}</div><div> else</div><div> {</div><div> // Read dicom image when dicom image folder path is specified</div><div> vtkSmartPointer<vtkDICOMImageReader> reader =</div><div> vtkSmartPointer<vtkDICOMImageReader>::New();</div>
<div> reader->SetDirectoryName(argv[1]);</div><div> reader->Update();</div><div> volume->DeepCopy(reader->GetOutput());</div><div> isoValue = atof(argv[2]);</div><div> }</div><div><br></div><div>
// Use MarchingCubes to generate iso-surface</div><div> vtkSmartPointer<vtkMarchingCubes> surface
= </div><div> vtkSmartPointer<vtkMarchingCubes>::New();</div><div><br></div><div> #if VTK_MAJOR_VERSION <= 5</div><div> surface->SetInput(volume);</div><div> #else</div><div> surface->SetInputData(volume);</div>
<div> #endif</div><div> surface->ComputeNormalsOn();</div><div> surface->ComputeScalarsOn();</div><div> surface->SetValue(0, isoValue);</div><div> std::cout << "Marching Cube finished..." << std::endl;</div>
<div><br></div><div> // Create polydata from iso-surface</div><div> vtkSmartPointer<vtkPolyData> marched =</div><div> vtkSmartPointer<vtkPolyData>::New();</div><div> surface->Update();</div><div> marched->DeepCopy(surface->GetOutput());</div>
<div> std::cout<<"Number of points: " << marched->GetNumberOfPoints() <<
std::endl;</div><div> </div><div> // Decimation to reduce the number of triangles</div><div> vtkSmartPointer<vtkDecimatePro> decimator = </div><div> vtkDecimatePro::New();</div><div> decimator->SetInputData(marched);</div>
<div> decimator->SetTargetReduction(0.5);</div><div> decimator->SetPreserveTopology(1);</div><div> decimator->Update();</div><div> std::cout << "Decimation finished...." << std::endl;</div>
<div><span style="white-space:pre-wrap"> </span></div><div> // Smoothing</div><div> vtkSmartPointer<vtkSmoothPolyDataFilter> smoother = vtkSmartPointer<vtkSmoothPolyDataFilter>::New();</div><div> smoother->SetInputData(decimator->GetOutput());</div>
<div> smoother->SetNumberOfIterations(5);</div><div> smoother->SetFeatureAngle(60);</div><div>
smoother->SetRelaxationFactor(0.05);</div><div> smoother->FeatureEdgeSmoothingOff();</div><div> std::cout << "Smoothing mesh finished...." << std::endl;</div><div><br></div><div> // Select the largest region</div>
<div> vtkSmartPointer<vtkPolyDataConnectivityFilter> connectivityFilter =</div><div><span style="white-space:pre-wrap"> </span>vtkSmartPointer<vtkPolyDataConnectivityFilter>::New();</div><div> connectivityFilter->SetInputConnection(smoother->GetOutputPort());</div>
<div> connectivityFilter->ScalarConnectivityOff();</div><div> connectivityFilter->SetExtractionModeToLargestRegion();</div><div> connectivityFilter->Update();</div><div> </div><div> // Create final polygon mesh</div>
<div> vtkSmartPointer<vtkPolyData> mesh =</div><div> vtkSmartPointer<vtkPolyData>::New();</div><div>
mesh->ShallowCopy(connectivityFilter->GetOutput());</div><div> std::cout<<"Number of points in the final polygon: " << mesh->GetNumberOfPoints() << std::endl;</div><div> </div><div>
// Visualization</div><div> vtkSmartPointer<vtkRenderer> renderer = </div><div> vtkSmartPointer<vtkRenderer>::New();</div><div> renderer->SetBackground(.1, .2, .3);</div><div><br></div><div> vtkSmartPointer<vtkRenderWindow> renderWindow = </div>
<div> vtkSmartPointer<vtkRenderWindow>::New();</div><div> renderWindow->AddRenderer(renderer);</div><div> vtkSmartPointer<vtkRenderWindowInteractor> interactor = </div><div> vtkSmartPointer<vtkRenderWindowInteractor>::New();</div>
<div> interactor->SetRenderWindow(renderWindow);</div><div><br></div><div> vtkSmartPointer<vtkPolyDataMapper> mapper
= </div><div> vtkSmartPointer<vtkPolyDataMapper>::New();</div><div> mapper->SetInputData(mesh);</div><div> mapper->ScalarVisibilityOff();</div><div><br></div><div> vtkSmartPointer<vtkActor> actor = </div>
<div> vtkSmartPointer<vtkActor>::New();</div><div> actor->SetMapper(mapper);</div><div><br></div><div> renderer->AddActor(actor);</div><div><br></div><div> renderWindow->Render();</div><div> interactor->Start();</div>
<div><br></div><div> // Save the mesh as a .vtk file</div><div> vtkSmartPointer<vtkPolyDataWriter> writer = vtkPolyDataWriter::New();</div><div> writer->SetFileName("mesh.vtk");</div><div> writer->SetInputData(mesh);</div>
<div> writer->SetFileTypeToASCII();</div><div> writer->Write();</div><div> </div><div> return
EXIT_SUCCESS;</div><div>}</div><div><br></div><div><br></div><div>Best Regards,</div><div><span></span></div><div>Harold</div><div><div class="h5"><div style="display:block"> <br> <br> <div style="font-family:HelveticaNeue,"Helvetica Neue",Helvetica,Arial,"Lucida Grande",sans-serif;font-size:10pt">
<div style="font-family:HelveticaNeue,"Helvetica Neue",Helvetica,Arial,"Lucida Grande",sans-serif;font-size:12pt"> <div dir="ltr"> <font face="Arial"> On Wednesday, November 27, 2013 12:06 PM, Miro Drahos <<a href="mailto:mdrahos@robodoc.com" target="_blank">mdrahos@robodoc.com</a>> wrote:<br>
</font> </div> <div>This is basically a segmentation problem, which in general is a <br clear="none">difficult problem.<br clear="none">One such pipeline could be as follows:<br clear="none">1.) image preprocessing (cropping, smoothing, etc).<br clear="none">
2.) thresholding to get a filled binary mask of your object (vertebra)<br clear="none">3.) running Marching
Cubes to get a mesh.<br clear="none"><br clear="none">This is a very crude and basic approach, but something to get you <br clear="none">started. Note that marching cubes might give you a mesh with topological <br clear="none">
defects (holes, handles), from the noise in your data.<br clear="none"><br clear="none">There are a variety of other sophisticated algorithms, e.g. using <br clear="none">statistical methods if you have some already segmented data (like <br clear="none">
creating a mean shape from the previously segmented data, and using PCA <br clear="none">to do shape analysis, or some machine learning approach).<br clear="none"><br clear="none">Another approach would be to use deformable models, i.e. define forces <br clear="none">
that would (iteratively) drive an initial mesh to fit the image. The <br clear="none">forces would be based on image intensity, gradient and/or laplacian, and <br clear="none">internal forces to keep the mesh regular
and relatively smooth.<br clear="none"><br clear="none">There is a lot that has been done in this domain, so look at some <br clear="none">articles on segmentation.<br clear="none"><br clear="none">HTH,<br clear="none">
Miro<br clear="none"><br clear="none"><div><br clear="none">On 11/25/2013 06:20 AM, Boogiedoll wrote:<br clear="none">> Hi,<br clear="none">> I would like to ask the experts about common way in vtk to create a spine bone (e.g L1) mesh from a CT DICOM stack, including the pre-processing.<br clear="none">
><br clear="none">> Thanks & Regards,<br clear="none">> Harold<br clear="none">> _______________________________________________<br clear="none">> Powered by <a href="http://www.kitware.com" target="_blank">www.kitware.com</a><br clear="none">
><br clear="none">> Visit other Kitware open-source projects at <a href="http://www.kitware.com/opensource/opensource.html" shape="rect" target="_blank">http://www.kitware.com/opensource/opensource.html</a><br clear="none">
><br clear="none">> Please keep messages on-topic and check the VTK FAQ at: <a href="http://www.vtk.org/Wiki/VTK_FAQ" shape="rect" target="_blank">http://www.vtk.org/Wiki/VTK_FAQ</a><br clear="none">><br clear="none">
> Follow this link to subscribe/unsubscribe:<br clear="none">> <a href="http://www.vtk.org/mailman/listinfo/vtkusers" shape="rect" target="_blank">http://www.vtk.org/mailman/listinfo/vtkusers</a><br clear="none"><br clear="none">
</div><br><br></div> </div> </div> </div> </div></div></div></div><br>_______________________________________________<br>
Powered by <a href="http://www.kitware.com" target="_blank">www.kitware.com</a><br>
<br>
Visit other Kitware open-source projects at <a href="http://www.kitware.com/opensource/opensource.html" target="_blank">http://www.kitware.com/opensource/opensource.html</a><br>
<br>
Please keep messages on-topic and check the VTK FAQ at: <a href="http://www.vtk.org/Wiki/VTK_FAQ" target="_blank">http://www.vtk.org/Wiki/VTK_FAQ</a><br>
<br>
Follow this link to subscribe/unsubscribe:<br>
<a href="http://www.vtk.org/mailman/listinfo/vtkusers" target="_blank">http://www.vtk.org/mailman/listinfo/vtkusers</a><br>
<br></blockquote></div><br></div>