|
aegir
Phlex-based simulation framework for the SHiP experiment.
|
aegir uses a direct worker integration pattern where Phlex framework threads act as Geant4 worker threads directly, with no separate G4 thread pool or event queue. This is inspired by the CMSSW OscarMTProducer pattern used in CMS.
We evaluated five integration strategies in the ship-phlex-sim benchmarking repository:
G4RunManager::BeamOn(1) on a dedicated thread.std::future for synchronisation.G4TaskRunManager with per-worker rendezvous slots.G4MTRunManager on a dedicated thread; G4 workers pop from a concurrent queue.BeamOn.Strategies 1–4 all involve cross-thread handoff overhead (queues, promises, condition variables). The direct worker approach eliminates this entirely: the thread that receives the event from the framework is the same thread that processes it through Geant4.
In benchmarks with particle gun events on a 4-core machine, the direct pattern matched the best queue-based strategy (sync MT slots) while being simpler to configure and reason about. It also removes the num_events configuration footgun present in the queue-based MT module, where the configured event count had to exactly match the driver's total.
A dedicated std::thread creates a G4MTRunManager and runs the full Geant4 initialisation sequence (geometry, physics, run initialisation). It then stores the world physical volume, physics list, and detector construction pointers for workers to use. The master thread blocks on a shutdown future until the module is destroyed.
Each Phlex thread lazily initialises on its first call to simulate():
G4Threading::G4SetThreadId()G4WorkerThread::BuildGeometryAndPhysicsVector()G4WorkerRunManagerKernel and defines the world volumeInstead of G4RunManager::BeamOn(), the module:
G4Event directly from the input MCParticle vector — keeping only final-state entries, see Primary selection — creating G4PrimaryVertex and G4PrimaryParticle objects (no G4VUserPrimaryGeneratorAction involved)G4State_GeomClosedG4EventManager::ProcessOneEvent(event)This bypasses the G4 run loop entirely, giving the framework full control over event scheduling.
SHiP::MCParticle::status carries the HepMC status code — the same convention Pythia8's Particle::statusHepMC() produces and EDM4hep uses for its generatorStatus:
| Code | Meaning | Tracked? |
|---|---|---|
| 0 | Empty entry, no meaningful information | no |
| 1 | Final state — not decayed by the generator; may still be unstable | yes |
| 2 | Decayed Standard Model hadron, tau or muon | no |
| 3 | Documentation entry | no |
| 4 | Incoming beam particle | no |
| 11-200 | Intermediate entry, generator-dependent classification | no |
Only status-1 entries become primaries. Note that "final state" does not mean stable: a final-state K_S is status 1 and is expected to decay in the detector — that is exactly the work being handed over. Everything else is either the generator's own bookkeeping or a particle the generator has already decayed, and tracking it alongside the daughters it already produced would double-count the event.
The predicate lives in src/hepmc_status.hpp, next to the code set. Geant4's own HepMC example uses the same status == 1 test.
Skipped primaries are counted and reported once per event in an aggregated warning; they remain in the output mc_particles product, so the record is not lost.
This filter is a no-op for a generator running in the default ‘record: 'final_state’mode, where every emitted entry is status 1. It matters underrecord: 'full'` — see generator_record.md.
Without it the failure mode is silent and partial rather than loud: Geant4 declines to track short-lived definitions (quarks, gluons, diquarks, strings) on its own, but happily tracks decayed hadrons and beam particles. Setting track_all_primaries on a full record demonstrates this — the primary count rises to exactly the number of status 1, 2 and 4 entries.
The G4MTRunManager is intentionally leaked at shutdown. Its destructor accesses global singletons that may already be torn down during plugin unloading, causing crashes. This is a known Geant4 limitation.
| Parameter | Type | Default | Description |
|---|---|---|---|
physics_list | string | FTFP_BERT | Geant4 reference physics list name |
concurrency | int | 1 | Number of concurrent worker threads |
verbosity | int | 0 | Geant4 run verbosity level |
sd_mode | string | scoring | Sensitive detector mode: scoring or crossing |
ke_threshold | double | 0.0 | Kinetic energy threshold for crossing SD (GeV) |
energy_cut | bool | false | Enable stepping energy cut |
energy_cut_threshold | double | ke_threshold | KE below which tracks are killed (GeV) |
particle_ke_cut | double | 0.0 | KE below which secondary particles are not recorded (GeV) |
track_all_primaries | bool | false | Hand every input MCParticle to Geant4 regardless of its generator status. By default only HepMC status-1 (final-state) entries are tracked; set this only for inputs whose status does not follow the HepMC convention |
regions | map | {} | Volume name pattern to production cut (mm) mapping |
scoring_planes | list | [] | Scoring planes to add to the geometry, see Scoring planes |
export_gdml | string | *(unset)* | Write the constructed geometry to this GDML file after initialisation. Errors if the file exists. Lets external tools (e.g. the GENIE event generator) use exactly the geometry Geant4 tracks in |
progress_interval | int | 100 | Log a progress line (event count and average rate) every this many simulated events; 0 disables |
For consumers that read the exported file with ROOT's TGeo importer, run scripts/gdml_fix_element_names.py on it first: ROOT confuses materials and elements that share a name after pointer-suffix stripping (Geant4's NIST materials routinely do — material "Iron" made of element "Iron") and silently imports them as empty mixtures. The script renames the colliding elements; the physics is unchanged. scripts/gdml_target_nuclei.py then lists the nuclide PDG codes of the geometry — e.g. the target list neutrino cross-section splines must cover.
Example workflow configuration:
A scoring plane is a thin box or cylindrical shell that records the particles crossing it. It is aegir's counterpart of FairShip's exitHadronAbsorber and is meant for flux studies: the muons leaving the hadron absorber, the neutrinos coming out of the target, and so on.
The planes live in a Geant4 parallel world (aegir_scoring_planes), overlaid on whatever geometry the provider built — builtin, GDML or GeoModel alike. Geant4 tracks every particle through both worlds at once (G4ParallelWorldPhysics, added to the physics list whenever planes are configured). A plane can therefore be any size, anywhere: it can cut through or extend beyond the detector's volumes — a muon-shield container, a magnet, the target — and the parallel world carries no material of its own, so the detector's materials apply everywhere. Planes must not overlap one another (touching is fine) or extend beyond the world volume.
A plane records one SimHit for every track that enters it through its surface, whatever the particle — neutrinos and other neutrals included — and whether or not it deposits energy. It does this through its own sensitive detector, so it works the same with ‘sd_mode: 'scoring’and sd_mode: 'crossing', alongside the geometry's own sensitive volumes, and thesensitive_volumesand fieldvolume_patternpatterns never match the planes. A track entering a plane twice (a cylinder crossed in and out, a particle scattered back) is recorded each time; setkill` to stop a track at the first plane it reaches.
Planes do change how the detector geometry is stepped through. With a parallel world registered, Geant4 uses coupled transportation, and steps in the detector geometry also end at every plane face. Hence:
writes one hit per step, so a plane cutting through a sensitive volume splits some of its hits in two; the total energy deposit is unchanged. -sd_mode: 'crossing'` is unaffected, since its first-step-in-volume flag only follows the detector's own boundaries.The hits are told apart from detector hits by their ids:
SimHit field | Scoring-plane value |
|---|---|
detectorId | 100 (SHiP::g4::kScoringPlaneDetectorId), clear of the SHiP subsystem ids and of the sensitive_volumes indices (which limits a geometry to 100 sensitive_volumes patterns) |
geometryNodeId | the plane's id |
position, momentum, time | at the entry point |
energyDeposit | 0 |
pathLength | track length from the track's origin to the entry point. A ‘sd_mode: 'crossing’` hit gives it at the end of its first step in the volume instead, so the two don't mean quite the same thing |
In the output module's validation histograms, plane hits fill their own h_plane_hit_* histograms (hits per event, plane id, momentum, PDG code) and are left out of the h_hit_* ones, which keep describing the detector's hits: with detector id 100 and no energy deposit, plane hits would otherwise land in the overflow of h_hit_detector and pile up at zero in h_hit_edep. h_plane_hit_id covers plane ids 0–99; other ids land in its underflow or overflow bin.
Each entry of scoring_planes takes:
| Key | Type | Default | Description |
|---|---|---|---|
name | string | — | Volume name; must not already exist in the geometry |
id | int | — | Written to SimHit::geometryNodeId, telling the planes' hits apart; required and unique |
shape | string | box | box or tube (cylindrical shell, axis along z) |
half_x, half_y | double | — | Box half-widths (mm); required for a box, rejected for a tube |
r_min, r_max | double | — | Tube inner and outer radius (mm); required for a tube, rejected for a box |
half_z | double | — | Half-length along z (mm); required for both shapes |
position | [x, y, z] | — | Centre in the global frame (mm) |
after_volume | string | — | Place the plane downstream of this detector volume, centred transversely on its global bounding box |
gap | double | 0 | With after_volume: distance from the volume's downstream end to the plane's upstream face (mm) |
ke_threshold | double | 0 | Record only tracks entering with at least this kinetic energy (GeV). With kill, tracks below it are stopped too, unrecorded |
kill | bool | false | Stop and kill every track that steps into the plane, as FairShip's exitHadronAbsorber does: recorded or not, including secondaries born inside the plane |
Exactly one of position and after_volume must be given. after_volume is an exact logical-volume name, and the volume must be placed exactly once in the geometry — otherwise the plane's position would be ambiguous. Planes are always aligned with the global axes. Unknown keys and inconsistent settings fail the job at configuration time. A plane that overlaps another plane or extends beyond the world fails it at initialisation, before any plane is built. That test is exact for these axis-aligned shapes, so even a small corner overlap between thin planes is caught.
workflows/gun_st_scoring_planes.jsonnet places planes of each kind in the builtin geometry, one of them cutting through the tungsten target, and scripts/scoring_planes.sh (ctest scoring_planes) checks the recorded hits against the geometry for neutrinos and muons.