ITK 6.0.0
Insight Toolkit
 
Loading...
Searching...
No Matches
itkImageRegistrationMethodImageSource.h
Go to the documentation of this file.
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#ifndef itkImageRegistrationMethodImageSource_h
19#define itkImageRegistrationMethodImageSource_h
20#include "itkImage.h"
24
34
36{
37
38template <typename TFixedPixelType, typename TMovingPixelType, unsigned int VDimension>
40{
41public:
47
49 itkNewMacro(Self);
50
52 itkOverrideGetNameOfClassMacro(ImageRegistrationMethodImageSource);
53
54
57
58 const MovingImageType *
60 {
61 return m_MovingImage.GetPointer();
62 }
63
64 const FixedImageType *
66 {
67 return m_FixedImage.GetPointer();
68 }
69
70 const ParametersType &
72 {
73 return m_Parameters;
74 }
75
76
77 void
79 {
80 const typename MovingImageType::RegionType region(size);
81
82 m_MovingImage->SetRegions(region);
83 m_MovingImage->Allocate();
84
85 m_FixedImage->SetRegions(region);
86 m_FixedImage->Allocate();
87
88 /* Fill images with a 2D gaussian*/
89 using MovingImageIteratorType = itk::ImageRegionIteratorWithIndex<MovingImageType>;
90
91 using FixedImageIteratorType = itk::ImageRegionIteratorWithIndex<FixedImageType>;
92
93
95 center[0] = static_cast<double>(region.GetSize()[0]) / 2.0;
96 center[1] = static_cast<double>(region.GetSize()[1]) / 2.0;
97
98 const double s = static_cast<double>(region.GetSize()[0]) / 2.0;
99
102
103 /* Set the displacement */
104 itk::Vector<double, 2> displacement;
105 displacement[0] = m_Parameters[0];
106 displacement[1] = m_Parameters[1];
107
108
109 MovingImageIteratorType ri(m_MovingImage, region);
110 FixedImageIteratorType ti(m_FixedImage, region);
111 while (!ri.IsAtEnd())
112 {
113 p[0] = ri.GetIndex()[0];
114 p[1] = ri.GetIndex()[1];
115 d = p - center;
116 d += displacement;
117 const double x = d[0];
118 const double y = d[1];
119 const double value = 200.0 * std::exp(-(x * x + y * y) / (s * s));
120 ri.Set(static_cast<typename MovingImageType::PixelType>(value));
121 ++ri;
122 }
123
124
125 while (!ti.IsAtEnd())
126 {
127 p[0] = ti.GetIndex()[0];
128 p[1] = ti.GetIndex()[1];
129 d = p - center;
130 const double x = d[0];
131 const double y = d[1];
132 const double value = 200.0 * std::exp(-(x * x + y * y) / (s * s));
133 ti.Set(static_cast<typename FixedImageType::PixelType>(value));
134 ++ti;
135 }
136 }
137
138protected:
147
148private:
151
153};
154
155} // namespace itk::testhelper
156
157// end namespace itk
158#endif
A multi-dimensional iterator templated over image type that walks pixels within a region and is speci...
const SizeType & GetSize() const
Templated n-dimensional image class.
Definition itkImage.h:89
ImageRegion< VImageDimension > RegionType
Base class for most ITK classes.
Definition itkObject.h:62
Class to hold and manage different parameter types used during optimization.
A templated class holding a geometric point in n-Dimensional space.
Definition itkPoint.h:54
Implements transparent reference counting.
A templated class holding a n-Dimensional vector.
Definition itkVector.h:63
void GenerateImages(const typename MovingImageType::SizeType &size)