feat: introduce G4Trap converter

Plugins/Geant4/include/Acts/Plugins/Geant4/Geant4Converters.hpp

@@ -125,14 +125,22 @@ struct Geant4ShapeConverter {
   std::tuple<std::shared_ptr<RectangleBounds>, std::array<int, 2u>, double>
   rectangleBounds(const G4Box& g4Box);
 
-  /// @brief Convert to trapezoid bounds - from Trap
+  /// @brief Convert to trapezoid bounds - from Trd
   ///
   /// @param g4Trd a Geant4 trapezoid shape
   ///
   /// @return an ACTS Trapezoid bounds object, axis orientation, and thickness
   std::tuple<std::shared_ptr<TrapezoidBounds>, std::array<int, 2u>, double>
   trapezoidBounds(const G4Trd& g4Trd);
 
+  /// @brief Convert to trapezoid bounds - from Trap
+  ///
+  /// @param g4Trap a Geant4 trapezoid shape
+  ///
+  /// @return an ACTS Trapezoid bounds object, axis orientation, and thickness
+  std::tuple<std::shared_ptr<TrapezoidBounds>, std::array<int, 2u>, double>
+  trapezoidBounds(const G4Trap& g4Trap);
+
   /// @brief Convert to general solid into a planar shape
   ///
   /// @param g4Solid a Geant4 solid shape

Plugins/Geant4/src/Geant4Converters.cpp

@@ -38,6 +38,7 @@
 #include "G4RotationMatrix.hh"
 #include "G4ThreeVector.hh"
 #include "G4Transform3D.hh"
+#include "G4Trap.hh"
 #include "G4Trd.hh"
 #include "G4Tubs.hh"
 #include "G4VPhysicalVolume.hh"
@@ -229,6 +230,72 @@ Acts::Geant4ShapeConverter::trapezoidBounds(const G4Trd& g4Trd) {
   return {std::move(tBounds), rAxes, thickness};
 }
 
+std::tuple<std::shared_ptr<Acts::TrapezoidBounds>, std::array<int, 2u>, double>
+Acts::Geant4ShapeConverter::trapezoidBounds(const G4Trap& g4Trap) {
+  // primary parameters
+  auto y1 = static_cast<double>(g4Trap.GetYHalfLength1());
+  auto y2 = static_cast<double>(g4Trap.GetYHalfLength2());
+  auto x1 = static_cast<double>(g4Trap.GetXHalfLength1());
+  auto x2 = static_cast<double>(g4Trap.GetXHalfLength2());
+  auto x3 = static_cast<double>(g4Trap.GetXHalfLength3());
+  auto x4 = static_cast<double>(g4Trap.GetXHalfLength4());
+  auto phi = static_cast<double>(g4Trap.GetPhi());
+  auto theta = static_cast<double>(g4Trap.GetTheta());
+  auto z = static_cast<double>(g4Trap.GetZHalfLength());
+
+  double hlX0 = (x1 + x2) * 0.5;
+  double hlX1 = 2 * z * tan(theta) * cos(phi) + (x3 + x4) * 0.5;
+  double hlY0 = y1;
+  double hlY1 = y2 + 2 * z * tan(theta) * sin(phi);
+  double hlZ = z;
+
+  std::vector<double> dXYZ = {(hlX0 + hlX1) * 0.5, (hlY0 + hlY1) * 0.5, hlZ};
+
+  auto minAt = std::ranges::min_element(dXYZ);
+  std::size_t minPos = std::distance(dXYZ.begin(), minAt);
+  double thickness = 2. * dXYZ[minPos];
+
+  double halfLengthXminY = 0.;
+  double halfLengthXmaxY = 0.;
+  double halfLengthY = 0.;
+
+  std::array<int, 2u> rAxes = {};
+  switch (minPos) {
+    case 0: {
+      halfLengthXminY = std::min(hlY0, hlY1);
+      halfLengthXmaxY = std::max(hlY0, hlY1);
+      halfLengthY = hlZ;
+      rAxes = {1, 2};
+    } break;
+    case 1: {
+      halfLengthXminY = std::min(hlX0, hlX1);
+      halfLengthXmaxY = std::max(hlX0, hlX1);
+      halfLengthY = hlZ;
+      rAxes = {0, -2};
+    } break;
+    case 2: {
+      if (std::abs(hlY0 - hlY1) < std::abs(hlX0 - hlX1)) {
+        halfLengthXminY = std::min(hlX0, hlX1);
+        halfLengthXmaxY = std::max(hlX0, hlX1);
+        halfLengthY = (hlY0 + hlY1) * 0.5;
+        rAxes = {0, 1};
+      } else {
+        halfLengthXminY = std::min(hlY0, hlY1);
+        halfLengthXmaxY = std::max(hlY0, hlY1);
+        halfLengthY = (hlX0 + hlX1) * 0.5;
+        rAxes = {-1, 0};
+      }
+    } break;
+    default: {
+      throw std::runtime_error("Geant4Converters: could not convert G4Trap.");
+    }
+  }
+
+  auto tBounds = std::make_shared<TrapezoidBounds>(
+      halfLengthXminY, halfLengthXmaxY, halfLengthY);
+  return std::make_tuple(std::move(tBounds), rAxes, thickness);
+}
+
 std::tuple<std::shared_ptr<Acts::PlanarBounds>, std::array<int, 2u>, double>
 Acts::Geant4ShapeConverter::planarBounds(const G4VSolid& g4Solid) {
   const G4Box* box = dynamic_cast<const G4Box*>(&g4Solid);
@@ -309,7 +376,7 @@ std::shared_ptr<Acts::Surface> Acts::Geant4PhysicalVolumeConverter::surface(
     }
   }
 
-  // Into a Trapezoid
+  // Into a Trapezoid - from Trd
   auto g4Trd = dynamic_cast<const G4Trd*>(g4Solid);
   if (g4Trd != nullptr) {
     if (forcedType == Surface::SurfaceType::Other ||
@@ -326,6 +393,23 @@ std::shared_ptr<Acts::Surface> Acts::Geant4PhysicalVolumeConverter::surface(
     }
   }
 
+  // Into a Trapezoid - from Trap
+  auto g4Trap = dynamic_cast<const G4Trap*>(g4Solid);
+  if (g4Trap != nullptr) {
+    if (forcedType == Surface::SurfaceType::Other ||
+        forcedType == Surface::SurfaceType::Plane) {
+      auto [bounds, axes, original] =
+          Geant4ShapeConverter{}.trapezoidBounds(*g4Trap);
+      auto orientedToGlobal = axesOriented(toGlobal, axes);
+      surface = Acts::Surface::makeShared<PlaneSurface>(orientedToGlobal,
+                                                        std::move(bounds));
+      assignMaterial(*surface.get(), original, compressed);
+      return surface;
+    } else {
+      throw std::runtime_error("Can not convert 'G4Trap' into forced shape.");
+    }
+  }
+
   // Into a Cylinder, disc or line
   auto g4Tubs = dynamic_cast<const G4Tubs*>(g4Solid);
   if (g4Tubs != nullptr) {

Tests/UnitTests/Plugins/Geant4/Geant4ConvertersTests.cpp

@@ -35,6 +35,7 @@
 #include "G4RotationMatrix.hh"
 #include "G4SystemOfUnits.hh"
 #include "G4ThreeVector.hh"
+#include "G4Trap.hh"
 #include "G4Trd.hh"
 #include "G4Tubs.hh"
 #include "G4VPhysicalVolume.hh"
@@ -142,7 +143,7 @@ BOOST_AUTO_TEST_CASE(Geant4BoxConversion) {
   CHECK_CLOSE_ABS(thicknessY2, 4., 10e-10);
 }
 
-BOOST_AUTO_TEST_CASE(Geant4TrapzoidConversion) {
+BOOST_AUTO_TEST_CASE(Geant4TrapzoidConversionTrd) {
   // Standard TRD: XY are already well defined
   G4Trd trdXY("trdXY", 100, 150, 200, 200, 2);
   auto [boundsXY, axesXY, thicknessZ] =
@@ -212,6 +213,59 @@ BOOST_AUTO_TEST_CASE(Geant4TrapzoidConversion) {
   CHECK_CLOSE_ABS(thicknessY, 2., 10e-10);
 }
 
+BOOST_AUTO_TEST_CASE(Geant4TrapzoidConversionTrap) {
+  // x value changes, y is the symmetric side, z is the thickness
+  G4Trap trapXY("trapXY", 10., 15., 20., 20., 0.125);
+  auto [boundsXY, axesXY, thicknessZ] =
+      Acts::Geant4ShapeConverter{}.trapezoidBounds(trapXY);
+  CHECK_CLOSE_ABS(
+      boundsXY->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthXnegY), 10,
+      10e-10);
+  CHECK_CLOSE_ABS(
+      boundsXY->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthXposY), 15,
+      10e-10);
+  CHECK_CLOSE_ABS(
+      boundsXY->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthY), 20,
+      10e-10);
+  auto refXY = std::array<int, 2u>{0, 1};
+  BOOST_CHECK(axesXY == refXY);
+  CHECK_CLOSE_ABS(thicknessZ, 0.25, 10e-10);
+
+  // y value changes, x is the symmetric side, z is the thickness
+  G4Trap trapYx("trapYx", 22., 22., 11., 16., 0.250);
+  auto [boundsYx, axesYx, thicknessYx] =
+      Acts::Geant4ShapeConverter{}.trapezoidBounds(trapYx);
+  CHECK_CLOSE_ABS(
+      boundsYx->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthXnegY), 11,
+      10e-10);
+  CHECK_CLOSE_ABS(
+      boundsYx->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthXposY), 16,
+      10e-10);
+  CHECK_CLOSE_ABS(
+      boundsYx->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthY), 22,
+      10e-10);
+  auto refYx = std::array<int, 2u>{-1, 0};
+  BOOST_CHECK(axesYx == refYx);
+  CHECK_CLOSE_ABS(thicknessYx, 0.5, 10e-10);
+
+  // x is the thickness, y changes, z is symmetric
+  G4Trap trapXz("trapXz", 0.5, 0.5, 8., 16., 10.);
+  auto [boundsYZ, axesYZ, thicknessXZ] =
+      Acts::Geant4ShapeConverter{}.trapezoidBounds(trapXz);
+  CHECK_CLOSE_ABS(
+      boundsYZ->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthXnegY), 8,
+      10e-10);
+  CHECK_CLOSE_ABS(
+      boundsYZ->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthXposY), 16,
+      10e-10);
+  CHECK_CLOSE_ABS(
+      boundsYZ->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthY), 10.,
+      10e-10);
+  auto refYZ = std::array<int, 2u>{1, 2};
+  BOOST_CHECK(axesYZ == refYZ);
+  CHECK_CLOSE_ABS(thicknessXZ, 1., 10e-10);
+}
+
 BOOST_AUTO_TEST_CASE(Geant4PlanarConversion) {
   G4Box boxXY("boxXY", 23., 34., 1.);
   auto pBoundsBox =