ITK  6.0.0
Insight Toolkit
itkDisplacementFieldToBSplineImageFilter.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 itkDisplacementFieldToBSplineImageFilter_h
19#define itkDisplacementFieldToBSplineImageFilter_h
20
22
24#include "itkPointSet.h"
25#include "itkVector.h"
26
27namespace itk
28{
29
41template <typename TInputImage,
42 typename TInputPointSet = PointSet<typename TInputImage::PixelType, TInputImage::ImageDimension>,
43 typename TOutputImage = TInputImage>
44class ITK_TEMPLATE_EXPORT DisplacementFieldToBSplineImageFilter : public ImageToImageFilter<TInputImage, TOutputImage>
45{
46public:
47 ITK_DISALLOW_COPY_AND_MOVE(DisplacementFieldToBSplineImageFilter);
48
53
55 itkNewMacro(Self);
56
58 itkOverrideGetNameOfClassMacro(DisplacementFieldToBSplineImageFilter);
59
61 static constexpr unsigned int ImageDimension = TInputImage::ImageDimension;
62
63 using InputFieldType = TInputImage;
64 using InputPointSetType = TInputPointSet;
65 using OutputFieldType = TOutputImage;
66
70
72 using PixelType = typename OutputFieldType::PixelType;
73 using VectorType = typename OutputFieldType::PixelType;
76
77 using SpacingType = typename OutputFieldType::SpacingType;
81
82 using RealType = typename VectorType::RealValueType;
84
87 using PointDataType = typename InputPointSetType::PixelType;
88 using PointsContainerType = typename InputPointSetType::PointsContainer;
89 using PointDataContainerType = typename InputPointSetType::PointDataContainer;
90
96
98 void
100 {
101 this->SetInput(0, field);
102 }
103
105 const InputFieldType *
107 {
108 return this->GetInput(0);
109 }
110
116 void
118 {
119 this->SetNthInput(1, const_cast<RealImageType *>(image));
120 }
121 void
123 {
124 this->SetConfidenceImage(image);
125 }
129 const RealImageType *
131 {
132 return static_cast<const RealImageType *>(this->ProcessObject::GetInput(1));
133 }
134
136 void
138 {
139 this->SetNthInput(2, const_cast<InputPointSetType *>(points));
140 }
141 void
143 {
144 this->SetPointSet(points);
145 }
149 const InputPointSetType *
151 {
152 return static_cast<const InputPointSetType *>(this->ProcessObject::GetInput(2));
153 }
154
156 void
158
162 {
163 return static_cast<const DisplacementFieldControlPointLatticeType *>(this->GetOutput(1));
164 }
165
167 void
169
171 void
173 {
174 this->SetBSplineDomainFromImage(const_cast<RealImageType *>(image));
175 }
176
178 void
180
182 void
184 {
185 this->SetBSplineDomainFromImage(const_cast<InputFieldType *>(field));
186 }
187
190
191 /* Set/Get b-spline domain origin. */
192 itkGetConstMacro(BSplineDomainOrigin, OriginType);
193
194 /* Set/Get b-spline domain spacing. */
195 itkGetConstMacro(BSplineDomainSpacing, SpacingType);
196
197 /* Set/Get b-spline domain size. */
198 itkGetConstMacro(BSplineDomainSize, SizeType);
199
200 /* Set/Get b-spline domain direction. */
201 itkGetConstMacro(BSplineDomainDirection, DirectionType);
202
203 /* Use input field to define the B-spline domain. */
204 itkSetMacro(UseInputFieldToDefineTheBSplineDomain, bool);
205 itkGetConstMacro(UseInputFieldToDefineTheBSplineDomain, bool);
206 itkBooleanMacro(UseInputFieldToDefineTheBSplineDomain);
207
211 itkSetMacro(SplineOrder, unsigned int);
212
216 itkGetConstMacro(SplineOrder, unsigned int);
217
225 itkSetMacro(NumberOfControlPoints, ArrayType);
226
234 itkGetConstMacro(NumberOfControlPoints, ArrayType);
235
242 itkSetMacro(NumberOfFittingLevels, ArrayType);
243
250 void
252 {
253 auto nlevels = MakeFilled<ArrayType>(n);
254 this->SetNumberOfFittingLevels(nlevels);
255 }
264 itkGetConstMacro(NumberOfFittingLevels, ArrayType);
265
269 itkBooleanMacro(EstimateInverse);
270 itkSetMacro(EstimateInverse, bool);
271 itkGetConstMacro(EstimateInverse, bool);
277 itkBooleanMacro(EnforceStationaryBoundary);
278 itkSetMacro(EnforceStationaryBoundary, bool);
279 itkGetConstMacro(EnforceStationaryBoundary, bool);
282protected:
285
288
290 void
291 PrintSelf(std::ostream & os, Indent indent) const override;
292
294 void
295 GenerateData() override;
296
297private:
298 bool m_EstimateInverse{ false };
299 bool m_EnforceStationaryBoundary{ true };
300 unsigned int m_SplineOrder{ 3 };
301 ArrayType m_NumberOfControlPoints{};
302 ArrayType m_NumberOfFittingLevels{};
303
304 typename WeightsContainerType::Pointer m_PointWeights{};
305 bool m_UsePointWeights{ false };
306
307 OriginType m_BSplineDomainOrigin{};
308 SpacingType m_BSplineDomainSpacing{};
309 SizeType m_BSplineDomainSize{};
310 DirectionType m_BSplineDomainDirection{};
311
312 bool m_BSplineDomainIsDefined{ true };
313 bool m_UseInputFieldToDefineTheBSplineDomain{ false };
314};
315
316} // end namespace itk
317
318#ifndef ITK_MANUAL_INSTANTIATION
319# include "itkDisplacementFieldToBSplineImageFilter.hxx"
320#endif
321
322#endif
Image filter which provides a B-spline output approximation.
Class which takes a dense displacement field image and/or a set of points with associated displacemen...
typename InputPointSetType::PointsContainer PointsContainerType
void SetPointSetConfidenceWeights(WeightsContainerType *weights)
const DisplacementFieldControlPointLatticeType * GetDisplacementFieldControlPointLattice() const
typename BSplineFilterType::PointDataImageType DisplacementFieldControlPointLatticeType
void PrintSelf(std::ostream &os, Indent indent) const override
void SetBSplineDomain(OriginType, SpacingType, SizeType, DirectionType)
typename BSplineFilterType::WeightsContainerType WeightsContainerType
void SetBSplineDomainFromImage(InputFieldType *)
typename InputPointSetType::PointDataContainer PointDataContainerType
void SetBSplineDomainFromImage(RealImageType *)
~DisplacementFieldToBSplineImageFilter() override=default
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
Light weight base class for most itk classes.
DataObject * GetInput(const DataObjectIdentifierType &key)
Return an input.
Define a front-end to the STL "vector" container that conforms to the IndexedContainerInterface.
SmartPointer< Self > Pointer
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....