SimG4Core/Application/src/RunManagerMTWorker.cc - master - githubmirror/cmssw

SimG4Core/Application/src/RunManagerMTWorker.cc
783 строки30 KB
Vladimir
code format
23 апр 2026, 16:04
23 апр 2026, 16:044a762ed
#include "SimG4Core/Application/interface/RunManagerMTWorker.h"
#include "SimG4Core/Application/interface/RunManagerMT.h"
#include "SimG4Core/Application/interface/SimRunInterface.h"
#include "SimG4Core/Application/interface/RunAction.h"
#include "SimG4Core/Application/interface/EventAction.h"
#include "SimG4Core/Application/interface/Phase2EventAction.h"
#include "SimG4Core/Application/interface/StackingAction.h"
#include "SimG4Core/Application/interface/TrackingAction.h"
#include "SimG4Core/Application/interface/Phase2TrackingAction.h"
#include "SimG4Core/Application/interface/SteppingAction.h"
#include "SimG4Core/Application/interface/Phase2SteppingAction.h"
#include "SimG4Core/Application/interface/Phase2StackingAction.h"
#include "SimG4Core/Application/interface/CMSSimEventManager.h"
#include "SimG4Core/Application/interface/CustomUIsessionThreadPrefix.h"
#include "SimG4Core/Application/interface/CustomUIsessionToFile.h"
#include "SimG4Core/Application/interface/ExceptionHandler.h"
#include "SimG4Core/Application/interface/CMSGDMLWriteStructure.h"

#include "SimG4Core/Physics/interface/CMSG4TrackInterface.h"
#include "SimG4Core/Geometry/interface/CustomUIsession.h"

#include "FWCore/ParameterSet/interface/ParameterSet.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"
#include "FWCore/Framework/interface/Event.h"
#include "FWCore/Framework/interface/EventSetup.h"
#include "FWCore/Framework/interface/ConsumesCollector.h"
#include "FWCore/ServiceRegistry/interface/Service.h"
#include "FWCore/Utilities/interface/Exception.h"
#include "FWCore/Utilities/interface/thread_safety_macros.h"

#include "SimG4Core/Notification/interface/SimActivityRegistry.h"
#include "SimG4Core/Notification/interface/BeginOfJob.h"
#include "SimG4Core/Notification/interface/CMSSteppingVerbose.h"
#include "SimG4Core/Notification/interface/SimTrackManager.h"
#include "SimG4Core/Notification/interface/CurrentG4Track.h"
#include "SimG4Core/Watcher/interface/SimWatcherFactory.h"

#include "SimG4Core/Geometry/interface/DDDWorld.h"
#include "SimG4Core/MagneticField/interface/FieldBuilder.h"
#include "SimG4Core/MagneticField/interface/CMSFieldManager.h"

#include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
#include "SimDataFormats/GeneratorProducts/interface/HepMC3Product.h"
#include "DataFormats/GeometryVector/interface/GlobalPoint.h"
#include "HepMC3/Print.h"

#include "SimG4Core/Physics/interface/PhysicsList.h"
#include "SimG4Core/Physics/interface/CMSG4TrackInterface.h"

#include "SimG4Core/SensitiveDetector/interface/AttachSD.h"
#include "SimG4Core/SensitiveDetector/interface/SensitiveTkDetector.h"
#include "SimG4Core/SensitiveDetector/interface/SensitiveCaloDetector.h"
#include "SimG4Core/SensitiveDetector/interface/SensitiveDetectorMakerBase.h"
#include "SimG4Core/SensitiveDetector/interface/sensitiveDetectorMakers.h"

#include "G4Timer.hh"
#include "G4Event.hh"
#include "G4Run.hh"
#include <CLHEP/Units/SystemOfUnits.h>
#include "G4Threading.hh"
#include "G4UImanager.hh"
#include "G4WorkerThread.hh"
#include "G4WorkerRunManagerKernel.hh"
#include "G4StateManager.hh"
#include "G4TransportationManager.hh"
#include "G4LossTableManager.hh"
#include "G4PhysListUtil.hh"
#include "G4Field.hh"
#include "G4FieldManager.hh"
#include "G4ScoringManager.hh"
#include "G4UserSteppingAction.hh"
#include "G4GDMLParser.hh"
#include "G4Threading.hh"

#include <atomic>
#include <memory>
#include <thread>
#include <sstream>
#include <vector>

static std::once_flag applyOnce;
static std::once_flag applyOnceEnd;
static std::once_flag applyOnceGDML;

// from https://hypernews.cern.ch/HyperNews/CMS/get/edmFramework/3302/2.html
namespace {
  std::atomic<int> thread_counter{0};

  int get_new_thread_index() { return thread_counter++; }

  bool createWatchers(const edm::ParameterSet& iP,
                      SimActivityRegistry* iReg,
                      std::vector<std::shared_ptr<SimWatcher>>& oWatchers,
                      std::vector<std::shared_ptr<SimProducer>>& oProds,
                      int threadID) {
    std::vector<edm::ParameterSet> ws = iP.getParameter<std::vector<edm::ParameterSet>>("Watchers");
    if (ws.empty()) {
      return false;
    }

    for (auto& watcher : ws) {
      std::unique_ptr<SimWatcherMakerBase> maker(
          SimWatcherFactory::get()->create(watcher.getParameter<std::string>("type")));
      if (maker == nullptr) {
        throw cms::Exception("Configuration")
            << "RunManagerMTWorker::createWatchers: "
            << "Unable to find the requested Watcher " << watcher.getParameter<std::string>("type");
      } else {
        std::shared_ptr<SimWatcher> watcherTemp;
        std::shared_ptr<SimProducer> producerTemp;
        maker->make(watcher, *(iReg), watcherTemp, producerTemp);
        if (nullptr != watcherTemp) {
          if (!watcherTemp->isMT() && 0 < threadID) {
            throw cms::Exception("Configuration")
                << "RunManagerMTWorker::createWatchers: "
                << "Unable to use Watcher " << watcher.getParameter<std::string>("type") << " if number of threads > 1";
          } else {
            oWatchers.push_back(watcherTemp);
            if (nullptr != producerTemp) {
              oProds.push_back(producerTemp);
            }
          }
        }
      }
    }
    return (!oWatchers.empty());
  }
}  // namespace

struct RunManagerMTWorker::TLSData {
  std::unique_ptr<G4RunManagerKernel> kernel;  //must be deleted last
  std::unique_ptr<RunAction> userRunAction;
  std::unique_ptr<SimRunInterface> runInterface;
  std::unique_ptr<SimActivityRegistry> registry;
  std::unique_ptr<SimTrackManager> trackManager;
  std::vector<SensitiveTkDetector*> sensTkDets;
  std::vector<SensitiveCaloDetector*> sensCaloDets;
  std::vector<std::shared_ptr<SimWatcher>> watchers;
  std::vector<std::shared_ptr<SimProducer>> producers;
  // G4Run can only be deleted if there is a G4RunManager
  // on the thread where the G4Run is being deleted,
  // else it causes a segmentation fault
  G4Run* currentRun = nullptr;
  std::unique_ptr<G4Event> currentEvent;
  edm::RunNumber_t currentRunNumber = 0;
  bool threadInitialized = false;
  bool runTerminated = false;

  TLSData() {}

  ~TLSData() = default;
};

// This can not be a smart pointer since we must delete some of the members
// before leaving main() else we get a segmentation fault caused by accessing
// other 'singletons' after those singletons have been deleted. Instead we
// atempt to delete all TLS at RunManagerMTWorker destructor. If that fails for
// some reason, it is better to leak than cause a crash.
// thread_local RunManagerMTWorker::TLSData* RunManagerMTWorker::m_tls{nullptr};

RunManagerMTWorker::RunManagerMTWorker(const edm::ParameterSet& p, edm::ConsumesCollector&& iC)
    : m_generator(p.getParameter<edm::ParameterSet>("Generator")),
      m_generator3(p.getParameter<edm::ParameterSet>("Generator")),
      m_InToken(iC.consumes<edm::HepMCProduct>(
          p.getParameter<edm::ParameterSet>("Generator").getParameter<edm::InputTag>("HepMCProductLabel"))),
      m_InToken3(iC.consumes<edm::HepMC3Product>(
          p.getParameter<edm::ParameterSet>("Generator").getParameter<edm::InputTag>("HepMCProductLabel"))),
      m_theLHCTlinkToken(iC.consumes<edm::LHCTransportLinkContainer>(p.getParameter<edm::InputTag>("theLHCTlinkTag"))),
      m_nonBeam(p.getParameter<bool>("NonBeamEvent")),
      m_UseG4EventManager(p.getParameter<bool>("UseG4EventManager")),
      m_pUseMagneticField(p.getParameter<bool>("UseMagneticField")),
      m_LHCTransport(p.getParameter<bool>("LHCTransport")),
      m_thread_index{(get_new_thread_index())},
      m_pField(p.getParameter<edm::ParameterSet>("MagneticField")),
      m_pRunAction(p.getParameter<edm::ParameterSet>("RunAction")),
      m_pEventAction(p.getParameter<edm::ParameterSet>("EventAction")),
      m_pStackingAction(p.getParameter<edm::ParameterSet>("StackingAction")),
      m_pTrackingAction(p.getParameter<edm::ParameterSet>("TrackingAction")),
      m_pSteppingAction(p.getParameter<edm::ParameterSet>("SteppingAction")),
      m_G4Commands(p.getParameter<std::vector<std::string>>("G4Commands")),
      m_G4CommandsEndRun(p.getParameter<std::vector<std::string>>("G4CommandsEndRun")),
      m_p(p) {
  int id = getThreadIndex();
  if (id + 1 > CurrentG4Track::numberOfWorkers()) {
    CurrentG4Track::setNumberOfWorkers(id + 1);
  }
  edm::LogVerbatim("SimG4CoreApplication") << "RunManagerMTWorker for the thread " << id;

  // Initialize per-thread output
  CMSG4TrackInterface::instance()->setThreadID(id);
  G4Threading::G4SetThreadId(id);
  G4UImanager::GetUIpointer()->SetUpForAThread(id);
  auto iPset = p.getUntrackedParameter<edm::ParameterSet>("CustomUIsession");
  const std::string& uitype = iPset.getUntrackedParameter<std::string>("Type");
  if (uitype == "MessageLogger") {
    m_UIsession = new CustomUIsession();
  } else if (uitype == "MessageLoggerThreadPrefix") {
    m_UIsession = new CustomUIsessionThreadPrefix(iPset.getUntrackedParameter<std::string>("ThreadPrefix"), id);
  } else if (uitype == "FilePerThread") {
    m_UIsession = new CustomUIsessionToFile(iPset.getUntrackedParameter<std::string>("ThreadFile"), id);
  } else {
    throw cms::Exception("Configuration")
        << "RunManagerMTWorker::initializeG4: Invalid value of CustomUIsession.Type '" << uitype
        << "', valid are MessageLogger, MessageLoggerThreadPrefix, FilePerThread";
  }
  G4UImanager::GetUIpointer()->SetCoutDestination(m_UIsession);

  // sensitive detectors
  std::vector<std::string> onlySDs = p.getParameter<std::vector<std::string>>("OnlySDs");
  m_sdMakers = sim::sensitiveDetectorMakers(p, iC, onlySDs);

  // TLS and watchers
  initializeTLS();
  if (m_hasWatchers) {
    for (auto& watcher : m_tls->watchers) {
      watcher->registerConsumes(iC);
    }
  }
  if (m_LHCTransport) {
    m_LHCToken = iC.consumes<edm::HepMCProduct>(edm::InputTag("LHCTransport"));
  }
  if (m_pUseMagneticField) {
    m_MagField = iC.esConsumes<MagneticField, IdealMagneticFieldRecord, edm::Transition::BeginRun>();
  }
  edm::LogVerbatim("SimG4CoreApplication") << "RunManagerMTWorker is constructed for the thread " << id;
  unsigned int k = 0;
  for (std::unordered_map<std::string, std::unique_ptr<SensitiveDetectorMakerBase>>::const_iterator itr =
           m_sdMakers.begin();
       itr != m_sdMakers.end();
       ++itr, ++k) {
    edm::LogVerbatim("SimG4CoreApplication") << "SD[" << k << "] " << itr->first;
  }
}

RunManagerMTWorker::~RunManagerMTWorker() {
  m_tls = nullptr;
  delete m_UIsession;
}

void RunManagerMTWorker::beginRun(edm::EventSetup const& es) {
  int id = getThreadIndex();
  edm::LogVerbatim("SimG4CoreApplication") << "RunManagerMTWorker::beginRun for the thread " << id;
  for (auto& maker : m_sdMakers) {
    maker.second->beginRun(es);
  }
  if (m_pUseMagneticField) {
    m_pMagField = &es.getData(m_MagField);
  }
  if (m_hasWatchers) {
    for (auto& watcher : m_tls->watchers) {
      watcher->beginRun(es);
    }
  }
  edm::LogVerbatim("SimG4CoreApplication") << "RunManagerMTWorker::beginRun done for the thread " << id;
}

void RunManagerMTWorker::endRun() {
  int id = getThreadIndex();
  edm::LogVerbatim("SimG4CoreApplication") << "RunManagerMTWorker::endRun for the thread " << id;
  terminateRun();
}

void RunManagerMTWorker::initializeTLS() {
  if (nullptr != m_tls) {
    return;
  }

  m_tls = new TLSData();
  m_tls->registry = std::make_unique<SimActivityRegistry>();

  edm::Service<SimActivityRegistry> otherRegistry;
  //Look for an outside SimActivityRegistry
  // this is used by the visualization code
  int thisID = getThreadIndex();
  if (otherRegistry) {
    m_tls->registry->connect(*otherRegistry);
    if (thisID > 0) {
      throw cms::Exception("Configuration")
          << "RunManagerMTWorker::initializeTLS : "
          << "SimActivityRegistry service (i.e. visualization) is not supported for more than 1 thread. "
          << " \n If this use case is needed, RunManagerMTWorker has to be updated.";
    }
  }
  m_hasWatchers = createWatchers(m_p, m_tls->registry.get(), m_tls->watchers, m_tls->producers, thisID);
}

void RunManagerMTWorker::initializeG4(RunManagerMT* runManagerMaster, const edm::EventSetup& es) {
  if (m_tls->threadInitialized)
    return;

  G4Timer timer;
  timer.Start();

  // I guess everything initialized here should be in thread_local storage
  initializeTLS();

  int thisID = getThreadIndex();
  edm::LogVerbatim("SimG4CoreApplication") << "RunManagerMTWorker::initializeG4 in thread " << thisID << " is started";

  // Initialize worker part of shared resources (geometry, physics)
  G4WorkerThread::BuildGeometryAndPhysicsVector();

  // Create worker run manager
  m_tls->kernel.reset(G4WorkerRunManagerKernel::GetRunManagerKernel());
  if (nullptr == m_tls->kernel) {
    m_tls->kernel = std::make_unique<G4WorkerRunManagerKernel>();
  }

  // Define G4 exception handler
  double th = m_p.getParameter<double>("ThresholdForGeometryExceptions") * CLHEP::GeV;
  bool tr = m_p.getParameter<bool>("TraceExceptions");
  G4StateManager::GetStateManager()->SetExceptionHandler(new ExceptionHandler(th, tr));

  // Set the geometry for the worker, share from master
  auto worldPV = runManagerMaster->world().GetWorldVolume();
  m_tls->kernel->WorkerDefineWorldVolume(worldPV);
  G4TransportationManager* tM = G4TransportationManager::GetTransportationManager();
  tM->SetWorldForTracking(worldPV);

  // flag of Phase2
  m_isPhase2 = runManagerMaster->isPhase2();

  // we need the track manager now
  int verbose = m_p.getParameter<int>("EventVerbose");
  m_tls->trackManager = std::make_unique<SimTrackManager>(&m_simEvent, verbose);

  // setup the magnetic field
  if (m_pUseMagneticField) {
    const GlobalPoint gg(0.f, 0.f, 0.f);

    sim::FieldBuilder fieldBuilder(m_pMagField, m_pField);

    CMSFieldManager* fieldManager = new CMSFieldManager();
    tM->SetFieldManager(fieldManager);
    fieldBuilder.build(fieldManager, tM->GetPropagatorInField());

    std::string fieldFile = m_p.getUntrackedParameter<std::string>("FileNameField", "");
    if (!fieldFile.empty()) {
      std::call_once(applyOnce, [this]() { m_dumpMF = true; });
      if (m_dumpMF) {
        edm::LogVerbatim("SimG4CoreApplication")
            << "RunManagerMTWorker::InitializeG4: Dump magnetic field to file " << fieldFile;
        DumpMagneticField(tM->GetFieldManager()->GetDetectorField(), fieldFile);
      }
    }
  }

  // attach sensitive detector
  auto sensDets = sim::attachSD(
      m_sdMakers, es, runManagerMaster->catalog(), m_p, m_tls->trackManager.get(), *(m_tls->registry.get()));

  m_tls->sensTkDets = sensDets.first;
  m_tls->sensCaloDets = sensDets.second;

  edm::LogVerbatim("SimG4CoreApplication")
      << "RunManagerMTWorker::InitializeG4: Sensitive Detectors are built in thread " << thisID << " found "
      << m_tls->sensTkDets.size() << " Tk type SD, and " << m_tls->sensCaloDets.size() << " Calo type SD";

  // geometry dump
  G4String writeFile = (G4String)m_p.getUntrackedParameter<std::string>("FileNameGDML");
  if (!writeFile.empty()) {
    std::call_once(applyOnceGDML, [this]() { m_dumpGDML = true; });
    edm::LogVerbatim("SimG4CoreApplication") << "DumpGDML:" << m_dumpGDML;
    if (m_dumpGDML) {
      G4int thID = G4Threading::G4GetThreadId();
      edm::LogVerbatim("SimG4CoreApplication") << "ThreadID=" << thID;
      G4Threading::G4SetThreadId(-1);
      G4GDMLParser gdml;
      gdml.SetRegionExport(true);
      gdml.SetEnergyCutsExport(true);
      gdml.SetSDExport(true);
      gdml.Write(writeFile, worldPV, true);
      G4Threading::G4SetThreadId(thID);
    }
  }

  // Enable couple transportation
  bool scorer = m_p.getParameter<bool>("UseCommandBaseScorer");
  if (scorer) {
    G4ScoringManager* scManager = G4ScoringManager::GetScoringManager();
    scManager->SetVerboseLevel(1);
  }

  // Set the physics list for the worker, share from master
  PhysicsList* physicsList = runManagerMaster->physicsListForWorker();

  edm::LogVerbatim("SimG4CoreApplication")
      << "RunManagerMTWorker::InitializeG4: start initialisation of PhysicsList for the thread " << thisID;

  // Geant4 UI commands in PreInit state
  if (!runManagerMaster->G4Commands().empty()) {
    G4cout << "RunManagerMTWorker::InitializeG4: Requested UI commands: " << G4endl;
    for (const std::string& command : runManagerMaster->G4Commands()) {
      G4cout << "          " << command << G4endl;
      G4UImanager::GetUIpointer()->ApplyCommand(command);
    }
  }
  G4StateManager::GetStateManager()->SetNewState(G4State_Init);

  physicsList->InitializeWorker();
  m_tls->kernel->SetPhysics(physicsList);
  m_tls->kernel->InitializePhysics();

  if (!m_tls->kernel->RunInitialization()) {
    throw cms::Exception("Configuration")
        << "RunManagerMTWorker::InitializeG4: Geant4 kernel initialization failed in thread " << thisID;
  }

  //tell all interesting parties that we are beginning the job
  BeginOfJob aBeginOfJob(&es);
  m_tls->registry->beginOfJobSignal_(&aBeginOfJob);

  G4int sv = m_p.getUntrackedParameter<int>("SteppingVerbosity");

  if (sv > 0) {
    G4double elim = m_p.getUntrackedParameter<double>("StepVerboseThreshold", 0.1) * CLHEP::GeV;
    std::vector<int> ve = m_p.getUntrackedParameter<std::vector<int>>("VerboseEvents");
    std::vector<int> vn = m_p.getUntrackedParameter<std::vector<int>>("VertexNumber");
    std::vector<int> vt = m_p.getUntrackedParameter<std::vector<int>>("VerboseTracks");
    m_sVerbose = std::make_unique<CMSSteppingVerbose>(sv, elim, ve, vn, vt);
  }
  if (!m_UseG4EventManager)
    m_evtManager = std::make_unique<CMSSimEventManager>(m_p);
  initializeUserActions();

  G4StateManager::GetStateManager()->SetNewState(G4State_Idle);

  timer.Stop();
  edm::LogVerbatim("SimG4CoreApplication")
      << "RunManagerMTWorker::initializeG4 done for the thread " << thisID << "  " << timer;
  m_tls->threadInitialized = true;
}

void RunManagerMTWorker::initializeUserActions() {
  m_tls->runInterface = std::make_unique<SimRunInterface>(this, false);
  m_tls->userRunAction = std::make_unique<RunAction>(m_pRunAction, m_tls->runInterface.get(), false);
  m_tls->userRunAction->SetMaster(false);
  Connect(m_tls->userRunAction.get());

  G4int ver = m_p.getParameter<int>("EventVerbose");
  G4EventManager* eventManager = m_tls->kernel->GetEventManager();
  eventManager->SetVerboseLevel(ver);

  // different event actions for Run2,3 and Phase2
  G4UserEventAction* userEventAction;
  if (m_isPhase2) {
    auto ptr =
        new Phase2EventAction(m_pEventAction, m_tls->runInterface.get(), m_tls->trackManager.get(), m_sVerbose.get());
    Connect(ptr);
    userEventAction = (G4UserEventAction*)ptr;
  } else {
    auto ptr = new EventAction(m_pEventAction, m_tls->runInterface.get(), m_tls->trackManager.get(), m_sVerbose.get());
    Connect(ptr);
    userEventAction = (G4UserEventAction*)ptr;
  }

  // different tracking actions for Run2,3 and Phase2
  G4UserTrackingAction* userTrackingAction;
  if (m_isPhase2) {
    auto ptr = new Phase2TrackingAction(m_tls->trackManager.get(), m_sVerbose.get(), m_pTrackingAction);
    Connect(ptr);
    userTrackingAction = (G4UserTrackingAction*)ptr;
  } else {
    auto ptr = new TrackingAction(m_tls->trackManager.get(), m_sVerbose.get(), m_pTrackingAction);
    Connect(ptr);
    userTrackingAction = (G4UserTrackingAction*)ptr;
  }

  // different stepping actions for Run2,3 and Phase2
  G4UserSteppingAction* userSteppingAction;
  bool dd4hep = m_p.getParameter<bool>("g4GeometryDD4hepSource");
  if (m_isPhase2) {
    auto ptr = new Phase2SteppingAction(m_sVerbose.get(), m_pSteppingAction, m_pStackingAction, m_hasWatchers, dd4hep);
    Connect(ptr);
    userSteppingAction = (G4UserSteppingAction*)ptr;
  } else {
    auto ptr = new SteppingAction(m_sVerbose.get(), m_pSteppingAction, m_hasWatchers, dd4hep);
    Connect(ptr);
    userSteppingAction = (G4UserSteppingAction*)ptr;
  }

  // stacking actions and event manager
  G4UserStackingAction* userStackingAction{nullptr};
  if (m_isPhase2) {
    userStackingAction = new Phase2StackingAction(m_pStackingAction, m_sVerbose.get());
  } else {
    userStackingAction = new StackingAction(m_pStackingAction, m_sVerbose.get());
  }
  if (m_UseG4EventManager) {
    eventManager->SetUserAction(userEventAction);
    eventManager->SetUserAction(userTrackingAction);
    eventManager->SetUserAction(userSteppingAction);
    eventManager->SetUserAction(userStackingAction);
  } else {
    m_evtManager->SetUserAction(userEventAction);
    m_evtManager->SetUserAction(userTrackingAction);
    m_evtManager->SetUserAction(userSteppingAction);
    m_evtManager->SetUserAction(userStackingAction);
  }
}

void RunManagerMTWorker::Connect(RunAction* runAction) {
  runAction->m_beginOfRunSignal.connect(m_tls->registry->beginOfRunSignal_);
  runAction->m_endOfRunSignal.connect(m_tls->registry->endOfRunSignal_);
}

void RunManagerMTWorker::Connect(EventAction* eventAction) {
  eventAction->m_beginOfEventSignal.connect(m_tls->registry->beginOfEventSignal_);
  eventAction->m_endOfEventSignal.connect(m_tls->registry->endOfEventSignal_);
}

void RunManagerMTWorker::Connect(Phase2EventAction* eventAction) {
  eventAction->m_beginOfEventSignal.connect(m_tls->registry->beginOfEventSignal_);
  eventAction->m_endOfEventSignal.connect(m_tls->registry->endOfEventSignal_);
}

void RunManagerMTWorker::Connect(TrackingAction* trackingAction) {
  trackingAction->m_beginOfTrackSignal.connect(m_tls->registry->beginOfTrackSignal_);
  trackingAction->m_endOfTrackSignal.connect(m_tls->registry->endOfTrackSignal_);
}

void RunManagerMTWorker::Connect(Phase2TrackingAction* trackingAction) {
  trackingAction->m_beginOfTrackSignal.connect(m_tls->registry->beginOfTrackSignal_);
  trackingAction->m_endOfTrackSignal.connect(m_tls->registry->endOfTrackSignal_);
}

void RunManagerMTWorker::Connect(SteppingAction* steppingAction) {
  steppingAction->m_g4StepSignal.connect(m_tls->registry->g4StepSignal_);
}

void RunManagerMTWorker::Connect(Phase2SteppingAction* steppingAction) {
  steppingAction->m_g4StepSignal.connect(m_tls->registry->g4StepSignal_);
}

SimTrackManager* RunManagerMTWorker::getSimTrackManager() {
  initializeTLS();
  return m_tls->trackManager.get();
}
std::vector<SensitiveTkDetector*>& RunManagerMTWorker::sensTkDetectors() {
  initializeTLS();
  return m_tls->sensTkDets;
}
std::vector<SensitiveCaloDetector*>& RunManagerMTWorker::sensCaloDetectors() {
  initializeTLS();
  return m_tls->sensCaloDets;
}
std::vector<std::shared_ptr<SimProducer>>& RunManagerMTWorker::producers() {
  initializeTLS();
  return m_tls->producers;
}

void RunManagerMTWorker::initializeRun() {
  m_tls->currentRun = new G4Run();
  G4StateManager::GetStateManager()->SetNewState(G4State_GeomClosed);
  if (nullptr != m_tls->userRunAction) {
    m_tls->userRunAction->BeginOfRunAction(m_tls->currentRun);
  }
  int id = getThreadIndex();
  edm::LogVerbatim("SimG4CoreApplication") << "RunManagerMTWorker::initializeRun done for thread " << id;
}

void RunManagerMTWorker::terminateRun() {
  int id = getThreadIndex();
  edm::LogVerbatim("SimG4CoreApplication") << "RunManagerMTWorker::terminateRun for thread " << id;
  if (nullptr == m_tls || m_tls->runTerminated) {
    return;
  }

  // Geant4 UI commands after the run
  if (!m_G4CommandsEndRun.empty()) {
    std::call_once(applyOnceEnd, [this]() { m_endOfRun = true; });
    if (m_endOfRun) {
      edm::LogVerbatim("SimG4CoreApplication") << "RunManagerMTWorker: Requested end of run UI commands: ";
      for (const std::string& command : m_G4CommandsEndRun) {
        edm::LogVerbatim("SimG4CoreApplication") << "    " << command;
        G4UImanager::GetUIpointer()->ApplyCommand(command);
      }
    }
  }
  if (m_tls->userRunAction) {
    m_tls->userRunAction->EndOfRunAction(m_tls->currentRun);
    m_tls->userRunAction.reset();
  }
  m_tls->currentEvent.reset();
  if (m_tls->kernel) {
    m_tls->kernel->RunTermination();
  }
  m_tls->runTerminated = true;
  edm::LogVerbatim("SimG4CoreApplication") << "RunManagerMTWorker::terminateRun done for thread " << id;
}

TmpSimEvent* RunManagerMTWorker::produce(const edm::Event& inpevt,
                                         const edm::EventSetup& es,
                                         RunManagerMT& runManagerMaster) {
  // We have to do the per-thread initialization, and per-thread
  // per-run initialization here by ourselves.

  edm::LogVerbatim("SimG4CoreApplication") << "RunManagerMTWorker::produce: start EventID=" << inpevt.id().event();

  assert(m_tls != nullptr and m_tls->threadInitialized);
  // Initialize run
  if (inpevt.id().run() != m_tls->currentRunNumber) {
    edm::LogVerbatim("SimG4CoreApplication")
        << "RunManagerMTWorker::produce: RunID= " << inpevt.id().run() << "  TLS RunID= " << m_tls->currentRunNumber;
    if (m_tls->currentRunNumber != 0 && !m_tls->runTerminated) {
      // If previous run in this thread was not terminated via endRun() call,
      // then terminate it now
      terminateRun();
    }
    initializeRun();
    m_tls->currentRunNumber = inpevt.id().run();
  }

  // event and primary

  m_tls->currentEvent.reset(generateEvent(inpevt));
  m_simEvent.clear();
  m_simEvent.setHepEvent(m_generator.genEvent());
  m_simEvent.setWeight(m_generator.eventWeight());

  if (m_generator.genVertex() != nullptr) {
    auto genVertex = m_generator.genVertex();
    m_simEvent.collisionPoint(math::XYZTLorentzVectorD(genVertex->x() / CLHEP::cm,
                                                       genVertex->y() / CLHEP::cm,
                                                       genVertex->z() / CLHEP::cm,
                                                       genVertex->t() / CLHEP::second));
  }
  if (m_tls->currentEvent->GetNumberOfPrimaryVertex() == 0) {
    throw cms::Exception("EventCorruption")
        << "RunManagerMTWorker::produce: event " << inpevt.id().event() << " with no G4PrimaryVertices"
        << " StreamID=" << inpevt.streamID() << " threadIndex=" << getThreadIndex();

  } else {
    edm::LogVerbatim("SimG4CoreApplication")
        << "RunManagerMTWorker::produce: start EventID=" << inpevt.id().event() << " StreamID=" << inpevt.streamID()
        << " threadIndex=" << getThreadIndex() << " weight=" << m_simEvent.weight()
        << " Nprimary: " << m_tls->currentEvent->GetNumberOfPrimaryVertex() << " vertices and " << m_simEvent.nTracks()
        << " particles";
    // process event
    if (m_UseG4EventManager) {
      m_tls->kernel->GetEventManager()->ProcessOneEvent(m_tls->currentEvent.get());
    } else {
      m_evtManager->ProcessOneEvent(m_tls->currentEvent.get());
    }
  }

  // remove memory only needed during event processing
  m_tls->currentEvent.reset();
  for (auto& sd : m_tls->sensCaloDets) {
    sd->reset();
  }

  edm::LogVerbatim("SimG4CoreApplication") << "RunManagerMTWorker::produce: ended Event " << inpevt.id().event();
  return &m_simEvent;
}

void RunManagerMTWorker::abortEvent() {
  if (m_tls->runTerminated) {
    return;
  }
  G4Track* t = m_tls->kernel->GetEventManager()->GetTrackingManager()->GetTrack();
  t->SetTrackStatus(fStopAndKill);

  // CMS-specific act
  //
  TrackingAction* uta = static_cast<TrackingAction*>(m_tls->kernel->GetEventManager()->GetUserTrackingAction());
  uta->PostUserTrackingAction(t);

  m_tls->currentEvent->SetEventAborted();
  m_tls->kernel->GetEventManager()->GetStackManager()->clear();
  m_tls->kernel->GetEventManager()->GetTrackingManager()->EventAborted();
}

void RunManagerMTWorker::abortRun(bool softAbort) {
  if (!softAbort) {
    abortEvent();
  }
  m_tls->currentRun = nullptr;
  terminateRun();
}

G4Event* RunManagerMTWorker::generateEvent(const edm::Event& inpevt) {
  m_tls->currentEvent.reset();

  // 64 bits event ID in CMSSW converted into Geant4 event ID
  G4int evtid = (G4int)inpevt.id().event();
  G4Event* evt = new G4Event(evtid);

  // required to reset the GenParticle Id for particles transported
  // along the beam pipe to their original value for SimTrack creation
  resetGenParticleId(inpevt);

  edm::Handle<edm::HepMCProduct> HepMCEvt;
  bool found = inpevt.getByToken(m_InToken, HepMCEvt);

  if (found) {  // HepMC event exists

    m_generator.setGenEvent(HepMCEvt->GetEvent());

    if (!m_nonBeam) {
      m_generator.HepMC2G4(HepMCEvt->GetEvent(), evt);
      if (m_LHCTransport) {
        edm::Handle<edm::HepMCProduct> LHCMCEvt;
        inpevt.getByToken(m_LHCToken, LHCMCEvt);
        m_generator.nonCentralEvent2G4(LHCMCEvt->GetEvent(), evt);
      }
    } else {
      m_generator.nonCentralEvent2G4(HepMCEvt->GetEvent(), evt);
    }

  } else {  // no HepMC event, try to get HepMC3 event

    edm::Handle<edm::HepMC3Product> HepMCEvt3;
    inpevt.getByToken(m_InToken3, HepMCEvt3);

    HepMC3::GenEvent genevt3;
    genevt3.read_data(*HepMCEvt3->GetEvent());
    m_generator3.setGenEvent(&genevt3);

    if (!m_nonBeam) {
      m_generator3.HepMC2G4(&genevt3, evt);
      if (m_LHCTransport) {
        edm::Handle<edm::HepMC3Product> LHCMCEvt;
        inpevt.getByToken(m_LHCToken, LHCMCEvt);
        //m_generator3.nonCentralEvent2G4(LHCMCEvt->GetEvent(), evt);
      }
    } else {
      //m_generator3.nonCentralEvent2G4(HepMCEvt->GetEvent(), evt);
    }
  }
  return evt;
}

void RunManagerMTWorker::resetGenParticleId(const edm::Event& inpevt) {
  edm::Handle<edm::LHCTransportLinkContainer> theLHCTlink;
  inpevt.getByToken(m_theLHCTlinkToken, theLHCTlink);
  if (theLHCTlink.isValid()) {
    m_tls->trackManager->setLHCTransportLink(theLHCTlink.product());
  }
}

void RunManagerMTWorker::DumpMagneticField(const G4Field* field, const std::string& file) const {
  std::ofstream fout(file.c_str(), std::ios::out);
  if (fout.fail()) {
    edm::LogWarning("SimG4CoreApplication")
        << "MTWorker::DumpMagneticField: error opening file <" << file << "> for magnetic field";
  } else {
    // CMS magnetic field volume
    double rmax = 9000 * CLHEP::mm;
    double zmax = 24000 * CLHEP::mm;

    double dx = 10 * CLHEP::cm;
    double dz = 10 * CLHEP::cm;

    int nx = G4lrint(2 * rmax / dx);
    int nz = G4lrint(2 * zmax / dz);

    double point[4] = {0.0, 0.0, 0.0, 0.0};
    double bfield[3] = {0.0, 0.0, 0.0};

    double d2 = 1. / CLHEP::mm;
    double d3 = 1. / CLHEP::tesla;

    fout << std::setprecision(6);
    fout << "### " << file << " CMS magnetic field: X(mm) Y(mm) Z(mm) Bx(tesla) By(tesla) Bz(tesla)  ###" << G4endl;
    for (int k = 0; k <= nx; ++k) {
      double x = k * dx - rmax;
      for (int i = 0; i <= nx; ++i) {
        double y = i * dx - rmax;
        for (int j = 0; j <= nz; ++j) {
          double z = j * dz - zmax;
          point[0] = x;
          point[1] = y;
          point[2] = z;
          field->GetFieldValue(point, bfield);
          fout << x * d2 << " " << y * d2 << " " << z * d2 << " " << bfield[0] * d3 << " " << bfield[1] * d3 << " "
               << bfield[2] * d3 << G4endl;
        }
      }
    }
    fout.close();
  }
}