ITK  5.4.0
Insight Toolkit
Examples/RegistrationITKv4/ImageRegistration4.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# https://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
19import itk
20from sys import argv
21
22
23#
24# Check input parameters
25# INPUTS(fixedImage): {BrainProtonDensitySliceBorder20.png}
26# INPUTS(movingImage): {BrainProtonDensitySliceShifted13x17y.png}
27#
28if len(argv) < 4:
29 print("Missing Parameters")
30 print(
31 "Usage: ImageRegistration4.py fixedImageFile movingImageFile outputImagefile"
32 )
33 exit()
34
35
36#
37# Define data types
38#
39FixedImageType = itk.Image[itk.F, 2]
40MovingImageType = itk.Image[itk.F, 2]
41TransformType = itk.TranslationTransform[itk.D, 2]
43RegistrationType = itk.ImageRegistrationMethodv4[FixedImageType, MovingImageType]
45 FixedImageType, MovingImageType
46]
47
48
49#
50# Read the fixed and moving images using filenames
51# from the command line arguments
52#
53fixedImageReader = itk.ImageFileReader[FixedImageType].New()
54movingImageReader = itk.ImageFileReader[MovingImageType].New()
55
56fixedImageReader.SetFileName(argv[1])
57movingImageReader.SetFileName(argv[2])
58
59fixedImageReader.Update()
60movingImageReader.Update()
61
62fixedImage = fixedImageReader.GetOutput()
63movingImage = movingImageReader.GetOutput()
64
65
66#
67# Instantiate the classes for the registration framework
68#
69registration = RegistrationType.New()
70imageMetric = MetricType.New()
71transform = TransformType.New()
72optimizer = OptimizerType.New()
73
74registration.SetOptimizer(optimizer)
75registration.SetMetric(imageMetric)
76
77numberOfBins = 24
78
79imageMetric.SetNumberOfHistogramBins(numberOfBins)
80imageMetric.SetUseMovingImageGradientFilter(False)
81imageMetric.SetUseFixedImageGradientFilter(False)
82
83registration.SetFixedImage(fixedImage)
84registration.SetMovingImage(movingImage)
85
86registration.SetInitialTransform(transform)
87
88
89#
90# Define optimizer parameters
91#
92optimizer.SetLearningRate(8.00)
93optimizer.SetMinimumStepLength(0.001)
94optimizer.SetNumberOfIterations(100)
95optimizer.ReturnBestParametersAndValueOn()
96optimizer.SetRelaxationFactor(0.8)
97
98
99#
100# One level registration process without shrinking and smoothing.
101#
102registration.SetNumberOfLevels(1)
103registration.SetSmoothingSigmasPerLevel([0])
104registration.SetShrinkFactorsPerLevel([1])
105
106registration.SetMetricSamplingStrategy(RegistrationType.RANDOM)
107registration.SetMetricSamplingPercentage(0.20)
108
109
110#
111# Iteration Observer
112#
113def iterationUpdate():
114 currentParameter = registration.GetOutput().Get().GetParameters()
115 print(
116 "M: %f P: %f %f "
117 % (
118 optimizer.GetValue(),
119 currentParameter.GetElement(0),
120 currentParameter.GetElement(1),
121 )
122 )
123
124
125iterationCommand = itk.PyCommand.New()
126iterationCommand.SetCommandCallable(iterationUpdate)
127optimizer.AddObserver(itk.IterationEvent(), iterationCommand)
128
129print("Starting registration")
130
131
132#
133# Start the registration process
134#
135registration.Update()
136
137
138#
139# Get the final parameters of the transformation
140#
141finalParameters = registration.GetOutput().Get().GetParameters()
142
143print("Final Registration Parameters ")
144print(f"Translation X = {finalParameters.GetElement(0):f}")
145print(f"Translation Y = {finalParameters.GetElement(1):f}")
146
147
148#
149# Now, we use the final transform for resampling the
150# moving image.
151#
152resampler = itk.ResampleImageFilter[MovingImageType, FixedImageType].New()
153resampler.SetTransform(registration.GetTransform())
154resampler.SetInput(movingImageReader.GetOutput())
155
156region = fixedImage.GetLargestPossibleRegion()
157
158resampler.SetSize(region.GetSize())
159resampler.SetOutputOrigin(fixedImage.GetOrigin())
160resampler.SetOutputSpacing(fixedImage.GetSpacing())
161resampler.SetOutputDirection(fixedImage.GetDirection())
162resampler.SetDefaultPixelValue(100)
163
164OutputImageType = itk.Image[itk.UC, 2]
165outputCast = itk.CastImageFilter[FixedImageType, OutputImageType].New()
166outputCast.SetInput(resampler.GetOutput())
167
168
169#
170# Write the resampled image
171#
172writer = itk.ImageFileWriter[OutputImageType].New()
173writer.SetFileName(argv[3])
174writer.SetInput(outputCast.GetOutput())
175writer.Update()
Casts input pixels to output pixel type.
Data source that reads image data from a single file.
Writes image data to a single file.
Interface method for the current registration framework.
Templated n-dimensional image class.
Definition: itkImage.h:89
Computes the mutual information between two images to be registered using the method of Mattes et al.
Resample an image via a coordinate transform.
Translation transformation of a vector space (e.g. space coordinates)
static Pointer New()