Geant4  10.02.p01
G4RDAtomicDeexcitation.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 //
27 // $Id: G4RDAtomicDeexcitation.cc,v 1.11
28 // GEANT4 tag $Name: geant4-09-01-ref-00 $
29 //
30 // Authors: Elena Guardincerri (Elena.Guardincerri@ge.infn.it)
31 // Alfonso Mantero (Alfonso.Mantero@ge.infn.it)
32 //
33 // History:
34 // -----------
35 //
36 // 16 Sept 2001 First committed to cvs
37 // 12 Sep 2003 Bug in auger production fixed
38 //
39 // -------------------------------------------------------------------
40 
42 #include "Randomize.hh"
43 #include "G4PhysicalConstants.hh"
44 #include "G4SystemOfUnits.hh"
45 #include "G4Gamma.hh"
46 #include "G4Electron.hh"
48 #include "G4RDFluoTransition.hh"
49 
51  minGammaEnergy(100.*eV),
52  minElectronEnergy(100.*eV),
53  fAuger(false)
54 {}
55 
57 {}
58 
59 std::vector<G4DynamicParticle*>* G4RDAtomicDeexcitation::GenerateParticles(G4int Z,G4int givenShellId)
60 {
61 
62  std::vector<G4DynamicParticle*>* vectorOfParticles;
63 
64  vectorOfParticles = new std::vector<G4DynamicParticle*>;
65  G4DynamicParticle* aParticle;
66  G4int provShellId = 0;
67  G4int counter = 0;
68 
69  // The aim of this loop is to generate more than one fluorecence photon
70  // from the same ionizing event
71  do
72  {
73  if (counter == 0)
74  // First call to GenerateParticles(...):
75  // givenShellId is given by the process
76  {
77  provShellId = SelectTypeOfTransition(Z, givenShellId);
78  //std::cout << "AtomicDeexcitation::Generate counter 0 - provShellId = "
79  //<< provShellId << std::endl;
80 
81  if ( provShellId >0)
82  {
83  aParticle = GenerateFluorescence(Z,givenShellId,provShellId);
84  //std::cout << "AtomicDeexcitation::Generate Fluo counter 0 " << std::endl;
85  }
86  else if ( provShellId == -1)
87  {
88  aParticle = GenerateAuger(Z, givenShellId);
89  //std::cout << "AtomicDeexcitation::Generate Auger counter 0 " << std::endl;
90  }
91  else
92  {
93  G4Exception("G4RDAtomicDeexcitation::GenerateParticles()",
94  "InvalidSetup", FatalException,
95  "Starting shell uncorrect: check it!");
96  }
97  }
98  else
99  // Following calls to GenerateParticles(...):
100  // newShellId is given by GenerateFluorescence(...)
101  {
102  provShellId = SelectTypeOfTransition(Z,newShellId);
103  //std::cout << "AtomicDeexcitation::Generate counter 0 - provShellId = "
104  //<< provShellId << ", new ShellId = "<< newShellId
105  //<< std::endl;
106 
107 
108  if (provShellId >0)
109  {
110  aParticle = GenerateFluorescence(Z,newShellId,provShellId);
111  //std::cout << "AtomicDeexcitation::Generate Fluo " << std::endl;
112  }
113  else if ( provShellId == -1)
114  {
115  aParticle = GenerateAuger(Z, newShellId);
116  //std::cout << "AtomicDeexcitation::Generate Auger " << std::endl;
117  }
118  else
119  {
120  G4Exception("G4RDAtomicDeexcitation::GenerateParticles()",
121  "InvalidSetup", FatalException,
122  "Starting shell uncorrect: check it!");
123  }
124  }
125  counter++;
126  if (aParticle != 0) {vectorOfParticles->push_back(aParticle);}
127  else {provShellId = -2;}
128  }
129 
130  // Look this in a particular way: only one auger emitted! //
131  while (provShellId > -2);
132 
133  return vectorOfParticles;
134 }
135 
137 {
138  if (shellId <=0 )
139  {
140  G4Exception("G4RDAtomicDeexcitation::SelectTypeOfTransition()",
141  "InvalidCondition", FatalException,
142  "Zero or negative shellId!");
143  }
144 
145  const G4RDAtomicTransitionManager* transitionManager =
147  G4int provShellId = -1;
148  G4int shellNum = 0;
149  G4int maxNumOfShells = transitionManager->NumberOfReachableShells(Z);
150 
151  //std::cout << "AtomicDeexcitation::SelectType - NumberOfReachableShells = "
152  //<< maxNumOfShells<< std::endl;
153 
154  const G4RDFluoTransition* refShell = transitionManager->ReachableShell(Z,maxNumOfShells-1);
155 
156  // This loop gives shellNum the value of the index of shellId
157  // in the vector storing the list of the shells reachable through
158  // a radiative transition
159  if ( shellId <= refShell->FinalShellId())
160  {
161  while (shellId != transitionManager->ReachableShell(Z,shellNum)->FinalShellId())
162  {
163  if(shellNum ==maxNumOfShells-1)
164  {
165  break;
166  }
167  shellNum++;
168  }
169  G4int transProb = 0; //AM change 29/6/07 was 1
170 
171  G4double partialProb = G4UniformRand();
172  G4double partSum = 0;
173  const G4RDFluoTransition* aShell = transitionManager->ReachableShell(Z,shellNum);
174  G4int trSize = (aShell->TransitionProbabilities()).size();
175 
176  // Loop over the shells wich can provide an electron for a
177  // radiative transition towards shellId:
178  // in every loop the partial sum of the first transProb shells
179  // is calculated and compared with a random number [0,1].
180  // If the partial sum is greater, the shell whose index is transProb
181  // is chosen as the starting shell for a radiative transition
182  // and its identity is returned
183  // Else, terminateded the loop, -1 is returned
184  while(transProb < trSize){
185 
186  partSum += aShell->TransitionProbability(transProb);
187 
188  if(partialProb <= partSum)
189  {
190  provShellId = aShell->OriginatingShellId(transProb);
191  break;
192  }
193  transProb++;
194  }
195 
196  // here provShellId is the right one or is -1.
197  // if -1, the control is passed to the Auger generation part of the package
198  }
199 
200 
201 
202  else
203  {
204 
205  provShellId = -1;
206 
207  }
208  return provShellId;
209 }
210 
212  G4int shellId,
213  G4int provShellId )
214 {
215 
216 
218  // G4int provenienceShell = provShellId;
219 
220  //isotropic angular distribution for the outcoming photon
221  G4double newcosTh = 1.-2.*G4UniformRand();
222  G4double newsinTh = std::sqrt(1.-newcosTh*newcosTh);
223  G4double newPhi = twopi*G4UniformRand();
224 
225  G4double xDir = newsinTh*std::sin(newPhi);
226  G4double yDir = newsinTh*std::cos(newPhi);
227  G4double zDir = newcosTh;
228 
229  G4ThreeVector newGammaDirection(xDir,yDir,zDir);
230 
231  G4int shellNum = 0;
232  G4int maxNumOfShells = transitionManager->NumberOfReachableShells(Z);
233 
234  // find the index of the shell named shellId
235  while (shellId != transitionManager->
236  ReachableShell(Z,shellNum)->FinalShellId())
237  {
238  if(shellNum == maxNumOfShells-1)
239  {
240  break;
241  }
242  shellNum++;
243  }
244  // number of shell from wich an electron can reach shellId
245  size_t transitionSize = transitionManager->
246  ReachableShell(Z,shellNum)->OriginatingShellIds().size();
247 
248  size_t index = 0;
249 
250  // find the index of the shell named provShellId in the vector
251  // storing the shells from which shellId can be reached
252  while (provShellId != transitionManager->
253  ReachableShell(Z,shellNum)->OriginatingShellId(index))
254  {
255  if(index == transitionSize-1)
256  {
257  break;
258  }
259  index++;
260  }
261  // energy of the gamma leaving provShellId for shellId
262  G4double transitionEnergy = transitionManager->
263  ReachableShell(Z,shellNum)->TransitionEnergy(index);
264 
265  // This is the shell where the new vacancy is: it is the same
266  // shell where the electron came from
267  newShellId = transitionManager->
268  ReachableShell(Z,shellNum)->OriginatingShellId(index);
269 
270 
272  newGammaDirection,
273  transitionEnergy);
274  return newPart;
275 }
276 
278 {
279  if(!fAuger) return 0;
280 
281 
282  const G4RDAtomicTransitionManager* transitionManager =
284 
285 
286 
287  if (shellId <=0 )
288  {
289  G4Exception("G4RDAtomicDeexcitation::GenerateAuger()",
290  "InvalidCondition", FatalException,
291  "Zero or negative shellId!");
292  }
293 
294  // G4int provShellId = -1;
295  G4int maxNumOfShells = transitionManager->NumberOfReachableAugerShells(Z);
296 
297  const G4RDAugerTransition* refAugerTransition =
298  transitionManager->ReachableAugerShell(Z,maxNumOfShells-1);
299 
300 
301  // This loop gives to shellNum the value of the index of shellId
302  // in the vector storing the list of the vacancies in the variuos shells
303  // that can originate a NON-radiative transition
304 
305  // ---- MGP ---- Next line commented out to remove compilation warning
306  // G4int p = refAugerTransition->FinalShellId();
307 
308  G4int shellNum = 0;
309 
310 
311  if ( shellId <= refAugerTransition->FinalShellId() )
312  //"FinalShellId" is final from the point of view of the elctron who makes the transition,
313  // being the Id of the shell in which there is a vacancy
314  {
315  G4int pippo = transitionManager->ReachableAugerShell(Z,shellNum)->FinalShellId();
316  if (shellId != pippo ) {
317  do {
318  shellNum++;
319  if(shellNum == maxNumOfShells)
320  {
321 // G4cout << "G4RDAtomicDeexcitation warning: No Auger transition found" << G4endl;
322 // G4cout << "Absorbed enrgy deposited locally" << G4endl;
323  return 0;
324 // // G4Exception("G4RDAtomicDeexcitation: No Auger transition found");
325  }
326  }
327  while (shellId != (transitionManager->ReachableAugerShell(Z,shellNum)->FinalShellId()) ) ;
328  }
329  /* {
330 
331  if(shellNum == maxNumOfShells-1)
332  {
333  G4Exception("G4RDAtomicDeexcitation: No Auger tramsition found");
334  }
335  shellNum++;
336  }*/
337 
338 
339 
340 
341  // Now we have that shellnum is the shellIndex of the shell named ShellId
342 
343  // G4cout << " the index of the shell is: "<<shellNum<<G4endl;
344 
345  // But we have now to select two shells: one for the transition,
346  // and another for the auger emission.
347 
348  G4int transitionLoopShellIndex = 0;
349  G4double partSum = 0;
350  const G4RDAugerTransition* anAugerTransition =
351  transitionManager->ReachableAugerShell(Z,shellNum);
352 
353  // G4cout << " corresponding to the ID: "<< anAugerTransition->FinalShellId() << G4endl;
354 
355 
356  G4int transitionSize =
357  (anAugerTransition->TransitionOriginatingShellIds())->size();
358  while (transitionLoopShellIndex < transitionSize) {
359 
360  std::vector<G4int>::const_iterator pos =
361  anAugerTransition->TransitionOriginatingShellIds()->begin();
362 
363  G4int transitionLoopShellId = *(pos+transitionLoopShellIndex);
364  G4int numberOfPossibleAuger =
365  (anAugerTransition->AugerTransitionProbabilities(transitionLoopShellId))->size();
366  G4int augerIndex = 0;
367  // G4int partSum2 = 0;
368 
369 
370  if (augerIndex < numberOfPossibleAuger) {
371 
372  do
373  {
374  G4double thisProb = anAugerTransition->AugerTransitionProbability(augerIndex,
375  transitionLoopShellId);
376  partSum += thisProb;
377  augerIndex++;
378 
379  } while (augerIndex < numberOfPossibleAuger);
380  }
381  transitionLoopShellIndex++;
382  }
383 
384 
385 
386  // Now we have the entire probability of an auger transition for the vacancy
387  // located in shellNum (index of shellId)
388 
389  // AM *********************** F I X E D **************************** AM
390  // Here we duplicate the previous loop, this time looking to the sum of the probabilities
391  // to be under the random number shoot by G4 UniformRdandom. This could have been done in the
392  // previuos loop, while integrating the probabilities. There is a bug that will be fixed
393  // 5 minutes from now: a line:
394  // G4int numberOfPossibleAuger = (anAugerTransition->
395  // AugerTransitionProbabilities(transitionLoopShellId))->size();
396  // to be inserted.
397  // AM *********************** F I X E D **************************** AM
398 
399  // Remains to get the same result with a single loop.
400 
401  // AM *********************** F I X E D **************************** AM
402  // Another Bug: in EADL Auger Transition are normalized to all the transitions deriving from
403  // a vacancy in one shell, but not all of these are present in data tables. So if a transition
404  // doesn't occur in the main one a local energy deposition must occur, instead of (like now)
405  // generating the last transition present in EADL data.
406  // AM *********************** F I X E D **************************** AM
407 
408 
409  G4double totalVacancyAugerProbability = partSum;
410 
411 
412  //And now we start to select the right auger transition and emission
413  G4int transitionRandomShellIndex = 0;
414  G4int transitionRandomShellId = 1;
415  G4int augerIndex = 0;
416  partSum = 0;
417  G4double partialProb = G4UniformRand();
418  // G4int augerOriginatingShellId = 0;
419 
420  G4int numberOfPossibleAuger =
421  (anAugerTransition->AugerTransitionProbabilities(transitionRandomShellId))->size();
422  G4bool foundFlag = false;
423 
424  while (transitionRandomShellIndex < transitionSize) {
425 
426  std::vector<G4int>::const_iterator pos =
427  anAugerTransition->TransitionOriginatingShellIds()->begin();
428 
429  transitionRandomShellId = *(pos+transitionRandomShellIndex);
430 
431  augerIndex = 0;
432  numberOfPossibleAuger = (anAugerTransition->
433  AugerTransitionProbabilities(transitionRandomShellId))->size();
434 
435  while (augerIndex < numberOfPossibleAuger) {
436  G4double thisProb =anAugerTransition->AugerTransitionProbability(augerIndex,
437  transitionRandomShellId);
438 
439  partSum += thisProb;
440 
441  if (partSum >= (partialProb*totalVacancyAugerProbability) ) { // was /
442  foundFlag = true;
443  break;
444  }
445  augerIndex++;
446  }
447  if (partSum >= (partialProb*totalVacancyAugerProbability) ) {break;} // was /
448  transitionRandomShellIndex++;
449  }
450 
451  // Now we have the index of the shell from wich comes the auger electron (augerIndex),
452  // and the id of the shell, from which the transition e- come (transitionRandomShellid)
453  // If no Transition has been found, 0 is returned.
454 
455  if (!foundFlag) {return 0;}
456 
457  // Isotropic angular distribution for the outcoming e-
458  G4double newcosTh = 1.-2.*G4UniformRand();
459  G4double newsinTh = std::sqrt(1.-newcosTh*newcosTh);
460  G4double newPhi = twopi*G4UniformRand();
461 
462  G4double xDir = newsinTh*std::sin(newPhi);
463  G4double yDir = newsinTh*std::cos(newPhi);
464  G4double zDir = newcosTh;
465 
466  G4ThreeVector newElectronDirection(xDir,yDir,zDir);
467 
468  // energy of the auger electron emitted
469 
470 
471  G4double transitionEnergy = anAugerTransition->AugerTransitionEnergy(augerIndex, transitionRandomShellId);
472  /*
473  G4cout << "AUger TransitionId " << anAugerTransition->FinalShellId() << G4endl;
474  G4cout << "augerIndex: " << augerIndex << G4endl;
475  G4cout << "transitionShellId: " << transitionRandomShellId << G4endl;
476  */
477 
478  // This is the shell where the new vacancy is: it is the same
479  // shell where the electron came from
480  newShellId = transitionRandomShellId;
481 
482 
484  newElectronDirection,
485  transitionEnergy);
486  return newPart;
487 
488  }
489  else
490  {
491  //G4Exception("G4RDAtomicDeexcitation: no auger transition found");
492  return 0;
493  }
494 
495 }
496 
498 {
499  minGammaEnergy = cut;
500 }
501 
503 {
504  minElectronEnergy = cut;
505 }
506 
508 {
509  fAuger = val;
510 }
511 
512 
513 
514 
515 
516 
517 
const G4RDFluoTransition * ReachableShell(G4int Z, size_t shellIndex) const
G4double TransitionProbability(G4int index) const
G4double AugerTransitionProbability(G4int index, G4int startShellId) const
G4DynamicParticle * GenerateFluorescence(G4int Z, G4int shellId, G4int provShellId)
G4double AugerTransitionEnergy(G4int index, G4int startShellId) const
void SetCutForSecondaryPhotons(G4double cut)
std::vector< G4DynamicParticle * > * GenerateParticles(G4int Z, G4int shellId)
const G4DataVector & TransitionProbabilities() const
void ActivateAugerElectronProduction(G4bool val)
CLHEP::Hep3Vector G4ThreeVector
const std::vector< G4int > * TransitionOriginatingShellIds() const
G4int NumberOfReachableAugerShells(G4int Z) const
int G4int
Definition: G4Types.hh:78
static G4RDAtomicTransitionManager * Instance()
#define G4UniformRand()
Definition: Randomize.hh:97
bool G4bool
Definition: G4Types.hh:79
static const double twopi
Definition: G4SIunits.hh:75
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
const G4RDAugerTransition * ReachableAugerShell(G4int Z, G4int shellIndex) const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4DynamicParticle * GenerateAuger(G4int Z, G4int shellId)
static const double eV
Definition: G4SIunits.hh:212
G4int OriginatingShellId(G4int index) const
void SetCutForAugerElectrons(G4double cut)
G4int SelectTypeOfTransition(G4int Z, G4int shellId)
static G4Electron * Electron()
Definition: G4Electron.cc:94
double G4double
Definition: G4Types.hh:76
G4int FinalShellId() const
static const G4double pos
const G4DataVector * AugerTransitionProbabilities(G4int startShellId) const