source: trunk/source/particles/management/src/G4DynamicParticle.cc @ 1340

Last change on this file since 1340 was 1340, checked in by garnier, 14 years ago

update ti head
File size: 15.0 KB
Line 
1//
2// ********************************************************************
3// * License and Disclaimer                                           *
4// *                                                                  *
5// * The  Geant4 software  is  copyright of the Copyright Holders  of *
6// * the Geant4 Collaboration.  It is provided  under  the terms  and *
7// * conditions of the Geant4 Software License,  included in the file *
8// * LICENSE and available at  http://cern.ch/geant4/license .  These *
9// * include a list of copyright holders.                             *
10// *                                                                  *
11// * Neither the authors of this software system, nor their employing *
12// * institutes,nor the agencies providing financial support for this *
13// * work  make  any representation or  warranty, express or implied, *
14// * regarding  this  software system or assume any liability for its *
15// * use.  Please see the license in the file  LICENSE  and URL above *
16// * for the full disclaimer and the limitation of liability.         *
17// *                                                                  *
18// * This  code  implementation is the result of  the  scientific and *
19// * technical work of the GEANT4 collaboration.                      *
20// * By using,  copying,  modifying or  distributing the software (or *
21// * any work based  on the software)  you  agree  to acknowledge its *
22// * use  in  resulting  scientific  publications,  and indicate your *
23// * acceptance of all terms of the Geant4 Software license.          *
24// ********************************************************************
25//
26//
27// $Id: G4DynamicParticle.cc,v 1.30 2010/08/10 15:47:42 kurasige Exp $
28// GEANT4 tag $Name: particles-V09-03-15 $
29//
30//
31// --------------------------------------------------------------
32//      GEANT 4 class implementation file
33//
34//      History: first implementation, based on object model of
35//      2nd December 1995, G.Cosmo
36//      ---------------- G4DynamicParticle  ----------------
37//      first implementation by Makoto Asai, 29 January 1996
38//      revised by G.Cosmo, 29 February 1996
39//      revised by H.Kurashige 06 May 1996
40//      revised by Hisaya Kurashige, 27 July 1996
41//         modify thePreAssignedDecayProducts
42//         add   void SetMomentum(G4ThreeVector &momentum)
43//         add   void Set4Momentum(G4LorentzVector &momentum)
44//         add   G4DynamicParticle(G4ParticleDefinition * aParticleDefinition,
45//                                 G4LorentzVector &p4vector)
46//      revised by Hisaya Kurashige, 19 Oct 1996
47//         add    theKillProcess
48//         add    ProperTime
49//      revised by Hisaya Kurashige, 26 Mar 1997
50//         modify destructor
51//      revised by Hisaya Kurashige, 05 June 1997
52//         modify DumpInfo()
53//      revised by Hisaya Kurashige, 5  June 1998
54//         remove    theKillProcess
55//      revised by Hisaya Kurashige, 5  Mar 2001
56//         fixed  SetDefinition()
57//      revised by V.Ivanchenko,    12 June 2003
58//         fixed problem of massless particles
59//      revised by V.Ivanchenko,    18 June 2003
60//         take into account the case of virtual photons
61//      revised by M.Kelsey         12 May 2010
62//         ensure that all constructors initialize all data members
63//--------------------------------------------------------------
64
65#include "G4DynamicParticle.hh"
66#include "G4DecayProducts.hh"
67#include "G4LorentzVector.hh"
68#include "G4ParticleDefinition.hh"
69#include "G4ParticleTable.hh"
70#include "G4IonTable.hh"
71
72G4Allocator<G4DynamicParticle> aDynamicParticleAllocator;
73
74static const G4double EnergyMomentumRelationAllowance = keV;
75
76////////////////////
77G4DynamicParticle::G4DynamicParticle():
78                   theMomentumDirection(0.0,0.0,1.0),
79                   theParticleDefinition(0),
80                   theKineticEnergy(0.0),
81                   theProperTime(0.0),
82                   theDynamicalMass(0.0),
83                   theDynamicalCharge(0.0),
84                   theDynamicalSpin(0.0),
85                   theDynamicalMagneticMoment(0.0),
86                   theElectronOccupancy(0),
87                   thePreAssignedDecayProducts(0),
88                   thePreAssignedDecayTime(-1.0),
89                   verboseLevel(1),
90                   primaryParticle(0),
91                   thePDGcode(0) {}
92
93////////////////////
94// -- constructors ----
95////////////////////
96G4DynamicParticle::G4DynamicParticle(const G4ParticleDefinition * aParticleDefinition,
97                                     const G4ThreeVector& aMomentumDirection,
98                                     G4double aKineticEnergy):
99                   theMomentumDirection(aMomentumDirection),
100                   theParticleDefinition(aParticleDefinition),
101                   theKineticEnergy(aKineticEnergy),
102                   theProperTime(0.0),
103                   theDynamicalMass(aParticleDefinition->GetPDGMass()),
104                   theDynamicalCharge(aParticleDefinition->GetPDGCharge()),
105                   theDynamicalSpin(aParticleDefinition->GetPDGSpin()),
106                   theDynamicalMagneticMoment(aParticleDefinition->GetPDGMagneticMoment()),
107                   theElectronOccupancy(0),
108                   thePreAssignedDecayProducts(0),
109                   thePreAssignedDecayTime(-1.0),
110                   verboseLevel(1),
111                   primaryParticle(0),
112                   thePDGcode(0) {}
113
114////////////////////
115G4DynamicParticle::G4DynamicParticle(const G4ParticleDefinition * aParticleDefinition,
116                                     const G4ThreeVector& aParticleMomentum):
117                   theParticleDefinition(aParticleDefinition),
118                   theKineticEnergy(0.0),
119                   theProperTime(0.0),
120                   theDynamicalMass(aParticleDefinition->GetPDGMass()),
121                   theDynamicalCharge(aParticleDefinition->GetPDGCharge()),
122                   theDynamicalSpin(aParticleDefinition->GetPDGSpin()),
123                   theDynamicalMagneticMoment(aParticleDefinition->GetPDGMagneticMoment()),
124                   theElectronOccupancy(0),
125                   thePreAssignedDecayProducts(0),
126                   thePreAssignedDecayTime(-1.0),
127                   verboseLevel(1),
128                   primaryParticle(0),
129                   thePDGcode(0)
130{
131  // 3-dim momentum is given
132  SetMomentum(aParticleMomentum);
133}
134
135////////////////////
136G4DynamicParticle::G4DynamicParticle(const G4ParticleDefinition * aParticleDefinition,
137                                     const G4LorentzVector   &aParticleMomentum):
138                   theParticleDefinition(aParticleDefinition),
139                   theKineticEnergy(0.0),
140                   theProperTime(0.0),
141                   theDynamicalMass(aParticleDefinition->GetPDGMass()),
142                   theDynamicalCharge(aParticleDefinition->GetPDGCharge()),
143                   theDynamicalSpin(aParticleDefinition->GetPDGSpin()),
144                   theDynamicalMagneticMoment(aParticleDefinition->GetPDGMagneticMoment()),
145                   theElectronOccupancy(0),
146                   thePreAssignedDecayProducts(0),
147                   thePreAssignedDecayTime(-1.0),
148                   verboseLevel(1),
149                   primaryParticle(0),
150                   thePDGcode(0)
151{
152  // 4-momentum vector (Lorentz vecotr) is given
153  Set4Momentum(aParticleMomentum);
154}
155
156G4DynamicParticle::G4DynamicParticle(const G4ParticleDefinition * aParticleDefinition,
157                                     G4double totalEnergy,
158                                     const G4ThreeVector &aParticleMomentum):
159                   theParticleDefinition(aParticleDefinition),
160                   theKineticEnergy(0.0),
161                   theProperTime(0.0),
162                   theDynamicalMass(aParticleDefinition->GetPDGMass()),
163                   theDynamicalCharge(aParticleDefinition->GetPDGCharge()),
164                   theDynamicalSpin(aParticleDefinition->GetPDGSpin()),
165                   theDynamicalMagneticMoment(aParticleDefinition->GetPDGMagneticMoment()),
166                   theElectronOccupancy(0),
167                   thePreAssignedDecayProducts(0),
168                   thePreAssignedDecayTime(-1.0),
169                   verboseLevel(1),
170                   primaryParticle(0),
171                   thePDGcode(0)
172{
173  // total energy and 3-dim momentum are given
174  G4double pModule2 = aParticleMomentum.mag2();
175  if (pModule2>0.0) {
176    G4double mass2 = totalEnergy*totalEnergy - pModule2;
177    if(mass2 < EnergyMomentumRelationAllowance*EnergyMomentumRelationAllowance) {
178      theDynamicalMass = 0.;
179      SetMomentumDirection(aParticleMomentum.unit());
180      SetKineticEnergy(totalEnergy);
181    } else {
182      theDynamicalMass = std::sqrt(mass2);
183      SetMomentum(aParticleMomentum);
184    }
185  } else {
186    SetMomentumDirection(1.0,0.0,0.0);
187    SetKineticEnergy(0.0);
188  }
189}
190
191////////////////////
192G4DynamicParticle::G4DynamicParticle(const G4DynamicParticle &right):
193  theMomentumDirection(right.theMomentumDirection),
194  theParticleDefinition(right.theParticleDefinition),
195  thePolarization(right.thePolarization),
196  theKineticEnergy(right.theKineticEnergy),
197  theProperTime(0.0),
198  theDynamicalMass(right.theDynamicalMass),
199  theDynamicalCharge(right.theDynamicalCharge),
200  theDynamicalSpin(right.theDynamicalSpin),
201  theDynamicalMagneticMoment(right.theDynamicalMagneticMoment),
202  theElectronOccupancy(0),
203  thePreAssignedDecayProducts(0),       // Do not copy preassignedDecayProducts
204  thePreAssignedDecayTime(-1.0),
205  verboseLevel(right.verboseLevel),
206  primaryParticle(right.primaryParticle),
207  thePDGcode(right.thePDGcode)
208{
209  if (right.theElectronOccupancy != 0) {
210      theElectronOccupancy =
211        new G4ElectronOccupancy(*right.theElectronOccupancy);
212  }
213}
214
215////////////////////
216// -- destructor ----
217////////////////////
218G4DynamicParticle::~G4DynamicParticle() {
219
220  //  delete thePreAssignedDecayProducts
221  if (thePreAssignedDecayProducts != 0) delete thePreAssignedDecayProducts;
222  thePreAssignedDecayProducts = 0;
223
224  if (theElectronOccupancy != 0) delete theElectronOccupancy;
225  theElectronOccupancy =0;
226}
227
228
229////////////////////
230// -- operators ----
231////////////////////
232G4DynamicParticle & G4DynamicParticle::operator=(const G4DynamicParticle &right)
233{
234  if (this != &right) {
235    theMomentumDirection = right.theMomentumDirection;
236    theParticleDefinition = right.theParticleDefinition;
237    thePolarization = right.thePolarization;
238    theKineticEnergy = right.theKineticEnergy;
239    theProperTime = right.theProperTime;
240
241    theDynamicalMass = right.theDynamicalMass;
242    theDynamicalCharge = right.theDynamicalCharge;
243    theDynamicalSpin = right.theDynamicalSpin;
244    theDynamicalMagneticMoment = right.theDynamicalMagneticMoment;
245
246    if (theElectronOccupancy != 0) delete theElectronOccupancy;
247    if (right.theElectronOccupancy != 0){
248      theElectronOccupancy =
249             new G4ElectronOccupancy(*right.theElectronOccupancy);
250    } else {
251      theElectronOccupancy = 0;
252    }
253
254    // thePreAssignedDecayProducts must not be copied.
255    thePreAssignedDecayProducts = 0;
256    thePreAssignedDecayTime = -1.0;
257
258    verboseLevel = right.verboseLevel;
259
260    // Primary particle information must be preserved
261    //*** primaryParticle = right.primaryParticle;
262
263    thePDGcode = right.thePDGcode;
264  }
265  return *this;
266}
267
268////////////////////
269void G4DynamicParticle::SetDefinition(const G4ParticleDefinition * aParticleDefinition)
270{
271  // remove preassigned decay
272  if (thePreAssignedDecayProducts != 0) {
273#ifdef G4VERBOSE
274    if (verboseLevel>0) {
275      G4cerr << " G4DynamicParticle::SetDefinition()::"
276             << "!!! Pre-assigned decay products is attached !!!! " << G4endl;
277      G4cerr << "!!! New Definition is " << aParticleDefinition->GetParticleName() 
278             << " !!! " << G4endl;
279      G4cerr << "!!! Pre-assigned decay products will be deleted !!!! " << G4endl;
280    }
281#endif
282    delete thePreAssignedDecayProducts;
283  }
284  thePreAssignedDecayProducts = 0;
285
286  theParticleDefinition = aParticleDefinition;
287
288  // set Dynamic mass/chrge
289  theDynamicalMass = theParticleDefinition->GetPDGMass();
290  theDynamicalCharge = theParticleDefinition->GetPDGCharge();
291  theDynamicalSpin = theParticleDefinition->GetPDGSpin();
292  theDynamicalMagneticMoment = theParticleDefinition->GetPDGMagneticMoment();
293
294  // Set electron orbits
295  if (theElectronOccupancy != 0) delete theElectronOccupancy;
296  theElectronOccupancy =0;
297  //AllocateElectronOccupancy();
298
299}
300
301////////////////////
302G4int G4DynamicParticle::operator==(const G4DynamicParticle &right) const
303{
304  return (this == (G4DynamicParticle *) &right);
305}
306
307////////////////////
308G4int G4DynamicParticle::operator!=(const G4DynamicParticle &right) const
309{
310  return (this != (G4DynamicParticle *) &right);
311}
312
313
314
315////////////////////
316// -- AllocateElectronOccupancy --
317////////////////////
318void  G4DynamicParticle::AllocateElectronOccupancy()
319{
320  const G4ParticleDefinition* particle = GetDefinition();
321
322  if (G4IonTable::IsIon(particle)) {
323    // Only ions can have ElectronOccupancy
324    theElectronOccupancy = new G4ElectronOccupancy();
325
326  } else {
327    theElectronOccupancy = 0;
328
329  }
330}
331
332////////////////////
333// -- methods for setting Energy/Momentum  --
334////////////////////
335void G4DynamicParticle::SetMomentum(const G4ThreeVector &momentum)
336{
337  G4double pModule2 = momentum.mag2();
338  if (pModule2>0.0) {
339    G4double mass = theDynamicalMass;
340    SetMomentumDirection(momentum.unit());
341    SetKineticEnergy(std::sqrt(pModule2 + mass*mass)-mass);
342  } else {
343    SetMomentumDirection(1.0,0.0,0.0);
344    SetKineticEnergy(0.0);
345  }
346}
347
348////////////////////
349void G4DynamicParticle::Set4Momentum(const G4LorentzVector &momentum )
350{
351  G4double pModule2 = momentum.vect().mag2();
352  if (pModule2>0.0) {
353    SetMomentumDirection(momentum.vect().unit());
354    G4double totalenergy = momentum.t();
355    G4double mass2 = totalenergy*totalenergy - pModule2;
356    if(mass2 < EnergyMomentumRelationAllowance*EnergyMomentumRelationAllowance) {
357      theDynamicalMass = 0.;
358      SetKineticEnergy(totalenergy);
359    } else {
360      theDynamicalMass = std::sqrt(mass2);
361      SetKineticEnergy(totalenergy-theDynamicalMass);
362    }
363  } else {
364    SetMomentumDirection(1.0,0.0,0.0);
365    SetKineticEnergy(0.0);
366  }
367}
368
369
370////////////////////
371//  --- Dump Information --
372////////////////////
373void G4DynamicParticle::DumpInfo(G4int mode) const
374{
375#ifdef G4VERBOSE
376  if (theParticleDefinition == 0) {
377    G4cout << " G4DynamicParticle::DumpInfo():: !!!Particle type not defined !!!! " << G4endl;
378  } else {
379    G4cout << " Particle type - " << theParticleDefinition->GetParticleName() << G4endl
380         << "   mass:        " << GetMass()/GeV <<  "[GeV]" <<G4endl
381         << "   charge:      " << GetCharge()/eplus <<  "[e]" <<G4endl
382         << "   Direction x: " << GetMomentumDirection().x() << ", y: "
383         << GetMomentumDirection().y() << ", z: "
384                             << GetMomentumDirection().z() << G4endl
385         << "   Total Momentum = " << GetTotalMomentum() /GeV << "[GeV]" << G4endl
386         << "   Momentum: "    << GetMomentum().x() /GeV << "[GeV]" << ", y: "
387                               << GetMomentum().y() /GeV << "[GeV]" << ", z: "
388                               << GetMomentum().z() /GeV << "[GeV]" << G4endl
389         << "   Total Energy   = " << GetTotalEnergy()/GeV << "[GeV]"  << G4endl
390         << "   Kinetic Energy = " << GetKineticEnergy() /GeV << "[GeV]" << G4endl
391         << " MagneticMoment  [MeV/T]: " << GetMagneticMoment()/MeV*tesla << G4endl
392         << "   ProperTime     = " << GetProperTime() /ns <<  "[ns]" << G4endl;
393
394    if (mode>0) {
395      if( theElectronOccupancy != 0) {
396        theElectronOccupancy->DumpInfo();
397      }
398    }
399  }
400#endif
401}
402
403////////////////////////
404G4double  G4DynamicParticle::GetElectronMass() const
405{
406  static G4double electronMass = 0.0;
407
408  // check if electron exits and get the mass
409  if (electronMass<=0.0) {
410    G4ParticleDefinition* electron = G4ParticleTable::GetParticleTable()->FindParticle("e-");
411    if (electron == 0) {
412      G4Exception("G4DynamicParticle: G4Electron is not defined !!");
413    }
414    electronMass = electron->GetPDGMass();
415  }
416
417  return electronMass;
418}
Note: See TracBrowser for help on using the repository browser.