Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4UAtomicDeexcitation Class Reference

#include <G4UAtomicDeexcitation.hh>

Inheritance diagram for G4UAtomicDeexcitation:
Collaboration diagram for G4UAtomicDeexcitation:

Public Member Functions

 G4UAtomicDeexcitation ()
 
virtual ~G4UAtomicDeexcitation ()
 
virtual void InitialiseForNewRun ()
 
virtual void InitialiseForExtraAtom (G4int Z)
 
void SetCutForSecondaryPhotons (G4double cut)
 
void SetCutForAugerElectrons (G4double cut)
 
virtual const G4AtomicShellGetAtomicShell (G4int Z, G4AtomicShellEnumerator shell)
 
virtual void GenerateParticles (std::vector< G4DynamicParticle * > *secVect, const G4AtomicShell *, G4int Z, G4double gammaCut, G4double eCut)
 
virtual G4double GetShellIonisationCrossSectionPerAtom (const G4ParticleDefinition *, G4int Z, G4AtomicShellEnumerator shell, G4double kinE, const G4Material *mat=0)
 
virtual G4double ComputeShellIonisationCrossSectionPerAtom (const G4ParticleDefinition *, G4int Z, G4AtomicShellEnumerator shell, G4double kinE, const G4Material *mat=0)
 
- Public Member Functions inherited from G4VAtomDeexcitation
 G4VAtomDeexcitation (const G4String &modname="Deexcitation")
 
virtual ~G4VAtomDeexcitation ()
 
void InitialiseAtomicDeexcitation ()
 
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)
 
void GenerateParticles (std::vector< G4DynamicParticle * > *secVect, const G4AtomicShell *, G4int Z, G4int coupleIndex)
 
void AlongStepDeexcitation (std::vector< G4Track * > &tracks, const G4Step &step, G4double &eLoss, G4int coupleIndex)
 

Detailed Description

Definition at line 61 of file G4UAtomicDeexcitation.hh.

Constructor & Destructor Documentation

G4UAtomicDeexcitation::G4UAtomicDeexcitation ( )

Definition at line 76 of file G4UAtomicDeexcitation.cc.

76  :
77  G4VAtomDeexcitation("UAtomDeexcitation"),
78  minGammaEnergy(DBL_MAX),
79  minElectronEnergy(DBL_MAX)
80 {
81  anaPIXEshellCS = nullptr;
82  PIXEshellCS = nullptr;
83  ePIXEshellCS = nullptr;
85  theElectron = G4Electron::Electron();
86  thePositron = G4Positron::Positron();
87  transitionManager = G4AtomicTransitionManager::Instance();
88 }
static G4LossTableManager * Instance()
G4EmCorrections * EmCorrections()
static G4Positron * Positron()
Definition: G4Positron.cc:94
static G4Electron * Electron()
Definition: G4Electron.cc:94
static G4AtomicTransitionManager * Instance()
#define DBL_MAX
Definition: templates.hh:83
G4VAtomDeexcitation(const G4String &modname="Deexcitation")

Here is the call graph for this function:

G4UAtomicDeexcitation::~G4UAtomicDeexcitation ( )
virtual

Definition at line 90 of file G4UAtomicDeexcitation.cc.

91 {
92  delete anaPIXEshellCS;
93  delete PIXEshellCS;
94  delete ePIXEshellCS;
95 }

Member Function Documentation

G4double G4UAtomicDeexcitation::ComputeShellIonisationCrossSectionPerAtom ( const G4ParticleDefinition p,
G4int  Z,
G4AtomicShellEnumerator  shell,
G4double  kinE,
const G4Material mat = 0 
)
virtual

Implements G4VAtomDeexcitation.

Definition at line 404 of file G4UAtomicDeexcitation.cc.

410 {
411  return GetShellIonisationCrossSectionPerAtom(p,Z,shell,kinE,mat);
412 }
virtual G4double GetShellIonisationCrossSectionPerAtom(const G4ParticleDefinition *, G4int Z, G4AtomicShellEnumerator shell, G4double kinE, const G4Material *mat=0)

Here is the call graph for this function:

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

Implements G4VAtomDeexcitation.

Definition at line 176 of file G4UAtomicDeexcitation.cc.

182 {
183 
184  // Defined initial conditions
185  G4int givenShellId = atomicShell->ShellId();
186  //G4cout << "generating particles for vacancy in shellId: "
187  // << givenShellId << G4endl; // debug
188  minGammaEnergy = gammaCut;
189  minElectronEnergy = eCut;
190 
191  // generation secondaries
192  G4DynamicParticle* aParticle=0;
193  G4int provShellId = 0;
194 
195 //ORIGINAL METHOD BY ALFONSO MANTERO
196 if (!IsAugerCascadeActive())
197 {
198 //----------------------------
199  G4int counter = 0;
200 
201  // let's check that 5<Z<100
202 
203  if (Z>5 && Z<100) {
204 
205  // The aim of this loop is to generate more than one fluorecence photon
206  // from the same ionizing event
207  do
208  {
209  if (counter == 0)
210  // First call to GenerateParticles(...):
211  // givenShellId is given by the process
212  {
213  provShellId = SelectTypeOfTransition(Z, givenShellId);
214 
215  if ( provShellId >0)
216  {
217  aParticle = GenerateFluorescence(Z,givenShellId,provShellId);
218  //if (aParticle != 0) {
219  // G4cout << "****FLUO!_1**** "
220  // << aParticle->GetParticleDefinition()->GetParticleType()
221  // << " " << aParticle->GetKineticEnergy()/keV << G4endl ;}
222  }
223  else if ( provShellId == -1)
224  {
225  // G4cout << "Try to generate Auger 1" << G4endl;
226  aParticle = GenerateAuger(Z, givenShellId);
227  // if (aParticle != 0) { G4cout << "****AUGER!****" << G4endl;}
228  }
229  else
230  {
231  //G4Exception("G4UAtomicDeexcitation::GenerateParticles()",
232  // "de0002",JustWarning, "Energy deposited locally");
233  }
234  }
235  else
236  // Following calls to GenerateParticles(...):
237  // newShellId is given by GenerateFluorescence(...)
238  {
239  provShellId = SelectTypeOfTransition(Z,newShellId);
240  if (provShellId >0)
241  {
242  aParticle = GenerateFluorescence(Z,newShellId,provShellId);
243  //if (aParticle != 0) { G4cout << "****FLUO!_2****" << aParticle->GetParticleDefinition()->GetParticleType() << " " << aParticle->GetKineticEnergy()/keV << G4endl;} //debug
244  }
245  else if ( provShellId == -1)
246  {
247  // G4cout << "Try to generate Auger 2" << G4endl; //debug
248  aParticle = GenerateAuger(Z, newShellId);
249  // if (aParticle != 0) { G4cout << "****AUGER!****" << G4endl;} //debug
250  }
251  else
252  {
253  //G4Exception("G4UAtomicDeexcitation::GenerateParticles()","de0002",JustWarning, "Energy deposited locally");
254  }
255  }
256  counter++;
257  if (aParticle != 0)
258  {
259  vectorOfParticles->push_back(aParticle);
260  //G4cout << "Deexcitation Occurred!" << G4endl; //debug
261  }
262  else {provShellId = -2;}
263  }
264  while (provShellId > -2);
265  }
266  else
267  {
268  //G4Exception("G4UAtomicDeexcitation::GenerateParticles()","de0001",JustWarning, "Energy deposited locally");
269  }
270 
271  //G4cout << "---------FATTO!---------" << G4endl; //debug
272 
273 } // Auger cascade is not active
274 
275 //END OF ORIGINAL METHOD BY ALFONSO MANTERO
276 //----------------------
277 
278 // NEW METHOD
279 // Auger cascade by Burkhant Suerfu on March 24 2015 (Bugzilla 1727)
280 
282 {
283 //----------------------
284 
285  vacancyArray.push_back(givenShellId);
286 
287  // let's check that 5<Z<100
288  if (Z<6 || Z>99){
289  //G4Exception("G4UAtomicDeexcitation::GenerateParticles()","de0001",JustWarning, "Energy deposited locally");
290  return;
291  }
292 
293  // as long as there is vacancy to be filled by either fluo or auger, stay in the loop.
294  while(!vacancyArray.empty()){
295 
296 // prepare to process the last element, and then delete it from the vector.
297  givenShellId = vacancyArray[0];
298  provShellId = SelectTypeOfTransition(Z,givenShellId);
299 
300  //G4cout<<"\n------ Atom Transition with Z: "<<Z<<"\tbetween current:"
301  // <<givenShellId<<" & target:"<<provShellId<<G4endl;
302  if(provShellId>0){
303  aParticle = GenerateFluorescence(Z,givenShellId,provShellId);
304 // if (aParticle != 0) {
305 // G4cout << "****FLUO!_1**** "
306 // << aParticle->GetParticleDefinition()->GetParticleType()
307 // << " " << aParticle->GetKineticEnergy()/keV << G4endl ;}
308  }
309  else if(provShellId == -1){
310  aParticle = GenerateAuger(Z, givenShellId);
311 // if (aParticle != 0) { G4cout << "****AUGER!****" <<
312 // aParticle->GetParticleDefinition()->GetParticleType()
313 // << " " << aParticle->GetKineticEnergy()/keV << G4endl ; }
314 // else G4cout<<G4endl;
315  }
316  //else
317  // G4Exception("G4UAtomicDeexcitation::GenerateParticles()","de0002",JustWarning, "Energy deposited locally");
318 
319 // if a particle is created, put it in the vector of new particles
320  if(aParticle!=0)
321  vectorOfParticles->push_back(aParticle);
322  else{;}
323 // one vacancy has been processed. Erase it.
324  vacancyArray.erase(vacancyArray.begin());
325  }
326 
327 
328 //----------------------
329 //End of Auger cascade by Burkhant Suerfu on March 24 2015 (Bugzilla 1727)
330 
331 } // Auger cascade is active
332 
333 //ENDSI
334 }
G4int ShellId() const
int G4int
Definition: G4Types.hh:78
G4bool IsAugerCascadeActive() const

Here is the call graph for this function:

const G4AtomicShell * G4UAtomicDeexcitation::GetAtomicShell ( G4int  Z,
G4AtomicShellEnumerator  shell 
)
virtual

Implements G4VAtomDeexcitation.

Definition at line 171 of file G4UAtomicDeexcitation.cc.

172 {
173  return transitionManager->Shell(Z, size_t(shell));
174 }
G4AtomicShell * Shell(G4int Z, size_t shellIndex) const

Here is the call graph for this function:

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

Implements G4VAtomDeexcitation.

Definition at line 337 of file G4UAtomicDeexcitation.cc.

343 {
344  // we must put a control on the shell that are passed:
345  // some shells should not pass (line "0" or "2")
346  //G4cout << pdef->GetParticleName() << " Z= " << Z << " Shell= " << shellEnum
347  // << " E= " << kineticEnergy << G4endl;
348 
349  // check atomic number
350  G4double xsec = 0.0;
351  if(Z > 93 || Z < 6 ) { return xsec; } //corrected by alf - Z<6 missing
352  G4int idx = G4int(shellEnum);
353  if(idx >= G4AtomicShells::GetNumberOfShells(Z)) { return xsec; }
354 
355  //G4cout << pdef->GetParticleName() << " Z= " << Z << " " << PIXEshellCS
356  // << " " << ePIXEshellCS << G4endl;
357  //
358  if(pdef == theElectron || pdef == thePositron) {
359  xsec = ePIXEshellCS->CrossSection(Z,shellEnum,kineticEnergy,0.0,mat);
360  return xsec;
361  }
362 
363  G4double mass = pdef->GetPDGMass();
364  G4double escaled = kineticEnergy;
365  G4double q2 = 0.0;
366 
367  // scaling to protons
368  if ((pdef->GetParticleName() != "proton" && pdef->GetParticleName() != "alpha" ) )
369  {
370  mass = proton_mass_c2;
371  escaled = kineticEnergy*mass/(pdef->GetPDGMass());
372 
373  if(mat) {
374  q2 = emcorr->EffectiveChargeSquareRatio(pdef,mat,kineticEnergy);
375  } else {
376  G4double q = pdef->GetPDGCharge()/eplus;
377  q2 = q*q;
378  }
379  }
380 
381  if(PIXEshellCS) { xsec = PIXEshellCS->CrossSection(Z,shellEnum,escaled,mass,mat); }
382  if(xsec < 1e-100) {
383 
384  xsec = anaPIXEshellCS->CrossSection(Z,shellEnum,escaled,mass,mat);
385 
386  }
387 
388  if (q2) {xsec *= q2;}
389 
390  return xsec;
391 }
G4double EffectiveChargeSquareRatio(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
static constexpr double eplus
Definition: G4SIunits.hh:199
float proton_mass_c2
Definition: hepunit.py:275
G4double GetPDGMass() const
double G4double
Definition: G4Types.hh:76
G4double GetPDGCharge() const
virtual G4double CrossSection(G4int Z, G4AtomicShellEnumerator shell, G4double incidentEnergy, G4double mass, const G4Material *mat)=0
static G4int GetNumberOfShells(G4int Z)

Here is the call graph for this function:

Here is the caller graph for this function:

void G4UAtomicDeexcitation::InitialiseForExtraAtom ( G4int  Z)
virtual

Implements G4VAtomDeexcitation.

Definition at line 167 of file G4UAtomicDeexcitation.cc.

168 {}
void G4UAtomicDeexcitation::InitialiseForNewRun ( )
virtual

Implements G4VAtomDeexcitation.

Definition at line 97 of file G4UAtomicDeexcitation.cc.

98 {
99  if(!IsFluoActive()) { return; }
100 
101  transitionManager->Initialise();
102 
103  if(!IsPIXEActive()) { return; }
104 
105  if(!anaPIXEshellCS) {
106  anaPIXEshellCS = new G4teoCrossSection("ECPSSR_Analytical");
107  }
108  G4cout << G4endl;
109  G4cout << "### === G4UAtomicDeexcitation::InitialiseForNewRun()" << G4endl;
110 
112  G4String namePIXExsModel = param->PIXECrossSectionModel();
113  G4String namePIXExsElectronModel = param->PIXEElectronCrossSectionModel();
114 
115  //G4cout << namePIXExsModel << " " << namePIXExsElectronModel << G4endl;
116 
117  // Check if old cross section for p/ion should be deleted
118  if(PIXEshellCS && namePIXExsModel != PIXEshellCS->GetName())
119  {
120  delete PIXEshellCS;
121  PIXEshellCS = nullptr;
122  }
123 
124  // Instantiate new proton/ion cross section
125  if(!PIXEshellCS) {
126  if (namePIXExsModel == "ECPSSR_FormFactor")
127  {
128  PIXEshellCS = new G4teoCrossSection(namePIXExsModel);
129  }
130  else if(namePIXExsModel == "Empirical")
131  {
132  PIXEshellCS = new G4empCrossSection(namePIXExsModel);
133  }
134  }
135  //G4cout << "PIXE is initialised" << G4endl;
136 
137  // Check if old cross section for e+- should be deleted
138  if(ePIXEshellCS && namePIXExsElectronModel != ePIXEshellCS->GetName())
139  {
140  delete ePIXEshellCS;
141  ePIXEshellCS = nullptr;
142  }
143 
144  // Instantiate new e+- cross section
145  if(!ePIXEshellCS)
146  {
147  if(namePIXExsElectronModel == "Empirical")
148  {
149  ePIXEshellCS = new G4empCrossSection("Empirical");
150  }
151  else if(namePIXExsElectronModel == "ECPSSR_Analytical")
152  {
153  ePIXEshellCS = new G4teoCrossSection("ECPSSR_Analytical");
154  }
155  else if (namePIXExsElectronModel == "Penelope")
156  {
157  ePIXEshellCS = new G4PenelopeIonisationCrossSection();
158  }
159  else
160  {
161  ePIXEshellCS = new G4LivermoreIonisationCrossSection();
162  }
163  }
164  //G4cout << "ePIXE is initialised" << G4endl;
165 }
G4bool IsFluoActive() const
G4bool IsPIXEActive() const
const G4String & GetName() const
G4GLOB_DLL std::ostream G4cout
const G4String & PIXECrossSectionModel()
static G4EmParameters * Instance()
#define G4endl
Definition: G4ios.hh:61
const G4String & PIXEElectronCrossSectionModel()

Here is the call graph for this function:

void G4UAtomicDeexcitation::SetCutForAugerElectrons ( G4double  cut)

Definition at line 398 of file G4UAtomicDeexcitation.cc.

399 {
400  minElectronEnergy = cut;
401 }
void G4UAtomicDeexcitation::SetCutForSecondaryPhotons ( G4double  cut)

Definition at line 393 of file G4UAtomicDeexcitation.cc.

394 {
395  minGammaEnergy = cut;
396 }

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