ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4RToEConvForPositron.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4RToEConvForPositron.cc
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 //
29 // --------------------------------------------------------------
30 // GEANT 4 class implementation file/ History:
31 // 5 Oct. 2002, H.Kuirashige : Structure created based on object model
32 // --------------------------------------------------------------
33 
34 #include "G4RToEConvForPositron.hh"
35 #include "G4ParticleDefinition.hh"
36 #include "G4ParticleTable.hh"
37 #include "G4Material.hh"
38 #include "G4PhysicsLogVector.hh"
39 
40 #include "G4ios.hh"
41 #include "G4PhysicalConstants.hh"
42 #include "G4SystemOfUnits.hh"
43 
46  Mass(0.0),
47  Z(-1.),
48  taul(0.0),
49  ionpot(0.0),
50  ionpotlog(-1.0e-10),
51  bremfactor(0.1)
52 {
54  if (theParticle ==0) {
55 #ifdef G4VERBOSE
56  if (GetVerboseLevel()>0) {
57  G4cout << " G4RToEConvForPositron::G4RToEConvForPositron() ";
58  G4cout << " Positron is not defined !!" << G4endl;
59  }
60 #endif
61  } else {
63  }
64 }
65 
67 {
68 }
69 
70 
71 
72 
73 // **********************************************************************
74 // ************************* ComputeLoss ********************************
75 // **********************************************************************
77  G4double KineticEnergy)
78 {
79  const G4double cbr1=0.02, cbr2=-5.7e-5, cbr3=1., cbr4=0.072;
80  const G4double Tlow=10.*keV, Thigh=1.*GeV;
81 
82  // calculate dE/dx for electrons
83  if( std::fabs(AtomicNumber-Z)>0.1 ) {
84  Z = AtomicNumber;
85  taul = Tlow/Mass;
86  ionpot = 1.6e-5*MeV*std::exp(0.9*std::log(Z))/Mass;
87  ionpotlog = std::log(ionpot);
88  }
89 
90  G4double tau = KineticEnergy/Mass;
91  G4double dEdx;
92 
93  if(tau<taul)
94  {
95  G4double t1 = taul+1.;
96  G4double t2 = taul+2.;
97  G4double tsq = taul*taul;
98  G4double beta2 = taul*t2/(t1*t1);
99  G4double f = 2.*std::log(taul)
100  -(6.*taul+1.5*tsq-taul*(1.-tsq/3.)/t2-tsq*(0.5-tsq/12.)/
101  (t2*t2))/(t1*t1);
102  dEdx = (std::log(2.*taul+4.)-2.*ionpotlog+f)/beta2;
103  dEdx = twopi_mc2_rcl2*Z*dEdx;
104  G4double clow = dEdx*std::sqrt(taul);
105  dEdx = clow/std::sqrt(KineticEnergy/Mass);
106 
107  } else {
108  G4double t1 = tau+1.;
109  G4double t2 = tau+2.;
110  G4double tsq = tau*tau;
111  G4double beta2 = tau*t2/(t1*t1);
112  G4double f = 2.*std::log(tau)
113  - (6.*tau+1.5*tsq-tau*(1.-tsq/3.)/t2-tsq*(0.5-tsq/12.)/
114  (t2*t2))/(t1*t1);
115  dEdx = (std::log(2.*tau+4.)-2.*ionpotlog+f)/beta2;
116  dEdx = twopi_mc2_rcl2*Z*dEdx;
117 
118  // loss from bremsstrahlung follows
119  G4double cbrem = (cbr1+cbr2*Z)
120  *(cbr3+cbr4*std::log(KineticEnergy/Thigh));
121  cbrem = Z*(Z+1.)*cbrem*tau/beta2;
122  cbrem *= bremfactor ;
123  dEdx += twopi_mc2_rcl2*cbrem;
124  }
125  return dEdx;
126 }
127 
128