Geant4 Cross Reference |
1 // 1 2 // ******************************************* 3 // * License and Disclaimer 4 // * 5 // * The Geant4 software is copyright of th 6 // * the Geant4 Collaboration. It is provided 7 // * conditions of the Geant4 Software License 8 // * LICENSE and available at http://cern.ch/ 9 // * include a list of copyright holders. 10 // * 11 // * Neither the authors of this software syst 12 // * institutes,nor the agencies providing fin 13 // * work make any representation or warran 14 // * regarding this software system or assum 15 // * use. Please see the license in the file 16 // * for the full disclaimer and the limitatio 17 // * 18 // * This code implementation is the result 19 // * technical work of the GEANT4 collaboratio 20 // * By using, copying, modifying or distri 21 // * any work based on the software) you ag 22 // * use in resulting scientific publicati 23 // * acceptance of all terms of the Geant4 Sof 24 // ******************************************* 25 // 26 // 16.05.17 V. Grichine first implementation 27 28 29 #include "G4NeutronElectronElXsc.hh" 30 #include "G4PhysicalConstants.hh" 31 #include "G4SystemOfUnits.hh" 32 #include "G4DynamicParticle.hh" 33 #include "G4ParticleDefinition.hh" 34 #include "G4ParticleTable.hh" 35 #include "G4IonTable.hh" 36 #include "G4HadTmpUtil.hh" 37 #include "G4NistManager.hh" 38 // #include "G4Integrator.hh" 39 40 #include "G4PhysicsLogVector.hh" 41 #include "G4PhysicsTable.hh" 42 43 #include "G4Neutron.hh" 44 #include "G4Electron.hh" 45 46 47 using namespace std; 48 using namespace CLHEP; 49 50 G4NeutronElectronElXsc::G4NeutronElectronElXsc 51 : G4VCrossSectionDataSet("NuElectronCcXsc") 52 { 53 // neutron magneton squared 54 55 fM = neutron_mass_c2; // neutron mass 56 fM2 = fM*fM; 57 fme = electron_mass_c2; 58 fme2 = fme*fme; 59 fMv2 = 0.7056*GeV*GeV; 60 fee = fme; 61 fee2 = fee*fee; 62 fAm = 0.001; 63 64 fCofXsc = pi*fine_structure_const*fine_stru 65 fCofXsc *= 3.6481; // neutron Fm^2(0) 66 fCofXsc /= fM*fM; 67 68 // G4cout<<"fCofXsc = "<<fCofXsc*GeV/cm2<<" 69 70 // G4cout<<"hbarc = "<<hbarc/MeV/fermi<<" Me 71 72 fCutEnergy = 0.; // default value 73 74 fEnergyBin = 200; 75 fMinEnergy = 1.*MeV; 76 fMaxEnergy = 10000.*GeV; 77 78 fEnergyXscVector = new G4PhysicsLogVector(fM 79 80 for( G4int iTkin = 0; iTkin < fEnergyBin; iT 81 82 fBiasingFactor = 1.; 83 84 // Initialise(); 85 } 86 87 G4NeutronElectronElXsc::~G4NeutronElectronElXs 88 { 89 if( fEnergyXscVector ) 90 { 91 delete fEnergyXscVector; 92 fEnergyXscVector = 0; 93 } 94 } 95 96 ////////////////////////////////////////////// 97 // 98 // For neutrons in the precalculated energy in 99 100 G4bool 101 G4NeutronElectronElXsc::IsElementApplicable( c 102 { 103 G4bool result = false; 104 G4String pName = aPart->GetDefinition()->Get 105 G4double Tkin = aPart->GetKineticEnergy(); 106 107 if( pName == "neutron" && 108 Tkin >= fMinEnergy && 109 Tkin <= fMaxEnergy ) result = true; 110 111 return result; 112 } 113 114 ////////////////////////////////////////////// 115 116 void G4NeutronElectronElXsc::Initialise() 117 { 118 G4int iTkin; 119 G4double Tkin, rosxsc, xsc, delta, err=1.e-5 120 const G4ThreeVector mDir = G4ThreeVector(0., 121 const G4ParticleDefinition* pD = G4Neutron:: 122 G4Material* mat = G4NistManager::Instance()- 123 124 G4DynamicParticle dP; 125 126 for( iTkin = 0; iTkin < fEnergyBin; iTkin++) 127 { 128 Tkin = fEnergyXscVector->GetLowEdgeEnergy 129 dP = G4DynamicParticle(pD, mDir, Tkin); 130 rosxsc = GetRosenbluthXsc(&dP, 1, mat); 131 fEnergyXscVector->PutValue(iTkin, rosxsc); 132 xsc= fEnergyXscVector->Value(Tkin); // xsc 133 delta = 0.5*std::abs( (rosxsc-xsc) )/(rosx 134 if(delta > err) G4cout<<Tkin/GeV<<" GeV, r 135 } 136 return; 137 } 138 139 ////////////////////////////////////////////// 140 141 G4double G4NeutronElectronElXsc:: 142 GetElementCrossSection(const G4DynamicParticle 143 const G4Material*) 144 { 145 G4double result = 0., Tkin; 146 147 Tkin = aPart->GetKineticEnergy(); 148 149 result = fEnergyXscVector->Value(Tkin); 150 151 result *= ZZ; // incoherent sum over all e 152 153 result *= fBiasingFactor; 154 155 return result; 156 } 157 158 ////////////////////////////////////////////// 159 // 160 // Integration of the Rosenbluth differential 161 162 G4double G4NeutronElectronElXsc:: 163 GetRosenbluthXsc(const G4DynamicParticle* aPar 164 const G4Material*) 165 { 166 G4double result = 0., momentum; 167 168 fee = aPart->GetTotalEnergy()*fme/fM; 169 fee2 = fee*fee; 170 momentum = sqrt( fee2 - fme2 ); 171 fAm = CalculateAm(momentum); 172 173 // G4Integrator<G4NeutronElectronElXsc, G4do 174 175 // result = integral.Legendre96( this, &G4Ne 176 177 result *= fCofXsc; 178 179 result *= ZZ; // incoherent sum over all e 180 181 return result; 182 } 183 184 ///////////////////////////////////////// 185 // 186 // Rosenbluth relation in the neutron rest fra 187 // x = sin^2(theta/2), theta is the electron s 188 // Magnetic form factor in the dipole approxim 189 190 G4double G4NeutronElectronElXsc::XscIntegrand( 191 { 192 G4double result = 1., q2, znq2, znf, znf2, z 193 194 znq2 = 1. + 2.*fee*x/fM; 195 196 q2 = 4.*fee2*x/znq2; 197 198 znf = 1 + q2/fMv2; 199 znf2 = znf*znf; 200 znf4 = znf2*znf2; 201 202 result /= ( x + fAm )*znq2*znq2*znf4; 203 204 result *= ( 1 - x )/( 1 + q2/4./fM2 ) + 2.*x 205 206 return result; 207 } 208 209 ////////////////////////////////////////////// 210 // 211 // Rosenbluth xsc in microbarn from 1*MeV to 1 212 213 const G4double G4NeutronElectronElXsc::fXscArr 214 1.52681, 1.54903, 1.57123, 1.59341, 1.61556, 1 215 1.72601, 1.74805, 1.77007, 1.79208, 1.81407, 1 216 1.92385, 1.94578, 1.96771, 1.98962, 2.01154, 2 217 2.12105, 2.14295, 2.16485, 2.18675, 2.20865, 2 218 2.31807, 2.33992, 2.36173, 2.38351, 2.40524, 2 219 2.51262, 2.53369, 2.55457, 2.57524, 2.59565, 2 220 2.69281, 2.71104, 2.72881, 2.74607, 2.76282, 2 221 2.83811, 2.85139, 2.86408, 2.87616, 2.88764, 2 222 2.93641, 2.94453, 2.95213, 2.95924, 2.96588, 2 223 2.99268, 2.9969, 3.00078, 3.00435, 3.00761, 3. 224 3.02003, 3.02185, 3.02347, 3.02491, 3.02619, 3 225 3.03049, 3.03099, 3.03139, 3.03169, 3.03191, 3 226 3.03177, 3.03152, 3.0312, 3.03081, 3.03034, 3. 227 3.02684, 3.02588, 3.02482, 3.02365, 3.02237, 3 228 3.01389, 3.01169, 3.00929, 3.00666, 3.00379, 3 229 2.98487, 2.97996, 2.97459, 2.9687, 2.96226, 2. 230 2.91965, 2.90858, 2.89649, 2.88329, 2.86889, 2 231 2.77594, 2.7526, 2.72754, 2.70071, 2.67209, 2. 232 2.50336, 2.46504, 2.42548, 2.38484, 2.34328, 2 233 2.12735, 2.08345, 2.03954, 1.99569, 1.95191, 1 234 1.7348, 1.69171, 1.64869, 1.60575, 1.56286, 1. 235 1.34877, 1.30596, 1.26314, 1.22031, 1.17746, 1 236 0.962892, 0.919908 }; 237