CMS 3D CMS Logo

CustomPhysicsList.cc
Go to the documentation of this file.
7 
11 
12 #include "G4hMultipleScattering.hh"
13 #include "G4hIonisation.hh"
14 #include "G4ProcessManager.hh"
15 
18 
19 using namespace CLHEP;
20 
21 G4ThreadLocal std::unique_ptr<G4ProcessHelper> CustomPhysicsList::myHelper;
22 
24  const edm::ParameterSet & p, bool apinew)
25  : G4VPhysicsConstructor(name)
26 {
27  myConfig = p;
28  if(apinew) {
29  dfactor = p.getParameter<double>("DarkMPFactor");
30  fHadronicInteraction = p.getParameter<bool>("RhadronPhysics");
31  } else {
32  // this is left for backward compatibility
33  dfactor = p.getParameter<double>("dark_factor");
34  fHadronicInteraction = p.getParameter<bool>("rhadronPhysics");
35  }
36  edm::FileInPath fp = p.getParameter<edm::FileInPath>("particlesDef");
37  particleDefFilePath = fp.fullPath();
39  myHelper.reset(nullptr);
40 
41  edm::LogInfo("SimG4CoreCustomPhysics")
42  << "CustomPhysicsList: Path for custom particle definition file: \n"
43  <<particleDefFilePath << "\n" << " dark_factor= " << dfactor;
44 }
45 
47 }
48 
50  edm::LogInfo("SimG4CoreCustomPhysics")
51  << "===== CustomPhysicsList::ConstructParticle ";
52  fParticleFactory.get()->loadCustomParticles(particleDefFilePath);
53 }
54 
56 
57  edm::LogInfo("SimG4CoreCustomPhysics")
58  <<"CustomPhysicsList: adding CustomPhysics processes "
59  << "for the list of particles";
60 
61  G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();
62 
63  for(auto particle : fParticleFactory.get()->GetCustomParticles()) {
64  CustomParticle* cp = dynamic_cast<CustomParticle*>(particle);
65  if(cp) {
66  G4ProcessManager* pmanager = particle->GetProcessManager();
67  edm::LogInfo("SimG4CoreCustomPhysics")
68  <<"CustomPhysicsList: " << particle->GetParticleName()
69  <<" PDGcode= " << particle->GetPDGEncoding()
70  << " Mass= " << particle->GetPDGMass()/GeV <<" GeV.";
71  if(pmanager) {
72  if(particle->GetPDGCharge() != 0.0) {
73  ph->RegisterProcess(new G4hMultipleScattering, particle);
74  ph->RegisterProcess(new G4hIonisation, particle);
75  }
76  if(cp->GetCloud() && fHadronicInteraction &&
77  CustomPDGParser::s_isRHadron(particle->GetPDGEncoding())) {
78  edm::LogInfo("SimG4CoreCustomPhysics")
79  <<"CustomPhysicsList: " << particle->GetParticleName()
80  <<" CloudMass= " <<cp->GetCloud()->GetPDGMass()/GeV
81  <<" GeV; SpectatorMass= " << cp->GetSpectator()->GetPDGMass()/GeV<<" GeV.";
82 
83  if(!myHelper.get()) { myHelper.reset(new G4ProcessHelper(myConfig, fParticleFactory.get())); }
84  pmanager->AddDiscreteProcess(new FullModelHadronicProcess(myHelper.get()));
85  }
86  if(particle->GetParticleType()=="darkpho"){
88  pmanager->AddDiscreteProcess(darkGamma);
89  }
90  }
91  }
92  }
93 }
T getParameter(std::string const &) const
CustomPhysicsList(const std::string &name, const edm::ParameterSet &p, bool useuni=false)
const double GeV
Definition: MathUtil.h:16
void ConstructProcess() override
edm::ParameterSet myConfig
std::string particleDefFilePath
G4ParticleDefinition * GetCloud()
static G4ThreadLocal std::unique_ptr< G4ProcessHelper > myHelper
void ConstructParticle() override
G4ParticleDefinition * GetSpectator()
static bool s_isRHadron(int pdg)
std::unique_ptr< CustomParticleFactory > fParticleFactory
~CustomPhysicsList() override