RTK  2.0.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  * 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 rtkConjugateGradientConeBeamReconstructionFilter_h
20 #define rtkConjugateGradientConeBeamReconstructionFilter_h
21 
22 #include <itkMultiplyImageFilter.h>
24 
30 #include "rtkConstantImageSource.h"
33 
34 #ifdef RTK_USE_CUDA
38 #endif
39 
40 namespace rtk
41 {
105 template< typename TOutputImage,
106  typename TSingleComponentImage = TOutputImage,
107  typename TWeightsImage = TOutputImage>
109 {
110 public:
111  ITK_DISALLOW_COPY_AND_ASSIGN(ConjugateGradientConeBeamReconstructionFilter);
112 
117 
119  itkNewMacro(Self)
120 
121 
123 
125  void SetInputVolume(const TOutputImage* vol);
126  void SetInputProjectionStack(const TOutputImage* projs);
127  void SetInputWeights(const TWeightsImage* weights);
129 
130  using ForwardProjectionFilterType = ForwardProjectionImageFilter< TOutputImage, TOutputImage >;
132  using BackProjectionFilterType = BackProjectionImageFilter< TOutputImage, TOutputImage >;
135  using MultiplyFilterType = itk::MultiplyImageFilter<TOutputImage, TSingleComponentImage, TOutputImage>;
137  TSingleComponentImage,
138  TWeightsImage>;
141  using OutputImagePointer = typename TOutputImage::Pointer;
142  using StatisticsFilterType = itk::StatisticsImageFilter<TSingleComponentImage>;
143 
146 
147  // If TOutputImage is an itk::Image of floats or double, so are the weights, and a simple Multiply filter is required
148  // If TOutputImage is an itk::Image of itk::Vector<float (or double)>, a BlockDiagonalMatrixVectorMultiply filter
149  // is needed. Thus the meta-programming construct
151  using PlainMultiplyFilterType = itk::MultiplyImageFilter<TOutputImage, TOutputImage, TOutputImage>;
152  typedef typename std::conditional<std::is_same< TSingleComponentImage, TOutputImage>::value,
155  using CPUOutputImageType = typename itk::Image< typename TOutputImage::PixelType,
156  TOutputImage::ImageDimension>;
157 #ifdef RTK_USE_CUDA
158  typedef typename std::conditional<!std::is_same< TOutputImage, CPUOutputImageType >::value &&
159  std::is_same< TSingleComponentImage, TOutputImage>::value,
160  CudaDisplacedDetectorImageFilter,
162  typedef typename std::conditional<!std::is_same< TOutputImage, CPUOutputImageType >::value &&
163  std::is_same< TSingleComponentImage, TOutputImage>::value,
164  CudaConstantVolumeSource,
166 #else
167  using DisplacedDetectorFilterType = DisplacedDetectorImageFilter<TWeightsImage>;
168  using ConstantImageSourceType = ConstantImageSource<TOutputImage>;
169 #endif
170 
172  void SetForwardProjectionFilter (ForwardProjectionType _arg) override;
173 
175  void SetBackProjectionFilter (BackProjectionType _arg) override;
176 
178  void SetSupportMask(const TSingleComponentImage *SupportMask);
179  typename TSingleComponentImage::ConstPointer GetSupportMask();
181 
183  itkSetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry)
184 
185  itkSetMacro(NumberOfIterations, int)
186  itkGetMacro(NumberOfIterations, int)
187 
188  itkSetMacro(IterationCosts, bool)
189  itkGetMacro(IterationCosts, bool)
190 
192  itkSetMacro(DisableDisplacedDetectorFilter, bool)
193  itkGetMacro(DisableDisplacedDetectorFilter, bool)
194 
197  itkSetMacro(Tikhonov, float)
198  itkGetMacro(Tikhonov, float)
199  itkSetMacro(Gamma, float)
200  itkGetMacro(Gamma, float)
201 
203  itkGetMacro(CudaConjugateGradient, bool)
204  itkSetMacro(CudaConjugateGradient, bool)
205 
207  const std::vector<double> &GetResidualCosts();
208 
209 protected:
210  ConjugateGradientConeBeamReconstructionFilter();
211  ~ConjugateGradientConeBeamReconstructionFilter() override = default;
212 
214  void GenerateData() override;
215 
222  typename ForwardProjectionImageFilter<TOutputImage, TOutputImage>::Pointer m_ForwardProjectionFilter;
223  typename BackProjectionImageFilter<TOutputImage, TOutputImage>::Pointer m_BackProjectionFilter;
224  typename BackProjectionImageFilter<TOutputImage, TOutputImage>::Pointer m_BackProjectionFilterForB;
225  typename DisplacedDetectorFilterType::Pointer m_DisplacedDetectorFilter;
226  typename ConstantImageSourceType::Pointer m_ConstantVolumeSource;
227  typename MultiplyWithWeightsFilterType::Pointer m_MultiplyWithWeightsFilter;
228 
232 #if ITK_VERSION_MAJOR<5
233  void VerifyInputInformation() override {}
234 #else
235  void VerifyInputInformation() const override {}
236 #endif
237 
238 
241  void GenerateInputRequestedRegion() override;
242  void GenerateOutputInformation() override;
244 
246  typename TOutputImage::ConstPointer GetInputVolume();
247  typename TOutputImage::ConstPointer GetInputProjectionStack();
248  typename TWeightsImage::ConstPointer GetInputWeights();
250 
251  template < typename ImageType, typename IterativeConeBeamReconstructionFilter<TOutputImage>::template EnableCudaScalarAndVectorType<ImageType>* = nullptr >
253 
254  template < typename ImageType, typename IterativeConeBeamReconstructionFilter<TOutputImage>::template DisableCudaScalarAndVectorType<ImageType>* = nullptr >
256 
257 private:
259 
261  float m_Gamma;
262  float m_Tikhonov;
266 };
267 } //namespace RTK
268 
269 
270 #ifndef ITK_MANUAL_INSTANTIATION
271 #include "rtkConjugateGradientConeBeamReconstructionFilter.hxx"
272 #endif
273 
274 #endif
Implements the operator A used in conjugate gradient reconstruction.
BackProjectionImageFilter< TOutputImage, TOutputImage >::Pointer m_BackProjectionFilter
TOutputImage::ConstPointer GetInputProjectionStack()
Generate an n-dimensional image with constant pixel values.
BackProjectionImageFilter< TOutputImage, TOutputImage >::Pointer m_BackProjectionFilterForB
TWeightsImage::ConstPointer GetInputWeights()
Weigting for displaced detectors.
STL namespace.
std::conditional< std::is_same< TSingleComponentImage, TOutputImage >::value, PlainMultiplyFilterType, MatrixVectorMultiplyFilterType >::type MultiplyWithWeightsFilterType
typename itk::Image< typename TOutputImage::PixelType, TOutputImage::ImageDimension > CPUOutputImageType
typename std::enable_if< !std::is_same< CPUImageType, ImageType >::value &&std::is_same< typename itk::PixelTraits< typename ImageType::PixelType >::ValueType, float >::value &&(itk::PixelTraits< typename ImageType::PixelType >::Dimension==1||itk::PixelTraits< typename ImageType::PixelType >::Dimension==3) >::type EnableCudaScalarAndVectorType
const std::vector< double > & GetResidualCosts()
ConjugateGradientFilterPointer InstantiateCudaConjugateGradientImageFilter()
Projection geometry for a source and a 2-D flat panel.
void SetInputProjectionStack(const TOutputImage *projs)
TSingleComponentImage::ConstPointer GetSupportMask()
void SetInputWeights(const TWeightsImage *weights)
void SetInputVolume(const TOutputImage *vol)
Mother class for cone beam reconstruction filters which need runtime selection of their forward and b...
Solves AX = B by conjugate gradient.
void SetBackProjectionFilter(BackProjectionType _arg) override
void SetForwardProjectionFilter(ForwardProjectionType _arg) override
typename std::enable_if< std::is_same< CPUImageType, ImageType >::value||!std::is_same< typename itk::PixelTraits< typename ImageType::PixelType >::ValueType, float >::value||(itk::PixelTraits< typename ImageType::PixelType >::Dimension!=1 &&itk::PixelTraits< typename ImageType::PixelType >::Dimension!=3) >::type DisableCudaScalarAndVectorType
ForwardProjectionImageFilter< TOutputImage, TOutputImage >::Pointer m_ForwardProjectionFilter
#define itkSetMacro(name, type)
void SetSupportMask(const TSingleComponentImage *SupportMask)