RTK  1.4.0
Reconstruction Toolkit
rtkFourDSARTConeBeamReconstructionFilter.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 rtkFourDSARTConeBeamReconstructionFilter_h
20 #define rtkFourDSARTConeBeamReconstructionFilter_h
21 
24 
25 #include <itkExtractImageFilter.h>
26 #include <itkMultiplyImageFilter.h>
27 #include <itkSubtractImageFilter.h>
28 #include <itkAddImageAdaptor.h>
29 #include <itkAddImageFilter.h>
32 
34 #include "rtkConstantImageSource.h"
39 
40 namespace rtk
41 {
42 
119 template<class VolumeSeriesType, class ProjectionStackType>
121  public rtk::IterativeConeBeamReconstructionFilter<VolumeSeriesType, ProjectionStackType>
122 {
123 public:
129 
131  typedef VolumeSeriesType InputImageType;
132  typedef VolumeSeriesType OutputImageType;
133  typedef ProjectionStackType VolumeType;
134 
135  typedef typename Superclass::ForwardProjectionType ForwardProjectionType;
136  typedef typename Superclass::BackProjectionType BackProjectionType;
137 
153 
155  itkNewMacro(Self);
156 
159 
161  void SetInputVolumeSeries(const VolumeSeriesType* VolumeSeries);
162  typename VolumeSeriesType::ConstPointer GetInputVolumeSeries();
164 
166  void SetInputProjectionStack(const ProjectionStackType* Projection);
167  typename ProjectionStackType::Pointer GetInputProjectionStack();
169 
171  itkGetModifiableObjectMacro(Geometry, ThreeDCircularProjectionGeometry);
172  itkSetObjectMacro(Geometry, ThreeDCircularProjectionGeometry);
174 
176  itkGetMacro(NumberOfIterations, unsigned int);
177  itkSetMacro(NumberOfIterations, unsigned int);
179 
181  itkGetMacro(NumberOfProjectionsPerSubset, unsigned int);
182  itkSetMacro(NumberOfProjectionsPerSubset, unsigned int);
184 
186  itkGetMacro(Lambda, double);
187  itkSetMacro(Lambda, double);
189 
191  itkGetMacro(EnforcePositivity, bool);
192  itkSetMacro(EnforcePositivity, bool);
194 
196  void SetForwardProjectionFilter (ForwardProjectionType _arg) ITK_OVERRIDE;
197 
199  void SetBackProjectionFilter (BackProjectionType _arg) ITK_OVERRIDE;
200 
202  void SetWeights(const itk::Array2D<float> _arg);
203 
205  virtual void SetSignal(const std::vector<double> signal);
206 
208  itkSetMacro(DisableDisplacedDetectorFilter, bool)
209  itkGetMacro(DisableDisplacedDetectorFilter, bool)
210 
211 protected:
213  virtual ~FourDSARTConeBeamReconstructionFilter() ITK_OVERRIDE {}
214 
215  void GenerateInputRequestedRegion() ITK_OVERRIDE;
216 
217  void GenerateOutputInformation() ITK_OVERRIDE;
218 
219  void GenerateData() ITK_OVERRIDE;
220 
223  void VerifyInputInformation() ITK_OVERRIDE {}
224 
226  typename ExtractFilterType::Pointer m_ExtractFilter;
227  typename ExtractFilterType::Pointer m_ExtractFilterRayBox;
243 
245  std::vector< unsigned int > m_ProjectionsOrder;
248  std::vector<double> m_Signal;
250 
251 private:
255 
256  //purposely not implemented
258  void operator=(const Self&);
259 
262 
264  unsigned int m_NumberOfIterations;
265 
268  double m_Lambda;
269 }; // end of class
270 
271 } // end namespace rtk
272 
273 #ifndef ITK_MANUAL_INSTANTIATION
274 #include "rtkFourDSARTConeBeamReconstructionFilter.hxx"
275 #endif
276 
277 #endif
itk::ExtractImageFilter< ProjectionStackType, ProjectionStackType > ExtractFilterType
rtk::BackProjectionImageFilter< VolumeType, VolumeType > BackProjectionFilterType
Base class for forward projection, i.e. accumulation along x-ray lines.
FourDToProjectionStackFilterType::Pointer m_FourDToProjectionStackFilter
rtk::FourDToProjectionStackImageFilter< ProjectionStackType, VolumeSeriesType > FourDToProjectionStackFilterType
ProjectionStackToFourDFilterType::Pointer m_ProjectionStackToFourDFilter
Generate an n-dimensional image with constant pixel values.
IterativeConeBeamReconstructionFilter< VolumeSeriesType, ProjectionStackType > Superclass
Weigting for displaced detectors.
Analytical projection of a BoxShape.
itk::DivideOrZeroOutImageFilter< ProjectionStackType, ProjectionStackType, ProjectionStackType > DivideFilterType
rtk::ProjectionStackToFourDImageFilter< VolumeSeriesType, ProjectionStackType > ProjectionStackToFourDFilterType
itk::MultiplyImageFilter< ProjectionStackType, ProjectionStackType, ProjectionStackType > MultiplyFilterType
rtk::DisplacedDetectorImageFilter< ProjectionStackType > DisplacedDetectorFilterType
itk::ThresholdImageFilter< VolumeSeriesType > ThresholdFilterType
Projection geometry for a source and a 2-D flat panel.
rtk::ForwardProjectionImageFilter< ProjectionStackType, ProjectionStackType > ForwardProjectionFilterType
itk::SubtractImageFilter< ProjectionStackType, ProjectionStackType > SubtractFilterType
Implements the 4D Simultaneous Algebraic Reconstruction Technique.
ConstantVolumeSeriesSourceType::Pointer m_ConstantVolumeSeriesSource
rtk::RayBoxIntersectionImageFilter< ProjectionStackType, ProjectionStackType > RayBoxIntersectionFilterType
ConstantProjectionStackSourceType::Pointer m_ConstantProjectionStackSource
Mother class for cone beam reconstruction filters which need runtime selection of their forward and b...
itk::AddImageFilter< VolumeSeriesType, VolumeSeriesType > AddFilterType
rtk::ConstantImageSource< VolumeSeriesType > ConstantVolumeSeriesSourceType
rtk::ConstantImageSource< ProjectionStackType > ConstantProjectionStackSourceType
Implements part of the 4D reconstruction by conjugate gradient.
#define itkSetMacro(name, type)