1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20""" A SimpleITK example demonstrating image registration with histogram
21 mutual information as a similarity measure. """
22
23import sys
24import os
25import SimpleITK as sitk
26
27
28def command_iteration(method):
29 """ Callback invoked when the optimization process is running. """
30 print(
31 f"{method.GetOptimizerIteration():3} "
32 + f" = {method.GetMetricValue():7.5f} "
33 + f" : {method.GetOptimizerPosition()}"
34 )
35
36
37def main(args):
38 """ A SimpleITK example demonstrating image registration with histogram
39 mutual information as a similarity measure. """
40 if len(args) < 3:
41 print(
42 "Usage:",
43 "ImageRegistrationMethod2",
44 "<fixedImageFilter> <movingImageFile> <outputTransformFile>",
45 )
46 sys.exit(1)
47
51
55
57
58 R.SetMetricAsJointHistogramMutualInformation()
59
60 R.SetOptimizerAsGradientDescentLineSearch(
61 learningRate=1.0,
62 numberOfIterations=200,
63 convergenceMinimumValue=1e-5,
64 convergenceWindowSize=5,
65 )
66
68
69 R.SetInterpolator(sitk.sitkLinear)
70
71 R.AddCommand(sitk.sitkIterationEvent, lambda: command_iteration(R))
72
73 outTx = R.Execute(fixed, moving)
74
75 print("-------")
76 print(outTx)
77 print(f"Optimizer stop condition: {R.GetOptimizerStopConditionDescription()}")
78 print(f" Iteration: {R.GetOptimizerIteration()}")
79 print(f" Metric value: {R.GetMetricValue()}")
80
82
84 resampler.SetReferenceImage(fixed)
85 resampler.SetInterpolator(sitk.sitkLinear)
86 resampler.SetDefaultPixelValue(1)
87 resampler.SetTransform(outTx)
88
89 out = resampler.Execute(moving)
90
93 cimg =
sitk.Compose(simg1, simg2, simg1 // 2.0 + simg2 // 2.0)
94
95 return {"fixed": fixed, "moving": moving, "composition": cimg}
96
97
98if __name__ == "__main__":
99 return_dict = main(sys.argv)
100 if "SITK_NOSHOW" not in os.environ:
101 sitk.Show(return_dict[
"composition"],
"ImageRegistration2 Composition")
An interface method to the modular ITKv4 registration framework.
Resample an image via a coordinate transform.
SITKBasicFilters_EXPORT Image DiscreteGaussian(const Image &image1, double variance, unsigned int maximumKernelWidth=32u, double maximumError=0.01, bool useImageSpacing=true)
Blurs an image by separable convolution with discrete gaussian kernels. This filter performs Gaussian...
SITKIO_EXPORT Image ReadImage(const PathType &filename, PixelIDValueEnum outputPixelType=sitkUnknown, const std::string &imageIO="")
ReadImage is a procedural interface to the ImageFileReader class which is convenient for most image r...
Image Normalize(const Image &image1)
Normalize an image by setting its mean to zero and variance to one.
Image Compose(const std::vector< Image > &images)
ComposeImageFilter combine several scalar images into a multicomponent image.
void SITKIO_EXPORT Show(const Image &image, const std::string &title="", const bool debugOn=ProcessObject::GetGlobalDefaultDebug())
SITKCommon_EXPORT void WriteTransform(const Transform &transform, const PathType &filename)
Image RescaleIntensity(Image &&image1, double outputMinimum=0, double outputMaximum=255)
Applies a linear transformation to the intensity levels of the input Image .
Image Cast(const Image &image, PixelIDValueEnum pixelID)