ITK  5.2.0
Insight Toolkit
Examples/RegistrationITKv4/ImageRegistration3.py
1 # ==========================================================================
2 #
3 # Copyright NumFOCUS
4 #
5 # Licensed under the Apache License, Version 2.0 (the "License");
6 # you may not use this file except in compliance with the License.
7 # You may obtain a copy of the License at
8 #
9 # http://www.apache.org/licenses/LICENSE-2.0.txt
10 #
11 # Unless required by applicable law or agreed to in writing, software
12 # distributed under the License is distributed on an "AS IS" BASIS,
13 # WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 # See the License for the specific language governing permissions and
15 # limitations under the License.
16 #
17 # ==========================================================================*/
18 
19 import itk
20 from sys import argv
21 
22 
23 #
24 # Check input parameters
25 # INPUTS(fixedImage): {BrainProtonDensitySliceBorder20.png}
26 # INPUTS(movingImage): {BrainProtonDensitySliceShifted13x17y.png}
27 #
28 if len(argv) < 4:
29  print("Missing Parameters")
30  print(
31  "Usage: ImageRegistration3.py fixedImageFile movingImageFile outputImagefile"
32  )
33  exit()
34 
35 
36 #
37 # Define data types
38 #
39 FixedImageType = itk.Image[itk.F, 2]
40 MovingImageType = itk.Image[itk.F, 2]
41 TransformType = itk.TranslationTransform[itk.D, 2]
42 OptimizerType = itk.RegularStepGradientDescentOptimizerv4[itk.D]
43 RegistrationType = itk.ImageRegistrationMethodv4[FixedImageType, MovingImageType]
44 MetricType = itk.MeanSquaresImageToImageMetricv4[FixedImageType, MovingImageType]
45 
46 
47 #
48 # Read the fixed and moving images using filenames
49 # from the command line arguments
50 #
51 fixedImageReader = itk.ImageFileReader[FixedImageType].New()
52 movingImageReader = itk.ImageFileReader[MovingImageType].New()
53 
54 fixedImageReader.SetFileName(argv[1])
55 movingImageReader.SetFileName(argv[2])
56 
57 fixedImageReader.Update()
58 movingImageReader.Update()
59 
60 fixedImage = fixedImageReader.GetOutput()
61 movingImage = movingImageReader.GetOutput()
62 
63 
64 #
65 # Instantiate the classes for the registration framework
66 #
67 registration = RegistrationType.New()
68 imageMetric = MetricType.New()
69 transform = TransformType.New()
70 optimizer = OptimizerType.New()
71 
72 registration.SetOptimizer(optimizer)
73 registration.SetMetric(imageMetric)
74 
75 registration.SetFixedImage(fixedImage)
76 registration.SetMovingImage(movingImage)
77 
78 registration.SetInitialTransform(transform)
79 
80 
81 #
82 # Define optimizer parameters
83 #
84 optimizer.SetLearningRate(4)
85 optimizer.SetMinimumStepLength(0.001)
86 optimizer.SetRelaxationFactor(0.5)
87 optimizer.SetNumberOfIterations(100)
88 
89 
90 #
91 # One level registration process without shrinking and smoothing.
92 #
93 registration.SetNumberOfLevels(1)
94 registration.SetSmoothingSigmasPerLevel([0])
95 registration.SetShrinkFactorsPerLevel([1])
96 
97 
98 #
99 # Iteration Observer
100 #
101 def iterationUpdate():
102  currentParameter = registration.GetOutput().Get().GetParameters()
103  print(
104  "M: %f P: %f %f "
105  % (
106  optimizer.GetValue(),
107  currentParameter.GetElement(0),
108  currentParameter.GetElement(1),
109  )
110  )
111 
112 
113 iterationCommand = itk.PyCommand.New()
114 iterationCommand.SetCommandCallable(iterationUpdate)
115 optimizer.AddObserver(itk.IterationEvent(), iterationCommand)
116 
117 print("Starting registration")
118 
119 
120 #
121 # Start the registration process
122 #
123 registration.Update()
124 
125 
126 #
127 # Get the final parameters of the transformation
128 #
129 finalParameters = registration.GetOutput().Get().GetParameters()
130 
131 print("Final Registration Parameters ")
132 print(f"Translation X = {finalParameters.GetElement(0):f}")
133 print(f"Translation Y = {finalParameters.GetElement(1):f}")
134 
135 
136 #
137 # Now, we use the final transform for resampling the
138 # moving image.
139 #
140 resampler = itk.ResampleImageFilter[MovingImageType, FixedImageType].New()
141 resampler.SetTransform(registration.GetTransform())
142 resampler.SetInput(movingImageReader.GetOutput())
143 
144 region = fixedImage.GetLargestPossibleRegion()
145 
146 resampler.SetSize(region.GetSize())
147 resampler.SetOutputOrigin(fixedImage.GetOrigin())
148 resampler.SetOutputSpacing(fixedImage.GetSpacing())
149 resampler.SetOutputDirection(fixedImage.GetDirection())
150 resampler.SetDefaultPixelValue(100)
151 
152 OutputImageType = itk.Image[itk.UC, 2]
153 outputCast = itk.CastImageFilter[FixedImageType, OutputImageType].New()
154 outputCast.SetInput(resampler.GetOutput())
155 
156 
157 #
158 # Write the resampled image
159 #
160 writer = itk.ImageFileWriter[OutputImageType].New()
161 writer.SetFileName(argv[3])
162 writer.SetInput(outputCast.GetOutput())
163 writer.Update()
itk::CastImageFilter
Casts input pixels to output pixel type.
Definition: itkCastImageFilter.h:104
itk::ImageFileReader
Data source that reads image data from a single file.
Definition: itkImageFileReader.h:75
itk::RegularStepGradientDescentOptimizerv4
Regular Step Gradient descent optimizer.
Definition: itkRegularStepGradientDescentOptimizerv4.h:47
itk::ImageFileWriter
Writes image data to a single file.
Definition: itkImageFileWriter.h:87
itk::TranslationTransform
Translation transformation of a vector space (e.g. space coordinates)
Definition: itkTranslationTransform.h:43
itk::ImageRegistrationMethodv4
Interface method for the current registration framework.
Definition: itkImageRegistrationMethodv4.h:117
itk::MeanSquaresImageToImageMetricv4
Class implementing a mean squares metric.
Definition: itkMeanSquaresImageToImageMetricv4.h:46
itk::ResampleImageFilter
Resample an image via a coordinate transform.
Definition: itkResampleImageFilter.h:90
itk::Image
Templated n-dimensional image class.
Definition: itkImage.h:86