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