RTK  2.6.0
Reconstruction Toolkit
rtkFFTProjectionsConvolutionImageFilter.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  * 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 
19 #ifndef rtkFFTProjectionsConvolutionImageFilter_h
20 #define rtkFFTProjectionsConvolutionImageFilter_h
21 
22 #include <itkImageToImageFilter.h>
23 #include <itkConceptChecking.h>
24 
25 #include "rtkConfiguration.h"
26 #include "rtkMacro.h"
27 
28 namespace rtk
29 {
30 
44 template <class TInputImage, class TOutputImage, class TFFTPrecision>
45 class ITK_TEMPLATE_EXPORT FFTProjectionsConvolutionImageFilter
46  : public itk::ImageToImageFilter<TInputImage, TOutputImage>
47 {
48 public:
49  ITK_DISALLOW_COPY_AND_MOVE(FFTProjectionsConvolutionImageFilter);
50 
56 
58  using InputImageType = TInputImage;
59  using OutputImageType = TOutputImage;
60  using RegionType = typename InputImageType::RegionType;
61  using IndexType = typename InputImageType::IndexType;
62  using SizeType = typename InputImageType::SizeType;
63 
65  using FFTInputImagePointer = typename FFTInputImageType::Pointer;
66  using FFTOutputImageType = typename itk::Image<std::complex<TFFTPrecision>, TInputImage::ImageDimension>;
67  using FFTOutputImagePointer = typename FFTOutputImageType::Pointer;
69 
71  static constexpr unsigned int ImageDimension = TOutputImage::ImageDimension;
72 
74 #ifdef itkOverrideGetNameOfClassMacro
75  itkOverrideGetNameOfClassMacro(FFTProjectionsConvolutionImageFilter);
76 #else
78 #endif
79 
80 
94  itkGetConstMacro(GreatestPrimeFactor, int);
95  itkSetMacro(GreatestPrimeFactor, int);
97 
101  itkGetConstMacro(TruncationCorrection, double);
102  itkSetMacro(TruncationCorrection, double);
104 
108  itkGetConstMacro(ZeroPadFactors, ZeroPadFactorsType);
109  virtual void
111  {
112  if (m_ZeroPadFactors != _arg)
113  {
114  m_ZeroPadFactors = _arg;
115  m_ZeroPadFactors[0] = std::max(m_ZeroPadFactors[0], 1);
116  m_ZeroPadFactors[1] = std::max(m_ZeroPadFactors[1], 1);
117  m_ZeroPadFactors[0] = std::min(m_ZeroPadFactors[0], 2);
118  m_ZeroPadFactors[1] = std::min(m_ZeroPadFactors[1], 2);
119  this->Modified();
120  }
121  }
123 
124 #ifdef ITK_USE_CONCEPT_CHECKING
126 #endif
127 
128 protected:
130  ~FFTProjectionsConvolutionImageFilter() override = default;
131 
132  void
133  GenerateInputRequestedRegion() override;
134 
135  void
136  BeforeThreadedGenerateData() override;
137 
138  void
139  AfterThreadedGenerateData() override;
140 
141  void
142  ThreadedGenerateData(const RegionType & outputRegionForThread, ThreadIdType threadId) override;
143 
148  virtual FFTInputImagePointer
149  PadInputImageRegion(const RegionType & inputRegion);
150  RegionType
151  GetPaddedImageRegion(const RegionType & inputRegion);
153 
154  void
155  PrintSelf(std::ostream & os, itk::Indent indent) const override;
156 
157  bool
158  IsPrime(int n) const;
159 
160  int
161  GreatestPrimeFactor(int n) const;
162 
165  virtual void
166  UpdateFFTProjectionsConvolutionKernel(const SizeType size) = 0;
167 
173  virtual void
174  UpdateTruncationMirrorWeights();
175  typename std::vector<TFFTPrecision> m_TruncationMirrorWeights;
176 
179  int m_KernelDimension{ 1 };
180 
185 
186 private:
190  double m_TruncationCorrection{ 0. };
191  int
192  GetTruncationCorrectionExtent();
194 
199 
203  int m_GreatestPrimeFactor{ 2 };
204  int m_BackupNumberOfThreads{ 1 };
205 }; // end of class
206 
207 } // end namespace rtk
208 
209 #ifndef ITK_MANUAL_INSTANTIATION
210 # include "rtkFFTProjectionsConvolutionImageFilter.hxx"
211 #endif
212 
213 #endif
Base class for 1D or 2D FFT based convolution of projections.
ITKCommon_EXPORT bool IsPrime(unsigned int n)
TInputImage InputImageType
#define itkSetMacro(name, type)
TOutputImage OutputImageType
typename itk::Image< TFFTPrecision, TInputImage::ImageDimension > FFTInputImageType
ITKCommon_EXPORT unsigned int GreatestPrimeFactor(unsigned int n)
unsigned int ThreadIdType
#define itkConceptMacro(name, concept)
typename itk::Image< std::complex< TFFTPrecision >, TInputImage::ImageDimension > FFTOutputImageType