major upgrade: separated TG4ParticlesManager, TG4G3PhysicsManager (now aggregated)
[u/mrichter/AliRoot.git] / TGeant4 / TG4PhysicsManager.cxx
CommitLineData
2817d3e2 1// $Id$
2// Category: physics
3//
4// See the class description in the header file.
5
6#include "TG4PhysicsManager.h"
b2030327 7#include "TG4ParticlesManager.h"
8#include "TG4G3PhysicsManager.h"
9#include "TG4PhysicsConstructorEM.h"
10#include "TG4PhysicsConstructorOptical.h"
11#include "TG4PhysicsConstructorHadron.h"
12#include "TG4PhysicsConstructorSpecialCuts.h"
13#include "TG4PhysicsConstructorSpecialControls.h"
14#include "TG4G3Cut.h"
15#include "TG4G3Control.h"
2817d3e2 16
17#include <G4ParticleDefinition.hh>
094a5851 18#include <G4VProcess.hh>
b2030327 19#include <G4VModularPhysicsList.hh>
2817d3e2 20
21#include <TDatabasePDG.h>
22
23TG4PhysicsManager* TG4PhysicsManager::fgInstance = 0;
24
b2030327 25TG4PhysicsManager::TG4PhysicsManager(G4VModularPhysicsList* physicsList)
26 : fPhysicsList(physicsList),
27 fSetEMPhysics(true),
28 fSetOpticalPhysics(false),
29 fSetHadronPhysics(false),
30 fSetSpecialCutsPhysics(false),
31 fSetSpecialControlsPhysics(false)
32
2817d3e2 33{
34//
35 if (fgInstance) {
36 TG4Globals::Exception(
37 "TG4PhysicsManager: attempt to create two instances of singleton.");
38 }
39
40 fgInstance = this;
b2030327 41
42 // create particles manager
43 fParticlesManager = new TG4ParticlesManager();
2817d3e2 44
b2030327 45 // create G3 physics manager
46 fG3PhysicsManager = new TG4G3PhysicsManager();
2817d3e2 47
0453c41f 48 // fill process name map
49 FillProcessMap();
2817d3e2 50}
51
b2030327 52TG4PhysicsManager::TG4PhysicsManager(){
53//
54 delete fParticlesManager;
55 delete fG3PhysicsManager;
56}
57
2817d3e2 58TG4PhysicsManager::TG4PhysicsManager(const TG4PhysicsManager& right) {
59//
60 TG4Globals::Exception(
61 "Attempt to copy TG4PhysicsManager singleton.");
62}
63
64TG4PhysicsManager::~TG4PhysicsManager() {
65//
2817d3e2 66}
67
68// operators
69
70TG4PhysicsManager&
71TG4PhysicsManager::operator=(const TG4PhysicsManager& right)
72{
73 // check assignement to self
74 if (this == &right) return *this;
75
76 TG4Globals::Exception(
77 "Attempt to assign TG4PhysicsManager singleton.");
78
79 return *this;
80}
81
82// private methods
83
0453c41f 84void TG4PhysicsManager::FillProcessMap()
85{
86// Fills fProcessMap.
87// The default G4 process names are used in the map.
88// ---
89
90 // multiple scattering
91 fProcessMap.Add("msc", kPMultipleScattering);
92 fProcessMap.Add("Imsc", kPMultipleScattering);
93
94 // continuous energy loss
95 // !! including delta rays
96 fProcessMap.Add("eIoni", kPEnergyLoss);
97 fProcessMap.Add("IeIoni", kPEnergyLoss);
98 fProcessMap.Add("LowEnergyIoni", kPEnergyLoss);
99 fProcessMap.Add("hIoni", kPEnergyLoss);
100 fProcessMap.Add("IhIoni", kPEnergyLoss);
101 fProcessMap.Add("hLowEIoni", kPEnergyLoss);
102 fProcessMap.Add("MuIoni", kPEnergyLoss);
103 fProcessMap.Add("IMuIonisation", kPEnergyLoss);
104 fProcessMap.Add("ionIoni", kPEnergyLoss);
105 fProcessMap.Add("ionLowEIoni", kPEnergyLoss);
106 fProcessMap.Add("PAIonisation", kPEnergyLoss);
107
108 // bending in mag. field
109 // kPMagneticFieldL
110
111 // particle decay
112 fProcessMap.Add("Decay", kPDecay);
113
114 // photon pair production or
115 // muon direct pair production
116 fProcessMap.Add("conv", kPPair);
117 fProcessMap.Add("LowEnConversion", kPPair);
118 fProcessMap.Add("MuPairProd", kPPair);
119 fProcessMap.Add("IMuPairProduction", kPPair);
120
121 // Compton scattering
122 fProcessMap.Add("compt", kPCompton);
123 fProcessMap.Add("LowEnCompton", kPCompton);
124 fProcessMap.Add("polarCompt", kPCompton);
125
126 // photoelectric effect
127 fProcessMap.Add("phot", kPPhotoelectric);
128 fProcessMap.Add("LowEnPhotoElec", kPPhotoelectric);
129
130 // bremsstrahlung
131 fProcessMap.Add("eBrem", kPBrem);
132 fProcessMap.Add("IeBrem", kPBrem);
133 fProcessMap.Add("MuBrem", kPBrem);
134 fProcessMap.Add("IMuBremsstrahlung", kPBrem);
135 fProcessMap.Add("LowEnBrem", kPBrem);
136
137 // delta-ray production
138 // kPDeltaRay
139 // has to be distinguished from kPEnergyLoss on flight
140
141 // positron annihilation
142 fProcessMap.Add("annihil", kPAnnihilation);
143 fProcessMap.Add("Iannihil", kPAnnihilation);
144
145 // hadronic interaction
146 // kPHadronic
147
148 // nuclear evaporation
149 // kPEvaporation
150
151 // nuclear fission
152 // kPNuclearFission
153
154 // nuclear absorption
155 fProcessMap.Add("PionMinusAbsorptionAtRest", kPNuclearAbsorption);
156 fProcessMap.Add("PiMinusAbsorptionAtRest", kPNuclearAbsorption);
157 fProcessMap.Add("KaonMinusAbsorption", kPNuclearAbsorption);
158 fProcessMap.Add("KaonMinusAbsorptionAtRest", kPNuclearAbsorption);
159
160 // antiproton annihilation
161 fProcessMap.Add("AntiProtonAnnihilationAtRest", kPPbarAnnihilation);
162 // fProcessMap.Add("AntiNeutronAnnihilationAtRest", not defined);
163
164 // neutron capture
165 fProcessMap.Add("NeutronCaptureAtRest", kPNCapture);
166 // fProcessMap.Add("LCapture", hadron capture not defined);
167
168 // hadronic elastic incoherent scattering
169 fProcessMap.Add("LElastic", kPHElastic);
170
171 // hadronic inelastic scattering
172 fProcessMap.Add("inelastic", kPHInhelastic);
173
174 // muon nuclear interaction
175 fProcessMap.Add("MuNucl", kPMuonNuclear);
176
177 // exceeded time of flight cut
178 // kPTOFlimit
179
180 // nuclear photofission
181 // kPPhotoFission
182
183 // Rayleigh scattering
184 fProcessMap.Add("Rayleigh Scattering", kPRayleigh);
185
186 // no mechanism is active, usually at the entrance of a new volume
094a5851 187 fProcessMap.Add("Transportation", kPNull);
0453c41f 188
189 // particle has fallen below energy threshold and tracking stops
190 // kPStop
191
192 // Cerenkov photon absorption
193 fProcessMap.Add("Absorption", kPLightAbsorption);
194
195 // Cerenkov photon reflection/refraction
196 // kPLightScattering, kPLightReflection, kPLightRefraction
197 // has to be inquired from the G4OpBoundary process
198
199 // synchrotron radiation
200 fProcessMap.Add("SynchrotronRadiation", kPSynchrotron);
201}
202
203
b2030327 204// public methods
205
206void TG4PhysicsManager::BuildPhysics()
2817d3e2 207{
b2030327 208// Empty function - not needed in G4.
209// (Physics is built within /run/initialize.)
2817d3e2 210// ---
211
b2030327 212 TG4Globals::Warning(
213 "TG4PhysicsManager::BuildPhysics: is empty function in G4 MC.");
214}
2817d3e2 215
b2030327 216void TG4PhysicsManager::CreatePhysicsConstructors()
0453c41f 217{
b2030327 218// Creates the selected physics constructors
219// and registeres them in the modular physics list.
0453c41f 220// ---
221
b2030327 222 // electromagnetic physics
223 if (fSetEMPhysics)
224 fPhysicsList->RegisterPhysics(new TG4PhysicsConstructorEM());
0453c41f 225
b2030327 226 // optical physics
227 if (fSetOpticalPhysics)
228 fPhysicsList->RegisterPhysics(new TG4PhysicsConstructorOptical());
2817d3e2 229
b2030327 230 // hadron physics
231 if (fSetHadronPhysics)
232 fPhysicsList->RegisterPhysics(new TG4PhysicsConstructorHadron());
2817d3e2 233
b2030327 234 if (fSetSpecialCutsPhysics)
235 fPhysicsList->RegisterPhysics(new TG4PhysicsConstructorSpecialCuts());
2817d3e2 236
b2030327 237 if (fSetSpecialControlsPhysics)
238 fPhysicsList->RegisterPhysics(new TG4PhysicsConstructorSpecialControls());
239
240 // all created physics constructors are deleted
241 // in the G4VModularPhysicsList destructor
242}
243
244void TG4PhysicsManager::SetCut(const char* cutName, Float_t cutValue)
245{
246// Sets the specified cut.
2817d3e2 247// ---
248
b2030327 249 fG3PhysicsManager->CheckLock();
250 TG4G3Cut g3Cut = fG3PhysicsManager->GetG3Cut(cutName);
251 if (g3Cut != kNoG3Cuts)
252 fG3PhysicsManager->SetCut(g3Cut, cutValue);
253 else {
254 G4String text = "TG4PhysicsManager::SetCut:\n";
255 text = text + " Parameter " + cutName;
256 text = text + " is not implemented.";
257 TG4Globals::Warning(text);
2817d3e2 258 }
2817d3e2 259}
b2030327 260
261void TG4PhysicsManager::SetProcess(const char* controlName, Int_t controlValue)
2817d3e2 262{
b2030327 263// Sets the specified process control.
2817d3e2 264// ---
265
b2030327 266 fG3PhysicsManager->CheckLock();
267 TG4G3Control control = fG3PhysicsManager->GetG3Control(controlName);
268 if (control != kNoG3Controls)
269 fG3PhysicsManager->SetProcess(control, controlValue);
270 else {
271 G4String text = "TG4PhysicsManager::SetProcess:\n";
272 text = text + " Parameter " + controlName;
273 text = text + " is not implemented.";
274 TG4Globals::Warning(text);
2817d3e2 275 }
b2030327 276}
2817d3e2 277
278Float_t TG4PhysicsManager::Xsec(char* ch, Float_t p1, Int_t i1, Int_t i2)
279{
280// Not yet implemented -> gives exception.
281// ---
282
283 TG4Globals::Exception(
284 "TG4PhysicsManager::Xsec: not yet implemented.");
285
286 return 0.;
287}
288
289Int_t TG4PhysicsManager::IdFromPDG(Int_t pdgID) const
290{
291// G4 does not use the integer particle identifiers
292// Id <-> PDG is identity.
293// ---
294
295 return pdgID;
296}
297
298Int_t TG4PhysicsManager::PDGFromId(Int_t mcID) const
299{
300// G4 does not use integer particle identifiers
301// Id <-> PDG is identity.
302// ---
303
304 return mcID;
305}
306
307void TG4PhysicsManager::DefineParticles()
308{
309 // ======
310 // Taken from TGeant3
311 //
312 // Use ENDF-6 mapping for ions = 10000*z+10*a+iso
313 // and add 1 000 000
314 // and numbers above 5 000 000 for special applications
315 //
316
317 const Int_t kion=10000000;
318 const Int_t kspe=50000000;
319
320 const Double_t kGeV=0.9314943228;
321 const Double_t kHslash = 1.0545726663e-27;
322 const Double_t kErgGeV = 1/1.6021773349e-3;
323 const Double_t kHshGeV = kHslash*kErgGeV;
324 const Double_t kYearsToSec = 3600*24*365.25;
325
326 TDatabasePDG *pdgDB = TDatabasePDG::Instance();
327
328 pdgDB->AddParticle("Deuteron","Deuteron",2*kGeV+8.071e-3,kTRUE,
329 0,1,"Ion",kion+10020);
330
331 pdgDB->AddParticle("Triton","Triton",3*kGeV+14.931e-3,kFALSE,
332 kHshGeV/(12.33*kYearsToSec),1,"Ion",kion+10030);
333
334 pdgDB->AddParticle("Alpha","Alpha",4*kGeV+2.424e-3,kTRUE,
335 kHshGeV/(12.33*kYearsToSec),2,"Ion",kion+20040);
336
337 pdgDB->AddParticle("HE3","HE3",3*kGeV+14.931e-3,kFALSE,
338 0,2,"Ion",kion+20030);
339
340 pdgDB->AddParticle("Cherenkov","Cherenkov",0,kFALSE,
341 0,0,"Special",kspe+50);
342
343 pdgDB->AddParticle("FeedbackPhoton","FeedbackPhoton",0,kFALSE,
344 0,0,"Special",kspe+51);
345
346
347 // To do: define the PDG database extension
348 // in a common part.
349 //
350 // AliMC::ExtendPDGDatabase();
351 //
352 // end of "common" implementation
353 // ======
354
b2030327 355 fParticlesManager->MapParticles();
2817d3e2 356}
357
2817d3e2 358
359void TG4PhysicsManager::SetProcessActivation()
360{
361// (In)Activates built processes according
b2030327 362// to the setup in fControlVector.
2817d3e2 363// ---
364
365 if (fPhysicsList) {
366 // temporarily excluded
367 // fPhysicsList->SetProcessActivation();
368 }
369 else {
370 G4String text = "TG4PhysicsManager::SetProcessActivation:\n";
371 text = text + " There is no physics list set.";
372 TG4Globals::Exception(text);
373 }
374}
375
2817d3e2 376
b2030327 377AliMCProcess TG4PhysicsManager::GetMCProcess(const G4VProcess* process)
2817d3e2 378{
b2030327 379// Returns the AliMCProcess code of the specified G4 process.
2817d3e2 380// ---
b2030327 381
382 if (!process) return kPNoProcess;
2817d3e2 383
b2030327 384 G4String name = process->GetProcessName();
385 G4int code = fProcessMap.GetSecond(name);
386
387 if (code == 0) return kPNoProcess;
388
389 return (AliMCProcess)code;
2817d3e2 390}
391