RTK  2.6.0
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  * 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 rtkFourDConjugateGradientConeBeamReconstructionFilter_h
20 #define rtkFourDConjugateGradientConeBeamReconstructionFilter_h
21 
29 
30 #include <itkExtractImageFilter.h>
31 #include <itkSubtractImageFilter.h>
32 #include <itkMultiplyImageFilter.h>
33 #include <itkIterationReporter.h>
34 #ifdef RTK_USE_CUDA
36 #endif
37 
38 namespace rtk
39 {
40 
95 template <typename VolumeSeriesType, typename ProjectionStackType>
97  : public rtk::IterativeConeBeamReconstructionFilter<VolumeSeriesType, ProjectionStackType>
98 {
99 public:
100  ITK_DISALLOW_COPY_AND_MOVE(FourDConjugateGradientConeBeamReconstructionFilter);
101 
107 
109  using InputImageType = VolumeSeriesType;
110  using OutputImageType = VolumeSeriesType;
111  using VolumeType = ProjectionStackType;
112 
120 
121  using ForwardProjectionType = typename Superclass::ForwardProjectionType;
122  using BackProjectionType = typename Superclass::BackProjectionType;
123 
125  using CPUVolumeSeriesType =
127 #ifdef RTK_USE_CUDA
128  typedef typename std::conditional<std::is_same<VolumeSeriesType, CPUVolumeSeriesType>::value,
130  CudaConjugateGradientImageFilter<VolumeSeriesType>>::type
132 #else
134 #endif
135 
137  itkNewMacro(Self);
138 
140 #ifdef itkOverrideGetNameOfClassMacro
141  itkOverrideGetNameOfClassMacro(FourDConjugateGradientConeBeamReconstructionFilter);
142 #else
144 #endif
145 
146 
148  itkGetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry);
149  itkSetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry);
151 
153  itkGetMacro(NumberOfIterations, unsigned int);
154  itkSetMacro(NumberOfIterations, unsigned int);
156 
158  itkGetMacro(CudaConjugateGradient, bool);
159  itkSetMacro(CudaConjugateGradient, bool);
161 
163  void
164  SetInputVolumeSeries(const VolumeSeriesType * VolumeSeries);
165  typename VolumeSeriesType::ConstPointer
166  GetInputVolumeSeries();
168 
170  void
171  SetInputProjectionStack(const ProjectionStackType * Projections);
172  typename ProjectionStackType::ConstPointer
173  GetInputProjectionStack();
175 
177  void
178  SetWeights(const itk::Array2D<float> _arg);
179 
181  virtual void
182  SetSignal(const std::vector<double> signal);
183 
185  itkSetMacro(DisableDisplacedDetectorFilter, bool);
186  itkGetMacro(DisableDisplacedDetectorFilter, bool);
188 
189 protected:
192 
194  void
195  VerifyPreconditions() const override;
196 
197  void
198  GenerateOutputInformation() override;
199 
200  void
201  GenerateInputRequestedRegion() override;
202 
203  void
204  GenerateData() override;
205 
208  void
209  VerifyInputInformation() const override
210  {}
211 
220 
222  std::vector<double> m_Signal;
224 
225  // Iteration reporting
227 
228 private:
231 
233  unsigned int m_NumberOfIterations;
234 
236  void
237  ReportProgress(itk::Object *, const itk::EventObject &);
238 
239 }; // end of class
240 
241 } // end namespace rtk
242 
243 #ifndef ITK_MANUAL_INSTANTIATION
244 # include "rtkFourDConjugateGradientConeBeamReconstructionFilter.hxx"
245 #endif
246 
247 #endif
Base class for forward projection, i.e. accumulation along x-ray lines.
Weigting for displaced detectors.
typename itk::Image< typename VolumeSeriesType::PixelType, VolumeSeriesType::ImageDimension > CPUVolumeSeriesType
Projection geometry for a source and a 2-D flat panel.
#define itkSetMacro(name, type)
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.