RTK  2.7.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  * 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 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, typename TSingleComponentImage = TOutputImage, typename TWeightsImage = TOutputImage>
118 class ITK_TEMPLATE_EXPORT ReconstructionConjugateGradientOperator : public ConjugateGradientOperator<TOutputImage>
119 {
120 public:
121  ITK_DISALLOW_COPY_AND_MOVE(ReconstructionConjugateGradientOperator);
122 
128  using GradientImageType =
129  typename TOutputImage::template RebindImageType<VectorPixelType, TOutputImage::ImageDimension>;
130 
132  itkNewMacro(Self);
133 
135  void
136  SetInputVolume(const TOutputImage * vol);
137  void
138  SetInputProjectionStack(const TOutputImage * projs);
139  void
140  SetInputWeights(const TWeightsImage * weights);
142 
144  itkOverrideGetNameOfClassMacro(ReconstructionConjugateGradientOperator);
145 
148 
151 
155 
156  // If TOutputImage is an itk::Image of floats or double, so are the weights, and a simple Multiply filter is required
157  // If TOutputImage is an itk::Image of itk::Vector<float (or double)>, a BlockDiagonalMatrixVectorMultiply filter
158  // is needed. Thus the meta-programming construct
161  using MultiplyWithWeightsFilterType = typename std::conditional_t<std::is_same_v<TSingleComponentImage, TOutputImage>,
164 
165  using OutputImagePointer = typename TOutputImage::Pointer;
166 
168  void
169  SetBackProjectionFilter(BackProjectionFilterType * _arg);
170 
172  void
173  SetForwardProjectionFilter(ForwardProjectionFilterType * _arg);
174 
176  void
177  SetSupportMask(const TSingleComponentImage * SupportMask);
178  typename TSingleComponentImage::ConstPointer
179  GetSupportMask();
181 
185  void
186  SetLocalRegularizationWeights(const TSingleComponentImage * localRegularizationWeights);
187  typename TSingleComponentImage::ConstPointer
188  GetLocalRegularizationWeights();
190 
192  itkSetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry);
193 
196  itkSetMacro(Gamma, float);
197  itkGetMacro(Gamma, float);
198  itkSetMacro(Tikhonov, float);
199  itkGetMacro(Tikhonov, float);
201 
202 protected:
204  ~ReconstructionConjugateGradientOperator() override = default;
205 
207  void
208  VerifyPreconditions() const override;
209 
211  void
212  GenerateData() override;
213 
214  template <typename ImageType>
215  typename std::enable_if<std::is_same_v<TSingleComponentImage, ImageType>, ImageType>::type::Pointer
216  ConnectGradientRegularization();
217 
218  template <typename ImageType>
219  typename std::enable_if<!std::is_same_v<TSingleComponentImage, ImageType>, ImageType>::type::Pointer
220  ConnectGradientRegularization();
221 
225 
236  typename MultiplyWithWeightsFilterType::Pointer m_MultiplyWithWeightsFilter;
237 
240  float m_Gamma{ 0 }; // Strength of the laplacian regularization
241  float m_Tikhonov{ 0 }; // Strength of the Tikhonov regularization
242 
244  typename TOutputImage::Pointer m_FloatingInputPointer, m_FloatingOutputPointer;
245 
248  void
249  VerifyInputInformation() const override
250  {}
251 
253  void
254  GenerateInputRequestedRegion() override;
255  void
256  GenerateOutputInformation() override;
258 
260  typename TOutputImage::ConstPointer
261  GetInputVolume();
262  typename TOutputImage::ConstPointer
263  GetInputProjectionStack();
264  typename TWeightsImage::ConstPointer
265  GetInputWeights();
266 };
267 } // namespace rtk
269 
270 
271 #ifndef ITK_MANUAL_INSTANTIATION
272 # include "rtkReconstructionConjugateGradientOperator.hxx"
273 #endif
274 
275 #endif
Implements the operator A used in conjugate gradient reconstruction.
Generate an n-dimensional image with constant pixel values.
typename BackProjectionFilterType::Pointer BackProjectionFilterPointer
Projection geometry for a source and a 2-D flat panel.
#define itkSetMacro(name, type)
typename ForwardProjectionFilterType::Pointer ForwardProjectionFilterPointer
typename TOutputImage::template RebindImageType< VectorPixelType, TOutputImage::ImageDimension > GradientImageType
itk::ImageToImageFilter< TOutputImage, TOutputImage >::Pointer m_LaplacianFilter
typename std::conditional_t< std::is_same_v< TSingleComponentImage, TOutputImage >, PlainMultiplyFilterType, MatrixVectorMultiplyFilterType > MultiplyWithWeightsFilterType