diff --git a/DDCore/include/DD4hep/Shapes.h b/DDCore/include/DD4hep/Shapes.h index 9907806885db826505531252514d936ca6484367..e73776517e6e0e54705acb2c106f3c99aee777db 100644 --- a/DDCore/include/DD4hep/Shapes.h +++ b/DDCore/include/DD4hep/Shapes.h @@ -85,7 +85,7 @@ namespace dd4hep { */ template <typename T> class Solid_type: public Handle<T> { protected: - void _setDimensions(double* param); + void _setDimensions(double* param) const; /// Assign pointrs and register solid to geometry void _assign(T* n, const std::string& nam, const std::string& tit, bool cbbox); @@ -170,6 +170,8 @@ namespace dd4hep { ShapelessSolid& operator=(ShapelessSolid&& copy) = default; /// Copy Assignment operator ShapelessSolid& operator=(const ShapelessSolid& copy) = default; + /// Set the assembly dimensions (noop, prints warning) + ShapelessSolid& setDimensions(const std::vector<double>& params); }; /// Class describing a box shape @@ -221,6 +223,8 @@ namespace dd4hep { Box& operator=(const Box& copy) = default; /// Set the box dimensions Box& setDimensions(double x_val, double y_val, double z_val); + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + Box& setDimensions(const std::vector<double>& params); /// Access half "length" of the box double x() const; /// Access half "width" of the box @@ -268,7 +272,9 @@ namespace dd4hep { HalfSpace& operator=(HalfSpace&& copy) = default; /// Copy Assignment operator HalfSpace& operator=(const HalfSpace& copy) = default; - }; + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + HalfSpace& setDimensions(const std::vector<double>& params); + }; /// Class describing a Polycone shape /** @@ -326,6 +332,8 @@ namespace dd4hep { /// Add Z-planes to the Polycone void addZPlanes(const std::vector<double>& rmin, const std::vector<double>& rmax, const std::vector<double>& z); + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + Polycone& setDimensions(const std::vector<double>& params); }; /// Class describing a cone segment shape @@ -369,6 +377,8 @@ namespace dd4hep { ConeSegment& setDimensions(double dz, double rmin1, double rmax1, double rmin2, double rmax2, double startPhi = 0.0, double endPhi = 2.0 * M_PI); + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + ConeSegment& setDimensions(const std::vector<double>& params); }; #if 0 /// Intermediate class to overcome drawing probles with the TGeoTubeSeg @@ -439,6 +449,8 @@ namespace dd4hep { Tube& operator=(const Tube& copy) = default; /// Set the tube dimensions Tube& setDimensions(double rmin, double rmax, double dz, double startPhi=0.0, double endPhi=2*M_PI); + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + Tube& setDimensions(const std::vector<double>& params); }; /// Class describing a tube shape of a section of a cut tube segment @@ -482,6 +494,8 @@ namespace dd4hep { CutTube& operator=(CutTube&& copy) = default; /// Copy Assignment operator CutTube& operator=(const CutTube& copy) = default; + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + CutTube& setDimensions(const std::vector<double>& params); }; @@ -526,6 +540,8 @@ namespace dd4hep { TruncatedTube& operator=(TruncatedTube&& copy) = default; /// Copy Assignment operator TruncatedTube& operator=(const TruncatedTube& copy) = default; + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + TruncatedTube& setDimensions(const std::vector<double>& params); }; /// Class describing a elliptical tube shape @@ -578,6 +594,8 @@ namespace dd4hep { EllipticalTube& operator=(const EllipticalTube& copy) = default; /// Set the tube dimensions EllipticalTube& setDimensions(double a, double b, double dz); + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + EllipticalTube& setDimensions(const std::vector<double>& params); }; /// Class describing a cone shape @@ -627,6 +645,8 @@ namespace dd4hep { Cone& operator=(const Cone& copy) = default; /// Set the box dimensions Cone& setDimensions(double z, double rmin1, double rmax1, double rmin2, double rmax2); + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + Cone& setDimensions(const std::vector<double>& params); }; /// Class describing a trap shape @@ -686,6 +706,8 @@ namespace dd4hep { Trap& setDimensions(double z, double theta, double phi, double h1, double bl1, double tl1, double alpha1, double h2, double bl2, double tl2, double alpha2); + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + Trap& setDimensions(const std::vector<double>& params); }; /// Class describing a pseudo trap shape (CMS'ism) @@ -728,6 +750,8 @@ namespace dd4hep { PseudoTrap& operator=(PseudoTrap&& copy) = default; /// Copy Assignment operator PseudoTrap& operator=(const PseudoTrap& copy) = default; + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + PseudoTrap& setDimensions(const std::vector<double>& params); }; /// Class describing a Trd1 shape @@ -779,6 +803,8 @@ namespace dd4hep { Trd1& operator=(const Trd1& copy) = default; /// Set the Trd1 dimensions Trd1& setDimensions(double x1, double x2, double y, double z); + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + Trd1& setDimensions(const std::vector<double>& params); }; /// Class describing a Trd2 shape @@ -830,6 +856,8 @@ namespace dd4hep { Trd2& operator=(const Trd2& copy) = default; /// Set the Trd2 dimensions Trd2& setDimensions(double x1, double x2, double y1, double y2, double z); + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + Trd2& setDimensions(const std::vector<double>& params); }; /// Shortcut name definition typedef Trd2 Trapezoid; @@ -881,6 +909,8 @@ namespace dd4hep { Torus& operator=(const Torus& copy) = default; /// Set the Torus dimensions Torus& setDimensions(double r, double rmin, double rmax, double startPhi, double deltaPhi); + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + Torus& setDimensions(const std::vector<double>& params); }; /// Class describing a sphere shape @@ -958,6 +988,8 @@ namespace dd4hep { Sphere& setDimensions(double rmin, double rmax, double startTheta, double endTheta, double startPhi, double endPhi); + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + Sphere& setDimensions(const std::vector<double>& params); }; /// Class describing a Paraboloid shape @@ -1007,6 +1039,8 @@ namespace dd4hep { Paraboloid& operator=(const Paraboloid& copy) = default; /// Set the Paraboloid dimensions Paraboloid& setDimensions(double r_low, double r_high, double delta_z); + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + Paraboloid& setDimensions(const std::vector<double>& params); }; /// Class describing a Hyperboloid shape @@ -1056,6 +1090,8 @@ namespace dd4hep { Hyperboloid& operator=(const Hyperboloid& copy) = default; /// Set the Hyperboloid dimensions Hyperboloid& setDimensions(double rin, double stin, double rout, double stout, double dz); + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + Hyperboloid& setDimensions(const std::vector<double>& params); }; /// Class describing a regular polyhedron shape @@ -1106,6 +1142,8 @@ namespace dd4hep { PolyhedraRegular& operator=(PolyhedraRegular&& copy) = default; /// Copy Assignment operator PolyhedraRegular& operator=(const PolyhedraRegular& copy) = default; + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + PolyhedraRegular& setDimensions(const std::vector<double>& params); }; /// Class describing a regular polyhedron shape @@ -1161,6 +1199,8 @@ namespace dd4hep { Polyhedra& operator=(Polyhedra&& copy) = default; /// Copy Assignment operator Polyhedra& operator=(const Polyhedra& copy) = default; + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + Polyhedra& setDimensions(const std::vector<double>& params); }; /// Class describing a extruded polygon shape @@ -1216,6 +1256,8 @@ namespace dd4hep { ExtrudedPolygon& operator=(ExtrudedPolygon&& copy) = default; /// Copy Assignment operator ExtrudedPolygon& operator=(const ExtrudedPolygon& copy) = default; + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + ExtrudedPolygon& setDimensions(const std::vector<double>& params); }; /// Class describing an arbitray solid defined by 8 vertices. @@ -1255,6 +1297,8 @@ namespace dd4hep { EightPointSolid& operator=(EightPointSolid&& copy) = default; /// Copy Assignment operator EightPointSolid& operator=(const EightPointSolid& copy) = default; + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + EightPointSolid& setDimensions(const std::vector<double>& params); }; /// Base class describing boolean (=union,intersection,subtraction) solids @@ -1331,6 +1375,8 @@ namespace dd4hep { SubtractionSolid& operator=(SubtractionSolid&& copy) = default; /// Copy Assignment operator SubtractionSolid& operator=(const SubtractionSolid& copy) = default; + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + SubtractionSolid& setDimensions(const std::vector<double>& params); }; /// Class describing boolean union solid @@ -1379,6 +1425,8 @@ namespace dd4hep { UnionSolid& operator=(UnionSolid&& copy) = default; /// Copy Assignment operator UnionSolid& operator=(const UnionSolid& copy) = default; + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + UnionSolid& setDimensions(const std::vector<double>& params); }; /// Class describing boolean intersection solid @@ -1427,6 +1475,8 @@ namespace dd4hep { IntersectionSolid& operator=(IntersectionSolid&& copy) = default; /// Copy Assignment operator IntersectionSolid& operator=(const IntersectionSolid& copy) = default; + /// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. + IntersectionSolid& setDimensions(const std::vector<double>& params); }; } /* End namespace dd4hep */ diff --git a/DDCore/src/Shapes.cpp b/DDCore/src/Shapes.cpp index 6b8afd02a83a3540edd0a49ed4ff686864b9f107..ca6f731f638890742e4db06284882ea970ab325b 100644 --- a/DDCore/src/Shapes.cpp +++ b/DDCore/src/Shapes.cpp @@ -22,6 +22,7 @@ // C/C++ include files #include <stdexcept> #include <iomanip> +#include <sstream> // ROOT includes #include "TClass.h" @@ -34,53 +35,95 @@ using namespace dd4hep; namespace units = dd4hep; namespace { + template <typename SHAPE> void invalidSetDimensionCall(SHAPE sh, const vector<double>& params) { + stringstream str; + str << "Shape: " << typeName(typeid(sh)) << "::setDimension: " + << "Invalid number of parameters: " << params.size(); + throw runtime_error(str.str()); + } + std::vector<double> get_shape_dimensions(TGeoShape* shape) { if (shape) { - if (shape->IsA() == TGeoShapeAssembly::Class()) { + TClass* cl = shape->IsA(); + if (cl == TGeoShapeAssembly::Class()) { const TGeoShapeAssembly* sh = (const TGeoShapeAssembly*) shape; return { sh->GetDX(), sh->GetDY(), sh->GetDZ() }; } - else if (shape->IsA() == TGeoBBox::Class()) { + else if (cl == TGeoBBox::Class()) { const TGeoBBox* sh = (const TGeoBBox*) shape; return { sh->GetDX(), sh->GetDY(), sh->GetDZ() }; } - else if (shape->IsA() == TGeoTube::Class()) { + else if (cl == TGeoHalfSpace::Class()) { + TGeoHalfSpace* sh = (TGeoHalfSpace*)(const_cast<TGeoShape*>(shape)); + return { sh->GetPoint()[0], sh->GetPoint()[1], sh->GetPoint()[2], + sh->GetNorm()[0], sh->GetNorm()[1], sh->GetNorm()[2] }; + } + else if (cl == TGeoPcon::Class()) { + const TGeoPcon* sh = (const TGeoPcon*) shape; + vector<double> pars { sh->GetPhi1()*units::deg, sh->GetDphi()*units::deg, double(sh->GetNz()) }; + for (Int_t i = 0; i < sh->GetNz(); ++i) { + pars.emplace_back(sh->GetZ(i)); + pars.emplace_back(sh->GetRmin(i)); + pars.emplace_back(sh->GetRmax(i)); + } + return pars; + } + else if (cl == TGeoConeSeg::Class()) { + const TGeoConeSeg* sh = (const TGeoConeSeg*) shape; + return { sh->GetDz(), sh->GetRmin1(), sh->GetRmax1(), sh->GetRmin2(), sh->GetRmax2(), + sh->GetPhi1()*units::deg, sh->GetPhi2()*units::deg }; + } + else if (cl == TGeoCone::Class()) { + const TGeoCone* sh = (const TGeoCone*) shape; + return { sh->GetDz(), sh->GetRmin1(), sh->GetRmax1(), sh->GetRmin2(), sh->GetRmax2(), sh->GetDz() }; + } + else if (cl == TGeoTube::Class()) { const TGeoTube* sh = (const TGeoTube*) shape; return { sh->GetRmin(), sh->GetRmax(), sh->GetDz() }; } - else if (shape->IsA() == TGeoTubeSeg::Class()) { + else if (cl == TGeoTubeSeg::Class()) { const TGeoTubeSeg* sh = (const TGeoTubeSeg*) shape; return { sh->GetRmin(), sh->GetRmax(), sh->GetDz(), sh->GetPhi1()*units::deg, sh->GetPhi2()*units::deg }; } - else if (shape->IsA() == TGeoEltu::Class()) { + else if (cl == TGeoEltu::Class()) { const TGeoEltu* sh = (const TGeoEltu*) shape; return { sh->GetA(), sh->GetB(), sh->GetDz() }; } - else if (shape->IsA() == TGeoTrd1::Class()) { + else if (cl == TGeoTrd1::Class()) { const TGeoTrd1* sh = (const TGeoTrd1*) shape; return { sh->GetDx1(), sh->GetDx2(), sh->GetDy(), sh->GetDz() }; } - else if (shape->IsA() == TGeoTrd2::Class()) { + else if (cl == TGeoTrd2::Class()) { const TGeoTrd2* sh = (const TGeoTrd2*) shape; return { sh->GetDx1(), sh->GetDx2(), sh->GetDy1(), sh->GetDy2(), sh->GetDz() }; } - else if (shape->IsA() == TGeoHype::Class()) { + else if (cl == TGeoParaboloid::Class()) { + const TGeoParaboloid* sh = (const TGeoParaboloid*) shape; + return { sh->GetRlo(), sh->GetRhi(), sh->GetDz() }; + } + else if (cl == TGeoHype::Class()) { const TGeoHype* sh = (const TGeoHype*) shape; - return { sh->GetDz(), sh->GetRmin(), sh->GetStIn(), sh->GetRmax(), sh->GetStOut(), sh->GetDz() }; + return { sh->GetRmin(), sh->GetStIn()*units::deg, + sh->GetRmax(), sh->GetStOut()*units::deg, sh->GetDz() }; } - else if (shape->IsA() == TGeoXtru::Class()) { - const TGeoXtru* sh = (const TGeoXtru*) shape; - Int_t nz = sh->GetNz(); - vector<double> pars { double(nz) }; - for(Int_t i=0; i<nz; ++i) { - pars.emplace_back(sh->GetZ(i)); - pars.emplace_back(sh->GetXOffset(i)); - pars.emplace_back(sh->GetYOffset(i)); - pars.emplace_back(sh->GetScale(i)); - } - return pars; + else if (cl == TGeoSphere::Class()) { + const TGeoSphere* sh = (const TGeoSphere*) shape; + return { sh->GetRmin(), sh->GetRmax(), + sh->GetTheta1()*units::deg, sh->GetTheta2()*units::deg, + sh->GetPhi1()*units::deg, sh->GetPhi2()*units::deg }; + } + else if (cl == TGeoTorus::Class()) { + const TGeoTorus* sh = (const TGeoTorus*) shape; + return { sh->GetR(), sh->GetRmin(), sh->GetRmax(), + sh->GetPhi1()*units::deg, sh->GetDphi()*units::deg }; } - else if (shape->IsA() == TGeoPgon::Class()) { + else if (cl == TGeoTrap::Class()) { + const TGeoTrap* sh = (const TGeoTrap*) shape; + return { sh->GetDz(), sh->GetTheta()*units::deg, sh->GetPhi()*units::deg, + sh->GetH1(), sh->GetBl1(), sh->GetBl2(), sh->GetAlpha1()*units::deg, + sh->GetH2(), sh->GetTl1(), sh->GetTl2(), sh->GetAlpha2()*units::deg }; + } + else if (cl == TGeoPgon::Class()) { const TGeoPgon* sh = (const TGeoPgon*) shape; vector<double> pars { sh->GetPhi1()*units::deg, sh->GetDphi()*units::deg, double(sh->GetNedges()), double(sh->GetNz()) }; for (Int_t i = 0; i < sh->GetNz(); ++i) { @@ -90,42 +133,19 @@ namespace { } return pars; } - else if (shape->IsA() == TGeoPcon::Class()) { - const TGeoPcon* sh = (const TGeoPcon*) shape; - vector<double> pars { sh->GetPhi1()*units::deg, sh->GetDphi()*units::deg, double(sh->GetNz()) }; - for (Int_t i = 0; i < sh->GetNz(); ++i) { + else if (cl == TGeoXtru::Class()) { + const TGeoXtru* sh = (const TGeoXtru*) shape; + Int_t nz = sh->GetNz(); + vector<double> pars { double(nz) }; + for(Int_t i=0; i<nz; ++i) { pars.emplace_back(sh->GetZ(i)); - pars.emplace_back(sh->GetRmin(i)); - pars.emplace_back(sh->GetRmax(i)); + pars.emplace_back(sh->GetXOffset(i)); + pars.emplace_back(sh->GetYOffset(i)); + pars.emplace_back(sh->GetScale(i)); } return pars; } - else if (shape->IsA() == TGeoCone::Class()) { - const TGeoCone* sh = (const TGeoCone*) shape; - return { sh->GetDz(), sh->GetRmin1(), sh->GetRmax1(), sh->GetRmin2(), sh->GetRmax2(), sh->GetDz() }; - } - else if (shape->IsA() == TGeoConeSeg::Class()) { - const TGeoConeSeg* sh = (const TGeoConeSeg*) shape; - return { sh->GetDz(), sh->GetRmin1(), sh->GetRmax1(), sh->GetRmin2(), sh->GetRmax2(), sh->GetPhi1()*units::deg, sh->GetPhi2()*units::deg }; - } - else if (shape->IsA() == TGeoParaboloid::Class()) { - const TGeoParaboloid* sh = (const TGeoParaboloid*) shape; - return { sh->GetRlo(), sh->GetRhi(), sh->GetDz() }; - } - else if (shape->IsA() == TGeoSphere::Class()) { - const TGeoSphere* sh = (const TGeoSphere*) shape; - return { sh->GetRmin(), sh->GetRmax(), sh->GetTheta1()*units::deg, sh->GetTheta2()*units::deg, sh->GetPhi1()*units::deg, sh->GetPhi2()*units::deg }; - } - else if (shape->IsA() == TGeoTorus::Class()) { - const TGeoTorus* sh = (const TGeoTorus*) shape; - return { sh->GetR(), sh->GetRmin(), sh->GetRmax(), sh->GetPhi1()*units::deg, sh->GetDphi()*units::deg }; - } - else if (shape->IsA() == TGeoTrap::Class()) { - const TGeoTrap* sh = (const TGeoTrap*) shape; - return { sh->GetDz(), sh->GetTheta()*units::deg, sh->GetPhi()*units::deg, sh->GetH1(), sh->GetH2(), - sh->GetBl1(), sh->GetBl2(), sh->GetTl1(), sh->GetTl2(), sh->GetAlpha1()*units::deg, sh->GetAlpha2()*units::deg }; - } - else if (shape->IsA() == TGeoArb8::Class()) { + else if (cl == TGeoArb8::Class()) { TGeoTrap* sh = (TGeoTrap*) shape; const Double_t* vertices = sh->GetVertices(); vector<double> pars; @@ -135,7 +155,7 @@ namespace { } return pars; } - else if (shape->IsA() == TGeoCompositeShape::Class()) { + else if (cl == TGeoCompositeShape::Class()) { const TGeoCompositeShape* sh = (const TGeoCompositeShape*) shape; const TGeoBoolNode* boolean = sh->GetBoolNode(); TGeoBoolNode::EGeoBoolType oper = boolean->GetBooleanOperator(); @@ -163,7 +183,7 @@ namespace { } else { printout(ERROR,"Solid","Failed to access dimensions for shape of type:%s.", - shape->IsA()->GetName()); + cl->GetName()); } return {}; } @@ -438,7 +458,7 @@ namespace dd4hep { } } -template <typename T> void Solid_type<T>::_setDimensions(double* param) { +template <typename T> void Solid_type<T>::_setDimensions(double* param) const { this->ptr()->SetDimensions(param); this->ptr()->ComputeBBox(); } @@ -499,6 +519,12 @@ ShapelessSolid::ShapelessSolid(const string& nam) { _assign(new TGeoShapeAssembly(), nam, "Assembly", true); } +/// Set the assembly dimensions (noop, prints warning) +ShapelessSolid& ShapelessSolid::setDimensions(const vector<double>& /* params */) { + printout(WARNING,"ShapelessSolid","+++ setDimensions is a compatibility call. Nothing implemented."); + return *this; +} + void Box::make(const std::string& nam, double x_val, double y_val, double z_val) { _assign(new TGeoBBox(nam.c_str(), x_val, y_val, z_val), "", "Box", true); } @@ -510,6 +536,17 @@ Box& Box::setDimensions(double x_val, double y_val, double z_val) { return *this; } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +Box& Box::setDimensions(const vector<double>& params) { + // No angles in this shape. Just pass on the parameters + if ( params.size() != 3 ) { + invalidSetDimensionCall(*this,params); + } + auto pars = params; + _setDimensions(&pars[0]); + return *this; +} + /// Access half "length" of the box double Box::x() const { return this->ptr()->GetDX(); @@ -530,6 +567,17 @@ void HalfSpace::make(const std::string& nam, const double* const point, const do _assign(new TGeoHalfSpace(nam.c_str(),(Double_t*)point, (Double_t*)normal), "", "HalfSpace",true); } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +HalfSpace& HalfSpace::setDimensions(const vector<double>& params) { + // No angles in this shape. Just pass on the parameters + if ( params.size() != 6 ) { + invalidSetDimensionCall(*this,params); + } + auto pars = params; + _setDimensions(&pars[0]); + return *this; +} + /// Constructor to be used when creating a new object Polycone::Polycone(double startPhi, double deltaPhi) { _assign(new TGeoPcon(startPhi/units::deg, deltaPhi/units::deg, 0), "", "Polycone", false); @@ -627,7 +675,7 @@ void Polycone::addZPlanes(const vector<double>& rmin, const vector<double>& rmax TGeoPcon* sh = *this; vector<double> params; size_t num = sh->GetNz(); - if (rmin.size() < 2) { + if (rmin.size() < 2) { throw runtime_error("dd4hep: PolyCone::addZPlanes> Not enough Z planes. minimum is 2!"); } params.emplace_back(sh->GetPhi1()); @@ -646,6 +694,25 @@ void Polycone::addZPlanes(const vector<double>& rmin, const vector<double>& rmax _setDimensions(¶ms[0]); } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +Polycone& Polycone::setDimensions(const vector<double>& params) { + if ( params.size() < 3 ) { + invalidSetDimensionCall(*this,params); + } + size_t nz = size_t(params[2]); + if ( params.size() != 3 + 3*nz ) { + invalidSetDimensionCall(*this,params); + } + vector<double> pars; + pars.emplace_back(params[0]/units::deg); + pars.emplace_back(params[1]/units::deg); + pars.emplace_back(params[2]); + for (size_t i = 3, n=3+3*nz; i < n; ++i) + pars.emplace_back(params[i]); + _setDimensions(&pars[0]); + return *this; +} + /// Constructor to be used when creating a new cone segment object ConeSegment::ConeSegment(double dz, double rmin1, double rmax1, @@ -677,6 +744,43 @@ ConeSegment& ConeSegment::setDimensions(double dz, return *this; } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +ConeSegment& ConeSegment::setDimensions(const vector<double>& params) { + if ( params.size() != 7 ) { + invalidSetDimensionCall(*this,params); + } + vector<double> pars; + for (size_t i = 0; i < 5; ++i) + pars.emplace_back(params[i]); + pars.emplace_back(params[5]/units::deg); + pars.emplace_back(params[6]/units::deg); + _setDimensions(&pars[0]); + return *this; +} + +/// Constructor to be used when creating a new object with attribute initialization +void Cone::make(const std::string& nam, double z, double rmin1, double rmax1, double rmin2, double rmax2) { + _assign(new TGeoCone(nam.c_str(), z, rmin1, rmax1, rmin2, rmax2 ), "", "Cone", true); +} + +/// Set the box dimensions (startPhi=0.0, endPhi=2*pi) +Cone& Cone::setDimensions(double z, double rmin1, double rmax1, double rmin2, double rmax2) { + double params[] = { z, rmin1, rmax1, rmin2, rmax2 }; + _setDimensions(params); + return *this; +} + +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +Cone& Cone::setDimensions(const vector<double>& params) { + // No angles in this shape. Just pass on the parameters + if ( params.size() != 5 ) { + invalidSetDimensionCall(*this,params); + } + auto pars = params; + _setDimensions(&pars[0]); + return *this; +} + /// Constructor to be used when creating a new object with attribute initialization void Tube::make(const string& nam, double rmin, double rmax, double z, double startPhi, double endPhi) { _assign(new TGeoTubeSeg(nam.c_str(), rmin, rmax, z, startPhi/units::deg, endPhi/units::deg),nam,"Tube",true); @@ -689,6 +793,16 @@ Tube& Tube::setDimensions(double rmin, double rmax, double z, double startPhi, d return *this; } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +Tube& Tube::setDimensions(const vector<double>& params) { + if ( params.size() != 3 ) { + invalidSetDimensionCall(*this,params); + } + auto pars = params; + _setDimensions(&pars[0]); + return *this; +} + /// Constructor to be used when creating a new object with attribute initialization CutTube::CutTube(double rmin, double rmax, double dz, double startPhi, double endPhi, double lx, double ly, double lz, double tx, double ty, double tz) { @@ -708,6 +822,12 @@ void CutTube::make(const std::string& nam, double rmin, double rmax, double dz, _assign(new TGeoCtub(nam.c_str(), rmin,rmax,dz,startPhi,endPhi,lx,ly,lz,tx,ty,tz),"","CutTube",true); } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +CutTube& CutTube::setDimensions(const vector<double>& params) { + if ( params.size() ) {} + throw runtime_error("CutTube: Special Boolean shape. setDimensions is not implemented!"); +} + /// Constructor to create a truncated tube object with attribute initialization TruncatedTube::TruncatedTube(double zhalf, double rmin, double rmax, double startPhi, double deltaPhi, double cutAtStart, double cutAtDelta, bool cutInside) @@ -790,20 +910,25 @@ void TruncatedTube::make(const std::string& nam, #endif } -/// Constructor to be used when creating a new object with attribute initialization -void EllipticalTube::make(const std::string& nam, double a, double b, double dz) { - _assign(new TGeoEltu(nam.c_str(), a, b, dz), "", "EllipticalTube", true); +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +TruncatedTube& TruncatedTube::setDimensions(const vector<double>& params) { + if ( params.size() ) {} + throw runtime_error("TruncatedTube: Special Boolean shape. setDimensions is not implemented!"); } /// Constructor to be used when creating a new object with attribute initialization -void Cone::make(const std::string& nam, double z, double rmin1, double rmax1, double rmin2, double rmax2) { - _assign(new TGeoCone(nam.c_str(), z, rmin1, rmax1, rmin2, rmax2 ), "", "Cone", true); +void EllipticalTube::make(const std::string& nam, double a, double b, double dz) { + _assign(new TGeoEltu(nam.c_str(), a, b, dz), "", "EllipticalTube", true); } -/// Set the box dimensions (startPhi=0.0, endPhi=2*pi) -Cone& Cone::setDimensions(double z, double rmin1, double rmax1, double rmin2, double rmax2) { - double params[] = { z, rmin1, rmax1, rmin2, rmax2 }; - _setDimensions(params); +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +EllipticalTube& EllipticalTube::setDimensions(const vector<double>& params) { + // No angles in this shape. Just pass on the parameters + if ( params.size() != 3 ) { + invalidSetDimensionCall(*this,params); + } + auto pars = params; + _setDimensions(&pars[0]); return *this; } @@ -819,6 +944,17 @@ Trd1& Trd1::setDimensions(double x1, double x2, double y, double z) { return *this; } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +Trd1& Trd1::setDimensions(const vector<double>& params) { + // No angles in this shape. Just pass on the parameters + if ( params.size() != 4 ) { + invalidSetDimensionCall(*this,params); + } + auto pars = params; + _setDimensions(&pars[0]); + return *this; +} + /// Constructor to be used when creating a new object with attribute initialization void Trd2::make(const std::string& nam, double x1, double x2, double y1, double y2, double z) { _assign(new TGeoTrd2(nam.c_str(), x1, x2, y1, y2, z ), "", "Trd2", true); @@ -831,6 +967,17 @@ Trd2& Trd2::setDimensions(double x1, double x2, double y1, double y2, double z) return *this; } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +Trd2& Trd2::setDimensions(const vector<double>& params) { + // No angles in this shape. Just pass on the parameters + if ( params.size() != 5 ) { + invalidSetDimensionCall(*this,params); + } + auto pars = params; + _setDimensions(&pars[0]); + return *this; +} + /// Constructor to be used when creating a new object with attribute initialization void Paraboloid::make(const std::string& nam, double r_low, double r_high, double delta_z) { _assign(new TGeoParaboloid(nam.c_str(), r_low, r_high, delta_z ), "", "Paraboloid", true); @@ -843,6 +990,17 @@ Paraboloid& Paraboloid::setDimensions(double r_low, double r_high, double delta_ return *this; } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +Paraboloid& Paraboloid::setDimensions(const vector<double>& params) { + // No angles in this shape. Just pass on the parameters + if ( params.size() != 3 ) { + invalidSetDimensionCall(*this,params); + } + auto pars = params; + _setDimensions(&pars[0]); + return *this; +} + /// Constructor to create a new anonymous object with attribute initialization void Hyperboloid::make(const std::string& nam, double rin, double stin, double rout, double stout, double dz) { _assign(new TGeoHype(nam.c_str(), rin, stin/units::deg, rout, stout/units::deg, dz), "", "Hyperboloid", true); @@ -855,6 +1013,16 @@ Hyperboloid& Hyperboloid::setDimensions(double rin, double stin, double rout, do return *this; } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +Hyperboloid& Hyperboloid::setDimensions(const vector<double>& params) { + if ( params.size() != 5 ) { + invalidSetDimensionCall(*this,params); + } + double pars[] = { params[0], params[1]/units::deg, params[2], params[3]/units::deg, params[4]}; + _setDimensions(pars); + return *this; +} + /// Constructor function to be used when creating a new object with attribute initialization void Sphere::make(const std::string& nam, double rmin, double rmax, double startTheta, double endTheta, double startPhi, double endPhi) { _assign(new TGeoSphere(nam.c_str(), rmin, rmax, @@ -869,6 +1037,20 @@ Sphere& Sphere::setDimensions(double rmin, double rmax, double startTheta, doubl return *this; } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +Sphere& Sphere::setDimensions(const vector<double>& params) { + if ( params.size() < 2 ) { + invalidSetDimensionCall(*this,params); + } + double pars[] = { params[0], params[1], 0., M_PI, 0., 2*M_PI}; + if (params.size() > 2) pars[2] = params[2]/units::deg; + if (params.size() > 3) pars[3] = params[3]/units::deg; + if (params.size() > 4) pars[4] = params[4]/units::deg; + if (params.size() > 5) pars[5] = params[5]/units::deg; + _setDimensions(pars); + return *this; +} + /// Constructor to be used when creating a new object with attribute initialization void Torus::make(const std::string& nam, double r, double rmin, double rmax, double startPhi, double deltaPhi) { _assign(new TGeoTorus(nam.c_str(), r, rmin, rmax, startPhi/units::deg, deltaPhi/units::deg), "", "Torus", true); @@ -881,6 +1063,17 @@ Torus& Torus::setDimensions(double r, double rmin, double rmax, double startPhi, return *this; } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +Torus& Torus::setDimensions(const vector<double>& params) { + // No angles in this shape. Just pass on the parameters + if ( params.size() != 5 ) { + invalidSetDimensionCall(*this,params); + } + double pars[] = { params[0], params[1], params[2], params[3]/units::deg, params[4]/units::deg }; + _setDimensions(pars); + return *this; +} + /// Constructor to be used when creating a new anonymous object with attribute initialization Trap::Trap(double z, double theta, double phi, double h1, double bl1, double tl1, double alpha1, @@ -925,6 +1118,20 @@ Trap& Trap::setDimensions(double z, double theta, double phi, return *this; } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +Trap& Trap::setDimensions(const vector<double>& params) { + auto pars = params; + if ( params.size() != 10 ) { + invalidSetDimensionCall(*this,params); + } + pars[1] /= units::deg; + pars[2] /= units::deg; + pars[6] /= units::deg; + pars[10] /= units::deg; + _setDimensions(&pars[0]); + return *this; +} + /// Internal helper method to support object construction void PseudoTrap::make(const std::string& nam, double x1, double x2, double y1, double y2, double z, double r, bool atMinusZ) { double x = atMinusZ ? x1 : x2; @@ -1017,6 +1224,12 @@ void PseudoTrap::make(const std::string& nam, double x1, double x2, double y1, d _assign(solid,"","PseudoTrap", true); } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +PseudoTrap& PseudoTrap::setDimensions(const vector<double>& params) { + if ( params.size() ) {} + throw runtime_error("PseudoTrap: Special Boolean shape. setDimensions is not implemented!"); +} + /// Helper function to create poly hedron void PolyhedraRegular::make(const std::string& nam, int nsides, double rmin, double rmax, double zpos, double zneg, double start, double delta) { @@ -1029,6 +1242,18 @@ void PolyhedraRegular::make(const std::string& nam, int nsides, double rmin, dou //_setDimensions(¶ms[0]); } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +PolyhedraRegular& PolyhedraRegular::setDimensions(const vector<double>& params) { + auto pars = params; + if ( params.size() != 3 + 3*size_t(params[2]) ) { + invalidSetDimensionCall(*this,params); + } + pars[0] /= units::deg; + pars[1] /= units::deg; + _setDimensions(&pars[0]); + return *this; +} + /// Helper function to create poly hedron void Polyhedra::make(const std::string& nam, int nsides, double start, double delta, const vector<double>& z, const vector<double>& rmin, const vector<double>& rmax) { @@ -1052,6 +1277,18 @@ void Polyhedra::make(const std::string& nam, int nsides, double start, double de _assign(new TGeoPgon(&temp[0]), nam, "Polyhedra", false); } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +Polyhedra& Polyhedra::setDimensions(const vector<double>& params) { + auto pars = params; + if ( params.size() != 3 + 3*size_t(params[2]) ) { + invalidSetDimensionCall(*this,params); + } + pars[0] /= units::deg; + pars[1] /= units::deg; + _setDimensions(&pars[0]); + return *this; +} + /// Helper function to create the polyhedron void ExtrudedPolygon::make(const std::string& nam, const vector<double>& pt_x, @@ -1069,11 +1306,31 @@ void ExtrudedPolygon::make(const std::string& nam, solid->DefineSection(i, sec_z[i], sec_x[i], sec_y[i], sec_scale[i]); } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +ExtrudedPolygon& ExtrudedPolygon::setDimensions(const vector<double>& params) { + if ( params.size() != 1 + 4*size_t(params[0]) ) { + invalidSetDimensionCall(*this,params); + } + auto pars = params; + _setDimensions(&pars[0]); + return *this; +} + /// Creator method void EightPointSolid::make(const std::string& nam, double dz, const double* vtx) { _assign(new TGeoArb8(nam.c_str(), dz, (double*)vtx), "", "EightPointSolid", true); } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +EightPointSolid& EightPointSolid::setDimensions(const vector<double>& params) { + if ( params.size() != 17 ) { + invalidSetDimensionCall(*this,params); + } + auto pars = params; + _setDimensions(&pars[0]); + return *this; +} + /// Constructor to be used when creating a new object. Position is identity, Rotation is the identity rotation SubtractionSolid::SubtractionSolid(const Solid& shape1, const Solid& shape2) { TGeoSubtraction* sub = new TGeoSubtraction(shape1, shape2, detail::matrix::_identity(), detail::matrix::_identity()); @@ -1134,6 +1391,11 @@ SubtractionSolid::SubtractionSolid(const std::string& nam, const Solid& shape1, _assign(new TGeoCompositeShape(nam.c_str(), sub), "", "subtraction", true); } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +SubtractionSolid& SubtractionSolid::setDimensions(const vector<double>& ) { + throw runtime_error("SubtractionSolid: Boolean shape. setDimensions is not implemented!"); +} + /// Constructor to be used when creating a new object. Position is identity, Rotation is identity rotation UnionSolid::UnionSolid(const Solid& shape1, const Solid& shape2) { TGeoUnion* uni = new TGeoUnion(shape1, shape2, detail::matrix::_identity(), detail::matrix::_identity()); @@ -1194,6 +1456,11 @@ UnionSolid::UnionSolid(const std::string& nam, const Solid& shape1, const Solid& _assign(new TGeoCompositeShape(nam.c_str(), uni), "", "Union", true); } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +UnionSolid& UnionSolid::setDimensions(const vector<double>& ) { + throw runtime_error("UnionSolid: Boolean shape. setDimensions is not implemented!"); +} + /// Constructor to be used when creating a new object. Position is identity, Rotation is identity rotation IntersectionSolid::IntersectionSolid(const Solid& shape1, const Solid& shape2) { TGeoIntersection* inter = new TGeoIntersection(shape1, shape2, detail::matrix::_identity(), detail::matrix::_identity()); @@ -1254,6 +1521,12 @@ IntersectionSolid::IntersectionSolid(const std::string& nam, const Solid& shape1 _assign(new TGeoCompositeShape(nam.c_str(), inter), "", "Intersection", true); } +/// Set the shape dimensions. As for the TGeo shape, but angles in rad rather than degrees. +IntersectionSolid& IntersectionSolid::setDimensions(const vector<double>& ) { + throw runtime_error("IntersectionSolid: Boolean shape. setDimensions is not implemented!"); +} + + #define INSTANTIATE(X) template class dd4hep::Solid_type<X> INSTANTIATE(TGeoShape); diff --git a/DDG4/plugins/Geant4FieldTrackingSetup.cpp b/DDG4/plugins/Geant4FieldTrackingSetup.cpp index 4f568335d581b4839b6d65cb811eff24a09fbbce..d6a129e2d04ef1f2e753b641581bfb1a6359765b 100644 --- a/DDG4/plugins/Geant4FieldTrackingSetup.cpp +++ b/DDG4/plugins/Geant4FieldTrackingSetup.cpp @@ -199,11 +199,17 @@ int Geant4FieldTrackingSetup::execute(Detector& description) { G4MagneticField* mag_field = new sim::Geant4Field(fld); G4Mag_EqRhs* mag_equation = PluginService::Create<G4Mag_EqRhs*>(eq_typ,mag_field); G4EquationOfMotion* mag_eq = mag_equation; - G4MagIntegratorStepper* fld_stepper = PluginService::Create<G4MagIntegratorStepper*>(stepper_typ,mag_eq); + if ( nullptr == mag_eq ) { + mag_eq = PluginService::Create<G4EquationOfMotion*>(eq_typ,mag_field); + if ( nullptr == mag_eq ) { + except("FieldSetup", "Cannot create G4EquationOfMotion of type: %s.",eq_typ.c_str()); + } + } + G4MagIntegratorStepper* fld_stepper = PluginService::Create<G4MagIntegratorStepper*>(stepper_typ,mag_eq); if ( nullptr == fld_stepper ) { fld_stepper = PluginService::Create<G4MagIntegratorStepper*>(stepper_typ,mag_equation); if ( nullptr == fld_stepper ) { - printout(WARNING,"FieldSetup", "Cannot create stepper of type: %s. Taking Geant4 defaults.",stepper_typ.c_str()); + except("FieldSetup", "Cannot create stepper of type: %s.",stepper_typ.c_str()); } } chordFinder = new G4ChordFinder(mag_field,min_chord_step,fld_stepper);