ITK  6.0.0
Insight Toolkit
itkInvertDisplacementFieldImageFilter.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 itkInvertDisplacementFieldImageFilter_h
19#define itkInvertDisplacementFieldImageFilter_h
20
24#include <mutex>
25
26namespace itk
27{
28
40template <typename TInputImage, typename TOutputImage = TInputImage>
41class ITK_TEMPLATE_EXPORT InvertDisplacementFieldImageFilter : public ImageToImageFilter<TInputImage, TOutputImage>
42{
43public:
44 ITK_DISALLOW_COPY_AND_MOVE(InvertDisplacementFieldImageFilter);
45
50
52 itkNewMacro(Self);
53
55 itkOverrideGetNameOfClassMacro(InvertDisplacementFieldImageFilter);
56
58 static constexpr unsigned int ImageDimension = TInputImage::ImageDimension;
59
60 using InputFieldType = TInputImage;
61 using OutputFieldType = TOutputImage;
62
65
67 using PixelType = typename OutputFieldType::PixelType;
68 using VectorType = typename OutputFieldType::PixelType;
72
74 using SpacingType = typename OutputFieldType::SpacingType;
78
80 using RealType = typename VectorType::ComponentType;
84
86 itkGetModifiableObjectMacro(Interpolator, InterpolatorType);
87
89 void
91 {
92 itkDebugMacro("setting deformation field to " << field);
93 if (field != this->GetInput(0))
94 {
95 this->SetInput(0, field);
96 this->Modified();
97 if (!this->m_Interpolator.IsNull())
98 {
99 this->m_Interpolator->SetInputImage(field);
100 }
101 }
102 }
108 const InputFieldType *
110 {
111 return this->GetInput(0);
112 }
113
115 itkSetInputMacro(InverseFieldInitialEstimate, InverseDisplacementFieldType);
116 itkGetInputMacro(InverseFieldInitialEstimate, InverseDisplacementFieldType);
119 /* Set the interpolator. */
120 virtual void
122
123 /* Set/Get the number of iterations */
124 itkSetMacro(MaximumNumberOfIterations, unsigned int);
125 itkGetConstMacro(MaximumNumberOfIterations, unsigned int);
126
127 /* Set/Get the mean stopping criterion */
128 itkSetMacro(MeanErrorToleranceThreshold, RealType);
129 itkGetConstMacro(MeanErrorToleranceThreshold, RealType);
130
131 /* Set/Get the max stopping criterion */
132 itkSetMacro(MaxErrorToleranceThreshold, RealType);
133 itkGetConstMacro(MaxErrorToleranceThreshold, RealType);
134
135 /* Get the max norm */
136 itkGetConstMacro(MaxErrorNorm, RealType);
137
138 /* Get the mean norm */
139 itkGetConstMacro(MeanErrorNorm, RealType);
140
141 /* Should we force the boundary to have zero displacement? */
142 itkSetMacro(EnforceBoundaryCondition, bool);
143 itkGetMacro(EnforceBoundaryCondition, bool);
144 itkBooleanMacro(EnforceBoundaryCondition);
145
146protected:
149
152
154 void
155 PrintSelf(std::ostream & os, Indent indent) const override;
156
158 void
159 GenerateData() override;
160
162 void
164
165private:
167 typename InterpolatorType::Pointer m_Interpolator{};
168
169 unsigned int m_MaximumNumberOfIterations{ 20 };
170
171 RealType m_MaxErrorToleranceThreshold{};
172 RealType m_MeanErrorToleranceThreshold{};
173
174 // internal ivars necessary for multithreading basic operations
175
176 typename DisplacementFieldType::Pointer m_ComposedField{};
177 typename RealImageType::Pointer m_ScaledNormImage{};
178
179 RealType m_MaxErrorNorm{};
180 RealType m_MeanErrorNorm{};
181 RealType m_Epsilon{};
182 SpacingType m_DisplacementFieldSpacing{};
183 bool m_DoThreadedEstimateInverse{ false };
184 bool m_EnforceBoundaryCondition{ true };
185 std::mutex m_Mutex{};
186};
187
188} // end namespace itk
189
190#ifndef ITK_MANUAL_INSTANTIATION
191# include "itkInvertDisplacementFieldImageFilter.hxx"
192#endif
193
194#endif
Base class for filters that take an image as input and produce an image as output.
Templated n-dimensional image class.
Definition: itkImage.h:89
Control indentation during Print() invocation.
Definition: itkIndent.h:50
Iteratively estimate the inverse field of a displacement field.
void PrintSelf(std::ostream &os, Indent indent) const override
itkSetInputMacro(InverseFieldInitialEstimate, InverseDisplacementFieldType)
virtual void SetInterpolator(InterpolatorType *interpolator)
itkGetInputMacro(InverseFieldInitialEstimate, InverseDisplacementFieldType)
typename OutputFieldType::DirectionType DirectionType
~InvertDisplacementFieldImageFilter() override=default
void DynamicThreadedGenerateData(const RegionType &) override
Light weight base class for most itk classes.
Base class for all vector image interpolators.
Linearly interpolate a vector image at specified positions.
SmartPointer< Self > Pointer
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....