RTK  2.0.1
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:
124  ITK_DISALLOW_COPY_AND_ASSIGN(FourDSARTConeBeamReconstructionFilter);
125 
131 
133  using InputImageType = VolumeSeriesType;
134  using OutputImageType = VolumeSeriesType;
135  using VolumeType = ProjectionStackType;
136 
137  using ForwardProjectionType = typename Superclass::ForwardProjectionType;
138  using BackProjectionType = typename Superclass::BackProjectionType;
139 
155 
157  itkNewMacro(Self);
158 
161 
163  void SetInputVolumeSeries(const VolumeSeriesType* VolumeSeries);
164  typename VolumeSeriesType::ConstPointer GetInputVolumeSeries();
166 
168  void SetInputProjectionStack(const ProjectionStackType* Projection);
169  typename ProjectionStackType::Pointer GetInputProjectionStack();
171 
173  itkGetModifiableObjectMacro(Geometry, ThreeDCircularProjectionGeometry);
174  itkSetObjectMacro(Geometry, ThreeDCircularProjectionGeometry);
176 
178  itkGetMacro(NumberOfIterations, unsigned int);
179  itkSetMacro(NumberOfIterations, unsigned int);
181 
183  itkGetMacro(NumberOfProjectionsPerSubset, unsigned int);
184  itkSetMacro(NumberOfProjectionsPerSubset, unsigned int);
186 
188  itkGetMacro(Lambda, double);
189  itkSetMacro(Lambda, double);
191 
193  itkGetMacro(EnforcePositivity, bool);
194  itkSetMacro(EnforcePositivity, bool);
196 
198  void SetForwardProjectionFilter (ForwardProjectionType _arg) override;
199 
201  void SetBackProjectionFilter (BackProjectionType _arg) override;
202 
204  void SetWeights(const itk::Array2D<float> _arg);
205 
207  virtual void SetSignal(const std::vector<double> signal);
208 
210  itkSetMacro(DisableDisplacedDetectorFilter, bool)
211  itkGetMacro(DisableDisplacedDetectorFilter, bool)
212 
213 protected:
215  ~FourDSARTConeBeamReconstructionFilter() override = default;
216 
217  void GenerateInputRequestedRegion() override;
218 
219  void GenerateOutputInformation() override;
220 
221  void GenerateData() override;
222 
225 #if ITK_VERSION_MAJOR<5
226  void VerifyInputInformation() override {}
227 #else
228  void VerifyInputInformation() const override {}
229 #endif
230 
231 
250 
252  std::vector< unsigned int > m_ProjectionsOrder;
255  std::vector<double> m_Signal;
257 
258 private:
262 
265 
267  unsigned int m_NumberOfIterations;
268 
271  double m_Lambda;
272 }; // end of class
273 
274 } // end namespace rtk
275 
276 #ifndef ITK_MANUAL_INSTANTIATION
277 #include "rtkFourDSARTConeBeamReconstructionFilter.hxx"
278 #endif
279 
280 #endif
Base class for forward projection, i.e. accumulation along x-ray lines.
FourDToProjectionStackFilterType::Pointer m_FourDToProjectionStackFilter
ProjectionStackToFourDFilterType::Pointer m_ProjectionStackToFourDFilter
Generate an n-dimensional image with constant pixel values.
Weigting for displaced detectors.
Analytical projection of a BoxShape.
Projection geometry for a source and a 2-D flat panel.
TOutputImage OutputImageType
Implements the 4D Simultaneous Algebraic Reconstruction Technique.
ConstantVolumeSeriesSourceType::Pointer m_ConstantVolumeSeriesSource
ConstantProjectionStackSourceType::Pointer m_ConstantProjectionStackSource
Mother class for cone beam reconstruction filters which need runtime selection of their forward and b...
typename Superclass::ForwardProjectionType ForwardProjectionType
Implements part of the 4D reconstruction by conjugate gradient.
#define itkSetMacro(name, type)