ITK  6.0.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 # 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 
19 import sys
20 
21 import itk
22 
23 #
24 # Check input parameters
25 # INPUTS(fixedImage): {BrainProtonDensitySliceBorder20.png}
26 # INPUTS(movingImage): {BrainProtonDensitySliceShifted13x17y.png}
27 #
28 if len(sys.argv) < 4:
29  print("Missing Parameters")
30  print(
31  "Usage: ImageRegistration3.py fixed_image_file moving_image_file output_image_file"
32  )
33  sys.exit(1)
34 fixed_image_file = sys.argv[1]
35 moving_image_file = sys.argv[2]
36 output_image_file = sys.argv[3]
37 
38 #
39 # Define data types
40 #
41 PixelType = itk.ctype("float")
42 
43 #
44 # Read the fixed and moving images using filenames
45 # from the command line arguments
46 #
47 fixed_image = itk.imread(fixed_image_file, PixelType)
48 moving_image = itk.imread(moving_image_file, PixelType)
49 Dimension = fixed_image.GetImageDimension()
50 
51 #
52 # Create the spatial transform we will optimize.
53 #
54 initial_transform = itk.TranslationTransform[itk.D, Dimension].New()
55 
56 #
57 # Create the matching metric we will use to compare the images during
58 # optimization.
59 matching_metric = itk.MeanSquaresImageToImageMetricv4[
60  type(fixed_image), type(moving_image)
61 ].New()
62 
63 #
64 # Create the optimizer we will use for optimization.
65 #
67 optimizer.SetLearningRate(4)
68 optimizer.SetMinimumStepLength(0.001)
69 optimizer.SetRelaxationFactor(0.5)
70 optimizer.SetNumberOfIterations(100)
71 #
72 # Iteration Observer
73 #
74 def iteration_update():
75  metric_value = optimizer.GetValue()
76  current_parameters = optimizer.GetCurrentPosition()
77  parameter_list = [current_parameters[i] for i in range(len(current_parameters))]
78  print(f"Metric: {metric_value:.8g} \tParameters: {parameter_list}")
79 
80 
81 iteration_command = itk.PyCommand.New()
82 iteration_command.SetCommandCallable(iteration_update)
83 optimizer.AddObserver(itk.IterationEvent(), iteration_command)
84 
85 #
86 # One level registration process without shrinking and smoothing.
87 #
89  fixed_image=fixed_image,
90  moving_image=moving_image,
91  optimizer=optimizer,
92  metric=matching_metric,
93  initial_transform=initial_transform,
94 )
95 registration.SetNumberOfLevels(1)
96 registration.SetSmoothingSigmasPerLevel([0])
97 registration.SetShrinkFactorsPerLevel([1])
98 
99 #
100 # Execute the registration
101 #
102 registration.Update()
103 
104 
105 #
106 # Get the final parameters of the transformation
107 #
108 final_parameters = registration.GetOutput().Get().GetParameters()
109 
110 print("\nFinal Registration Parameters: ")
111 print(f"Translation X = {final_parameters[0]:f}")
112 print(f"Translation Y = {final_parameters[1]:f}")
113 
114 
115 #
116 # Now, we use the final transform for resampling the
117 # moving image.
118 #
119 resampled_moving = itk.resample_image_filter(
120  moving_image,
121  transform=registration.GetTransform(),
122  use_reference_image=True,
123  reference_image=fixed_image,
124  default_pixel_value=100,
125 )
126 
127 #
128 # Cast to 8-bit unsigned integer pixels, supported by the PNG file format
129 #
130 OutputPixelType = itk.ctype("unsigned char")
131 resampled_cast = resampled_moving.astype(OutputPixelType)
132 
133 #
134 # Write the resampled image
135 #
136 itk.imwrite(resampled_cast, output_image_file)
itk::ImageRegistrationMethodv4::New
static Pointer New()
itk::RegularStepGradientDescentOptimizerv4
Regular Step Gradient descent optimizer.
Definition: itkRegularStepGradientDescentOptimizerv4.h:47
itk::TranslationTransform
Translation transformation of a vector space (e.g. space coordinates)
Definition: itkTranslationTransform.h:43
itk::MeanSquaresImageToImageMetricv4
Class implementing a mean squares metric.
Definition: itkMeanSquaresImageToImageMetricv4.h:46
New
static Pointer New()