<html><body><div style="color:#000; background-color:#fff; 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 http://www.vtk.org/Wiki/VTK/Examples/Cxx/Modelling/MarchingCubes". 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 class="Apple-tab-span" style="white-space:pre">        </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 class="Apple-tab-span" style="white-space:pre">                </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 class="Apple-tab-span" style="white-space:pre">        </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 class="yahoo_quoted" 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 size="2" face="Arial"> On Wednesday, November 27, 2013 12:06 PM, Miro Drahos <mdrahos@robodoc.com> wrote:<br> </font> </div> <div class="y_msg_container">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 class="yqt5310183674" id="yqtfd95729"><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 www.kitware.com<br clear="none">><br clear="none">> Visit other Kitware open-source projects at <a shape="rect" href="http://www.kitware.com/opensource/opensource.html"
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 shape="rect" href="http://www.vtk.org/Wiki/VTK_FAQ" 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 shape="rect" href="http://www.vtk.org/mailman/listinfo/vtkusers" target="_blank">http://www.vtk.org/mailman/listinfo/vtkusers</a><br clear="none"><br clear="none"></div><br><br></div> </div> </div> </div> </div></body></html>