RTK  2.0.0
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:
115 
117  itkNewMacro(Self)
118 
119 
121 
123  void SetInputVolume(const TOutputImage* vol);
124  void SetInputProjectionStack(const TOutputImage* projs);
125  void SetInputWeights(const TWeightsImage* weights);
127 
128  typedef ForwardProjectionImageFilter< TOutputImage, TOutputImage > ForwardProjectionFilterType;
129  typedef typename ForwardProjectionFilterType::Pointer ForwardProjectionFilterPointer;
130  typedef BackProjectionImageFilter< TOutputImage, TOutputImage > BackProjectionFilterType;
132  typedef typename ConjugateGradientFilterType::Pointer ConjugateGradientFilterPointer;
133  typedef itk::MultiplyImageFilter<TOutputImage, TSingleComponentImage, TOutputImage> MultiplyFilterType;
134  typedef ReconstructionConjugateGradientOperator<TOutputImage,
135  TSingleComponentImage,
136  TWeightsImage> CGOperatorFilterType;
139  typedef typename TOutputImage::Pointer OutputImagePointer;
140  typedef itk::StatisticsImageFilter<TSingleComponentImage> StatisticsFilterType;
141 
142  typedef typename Superclass::ForwardProjectionType ForwardProjectionType;
143  typedef typename Superclass::BackProjectionType BackProjectionType;
144 
145  // If TOutputImage is an itk::Image of floats or double, so are the weights, and a simple Multiply filter is required
146  // If TOutputImage is an itk::Image of itk::Vector<float (or double)>, a BlockDiagonalMatrixVectorMultiply filter
147  // is needed. Thus the meta-programming construct
149  typedef itk::MultiplyImageFilter<TOutputImage, TOutputImage, TOutputImage> PlainMultiplyFilterType;
150  typedef typename std::conditional<std::is_same< TSingleComponentImage, TOutputImage>::value,
151  PlainMultiplyFilterType,
152  MatrixVectorMultiplyFilterType>::type MultiplyWithWeightsFilterType;
153  typedef typename itk::Image< typename TOutputImage::PixelType,
154  TOutputImage::ImageDimension> CPUOutputImageType;
155 #ifdef RTK_USE_CUDA
156  typedef typename std::conditional<!std::is_same< TOutputImage, CPUOutputImageType >::value &&
157  std::is_same< TSingleComponentImage, TOutputImage>::value,
158  CudaDisplacedDetectorImageFilter,
160  typedef typename std::conditional<!std::is_same< TOutputImage, CPUOutputImageType >::value &&
161  std::is_same< TSingleComponentImage, TOutputImage>::value,
162  CudaConstantVolumeSource,
164 #else
167 #endif
168 
170  void SetForwardProjectionFilter (ForwardProjectionType _arg) ITK_OVERRIDE;
171 
173  void SetBackProjectionFilter (BackProjectionType _arg) ITK_OVERRIDE;
174 
176  void SetSupportMask(const TSingleComponentImage *SupportMask);
177  typename TSingleComponentImage::ConstPointer GetSupportMask();
179 
181  itkSetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry)
182 
183  itkSetMacro(NumberOfIterations, int)
184  itkGetMacro(NumberOfIterations, int)
185 
186  itkSetMacro(IterationCosts, bool)
187  itkGetMacro(IterationCosts, bool)
188 
190  itkSetMacro(DisableDisplacedDetectorFilter, bool)
191  itkGetMacro(DisableDisplacedDetectorFilter, bool)
192 
195  itkSetMacro(Tikhonov, float)
196  itkGetMacro(Tikhonov, float)
197  itkSetMacro(Gamma, float)
198  itkGetMacro(Gamma, float)
199 
201  itkGetMacro(CudaConjugateGradient, bool)
202  itkSetMacro(CudaConjugateGradient, bool)
203 
205  const std::vector<double> &GetResidualCosts();
206 
207 protected:
208  ConjugateGradientConeBeamReconstructionFilter();
209  virtual ~ConjugateGradientConeBeamReconstructionFilter() ITK_OVERRIDE {}
210 
212  void GenerateData() ITK_OVERRIDE;
213 
218  ConjugateGradientFilterPointer m_ConjugateGradientFilter;
225  typename MultiplyWithWeightsFilterType::Pointer m_MultiplyWithWeightsFilter;
226 
230 #if ITK_VERSION_MAJOR<5
231  void VerifyInputInformation() ITK_OVERRIDE {}
232 #else
233  void VerifyInputInformation() const ITK_OVERRIDE {}
234 #endif
235 
236 
239  void GenerateInputRequestedRegion() ITK_OVERRIDE;
240  void GenerateOutputInformation() ITK_OVERRIDE;
242 
244  typename TOutputImage::ConstPointer GetInputVolume();
245  typename TOutputImage::ConstPointer GetInputProjectionStack();
246  typename TWeightsImage::ConstPointer GetInputWeights();
248 
249  template < typename ImageType, typename IterativeConeBeamReconstructionFilter<TOutputImage>::template EnableCudaScalarAndVectorType<ImageType>* = ITK_NULLPTR >
250  ConjugateGradientFilterPointer InstantiateCudaConjugateGradientImageFilter();
251 
252  template < typename ImageType, typename IterativeConeBeamReconstructionFilter<TOutputImage>::template DisableCudaScalarAndVectorType<ImageType>* = ITK_NULLPTR >
253  ConjugateGradientFilterPointer InstantiateCudaConjugateGradientImageFilter();
254 
255 private:
256  ConjugateGradientConeBeamReconstructionFilter(const Self &); //purposely not implemented
257  void operator=(const Self &); //purposely not implemented
258 
260 
262  float m_Gamma;
263  float m_Tikhonov;
267 };
268 } //namespace RTK
269 
270 
271 #ifndef ITK_MANUAL_INSTANTIATION
272 #include "rtkConjugateGradientConeBeamReconstructionFilter.hxx"
273 #endif
274 
275 #endif
Implements the operator A used in conjugate gradient reconstruction.
Base class for forward projection, i.e. accumulation along x-ray lines.
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 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)