#include <G4AnnihiToMuPair.hh>
Inheritance diagram for G4AnnihiToMuPair:
Public Member Functions | |
G4AnnihiToMuPair (const G4String &processName="AnnihiToMuPair", G4ProcessType type=fElectromagnetic) | |
~G4AnnihiToMuPair () | |
G4bool | IsApplicable (const G4ParticleDefinition &) |
void | BuildPhysicsTable (const G4ParticleDefinition &) |
void | PrintInfoDefinition () |
void | SetCrossSecFactor (G4double fac) |
G4double | GetCrossSecFactor () |
G4double | CrossSectionPerVolume (G4double PositronEnergy, const G4Material *) |
G4double | ComputeCrossSectionPerAtom (G4double PositronEnergy, G4double AtomicZ) |
G4double | GetMeanFreePath (const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *) |
G4VParticleChange * | PostStepDoIt (const G4Track &aTrack, const G4Step &aStep) |
Definition at line 57 of file G4AnnihiToMuPair.hh.
G4AnnihiToMuPair::G4AnnihiToMuPair | ( | const G4String & | processName = "AnnihiToMuPair" , |
|
G4ProcessType | type = fElectromagnetic | |||
) |
Definition at line 58 of file G4AnnihiToMuPair.cc.
References G4ParticleDefinition::GetPDGMass(), G4MuonPlus::MuonPlus(), and G4VProcess::SetProcessSubType().
00059 :G4VDiscreteProcess (processName, type) 00060 { 00061 //e+ Energy threshold 00062 const G4double Mu_massc2 = G4MuonPlus::MuonPlus()->GetPDGMass(); 00063 LowestEnergyLimit = 2*Mu_massc2*Mu_massc2/electron_mass_c2 - electron_mass_c2; 00064 00065 //modele ok up to 1000 TeV due to neglected Z-interference 00066 HighestEnergyLimit = 1000*TeV; 00067 00068 CurrentSigma = 0.0; 00069 CrossSecFactor = 1.; 00070 SetProcessSubType(6); 00071 00072 }
G4AnnihiToMuPair::~G4AnnihiToMuPair | ( | ) |
void G4AnnihiToMuPair::BuildPhysicsTable | ( | const G4ParticleDefinition & | ) | [virtual] |
Reimplemented from G4VProcess.
Definition at line 88 of file G4AnnihiToMuPair.cc.
References PrintInfoDefinition().
00091 { 00092 CurrentSigma = 0.0; 00093 PrintInfoDefinition(); 00094 }
Definition at line 108 of file G4AnnihiToMuPair.cc.
References G4ParticleDefinition::GetPDGMass(), G4MuonPlus::MuonPlus(), and G4INCL::Math::pi.
Referenced by CrossSectionPerVolume().
00111 { 00112 static const G4double Mmuon = G4MuonPlus::MuonPlus()->GetPDGMass(); 00113 static const G4double Rmuon = elm_coupling/Mmuon; //classical particle radius 00114 static const G4double Sig0 = pi*Rmuon*Rmuon/3.; //constant in crossSection 00115 00116 G4double CrossSection = 0.; 00117 if (Epos < LowestEnergyLimit) return CrossSection; 00118 00119 G4double xi = LowestEnergyLimit/Epos; 00120 G4double SigmaEl = Sig0*xi*(1.+xi/2.)*sqrt(1.-xi); // per electron 00121 CrossSection = SigmaEl*Z; // number of electrons per atom 00122 return CrossSection; 00123 }
G4double G4AnnihiToMuPair::CrossSectionPerVolume | ( | G4double | PositronEnergy, | |
const G4Material * | ||||
) |
Definition at line 127 of file G4AnnihiToMuPair.cc.
References ComputeCrossSectionPerAtom(), G4Material::GetElementVector(), G4Material::GetNumberOfElements(), and G4Material::GetVecNbOfAtomsPerVolume().
Referenced by GetMeanFreePath(), and PostStepDoIt().
00129 { 00130 const G4ElementVector* theElementVector = aMaterial->GetElementVector(); 00131 const G4double* NbOfAtomsPerVolume = aMaterial->GetVecNbOfAtomsPerVolume(); 00132 00133 G4double SIGMA = 0.0; 00134 00135 for ( size_t i=0 ; i < aMaterial->GetNumberOfElements() ; ++i ) 00136 { 00137 G4double AtomicZ = (*theElementVector)[i]->GetZ(); 00138 SIGMA += NbOfAtomsPerVolume[i] * 00139 ComputeCrossSectionPerAtom(PositronEnergy,AtomicZ); 00140 } 00141 return SIGMA; 00142 }
G4double G4AnnihiToMuPair::GetCrossSecFactor | ( | ) | [inline] |
G4double G4AnnihiToMuPair::GetMeanFreePath | ( | const G4Track & | aTrack, | |
G4double | previousStepSize, | |||
G4ForceCondition * | ||||
) | [virtual] |
Implements G4VDiscreteProcess.
Definition at line 146 of file G4AnnihiToMuPair.cc.
References CrossSectionPerVolume(), DBL_MAX, DBL_MIN, G4Track::GetDynamicParticle(), G4DynamicParticle::GetKineticEnergy(), and G4Track::GetMaterial().
00151 { 00152 const G4DynamicParticle* aDynamicPositron = aTrack.GetDynamicParticle(); 00153 G4double PositronEnergy = aDynamicPositron->GetKineticEnergy() 00154 +electron_mass_c2; 00155 G4Material* aMaterial = aTrack.GetMaterial(); 00156 CurrentSigma = CrossSectionPerVolume(PositronEnergy, aMaterial); 00157 00158 // increase the CrossSection by CrossSecFactor (default 1) 00159 G4double mfp = DBL_MAX; 00160 if(CurrentSigma > DBL_MIN) mfp = 1.0/(CurrentSigma*CrossSecFactor); 00161 00162 return mfp; 00163 }
G4bool G4AnnihiToMuPair::IsApplicable | ( | const G4ParticleDefinition & | ) | [virtual] |
Reimplemented from G4VProcess.
Definition at line 81 of file G4AnnihiToMuPair.cc.
References G4Positron::Positron().
00082 { 00083 return ( &particle == G4Positron::Positron() ); 00084 }
G4VParticleChange * G4AnnihiToMuPair::PostStepDoIt | ( | const G4Track & | aTrack, | |
const G4Step & | aStep | |||
) | [virtual] |
Reimplemented from G4VDiscreteProcess.
Definition at line 167 of file G4AnnihiToMuPair.cc.
References G4ParticleChange::AddSecondary(), G4VProcess::aParticleChange, CrossSectionPerVolume(), fStopAndKill, G4UniformRand, G4Track::GetDynamicParticle(), G4DynamicParticle::GetKineticEnergy(), G4Track::GetMaterial(), G4DynamicParticle::GetMomentumDirection(), G4ParticleDefinition::GetPDGMass(), G4ParticleChange::Initialize(), G4MuonMinus::MuonMinus(), G4MuonPlus::MuonPlus(), G4INCL::Math::pi, G4VDiscreteProcess::PostStepDoIt(), G4ParticleChange::ProposeEnergy(), G4VParticleChange::ProposeTrackStatus(), and G4VParticleChange::SetNumberOfSecondaries().
00172 { 00173 00174 aParticleChange.Initialize(aTrack); 00175 static const G4double Mele=electron_mass_c2; 00176 static const G4double Mmuon=G4MuonPlus::MuonPlus()->GetPDGMass(); 00177 00178 // current Positron energy and direction, return if energy too low 00179 const G4DynamicParticle *aDynamicPositron = aTrack.GetDynamicParticle(); 00180 G4double Epos = aDynamicPositron->GetKineticEnergy() + Mele; 00181 00182 // test of cross section 00183 if(CurrentSigma*G4UniformRand() > 00184 CrossSectionPerVolume(Epos, aTrack.GetMaterial())) 00185 { 00186 return G4VDiscreteProcess::PostStepDoIt(aTrack,aStep); 00187 } 00188 00189 if (Epos < LowestEnergyLimit) { 00190 return G4VDiscreteProcess::PostStepDoIt(aTrack,aStep); 00191 } 00192 00193 G4ParticleMomentum PositronDirection = 00194 aDynamicPositron->GetMomentumDirection(); 00195 G4double xi = LowestEnergyLimit/Epos; // xi is always less than 1, 00196 // goes to 0 at high Epos 00197 00198 // generate cost 00199 // 00200 G4double cost; 00201 do cost = 2.*G4UniformRand()-1.; 00202 while (2.*G4UniformRand() > 1.+xi+cost*cost*(1.-xi) ); 00203 //1+cost**2 at high Epos 00204 G4double sint = sqrt(1.-cost*cost); 00205 00206 // generate phi 00207 // 00208 G4double phi=2.*pi*G4UniformRand(); 00209 00210 G4double Ecm = sqrt(0.5*Mele*(Epos+Mele)); 00211 G4double Pcm = sqrt(Ecm*Ecm-Mmuon*Mmuon); 00212 G4double beta = sqrt((Epos-Mele)/(Epos+Mele)); 00213 G4double gamma = Ecm/Mele; // =sqrt((Epos+Mele)/(2.*Mele)); 00214 G4double Pt = Pcm*sint; 00215 00216 // energy and momentum of the muons in the Lab 00217 // 00218 G4double EmuPlus = gamma*( Ecm+cost*beta*Pcm); 00219 G4double EmuMinus = gamma*( Ecm-cost*beta*Pcm); 00220 G4double PmuPlusZ = gamma*(beta*Ecm+cost* Pcm); 00221 G4double PmuMinusZ = gamma*(beta*Ecm-cost* Pcm); 00222 G4double PmuPlusX = Pt*cos(phi); 00223 G4double PmuPlusY = Pt*sin(phi); 00224 G4double PmuMinusX =-Pt*cos(phi); 00225 G4double PmuMinusY =-Pt*sin(phi); 00226 // absolute momenta 00227 G4double PmuPlus = sqrt(Pt*Pt+PmuPlusZ *PmuPlusZ ); 00228 G4double PmuMinus = sqrt(Pt*Pt+PmuMinusZ*PmuMinusZ); 00229 00230 // mu+ mu- directions for Positron in z-direction 00231 // 00232 G4ThreeVector 00233 MuPlusDirection ( PmuPlusX/PmuPlus, PmuPlusY/PmuPlus, PmuPlusZ/PmuPlus ); 00234 G4ThreeVector 00235 MuMinusDirection(PmuMinusX/PmuMinus,PmuMinusY/PmuMinus,PmuMinusZ/PmuMinus); 00236 00237 // rotate to actual Positron direction 00238 // 00239 MuPlusDirection.rotateUz(PositronDirection); 00240 MuMinusDirection.rotateUz(PositronDirection); 00241 00242 aParticleChange.SetNumberOfSecondaries(2); 00243 // create G4DynamicParticle object for the particle1 00244 G4DynamicParticle* aParticle1= new G4DynamicParticle( 00245 G4MuonPlus::MuonPlus(),MuPlusDirection,EmuPlus-Mmuon); 00246 aParticleChange.AddSecondary(aParticle1); 00247 // create G4DynamicParticle object for the particle2 00248 G4DynamicParticle* aParticle2= new G4DynamicParticle( 00249 G4MuonMinus::MuonMinus(),MuMinusDirection,EmuMinus-Mmuon); 00250 aParticleChange.AddSecondary(aParticle2); 00251 00252 // Kill the incident positron 00253 // 00254 aParticleChange.ProposeEnergy(0.); 00255 aParticleChange.ProposeTrackStatus(fStopAndKill); 00256 00257 return &aParticleChange; 00258 }
void G4AnnihiToMuPair::PrintInfoDefinition | ( | ) |
Definition at line 262 of file G4AnnihiToMuPair.cc.
References G4cout, G4endl, G4VProcess::GetProcessName(), and G4VProcess::GetProcessSubType().
Referenced by BuildPhysicsTable().
00263 { 00264 G4String comments ="e+e->mu+mu- annihilation, atomic e- at rest, SubType=."; 00265 G4cout << G4endl << GetProcessName() << ": " << comments 00266 << GetProcessSubType() << G4endl; 00267 G4cout << " threshold at " << LowestEnergyLimit/GeV << " GeV" 00268 << " good description up to " 00269 << HighestEnergyLimit/TeV << " TeV for all Z." << G4endl; 00270 }
void G4AnnihiToMuPair::SetCrossSecFactor | ( | G4double | fac | ) |