#include <G4eCoulombScatteringModel.hh>
Inheritance diagram for G4eCoulombScatteringModel:
Definition at line 78 of file G4eCoulombScatteringModel.hh.
G4eCoulombScatteringModel::G4eCoulombScatteringModel | ( | const G4String & | nam = "eCoulombScattering" |
) |
Definition at line 79 of file G4eCoulombScatteringModel.cc.
References cosTetMaxNuc, cosTetMinNuc, currentCouple, currentMaterial, currentMaterialIndex, elecRatio, fNistManager, fParticleChange, G4ParticleTable::GetParticleTable(), G4NistManager::Instance(), lowEnergyThreshold, mass, particle, pCuts, G4Proton::Proton(), recoilThreshold, theParticleTable, theProton, and wokvi.
00080 : G4VEmModel(nam), 00081 cosThetaMin(1.0), 00082 cosThetaMax(-1.0), 00083 isInitialised(false) 00084 { 00085 fParticleChange = 0; 00086 fNistManager = G4NistManager::Instance(); 00087 theParticleTable = G4ParticleTable::GetParticleTable(); 00088 theProton = G4Proton::Proton(); 00089 currentMaterial = 0; 00090 00091 pCuts = 0; 00092 00093 lowEnergyThreshold = 1*keV; // particle will be killed for lower energy 00094 recoilThreshold = 0.*keV; // by default does not work 00095 00096 particle = 0; 00097 currentCouple = 0; 00098 wokvi = new G4WentzelOKandVIxSection(); 00099 00100 currentMaterialIndex = 0; 00101 00102 cosTetMinNuc = 1.0; 00103 cosTetMaxNuc = -1.0; 00104 elecRatio = 0.0; 00105 mass = proton_mass_c2; 00106 }
G4eCoulombScatteringModel::~G4eCoulombScatteringModel | ( | ) | [virtual] |
Definition at line 110 of file G4eCoulombScatteringModel.cc.
References wokvi.
00111 { 00112 delete wokvi; 00113 }
G4double G4eCoulombScatteringModel::ComputeCrossSectionPerAtom | ( | const G4ParticleDefinition * | , | |
G4double | kinEnergy, | |||
G4double | Z, | |||
G4double | A, | |||
G4double | cut, | |||
G4double | emax | |||
) | [virtual] |
Reimplemented from G4VEmModel.
Definition at line 148 of file G4eCoulombScatteringModel.cc.
References G4WentzelOKandVIxSection::ComputeElectronCrossSection(), G4WentzelOKandVIxSection::ComputeNuclearCrossSection(), cosTetMaxNuc, cosTetMinNuc, cosThetaMax, G4VEmModel::CurrentCouple(), currentMaterial, DefineMaterial(), elecRatio, particle, G4WentzelOKandVIxSection::SetupKinematic(), SetupParticle(), G4WentzelOKandVIxSection::SetupTarget(), theProton, and wokvi.
Referenced by SampleSecondaries().
00153 { 00154 //G4cout << "### G4eCoulombScatteringModel::ComputeCrossSectionPerAtom for " 00155 // << p->GetParticleName()<<" Z= "<<Z<<" e(MeV)= "<< kinEnergy/MeV << G4endl; 00156 G4double cross = 0.0; 00157 if(p != particle) { SetupParticle(p); } 00158 00159 // cross section is set to zero to avoid problems in sample secondary 00160 if(kinEnergy <= 0.0) { return cross; } 00161 DefineMaterial(CurrentCouple()); 00162 cosTetMinNuc = wokvi->SetupKinematic(kinEnergy, currentMaterial); 00163 if(cosThetaMax < cosTetMinNuc) { 00164 G4int iz = G4int(Z); 00165 cosTetMinNuc = wokvi->SetupTarget(iz, cutEnergy); 00166 cosTetMaxNuc = cosThetaMax; 00167 if(iz == 1 && cosTetMaxNuc < 0.0 && particle == theProton) { 00168 cosTetMaxNuc = 0.0; 00169 } 00170 cross = wokvi->ComputeNuclearCrossSection(cosTetMinNuc, cosTetMaxNuc); 00171 elecRatio = wokvi->ComputeElectronCrossSection(cosTetMinNuc, cosThetaMax); 00172 cross += elecRatio; 00173 if(cross > 0.0) { elecRatio /= cross; } 00174 } 00175 /* 00176 if(p->GetParticleName() == "e-") 00177 G4cout << "e(MeV)= " << kinEnergy/MeV << " cross(b)= " << cross/barn 00178 << " 1-cosTetMinNuc= " << 1-cosTetMinNuc 00179 << " 1-cosTetMaxNuc= " << 1-cosTetMaxNuc 00180 << G4endl; 00181 */ 00182 return cross; 00183 }
void G4eCoulombScatteringModel::DefineMaterial | ( | const G4MaterialCutsCouple * | ) | [inline, protected] |
Definition at line 156 of file G4eCoulombScatteringModel.hh.
References currentCouple, currentMaterial, currentMaterialIndex, G4MaterialCutsCouple::GetIndex(), and G4MaterialCutsCouple::GetMaterial().
Referenced by ComputeCrossSectionPerAtom(), and SampleSecondaries().
00157 { 00158 if(cup != currentCouple) { 00159 currentCouple = cup; 00160 currentMaterial = cup->GetMaterial(); 00161 currentMaterialIndex = currentCouple->GetIndex(); 00162 } 00163 }
void G4eCoulombScatteringModel::Initialise | ( | const G4ParticleDefinition * | , | |
const G4DataVector & | ||||
) | [virtual] |
Implements G4VEmModel.
Definition at line 117 of file G4eCoulombScatteringModel.cc.
References cosThetaMin, currentCouple, fParticleChange, G4ProductionCutsTable::GetEnergyCutsVector(), G4VEmModel::GetParticleChangeForGamma(), G4ProductionCutsTable::GetProductionCutsTable(), G4WentzelOKandVIxSection::Initialise(), G4VEmModel::InitialiseElementSelectors(), mass, pCuts, G4VEmModel::PolarAngleLimit(), SetupParticle(), and wokvi.
00119 { 00120 SetupParticle(p); 00121 currentCouple = 0; 00122 cosThetaMin = cos(PolarAngleLimit()); 00123 wokvi->Initialise(p, cosThetaMin); 00124 /* 00125 G4cout << "G4eCoulombScatteringModel: " << particle->GetParticleName() 00126 << " 1-cos(ThetaLimit)= " << 1 - cosThetaMin 00127 << " cos(thetaMax)= " << cosThetaMax 00128 << G4endl; 00129 */ 00130 pCuts = G4ProductionCutsTable::GetProductionCutsTable()->GetEnergyCutsVector(3); 00131 /* 00132 G4cout << "!!! G4eCoulombScatteringModel::Initialise for " 00133 << p->GetParticleName() << " cos(TetMin)= " << cosThetaMin 00134 << " cos(TetMax)= " << cosThetaMax <<G4endl; 00135 G4cout << "cut0= " << cuts[0] << " cut1= " << cuts[1] << G4endl; 00136 */ 00137 if(!isInitialised) { 00138 isInitialised = true; 00139 fParticleChange = GetParticleChangeForGamma(); 00140 } 00141 if(mass < GeV) { 00142 InitialiseElementSelectors(p,cuts); 00143 } 00144 }
void G4eCoulombScatteringModel::SampleSecondaries | ( | std::vector< G4DynamicParticle * > * | , | |
const G4MaterialCutsCouple * | , | |||
const G4DynamicParticle * | , | |||
G4double | tmin, | |||
G4double | maxEnergy | |||
) | [virtual] |
Implements G4VEmModel.
Definition at line 187 of file G4eCoulombScatteringModel.cc.
References ComputeCrossSectionPerAtom(), cosTetMinNuc, cosThetaMax, currentMaterialIndex, DefineMaterial(), elecRatio, fParticleChange, G4DynamicParticle::GetDefinition(), G4ParticleTable::GetIon(), G4DynamicParticle::GetKineticEnergy(), G4DynamicParticle::GetMomentumDirection(), G4WentzelOKandVIxSection::GetMomentumSquare(), G4NucleiProperties::GetNuclearMass(), G4Element::GetZ(), lowEnergyThreshold, mass, particle, pCuts, G4VParticleChange::ProposeLocalEnergyDeposit(), G4ParticleChangeForGamma::ProposeMomentumDirection(), G4VParticleChange::ProposeNonIonizingEnergyDeposit(), recoilThreshold, G4WentzelOKandVIxSection::SampleSingleScattering(), G4VEmModel::SelectIsotopeNumber(), G4VEmModel::SelectRandomAtom(), G4ParticleChangeForGamma::SetProposedKineticEnergy(), G4WentzelOKandVIxSection::SetTargetMass(), SetupParticle(), theParticleTable, and wokvi.
00193 { 00194 G4double kinEnergy = dp->GetKineticEnergy(); 00195 00196 // absorb particle below low-energy limit to avoid situation 00197 // when a particle has no energy loss 00198 if(kinEnergy < lowEnergyThreshold) { 00199 fParticleChange->SetProposedKineticEnergy(0.0); 00200 fParticleChange->ProposeLocalEnergyDeposit(kinEnergy); 00201 fParticleChange->ProposeNonIonizingEnergyDeposit(kinEnergy); 00202 return; 00203 } 00204 SetupParticle(dp->GetDefinition()); 00205 DefineMaterial(couple); 00206 /* 00207 G4cout << "G4eCoulombScatteringModel::SampleSecondaries e(MeV)= " 00208 << kinEnergy << " " << particle->GetParticleName() 00209 << " cut= " << cutEnergy<< G4endl; 00210 */ 00211 // Choose nucleus 00212 const G4Element* currentElement = 00213 SelectRandomAtom(couple,particle,kinEnergy,cutEnergy,kinEnergy); 00214 00215 G4double Z = currentElement->GetZ(); 00216 00217 if(ComputeCrossSectionPerAtom(particle,kinEnergy, Z, 00218 kinEnergy, cutEnergy, kinEnergy) == 0.0) 00219 { return; } 00220 00221 G4int iz = G4int(Z); 00222 G4int ia = SelectIsotopeNumber(currentElement); 00223 G4double targetMass = G4NucleiProperties::GetNuclearMass(ia, iz); 00224 wokvi->SetTargetMass(targetMass); 00225 00226 G4ThreeVector newDirection = 00227 wokvi->SampleSingleScattering(cosTetMinNuc, cosThetaMax, elecRatio); 00228 G4double cost = newDirection.z(); 00229 00230 G4ThreeVector direction = dp->GetMomentumDirection(); 00231 newDirection.rotateUz(direction); 00232 00233 fParticleChange->ProposeMomentumDirection(newDirection); 00234 00235 // recoil sampling assuming a small recoil 00236 // and first order correction to primary 4-momentum 00237 G4double mom2 = wokvi->GetMomentumSquare(); 00238 G4double trec = mom2*(1.0 - cost)/(targetMass + (mass + kinEnergy)*(1.0 - cost)); 00239 G4double finalT = kinEnergy - trec; 00240 //G4cout<<"G4eCoulombScatteringModel: finalT= "<<finalT<<" Trec= "<<trec<<G4endl; 00241 if(finalT <= lowEnergyThreshold) { 00242 trec = kinEnergy; 00243 finalT = 0.0; 00244 } 00245 00246 fParticleChange->SetProposedKineticEnergy(finalT); 00247 G4double tcut = recoilThreshold; 00248 if(pCuts) { tcut= std::max(tcut,(*pCuts)[currentMaterialIndex]); } 00249 00250 if(trec > tcut) { 00251 G4ParticleDefinition* ion = theParticleTable->GetIon(iz, ia, 0.0); 00252 G4ThreeVector dir = (direction*sqrt(mom2) - 00253 newDirection*sqrt(finalT*(2*mass + finalT))).unit(); 00254 G4DynamicParticle* newdp = new G4DynamicParticle(ion, dir, trec); 00255 fvect->push_back(newdp); 00256 } else { 00257 fParticleChange->ProposeLocalEnergyDeposit(trec); 00258 fParticleChange->ProposeNonIonizingEnergyDeposit(trec); 00259 } 00260 00261 return; 00262 }
void G4eCoulombScatteringModel::SetLowEnergyThreshold | ( | G4double | val | ) | [inline] |
Definition at line 180 of file G4eCoulombScatteringModel.hh.
References lowEnergyThreshold.
00181 { 00182 lowEnergyThreshold = val; 00183 }
void G4eCoulombScatteringModel::SetRecoilThreshold | ( | G4double | eth | ) | [inline] |
Definition at line 187 of file G4eCoulombScatteringModel.hh.
References recoilThreshold.
00188 { 00189 recoilThreshold = eth; 00190 }
void G4eCoulombScatteringModel::SetupParticle | ( | const G4ParticleDefinition * | ) | [inline, protected] |
Definition at line 168 of file G4eCoulombScatteringModel.hh.
References G4ParticleDefinition::GetPDGMass(), mass, particle, G4WentzelOKandVIxSection::SetupParticle(), and wokvi.
Referenced by ComputeCrossSectionPerAtom(), Initialise(), and SampleSecondaries().
00169 { 00170 // Initialise mass and charge 00171 if(p != particle) { 00172 particle = p; 00173 mass = particle->GetPDGMass(); 00174 wokvi->SetupParticle(p); 00175 } 00176 }
G4double G4eCoulombScatteringModel::cosTetMaxNuc [protected] |
Definition at line 137 of file G4eCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), and G4eCoulombScatteringModel().
G4double G4eCoulombScatteringModel::cosTetMinNuc [protected] |
Definition at line 136 of file G4eCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), G4eCoulombScatteringModel(), and SampleSecondaries().
G4double G4eCoulombScatteringModel::cosThetaMax [protected] |
Definition at line 135 of file G4eCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), and SampleSecondaries().
G4double G4eCoulombScatteringModel::cosThetaMin [protected] |
const G4MaterialCutsCouple* G4eCoulombScatteringModel::currentCouple [protected] |
Definition at line 130 of file G4eCoulombScatteringModel.hh.
Referenced by DefineMaterial(), G4eCoulombScatteringModel(), and Initialise().
const G4Material* G4eCoulombScatteringModel::currentMaterial [protected] |
Definition at line 131 of file G4eCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), DefineMaterial(), and G4eCoulombScatteringModel().
G4int G4eCoulombScatteringModel::currentMaterialIndex [protected] |
Definition at line 132 of file G4eCoulombScatteringModel.hh.
Referenced by DefineMaterial(), G4eCoulombScatteringModel(), and SampleSecondaries().
G4double G4eCoulombScatteringModel::elecRatio [protected] |
Definition at line 139 of file G4eCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), G4eCoulombScatteringModel(), and SampleSecondaries().
G4NistManager* G4eCoulombScatteringModel::fNistManager [protected] |
Definition at line 126 of file G4eCoulombScatteringModel.hh.
Referenced by G4eCoulombScatteringModel().
Definition at line 124 of file G4eCoulombScatteringModel.hh.
Referenced by G4eCoulombScatteringModel(), Initialise(), and SampleSecondaries().
Definition at line 146 of file G4eCoulombScatteringModel.hh.
Referenced by G4eCoulombScatteringModel(), SampleSecondaries(), and SetLowEnergyThreshold().
G4double G4eCoulombScatteringModel::mass [protected] |
Definition at line 140 of file G4eCoulombScatteringModel.hh.
Referenced by G4eCoulombScatteringModel(), Initialise(), SampleSecondaries(), and SetupParticle().
const G4ParticleDefinition* G4eCoulombScatteringModel::particle [protected] |
Definition at line 143 of file G4eCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), G4eCoulombScatteringModel(), SampleSecondaries(), and SetupParticle().
const std::vector<G4double>* G4eCoulombScatteringModel::pCuts [protected] |
Definition at line 128 of file G4eCoulombScatteringModel.hh.
Referenced by G4eCoulombScatteringModel(), Initialise(), and SampleSecondaries().
G4double G4eCoulombScatteringModel::recoilThreshold [protected] |
Definition at line 138 of file G4eCoulombScatteringModel.hh.
Referenced by G4eCoulombScatteringModel(), SampleSecondaries(), and SetRecoilThreshold().
Definition at line 123 of file G4eCoulombScatteringModel.hh.
Referenced by G4eCoulombScatteringModel(), and SampleSecondaries().
const G4ParticleDefinition* G4eCoulombScatteringModel::theProton [protected] |
Definition at line 144 of file G4eCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), and G4eCoulombScatteringModel().
G4WentzelOKandVIxSection* G4eCoulombScatteringModel::wokvi [protected] |
Definition at line 125 of file G4eCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), G4eCoulombScatteringModel(), Initialise(), SampleSecondaries(), SetupParticle(), and ~G4eCoulombScatteringModel().