Geant4 11.4.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4LivermorePhotoElectricModel Class Reference

#include <G4LivermorePhotoElectricModel.hh>

Inheritance diagram for G4LivermorePhotoElectricModel:

Public Member Functions

 G4LivermorePhotoElectricModel (const G4String &nam="LivermorePhElectric")
 ~G4LivermorePhotoElectricModel () override
void Initialise (const G4ParticleDefinition *, const G4DataVector &) override
G4double CrossSectionPerVolume (const G4Material *, const G4ParticleDefinition *, G4double energy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX) override
G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double energy, G4double Z, G4double A=0, G4double cut=0, G4double emax=DBL_MAX) override
void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy) override
void SetLimitNumberOfShells (G4int n)
G4double GetBindingEnergy (G4int Z, G4int shell)
G4LivermorePhotoElectricModel & operator= (const G4LivermorePhotoElectricModel &right)=delete
 G4LivermorePhotoElectricModel (const G4LivermorePhotoElectricModel &)=delete
Public Member Functions inherited from G4VEmModel
 G4VEmModel (const G4String &nam)
virtual ~G4VEmModel ()
virtual void InitialiseLocal (const G4ParticleDefinition *, G4VEmModel *masterModel)
virtual void InitialiseForMaterial (const G4ParticleDefinition *, const G4Material *)
virtual void InitialiseForElement (const G4ParticleDefinition *, G4int Z)
virtual G4double ComputeDEDXPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
virtual G4double GetPartialCrossSection (const G4Material *, G4int level, const G4ParticleDefinition *, G4double kineticEnergy)
virtual G4double ComputeCrossSectionPerShell (const G4ParticleDefinition *, G4int Z, G4int shellIdx, G4double kinEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
virtual G4double ChargeSquareRatio (const G4Track &)
virtual G4double GetChargeSquareRatio (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
virtual G4double GetParticleCharge (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
virtual void StartTracking (G4Track *)
virtual void CorrectionsAlongStep (const G4Material *, const G4ParticleDefinition *, const G4double kinEnergy, const G4double cutEnergy, const G4double &length, G4double &eloss)
virtual G4double Value (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy)
virtual G4double MinPrimaryEnergy (const G4Material *, const G4ParticleDefinition *, G4double cut=0.0)
virtual G4double MinEnergyCut (const G4ParticleDefinition *, const G4MaterialCutsCouple *)
virtual void SetupForMaterial (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
virtual void DefineForRegion (const G4Region *)
virtual void FillNumberOfSecondaries (G4int &numberOfTriplets, G4int &numberOfRecoil)
virtual void ModelDescription (std::ostream &outFile) const
void InitialiseElementSelectors (const G4ParticleDefinition *, const G4DataVector &)
std::vector< G4EmElementSelector * > * GetElementSelectors ()
void SetElementSelectors (std::vector< G4EmElementSelector * > *)
G4double ComputeDEDX (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
G4double CrossSection (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
G4double ComputeMeanFreePath (const G4ParticleDefinition *, G4double kineticEnergy, const G4Material *, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, const G4Element *, G4double kinEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
const G4Element * SelectRandomAtom (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
const G4Element * SelectTargetAtom (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double logKineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
const G4Element * SelectRandomAtom (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
const G4Element * GetCurrentElement (const G4Material *mat=nullptr) const
G4int SelectRandomAtomNumber (const G4Material *) const
const G4Isotope * GetCurrentIsotope (const G4Element *elm=nullptr) const
G4int SelectIsotopeNumber (const G4Element *) const
void SetParticleChange (G4VParticleChange *, G4VEmFluctuationModel *f=nullptr)
void SetCrossSectionTable (G4PhysicsTable *, G4bool isLocal)
G4ElementData * GetElementData ()
G4PhysicsTable * GetCrossSectionTable ()
G4VEmFluctuationModel * GetModelOfFluctuations ()
G4VEmAngularDistribution * GetAngularDistribution ()
G4VEmModel * GetTripletModel ()
void SetTripletModel (G4VEmModel *)
void SetAngularDistribution (G4VEmAngularDistribution *)
G4double HighEnergyLimit () const
G4double LowEnergyLimit () const
G4double HighEnergyActivationLimit () const
G4double LowEnergyActivationLimit () const
G4double PolarAngleLimit () const
G4double SecondaryThreshold () const
G4bool DeexcitationFlag () const
G4bool ForceBuildTableFlag () const
G4bool UseAngularGeneratorFlag () const
void SetAngularGeneratorFlag (G4bool)
void SetHighEnergyLimit (G4double)
void SetLowEnergyLimit (G4double)
void SetActivationHighEnergyLimit (G4double)
void SetActivationLowEnergyLimit (G4double)
G4bool IsActive (G4double kinEnergy) const
void SetPolarAngleLimit (G4double)
void SetSecondaryThreshold (G4double)
void SetDeexcitationFlag (G4bool val)
void SetForceBuildTable (G4bool val)
void SetFluctuationFlag (G4bool val)
G4bool IsMaster () const
void SetUseBaseMaterials (G4bool val)
G4bool UseBaseMaterials () const
G4double MaxSecondaryKinEnergy (const G4DynamicParticle *dynParticle)
const G4String & GetName () const
void SetCurrentCouple (const G4MaterialCutsCouple *)
G4bool IsLocked () const
void SetLocked (G4bool)
void SetLPMFlag (G4bool)
void SetMasterThread (G4bool)
G4VEmModel & operator= (const G4VEmModel &right)=delete
 G4VEmModel (const G4VEmModel &)=delete

Protected Attributes

G4ParticleChangeForGamma * fParticleChange = nullptr
Protected Attributes inherited from G4VEmModel
G4ElementData * fElementData = nullptr
G4VParticleChange * pParticleChange = nullptr
G4PhysicsTable * xSectionTable = nullptr
const G4Material * pBaseMaterial = nullptr
const std::vector< G4double > * theDensityFactor = nullptr
const std::vector< G4int > * theDensityIdx = nullptr
G4double inveplus
G4double pFactor = 1.0
std::size_t currentCoupleIndex = 0
std::size_t basedCoupleIndex = 0
G4bool lossFlucFlag = true

Additional Inherited Members

Protected Member Functions inherited from G4VEmModel
G4ParticleChangeForLoss * GetParticleChangeForLoss ()
G4ParticleChangeForGamma * GetParticleChangeForGamma ()
virtual G4double MaxSecondaryEnergy (const G4ParticleDefinition *, G4double kineticEnergy)
const G4MaterialCutsCouple * CurrentCouple () const
void SetCurrentElement (const G4Element *)

Detailed Description

Definition at line 49 of file G4LivermorePhotoElectricModel.hh.

Constructor & Destructor Documentation

◆ G4LivermorePhotoElectricModel() [1/2]

G4LivermorePhotoElectricModel::G4LivermorePhotoElectricModel ( const G4String & nam = "LivermorePhElectric")
explicit

Definition at line 71 of file G4LivermorePhotoElectricModel.cc.

72 : G4VEmModel(nam)
73{
74 verboseLevel = 0;
75 // Verbosity scale:
76 // 0 = nothing
77 // 1 = warning for energy non-conservation
78 // 2 = details of energy budget
79 // 3 = calculation of cross sections, file openings, sampling of atoms
80 // 4 = entering in methods
81
82 theGamma = G4Gamma::Gamma();
83 theElectron = G4Electron::Electron();
84
85 // default generator
86 SetAngularDistribution(new G4SauterGavrilaAngularDistribution());
87
88 if (verboseLevel > 0) {
89 G4cout << "Livermore PhotoElectric is constructed "
90 << " nShellLimit= " << nShellLimit << G4endl;
91 }
92
93 // Mark this model as "applicable" for atomic deexcitation
94 SetDeexcitationFlag(true);
95
96 // For water
97 fSandiaCof.resize(4, 0.0);
98
99 FindDirectoryPath();
100
101 if (fCrossSection == nullptr) {
102 fCrossSection = new G4ElementData(ZMAXPE);
103 fCrossSection->SetName("PhotoEffXS");
104 fCrossSectionLE = new G4ElementData(ZMAXPE);
105 fCrossSectionLE->SetName("PhotoEffLowXS");
106 for (G4int i = 0; i < ZMAXPE; ++i) {
107 fParamHigh[i] = nullptr;
108 fParamLow[i] = nullptr;
109 fNShells[i] = 0;
110 fNShellsUsed[i] = 0;
111 }
112 }
113}
G4int
int G4int
Definition G4Types.hh:85
G4endl
#define G4endl
Definition G4ios.hh:67
G4cout
G4GLOB_DLL std::ostream G4cout
G4Electron::Electron
static G4Electron * Electron()
Definition G4Electron.cc:91
G4Gamma::Gamma
static G4Gamma * Gamma()
Definition G4Gamma.cc:81
G4VEmModel::SetDeexcitationFlag
void SetDeexcitationFlag(G4bool val)
Definition G4VEmModel.hh:783
G4VEmModel::G4VEmModel
G4VEmModel(const G4String &nam)
Definition G4VEmModel.cc:67
G4VEmModel::SetAngularDistribution
void SetAngularDistribution(G4VEmAngularDistribution *)
Definition G4VEmModel.hh:609

Referenced by G4LivermorePhotoElectricModel(), and operator=().

◆ ~G4LivermorePhotoElectricModel()

G4LivermorePhotoElectricModel::~G4LivermorePhotoElectricModel ( )
override

Definition at line 117 of file G4LivermorePhotoElectricModel.cc.

118{
119 if (isInitializer) {
120 for (G4int i = 0; i < ZMAXPE; ++i) {
121 delete fParamHigh[i];
122 fParamHigh[i] = nullptr;
123 delete fParamLow[i];
124 fParamLow[i] = nullptr;
125 }
126 }
127}

◆ G4LivermorePhotoElectricModel() [2/2]

G4LivermorePhotoElectricModel::G4LivermorePhotoElectricModel ( const G4LivermorePhotoElectricModel & )
delete

Member Function Documentation

◆ ComputeCrossSectionPerAtom()

G4double G4LivermorePhotoElectricModel::ComputeCrossSectionPerAtom ( const G4ParticleDefinition * ,
G4double energy,
G4double Z,
G4double A = 0,
G4double cut = 0,
G4double emax = DBL_MAX )
overridevirtual

Reimplemented from G4VEmModel.

Definition at line 217 of file G4LivermorePhotoElectricModel.cc.

221{
222 if (verboseLevel > 3) {
223 G4cout << "\n G4LivermorePhotoElectricModel::ComputeCrossSectionPerAtom():"
224 << " Z= " << ZZ << " R(keV)= " << energy / keV << G4endl;
225 }
226 G4double cs = 0.0;
227 G4int Z = G4lrint(ZZ);
228 if (Z >= ZMAXPE || Z <= 0) {
229 return cs;
230 }
231 // if element was not initialised
232
233 // do initialisation safely for MT mode
234 if (fCrossSection->GetElementData(Z) == nullptr) {
235 InitialiseOnFly(Z);
236 if (fCrossSection->GetElementData(Z) == nullptr) { return cs; }
237 }
238
239 // 7: rows in the parameterization file; 5: number of parameters
240 G4int idx = fNShells[Z] * 7 - 5;
241
242 energy = std::max(energy, (*(fParamHigh[Z]))[idx - 1]);
243
244 G4double x1 = 1.0 / energy;
245 G4double x2 = x1 * x1;
246 G4double x3 = x2 * x1;
247
248 // high energy parameterisation
249 if (energy >= (*(fParamHigh[Z]))[0]) {
250 G4double x4 = x2 * x2;
251 G4double x5 = x4 * x1;
252
253 cs = x1 *
254 ((*(fParamHigh[Z]))[idx] + x1 * (*(fParamHigh[Z]))[idx + 1]
255 + x2 * (*(fParamHigh[Z]))[idx + 2] + x3 * (*(fParamHigh[Z]))[idx + 3]
256 + x4 * (*(fParamHigh[Z]))[idx + 4] + x5 * (*(fParamHigh[Z]))[idx + 5]);
257 }
258 // low energy parameterisation
259 else if (energy >= (*(fParamLow[Z]))[0]) {
260 G4double x4 = x2 * x2;
261 G4double x5 = x4 * x1; // this variable usage can probably be optimized
262 cs = x1 *
263 ((*(fParamLow[Z]))[idx] + x1 * (*(fParamLow[Z]))[idx + 1]
264 + x2 * (*(fParamLow[Z]))[idx + 2] + x3 * (*(fParamLow[Z]))[idx + 3]
265 + x4 * (*(fParamLow[Z]))[idx + 4] + x5 * (*(fParamLow[Z]))[idx + 5]);
266 }
267 // Tabulated values above k-shell ionization energy
268 else if (energy >= (*(fParamHigh[Z]))[1]) {
269 cs = x3 * (fCrossSection->GetElementData(Z))->Value(energy);
270 }
271 // Tabulated values below k-shell ionization energy
272 else {
273 cs = x3 * (fCrossSectionLE->GetElementData(Z))->Value(energy);
274 }
275 if (verboseLevel > 1) {
276 G4cout << "G4LivermorePhotoElectricModel: E(keV)= " << energy / keV << " Z= " << Z
277 << " cross(barn)= " << cs / barn << G4endl;
278 }
279 return cs;
280}
G4double
double G4double
Definition G4Types.hh:83
G4INCL::KinematicsUtils::energy
G4double energy(const ThreeVector &p, const G4double m)
Definition G4INCLKinematicsUtils.cc:185
G4lrint
int G4lrint(double ad)
Definition templates.hh:134

◆ CrossSectionPerVolume()

G4double G4LivermorePhotoElectricModel::CrossSectionPerVolume ( const G4Material * material,
const G4ParticleDefinition * p,
G4double energy,
G4double cutEnergy = 0.0,
G4double maxEnergy = DBL_MAX )
overridevirtual

Reimplemented from G4VEmModel.

Definition at line 188 of file G4LivermorePhotoElectricModel.cc.

192{
193 fCurrSection = 0.0;
194 if (nullptr != fWater &&
195 (material == fWater || material->GetBaseMaterial() == fWater)) {
196 if (energy <= fWaterEnergyLimit) {
197 fWater->GetSandiaTable()->GetSandiaCofWater(energy, fSandiaCof);
198
199 G4double energy2 = energy * energy;
200 G4double energy3 = energy * energy2;
201 G4double energy4 = energy2 * energy2;
202
203 fCurrSection = material->GetDensity()
204 * (fSandiaCof[0] / energy + fSandiaCof[1] / energy2 +
205 fSandiaCof[2] / energy3 + fSandiaCof[3] / energy4);
206 }
207 }
208 if (0.0 == fCurrSection) {
209 fCurrSection = G4VEmModel::CrossSectionPerVolume(material, p, energy);
210 }
211 return fCurrSection;
212}
G4Material::GetDensity
G4double GetDensity() const
Definition G4Material.hh:174
G4Material::GetBaseMaterial
const G4Material * GetBaseMaterial() const
Definition G4Material.hh:218
G4VEmModel::CrossSectionPerVolume
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition G4VEmModel.cc:178

◆ GetBindingEnergy()

G4double G4LivermorePhotoElectricModel::GetBindingEnergy ( G4int Z,
G4int shell )

Definition at line 730 of file G4LivermorePhotoElectricModel.cc.

731{
732 if (Z < 1 || Z >= ZMAXPE) {
733 return -1;
734 } // If Z is out of the supported return 0
735
736 InitialiseOnFly(Z);
737 if (fCrossSection->GetElementData(Z) == nullptr ||
738 shell < 0 || shell >= fNShellsUsed[Z]) {
739 return -1;
740 }
741
742 if (Z > 2) {
743 return fCrossSection->GetComponentDataByIndex(Z, shell)->Energy(0);
744 }
745 else {
746 return fCrossSection->GetElementData(Z)->Energy(0);
747 }
748}

◆ Initialise()

void G4LivermorePhotoElectricModel::Initialise ( const G4ParticleDefinition * ,
const G4DataVector &  )
overridevirtual

Implements G4VEmModel.

Definition at line 131 of file G4LivermorePhotoElectricModel.cc.

133{
134 if (verboseLevel > 1) {
135 G4cout << "Calling G4LivermorePhotoElectricModel::Initialise() " << G4endl;
136 }
137
138 // initialise static tables only once
139 std::call_once(applyOnce, [this]() { isInitializer = true; });
140
141 if (isInitializer) {
142 G4AutoLock l(&livPhotoeffMutex);
143 if (fWater == nullptr) {
144 fWater = G4Material::GetMaterial("G4_WATER", false);
145 if (fWater == nullptr) {
146 fWater = G4Material::GetMaterial("Water", false);
147 }
148 if (fWater != nullptr) {
149 fWaterEnergyLimit = 13.6 * CLHEP::eV;
150 }
151 }
152
153 const G4ElementTable* elemTable = G4Element::GetElementTable();
154 std::size_t numElems = (*elemTable).size();
155 for (std::size_t ie = 0; ie < numElems; ++ie) {
156 const G4Element* elem = (*elemTable)[ie];
157 G4int Z = elem->GetZasInt();
158 if (Z < ZMAXPE) {
159 if (fCrossSection->GetElementData(Z) == nullptr) {
160 ReadData(Z);
161 }
162 }
163 }
164 l.unlock();
165 }
166
167 if (verboseLevel > 1) {
168 G4cout << "Loaded cross section files for new LivermorePhotoElectric model" << G4endl;
169 }
170 if (nullptr == fParticleChange) {
171 fParticleChange = GetParticleChangeForGamma();
172 fAtomDeexcitation = G4LossTableManager::Instance()->AtomDeexcitation();
173 }
174
175 fDeexcitationActive = false;
176 if (nullptr != fAtomDeexcitation) {
177 fDeexcitationActive = fAtomDeexcitation->IsFluoActive();
178 }
179
180 if (verboseLevel > 1) {
181 G4cout << "LivermorePhotoElectric model is initialized " << G4endl;
182 }
183}
G4AutoLock
G4TemplateAutoLock< G4Mutex > G4AutoLock
Definition G4AutoLock.hh:582
G4ElementTable
std::vector< G4Element * > G4ElementTable
Definition G4ElementTable.hh:33
elem
#define elem(i, j)
G4Element::GetElementTable
static const G4ElementTable * GetElementTable()
Definition G4Element.cc:401
G4LossTableManager::Instance
static G4LossTableManager * Instance()
Definition G4LossTableManager.cc:90
G4LossTableManager::AtomDeexcitation
G4VAtomDeexcitation * AtomDeexcitation()
Definition G4LossTableManager.hh:396
G4Material::GetMaterial
static G4Material * GetMaterial(const G4String &name, G4bool warning=true)
Definition G4Material.cc:663
G4VEmModel::GetParticleChangeForGamma
G4ParticleChangeForGamma * GetParticleChangeForGamma()
Definition G4VEmModel.cc:121

◆ operator=()

G4LivermorePhotoElectricModel & G4LivermorePhotoElectricModel::operator= ( const G4LivermorePhotoElectricModel & right)
delete

◆ SampleSecondaries()

void G4LivermorePhotoElectricModel::SampleSecondaries ( std::vector< G4DynamicParticle * > * fvect,
const G4MaterialCutsCouple * couple,
const G4DynamicParticle * aDynamicGamma,
G4double tmin,
G4double maxEnergy )
overridevirtual

Implements G4VEmModel.

Definition at line 284 of file G4LivermorePhotoElectricModel.cc.

288{
289 G4double gammaEnergy = aDynamicGamma->GetKineticEnergy();
290 if (verboseLevel > 3) {
291 G4cout << "G4LivermorePhotoElectricModel::SampleSecondaries() Egamma(keV)= "
292 << gammaEnergy / keV << G4endl;
293 }
294
295 // kill incident photon
296 fParticleChange->ProposeTrackStatus(fStopAndKill);
297 fParticleChange->SetProposedKineticEnergy(0.);
298
299 // low-energy photo-effect in water - full absorption
300 const G4Material* material = couple->GetMaterial();
301 if (fWater && (material == fWater || material->GetBaseMaterial() == fWater)) {
302 if (gammaEnergy <= fWaterEnergyLimit) {
303 fParticleChange->ProposeLocalEnergyDeposit(gammaEnergy);
304 return;
305 }
306 }
307
308 // Returns the normalized direction of the momentum
309 G4ThreeVector photonDirection = aDynamicGamma->GetMomentumDirection();
310
311 // Select randomly one element in the current material
312 const G4Element* elm = SelectRandomAtom(material, theGamma, gammaEnergy);
313 G4int Z = elm->GetZasInt();
314
315 // Select the ionised shell in the current atom according to shell
316 // cross sections
317
318 // If element was not initialised gamma should be absorbed
319 if (Z >= ZMAXPE || fCrossSection->GetElementData(Z) == nullptr) {
320 fParticleChange->ProposeLocalEnergyDeposit(gammaEnergy);
321 return;
322 }
323
324 // SAMPLING OF THE SHELL INDEX
325 std::size_t shellIdx = 0;
326 std::size_t nn = fNShellsUsed[Z];
327 if (nn > 1) {
328 if (gammaEnergy >= (*(fParamHigh[Z]))[0]) {
329 G4double x1 = 1.0 / gammaEnergy;
330 G4double x2 = x1 * x1;
331 G4double x3 = x2 * x1;
332 G4double x4 = x3 * x1;
333 G4double x5 = x4 * x1;
334 std::size_t idx = nn * 7 - 5;
335 // when do sampling common factors are not taken into account
336 // so cross section is not real
337
338 G4double rand = G4UniformRand();
339 G4double cs0 = rand *
340 ((*(fParamHigh[Z]))[idx] + x1 * (*(fParamHigh[Z]))[idx + 1]
341 + x2 * (*(fParamHigh[Z]))[idx + 2] + x3 * (*(fParamHigh[Z]))[idx + 3]
342 + x4 * (*(fParamHigh[Z]))[idx + 4] + x5 * (*(fParamHigh[Z]))[idx + 5]);
343
344 for (shellIdx = 0; shellIdx < nn; ++shellIdx) {
345 idx = shellIdx * 7 + 2;
346 if (gammaEnergy > (*(fParamHigh[Z]))[idx - 1]) {
347 G4double cs = (*(fParamHigh[Z]))[idx] + x1 * (*(fParamHigh[Z]))[idx + 1]
348 + x2 * (*(fParamHigh[Z]))[idx + 2] + x3 * (*(fParamHigh[Z]))[idx + 3]
349 + x4 * (*(fParamHigh[Z]))[idx + 4] + x5 * (*(fParamHigh[Z]))[idx + 5];
350
351 if (cs >= cs0) {
352 break;
353 }
354 }
355 }
356 if (shellIdx >= nn) {
357 shellIdx = nn - 1;
358 }
359 }
360 else if (gammaEnergy >= (*(fParamLow[Z]))[0]) {
361 G4double x1 = 1.0 / gammaEnergy;
362 G4double x2 = x1 * x1;
363 G4double x3 = x2 * x1;
364 G4double x4 = x3 * x1;
365 G4double x5 = x4 * x1;
366 std::size_t idx = nn * 7 - 5;
367 // when do sampling common factors are not taken into account
368 // so cross section is not real
369 G4double cs0 = G4UniformRand() *
370 ((*(fParamLow[Z]))[idx] + x1 * (*(fParamLow[Z]))[idx + 1]
371 + x2 * (*(fParamLow[Z]))[idx + 2] + x3 * (*(fParamLow[Z]))[idx + 3]
372 + x4 * (*(fParamLow[Z]))[idx + 4] + x5 * (*(fParamLow[Z]))[idx + 5]);
373 for (shellIdx = 0; shellIdx < nn; ++shellIdx) {
374 idx = shellIdx * 7 + 2;
375 if (gammaEnergy > (*(fParamLow[Z]))[idx - 1]) {
376 G4double cs = (*(fParamLow[Z]))[idx] + x1 * (*(fParamLow[Z]))[idx + 1]
377 + x2 * (*(fParamLow[Z]))[idx + 2] + x3 * (*(fParamLow[Z]))[idx + 3]
378 + x4 * (*(fParamLow[Z]))[idx + 4] + x5 * (*(fParamLow[Z]))[idx + 5];
379 if (cs >= cs0) {
380 break;
381 }
382 }
383 }
384 if (shellIdx >= nn) {
385 shellIdx = nn - 1;
386 }
387 }
388 else {
389 // when do sampling common factors are not taken into account
390 // so cross section is not real
391 G4double cs = G4UniformRand();
392
393 if (gammaEnergy >= (*(fParamHigh[Z]))[1]) {
394 // above K-shell binding energy
395 cs *= fCrossSection->GetElementData(Z)->Value(gammaEnergy);
396 }
397 else {
398 // below K-shell binding energy
399 cs *= fCrossSectionLE->GetElementData(Z)->Value(gammaEnergy);
400 }
401
402 for (G4int j = 0; j < (G4int)nn; ++j) {
403 shellIdx = (std::size_t)fCrossSection->GetComponentID(Z, j);
404 if (gammaEnergy > (*(fParamLow[Z]))[7 * shellIdx + 1]) {
405 cs -= fCrossSection->GetValueForComponent(Z, j, gammaEnergy);
406 }
407 if (cs <= 0.0 || j + 1 == (G4int)nn) {
408 break;
409 }
410 }
411 }
412 }
413 // END: SAMPLING OF THE SHELL
414
415 G4double bindingEnergy = (*(fParamHigh[Z]))[shellIdx * 7 + 1];
416 const G4AtomicShell* shell = nullptr;
417
418 // no de-excitation from the last shell
419 if (fDeexcitationActive && shellIdx + 1 < nn) {
420 auto as = G4AtomicShellEnumerator(shellIdx);
421 shell = fAtomDeexcitation->GetAtomicShell(Z, as);
422 }
423
424 // If binding energy of the selected shell is larger than photon energy
425 // do not generate secondaries
426 if (gammaEnergy < bindingEnergy) {
427 fParticleChange->ProposeLocalEnergyDeposit(gammaEnergy);
428 return;
429 }
430
431 // Primary outcoming electron
432 G4double eKineticEnergy = gammaEnergy - bindingEnergy;
433 G4double edep = bindingEnergy;
434
435 // Calculate direction of the photoelectron
436 G4ThreeVector electronDirection = GetAngularDistribution()->SampleDirection(
437 aDynamicGamma, eKineticEnergy, (G4int)shellIdx, couple->GetMaterial());
438
439 // The electron is created
440 auto electron =
441 new G4DynamicParticle(theElectron, electronDirection, eKineticEnergy);
442 fvect->push_back(electron);
443
444 // Sample deexcitation
445 if (nullptr != shell) {
446 G4int index = couple->GetIndex();
447 if (fAtomDeexcitation->CheckDeexcitationActiveRegion(index)) {
448 std::size_t nbefore = fvect->size();
449
450 fAtomDeexcitation->GenerateParticles(fvect, shell, Z, index);
451 std::size_t nafter = fvect->size();
452 if (nafter > nbefore) {
453 G4double esec = 0.0;
454 for (std::size_t j = nbefore; j < nafter; ++j) {
455 G4double e = ((*fvect)[j])->GetKineticEnergy();
456 if (esec + e > edep) {
457 // correct energy in order to have energy balance
458 e = edep - esec;
459 ((*fvect)[j])->SetKineticEnergy(e);
460 esec += e;
461 // delete the rest of secondaries (should not happens)
462 for (std::size_t jj = nafter - 1; jj > j; --jj) {
463 delete (*fvect)[jj];
464 fvect->pop_back();
465 }
466 break;
467 }
468 esec += e;
469 }
470 edep -= esec;
471 }
472 }
473 }
474 // energy balance - excitation energy left
475 if (edep > 0.0) {
476 fParticleChange->ProposeLocalEnergyDeposit(edep);
477 }
478}
G4ThreeVector
CLHEP::Hep3Vector G4ThreeVector
Definition G4ThreeVector.hh:36
fStopAndKill
@ fStopAndKill
Definition G4TrackStatus.hh:46
G4UniformRand
#define G4UniformRand()
Definition Randomize.hh:52
G4DynamicParticle::GetMomentumDirection
const G4ThreeVector & GetMomentumDirection() const
G4DynamicParticle::GetKineticEnergy
G4double GetKineticEnergy() const
G4Element::GetZasInt
G4int GetZasInt() const
Definition G4Element.hh:120
G4MaterialCutsCouple::GetMaterial
const G4Material * GetMaterial() const
Definition G4MaterialCutsCouple.hh:166
G4VEmAngularDistribution::SampleDirection
virtual G4ThreeVector & SampleDirection(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, const G4Material *)=0
G4VEmModel::GetAngularDistribution
G4VEmAngularDistribution * GetAngularDistribution()
Definition G4VEmModel.hh:602
G4VEmModel::SelectRandomAtom
const G4Element * SelectRandomAtom(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition G4VEmModel.hh:563
G4InuclSpecialFunctions::bindingEnergy
G4double bindingEnergy(G4int A, G4int Z)
Definition bindingEnergy.cc:36

◆ SetLimitNumberOfShells()

void G4LivermorePhotoElectricModel::SetLimitNumberOfShells ( G4int n)
inline

Definition at line 70 of file G4LivermorePhotoElectricModel.hh.

70{ nShellLimit = n; };

Member Data Documentation

◆ fParticleChange

G4ParticleChangeForGamma* G4LivermorePhotoElectricModel::fParticleChange = nullptr
protected

Definition at line 78 of file G4LivermorePhotoElectricModel.hh.

Referenced by Initialise(), and SampleSecondaries().

The documentation for this class was generated from the following files: