RTK  2.7.0
Reconstruction Toolkit
rtkADMMTotalVariationConjugateGradientOperator.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 rtkADMMTotalVariationConjugateGradientOperator_h
20 #define rtkADMMTotalVariationConjugateGradientOperator_h
21 
22 #include <itkMultiplyImageFilter.h>
23 #include <itkSubtractImageFilter.h>
24 
33 
34 #ifdef RTK_USE_CUDA
35 # include <itkCudaImage.h>
36 #endif
37 
38 namespace rtk
39 {
40 
107 template <typename TOutputImage>
108 class ITK_TEMPLATE_EXPORT ADMMTotalVariationConjugateGradientOperator : public ConjugateGradientOperator<TOutputImage>
109 {
110 public:
111  ITK_DISALLOW_COPY_AND_MOVE(ADMMTotalVariationConjugateGradientOperator);
112 
117 
119  itkNewMacro(Self);
120 
122  itkOverrideGetNameOfClassMacro(ADMMTotalVariationConjugateGradientOperator);
123 
126 
129 
132 
135 #ifdef RTK_USE_CUDA
136  typedef
137  typename std::conditional<std::is_same<TOutputImage, CPUImageType>::value,
139  itk::CudaImage<VectorPixelType, TOutputImage::ImageDimension>>::type GradientImageType;
140 #else
142 #endif
143 
145  typename TOutputImage::ValueType,
146  typename TOutputImage::ValueType,
151 
153  void
154  SetBackProjectionFilter(BackProjectionFilterType * _arg);
155 
157  void
158  SetForwardProjectionFilter(ForwardProjectionFilterType * _arg);
159 
161  itkSetObjectMacro(Geometry, ThreeDCircularProjectionGeometry);
162 
164  itkSetMacro(Beta, float);
165 
167  void
168  SetGatingWeights(std::vector<float> weights);
169 
171  itkSetMacro(DisableDisplacedDetectorFilter, bool);
172  itkGetMacro(DisableDisplacedDetectorFilter, bool);
174 
175 protected:
177  ~ADMMTotalVariationConjugateGradientOperator() override = default;
178 
180  void
181  VerifyPreconditions() const override;
182 
184  void
185  GenerateData() override;
186 
190 
199 
201  float m_Beta;
203 
206  bool m_IsGated;
207  std::vector<float> m_GatingWeights;
208 
212  void
213  VerifyInputInformation() const override
214  {}
215 
218  void
219  GenerateInputRequestedRegion() override;
220  void
221  GenerateOutputInformation() override;
222 };
223 } // namespace rtk
225 
226 
227 #ifndef ITK_MANUAL_INSTANTIATION
228 # include "rtkADMMTotalVariationConjugateGradientOperator.hxx"
229 #endif
230 
231 #endif
Weigting for displaced detectors.
typename ForwardProjectionFilterType::Pointer ForwardProjectionFilterPointer
Projection geometry for a source and a 2-D flat panel.
#define itkSetMacro(name, type)
Computes the gradient of an image using forward difference.
Implements the operator A used in the conjugate gradient step of ADMM reconstruction with total varia...
Multiplies each (n-1) dimension image by the corresponding element in a vector.
std::conditional< std::is_same< TOutputImage, CPUImageType >::value, itk::Image< VectorPixelType, TOutputImage::ImageDimension >, itk::CudaImage< VectorPixelType, TOutputImage::ImageDimension > >::type GradientImageType
Computes the backward differences divergence (adjoint of the forward differences gradient) of the inp...