Line data    Source code

       1             : #include "LimitedAngleTrajectoryGenerator.h"
       2             : #include "VolumeDescriptor.h"
       3             : #include "PlanarDetectorDescriptor.h"
       4             : #include "Logger.h"
       5             : 
       6             : #include <stdexcept>
       7             : 
       8             : namespace elsa
       9             : {
      10             :     std::unique_ptr<DetectorDescriptor> LimitedAngleTrajectoryGenerator::createTrajectory(
      11             :         index_t numberOfPoses,
      12             :         std::pair<elsa::geometry::Degree, elsa::geometry::Degree> missingWedgeAngles,
      13             :         const DataDescriptor& volumeDescriptor, index_t arcDegrees, real_t sourceToCenter,
      14             :         real_t centerToDetector, bool mirrored)
      15           4 :     {
      16             :         // pull in geometry namespace, to reduce cluttering
      17           4 :         using namespace geometry;
      18             : 
      19             :         // sanity check
      20           4 :         const auto dim = volumeDescriptor.getNumberOfDimensions();
      21           4 :         if (dim != 2) {
      22           0 :             throw InvalidArgumentError("LimitedAngleTrajectoryGenerator: can only handle 2D");
      23           0 :         }
      24             : 
      25           4 :         Logger::get("LimitedAngleTrajectoryGenerator")
      26           4 :             ->info("creating 2D trajectory with {} poses in an {} degree arc", numberOfPoses,
      27           4 :                    arcDegrees);
      28             : 
      29           4 :         const auto [coeffs, spacing] =
      30           4 :             calculateSizeAndSpacingPerGeometry(volumeDescriptor, numberOfPoses);
      31             : 
      32             :         // create vector and reserve the necessary size, minor optimization such that no new
      33             :         // allocations are necessary in the loop
      34           4 :         std::vector<Geometry> geometryList;
      35           4 :         geometryList.reserve(static_cast<std::size_t>(numberOfPoses));
      36             : 
      37           4 :         real_t wedgeArc = mirrored ? 2 * (missingWedgeAngles.second - missingWedgeAngles.first)
      38           4 :                                    : missingWedgeAngles.second - missingWedgeAngles.first;
      39             : 
      40           4 :         const real_t angleIncrement =
      41           4 :             (static_cast<real_t>(arcDegrees) - wedgeArc) / (static_cast<real_t>(numberOfPoses));
      42             : 
      43        1679 :         for (index_t i = 0;; ++i) {
      44        1679 :             Radian angle = Degree{static_cast<real_t>(i) * angleIncrement};
      45             : 
      46        1679 :             if (angle.to_degree() >= static_cast<real_t>(arcDegrees)) {
      47           4 :                 break;
      48           4 :             }
      49             : 
      50        1675 :             if (notInMissingWedge(angle, missingWedgeAngles, mirrored)) {
      51             :                 // use emplace_back, then no copy is created
      52        1281 :                 geometryList.emplace_back(SourceToCenterOfRotation{sourceToCenter},
      53        1281 :                                           CenterOfRotationToDetector{centerToDetector},
      54        1281 :                                           Radian{angle},
      55        1281 :                                           VolumeData2D{volumeDescriptor.getSpacingPerDimension(),
      56        1281 :                                                        volumeDescriptor.getLocationOfOrigin()},
      57        1281 :                                           SinogramData2D{Size2D{coeffs}, Spacing2D{spacing}});
      58        1281 :             }
      59        1675 :         }
      60             : 
      61           4 :         return std::make_unique<PlanarDetectorDescriptor>(coeffs, spacing, geometryList);
      62           4 :     }
      63             : 
      64             :     bool LimitedAngleTrajectoryGenerator::notInMissingWedge(
      65             :         elsa::geometry::Radian angle,
      66             :         std::pair<elsa::geometry::Degree, elsa::geometry::Degree> missingWedgeAngles, bool mirrored)
      67        1675 :     {
      68        1675 :         if (!mirrored) {
      69         577 :             return !(angle.to_degree() >= missingWedgeAngles.first
      70         577 :                      && angle.to_degree() <= missingWedgeAngles.second);
      71        1098 :         } else {
      72        1098 :             return !((angle.to_degree() >= missingWedgeAngles.first
      73        1098 :                       && angle.to_degree() <= missingWedgeAngles.second)
      74        1098 :                      || (angle.to_degree() >= (missingWedgeAngles.first + 180)
      75         805 :                          && angle.to_degree() <= (missingWedgeAngles.second + 180)));
      76        1098 :         }
      77        1675 :     }
      78             : } // namespace elsa