ITK  4.6.0
Insight Segmentation and Registration Toolkit
RegistrationITKv3/ImageRegistration3.py
1 #==========================================================================
2 #
3 # Copyright Insight Software Consortium
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 from InsightToolkit import *
20 
21 from sys import argv
22 
23 #
24 # Read the fixed and moving images using filenames
25 # from the command line arguments
26 #
27 fixedImageReader = itkImageFileReaderF2_New()
28 movingImageReader = itkImageFileReaderF2_New()
29 
30 fixedImageReader.SetFileName( argv[1] )
31 movingImageReader.SetFileName( argv[2] )
32 
33 fixedImageReader.Update()
34 movingImageReader.Update()
35 
36 fixedImage = fixedImageReader.GetOutput()
37 movingImage = movingImageReader.GetOutput()
38 
39 
40 
41 
42 #
43 # Instantiate the classes for the registration framework
44 #
45 registration = itkImageRegistrationMethodF2F2_New()
46 imageMetric = itkMeanSquaresImageToImageMetricF2F2_New()
47 transform = itkTranslationTransform2_New()
48 optimizer = itkRegularStepGradientDescentOptimizer_New()
49 interpolator = itkLinearInterpolateImageFunctionF2D_New()
50 
51 
52 registration.SetOptimizer( optimizer.GetPointer() )
53 registration.SetTransform( transform.GetPointer() )
54 registration.SetInterpolator( interpolator.GetPointer() )
55 registration.SetMetric( imageMetric.GetPointer() )
56 registration.SetFixedImage( fixedImage )
57 registration.SetMovingImage( movingImage )
58 
59 registration.SetFixedImageRegion( fixedImage.GetBufferedRegion() )
60 
61 transform.SetIdentity()
62 initialParameters = transform.GetParameters()
63 
64 registration.SetInitialTransformParameters( initialParameters )
65 
66 #
67 # Iteration Observer
68 #
69 def iterationUpdate():
70  currentParameter = transform.GetParameters()
71  print "M: %f P: %f %f " % ( optimizer.GetValue(),
72  currentParameter.GetElement(0),
73  currentParameter.GetElement(1) )
74 
75 iterationCommand = itkPyCommand_New()
76 iterationCommand.SetCommandCallable( iterationUpdate )
77 optimizer.AddObserver( itkIterationEvent(), iterationCommand.GetPointer() )
78 
79 
80 
81 #
82 # Define optimizer parameters
83 #
84 optimizer.SetMaximumStepLength( 4.00 )
85 optimizer.SetMinimumStepLength( 0.01 )
86 optimizer.SetNumberOfIterations( 200 )
87 
88 
89 print "Starting registration"
90 
91 #
92 # Start the registration process
93 #
94 
95 registration.Update()
96 
97 
98 #
99 # Get the final parameters of the transformation
100 #
101 finalParameters = registration.GetLastTransformParameters()
102 
103 print "Final Registration Parameters "
104 print "Translation X = %f" % (finalParameters.GetElement(0),)
105 print "Translation Y = %f" % (finalParameters.GetElement(1),)
106 
107 
108 
109 
110 #
111 # Now, we use the final transform for resampling the
112 # moving image.
113 #
114 resampler = itkResampleImageFilterF2F2_New()
115 resampler.SetTransform( transform.GetPointer() )
116 resampler.SetInput( movingImage )
117 
118 region = fixedImage.GetLargestPossibleRegion()
119 
120 resampler.SetSize( region.GetSize() )
121 
122 resampler.SetOutputSpacing( fixedImage.GetSpacing() )
123 resampler.SetOutputOrigin( fixedImage.GetOrigin() )
124 resampler.SetOutputDirection( fixedImage.GetDirection() )
125 resampler.SetDefaultPixelValue( 100 )
126 
127 outputCast = itkRescaleIntensityImageFilterF2US2_New()
128 outputCast.SetOutputMinimum( 0 )
129 outputCast.SetOutputMaximum( 65535 )
130 outputCast.SetInput(resampler.GetOutput())
131 
132 #
133 # Write the resampled image
134 #
135 writer = itkImageFileWriterUS2_New()
136 
137 writer.SetFileName( argv[3] )
138 writer.SetInput( outputCast.GetOutput() )
139 writer.Update()