Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2026-05-21 07:47:42

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 CERN for the benefit of the ACTS project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at https://mozilla.org/MPL/2.0/.
0008 
0009 #include "ActsExamples/Geant4/Geant4Simulation.hpp"
0010 
0011 #include "Acts/Definitions/Algebra.hpp"
0012 #include "Acts/Utilities/Logger.hpp"
0013 #include "ActsExamples/Framework/AlgorithmContext.hpp"
0014 #include "ActsExamples/Framework/IAlgorithm.hpp"
0015 #include "ActsExamples/Framework/RandomNumbers.hpp"
0016 #include "ActsExamples/Framework/WhiteBoard.hpp"
0017 #include "ActsExamples/Geant4/EventStore.hpp"
0018 #include "ActsExamples/Geant4/Geant4Manager.hpp"
0019 #include "ActsExamples/Geant4/MagneticFieldWrapper.hpp"
0020 #include "ActsExamples/Geant4/MaterialPhysicsList.hpp"
0021 #include "ActsExamples/Geant4/MaterialSteppingAction.hpp"
0022 #include "ActsExamples/Geant4/ParticleKillAction.hpp"
0023 #include "ActsExamples/Geant4/ParticleTrackingAction.hpp"
0024 #include "ActsExamples/Geant4/SensitiveSteppingAction.hpp"
0025 #include "ActsExamples/Geant4/SensitiveSurfaceMapper.hpp"
0026 #include "ActsExamples/Geant4/SimParticleTranslation.hpp"
0027 #include "ActsExamples/Geant4/SteppingActionList.hpp"
0028 #include "ActsPlugins/FpeMonitoring/FpeMonitor.hpp"
0029 
0030 #include <stdexcept>
0031 #include <utility>
0032 
0033 #include <G4FieldManager.hh>
0034 #include <G4RunManager.hh>
0035 #include <G4TransportationManager.hh>
0036 #include <G4UniformMagField.hh>
0037 #include <G4UserEventAction.hh>
0038 #include <G4UserLimits.hh>
0039 #include <G4UserRunAction.hh>
0040 #include <G4UserSteppingAction.hh>
0041 #include <G4UserTrackingAction.hh>
0042 #include <G4VUserDetectorConstruction.hh>
0043 #include <G4VUserPhysicsList.hh>
0044 #include <G4Version.hh>
0045 #include <Randomize.hh>
0046 
0047 namespace ActsExamples {
0048 
0049 Geant4SimulationBase::Geant4SimulationBase(
0050     const Config& cfg, const std::string& name,
0051     std::unique_ptr<const Acts::Logger> logger)
0052     : IAlgorithm(name, std::move(logger)) {
0053   if (cfg.inputParticles.empty()) {
0054     throw std::invalid_argument("Missing input particle collection");
0055   }
0056   if (cfg.detector == nullptr) {
0057     throw std::invalid_argument("Missing detector construction factory");
0058   }
0059   if (cfg.randomNumbers == nullptr) {
0060     throw std::invalid_argument("Missing random numbers");
0061   }
0062 
0063   m_eventStore = std::make_shared<Geant4::EventStore>();
0064 
0065   // tweak logging
0066   // If we are in VERBOSE mode, set the verbose level in Geant4 to 2.
0067   // 3 would be also possible, but that produces infinite amount of output.
0068   m_geant4Level = this->logger().level() == Acts::Logging::VERBOSE ? 2 : 0;
0069 }
0070 
0071 Geant4SimulationBase::~Geant4SimulationBase() = default;
0072 
0073 void Geant4SimulationBase::commonInitialization() {
0074   // Set the detector construction
0075   {
0076     // Clear detector construction if it exists
0077     if (runManager().GetUserDetectorConstruction() != nullptr) {
0078       delete runManager().GetUserDetectorConstruction();
0079     }
0080     // G4RunManager will take care of deletion
0081     m_detectorConstruction =
0082         config()
0083             .detector
0084             ->buildGeant4DetectorConstruction(config().constructionOptions)
0085             .release();
0086     runManager().SetUserInitialization(m_detectorConstruction);
0087     runManager().InitializeGeometry();
0088   }
0089 
0090   m_geant4Instance->tweakLogging(m_geant4Level);
0091 }
0092 
0093 G4RunManager& Geant4SimulationBase::runManager() const {
0094   return *m_geant4Instance->runManager;
0095 }
0096 
0097 Geant4::EventStore& Geant4SimulationBase::eventStore() const {
0098   return *m_eventStore;
0099 }
0100 
0101 ProcessCode Geant4SimulationBase::initialize() {
0102   // Initialize the Geant4 run manager
0103   runManager().Initialize();
0104 
0105   return ProcessCode::SUCCESS;
0106 }
0107 
0108 ProcessCode Geant4SimulationBase::execute(const AlgorithmContext& ctx) const {
0109   // Ensure exclusive access to the Geant4 run manager
0110   std::lock_guard<std::mutex> guard(m_geant4Instance->mutex);
0111 
0112   // Set the seed new per event, so that we get reproducible results
0113   G4Random::setTheSeed(config().randomNumbers->generateSeed(ctx));
0114 
0115   // Get and reset event registry state
0116   eventStore() = Geant4::EventStore{};
0117 
0118   // Register the current event store to the registry
0119   // this will allow access from the User*Actions
0120   eventStore().store = &(ctx.eventStore);
0121 
0122   // Register the input particle read handle
0123   eventStore().inputParticles = &m_inputParticles;
0124 
0125   eventStore().geoContext = ctx.geoContext;
0126 
0127   ACTS_DEBUG("Sending Geant RunManager the BeamOn() command.");
0128   {
0129     ActsPlugins::FpeMonitor mon{0};  // disable all FPEs while we're in Geant4
0130     // Start simulation. each track is simulated as a separate Geant4 event.
0131     runManager().BeamOn(1);
0132   }
0133 
0134   // Print out warnings about possible particle collision if happened
0135   if (eventStore().particleIdCollisionsInitial > 0 ||
0136       eventStore().particleIdCollisionsFinal > 0 ||
0137       eventStore().parentIdNotFound > 0) {
0138     ACTS_WARNING(
0139         "Particle ID collisions detected, don't trust the particle "
0140         "identification!");
0141     ACTS_WARNING(
0142         "- initial states: " << eventStore().particleIdCollisionsInitial);
0143     ACTS_WARNING("- final states: " << eventStore().particleIdCollisionsFinal);
0144     ACTS_WARNING("- parent ID not found: " << eventStore().parentIdNotFound);
0145   }
0146 
0147   if (eventStore().hits.empty()) {
0148     ACTS_DEBUG("Step merging: No steps recorded");
0149   } else {
0150     ACTS_DEBUG("Step merging: mean hits per hit: "
0151                << static_cast<double>(eventStore().numberGeantSteps) /
0152                       eventStore().hits.size());
0153     ACTS_DEBUG(
0154         "Step merging: max hits per hit: " << eventStore().maxStepsForHit);
0155   }
0156 
0157   return ProcessCode::SUCCESS;
0158 }
0159 
0160 std::shared_ptr<Geant4Handle> Geant4SimulationBase::geant4Handle() const {
0161   return m_geant4Instance;
0162 }
0163 
0164 Geant4Simulation::Geant4Simulation(const Config& cfg,
0165                                    std::unique_ptr<const Acts::Logger> logger)
0166     : Geant4SimulationBase(cfg, "Geant4Simulation", std::move(logger)),
0167       m_cfg(cfg) {
0168   m_geant4Instance =
0169       m_cfg.geant4Handle
0170           ? m_cfg.geant4Handle
0171           : Geant4Manager::instance().createHandle(m_cfg.physicsList);
0172   if (m_geant4Instance->physicsListName != m_cfg.physicsList) {
0173     throw std::runtime_error("inconsistent physics list");
0174   }
0175 
0176   commonInitialization();
0177 
0178   // Set the primarty generator
0179   {
0180     // Clear primary generation action if it exists
0181     if (runManager().GetUserPrimaryGeneratorAction() != nullptr) {
0182       delete runManager().GetUserPrimaryGeneratorAction();
0183     }
0184     Geant4::SimParticleTranslation::Config prCfg;
0185     prCfg.eventStore = m_eventStore;
0186     // G4RunManager will take care of deletion
0187     auto primaryGeneratorAction = new Geant4::SimParticleTranslation(
0188         prCfg, this->logger().cloneWithSuffix("SimParticleTranslation"));
0189     // Set the primary generator action
0190     runManager().SetUserAction(primaryGeneratorAction);
0191   }
0192 
0193   // Particle action
0194   {
0195     // Clear tracking action if it exists
0196     if (runManager().GetUserTrackingAction() != nullptr) {
0197       delete runManager().GetUserTrackingAction();
0198     }
0199     Geant4::ParticleTrackingAction::Config trackingCfg;
0200     trackingCfg.eventStore = m_eventStore;
0201     trackingCfg.keepParticlesWithoutHits = cfg.keepParticlesWithoutHits;
0202     // G4RunManager will take care of deletion
0203     auto trackingAction = new Geant4::ParticleTrackingAction(
0204         trackingCfg, this->logger().cloneWithSuffix("ParticleTracking"));
0205     runManager().SetUserAction(trackingAction);
0206   }
0207 
0208   // Stepping actions
0209   Geant4::SensitiveSteppingAction* sensitiveSteppingActionAccess = nullptr;
0210   {
0211     // Clear stepping action if it exists
0212     if (runManager().GetUserSteppingAction() != nullptr) {
0213       delete runManager().GetUserSteppingAction();
0214     }
0215 
0216     Geant4::ParticleKillAction::Config particleKillCfg;
0217     particleKillCfg.eventStore = m_eventStore;
0218     particleKillCfg.volume = cfg.killVolume;
0219     particleKillCfg.maxTime = cfg.killAfterTime;
0220     particleKillCfg.secondaries = cfg.killSecondaries;
0221 
0222     Geant4::SensitiveSteppingAction::Config stepCfg;
0223     stepCfg.eventStore = m_eventStore;
0224     stepCfg.charged = cfg.recordHitsOfCharged;
0225     stepCfg.neutral = cfg.recordHitsOfNeutrals;
0226     stepCfg.primary = cfg.recordHitsOfPrimaries;
0227     stepCfg.secondary = cfg.recordHitsOfSecondaries;
0228     stepCfg.stepLogging = cfg.recordPropagationSummaries;
0229 
0230     Geant4::SteppingActionList::Config steppingCfg;
0231     steppingCfg.actions.push_back(std::make_unique<Geant4::ParticleKillAction>(
0232         particleKillCfg, this->logger().cloneWithSuffix("Killer")));
0233 
0234     auto sensitiveSteppingAction =
0235         std::make_unique<Geant4::SensitiveSteppingAction>(
0236             stepCfg, this->logger().cloneWithSuffix("SensitiveStepping"));
0237     sensitiveSteppingActionAccess = sensitiveSteppingAction.get();
0238 
0239     steppingCfg.actions.push_back(std::move(sensitiveSteppingAction));
0240 
0241     // G4RunManager will take care of deletion
0242     auto steppingAction = new Geant4::SteppingActionList(steppingCfg);
0243     runManager().SetUserAction(steppingAction);
0244   }
0245 
0246   // Get the g4World cache
0247   G4VPhysicalVolume* g4World = m_detectorConstruction->Construct();
0248 
0249   // Please note:
0250   // The following two blocks rely on the fact that the Acts
0251   // detector constructions cache the world volume
0252 
0253   // Set the magnetic field
0254   if (cfg.magneticField) {
0255     ACTS_LOG_WITH_LOGGER(this->logger(), Acts::Logging::INFO,
0256                          "Setting ACTS configured field to Geant4.");
0257 
0258     Geant4::MagneticFieldWrapper::Config g4FieldCfg;
0259     g4FieldCfg.magneticField = cfg.magneticField;
0260     m_magneticField =
0261         std::make_unique<Geant4::MagneticFieldWrapper>(g4FieldCfg);
0262 
0263     // Set the field or the G4Field manager
0264     m_fieldManager = std::make_unique<G4FieldManager>();
0265     m_fieldManager->SetDetectorField(m_magneticField.get());
0266     m_fieldManager->CreateChordFinder(m_magneticField.get());
0267 
0268     // Propagate down to all childrend
0269     g4World->GetLogicalVolume()->SetFieldManager(m_fieldManager.get(), true);
0270   }
0271 
0272   // ACTS sensitive surfaces are provided, so hit creation is turned on
0273   if (cfg.sensitiveSurfaceMapper != nullptr) {
0274     Geant4::SensitiveSurfaceMapper::State sState;
0275     ACTS_LOG_WITH_LOGGER(this->logger(), Acts::Logging::INFO,
0276                          "Remapping selected volumes from Geant4 to "
0277                          "Acts::Surface::GeometryID");
0278     cfg.sensitiveSurfaceMapper->remapSensitiveNames(
0279         sState, Acts::GeometryContext::dangerouslyDefaultConstruct(), g4World,
0280         Acts::Transform3::Identity());
0281 
0282     auto allSurfacesMapped = cfg.sensitiveSurfaceMapper->checkMapping(
0283         sState, Acts::GeometryContext::dangerouslyDefaultConstruct(), false,
0284         false);
0285     if (!allSurfacesMapped) {
0286       ACTS_LOG_WITH_LOGGER(this->logger(), Acts::Logging::WARNING,
0287                            "Not all sensitive surfaces have been mapped to "
0288                            "Geant4 volumes!");
0289     }
0290 
0291     sensitiveSteppingActionAccess->assignSurfaceMapping(
0292         sState.g4VolumeToSurfaces);
0293   }
0294 
0295   m_inputParticles.initialize(cfg.inputParticles);
0296   m_outputSimHits.initialize(cfg.outputSimHits);
0297   m_outputParticles.initialize(cfg.outputParticles);
0298 
0299   if (cfg.recordPropagationSummaries) {
0300     m_outputPropagationSummaries.initialize(cfg.outputPropagationSummaries);
0301   }
0302 }
0303 
0304 Geant4Simulation::~Geant4Simulation() = default;
0305 
0306 ProcessCode Geant4Simulation::execute(const AlgorithmContext& ctx) const {
0307   auto ret = Geant4SimulationBase::execute(ctx);
0308   if (ret != ProcessCode::SUCCESS) {
0309     return ret;
0310   }
0311 
0312   // Output handling: Simulation
0313   m_outputParticles(
0314       ctx, SimParticleContainer(eventStore().particlesSimulated.begin(),
0315                                 eventStore().particlesSimulated.end()));
0316 
0317   m_outputSimHits(
0318       ctx, SimHitContainer(eventStore().hits.begin(), eventStore().hits.end()));
0319 
0320   // Output the propagation summaries if requested
0321   if (m_cfg.recordPropagationSummaries) {
0322     PropagationSummaries summaries;
0323     summaries.reserve(eventStore().propagationRecords.size());
0324     for (auto& [trackId, summary] : eventStore().propagationRecords) {
0325       summaries.push_back(std::move(summary));
0326     }
0327     m_outputPropagationSummaries(ctx, std::move(summaries));
0328   }
0329 
0330   return ProcessCode::SUCCESS;
0331 }
0332 
0333 Geant4MaterialRecording::Geant4MaterialRecording(
0334     const Config& cfg, std::unique_ptr<const Acts::Logger> logger)
0335     : Geant4SimulationBase(cfg, "Geant4Simulation", std::move(logger)),
0336       m_cfg(cfg) {
0337   auto physicsListName = "MaterialPhysicsList";
0338   m_geant4Instance =
0339       m_cfg.geant4Handle
0340           ? m_cfg.geant4Handle
0341           : Geant4Manager::instance().createHandle(
0342                 std::make_unique<Geant4::MaterialPhysicsList>(
0343                     this->logger().cloneWithSuffix("MaterialPhysicsList")),
0344                 physicsListName);
0345   if (m_geant4Instance->physicsListName != physicsListName) {
0346     throw std::runtime_error("inconsistent physics list");
0347   }
0348 
0349   commonInitialization();
0350 
0351   // Set the primarty generator
0352   {
0353     // Clear primary generation action if it exists
0354     if (runManager().GetUserPrimaryGeneratorAction() != nullptr) {
0355       delete runManager().GetUserPrimaryGeneratorAction();
0356     }
0357 
0358     Geant4::SimParticleTranslation::Config prCfg;
0359     prCfg.eventStore = m_eventStore;
0360     prCfg.forcedPdgCode = 0;
0361     prCfg.forcedCharge = 0.;
0362     prCfg.forcedMass = 0.;
0363 
0364     // G4RunManager will take care of deletion
0365     auto primaryGeneratorAction = new Geant4::SimParticleTranslation(
0366         prCfg, this->logger().cloneWithSuffix("SimParticleTranslation"));
0367     // Set the primary generator action
0368     runManager().SetUserAction(primaryGeneratorAction);
0369   }
0370 
0371   // Particle action
0372   {
0373     // Clear tracking action if it exists
0374     if (runManager().GetUserTrackingAction() != nullptr) {
0375       delete runManager().GetUserTrackingAction();
0376     }
0377     Geant4::ParticleTrackingAction::Config trackingCfg;
0378     trackingCfg.eventStore = m_eventStore;
0379     trackingCfg.keepParticlesWithoutHits = true;
0380     // G4RunManager will take care of deletion
0381     auto trackingAction = new Geant4::ParticleTrackingAction(
0382         trackingCfg, this->logger().cloneWithSuffix("ParticleTracking"));
0383     runManager().SetUserAction(trackingAction);
0384   }
0385 
0386   // Stepping action
0387   {
0388     // Clear stepping action if it exists
0389     if (runManager().GetUserSteppingAction() != nullptr) {
0390       delete runManager().GetUserSteppingAction();
0391     }
0392     Geant4::MaterialSteppingAction::Config steppingCfg;
0393     steppingCfg.eventStore = m_eventStore;
0394     steppingCfg.excludeMaterials = m_cfg.excludeMaterials;
0395     steppingCfg.recordElementFractions = m_cfg.recordElementFractions;
0396     // G4RunManager will take care of deletion
0397     auto steppingAction = new Geant4::MaterialSteppingAction(
0398         steppingCfg, this->logger().cloneWithSuffix("MaterialSteppingAction"));
0399     runManager().SetUserAction(steppingAction);
0400   }
0401 
0402   runManager().Initialize();
0403 
0404   m_inputParticles.initialize(cfg.inputParticles);
0405   m_outputMaterialTracks.initialize(cfg.outputMaterialTracks);
0406 }
0407 
0408 Geant4MaterialRecording::~Geant4MaterialRecording() = default;
0409 
0410 ProcessCode Geant4MaterialRecording::execute(
0411     const AlgorithmContext& ctx) const {
0412   const auto ret = Geant4SimulationBase::execute(ctx);
0413   if (ret != ProcessCode::SUCCESS) {
0414     return ret;
0415   }
0416 
0417   // Output handling: Material tracks
0418   m_outputMaterialTracks(
0419       ctx, decltype(eventStore().materialTracks)(eventStore().materialTracks));
0420 
0421   return ProcessCode::SUCCESS;
0422 }
0423 
0424 }  // namespace ActsExamples