[ITK-users] Simplex Mesh Deformable Model doesn't change

[bsouza]

Dear itk-developers,

I'm trying to use the itkDeformableSimplexMesh3D filter to make a local
adjustment of a simplex mesh. The mesh that I want to deform using this
filter is already initialized close to the target structure in the image.
However, after playing with all parameter values of the filter, it seems
that there is not significant differences between the input and output
meshes.

I'm reading the mesh from a file (stl or vtk extension) and I want to deform
it using a brain MRI (nifiti format) image as a reference. After the filter
application, the result seems to not approximate the mesh to the local
structure in the image (only occurs a small smoothness of the model). I've
tried to change the parameters (alpha, beta, gamma, rigidity, number of
iterations) values to check if some change would occur, but the result
always seems the same. I would like that deformation to be able to move the
mesh to closer to the structure in the image. My code is below. Is there
something wrong with my code? Why it seems that the external force has no
influence in the mesh deform result?

The image and the mesh I'm using in this test is available in the link below
(I don't know how much time until the link expires):

http://www73.zippyshare.com/v/sqQ6cgjN/file.html

Here is my code:

=================================================================
enum PARAMS {
    PROGRAM_NAME = 0,
    IMAGE_FILENAME,
    INPUT_MESH_FILENAME,
    OUTPUT_MESH_FILENAME,
    NUM_ARGS
};


using namespace std;

int main( int argc, char *argv[] )
{

    if( argc < NUM_ARGS )
    {
        std::cerr << "\n\n";
        std::cerr << "Usage: " << argv[PROGRAM_NAME] << " imageFileName
inputMeshFileName ouputMeshFileName";
        std::cerr << "\n\n";
        return 1;
    }

    itk::STLMeshIOFactory::RegisterOneFactory();

  const int dimension = 3;

  typedef float PixelType;
  typedef itk::Image<PixelType, dimension> ImageType;

  typedef itk::DefaultDynamicMeshTraits<PixelType, dimension>      
TriangleMeshTraits;
  typedef itk::DefaultDynamicMeshTraits<PixelType, dimension>      
SimplexMeshTraits;
  typedef itk::Mesh<PixelType, dimension, TriangleMeshTraits>      
TriangleMeshType;
  typedef itk::SimplexMesh<PixelType, dimension, SimplexMeshTraits>
SimplexMeshType;

  typedef itk::TriangleMeshToSimplexMeshFilter<TriangleMeshType,
SimplexMeshType> TriangleToSimplexFilterType;
  typedef itk::SimplexMeshToTriangleMeshFilter<SimplexMeshType,
TriangleMeshType> SimplexToTriangleFilterType;

  typedef itk::DeformableSimplexMesh3DFilter<SimplexMeshType,
SimplexMeshType> DeformFilterType;
  typedef DeformFilterType::GradientImageType                                 
GradientType;

  typedef itk::GradientAnisotropicDiffusionImageFilter <ImageType,
ImageType> DiffusionFilterType;

  typedef itk::GradientAnisotropicDiffusionImageFilter < ImageType,
ImageType >  GradientAnisotropicImageType;
  typedef itk::GradientMagnitudeRecursiveGaussianImageFilter < ImageType,
ImageType >  GradientMagnitudeType;
  typedef itk::SigmoidImageFilter< ImageType, ImageType > SigmoidImageType;
  typedef itk::GradientRecursiveGaussianImageFilter<ImageType,GradientType>
GradientFilterType;

  typedef itk::ImageFileReader<ImageType>       ImageReaderType;
  typedef itk::ImageFileWriter<ImageType>       ImageWriterType;
  typedef itk::MeshFileReader<TriangleMeshType> MeshReaderType;
  typedef itk::MeshFileWriter<TriangleMeshType>  MeshWriterType;


  ImageReaderType::Pointer imageReader = ImageReaderType::New();
  imageReader->SetFileName(argv[IMAGE_FILENAME]);
  imageReader->Update();

  MeshReaderType::Pointer meshReader = MeshReaderType::New();
  meshReader->SetFileName(argv[INPUT_MESH_FILENAME]);
  meshReader->Update();

  std::cout << " starting to Filter Image" << std::endl;
  GradientAnisotropicImageType::Pointer gradientanisotropicfilter =
GradientAnisotropicImageType::New();
  gradientanisotropicfilter->SetInput(imageReader->GetOutput());
  gradientanisotropicfilter->SetNumberOfIterations(5);
  gradientanisotropicfilter->SetTimeStep(0.0625);
  gradientanisotropicfilter->SetConductanceParameter(3);
  gradientanisotropicfilter->Update();
  std::cout << "GradientAnisotropicDiffusion is DONE!" << std::endl;

  GradientMagnitudeType::Pointer gradientmagnitudefilter =
GradientMagnitudeType::New();
  gradientmagnitudefilter->SetInput( gradientanisotropicfilter->GetOutput()
);
  gradientmagnitudefilter->SetSigma(1.0);
  gradientmagnitudefilter->Update();
  std::cout << "GradientMagnitude is DONE!" << std::endl;

  SigmoidImageType::Pointer sigmoidimagefilter = SigmoidImageType::New();
  sigmoidimagefilter->SetInput( gradientmagnitudefilter->GetOutput());
  sigmoidimagefilter->SetOutputMinimum(0);
  sigmoidimagefilter->SetOutputMaximum(1);
  sigmoidimagefilter->SetAlpha(10);
  sigmoidimagefilter->SetBeta(100);
  sigmoidimagefilter->Update();
  std::cout << "Sigmoid is DONE!" << std::endl;

  GradientFilterType::Pointer gradientFilter = GradientFilterType::New();
  gradientFilter->SetInput( sigmoidimagefilter->GetOutput() );
  gradientFilter->SetSigma(1.0);
  gradientFilter->Update();
  std::cout << "GradientMagnitude is DONE!" << std::endl;

  std::cout << "Gradient filter ok\n";

  TriangleToSimplexFilterType::Pointer triangleToSimplexFilter =
TriangleToSimplexFilterType::New();
  triangleToSimplexFilter->SetInput(meshReader->GetOutput());
  triangleToSimplexFilter->Update();

  GradientType::Pointer gradientImage    = gradientFilter->GetOutput();
  SimplexMeshType::Pointer simplexMesh   =
triangleToSimplexFilter->GetOutput();
  DeformFilterType::Pointer deformFilter = DeformFilterType::New();

  const unsigned int numberOfCycles = 100;

  for (unsigned int i = 0; i < numberOfCycles; i++)
    {
    // must disconnect the pipeline
    simplexMesh->DisconnectPipeline();
    deformFilter->SetInput( simplexMesh );
    deformFilter->SetGradient( gradientImage );
    deformFilter->SetAlpha(0.1);
    deformFilter->SetBeta(-0.1);
    deformFilter->SetIterations(5);
    deformFilter->SetRigidity(1);
    deformFilter->Update();
    }
  SimplexMeshType::Pointer deformResult =  deformFilter->GetOutput();


   SimplexToTriangleFilterType::Pointer simplexToTriangleFilter =
SimplexToTriangleFilterType::New();
   simplexToTriangleFilter->SetInput(deformResult);
   simplexToTriangleFilter->Update();
   TriangleMeshType::Pointer conversionResult =
simplexToTriangleFilter->GetOutput();

  conversionResult->DisconnectPipeline();

  MeshWriterType::Pointer meshWriter = MeshWriterType::New();
  meshWriter->SetFileName(argv[OUTPUT_MESH_FILENAME]);
  meshWriter->SetInput(conversionResult);
  meshWriter->Update();

  return EXIT_SUCCESS;
}

=================================================================

I really appreciate any help you can provide me. 

Thank you very much,
Breno




--
Sent from: http://itk-insight-users.2283740.n2.nabble.com/