CMS 3D CMS Logo

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