ITK  6.0.0
Insight Toolkit
itkBSplineControlPointImageFilter.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 itkBSplineControlPointImageFilter_h
19#define itkBSplineControlPointImageFilter_h
20
24#include "itkFixedArray.h"
25#include "itkPointSet.h"
27#include "itkVector.h"
28#include "itkVectorContainer.h"
29
30#include "vnl/vnl_matrix.h"
31#include "vnl/vnl_vector.h"
32
33namespace itk
34{
35
58template <typename TInputImage, typename TOutputImage = TInputImage>
59class ITK_TEMPLATE_EXPORT BSplineControlPointImageFilter : public ImageToImageFilter<TInputImage, TOutputImage>
60{
61public:
62 ITK_DISALLOW_COPY_AND_MOVE(BSplineControlPointImageFilter);
63
68
70 itkNewMacro(Self);
71
73 itkOverrideGetNameOfClassMacro(BSplineControlPointImageFilter);
74
76 static constexpr unsigned int ImageDimension = TInputImage::ImageDimension;
77
78 using ControlPointLatticeType = TInputImage;
79 using OutputImageType = TOutputImage;
80
82 using PixelType = typename OutputImageType::PixelType;
87
88 using SpacingType = typename OutputImageType::SpacingType;
92
94 using RealType = float;
97
100
107
114
119 void
120 SetSplineOrder(unsigned int);
121
127
131 itkGetConstReferenceMacro(SplineOrder, ArrayType);
132
149 itkSetMacro(CloseDimension, ArrayType);
150 itkGetConstReferenceMacro(CloseDimension, ArrayType);
156 itkSetMacro(Spacing, SpacingType);
157 itkGetConstMacro(Spacing, SpacingType);
163 itkSetMacro(Origin, OriginType);
164 itkGetConstMacro(Origin, OriginType);
170 itkSetMacro(Size, SizeType);
171 itkGetConstMacro(Size, SizeType);
188 itkSetMacro(Direction, DirectionType);
189
193 itkGetConstMacro(Direction, DirectionType);
194
203
204protected:
207 void
208 PrintSelf(std::ostream & os, Indent indent) const override;
209
211 void
213
214
215private:
220 void
222
227 unsigned int
228 SplitRequestedRegion(unsigned int, unsigned int, OutputImageRegionType &) override;
229
234 void
236
241
243 SizeType m_Size{};
244 SpacingType m_Spacing{};
245 OriginType m_Origin{};
246 DirectionType m_Direction{};
247
248 bool m_DoMultilevel{ false };
249 unsigned int m_MaximumNumberOfLevels{ 1 };
250 ArrayType m_NumberOfControlPoints{};
251 ArrayType m_CloseDimension{};
252 ArrayType m_SplineOrder{};
253 ArrayType m_NumberOfLevels{};
254
255 vnl_matrix<RealType> m_RefinedLatticeCoefficients[ImageDimension]{};
256
257 typename KernelType::Pointer m_Kernel[ImageDimension]{};
258 typename KernelOrder0Type::Pointer m_KernelOrder0{};
259 typename KernelOrder1Type::Pointer m_KernelOrder1{};
260 typename KernelOrder2Type::Pointer m_KernelOrder2{};
261 typename KernelOrder3Type::Pointer m_KernelOrder3{};
262
263 RealType m_BSplineEpsilon{ static_cast<RealType>(1e-3) };
264
265 inline typename RealImageType::IndexType
266 NumberToIndex(unsigned int number, typename RealImageType::SizeType size)
267 {
268 typename RealImageType::IndexType k;
269 k[0] = 1;
270
271 for (unsigned int i = 1; i < ImageDimension; ++i)
272 {
273 k[i] = size[ImageDimension - i - 1] * k[i - 1];
274 }
275 typename RealImageType::IndexType index;
276 for (unsigned int i = 0; i < ImageDimension; ++i)
277 {
278 index[ImageDimension - i - 1] = static_cast<unsigned int>(number / k[ImageDimension - i - 1]);
279 number %= k[ImageDimension - i - 1];
280 }
281 return index;
282 }
283};
284
285} // end namespace itk
286
287#ifndef ITK_MANUAL_INSTANTIATION
288# include "itkBSplineControlPointImageFilter.hxx"
289#endif
290
291#endif
Process a given a B-spline grid of control points.
void BeforeThreadedGenerateData() override
void PrintSelf(std::ostream &os, Indent indent) const override
typename OutputImageType::SpacingType SpacingType
unsigned int SplitRequestedRegion(unsigned int, unsigned int, OutputImageRegionType &) override
void CollapsePhiLattice(PointDataImageType *, PointDataImageType *, const RealType, const unsigned int)
typename OutputImageType::DirectionType DirectionType
typename OutputImageType::PointType PointType
~BSplineControlPointImageFilter() override=default
RealImageType::IndexType NumberToIndex(unsigned int number, typename RealImageType::SizeType size)
typename OutputImageType::RegionType OutputImageRegionType
typename OutputImageType::PointType OriginType
ControlPointLatticeType::Pointer RefineControlPointLattice(ArrayType)
typename PointSetType::PointDataContainer PointDataContainerType
typename OutputImageType::PixelType PixelType
void DynamicThreadedGenerateData(const OutputImageRegionType &) override
typename OutputImageType::RegionType RegionType
typename OutputImageType::IndexType IndexType
typename PointDataImageType::Pointer PointDataImagePointer
BSpline kernel used for density estimation and nonparametric regression.
BSpline kernel used for density estimation and nonparametric regression.
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.
A superclass of the N-dimensional mesh structure; supports point (geometric coordinate and attribute)...
Definition: itkPointSet.h:82
typename MeshTraits::PixelType PixelType
Definition: itkPointSet.h:100
typename MeshTraits::PointDataContainer PointDataContainer
Definition: itkPointSet.h:105
SmartPointer< Self > Pointer
static constexpr double e
Definition: itkMath.h:56
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
Represent a n-dimensional size (bounds) of a n-dimensional image.
Definition: itkSize.h:70