Line data    Source code

       1             : #pragma once
       2             : 
       3             : #include "LinearOperator.h"
       4             : #include "Geometry.h"
       5             : #include "BoundingBox.h"
       6             : #include "VolumeDescriptor.h"
       7             : #include "DetectorDescriptor.h"
       8             : 
       9             : #include "XrayProjector.h"
      10             : 
      11             : #include <vector>
      12             : #include <utility>
      13             : 
      14             : #include <Eigen/Geometry>
      15             : 
      16             : namespace elsa
      17             : {
      18             :     template <typename data_t = real_t>
      19             :     class SiddonsMethod;
      20             : 
      21             :     template <typename data_t>
      22             :     struct XrayProjectorInnerTypes<SiddonsMethod<data_t>> {
      23             :         using value_type = data_t;
      24             :         using forward_tag = ray_driven_tag;
      25             :         using backward_tag = ray_driven_tag;
      26             :     };
      27             :     /**
      28             :      * @brief Operator representing the discretized X-ray transform in 2d/3d using Siddon's method.
      29             :      *
      30             :      * The volume is traversed along the rays as specified by the Geometry. Each ray is traversed in
      31             :      * a contiguous fashion (i.e. along long voxel borders, not diagonally) and each traversed
      32             :      * voxel is counted as a hit with weight according to the length of the path of the ray through
      33             :      * the voxel.
      34             :      *
      35             :      * The geometry is represented as a list of projection matrices (see class Geometry), one for
      36             :      * each acquisition pose.
      37             :      *
      38             :      * Forward projection is accomplished using apply(), backward projection using applyAdjoint().
      39             :      * This projector is matched.
      40             :      *
      41             :      * @author David Frank - initial code, refactor to XrayProjector
      42             :      * @author Nikola Dinev - modularization, fixes
      43             :      *
      44             :      * @tparam data_t data type for the domain and range of the operator, defaulting to real_t
      45             :      *
      46             :      */
      47             :     template <typename data_t>
      48             :     class SiddonsMethod : public XrayProjector<SiddonsMethod<data_t>>
      49             :     {
      50             :     public:
      51             :         using self_type = SiddonsMethod<data_t>;
      52             :         using base_type = XrayProjector<self_type>;
      53             :         using value_type = typename base_type::value_type;
      54             :         using forward_tag = typename base_type::forward_tag;
      55             :         using backward_tag = typename base_type::backward_tag;
      56             : 
      57             :         /**
      58             :          * @brief Constructor for Siddon's method traversal.
      59             :          *
      60             :          * @param[in] domainDescriptor describing the domain of the operator (the volume)
      61             :          * @param[in] rangeDescriptor describing the range of the operator (the sinogram)
      62             :          *
      63             :          * The domain is expected to be 2 or 3 dimensional (volSizeX, volSizeY, [volSizeZ]),
      64             :          * the range is expected to be matching the domain (detSizeX, [detSizeY], acqPoses).
      65             :          */
      66             :         SiddonsMethod(const VolumeDescriptor& domainDescriptor,
      67             :                       const DetectorDescriptor& rangeDescriptor);
      68             : 
      69             :         /// default destructor
      70          69 :         ~SiddonsMethod() override = default;
      71             : 
      72             :     protected:
      73             :         /// default copy constructor, hidden from non-derived classes to prevent potential slicing
      74             :         SiddonsMethod(const SiddonsMethod<data_t>&) = default;
      75             : 
      76             :     private:
      77             :         /// implement the polymorphic clone operation
      78             :         SiddonsMethod<data_t>* _cloneImpl() const;
      79             : 
      80             :         /// implement the polymorphic comparison operation
      81             :         bool _isEqual(const LinearOperator<data_t>& other) const;
      82             : 
      83             :         data_t traverseRayForward(BoundingBox aabb, const RealRay_t& ray,
      84             :                                   const DataContainer<data_t>& x) const;
      85             : 
      86             :         void traverseRayBackward(BoundingBox aabb, const RealRay_t& ray,
      87             :                                  const value_type& detectorValue, DataContainer<data_t>& Aty) const;
      88             : 
      89             :         template <int dim>
      90             :         data_t doTraverseRayForward(BoundingBox aabb, const RealRay_t& ray,
      91             :                                     const DataContainer<data_t>& x,
      92             :                                     const DataDescriptor& domain) const;
      93             :         template <int dim>
      94             :         void doTraverseRayBackward(BoundingBox aabb, const RealRay_t& ray,
      95             :                                    const value_type& detectorValue,
      96             :                                    DataContainer<data_t>& Aty) const;
      97             :         friend class XrayProjector<self_type>;
      98             :     };
      99             : 
     100             : } // namespace elsa