/
githubmirror
/
cmssw
Обзор
Документация
Войти
/
githubmirror
/
cmssw
Код
Запросы
0
Пакеты
0
Релизы
0
Аналитика
Безопасность
master
SimG4Core/Geometry/src/DDG4SolidConverter.cc
509 строк
21 KB
Shahzad Malik Muzaffar
[SIMULATION] [LLVM21] Apply clang-tidy changes
26 ноя 2025, 17:10
26 ноя 2025, 17:10
9591585
Код
Авторство
О чём код?
#include "G4VSolid.hh" #include "SimG4Core/Geometry/interface/DDG4SolidConverter.h" #include "DetectorDescription/Core/interface/DDSolid.h" #include "FWCore/MessageLogger/interface/MessageLogger.h" #include "FWCore/Utilities/interface/Exception.h" #include "G4ios.hh" #include <CLHEP/Units/SystemOfUnits.h> //#define EDM_ML_DEBUG const std::vector<double> *DDG4SolidConverter::par_ = nullptr; G4RotationMatrix *DDG4SolidConverter::rot = nullptr; DDG4SolidConverter::DDG4SolidConverter() { // could also be done 'dynamically' from outside // would then need to have a 'register' method ... convDispatch_[DDSolidShape::ddbox] = DDG4SolidConverter::box; convDispatch_[DDSolidShape::ddtubs] = DDG4SolidConverter::tubs; convDispatch_[DDSolidShape::ddcuttubs] = DDG4SolidConverter::cuttubs; convDispatch_[DDSolidShape::ddtrap] = DDG4SolidConverter::trap; convDispatch_[DDSolidShape::ddcons] = DDG4SolidConverter::cons; convDispatch_[DDSolidShape::ddpolycone_rrz] = DDG4SolidConverter::polycone_rrz; convDispatch_[DDSolidShape::ddpolycone_rz] = DDG4SolidConverter::polycone_rz; convDispatch_[DDSolidShape::ddpolyhedra_rrz] = DDG4SolidConverter::polyhedra_rrz; convDispatch_[DDSolidShape::ddpolyhedra_rz] = DDG4SolidConverter::polyhedra_rz; convDispatch_[DDSolidShape::ddextrudedpolygon] = DDG4SolidConverter::extrudedpolygon; convDispatch_[DDSolidShape::ddtorus] = DDG4SolidConverter::torus; convDispatch_[DDSolidShape::ddunion] = DDG4SolidConverter::unionsolid; convDispatch_[DDSolidShape::ddintersection] = DDG4SolidConverter::intersection; convDispatch_[DDSolidShape::ddsubtraction] = DDG4SolidConverter::subtraction; convDispatch_[DDSolidShape::ddpseudotrap] = DDG4SolidConverter::pseudotrap; convDispatch_[DDSolidShape::ddtrunctubs] = DDG4SolidConverter::trunctubs; convDispatch_[DDSolidShape::ddsphere] = DDG4SolidConverter::sphere; convDispatch_[DDSolidShape::ddellipticaltube] = DDG4SolidConverter::ellipticaltube; } DDG4SolidConverter::~DDG4SolidConverter() {} G4VSolid *DDG4SolidConverter::convert(const DDSolid &solid) { if (!solid) { edm::LogError("SimG4CoreGeometry") << " DDG4SolidConverter::convert(..) found an undefined DDSolid " << solid.toString(); throw cms::Exception("SimG4CoreGeometry", "DDG4SolidConverter::convert(..) found an undefined DDSolid " + solid.toString()); } G4VSolid *result = nullptr; par_ = &(solid.parameters()); std::map<DDSolidShape, FNPTR>::iterator it = convDispatch_.find(solid.shape()); if (it != convDispatch_.end()) { result = it->second(solid); } else { throw cms::Exception("DetectorDescriptionFault") << "DDG4SolidConverter::convert: conversion failed for s=" << solid << "\n solid.shape()=" << DDSolidShapesName::name(solid.shape()) << G4endl; } return result; } #include "G4Box.hh" G4VSolid *DDG4SolidConverter::box(const DDSolid &solid) { edm::LogVerbatim("SimG4CoreGeometry") << "DDG4SolidConverter: box = " << solid; return new G4Box(solid.name().name(), (*par_)[0], (*par_)[1], (*par_)[2]); } #include "G4Tubs.hh" G4VSolid *DDG4SolidConverter::tubs(const DDSolid &solid) { edm::LogVerbatim("SimG4CoreGeometry") << "DDG4SolidConverter: tubs = " << solid; return new G4Tubs(solid.name().name(), (*par_)[1], // rmin (*par_)[2], // rmax (*par_)[0], // dzHalf (*par_)[3], // phiStart (*par_)[4]); // deltaPhi } #include "G4CutTubs.hh" G4VSolid *DDG4SolidConverter::cuttubs(const DDSolid &solid) { edm::LogVerbatim("SimG4CoreGeometry") << "DDG4SolidConverter: tubs = " << solid; return new G4CutTubs(solid.name().name(), (*par_)[1], // rmin (*par_)[2], // rmax (*par_)[0], // dzHalf (*par_)[3], // phiStart (*par_)[4], // deltaPhi G4ThreeVector((*par_)[5], (*par_)[6], (*par_)[7]), G4ThreeVector((*par_)[8], (*par_)[9], (*par_)[10])); } #include "G4Trap.hh" G4VSolid *DDG4SolidConverter::trap(const DDSolid &solid) { edm::LogVerbatim("SimG4CoreGeometry") << "DDG4SolidConverter: trap = " << solid; return new G4Trap(solid.name().name(), (*par_)[0], // pDz (*par_)[1], // theta (*par_)[2], // phi (*par_)[3], // y1 (*par_)[4], // x1 (*par_)[5], // x2 (*par_)[6], // alpha1 (*par_)[7], // y2 (*par_)[8], // x3 (*par_)[9], // x4 (*par_)[10]); // alpha2 } #include "G4Cons.hh" G4VSolid *DDG4SolidConverter::cons(const DDSolid &solid) { edm::LogVerbatim("SimG4CoreGeometry") << "DDG4SolidConverter: cons = " << solid; return new G4Cons(solid.name().name(), (*par_)[1], // rmin -z (*par_)[2], // rmax -z (*par_)[3], // rmin +z (*par_)[4], // rmax +z (*par_)[0], // zHalf (*par_)[5], // phistart (*par_)[6]); // deltaphi } #include "G4Polycone.hh" G4VSolid *DDG4SolidConverter::polycone_rz(const DDSolid &solid) { edm::LogVerbatim("SimG4CoreGeometry") << "DDG4SolidConverter: pcon_rz = " << solid; std::vector<double> r; std::vector<double> z; std::vector<double>::const_iterator i = (*par_).begin() + 2; int count = 0; for (; i != (*par_).end(); ++i) { edm::LogVerbatim("SimG4CoreGeometry") << " z=" << *i / CLHEP::cm; z.push_back(*i); ++i; edm::LogVerbatim("SimG4CoreGeometry") << " r=" << *i / CLHEP::cm; r.push_back(*i); count++; } edm::LogVerbatim("SimG4CoreGeometry") << "sp=" << (*par_)[0] / CLHEP::deg << " ep=" << (*par_)[1] / CLHEP::deg; #ifdef EDM_ML_DEBUG G4cout << "### Polycone_RZ: " << "sp=" << (*par_)[0] / CLHEP::deg << " ep=" << (*par_)[1] / CLHEP::deg << " N= " << count << G4endl; for (int i = 0; i < count; ++i) { G4cout << " R= " << r[i] << " Z= " << z[i] << G4endl; } #endif return new G4Polycone(solid.name().name(), (*par_)[0], (*par_)[1], // start,delta-phi count, // numRZ &(r[0]), &(z[0])); } G4VSolid *DDG4SolidConverter::polycone_rrz(const DDSolid &solid) { edm::LogVerbatim("SimG4CoreGeometry") << "DDG4SolidConverter: pcon_rrz = " << solid; std::vector<double> z_p; std::vector<double> rmin_p; std::vector<double> rmax_p; std::vector<double>::const_iterator i = par_->begin() + 2; int count = 0; for (; i != par_->end(); ++i) { edm::LogVerbatim("SimG4CoreGeometry") << "z=" << *i / CLHEP::cm; z_p.push_back(*i); ++i; edm::LogVerbatim("SimG4CoreGeometry") << "rmin=" << *i / CLHEP::cm; rmin_p.push_back(*i); ++i; edm::LogVerbatim("SimG4CoreGeometry") << "rmax=" << *i / CLHEP::cm; rmax_p.push_back(*i); count++; } edm::LogVerbatim("SimG4CoreGeometry") << "sp=" << (*par_)[0] / CLHEP::deg << " ep=" << (*par_)[1] / CLHEP::deg; #ifdef EDM_ML_DEBUG G4cout << "### Polycone_RRZ: " << "sp=" << (*par_)[0] / CLHEP::deg << " ep=" << (*par_)[1] / CLHEP::deg << " N= " << count << G4endl; for (int i = 0; i < count; ++i) { G4cout << " R1= " << rmin_p[i] << " R1= " << rmax_p[i] << " Z= " << z_p[i] << G4endl; } #endif return new G4Polycone(solid.name().name(), (*par_)[0], (*par_)[1], // start,delta-phi count, // sections &(z_p[0]), &(rmin_p[0]), &(rmax_p[0])); } #include "G4Polyhedra.hh" G4VSolid *DDG4SolidConverter::polyhedra_rz(const DDSolid &solid) { edm::LogVerbatim("SimG4CoreGeometry") << "DDG4SolidConverter: phed_rz = " << solid; std::vector<double> r; std::vector<double> z; std::vector<double>::const_iterator i = par_->begin() + 3; int count = 0; for (; i != par_->end(); ++i) { z.push_back(*i); ++i; r.push_back(*i); count++; } return new G4Polyhedra(solid.name().name(), (*par_)[1], (*par_)[2], int((*par_)[0]), // start,delta-phi;sides count, // numRZ &(r[0]), &(z[0])); } G4VSolid *DDG4SolidConverter::polyhedra_rrz(const DDSolid &solid) { edm::LogVerbatim("SimG4CoreGeometry") << "DDG4SolidConverter: phed_rrz = " << solid; std::vector<double> z_p; std::vector<double> rmin_p; std::vector<double> rmax_p; std::vector<double>::const_iterator i = par_->begin() + 3; int count = 0; for (; i != par_->end(); ++i) { edm::LogVerbatim("SimG4CoreGeometry") << "z=" << *i / CLHEP::cm; z_p.push_back(*i); ++i; edm::LogVerbatim("SimG4CoreGeometry") << "rmin=" << *i / CLHEP::cm; rmin_p.push_back(*i); ++i; edm::LogVerbatim("SimG4CoreGeometry") << "rmax=" << *i / CLHEP::cm; rmax_p.push_back(*i); count++; } edm::LogVerbatim("SimG4CoreGeometry") << "sp=" << (*par_)[0] / CLHEP::deg << " ep=" << (*par_)[1] / CLHEP::deg; return new G4Polyhedra(solid.name().name(), (*par_)[1], (*par_)[2], int((*par_)[0]), // start,delta-phi,sides count, // sections &(z_p[0]), &(rmin_p[0]), &(rmax_p[0])); } #include "G4ExtrudedSolid.hh" G4VSolid *DDG4SolidConverter::extrudedpolygon(const DDSolid &solid) { edm::LogVerbatim("SimG4CoreGeometry") << "DDG4SolidConverter: extr_pgon = " << solid; std::vector<double> x = static_cast<DDExtrudedPolygon>(solid).xVec(); std::vector<double> y = static_cast<DDExtrudedPolygon>(solid).yVec(); std::vector<double> z = static_cast<DDExtrudedPolygon>(solid).zVec(); std::vector<double> zx = static_cast<DDExtrudedPolygon>(solid).zxVec(); std::vector<double> zy = static_cast<DDExtrudedPolygon>(solid).zyVec(); std::vector<double> zs = static_cast<DDExtrudedPolygon>(solid).zscaleVec(); std::vector<G4TwoVector> polygon; std::vector<G4ExtrudedSolid::ZSection> zsections; polygon.reserve(x.size()); for (unsigned int it = 0; it < x.size(); ++it) polygon.emplace_back(x[it], y[it]); zsections.reserve(z.size()); for (unsigned int it = 0; it < z.size(); ++it) zsections.emplace_back(z[it], G4TwoVector(zx[it], zy[it]), zs[it]); return new G4ExtrudedSolid(solid.name().name(), polygon, zsections); } #include "G4Torus.hh" G4VSolid *DDG4SolidConverter::torus(const DDSolid &solid) { edm::LogVerbatim("SimG4CoreGeometry") << "DDG4SolidConverter: torus = " << solid; return new G4Torus(solid.name().name(), (*par_)[0], // rmin (*par_)[1], // rmax (*par_)[2], // Rtor (*par_)[3], // phiStart (*par_)[4]); // deltaPhi } #include "G4UnionSolid.hh" G4VSolid *DDG4SolidConverter::unionsolid(const DDSolid &solid) { edm::LogVerbatim("SimG4CoreGeometry") << "DDG4SolidConverter: unionsolid = " << solid.name(); G4UnionSolid *us = nullptr; DDBooleanSolid bs(solid); if (bs) { edm::LogVerbatim("SimG4CoreGeometry") << "SolidA=" << bs.solidA(); G4VSolid *sa = DDG4SolidConverter().convert(bs.solidA()); edm::LogVerbatim("SimG4CoreGeometry") << "SolidB=" << bs.solidB(); G4VSolid *sb = DDG4SolidConverter().convert(bs.solidB()); edm::LogVerbatim("SimG4CoreGeometry") << " name:" << solid.name() << " t=" << bs.translation() << std::flush; edm::LogVerbatim("SimG4CoreGeometry") << " " << bs.rotation().rotation().Inverse() << std::flush; std::vector<double> tdbl(9); bs.rotation().rotation().Inverse().GetComponents(tdbl.begin(), tdbl.end()); CLHEP::HepRep3x3 temprep(tdbl[0], tdbl[1], tdbl[2], tdbl[3], tdbl[4], tdbl[5], tdbl[6], tdbl[7], tdbl[8]); CLHEP::Hep3Vector temphvec(bs.translation().X(), bs.translation().Y(), bs.translation().Z()); us = new G4UnionSolid(solid.name().name(), sa, sb, new CLHEP::HepRotation(temprep), temphvec); } // else? return us; } #include "G4SubtractionSolid.hh" #include <sstream> G4VSolid *DDG4SolidConverter::subtraction(const DDSolid &solid) { edm::LogVerbatim("SimG4CoreGeometry") << "DDG4SolidConverter: subtraction = " << solid; G4SubtractionSolid *us = nullptr; DDBooleanSolid bs(solid); if (bs) { G4VSolid *sa = DDG4SolidConverter().convert(bs.solidA()); G4VSolid *sb = DDG4SolidConverter().convert(bs.solidB()); edm::LogVerbatim("SimG4CoreGeometry") << " name:" << solid.name() << " t=" << bs.translation() << std::flush; edm::LogVerbatim("SimG4CoreGeometry") << " " << bs.rotation().rotation().Inverse() << std::flush; std::vector<double> tdbl(9); bs.rotation().rotation().Inverse().GetComponents(tdbl.begin(), tdbl.end()); CLHEP::HepRep3x3 temprep(tdbl[0], tdbl[1], tdbl[2], tdbl[3], tdbl[4], tdbl[5], tdbl[6], tdbl[7], tdbl[8]); CLHEP::Hep3Vector temphvec(bs.translation().X(), bs.translation().Y(), bs.translation().Z()); us = new G4SubtractionSolid(solid.name().name(), sa, sb, new CLHEP::HepRotation(temprep), temphvec); } return us; } #include "G4IntersectionSolid.hh" G4VSolid *DDG4SolidConverter::intersection(const DDSolid &solid) { edm::LogVerbatim("SimG4CoreGeometry") << "DDG4SolidConverter: intersection = " << solid; G4IntersectionSolid *us = nullptr; DDBooleanSolid bs(solid); if (bs) { G4VSolid *sa = DDG4SolidConverter().convert(bs.solidA()); G4VSolid *sb = DDG4SolidConverter().convert(bs.solidB()); edm::LogVerbatim("SimG4CoreGeometry") << " name:" << solid.name() << " t=" << bs.translation() << std::flush; edm::LogVerbatim("SimG4CoreGeometry") << " " << bs.rotation().rotation().Inverse() << std::flush; std::vector<double> tdbl(9); bs.rotation().rotation().Inverse().GetComponents(tdbl.begin(), tdbl.end()); CLHEP::HepRep3x3 temprep(tdbl[0], tdbl[1], tdbl[2], tdbl[3], tdbl[4], tdbl[5], tdbl[6], tdbl[7], tdbl[8]); CLHEP::Hep3Vector temphvec(bs.translation().X(), bs.translation().Y(), bs.translation().Z()); us = new G4IntersectionSolid(solid.name().name(), sa, sb, new CLHEP::HepRotation(temprep), temphvec); } return us; } #include "G4Trd.hh" G4VSolid *DDG4SolidConverter::pseudotrap(const DDSolid &solid) { if (nullptr == rot) { rot = new G4RotationMatrix; rot->rotateX(90. * CLHEP::deg); } edm::LogVerbatim("SimG4CoreGeometry") << "DDG4SolidConverter: pseudoTrap = " << solid; G4Trd *trap = nullptr; G4Tubs *tubs = nullptr; G4VSolid *result = nullptr; DDPseudoTrap pt(solid); // pt...PseudoTrap double r = pt.radius(); bool atMinusZ = pt.atMinusZ(); double x = 0; double h = 0; bool intersec = false; // union or intersection solid if (pt.atMinusZ()) { x = pt.x1(); // tubs radius } else { x = pt.x2(); // tubs radius } double openingAngle = 2. * asin(x / std::abs(r)); // trap = new G4Trd(solid.name().name(), double displacement = 0; double startPhi = 0; /* calculate the displacement of the tubs w.r.t. to the trap, determine the opening angle of the tubs */ double delta = sqrt(r * r - x * x); if (r < 0 && std::abs(r) >= x) { intersec = true; // intersection solid h = pt.y1() < pt.y2() ? pt.y2() : pt.y1(); // tubs half height h += h / 20.; // enlarge a bit - for subtraction solid if (atMinusZ) { displacement = -pt.halfZ() - delta; startPhi = 270. * CLHEP::deg - openingAngle / 2.; } else { displacement = pt.halfZ() + delta; startPhi = 90. * CLHEP::deg - openingAngle / 2.; } } else if (r > 0 && std::abs(r) >= x) { if (atMinusZ) { displacement = -pt.halfZ() + delta; startPhi = 90. * CLHEP::deg - openingAngle / 2.; h = pt.y1(); } else { displacement = pt.halfZ() - delta; startPhi = 270. * CLHEP::deg - openingAngle / 2.; h = pt.y2(); } } else { throw cms::Exception("DetectorDescriptionFault", "Check parameters of the PseudoTrap! name=" + pt.name().name()); } G4ThreeVector displ(0., 0., displacement); // displacement of the tubs w.r.t. trap edm::LogVerbatim("SimG4CoreGeometry") << "DDSolidConverter::pseudotrap(): displacement=" << displacement << " openingAngle=" << openingAngle / CLHEP::deg << " x=" << x << " h=" << h; // Now create two solids (trd & tubs), and a boolean solid out of them std::string name = pt.name().name(); trap = new G4Trd(name, pt.x1(), pt.x2(), pt.y1(), pt.y2(), pt.halfZ()); tubs = new G4Tubs(name, 0., // rMin std::abs(r), // rMax h, // half height startPhi, // start angle openingAngle); if (intersec) { result = new G4SubtractionSolid(name, trap, tubs, rot, displ); } else { G4VSolid *tubicCap = new G4SubtractionSolid( name, tubs, new G4Box(name, 1.1 * x, sqrt(r * r - x * x), 1.1 * h), nullptr, G4ThreeVector()); result = new G4UnionSolid(name, trap, tubicCap, rot, displ); } return result; } G4VSolid *DDG4SolidConverter::trunctubs(const DDSolid &solid) { // truncated tube-section: a boolean subtraction solid: // from a tube-section a box is subtracted according // to the given parameters edm::LogVerbatim("SimG4CoreGeometry") << "MantisConverter: solidshape=" << DDSolidShapesName::name(solid.shape()) << " " << solid; edm::LogVerbatim("SimG4CoreGeometry") << "before"; DDTruncTubs tt(solid); edm::LogVerbatim("SimG4CoreGeometry") << "after"; double rIn(tt.rIn()), rOut(tt.rOut()), zHalf(tt.zHalf()), startPhi(tt.startPhi()), deltaPhi(tt.deltaPhi()), cutAtStart(tt.cutAtStart()), cutAtDelta(tt.cutAtDelta()); bool cutInside(bool(tt.cutInside())); std::string name = tt.name().name(); // check the parameters if (rIn <= 0 || rOut <= 0 || cutAtStart <= 0 || cutAtDelta <= 0) { throw cms::Exception( "DetectorDescriptionFault", "TruncTubs " + std::string(tt.name().fullname()) + ": 0 <= rIn,cutAtStart,rOut,cutAtDelta,rOut violated!"); } if (rIn >= rOut) { throw cms::Exception("DetectorDescriptionFault", "TruncTubs " + std::string(tt.name().fullname()) + ": rIn<rOut violated!"); } if (startPhi != 0.) { throw cms::Exception("DetectorDescriptionFault", "TruncTubs " + std::string(tt.name().fullname()) + ": startPhi != 0 not supported!"); } startPhi = 0.; double r(cutAtStart), R(cutAtDelta); G4VSolid *result(nullptr); G4VSolid *tubs = new G4Tubs(name, rIn, rOut, zHalf, startPhi, deltaPhi); edm::LogVerbatim("SimG4CoreGeometry") << "G4Tubs: " << rIn / CLHEP::cm << ' ' << rOut / CLHEP::cm << ' ' << zHalf / CLHEP::cm << ' ' << startPhi / CLHEP::deg << ' ' << deltaPhi / CLHEP::deg; edm::LogVerbatim("SimG4CoreGeometry") << solid; // length & hight of the box double boxX(rOut * std::sqrt(2) * 1.1), boxY(rOut * std::sqrt(2) * 1.1); // exaggerate dimensions - does not matter, it's subtracted! // width of the box > width of the tubs double boxZ(1.1 * zHalf); // angle of the box w.r.t. tubs double cath = r - R * cos(deltaPhi); double hypo = sqrt(r * r + R * R - 2. * r * R * cos(deltaPhi)); double cos_alpha = cath / hypo; double alpha = -acos(cos_alpha); edm::LogVerbatim("SimG4CoreGeometry") << "cath=" << cath / CLHEP::cm; edm::LogVerbatim("SimG4CoreGeometry") << "hypo=" << hypo / CLHEP::cm; edm::LogVerbatim("SimG4CoreGeometry") << "al=" << acos(cath / hypo) / CLHEP::deg; edm::LogVerbatim("SimG4CoreGeometry") << "deltaPhi=" << deltaPhi / CLHEP::deg << "\n" << "r=" << r / CLHEP::cm << "\n" << "R=" << R / CLHEP::cm; edm::LogVerbatim("SimG4CoreGeometry") << "alpha=" << alpha / CLHEP::deg; // rotationmatrix of box w.r.t. tubs G4RotationMatrix *rot = new G4RotationMatrix; rot->rotateZ(-alpha); edm::LogVerbatim("SimG4CoreGeometry") << (*rot); // center point of the box double xBox; if (!cutInside) { xBox = r + boxY / sin(std::abs(alpha)); } else { xBox = -(boxY / sin(std::abs(alpha)) - r); } G4ThreeVector trans(xBox, 0., 0.); edm::LogVerbatim("SimG4CoreGeometry") << "trans=" << trans; G4VSolid *box = new G4Box(name, boxX, boxY, boxZ); result = new G4SubtractionSolid(name, tubs, box, rot, trans); return result; } #include "G4Sphere.hh" G4VSolid *DDG4SolidConverter::sphere(const DDSolid &solid) { edm::LogVerbatim("SimG4CoreGeometry") << "DDG4SolidConverter: sphere = " << solid; DDSphere sp(solid); return new G4Sphere(solid.name().name(), sp.innerRadius(), sp.outerRadius(), sp.startPhi(), sp.deltaPhi(), sp.startTheta(), sp.deltaTheta()); } #include "G4EllipticalTube.hh" G4VSolid *DDG4SolidConverter::ellipticaltube(const DDSolid &solid) { edm::LogVerbatim("SimG4CoreGeometry") << "DDG4SolidConverter: ellipticaltube = " << solid; DDEllipticalTube sp(solid); return new G4EllipticalTube(solid.name().name(), sp.xSemiAxis(), sp.ySemiAxis(), sp.zHeight()); }