[Insight-users] DisplacementFieldTransform on CT 3D images

[agatte]

Hi ;)

I am fresh itk user. I am trying to do on 3D images  example
DeformationFieldTransform from itk wiki:
http://www.itk.org/Wiki/ITK/Examples/Registration/DeformationFieldTransform
I have already  made some changes for 3D data in this example. But I have a
problem with debuging.
I received an exception.

Unhandled exception at 0x758cd36f in DeformationFieldTransform.exe:
Microsoft C++ exception: itk::MemoryAllocationError at memory location
0x0199f5b4..

 if ( !data )
    {
    // We cannot construct an error string here because we may be out
    // of memory.  Do not use the exception macro.
    throw MemoryAllocationError(__FILE__, __LINE__,
                                "Failed to allocate memory for image.",
                                ITK_LOCATION);
    }
  return data;
}


Could anyone look at this code ?
I would appreciate for any help.
I don't have any idea. I don't know what is it wrong ?

Here my code :


#include "itkImageFileReader.h"
#include "itkImageFileWriter.h"
#include "itkImage.h"
#include "itkVector.h"
#if ITK_VERSION_MAJOR < 4
#include "itkDeformationFieldTransform.h"
#include "itkDeformationFieldSource.h"
#else
#include "itkVectorLinearInterpolateImageFunction.h"
#include "itkDisplacementFieldTransform.h"
#include "itkLandmarkDisplacementFieldSource.h"
#endif
#include "itkResampleImageFilter.h"
 
const     unsigned int   Dimension = 3;
typedef   signed short PixelType;
typedef   itk::Image< PixelType, Dimension > ImageType;
typedef itk::ImageFileReader<ImageType> ReaderType;
typedef itk::ImageFileWriter<ImageType> WriterType;



int main()
{
	
	double tab[24] = {1,1,1,1,1,1, 0,0,0,0,0,0, 0,0,0,0,0,0, 1,1,1,1,1,1};


#if ITK_VERSION_MAJOR < 4
  typedef   float           VectorComponentType;
#else
  typedef   double          VectorComponentType;
#endif
  typedef   itk::Vector< VectorComponentType, Dimension >    VectorType;
  typedef   itk::Image< VectorType,  Dimension >   DeformationFieldType;

  ReaderType::Pointer readerFixed = ReaderType::New();
  ReaderType::Pointer readerMoving = ReaderType::New();

  readerFixed->SetFileName("fantom2V.vtk");
  readerMoving->SetFileName("fantom3V.vtk");
  readerFixed->Update();
  readerMoving->Update();

  std::cout<<readerMoving<<std::endl;
  std::cout<<readerFixed<<std::endl;

  ImageType::Pointer fixedImage = ImageType::New();
  ImageType::Pointer movingImage = ImageType::New();
  fixedImage = readerFixed->GetOutput();
  movingImage = readerMoving->GetOutput();

#if ITK_VERSION_MAJOR < 4
  typedef itk::DeformationFieldSource<DeformationFieldType> 
DeformationFieldSourceType;
#else
  typedef itk::LandmarkDisplacementFieldSource<DeformationFieldType> 
DeformationFieldSourceType;
#endif
  DeformationFieldSourceType::Pointer deformationFieldSource =
DeformationFieldSourceType::New();
  deformationFieldSource->SetOutputSpacing( fixedImage->GetSpacing() );
  deformationFieldSource->SetOutputOrigin( fixedImage->GetOrigin() );
  deformationFieldSource->SetOutputRegion(
fixedImage->GetLargestPossibleRegion() );
  deformationFieldSource->SetOutputDirection( fixedImage->GetDirection() );
 
  //  Create source and target landmarks.
  typedef DeformationFieldSourceType::LandmarkContainerPointer  
LandmarkContainerPointer;
  typedef DeformationFieldSourceType::LandmarkContainer         
LandmarkContainerType;
  typedef DeformationFieldSourceType::LandmarkPointType         
LandmarkPointType;
 
  LandmarkContainerType::Pointer sourceLandmarks =
LandmarkContainerType::New();
  LandmarkContainerType::Pointer targetLandmarks =
LandmarkContainerType::New();
 
  LandmarkPointType sourcePoint;
  LandmarkPointType targetPoint;
 
  sourcePoint[0] = tab[0];
  sourcePoint[1] = tab[1];
  sourcePoint[2] = tab[2];
  targetPoint[0] = tab[3];
  targetPoint[1] = tab[4];
  targetPoint[2] = tab[5];
  sourceLandmarks->InsertElement( 0, sourcePoint );
  targetLandmarks->InsertElement( 0, targetPoint );
 
  sourcePoint[0] = tab[6];
  sourcePoint[1] = tab[7];
  sourcePoint[2] = tab[8];
  targetPoint[0] = tab[9];
  targetPoint[1] = tab[10];
  targetPoint[2] = tab[11];
  sourceLandmarks->InsertElement( 1, sourcePoint );
  targetLandmarks->InsertElement( 1, targetPoint );
 
  sourcePoint[0] = tab[12];
  sourcePoint[1] = tab[13];
  sourcePoint[2] = tab[14];
  targetPoint[0] = tab[15];
  targetPoint[1] = tab[16];
  targetPoint[2] = tab[17];
  sourceLandmarks->InsertElement( 2, sourcePoint );
  targetLandmarks->InsertElement( 2, targetPoint );
 
  sourcePoint[0] = tab[18];
  sourcePoint[1] = tab[19];
  sourcePoint[2] = tab[20];
  targetPoint[0] = tab[21];
  targetPoint[1] = tab[22];
  targetPoint[2] = tab[23];
  sourceLandmarks->InsertElement( 3, sourcePoint );
  targetLandmarks->InsertElement( 3, targetPoint );

  //sourcePoint[0] = tab[24];
  //sourcePoint[1] = tab[25];
  //sourcePoint[2] = tab[26];
  //targetPoint[0] = tab[27];
  //targetPoint[1] = tab[28];
  //targetPoint[2] = tab[29];
  //sourceLandmarks->InsertElement( 4, sourcePoint );
  //targetLandmarks->InsertElement( 4, targetPoint );
 
  deformationFieldSource->SetSourceLandmarks( sourceLandmarks.GetPointer()
);
  deformationFieldSource->SetTargetLandmarks( targetLandmarks.GetPointer()
);
  deformationFieldSource->UpdateLargestPossibleRegion();
 
  //// Write the deformation field
 // {
  typedef itk::ImageFileWriter<  DeformationFieldType  > WriterType;
  WriterType::Pointer writer = WriterType::New();
  writer->SetInput (  deformationFieldSource->GetOutput() );
  writer->SetFileName( "deformationFieldfantomV.vtk" );
  writer->Update();
 // }
 
#if ITK_VERSION_MAJOR < 4
  typedef itk::DeformationFieldTransform<VectorComponentType, 2> 
DeformationFieldTransformType;
#else
 // typedef itk::DisplacementFieldTransform<VectorComponentType, 3> 
DeformationFieldTransformType;
 typedef itk::DisplacementFieldTransform<VectorComponentType, 3> 
DeformationFieldTransformType;
#endif
  DeformationFieldTransformType::Pointer deformationFieldTransform =
DeformationFieldTransformType::New();
 
#if ITK_VERSION_MAJOR < 4
  deformationFieldTransform->SetDeformationField(
deformationFieldSource->GetOutput() );
#else
  deformationFieldTransform->SetDisplacementField(
deformationFieldSource->GetOutput() );
#endif  

  std::cout<<" "<<deformationFieldTransform<<std::endl;

  system("Pause");


  typedef itk::ResampleImageFilter<ImageType, ImageType, VectorComponentType
>    ResampleFilterType;
  ResampleFilterType::Pointer resampleFilter = ResampleFilterType::New();
  resampleFilter->SetInput( movingImage );
  resampleFilter->SetTransform( deformationFieldTransform );
  resampleFilter->SetSize( fixedImage->GetLargestPossibleRegion().GetSize()
);
  resampleFilter->SetOutputOrigin(  fixedImage->GetOrigin() );
  resampleFilter->SetOutputSpacing( fixedImage->GetSpacing() );
  resampleFilter->SetOutputDirection( fixedImage->GetDirection() );
  resampleFilter->SetDefaultPixelValue( 200 );
  resampleFilter->GetOutput();
 
  // Write the output
  typedef itk::ImageFileWriter<  ImageType  > WriterType2;
  WriterType2::Pointer writerOutput = WriterType2::New();
  writerOutput ->SetInput (  resampleFilter->GetOutput() );
  writerOutput ->SetFileName( "outputfantomV.dcm" );
  writerOutput ->Update();
 
  return EXIT_SUCCESS;
}
 
-- 
View this message in context: http://old.nabble.com/DisplacementFieldTransform-on--CT--3D-images-tp33760464p33760464.html
Sent from the ITK - Users mailing list archive at Nabble.com.