ITK Registration Methods

1. ITK Registration Methods Kitware Inc.

2. Overview • Image Resampling • Registration Framework • Multi-Modality • Multi-Resolution • Deformable registration

3. Image Resampling

4. Why Resampling ? Resamplingis the Essence of Intensity-BasedImage Registration

5. What is an Image ? An Image is a sampling of a continuous fieldusing a discrete grid

6. y x Image Origin & Spacing Spacing (Sx) Spacing (Sy) Origin (Ox,Oy)

7. Image Sampling Grid Spacing (Sx) Spacing (Sy) Origin (Ox,Oy)

8. Image Pixel Spacing (Sx) Pixel Value Pixel Region Spacing (Sy) Origin (Ox,Oy)

9. Image Indices Spacing (Sx) Pixel Index [0,7] [4,7] [0,6] [0,5] [0,4] [0,3] [0,2] Spacing (Sy) [0,1] [0,0] [1,0] [2,0] [3,0] [4,0] [5,0] Origin (Ox,Oy)

10. Index to Physical Coordinates Spacing (Sx) Pixel Index [0,7] [4,7] [0,6] [0,5] P[0] = Index[0] x Spacing[0] + Origin[0] P[1] = Index[1] x Spacing[1] + Origin[1] [0,4] Index[0] = floor( ( P[0] - Origin[0] ) / Spacing[0] + 0.5 ) [0,3] Index[1] = floor( ( P[1] - Origin[1] ) / Spacing[1] + 0.5 ) [0,2] Spacing (Sy) [0,1] [0,0] [1,0] [2,0] [3,0] [4,0] [5,0] Origin (Ox,Oy)

11. Image Region Spacing (Sx) Pixel Value Image Region Pixel Region Spacing (Sy) Origin (Ox,Oy)

12. Image Region Spacing (Sx) Region Size Image Region [3,5] Starting Index [2,3] Spacing (Sy) Origin (Ox,Oy)

13. Basic Resampling ResamplingTrivial Cases

14. Spacing ( 2 x Sx ) Spacing ( 2 x Sy ) Sub-Sampling by Half Spacing (Sx) Image Region Spacing (Sy) Origin (Ox,Oy)

15. Sub-Sampling by Half New Spacing S’x New Spacing S’y New Origin (O’x,O’y) Origin (Ox,Oy)

16. Super-Sampling by Double Spacing ( Sx/2 ) Spacing (Sx) Spacing ( Sy/2 ) Image Region Spacing (Sy) Origin (Ox,Oy)

17. New Spacing S’y Super-Sampling by Double New Spacing S’x New Origin (O’x,O’y) Origin (Ox,Oy)

18. Resampling in ITK itk::ResampleImageFilter

19. Input Image Output Image Resampling in ITK Origin Spacing Region Size Resample Filter Region Start Transform Interpolator

20. Resample Image Filter #include "itkImage.h" #include "itkResampleImageFilter.h" #include "itkIdentityTransform.h“ #include "itkLinearInterpolateImageFunction.h" typedef itk::Image< char, 2 > ImageType; ImageType::Pointer inputImage = GetImageSomeHow(); typedef itk::ResampleImageFilter< ImageType > FilterType; FilterType::Pointerresampler = FilterType::New(); ImageType::SizeTypesize; size[0] = 200; size[1] = 300; ImageType::IndexTypestart; start[0] = 0; start[1] = 0;

21. Resample Image Filter ImageType::PointTypeorigin; origin[0] = 10.0; // millimeters origin[1] = 25.5; // millimeters ImageType::SpacingTypespacing; spacing[0] = 2.0; // millimeters spacing[1] = 1.5; // millimeters resampler->SetOutputSpacing( spacing ); resampler->SetOutputOrigin( origin ); resampler->SetSize( size ); resampler->SetOutputStartIndex( start ); resampler->SetDefaultPixelValue( 100 ); resampler->SetInput( inputImage );

22. Resample Image Filter typedef itk::LinearInterpolateImageFunction< ImageType, double > InterpolatorType; InterpolatorType::Pointer interpolator = InterpolatorType::New(); typedef itk::TranslationTransform< double, 2 > TransformType; TransformType::Pointertransform = TransformType::New(); transform->SetIdentity(); resampler->SetInterpolator( interpolator ); resampler->SetTransform( transform ); resampler->Update(); const ImageType * outputImage = resampler->GetOutput();

23. Basic Registration

24. j j i i Image Grid Image Grid y’ y y x Space Transform x’ x Physical Coordinates Physical Coordinates Coordinate System Conversions

25. I will not register images in pixel space I will not register images in pixel space I will not register images in pixel space I will not register images in pixel space I will not register images in pixel space I will not register images in pix Things I will not do…

26. j j i i Moving Image Grid Fixed Image Grid y’ y Space Transform x’ x Moving ImagePhysical Coordinates Fixed ImagePhysical Coordinates Fixed Image & Moving Image

27. Selecting Moving & Fixed Images In principle the denomination of Fixed Image & Moving Image is arbitrary In practice the moving image is the one that will be resampled into the fixed image coordinate system

28. MRI-T2 PET 128 x 128 pixels 256 x 256 pixels Quiz #1 Images from the same patient Moving Image ? Fixed Image ? Images provided as part of the project: “Retrospective Image Registration Evaluation”, NIH, Project No. 8R01EB002124-03, Principal Investigator, J. Michael Fitzpatrick, Vanderbilt University, Nashville, TN.

29. MRI-T2 PET 128 x 128 pixels 256 x 256 pixels Scaling Transform Quiz #2 Images from the same patient What scale factor ? • 2.0 • 1.0 • 0.5 Images provided as part of the project: “Retrospective Image Registration Evaluation”, NIH, Project No. 8R01EB002124-03, Principal Investigator, J. Michael Fitzpatrick, Vanderbilt University, Nashville, TN.

30. I will not register images in pixel space I will not register images in pixel space I will not register images in pixel space I will not register images in pixel space I will not register images in pixel space I will not register images in pix Things I will not do…

31. Components Registration in ITK MultiResolutionRegistrationFramework ImageRegistrationFramework PDEBasedRegistration FEMBasedRegistration

32. Registration Framework

33. Components Registration Method FixedImage Metric Optimizer Interpolator MovingImage Transform

34. Image Metrics • Mean Squares • Normalized Correlation • Mean Reciprocal Square Difference • Mutual Information- Viola-Wells- Mattes- Histogram based- Histogram normalized

35. Transforms • Translation • Scaling • Rotation • Rigid3D • Rigid2D • Affine • BSplines • Splines: TPS, EBS, VS

36. Optimizers • Gradient Descent • Regular Step Gradient Descent • Conjugate Gradient • Levenberg-Marquardt • One plus One Evolutionary Algorithm

37. Interpolators • Nearest Neighbor • Linear • BSpline

38. Exercise 20

39. Registered Moving Image Fixed Image Moving Image Image Registration

40. Image Registration #include ”itkImageRegistrationMethod.h” #include ”itkTranslationTransform.h” #include ”itkMeanSquaresImageToImageMetric.h” #include ”itkLinearInterpolateImageFunction.h” #include ”itkRegularStepGradientDescentOptimizer.h” #include ”itkImage.h” #include ”itkImageFileReader.h” #include ”itkImageFileWriter.h” #include ”itkResampleImageFilter.h”

41. Image Registration const unsigned int Dimension = 2; typedef unsigned char PixelType; typedef itk::Image< PixelType , Dimension > FixedImageType; typedef itk::Image< PixelType , Dimension > MovingImageType; typedef itk::TranslationTransform< double, Dimension > TransformType; typedef itk::RegularStepGradientDescentOptimizer OptimizerType; typedef itk::LinearInterpolateImageFunction< MovingImageType , double > InterpolatorType; typedef itk::MeanSquaresImageToImageMetric< FixedImageType , MovingImageType> MetricType; typedef itk::ImageRegistrationMethod< FixedImageType , MovingImageType> RegistrationType;

42. Image Registration TransformType::Pointertransform = TransformType::New(); OptimizerType::Pointeroptimizer = OptimizerType::New(); InterpolatorType::Pointerinterpolator = InterpolatorType::New(); MetricType::Pointermetric = MetricType::New(); RegistrationType::Pointerregistrator = RegistrationType::New(); registrator->SetTransform( transform ); registrator->SetOptimizer( optimizer ); registrator->SetInterpolator( interpolator ); registrator->SetMetric( metric ); registrator->SetFixedImage( fixedImageReader->GetOutput() ); registrator->SetMovingImage( movingImageReader->GetOutput() );

43. Image Registration registrator->SetFixedImageRegion( fixedImageReader->GetOutput()->GetBufferedRegion() ); typedef RegistrationType::ParametersType ParametersType; transform->SetIdentity(); registrator->SetInitialTransformParameters(transform->GetParameters() ); optimizer->SetMaximumStepLength( 4.00 ); optimizer->SetMinimumStepLength( 0.01 ); optimizer->SetNumberOfIterations( 100 ); optimizer->MaximizeOff();

44. Image Registration try { registrator->StartRegistration (); } catch(itk::ExceptionObject & excp ) { std::cerr << “Error in registration” << std::endl;std::cerr << excp<< std::endl; } transform->SetParameters( registrator->GetLastTransformParameters() );

45. Image Registration typedef itk::ResampleImageFilter<FixedImageType , MovingImageType> ResamplerType; ResamplerType::Pointerresampler = ResamplerType::New(); resampler->SetTransform ( transform ); resampler->SetInput( movingImageReader->GetOutput()); FixedImageType::PointerfixedImage = fixedImageReader->GetOutput(); resampler->SetOrigin( fixedImage->GetOrigin()); resampler->SetSpacing( fixedImage->GetSpacing()); resampler->SetSize( fixedImage->GetLargestPossibleRegion()->GetSize()); resampler->Update();

46. Tracking theRegistration Process

47. Observing Registration #include ”itkCommand.h” class CommandIteration : public itk::Command { public:typedef CommandIterationSelf; typedef itk::CommandSuperClass; typedef itk::SmartPointer< Self>Pointer; itkNewMacro( Self ); protected: CommandIteration() {}; public: typedef itk::RegularStepGradientDescentOptimizer OptimizerType; typedef const OptimizerType * OptimizerPointer;

48. Observing Registration voidExecute( itk::Object *caller, const itk::EventObject &event) {this-> Execute((const itk::Object *)caller, event ); } voidExecute( const itk::Object *caller, const itk::EventObject &event) {OptimizerPointeroptimizer = dynamic_cast<OptimizerPointer>(caller); if( typeid(event ) == typeid(itk::IterationEvent) ) { std::cout << optimizer->GetCurrentIteration() << “ : “; std::cout << optimizer->GetValue() << “ : “; std::cout << optimizer->GetCurrentPosition() << std::endl; } }

49. Image Registration CommandIteration::Pointerobserver = CommandIteration::New(); optimizer->AddObserver( itk::IterationEvent(),observer) try { registrator->StartRegistration (); } catch(itk::ExceptionObject & excp ) { std::cerr << “Error in registration” << std::endl;std::cerr << excp<< std::endl; }

50. Image Similarity Metrics