ITK  6.0.0
Insight Toolkit
itkWindowedSincInterpolateImageFunction.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 itkWindowedSincInterpolateImageFunction_h
19#define itkWindowedSincInterpolateImageFunction_h
20
24#include "itkMath.h"
25
26namespace itk
27{
28namespace Function
29{
37template <unsigned int VRadius, typename TInput = double, typename TOutput = double>
38class ITK_TEMPLATE_EXPORT CosineWindowFunction
39{
40public:
41 inline TOutput
42 operator()(const TInput & A) const
43 {
44 return static_cast<TOutput>(std::cos(A * m_Factor));
45 }
48private:
50 static const double m_Factor;
51};
52
60template <unsigned int VRadius, typename TInput = double, typename TOutput = double>
61class ITK_TEMPLATE_EXPORT HammingWindowFunction
62{
63public:
64 inline TOutput
65 operator()(const TInput & A) const
66 {
67 return static_cast<TOutput>(0.54 + 0.46 * std::cos(A * m_Factor));
68 }
71private:
73 static const double m_Factor;
74};
75
83template <unsigned int VRadius, typename TInput = double, typename TOutput = double>
84class ITK_TEMPLATE_EXPORT WelchWindowFunction
85{
86public:
87 inline TOutput
88 operator()(const TInput & A) const
89 {
90 return static_cast<TOutput>(1.0 - A * m_Factor * A);
91 }
94private:
96 static const double m_Factor;
97};
98
108template <unsigned int VRadius, typename TInput = double, typename TOutput = double>
109class ITK_TEMPLATE_EXPORT LanczosWindowFunction
110{
111public:
112 inline TOutput
113 operator()(const TInput & A) const
114 {
115 if (A == 0.0)
116 {
117 return static_cast<TOutput>(1.0);
118 }
119 double z = m_Factor * A;
120 return static_cast<TOutput>(std::sin(z) / z);
121 }
124private:
126 static const double m_Factor;
127};
128
136template <unsigned int VRadius, typename TInput = double, typename TOutput = double>
137class ITK_TEMPLATE_EXPORT BlackmanWindowFunction
138{
139public:
140 inline TOutput
141 operator()(const TInput & A) const
142 {
143 return static_cast<TOutput>(0.42 + 0.5 * std::cos(A * m_Factor1) + 0.08 * std::cos(A * m_Factor2));
144 }
147private:
149 static const double m_Factor1;
150
152 static const double m_Factor2;
153};
154} // namespace Function
155
266template <typename TInputImage,
267 unsigned int VRadius,
268 typename TWindowFunction = Function::HammingWindowFunction<VRadius>,
270 class TCoordRep = double>
271class ITK_TEMPLATE_EXPORT WindowedSincInterpolateImageFunction : public InterpolateImageFunction<TInputImage, TCoordRep>
272{
273public:
274 ITK_DISALLOW_COPY_AND_MOVE(WindowedSincInterpolateImageFunction);
280
283
285 itkOverrideGetNameOfClassMacro(WindowedSincInterpolateImageFunction);
286
288 itkNewMacro(Self);
289
291 using typename Superclass::OutputType;
292
294 using typename Superclass::InputImageType;
295
297 using typename Superclass::RealType;
298
300 static constexpr unsigned int ImageDimension = Superclass::ImageDimension;
301
303 using typename Superclass::IndexType;
304 using typename Superclass::IndexValueType;
305
307 using typename Superclass::SizeType;
308
310 using ImageType = TInputImage;
311
313 using typename Superclass::ContinuousIndexType;
314
315 void
316 SetInputImage(const ImageType * image) override;
317
325 EvaluateAtContinuousIndex(const ContinuousIndexType & index) const override;
326
328 GetRadius() const override
329 {
330 auto radius = SizeType::Filled(VRadius);
331 return radius;
332 }
333
334protected:
337 void
338 PrintSelf(std::ostream & os, Indent indent) const override;
339
340private:
341 // Internal type alias
343
344 // Constant to store twice the radius
345 static constexpr unsigned int m_WindowSize{ 2 * VRadius };
346
348 TWindowFunction m_WindowFunction{};
349
351 static constexpr unsigned int m_OffsetTableSize = Math::UnsignedPower(m_WindowSize, ImageDimension);
352
355 unsigned int m_OffsetTable[m_OffsetTableSize]{};
356
358 unsigned int m_WeightOffsetTable[m_OffsetTableSize][ImageDimension]{};
359
361 inline double
362 Sinc(double x) const
363 {
364 const double px = itk::Math::pi * x;
365
366 return (x == 0.0) ? 1.0 : std::sin(px) / px;
367 }
368};
369} // namespace itk
370
371#ifndef ITK_MANUAL_INSTANTIATION
372# include "itkWindowedSincInterpolateImageFunction.hxx"
373#endif
374
375#endif // _itkWindowedSincInterpolateImageFunction_h
Const version of NeighborhoodIterator, defining iteration of a local N-dimensional neighborhood of pi...
Control indentation during Print() invocation.
Definition: itkIndent.h:50
Base class for all image interpolators.
Light weight base class for most itk classes.
Use the windowed sinc function to interpolate.
void SetInputImage(const ImageType *image) override
void PrintSelf(std::ostream &os, Indent indent) const override
~WindowedSincInterpolateImageFunction() override=default
OutputType EvaluateAtContinuousIndex(const ContinuousIndexType &index) const override
constexpr TReturnType UnsignedPower(const uintmax_t base, const uintmax_t exponent) noexcept
Definition: itkMath.h:797
static constexpr double pi
Definition: itkMath.h:66
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
long IndexValueType
Definition: itkIntTypes.h:93