Geant4_10
G4LivermorePolarizedPhotoElectricModel.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 // $Id: G4LivermorePolarizedPhotoElectricModel.hh,v 1.2 2010-11-23 16:42:15 flongo Exp $
27 //
28 // Authors: G.Depaola & F.Longo
29 //
30 
32 #include "G4PhysicalConstants.hh"
33 #include "G4SystemOfUnits.hh"
34 #include "G4LossTableManager.hh"
35 #include "G4VAtomDeexcitation.hh"
38 #include "G4AtomicShell.hh"
39 #include "G4ProductionCutsTable.hh"
40 #include "G4Gamma.hh"
41 
42 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
43 
44 using namespace std;
45 
46 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
47 
49  const G4String& nam)
50  :G4VEmModel(nam),fParticleChange(0),
51  crossSectionHandler(0), shellCrossSectionHandler(0)
52 {
53  lowEnergyLimit = 10 * eV; // SI - Could be 10 eV ?
54  highEnergyLimit = 100 * GeV;
55 
56  verboseLevel= 0;
57  // Verbosity scale:
58  // 0 = nothing
59  // 1 = warning for energy non-conservation
60  // 2 = details of energy budget
61  // 3 = calculation of cross sections, file openings, sampling of atoms
62  // 4 = entering in methods
63 
64  theGamma = G4Gamma::Gamma();
65  theElectron = G4Electron::Electron();
66 
67  // use default generator
69 
70  // polarized generator needs fix
71  //SetAngularDistribution(new G4PhotoElectricAngularGeneratorPolarized());
72  SetDeexcitationFlag(true);
73  fAtomDeexcitation = 0;
74  fDeexcitationActive = false;
75 
76  if (verboseLevel > 1) {
77  G4cout << "Livermore Polarized PhotoElectric is constructed " << G4endl
78  << "Energy range: "
79  << lowEnergyLimit / keV << " keV - "
80  << highEnergyLimit / GeV << " GeV"
81  << G4endl;
82  }
83 }
84 
85 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
86 
88 {
89  delete crossSectionHandler;
90  delete shellCrossSectionHandler;
91 }
92 
93 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
94 
95 void
97  const G4ParticleDefinition*,
98  const G4DataVector&)
99 {
100  if (verboseLevel > 3) {
101  G4cout << "Calling G4LivermorePolarizedPhotoElectricModel::Initialise()" << G4endl;
102  }
103  if (crossSectionHandler)
104  {
105  crossSectionHandler->Clear();
106  delete crossSectionHandler;
107  }
108 
109  if (shellCrossSectionHandler)
110  {
111  shellCrossSectionHandler->Clear();
112  delete shellCrossSectionHandler;
113  }
114 
115 
116  // Reading of data files - all materials are read
117  crossSectionHandler = new G4CrossSectionHandler;
118  crossSectionHandler->Clear();
119  G4String crossSectionFile = "phot/pe-cs-";
120  crossSectionHandler->LoadData(crossSectionFile);
121 
122  shellCrossSectionHandler = new G4CrossSectionHandler();
123  shellCrossSectionHandler->Clear();
124  G4String shellCrossSectionFile = "phot/pe-ss-cs-";
125  shellCrossSectionHandler->LoadShellData(shellCrossSectionFile);
126 
128 
129  fAtomDeexcitation = G4LossTableManager::Instance()->AtomDeexcitation();
130  if(fAtomDeexcitation) {
131  fDeexcitationActive = fAtomDeexcitation->IsFluoActive();
132  }
133  if (verboseLevel > 1) {
134  G4cout << "Livermore Polarized PhotoElectric model is initialized "
135  << G4endl
136  << "Energy range: "
137  << LowEnergyLimit() / keV << " keV - "
138  << HighEnergyLimit() / GeV << " GeV"
139  << G4endl;
140  }
141 }
142 
143 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
144 
146  const G4ParticleDefinition*,
147  G4double GammaEnergy,
150 {
151  G4double cs = crossSectionHandler->FindValue(G4lrint(Z), GammaEnergy);
152 
153  if (verboseLevel > 3) {
154  G4cout << "G4LivermorePolarizedPhotoElectricModel::ComputeCrossSectionPerAtom()"
155  << G4endl;
156  G4cout << " E(keV)= " << GammaEnergy/keV << " Z= " << Z
157  << " CrossSection(barn)= "
158  << cs/barn << G4endl;
159  }
160  return cs;
161 }
162 
163 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
164 
166  std::vector<G4DynamicParticle*>* fvect,
167  const G4MaterialCutsCouple* couple,
168  const G4DynamicParticle* aDynamicGamma,
169  G4double,
170  G4double)
171 {
172  if (verboseLevel > 3) {
173  G4cout << "Calling SampleSecondaries() of G4LivermorePolarizedPhotoElectricModel"
174  << G4endl;
175  }
176 
177  G4double photonEnergy = aDynamicGamma->GetKineticEnergy();
178 
179  // kill incident photon
182 
183  // low-energy gamma is absorpted by this process
184 
185  if (photonEnergy <= lowEnergyLimit)
186  {
188  return;
189  }
190 
191  // Select randomly one element in the current material
192  //G4cout << "Select random atom Egamma(keV)= " << photonEnergy/keV << G4endl;
193  const G4Element* elm =
194  SelectRandomAtom(couple->GetMaterial(),theGamma,photonEnergy);
195  G4int Z = G4lrint(elm->GetZ());
196 
197  // Select the ionised shell in the current atom according to shell cross sections
198  //G4cout << "Select random shell Z= " << Z << G4endl;
199  size_t shellIndex = shellCrossSectionHandler->SelectRandomShell(Z,photonEnergy);
200  //G4cout << "Shell index= " << shellIndex
201  // << " nShells= " << elm->GetNbOfAtomicShells() << G4endl;
203  const G4AtomicShell* shell = 0;
204  if(fDeexcitationActive) {
206  shell = fAtomDeexcitation->GetAtomicShell(Z, as);
207  bindingEnergy = shell->BindingEnergy();
208  } else {
209  G4int nshells = elm->GetNbOfAtomicShells() - 1;
210  if(G4int(shellIndex) > nshells) { shellIndex = std::max(0, nshells); }
211  bindingEnergy = elm->GetAtomicShell(shellIndex);
212  }
213 
214  // There may be cases where the binding energy of the selected
215  // shell is > photon energy
216  // In such cases do not generate secondaries
217  if(photonEnergy < bindingEnergy) {
219  return;
220  }
221  //G4cout << "Z= " << Z << " shellIndex= " << shellIndex
222  // << " Ebind(keV)= " << bindingEnergy/keV << G4endl;
223 
224 
225  // Primary outcoming electron
226  G4double eKineticEnergy = photonEnergy - bindingEnergy;
228 
229  // Calculate direction of the photoelectron
230  G4ThreeVector electronDirection =
231  GetAngularDistribution()->SampleDirection(aDynamicGamma,
232  eKineticEnergy,
233  shellIndex,
234  couple->GetMaterial());
235 
236  // The electron is created ...
237  G4DynamicParticle* electron = new G4DynamicParticle (theElectron,
238  electronDirection,
239  eKineticEnergy);
240  fvect->push_back(electron);
241 
242  // Sample deexcitation
243  if(fDeexcitationActive) {
244  G4int index = couple->GetIndex();
245  if(fAtomDeexcitation->CheckDeexcitationActiveRegion(index)) {
246  size_t nbefore = fvect->size();
247 
248  fAtomDeexcitation->GenerateParticles(fvect, shell, Z, index);
249  size_t nafter = fvect->size();
250  if(nafter > nbefore) {
251  for (size_t j=nbefore; j<nafter; ++j) {
252  edep -= ((*fvect)[j])->GetKineticEnergy();
253  }
254  }
255  }
256  }
257 
258  // energy balance - excitation energy left
259  if(edep > 0.0) {
261  }
262 }
263 
264 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
G4double LowEnergyLimit() const
Definition: G4VEmModel.hh:599
G4bool CheckDeexcitationActiveRegion(G4int coupleIndex)
static G4LossTableManager * Instance()
G4int GetNbOfAtomicShells() const
Definition: G4Element.hh:146
G4bool IsFluoActive() const
G4double GetKineticEnergy() const
Int_t index
Definition: macro.C:9
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:592
G4int SelectRandomShell(G4int Z, G4double e) const
G4double GetZ() const
Definition: G4Element.hh:131
G4VEmAngularDistribution * GetAngularDistribution()
Definition: G4VEmModel.hh:578
G4double BindingEnergy() const
int G4int
Definition: G4Types.hh:78
Double_t edep
Definition: macro.C:13
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
virtual void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)
virtual const G4AtomicShell * GetAtomicShell(G4int Z, G4AtomicShellEnumerator shell)=0
G4double FindValue(G4int Z, G4double e) const
G4LivermorePolarizedPhotoElectricModel(const G4String &nam="LivermorePolarizedPhotoElectric")
G4GLOB_DLL std::ostream G4cout
virtual G4ThreeVector & SampleDirection(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, const G4Material *)=0
Float_t Z
Definition: plot.C:39
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0, G4double cut=0, G4double emax=DBL_MAX)
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
void LoadShellData(const G4String &dataFile)
int G4lrint(double ad)
Definition: templates.hh:163
T max(const T t1, const T t2)
brief Return the largest of the two arguments
void SetAngularDistribution(G4VEmAngularDistribution *)
Definition: G4VEmModel.hh:585
void LoadData(const G4String &dataFile)
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &)
static G4Electron * Electron()
Definition: G4Electron.cc:94
void SetProposedKineticEnergy(G4double proposedKinEnergy)
#define G4endl
Definition: G4ios.hh:61
G4VAtomDeexcitation * AtomDeexcitation()
G4double GetAtomicShell(G4int index) const
Definition: G4Element.cc:365
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
void GenerateParticles(std::vector< G4DynamicParticle * > *secVect, const G4AtomicShell *, G4int Z, G4int coupleIndex)
void SetDeexcitationFlag(G4bool val)
Definition: G4VEmModel.hh:739
G4double bindingEnergy(G4int A, G4int Z)
G4AtomicShellEnumerator
const G4Element * SelectRandomAtom(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.hh:510
G4ParticleChangeForGamma * GetParticleChangeForGamma()
Definition: G4VEmModel.cc:121
const G4Material * GetMaterial() const