RTK  1.4.0
Reconstruction Toolkit
rtkADMMTotalVariationConeBeamReconstructionFilter.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 rtkADMMTotalVariationConeBeamReconstructionFilter_h
20 #define rtkADMMTotalVariationConeBeamReconstructionFilter_h
21 
22 #include <itkImageToImageFilter.h>
23 #include <itkAddImageFilter.h>
24 #include <itkSubtractImageFilter.h>
25 #include <itkMultiplyImageFilter.h>
26 
34 
35 namespace rtk
36 {
131 template< typename TOutputImage, typename TGradientOutputImage =
133  TOutputImage::ImageDimension > >
135 {
136 public:
141 
144 
146  itkNewMacro(Self)
147 
148 
150 
151  typedef rtk::ForwardProjectionImageFilter< TOutputImage, TOutputImage > ForwardProjectionFilterType;
152  typedef typename ForwardProjectionFilterType::Pointer ForwardProjectionFilterPointer;
153  typedef rtk::BackProjectionImageFilter< TOutputImage, TOutputImage > BackProjectionFilterType;
154  typedef typename BackProjectionFilterType::Pointer BackProjectionFilterPointer;
156  typedef ForwardDifferenceGradientImageFilter<TOutputImage,
157  typename TOutputImage::ValueType,
158  typename TOutputImage::ValueType,
159  TGradientOutputImage> ImageGradientFilterType;
161  <TGradientOutputImage, TOutputImage> ImageDivergenceFilterType;
163  <TGradientOutputImage> SoftThresholdTVFilterType;
165  typedef itk::AddImageFilter<TGradientOutputImage> AddGradientsFilterType;
167  typedef itk::MultiplyImageFilter<TGradientOutputImage> MultiplyGradientFilterType;
168  typedef itk::SubtractImageFilter<TGradientOutputImage> SubtractGradientsFilterType;
172 
174  void SetForwardProjectionFilter (ForwardProjectionType _arg) ITK_OVERRIDE;
175 
177  void SetBackProjectionFilter (BackProjectionType _arg) ITK_OVERRIDE;
178 
180  itkSetObjectMacro(Geometry, ThreeDCircularProjectionGeometry)
181 
183  void SetBetaForCurrentIteration(int iter);
184 
186  void SetGatingWeights(std::vector<float> weights);
187 
188  itkSetMacro(Alpha, float)
189  itkGetMacro(Alpha, float)
190 
191  itkSetMacro(Beta, float)
192  itkGetMacro(Beta, float)
193 
194  itkSetMacro(AL_iterations, float)
195  itkGetMacro(AL_iterations, float)
196 
197  itkSetMacro(CG_iterations, float)
198  itkGetMacro(CG_iterations, float)
199 
201  itkSetMacro(DisableDisplacedDetectorFilter, bool)
202  itkGetMacro(DisableDisplacedDetectorFilter, bool)
203 
204 protected:
205  ADMMTotalVariationConeBeamReconstructionFilter();
206  virtual ~ADMMTotalVariationConeBeamReconstructionFilter() ITK_OVERRIDE {}
207 
209  void GenerateData() ITK_OVERRIDE;
210 
230 
234 #if ITK_VERSION_MAJOR<5
235  void VerifyInputInformation() ITK_OVERRIDE {}
236 #else
237  void VerifyInputInformation() const ITK_OVERRIDE {}
238 #endif
239 
240 
243  void GenerateInputRequestedRegion() ITK_OVERRIDE;
244  void GenerateOutputInformation() ITK_OVERRIDE;
246 
249  bool m_IsGated;
250  std::vector<float> m_GatingWeights;
252 
253 private:
254  ADMMTotalVariationConeBeamReconstructionFilter(const Self &); //purposely not implemented
255  void operator=(const Self &); //purposely not implemented
256 
257  float m_Alpha;
258  float m_Beta;
259  unsigned int m_AL_iterations;
260  unsigned int m_CG_iterations;
261 
263 };
264 } //namespace ITK
265 
266 
267 #ifndef ITK_MANUAL_INSTANTIATION
268 #include "rtkADMMTotalVariationConeBeamReconstructionFilter.hxx"
269 #endif
270 
271 #endif
void SetGatingWeights(std::vector< float > weights)
Base class for forward projection, i.e. accumulation along x-ray lines.
Weigting for displaced detectors.
STL namespace.
void SetBackProjectionFilter(BackProjectionType _arg) override
Projection geometry for a source and a 2-D flat panel.
BackProjectionImageFilter< TOutputImage, TOutputImage >::Pointer m_BackProjectionFilter
ForwardProjectionImageFilter< TOutputImage, TOutputImage >::Pointer m_ForwardProjectionFilter
Computes the gradient of an image using forward difference.
void SetForwardProjectionFilter(ForwardProjectionType _arg) override
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.
IterativeConeBeamReconstructionFilter< TOutputImage, TOutputImage > Superclass
Mother class for cone beam reconstruction filters which need runtime selection of their forward and b...
Computes the Total Variation from a gradient input image (pixels are vectors), soft thresholds it...
Computes the backward differences divergence (adjoint of the forward differences gradient) of the inp...
Implements the ADMM reconstruction with total variation regularization.
Solves AX = B by conjugate gradient.
BackProjectionImageFilter< TOutputImage, TOutputImage >::Pointer m_BackProjectionFilterForConjugateGradient
#define itkSetMacro(name, type)