#include <G4hCoulombScatteringModel.hh>
Inheritance diagram for G4hCoulombScatteringModel:
Definition at line 67 of file G4hCoulombScatteringModel.hh.
G4hCoulombScatteringModel::G4hCoulombScatteringModel | ( | const G4String & | nam = "eCoulombScattering" |
) |
Definition at line 69 of file G4hCoulombScatteringModel.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.
00070 : G4VEmModel(nam), 00071 cosThetaMin(1.0), 00072 cosThetaMax(-1.0), 00073 isInitialised(false) 00074 { 00075 fParticleChange = 0; 00076 fNistManager = G4NistManager::Instance(); 00077 theParticleTable = G4ParticleTable::GetParticleTable(); 00078 theProton = G4Proton::Proton(); 00079 currentMaterial = 0; 00080 00081 pCuts = 0; 00082 00083 lowEnergyThreshold = 1*keV; // particle will be killed for lower energy 00084 recoilThreshold = 0.*keV; // by default does not work 00085 00086 particle = 0; 00087 currentCouple = 0; 00088 wokvi = new G4WentzelVIRelXSection(); 00089 00090 currentMaterialIndex = 0; 00091 00092 cosTetMinNuc = 1.0; 00093 cosTetMaxNuc = -1.0; 00094 elecRatio = 0.0; 00095 mass = proton_mass_c2; 00096 }
G4hCoulombScatteringModel::~G4hCoulombScatteringModel | ( | ) | [virtual] |
Definition at line 100 of file G4hCoulombScatteringModel.cc.
References wokvi.
00101 { 00102 delete wokvi; 00103 }
G4double G4hCoulombScatteringModel::ComputeCrossSectionPerAtom | ( | const G4ParticleDefinition * | , | |
G4double | kinEnergy, | |||
G4double | Z, | |||
G4double | A, | |||
G4double | cut, | |||
G4double | emax | |||
) | [virtual] |
Reimplemented from G4VEmModel.
Definition at line 136 of file G4hCoulombScatteringModel.cc.
References G4WentzelVIRelXSection::ComputeElectronCrossSection(), G4WentzelVIRelXSection::ComputeNuclearCrossSection(), cosTetMaxNuc, cosTetMinNuc, cosThetaMax, G4VEmModel::CurrentCouple(), currentMaterial, DefineMaterial(), elecRatio, particle, G4WentzelVIRelXSection::SetupKinematic(), SetupParticle(), G4WentzelVIRelXSection::SetupTarget(), theProton, and wokvi.
Referenced by SampleSecondaries().
00141 { 00142 //G4cout << "### G4hCoulombScatteringModel::ComputeCrossSectionPerAtom for " 00143 // << p->GetParticleName()<<" Z= "<<Z<<" e(MeV)= "<< kinEnergy/MeV << G4endl; 00144 G4double cross = 0.0; 00145 if(p != particle) { SetupParticle(p); } 00146 00147 // cross section is set to zero to avoid problems in sample secondary 00148 if(kinEnergy <= 0.0) { return cross; } 00149 DefineMaterial(CurrentCouple()); 00150 cosTetMinNuc = wokvi->SetupKinematic(kinEnergy, currentMaterial); 00151 if(cosThetaMax < cosTetMinNuc) { 00152 G4int iz = G4int(Z); 00153 cosTetMinNuc = wokvi->SetupTarget(iz, cutEnergy); 00154 cosTetMaxNuc = cosThetaMax; 00155 if(iz == 1 && cosTetMaxNuc < 0.0 && particle == theProton) { 00156 cosTetMaxNuc = 0.0; 00157 } 00158 cross = wokvi->ComputeNuclearCrossSection(cosTetMinNuc, cosTetMaxNuc); 00159 elecRatio = wokvi->ComputeElectronCrossSection(cosTetMinNuc, cosThetaMax); 00160 cross += elecRatio; 00161 if(cross > 0.0) { elecRatio /= cross; } 00162 } 00163 /* 00164 if(p->GetParticleName() == "e-") 00165 G4cout << "e(MeV)= " << kinEnergy/MeV << " cross(b)= " << cross/barn 00166 << " 1-cosTetMinNuc= " << 1-cosTetMinNuc 00167 << " 1-cosTetMaxNuc= " << 1-cosTetMaxNuc 00168 << G4endl; 00169 */ 00170 return cross; 00171 }
void G4hCoulombScatteringModel::DefineMaterial | ( | const G4MaterialCutsCouple * | ) | [inline, protected] |
Definition at line 145 of file G4hCoulombScatteringModel.hh.
References currentCouple, currentMaterial, currentMaterialIndex, G4MaterialCutsCouple::GetIndex(), and G4MaterialCutsCouple::GetMaterial().
Referenced by ComputeCrossSectionPerAtom(), and SampleSecondaries().
00146 { 00147 if(cup != currentCouple) { 00148 currentCouple = cup; 00149 currentMaterial = cup->GetMaterial(); 00150 currentMaterialIndex = currentCouple->GetIndex(); 00151 } 00152 }
void G4hCoulombScatteringModel::Initialise | ( | const G4ParticleDefinition * | , | |
const G4DataVector & | ||||
) | [virtual] |
Implements G4VEmModel.
Definition at line 107 of file G4hCoulombScatteringModel.cc.
References cosThetaMin, currentCouple, fParticleChange, G4ProductionCutsTable::GetEnergyCutsVector(), G4VEmModel::GetParticleChangeForGamma(), G4ProductionCutsTable::GetProductionCutsTable(), G4WentzelVIRelXSection::Initialise(), G4VEmModel::InitialiseElementSelectors(), mass, pCuts, G4VEmModel::PolarAngleLimit(), SetupParticle(), and wokvi.
00109 { 00110 SetupParticle(p); 00111 currentCouple = 0; 00112 cosThetaMin = cos(PolarAngleLimit()); 00113 wokvi->Initialise(p, cosThetaMin); 00114 /* 00115 G4cout << "G4hCoulombScatteringModel: " << particle->GetParticleName() 00116 << " 1-cos(ThetaLimit)= " << 1 - cosThetaMin 00117 << " cos(thetaMax)= " << cosThetaMax 00118 << G4endl; 00119 */ 00120 pCuts = G4ProductionCutsTable::GetProductionCutsTable()->GetEnergyCutsVector(3); 00121 //G4cout << "!!! G4hCoulombScatteringModel::Initialise for " 00122 // << p->GetParticleName() << " cos(TetMin)= " << cosThetaMin 00123 // << " cos(TetMax)= " << cosThetaMax <<G4endl; 00124 // G4cout << "cut0= " << cuts[0] << " cut1= " << cuts[1] << G4endl; 00125 if(!isInitialised) { 00126 isInitialised = true; 00127 fParticleChange = GetParticleChangeForGamma(); 00128 } 00129 if(mass < GeV) { 00130 InitialiseElementSelectors(p,cuts); 00131 } 00132 }
void G4hCoulombScatteringModel::SampleSecondaries | ( | std::vector< G4DynamicParticle * > * | , | |
const G4MaterialCutsCouple * | , | |||
const G4DynamicParticle * | , | |||
G4double | tmin, | |||
G4double | maxEnergy | |||
) | [virtual] |
Implements G4VEmModel.
Definition at line 175 of file G4hCoulombScatteringModel.cc.
References ComputeCrossSectionPerAtom(), cosTetMinNuc, cosThetaMax, currentMaterialIndex, DefineMaterial(), elecRatio, fParticleChange, G4DynamicParticle::GetDefinition(), G4ParticleTable::GetIon(), G4DynamicParticle::GetKineticEnergy(), G4DynamicParticle::GetMomentumDirection(), G4WentzelVIRelXSection::GetMomentumSquare(), G4NucleiProperties::GetNuclearMass(), G4Element::GetZ(), lowEnergyThreshold, mass, particle, pCuts, G4VParticleChange::ProposeLocalEnergyDeposit(), G4ParticleChangeForGamma::ProposeMomentumDirection(), G4VParticleChange::ProposeNonIonizingEnergyDeposit(), recoilThreshold, G4WentzelVIRelXSection::SampleSingleScattering(), G4VEmModel::SelectIsotopeNumber(), G4VEmModel::SelectRandomAtom(), G4ParticleChangeForGamma::SetProposedKineticEnergy(), G4WentzelVIRelXSection::SetTargetMass(), SetupParticle(), theParticleTable, and wokvi.
00181 { 00182 G4double kinEnergy = dp->GetKineticEnergy(); 00183 00184 // absorb particle below low-energy limit to avoid situation 00185 // when a particle has no energy loss 00186 if(kinEnergy < lowEnergyThreshold) { 00187 fParticleChange->SetProposedKineticEnergy(0.0); 00188 fParticleChange->ProposeLocalEnergyDeposit(kinEnergy); 00189 fParticleChange->ProposeNonIonizingEnergyDeposit(kinEnergy); 00190 return; 00191 } 00192 SetupParticle(dp->GetDefinition()); 00193 DefineMaterial(couple); 00194 00195 //G4cout << "G4hCoulombScatteringModel::SampleSecondaries e(MeV)= " 00196 // << kinEnergy << " " << particle->GetParticleName() 00197 // << " cut= " << cutEnergy<< G4endl; 00198 00199 // Choose nucleus 00200 const G4Element* currentElement = 00201 SelectRandomAtom(couple,particle,kinEnergy,cutEnergy,kinEnergy); 00202 00203 G4double Z = currentElement->GetZ(); 00204 00205 if(ComputeCrossSectionPerAtom(particle,kinEnergy, Z, 00206 kinEnergy, cutEnergy, kinEnergy) == 0.0) 00207 { return; } 00208 00209 G4int iz = G4int(Z); 00210 G4int ia = SelectIsotopeNumber(currentElement); 00211 G4double targetMass = G4NucleiProperties::GetNuclearMass(ia, iz); 00212 wokvi->SetTargetMass(targetMass); 00213 00214 G4ThreeVector newDirection = 00215 wokvi->SampleSingleScattering(cosTetMinNuc, cosThetaMax, elecRatio); 00216 G4double cost = newDirection.z(); 00217 00218 G4ThreeVector direction = dp->GetMomentumDirection(); 00219 newDirection.rotateUz(direction); 00220 00221 fParticleChange->ProposeMomentumDirection(newDirection); 00222 00223 // recoil sampling assuming a small recoil 00224 // and first order correction to primary 4-momentum 00225 G4double mom2 = wokvi->GetMomentumSquare(); 00226 G4double trec = mom2*(1.0 - cost)/(targetMass + (mass + kinEnergy)*(1.0 - cost)); 00227 G4double finalT = kinEnergy - trec; 00228 //G4cout<<"G4hCoulombScatteringModel: finalT= "<<finalT<<" Trec= "<<trec<<G4endl; 00229 if(finalT <= lowEnergyThreshold) { 00230 trec = kinEnergy; 00231 finalT = 0.0; 00232 } 00233 00234 fParticleChange->SetProposedKineticEnergy(finalT); 00235 G4double tcut = recoilThreshold; 00236 if(pCuts) { tcut= std::max(tcut,(*pCuts)[currentMaterialIndex]); } 00237 00238 if(trec > tcut) { 00239 G4ParticleDefinition* ion = theParticleTable->GetIon(iz, ia, 0.0); 00240 G4ThreeVector dir = (direction*sqrt(mom2) - 00241 newDirection*sqrt(finalT*(2*mass + finalT))).unit(); 00242 G4DynamicParticle* newdp = new G4DynamicParticle(ion, dir, trec); 00243 fvect->push_back(newdp); 00244 } else { 00245 fParticleChange->ProposeLocalEnergyDeposit(trec); 00246 fParticleChange->ProposeNonIonizingEnergyDeposit(trec); 00247 } 00248 00249 return; 00250 }
void G4hCoulombScatteringModel::SetLowEnergyThreshold | ( | G4double | val | ) | [inline] |
Definition at line 169 of file G4hCoulombScatteringModel.hh.
References lowEnergyThreshold.
00170 { 00171 lowEnergyThreshold = val; 00172 }
void G4hCoulombScatteringModel::SetRecoilThreshold | ( | G4double | eth | ) | [inline] |
Definition at line 176 of file G4hCoulombScatteringModel.hh.
References recoilThreshold.
00177 { 00178 recoilThreshold = eth; 00179 }
void G4hCoulombScatteringModel::SetupParticle | ( | const G4ParticleDefinition * | ) | [inline, protected] |
Definition at line 157 of file G4hCoulombScatteringModel.hh.
References G4ParticleDefinition::GetPDGMass(), mass, particle, G4WentzelVIRelXSection::SetupParticle(), and wokvi.
Referenced by ComputeCrossSectionPerAtom(), Initialise(), and SampleSecondaries().
00158 { 00159 // Initialise mass and charge 00160 if(p != particle) { 00161 particle = p; 00162 mass = particle->GetPDGMass(); 00163 wokvi->SetupParticle(p); 00164 } 00165 }
G4double G4hCoulombScatteringModel::cosTetMaxNuc [protected] |
Definition at line 126 of file G4hCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), and G4hCoulombScatteringModel().
G4double G4hCoulombScatteringModel::cosTetMinNuc [protected] |
Definition at line 125 of file G4hCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), G4hCoulombScatteringModel(), and SampleSecondaries().
G4double G4hCoulombScatteringModel::cosThetaMax [protected] |
Definition at line 124 of file G4hCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), and SampleSecondaries().
G4double G4hCoulombScatteringModel::cosThetaMin [protected] |
const G4MaterialCutsCouple* G4hCoulombScatteringModel::currentCouple [protected] |
Definition at line 119 of file G4hCoulombScatteringModel.hh.
Referenced by DefineMaterial(), G4hCoulombScatteringModel(), and Initialise().
const G4Material* G4hCoulombScatteringModel::currentMaterial [protected] |
Definition at line 120 of file G4hCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), DefineMaterial(), and G4hCoulombScatteringModel().
G4int G4hCoulombScatteringModel::currentMaterialIndex [protected] |
Definition at line 121 of file G4hCoulombScatteringModel.hh.
Referenced by DefineMaterial(), G4hCoulombScatteringModel(), and SampleSecondaries().
G4double G4hCoulombScatteringModel::elecRatio [protected] |
Definition at line 128 of file G4hCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), G4hCoulombScatteringModel(), and SampleSecondaries().
G4NistManager* G4hCoulombScatteringModel::fNistManager [protected] |
Definition at line 115 of file G4hCoulombScatteringModel.hh.
Referenced by G4hCoulombScatteringModel().
Definition at line 113 of file G4hCoulombScatteringModel.hh.
Referenced by G4hCoulombScatteringModel(), Initialise(), and SampleSecondaries().
Definition at line 135 of file G4hCoulombScatteringModel.hh.
Referenced by G4hCoulombScatteringModel(), SampleSecondaries(), and SetLowEnergyThreshold().
G4double G4hCoulombScatteringModel::mass [protected] |
Definition at line 129 of file G4hCoulombScatteringModel.hh.
Referenced by G4hCoulombScatteringModel(), Initialise(), SampleSecondaries(), and SetupParticle().
const G4ParticleDefinition* G4hCoulombScatteringModel::particle [protected] |
Definition at line 132 of file G4hCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), G4hCoulombScatteringModel(), SampleSecondaries(), and SetupParticle().
const std::vector<G4double>* G4hCoulombScatteringModel::pCuts [protected] |
Definition at line 117 of file G4hCoulombScatteringModel.hh.
Referenced by G4hCoulombScatteringModel(), Initialise(), and SampleSecondaries().
G4double G4hCoulombScatteringModel::recoilThreshold [protected] |
Definition at line 127 of file G4hCoulombScatteringModel.hh.
Referenced by G4hCoulombScatteringModel(), SampleSecondaries(), and SetRecoilThreshold().
Definition at line 112 of file G4hCoulombScatteringModel.hh.
Referenced by G4hCoulombScatteringModel(), and SampleSecondaries().
const G4ParticleDefinition* G4hCoulombScatteringModel::theProton [protected] |
Definition at line 133 of file G4hCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), and G4hCoulombScatteringModel().
G4WentzelVIRelXSection* G4hCoulombScatteringModel::wokvi [protected] |
Definition at line 114 of file G4hCoulombScatteringModel.hh.
Referenced by ComputeCrossSectionPerAtom(), G4hCoulombScatteringModel(), Initialise(), SampleSecondaries(), SetupParticle(), and ~G4hCoulombScatteringModel().