Geant4 9.6.0
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", const G4String &pixename="")
 
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 SetPIXE (G4bool)
 
G4bool IsPIXEActive () const
 
const G4StringGetName () const
 
void SetPIXECrossSectionModel (const G4String &)
 
const G4StringPIXECrossSectionModel () const
 
void SetPIXEElectronCrossSectionModel (const G4String &)
 
const G4StringPIXEElectronCrossSectionModel () 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=0)=0
 
virtual G4double ComputeShellIonisationCrossSectionPerAtom (const G4ParticleDefinition *, G4int Z, G4AtomicShellEnumerator shell, G4double kinE, const G4Material *mat=0)=0
 
void AlongStepDeexcitation (std::vector< G4Track * > &tracks, const G4Step &step, G4double &eLoss, G4int coupleIndex)
 

Detailed Description

Definition at line 63 of file G4VAtomDeexcitation.hh.

Constructor & Destructor Documentation

◆ G4VAtomDeexcitation()

G4VAtomDeexcitation::G4VAtomDeexcitation ( const G4String modname = "Deexcitation",
const G4String pixename = "" 
)

Definition at line 65 of file G4VAtomDeexcitation.cc.

67 : lowestKinEnergy(keV), verbose(1), name(modname), namePIXE(pname),
68 nameElectronPIXE(""), isActive(false), flagAuger(false), flagPIXE(false)
69{
70 vdyn.reserve(5);
71 theCoupleTable = 0;
72}

◆ ~G4VAtomDeexcitation()

G4VAtomDeexcitation::~G4VAtomDeexcitation ( )
virtual

Definition at line 76 of file G4VAtomDeexcitation.cc.

77{}

Member Function Documentation

◆ AlongStepDeexcitation()

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

Definition at line 198 of file G4VAtomDeexcitation.cc.

202{
203 G4double truelength = step.GetStepLength();
204 if(!flagPIXE && !activePIXEMedia[coupleIndex]) { return; }
205 if(eLossMax <= 0.0 || truelength <= 0.0) { return; }
206
207 // step parameters
208 const G4StepPoint* preStep = step.GetPreStepPoint();
209 G4ThreeVector prePos = preStep->GetPosition();
210 G4ThreeVector delta = step.GetPostStepPoint()->GetPosition() - prePos;
211 G4double preTime = preStep->GetGlobalTime();
212 G4double dt = step.GetPostStepPoint()->GetGlobalTime() - preTime;
213
214 // particle parameters
215 const G4Track* track = step.GetTrack();
216 const G4ParticleDefinition* part = track->GetDefinition();
217 G4double ekin = preStep->GetKineticEnergy();
218
219 // media parameters
220 G4double gCut = (*theCoupleTable->GetEnergyCutsVector(0))[coupleIndex];
221 G4double eCut = DBL_MAX;
222 if(CheckAugerActiveRegion(coupleIndex)) {
223 eCut = (*theCoupleTable->GetEnergyCutsVector(1))[coupleIndex];
224 }
225
226 //G4cout<<"!Sample PIXE gCut(MeV)= "<<gCut<<" eCut(MeV)= "<<eCut
227 // <<" Ekin(MeV)= " << ekin/MeV << G4endl;
228
229 const G4Material* material = preStep->GetMaterial();
230 const G4ElementVector* theElementVector = material->GetElementVector();
231 const G4double* theAtomNumDensityVector = material->GetVecNbOfAtomsPerVolume();
232 G4int nelm = material->GetNumberOfElements();
233
234 // loop over deexcitations
235 for(G4int i=0; i<nelm; ++i) {
236 G4int Z = G4lrint((*theElementVector)[i]->GetZ());
237 if(activeZ[Z] && Z < 93) {
238 G4int nshells = std::min(9,(*theElementVector)[i]->GetNbOfAtomicShells());
239 G4double rho = truelength*theAtomNumDensityVector[i];
240 //G4cout << " Z " << Z <<" is active x(mm)= " << truelength/mm << G4endl;
241 for(G4int ii=0; ii<nshells; ++ii) {
243 const G4AtomicShell* shell = GetAtomicShell(Z, as);
245
246 if(gCut > bindingEnergy) { break; }
247
248 if(eLossMax > bindingEnergy) {
249 G4double sig = rho*
250 GetShellIonisationCrossSectionPerAtom(part, Z, as, ekin, material);
251
252 // mfp is mean free path in units of step size
253 if(sig > 0.0) {
254 G4double mfp = 1.0/sig;
255 G4double stot = 0.0;
256 //G4cout << " Shell " << ii << " mfp(mm)= " << mfp/mm << G4endl;
257 // sample ionisation points
258 do {
259 stot -= mfp*std::log(G4UniformRand());
260 if( stot > 1.0 || eLossMax < bindingEnergy) { break; }
261 // sample deexcitation
262 vdyn.clear();
263 GenerateParticles(&vdyn, shell, Z, gCut, eCut);
264 G4int nsec = vdyn.size();
265 if(nsec > 0) {
266 G4ThreeVector r = prePos + stot*delta;
267 G4double time = preTime + stot*dt;
268 for(G4int j=0; j<nsec; ++j) {
269 G4DynamicParticle* dp = vdyn[j];
270 G4double e = dp->GetKineticEnergy();
271
272 // save new secondary if there is enough energy
273 if(eLossMax >= e) {
274 eLossMax -= e;
275 G4Track* t = new G4Track(dp, time, r);
276 tracks.push_back(t);
277 } else {
278 delete dp;
279 }
280 }
281 }
282 } while (stot < 1.0);
283 }
284 }
285 }
286 }
287 }
288 return;
289}
G4AtomicShellEnumerator
std::vector< G4Element * > G4ElementVector
double G4double
Definition: G4Types.hh:64
int G4int
Definition: G4Types.hh:66
#define G4UniformRand()
Definition: Randomize.hh:53
G4double BindingEnergy() const
G4double GetKineticEnergy() const
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:189
size_t GetNumberOfElements() const
Definition: G4Material.hh:185
const G4double * GetVecNbOfAtomsPerVolume() const
Definition: G4Material.hh:205
const std::vector< G4double > * GetEnergyCutsVector(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
virtual G4double GetShellIonisationCrossSectionPerAtom(const G4ParticleDefinition *, G4int Z, G4AtomicShellEnumerator shell, G4double kinE, const G4Material *mat=0)=0
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)
G4double bindingEnergy(G4int A, G4int Z)
int G4lrint(double ad)
Definition: templates.hh:163
#define DBL_MAX
Definition: templates.hh:83

Referenced by G4VEnergyLossProcess::AlongStepDoIt().

◆ CheckAugerActiveRegion()

G4bool G4VAtomDeexcitation::CheckAugerActiveRegion ( G4int  coupleIndex)
inline

Definition at line 277 of file G4VAtomDeexcitation.hh.

278{
279
280 return (flagAuger || activeAugerMedia[coupleIndex]);
281}

Referenced by AlongStepDeexcitation(), and GenerateParticles().

◆ CheckDeexcitationActiveRegion()

◆ ComputeShellIonisationCrossSectionPerAtom()

virtual G4double G4VAtomDeexcitation::ComputeShellIonisationCrossSectionPerAtom ( const G4ParticleDefinition ,
G4int  Z,
G4AtomicShellEnumerator  shell,
G4double  kinE,
const G4Material mat = 0 
)
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]

◆ GetAtomicShell()

◆ GetListOfActiveAtoms()

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

Definition at line 255 of file G4VAtomDeexcitation.hh.

256{
257 return activeZ;
258}

◆ GetName()

const G4String & G4VAtomDeexcitation::GetName ( ) const
inline

Definition at line 225 of file G4VAtomDeexcitation.hh.

226{
227 return name;
228}

◆ GetShellIonisationCrossSectionPerAtom()

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

◆ GetVerboseLevel()

G4int G4VAtomDeexcitation::GetVerboseLevel ( ) const
inline

Definition at line 265 of file G4VAtomDeexcitation.hh.

266{
267 return verbose;
268}

◆ InitialiseAtomicDeexcitation()

void G4VAtomDeexcitation::InitialiseAtomicDeexcitation ( )

Definition at line 81 of file G4VAtomDeexcitation.cc.

82{
83 // Define list of couples
85 size_t numOfCouples = theCoupleTable->GetTableSize();
86
87 // needed for unit tests
88 if(0 == numOfCouples) { numOfCouples = 1; }
89
90 activeDeexcitationMedia.resize(numOfCouples, false);
91 activeAugerMedia.resize(numOfCouples, false);
92 activePIXEMedia.resize(numOfCouples, false);
93 activeZ.resize(93, false);
94
95 // check if deexcitation is active for the given run
96 if( !isActive ) { return; }
97
98 // Define list of regions
99 size_t nRegions = deRegions.size();
100
101 if(0 == nRegions) {
102 SetDeexcitationActiveRegion("World",isActive,flagAuger,flagPIXE);
103 nRegions = 1;
104 }
105
106 if(0 < verbose) {
107 G4cout << G4endl;
108 G4cout << "### === Deexcitation model " << name
109 << " is activated for " << nRegions;
110 if(1 == nRegions) { G4cout << " region:" << G4endl; }
111 else { G4cout << " regions:" << G4endl;}
112 }
113
114 // Identify active media
116 for(size_t j=0; j<nRegions; ++j) {
117 const G4Region* reg = regionStore->GetRegion(activeRegions[j], false);
118 const G4ProductionCuts* rpcuts = reg->GetProductionCuts();
119 if(0 < verbose) {
120 G4cout << " " << activeRegions[j] << G4endl;
121 }
122
123 for(size_t i=0; i<numOfCouples; ++i) {
124 const G4MaterialCutsCouple* couple =
125 theCoupleTable->GetMaterialCutsCouple(i);
126 if (couple->GetProductionCuts() == rpcuts) {
127 activeDeexcitationMedia[i] = deRegions[j];
128 activeAugerMedia[i] = AugerRegions[j];
129 activePIXEMedia[i] = PIXERegions[j];
130 const G4Material* mat = couple->GetMaterial();
131 const G4ElementVector* theElementVector =
132 mat->GetElementVector();
133 G4int nelm = mat->GetNumberOfElements();
134 if(deRegions[j]) {
135 for(G4int k=0; k<nelm; ++k) {
136 G4int Z = G4lrint(((*theElementVector)[k])->GetZ());
137 if(Z > 5 && Z < 93) {
138 activeZ[Z] = true;
139 //G4cout << "!!! Active de-excitation Z= " << Z << G4endl;
140 }
141 }
142 }
143 }
144 }
145 }
146
147 // Initialise derived class
149
150 if(0 < verbose && flagPIXE) {
151 G4cout << "### === PIXE model for hadrons: " << namePIXE
152 << " " << IsPIXEActive()
153 << G4endl;
154 G4cout << "### === PIXE model for e+-: " << nameElectronPIXE
155 << " " << IsPIXEActive()
156 << G4endl;
157 }
158}
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cout
const G4Material * GetMaterial() const
G4ProductionCuts * GetProductionCuts() const
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) 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 G4RadioactiveDecay::BuildPhysicsTable(), and G4LossTableManager::BuildPhysicsTable().

◆ 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 209 of file G4VAtomDeexcitation.hh.

210{
211 return flagAuger;
212}

Referenced by G4UAtomicDeexcitation::GenerateParticles().

◆ IsFluoActive()

G4bool G4VAtomDeexcitation::IsFluoActive ( ) const
inline

◆ IsPIXEActive()

G4bool G4VAtomDeexcitation::IsPIXEActive ( ) const
inline

◆ PIXECrossSectionModel()

const G4String & G4VAtomDeexcitation::PIXECrossSectionModel ( ) const
inline

Definition at line 243 of file G4VAtomDeexcitation.hh.

244{
245 return namePIXE;
246}

Referenced by G4UAtomicDeexcitation::InitialiseForNewRun().

◆ PIXEElectronCrossSectionModel()

const G4String & G4VAtomDeexcitation::PIXEElectronCrossSectionModel ( ) const
inline

Definition at line 249 of file G4VAtomDeexcitation.hh.

250{
251 return nameElectronPIXE;
252}

Referenced by G4UAtomicDeexcitation::InitialiseForNewRun().

◆ SetAuger()

void G4VAtomDeexcitation::SetAuger ( G4bool  val)
inline

Definition at line 203 of file G4VAtomDeexcitation.hh.

204{
205 flagAuger = val;
206 if(val) { isActive = true; }
207}

Referenced by G4EmProcessOptions::SetAuger().

◆ SetDeexcitationActiveRegion()

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

Definition at line 163 of file G4VAtomDeexcitation.cc.

167{
168 G4String ss = rname;
169 //G4cout << "### G4VAtomDeexcitation::SetDeexcitationActiveRegion " << ss
170 // << " " << valDeexcitation << " " << valAuger
171 // << " " << valPIXE << G4endl;
172 if(ss == "world" || ss == "World" || ss == "WORLD") {
173 ss = "DefaultRegionForTheWorld";
174 }
175 size_t n = deRegions.size();
176 if(n > 0) {
177 for(size_t i=0; i<n; ++i) {
178
179 // Region already exist
180 if(ss == activeRegions[i]) {
181 deRegions[i] = valDeexcitation;
182 AugerRegions[i] = valAuger;
183 PIXERegions[i] = valPIXE;
184 return;
185 }
186 }
187 }
188 // New region
189 activeRegions.push_back(ss);
190 deRegions.push_back(valDeexcitation);
191 AugerRegions.push_back(valAuger);
192 PIXERegions.push_back(valPIXE);
193}

Referenced by InitialiseAtomicDeexcitation(), and G4EmProcessOptions::SetDeexcitationActiveRegion().

◆ SetFluo()

◆ SetPIXE()

void G4VAtomDeexcitation::SetPIXE ( G4bool  val)
inline

Definition at line 214 of file G4VAtomDeexcitation.hh.

215{
216 flagPIXE = val;
217 if(val) { isActive = true; }
218}

Referenced by G4EmProcessOptions::SetPIXE().

◆ SetPIXECrossSectionModel()

void G4VAtomDeexcitation::SetPIXECrossSectionModel ( const G4String n)
inline

Definition at line 231 of file G4VAtomDeexcitation.hh.

232{
233 namePIXE = n;
234}

Referenced by G4UAtomicDeexcitation::InitialiseForNewRun(), and G4EmProcessOptions::SetPIXECrossSectionModel().

◆ SetPIXEElectronCrossSectionModel()

void G4VAtomDeexcitation::SetPIXEElectronCrossSectionModel ( const G4String n)
inline

Definition at line 237 of file G4VAtomDeexcitation.hh.

238{
239 nameElectronPIXE = n;
240}

Referenced by G4UAtomicDeexcitation::InitialiseForNewRun(), and G4EmProcessOptions::SetPIXEElectronCrossSectionModel().

◆ SetVerboseLevel()

void G4VAtomDeexcitation::SetVerboseLevel ( G4int  val)
inline

Definition at line 260 of file G4VAtomDeexcitation.hh.

261{
262 verbose = val;
263}

Referenced by G4LossTableManager::SetVerbose().


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