RTK  2.0.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  * 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 
19 #ifndef rtkFourDConjugateGradientConeBeamReconstructionFilter_h
20 #define rtkFourDConjugateGradientConeBeamReconstructionFilter_h
21 
29 
30 #include <itkExtractImageFilter.h>
31 #include <itkSubtractImageFilter.h>
32 #include <itkMultiplyImageFilter.h>
33 #ifdef RTK_USE_CUDA
35 #endif
36 
37 namespace rtk
38 {
39 
94 template<typename VolumeSeriesType, typename ProjectionStackType>
96  public rtk::IterativeConeBeamReconstructionFilter<VolumeSeriesType, ProjectionStackType>
97 {
98 public:
104 
106  typedef VolumeSeriesType InputImageType;
107  typedef VolumeSeriesType OutputImageType;
108  typedef ProjectionStackType VolumeType;
109 
117 
120 
122  typedef typename itk::Image< typename VolumeSeriesType::PixelType,
123  VolumeSeriesType::ImageDimension> CPUVolumeSeriesType;
124 #ifdef RTK_USE_CUDA
125  typedef typename std::conditional<std::is_same< VolumeSeriesType, CPUVolumeSeriesType >::value,
127  CudaConjugateGradientImageFilter<VolumeSeriesType> >::type CudaConjugateGradientImageFilterType;
128 #else
129  typedef ConjugateGradientImageFilter<VolumeSeriesType> CudaConjugateGradientImageFilterType;
130 #endif
131 
133  itkNewMacro(Self)
134 
135 
137 
139  itkGetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry)
140  itkSetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry)
141 
143  itkGetMacro(NumberOfIterations, unsigned int)
144  itkSetMacro(NumberOfIterations, unsigned int)
145 
147  itkGetMacro(CudaConjugateGradient, bool)
148  itkSetMacro(CudaConjugateGradient, bool)
149 
151  void SetInputVolumeSeries(const VolumeSeriesType* VolumeSeries);
152  typename VolumeSeriesType::ConstPointer GetInputVolumeSeries();
154 
156  void SetInputProjectionStack(const ProjectionStackType* Projections);
157  typename ProjectionStackType::ConstPointer GetInputProjectionStack();
159 
161  virtual void SetForwardProjectionFilter (ForwardProjectionType _arg) ITK_OVERRIDE;
162 
164  virtual void SetBackProjectionFilter (BackProjectionType _arg) ITK_OVERRIDE;
165 
167  void SetWeights(const itk::Array2D<float> _arg);
168 
170  virtual void SetSignal(const std::vector<double> signal);
171 
173  itkSetMacro(DisableDisplacedDetectorFilter, bool)
174  itkGetMacro(DisableDisplacedDetectorFilter, bool)
175 protected:
176  FourDConjugateGradientConeBeamReconstructionFilter();
177  virtual ~FourDConjugateGradientConeBeamReconstructionFilter() ITK_OVERRIDE {}
178 
179  void GenerateOutputInformation() ITK_OVERRIDE;
180 
181  void GenerateInputRequestedRegion() ITK_OVERRIDE;
182 
183  void GenerateData() ITK_OVERRIDE;
184 
187 #if ITK_VERSION_MAJOR<5
188  void VerifyInputInformation() ITK_OVERRIDE {}
189 #else
190  void VerifyInputInformation() const ITK_OVERRIDE {}
191 #endif
192 
193 
202 
204  std::vector<double> m_Signal;
206 
207 private:
208  //purposely not implemented
210  void operator=(const Self&);
211 
214 
216  unsigned int m_NumberOfIterations;
217 
218 }; // end of class
219 
220 } // end namespace rtk
221 
222 #ifndef ITK_MANUAL_INSTANTIATION
223 #include "rtkFourDConjugateGradientConeBeamReconstructionFilter.hxx"
224 #endif
225 
226 #endif
FourDReconstructionConjugateGradientOperator< VolumeSeriesType, ProjectionStackType > CGOperatorFilterType
Base class for forward projection, i.e. accumulation along x-ray lines.
Weigting for displaced detectors.
STL namespace.
IterativeConeBeamReconstructionFilter< VolumeSeriesType, ProjectionStackType > Superclass
Projection geometry for a source and a 2-D flat panel.
Implements part of the 4D reconstruction by conjugate gradient.
itk::Image< typename VolumeSeriesType::PixelType, VolumeSeriesType::ImageDimension > CPUVolumeSeriesType
Mother class for cone beam reconstruction filters which need runtime selection of their forward and b...
BackProjectionImageFilter< ProjectionStackType, VolumeType > BackProjectionFilterType
ForwardProjectionImageFilter< VolumeType, ProjectionStackType > ForwardProjectionFilterType
Solves AX = B by conjugate gradient.
Implements part of the 4D reconstruction by conjugate gradient.
ProjectionStackToFourDImageFilter< VolumeSeriesType, ProjectionStackType > ProjStackToFourDFilterType
#define itkSetMacro(name, type)