RTK  2.0.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  * 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:
136 
138  typedef ProjectionStackType VolumeType;
139 
141  itkNewMacro(Self)
142 
143 
145 
147  void SetInputVolumeSeries(const VolumeSeriesType* VolumeSeries);
148  typename VolumeSeriesType::ConstPointer GetInputVolumeSeries();
150 
152  void SetInputProjectionStack(const ProjectionStackType* Projections);
153  typename ProjectionStackType::ConstPointer GetInputProjectionStack();
155 
156  typedef BackProjectionImageFilter< ProjectionStackType, ProjectionStackType > BackProjectionFilterType;
157  typedef ForwardProjectionImageFilter< ProjectionStackType, ProjectionStackType > ForwardProjectionFilterType;
159  typedef SplatWithKnownWeightsImageFilter<VolumeSeriesType, VolumeType> SplatFilterType;
163 
165  typedef typename itk::Image< typename ProjectionStackType::PixelType,
166  ProjectionStackType::ImageDimension> CPUProjectionStackType;
167 #ifdef RTK_USE_CUDA
168  typedef typename std::conditional< std::is_same< ProjectionStackType, CPUProjectionStackType >::value,
170  CudaDisplacedDetectorImageFilter >::type DisplacedDetectorFilterType;
171  typedef typename std::conditional< std::is_same< ProjectionStackType, CPUProjectionStackType >::value,
173  CudaInterpolateImageFilter >::type CudaInterpolateImageFilterType;
174  typedef typename std::conditional< std::is_same< ProjectionStackType, CPUProjectionStackType >::value,
176  CudaSplatImageFilter >::type CudaSplatImageFilterType;
177  typedef typename std::conditional< std::is_same< ProjectionStackType, CPUProjectionStackType >::value,
179  CudaConstantVolumeSource >::type CudaConstantVolumeSourceType;
180  typedef typename std::conditional< std::is_same< ProjectionStackType, CPUProjectionStackType >::value,
182  CudaConstantVolumeSeriesSource >::type CudaConstantVolumeSeriesSourceType;
183 #else
184  typedef DisplacedDetectorImageFilter<ProjectionStackType> DisplacedDetectorFilterType;
185  typedef InterpolationFilterType CudaInterpolateImageFilterType;
186  typedef SplatFilterType CudaSplatImageFilterType;
187  typedef ConstantVolumeSourceType CudaConstantVolumeSourceType;
188  typedef ConstantVolumeSeriesSourceType CudaConstantVolumeSeriesSourceType;
189 #endif
190 
193 
196 
198  itkSetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry)
199 
200 
201  itkSetMacro(UseCudaInterpolation, bool)
202  itkGetMacro(UseCudaInterpolation, bool)
203  itkSetMacro(UseCudaSplat, bool)
204  itkGetMacro(UseCudaSplat, bool)
205  itkSetMacro(UseCudaSources, bool)
206  itkGetMacro(UseCudaSources, bool)
207 
209  itkGetMacro(Weights, itk::Array2D<float>)
210  itkSetMacro(Weights, itk::Array2D<float>)
211 
213  virtual void SetSignal(const std::vector<double> signal);
214 
216  itkSetMacro(DisableDisplacedDetectorFilter, bool)
217  itkGetMacro(DisableDisplacedDetectorFilter, bool)
218 
219 protected:
220  FourDReconstructionConjugateGradientOperator();
221  virtual ~FourDReconstructionConjugateGradientOperator() ITK_OVERRIDE {}
222 
224  void GenerateOutputInformation() ITK_OVERRIDE;
225 
227  void GenerateInputRequestedRegion() ITK_OVERRIDE;
228 
230  void GenerateData() ITK_OVERRIDE;
231 
234 
245 
251  std::vector<double> m_Signal;
253 
254 private:
255  FourDReconstructionConjugateGradientOperator(const Self &); //purposely not implemented
256  void operator=(const Self &); //purposely not implemented
257 };
258 } //namespace ITK
259 
260 
261 #ifndef ITK_MANUAL_INSTANTIATION
262 #include "rtkFourDReconstructionConjugateGradientOperator.hxx"
263 #endif
264 
265 #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.
STL namespace.
void SetInputVolumeSeries(const VolumeSeriesType *VolumeSeries)
ProjectionStackType::ConstPointer GetInputProjectionStack()
SplatWithKnownWeightsImageFilter< VolumeSeriesType, VolumeType > SplatFilterType
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.
Implements part of the 4D reconstruction by conjugate gradient.
InterpolatorWithKnownWeightsImageFilter< VolumeType, VolumeSeriesType > InterpolationFilterType
void SetBackProjectionFilter(const typename BackProjectionFilterType::Pointer _arg)
Interpolates (linearly) in a 3D+t sequence of volumes to get a 3D volume.
DisplacedDetectorImageFilter< ProjectionStackType > DisplacedDetectorFilterType
#define itkSetMacro(name, type)