00001 // 00002 // ******************************************************************** 00003 // * License and Disclaimer * 00004 // * * 00005 // * The Geant4 software is copyright of the Copyright Holders of * 00006 // * the Geant4 Collaboration. It is provided under the terms and * 00007 // * conditions of the Geant4 Software License, included in the file * 00008 // * LICENSE and available at http://cern.ch/geant4/license . These * 00009 // * include a list of copyright holders. * 00010 // * * 00011 // * Neither the authors of this software system, nor their employing * 00012 // * institutes,nor the agencies providing financial support for this * 00013 // * work make any representation or warranty, express or implied, * 00014 // * regarding this software system or assume any liability for its * 00015 // * use. Please see the license in the file LICENSE and URL above * 00016 // * for the full disclaimer and the limitation of liability. * 00017 // * * 00018 // * This code implementation is the result of the scientific and * 00019 // * technical work of the GEANT4 collaboration. * 00020 // * By using, copying, modifying or distributing the software (or * 00021 // * any work based on the software) you agree to acknowledge its * 00022 // * use in resulting scientific publications, and indicate your * 00023 // * acceptance of all terms of the Geant4 Software license. * 00024 // ******************************************************************** 00025 // 00026 // $Id: G4hPairProductionModel.cc 66996 2013-01-29 14:50:52Z gcosmo $ 00027 // 00028 // ------------------------------------------------------------------- 00029 // 00030 // GEANT4 Class file 00031 // 00032 // 00033 // File name: G4hPairProductionModel 00034 // 00035 // Author: Vladimir Ivanchenko on base of G4MuPairProductionModel 00036 // 00037 // Creation date: 28.02.2008 00038 // 00039 // Modifications: 00040 // 00041 00042 // 00043 // Class Description: 00044 // 00045 // 00046 // ------------------------------------------------------------------- 00047 // 00048 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 00049 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 00050 00051 #include "G4hPairProductionModel.hh" 00052 #include "G4PhysicalConstants.hh" 00053 00054 using namespace std; 00055 00056 G4hPairProductionModel::G4hPairProductionModel(const G4ParticleDefinition* p, 00057 const G4String& nam) 00058 : G4MuPairProductionModel(p, nam) 00059 {} 00060 00061 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 00062 00063 G4hPairProductionModel::~G4hPairProductionModel() 00064 {} 00065 00066 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 00067 00068 G4double G4hPairProductionModel::ComputeDMicroscopicCrossSection( 00069 G4double tkin, 00070 G4double Z, 00071 G4double pairEnergy) 00072 // differential cross section 00073 { 00074 G4double bbbtf= 183. ; 00075 G4double bbbh = 202.4 ; 00076 G4double g1tf = 1.95e-5 ; 00077 G4double g2tf = 5.3e-5 ; 00078 G4double g1h = 4.4e-5 ; 00079 G4double g2h = 4.8e-5 ; 00080 00081 G4double totalEnergy = tkin + particleMass; 00082 G4double residEnergy = totalEnergy - pairEnergy; 00083 G4double massratio = particleMass/electron_mass_c2 ; 00084 G4double massratio2 = massratio*massratio ; 00085 G4double cross = 0.; 00086 00087 SetCurrentElement(Z); 00088 00089 G4double c3 = 0.75*sqrte*particleMass; 00090 if (residEnergy <= c3*z13) return cross; 00091 00092 G4double c7 = 4.*electron_mass_c2; 00093 G4double c8 = 6.*particleMass*particleMass; 00094 G4double alf = c7/pairEnergy; 00095 G4double a3 = 1. - alf; 00096 if (a3 <= 0.) return cross; 00097 00098 // zeta calculation 00099 G4double bbb,g1,g2; 00100 if( Z < 1.5 ) { bbb = bbbh ; g1 = g1h ; g2 = g2h ; } 00101 else { bbb = bbbtf; g1 = g1tf; g2 = g2tf; } 00102 00103 G4double zeta = 0; 00104 G4double zeta1 = 0.073*log(totalEnergy/(particleMass+g1*z23*totalEnergy))-0.26; 00105 if ( zeta1 > 0.) 00106 { 00107 G4double zeta2 = 0.058*log(totalEnergy/(particleMass+g2*z13*totalEnergy))-0.14; 00108 zeta = zeta1/zeta2 ; 00109 } 00110 00111 G4double z2 = Z*(Z+zeta); 00112 G4double screen0 = 2.*electron_mass_c2*sqrte*bbb/(z13*pairEnergy); 00113 G4double a0 = totalEnergy*residEnergy; 00114 G4double a1 = pairEnergy*pairEnergy/a0; 00115 G4double bet = 0.5*a1; 00116 G4double xi0 = 0.25*massratio2*a1; 00117 G4double del = c8/a0; 00118 00119 G4double rta3 = sqrt(a3); 00120 G4double tmnexp = alf/(1. + rta3) + del*rta3; 00121 if(tmnexp >= 1.0) return cross; 00122 00123 G4double tmn = log(tmnexp); 00124 G4double sum = 0.; 00125 00126 // Gaussian integration in ln(1-ro) ( with 8 points) 00127 for (G4int i=0; i<8; i++) 00128 { 00129 G4double a4 = exp(tmn*xgi[i]); // a4 = (1.-asymmetry) 00130 G4double a5 = a4*(2.-a4) ; 00131 G4double a6 = 1.-a5 ; 00132 G4double a7 = 1.+a6 ; 00133 G4double a9 = 3.+a6 ; 00134 G4double xi = xi0*a5 ; 00135 G4double xii = 1./xi ; 00136 G4double xi1 = 1.+xi ; 00137 G4double screen = screen0*xi1/a5 ; 00138 G4double yeu = 5.-a6+4.*bet*a7 ; 00139 G4double yed = 2.*(1.+3.*bet)*log(3.+xii)-a6-a1*(2.-a6) ; 00140 G4double ye1 = 1.+yeu/yed ; 00141 G4double ale=log(bbb/z13*sqrt(xi1*ye1)/(1.+screen*ye1)) ; 00142 G4double cre = 0.5*log(1.+2.25*z23*xi1*ye1/massratio2) ; 00143 G4double be; 00144 00145 if (xi <= 1.e3) be = ((2.+a6)*(1.+bet)+xi*a9)*log(1.+xii)+(a5-bet)/xi1-a9; 00146 else be = (3.-a6+a1*a7)/(2.*xi); 00147 00148 G4double fe = (ale-cre)*be; 00149 if ( fe < 0.) fe = 0. ; 00150 00151 G4double ymu = 4.+a6 +3.*bet*a7 ; 00152 G4double ymd = a7*(1.5+a1)*log(3.+xi)+1.-1.5*a6 ; 00153 G4double ym1 = 1.+ymu/ymd ; 00154 G4double alm_crm = log(bbb*massratio/(1.5*z23*(1.+screen*ym1))); 00155 G4double a10,bm; 00156 if ( xi >= 1.e-3) 00157 { 00158 a10 = (1.+a1)*a5 ; 00159 bm = (a7*(1.+1.5*bet)-a10*xii)*log(xi1)+xi*(a5-bet)/xi1+a10; 00160 } else { 00161 bm = (5.-a6+bet*a9)*(xi/2.); 00162 } 00163 00164 G4double fm = alm_crm*bm; 00165 if ( fm < 0.) fm = 0. ; 00166 00167 sum += wgi[i]*a4*(fe+fm/massratio2); 00168 } 00169 00170 cross = -tmn*sum*factorForCross*z2*residEnergy/(totalEnergy*pairEnergy); 00171 00172 return cross; 00173 } 00174 00175 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......