#include <G4IonCoulombScatteringModel.hh>
Inheritance diagram for G4IonCoulombScatteringModel:
Definition at line 71 of file G4IonCoulombScatteringModel.hh.
G4IonCoulombScatteringModel::G4IonCoulombScatteringModel | ( | const G4String & | nam = "IonCoulombScattering" |
) |
Definition at line 73 of file G4IonCoulombScatteringModel.cc.
References currentCouple, currentElement, currentMaterial, currentMaterialIndex, fNistManager, G4ParticleTable::GetParticleTable(), heavycorr, G4NistManager::Instance(), ioncross, lowEnergyLimit, mass, particle, pCuts, G4Proton::Proton(), recoilThreshold, theParticleTable, and theProton.
00074 : G4VEmModel(nam), 00075 00076 cosThetaMin(1.0), 00077 isInitialised(false) 00078 { 00079 fNistManager = G4NistManager::Instance(); 00080 theParticleTable = G4ParticleTable::GetParticleTable(); 00081 theProton = G4Proton::Proton(); 00082 00083 pCuts=0; 00084 currentMaterial = 0; 00085 currentElement = 0; 00086 currentCouple = 0; 00087 00088 lowEnergyLimit = 100*eV; 00089 recoilThreshold = 0.*eV; 00090 heavycorr =0; 00091 particle = 0; 00092 mass=0; 00093 currentMaterialIndex = -1; 00094 00095 ioncross = new G4IonCoulombCrossSection(); 00096 00097 }
G4IonCoulombScatteringModel::~G4IonCoulombScatteringModel | ( | ) | [virtual] |
Definition at line 102 of file G4IonCoulombScatteringModel.cc.
References ioncross.
00103 { delete ioncross;}
G4double G4IonCoulombScatteringModel::ComputeCrossSectionPerAtom | ( | const G4ParticleDefinition * | , | |
G4double | kinEnergy, | |||
G4double | Z, | |||
G4double | A, | |||
G4double | cut, | |||
G4double | emax | |||
) | [virtual] |
Reimplemented from G4VEmModel.
Definition at line 127 of file G4IonCoulombScatteringModel.cc.
References G4VEmModel::CurrentCouple(), DefineMaterial(), heavycorr, ioncross, lowEnergyLimit, G4IonCoulombCrossSection::NuclearCrossSection(), G4IonCoulombCrossSection::SetupKinematic(), SetupParticle(), and G4IonCoulombCrossSection::SetupTarget().
Referenced by SampleSecondaries().
00134 { 00135 00136 SetupParticle(p); 00137 00138 G4double cross =0.0; 00139 if(kinEnergy < lowEnergyLimit) return cross; 00140 00141 DefineMaterial(CurrentCouple()); 00142 00143 G4int iz = G4int(Z); 00144 00145 //from lab to pCM & mu_rel of effective particle 00146 ioncross->SetupKinematic(kinEnergy, cutEnergy,iz); 00147 00148 00149 ioncross->SetupTarget(Z, kinEnergy, heavycorr); 00150 00151 cross = ioncross->NuclearCrossSection(); 00152 00153 //cout<< "..........cross "<<G4BestUnit(cross,"Surface") <<endl; 00154 return cross; 00155 }
void G4IonCoulombScatteringModel::DefineMaterial | ( | const G4MaterialCutsCouple * | ) | [inline, protected] |
Definition at line 157 of file G4IonCoulombScatteringModel.hh.
References currentCouple, currentMaterial, currentMaterialIndex, G4MaterialCutsCouple::GetIndex(), and G4MaterialCutsCouple::GetMaterial().
Referenced by ComputeCrossSectionPerAtom(), and SampleSecondaries().
00158 { 00159 if(cup != currentCouple) { 00160 currentCouple = cup; 00161 currentMaterial = cup->GetMaterial(); 00162 currentMaterialIndex = currentCouple->GetIndex(); 00163 00164 } 00165 }
G4int G4IonCoulombScatteringModel::GetHeavyIonCorr | ( | ) | [inline] |
Definition at line 100 of file G4IonCoulombScatteringModel.hh.
References heavycorr.
00100 {return heavycorr; };
void G4IonCoulombScatteringModel::Initialise | ( | const G4ParticleDefinition * | , | |
const G4DataVector & | ||||
) | [virtual] |
Implements G4VEmModel.
Definition at line 107 of file G4IonCoulombScatteringModel.cc.
References cosThetaMin, currentCouple, currentMaterialIndex, fParticleChange, G4ProductionCutsTable::GetEnergyCutsVector(), G4VEmModel::GetParticleChangeForGamma(), G4ProductionCutsTable::GetProductionCutsTable(), G4IonCoulombCrossSection::Initialise(), ioncross, pCuts, G4VEmModel::PolarAngleLimit(), and SetupParticle().
00109 { 00110 SetupParticle(p); 00111 currentCouple = 0; 00112 currentMaterialIndex = -1; 00113 cosThetaMin = cos(PolarAngleLimit()); 00114 ioncross->Initialise(p,cosThetaMin); 00115 00116 pCuts = G4ProductionCutsTable::GetProductionCutsTable()->GetEnergyCutsVector(3); 00117 00118 00119 if(!isInitialised) { 00120 isInitialised = true; 00121 fParticleChange = GetParticleChangeForGamma(); 00122 } 00123 }
void G4IonCoulombScatteringModel::SampleSecondaries | ( | std::vector< G4DynamicParticle * > * | , | |
const G4MaterialCutsCouple * | , | |||
const G4DynamicParticle * | , | |||
G4double | tmin, | |||
G4double | maxEnergy | |||
) | [virtual] |
Implements G4VEmModel.
Definition at line 159 of file G4IonCoulombScatteringModel.cc.
References ComputeCrossSectionPerAtom(), currentElement, currentMaterialIndex, DefineMaterial(), fParticleChange, G4UniformRand, G4InuclParticleNames::gam, G4DynamicParticle::GetDefinition(), G4ParticleTable::GetIon(), G4DynamicParticle::GetKineticEnergy(), G4IonCoulombCrossSection::GetMomentum2(), G4DynamicParticle::GetMomentumDirection(), G4NucleiProperties::GetNuclearMass(), ioncross, lowEnergyLimit, mass, particle, pCuts, G4VParticleChange::ProposeLocalEnergyDeposit(), G4ParticleChangeForGamma::ProposeMomentumDirection(), G4VParticleChange::ProposeNonIonizingEnergyDeposit(), recoilThreshold, G4IonCoulombCrossSection::SampleCosineTheta(), G4VEmModel::SelectIsotopeNumber(), G4VEmModel::SelectRandomAtom(), G4ParticleChangeForGamma::SetProposedKineticEnergy(), SetupParticle(), and theParticleTable.
00165 { 00166 G4double kinEnergy = dp->GetKineticEnergy(); 00167 00168 if(kinEnergy < lowEnergyLimit) return; 00169 00170 DefineMaterial(couple); 00171 00172 SetupParticle(dp->GetDefinition()); 00173 00174 // Choose nucleus 00175 currentElement = SelectRandomAtom(couple,particle, 00176 kinEnergy,cutEnergy,kinEnergy); 00177 00178 G4double Z = currentElement->GetZ(); 00179 G4int iz = G4int(Z); 00180 G4int ia = SelectIsotopeNumber(currentElement); 00181 G4double mass2 = G4NucleiProperties::GetNuclearMass(ia, iz); 00182 00183 00184 00185 G4double cross= ComputeCrossSectionPerAtom(particle,kinEnergy, Z, 00186 kinEnergy, cutEnergy, kinEnergy) ; 00187 if(cross == 0.0) { return; } 00188 00189 //scattering angle, z1 == (1-cost) 00190 G4double z1 = ioncross->SampleCosineTheta(); 00191 if(z1 > 2.0) { z1 = 2.0; } 00192 else if(z1 < 0.0) { z1 = 0.0; } 00193 00194 G4double cost = 1.0 - z1; 00195 G4double sint = sqrt(z1*(1.0 + cost)); 00196 G4double phi = twopi * G4UniformRand(); 00197 00198 00199 // kinematics in the Lab system 00200 G4double etot = kinEnergy + mass; 00201 G4double mom2= kinEnergy*(kinEnergy+2.0*mass); 00202 G4double ptot = sqrt(mom2); 00203 00204 //CM particle 1 00205 G4double bet = ptot/(etot + mass2); 00206 G4double gam = 1.0/sqrt((1.0 - bet)*(1.0 + bet)); 00207 00208 //CM 00209 G4double momCM2= ioncross->GetMomentum2(); 00210 G4double momCM =std::sqrt(momCM2); 00211 //energy & momentum after scattering of incident particle 00212 G4double pxCM = momCM*sint*cos(phi); 00213 G4double pyCM = momCM*sint*sin(phi); 00214 G4double pzCM = momCM*cost; 00215 G4double eCM = sqrt(momCM2 + mass*mass); 00216 00217 //CM--->Lab 00218 G4ThreeVector v1(pxCM , pyCM, gam*(pzCM + bet*eCM)); 00219 G4ThreeVector dir = dp->GetMomentumDirection(); 00220 00221 G4ThreeVector newDirection = v1.unit(); 00222 newDirection.rotateUz(dir); 00223 00224 fParticleChange->ProposeMomentumDirection(newDirection); 00225 00226 // V.Ivanchenko fix of final energies after scattering 00227 // recoil....................................... 00228 //G4double trec =(1.0 - cost)* mass2*(etot*etot - mass*mass )/ 00229 // (mass*mass + mass2*mass2+ 2.*mass2*etot); 00230 //G4double finalT = kinEnergy - trec; 00231 00232 // new computation 00233 G4double finalT = gam*(eCM + bet*pzCM) - mass; 00234 G4double trec = kinEnergy - finalT; 00235 00236 if(finalT <= lowEnergyLimit) { 00237 trec = kinEnergy; 00238 finalT = 0.0; 00239 } else if(trec < 0.0) { 00240 trec = 0.0; 00241 finalT = kinEnergy; 00242 } 00243 00244 fParticleChange->SetProposedKineticEnergy(finalT); 00245 00246 G4double tcut = recoilThreshold; 00247 if(pCuts) { tcut= std::max(tcut,(*pCuts)[currentMaterialIndex]); 00248 00249 //G4cout<<" tcut eV "<<tcut/eV<<endl; 00250 } 00251 00252 if(trec > tcut) { 00253 G4ParticleDefinition* ion = theParticleTable->GetIon(iz, ia, 0.0); 00254 G4double plab = sqrt(finalT*(finalT + 2.0*mass)); 00255 G4ThreeVector p2 = (ptot*dir - plab*newDirection).unit(); 00256 G4DynamicParticle* newdp = new G4DynamicParticle(ion, p2, trec); 00257 fvect->push_back(newdp); 00258 } else if(trec > 0.0) { 00259 fParticleChange->ProposeLocalEnergyDeposit(trec); 00260 fParticleChange->ProposeNonIonizingEnergyDeposit(trec); 00261 } 00262 00263 00264 }
void G4IonCoulombScatteringModel::SetHeavyIonCorr | ( | G4int | b | ) | [inline] |
Definition at line 99 of file G4IonCoulombScatteringModel.hh.
References heavycorr.
00099 {heavycorr=b; };
void G4IonCoulombScatteringModel::SetRecoilThreshold | ( | G4double | eth | ) | [inline] |
Definition at line 180 of file G4IonCoulombScatteringModel.hh.
References recoilThreshold.
00181 { 00182 recoilThreshold = eth; 00183 }
void G4IonCoulombScatteringModel::SetupParticle | ( | const G4ParticleDefinition * | ) | [inline, protected] |
Definition at line 170 of file G4IonCoulombScatteringModel.hh.
References G4ParticleDefinition::GetPDGMass(), ioncross, mass, particle, and G4IonCoulombCrossSection::SetupParticle().
Referenced by ComputeCrossSectionPerAtom(), Initialise(), and SampleSecondaries().
00171 { 00172 if(p != particle) { 00173 particle = p; 00174 mass = particle->GetPDGMass(); 00175 ioncross->SetupParticle(p); 00176 } 00177 }
G4double G4IonCoulombScatteringModel::cosThetaMin [protected] |
const G4MaterialCutsCouple* G4IonCoulombScatteringModel::currentCouple [protected] |
Definition at line 131 of file G4IonCoulombScatteringModel.hh.
Referenced by DefineMaterial(), G4IonCoulombScatteringModel(), and Initialise().
const G4Element* G4IonCoulombScatteringModel::currentElement [protected] |
Definition at line 133 of file G4IonCoulombScatteringModel.hh.
Referenced by G4IonCoulombScatteringModel(), and SampleSecondaries().
const G4Material* G4IonCoulombScatteringModel::currentMaterial [protected] |
Definition at line 132 of file G4IonCoulombScatteringModel.hh.
Referenced by DefineMaterial(), and G4IonCoulombScatteringModel().
Definition at line 134 of file G4IonCoulombScatteringModel.hh.
Referenced by DefineMaterial(), G4IonCoulombScatteringModel(), Initialise(), and SampleSecondaries().
G4NistManager* G4IonCoulombScatteringModel::fNistManager [protected] |
Definition at line 127 of file G4IonCoulombScatteringModel.hh.
Referenced by G4IonCoulombScatteringModel().
Definition at line 126 of file G4IonCoulombScatteringModel.hh.
Referenced by Initialise(), and SampleSecondaries().
G4int G4IonCoulombScatteringModel::heavycorr [protected] |
Definition at line 137 of file G4IonCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), G4IonCoulombScatteringModel(), GetHeavyIonCorr(), and SetHeavyIonCorr().
Definition at line 128 of file G4IonCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), G4IonCoulombScatteringModel(), Initialise(), SampleSecondaries(), SetupParticle(), and ~G4IonCoulombScatteringModel().
G4double G4IonCoulombScatteringModel::lowEnergyLimit [protected] |
Definition at line 147 of file G4IonCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), G4IonCoulombScatteringModel(), and SampleSecondaries().
G4double G4IonCoulombScatteringModel::mass [protected] |
Definition at line 146 of file G4IonCoulombScatteringModel.hh.
Referenced by G4IonCoulombScatteringModel(), SampleSecondaries(), and SetupParticle().
const G4ParticleDefinition* G4IonCoulombScatteringModel::particle [protected] |
Definition at line 144 of file G4IonCoulombScatteringModel.hh.
Referenced by G4IonCoulombScatteringModel(), SampleSecondaries(), and SetupParticle().
const std::vector<G4double>* G4IonCoulombScatteringModel::pCuts [protected] |
Definition at line 130 of file G4IonCoulombScatteringModel.hh.
Referenced by G4IonCoulombScatteringModel(), Initialise(), and SampleSecondaries().
G4double G4IonCoulombScatteringModel::recoilThreshold [protected] |
Definition at line 140 of file G4IonCoulombScatteringModel.hh.
Referenced by G4IonCoulombScatteringModel(), SampleSecondaries(), and SetRecoilThreshold().
Definition at line 125 of file G4IonCoulombScatteringModel.hh.
Referenced by G4IonCoulombScatteringModel(), and SampleSecondaries().
const G4ParticleDefinition* G4IonCoulombScatteringModel::theProton [protected] |
Definition at line 145 of file G4IonCoulombScatteringModel.hh.
Referenced by G4IonCoulombScatteringModel().