Belle II Software  release-05-01-25
Belle2PhysicsList.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(c) 2018 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributor: Dennis Wright (SLAC) *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <simulation/physicslist/Belle2PhysicsList.h>
12 #include "G4SystemOfUnits.hh"
13 #include "G4RegionStore.hh"
14 #include "G4ProductionCuts.hh"
15 
16 // EM and decay physics
17 #include "G4EmStandardPhysics.hh"
18 #include "G4EmStandardPhysics_option1.hh"
19 #include "G4OpticalPhysics.hh"
20 #include "G4DecayPhysics.hh"
21 #include <simulation/physicslist/Geant4ePhysics.h>
22 
23 // Hadronic physics
24 #include <simulation/physicslist/ProtonPhysics.h>
25 #include <simulation/physicslist/NeutronPhysics.h>
26 #include <simulation/physicslist/PionPhysics.h>
27 #include <simulation/physicslist/KaonPhysics.h>
28 #include <simulation/physicslist/HyperonPhysics.h>
29 #include <simulation/physicslist/AntiBaryonPhysics.h>
30 #include <simulation/physicslist/IonPhysics.h>
31 #include <simulation/physicslist/GammaLeptoNuclearPhysics.h>
32 
33 // Particles
34 #include "G4BosonConstructor.hh"
35 #include "G4LeptonConstructor.hh"
36 #include "G4MesonConstructor.hh"
37 #include "G4BaryonConstructor.hh"
38 #include "G4IonConstructor.hh"
39 #include "G4ShortLivedConstructor.hh"
40 
41 // Long lived
42 #include "simulation/longlivedneutral/G4LongLivedNeutralPhysics.h"
43 #include <framework/logging/Logger.h>
44 
45 #define g4ePDGcode 0
46 
47 using namespace Belle2;
48 using namespace Simulation;
49 
50 
51 Belle2PhysicsList::Belle2PhysicsList(const G4String& physicsListName)
52  : G4VModularPhysicsList(), m_globalCutValue(0.07)
53 {
54  B2INFO("Using " << physicsListName << " physics list");
55 
56  m_pxdCutValue = m_globalCutValue;
57  m_svdCutValue = m_globalCutValue;
58  m_cdcCutValue = m_globalCutValue;
59  m_arichtopCutValue = m_globalCutValue;
60  m_eclCutValue = m_globalCutValue;
61  m_klmCutValue = m_globalCutValue;
62 
63  // Decay
64  RegisterPhysics(new G4DecayPhysics());
65 
66  // Hadronic physics
67  RegisterPhysics(new ProtonPhysics());
68  RegisterPhysics(new NeutronPhysics());
69  RegisterPhysics(new PionPhysics());
70  RegisterPhysics(new KaonPhysics());
71  RegisterPhysics(new HyperonPhysics());
72  RegisterPhysics(new AntiBaryonPhysics());
73  RegisterPhysics(new IonPhysics());
74  RegisterPhysics(new GammaLeptoNuclearPhysics());
75  RegisterPhysics(new Geant4ePhysics());
76 }
77 
78 
79 Belle2PhysicsList::~Belle2PhysicsList()
80 {}
81 
82 
83 void Belle2PhysicsList::ConstructParticle()
84 {
85  G4BosonConstructor pBosonConstructor;
86  pBosonConstructor.ConstructParticle();
87 
88  G4LeptonConstructor pLeptonConstructor;
89  pLeptonConstructor.ConstructParticle();
90 
91  G4MesonConstructor pMesonConstructor;
92  pMesonConstructor.ConstructParticle();
93 
94  G4BaryonConstructor pBaryonConstructor;
95  pBaryonConstructor.ConstructParticle();
96 
97  G4IonConstructor pIonConstructor;
98  pIonConstructor.ConstructParticle();
99 
100  G4ShortLivedConstructor pShortLivedConstructor;
101  pShortLivedConstructor.ConstructParticle();
102 
103  ConstructG4eParticles();
104 }
105 
106 
107 void Belle2PhysicsList::ConstructG4eParticles()
108 {
109  static G4ParticleDefinition* g4eParticle = NULL;
110 
111  if (g4eParticle == NULL) {
112  // Bohr Magneton for positron and positive muon
113  G4double muBpositron = 0.5 * CLHEP::eplus * CLHEP::hbar_Planck / (0.51099906 * CLHEP::MeV / CLHEP::c_squared);
114  G4double muBmuon = 0.5 * CLHEP::eplus * CLHEP::hbar_Planck / (0.1056584 * CLHEP::GeV / CLHEP::c_squared);
115  G4double muNucleon = CLHEP::eplus * CLHEP::hbar_Planck / 2. / (CLHEP::proton_mass_c2 / CLHEP::c_squared);
116  // Copied from G4Gamma.cc
117  g4eParticle = new G4ParticleDefinition(
118  "g4e_gamma", 0.0 * CLHEP::MeV, 0.0 * CLHEP::MeV, 0.0,
119  2, -1, -1,
120  0, 0, 0,
121  "gamma", 0, 0, g4ePDGcode,
122  true, 0.0, NULL,
123  false, "photon", g4ePDGcode
124  );
125  // Copied from G4Electron.cc
126  new G4ParticleDefinition(
127  "g4e_e-", 0.51099906 * CLHEP::MeV, 0.0 * CLHEP::MeV, -1.0 * CLHEP::eplus,
128  1, 0, 0,
129  0, 0, 0,
130  "lepton", 1, 0, g4ePDGcode,
131  true, -1.0, NULL,
132  false, "e", g4ePDGcode, -1.0011596521859 * muBpositron
133  );
134  // Copied from G4Positron.cc
135  new G4ParticleDefinition(
136  "g4e_e+", 0.51099906 * CLHEP::MeV, 0.0 * CLHEP::MeV, +1.0 * CLHEP::eplus,
137  1, 0, 0,
138  0, 0, 0,
139  "lepton", -1, 0, g4ePDGcode,
140  true, -1.0, NULL,
141  false, "e", g4ePDGcode, 1.0011596521859 * muBpositron
142  );
143  // Copied from G4MuonPlus.cc
144  new G4ParticleDefinition(
145  "g4e_mu+", 0.1056584 * CLHEP::GeV, 2.99591e-16 * CLHEP::MeV, +1.0 * CLHEP::eplus,
146  1, 0, 0,
147  0, 0, 0,
148  "lepton", -1, 0, g4ePDGcode,
149  true, -1.0, NULL,
150  false, "mu", g4ePDGcode, 1.0011659208 * muBmuon
151  );
152  // Copied from G4MuonMinus.cc
153  new G4ParticleDefinition(
154  "g4e_mu-", 0.1056584 * CLHEP::GeV, 2.99591e-16 * CLHEP::MeV, -1.0 * CLHEP::eplus,
155  1, 0, 0,
156  0, 0, 0,
157  "lepton", 1, 0, g4ePDGcode,
158  true, -1.0, NULL,
159  false, "mu", g4ePDGcode, -1.0011659208 * muBmuon
160  );
161  // Copied from G4PionPlus.cc
162  new G4ParticleDefinition(
163  "g4e_pi+", 0.1395700 * CLHEP::GeV, 2.5284e-14 * CLHEP::MeV, +1.0 * CLHEP::eplus,
164  0, -1, 0,
165  2, +2, -1,
166  "meson", 0, 0, g4ePDGcode,
167  true, -1.0, NULL,
168  false, "pi", g4ePDGcode
169  );
170  // Copied from G4PionMinus.cc
171  new G4ParticleDefinition(
172  "g4e_pi-", 0.1395700 * CLHEP::GeV, 2.5284e-14 * CLHEP::MeV, -1.0 * CLHEP::eplus,
173  0, -1, 0,
174  2, -2, -1,
175  "meson", 0, 0, g4ePDGcode,
176  true, -1.0, NULL,
177  false, "pi", g4ePDGcode
178  );
179  // Copied from G4KaonPlus.cc
180  new G4ParticleDefinition(
181  "g4e_kaon+", 0.493677 * CLHEP::GeV, 5.315e-14 * CLHEP::MeV, +1.0 * CLHEP::eplus,
182  0, -1, 0,
183  1, +1, 0,
184  "meson", 0, 0, g4ePDGcode,
185  true, -1.0, NULL,
186  false, "kaon", g4ePDGcode
187  );
188  // Copied from G4KaonMinus.cc
189  new G4ParticleDefinition(
190  "g4e_kaon-", 0.493677 * CLHEP::GeV, 5.315e-14 * CLHEP::MeV, -1.0 * CLHEP::eplus,
191  0, -1, 0,
192  1, -1, 0,
193  "meson", 0, 0, g4ePDGcode,
194  true, -1.0, NULL,
195  false, "kaon", g4ePDGcode
196  );
197  // Copied from G4Proton.cc except use G4ParticleDefinition instead of G4Ions
198  new G4ParticleDefinition(
199  "g4e_proton", 0.9382723 * CLHEP::GeV, 0.0 * CLHEP::MeV, +1.0 * CLHEP::eplus,
200  1, +1, 0,
201  1, +1, 0,
202  "baryon", 0, +1, g4ePDGcode,
203  true, -1.0, NULL,
204  false, "nucleon", g4ePDGcode, 2.792847351 * muNucleon
205  );
206  // Copied from G4AntiProton.cc except use G4ParticleDefinition instead of G4Ions
207  new G4ParticleDefinition(
208  "g4e_anti_proton", 0.9382723 * CLHEP::GeV, 0.0 * CLHEP::MeV, -1.0 * CLHEP::eplus,
209  1, +1, 0,
210  1, -1, 0,
211  "baryon", 0, -1, g4ePDGcode,
212  true, -1.0, NULL,
213  false, "nucleon", g4ePDGcode, -2.792847351 * muNucleon
214  );
215  // copied from G4Deuteron.hh except use G4ParticleDefinition instead of G4Ions
216  new G4ParticleDefinition(
217  "g4e_deuteron", 1.875613 * CLHEP::GeV, 0.0 * CLHEP::MeV, +1.0 * CLHEP::eplus,
218  2, +1, 0,
219  0, 0, 0,
220  "nucleus", 0, +2, g4ePDGcode,
221  true, -1.0, NULL,
222  false, "static", g4ePDGcode, 0.857438230 * muNucleon
223  );
224  // copied from G4AntiDeuteron.hh except use G4ParticleDefinition instead of G4Ions
225  new G4ParticleDefinition(
226  "g4e_anti_deuteron", 1.875613 * CLHEP::GeV, 0.0 * CLHEP::MeV, -1.0 * CLHEP::eplus,
227  2, +1, 0,
228  0, 0, 0,
229  "anti_nucleus", 0, -2, g4ePDGcode,
230  true, -1.0, NULL,
231  false, "static", g4ePDGcode, -0.857438230 * muNucleon
232  );
233  }
234 }
235 
236 void Belle2PhysicsList::setRegionCuts(const std::string& name, const std::vector<std::string>& regions, double cutValue)
237 {
238  G4RegionStore* theRegionStore = G4RegionStore::GetInstance();
239  if (cutValue == 0) {
240  cutValue = m_globalCutValue;
241  } else {
242  B2INFO("Set production cut for detector region" << LogVar("detector", name) << LogVar("production_cut", cutValue));
243  }
244  auto* regionCuts = new G4ProductionCuts;
245  regionCuts->SetProductionCut(cutValue * cm);
246  bool foundOne{false};
247  for (const auto& regionName : regions) {
248  auto* region = theRegionStore->GetRegion(regionName, false);
249  if (!region) {
250  B2WARNING("Cannot find Geant4 region for sub detector. Probably detector not present?"
251  << LogVar("detector", name) << LogVar("region", regionName));
252  continue;
253  }
254  region->SetProductionCuts(regionCuts);
255  foundOne = true;
256  }
257  if (!foundOne) delete regionCuts;
258 }
259 
260 
261 void Belle2PhysicsList::SetCuts()
262 {
263  // Belle2 assumes input units are cm
264  B2INFO("Setting global Geant4 production cut" << LogVar("cutValue", m_globalCutValue));
265 
266  SetCutValue(m_globalCutValue * cm, "proton");
267  SetCutValue(m_globalCutValue * cm, "e-");
268  SetCutValue(m_globalCutValue * cm, "e+");
269  SetCutValue(m_globalCutValue * cm, "gamma");
270 
271  setRegionCuts("PXD", {"PXDEnvelope"}, m_pxdCutValue);
272  setRegionCuts("SVD", {"SVDEnvelope"}, m_svdCutValue);
273  setRegionCuts("CDC", {"CDCEnvelope"}, m_cdcCutValue);
274  setRegionCuts("ARICH", {"ARICHEnvelope"}, m_arichtopCutValue);
275  setRegionCuts("TOP", {"TOPEnvelope"}, m_arichtopCutValue);
276  setRegionCuts("ECL", {"ECLForwardEnvelope", "ECLBarrelSector", "ECLBackwardEnvelope"}, m_eclCutValue);
277  setRegionCuts("KLM", {"BKLMEnvelope", "EKLMEnvelope"}, m_klmCutValue);
278 }
279 
280 
281 void Belle2PhysicsList::SetVerbosity(G4int verb)
282 {
283  SetVerboseLevel(verb);
284 }
285 
286 
287 void Belle2PhysicsList::SetProductionCutValue(G4double value)
288 {
289  m_globalCutValue = value;
290 }
291 
292 
293 void Belle2PhysicsList::SetPXDProductionCutValue(G4double value)
294 {
295  m_pxdCutValue = value;
296 }
297 
298 
299 void Belle2PhysicsList::SetSVDProductionCutValue(G4double value)
300 {
301  m_svdCutValue = value;
302 }
303 
304 
305 void Belle2PhysicsList::SetCDCProductionCutValue(G4double value)
306 {
307  m_cdcCutValue = value;
308 }
309 
310 
311 void Belle2PhysicsList::SetARICHTOPProductionCutValue(G4double value)
312 {
313  m_arichtopCutValue = value;
314 }
315 
316 
317 void Belle2PhysicsList::SetECLProductionCutValue(G4double value)
318 {
319  m_eclCutValue = value;
320 }
321 
322 
323 void Belle2PhysicsList::SetKLMProductionCutValue(G4double value)
324 {
325  m_klmCutValue = value;
326 }
327 
328 
329 void Belle2PhysicsList::UseStandardEMPhysics(G4bool yesno)
330 {
331  if (yesno) {
332  RegisterPhysics(new G4EmStandardPhysics());
333  } else {
334  RegisterPhysics(new G4EmStandardPhysics_option1());
335  }
336 }
337 
338 
339 void Belle2PhysicsList::UseOpticalPhysics(G4bool yesno)
340 {
341  if (yesno) RegisterPhysics(new G4OpticalPhysics());
342 }
343 
344 
345 void Belle2PhysicsList::UseHighPrecisionNeutrons(G4bool yesno)
346 {
347  if (yesno) G4cout << " High precision neutron option not yet ready " << G4endl;
348 }
349 
350 void Belle2PhysicsList::UseLongLivedNeutralParticles()
351 {
352  G4LongLivedNeutralPhysics* pLongLivedNeutral = new G4LongLivedNeutralPhysics();
353  RegisterPhysics(pLongLivedNeutral);
354  pLongLivedNeutral->ConstructParticle();
355 }
Belle2::Simulation::NeutronPhysics
Neutron hadronic physics constructor for Belle II physics list.
Definition: NeutronPhysics.h:33
Belle2::Simulation::Belle2PhysicsList::m_klmCutValue
G4double m_klmCutValue
threshold for BKLM and EKLM
Definition: Belle2PhysicsList.h:106
Belle2::Simulation::Belle2PhysicsList::SetKLMProductionCutValue
void SetKLMProductionCutValue(G4double)
Set cut value for BKLM and EKLM envelopes.
Definition: Belle2PhysicsList.cc:323
Belle2::Simulation::Belle2PhysicsList::SetECLProductionCutValue
void SetECLProductionCutValue(G4double)
Set cut value for ECL barrel, forward and backward envelopes.
Definition: Belle2PhysicsList.cc:317
Belle2::Simulation::Belle2PhysicsList::m_svdCutValue
G4double m_svdCutValue
threshold for SVD
Definition: Belle2PhysicsList.h:102
Belle2::Simulation::IonPhysics
Ion hadronic physics constructor for Belle II physics list.
Definition: IonPhysics.h:35
Belle2::Simulation::Belle2PhysicsList::UseHighPrecisionNeutrons
void UseHighPrecisionNeutrons(G4bool)
Use high precision neutron models below 20 MeV.
Definition: Belle2PhysicsList.cc:345
Belle2::Simulation::Belle2PhysicsList::m_pxdCutValue
G4double m_pxdCutValue
threshold for PXD
Definition: Belle2PhysicsList.h:101
Belle2::Simulation::Belle2PhysicsList::SetProductionCutValue
void SetProductionCutValue(G4double)
Use parameter to set global cut value.
Definition: Belle2PhysicsList.cc:287
Belle2::Simulation::KaonPhysics
Kaon hadronic physics constructor for Belle II physics list.
Definition: KaonPhysics.h:33
Belle2::Simulation::Belle2PhysicsList::UseStandardEMPhysics
void UseStandardEMPhysics(G4bool)
Use standard EM physics instead of EM option1.
Definition: Belle2PhysicsList.cc:329
Belle2::Simulation::Belle2PhysicsList::UseOpticalPhysics
void UseOpticalPhysics(G4bool)
Add optical photon physics.
Definition: Belle2PhysicsList.cc:339
Belle2::Simulation::Belle2PhysicsList::setRegionCuts
void setRegionCuts(const std::string &name, const std::vector< std::string > &regions, double cutValue)
Set the produciton cuts to the given value for a list of regions belonging to a sub detector.
Definition: Belle2PhysicsList.cc:236
Belle2::Simulation::Belle2PhysicsList::Belle2PhysicsList
Belle2PhysicsList(const G4String &physicsListName)
Constructor.
Definition: Belle2PhysicsList.cc:51
Belle2::Simulation::Belle2PhysicsList::m_cdcCutValue
G4double m_cdcCutValue
threshold for CDC
Definition: Belle2PhysicsList.h:103
Belle2::Simulation::Belle2PhysicsList::m_globalCutValue
G4double m_globalCutValue
Secondary production thresholds.
Definition: Belle2PhysicsList.h:99
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::Simulation::Belle2PhysicsList::ConstructG4eParticles
void ConstructG4eParticles()
Construct parallel particle types needed for reco.
Definition: Belle2PhysicsList.cc:107
Belle2::Simulation::Belle2PhysicsList::m_arichtopCutValue
G4double m_arichtopCutValue
threshold for ARICH and TOP
Definition: Belle2PhysicsList.h:104
Belle2::Simulation::Belle2PhysicsList::m_eclCutValue
G4double m_eclCutValue
threshold for ECL
Definition: Belle2PhysicsList.h:105
LogVar
Class to store variables with their name which were sent to the logging service.
Definition: LogVariableStream.h:24
Belle2::Simulation::PionPhysics
Pion hadronic physics constructor for Belle II physics list.
Definition: PionPhysics.h:33
Belle2::Simulation::Geant4ePhysics
Define geant4e-specific physics.
Definition: Geant4ePhysics.h:41
Belle2::Simulation::Belle2PhysicsList::SetARICHTOPProductionCutValue
void SetARICHTOPProductionCutValue(G4double)
Set cut value for ARICH and TOP envelopes.
Definition: Belle2PhysicsList.cc:311
Belle2::Simulation::Belle2PhysicsList::ConstructParticle
virtual void ConstructParticle()
Build all particle types used in physics list.
Definition: Belle2PhysicsList.cc:83
Belle2::Simulation::Belle2PhysicsList::UseLongLivedNeutralParticles
void UseLongLivedNeutralParticles()
Simulate neutral long-lived particles with given pdg and mass value.
Definition: Belle2PhysicsList.cc:350
Belle2::Simulation::GammaLeptoNuclearPhysics
Gamma-nuclear, electro-nuclear and muon-nuclear physics constructor for Belle II physics list.
Definition: GammaLeptoNuclearPhysics.h:35
Belle2::Simulation::Belle2PhysicsList::SetCuts
virtual void SetCuts()
Set the secondary particle production thresholds.
Definition: Belle2PhysicsList.cc:261
Belle2::Simulation::HyperonPhysics
Hyperon hadronic physics constructor for Belle II physics list.
Definition: HyperonPhysics.h:33
Belle2::Simulation::Belle2PhysicsList::SetVerbosity
void SetVerbosity(G4int verb)
Run/event verbosity level.
Definition: Belle2PhysicsList.cc:281
Belle2::G4LongLivedNeutralPhysics::ConstructParticle
virtual void ConstructParticle()
Adds monopole and anti-monopole to GEANT4 with a pdg of +/-99666 and parameters taken from current cl...
Definition: G4LongLivedNeutralPhysics.cc:37
Belle2::Simulation::Belle2PhysicsList::SetPXDProductionCutValue
void SetPXDProductionCutValue(G4double)
Set cut value for PXD envelope.
Definition: Belle2PhysicsList.cc:293
Belle2::Simulation::AntiBaryonPhysics
Anti-baryon hadronic physics constructor for Belle II physics list.
Definition: AntiBaryonPhysics.h:34
Belle2::Simulation::Belle2PhysicsList::SetCDCProductionCutValue
void SetCDCProductionCutValue(G4double)
Set cut value for CDC envelope.
Definition: Belle2PhysicsList.cc:305
Belle2::G4LongLivedNeutralPhysics
LongLivedNeutral physics Class – to be registered in the physics list.
Definition: G4LongLivedNeutralPhysics.h:37
Belle2::Simulation::ProtonPhysics
Proton hadronic physics constructor for Belle II physics list.
Definition: ProtonPhysics.h:33
Belle2::Simulation::Belle2PhysicsList::~Belle2PhysicsList
~Belle2PhysicsList()
destructor
Definition: Belle2PhysicsList.cc:79
Belle2::Simulation::Belle2PhysicsList::SetSVDProductionCutValue
void SetSVDProductionCutValue(G4double)
Set cut value for SVD envelope.
Definition: Belle2PhysicsList.cc:299