#include <G4KokoulinMuonNuclearXS.hh>
Inheritance diagram for G4KokoulinMuonNuclearXS:
Public Member Functions | |
G4KokoulinMuonNuclearXS () | |
virtual | ~G4KokoulinMuonNuclearXS () |
G4bool | IsElementApplicable (const G4DynamicParticle *particle, G4int Z, const G4Material *) |
G4double | GetElementCrossSection (const G4DynamicParticle *particle, G4int Z, const G4Material *) |
void | BuildCrossSectionTable () |
G4double | ComputeDDMicroscopicCrossSection (G4double incidentKE, G4double, G4double AtomicWeight, G4double epsilon) |
Definition at line 55 of file G4KokoulinMuonNuclearXS.hh.
G4KokoulinMuonNuclearXS::G4KokoulinMuonNuclearXS | ( | ) |
Definition at line 47 of file G4KokoulinMuonNuclearXS.cc.
References BuildCrossSectionTable().
00048 :G4VCrossSectionDataSet("KokoulinMuonNuclearXS"), theCrossSectionTable(0), 00049 LowestKineticEnergy(1*GeV), HighestKineticEnergy(1*PeV), 00050 TotBin(60), CutFixed(0.2*GeV) 00051 { 00052 BuildCrossSectionTable(); 00053 }
G4KokoulinMuonNuclearXS::~G4KokoulinMuonNuclearXS | ( | ) | [virtual] |
Definition at line 55 of file G4KokoulinMuonNuclearXS.cc.
References G4PhysicsTable::clearAndDestroy().
00056 { 00057 if (theCrossSectionTable) { 00058 theCrossSectionTable->clearAndDestroy(); 00059 delete theCrossSectionTable; 00060 } 00061 }
void G4KokoulinMuonNuclearXS::BuildCrossSectionTable | ( | ) |
Definition at line 71 of file G4KokoulinMuonNuclearXS.cc.
References G4PhysicsVector::Energy(), G4Element::GetElementTable(), G4Element::GetNumberOfElements(), G4PhysicsTable::insertAt(), and G4PhysicsVector::PutValue().
Referenced by G4KokoulinMuonNuclearXS().
00072 { 00073 G4double lowEdgeEnergy; 00074 G4PhysicsLogVector* ptrVector; 00075 00076 G4int nEl = G4Element::GetNumberOfElements(); 00077 const G4ElementTable* theElementTable = G4Element::GetElementTable(); 00078 00079 theCrossSectionTable = new G4PhysicsTable(nEl); 00080 00081 G4double AtomicNumber; 00082 G4double AtomicWeight; 00083 G4double Value; 00084 00085 for (G4int j = 0; j < nEl; j++) { 00086 ptrVector = new G4PhysicsLogVector(LowestKineticEnergy, 00087 HighestKineticEnergy, TotBin); 00088 AtomicNumber = (*theElementTable)[j]->GetZ(); 00089 AtomicWeight = (*theElementTable)[j]->GetA(); 00090 00091 for (G4int i = 0; i <= TotBin; ++i) { 00092 lowEdgeEnergy = ptrVector->Energy(i); 00093 Value = ComputeMicroscopicCrossSection(lowEdgeEnergy, 00094 AtomicNumber, AtomicWeight); 00095 ptrVector->PutValue(i,Value); 00096 } 00097 00098 theCrossSectionTable->insertAt(j, ptrVector); 00099 } 00100 00101 // Build (Z,element) look-up table (for use in GetZandACrossSection) 00102 for (G4int j = 0; j < nEl; j++) { 00103 G4int ZZ = G4int((*theElementTable)[j]->GetZ()); 00104 zelMap[ZZ] = j; 00105 } 00106 }
G4double G4KokoulinMuonNuclearXS::ComputeDDMicroscopicCrossSection | ( | G4double | incidentKE, | |
G4double | , | |||
G4double | AtomicWeight, | |||
G4double | epsilon | |||
) |
Definition at line 154 of file G4KokoulinMuonNuclearXS.cc.
References G4ParticleDefinition::GetPDGMass(), G4MuonMinus::MuonMinus(), and G4INCL::Math::pi.
00157 { 00158 // Calculate the double-differential microscopic cross section (in muon 00159 // incident kinetic energy and energy loss) using the cross section formula 00160 // of R.P. Kokoulin (18/01/98) 00161 00162 const G4double alam2 = 0.400*GeV*GeV; 00163 const G4double alam = 0.632456*GeV; 00164 const G4double coeffn = fine_structure_const/pi; 00165 00166 G4double ParticleMass = G4MuonMinus::MuonMinus()->GetPDGMass(); 00167 G4double TotalEnergy = KineticEnergy + ParticleMass; 00168 00169 G4double DCrossSection = 0.; 00170 00171 if ((epsilon >= TotalEnergy - 0.5*proton_mass_c2) || 00172 (epsilon <= CutFixed) ) return DCrossSection; 00173 00174 G4double ep = epsilon/GeV; 00175 G4double a = AtomicWeight/(g/mole); 00176 G4double aeff = 0.22*a+0.78*std::exp(0.89*std::log(a)); //shadowing 00177 G4double sigph = (49.2+11.1*std::log(ep)+151.8/std::sqrt(ep))*microbarn; 00178 00179 G4double v = epsilon/TotalEnergy; 00180 G4double v1 = 1.-v; 00181 G4double v2 = v*v; 00182 G4double mass2 = ParticleMass*ParticleMass; 00183 00184 G4double up = TotalEnergy*TotalEnergy*v1/mass2*(1.+mass2*v2/(alam2*v1)); 00185 G4double down = 1.+epsilon/alam*(1.+alam/(2.*proton_mass_c2)+epsilon/alam); 00186 00187 DCrossSection = coeffn*aeff*sigph/epsilon* 00188 (-v1+(v1+0.5*v2*(1.+2.*mass2/alam2))*std::log(up/down)); 00189 00190 if (DCrossSection < 0.) DCrossSection = 0.; 00191 return DCrossSection; 00192 }
G4double G4KokoulinMuonNuclearXS::GetElementCrossSection | ( | const G4DynamicParticle * | particle, | |
G4int | Z, | |||
const G4Material * | ||||
) | [virtual] |
Reimplemented from G4VCrossSectionDataSet.
Definition at line 195 of file G4KokoulinMuonNuclearXS.cc.
References G4DynamicParticle::GetKineticEnergy().
00197 { 00198 G4int j = zelMap[ZZ]; 00199 return (*theCrossSectionTable)[j]->Value(aPart->GetKineticEnergy()); 00200 }
G4bool G4KokoulinMuonNuclearXS::IsElementApplicable | ( | const G4DynamicParticle * | particle, | |
G4int | Z, | |||
const G4Material * | ||||
) | [virtual] |
Reimplemented from G4VCrossSectionDataSet.
Definition at line 64 of file G4KokoulinMuonNuclearXS.cc.