#include <G4LivermorePhotoElectricModel.hh>
Inheritance diagram for G4LivermorePhotoElectricModel:
Public Member Functions | |
G4LivermorePhotoElectricModel (const G4String &nam="LivermorePhElectric") | |
virtual | ~G4LivermorePhotoElectricModel () |
virtual void | Initialise (const G4ParticleDefinition *, const G4DataVector &) |
virtual G4double | ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0, G4double cut=0, G4double emax=DBL_MAX) |
virtual void | SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy) |
void | SetLimitNumberOfShells (G4int) |
Protected Attributes | |
G4ParticleChangeForGamma * | fParticleChange |
Definition at line 50 of file G4LivermorePhotoElectricModel.hh.
G4LivermorePhotoElectricModel::G4LivermorePhotoElectricModel | ( | const G4String & | nam = "LivermorePhElectric" |
) |
Definition at line 54 of file G4LivermorePhotoElectricModel.cc.
References G4Electron::Electron(), G4cout, G4endl, G4Gamma::Gamma(), G4VEmModel::SetAngularDistribution(), and G4VEmModel::SetDeexcitationFlag().
00056 : G4VEmModel(nam),fParticleChange(0),maxZ(99), 00057 nShellLimit(100),fDeexcitationActive(false),isInitialised(false), 00058 fAtomDeexcitation(0) 00059 { 00060 verboseLevel= 0; 00061 // Verbosity scale: 00062 // 0 = nothing 00063 // 1 = warning for energy non-conservation 00064 // 2 = details of energy budget 00065 // 3 = calculation of cross sections, file openings, sampling of atoms 00066 // 4 = entering in methods 00067 00068 theGamma = G4Gamma::Gamma(); 00069 theElectron = G4Electron::Electron(); 00070 00071 // default generator 00072 SetAngularDistribution(new G4SauterGavrilaAngularDistribution()); 00073 00074 for(G4int i=0; i<maxZ; ++i) { 00075 fCrossSection[i] = 0; 00076 fCrossSectionLE[i] = 0; 00077 fNShells[i] = 0; 00078 fNShellsUsed[i] = 0; 00079 } 00080 00081 if(verboseLevel>0) { 00082 G4cout << "Livermore PhotoElectric is constructed " 00083 << " nShellLimit= " << nShellLimit << G4endl; 00084 } 00085 00086 //Mark this model as "applicable" for atomic deexcitation 00087 SetDeexcitationFlag(true); 00088 }
G4LivermorePhotoElectricModel::~G4LivermorePhotoElectricModel | ( | ) | [virtual] |
Definition at line 92 of file G4LivermorePhotoElectricModel.cc.
00093 { 00094 for(G4int i=0; i<maxZ; ++i) { 00095 delete fCrossSection[i]; 00096 delete fCrossSectionLE[i]; 00097 } 00098 }
G4double G4LivermorePhotoElectricModel::ComputeCrossSectionPerAtom | ( | const G4ParticleDefinition * | , | |
G4double | kinEnergy, | |||
G4double | Z, | |||
G4double | A = 0 , |
|||
G4double | cut = 0 , |
|||
G4double | emax = DBL_MAX | |||
) | [virtual] |
Reimplemented from G4VEmModel.
Definition at line 155 of file G4LivermorePhotoElectricModel.cc.
References G4cout, G4endl, G4lrint(), and G4VEmModel::Value().
00160 { 00161 if (verboseLevel > 3) { 00162 G4cout << "G4LivermorePhotoElectricModel::Calling ComputeCrossSectionPerAtom()" 00163 << " Z= " << ZZ << " R(keV)= " << energy/keV << G4endl; 00164 } 00165 G4double cs = 0.0; 00166 G4double gammaEnergy = energy; 00167 00168 G4int Z = G4lrint(ZZ); 00169 if(Z < 1 || Z >= maxZ) { return cs; } 00170 00171 // element was not initialised 00172 if(!fCrossSection[Z]) { 00173 char* path = getenv("G4LEDATA"); 00174 ReadData(Z, path); 00175 if(!fCrossSection[Z]) { return cs; } 00176 } 00177 00178 G4int idx = fNShells[Z]*6 - 4; 00179 if (gammaEnergy <= (fParam[Z])[idx-1]) { return cs; } 00180 00181 G4double x1 = 1.0/gammaEnergy; 00182 G4double x2 = x1*x1; 00183 G4double x3 = x2*x1; 00184 00185 // parameterisation 00186 if(gammaEnergy >= (fParam[Z])[0]) { 00187 G4double x4 = x2*x2; 00188 cs = x1*((fParam[Z])[idx] + x1*(fParam[Z])[idx+1] 00189 + x2*(fParam[Z])[idx+2] + x3*(fParam[Z])[idx+3] 00190 + x4*(fParam[Z])[idx+4]); 00191 // high energy part 00192 } else if(gammaEnergy >= (fParam[Z])[1]) { 00193 cs = x3*(fCrossSection[Z])->Value(gammaEnergy); 00194 00195 // low energy part 00196 } else { 00197 cs = x3*(fCrossSectionLE[Z])->Value(gammaEnergy); 00198 } 00199 if (verboseLevel > 1) { 00200 G4cout << "LivermorePhotoElectricModel: E(keV)= " << gammaEnergy/keV 00201 << " Z= " << Z << " cross(barn)= " << cs/barn << G4endl; 00202 } 00203 return cs; 00204 }
void G4LivermorePhotoElectricModel::Initialise | ( | const G4ParticleDefinition * | , | |
const G4DataVector & | ||||
) | [virtual] |
Implements G4VEmModel.
Definition at line 103 of file G4LivermorePhotoElectricModel.cc.
References G4LossTableManager::AtomDeexcitation(), fParticleChange, G4cout, G4endl, G4Material::GetElementVector(), G4MaterialCutsCouple::GetMaterial(), G4ProductionCutsTable::GetMaterialCutsCouple(), G4Material::GetNumberOfElements(), G4VEmModel::GetParticleChangeForGamma(), G4ProductionCutsTable::GetProductionCutsTable(), G4ProductionCutsTable::GetTableSize(), G4LossTableManager::Instance(), and G4VAtomDeexcitation::IsFluoActive().
00105 { 00106 if (verboseLevel > 2) { 00107 G4cout << "Calling G4LivermorePhotoElectricModel::Initialise()" << G4endl; 00108 } 00109 00110 char* path = getenv("G4LEDATA"); 00111 00112 G4ProductionCutsTable* theCoupleTable = 00113 G4ProductionCutsTable::GetProductionCutsTable(); 00114 G4int numOfCouples = theCoupleTable->GetTableSize(); 00115 00116 for(G4int i=0; i<numOfCouples; ++i) 00117 { 00118 const G4MaterialCutsCouple* couple = theCoupleTable->GetMaterialCutsCouple(i); 00119 const G4Material* material = couple->GetMaterial(); 00120 const G4ElementVector* theElementVector = material->GetElementVector(); 00121 G4int nelm = material->GetNumberOfElements(); 00122 00123 for (G4int j=0; j<nelm; ++j) 00124 { 00125 G4int Z = (G4int)(*theElementVector)[j]->GetZ(); 00126 if(Z < 1) { Z = 1; } 00127 else if(Z > maxZ) { Z = maxZ; } 00128 if(!fCrossSection[Z]) { ReadData(Z, path); } 00129 } 00130 } 00131 // 00132 if (verboseLevel > 2) { 00133 G4cout << "Loaded cross section files for LivermorePhotoElectric model" 00134 << G4endl; 00135 } 00136 if(!isInitialised) { 00137 isInitialised = true; 00138 fParticleChange = GetParticleChangeForGamma(); 00139 00140 fAtomDeexcitation = G4LossTableManager::Instance()->AtomDeexcitation(); 00141 } 00142 fDeexcitationActive = false; 00143 if(fAtomDeexcitation) { 00144 fDeexcitationActive = fAtomDeexcitation->IsFluoActive(); 00145 } 00146 00147 if (verboseLevel > 0) { 00148 G4cout << "LivermorePhotoElectric model is initialized " << G4endl 00149 << G4endl; 00150 } 00151 }
void G4LivermorePhotoElectricModel::SampleSecondaries | ( | std::vector< G4DynamicParticle * > * | , | |
const G4MaterialCutsCouple * | , | |||
const G4DynamicParticle * | , | |||
G4double | tmin, | |||
G4double | maxEnergy | |||
) | [virtual] |
Implements G4VEmModel.
Definition at line 209 of file G4LivermorePhotoElectricModel.cc.
References G4InuclSpecialFunctions::bindingEnergy(), G4VAtomDeexcitation::CheckDeexcitationActiveRegion(), fParticleChange, fStopAndKill, G4cout, G4endl, G4lrint(), G4UniformRand, G4VAtomDeexcitation::GenerateParticles(), G4VEmModel::GetAngularDistribution(), G4VAtomDeexcitation::GetAtomicShell(), G4ElementData::GetComponentID(), G4MaterialCutsCouple::GetIndex(), G4DynamicParticle::GetKineticEnergy(), G4MaterialCutsCouple::GetMaterial(), G4DynamicParticle::GetMomentumDirection(), G4ElementData::GetValueForComponent(), G4Element::GetZ(), G4InuclParticleNames::nn, G4VParticleChange::ProposeLocalEnergyDeposit(), G4VParticleChange::ProposeTrackStatus(), G4VEmAngularDistribution::SampleDirection(), G4VEmModel::SelectRandomAtom(), G4ParticleChangeForGamma::SetProposedKineticEnergy(), and G4VEmModel::Value().
00215 { 00216 G4double gammaEnergy = aDynamicGamma->GetKineticEnergy(); 00217 if (verboseLevel > 3) { 00218 G4cout << "G4LivermorePhotoElectricModel::SampleSecondaries() Egamma(keV)= " 00219 << gammaEnergy/keV << G4endl; 00220 } 00221 00222 // kill incident photon 00223 fParticleChange->SetProposedKineticEnergy(0.); 00224 fParticleChange->ProposeTrackStatus(fStopAndKill); 00225 00226 // Returns the normalized direction of the momentum 00227 G4ThreeVector photonDirection = aDynamicGamma->GetMomentumDirection(); 00228 00229 // Select randomly one element in the current material 00230 //G4cout << "Select random atom Egamma(keV)= " << gammaEnergy/keV << G4endl; 00231 const G4Element* elm = SelectRandomAtom(couple->GetMaterial(),theGamma, 00232 gammaEnergy); 00233 G4int Z = G4lrint(elm->GetZ()); 00234 00235 // Select the ionised shell in the current atom according to shell 00236 // cross sections 00237 // G4cout << "Select random shell Z= " << Z << G4endl; 00238 00239 if(Z >= maxZ) { Z = maxZ-1; } 00240 00241 // element was not initialised 00242 if(!fCrossSection[Z]) { 00243 char* path = getenv("G4LEDATA"); 00244 ReadData(Z, path); 00245 if(!fCrossSection[Z]) { 00246 fParticleChange->ProposeLocalEnergyDeposit(gammaEnergy); 00247 return; 00248 } 00249 } 00250 00251 // shell index 00252 size_t shellIdx = 0; 00253 size_t nn = fNShellsUsed[Z]; 00254 00255 if(nn > 1) { 00256 if(gammaEnergy >= (fParam[Z])[0]) { 00257 G4double x1 = 1.0/gammaEnergy; 00258 G4double x2 = x1*x1; 00259 G4double x3 = x2*x1; 00260 G4double x4 = x3*x1; 00261 G4int idx = nn*6 - 4; 00262 // when do sampling common factors are not taken into account 00263 // so cross section is not real 00264 G4double cs0 = G4UniformRand()*((fParam[Z])[idx] + x1*(fParam[Z])[idx+1] 00265 + x2*(fParam[Z])[idx+2] 00266 + x3*(fParam[Z])[idx+3] 00267 + x4*(fParam[Z])[idx+4]); 00268 for(shellIdx=0; shellIdx<nn; ++shellIdx) { 00269 idx = shellIdx*6 + 2; 00270 if(gammaEnergy > (fParam[Z])[idx-1]) { 00271 G4double cs = (fParam[Z])[idx] + x1*(fParam[Z])[idx+1] 00272 + x2*(fParam[Z])[idx+2] + x3*(fParam[Z])[idx+3] 00273 + x4*(fParam[Z])[idx+4]; 00274 if(cs >= cs0) { break; } 00275 } 00276 } 00277 if(shellIdx >= nn) { shellIdx = nn-1; } 00278 00279 } else { 00280 00281 // when do sampling common factors are not taken into account 00282 // so cross section is not real 00283 G4double cs = G4UniformRand(); 00284 00285 if(gammaEnergy >= (fParam[Z])[1]) { 00286 cs *= (fCrossSection[Z])->Value(gammaEnergy); 00287 } else { 00288 cs *= (fCrossSectionLE[Z])->Value(gammaEnergy); 00289 } 00290 00291 for(size_t j=0; j<nn; ++j) { 00292 shellIdx = (size_t)fShellCrossSection.GetComponentID(Z, j); 00293 if(gammaEnergy > (fParam[Z])[6*shellIdx+1]) { 00294 cs -= fShellCrossSection.GetValueForComponent(Z, j, gammaEnergy); 00295 } 00296 if(cs <= 0.0 || j+1 == nn) { break; } 00297 } 00298 } 00299 } 00300 00301 G4double bindingEnergy = (fParam[Z])[shellIdx*6 + 1]; 00302 //G4cout << "Z= " << Z << " shellIdx= " << shellIdx 00303 // << " nShells= " << fNShells[Z] 00304 // << " Ebind(keV)= " << bindingEnergy/keV 00305 // << " Egamma(keV)= " << gammaEnergy/keV << G4endl; 00306 00307 const G4AtomicShell* shell = 0; 00308 00309 // no de-excitation from the last shell 00310 if(fDeexcitationActive && shellIdx + 1 < nn) { 00311 G4AtomicShellEnumerator as = G4AtomicShellEnumerator(shellIdx); 00312 shell = fAtomDeexcitation->GetAtomicShell(Z, as); 00313 } 00314 00315 // If binding energy of the selected shell is larger than photon energy 00316 // do not generate secondaries 00317 if(gammaEnergy < bindingEnergy) { 00318 fParticleChange->ProposeLocalEnergyDeposit(gammaEnergy); 00319 return; 00320 } 00321 00322 // Primary outcoming electron 00323 G4double eKineticEnergy = gammaEnergy - bindingEnergy; 00324 G4double edep = bindingEnergy; 00325 00326 // Calculate direction of the photoelectron 00327 G4ThreeVector electronDirection = 00328 GetAngularDistribution()->SampleDirection(aDynamicGamma, 00329 eKineticEnergy, 00330 shellIdx, 00331 couple->GetMaterial()); 00332 00333 // The electron is created 00334 G4DynamicParticle* electron = new G4DynamicParticle (theElectron, 00335 electronDirection, 00336 eKineticEnergy); 00337 fvect->push_back(electron); 00338 00339 // Sample deexcitation 00340 if(shell) { 00341 G4int index = couple->GetIndex(); 00342 if(fAtomDeexcitation->CheckDeexcitationActiveRegion(index)) { 00343 size_t nbefore = fvect->size(); 00344 00345 fAtomDeexcitation->GenerateParticles(fvect, shell, Z, index); 00346 size_t nafter = fvect->size(); 00347 if(nafter > nbefore) { 00348 for (size_t j=nbefore; j<nafter; ++j) { 00349 edep -= ((*fvect)[j])->GetKineticEnergy(); 00350 } 00351 } 00352 } 00353 } 00354 // energy balance - excitation energy left 00355 if(edep > 0.0) { 00356 fParticleChange->ProposeLocalEnergyDeposit(edep); 00357 } 00358 }
void G4LivermorePhotoElectricModel::SetLimitNumberOfShells | ( | G4int | ) | [inline] |
Definition at line 110 of file G4LivermorePhotoElectricModel.hh.
00111 { 00112 nShellLimit = n; 00113 }
Definition at line 80 of file G4LivermorePhotoElectricModel.hh.
Referenced by Initialise(), and SampleSecondaries().