RTK  1.4.0
Reconstruction Toolkit
rtkReconstructionConjugateGradientOperator.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 rtkReconstructionConjugateGradientOperator_h
20 #define rtkReconstructionConjugateGradientOperator_h
21 
22 #include <itkMultiplyImageFilter.h>
23 #include <itkAddImageFilter.h>
24 
25 #include "rtkConstantImageSource.h"
26 
30 
34 
35 #ifdef RTK_USE_CUDA
38 #endif
39 
40 namespace rtk
41 {
42 
117 template< typename TOutputImage,
118  typename TSingleComponentImage = TOutputImage,
119  typename TWeightsImage = TOutputImage>
121 {
122 public:
127 #ifdef RTK_USE_CUDA
128  typedef itk::CudaImage<itk::CovariantVector<typename TOutputImage::PixelType, TOutputImage::ImageDimension>, TOutputImage::ImageDimension > GradientImageType;
129 #else
131 #endif
132 
134  itkNewMacro(Self)
135 
136 
137  void SetInputVolume(const TOutputImage* vol);
138  void SetInputProjectionStack(const TOutputImage* projs);
139  void SetInputWeights(const TWeightsImage* weights);
141 
143  itkTypeMacro(rtkReconstructionConjugateGradientOperator, ConjugateGradientOperator)
144 
145  typedef rtk::BackProjectionImageFilter< TOutputImage, TOutputImage > BackProjectionFilterType;
146  typedef typename BackProjectionFilterType::Pointer BackProjectionFilterPointer;
147 
148  typedef rtk::ForwardProjectionImageFilter< TOutputImage, TOutputImage > ForwardProjectionFilterType;
149  typedef typename ForwardProjectionFilterType::Pointer ForwardProjectionFilterPointer;
150 
152  typedef itk::MultiplyImageFilter<TOutputImage, TSingleComponentImage> MultiplyFilterType;
153  typedef itk::AddImageFilter<TOutputImage> AddFilterType;
154 
155  // If TOutputImage is an itk::Image of floats or double, so are the weights, and a simple Multiply filter is required
156  // If TOutputImage is an itk::Image of itk::Vector<float (or double)>, a BlockDiagonalMatrixVectorMultiply filter
157  // is needed. Thus the meta-programming construct
159  typedef itk::MultiplyImageFilter<TOutputImage, TOutputImage, TOutputImage> PlainMultiplyFilterType;
160  typedef typename std::conditional<std::is_same< TSingleComponentImage, TOutputImage>::value,
161  PlainMultiplyFilterType,
162  MatrixVectorMultiplyFilterType>::type MultiplyWithWeightsFilterType;
163 
164 // typedef rtk::LaplacianImageFilter<TOutputImage, GradientImageType> LaplacianFilterType;
165 
166  typedef typename TOutputImage::Pointer OutputImagePointer;
167 
169  void SetBackProjectionFilter (const BackProjectionFilterPointer _arg);
170 
172  void SetForwardProjectionFilter (const ForwardProjectionFilterPointer _arg);
173 
175  void SetSupportMask(const TSingleComponentImage *SupportMask);
176  typename TSingleComponentImage::ConstPointer GetSupportMask();
178 
180  itkSetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry)
181 
184  itkSetMacro(Gamma, float)
185  itkGetMacro(Gamma, float)
186  itkSetMacro(Tikhonov, float)
187  itkGetMacro(Tikhonov, float)
188 
189 protected:
191  virtual ~ReconstructionConjugateGradientOperator() ITK_OVERRIDE {}
192 
194  void GenerateData() ITK_OVERRIDE;
195 
199 
204 // typename MultiplyFilterType::Pointer m_MultiplyLaplacianFilter;
208 // typename LaplacianFilterType::Pointer m_LaplacianFilter;
209  typename MultiplyWithWeightsFilterType::Pointer m_MultiplyWithWeightsFilter;
210 
213  float m_Gamma; //Strength of the laplacian regularization
214  float m_Tikhonov; //Strength of the Tikhonov regularization
215 
217  typename TOutputImage::Pointer m_FloatingInputPointer, m_FloatingOutputPointer;
218 
221  void VerifyInputInformation() ITK_OVERRIDE {}
222 
224  void GenerateInputRequestedRegion() ITK_OVERRIDE;
225  void GenerateOutputInformation() ITK_OVERRIDE;
227 
229  typename TOutputImage::ConstPointer GetInputVolume();
230  typename TOutputImage::ConstPointer GetInputProjectionStack();
231  typename TWeightsImage::ConstPointer GetInputWeights();
233 
234 private:
235  ReconstructionConjugateGradientOperator(const Self &); //purposely not implemented
236  void operator=(const Self &); //purposely not implemented
237 
238 };
239 } //namespace RTK
240 
241 
242 #ifndef ITK_MANUAL_INSTANTIATION
243 #include "rtkReconstructionConjugateGradientOperator.hxx"
244 #endif
245 
246 #endif
itk::Image< itk::CovariantVector< typename TOutputImage::PixelType, TOutputImage::ImageDimension >, TOutputImage::ImageDimension > GradientImageType
TOutputImage::ConstPointer GetInputVolume()
Implements the operator A used in conjugate gradient reconstruction.
Base class for forward projection, i.e. accumulation along x-ray lines.
void SetBackProjectionFilter(const BackProjectionFilterPointer _arg)
rtk::ThreeDCircularProjectionGeometry::ConstPointer m_Geometry
Generate an n-dimensional image with constant pixel values.
TWeightsImage::ConstPointer GetInputWeights()
void SetForwardProjectionFilter(const ForwardProjectionFilterPointer _arg)
STL namespace.
TSingleComponentImage::ConstPointer GetSupportMask()
void SetInputWeights(const TWeightsImage *weights)
TOutputImage::ConstPointer GetInputProjectionStack()
Projection geometry for a source and a 2-D flat panel.
void SetSupportMask(const TSingleComponentImage *SupportMask)
void SetInputProjectionStack(const TOutputImage *projs)
void SetInputVolume(const TOutputImage *vol)
std::conditional< std::is_same< TSingleComponentImage, TOutputImage >::value, PlainMultiplyFilterType, MatrixVectorMultiplyFilterType >::type MultiplyWithWeightsFilterType
#define itkSetMacro(name, type)