]>
Commit | Line | Data |
---|---|---|
f04eb3a5 | 1 | // $Id$ |
2 | // Category: physics | |
3 | // | |
e5967ab3 | 4 | // Author: I. Hrivnacova |
5 | // | |
6 | // Class TG4ModularPhysicsList | |
7 | // --------------------------- | |
f04eb3a5 | 8 | // See the class description in the header file. |
9 | ||
10 | #include "TG4ModularPhysicsList.h" | |
e5967ab3 | 11 | #include "TG4GeometryServices.h" |
f04eb3a5 | 12 | #include "TG4G3PhysicsManager.h" |
13 | #include "TG4G3ControlVector.h" | |
e5967ab3 | 14 | #include "TG4ProcessControlMap.h" |
f04eb3a5 | 15 | |
16 | #include <G4ParticleDefinition.hh> | |
17 | #include <G4ProcessManager.hh> | |
18 | #include <G4BosonConstructor.hh> | |
19 | #include <G4LeptonConstructor.hh> | |
20 | #include <G4MesonConstructor.hh> | |
21 | #include <G4BaryonConstructor.hh> | |
22 | #include <G4IonConstructor.hh> | |
23 | #include <G4ShortLivedConstructor.hh> | |
24 | #include <G4ProcessTable.hh> | |
f04eb3a5 | 25 | |
696e37fa | 26 | const G4double TG4ModularPhysicsList::fgkDefaultCutValue = 1.0 * mm; |
f04eb3a5 | 27 | |
72095f7c | 28 | //_____________________________________________________________________________ |
f04eb3a5 | 29 | TG4ModularPhysicsList::TG4ModularPhysicsList() |
696e37fa | 30 | : G4VModularPhysicsList() { |
8496e315 | 31 | // |
e5967ab3 | 32 | defaultCutValue = fgkDefaultCutValue; |
f04eb3a5 | 33 | |
f179c7fc | 34 | SetVerboseLevel(1); |
f04eb3a5 | 35 | } |
36 | ||
72095f7c | 37 | //_____________________________________________________________________________ |
8496e315 | 38 | TG4ModularPhysicsList::TG4ModularPhysicsList(const TG4ModularPhysicsList& right) |
39 | { | |
40 | // | |
41 | TG4Globals::Exception("TG4ModularPhysicsList is protected from copying."); | |
42 | } | |
43 | ||
72095f7c | 44 | //_____________________________________________________________________________ |
f04eb3a5 | 45 | TG4ModularPhysicsList::~TG4ModularPhysicsList() { |
46 | // | |
8496e315 | 47 | //delete fExtDecayer; |
48 | // fExtDecayer is deleted in G4Decay destructor | |
49 | } | |
50 | ||
51 | // operators | |
52 | ||
72095f7c | 53 | //_____________________________________________________________________________ |
8496e315 | 54 | TG4ModularPhysicsList& |
55 | TG4ModularPhysicsList::operator=(const TG4ModularPhysicsList &right) | |
56 | { | |
57 | // check assignement to self | |
58 | if (this == &right) return *this; | |
59 | ||
60 | TG4Globals::Exception("TG4ModularPhysicsList is protected from assigning."); | |
61 | ||
62 | return *this; | |
f04eb3a5 | 63 | } |
64 | ||
e5967ab3 | 65 | // private methods |
66 | ||
67 | //_____________________________________________________________________________ | |
68 | void TG4ModularPhysicsList::SetProcessActivation( | |
69 | G4ProcessManager* processManager, | |
70 | G4int processId, G4bool activation) | |
71 | { | |
72 | // Sets process activation for the given process. | |
73 | // --- | |
74 | ||
75 | G4String strActivation = "activation"; | |
76 | if (!activation) strActivation = "inactivation"; | |
77 | ||
78 | if (verboseLevel>1) { | |
79 | G4cout << "Set process " << strActivation << " for " | |
80 | << (*processManager->GetProcessList())[processId]->GetProcessName() | |
81 | << G4endl; | |
82 | } | |
83 | ||
84 | processManager->SetProcessActivation(processId, activation); | |
85 | } | |
86 | ||
f04eb3a5 | 87 | // protected methods |
88 | ||
72095f7c | 89 | //_____________________________________________________________________________ |
f04eb3a5 | 90 | void TG4ModularPhysicsList::ConstructParticle() |
91 | { | |
92 | // In this method, static member functions should be called | |
93 | // for all particles which you want to use. | |
94 | // This ensures that objects of these particle types will be | |
95 | // created in the program. | |
96 | // --- | |
97 | ||
f04eb3a5 | 98 | // lock physics manager |
99 | TG4G3PhysicsManager* g3PhysicsManager = TG4G3PhysicsManager::Instance(); | |
100 | g3PhysicsManager->Lock(); | |
f04eb3a5 | 101 | |
102 | // create particles for registered physics | |
103 | G4VModularPhysicsList::ConstructParticle(); | |
104 | } | |
105 | ||
72095f7c | 106 | //_____________________________________________________________________________ |
f04eb3a5 | 107 | void TG4ModularPhysicsList::ConstructProcess() |
108 | { | |
109 | // Constructs all processes. | |
110 | // --- | |
111 | ||
f04eb3a5 | 112 | // create processes for registered physics |
113 | G4VModularPhysicsList::ConstructProcess(); | |
114 | ||
f04eb3a5 | 115 | // verbose |
696e37fa | 116 | if (verboseLevel>1) DumpAllProcesses(); |
f04eb3a5 | 117 | } |
118 | ||
f04eb3a5 | 119 | |
120 | // public methods | |
121 | ||
72095f7c | 122 | //_____________________________________________________________________________ |
f04eb3a5 | 123 | void TG4ModularPhysicsList::SetCuts() |
124 | { | |
125 | // Sets the default cut value for all particle types | |
126 | // other then e-/e+. | |
127 | // The cut value for e-/e+ is high in oredr to supress | |
128 | // tracking of delta electrons. | |
129 | // --- | |
130 | ||
131 | // SetCutsWithDefault(); | |
132 | // "G4VUserPhysicsList::SetCutsWithDefault" method sets | |
133 | // the default cut value for all particle types. | |
134 | ||
135 | // default cut value | |
136 | G4double cut = defaultCutValue; | |
8496e315 | 137 | //G4double ecut = 10.*m; |
138 | G4double ecut = cut; | |
f04eb3a5 | 139 | |
140 | #ifdef G4VERBOSE | |
141 | if (verboseLevel >1){ | |
696e37fa | 142 | G4cout << "TG4ModularPhysicsList::SetCutsWithDefault:"; |
f04eb3a5 | 143 | G4cout << "CutLength : " << cut/mm << " (mm)" << G4endl; |
144 | } | |
145 | #endif | |
146 | ||
147 | // set cut values for gamma at first and for e- second and next for e+, | |
148 | // because some processes for e+/e- need cut values for gamma | |
149 | SetCutValue(cut, "gamma"); | |
150 | SetCutValue(ecut, "e-"); | |
151 | SetCutValue(ecut, "e+"); | |
152 | ||
153 | // set cut values for proton and anti_proton before all other hadrons | |
154 | // because some processes for hadrons need cut values for proton/anti_proton | |
155 | SetCutValue(cut, "proton"); | |
156 | SetCutValue(cut, "anti_proton"); | |
157 | ||
158 | SetCutValueForOthers(cut); | |
159 | ||
160 | if (verboseLevel>1) { | |
161 | DumpCutValuesTable(); | |
162 | } | |
163 | } | |
164 | ||
72095f7c | 165 | //_____________________________________________________________________________ |
f04eb3a5 | 166 | void TG4ModularPhysicsList::SetProcessActivation() |
167 | { | |
168 | // (In)Activates built processes according | |
169 | // to the setup in TG4G3PhysicsManager::fControlVector. | |
170 | // --- | |
171 | ||
e5967ab3 | 172 | TG4G3ControlVector* controlVector |
173 | = TG4G3PhysicsManager::Instance()->GetControlVector(); | |
174 | ||
175 | G4bool specialControls | |
176 | = TG4GeometryServices::Instance()->IsSpecialControls(); | |
177 | ||
178 | if (!specialControls) | |
179 | G4cout << G4endl | |
180 | << "### No special controls in user limits are set." << G4endl; | |
181 | ||
182 | if (!controlVector) { | |
f04eb3a5 | 183 | G4String text = "TG4ModularPhysicsList::SetProcessActivation: \n"; |
184 | text = text + " Vector of processes controls is not set."; | |
185 | TG4Globals::Warning(text); | |
e5967ab3 | 186 | return; |
f04eb3a5 | 187 | } |
e5967ab3 | 188 | |
189 | theParticleIterator->reset(); | |
190 | while ((*theParticleIterator)()) | |
191 | { | |
192 | G4ProcessManager* processManager | |
193 | = theParticleIterator->value()->GetProcessManager(); | |
194 | ||
195 | G4ProcessVector* processVector = processManager->GetProcessList(); | |
196 | ||
197 | // set processes controls | |
198 | for (G4int i=0; i<processVector->length(); i++) { | |
199 | ||
200 | TG4G3ControlValue control | |
201 | = controlVector->GetControlValue((*processVector)[i]); | |
202 | G4bool activation = processManager->GetProcessActivation(i); | |
203 | ||
204 | if (control != kUnset) { | |
205 | if (!TG4Globals::Compare(activation, control)) { | |
206 | ||
207 | // set new process activation | |
208 | G4bool activate; | |
209 | if (control == kInActivate) activate = false; | |
210 | else activate = true; | |
211 | ||
212 | SetProcessActivation(processManager, i, activate); | |
213 | } | |
214 | } | |
215 | else { | |
216 | // control == kUnset | |
217 | if ((*processVector)[i]->GetProcessName().find("specialControl") | |
218 | != G4String::npos) { | |
219 | ||
220 | SetProcessActivation(processManager, i, specialControls); | |
221 | } | |
222 | } | |
223 | } | |
224 | } | |
f04eb3a5 | 225 | } |
226 | ||
72095f7c | 227 | //_____________________________________________________________________________ |
f04eb3a5 | 228 | void TG4ModularPhysicsList::PrintAllProcesses() const |
229 | { | |
230 | // Prints all processes. | |
231 | // --- | |
232 | ||
233 | G4cout << "TG4ModularPhysicsList processes: " << G4endl; | |
696e37fa | 234 | G4cout << "================================ " << G4endl; |
f04eb3a5 | 235 | |
236 | G4ProcessTable* processTable = G4ProcessTable::GetProcessTable(); | |
237 | G4ProcessTable::G4ProcNameVector* processNameList | |
238 | = processTable->GetNameList(); | |
239 | ||
240 | for (G4int i=0; i <processNameList->size(); i++){ | |
241 | G4cout << " " << (*processNameList)[i] << G4endl; | |
242 | } | |
243 | } | |
244 | ||
696e37fa | 245 | //_____________________________________________________________________________ |
246 | void TG4ModularPhysicsList::DumpAllProcesses() const | |
247 | { | |
248 | // Dumps all particles and their processes. | |
249 | // --- | |
250 | ||
251 | G4cout << "TG4ModularPhysicsList particles and processes: " << G4endl; | |
252 | G4cout << "============================================== " << G4endl; | |
253 | ||
254 | theParticleIterator->reset(); | |
255 | while ((*theParticleIterator)()) | |
256 | { | |
257 | // print particle name | |
258 | G4cout << "Particle: " | |
259 | << theParticleIterator->value()->GetParticleName() | |
260 | << G4endl; | |
261 | ||
262 | // dump particle processes | |
263 | G4ProcessVector* processVector | |
264 | = theParticleIterator->value()->GetProcessManager()->GetProcessList(); | |
265 | for (G4int i=0; i<processVector->length(); i++) | |
266 | (*processVector)[i]->DumpInfo(); | |
267 | ||
268 | G4cout << G4endl; | |
269 | } | |
270 | } | |
271 |