source: trunk/source/processes/electromagnetic/lowenergy/src/G4AtomicDeexcitation.cc@ 1058

Last change on this file since 1058 was 1055, checked in by garnier, 17 years ago

maj sur la beta de geant 4.9.3
File size: 15.8 KB
RevLine 
[819]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: G4AtomicDeexcitation.cc,v 1.11
[1055]28// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
[819]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
41#include "G4AtomicDeexcitation.hh"
42#include "Randomize.hh"
43#include "G4Gamma.hh"
44#include "G4Electron.hh"
45#include "G4AtomicTransitionManager.hh"
46#include "G4FluoTransition.hh"
47
48G4AtomicDeexcitation::G4AtomicDeexcitation():
49 minGammaEnergy(100.*eV),
50 minElectronEnergy(100.*eV),
51 fAuger(false)
52{}
53
54G4AtomicDeexcitation::~G4AtomicDeexcitation()
55{}
56
57std::vector<G4DynamicParticle*>* G4AtomicDeexcitation::GenerateParticles(G4int Z,G4int givenShellId)
58{
59
60 std::vector<G4DynamicParticle*>* vectorOfParticles;
61
62 vectorOfParticles = new std::vector<G4DynamicParticle*>;
63 G4DynamicParticle* aParticle;
64 G4int provShellId = 0;
65 G4int counter = 0;
66
67 // The aim of this loop is to generate more than one fluorecence photon
68 // from the same ionizing event
69 do
70 {
71 if (counter == 0)
72 // First call to GenerateParticles(...):
73 // givenShellId is given by the process
74 {
75 provShellId = SelectTypeOfTransition(Z, givenShellId);
76
77 if ( provShellId >0)
78 {
79 aParticle = GenerateFluorescence(Z,givenShellId,provShellId);
80 }
81 else if ( provShellId == -1)
82 {
83 aParticle = GenerateAuger(Z, givenShellId);
84 }
85 else
86 {
87 G4Exception("G4AtomicDeexcitation: starting shell uncorrect: check it");
88 }
89 }
90 else
91 // Following calls to GenerateParticles(...):
92 // newShellId is given by GenerateFluorescence(...)
93 {
94 provShellId = SelectTypeOfTransition(Z,newShellId);
95 if (provShellId >0)
96 {
97 aParticle = GenerateFluorescence(Z,newShellId,provShellId);
98 }
99 else if ( provShellId == -1)
100 {
101 aParticle = GenerateAuger(Z, newShellId);
102 }
103 else
104 {
105 G4Exception("G4AtomicDeexcitation: starting shell uncorrect: check it");
106 }
107 }
108 counter++;
109 if (aParticle != 0) {vectorOfParticles->push_back(aParticle);}
110 else {provShellId = -2;}
111 }
112
[1055]113 // Look this in a particular way: only one auger emitted! // ????
[819]114 while (provShellId > -2);
115
116 return vectorOfParticles;
117}
118
119G4int G4AtomicDeexcitation::SelectTypeOfTransition(G4int Z, G4int shellId)
120{
121 if (shellId <=0 )
122 {G4Exception("G4AtomicDeexcitation: zero or negative shellId");}
123
124 G4bool fluoTransitionFoundFlag = false;
125
126 const G4AtomicTransitionManager* transitionManager =
127 G4AtomicTransitionManager::Instance();
128 G4int provShellId = -1;
129 G4int shellNum = 0;
130 G4int maxNumOfShells = transitionManager->NumberOfReachableShells(Z);
131
132 const G4FluoTransition* refShell = transitionManager->ReachableShell(Z,maxNumOfShells-1);
133
134 // This loop gives shellNum the value of the index of shellId
135 // in the vector storing the list of the shells reachable through
136 // a radiative transition
137 if ( shellId <= refShell->FinalShellId())
138 {
139 while (shellId != transitionManager->ReachableShell(Z,shellNum)->FinalShellId())
140 {
141 if(shellNum ==maxNumOfShells-1)
142 {
143 break;
144 }
145 shellNum++;
146 }
147 G4int transProb = 0; //AM change 29/6/07 was 1
148
149 G4double partialProb = G4UniformRand();
150 G4double partSum = 0;
151 const G4FluoTransition* aShell = transitionManager->ReachableShell(Z,shellNum);
152 G4int trSize = (aShell->TransitionProbabilities()).size();
153
154 // Loop over the shells wich can provide an electron for a
155 // radiative transition towards shellId:
156 // in every loop the partial sum of the first transProb shells
157 // is calculated and compared with a random number [0,1].
158 // If the partial sum is greater, the shell whose index is transProb
159 // is chosen as the starting shell for a radiative transition
160 // and its identity is returned
161 // Else, terminateded the loop, -1 is returned
162 while(transProb < trSize){
163
164 partSum += aShell->TransitionProbability(transProb);
165
166 if(partialProb <= partSum)
167 {
168 provShellId = aShell->OriginatingShellId(transProb);
169 fluoTransitionFoundFlag = true;
170
171 break;
172 }
173 transProb++;
174 }
175
176 // here provShellId is the right one or is -1.
177 // if -1, the control is passed to the Auger generation part of the package
178 }
179
180
181
182 else
183 {
184
185 provShellId = -1;
186
187 }
188 return provShellId;
189}
190
191G4DynamicParticle* G4AtomicDeexcitation::GenerateFluorescence(G4int Z,
192 G4int shellId,
193 G4int provShellId )
194{
195
196
197 const G4AtomicTransitionManager* transitionManager = G4AtomicTransitionManager::Instance();
198 // G4int provenienceShell = provShellId;
199
200 //isotropic angular distribution for the outcoming photon
201 G4double newcosTh = 1.-2.*G4UniformRand();
202 G4double newsinTh = std::sqrt(1.-newcosTh*newcosTh);
203 G4double newPhi = twopi*G4UniformRand();
204
205 G4double xDir = newsinTh*std::sin(newPhi);
206 G4double yDir = newsinTh*std::cos(newPhi);
207 G4double zDir = newcosTh;
208
209 G4ThreeVector newGammaDirection(xDir,yDir,zDir);
210
211 G4int shellNum = 0;
212 G4int maxNumOfShells = transitionManager->NumberOfReachableShells(Z);
213
214 // find the index of the shell named shellId
215 while (shellId != transitionManager->
216 ReachableShell(Z,shellNum)->FinalShellId())
217 {
218 if(shellNum == maxNumOfShells-1)
219 {
220 break;
221 }
222 shellNum++;
223 }
224 // number of shell from wich an electron can reach shellId
225 size_t transitionSize = transitionManager->
226 ReachableShell(Z,shellNum)->OriginatingShellIds().size();
227
228 size_t index = 0;
229
230 // find the index of the shell named provShellId in the vector
231 // storing the shells from which shellId can be reached
232 while (provShellId != transitionManager->
233 ReachableShell(Z,shellNum)->OriginatingShellId(index))
234 {
235 if(index == transitionSize-1)
236 {
237 break;
238 }
239 index++;
240 }
241 // energy of the gamma leaving provShellId for shellId
242 G4double transitionEnergy = transitionManager->
243 ReachableShell(Z,shellNum)->TransitionEnergy(index);
244
245 // This is the shell where the new vacancy is: it is the same
246 // shell where the electron came from
247 newShellId = transitionManager->
248 ReachableShell(Z,shellNum)->OriginatingShellId(index);
249
250
251 G4DynamicParticle* newPart = new G4DynamicParticle(G4Gamma::Gamma(),
252 newGammaDirection,
253 transitionEnergy);
254 return newPart;
255}
256
257G4DynamicParticle* G4AtomicDeexcitation::GenerateAuger(G4int Z, G4int shellId)
258{
259 if(!fAuger) return 0;
260
261
262 const G4AtomicTransitionManager* transitionManager =
263 G4AtomicTransitionManager::Instance();
264
265
266
267 if (shellId <=0 )
268 {G4Exception("G4AtomicDeexcitation: zero or negative shellId");}
269
270 // G4int provShellId = -1;
271 G4int maxNumOfShells = transitionManager->NumberOfReachableAugerShells(Z);
272
273 const G4AugerTransition* refAugerTransition =
274 transitionManager->ReachableAugerShell(Z,maxNumOfShells-1);
275
276
277 // This loop gives to shellNum the value of the index of shellId
278 // in the vector storing the list of the vacancies in the variuos shells
279 // that can originate a NON-radiative transition
280
281 // ---- MGP ---- Next line commented out to remove compilation warning
282 // G4int p = refAugerTransition->FinalShellId();
283
284 G4int shellNum = 0;
285
286
287 if ( shellId <= refAugerTransition->FinalShellId() )
288 //"FinalShellId" is final from the point of view of the elctron who makes the transition,
289 // being the Id of the shell in which there is a vacancy
290 {
291 G4int pippo = transitionManager->ReachableAugerShell(Z,shellNum)->FinalShellId();
292 if (shellId != pippo ) {
293 do {
294 shellNum++;
295 if(shellNum == maxNumOfShells)
296 {
297// G4cout << "G4AtomicDeexcitation warning: No Auger transition found" << G4endl;
298// G4cout << "Absorbed enrgy deposited locally" << G4endl;
[1055]299 G4Exception("G4AtomicDeexcitation: No Auger transition found");
300 return 0;
[819]301 }
302 }
303 while (shellId != (transitionManager->ReachableAugerShell(Z,shellNum)->FinalShellId()) ) ;
304 }
305 /* {
306
307 if(shellNum == maxNumOfShells-1)
308 {
309 G4Exception("G4AtomicDeexcitation: No Auger tramsition found");
310 }
311 shellNum++;
312 }*/
313
314
315
316
317 // Now we have that shellnum is the shellIndex of the shell named ShellId
318
319 // G4cout << " the index of the shell is: "<<shellNum<<G4endl;
320
321 // But we have now to select two shells: one for the transition,
322 // and another for the auger emission.
323
324 G4int transitionLoopShellIndex = 0;
325 G4double partSum = 0;
326 const G4AugerTransition* anAugerTransition =
327 transitionManager->ReachableAugerShell(Z,shellNum);
328
329 // G4cout << " corresponding to the ID: "<< anAugerTransition->FinalShellId() << G4endl;
330
331
332 G4int transitionSize =
333 (anAugerTransition->TransitionOriginatingShellIds())->size();
334 while (transitionLoopShellIndex < transitionSize) {
335
336 std::vector<G4int>::const_iterator pos =
337 anAugerTransition->TransitionOriginatingShellIds()->begin();
338
339 G4int transitionLoopShellId = *(pos+transitionLoopShellIndex);
340 G4int numberOfPossibleAuger =
341 (anAugerTransition->AugerTransitionProbabilities(transitionLoopShellId))->size();
342 G4int augerIndex = 0;
343 // G4int partSum2 = 0;
344
345
346 if (augerIndex < numberOfPossibleAuger) {
347
348 do
349 {
350 G4double thisProb = anAugerTransition->AugerTransitionProbability(augerIndex,
351 transitionLoopShellId);
352 partSum += thisProb;
353 augerIndex++;
354
355 } while (augerIndex < numberOfPossibleAuger);
356 }
357 transitionLoopShellIndex++;
358 }
359
360
361
362 // Now we have the entire probability of an auger transition for the vacancy
363 // located in shellNum (index of shellId)
364
365 // AM *********************** F I X E D **************************** AM
366 // Here we duplicate the previous loop, this time looking to the sum of the probabilities
367 // to be under the random number shoot by G4 UniformRdandom. This could have been done in the
368 // previuos loop, while integrating the probabilities. There is a bug that will be fixed
369 // 5 minutes from now: a line:
370 // G4int numberOfPossibleAuger = (anAugerTransition->
371 // AugerTransitionProbabilities(transitionLoopShellId))->size();
372 // to be inserted.
373 // AM *********************** F I X E D **************************** AM
374
375 // Remains to get the same result with a single loop.
376
377 // AM *********************** F I X E D **************************** AM
378 // Another Bug: in EADL Auger Transition are normalized to all the transitions deriving from
379 // a vacancy in one shell, but not all of these are present in data tables. So if a transition
380 // doesn't occur in the main one a local energy deposition must occur, instead of (like now)
381 // generating the last transition present in EADL data.
382 // AM *********************** F I X E D **************************** AM
383
384
385 G4double totalVacancyAugerProbability = partSum;
386
387
388 //And now we start to select the right auger transition and emission
389 G4int transitionRandomShellIndex = 0;
390 G4int transitionRandomShellId = 1;
391 G4int augerIndex = 0;
392 partSum = 0;
393 G4double partialProb = G4UniformRand();
394 // G4int augerOriginatingShellId = 0;
395
[1055]396 G4int numberOfPossibleAuger = 0;
397 numberOfPossibleAuger = anAugerTransition->AugerTransitionProbabilities(transitionRandomShellId)->size();
398
399
[819]400 G4bool foundFlag = false;
401
402 while (transitionRandomShellIndex < transitionSize) {
403
404 std::vector<G4int>::const_iterator pos =
405 anAugerTransition->TransitionOriginatingShellIds()->begin();
406
407 transitionRandomShellId = *(pos+transitionRandomShellIndex);
408
409 augerIndex = 0;
410 numberOfPossibleAuger = (anAugerTransition->
411 AugerTransitionProbabilities(transitionRandomShellId))->size();
412
413 while (augerIndex < numberOfPossibleAuger) {
414 G4double thisProb =anAugerTransition->AugerTransitionProbability(augerIndex,
415 transitionRandomShellId);
416
417 partSum += thisProb;
418
419 if (partSum >= (partialProb*totalVacancyAugerProbability) ) { // was /
420 foundFlag = true;
421 break;
422 }
423 augerIndex++;
424 }
425 if (partSum >= (partialProb*totalVacancyAugerProbability) ) {break;} // was /
426 transitionRandomShellIndex++;
427 }
428
429 // Now we have the index of the shell from wich comes the auger electron (augerIndex),
430 // and the id of the shell, from which the transition e- come (transitionRandomShellid)
431 // If no Transition has been found, 0 is returned.
432
433 if (!foundFlag) {return 0;}
434
435 // Isotropic angular distribution for the outcoming e-
436 G4double newcosTh = 1.-2.*G4UniformRand();
437 G4double newsinTh = std::sqrt(1.-newcosTh*newcosTh);
438 G4double newPhi = twopi*G4UniformRand();
439
440 G4double xDir = newsinTh*std::sin(newPhi);
441 G4double yDir = newsinTh*std::cos(newPhi);
442 G4double zDir = newcosTh;
443
444 G4ThreeVector newElectronDirection(xDir,yDir,zDir);
445
446 // energy of the auger electron emitted
447
448
449 G4double transitionEnergy = anAugerTransition->AugerTransitionEnergy(augerIndex, transitionRandomShellId);
450 /*
451 G4cout << "AUger TransitionId " << anAugerTransition->FinalShellId() << G4endl;
452 G4cout << "augerIndex: " << augerIndex << G4endl;
453 G4cout << "transitionShellId: " << transitionRandomShellId << G4endl;
454 */
455
456 // This is the shell where the new vacancy is: it is the same
457 // shell where the electron came from
458 newShellId = transitionRandomShellId;
459
460
461 G4DynamicParticle* newPart = new G4DynamicParticle(G4Electron::Electron(),
462 newElectronDirection,
463 transitionEnergy);
464 return newPart;
465
466 }
467 else
468 {
469 //G4Exception("G4AtomicDeexcitation: no auger transition found");
470 return 0;
471 }
472
473}
474
475void G4AtomicDeexcitation::SetCutForSecondaryPhotons(G4double cut)
476{
477 minGammaEnergy = cut;
478}
479
480void G4AtomicDeexcitation::SetCutForAugerElectrons(G4double cut)
481{
482 minElectronEnergy = cut;
483}
484
485void G4AtomicDeexcitation::ActivateAugerElectronProduction(G4bool val)
486{
487 fAuger = val;
488}
489
490
491
492
493
494
495
Note: See TracBrowser for help on using the repository browser.