ITK  6.0.0
Insight Toolkit
itkCurvatureRegistrationFilter.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 itkCurvatureRegistrationFilter_h
19#define itkCurvatureRegistrationFilter_h
20
23
24#if !defined(ITK_USE_CUFFTW) && (defined(ITK_USE_FFTWF) || defined(ITK_USE_FFTWD))
25# include "fftw3.h"
26
27namespace itk
28{
97template <typename TFixedImage,
98 typename TMovingImage,
99 typename TDisplacementField,
100 typename TImageForceFunction = MeanSquareRegistrationFunction<TFixedImage, TMovingImage, TDisplacementField>>
101class ITK_TEMPLATE_EXPORT CurvatureRegistrationFilter
102 : public PDEDeformableRegistrationFilter<TFixedImage, TMovingImage, TDisplacementField>
103{
104public:
105 ITK_DISALLOW_COPY_AND_MOVE(CurvatureRegistrationFilter);
106
112
114 itkNewMacro(Self);
115
117 itkOverrideGetNameOfClassMacro(CurvatureRegistrationFilter);
118
120 using typename Superclass::TimeStepType;
121
123 using typename Superclass::FixedImageType;
124 using typename Superclass::FixedImagePointer;
125 static constexpr unsigned int ImageDimension = FixedImageType::ImageDimension;
126
128 using typename Superclass::MovingImageType;
129 using typename Superclass::MovingImagePointer;
130
132 using typename Superclass::DisplacementFieldType;
133 using typename Superclass::DisplacementFieldPointer;
134
135 using DisplacementFieldPixelType = typename TDisplacementField::PixelType;
136 using DisplacementFieldComponentType = typename DisplacementFieldPixelType::ValueType;
137 static constexpr unsigned int DeformationVectorDimension = DisplacementFieldPixelType::Dimension;
138
139# if defined(ITK_USE_FFTWD)
140 // Prefer to use double precision
141 using RealTypeDFT = double;
142# else
143# if defined(ITK_USE_FFTWF)
144 // Allow to use single precision
145# warning "Using single precision for FFT computations!"
146 using RealTypeDFT = double;
147# endif
148# endif
149
152
154 using typename Superclass::FiniteDifferenceFunctionType;
155
157 using RegistrationFunctionType = TImageForceFunction;
158
160 void
161 SetConstraintWeight(const float w)
162 {
163 m_ConstraintWeight = w;
164 }
165
167 void
169 {
170 m_TimeStep = ts;
171 }
172
178 virtual double
179 GetMetric() const;
180
181protected:
184 void
185 PrintSelf(std::ostream & os, Indent indent) const override;
186
188 void
189 Initialize() override;
190
192 void
193 ApplyUpdate(const TimeStepType & dt) override;
194
195private:
196 unsigned int m_FixedImageDimensions[ImageDimension]{};
197
198 RealTypeDFT * m_DisplacementFieldComponentImage{};
199 RealTypeDFT * m_DisplacementFieldComponentImageDCT{};
200
201 float m_ConstraintWeight{};
202
203 fftw_plan m_PlanForwardDCT{};
204 fftw_plan m_PlanBackwardDCT{};
205
206 TimeStepType m_TimeStep{};
207
208 RealTypeDFT * m_DiagonalElements[ImageDimension]{};
209};
210} // end namespace itk
211
212# ifndef ITK_MANUAL_INSTANTIATION
213# include "itkCurvatureRegistrationFilter.hxx"
214# endif
215
216#endif // defined(ITK_USE_FFTWF) || defined(ITK_USE_FFTWD)
217
218#endif
Deformably register two images using the fast curvature algorithm.
void PrintSelf(std::ostream &os, Indent indent) const override
typename DisplacementFieldPixelType::ValueType DisplacementFieldComponentType
typename DisplacementFieldComponentImageType::Pointer DisplacementFieldComponentImagePointer
virtual double GetMetric() const
typename TDisplacementField::PixelType DisplacementFieldPixelType
void ApplyUpdate(const TimeStepType &dt) override
Templated n-dimensional image class.
Definition: itkImage.h:89
Control indentation during Print() invocation.
Definition: itkIndent.h:50
Light weight base class for most itk classes.
Deformably register two images using a PDE algorithm.
SmartPointer< Self > Pointer
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....