RTK  2.0.1
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  * 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 #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 {
128 template< typename VolumeSeriesType, typename ProjectionStackType>
130 {
131 public:
132  ITK_DISALLOW_COPY_AND_ASSIGN(FourDReconstructionConjugateGradientOperator);
133 
138 
140  using VolumeType = ProjectionStackType;
141 
143  itkNewMacro(Self)
144 
145 
147 
149  void SetInputVolumeSeries(const VolumeSeriesType* VolumeSeries);
150  typename VolumeSeriesType::ConstPointer GetInputVolumeSeries();
152 
154  void SetInputProjectionStack(const ProjectionStackType* Projections);
155  typename ProjectionStackType::ConstPointer GetInputProjectionStack();
157 
158  using BackProjectionFilterType = BackProjectionImageFilter< ProjectionStackType, ProjectionStackType >;
159  using ForwardProjectionFilterType = ForwardProjectionImageFilter< ProjectionStackType, ProjectionStackType >;
161  using SplatFilterType = SplatWithKnownWeightsImageFilter<VolumeSeriesType, VolumeType>;
165 
167  using CPUProjectionStackType = typename itk::Image< typename ProjectionStackType::PixelType,
168  ProjectionStackType::ImageDimension>;
169 #ifdef RTK_USE_CUDA
170  typedef typename std::conditional< std::is_same< ProjectionStackType, CPUProjectionStackType >::value,
172  CudaDisplacedDetectorImageFilter >::type DisplacedDetectorFilterType;
173  typedef typename std::conditional< std::is_same< ProjectionStackType, CPUProjectionStackType >::value,
175  CudaInterpolateImageFilter >::type CudaInterpolateImageFilterType;
176  typedef typename std::conditional< std::is_same< ProjectionStackType, CPUProjectionStackType >::value,
178  CudaSplatImageFilter >::type CudaSplatImageFilterType;
179  typedef typename std::conditional< std::is_same< ProjectionStackType, CPUProjectionStackType >::value,
181  CudaConstantVolumeSource >::type CudaConstantVolumeSourceType;
182  typedef typename std::conditional< std::is_same< ProjectionStackType, CPUProjectionStackType >::value,
184  CudaConstantVolumeSeriesSource >::type CudaConstantVolumeSeriesSourceType;
185 #else
186  using DisplacedDetectorFilterType = DisplacedDetectorImageFilter<ProjectionStackType>;
187  using CudaInterpolateImageFilterType = InterpolationFilterType;
188  using CudaSplatImageFilterType = SplatFilterType;
189  using CudaConstantVolumeSourceType = ConstantVolumeSourceType;
190  using CudaConstantVolumeSeriesSourceType = ConstantVolumeSeriesSourceType;
191 #endif
192 
195 
198 
200  itkSetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry)
201 
202 
203  itkSetMacro(UseCudaInterpolation, bool)
204  itkGetMacro(UseCudaInterpolation, bool)
205  itkSetMacro(UseCudaSplat, bool)
206  itkGetMacro(UseCudaSplat, bool)
207  itkSetMacro(UseCudaSources, bool)
208  itkGetMacro(UseCudaSources, bool)
209 
211  itkGetMacro(Weights, itk::Array2D<float>)
212  itkSetMacro(Weights, itk::Array2D<float>)
213 
215  virtual void SetSignal(const std::vector<double> signal);
216 
218  itkSetMacro(DisableDisplacedDetectorFilter, bool)
219  itkGetMacro(DisableDisplacedDetectorFilter, bool)
220 
221 protected:
222  FourDReconstructionConjugateGradientOperator();
223  ~FourDReconstructionConjugateGradientOperator() override = default;
224 
226  void GenerateOutputInformation() override;
227 
229  void GenerateInputRequestedRegion() override;
230 
232  void GenerateData() override;
233 
236 
240  typename InterpolationFilterType::Pointer m_InterpolationFilter;
241  typename SplatFilterType::Pointer m_SplatFilter;
242  typename ConstantVolumeSourceType::Pointer m_ConstantVolumeSource1;
243  typename ConstantVolumeSourceType::Pointer m_ConstantVolumeSource2;
245  typename ConstantVolumeSeriesSourceType::Pointer m_ConstantVolumeSeriesSource;
246  typename DisplacedDetectorFilterType::Pointer m_DisplacedDetectorFilter;
247 
253  std::vector<double> m_Signal;
255 };
256 } //namespace ITK
257 
258 
259 #ifndef ITK_MANUAL_INSTANTIATION
260 #include "rtkFourDReconstructionConjugateGradientOperator.hxx"
261 #endif
262 
263 #endif
Base class for forward projection, i.e. accumulation along x-ray lines.
VolumeSeriesType::ConstPointer GetInputVolumeSeries()
Generate an n-dimensional image with constant pixel values.
Weigting for displaced detectors.
SplatWithKnownWeightsImageFilter< VolumeSeriesType, VolumeType > SplatFilterType
STL namespace.
void SetInputVolumeSeries(const VolumeSeriesType *VolumeSeries)
ProjectionStackType::ConstPointer GetInputProjectionStack()
Splats (linearly) a 3D volume into a 3D+t sequence of volumes.
void SetForwardProjectionFilter(const typename ForwardProjectionFilterType::Pointer _arg)
void SetInputProjectionStack(const ProjectionStackType *Projections)
virtual void SetSignal(const std::vector< double > signal)
Projection geometry for a source and a 2-D flat panel.
InterpolatorWithKnownWeightsImageFilter< VolumeType, VolumeSeriesType > InterpolationFilterType
Implements part of the 4D reconstruction by conjugate gradient.
void SetBackProjectionFilter(const typename BackProjectionFilterType::Pointer _arg)
DisplacedDetectorImageFilter< ProjectionStackType > DisplacedDetectorFilterType
typename itk::Image< typename ProjectionStackType::PixelType, ProjectionStackType::ImageDimension > CPUProjectionStackType
Interpolates (linearly) in a 3D+t sequence of volumes to get a 3D volume.
#define itkSetMacro(name, type)