Skip to content

Commit 26ff503

Browse files
MarkusFrankATcernchMarkusFrankATcernch
committedMar 3, 2023
Propagation of World settings from dd4hep to Geant4
1 parent 6166ded commit 26ff503

File tree

4 files changed

+151
-101
lines changed

4 files changed

+151
-101
lines changed
 

‎DDG4/include/DDG4/Geant4Converter.h

+2
Original file line numberDiff line numberDiff line change
@@ -45,6 +45,8 @@ namespace dd4hep {
4545
bool debugReflections = false;
4646
/// Property: Flag to debug regions during conversion mechanism
4747
bool debugRegions = false;
48+
/// Property: Flag to debug LimitSets during conversion mechanism
49+
bool debugLimits = false;
4850
/// Property: Flag to debug surfaces during conversion mechanism
4951
bool debugSurfaces = false;
5052

‎DDG4/include/DDG4/Geant4UserLimits.h

+3-1
Original file line numberDiff line numberDiff line change
@@ -72,9 +72,11 @@ namespace dd4hep {
7272

7373
public:
7474
/// Initializing Constructor
75-
Geant4UserLimits(LimitSet ls);
75+
Geant4UserLimits(LimitSet limitset);
7676
/// Standard destructor
7777
virtual ~Geant4UserLimits();
78+
/// Update the object
79+
virtual void update(LimitSet limitset);
7880
/// Access the user tracklength for a G4 track object
7981
virtual G4double GetMaxAllowedStep(const G4Track& track)
8082
{ return maxStepLength.value(track); }

‎DDG4/src/Geant4Converter.cpp

+124-90
Original file line numberDiff line numberDiff line change
@@ -55,10 +55,11 @@
5555
#include <G4Isotope.hh>
5656
#include <G4Material.hh>
5757
#include <G4UserLimits.hh>
58+
#include <G4RegionStore.hh>
5859
#include <G4FieldManager.hh>
5960
#include <G4LogicalVolume.hh>
60-
#include <G4ReflectionFactory.hh>
6161
#include <G4OpticalSurface.hh>
62+
#include <G4ReflectionFactory.hh>
6263
#include <G4LogicalSkinSurface.hh>
6364
#include <G4ElectroMagneticField.hh>
6465
#include <G4LogicalBorderSurface.hh>
@@ -73,6 +74,7 @@
7374
#include <iostream>
7475
#include <iomanip>
7576
#include <sstream>
77+
#include <limits>
7678

7779
namespace units = dd4hep;
7880
using namespace dd4hep::detail;
@@ -711,70 +713,94 @@ void* Geant4Converter::handleVolume(const string& name, const TGeoVolume* volume
711713
return nullptr;
712714
}
713715
else if (volIt == info.g4Volumes.end() ) {
714-
string n = volume->GetName();
715-
TGeoMedium* med = volume->GetMedium();
716-
TGeoShape* sh = volume->GetShape();
717-
G4VSolid* solid = (G4VSolid*) handleSolid(sh->GetName(), sh);
718-
bool is_assembly = sh->IsA() == TGeoShapeAssembly::Class() || volume->IsA() == TGeoVolumeAssembly::Class();
716+
const char* vnam = volume->GetName();
717+
TGeoMedium* med = volume->GetMedium();
718+
Solid sh = volume->GetShape();
719+
bool is_assembly = sh->IsA() == TGeoShapeAssembly::Class() || volume->IsAssembly();
719720

720721
printout(lvl, "Geant4Converter", "++ Convert Volume %-32s: %p %s/%s assembly:%s",
721-
n.c_str(), volume, sh->IsA()->GetName(), volume->IsA()->GetName(), yes_no(is_assembly));
722-
722+
vnam, volume, sh.type(), _v.type(), yes_no(is_assembly));
723723
if ( is_assembly ) {
724724
return nullptr;
725725
}
726-
Region reg = _v.region();
727-
LimitSet lim = _v.limitSet();
728-
VisAttr vis = _v.visAttributes();
729-
G4Region* region = reg.isValid() ? info.g4Regions[reg] : nullptr;
730-
G4UserLimits* limits = lim.isValid() ? info.g4Limits[lim] : nullptr;
731-
G4Material* medium = (G4Material*) handleMaterial(med->GetName(), Material(med));
726+
Region reg = _v.region();
727+
LimitSet lim = _v.limitSet();
728+
VisAttr vis = _v.visAttributes();
729+
G4Region* g4region = reg.isValid() ? info.g4Regions[reg] : nullptr;
730+
G4UserLimits* g4limits = lim.isValid() ? info.g4Limits[lim] : nullptr;
731+
G4VSolid* g4solid = (G4VSolid*) handleSolid(sh->GetName(), sh);
732+
G4Material* g4medium = (G4Material*) handleMaterial(med->GetName(), Material(med));
732733
/// Check all pre-conditions
733-
if ( !solid ) {
734-
except("G4Converter","++ No Geant4 Solid present for volume:" + n);
734+
if ( !g4solid ) {
735+
except("G4Converter","++ No Geant4 Solid present for volume: %s", vnam);
735736
}
736-
else if ( !medium ) {
737-
except("G4Converter","++ No Geant4 material present for volume:" + n);
737+
else if ( !g4medium ) {
738+
except("G4Converter","++ No Geant4 material present for volume: %s", vnam);
738739
}
739-
else if ( reg.isValid() && !region ) {
740-
except("G4Cnv::volume["+name+"]"," ++ Failed to access Geant4 region %s.",reg.name());
740+
else if ( reg.isValid() && !g4region ) {
741+
except("G4Cnv::volume["+name+"]"," ++ Failed to access Geant4 region %s.", reg.name());
741742
}
742-
else if ( lim.isValid() && !limits ) {
743-
except("G4Cnv::volume["+name+"]","++ FATAL Failed to access Geant4 user limits %s.",lim.name());
743+
else if ( lim.isValid() && !g4limits ) {
744+
except("G4Cnv::volume["+name+"]","++ FATAL Failed to access Geant4 user limits %s.", lim.name());
744745
}
745-
else if ( limits ) {
746-
printout(lvl, "Geant4Converter", "++ Volume + Apply LIMITS settings:%-24s to volume %s.",
747-
lim.name(), _v.name());
746+
else if ( g4limits ) {
747+
printout(lvl, "Geant4Converter", "++ Volume + Apply LIMITS settings: %-24s to volume %s.",
748+
lim.name(), vnam);
748749
}
749750

750751
G4LogicalVolume* g4vol = nullptr;
751752
if ( _v.hasProperties() && !_v.getProperty(GEANT4_TAG_PLUGIN,"").empty() ) {
752753
Detector* det = const_cast<Detector*>(&m_detDesc);
753754
string plugin = _v.getProperty(GEANT4_TAG_PLUGIN,"");
754-
g4vol = PluginService::Create<G4LogicalVolume*>(plugin, det, _v, solid, medium);
755+
g4vol = PluginService::Create<G4LogicalVolume*>(plugin, det, _v, g4solid, g4medium);
755756
if ( !g4vol ) {
756757
except("G4Cnv::volume["+name+"]","++ FATAL Failed to call plugin to create logical volume.");
757758
}
758759
}
759760
else {
760-
g4vol = new G4LogicalVolume(solid, medium, n, nullptr, nullptr, nullptr);
761+
g4vol = new G4LogicalVolume(g4solid, g4medium, vnam, nullptr, nullptr, nullptr);
761762
}
762763

763-
if ( limits ) {
764-
g4vol->SetUserLimits(limits);
764+
if ( g4limits ) {
765+
g4vol->SetUserLimits(g4limits);
765766
}
766-
if ( region ) {
767-
printout(lvl, "Geant4Converter", "++ Volume + Apply REGION settings: %s to volume %s.",
768-
reg.name(), _v.name());
769-
g4vol->SetRegion(region);
770-
region->AddRootLogicalVolume(g4vol);
767+
if ( g4region ) {
768+
PrintLevel plevel = (debugVolumes||debugRegions||debugLimits) ? ALWAYS : outputLevel;
769+
printout(plevel, "Geant4Converter", "++ Volume + Apply REGION settings: %-24s to volume %s.",
770+
reg.name(), vnam);
771+
// Handle the region settings for the world volume seperately.
772+
// Geant4 does NOT WANT any regions assigned to the workd volume.
773+
// The world's region is created in the G4RunManagerKernel!
774+
if ( _v == m_detDesc.worldVolume() ) {
775+
const char* wrd_nam = "DefaultRegionForTheWorld";
776+
const char* src_nam = g4region->GetName().c_str();
777+
auto* world_region = G4RegionStore::GetInstance()->GetRegion(wrd_nam, false);
778+
if ( auto* cuts = g4region->GetProductionCuts() ) {
779+
world_region->SetProductionCuts(cuts);
780+
printout(plevel, "Geant4Converter",
781+
"++ Volume %s Region: %s. Apply production cuts from %s",
782+
vnam, wrd_nam, src_nam);
783+
}
784+
if ( auto* lims = g4region->GetUserLimits() ) {
785+
world_region->SetUserLimits(lims);
786+
printout(plevel, "Geant4Converter",
787+
"++ Volume %s Region: %s. Apply user limits from %s",
788+
vnam, wrd_nam, src_nam);
789+
}
790+
}
791+
else {
792+
g4vol->SetRegion(g4region);
793+
g4region->AddRootLogicalVolume(g4vol);
794+
}
771795
}
772-
G4VisAttributes* vattr = vis.isValid() ? (G4VisAttributes*)handleVis(vis.name(), vis) : nullptr;
773-
if ( vattr ) {
774-
g4vol->SetVisAttributes(vattr);
796+
G4VisAttributes* g4vattr = vis.isValid()
797+
? (G4VisAttributes*)handleVis(vis.name(), vis) : nullptr;
798+
if ( g4vattr ) {
799+
g4vol->SetVisAttributes(g4vattr);
775800
}
776801
info.g4Volumes[volume] = g4vol;
777-
printout(lvl, "Geant4Converter", "++ Volume + %s converted: %p ---> G4: %p", n.c_str(), volume, g4vol);
802+
printout(lvl, "Geant4Converter",
803+
"++ Volume + %s converted: %p ---> G4: %p", vnam, volume, g4vol);
778804
}
779805
return nullptr;
780806
}
@@ -786,13 +812,16 @@ void* Geant4Converter::collectVolume(const string& /* name */, const TGeoVolume*
786812
Region reg = _v.region();
787813
LimitSet lim = _v.limitSet();
788814
SensitiveDetector det = _v.sensitiveDetector();
789-
790-
if ( lim.isValid() )
791-
info.limits[lim].insert(volume);
792-
if ( reg.isValid() )
793-
info.regions[reg].insert(volume);
794-
if ( det.isValid() )
795-
info.sensitives[det].insert(volume);
815+
bool world = (volume == m_detDesc.worldVolume().ptr());
816+
817+
if ( !world ) {
818+
if ( lim.isValid() )
819+
info.limits[lim].insert(volume);
820+
if ( reg.isValid() )
821+
info.regions[reg].insert(volume);
822+
if ( det.isValid() )
823+
info.sensitives[det].insert(volume);
824+
}
796825
return (void*)volume;
797826
}
798827

@@ -1043,8 +1072,8 @@ void* Geant4Converter::handleRegion(Region region, const set<const TGeoVolume*>&
10431072
G4Region* g4 = data().g4Regions[region];
10441073
if ( !g4 ) {
10451074
PrintLevel lvl = debugRegions ? ALWAYS : outputLevel;
1046-
Region r = region;
1047-
g4 = new G4Region(r.name());
1075+
Region r = region;
1076+
g4 = new G4Region(region.name());
10481077

10491078
// create region info with storeSecondaries flag
10501079
if( not r.wasThresholdSet() and r.storeSecondaries() ) {
@@ -1065,49 +1094,49 @@ void* Geant4Converter::handleRegion(Region region, const set<const TGeoVolume*>&
10651094
cuts = new G4ProductionCuts();
10661095
cuts->SetProductionCut(r.cut()*CLHEP::mm/units::mm);
10671096
printout(lvl, "Geant4Converter", "++ %s: Using default cut: %f [mm]",
1068-
r.name(), r.cut()*CLHEP::mm/units::mm);
1097+
r.name(), r.cut()*CLHEP::mm/units::mm);
10691098
}
10701099
for( const auto& nam : limits ) {
10711100
LimitSet ls = m_detDesc.limitSet(nam);
1072-
if (ls.isValid()) {
1073-
const LimitSet::Set& cts = ls.cuts();
1074-
for (const auto& c : cts ) {
1075-
int pid = 0;
1076-
if ( c.particles == "*" ) pid = -1;
1077-
else if ( c.particles == "e-" ) pid = idxG4ElectronCut;
1078-
else if ( c.particles == "e+" ) pid = idxG4PositronCut;
1079-
else if ( c.particles == "e[+-]" ) pid = -idxG4PositronCut-idxG4ElectronCut;
1080-
else if ( c.particles == "e[-+]" ) pid = -idxG4PositronCut-idxG4ElectronCut;
1081-
else if ( c.particles == "gamma" ) pid = idxG4GammaCut;
1082-
else if ( c.particles == "proton" ) pid = idxG4ProtonCut;
1083-
else throw runtime_error("G4Region: Invalid production cut particle-type:" + c.particles);
1084-
if ( !cuts ) cuts = new G4ProductionCuts();
1085-
if ( pid == -(idxG4PositronCut+idxG4ElectronCut) ) {
1086-
cuts->SetProductionCut(c.value*CLHEP::mm/units::mm, idxG4PositronCut);
1087-
cuts->SetProductionCut(c.value*CLHEP::mm/units::mm, idxG4ElectronCut);
1088-
}
1089-
else {
1090-
cuts->SetProductionCut(c.value*CLHEP::mm/units::mm, pid);
1091-
}
1092-
printout(lvl, "Geant4Converter", "++ %s: Set cut [%s/%d] = %f [mm]",
1093-
r.name(), c.particles.c_str(), pid, c.value*CLHEP::mm/units::mm);
1094-
}
1095-
bool found = false;
1096-
const auto& lm = data().g4Limits;
1097-
for (const auto& j : lm ) {
1098-
if (nam == j.first->GetName()) {
1099-
g4->SetUserLimits(j.second);
1100-
printout(lvl, "Geant4Converter", "++ %s: Set limits %s to region type %s",
1101-
r.name(), nam.c_str(), j.second->GetType().c_str());
1102-
found = true;
1103-
break;
1104-
}
1105-
}
1106-
if ( found ) {
1107-
continue;
1108-
}
1101+
if ( ls.isValid() ) {
1102+
const LimitSet::Set& cts = ls.cuts();
1103+
for (const auto& c : cts ) {
1104+
int pid = 0;
1105+
if ( c.particles == "*" ) pid = -1;
1106+
else if ( c.particles == "e-" ) pid = idxG4ElectronCut;
1107+
else if ( c.particles == "e+" ) pid = idxG4PositronCut;
1108+
else if ( c.particles == "e[+-]" ) pid = -idxG4PositronCut-idxG4ElectronCut;
1109+
else if ( c.particles == "e[-+]" ) pid = -idxG4PositronCut-idxG4ElectronCut;
1110+
else if ( c.particles == "gamma" ) pid = idxG4GammaCut;
1111+
else if ( c.particles == "proton" ) pid = idxG4ProtonCut;
1112+
else throw runtime_error("G4Region: Invalid production cut particle-type:" + c.particles);
1113+
if ( !cuts ) cuts = new G4ProductionCuts();
1114+
if ( pid == -(idxG4PositronCut+idxG4ElectronCut) ) {
1115+
cuts->SetProductionCut(c.value*CLHEP::mm/units::mm, idxG4PositronCut);
1116+
cuts->SetProductionCut(c.value*CLHEP::mm/units::mm, idxG4ElectronCut);
1117+
}
1118+
else {
1119+
cuts->SetProductionCut(c.value*CLHEP::mm/units::mm, pid);
1120+
}
1121+
printout(lvl, "Geant4Converter", "++ %s: Set cut [%s/%d] = %f [mm]",
1122+
r.name(), c.particles.c_str(), pid, c.value*CLHEP::mm/units::mm);
1123+
}
1124+
bool found = false;
1125+
const auto& lm = data().g4Limits;
1126+
for (const auto& j : lm ) {
1127+
if (nam == j.first->GetName()) {
1128+
g4->SetUserLimits(j.second);
1129+
printout(lvl, "Geant4Converter", "++ %s: Set limits %s to region type %s",
1130+
r.name(), nam.c_str(), j.second->GetType().c_str());
1131+
found = true;
1132+
break;
1133+
}
1134+
}
1135+
if ( found ) {
1136+
continue;
1137+
}
11091138
}
1110-
throw runtime_error("G4Region: Failed to resolve user limitset:" + nam);
1139+
except("Geant4Converter", "++ G4Region: Failed to resolve limitset: " + nam);
11111140
}
11121141
/// Assign cuts to region if they were created
11131142
if ( cuts ) g4->SetProductionCuts(cuts);
@@ -1120,27 +1149,32 @@ void* Geant4Converter::handleRegion(Region region, const set<const TGeoVolume*>&
11201149
void* Geant4Converter::handleLimitSet(LimitSet limitset, const set<const TGeoVolume*>& /* volumes */) const {
11211150
G4UserLimits* g4 = data().g4Limits[limitset];
11221151
if ( !g4 ) {
1152+
PrintLevel lvl = debugLimits || debugRegions ? ALWAYS : outputLevel;
11231153
struct LimitPrint {
11241154
const LimitSet& ls;
11251155
LimitPrint(const LimitSet& lset) : ls(lset) {}
11261156
const LimitPrint& operator()(const std::string& pref, const Geant4UserLimits::Handler& h) const {
11271157
if ( !h.particleLimits.empty() ) {
11281158
printout(ALWAYS,"Geant4Converter",
1129-
"+++ LimitSet: Limit %s.%s applied for particles:",ls.name(), pref.c_str());
1159+
"+++ LimitSet: Explicit Limit %s.%s applied for particles:",ls.name(), pref.c_str());
11301160
for(const auto& p : h.particleLimits)
11311161
printout(ALWAYS,"Geant4Converter","+++ LimitSet: Particle type: %-18s PDG: %-6d : %f",
11321162
p.first->GetParticleName().c_str(), p.first->GetPDGEncoding(), p.second);
11331163
}
1134-
else {
1164+
else if ( h.defaultValue > std::numeric_limits<double>::epsilon() ) {
11351165
printout(ALWAYS,"Geant4Converter",
1136-
"+++ LimitSet: Limit %s.%s NOT APPLIED.",ls.name(), pref.c_str());
1166+
"+++ LimitSet: Implicit Limit %s.%s for wildcard particles: %f",
1167+
ls.name(), pref.c_str(), float(h.defaultValue));
11371168
}
11381169
return *this;
11391170
}
11401171
};
11411172
Geant4UserLimits* limits = new Geant4UserLimits(limitset);
11421173
g4 = limits;
1143-
if ( debugRegions ) {
1174+
printout(lvl, "Geant4Converter",
1175+
"++ Successfully converted LimitSet: %s [%ld cuts, %ld limits]",
1176+
limitset.name(), limitset.cuts().size(), limitset.limits().size());
1177+
if ( debugRegions || debugLimits ) {
11441178
LimitPrint print(limitset);
11451179
print("maxTime", limits->maxTime)
11461180
("minEKine", limits->minEKine)

0 commit comments

Comments
 (0)