Geant4 11.1.1
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4VAtomDeexcitation Class Referenceabstract

#include <G4VAtomDeexcitation.hh>

+ Inheritance diagram for G4VAtomDeexcitation:

Public Member Functions

 G4VAtomDeexcitation (const G4String &modname="Deexcitation")
 
virtual ~G4VAtomDeexcitation ()
 
void InitialiseAtomicDeexcitation ()
 
virtual void InitialiseForNewRun ()=0
 
virtual void InitialiseForExtraAtom (G4int Z)=0
 
void SetDeexcitationActiveRegion (const G4String &rname, G4bool valDeexcitation, G4bool valAuger, G4bool valPIXE)
 
void SetFluo (G4bool)
 
G4bool IsFluoActive () const
 
void SetAuger (G4bool)
 
G4bool IsAugerActive () const
 
void SetAugerCascade (G4bool)
 
G4bool IsAugerCascadeActive () const
 
void SetPIXE (G4bool)
 
G4bool IsPIXEActive () const
 
const G4StringGetName () const
 
const std::vector< G4bool > & GetListOfActiveAtoms () const
 
void SetVerboseLevel (G4int)
 
G4int GetVerboseLevel () const
 
G4bool CheckDeexcitationActiveRegion (G4int coupleIndex)
 
G4bool CheckAugerActiveRegion (G4int coupleIndex)
 
virtual const G4AtomicShellGetAtomicShell (G4int Z, G4AtomicShellEnumerator shell)=0
 
void GenerateParticles (std::vector< G4DynamicParticle * > *secVect, const G4AtomicShell *, G4int Z, G4int coupleIndex)
 
virtual void GenerateParticles (std::vector< G4DynamicParticle * > *secVect, const G4AtomicShell *, G4int Z, G4double gammaCut, G4double eCut)=0
 
virtual G4double GetShellIonisationCrossSectionPerAtom (const G4ParticleDefinition *, G4int Z, G4AtomicShellEnumerator shell, G4double kinE, const G4Material *mat=nullptr)=0
 
virtual G4double ComputeShellIonisationCrossSectionPerAtom (const G4ParticleDefinition *, G4int Z, G4AtomicShellEnumerator shell, G4double kinE, const G4Material *mat=nullptr)=0
 
void AlongStepDeexcitation (std::vector< G4Track * > &tracks, const G4Step &step, G4double &eLoss, G4int coupleIndex)
 
 G4VAtomDeexcitation (G4VAtomDeexcitation &)=delete
 
G4VAtomDeexcitationoperator= (const G4VAtomDeexcitation &right)=delete
 

Detailed Description

Definition at line 64 of file G4VAtomDeexcitation.hh.

Constructor & Destructor Documentation

◆ G4VAtomDeexcitation() [1/2]

G4VAtomDeexcitation::G4VAtomDeexcitation ( const G4String modname = "Deexcitation")
explicit

Definition at line 72 of file G4VAtomDeexcitation.cc.

73 : name(modname)
74{
75 vdyn.reserve(5);
76 theCoupleTable = nullptr;
77 gamma = G4Gamma::Gamma();
78}
static G4Gamma * Gamma()
Definition: G4Gamma.cc:85

◆ ~G4VAtomDeexcitation()

G4VAtomDeexcitation::~G4VAtomDeexcitation ( )
virtualdefault

◆ G4VAtomDeexcitation() [2/2]

G4VAtomDeexcitation::G4VAtomDeexcitation ( G4VAtomDeexcitation )
delete

Member Function Documentation

◆ AlongStepDeexcitation()

void G4VAtomDeexcitation::AlongStepDeexcitation ( std::vector< G4Track * > &  tracks,
const G4Step step,
G4double eLoss,
G4int  coupleIndex 
)

Definition at line 260 of file G4VAtomDeexcitation.cc.

264{
265 G4double truelength = step.GetStepLength();
266 if(!flagPIXE && !activePIXEMedia[coupleIndex]) { return; }
267 if(eLossMax <= 0.0 || truelength <= 0.0) { return; }
268
269 // step parameters
270 const G4StepPoint* preStep = step.GetPreStepPoint();
271 const G4ThreeVector prePos = preStep->GetPosition();
272 const G4ThreeVector delta = step.GetPostStepPoint()->GetPosition() - prePos;
273 const G4double preTime = preStep->GetGlobalTime();
274 const G4double dt = step.GetPostStepPoint()->GetGlobalTime() - preTime;
275
276 // particle parameters
277 const G4Track* track = step.GetTrack();
278 const G4ParticleDefinition* part = track->GetDefinition();
279 G4double ekin = preStep->GetKineticEnergy();
280
281 // media parameters
282 G4double gCut = (*theCoupleTable->GetEnergyCutsVector(0))[coupleIndex];
283 if(ignoreCuts) { gCut = 0.0; }
284 G4double eCut = DBL_MAX;
285 if(CheckAugerActiveRegion(coupleIndex)) {
286 eCut = (*theCoupleTable->GetEnergyCutsVector(1))[coupleIndex];
287 if(ignoreCuts) { eCut = 0.0; }
288 }
289
290 //G4cout<<"!Sample PIXE gCut(MeV)= "<<gCut<<" eCut(MeV)= "<<eCut
291 // <<" Ekin(MeV)= " << ekin/MeV << G4endl;
292
293 const G4Material* material = preStep->GetMaterial();
294 const G4ElementVector* theElementVector = material->GetElementVector();
295 const G4double* theAtomNumDensityVector =
296 material->GetVecNbOfAtomsPerVolume();
297 const std::size_t nelm = material->GetNumberOfElements();
298
299 // loop over deexcitations
300 for(std::size_t i=0; i<nelm; ++i) {
301 G4int Z = (*theElementVector)[i]->GetZasInt();
302 if(activeZ[Z] && Z < 93) {
303 G4int nshells =
304 std::min(9,(*theElementVector)[i]->GetNbOfAtomicShells());
305 G4double rho = truelength*theAtomNumDensityVector[i];
306 //G4cout<<" Z "<< Z <<" is active x(mm)= " << truelength/mm << G4endl;
307 for(G4int ii=0; ii<nshells; ++ii) {
308 auto as = (G4AtomicShellEnumerator)(ii);
309 const G4AtomicShell* shell = GetAtomicShell(Z, as);
310 const G4double bindingEnergy = shell->BindingEnergy();
311
312 if(gCut > bindingEnergy) { break; }
313
314 if(eLossMax > bindingEnergy) {
315 G4double sig = rho*
316 GetShellIonisationCrossSectionPerAtom(part, Z, as, ekin, material);
317
318 // mfp is mean free path in units of step size
319 if(sig > 0.0) {
320 G4double mfp = 1.0/sig;
321 G4double stot = 0.0;
322 //G4cout << " Shell " << ii << " mfp(mm)= " << mfp/mm << G4endl;
323 // sample ionisation points
324 do {
325 stot -= mfp*G4Log(G4UniformRand());
326 if( stot > 1.0 || eLossMax < bindingEnergy) { break; }
327 // sample deexcitation
328 vdyn.clear();
329 GenerateParticles(&vdyn, shell, Z, gCut, eCut);
330 std::size_t nsec = vdyn.size();
331 if(nsec > 0) {
332 G4ThreeVector r = prePos + stot*delta;
333 G4double time = preTime + stot*dt;
334 for(std::size_t j=0; j<nsec; ++j) {
335 G4DynamicParticle* dp = vdyn[j];
336 G4double e = dp->GetKineticEnergy();
337
338 // save new secondary if there is enough energy
339 if(eLossMax >= e) {
340 eLossMax -= e;
341 G4Track* t = new G4Track(dp, time, r);
342
343 // defined secondary type
344 if(dp->GetDefinition() == gamma) {
346 } else {
348 }
349 tracks.push_back(t);
350 } else {
351 delete dp;
352 }
353 }
354 }
355 // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
356 } while (stot < 1.0);
357 }
358 }
359 }
360 }
361 }
362 return;
363}
G4AtomicShellEnumerator
std::vector< const G4Element * > G4ElementVector
G4double G4Log(G4double x)
Definition: G4Log.hh:227
double G4double
Definition: G4Types.hh:83
int G4int
Definition: G4Types.hh:85
const G4int Z[17]
#define G4UniformRand()
Definition: Randomize.hh:52
G4double BindingEnergy() const
G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:185
size_t GetNumberOfElements() const
Definition: G4Material.hh:181
const G4double * GetVecNbOfAtomsPerVolume() const
Definition: G4Material.hh:201
const std::vector< G4double > * GetEnergyCutsVector(std::size_t pcIdx) const
G4double GetGlobalTime() const
G4Material * GetMaterial() const
const G4ThreeVector & GetPosition() const
G4double GetKineticEnergy() const
G4Track * GetTrack() const
G4StepPoint * GetPreStepPoint() const
G4double GetStepLength() const
G4StepPoint * GetPostStepPoint() const
G4ParticleDefinition * GetDefinition() const
void SetCreatorModelID(const G4int id)
virtual const G4AtomicShell * GetAtomicShell(G4int Z, G4AtomicShellEnumerator shell)=0
G4bool CheckAugerActiveRegion(G4int coupleIndex)
void GenerateParticles(std::vector< G4DynamicParticle * > *secVect, const G4AtomicShell *, G4int Z, G4int coupleIndex)
virtual G4double GetShellIonisationCrossSectionPerAtom(const G4ParticleDefinition *, G4int Z, G4AtomicShellEnumerator shell, G4double kinE, const G4Material *mat=nullptr)=0
G4double bindingEnergy(G4int A, G4int Z)
#define DBL_MAX
Definition: templates.hh:62

Referenced by G4VEnergyLossProcess::AlongStepDoIt().

◆ CheckAugerActiveRegion()

G4bool G4VAtomDeexcitation::CheckAugerActiveRegion ( G4int  coupleIndex)
inline

Definition at line 266 of file G4VAtomDeexcitation.hh.

267{
268 return (idx < nCouples) ? activeAugerMedia[idx] : false;
269}

Referenced by AlongStepDeexcitation(), and GenerateParticles().

◆ CheckDeexcitationActiveRegion()

◆ ComputeShellIonisationCrossSectionPerAtom()

virtual G4double G4VAtomDeexcitation::ComputeShellIonisationCrossSectionPerAtom ( const G4ParticleDefinition ,
G4int  Z,
G4AtomicShellEnumerator  shell,
G4double  kinE,
const G4Material mat = nullptr 
)
pure virtual

◆ GenerateParticles() [1/2]

virtual void G4VAtomDeexcitation::GenerateParticles ( std::vector< G4DynamicParticle * > *  secVect,
const G4AtomicShell ,
G4int  Z,
G4double  gammaCut,
G4double  eCut 
)
pure virtual

Implemented in G4UAtomicDeexcitation.

◆ GenerateParticles() [2/2]

void G4VAtomDeexcitation::GenerateParticles ( std::vector< G4DynamicParticle * > *  secVect,
const G4AtomicShell as,
G4int  Z,
G4int  coupleIndex 
)

Definition at line 234 of file G4VAtomDeexcitation.cc.

237{
238 G4double gCut = DBL_MAX;
239 if(ignoreCuts) {
240 gCut = 0.0;
241 } else if (nullptr != theCoupleTable) {
242 gCut = (*(theCoupleTable->GetEnergyCutsVector(0)))[idx];
243 }
244 if(gCut < as->BindingEnergy()) {
245 G4double eCut = DBL_MAX;
246 if(CheckAugerActiveRegion(idx)) {
247 if(ignoreCuts) {
248 eCut = 0.0;
249 } else if (nullptr != theCoupleTable) {
250 eCut = (*(theCoupleTable->GetEnergyCutsVector(1)))[idx];
251 }
252 }
253 GenerateParticles(v, as, Z, gCut, eCut);
254 }
255}

Referenced by AlongStepDeexcitation(), G4ECDecay::DecayIt(), G4ITDecay::DecayIt(), GenerateParticles(), G4DNABornIonisationModel1::SampleSecondaries(), G4DNABornIonisationModel2::SampleSecondaries(), G4DNACPA100IonisationModel::SampleSecondaries(), G4DNAEmfietzoglouIonisationModel::SampleSecondaries(), G4DNARelativisticIonisationModel::SampleSecondaries(), G4DNARuddIonisationExtendedModel::SampleSecondaries(), G4DNARuddIonisationModel::SampleSecondaries(), G4DNARPWBAIonisationModel::SampleSecondaries(), G4LivermoreComptonModel::SampleSecondaries(), G4LivermorePhotoElectricModel::SampleSecondaries(), G4LivermorePolarizedComptonModel::SampleSecondaries(), G4LowEPComptonModel::SampleSecondaries(), G4LowEPPolarizedComptonModel::SampleSecondaries(), G4MicroElecInelasticModel::SampleSecondaries(), G4MicroElecInelasticModel_new::SampleSecondaries(), G4PenelopeComptonModel::SampleSecondaries(), G4PenelopeIonisationModel::SampleSecondaries(), G4PenelopePhotoElectricModel::SampleSecondaries(), G4KleinNishinaModel::SampleSecondaries(), and G4PEEffectFluoModel::SampleSecondaries().

◆ GetAtomicShell()

◆ GetListOfActiveAtoms()

const std::vector< G4bool > & G4VAtomDeexcitation::GetListOfActiveAtoms ( ) const
inline

Definition at line 244 of file G4VAtomDeexcitation.hh.

245{
246 return activeZ;
247}

◆ GetName()

const G4String & G4VAtomDeexcitation::GetName ( ) const
inline

Definition at line 238 of file G4VAtomDeexcitation.hh.

239{
240 return name;
241}

Referenced by SetDeexcitationActiveRegion().

◆ GetShellIonisationCrossSectionPerAtom()

virtual G4double G4VAtomDeexcitation::GetShellIonisationCrossSectionPerAtom ( const G4ParticleDefinition ,
G4int  Z,
G4AtomicShellEnumerator  shell,
G4double  kinE,
const G4Material mat = nullptr 
)
pure virtual

◆ GetVerboseLevel()

G4int G4VAtomDeexcitation::GetVerboseLevel ( ) const
inline

Definition at line 254 of file G4VAtomDeexcitation.hh.

255{
256 return verbose;
257}

◆ InitialiseAtomicDeexcitation()

void G4VAtomDeexcitation::InitialiseAtomicDeexcitation ( )

Definition at line 86 of file G4VAtomDeexcitation.cc.

87{
89 theParameters->DefineRegParamForDeex(this);
90
91 // Define list of couples
93 nCouples = (G4int)theCoupleTable->GetTableSize();
94
95 // needed for unit tests
96 std::size_t nn = std::max(nCouples, 1);
97 if(activeDeexcitationMedia.size() != nn) {
98 activeDeexcitationMedia.resize(nn, false);
99 activeAugerMedia.resize(nn, false);
100 activePIXEMedia.resize(nn, false);
101 }
102 if(activeZ.size() != 93) { activeZ.resize(93, false); }
103
104 // initialisation of flags and options
105 // normally there is no locksed flags
106 if(!isActiveLocked) { isActive = theParameters->Fluo(); }
107 if(!isAugerLocked) { flagAuger = theParameters->Auger(); }
108 if(!isPIXELocked) { flagPIXE = theParameters->Pixe(); }
109 ignoreCuts = theParameters->DeexcitationIgnoreCut();
110
111 // Define list of regions
112 std::size_t nRegions = deRegions.size();
113 // check if deexcitation is active for the given run
114 if(!isActive && 0 == nRegions) { return; }
115
116 // if no active regions add a world
117 if(0 == nRegions) {
118 SetDeexcitationActiveRegion("World",isActive,flagAuger,flagPIXE);
119 nRegions = deRegions.size();
120 }
121
122 if(0 < verbose) {
123 G4cout << G4endl;
124 G4cout << "### === Deexcitation model " << name
125 << " is activated for " << nRegions;
126 if(1 == nRegions) { G4cout << " region:" << G4endl; }
127 else { G4cout << " regions:" << G4endl;}
128 }
129
130 // Identify active media
131 const G4RegionStore* regionStore = G4RegionStore::GetInstance();
132 for(std::size_t j=0; j<nRegions; ++j) {
133 const G4Region* reg = regionStore->GetRegion(activeRegions[j], false);
134 if(nullptr != reg && 0 < nCouples) {
135 const G4ProductionCuts* rpcuts = reg->GetProductionCuts();
136 if(0 < verbose) {
137 G4cout << " " << activeRegions[j]
138 << " " << deRegions[j] << " " << AugerRegions[j]
139 << " " << PIXERegions[j] << G4endl;
140 }
141 for(G4int i=0; i<nCouples; ++i) {
142 const G4MaterialCutsCouple* couple =
143 theCoupleTable->GetMaterialCutsCouple(i);
144 if (couple->GetProductionCuts() == rpcuts) {
145 activeDeexcitationMedia[i] = deRegions[j];
146 activeAugerMedia[i] = AugerRegions[j];
147 activePIXEMedia[i] = PIXERegions[j];
148 }
149 }
150 }
151 }
152 std::size_t nelm = G4Element::GetNumberOfElements();
153 //G4cout << nelm << G4endl;
154 for(std::size_t k=0; k<nelm; ++k) {
155 G4int Z = (*(G4Element::GetElementTable()))[k]->GetZasInt();
156 if(Z > 5 && Z < 93) {
157 activeZ[Z] = true;
158 //G4cout << "!!! Active de-excitation Z= " << Z << G4endl;
159 }
160 }
161
162 // Initialise derived class
164
165 if(0 < verbose && flagAuger) {
166 G4cout << "### === Auger flag: " << flagAuger
167 << G4endl;
168 }
169 if(0 < verbose) {
170 G4cout << "### === Ignore cuts flag: " << ignoreCuts
171 << G4endl;
172 }
173 if(0 < verbose && flagPIXE) {
174 G4cout << "### === PIXE model for hadrons: "
175 << theParameters->PIXECrossSectionModel()
176 << G4endl;
177 G4cout << "### === PIXE model for e+-: "
178 << theParameters->PIXEElectronCrossSectionModel()
179 << G4endl;
180 }
181}
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
static G4ElementTable * GetElementTable()
Definition: G4Element.cc:403
static size_t GetNumberOfElements()
Definition: G4Element.cc:410
static G4EmParameters * Instance()
const G4String & PIXECrossSectionModel()
const G4String & PIXEElectronCrossSectionModel()
G4bool Fluo() const
G4bool Pixe() const
void DefineRegParamForDeex(G4VAtomDeexcitation *) const
G4bool DeexcitationIgnoreCut() const
G4bool Auger() const
G4ProductionCuts * GetProductionCuts() const
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
std::size_t GetTableSize() const
static G4ProductionCutsTable * GetProductionCutsTable()
static G4RegionStore * GetInstance()
G4Region * GetRegion(const G4String &name, G4bool verbose=true) const
G4ProductionCuts * GetProductionCuts() const
void SetDeexcitationActiveRegion(const G4String &rname, G4bool valDeexcitation, G4bool valAuger, G4bool valPIXE)
virtual void InitialiseForNewRun()=0

Referenced by LBE::ConstructGeneral(), and G4LossTableManager::ResetParameters().

◆ InitialiseForExtraAtom()

virtual void G4VAtomDeexcitation::InitialiseForExtraAtom ( G4int  Z)
pure virtual

Implemented in G4UAtomicDeexcitation.

◆ InitialiseForNewRun()

virtual void G4VAtomDeexcitation::InitialiseForNewRun ( )
pure virtual

Implemented in G4UAtomicDeexcitation.

Referenced by InitialiseAtomicDeexcitation().

◆ IsAugerActive()

G4bool G4VAtomDeexcitation::IsAugerActive ( ) const
inline

Definition at line 213 of file G4VAtomDeexcitation.hh.

214{
215 return flagAuger;
216}

◆ IsAugerCascadeActive()

G4bool G4VAtomDeexcitation::IsAugerCascadeActive ( ) const
inline

Definition at line 223 of file G4VAtomDeexcitation.hh.

224{
225 return flagAuger;
226}

Referenced by G4UAtomicDeexcitation::GenerateParticles().

◆ IsFluoActive()

G4bool G4VAtomDeexcitation::IsFluoActive ( ) const
inline

◆ IsPIXEActive()

G4bool G4VAtomDeexcitation::IsPIXEActive ( ) const
inline

◆ operator=()

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

◆ SetAuger()

void G4VAtomDeexcitation::SetAuger ( G4bool  val)
inline

Definition at line 208 of file G4VAtomDeexcitation.hh.

209{
210 if(!isAugerLocked) { flagAuger = val; isAugerLocked = true; }
211}

Referenced by SetAugerCascade().

◆ SetAugerCascade()

void G4VAtomDeexcitation::SetAugerCascade ( G4bool  val)
inline

Definition at line 218 of file G4VAtomDeexcitation.hh.

219{
220 SetAuger(val);
221}

◆ SetDeexcitationActiveRegion()

void G4VAtomDeexcitation::SetDeexcitationActiveRegion ( const G4String rname,
G4bool  valDeexcitation,
G4bool  valAuger,
G4bool  valPIXE 
)

Definition at line 186 of file G4VAtomDeexcitation.cc.

190{
191 // no PIXE in parallel world
192 if(rname == "DefaultRegionForParallelWorld") { return; }
193
194 G4String ss = rname;
195 /*
196 G4cout << "### G4VAtomDeexcitation::SetDeexcitationActiveRegion " << ss
197 << " " << valDeexcitation << " " << valAuger
198 << " " << valPIXE << G4endl;
199 */
200 if(ss == "world" || ss == "World" || ss == "WORLD") {
201 ss = "DefaultRegionForTheWorld";
202 }
203 std::size_t n = deRegions.size();
204 for(std::size_t i=0; i<n; ++i) {
205
206 // Region already exist
207 if(ss == activeRegions[i]) {
208 deRegions[i] = valDeexcitation;
209 AugerRegions[i] = valAuger;
210 PIXERegions[i] = valPIXE;
211 return;
212 }
213 }
214 // New region
215 activeRegions.push_back(ss);
216 deRegions.push_back(valDeexcitation);
217 AugerRegions.push_back(valAuger);
218 PIXERegions.push_back(valPIXE);
219
220 // if de-excitation defined for the world volume
221 // it should be active for all G4Regions
222 if(ss == "DefaultRegionForTheWorld") {
224 std::size_t nn = regions->size();
225 for(std::size_t i=0; i<nn; ++i) {
226 if(ss == (*regions)[i]->GetName()) { continue; }
227 SetDeexcitationActiveRegion((*regions)[i]->GetName(), valDeexcitation,
228 valAuger, valPIXE);
229
230 }
231 }
232}
const G4String & GetName() const

Referenced by G4EmLowEParameters::DefineRegParamForDeex(), InitialiseAtomicDeexcitation(), and SetDeexcitationActiveRegion().

◆ SetFluo()

void G4VAtomDeexcitation::SetFluo ( G4bool  val)
inline

Definition at line 198 of file G4VAtomDeexcitation.hh.

199{
200 if(!isActiveLocked) { isActive = val; isActiveLocked = true; }
201}

◆ SetPIXE()

void G4VAtomDeexcitation::SetPIXE ( G4bool  val)
inline

Definition at line 228 of file G4VAtomDeexcitation.hh.

229{
230 if(!isPIXELocked) { flagPIXE = val; isPIXELocked = true; }
231}

◆ SetVerboseLevel()

void G4VAtomDeexcitation::SetVerboseLevel ( G4int  val)
inline

Definition at line 249 of file G4VAtomDeexcitation.hh.

250{
251 verbose = val;
252}

Referenced by G4LossTableManager::ResetParameters().


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