RTK  1.4.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 rtk::BackProjectionImageFilter< ProjectionStackType, ProjectionStackType > BackProjectionFilterType;
157  typedef rtk::ForwardProjectionImageFilter< ProjectionStackType, ProjectionStackType > ForwardProjectionFilterType;
159  typedef rtk::SplatWithKnownWeightsImageFilter<VolumeSeriesType, VolumeType> SplatFilterType;
164 
166  void SetBackProjectionFilter (const typename BackProjectionFilterType::Pointer _arg);
167 
169  void SetForwardProjectionFilter (const typename ForwardProjectionFilterType::Pointer _arg);
170 
172  itkSetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry)
173 
175  itkSetMacro(UseCudaInterpolation, bool)
176  itkGetMacro(UseCudaInterpolation, bool)
177  itkSetMacro(UseCudaSplat, bool)
178  itkGetMacro(UseCudaSplat, bool)
179  itkSetMacro(UseCudaSources, bool)
180  itkGetMacro(UseCudaSources, bool)
181 
183  itkGetMacro(Weights, itk::Array2D<float>)
184  itkSetMacro(Weights, itk::Array2D<float>)
185 
187  virtual void SetSignal(const std::vector<double> signal);
188 
190  itkSetMacro(DisableDisplacedDetectorFilter, bool)
191  itkGetMacro(DisableDisplacedDetectorFilter, bool)
192 
193 protected:
194  FourDReconstructionConjugateGradientOperator();
195  virtual ~FourDReconstructionConjugateGradientOperator() ITK_OVERRIDE {}
196 
198  void GenerateOutputInformation() ITK_OVERRIDE;
199 
201  void GenerateInputRequestedRegion() ITK_OVERRIDE;
202 
204  void GenerateData() ITK_OVERRIDE;
205 
208 
219 
225  std::vector<double> m_Signal;
227 
228 private:
229  FourDReconstructionConjugateGradientOperator(const Self &); //purposely not implemented
230  void operator=(const Self &); //purposely not implemented
231 };
232 } //namespace ITK
233 
234 
235 #ifndef ITK_MANUAL_INSTANTIATION
236 #include "rtkFourDReconstructionConjugateGradientOperator.hxx"
237 #endif
238 
239 #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()
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.
void SetBackProjectionFilter(const typename BackProjectionFilterType::Pointer _arg)
Interpolates (linearly) in a 3D+t sequence of volumes to get a 3D volume.
#define itkSetMacro(name, type)