RTK  2.0.1
Reconstruction Toolkit
rtkFourDConjugateGradientConeBeamReconstructionFilter.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 rtkFourDConjugateGradientConeBeamReconstructionFilter_h
20 #define rtkFourDConjugateGradientConeBeamReconstructionFilter_h
21 
29 
30 #include <itkExtractImageFilter.h>
31 #include <itkSubtractImageFilter.h>
32 #include <itkMultiplyImageFilter.h>
33 #ifdef RTK_USE_CUDA
35 #endif
36 
37 namespace rtk
38 {
39 
94 template<typename VolumeSeriesType, typename ProjectionStackType>
96  public rtk::IterativeConeBeamReconstructionFilter<VolumeSeriesType, ProjectionStackType>
97 {
98 public:
99  ITK_DISALLOW_COPY_AND_ASSIGN(FourDConjugateGradientConeBeamReconstructionFilter);
100 
106 
108  using InputImageType = VolumeSeriesType;
109  using OutputImageType = VolumeSeriesType;
110  using VolumeType = ProjectionStackType;
111 
119 
120  using ForwardProjectionType = typename Superclass::ForwardProjectionType;
121  using BackProjectionType = typename Superclass::BackProjectionType;
122 
124  using CPUVolumeSeriesType = typename itk::Image< typename VolumeSeriesType::PixelType,
125  VolumeSeriesType::ImageDimension>;
126 #ifdef RTK_USE_CUDA
127  typedef typename std::conditional<std::is_same< VolumeSeriesType, CPUVolumeSeriesType >::value,
129  CudaConjugateGradientImageFilter<VolumeSeriesType> >::type CudaConjugateGradientImageFilterType;
130 #else
131  using CudaConjugateGradientImageFilterType = ConjugateGradientImageFilter<VolumeSeriesType>;
132 #endif
133 
135  itkNewMacro(Self)
136 
137 
139 
141  itkGetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry)
142  itkSetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry)
143 
145  itkGetMacro(NumberOfIterations, unsigned int)
146  itkSetMacro(NumberOfIterations, unsigned int)
147 
149  itkGetMacro(CudaConjugateGradient, bool)
150  itkSetMacro(CudaConjugateGradient, bool)
151 
153  void SetInputVolumeSeries(const VolumeSeriesType* VolumeSeries);
154  typename VolumeSeriesType::ConstPointer GetInputVolumeSeries();
156 
158  void SetInputProjectionStack(const ProjectionStackType* Projections);
159  typename ProjectionStackType::ConstPointer GetInputProjectionStack();
161 
163  void SetForwardProjectionFilter (ForwardProjectionType _arg) override;
164 
166  void SetBackProjectionFilter (BackProjectionType _arg) override;
167 
169  void SetWeights(const itk::Array2D<float> _arg);
170 
172  virtual void SetSignal(const std::vector<double> signal);
173 
175  itkSetMacro(DisableDisplacedDetectorFilter, bool)
176  itkGetMacro(DisableDisplacedDetectorFilter, bool)
177 protected:
178  FourDConjugateGradientConeBeamReconstructionFilter();
179  ~FourDConjugateGradientConeBeamReconstructionFilter() override = default;
181 
182  void GenerateOutputInformation() override;
183 
184  void GenerateInputRequestedRegion() override;
185 
186  void GenerateData() override;
187 
190 #if ITK_VERSION_MAJOR<5
191  void VerifyInputInformation() override {}
192 #else
193  void VerifyInputInformation() const override {}
194 #endif
195 
196 
205 
207  std::vector<double> m_Signal;
209 
210 private:
213 
215  unsigned int m_NumberOfIterations;
216 
217 }; // end of class
218 
219 } // end namespace rtk
220 
221 #ifndef ITK_MANUAL_INSTANTIATION
222 #include "rtkFourDConjugateGradientConeBeamReconstructionFilter.hxx"
223 #endif
224 
225 #endif
Base class for forward projection, i.e. accumulation along x-ray lines.
typename itk::Image< typename VolumeSeriesType::PixelType, VolumeSeriesType::ImageDimension > CPUVolumeSeriesType
Weigting for displaced detectors.
STL namespace.
Projection geometry for a source and a 2-D flat panel.
Implements part of the 4D reconstruction by conjugate gradient.
TOutputImage OutputImageType
Mother class for cone beam reconstruction filters which need runtime selection of their forward and b...
Solves AX = B by conjugate gradient.
Implements part of the 4D reconstruction by conjugate gradient.
#define itkSetMacro(name, type)