Merge pull request #27543 from civanch/g4_dd4hep_geometry

Geant4 initialisation clean-up for DD4Hep migration

SimG4Core/Application/interface/OscarMTMasterThread.h

@@ -46,7 +46,6 @@ class OscarMTMasterThread {
 
   enum class ThreadState { NotExist = 0, BeginRun = 1, EndRun = 2, Destruct = 3 };
 
-  const bool m_pUseMagneticField;
   const bool m_pGeoFromDD4hep;
 
   std::shared_ptr<RunManagerMT> m_runManagerMaster;
@@ -57,7 +56,6 @@ class OscarMTMasterThread {
   mutable edm::ESWatcher<IdealMagneticFieldRecord> idealMagRcdWatcher_;
   mutable const DDCompactView* m_pDD;
   mutable const cms::DDCompactView* m_pDD4hep;
-  mutable const MagneticField* m_pMF;
   mutable const HepPDT::ParticleDataTable* m_pTable;
 
   mutable std::mutex m_protectMutex;

SimG4Core/Application/interface/RunManagerMT.h

@@ -31,7 +31,6 @@ namespace cms {
 }
 
 class DDDWorld;
-class MagneticField;
 
 class G4MTRunManagerKernel;
 class G4Run;
@@ -60,7 +59,8 @@ class RunManagerMT {
   explicit RunManagerMT(edm::ParameterSet const&);
   ~RunManagerMT();
 
-  void initG4(const DDCompactView*, const cms::DDCompactView*, const MagneticField*, const HepPDT::ParticleDataTable*);
+  //  void initG4(const DDCompactView*, const cms::DDCompactView*, const MagneticField*, const HepPDT::ParticleDataTable*);
+  void initG4(const DDCompactView*, const cms::DDCompactView*, const HepPDT::ParticleDataTable*);
 
   void initializeUserActions();
 
@@ -92,7 +92,6 @@ class RunManagerMT {
   std::unique_ptr<PhysicsList> m_physicsList;
   bool m_managerInitialized;
   bool m_runTerminated;
-  bool m_pUseMagneticField;
   RunAction* m_userRunAction;
   G4Run* m_currentRun;
   G4StateManager* m_stateManager;

SimG4Core/Application/src/OscarMTMasterThread.cc

@@ -8,9 +8,6 @@
 #include "FWCore/Framework/interface/ESTransientHandle.h"
 #include "FWCore/Utilities/interface/EDMException.h"
 
-#include "MagneticField/Engine/interface/MagneticField.h"
-#include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
-
 #include "Geometry/Records/interface/IdealGeometryRecord.h"
 #include "DetectorDescription/Core/interface/DDCompactView.h"
 #include "DetectorDescription/DDCMS/interface/DDCompactView.h"
@@ -21,11 +18,9 @@
 #include "G4PhysicalVolumeStore.hh"
 
 OscarMTMasterThread::OscarMTMasterThread(const edm::ParameterSet& iConfig)
-    : m_pUseMagneticField(iConfig.getParameter<bool>("UseMagneticField")),
-      m_pGeoFromDD4hep(iConfig.getParameter<bool>("g4GeometryDD4hepSource")),
+    : m_pGeoFromDD4hep(iConfig.getParameter<bool>("g4GeometryDD4hepSource")),
       m_pDD(nullptr),
       m_pDD4hep(nullptr),
-      m_pMF(nullptr),
       m_pTable(nullptr),
       m_masterThreadState(ThreadState::NotExist),
       m_masterCanProceed(false),
@@ -79,7 +74,7 @@ OscarMTMasterThread::OscarMTMasterThread(const edm::ParameterSet& iConfig)
       if (m_masterThreadState == ThreadState::BeginRun) {
         // Initialize Geant4
         edm::LogInfo("OscarMTMasterThread") << "Master thread: Initializing Geant4";
-        runManagerMaster->initG4(m_pDD, m_pDD4hep, m_pMF, m_pTable);
+        runManagerMaster->initG4(m_pDD, m_pDD4hep, m_pTable);
         isG4Alive = true;
       } else if (m_masterThreadState == ThreadState::EndRun) {
         // Stop Geant4
@@ -195,15 +190,6 @@ void OscarMTMasterThread::readES(const edm::EventSetup& iSetup) const {
         << "The Geometry configuration is changed during the job execution\n"
         << "this is not allowed, the geometry must stay unchanged";
   }
-  if (m_pUseMagneticField) {
-    bool magChanged = idealMagRcdWatcher_.check(iSetup);
-    if (magChanged && (!m_firstRun)) {
-      throw edm::Exception(edm::errors::Configuration)
-          << "[SimG4Core OscarMTMasterThread]\n"
-          << "The MagneticField configuration is changed during the job execution\n"
-          << "this is not allowed, the MagneticField must stay unchanged";
-    }
-  }
   // Don't read from ES if not the first run, just as in
   if (!m_firstRun)
     return;
@@ -219,12 +205,6 @@ void OscarMTMasterThread::readES(const edm::EventSetup& iSetup) const {
     m_pDD = pDD.product();
   }
 
-  if (m_pUseMagneticField) {
-    edm::ESHandle<MagneticField> pMF;
-    iSetup.get<IdealMagneticFieldRecord>().get(pMF);
-    m_pMF = pMF.product();
-  }
-
   edm::ESHandle<HepPDT::ParticleDataTable> fTable;
   iSetup.get<PDTRecord>().get(fTable);
   m_pTable = fTable.product();

SimG4Core/Application/src/RunManager.cc

@@ -12,9 +12,7 @@
 #include "SimG4Core/Application/interface/ExceptionHandler.h"
 
 #include "SimG4Core/Geometry/interface/DDDWorld.h"
-#include "SimG4Core/Geometry/interface/G4LogicalVolumeToDDLogicalPartMap.h"
 #include "SimG4Core/Geometry/interface/SensitiveDetectorCatalog.h"
-#include "SimG4Core/Geometry/interface/DDG4ProductionCuts.h"
 
 #include "SimG4Core/SensitiveDetector/interface/AttachSD.h"
 
@@ -187,6 +185,16 @@ void RunManager::initG4(const edm::EventSetup& es) {
                                       << "The Geometry configuration is changed during the job execution\n"
                                       << "this is not allowed, the geometry must stay unchanged\n";
   }
+  bool geoFromDD4hep = m_p.getParameter<bool>("g4GeometryDD4hepSource");
+  bool cuts = m_pPhysics.getParameter<bool>("CutsPerRegion");
+  bool protonCut = m_pPhysics.getUntrackedParameter<bool>("CutsOnProton", true);
+  int verb =
+      std::max(m_pPhysics.getUntrackedParameter<int>("Verbosity", 0), m_p.getParameter<int>("SteppingVerbosity"));
+  edm::LogVerbatim("SimG4CoreApplication")
+      << "RunManagerMT: start initialising of geometry DD4Hep: " << geoFromDD4hep << "\n"
+      << "              cutsPerRegion: " << cuts << " cutForProton: " << protonCut << "\n"
+      << "              G4 verbosity: " << verb;
+
   if (m_pUseMagneticField) {
     bool magChanged = idealMagRcdWatcher_.check(es);
     if (magChanged && (!firstRun)) {
@@ -200,22 +208,22 @@ void RunManager::initG4(const edm::EventSetup& es) {
   if (m_managerInitialized)
     return;
 
-  // DDDWorld: get the DDCV from the ES and use it to build the World
-  edm::ESTransientHandle<DDCompactView> pDD;
-  es.get<IdealGeometryRecord>().get(pDD);
-
-  G4LogicalVolumeToDDLogicalPartMap map_;
-  SensitiveDetectorCatalog catalog_;
-  const DDDWorld* world = new DDDWorld(&(*pDD), map_, catalog_, false);
-  m_registry.dddWorldSignal_(world);
-
-  int verb =
-      std::max(m_pPhysics.getUntrackedParameter<int>("Verbosity", 0), m_p.getParameter<int>("SteppingVerbosity"));
-  m_kernel->SetVerboseLevel(verb);
-  auto cuts = m_pPhysics.getParameter<bool>("CutsPerRegion");
-  if (cuts) {
-    DDG4ProductionCuts pcuts(map_, verb, m_pPhysics);
+  // initialise geometry
+  const DDCompactView* pDD = nullptr;
+  const cms::DDCompactView* pDD4hep = nullptr;
+  if (geoFromDD4hep) {
+    edm::ESTransientHandle<cms::DDCompactView> pDD;
+    es.get<IdealGeometryRecord>().get(pDD);
+    pDD4hep = pDD.product();
+  } else {
+    edm::ESTransientHandle<DDCompactView> pDD;
+    es.get<IdealGeometryRecord>().get(pDD);
+    pDD = pDD.product();
   }
+  SensitiveDetectorCatalog catalog;
+  const DDDWorld* world = new DDDWorld(pDD, pDD4hep, catalog, verb, cuts, protonCut);
+  G4VPhysicalVolume* pworld = world->GetWorldVolume();
+  m_registry.dddWorldSignal_(world);
 
   if (m_pUseMagneticField) {
     // setup the magnetic field
@@ -237,20 +245,17 @@ void RunManager::initG4(const edm::EventSetup& es) {
   // we need the track manager now
   m_trackManager = std::unique_ptr<SimTrackManager>(new SimTrackManager);
 
-  {
-    // attach sensitive detector
-
-    AttachSD attach;
-    auto sensDets = attach.create(es, catalog_, m_p, m_trackManager.get(), m_registry);
+  // attach sensitive detector
+  AttachSD attach;
+  auto sensDets = attach.create(es, catalog, m_p, m_trackManager.get(), m_registry);
 
-    m_sensTkDets.swap(sensDets.first);
-    m_sensCaloDets.swap(sensDets.second);
-  }
+  m_sensTkDets.swap(sensDets.first);
+  m_sensCaloDets.swap(sensDets.second);
 
-  edm::LogInfo("SimG4CoreApplication") << " RunManager: Sensitive Detector "
-                                       << "building finished; found " << m_sensTkDets.size()
-                                       << " Tk type Producers, and " << m_sensCaloDets.size()
-                                       << " Calo type producers ";
+  edm::LogVerbatim("SimG4CoreApplication")
+      << " RunManager: Sensitive Detector "
+      << "building finished; found " << m_sensTkDets.size() << " Tk type Producers, and " << m_sensCaloDets.size()
+      << " Calo type producers ";
 
   edm::ESHandle<HepPDT::ParticleDataTable> fTable;
   es.get<PDTRecord>().get(fTable);
@@ -337,7 +342,7 @@ void RunManager::initG4(const edm::EventSetup& es) {
     G4GDMLParser gdml;
     gdml.SetRegionExport(true);
     gdml.SetEnergyCutsExport(true);
-    gdml.Write(m_WriteFile, world->GetWorldVolume(), true);
+    gdml.Write(m_WriteFile, pworld, true);
   }
 
   if (!m_RegionFile.empty()) {

SimG4Core/Application/src/RunManagerMT.cc

@@ -1,3 +1,4 @@
+
 #include "SimG4Core/Application/interface/RunManagerMT.h"
 #include "SimG4Core/Application/interface/PrimaryTransformer.h"
 #include "SimG4Core/Application/interface/SimRunInterface.h"
@@ -7,21 +8,12 @@
 #include "SimG4Core/Application/interface/ExceptionHandler.h"
 
 #include "SimG4Core/Geometry/interface/DDDWorld.h"
-#include "SimG4Core/Geometry/interface/G4LogicalVolumeToDDLogicalPartMap.h"
-#include "SimG4Core/Geometry/interface/DDG4ProductionCuts.h"
-
-#include "SimG4Core/SensitiveDetector/interface/AttachSD.h"
 
 #include "SimG4Core/Physics/interface/PhysicsListFactory.h"
 #include "SimG4Core/PhysicsLists/interface/CMSMonopolePhysics.h"
 #include "SimG4Core/CustomPhysics/interface/CMSExoticaPhysics.h"
 
 #include "SimG4Core/Watcher/interface/SimWatcherFactory.h"
-#include "SimG4Core/MagneticField/interface/FieldBuilder.h"
-#include "SimG4Core/MagneticField/interface/Field.h"
-#include "SimG4Core/MagneticField/interface/CMSFieldManager.h"
-
-#include "MagneticField/Engine/interface/MagneticField.h"
 
 #include "SimG4Core/Notification/interface/G4SimEvent.h"
 #include "SimG4Core/Notification/interface/SimTrackManager.h"
@@ -60,8 +52,6 @@
 #include "G4Region.hh"
 #include "G4RegionStore.hh"
 
-#include "DDG4/Geant4Mapping.h"
-
 #include <iostream>
 #include <sstream>
 #include <fstream>
@@ -72,7 +62,6 @@
 RunManagerMT::RunManagerMT(edm::ParameterSet const& p)
     : m_managerInitialized(false),
       m_runTerminated(false),
-      m_pUseMagneticField(p.getParameter<bool>("UseMagneticField")),
       m_PhysicsTablesDir(p.getParameter<std::string>("PhysicsTablesDirectory")),
       m_StorePhysicsTables(p.getParameter<bool>("StorePhysicsTables")),
       m_RestorePhysicsTables(p.getParameter<bool>("RestorePhysicsTables")),
@@ -102,35 +91,28 @@ RunManagerMT::~RunManagerMT() { stopG4(); }
 
 void RunManagerMT::initG4(const DDCompactView* pDD,
                           const cms::DDCompactView* pDD4hep,
-                          const MagneticField* pMF,
                           const HepPDT::ParticleDataTable* fPDGTable) {
   if (m_managerInitialized) {
     edm::LogWarning("SimG4CoreApplication") << "RunManagerMT::initG4 was already done - exit";
     return;
   }
-  auto geoFromDD4hep = m_p.getParameter<bool>("g4GeometryDD4hepSource");
-  edm::LogVerbatim("SimG4CoreApplication") << "RunManagerMT: start initialising of geometry DD4Hep: " << geoFromDD4hep;
-
-  // DDDWorld: get the DDCV from the ES and use it to build the World
-  G4LogicalVolumeToDDLogicalPartMap map_lv;
-  dd4hep::sim::Geant4GeometryMaps::VolumeMap lvMap;
-  if (geoFromDD4hep) {
-    m_world.reset(new DDDWorld(pDD4hep->detector(), lvMap));
-  } else {
-    m_world.reset(new DDDWorld(pDD, map_lv, m_catalog, false));
-  }
-  G4VPhysicalVolume* world = m_world.get()->GetWorldVolume();
-
+  bool geoFromDD4hep = m_p.getParameter<bool>("g4GeometryDD4hepSource");
   bool cuts = m_pPhysics.getParameter<bool>("CutsPerRegion");
+  bool protonCut = m_pPhysics.getUntrackedParameter<bool>("CutsOnProton", true);
   int verb =
       std::max(m_pPhysics.getUntrackedParameter<int>("Verbosity", 0), m_p.getParameter<int>("SteppingVerbosity"));
+  edm::LogVerbatim("SimG4CoreApplication")
+      << "RunManagerMT: start initialising of geometry DD4Hep: " << geoFromDD4hep << "\n"
+      << "              cutsPerRegion: " << cuts << " cutForProton: " << protonCut << "\n"
+      << "              G4 verbosity: " << verb;
+
+  m_world.reset(new DDDWorld(pDD, pDD4hep, m_catalog, verb, cuts, protonCut));
+  G4VPhysicalVolume* world = m_world.get()->GetWorldVolume();
+
   m_kernel->SetVerboseLevel(verb);
   edm::LogVerbatim("SimG4CoreApplication")
       << "RunManagerMT: Define cuts: " << cuts << " Geant4 run manager verbosity: " << verb;
 
-  if (cuts) {
-    DDG4ProductionCuts pcuts(map_lv, verb, m_pPhysics);
-  }
   const G4RegionStore* regStore = G4RegionStore::GetInstance();
   const G4PhysicalVolumeStore* pvs = G4PhysicalVolumeStore::GetInstance();
   const G4LogicalVolumeStore* lvs = G4LogicalVolumeStore::GetInstance();

SimG4Core/Application/src/RunManagerMTWorker.cc

@@ -32,8 +32,6 @@
 #include "FWCore/Framework/interface/ESHandle.h"
 #include "FWCore/Framework/interface/ESTransientHandle.h"
 #include "Geometry/Records/interface/IdealGeometryRecord.h"
-#include "DetectorDescription/Core/interface/DDCompactView.h"
-#include "DetectorDescription/DDCMS/interface/DDCompactView.h"
 
 #include "SimG4Core/Geometry/interface/DDDWorld.h"
 #include "SimG4Core/MagneticField/interface/FieldBuilder.h"
@@ -66,7 +64,6 @@
 #include <thread>
 #include <sstream>
 #include <vector>
-//#include <mutex>
 
 static std::once_flag applyOnce;
 thread_local bool RunManagerMTWorker::dumpMF = false;
@@ -289,17 +286,6 @@ void RunManagerMTWorker::initializeThread(RunManagerMT& runManagerMaster, const
   // we need the track manager now
   m_tls->trackManager.reset(new SimTrackManager());
 
-  // Get DDCompactView, or would it be better to get the object from
-  // runManagerMaster instead of EventSetup in here?
-  auto geoFromDD4hep = m_p.getParameter<bool>("g4GeometryDD4hepSource");
-  edm::ESTransientHandle<cms::DDCompactView> pDD4hep;
-  edm::ESTransientHandle<DDCompactView> pDD;
-  if (geoFromDD4hep) {
-    es.get<IdealGeometryRecord>().get(pDD4hep);
-  } else {
-    es.get<IdealGeometryRecord>().get(pDD);
-  }
-
   // setup the magnetic field
   if (m_pUseMagneticField) {
     const GlobalPoint g(0., 0., 0.);

SimG4Core/Geometry/interface/DDDWorld.h

@@ -4,19 +4,22 @@
 #include "SimG4Core/Geometry/interface/G4LogicalVolumeToDDLogicalPartMap.h"
 #include "SimG4Core/Geometry/interface/SensitiveDetectorCatalog.h"
 #include "G4VPhysicalVolume.hh"
-#include "DDG4/Geant4GeometryInfo.h"
 
-class DDG4Builder;
 class DDCompactView;
 
 namespace cms {
-  class DDDetector;
-}  // namespace cms
+  class DDCompactView;
+}
 
 class DDDWorld {
 public:
+  DDDWorld(const DDCompactView *pDD,
+           const cms::DDCompactView *pDD4hep,
+           SensitiveDetectorCatalog &,
+           int verb,
+           bool cuts,
+           bool pcut);
   DDDWorld(const DDCompactView *, G4LogicalVolumeToDDLogicalPartMap &, SensitiveDetectorCatalog &, bool check);
-  DDDWorld(const cms::DDDetector *, dd4hep::sim::Geant4GeometryMaps::VolumeMap &);
   ~DDDWorld();
   G4VPhysicalVolume *GetWorldVolume() const { return m_world; }
 

SimG4Core/Geometry/interface/DDG4ProductionCuts.h

@@ -1,7 +1,6 @@
 #ifndef SimG4Core_DDG4ProductionCuts_H
 #define SimG4Core_DDG4ProductionCuts_H
 
-#include "FWCore/ParameterSet/interface/ParameterSet.h"
 #include "SimG4Core/Geometry/interface/G4LogicalVolumeToDDLogicalPartMap.h"
 
 #include <string>
@@ -10,23 +9,21 @@
 class DDLogicalPart;
 class G4Region;
 class G4LogicalVolume;
-class G4ProductionCuts;
 
 class DDG4ProductionCuts {
 public:
-  DDG4ProductionCuts(const G4LogicalVolumeToDDLogicalPartMap&, int, const edm::ParameterSet&);
+  explicit DDG4ProductionCuts(const G4LogicalVolumeToDDLogicalPartMap&, int, bool);
   ~DDG4ProductionCuts();
 
 private:
   void initialize();
   void setProdCuts(const DDLogicalPart, G4Region*);
 
   const G4LogicalVolumeToDDLogicalPartMap& map_;
-  // Legacy flag
-  bool protonCut_;
-  std::string keywordRegion_;
-  int verbosity_;
   G4LogicalVolumeToDDLogicalPartMap::Vector vec_;
+  const std::string keywordRegion_;
+  const int verbosity_;
+  const bool protonCut_;
 };
 
 #endif