ITK/Examples/Registration/LandmarkBasedTransformInitializer

LandmarkBasedTransformInitializer.cxx

<source lang="cpp">

1. include "itkImageFileWriter.h"
2. include "itkImage.h"
3. include "itkVector.h"
4. include "itkResampleImageFilter.h"
5. include "itkLandmarkBasedTransformInitializer.h"
6. include "itkRigid2DTransform.h"

const unsigned int Dimension = 2; typedef unsigned char PixelType; typedef itk::Image< PixelType, Dimension > ImageType;

static void CreateFixedImage(ImageType::Pointer image); static void CreateMovingImage(ImageType::Pointer image);

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

 typedef   float          VectorComponentType;

 typedef   itk::Vector< VectorComponentType, Dimension >    VectorType;
 typedef   itk::Image< VectorType,  Dimension >   DeformationFieldType;

 ImageType::Pointer fixedImage = ImageType::New();
 CreateFixedImage(fixedImage);
 
 ImageType::Pointer movingImage = ImageType::New();
 CreateMovingImage(movingImage);

 typedef itk::Rigid2DTransform< double > Rigid2DTransformType;
 typedef itk::LandmarkBasedTransformInitializer< Rigid2DTransformType, ImageType, ImageType > 
     LandmarkBasedTransformInitializerType;
     
 LandmarkBasedTransformInitializerType::Pointer landmarkBasedTransformInitializer =
   LandmarkBasedTransformInitializerType::New();
 //  Create source and target landmarks.
 typedef LandmarkBasedTransformInitializerType::LandmarkPointContainer     LandmarkContainerType;
 typedef LandmarkBasedTransformInitializerType::LandmarkPointType          LandmarkPointType;

 LandmarkContainerType fixedLandmarks;
 LandmarkContainerType movingLandmarks;

 LandmarkPointType fixedPoint;
 LandmarkPointType movingPoint;

 fixedPoint[0] = 10;
 fixedPoint[1] = 10;
 movingPoint[0] = 50;
 movingPoint[1] = 50;
 fixedLandmarks.push_back( fixedPoint );
 movingLandmarks.push_back( movingPoint );

 fixedPoint[0] = 10;
 fixedPoint[1] = 20;
 movingPoint[0] = 50;
 movingPoint[1] = 60;
 fixedLandmarks.push_back( fixedPoint );
 movingLandmarks.push_back( movingPoint );

 fixedPoint[0] = 20;
 fixedPoint[1] = 10;
 movingPoint[0] = 60;
 movingPoint[1] = 50;
 fixedLandmarks.push_back( fixedPoint );
 movingLandmarks.push_back( movingPoint );

 fixedPoint[0] = 20;
 fixedPoint[1] = 20;
 movingPoint[0] = 60;
 movingPoint[1] = 60;
 fixedLandmarks.push_back( fixedPoint );
 movingLandmarks.push_back( movingPoint );

 landmarkBasedTransformInitializer->SetFixedLandmarks( fixedLandmarks );
 landmarkBasedTransformInitializer->SetMovingLandmarks( movingLandmarks );
 
 Rigid2DTransformType::Pointer transform = Rigid2DTransformType::New();
 
 transform->SetIdentity();
 landmarkBasedTransformInitializer->SetTransform(transform);
 landmarkBasedTransformInitializer->InitializeTransform();
 
 typedef itk::ResampleImageFilter<ImageType, ImageType, double >    ResampleFilterType;
 ResampleFilterType::Pointer resampleFilter = ResampleFilterType::New();
 resampleFilter->SetInput( movingImage );
 resampleFilter->SetTransform( transform );
 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  > WriterType;
 WriterType::Pointer writer = WriterType::New();
 writer->SetInput (  resampleFilter->GetOutput() );
 writer->SetFileName( "output.png" );
 writer->Update();
 
 return EXIT_SUCCESS;

}

void CreateFixedImage(ImageType::Pointer image) {

 // Create a black image with a white square
 ImageType::IndexType start;
 start.Fill(0);

 ImageType::SizeType size;
 size.Fill(100);

 ImageType::RegionType region;
 region.SetSize(size);
 region.SetIndex(start);

 image->SetRegions(region);
 image->Allocate();
 image->FillBuffer(0);

 itk::ImageRegionIterator<ImageType> imageIterator(image,region);

 while(!imageIterator.IsAtEnd())
   {
   if(imageIterator.GetIndex()[0] > 10 && imageIterator.GetIndex()[0] < 20 &&
     imageIterator.GetIndex()[1] > 10 && imageIterator.GetIndex()[1] < 20)
     {
     imageIterator.Set(255);
     }
   ++imageIterator;
   }

 // Write the deformation field
 typedef itk::ImageFileWriter<  ImageType  > WriterType;
 WriterType::Pointer writer = WriterType::New();
 writer->SetInput (  image );
 writer->SetFileName( "fixed.png" );
 writer->Update();

}

void CreateMovingImage(ImageType::Pointer image) {

 // Create a black image with a white square
 ImageType::IndexType start;
 start.Fill(0);

 ImageType::SizeType size;
 size.Fill(100);

 ImageType::RegionType region;
 region.SetSize(size);
 region.SetIndex(start);

 image->SetRegions(region);
 image->Allocate();
 image->FillBuffer(0);

 itk::ImageRegionIterator<ImageType> imageIterator(image,region);

 while(!imageIterator.IsAtEnd())
   {
   if(imageIterator.GetIndex()[0] > 50 && imageIterator.GetIndex()[0] < 60 &&
     imageIterator.GetIndex()[1] > 50 && imageIterator.GetIndex()[1] < 60)
     {
     imageIterator.Set(100);
     }
   ++imageIterator;
   }

 // Write the deformation field
 typedef itk::ImageFileWriter<  ImageType  > WriterType;
 WriterType::Pointer writer = WriterType::New();
 writer->SetInput (  image );
 writer->SetFileName( "moving.png" );
 writer->Update();

}

</source>

LandmarkBasedTransformInitializer.py

<source lang="python">

1. !/usr/bin/env python

import itk

Dimension = 2 PixelType = itk.ctype('unsigned char') ImageType = itk.Image[PixelType, Dimension]

def CreateFixedImage(image):

   start = itk.Index[Dimension]()
   start.Fill(0)

   size = itk.Size[Dimension]()
   size.Fill(100)

   region = itk.ImageRegion[Dimension]()
   region.SetSize(size)
   region.SetIndex(start)

   image.SetRegions(region)
   image.Allocate()
   image.FillBuffer(0)

   index = itk.Index[Dimension]()
   for ii in range(10, 20):
       for jj in range(10, 20):
           index[0] = ii
           index[1] = jj
           image.SetPixel(index, 255)

   writer = itk.ImageFileWriter.New(Input=image)
   writer.SetFileName("fixed.png")
   writer.Update()

def CreateMovingImage(image):

   start = itk.Index[Dimension]()
   start.Fill(0)

   size = itk.Size[Dimension]()
   size.Fill(100)

   region = itk.ImageRegion[Dimension]()
   region.SetSize(size)
   region.SetIndex(start)

   image.SetRegions(region)
   image.Allocate()
   image.FillBuffer(0)

   index = itk.Index[Dimension]()
   for ii in range(50, 60):
       for jj in range(50, 60):
           index[0] = ii
           index[1] = jj
           image.SetPixel(index, 100)

   writer = itk.ImageFileWriter.New(Input=image)
   writer.SetFileName("moving.png")
   writer.Update()

fixed_image = ImageType.New() CreateFixedImage(fixed_image)

moving_image = ImageType.New() CreateMovingImage(moving_image)

VectorComponentType = itk.ctype('float') VectorType = itk.Vector[VectorComponentType, Dimension] DisplacementFieldType = itk.Image[VectorType, Dimension]

Rigid2DTransformType = itk.Rigid2DTransform[itk.D] landmark_based_transform_initializer = \

       itk.LandmarkBasedTransformInitializer[itk.Transform[itk.D, Dimension,
           Dimension]].New()

LandmarkPointType = itk.Point[itk.D, Dimension] LandmarkContainerType = itk.vector[LandmarkPointType]

fixed_landmarks = LandmarkContainerType() moving_landmarks = LandmarkContainerType()

fixed_point = LandmarkPointType() moving_point = LandmarkPointType()

fixed_point[0] = 10 fixed_point[1] = 10 moving_point[0] = 50 moving_point[1] = 50 fixed_landmarks.push_back(fixed_point) moving_landmarks.push_back(moving_point)

fixed_point[0] = 10 fixed_point[1] = 20 moving_point[0] = 50 moving_point[1] = 60 fixed_landmarks.push_back(fixed_point) moving_landmarks.push_back(moving_point)

fixed_point[0] = 20 fixed_point[1] = 10 moving_point[0] = 60 moving_point[1] = 50 fixed_landmarks.push_back(fixed_point) moving_landmarks.push_back(moving_point)

fixed_point[0] = 20 fixed_point[1] = 20 moving_point[0] = 60 moving_point[1] = 60 fixed_landmarks.push_back(fixed_point) moving_landmarks.push_back(moving_point)

landmark_based_transform_initializer.SetFixedLandmarks(fixed_landmarks) landmark_based_transform_initializer.SetMovingLandmarks(moving_landmarks)

transform = Rigid2DTransformType.New() transform.SetIdentity() landmark_based_transform_initializer.SetTransform(transform) landmark_based_transform_initializer.InitializeTransform()

resampler = itk.ResampleImageFilter.New(Input=moving_image) resampler.SetTransform(transform)

1. resampler.SetReferenceImage(fixed_image)

resampler.SetSize(fixed_image.GetLargestPossibleRegion().GetSize()) resampler.SetOutputOrigin(fixed_image.GetOrigin()) resampler.SetOutputSpacing(fixed_image.GetSpacing()) resampler.SetOutputDirection(fixed_image.GetDirection()) resampler.SetDefaultPixelValue(200) resampler.UpdateLargestPossibleRegion()

writer = itk.ImageFileWriter.New(Input=resampler.GetOutput()) writer.SetFileName("output.png") writer.UpdateLargestPossibleRegion() </source>