diff --git a/sysrap/sn.h b/sysrap/sn.h index f32c0a5be..c9045274c 100644 --- a/sysrap/sn.h +++ b/sysrap/sn.h @@ -422,6 +422,7 @@ struct SYSRAP_API sn **/ static sn* Cylinder(double radius, double z1, double z2) ; + static sn* Trapezoid(double z, double y, double x, double ltx); static sn* Cone(double r1, double z1, double r2, double z2); static sn* Sphere(double radius); static sn* ZSphere(double radius, double z1, double z2); @@ -2802,6 +2803,24 @@ inline sn* sn::Cylinder(double radius, double z1, double z2) // static nd->setBB( -radius, -radius, z1, +radius, +radius, z2 ); return nd ; } + +/** + * Support right angular wedge from STEP. + * + * References: + * + * https://geant4-userdoc.web.cern.ch/UsersGuides/ForApplicationDeveloper/html/Detector/Geometry/geomSolids.html#constructed-solid-geometry-csg-solids + * https://github.com/Geant4/geant4/blob/master/source/geometry/solids/CSG/include/G4Trap.hh + */ +inline sn* sn::Trapezoid(double z, double y, double x, double ltx) +{ + assert( x > 0 && y > 0 && z > 0 && ltx > 0 && ltx < x ); + sn* nd = Create(CSG_TRAPEZOID); + nd->setPA(x, y, z, 0.f, ltx, 0.f); + nd->setBB(-0.5*x, -0.5*y, -0.5*z, +0.5*x, +0.5*y, +0.5*z); + return nd; +} + inline sn* sn::Cone(double r1, double z1, double r2, double z2) // static { assert( z2 > z1 ); @@ -4889,6 +4908,12 @@ inline void sn::setAABB_LeafFrame() assert( px == 0. && py == 0. && a == 0. ); setBB( px-radius, py-radius, z1, px+radius, py+radius, z2 ); } + else if( typecode == CSG_TRAPEZOID ) + { + double x, y, z, unused1, ltx, unused2; + getParam_(x, y, z, unused1, ltx, unused2); + setBB(-0.5*x, -0.5*y, -0.5*z, +0.5*x, +0.5*y, +0.5*z); + } else if( typecode == CSG_DISC ) { double px, py, ir, r, z1, z2 ; diff --git a/u4/U4Solid.h b/u4/U4Solid.h index ac4242096..0082a19a0 100644 --- a/u4/U4Solid.h +++ b/u4/U4Solid.h @@ -42,6 +42,7 @@ npy/NNodeUncoincide npy/NNodeNudger #include "G4Ellipsoid.hh" #include "G4Box.hh" #include "G4Tubs.hh" +#include "G4Trap.hh" #include "G4Polycone.hh" #include "G4Cons.hh" #include "G4Hype.hh" @@ -67,6 +68,7 @@ enum { _G4Ellipsoid, _G4Box, _G4Tubs, + _G4Trap, _G4Polycone, _G4Cons, _G4Hype, @@ -85,6 +87,7 @@ struct U4Solid static constexpr const char* G4Ellipsoid_ = "Ell" ; static constexpr const char* G4Box_ = "Box" ; static constexpr const char* G4Tubs_ = "Tub" ; + static constexpr const char* G4Trap_ = "Trp" ; static constexpr const char* G4Polycone_ = "Pol" ; static constexpr const char* G4Cons_ = "Con" ; static constexpr const char* G4Hype_ = "Hyp" ; @@ -132,6 +135,7 @@ struct U4Solid void init_Ellipsoid(); void init_Box(); void init_Tubs(); + void init_Trap(); void init_Polycone(); void init_Cons(); void init_Hype(); @@ -279,6 +283,7 @@ inline int U4Solid::Type(const char* name) // static if( strcmp(name, "G4Ellipsoid") == 0 ) type = _G4Ellipsoid ; if( strcmp(name, "G4Box") == 0 ) type = _G4Box ; if( strcmp(name, "G4Tubs") == 0 ) type = _G4Tubs ; + if( strcmp(name, "G4Trap") == 0 ) type = _G4Trap ; if( strcmp(name, "G4Polycone") == 0 ) type = _G4Polycone ; if( strcmp(name, "G4Cons") == 0 ) type = _G4Cons ; if( strcmp(name, "G4Hype") == 0 ) type = _G4Hype ; @@ -301,6 +306,7 @@ inline const char* U4Solid::Tag(int type) // static case _G4Ellipsoid: tag = G4Ellipsoid_ ; break ; case _G4Box: tag = G4Box_ ; break ; case _G4Tubs: tag = G4Tubs_ ; break ; + case _G4Trap: tag = G4Trap_ ; break ; case _G4Polycone: tag = G4Polycone_ ; break ; case _G4Cons: tag = G4Cons_ ; break ; case _G4Hype: tag = G4Hype_ ; break ; @@ -397,6 +403,7 @@ inline void U4Solid::init_Constituents() case _G4Ellipsoid : init_Ellipsoid() ; break ; case _G4Box : init_Box() ; break ; case _G4Tubs : init_Tubs() ; break ; + case _G4Trap : init_Trap() ; break ; case _G4Polycone : init_Polycone() ; break ; case _G4Cons : init_Cons() ; break ; case _G4Hype : init_Hype() ; break ; @@ -790,6 +797,20 @@ inline void U4Solid::init_Tubs() } +inline void U4Solid::init_Trap() +{ + const G4Trap* trap = dynamic_cast(solid); + assert(trap); + + double z = 2*trap->GetZHalfLength()/CLHEP::mm; + double y = 2*trap->GetYHalfLength1()/CLHEP::mm; + double x = 2*trap->GetXHalfLength1()/CLHEP::mm; + double ltx = 2*trap->GetXHalfLength2()/CLHEP::mm; + + root = sn::Trapezoid(z, y, x, ltx); +} + + inline void U4Solid::init_Polycone() { const G4Polycone* polycone = dynamic_cast(solid);