ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4LowEnergyPhotoElectric.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4LowEnergyPhotoElectric.cc
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 //
28 // Author: A. Forti
29 // Maria Grazia Pia (Maria.Grazia.Pia@cern.ch)
30 //
31 // History:
32 // --------
33 // October 1998 - low energy modifications by Alessandra Forti
34 // Added Livermore data table construction methods A. Forti
35 // Modified BuildMeanFreePath to read new data tables A. Forti
36 // Added EnergySampling method A. Forti
37 // Modified PostStepDoIt to insert sampling with EPDL97 data A. Forti
38 // Added SelectRandomAtom A. Forti
39 // Added map of the elements A. Forti
40 // 10.04.2000 VL
41 // - Correcting Fluorescence transition probabilities in order to take into account
42 // non-radiative transitions. No Auger electron simulated yet: energy is locally deposited.
43 // 17.02.2000 Veronique Lefebure
44 // - bugs corrected in fluorescence simulation:
45 // . when final use of binding energy: no photon was ever created
46 // . no Fluorescence was simulated when the photo-electron energy
47 // was below production threshold.
48 //
49 // 07-09-99, if no e- emitted: edep=photon energy, mma
50 // 24.04.01 V.Ivanchenko remove RogueWave
51 // 12.08.2001 MGP Revised according to a design iteration
52 // 16.09.2001 E. Guardincerri Added fluorescence generation
53 // 06.10.2001 MGP Added protection to avoid negative electron energies
54 // when binding energy of selected shell > photon energy
55 // 18.04.2001 V.Ivanchenko Fix problem with low energy gammas from fluorescence
56 // MeanFreePath is calculated by crosSectionHandler directly
57 // 31.05.2002 V.Ivanchenko Add path of Fluo + Auger cuts to AtomicDeexcitation
58 // 14.06.2002 V.Ivanchenko By default do not cheak range of e-
59 // 21.01.2003 V.Ivanchenko Cut per region
60 // 10.05.2004 P.Rodrigues Changes to accommodate new angular generators
61 // 20.01.2006 A.Trindade Changes to accommodate polarized angular generator
62 //
63 // --------------------------------------------------------------
64 
66 
71 
72 #include "G4PhysicalConstants.hh"
73 #include "G4SystemOfUnits.hh"
74 #include "G4ParticleDefinition.hh"
75 #include "G4Track.hh"
76 #include "G4Step.hh"
77 #include "G4ForceCondition.hh"
78 #include "G4Gamma.hh"
79 #include "G4Electron.hh"
80 #include "G4DynamicParticle.hh"
81 #include "G4VParticleChange.hh"
82 #include "G4ThreeVector.hh"
85 #include "G4RDVEMDataSet.hh"
87 #include "G4RDVDataSetAlgorithm.hh"
89 #include "G4RDVRangeTest.hh"
90 #include "G4RDRangeNoTest.hh"
92 #include "G4RDAtomicShell.hh"
93 #include "G4ProductionCutsTable.hh"
94 
96  : G4VDiscreteProcess(processName), lowEnergyLimit(250*eV), highEnergyLimit(100*GeV),
97  intrinsicLowEnergyLimit(10*eV),
98  intrinsicHighEnergyLimit(100*GeV),
99  cutForLowEnergySecondaryPhotons(250.*eV),
100  cutForLowEnergySecondaryElectrons(250.*eV)
101 {
104  {
105  G4Exception("G4LowEnergyPhotoElectric::G4LowEnergyPhotoElectric()",
106  "OutOfRange", FatalException,
107  "Energy limit outside intrinsic process validity range!");
108  }
109 
112  meanFreePathTable = 0;
114  generatorName = "geant4.6.2";
115  ElectronAngularGenerator = new G4RDPhotoElectricAngularGeneratorSimple("GEANTSimpleGenerator"); // default generator
116 
117 
118  if (verboseLevel > 0)
119  {
120  G4cout << GetProcessName() << " is created " << G4endl
121  << "Energy range: "
122  << lowEnergyLimit / keV << " keV - "
123  << highEnergyLimit / GeV << " GeV"
124  << G4endl;
125  }
126 }
127 
129 {
130  delete crossSectionHandler;
132  delete meanFreePathTable;
133  delete rangeTest;
135 }
136 
138 {
139 
141  G4String crossSectionFile = "phot/pe-cs-";
142  crossSectionHandler->LoadData(crossSectionFile);
143 
145  G4String shellCrossSectionFile = "phot/pe-ss-cs-";
146  shellCrossSectionHandler->LoadShellData(shellCrossSectionFile);
147 
148  delete meanFreePathTable;
150 }
151 
153  const G4Step& aStep)
154 {
155  // Fluorescence generated according to:
156  // J. Stepanek ,"A program to determine the radiation spectra due to a single atomic
157  // subshell ionisation by a particle or due to deexcitation or decay of radionuclides",
158  // Comp. Phys. Comm. 1206 pp 1-1-9 (1997)
159 
160  aParticleChange.Initialize(aTrack);
161 
162  const G4DynamicParticle* incidentPhoton = aTrack.GetDynamicParticle();
163  G4double photonEnergy = incidentPhoton->GetKineticEnergy();
164  if (photonEnergy <= lowEnergyLimit)
165  {
169  return G4VDiscreteProcess::PostStepDoIt(aTrack,aStep);
170  }
171 
172  G4ThreeVector photonDirection = incidentPhoton->GetMomentumDirection(); // Returns the normalized direction of the momentum
173 
174  // Select randomly one element in the current material
175  const G4MaterialCutsCouple* couple = aTrack.GetMaterialCutsCouple();
176  G4int Z = crossSectionHandler->SelectRandomAtom(couple,photonEnergy);
177 
178  // Select the ionised shell in the current atom according to shell cross sections
179  size_t shellIndex = shellCrossSectionHandler->SelectRandomShell(Z,photonEnergy);
180 
181  // Retrieve the corresponding identifier and binding energy of the selected shell
183  const G4RDAtomicShell* shell = transitionManager->Shell(Z,shellIndex);
185  G4int shellId = shell->ShellId();
186 
187  // Create lists of pointers to DynamicParticles (photons and electrons)
188  // (Is the electron vector necessary? To be checked)
189  std::vector<G4DynamicParticle*>* photonVector = 0;
190  std::vector<G4DynamicParticle*> electronVector;
191 
192  G4double energyDeposit = 0.0;
193 
194  // Primary outcoming electron
195  G4double eKineticEnergy = photonEnergy - bindingEnergy;
196 
197  // There may be cases where the binding energy of the selected shell is > photon energy
198  // In such cases do not generate secondaries
199  if (eKineticEnergy > 0.)
200  {
201  // Generate the electron only if with large enough range w.r.t. cuts and safety
202  G4double safety = aStep.GetPostStepPoint()->GetSafety();
203 
204  if (rangeTest->Escape(G4Electron::Electron(),couple,eKineticEnergy,safety))
205  {
206 
207  // Calculate direction of the photoelectron
208  G4ThreeVector gammaPolarization = incidentPhoton->GetPolarization();
209  G4ThreeVector electronDirection = ElectronAngularGenerator->GetPhotoElectronDirection(photonDirection,eKineticEnergy,gammaPolarization,shellId);
210 
211  // The electron is created ...
213  electronDirection,
214  eKineticEnergy);
215  electronVector.push_back(electron);
216  }
217  else
218  {
219  energyDeposit += eKineticEnergy;
220  }
221  }
222  else
223  {
224  bindingEnergy = photonEnergy;
225  }
226 
227  G4int nElectrons = electronVector.size();
228  size_t nTotPhotons = 0;
229  G4int nPhotons=0;
230  const G4ProductionCutsTable* theCoupleTable=
232 
233  size_t index = couple->GetIndex();
234  G4double cutg = (*(theCoupleTable->GetEnergyCutsVector(0)))[index];
236 
237  G4double cute = (*(theCoupleTable->GetEnergyCutsVector(1)))[index];
239 
240  G4DynamicParticle* aPhoton;
241 
242  // Generation of fluorescence
243  // Data in EADL are available only for Z > 5
244  // Protection to avoid generating photons in the unphysical case of
245  // shell binding energy > photon energy
246  if (Z > 5 && (bindingEnergy > cutg || bindingEnergy > cute))
247  {
248  photonVector = deexcitationManager.GenerateParticles(Z,shellId);
249  nTotPhotons = photonVector->size();
250  for (size_t k=0; k<nTotPhotons; k++)
251  {
252  aPhoton = (*photonVector)[k];
253  if (aPhoton)
254  {
255  G4double itsCut = cutg;
256  if(aPhoton->GetDefinition() == G4Electron::Electron()) itsCut = cute;
257  G4double itsEnergy = aPhoton->GetKineticEnergy();
258 
259  if (itsEnergy > itsCut && itsEnergy <= bindingEnergy)
260  {
261  nPhotons++;
262  // Local energy deposit is given as the sum of the
263  // energies of incident photons minus the energies
264  // of the outcoming fluorescence photons
265  bindingEnergy -= itsEnergy;
266 
267  }
268  else
269  {
270  delete aPhoton;
271  (*photonVector)[k] = 0;
272  }
273  }
274  }
275  }
276 
277  energyDeposit += bindingEnergy;
278 
279  G4int nSecondaries = nElectrons + nPhotons;
281 
282  for (G4int l = 0; l<nElectrons; l++ )
283  {
284  aPhoton = electronVector[l];
285  if(aPhoton) {
286  aParticleChange.AddSecondary(aPhoton);
287  }
288  }
289  for ( size_t ll = 0; ll < nTotPhotons; ll++)
290  {
291  aPhoton = (*photonVector)[ll];
292  if(aPhoton) {
293  aParticleChange.AddSecondary(aPhoton);
294  }
295  }
296 
297  delete photonVector;
298 
299  if (energyDeposit < 0)
300  {
301  G4cout << "WARNING - "
302  << "G4LowEnergyPhotoElectric::PostStepDoIt - Negative energy deposit"
303  << G4endl;
304  energyDeposit = 0;
305  }
306 
307  // Kill the incident photon
310 
313 
314  // Reset NbOfInteractionLengthLeft and return aParticleChange
315  return G4VDiscreteProcess::PostStepDoIt( aTrack, aStep );
316 }
317 
319 {
320  return ( &particle == G4Gamma::Gamma() );
321 }
322 
324  G4double, // previousStepSize
326 {
328  G4double energy = photon->GetKineticEnergy();
329  G4Material* material = track.GetMaterial();
330  // size_t materialIndex = material->GetIndex();
331 
332  G4double meanFreePath = DBL_MAX;
333 
334  // if (energy > highEnergyLimit)
335  // meanFreePath = meanFreePathTable->FindValue(highEnergyLimit,materialIndex);
336  // else if (energy < lowEnergyLimit) meanFreePath = DBL_MAX;
337  // else meanFreePath = meanFreePathTable->FindValue(energy,materialIndex);
338 
340  if(cross > 0.0) meanFreePath = 1.0/cross;
341 
342  return meanFreePath;
343 }
344 
346 {
349 }
350 
352 {
355 }
356 
358 {
360 }
361 
363 {
364  ElectronAngularGenerator = distribution;
366 }
367 
369 {
370  if (name == "default")
371  {
373  ElectronAngularGenerator = new G4RDPhotoElectricAngularGeneratorSimple("GEANT4LowEnergySimpleGenerator");
375  }
376  else if (name == "standard")
377  {
379  ElectronAngularGenerator = new G4RDPhotoElectricAngularGeneratorSauterGavrila("GEANT4SauterGavrilaGenerator");
381  }
382  else if (name == "polarized")
383  {
385  ElectronAngularGenerator = new G4RDPhotoElectricAngularGeneratorPolarized("GEANT4LowEnergyPolarizedGenerator");
387  }
388  else
389  {
390  G4Exception("G4LowEnergyPhotoElectric::SetAngularGenerator()",
391  "InvalidSetup", FatalException,
392  "Generator does not exist!");
393  }
394 
396 }