Hello,<br><br>I have a binary mask dicom image with orientation 0 0 -1 / 0 1 0<br><br>I ran it through the DeformableModel1.cxx program,<br>to create and deform a mesh over the mask,<br>but my outputted deformed dicom image has the orientation 1 0 0 / 0 1 0<br>
<br>I tried copying the MetaDataDictionary from the mask to the mesh,<br>but this did not help.<br>I tried setting the direction of the outputted dicom image,<br>but when viewed in 3D Slicer, there was no image (all black).<br>
I suppose the points on the mesh would have to be transformed to fit the 0 0 -1 / 0 1 0 orientation.<br><br>Can anyone help?<br><br>Below is my code: ( a slightly modified version of DeformableModel1.cxx)<br><br>/*=========================================================================<br>
<br> Program: Insight Segmentation & Registration Toolkit<br> Module: $RCSfile: DeformableModel1.cxx,v $<br> Language: C++<br> Date: $Date: 2009-03-17 21:44:42 $<br> Version: $Revision: 1.30 $<br><br>
Copyright (c) Insight Software Consortium. All rights reserved.<br> See ITKCopyright.txt or <a href="http://www.itk.org/HTML/Copyright.htm">http://www.itk.org/HTML/Copyright.htm</a> for details.<br><br> This software is distributed WITHOUT ANY WARRANTY; without even <br>
the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR <br> PURPOSE. See the above copyright notices for more information.<br><br> =========================================================================*/<br>
#if defined(_MSC_VER)<br>#pragma warning ( disable : 4786 )<br>#endif<br><br>#ifdef __BORLANDC__<br>#define ITK_LEAN_AND_MEAN<br>#endif<br><br>// Software Guide : BeginLatex<br>//<br>// This example illustrates the use of the \doxygen{DeformableMesh3DFilter}<br>
// and \doxygen{BinaryMask3DMeshSource} in the hybrid segmentation framework.<br>//<br>// \begin{figure} \center<br>// \includegraphics[width=\textwidth]{DeformableModelCollaborationDiagram.eps}<br>// \itkcaption[Deformable model collaboration diagram]{Collaboration<br>
// diagram for the DeformableMesh3DFilter applied to a segmentation task.}<br>// \label{fig:DeformableModelCollaborationDiagram}<br>// \end{figure}<br>//<br>// The purpose of the DeformableMesh3DFilter is to take an initial surface<br>
// described by an \doxygen{Mesh} and deform it in order to adapt it to the<br>// shape of an anatomical structure in an image.<br>// Figure~\ref{fig:DeformableModelCollaborationDiagram} illustrates a typical<br>// setup for a segmentation method based on deformable models. First, an<br>
// initial mesh is generated using a binary mask and an isocontouring<br>// algorithm (such as marching cubes) to produce an initial mesh. The binary<br>// mask used here contains a simple shape which vaguely resembles the<br>
// anatomical structure that we want to segment. The application of the<br>// isocontouring algorithm produces a $3D$ mesh that has the shape of this<br>// initial structure. This initial mesh is passed as input to the deformable<br>
// model which will apply forces to the mesh points in order to reshape the<br>// surface until make it fit to the anatomical structures in the image.<br>//<br>// The forces to be applied on the surface are computed from an approximate<br>
// physical model that simulates an elastic deformation. Among the forces to<br>// be applied we need one that will pull the surface to the position of the<br>// edges in the anatomical structure. This force component is represented<br>
// here in the form of a vector field and is computed as illustrated in the<br>// lower left of Figure~\ref{fig:DeformableModelCollaborationDiagram}. The<br>// input image is passed to a<br>// \doxygen{GradientMagnitudeRecursiveGaussianImageFilter}, which computes<br>
// the magnitude of the image gradient. This scalar image is then passed to<br>// another gradient filter<br>// (\doxygen{GradientRecursiveGaussianImageFilter}). The output of this<br>// second gradient filter is a vector field in which every vector points to<br>
// the closest edge in the image and has a magnitude proportional to the<br>// second derivative of the image intensity along the direction of the<br>// gradient. Since this vector field is computed using Gaussian derivatives,<br>
// it is possible to regulate the smoothness of the vector field by playing<br>// with the value of sigma assigned to the Gaussian. Large values of sigma<br>// will result in a large capture radius, but will have poor precision in the<br>
// location of the edges. A reasonable strategy may involve the use of large<br>// sigmas for the initial iterations of the model and small sigmas to refine<br>// the model when it is close to the edges. A similar effect could be<br>
// achieved using multiresolution and taking advantage of the image pyramid<br>// structures already illustrated in the registration framework.<br>//<br>// \index{Deformable Models}<br>// \index{DeformableMesh3DFilter}<br>
//<br>// Software Guide : EndLatex <br><br><br>#include <iostream><br><br><br>// Software Guide : BeginLatex<br>//<br>// We start by including the headers of the main classes required for this<br>// example. The BinaryMask3DMeshSource is used to produce an initial<br>
// approximation of the shape to be segmented. This filter takes a binary<br>// image as input and produces a Mesh as output using the marching cube<br>// isocontouring algorithm.<br>//<br>// \index{itk::BinaryMask3DMeshSource!Header}<br>
//<br>// Software Guide : EndLatex <br><br>// Software Guide : BeginCodeSnippet<br>#include "itkBinaryMask3DMeshSource.h"<br>// Software Guide : EndCodeSnippet<br><br><br>// Software Guide : BeginLatex<br>// <br>
// Then we include the header of the DeformableMesh3DFilter that<br>// implements the deformable model algorithm.<br>//<br>// \index{itk::DeformableMesh3DFilter!Header}<br>//<br>// Software Guide : EndLatex <br><br>// Software Guide : BeginCodeSnippet<br>
#include "itkDeformableMesh3DFilter.h"<br>// Software Guide : EndCodeSnippet <br><br><br>// Software Guide : BeginLatex<br>// <br>// We also need the headers of the gradient filters that will be used for<br>
// computing the vector field. In our case they are the<br>// GradientMagnitudeRecursiveGaussianImageFilter and<br>// GradientRecursiveGaussianImageFilter.<br>//<br>// Software Guide : EndLatex <br><br>// Software Guide : BeginCodeSnippet<br>
#include "itkGradientRecursiveGaussianImageFilter.h"<br>#include "itkGradientMagnitudeRecursiveGaussianImageFilter.h"<br>// Software Guide : EndCodeSnippet <br><br><br>// Software Guide : BeginLatex<br>
// <br>// The main data structures required in this example are the Image<br>// and the Mesh classes. The deformable model \emph{per se} is<br>// represented as a Mesh.<br>//<br>// Software Guide : EndLatex <br><br>// Software Guide : BeginCodeSnippet<br>
#include "itkImage.h"<br>#include "itkMesh.h"<br>// Software Guide : EndCodeSnippet <br><br><br>// Software Guide : BeginLatex<br>// <br>// The \code{PixelType} of the image derivatives is represented with a<br>
// \doxygen{CovariantVector}. We include its header in the following line.<br>//<br>// Software Guide : EndLatex <br><br>// Software Guide : BeginCodeSnippet<br>#include "itkCovariantVector.h"<br>// Software Guide : EndCodeSnippet <br>
<br>// Software Guide : BeginLatex<br>// <br>// The deformed mesh is converted into a binary image using the<br>// \doxygen{PointSetToImageFilter}.<br>//<br>// Software Guide : EndLatex <br><br>// Software Guide : BeginCodeSnippet<br>
#include "itkPointSetToImageFilter.h"<br>// Software Guide : EndCodeSnippet <br><br>// Software Guide : BeginLatex<br>// <br>// In order to read both the input image and the mask image, we need the<br>// \doxygen{ImageFileReader} class. We also need the \doxygen{ImageFileWriter}<br>
// to save the resulting deformed mask image.<br>//<br>// Software Guide : EndLatex <br><br>// Software Guide : BeginCodeSnippet<br>#include "itkImageFileReader.h"<br>#include "itkImageFileWriter.h"<br>
// Software Guide : EndCodeSnippet <br><br>int main( int argc, char *argv[] )<br>{<br><br> if( argc < 4 )<br> {<br> std::cerr << "Missing Parameters " << std::endl;<br> std::cerr << "Usage: " << argv[0];<br>
std::cerr << " InputImage BinaryImage DeformedMaskImage" << std::endl;<br> return 1;<br> }<br> <br><br> // Software Guide : BeginLatex<br> //<br> // Here we declare the type of the image to be processed. This implies a<br>
// decision about the \code{PixelType} and the dimension. The<br> // DeformableMesh3DFilter is specialized for $3D$, so the choice<br> // is clear in our case.<br> //<br> // Software Guide : EndLatex <br><br> // Software Guide : BeginCodeSnippet<br>
<br> const unsigned int Dimension = 3;<br> <br> //typedef double PixelType;<br> typedef signed short PixelType;<br><br> typedef itk::Image<PixelType, Dimension> ImageType;<br>
// Software Guide : EndCodeSnippet<br><br><br> // Software Guide : BeginLatex<br> // <br> // The input to BinaryMask3DMeshSource is a binary mask that we<br> // will read from a file. This mask could be the result of a rough<br>
// segmentation algorithm applied previously to the same anatomical<br> // structure. We declare below the type of the binary mask image.<br> //<br> // Software Guide : EndLatex <br><br> // Software Guide : BeginCodeSnippet<br>
//typedef itk::Image< unsigned char, Dimension > BinaryImageType;<br> typedef itk::Image< signed short, Dimension > BinaryImageType;<br> // Software Guide : EndCodeSnippet<br><br><br> // Software Guide : BeginLatex<br>
// <br> // Then we define the type of the deformable mesh. We represent the<br> // deformable model using the Mesh class. The \code{double} type used as<br> // template parameter here is to be used to assign values to every point<br>
// of the Mesh.<br> //<br> // Software Guide : EndLatex <br><br> // Software Guide : BeginCodeSnippet<br> //typedef itk::Mesh<double> MeshType;<br> typedef itk::Mesh<signed short> MeshType;<br>
// Software Guide : EndCodeSnippet<br><br><br> // Software Guide : BeginLatex<br> // <br> // The following lines declare the type of the gradient image:<br> //<br> // Software Guide : EndLatex <br><br> // Software Guide : BeginCodeSnippet<br>
//typedef itk::CovariantVector< double, Dimension > GradientPixelType;<br> typedef itk::CovariantVector< signed short, Dimension > GradientPixelType;<br> typedef itk::Image< GradientPixelType, Dimension > GradientImageType;<br>
// Software Guide : EndCodeSnippet<br><br><br> // Software Guide : BeginLatex<br> // <br> // With it we can declare the type of the gradient filter and the gradient<br> // magnitude filter:<br> //<br> // Software Guide : EndLatex <br>
<br> // Software Guide : BeginCodeSnippet<br> typedef itk::GradientRecursiveGaussianImageFilter<ImageType, GradientImageType><br> GradientFilterType;<br> typedef itk::GradientMagnitudeRecursiveGaussianImageFilter<ImageType,ImageType><br>
GradientMagnitudeFilterType;<br> // Software Guide : EndCodeSnippet <br><br><br> // Software Guide : BeginLatex<br> // <br> // The filter implementing the isocontouring algorithm is the<br> // BinaryMask3DMeshSource filter.<br>
// <br> // \index{itk::BinaryMask3DMeshSource!Instantiation} <br> //<br> // Software Guide : EndLatex <br><br> // Software Guide : BeginCodeSnippet<br> typedef itk::BinaryMask3DMeshSource< BinaryImageType, MeshType > MeshSourceType;<br>
// Software Guide : EndCodeSnippet <br> // typedef itk::BinaryMaskToNarrowBandPointSetFilter< <br> // BinaryImageType, MeshType > MeshSourceType;<br><br><br> // Software Guide : BeginLatex<br>
//<br> // Now we instantiate the type of the DeformableMesh3DFilter that<br> // implements the deformable model algorithm. Note that both the input and<br> // output types of this filter are \doxygen{Mesh} classes.<br>
//<br> // \index{DeformableMesh3DFilter!Instantiation}<br> //<br> // Software Guide : EndLatex <br><br> // Software Guide : BeginCodeSnippet<br> typedef itk::DeformableMesh3DFilter<MeshType,MeshType> DeformableFilterType;<br>
// Software Guide : EndCodeSnippet<br><br><br> // Software Guide : BeginLatex<br> //<br> // Let's declare two readers. The first will read the image to be<br> // segmented. The second will read the binary mask containing a first<br>
// approximation of the segmentation that will be used to initialize a<br> // mesh for the deformable model.<br> //<br> // Software Guide : EndLatex <br><br> // Software Guide : BeginCodeSnippet<br> typedef itk::ImageFileReader< ImageType > ReaderType;<br>
typedef itk::ImageFileReader< BinaryImageType > BinaryReaderType;<br> ReaderType::Pointer imageReader = ReaderType::New();<br> BinaryReaderType::Pointer maskReader = BinaryReaderType::New();<br> // Software Guide : EndCodeSnippet<br>
<br><br> // Software Guide : BeginLatex<br> // <br> // In this example we take the filenames of the input image and the binary<br> // mask from the command line arguments.<br> //<br> // Software Guide : EndLatex <br>
<br> // Software Guide : BeginCodeSnippet<br> imageReader->SetFileName( argv[1] );<br> maskReader->SetFileName( argv[2] );<br> // Software Guide : EndCodeSnippet <br><br><br> // Software Guide : BeginLatex<br>
//<br> // We create here the GradientMagnitudeRecursiveGaussianImageFilter that<br> // will be used to compute the magnitude of the input image gradient. As<br> // usual, we invoke its \code{New()} method and assign the result to a<br>
// \doxygen{SmartPointer}.<br> //<br> // Software Guide : EndLatex <br><br> // Software Guide : BeginCodeSnippet<br> GradientMagnitudeFilterType::Pointer gradientMagnitudeFilter<br> = GradientMagnitudeFilterType::New();<br>
// Software Guide : EndCodeSnippet<br><br><br> // Software Guide : BeginLatex<br> //<br> // The output of the image reader is connected as input to the gradient<br> // magnitude filter. Then the value of sigma used to blur the image is<br>
// selected using the method \code{SetSigma()}.<br> //<br> // \index{itk::Gradient\-Magnitude\-Recursive\-Gaussian\-Image\-Filter!SetInput()}<br> // \index{itk::Gradient\-Magnitude\-Recursive\-Gaussian\-Image\-Filter!SetSigma()}<br>
//<br> // Software Guide : EndLatex <br><br> // Software Guide : BeginCodeSnippet<br> gradientMagnitudeFilter->SetInput( imageReader->GetOutput() ); <br> gradientMagnitudeFilter->SetSigma( 1.0 );<br> // Software Guide : EndCodeSnippet<br>
<br><br> // Software Guide : BeginLatex<br> // <br> // In the following line, we construct the gradient filter that will take<br> // the gradient magnitude of the input image that will be passed to the<br> // deformable model algorithm.<br>
//<br> // Software Guide : EndLatex <br><br> // Software Guide : BeginCodeSnippet<br> GradientFilterType::Pointer gradientMapFilter = GradientFilterType::New();<br> // Software Guide : EndCodeSnippet <br><br><br>
// Software Guide : BeginLatex<br> //<br> // The magnitude of the gradient is now passed to the next step of<br> // gradient computation. This allows us to obtain a second derivative of<br> // the initial image with the gradient vector pointing to the maxima of<br>
// the input image gradient. This gradient map will have the properties<br> // desirable for attracting the deformable model to the edges of the<br> // anatomical structure on the image. Once again we must select the value<br>
// of sigma to be used in the blurring process.<br> //<br> // Software Guide : EndLatex <br><br> // Software Guide : BeginCodeSnippet<br> gradientMapFilter->SetInput( gradientMagnitudeFilter->GetOutput());<br>
gradientMapFilter->SetSigma( 1.0 );<br> // Software Guide : EndCodeSnippet<br><br><br> // Software Guide : BeginLatex<br> // <br> // At this point, we are ready to compute the vector field. This is done<br> // simply by invoking the \code{Update()} method on the second derivative<br>
// filter. This was illustrated in<br> // Figure~\ref{fig:DeformableModelCollaborationDiagram}.<br> //<br> // Software Guide : EndLatex <br><br> try<br> {<br> // Software Guide : BeginCodeSnippet<br> gradientMapFilter->Update();<br>
// Software Guide : EndCodeSnippet<br> }<br> catch( itk::ExceptionObject & e )<br> {<br> std::cerr << "Exception caught when updating gradientMapFilter " << std::endl;<br> std::cerr << e << std::endl;<br>
return -1;<br> }<br> <br><br><br> std::cout << "The gradient map created!" << std::endl;<br><br><br> // Software Guide : BeginLatex<br> // <br> // Now we can construct the mesh source filter that implements the<br>
// isocontouring algorithm.<br> //<br> // \index{BinaryMask3DMeshSource!New()}<br> // \index{BinaryMask3DMeshSource!Pointer}<br> //<br> // Software Guide : EndLatex <br><br> // Software Guide : BeginCodeSnippet<br>
MeshSourceType::Pointer meshSource = MeshSourceType::New();<br> // Software Guide : EndCodeSnippet<br><br><br> // Software Guide : BeginLatex<br> // <br> // Then we create the filter implementing the deformable model and set its<br>
// input to the output of the binary mask mesh source. We also set the<br> // vector field using the \code{SetGradient()} method.<br> // <br> // Software Guide : EndLatex <br><br> // Software Guide : BeginCodeSnippet<br>
DeformableFilterType::Pointer deformableModelFilter = <br> DeformableFilterType::New();<br> deformableModelFilter->SetGradient( gradientMapFilter->GetOutput() );<br> // Software Guide : EndCodeSnippet<br>
<br><br> // Software Guide : BeginLatex<br> // <br> // Here we connect the output of the binary mask reader to the input of<br> // the BinaryMask3DMeshSource that will apply the isocontouring algorithm<br> // and generate the initial mesh to be deformed. We must also select the<br>
// value to be used for representing the binary object in the image. In<br> // this case we select the value $200$ and pass it to the filter using its<br> // method \code{SetObjectValue()}.<br> //<br> // \index{itk::BinaryMask3DMeshSource!SetInput()}<br>
// \index{itk::BinaryMask3DMeshSource!SetObjectValue()}<br> //<br> // Software Guide : EndLatex <br><br> // Software Guide : BeginCodeSnippet<br> BinaryImageType::Pointer mask = maskReader->GetOutput();<br> meshSource->SetInput( mask );<br>
//meshSource->SetObjectValue( 200 );<br> meshSource->SetObjectValue( 1 );<br><br> std::cout << "Creating mesh..." << std::endl;<br> try <br> {<br> meshSource->Update();<br> }<br>
catch( itk::ExceptionObject & excep )<br> {<br> std::cerr << "Exception Caught !" << std::endl;<br> std::cerr << excep << std::endl;<br> }<br><br> deformableModelFilter->SetInput( meshSource->GetOutput() );<br>
// Software Guide : EndCodeSnippet<br><br> meshSource->GetOutput()->Print(std::cout);<br><br> std::cout << "Deformable mesh created using Marching Cube!" << std::endl;<br><br> // Software Guide : BeginLatex<br>
// <br> // Next, we set the parameters of the deformable model computation.<br> // \code{Stiffness} defines the model stiffness in the vertical and<br> // horizontal directions on the deformable surface. \code{Scale} helps to<br>
// accommodate the deformable mesh to gradient maps of different size.<br> //<br> // \index{itk::DeformableMesh3DFilter!SetStiffness()}<br> // \index{itk::DeformableMesh3DFilter!SetScale()}<br> //<br> // Software Guide : EndLatex <br>
<br> // Software Guide : BeginCodeSnippet<br> typedef itk::CovariantVector<double, 2> double2DVector;<br> typedef itk::CovariantVector<double, 3> double3DVector;<br><br> double2DVector stiffness;<br>
stiffness[0] = 0.0001;<br> stiffness[1] = 0.1;<br><br> double3DVector scale;<br> scale[0] = 1.0;<br> scale[1] = 1.0; <br> scale[2] = 1.0;<br><br> deformableModelFilter->SetStiffness( stiffness );<br> deformableModelFilter->SetScale( scale );<br>
// Software Guide : EndCodeSnippet<br><br><br> // Software Guide : BeginLatex<br> // <br> // Other parameters to be set are the gradient magnitude, the time step and<br> // the step threshold. The gradient magnitude controls the magnitude of the<br>
// external force. The time step controls the length of each step during<br> // deformation. Step threshold is the number of the steps the model will<br> // deform.<br> //<br> // \index{itk::DeformableMesh3DFilter!SetGradientMagnitude()}<br>
// \index{itk::DeformableMesh3DFilter!SetTimeStep()}<br> // \index{itk::DeformableMesh3DFilter!SetStepThreshold()}<br> //<br> // Software Guide : EndLatex <br><br> // Software Guide : BeginCodeSnippet<br> deformableModelFilter->SetGradientMagnitude( 0.8 );<br>
deformableModelFilter->SetTimeStep( 0.01 );<br> deformableModelFilter->SetStepThreshold( 60 );<br> // Software Guide : EndCodeSnippet <br><br><br> std::cout << "Deformable mesh fitting...";<br>
<br><br> // Software Guide : BeginLatex<br> // <br> // Finally, we trigger the execution of the deformable model computation<br> // using the \code{Update()} method of the DeformableMesh3DFilter. As<br> // usual, the call to \code{Update()} should be placed in a<br>
// \code{try/catch} block in case any exceptions are thrown.<br> // <br> // \index{DeformableMesh3DFilter!Update()}<br> // <br> // Software Guide : EndLatex <br><br> // Software Guide : BeginCodeSnippet<br> try <br>
{<br> deformableModelFilter->Update();<br> }<br> catch( itk::ExceptionObject & excep )<br> {<br> std::cerr << "Exception Caught !" << std::endl;<br> std::cerr << excep << std::endl;<br>
}<br> // Software Guide : EndCodeSnippet<br><br> std::cout << "Mesh Source: " << meshSource;<br><br> // Software Guide : BeginLatex<br> // <br> // The \doxygen{PointSetToImageFilter} takes the deformed<br>
// mesh and produce a binary image corresponding to the node<br> // of the mesh. Note that only the nodes are producing the image<br> // and not the cells. See the section on SpatialObjects to produce<br> // a complete binary image from cells using the \doxygen{MeshSpatialObject}<br>
// combined with the \doxygen{SpatialObjectToImageFilter}.<br> // However, using SpatialObjects is computationally more expensive.<br> // <br> // Software Guide : EndLatex <br> <br> // Software Guide : BeginCodeSnippet<br>
typedef itk::PointSetToImageFilter<MeshType,ImageType> MeshFilterType;<br> MeshFilterType::Pointer meshFilter = MeshFilterType::New();<br> meshFilter->SetOrigin(mask->GetOrigin());<br> meshFilter->SetSize(mask->GetLargestPossibleRegion().GetSize());<br>
meshFilter->SetSpacing(mask->GetSpacing());<br><br> //added<br> //meshFilter->SetDirection(mask->GetDirection());<br><br> meshFilter->SetInput(meshSource->GetOutput());<br> try <br> {<br> meshFilter->Update();<br>
}<br> catch( itk::ExceptionObject & excep )<br> {<br> std::cerr << "Exception Caught !" << std::endl;<br> std::cerr << excep << std::endl;<br> }<br> // Software Guide : EndCodeSnippet<br>
<br> // Software Guide : BeginLatex<br> // <br> // The resulting deformed binary mask can be written on disk<br> // using the \doxygen{ImageFileWriter}.<br> // <br> // Software Guide : EndLatex <br><br> // Software Guide : BeginCodeSnippet<br>
typedef itk::ImageFileWriter<ImageType> WriterType;<br> WriterType::Pointer writer = WriterType::New();<br> writer->SetInput(meshFilter->GetOutput());<br> writer->SetFileName(argv[3]);<br> writer->Update();<br>
// Software Guide : EndCodeSnippet<br><br> // Software Guide : BeginLatex<br> //<br> // Note that in order to successfully segment images, input<br> // parameters must be adjusted to reflect the characteristics of the<br>
// data. The output of the filter is an Mesh. Users can use<br> // their own visualization packages to see the segmentation results.<br> // <br> // Software Guide : EndLatex <br><br> return EXIT_SUCCESS;<br>}<br>
<br>