Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
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 // $Id: G4RToEConvForPositron.cc 70745 2013-06-05 10:54:00Z gcosmo $
28 //
29 //
30 // --------------------------------------------------------------
31 // GEANT 4 class implementation file/ History:
32 // 5 Oct. 2002, H.Kuirashige : Structure created based on object model
33 // --------------------------------------------------------------
34 
35 #include "G4RToEConvForPositron.hh"
36 #include "G4ParticleDefinition.hh"
37 #include "G4ParticleTable.hh"
38 #include "G4Material.hh"
39 #include "G4PhysicsLogVector.hh"
40 
41 #include "G4ios.hh"
42 #include "G4PhysicalConstants.hh"
43 #include "G4SystemOfUnits.hh"
44 
47  Mass(0.0),
48  Z(-1.),
49  taul(0.0),
50  ionpot(0.0),
51  ionpotlog(-1.0e-10),
52  bremfactor(0.1)
53 {
55  if (theParticle ==0) {
56 #ifdef G4VERBOSE
57  if (GetVerboseLevel()>0) {
58  G4cout << " G4RToEConvForPositron::G4RToEConvForPositron() ";
59  G4cout << " Positron is not defined !!" << G4endl;
60  }
61 #endif
62  } else {
64  }
65 }
66 
68 {
69 }
70 
71 
72 
73 
74 // **********************************************************************
75 // ************************* ComputeLoss ********************************
76 // **********************************************************************
78  G4double KineticEnergy)
79 {
80  const G4double cbr1=0.02, cbr2=-5.7e-5, cbr3=1., cbr4=0.072;
81  const G4double Tlow=10.*keV, Thigh=1.*GeV;
82 
83  // calculate dE/dx for electrons
84  if( std::fabs(AtomicNumber-Z)>0.1 ) {
85  Z = AtomicNumber;
86  taul = Tlow/Mass;
87  ionpot = 1.6e-5*MeV*std::exp(0.9*std::log(Z))/Mass;
88  ionpotlog = std::log(ionpot);
89  }
90 
91  G4double tau = KineticEnergy/Mass;
92  G4double dEdx;
93 
94  if(tau<taul)
95  {
96  G4double t1 = taul+1.;
97  G4double t2 = taul+2.;
98  G4double tsq = taul*taul;
99  G4double beta2 = taul*t2/(t1*t1);
100  G4double f = 2.*std::log(taul)
101  -(6.*taul+1.5*tsq-taul*(1.-tsq/3.)/t2-tsq*(0.5-tsq/12.)/
102  (t2*t2))/(t1*t1);
103  dEdx = (std::log(2.*taul+4.)-2.*ionpotlog+f)/beta2;
104  dEdx = twopi_mc2_rcl2*Z*dEdx;
105  G4double clow = dEdx*std::sqrt(taul);
106  dEdx = clow/std::sqrt(KineticEnergy/Mass);
107 
108  } else {
109  G4double t1 = tau+1.;
110  G4double t2 = tau+2.;
111  G4double tsq = tau*tau;
112  G4double beta2 = tau*t2/(t1*t1);
113  G4double f = 2.*std::log(tau)
114  - (6.*tau+1.5*tsq-tau*(1.-tsq/3.)/t2-tsq*(0.5-tsq/12.)/
115  (t2*t2))/(t1*t1);
116  dEdx = (std::log(2.*tau+4.)-2.*ionpotlog+f)/beta2;
117  dEdx = twopi_mc2_rcl2*Z*dEdx;
118 
119  // loss from bremsstrahlung follows
120  G4double cbrem = (cbr1+cbr2*Z)
121  *(cbr3+cbr4*std::log(KineticEnergy/Thigh));
122  cbrem = Z*(Z+1.)*cbrem*tau/beta2;
123  cbrem *= bremfactor ;
124  dEdx += twopi_mc2_rcl2*cbrem;
125  }
126  return dEdx;
127 }
128 
129 
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
TTree * t1
Definition: plottest35.C:26
static G4ParticleTable * GetParticleTable()
static constexpr double MeV
Definition: G4SIunits.hh:214
static constexpr double keV
Definition: G4SIunits.hh:216
static constexpr double twopi_mc2_rcl2
#define G4endl
Definition: G4ios.hh:61
G4double GetPDGMass() const
Float_t Z
double G4double
Definition: G4Types.hh:76
TTree * t2
Definition: plottest35.C:36
virtual G4double ComputeLoss(G4double AtomicNumber, G4double KineticEnergy)
G4GLOB_DLL std::ostream G4cout
const G4ParticleDefinition * theParticle
static constexpr double GeV
Definition: G4SIunits.hh:217