From 9edb03766b96659f6d9d7ba1153fa77d32d950e2 Mon Sep 17 00:00:00 2001 From: Jonas Hahnfeld Date: Fri, 14 Oct 2022 17:39:26 +0200 Subject: [PATCH] Add G4HepEm to CMSSW The new EMH physics list replaces the EM processes for electrons, positrons, and gammas with G4HepEm. Note that gamma-lepto-nuclear interactions are NOT implemented yet, and disabled for this physics list (by removing G4EmExtraPhysics). --- BigProducts/Simulation/BuildFile.xml | 2 + SimG4Core/PhysicsLists/BuildFile.xml | 1 + .../interface/CMSEmStandardPhysicsEMH.h | 26 +++++ .../PhysicsLists/plugins/FTFPCMS_BERT_EMH.cc | 44 +++++++++ .../PhysicsLists/plugins/FTFPCMS_BERT_EMH.h | 12 +++ SimG4Core/PhysicsLists/plugins/module.cc | 3 + .../src/CMSEmStandardPhysicsEMH.cc | 94 +++++++++++++++++++ 7 files changed, 182 insertions(+) create mode 100644 SimG4Core/PhysicsLists/interface/CMSEmStandardPhysicsEMH.h create mode 100644 SimG4Core/PhysicsLists/plugins/FTFPCMS_BERT_EMH.cc create mode 100644 SimG4Core/PhysicsLists/plugins/FTFPCMS_BERT_EMH.h create mode 100644 SimG4Core/PhysicsLists/src/CMSEmStandardPhysicsEMH.cc diff --git a/BigProducts/Simulation/BuildFile.xml b/BigProducts/Simulation/BuildFile.xml index fbc1d9055af1d..351008a7368d8 100644 --- a/BigProducts/Simulation/BuildFile.xml +++ b/BigProducts/Simulation/BuildFile.xml @@ -31,7 +31,9 @@ + + diff --git a/SimG4Core/PhysicsLists/BuildFile.xml b/SimG4Core/PhysicsLists/BuildFile.xml index 6594fc7881be7..57463ab9eb23f 100644 --- a/SimG4Core/PhysicsLists/BuildFile.xml +++ b/SimG4Core/PhysicsLists/BuildFile.xml @@ -1,6 +1,7 @@ + diff --git a/SimG4Core/PhysicsLists/interface/CMSEmStandardPhysicsEMH.h b/SimG4Core/PhysicsLists/interface/CMSEmStandardPhysicsEMH.h new file mode 100644 index 0000000000000..81b461da6da8a --- /dev/null +++ b/SimG4Core/PhysicsLists/interface/CMSEmStandardPhysicsEMH.h @@ -0,0 +1,26 @@ +#ifndef SimG4Core_PhysicsLists_CMSEmStandardPhysicsEMH_h +#define SimG4Core_PhysicsLists_CMSEmStandardPhysicsEMH_h + +#include "G4VPhysicsConstructor.hh" +#include "globals.hh" +#include "G4MscStepLimitType.hh" + +#include "FWCore/ParameterSet/interface/ParameterSet.h" + +class CMSEmStandardPhysicsEMH : public G4VPhysicsConstructor { +public: + CMSEmStandardPhysicsEMH(G4int ver, const edm::ParameterSet& p); + ~CMSEmStandardPhysicsEMH() override; + + void ConstructParticle() override; + void ConstructProcess() override; + +private: + G4double fRangeFactor; + G4double fGeomFactor; + G4double fSafetyFactor; + G4double fLambdaLimit; + G4MscStepLimitType fStepLimitType; +}; + +#endif diff --git a/SimG4Core/PhysicsLists/plugins/FTFPCMS_BERT_EMH.cc b/SimG4Core/PhysicsLists/plugins/FTFPCMS_BERT_EMH.cc new file mode 100644 index 0000000000000..5374d0d906cf7 --- /dev/null +++ b/SimG4Core/PhysicsLists/plugins/FTFPCMS_BERT_EMH.cc @@ -0,0 +1,44 @@ +#include "FTFPCMS_BERT_EMH.h" +#include "SimG4Core/PhysicsLists/interface/CMSEmStandardPhysicsEMH.h" +#include "SimG4Core/PhysicsLists/interface/CMSHadronPhysicsFTFP_BERT.h" +#include "FWCore/MessageLogger/interface/MessageLogger.h" + +#include "G4DecayPhysics.hh" +#include "G4IonPhysics.hh" +#include "G4StoppingPhysics.hh" +#include "G4HadronElasticPhysics.hh" + +FTFPCMS_BERT_EMH::FTFPCMS_BERT_EMH(const edm::ParameterSet& p) : PhysicsList(p) { + int ver = p.getUntrackedParameter("Verbosity", 0); + bool emPhys = p.getUntrackedParameter("EMPhysics", true); + bool hadPhys = p.getUntrackedParameter("HadPhysics", true); + double minFTFP = p.getParameter("EminFTFP") * CLHEP::GeV; + double maxBERT = p.getParameter("EmaxBERT") * CLHEP::GeV; + double maxBERTpi = p.getParameter("EmaxBERTpi") * CLHEP::GeV; + edm::LogVerbatim("PhysicsList") << "CMS Physics List FTFP_BERT_EMH: " + << "\n Flags for EM Physics: " << emPhys << "; Hadronic Physics: " << hadPhys + << "\n Transition energy Bertini/FTFP from " << minFTFP / CLHEP::GeV << " to " + << maxBERT / CLHEP::GeV << "; for pions to " << maxBERTpi / CLHEP::GeV << " GeV"; + + if (emPhys) { + // EM Physics + RegisterPhysics(new CMSEmStandardPhysicsEMH(ver, p)); + } + + // Decays + this->RegisterPhysics(new G4DecayPhysics(ver)); + + if (hadPhys) { + // Hadron Elastic scattering + RegisterPhysics(new G4HadronElasticPhysics(ver)); + + // Hadron Physics + RegisterPhysics(new CMSHadronPhysicsFTFP_BERT(minFTFP, maxBERT, maxBERTpi, minFTFP, maxBERT)); + + // Stopping Physics + RegisterPhysics(new G4StoppingPhysics(ver)); + + // Ion Physics + RegisterPhysics(new G4IonPhysics(ver)); + } +} diff --git a/SimG4Core/PhysicsLists/plugins/FTFPCMS_BERT_EMH.h b/SimG4Core/PhysicsLists/plugins/FTFPCMS_BERT_EMH.h new file mode 100644 index 0000000000000..51051b811b957 --- /dev/null +++ b/SimG4Core/PhysicsLists/plugins/FTFPCMS_BERT_EMH.h @@ -0,0 +1,12 @@ +#ifndef SimG4Core_PhysicsLists_FTFPCMS_BERT_EMH_H +#define SimG4Core_PhysicsLists_FTFPCMS_BERT_EMH_H + +#include "SimG4Core/Physics/interface/PhysicsList.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" + +class FTFPCMS_BERT_EMH : public PhysicsList { +public: + FTFPCMS_BERT_EMH(const edm::ParameterSet& p); +}; + +#endif diff --git a/SimG4Core/PhysicsLists/plugins/module.cc b/SimG4Core/PhysicsLists/plugins/module.cc index 35d6e61a184e9..c55d9ee5b440e 100644 --- a/SimG4Core/PhysicsLists/plugins/module.cc +++ b/SimG4Core/PhysicsLists/plugins/module.cc @@ -3,6 +3,7 @@ #include "DummyPhysics.h" #include "FTFCMS_BIC.h" #include "FTFPCMS_BERT.h" +#include "FTFPCMS_BERT_EMH.h" #include "FTFPCMS_BERT_EML.h" #include "FTFPCMS_BERT_EMM.h" #include "G4Version.hh" @@ -33,6 +34,8 @@ typedef FTFCMS_BIC FTF_BIC; DEFINE_PHYSICSLIST(FTF_BIC); typedef FTFPCMS_BERT FTFP_BERT; DEFINE_PHYSICSLIST(FTFP_BERT); +typedef FTFPCMS_BERT_EMH FTFP_BERT_EMH; +DEFINE_PHYSICSLIST(FTFP_BERT_EMH); typedef FTFPCMS_BERT_EML FTFP_BERT_EML; DEFINE_PHYSICSLIST(FTFP_BERT_EML); typedef FTFPCMS_BERT_EMM FTFP_BERT_EMM; diff --git a/SimG4Core/PhysicsLists/src/CMSEmStandardPhysicsEMH.cc b/SimG4Core/PhysicsLists/src/CMSEmStandardPhysicsEMH.cc new file mode 100644 index 0000000000000..6690c0a57bf24 --- /dev/null +++ b/SimG4Core/PhysicsLists/src/CMSEmStandardPhysicsEMH.cc @@ -0,0 +1,94 @@ +#include "SimG4Core/PhysicsLists/interface/CMSEmStandardPhysicsEMH.h" +#include "FWCore/MessageLogger/interface/MessageLogger.h" + +#include "G4SystemOfUnits.hh" +#include "G4ParticleDefinition.hh" +#include "G4EmParameters.hh" +#include "G4EmBuilder.hh" + +#include "G4MscStepLimitType.hh" + +#include "G4hIonisation.hh" +#include "G4hMultipleScattering.hh" +#include "G4ionIonisation.hh" + +#include "G4ParticleTable.hh" +#include "G4Gamma.hh" +#include "G4Electron.hh" +#include "G4Positron.hh" +#include "G4GenericIon.hh" + +#include "G4PhysicsListHelper.hh" +#include "G4BuilderType.hh" +#include "G4ProcessManager.hh" + +#include "G4HepEmProcess.hh" + +#include + +CMSEmStandardPhysicsEMH::CMSEmStandardPhysicsEMH(G4int ver, const edm::ParameterSet& p) + : G4VPhysicsConstructor("CMSEmStandard_emh") { + SetVerboseLevel(ver); + G4EmParameters* param = G4EmParameters::Instance(); + param->SetDefaults(); + param->SetVerbose(ver); + param->SetApplyCuts(true); + param->SetStepFunction(0.8, 1 * CLHEP::mm); + param->SetMscRangeFactor(0.2); + param->SetMscStepLimitType(fUseSafety); + param->SetFluo(false); + SetPhysicsType(bElectromagnetic); + fRangeFactor = p.getParameter("G4MscRangeFactor"); + fGeomFactor = p.getParameter("G4MscGeomFactor"); + fSafetyFactor = p.getParameter("G4MscSafetyFactor"); + fLambdaLimit = p.getParameter("G4MscLambdaLimit") * CLHEP::mm; + std::string msc = p.getParameter("G4MscStepLimit"); + fStepLimitType = fUseSafety; + if (msc == "UseSafetyPlus") { + fStepLimitType = fUseSafetyPlus; + } + if (msc == "Minimal") { + fStepLimitType = fMinimal; + } + double tcut = p.getParameter("G4TrackingCut") * CLHEP::MeV; + param->SetLowestElectronEnergy(tcut); + param->SetLowestMuHadEnergy(tcut); +} + +CMSEmStandardPhysicsEMH::~CMSEmStandardPhysicsEMH() {} + +void CMSEmStandardPhysicsEMH::ConstructParticle() { + // minimal set of particles for EM physics + G4EmBuilder::ConstructMinimalEmSet(); +} + +void CMSEmStandardPhysicsEMH::ConstructProcess() { + if (verboseLevel > 0) { + edm::LogVerbatim("PhysicsList") << "### " << GetPhysicsName() << " Construct EM Processes"; + } + + // This EM builder takes default models of Geant4 10 EMV. + // Multiple scattering by WentzelVI for all particles except: + // a) e+e- below 100 MeV for which the Urban model is used + // b) ions for which Urban model is used + G4EmBuilder::PrepareEMPhysics(); + + G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper(); + // processes used by several particles + G4hMultipleScattering* hmsc = new G4hMultipleScattering("ionmsc"); + G4NuclearStopping* pnuc(nullptr); + + G4HepEmProcess* hepEmProcess = new G4HepEmProcess(); + G4Electron::Electron()->GetProcessManager()->AddProcess(hepEmProcess, -1, -1, 1); + G4Positron::Positron()->GetProcessManager()->AddProcess(hepEmProcess, -1, -1, 1); + G4Gamma::Gamma()->GetProcessManager()->AddProcess(hepEmProcess, -1, -1, 1); + + // generic ion + G4ParticleDefinition* particle = G4GenericIon::GenericIon(); + G4ionIonisation* ionIoni = new G4ionIonisation(); + ph->RegisterProcess(hmsc, particle); + ph->RegisterProcess(ionIoni, particle); + + // muons, hadrons ions + G4EmBuilder::ConstructCharged(hmsc, pnuc); +}