extention of neutron background simulation

cms-sw · Mar 8, 2017 · 8f5ccbf · 8f5ccbf
1 parent 1cdb34b
commit 8f5ccbf
Show file tree

Hide file tree

Showing 6 changed files with 80 additions and 57 deletions.
diff --git a/SimG4Core/Application/interface/CMSEmParticleList.h b/SimG4Core/Application/interface/CMSEmParticleList.h
@@ -0,0 +1,24 @@
+#ifndef SimG4Core_Application_CMSEmParticleList_h
+#define SimG4Core_Application_CMSEmParticleList_h 1
+
+// V.Ivanchenko 6 March 2017
+// List of Geant4 basic particle names used in SIM step 
+
+#include "globals.hh"
+
+static const G4int nEmParticles = 39;
+static const G4String EmPartNames[nEmParticles] = 
+{ 
+        "gamma",            "e-",           "e+",           "mu+",        "mu-",
+          "pi+",           "pi-",        "kaon+",         "kaon-",     "proton",
+  "anti_proton",         "alpha",          "He3",    "GenericIon",         "B+",
+           "B-",            "D+",           "D-",           "Ds+",        "Ds-",
+     "anti_He3",    "anti_alpha","anti_deuteron","anti_lambda_c+","anti_omega-",
+"anti_sigma_c+","anti_sigma_c++",  "anti_sigma+",   "anti_sigma-","anti_triton",
+     "sigma_c+",     "sigma_c++",       "sigma+",        "sigma-",       "tau+",
+         "tau-",        "triton",        "xi_c+",           "xi-"
+};
+
+
+#endif
+
diff --git a/SimG4Core/Application/python/NeutronBGforMuonsHP_cff.py b/SimG4Core/Application/python/NeutronBGforMuonsHP_cff.py
@@ -2,35 +2,16 @@
 
 def customise(process):
 
-    # fragment allowing to simulate neutron background in muon system
+  # fragment allowing to simulate neutron background in muon system
+  # using HP neutron package and thermal neutron scattering
+
+  from SimG4Core.Application.NeutronBGforMuons_cff import neutronBG
+
+  process = neutronBG(process)
 
   if hasattr(process,'g4SimHits'):
-    # time window 100 millisecond
-    process.common_maximum_time.MaxTrackTime = cms.double(100000000.0)
-    process.common_maximum_time.DeadRegions = cms.vstring()
-    # Physics List HP
     process.g4SimHits.Physics.type = cms.string('SimG4Core/Physics/FTFP_BERT_HP_EML')
-    process.g4SimHits.Physics.CutsOnProton  = cms.untracked.bool(False)
-    process.g4SimHits.Physics.FlagFluo    = cms.bool(True)
-    process.g4SimHits.Physics.ThermalNeutrons = cms.untracked.bool(False)
-    # Eta cut
-    process.g4SimHits.Generator.MinEtaCut = cms.double(-7.0)
-    process.g4SimHits.Generator.MaxEtaCut = cms.double(7.0)
-    # stacking action
-    process.g4SimHits.StackingAction.MaxTrackTime = cms.double(100000000.0)
-    process.g4SimHits.StackingAction.DeadRegions = cms.vstring()
-    process.g4SimHits.StackingAction.KillHeavy = cms.bool(True)
-    process.g4SimHits.StackingAction.IonThreshold = cms.double(0.0)
-    process.g4SimHits.StackingAction.ProtonThreshold = cms.double(0.0)
-    process.g4SimHits.StackingAction.NeutronThreshold = cms.double(0.0)
-    process.g4SimHits.StackingAction.GammaThreshold = cms.double(0.0)
-    # stepping action
-    process.g4SimHits.SteppingAction.MaxTrackTime = cms.double(100000000.0)
-    process.g4SimHits.SteppingAction.DeadRegions = cms.vstring()
-    # Russian roulette disabled
-    process.g4SimHits.StackingAction.RusRoGammaEnergyLimit = cms.double(0.0)
-    process.g4SimHits.StackingAction.RusRoNeutronEnergyLimit = cms.double(0.0)
-    # full simulation of HF 
-    process.g4SimHits.HFShower.UseHFGflash = cms.bool(False)
+    process.g4SimHits.Physics.ThermalNeutrons = cms.untracked.bool(True)
 
     return(process)
+
diff --git a/SimG4Core/Application/python/NeutronBGforMuonsXS_cff.py b/SimG4Core/Application/python/NeutronBGforMuonsXS_cff.py
@@ -2,31 +2,12 @@
 
 def customise(process):
 
-    # fragment allowing to simulate neutron background in muon system
+  # fragment allowing to simulate neutron background in muon system
+  # using default neutron tracking 
 
-  if hasattr(process,'g4SimHits'):
-    # time window 100 millisecond
-    process.common_maximum_time.MaxTrackTime = cms.double(100000000.0)
-    process.common_maximum_time.DeadRegions = cms.vstring()
-    # Physics List XS
-    process.g4SimHits.Physics.type = cms.string('SimG4Core/Physics/FTFP_BERT_XS_EML')
-    process.g4SimHits.Physics.CutsOnProton  = cms.untracked.bool(True)
-    process.g4SimHits.Physics.FlagFluo    = cms.bool(True)
-    # Eta cut
-    process.g4SimHits.Generator.MinEtaCut = cms.double(-7.0)
-    process.g4SimHits.Generator.MaxEtaCut = cms.double(7.0)
-    # stacking action
-    process.g4SimHits.StackingAction.MaxTrackTime = cms.double(100000000.0)
-    process.g4SimHits.StackingAction.DeadRegions = cms.vstring()
-    process.g4SimHits.StackingAction.KillHeavy = cms.bool(False)
-    process.g4SimHits.StackingAction.GammaThreshold = cms.double(0.0)
-    # stepping action
-    process.g4SimHits.SteppingAction.MaxTrackTime = cms.double(100000000.0)
-    process.g4SimHits.SteppingAction.DeadRegions = cms.vstring()
-    # Russian roulette disabled
-    process.g4SimHits.StackingAction.RusRoGammaEnergyLimit = cms.double(0.0)
-    process.g4SimHits.StackingAction.RusRoNeutronEnergyLimit = cms.double(0.0)
-    # full simulation of HF 
-    process.g4SimHits.HFShower.UseHFGflash = cms.bool(False)
+  from SimG4Core.Application.NeutronBGforMuons_cff import neutronBG
+
+  process = neutronBG(process)
+
+  return(process)
 
-    return(process)
diff --git a/SimG4Core/Application/python/NeutronBGforMuons_cff.py b/SimG4Core/Application/python/NeutronBGforMuons_cff.py
@@ -0,0 +1,33 @@
+import FWCore.ParameterSet.Config as cms
+
+def neutronBG(process):
+
+  # common fragment allowing to simulate neutron background in muon system
+
+  if hasattr(process,'g4SimHits'):
+  # time window 100 millisecond
+    process.common_maximum_time.MaxTrackTime = cms.double(10000000000.0)
+    process.common_maximum_time.DeadRegions = cms.vstring()
+    # Physics List XS
+    process.g4SimHits.Physics.type = cms.string('SimG4Core/Physics/FTFP_BERT_XS_EML')
+    process.g4SimHits.Physics.CutsOnProton  = cms.untracked.bool(True)
+    process.g4SimHits.Physics.FlagFluo    = cms.bool(True)
+    process.g4SimHits.Physics.ThermalNeutrons = cms.untracked.bool(False)
+    # Eta cut
+    process.g4SimHits.Generator.MinEtaCut = cms.double(-7.0)
+    process.g4SimHits.Generator.MaxEtaCut = cms.double(7.0)
+    # stacking action
+    process.g4SimHits.StackingAction.MaxTrackTime = cms.double(10000000000.0)
+    process.g4SimHits.StackingAction.DeadRegions = cms.vstring()
+    process.g4SimHits.StackingAction.KillHeavy = cms.bool(False)
+    process.g4SimHits.StackingAction.GammaThreshold = cms.double(0.0)
+    # stepping action
+    process.g4SimHits.SteppingAction.MaxTrackTime = cms.double(10000000000.0)
+    process.g4SimHits.SteppingAction.DeadRegions = cms.vstring()
+    # Russian roulette disabled
+    process.g4SimHits.StackingAction.RusRoGammaEnergyLimit = cms.double(0.0)
+    process.g4SimHits.StackingAction.RusRoNeutronEnergyLimit = cms.double(0.0)
+    # full simulation of HF 
+    process.g4SimHits.HFShower.UseHFGflash = cms.bool(False)
+
+    return(process)
diff --git a/SimG4Core/PhysicsLists/plugins/FTFPCMS_BERT_XS_EML.cc b/SimG4Core/PhysicsLists/plugins/FTFPCMS_BERT_XS_EML.cc
@@ -1,6 +1,7 @@
 #include "FTFPCMS_BERT_XS_EML.hh"
 #include "SimG4Core/PhysicsLists/interface/CMSEmStandardPhysicsXS.h"
 #include "SimG4Core/PhysicsLists/interface/CMSMonopolePhysics.h"
+#include "SimG4Core/PhysicsLists/interface/CMSThermalNeutrons.h"
 #include "FWCore/MessageLogger/interface/MessageLogger.h"
 
 #include "G4DecayPhysics.hh"
@@ -26,12 +27,14 @@ FTFPCMS_BERT_XS_EML::FTFPCMS_BERT_XS_EML(G4LogicalVolumeToDDLogicalPartMap& map,
   bool emPhys  = p.getUntrackedParameter<bool>("EMPhysics",true);
   bool hadPhys = p.getUntrackedParameter<bool>("HadPhysics",true);
   bool tracking= p.getParameter<bool>("TrackingCut");
+  bool thermal = p.getUntrackedParameter<bool>("ThermalNeutrons");
   double timeLimit = p.getParameter<double>("MaxTrackTime")*ns;
   edm::LogInfo("PhysicsList") << "You are using the simulation engine: "
 			      << "FTFP_BERT_XS_EML \n Flags for EM Physics "
 			      << emPhys << ", for Hadronic Physics "
 			      << hadPhys << " and tracking cut " << tracking
-			      << "   t(ns)= " << timeLimit/ns;
+			      << "   t(ns)= " << timeLimit/ns
+			      << " ThermalNeutrons: " << thermal;
 
   if (emPhys) {
     // EM Physics
@@ -60,14 +63,15 @@ FTFPCMS_BERT_XS_EML::FTFPCMS_BERT_XS_EML(G4LogicalVolumeToDDLogicalPartMap& map,
     // Ion Physics
     RegisterPhysics( new G4IonPhysics(ver));
 
-    RegisterPhysics( new G4NeutronCrossSectionXS(ver));
-
     // Neutron tracking cut
     if (tracking) {
       G4NeutronTrackingCut* ncut= new G4NeutronTrackingCut(ver);
       ncut->SetTimeLimit(timeLimit);
       RegisterPhysics(ncut);
     }
+    if(thermal) {
+      RegisterPhysics(new CMSThermalNeutrons(ver));
+    }
   }
 
   // Monopoles

diff --git a/SimG4Core/PhysicsLists/src/CMSThermalNeutrons.cc b/SimG4Core/PhysicsLists/src/CMSThermalNeutrons.cc
@@ -31,7 +31,7 @@ void CMSThermalNeutrons::ConstructProcess() {
   }
 
   G4int ni = (hpel->GetHadronicInteractionList()).size();
-  if(ni < 2) {
+  if(ni < 1) {
     G4cout << "### " << GetPhysicsName() 
 	   << " WARNING: Fail to add thermal neutron scattering - Nint= " 
 	   << ni << G4endl;