RTK  2.4.1
Reconstruction Toolkit
rtkConjugateGradientConeBeamReconstructionFilter.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 rtkConjugateGradientConeBeamReconstructionFilter_h
20 #define rtkConjugateGradientConeBeamReconstructionFilter_h
21 
22 #include <itkMultiplyImageFilter.h>
24 #include <itkProcessObject.h>
25 #include <itkObject.h>
26 #include <itkCommand.h>
27 #include <itkIterationReporter.h>
28 
34 #include "rtkConstantImageSource.h"
37 
38 #ifdef RTK_USE_CUDA
42 #endif
43 
44 namespace rtk
45 {
109 template <typename TOutputImage, typename TSingleComponentImage = TOutputImage, typename TWeightsImage = TOutputImage>
111  : public IterativeConeBeamReconstructionFilter<TOutputImage>
112 {
113 public:
114  ITK_DISALLOW_COPY_AND_MOVE(ConjugateGradientConeBeamReconstructionFilter);
115 
120 
122  itkNewMacro(Self);
123 
126 
128  void
129  SetInputVolume(const TOutputImage * vol);
130  void
131  SetInputProjectionStack(const TOutputImage * projs);
132  void
133  SetInputWeights(const TWeightsImage * weights);
135 
142  using CGOperatorFilterType =
146  using OutputImagePointer = typename TOutputImage::Pointer;
148 
149  using ForwardProjectionType = typename Superclass::ForwardProjectionType;
150  using BackProjectionType = typename Superclass::BackProjectionType;
151 
152  // If TOutputImage is an itk::Image of floats or double, so are the weights, and a simple Multiply filter is required
153  // If TOutputImage is an itk::Image of itk::Vector<float (or double)>, a BlockDiagonalMatrixVectorMultiply filter
154  // is needed. Thus the meta-programming construct
157  typedef typename std::conditional<std::is_same<TSingleComponentImage, TOutputImage>::value,
161 #ifdef RTK_USE_CUDA
162  typedef typename std::conditional<!std::is_same<TOutputImage, CPUOutputImageType>::value &&
163  std::is_same<TSingleComponentImage, TOutputImage>::value,
164  CudaDisplacedDetectorImageFilter,
166  typedef typename std::conditional<!std::is_same<TOutputImage, CPUOutputImageType>::value &&
167  std::is_same<TSingleComponentImage, TOutputImage>::value,
168  CudaConstantVolumeSource,
170 #else
173 #endif
174 
176  void
177  SetSupportMask(const TSingleComponentImage * SupportMask);
178  typename TSingleComponentImage::ConstPointer
179  GetSupportMask();
181 
183  itkSetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry);
184 
185  itkSetMacro(NumberOfIterations, int);
186  itkGetMacro(NumberOfIterations, int);
187 
189  itkSetMacro(DisableDisplacedDetectorFilter, bool);
190  itkGetMacro(DisableDisplacedDetectorFilter, bool);
192 
195  itkSetMacro(Tikhonov, float);
196  itkGetMacro(Tikhonov, float);
197  itkSetMacro(Gamma, float);
198  itkGetMacro(Gamma, float);
200 
202  itkGetMacro(CudaConjugateGradient, bool);
203  itkSetMacro(CudaConjugateGradient, bool);
205 
206 protected:
208  ~ConjugateGradientConeBeamReconstructionFilter() override = default;
209 
211  void
212  VerifyPreconditions() ITKv5_CONST override;
213 
215  void
216  GenerateData() override;
217 
219  typename MultiplyFilterType::Pointer m_MultiplyProjectionsFilter;
220  typename MultiplyFilterType::Pointer m_MultiplyVolumeFilter;
221  typename MultiplyFilterType::Pointer m_MultiplyOutputFilter;
222  ConjugateGradientFilterPointer m_ConjugateGradientFilter;
223  typename CGOperatorFilterType::Pointer m_CGOperator;
224  typename ForwardProjectionImageFilter<TOutputImage, TOutputImage>::Pointer m_ForwardProjectionFilter;
225  typename BackProjectionImageFilter<TOutputImage, TOutputImage>::Pointer m_BackProjectionFilter;
226  typename BackProjectionImageFilter<TOutputImage, TOutputImage>::Pointer m_BackProjectionFilterForB;
227  typename DisplacedDetectorFilterType::Pointer m_DisplacedDetectorFilter;
228  typename ConstantImageSourceType::Pointer m_ConstantVolumeSource;
229  typename MultiplyWithWeightsFilterType::Pointer m_MultiplyWithWeightsFilter;
230 
234  void
235  VerifyInputInformation() const override
236  {}
237 
240  void
241  GenerateInputRequestedRegion() override;
242  void
243  GenerateOutputInformation() override;
245 
247  typename TOutputImage::ConstPointer
248  GetInputVolume();
249  typename TOutputImage::ConstPointer
250  GetInputProjectionStack();
251  typename TWeightsImage::ConstPointer
252  GetInputWeights();
254 
255  template <typename ImageType,
256  typename IterativeConeBeamReconstructionFilter<TOutputImage>::template EnableCudaScalarAndVectorType<
257  ImageType> * = nullptr>
258  ConjugateGradientFilterPointer
259  InstantiateCudaConjugateGradientImageFilter();
260 
261  template <typename ImageType,
262  typename IterativeConeBeamReconstructionFilter<TOutputImage>::template DisableCudaScalarAndVectorType<
263  ImageType> * = nullptr>
264  ConjugateGradientFilterPointer
265  InstantiateCudaConjugateGradientImageFilter();
266 
268  void
269  ReportProgress(itk::Object *, const itk::EventObject &);
270 
271 private:
273 
275  float m_Gamma;
276  float m_Tikhonov;
279 
280  // Iteration reporting
282 };
283 } // namespace rtk
284 
285 
286 #ifndef ITK_MANUAL_INSTANTIATION
287 # include "rtkConjugateGradientConeBeamReconstructionFilter.hxx"
288 #endif
289 
290 #endif
typename OutputImageType::Pointer OutputImagePointer
Implements the operator A used in conjugate gradient reconstruction.
typename itk::Image< typename TOutputImage::PixelType, TOutputImage::ImageDimension > CPUOutputImageType
Generate an n-dimensional image with constant pixel values.
std::conditional< std::is_same< TSingleComponentImage, TOutputImage >::value, PlainMultiplyFilterType, MatrixVectorMultiplyFilterType >::type MultiplyWithWeightsFilterType
Weigting for displaced detectors.
Projection geometry for a source and a 2-D flat panel.
#define itkSetMacro(name, type)
Mother class for cone beam reconstruction filters which need runtime selection of their forward and b...
Solves AX = B by conjugate gradient.