RTK  1.4.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 rtk::ForwardProjectionImageFilter< TOutputImage, TOutputImage > ForwardProjectionFilterType;
129  typedef typename ForwardProjectionFilterType::Pointer ForwardProjectionFilterPointer;
130  typedef rtk::BackProjectionImageFilter< TOutputImage, TOutputImage > BackProjectionFilterType;
132  typedef itk::MultiplyImageFilter<TOutputImage, TSingleComponentImage, TOutputImage> MultiplyFilterType;
133  typedef rtk::ReconstructionConjugateGradientOperator<TOutputImage,
134  TSingleComponentImage,
135  TWeightsImage> CGOperatorFilterType;
140  typedef typename TOutputImage::Pointer OutputImagePointer;
141  typedef itk::StatisticsImageFilter<TSingleComponentImage> StatisticsFilterType;
142 
143  typedef typename Superclass::ForwardProjectionType ForwardProjectionType;
144  typedef typename Superclass::BackProjectionType BackProjectionType;
145 
146  // If TOutputImage is an itk::Image of floats or double, so are the weights, and a simple Multiply filter is required
147  // If TOutputImage is an itk::Image of itk::Vector<float (or double)>, a BlockDiagonalMatrixVectorMultiply filter
148  // is needed. Thus the meta-programming construct
150  typedef itk::MultiplyImageFilter<TOutputImage, TOutputImage, TOutputImage> PlainMultiplyFilterType;
151  typedef typename std::conditional<std::is_same< TSingleComponentImage, TOutputImage>::value,
152  PlainMultiplyFilterType,
153  MatrixVectorMultiplyFilterType>::type MultiplyWithWeightsFilterType;
154 
156  void SetForwardProjectionFilter (ForwardProjectionType _arg) ITK_OVERRIDE;
157 
159  void SetBackProjectionFilter (BackProjectionType _arg) ITK_OVERRIDE;
160 
162  void SetSupportMask(const TSingleComponentImage *SupportMask);
163  typename TSingleComponentImage::ConstPointer GetSupportMask();
165 
167  itkSetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry)
168 
169  itkSetMacro(NumberOfIterations, int)
170  itkGetMacro(NumberOfIterations, int)
171 
172  itkSetMacro(IterationCosts, bool)
173  itkGetMacro(IterationCosts, bool)
174 
176  itkSetMacro(DisableDisplacedDetectorFilter, bool)
177  itkGetMacro(DisableDisplacedDetectorFilter, bool)
178 
181  itkSetMacro(Tikhonov, float)
182  itkGetMacro(Tikhonov, float)
183  itkSetMacro(Gamma, float)
184  itkGetMacro(Gamma, float)
185 
187  itkGetMacro(CudaConjugateGradient, bool)
188  itkSetMacro(CudaConjugateGradient, bool)
189 
191  const std::vector<double> &GetResidualCosts();
192 
193 protected:
194  ConjugateGradientConeBeamReconstructionFilter();
195  virtual ~ConjugateGradientConeBeamReconstructionFilter() ITK_OVERRIDE {}
196 
198  void GenerateData() ITK_OVERRIDE;
199 
211  typename MultiplyWithWeightsFilterType::Pointer m_MultiplyWithWeightsFilter;
212 
216 #if ITK_VERSION_MAJOR<5
217  void VerifyInputInformation() ITK_OVERRIDE {}
218 #else
219  void VerifyInputInformation() const ITK_OVERRIDE {}
220 #endif
221 
222 
225  void GenerateInputRequestedRegion() ITK_OVERRIDE;
226  void GenerateOutputInformation() ITK_OVERRIDE;
228 
230  typename TOutputImage::ConstPointer GetInputVolume();
231  typename TOutputImage::ConstPointer GetInputProjectionStack();
232  typename TWeightsImage::ConstPointer GetInputWeights();
234 
235 private:
236  ConjugateGradientConeBeamReconstructionFilter(const Self &); //purposely not implemented
237  void operator=(const Self &); //purposely not implemented
238 
240 
242  float m_Gamma;
243  float m_Tikhonov;
247 };
248 } //namespace RTK
249 
250 
251 #ifndef ITK_MANUAL_INSTANTIATION
252 #include "rtkConjugateGradientConeBeamReconstructionFilter.hxx"
253 #endif
254 
255 #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
const std::vector< double > & GetResidualCosts()
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
ForwardProjectionImageFilter< TOutputImage, TOutputImage >::Pointer m_ForwardProjectionFilter
#define itkSetMacro(name, type)
IterativeConeBeamReconstructionFilter< TOutputImage, TOutputImage > Superclass
void SetSupportMask(const TSingleComponentImage *SupportMask)