[Insight-users] Problems about 3D construction with tetrahedral mesh from DICOM series

[smallping]

I want to read the DICOM series, and then create tetrahedral mesh for volume
modeling.

I do like this:
1) Read the 3D image data from a series of DICOM files
2) Use CurvatureFlowImageFilter to smooth the image
3) Run an image segmentation method 
4) Use the vtkContour filter to generate the surface mesh
5) Save the surface mesh data to a .vtk file 
6) Using TetGen, to generate the tetrahedral mesh from the surface mesh, and
save into .vtk file

But I have problems at step 6, when I try to use TetGen to generate the
tetrahedral mesh from the surface mesh, it failed.
Then I use ParaView to visualize my surface mesh, I find maybe there are
some problems with the surface mesh.
Because in the surface mesh, there are many holes, and a surface circle
around the object, I do not know what is this.
The surface mesh is shown as the following picture. 

http://itk-insight-users.2283740.n2.nabble.com/file/n6653520/Surface_mesh.jpg 

Can somebody know how to fix this ?
Or give me an example to do this ?
I need help from you.
Thanks

The following is my codes for step 2 to step 6 to generate the surface mesh
/*=========================================================================
#include "itkCommand.h"
#include "itkImage.h"
#include "itkVTKImageExport.h"
#include "itkVTKImageImport.h"
#include "itkConfidenceConnectedImageFilter.h"
#include "itkCastImageFilter.h"
#include "itkRGBPixel.h"
#include "itkImageFileReader.h"
#include "itkImageFileWriter.h"

#include "vtkImageImport.h"
#include "vtkImageExport.h"
#include "vtkRenderer.h"
#include "vtkRenderWindow.h"
#include "vtkRenderWindowInteractor.h"
#include "vtkActor.h"
#include "vtkPolyData.h"
#include "vtkPolyDataMapper.h"
#include "vtkContourFilter.h"
#include "vtkImageData.h"
#include "vtkDataSet.h"
#include "vtkProperty.h"
#include "vtkImagePlaneWidget.h"
#include "vtkCellPicker.h"
#include "vtkPolyDataWriter.h"

#include "itkCurvatureFlowImageFilter.h"
#include "vtkSmartPointer.h"
#include "vtkSTLWriter.h"


template <typename ITK_Exporter, typename VTK_Importer>
void ConnectPipelines(ITK_Exporter exporter, VTK_Importer* importer)
{
 
importer->SetUpdateInformationCallback(exporter->GetUpdateInformationCallback());
 
importer->SetPipelineModifiedCallback(exporter->GetPipelineModifiedCallback());
  importer->SetWholeExtentCallback(exporter->GetWholeExtentCallback());
  importer->SetSpacingCallback(exporter->GetSpacingCallback());
  importer->SetOriginCallback(exporter->GetOriginCallback());
  importer->SetScalarTypeCallback(exporter->GetScalarTypeCallback());
 
importer->SetNumberOfComponentsCallback(exporter->GetNumberOfComponentsCallback());
 
importer->SetPropagateUpdateExtentCallback(exporter->GetPropagateUpdateExtentCallback());
  importer->SetUpdateDataCallback(exporter->GetUpdateDataCallback());
  importer->SetDataExtentCallback(exporter->GetDataExtentCallback());
  importer->SetBufferPointerCallback(exporter->GetBufferPointerCallback());
  importer->SetCallbackUserData(exporter->GetCallbackUserData());
}

template <typename VTK_Exporter, typename ITK_Importer>
void ConnectPipelines(VTK_Exporter* exporter, ITK_Importer importer)
{
 
importer->SetUpdateInformationCallback(exporter->GetUpdateInformationCallback());
 
importer->SetPipelineModifiedCallback(exporter->GetPipelineModifiedCallback());
  importer->SetWholeExtentCallback(exporter->GetWholeExtentCallback());
  importer->SetSpacingCallback(exporter->GetSpacingCallback());
  importer->SetOriginCallback(exporter->GetOriginCallback());
  importer->SetScalarTypeCallback(exporter->GetScalarTypeCallback());
 
importer->SetNumberOfComponentsCallback(exporter->GetNumberOfComponentsCallback());
 
importer->SetPropagateUpdateExtentCallback(exporter->GetPropagateUpdateExtentCallback());
  importer->SetUpdateDataCallback(exporter->GetUpdateDataCallback());
  importer->SetDataExtentCallback(exporter->GetDataExtentCallback());
  importer->SetBufferPointerCallback(exporter->GetBufferPointerCallback());
  importer->SetCallbackUserData(exporter->GetCallbackUserData());
}


int main(int argc, char * argv [] )
{ 
  //DICOM-IMG.vtk is read from DICOM series
  argv[1]="/home/smallping/Codes/ITK/Test/MeshSurface/DICOM-IMG.vtk"; 
  
  try
  {
    typedef float        InputPixelType;
    typedef unsigned char MaskPixelType;

    const unsigned int Dimension = 3;
    typedef itk::Image< InputPixelType, Dimension > InputImageType;
    typedef itk::Image< MaskPixelType,  Dimension > MaskImageType ;

    typedef itk::ImageFileReader< InputImageType >  ReaderType    ;

    ReaderType::Pointer reader  = ReaderType::New();
    reader->SetFileName( argv[1] );
    reader->Update();

    ////////////////////////////////////////////////////////////// Smooth
    typedef itk::CurvatureFlowImageFilter< InputImageType, InputImageType >
CurvatureFlowImageFilterType; 
    CurvatureFlowImageFilterType::Pointer smoothing =
CurvatureFlowImageFilterType::New();

    smoothing->SetInput( reader->GetOutput() );
    smoothing->SetNumberOfIterations( 10 );
    smoothing->SetTimeStep( 0.0625 );

   
//////////////////////////////////////////////////////////////Segmentation
    typedef itk::ConfidenceConnectedImageFilter&lt;InputImageType,
MaskImageType&gt; ConnectedFilterType;
    ConnectedFilterType::Pointer filter = ConnectedFilterType::New();

    filter->SetInput( smoothing->GetOutput() );

    filter->SetNumberOfIterations(2);
    filter->SetReplaceValue(255);
    filter->SetMultiplier(2.5);
     
    InputImageType::Pointer inputImage = reader->GetOutput();
    InputImageType::SizeType  size  =
inputImage->GetBufferedRegion().GetSize();
    InputImageType::IndexType start =
inputImage->GetBufferedRegion().GetIndex();

    InputImageType::IndexType seed;
    seed[0] = start[0] + size[0] / 2;
    seed[1] = start[1] + size[1] / 2;
    seed[2] = start[2] + size[2] / 2;

    filter->SetSeed( seed );

    //////////////////////////////////////////////////////////////
    typedef itk::VTKImageExport< InputImageType > ExportFilter1Type;
    typedef itk::VTKImageExport< MaskImageType  > ExportFilter2Type;

    ExportFilter1Type::Pointer itkExporter1 = ExportFilter1Type::New();
    ExportFilter2Type::Pointer itkExporter2 = ExportFilter2Type::New();

    itkExporter1->SetInput( reader->GetOutput() );
    itkExporter2->SetInput( filter->GetOutput() );

    vtkImageImport* vtkImporter1 = vtkImageImport::New();  
    ConnectPipelines(itkExporter1, vtkImporter1);
    vtkImageImport* vtkImporter2 = vtkImageImport::New();  
    ConnectPipelines(itkExporter2, vtkImporter2);
    
    vtkImporter1->Update();
     
    //////////////////////////////////////////////////////////////Surface
    vtkContourFilter * contour = vtkContourFilter::New();
    contour->SetInput( vtkImporter2->GetOutput() );
    contour->SetValue(0, 128);

    //////////////////////////////////////////////////////////////
    vtkPolyDataWriter * vtkWriter = vtkPolyDataWriter::New();
   
vtkWriter->SetFileName("/home/smallping/Codes/ITK/Test/MeshSurface/DICOM-SURF.vtk");
//save as vtk
    vtkWriter->SetInput( contour->GetOutput() );
    vtkWriter->Write();

    vtkSmartPointer< vtkSTLWriter > stlWriter = vtkSmartPointer<
vtkSTLWriter >::New();  //save as stl
    stlWriter->SetInput( contour->GetOutput() );
    stlWriter->SetFileName(
"/home/smallping/Codes/ITK/Test/MeshSurface/DICOM-SURF.stl" );
    stlWriter->Update();
 
    //////////////////////////////////////////////////////////////
    vtkImporter1->Delete();
    vtkImporter2->Delete();
    contour->Delete();
  }
  catch( itk::ExceptionObject & e )
  {
    std::cerr << "Exception catched !! " << e << std::endl;
  }

  return 0;
}

--
View this message in context: http://itk-insight-users.2283740.n2.nabble.com/Problems-about-3D-construction-with-tetrahedral-mesh-from-DICOM-series-tp6653520p6653520.html
Sent from the ITK Insight Users mailing list archive at Nabble.com.