| Geant4 Cross Reference |
1 // 1
2 // *******************************************
3 // * License and Disclaimer
4 // *
5 // * The Geant4 software is copyright of th
6 // * the Geant4 Collaboration. It is provided
7 // * conditions of the Geant4 Software License
8 // * LICENSE and available at http://cern.ch/
9 // * include a list of copyright holders.
10 // *
11 // * Neither the authors of this software syst
12 // * institutes,nor the agencies providing fin
13 // * work make any representation or warran
14 // * regarding this software system or assum
15 // * use. Please see the license in the file
16 // * for the full disclaimer and the limitatio
17 // *
18 // * This code implementation is the result
19 // * technical work of the GEANT4 collaboratio
20 // * By using, copying, modifying or distri
21 // * any work based on the software) you ag
22 // * use in resulting scientific publicati
23 // * acceptance of all terms of the Geant4 Sof
24 // *******************************************
25 //
26 // G4DynamicParticle class implementation
27 //
28 // History:
29 // - 2 December 1995, G.Cosmo - first design,
30 // - 29 January 1996, M.Asai - first implement
31 // - 1996 - 2007, H.Kurashige - revisions.
32 // - 15 March 2019, M.Novak - log-kinetic en
33 // on demand if its stored
34 //--------------------------------------------
35
36 #include "G4DynamicParticle.hh"
37
38 #include "G4DecayProducts.hh"
39 #include "G4IonTable.hh"
40 #include "G4LorentzVector.hh"
41 #include "G4ParticleDefinition.hh"
42 #include "G4PrimaryParticle.hh"
43
44 G4Allocator<G4DynamicParticle>*& pDynamicParti
45 {
46 G4ThreadLocalStatic G4Allocator<G4DynamicPar
47 return _instance;
48 }
49
50 static const G4double EnergyMomentumRelationAl
51 static const G4double EnergyMRA2 =
52 EnergyMomentumRelationAllowance * EnergyMome
53
54 G4DynamicParticle::G4DynamicParticle()
55 : theMomentumDirection(0.0, 0.0, 1.0), thePo
56 {}
57
58 G4DynamicParticle::G4DynamicParticle(const G4P
59 const G4T
60 G4double
61 : theMomentumDirection(aMomentumDirection),
62 thePolarization(0.0, 0.0, 0.0),
63 theParticleDefinition(aParticleDefinition)
64 theKineticEnergy(aKineticEnergy),
65 theDynamicalMass(aParticleDefinition->GetP
66 theDynamicalCharge(aParticleDefinition->Ge
67 theDynamicalSpin(aParticleDefinition->GetP
68 theDynamicalMagneticMoment(aParticleDefini
69 {}
70
71 G4DynamicParticle::G4DynamicParticle(const G4P
72 const G4T
73 G4double
74 : theMomentumDirection(aMomentumDirection),
75 thePolarization(0.0, 0.0, 0.0),
76 theParticleDefinition(aParticleDefinition)
77 theKineticEnergy(aKineticEnergy),
78 theDynamicalMass(aParticleDefinition->GetP
79 theDynamicalCharge(aParticleDefinition->Ge
80 theDynamicalSpin(aParticleDefinition->GetP
81 theDynamicalMagneticMoment(aParticleDefini
82 {
83 if (std::abs(theDynamicalMass - dynamicalMas
84 if (dynamicalMass > EnergyMomentumRelation
85 theDynamicalMass = dynamicalMass;
86 else
87 theDynamicalMass = 0.0;
88 }
89 }
90
91 G4DynamicParticle::G4DynamicParticle(const G4P
92 const G4T
93 : thePolarization(0.0, 0.0, 0.0),
94 theParticleDefinition(aParticleDefinition)
95 theDynamicalMass(aParticleDefinition->GetP
96 theDynamicalCharge(aParticleDefinition->Ge
97 theDynamicalSpin(aParticleDefinition->GetP
98 theDynamicalMagneticMoment(aParticleDefini
99 {
100 SetMomentum(aParticleMomentum); // 3-dim mo
101 }
102
103 G4DynamicParticle::G4DynamicParticle(const G4P
104 const G4L
105 : thePolarization(0.0, 0.0, 0.0),
106 theParticleDefinition(aParticleDefinition)
107 theDynamicalMass(aParticleDefinition->GetP
108 theDynamicalCharge(aParticleDefinition->Ge
109 theDynamicalSpin(aParticleDefinition->GetP
110 theDynamicalMagneticMoment(aParticleDefini
111 {
112 Set4Momentum(aParticleMomentum); // 4-momen
113 }
114
115 G4DynamicParticle::G4DynamicParticle(const G4P
116 G4double
117 : thePolarization(0.0, 0.0, 0.0),
118 theParticleDefinition(aParticleDefinition)
119 theDynamicalMass(aParticleDefinition->GetP
120 theDynamicalCharge(aParticleDefinition->Ge
121 theDynamicalSpin(aParticleDefinition->GetP
122 theDynamicalMagneticMoment(aParticleDefini
123 {
124 // total energy and 3-dim momentum are given
125 G4double pModule2 = aParticleMomentum.mag2()
126 if (pModule2 > 0.0) {
127 G4double mass2 = totalEnergy * totalEnergy
128 G4double PDGmass2 = (aParticleDefinition->
129 SetMomentumDirection(aParticleMomentum.uni
130 if (mass2 < EnergyMRA2) {
131 theDynamicalMass = 0.;
132 SetKineticEnergy(totalEnergy);
133 }
134 else {
135 if (std::abs(PDGmass2 - mass2) > EnergyM
136 theDynamicalMass = std::sqrt(mass2);
137 SetKineticEnergy(totalEnergy - theDyna
138 }
139 else {
140 SetKineticEnergy(totalEnergy - theDyna
141 }
142 }
143 }
144 else {
145 SetMomentumDirection(1.0, 0.0, 0.0);
146 SetKineticEnergy(0.0);
147 }
148 }
149
150 G4DynamicParticle::G4DynamicParticle(const G4D
151 : theMomentumDirection(right.theMomentumDire
152 thePolarization(right.thePolarization),
153 theParticleDefinition(right.theParticleDef
154 // Don't copy preassignedDecayProducts
155 primaryParticle(right.primaryParticle),
156 theKineticEnergy(right.theKineticEnergy),
157 theLogKineticEnergy(right.theLogKineticEne
158 theBeta(right.theBeta),
159 theProperTime(right.theProperTime),
160 theDynamicalMass(right.theDynamicalMass),
161 theDynamicalCharge(right.theDynamicalCharg
162 theDynamicalSpin(right.theDynamicalSpin),
163 theDynamicalMagneticMoment(right.theDynami
164
165 verboseLevel(right.verboseLevel),
166 thePDGcode(right.thePDGcode)
167 {
168 if (right.theElectronOccupancy != nullptr) {
169 theElectronOccupancy = new G4ElectronOccup
170 }
171 }
172
173 G4DynamicParticle::G4DynamicParticle(G4Dynamic
174 : theMomentumDirection(from.theMomentumDirec
175 thePolarization(from.thePolarization),
176 theParticleDefinition(from.theParticleDefi
177 theElectronOccupancy(from.theElectronOccup
178 // Don't move preassignedDecayProducts
179 primaryParticle(from.primaryParticle),
180 theKineticEnergy(from.theKineticEnergy),
181 theLogKineticEnergy(from.theLogKineticEner
182 theBeta(from.theBeta),
183 theProperTime(from.theProperTime),
184 theDynamicalMass(from.theDynamicalMass),
185 theDynamicalCharge(from.theDynamicalCharge
186 theDynamicalSpin(from.theDynamicalSpin),
187 theDynamicalMagneticMoment(from.theDynamic
188
189 verboseLevel(from.verboseLevel),
190 thePDGcode(from.thePDGcode)
191 {
192 // Release the data from the source object
193 from.theParticleDefinition = nullptr;
194 from.theElectronOccupancy = nullptr;
195 from.thePreAssignedDecayProducts = nullptr;
196 from.primaryParticle = nullptr;
197 }
198
199 G4DynamicParticle::~G4DynamicParticle()
200 {
201 delete thePreAssignedDecayProducts;
202 thePreAssignedDecayProducts = nullptr;
203
204 delete theElectronOccupancy;
205 theElectronOccupancy = nullptr;
206 }
207
208 G4DynamicParticle& G4DynamicParticle::operator
209 {
210 if (this != &right) {
211 theMomentumDirection = right.theMomentumDi
212 theParticleDefinition = right.theParticleD
213 thePolarization = right.thePolarization;
214 theKineticEnergy = right.theKineticEnergy;
215 theProperTime = right.theProperTime;
216
217 theDynamicalMass = right.theDynamicalMass;
218 theDynamicalCharge = right.theDynamicalCha
219 theDynamicalSpin = right.theDynamicalSpin;
220 theDynamicalMagneticMoment = right.theDyna
221
222 delete theElectronOccupancy;
223 if (right.theElectronOccupancy != nullptr)
224 theElectronOccupancy = new G4ElectronOcc
225 }
226 else {
227 theElectronOccupancy = nullptr;
228 }
229
230 // thePreAssignedDecayProducts must not be
231 thePreAssignedDecayProducts = nullptr;
232 thePreAssignedDecayTime = -1.0;
233
234 verboseLevel = right.verboseLevel;
235
236 // Primary particle information must be pr
237 //*** primaryParticle = right.primaryParti
238
239 thePDGcode = right.thePDGcode;
240 }
241 return *this;
242 }
243
244 G4DynamicParticle& G4DynamicParticle::operator
245 {
246 if (this != &from) {
247 theMomentumDirection = from.theMomentumDir
248 thePolarization = from.thePolarization;
249 theKineticEnergy = from.theKineticEnergy;
250 theProperTime = from.theProperTime;
251
252 theDynamicalMass = from.theDynamicalMass;
253 theDynamicalCharge = from.theDynamicalChar
254 theDynamicalSpin = from.theDynamicalSpin;
255 theDynamicalMagneticMoment = from.theDynam
256
257 delete theElectronOccupancy;
258 theElectronOccupancy = from.theElectronOcc
259 from.theElectronOccupancy = nullptr;
260
261 // thePreAssignedDecayProducts must not be
262 thePreAssignedDecayProducts = nullptr;
263 from.thePreAssignedDecayProducts = nullptr
264 thePreAssignedDecayTime = -1.0;
265
266 theParticleDefinition = from.theParticleDe
267 from.theParticleDefinition = nullptr;
268
269 verboseLevel = from.verboseLevel;
270
271 primaryParticle = from.primaryParticle;
272 from.primaryParticle = nullptr;
273
274 thePDGcode = from.thePDGcode;
275 }
276 return *this;
277 }
278
279 void G4DynamicParticle::SetDefinition(const G4
280 {
281 // remove preassigned decay
282 if (thePreAssignedDecayProducts != nullptr)
283 #ifdef G4VERBOSE
284 if (verboseLevel > 0) {
285 G4cout << " G4DynamicParticle::SetDefini
286 << "!!! Pre-assigned decay produc
287 G4cout << "!!! New Definition is " << aP
288 << G4endl;
289 G4cout << "!!! Pre-assigned decay produc
290 }
291 #endif
292 delete thePreAssignedDecayProducts;
293 }
294 thePreAssignedDecayProducts = nullptr;
295
296 theParticleDefinition = aParticleDefinition;
297
298 // set Dynamic mass/charge
299 SetMass(theParticleDefinition->GetPDGMass())
300 theDynamicalCharge = theParticleDefinition->
301 theDynamicalSpin = theParticleDefinition->Ge
302 theDynamicalMagneticMoment = theParticleDefi
303
304 // Set electron orbits
305 if (theElectronOccupancy != nullptr) {
306 delete theElectronOccupancy;
307 theElectronOccupancy = nullptr;
308 }
309 }
310
311 G4bool G4DynamicParticle::operator==(const G4D
312 {
313 return (this == (G4DynamicParticle*)&right);
314 }
315
316 G4bool G4DynamicParticle::operator!=(const G4D
317 {
318 return (this != (G4DynamicParticle*)&right);
319 }
320
321 void G4DynamicParticle::AllocateElectronOccupa
322 {
323 if (G4IonTable::IsIon(theParticleDefinition)
324 // Only ions can have ElectronOccupancy
325 theElectronOccupancy = new G4ElectronOccup
326 }
327 else {
328 theElectronOccupancy = nullptr;
329 }
330 }
331
332 void G4DynamicParticle::SetMomentum(const G4Th
333 {
334 G4double pModule2 = momentum.mag2();
335 if (pModule2 > 0.0) {
336 const G4double mass = theDynamicalMass;
337 SetMomentumDirection(momentum.unit());
338 SetKineticEnergy(pModule2 / (std::sqrt(pMo
339 }
340 else {
341 SetMomentumDirection(1.0, 0.0, 0.0);
342 SetKineticEnergy(0.0);
343 }
344 }
345
346 void G4DynamicParticle::Set4Momentum(const G4L
347 {
348 G4double pModule2 = momentum.vect().mag2();
349 if (pModule2 > 0.0) {
350 SetMomentumDirection(momentum.vect().unit(
351 const G4double totalenergy = momentum.t();
352 const G4double mass2 = totalenergy * total
353 const G4double PDGmass2 =
354 (theParticleDefinition->GetPDGMass()) *
355 if (mass2 < EnergyMRA2) {
356 theDynamicalMass = 0.;
357 }
358 else if (std::abs(PDGmass2 - mass2) > Ener
359 theDynamicalMass = std::sqrt(mass2);
360 }
361 SetKineticEnergy(totalenergy - theDynamica
362 }
363 else {
364 SetMomentumDirection(1.0, 0.0, 0.0);
365 SetKineticEnergy(0.0);
366 }
367 }
368
369 #ifdef G4VERBOSE
370 void G4DynamicParticle::DumpInfo(G4int mode) c
371 {
372 if (theParticleDefinition == nullptr) {
373 G4cout << " G4DynamicParticle::DumpInfo()
374 }
375 else {
376 G4cout << " Particle type - " << thePartic
377 << " mass: " << GetMass()
378 << " charge: " << GetCharge(
379 << " Direction x: " << GetMomentu
380 << ", y: " << GetMomentumDirection(
381 << G4endl << " Total Momentum = "
382 << G4endl << " Momentum: " << Get
383 << ", y: " << GetMomentum().y() / C
384 << ", z: " << GetMomentum().z() / C
385 << " Total Energy = " << GetTot
386 << " Kinetic Energy = " << GetKin
387 << " MagneticMoment [MeV/T]: " <<
388 << G4endl << " ProperTime = "
389
390 if (mode > 0) {
391 if (theElectronOccupancy != nullptr) {
392 theElectronOccupancy->DumpInfo();
393 }
394 }
395 }
396 }
397 #else
398 void G4DynamicParticle::DumpInfo(G4int) const
399 {
400 return;
401 }
402 #endif
403
404 G4double G4DynamicParticle::GetElectronMass()
405 {
406 return CLHEP::electron_mass_c2;
407 }
408