Geant4  10.01.p02
G4KineticTrack.hh
Go to the documentation of this file.
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 //
28 // $Id: G4KineticTrack.hh,v 1.0 1998/05/20
29 // -----------------------------------------------------------------------------
30 // GEANT 4 class header file
31 //
32 // History: first implementation, A. Feliciello, 20th May 1998
33 // -----------------------------------------------------------------------------
34 
35 #ifndef G4KineticTrack_h
36 #define G4KineticTrack_h 1
37 
38 #include <CLHEP/Units/PhysicalConstants.h>
39 
40 #include "globals.hh"
41 #include "G4ios.hh"
42 
43 
44 #include "Randomize.hh"
45 #include "G4ThreeVector.hh"
46 #include "G4LorentzVector.hh"
47 #include "G4VKineticNucleon.hh"
48 #include "G4Nucleon.hh"
49 #include "G4ParticleDefinition.hh"
50 #include "G4VDecayChannel.hh"
51 
52 // #include "G4Allocator.hh"
53 
55 
56 
57 
58 
59 
61 {
62  public:
63 
65 
67 
68  G4KineticTrack(const G4ParticleDefinition* aDefinition,
69  G4double aFormationTime,
70  const G4ThreeVector& aPosition,
71  const G4LorentzVector& a4Momentum);
73  const G4ThreeVector& aPosition,
74  const G4LorentzVector& a4Momentum);
75 
77 
79 
80  G4int operator==(const G4KineticTrack& right) const;
81 
82  G4int operator!=(const G4KineticTrack& right) const;
83 /*
84  inline void *operator new(size_t);
85  inline void operator delete(void *aTrack);
86 */
87  const G4ParticleDefinition* GetDefinition() const;
88  void SetDefinition(const G4ParticleDefinition* aDefinition);
89 
90  G4double GetFormationTime() const;
91  void SetFormationTime(G4double aFormationTime);
92 
93  const G4ThreeVector& GetPosition() const;
94  void SetPosition(const G4ThreeVector aPosition);
95 
96  const G4LorentzVector& Get4Momentum() const;
97  void Set4Momentum(const G4LorentzVector& a4Momentum);
98  void Update4Momentum(G4double aEnergy); // update E and p, not changing mass
99  void Update4Momentum(const G4ThreeVector & aMomentum); // idem
100  void SetTrackingMomentum(const G4LorentzVector& a4Momentum);
101  void UpdateTrackingMomentum(G4double aEnergy); // update E and p, not changing mass
102  void UpdateTrackingMomentum(const G4ThreeVector & aMomentum); // idem
103 
104  const G4LorentzVector& GetTrackingMomentum() const;
105 
107 
108  void Hit();
109  void SetNucleon(G4Nucleon * aN) {theNucleon = aN;}
110 
111  G4bool IsParticipant() const;
112 
114 
115  // LB move to public (before was private) LB
116  G4double* GetActualWidth() const;
117 
118  G4double GetActualMass() const;
119  G4int GetnChannels() const;
120 
121 // position relativ to nucleus "state"
124 
125  CascadeState SetState(const CascadeState new_state);
126  CascadeState GetState() const;
127  void SetProjectilePotential(const G4double aPotential);
129 
130 
131  private:
132 
133 
134  void SetnChannels(const G4int aChannel);
135 
136  void SetActualWidth(G4double* anActualWidth);
137 
139 
141  const G4double* m_ij) const;
142 
143  G4double IntegrateCMMomentum(const G4double lowerLimit) const;
144 
145  G4double IntegrateCMMomentum(const G4double lowerLimit ,const G4double polemass) const;
146 
148 
149  public:
150 
151  G4double BrWig(const G4double Gamma,
152  const G4double rmass,
153  const G4double mass) const;
154 
155 private:
156  G4double IntegrandFunction1 (G4double xmass) const;
157  G4double IntegrandFunction2 (G4double xmass) const;
158  G4double IntegrandFunction3 (G4double xmass) const;
159  G4double IntegrandFunction4 (G4double xmass) const;
160 public:
161  // friend G4double IntegrandFunction3 (G4double xmass);
162 
163  // friend G4double IntegrandFunction4 (G4double xmass);
164 
165  private:
166 
168 
170 
172 
176 
178 
180 
182 
184 
185  // Temporary storage for daughter masses and widths
186  // (needed because Integrand Function cannot take > 1 argument)
189 
191 
193 };
194 
195 // extern G4Allocator<G4KineticTrack> theKTAllocator;
196 
197 
198 // Class G4KineticTrack
199 /*
200 inline void * G4KineticTrack::operator new(size_t)
201 {
202  void * aT;
203  aT = (void *) theKTAllocator.MallocSingle();
204  return aT;
205 }
206 
207 inline void G4KineticTrack::operator delete(void * aT)
208 {
209  theKTAllocator.FreeSingle((G4KineticTrack *) aT);
210 }
211 */
212 
214 {
215  return theDefinition;
216 }
217 
219 {
220  theDefinition = aDefinition;
221 }
222 
223 
224 
226 {
227  return theFormationTime;
228 }
229 
230 inline void G4KineticTrack::SetFormationTime(G4double aFormationTime)
231 {
232  theFormationTime = aFormationTime;
233 }
234 
235 
236 
238 {
239  return thePosition;
240 }
241 
242 inline void G4KineticTrack::SetPosition(const G4ThreeVector aPosition)
243 {
244  thePosition = aPosition;
245 }
246 
247 
249 {
250  return theTotal4Momentum;
251 }
252 
254 {
255  return the4Momentum;
256 }
257 
258 inline void G4KineticTrack::Set4Momentum(const G4LorentzVector& a4Momentum)
259 {
260 // set the4Momentum and update theTotal4Momentum
261 
262  theTotal4Momentum=a4Momentum;
265 }
266 
268 {
269 // update the4Momentum with aEnergy at constant mass (the4Momentum.mag()
270 // updates theTotal4Momentum as well.
271  G4double newP(0);
272  G4double mass2=theTotal4Momentum.mag2();
273  if ( sqr(aEnergy) > mass2 )
274  {
275  newP = std::sqrt(sqr(aEnergy) - mass2 );
276  } else
277  {
278  aEnergy=std::sqrt(mass2);
279  }
280  Set4Momentum(G4LorentzVector(newP*the4Momentum.vect().unit(), aEnergy));
281 }
282 
283 inline void G4KineticTrack::Update4Momentum(const G4ThreeVector & aMomentum)
284 {
285 // update the4Momentum with aMomentum at constant mass (the4Momentum.mag()
286 // updates theTotal4Momentum as well.
287  G4double newE=std::sqrt(theTotal4Momentum.mag2() + aMomentum.mag2());
288  Set4Momentum(G4LorentzVector(aMomentum, newE));
289 }
290 
292 {
293 // set the4Momentum and update theTotal4Momentum, keep the mass of aMomentum
294 
295  the4Momentum = aMomentum;
297 // keep mass of aMomentum for the total momentum
298  G4double mass2 = aMomentum.mag2();
299  G4double p2=theTotal4Momentum.vect().mag2();
300  theTotal4Momentum.setE(std::sqrt(mass2+p2));
301 }
302 
304 {
305 // update the4Momentum with aEnergy at constant mass (the4Momentum.mag()
306 // updates theTotal4Momentum as well.
307  G4double newP(0);
308  G4double mass2=theTotal4Momentum.mag2();
309  if ( sqr(aEnergy) > mass2 )
310  {
311  newP = std::sqrt(sqr(aEnergy) - mass2 );
312  } else
313  {
314  aEnergy=std::sqrt(mass2);
315  }
316  SetTrackingMomentum(G4LorentzVector(newP*the4Momentum.vect().unit(), aEnergy));
317 }
318 
320 {
321 // update the4Momentum with aMomentum at constant mass (the4Momentum.mag()
322 // updates theTotal4Momentum as well.
323  G4double newE=std::sqrt(theTotal4Momentum.mag2() + aMomentum.mag2());
324  SetTrackingMomentum(G4LorentzVector(aMomentum, newE));
325 }
326 
327 
328 
329 
331 {
332  return std::sqrt(std::abs(the4Momentum.mag2()));
333 }
334 
335 
336 
338 {
339  return nChannels;
340 }
341 
342 inline void G4KineticTrack::SetnChannels(const G4int numberOfChannels)
343 {
344  nChannels = numberOfChannels;
345 }
346 
347 
348 
350 {
351  return theActualWidth;
352 }
353 
354 inline void G4KineticTrack::SetActualWidth(G4double* anActualWidth)
355 {
356  theActualWidth = anActualWidth;
357 }
358 
359 
360 
362 {
363  G4int index;
364  G4double theTotalActualWidth = 0.0;
365  for (index = nChannels - 1; index >= 0; index--)
366  {
367  theTotalActualWidth += theActualWidth[index];
368  }
369  return theTotalActualWidth;
370 }
371 
372 
373 
375 {
376  G4double theTotalActualWidth = this->EvaluateTotalActualWidth();
377  G4double tau = CLHEP::hbar_Planck * (-1.0 / theTotalActualWidth);
378  G4double theResidualLifetime = tau * std::log(G4UniformRand());
379  return theResidualLifetime*the4Momentum.gamma();
380 }
381 
382 
383 
385  const G4double* m_ij) const
386 {
387  G4double theCMMomentum;
388  if((m_ij[0]+m_ij[1])<mass)
389  theCMMomentum = 1 / (2 * mass) *
390  std::sqrt (((mass * mass) - (m_ij[0] + m_ij[1]) * (m_ij[0] + m_ij[1])) *
391  ((mass * mass) - (m_ij[0] - m_ij[1]) * (m_ij[0] - m_ij[1])));
392  else
393  theCMMomentum=0.;
394 
395  return theCMMomentum;
396 }
397 
398 inline G4double G4KineticTrack::BrWig(const G4double Gamma, const G4double rmass, const G4double mass) const
399 {
400  G4double Norm = CLHEP::twopi;
401  return (Gamma/((mass-rmass)*(mass-rmass)+Gamma*Gamma/4.))/Norm;
402 }
403 
404 inline
406 {
407  if(theNucleon)
408  {
409  theNucleon->Hit(1);
410  }
411 }
412 
413 inline
415 {
416  if(!theNucleon) return true;
417  return theNucleon->AreYouHit();
418 }
419 
420 inline
422 {
423  return theStateToNucleus;
424 }
425 
426 inline
428 {
430  theStateToNucleus=new_state;
431  return old_state;
432 }
433 
434 inline
436 {
437  theProjectilePotential = aPotential;
438 }
439 inline
441 {
442  return theProjectilePotential;
443 }
444 
445 #endif
446 
447 
448 
void SetDefinition(const G4ParticleDefinition *aDefinition)
G4bool IsParticipant() const
void Update4Momentum(G4double aEnergy)
G4double IntegrateCMMomentum2() const
G4int operator==(const G4KineticTrack &right) const
CascadeState GetState() const
G4double theFormationTime
G4int GetnChannels() const
G4LorentzVector theTotal4Momentum
CLHEP::Hep3Vector G4ThreeVector
G4double * theDaughterWidth
const G4ThreeVector & GetPosition() const
G4KineticTrackVector * Decay()
void SetFormationTime(G4double aFormationTime)
void UpdateTrackingMomentum(G4double aEnergy)
G4double GetActualMass() const
G4ThreeVector thePosition
void SetProjectilePotential(const G4double aPotential)
int G4int
Definition: G4Types.hh:78
G4double IntegrateCMMomentum(const G4double lowerLimit) const
G4double IntegrandFunction3(G4double xmass) const
G4double * theDaughterMass
G4LorentzVector theFermi3Momentum
G4double IntegrandFunction2(G4double xmass) const
G4double * GetActualWidth() const
G4int operator!=(const G4KineticTrack &right) const
#define G4UniformRand()
Definition: Randomize.hh:93
void SetActualWidth(G4double *anActualWidth)
CascadeState SetState(const CascadeState new_state)
G4double IntegrandFunction4(G4double xmass) const
G4bool nucleon(G4int ityp)
G4double GetFormationTime() const
bool G4bool
Definition: G4Types.hh:79
G4LorentzVector the4Momentum
void SetNucleon(G4Nucleon *aN)
G4bool AreYouHit() const
Definition: G4Nucleon.hh:97
void SetPosition(const G4ThreeVector aPosition)
const G4double p2
G4double EvaluateTotalActualWidth()
G4KineticTrack & operator=(const G4KineticTrack &right)
void Set4Momentum(const G4LorentzVector &a4Momentum)
void SetnChannels(const G4int aChannel)
G4double EvaluateCMMomentum(const G4double mass, const G4double *m_ij) const
const G4LorentzVector & GetTrackingMomentum() const
G4double * theActualWidth
G4double IntegrandFunction1(G4double xmass) const
void SetTrackingMomentum(const G4LorentzVector &a4Momentum)
void Hit(G4VSplitableHadron *aHit)
Definition: G4Nucleon.hh:90
G4Nucleon * theNucleon
G4double theActualMass
T sqr(const T &x)
Definition: templates.hh:145
double G4double
Definition: G4Types.hh:76
const G4LorentzVector & Get4Momentum() const
G4double theProjectilePotential
G4double SampleResidualLifetime()
const G4ParticleDefinition * GetDefinition() const
G4double BrWig(const G4double Gamma, const G4double rmass, const G4double mass) const
G4double GetProjectilePotential() const
CLHEP::HepLorentzVector G4LorentzVector
const G4ParticleDefinition * theDefinition
CascadeState theStateToNucleus