RTK  2.6.0
Reconstruction Toolkit
rtkFourDReconstructionConjugateGradientOperator.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 #ifndef rtkFourDReconstructionConjugateGradientOperator_h
19 #define rtkFourDReconstructionConjugateGradientOperator_h
20 
22 
23 #include <itkArray2D.h>
24 #include <itkMultiplyImageFilter.h>
25 
26 #include "rtkConstantImageSource.h"
33 
34 #ifdef RTK_USE_CUDA
36 # include "rtkCudaSplatImageFilter.h"
40 #endif
41 
42 namespace rtk
43 {
129 template <typename VolumeSeriesType, typename ProjectionStackType>
131  : public ConjugateGradientOperator<VolumeSeriesType>
132 {
133 public:
134  ITK_DISALLOW_COPY_AND_MOVE(FourDReconstructionConjugateGradientOperator);
135 
140 
142  using VolumeType = ProjectionStackType;
143 
145  itkNewMacro(Self);
146 
148 #ifdef itkOverrideGetNameOfClassMacro
149  itkOverrideGetNameOfClassMacro(FourDReconstructionConjugateGradientOperator);
150 #else
152 #endif
153 
154 
156  void
157  SetInputVolumeSeries(const VolumeSeriesType * VolumeSeries);
158  typename VolumeSeriesType::ConstPointer
159  GetInputVolumeSeries();
161 
163  void
164  SetInputProjectionStack(const ProjectionStackType * Projections);
165  typename ProjectionStackType::ConstPointer
166  GetInputProjectionStack();
168 
176 
178  using CPUProjectionStackType =
180 #ifdef RTK_USE_CUDA
181  typedef typename std::conditional<std::is_same<ProjectionStackType, CPUProjectionStackType>::value,
183  CudaDisplacedDetectorImageFilter>::type DisplacedDetectorFilterType;
184  typedef typename std::conditional<std::is_same<ProjectionStackType, CPUProjectionStackType>::value,
186  CudaInterpolateImageFilter>::type CudaInterpolateImageFilterType;
187  typedef typename std::conditional<std::is_same<ProjectionStackType, CPUProjectionStackType>::value,
189  CudaSplatImageFilter>::type CudaSplatImageFilterType;
190  typedef typename std::conditional<std::is_same<ProjectionStackType, CPUProjectionStackType>::value,
192  CudaConstantVolumeSource>::type CudaConstantVolumeSourceType;
193  typedef typename std::conditional<std::is_same<ProjectionStackType, CPUProjectionStackType>::value,
195  CudaConstantVolumeSeriesSource>::type CudaConstantVolumeSeriesSourceType;
196 #else
202 #endif
203 
205  void
206  SetBackProjectionFilter(const typename BackProjectionFilterType::Pointer _arg);
207 
209  void
210  SetForwardProjectionFilter(const typename ForwardProjectionFilterType::Pointer _arg);
211 
213  itkSetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry);
214 
216  itkSetMacro(UseCudaInterpolation, bool);
217  itkGetMacro(UseCudaInterpolation, bool);
218  itkSetMacro(UseCudaSplat, bool);
219  itkGetMacro(UseCudaSplat, bool);
220  itkSetMacro(UseCudaSources, bool);
221  itkGetMacro(UseCudaSources, bool);
223 
225  itkGetMacro(Weights, itk::Array2D<float>);
228 
230  virtual void
231  SetSignal(const std::vector<double> signal);
232 
234  itkSetMacro(DisableDisplacedDetectorFilter, bool);
235  itkGetMacro(DisableDisplacedDetectorFilter, bool);
237 
238 protected:
240  ~FourDReconstructionConjugateGradientOperator() override = default;
241 
243  void
244  VerifyPreconditions() const override;
245 
247  void
248  GenerateOutputInformation() override;
249 
251  void
252  GenerateInputRequestedRegion() override;
253 
255  void
256  GenerateData() override;
257 
259  void
260  InitializeConstantSources();
261 
272 
278  std::vector<double> m_Signal;
280 };
281 } // namespace rtk
282 
283 
284 #ifndef ITK_MANUAL_INSTANTIATION
285 # include "rtkFourDReconstructionConjugateGradientOperator.hxx"
286 #endif
287 
288 #endif
Base class for forward projection, i.e. accumulation along x-ray lines.
typename itk::Image< typename ProjectionStackType::PixelType, ProjectionStackType::ImageDimension > CPUProjectionStackType
Generate an n-dimensional image with constant pixel values.
Weigting for displaced detectors.
Splats (linearly) a 3D volume into a 3D+t sequence of volumes.
Projection geometry for a source and a 2-D flat panel.
#define itkSetMacro(name, type)
Implements part of the 4D reconstruction by conjugate gradient.
Interpolates (linearly) in a 3D+t sequence of volumes to get a 3D volume.