source: trunk/source/processes/electromagnetic/lowenergy/src/G4LivermorePhotoElectricModel.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: 12.5 KB
Line 
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// $Id: G4LivermorePhotoElectricModel.cc,v 1.7 2009/04/18 18:29:34 vnivanch Exp $
27// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
28//
29//
30// Author: Sebastien Inserti
31// 30 October 2008
32//
33// History:
34// --------
35// 15 Apr 2009 V Ivanchenko Cleanup initialisation and generation of secondaries:
36// - apply internal high-energy limit only in constructor
37// - do not apply low-energy limit (default is 0)
38// - remove GetMeanFreePath method and table
39// - simplify sampling of deexcitation by using cut in energy
40// - added protection against numerical problem in energy sampling
41// - use G4ElementSelector
42
43#include "G4LivermorePhotoElectricModel.hh"
44
45//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
46
47using namespace std;
48
49//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
50
51G4LivermorePhotoElectricModel::G4LivermorePhotoElectricModel(const G4ParticleDefinition*,
52 const G4String& nam)
53:G4VEmModel(nam),isInitialised(false),meanFreePathTable(0),
54 crossSectionHandler(0),shellCrossSectionHandler(0),ElectronAngularGenerator(0)
55{
56 lowEnergyLimit = 250 * eV;
57 highEnergyLimit = 100 * GeV;
58 // SetLowEnergyLimit(lowEnergyLimit);
59 SetHighEnergyLimit(highEnergyLimit);
60
61 ActivateAuger(false);
62
63 verboseLevel= 0;
64 // Verbosity scale:
65 // 0 = nothing
66 // 1 = warning for energy non-conservation
67 // 2 = details of energy budget
68 // 3 = calculation of cross sections, file openings, sampling of atoms
69 // 4 = entering in methods
70 if(verboseLevel>0) {
71 G4cout << "Livermore PhotoElectric is constructed " << G4endl
72 << "Energy range: "
73 << lowEnergyLimit / eV << " eV - "
74 << highEnergyLimit / GeV << " GeV"
75 << G4endl;
76 }
77}
78
79//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
80
81G4LivermorePhotoElectricModel::~G4LivermorePhotoElectricModel()
82{
83 if (crossSectionHandler) delete crossSectionHandler;
84 if (shellCrossSectionHandler) delete shellCrossSectionHandler;
85 if (ElectronAngularGenerator) delete ElectronAngularGenerator;
86}
87
88//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
89
90void
91G4LivermorePhotoElectricModel::Initialise(const G4ParticleDefinition*,
92 const G4DataVector&)
93{
94 if (verboseLevel > 3)
95 G4cout << "Calling G4LivermorePhotoElectricModel::Initialise()" << G4endl;
96
97 if (crossSectionHandler)
98 {
99 crossSectionHandler->Clear();
100 delete crossSectionHandler;
101 }
102
103 if (shellCrossSectionHandler)
104 {
105 shellCrossSectionHandler->Clear();
106 delete shellCrossSectionHandler;
107 }
108
109 // Read data tables for all materials
110
111 crossSectionHandler = new G4CrossSectionHandler();
112 crossSectionHandler->Clear();
113 G4String crossSectionFile = "phot/pe-cs-";
114 crossSectionHandler->LoadData(crossSectionFile);
115
116 shellCrossSectionHandler = new G4CrossSectionHandler();
117 shellCrossSectionHandler->Clear();
118 G4String shellCrossSectionFile = "phot/pe-ss-cs-";
119 shellCrossSectionHandler->LoadShellData(shellCrossSectionFile);
120
121 // SI - Simple generator is buggy
122 //generatorName = "geant4.6.2";
123 //ElectronAngularGenerator = new G4PhotoElectricAngularGeneratorSimple("GEANTSimpleGenerator"); // default generator
124 ElectronAngularGenerator =
125 new G4PhotoElectricAngularGeneratorSauterGavrila("GEANTSauterGavrilaGenerator");
126
127 //
128 if (verboseLevel > 2)
129 G4cout << "Loaded cross section files for Livermore PhotoElectric model" << G4endl;
130
131 // InitialiseElementSelectors(particle,cuts);
132
133 if (verboseLevel > 0) {
134 G4cout << "Livermore PhotoElectric model is initialized " << G4endl
135 << "Energy range: "
136 << LowEnergyLimit() / eV << " eV - "
137 << HighEnergyLimit() / GeV << " GeV"
138 << G4endl;
139 }
140
141 if(isInitialised) return;
142 fParticleChange = GetParticleChangeForGamma();
143 isInitialised = true;
144}
145
146//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
147
148G4double G4LivermorePhotoElectricModel::ComputeCrossSectionPerAtom(
149 const G4ParticleDefinition*,
150 G4double GammaEnergy,
151 G4double Z, G4double,
152 G4double, G4double)
153{
154 if (verboseLevel > 3)
155 G4cout << "Calling ComputeCrossSectionPerAtom() of G4LivermorePhotoElectricModel"
156 << G4endl;
157
158 if (GammaEnergy < lowEnergyLimit || GammaEnergy > highEnergyLimit)
159 return 0;
160
161 G4double cs = crossSectionHandler->FindValue(G4int(Z), GammaEnergy);
162 return cs;
163}
164
165//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
166
167void
168G4LivermorePhotoElectricModel::SampleSecondaries(std::vector<G4DynamicParticle*>* fvect,
169 const G4MaterialCutsCouple* couple,
170 const G4DynamicParticle* aDynamicGamma,
171 G4double,
172 G4double)
173{
174
175 // Fluorescence generated according to:
176 // J. Stepanek ,"A program to determine the radiation spectra due to a single atomic
177 // subshell ionisation by a particle or due to deexcitation or decay of radionuclides",
178 // Comp. Phys. Comm. 1206 pp 1-1-9 (1997)
179
180 if (verboseLevel > 3)
181 G4cout << "Calling SampleSecondaries() of G4LivermorePhotoElectricModel" << G4endl;
182
183 G4double photonEnergy = aDynamicGamma->GetKineticEnergy();
184
185 // kill incident photon
186 fParticleChange->SetProposedKineticEnergy(0.);
187 fParticleChange->ProposeTrackStatus(fStopAndKill);
188
189 // low-energy gamma is absorpted by this process
190 if (photonEnergy <= lowEnergyLimit)
191 {
192 fParticleChange->ProposeLocalEnergyDeposit(photonEnergy);
193 return;
194 }
195
196 // Returns the normalized direction of the momentum
197 G4ThreeVector photonDirection = aDynamicGamma->GetMomentumDirection();
198
199 // Select randomly one element in the current material
200 // G4int Z = crossSectionHandler->SelectRandomAtom(couple,photonEnergy);
201 const G4ParticleDefinition* particle = aDynamicGamma->GetDefinition();
202 const G4Element* elm = SelectRandomAtom(couple->GetMaterial(),particle,photonEnergy);
203 G4int Z = (G4int)elm->GetZ();
204
205 // Select the ionised shell in the current atom according to shell cross sections
206 size_t shellIndex = shellCrossSectionHandler->SelectRandomShell(Z,photonEnergy);
207
208 // Retrieve the corresponding identifier and binding energy of the selected shell
209 const G4AtomicTransitionManager* transitionManager = G4AtomicTransitionManager::Instance();
210 const G4AtomicShell* shell = transitionManager->Shell(Z,shellIndex);
211 G4double bindingEnergy = shell->BindingEnergy();
212 G4int shellId = shell->ShellId();
213
214 // Primary outcoming electron
215 G4double eKineticEnergy = photonEnergy - bindingEnergy;
216
217 // There may be cases where the binding energy of the selected shell is > photon energy
218 // In such cases do not generate secondaries
219 if (eKineticEnergy > 0.)
220 {
221 // Calculate direction of the photoelectron
222 G4ThreeVector gammaPolarization = aDynamicGamma->GetPolarization();
223 G4ThreeVector electronDirection =
224 ElectronAngularGenerator->GetPhotoElectronDirection(photonDirection,
225 eKineticEnergy,
226 gammaPolarization,
227 shellId);
228
229 // The electron is created ...
230 G4DynamicParticle* electron = new G4DynamicParticle (G4Electron::Electron(),
231 electronDirection,
232 eKineticEnergy);
233 fvect->push_back(electron);
234 }
235 else
236 {
237 bindingEnergy = photonEnergy;
238 }
239
240 // deexcitation
241 if(DeexcitationFlag() && Z > 5) {
242
243 const G4ProductionCutsTable* theCoupleTable=
244 G4ProductionCutsTable::GetProductionCutsTable();
245
246 size_t index = couple->GetIndex();
247 G4double cutg = (*(theCoupleTable->GetEnergyCutsVector(0)))[index];
248 G4double cute = (*(theCoupleTable->GetEnergyCutsVector(1)))[index];
249
250 // Generation of fluorescence
251 // Data in EADL are available only for Z > 5
252 // Protection to avoid generating photons in the unphysical case of
253 // shell binding energy > photon energy
254 if (bindingEnergy > cutg || bindingEnergy > cute)
255 {
256 G4DynamicParticle* aPhoton;
257 deexcitationManager.SetCutForSecondaryPhotons(cutg);
258 deexcitationManager.SetCutForAugerElectrons(cute);
259
260 std::vector<G4DynamicParticle*>* photonVector =
261 deexcitationManager.GenerateParticles(Z,shellId);
262 size_t nTotPhotons = photonVector->size();
263 for (size_t k=0; k<nTotPhotons; k++)
264 {
265 aPhoton = (*photonVector)[k];
266 if (aPhoton)
267 {
268 G4double itsEnergy = aPhoton->GetKineticEnergy();
269 if (itsEnergy <= bindingEnergy)
270 {
271 // Local energy deposit is given as the sum of the
272 // energies of incident photons minus the energies
273 // of the outcoming fluorescence photons
274 bindingEnergy -= itsEnergy;
275 fvect->push_back(aPhoton);
276 }
277 else
278 {
279 // abnormal case of energy non-conservation
280 delete aPhoton;
281 (*photonVector)[k] = 0;
282 }
283 }
284 }
285 delete photonVector;
286 }
287 }
288 // excitation energy left
289 fParticleChange->ProposeLocalEnergyDeposit(bindingEnergy);
290}
291
292//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
293
294void G4LivermorePhotoElectricModel::SetCutForLowEnSecPhotons(G4double cut)
295{
296 cutForLowEnergySecondaryPhotons = cut;
297 deexcitationManager.SetCutForSecondaryPhotons(cut);
298}
299
300//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
301
302void G4LivermorePhotoElectricModel::SetCutForLowEnSecElectrons(G4double cut)
303{
304 cutForLowEnergySecondaryElectrons = cut;
305 deexcitationManager.SetCutForAugerElectrons(cut);
306}
307
308//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
309
310void G4LivermorePhotoElectricModel::ActivateAuger(G4bool val)
311{
312 deexcitationManager.ActivateAugerElectronProduction(val);
313}
314
315//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
316
317void
318G4LivermorePhotoElectricModel::SetAngularGenerator(G4VPhotoElectricAngularDistribution* dist)
319{
320 ElectronAngularGenerator = dist;
321 ElectronAngularGenerator->PrintGeneratorInformation();
322}
323
324//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
325
326void G4LivermorePhotoElectricModel::SetAngularGenerator(const G4String& name)
327{
328 if (name == "default")
329 {
330 delete ElectronAngularGenerator;
331 ElectronAngularGenerator =
332 new G4PhotoElectricAngularGeneratorSimple("GEANT4LowEnergySimpleGenerator");
333 generatorName = name;
334 }
335 else if (name == "standard")
336 {
337 delete ElectronAngularGenerator;
338 ElectronAngularGenerator =
339 new G4PhotoElectricAngularGeneratorSauterGavrila("GEANT4SauterGavrilaGenerator");
340 generatorName = name;
341 }
342 else if (name == "polarized")
343 {
344 delete ElectronAngularGenerator;
345 ElectronAngularGenerator =
346 new G4PhotoElectricAngularGeneratorPolarized("GEANT4LowEnergyPolarizedGenerator");
347 generatorName = name;
348 }
349 else
350 {
351 G4Exception("G4LowEnergyPhotoElectric::SetAngularGenerator - generator does not exist");
352 }
353
354 ElectronAngularGenerator->PrintGeneratorInformation();
355}
356
Note: See TracBrowser for help on using the repository browser.