RTK  2.1.0
Reconstruction Toolkit
rtkFFTRampImageFilter.h
Go to the documentation of this file.
1 /*=========================================================================
2  *
3  * Copyright RTK Consortium
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  * http://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 
19 #ifndef rtkFFTRampImageFilter_h
20 #define rtkFFTRampImageFilter_h
21 
22 #include <itkConceptChecking.h>
23 #include "rtkConfiguration.h"
25 #include "rtkMacro.h"
26 
27 // The Set macro is redefined to clear the current FFT kernel when a parameter
28 // is modified.
29 #undef itkSetMacro
30 #define itkSetMacro(name, type) \
31  virtual void Set##name(const type _arg) \
32  { \
33  itkDebugMacro("setting " #name " to " << _arg); \
34  CLANG_PRAGMA_PUSH \
35  CLANG_SUPPRESS_Wfloat_equal if (this->m_##name != _arg) \
36  { \
37  this->m_##name = _arg; \
38  this->Modified(); \
39  this->m_KernelFFT = nullptr; \
40  } \
41  CLANG_PRAGMA_POP \
42  }
43 
44 namespace rtk
45 {
46 
60 template <class TInputImage, class TOutputImage = TInputImage, class TFFTPrecision = double>
61 class ITK_EXPORT FFTRampImageFilter
62  : public rtk::FFTProjectionsConvolutionImageFilter<TInputImage, TOutputImage, TFFTPrecision>
63 {
64 public:
65 #if ITK_VERSION_MAJOR == 5 && ITK_VERSION_MINOR == 1
66  ITK_DISALLOW_COPY_AND_ASSIGN(FFTRampImageFilter);
67 #else
68  ITK_DISALLOW_COPY_AND_MOVE(FFTRampImageFilter);
69 #endif
70 
76 
78  using InputImageType = TInputImage;
79  using OutputImageType = TOutputImage;
80  using FFTPrecisionType = TFFTPrecision;
81  using IndexType = typename InputImageType::IndexType;
82  using SizeType = typename InputImageType::SizeType;
83 
84  using FFTInputImageType = typename Superclass::FFTInputImageType;
85  using FFTInputImagePointer = typename FFTInputImageType::Pointer;
86  using FFTOutputImageType = typename Superclass::FFTOutputImageType;
87  using FFTOutputImagePointer = typename FFTOutputImageType::Pointer;
88 
90  itkNewMacro(Self);
91 
94 
96  itkGetConstMacro(HannCutFrequency, double);
97  itkSetMacro(HannCutFrequency, double);
99 
101  itkGetConstMacro(CosineCutFrequency, double);
102  itkSetMacro(CosineCutFrequency, double);
104 
106  itkGetConstMacro(HammingFrequency, double);
107  itkSetMacro(HammingFrequency, double);
109 
111  itkGetConstMacro(HannCutFrequencyY, double);
112  itkSetMacro(HannCutFrequencyY, double);
114 
122  itkGetConstMacro(RamLakCutFrequency, double);
123  itkSetMacro(RamLakCutFrequency, double);
125 
133  itkGetConstMacro(SheppLoganCutFrequency, double);
134  itkSetMacro(SheppLoganCutFrequency, double);
136 
137 protected:
139  ~FFTRampImageFilter() override = default;
140 
141  void
142  GenerateInputRequestedRegion() override;
143 
146  void
147  UpdateFFTProjectionsConvolutionKernel(const SizeType s) override;
148 
149 private:
154  double m_HannCutFrequency{ 0. };
155  double m_CosineCutFrequency{ 0. };
156  double m_HammingFrequency{ 0. };
157  double m_HannCutFrequencyY{ 0. };
158 
161  double m_RamLakCutFrequency{ 0. };
162  double m_SheppLoganCutFrequency{ 0. };
163 
165 }; // end of class
166 
167 } // end namespace rtk
168 
169 #ifndef ITK_MANUAL_INSTANTIATION
170 # include "rtkFFTRampImageFilter.hxx"
171 #endif
172 
173 // Rollback to the original definition of the Set macro
174 #undef itkSetMacro
175 #define itkSetMacro(name, type) \
176  virtual void Set##name(const type _arg) \
177  { \
178  itkDebugMacro("setting " #name " to " << _arg); \
179  CLANG_PRAGMA_PUSH \
180  CLANG_SUPPRESS_Wfloat_equal if (this->m_##name != _arg) \
181  { \
182  this->m_##name = _arg; \
183  this->Modified(); \
184  } \
185  CLANG_PRAGMA_POP \
186  }
187 
188 #endif
Base class for 1D or 2D FFT based convolution of projections.
TInputImage InputImageType
#define itkSetMacro(name, type)
TOutputImage OutputImageType
typename itk::Image< TFFTPrecision, TInputImage::ImageDimension > FFTInputImageType
Implements the ramp image filter of the filtered backprojection algorithm.
typename itk::Image< std::complex< TFFTPrecision >, TInputImage::ImageDimension > FFTOutputImageType