File indexing completed on 2025-01-18 09:14:19
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016 #include <DDG4/Geant4DetectorConstruction.h>
0017
0018
0019 namespace dd4hep {
0020
0021
0022 namespace sim {
0023
0024
0025
0026
0027
0028
0029
0030
0031
0032
0033
0034 class Geant4DetectorGeometryConstruction : public Geant4DetectorConstruction {
0035
0036 unsigned long m_dumpHierarchy = 0;
0037
0038 bool m_debugMaterials = false;
0039
0040 bool m_debugElements = false;
0041
0042 bool m_debugShapes = false;
0043
0044 bool m_debugVolumes = false;
0045
0046 bool m_debugPlacements = false;
0047
0048 bool m_debugReflections = false;
0049
0050 bool m_debugRegions = false;
0051
0052 bool m_debugLimits = false;
0053
0054 bool m_debugSurfaces = false;
0055
0056
0057 bool m_printPlacements = false;
0058
0059 bool m_printSensitives = false;
0060
0061
0062 int m_geoInfoPrintLevel;
0063
0064 std::string m_dumpGDML;
0065
0066
0067 int writeGDML(const char* gdml_output);
0068
0069 int printVolumeObj(const char* vol_path, PlacedVolume pv, int flg);
0070
0071 int printVolumeTree(const char* vol_path);
0072
0073 int printVolTree(const char* vol_path);
0074
0075 int printG4Tree(const char* vol_path);
0076
0077 int printVolume(const char* vol_path);
0078
0079 int checkVolume(const char* vol_path);
0080
0081 int printMaterial(const char* mat_name);
0082
0083 std::pair<std::string, PlacedVolume> resolve_path(const char* vol_path) const;
0084 void printG4(const std::string& prefix, const G4VPhysicalVolume* g4pv) const;
0085
0086 public:
0087
0088 Geant4DetectorGeometryConstruction(Geant4Context* ctxt, const std::string& nam);
0089
0090 virtual ~Geant4DetectorGeometryConstruction();
0091
0092 void constructGeo(Geant4DetectorConstructionContext* ctxt) override;
0093
0094 virtual void installCommandMessenger() override;
0095 };
0096 }
0097 }
0098
0099
0100
0101 #include <DD4hep/InstanceCount.h>
0102 #include <DD4hep/DetectorTools.h>
0103 #include <DD4hep/DD4hepUnits.h>
0104 #include <DD4hep/Printout.h>
0105 #include <DD4hep/Detector.h>
0106
0107 #include <DDG4/Geant4HierarchyDump.h>
0108 #include <DDG4/Geant4UIMessenger.h>
0109 #include <DDG4/Geant4Converter.h>
0110 #include <DDG4/Geant4Kernel.h>
0111 #include <DDG4/Factories.h>
0112
0113 #include <TGeoScaledShape.h>
0114
0115
0116 #include <G4LogicalVolume.hh>
0117 #include <G4PVPlacement.hh>
0118 #include <G4Material.hh>
0119 #include <G4Version.hh>
0120 #include <G4VSolid.hh>
0121 #include <CLHEP/Units/SystemOfUnits.h>
0122
0123 #ifndef GEANT4_NO_GDML
0124 #include <G4GDMLParser.hh>
0125 #endif
0126
0127 #include <cmath>
0128
0129 using namespace dd4hep::sim;
0130 DECLARE_GEANT4ACTION(Geant4DetectorGeometryConstruction)
0131
0132
0133 Geant4DetectorGeometryConstruction::Geant4DetectorGeometryConstruction(Geant4Context* ctxt, const std::string& nam)
0134 : Geant4DetectorConstruction(ctxt,nam)
0135 {
0136 declareProperty("DebugMaterials", m_debugMaterials);
0137 declareProperty("DebugElements", m_debugElements);
0138 declareProperty("DebugShapes", m_debugShapes);
0139 declareProperty("DebugVolumes", m_debugVolumes);
0140 declareProperty("DebugPlacements", m_debugPlacements);
0141 declareProperty("DebugReflections", m_debugReflections);
0142 declareProperty("DebugRegions", m_debugRegions);
0143 declareProperty("DebugLimits", m_debugLimits);
0144 declareProperty("DebugSurfaces", m_debugSurfaces);
0145
0146 declareProperty("PrintPlacements", m_printPlacements);
0147 declareProperty("PrintSensitives", m_printSensitives);
0148 declareProperty("GeoInfoPrintLevel", m_geoInfoPrintLevel = DEBUG);
0149
0150 declareProperty("DumpHierarchy", m_dumpHierarchy);
0151 declareProperty("DumpGDML", m_dumpGDML="");
0152 InstanceCount::increment(this);
0153 }
0154
0155
0156 Geant4DetectorGeometryConstruction::~Geant4DetectorGeometryConstruction() {
0157 InstanceCount::decrement(this);
0158 }
0159
0160
0161 void Geant4DetectorGeometryConstruction::constructGeo(Geant4DetectorConstructionContext* ctxt) {
0162 Geant4Mapping& g4map = Geant4Mapping::instance();
0163 DetElement world = ctxt->description.world();
0164 Geant4Converter conv(ctxt->description, outputLevel());
0165 conv.debugMaterials = m_debugMaterials;
0166 conv.debugElements = m_debugElements;
0167 conv.debugShapes = m_debugShapes;
0168 conv.debugVolumes = m_debugVolumes;
0169 conv.debugRegions = m_debugRegions;
0170 conv.debugSurfaces = m_debugSurfaces;
0171 conv.debugPlacements = m_debugPlacements;
0172 conv.debugReflections = m_debugReflections;
0173 conv.debugLimits = m_debugLimits;
0174 conv.printPlacements = m_printPlacements;
0175 conv.printSensitives = m_printSensitives;
0176
0177 ctxt->geometry = conv.create(world).detach();
0178 ctxt->geometry->printLevel = outputLevel();
0179 g4map.attach(ctxt->geometry);
0180 G4VPhysicalVolume* w = ctxt->geometry->world();
0181
0182 context()->kernel().setWorld(w);
0183
0184 g4map.volumeManager();
0185 if ( m_dumpHierarchy != 0 ) {
0186 Geant4HierarchyDump dmp(ctxt->description, m_dumpHierarchy);
0187 dmp.dump("",w);
0188 }
0189 ctxt->world = w;
0190 if ( !m_dumpGDML.empty() ) writeGDML(m_dumpGDML.c_str());
0191 else if ( ::getenv("DUMP_GDML") ) writeGDML(::getenv("DUMP_GDML"));
0192 enableUI();
0193 }
0194
0195 std::pair<std::string, dd4hep::PlacedVolume>
0196 Geant4DetectorGeometryConstruction::resolve_path(const char* vol_path) const {
0197 std::string p = vol_path;
0198 Detector& det = context()->kernel().detectorDescription();
0199 PlacedVolume top = det.world().placement();
0200 PlacedVolume pv = detail::tools::findNode(top, p);
0201 if ( !pv.isValid() ) {
0202 DetElement de = detail::tools::findElement(det, p);
0203 if ( de.isValid() ) {
0204 pv = de.placement();
0205 p = detail::tools::placementPath(de);
0206 }
0207 }
0208 return make_pair(p,pv);
0209 }
0210
0211
0212 int Geant4DetectorGeometryConstruction::printMaterial(const char* mat_name) {
0213 if ( mat_name ) {
0214 auto& g4map = Geant4Mapping::instance().data().g4Materials;
0215 for ( auto it = g4map.begin(); it != g4map.end(); ++it ) {
0216 const auto* mat = (*it).second;
0217 if ( mat->GetName() == mat_name ) {
0218 std::stringstream output;
0219 const auto* ion = mat->GetIonisation();
0220 printP2("+++ Dump of GEANT4 material: %s", mat_name);
0221 output << mat;
0222 if ( ion ) {
0223 output << " MEE: ";
0224 output << std::setprecision(12);
0225 output << ion->GetMeanExcitationEnergy()/CLHEP::eV;
0226 output << " [eV]";
0227 }
0228 else
0229 output << " MEE: UNKNOWN";
0230 always("+++ printMaterial: \n%s\n", output.str().c_str());
0231 return 1;
0232 }
0233 }
0234 warning("+++ printMaterial: FAILED to find the material %s", mat_name);
0235 }
0236 warning("+++ printMaterial: Property materialName not set!");
0237 return 0;
0238 }
0239
0240
0241 int Geant4DetectorGeometryConstruction::printVolumeObj(const char* vol_path, PlacedVolume pv, int flg) {
0242 if ( pv.isValid() ) {
0243 const G4LogicalVolume* vol = 0;
0244 auto& g4map = Geant4Mapping::instance().data();
0245 auto pit = g4map.g4Placements.find(pv.ptr());
0246 auto vit = g4map.g4Volumes.find(pv.volume());
0247 warning("+++ printVolume: %s", vol_path);
0248 if ( vit != g4map.g4Volumes.end() ) {
0249 vol = (*vit).second;
0250 auto* sol = vol->GetSolid();
0251 const auto* mat = vol->GetMaterial();
0252 const auto* ion = mat->GetIonisation();
0253 Solid sh = pv.volume().solid();
0254 if ( flg ) {
0255 printP2("++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++");
0256 printP2("+++ Dump of GEANT4 solid: %s", vol_path);
0257 }
0258 std::stringstream output;
0259 if ( flg ) {
0260 output << mat;
0261 if ( ion ) {
0262 output << " MEE: ";
0263 output << std::setprecision(12);
0264 output << ion->GetMeanExcitationEnergy()/CLHEP::eV;
0265 output << " [eV]";
0266 }
0267 else
0268 output << " MEE: UNKNOWN";
0269 }
0270 if ( flg ) {
0271 output << std::endl << *sol;
0272 printP2("%s", output.str().c_str());
0273 printP2("+++ Dump of ROOT solid: %s", vol_path);
0274 sh->InspectShape();
0275 if ( sh->IsA() == TGeoScaledShape::Class() ) {
0276 TGeoScaledShape* scaled = (TGeoScaledShape*)sh.ptr();
0277 const Double_t* scale = scaled->GetScale()->GetScale();
0278 double dot = scale[0]*scale[1]*scale[2];
0279 printP2("+++ TGeoScaledShape: %8.3g %8.3g %8.3g [%s]", scale[0], scale[1], scale[2],
0280 dot > 0e0 ? "RIGHT handed" : "LEFT handed");
0281 }
0282 else if ( pit != g4map.g4Placements.end() ) {
0283 const G4VPhysicalVolume* pl = (*pit).second;
0284 const G4RotationMatrix* rot = pl->GetRotation();
0285 const G4ThreeVector& tr = pl->GetTranslation();
0286 G4Transform3D transform(rot ? *rot : G4RotationMatrix(), tr);
0287 HepGeom::Scale3D sc;
0288 HepGeom::Rotate3D rr;
0289 G4Translate3D tt;
0290 transform.getDecomposition(sc,rr,tt);
0291 double dot = sc(0,0)*sc(1,1)*sc(2,2);
0292 printP2("+++ TGeoShape: %8.3g %8.3g %8.3g [%s]", sc(0,0), sc(1,1), sc(2,2),
0293 dot > 0e0 ? "RIGHT handed" : "LEFT handed");
0294 }
0295 const TGeoMatrix* matrix = pv->GetMatrix();
0296 printP2("+++ TGeoMatrix: %s",
0297 matrix->TestBit(TGeoMatrix::kGeoReflection) ? "LEFT handed" : "RIGHT handed");
0298 printP2("+++ Shape: %s cubic volume: %8.3g mm^3 area: %8.3g mm^2",
0299 sol->GetName().c_str(), sol->GetCubicVolume(), sol->GetSurfaceArea());
0300 }
0301 return 1;
0302 }
0303 else {
0304 auto ai = g4map.g4AssemblyVolumes.find(pv.ptr());
0305 if ( ai != g4map.g4AssemblyVolumes.end() ) {
0306 Volume v = pv.volume();
0307 warning("+++ printVolume: volume %s is an assembly...need to resolve imprint",vol_path);
0308 for(Int_t i=0; i < v->GetNdaughters(); ++i) {
0309 TGeoNode* dau_nod = v->GetNode(i);
0310 std::string p = vol_path + std::string("/") + dau_nod->GetName();
0311 printVolumeObj(p.c_str(), dau_nod, flg);
0312 }
0313 return 0;
0314 }
0315 }
0316 warning("+++ printVolume: FAILED to find the volume %s in geant4 mapping...",vol_path);
0317 return 0;
0318 }
0319 warning("+++ printVolume: FAILED to dump invalid volume",vol_path);
0320 return 0;
0321 }
0322
0323
0324 int Geant4DetectorGeometryConstruction::printVolume(const char* vol_path) {
0325 if ( vol_path ) {
0326 auto physVol = resolve_path(vol_path);
0327 if ( physVol.second.isValid() ) {
0328 return printVolumeObj(vol_path, physVol.second, ~0x0);
0329 }
0330 }
0331 warning("+++ printVolume: Property VolumePath not set. [Ignored]");
0332 return 0;
0333 }
0334
0335
0336 int Geant4DetectorGeometryConstruction::printVolumeTree(const char* vol_path) {
0337 if ( vol_path ) {
0338 auto [p, pv] = resolve_path(vol_path);
0339 if ( pv.isValid() ) {
0340 if ( printVolumeObj(p.c_str(), pv, ~0x0) ) {
0341 TGeoVolume* vol = pv->GetVolume();
0342 for(Int_t i=0; i < vol->GetNdaughters(); ++i) {
0343 PlacedVolume dau_pv(vol->GetNode(i));
0344 std::string path = (p + "/") + dau_pv.name();
0345 if ( printVolumeTree(path.c_str()) ) {
0346 }
0347 }
0348 }
0349 return 1;
0350 }
0351 }
0352 warning("+++ printVolume: Could not access Volume/DetElement '%s'", vol_path ? vol_path : "UNKNOWN");
0353 return 0;
0354 }
0355
0356 int Geant4DetectorGeometryConstruction::printVolTree(const char* vol_path) {
0357 if ( vol_path ) {
0358 auto [p, pv] = resolve_path(vol_path);
0359 if ( pv.isValid() ) {
0360 if ( printVolumeObj(p.c_str(), pv, 0) ) {
0361 TGeoVolume* vol = pv->GetVolume();
0362 for(Int_t i=0; i < vol->GetNdaughters(); ++i) {
0363 PlacedVolume dau_pv(vol->GetNode(i));
0364 std::string path = (p + "/") + dau_pv.name();
0365 if ( printVolTree(path.c_str()) ) {
0366 }
0367 }
0368 }
0369 return 1;
0370 }
0371 }
0372 warning("+++ printVolume: Could not access Volume/DetElement '%s'", vol_path ? vol_path : "UNKNOWN");
0373 return 0;
0374 }
0375
0376
0377 int Geant4DetectorGeometryConstruction::checkVolume(const char* vol_path) {
0378 if ( vol_path ) {
0379 auto physVol = resolve_path(vol_path);
0380 auto pv = physVol.second;
0381 if ( pv.isValid() ) {
0382 auto& g4map = Geant4Mapping::instance().data().g4Volumes;
0383 auto it = g4map.find(pv.volume());
0384 if ( it != g4map.end() ) {
0385 const G4LogicalVolume* vol = (*it).second;
0386 auto* g4_sol = vol->GetSolid();
0387 Box rt_sol = pv.volume().solid();
0388 printP2("Geant4 Shape: %s cubic volume: %8.3g mm^3 area: %8.3g mm^2",
0389 g4_sol->GetName().c_str(), g4_sol->GetCubicVolume(), g4_sol->GetSurfaceArea());
0390 #if G4VERSION_NUMBER>=1040
0391 G4ThreeVector pMin, pMax;
0392 double conv = (dd4hep::centimeter/CLHEP::centimeter)/2.0;
0393 g4_sol->BoundingLimits(pMin,pMax);
0394 printP2("Geant4 Bounding box extends: %8.3g %8.3g %8.3g",
0395 (pMax.x()-pMin.x())*conv, (pMax.y()-pMin.y())*conv, (pMax.z()-pMin.z())*conv);
0396 #endif
0397 printP2("ROOT Bounding box dimensions: %8.3g %8.3g %8.3g",
0398 rt_sol->GetDX(), rt_sol->GetDY(), rt_sol->GetDZ());
0399
0400 return 1;
0401 }
0402 }
0403 warning("+++ checkVolume: FAILED to find the volume %s from the top volume",vol_path);
0404 }
0405 warning("+++ checkVolume: Property VolumePath not set. [Ignored]");
0406 return 0;
0407 }
0408
0409
0410 int Geant4DetectorGeometryConstruction::writeGDML(const char* output) {
0411 #ifdef GEANT4_NO_GDML
0412 warning("+++ writeGDML: GDML not found in the present Geant4 build! Output: %s not written", output);
0413 #else
0414 G4VPhysicalVolume* w = context()->world();
0415 if ( output && ::strlen(output) > 0 && output != m_dumpGDML.c_str() )
0416 m_dumpGDML = output;
0417
0418 if ( !m_dumpGDML.empty() ) {
0419 G4GDMLParser parser;
0420 parser.Write(m_dumpGDML.c_str(), w);
0421 info("+++ writeGDML: Wrote GDML file: %s", m_dumpGDML.c_str());
0422 return 1;
0423 }
0424 else {
0425 const char* gdml_dmp = ::getenv("DUMP_GDML");
0426 if ( gdml_dmp ) {
0427 G4GDMLParser parser;
0428 parser.Write(gdml_dmp, w);
0429 info("+++ writeGDML: Wrote GDML file: %s", gdml_dmp);
0430 return 1;
0431 }
0432 }
0433 warning("+++ writeGDML: Neither property DumpGDML nor environment DUMP_GDML set. No file written!");
0434 #endif
0435 return 0;
0436 }
0437
0438 void Geant4DetectorGeometryConstruction::printG4(const std::string& prefix, const G4VPhysicalVolume* g4pv) const {
0439 std::string path = prefix + "/";
0440 printP2("+++ GEANT4 volume: %s", prefix.c_str());
0441 auto* g4v = g4pv->GetLogicalVolume();
0442 for(size_t i=0, n=g4v->GetNoDaughters(); i<n; ++i) {
0443 auto* dau = g4v->GetDaughter(i);
0444 printG4(path + dau->GetName(), dau);
0445 }
0446 }
0447
0448 int Geant4DetectorGeometryConstruction::printG4Tree(const char* vol_path) {
0449 if ( vol_path ) {
0450 auto [p, pv] = resolve_path(vol_path);
0451 if ( pv.isValid() ) {
0452 auto& g4map = Geant4Mapping::instance().data().g4Placements;
0453 auto it = g4map.find(pv);
0454 if ( it != g4map.end() ) {
0455 printG4(p, (*it).second);
0456 }
0457 return 1;
0458 }
0459 }
0460 warning("+++ printVolume: Could not access Volume/DetElement '%s'", vol_path ? vol_path : "UNKNOWN");
0461 return 0;
0462 }
0463
0464
0465 void Geant4DetectorGeometryConstruction::installCommandMessenger() {
0466 this->Geant4DetectorConstruction::installCommandMessenger();
0467 m_control->addCall("writeGDML", "GDML: write geometry to file: '"+m_dumpGDML+
0468 "' [uses argument - or - property DumpGDML]",
0469 Callback(this).make(&Geant4DetectorGeometryConstruction::writeGDML),1);
0470 m_control->addCall("printVolume", "Print Geant4 volume properties [uses argument]",
0471 Callback(this).make(&Geant4DetectorGeometryConstruction::printVolume),1);
0472 m_control->addCall("printTree", "Print volume tree WITHOUT properties [uses argument]",
0473 Callback(this).make(&Geant4DetectorGeometryConstruction::printVolTree),1);
0474 m_control->addCall("printG4Tree", "Print Geant4 volume tree [uses argument]",
0475 Callback(this).make(&Geant4DetectorGeometryConstruction::printG4Tree),1);
0476 m_control->addCall("printVolumeTree", "Print volume tree with properties [uses argument]",
0477 Callback(this).make(&Geant4DetectorGeometryConstruction::printVolumeTree),1);
0478 m_control->addCall("checkVolume", "Check Geant4 volume properties [uses argument]",
0479 Callback(this).make(&Geant4DetectorGeometryConstruction::checkVolume),1);
0480 m_control->addCall("printMaterial", "Print Geant4 material properties [uses argument]",
0481 Callback(this).make(&Geant4DetectorGeometryConstruction::printMaterial),1);
0482 }
