Geant4_10
G4PEEffectFluoModel.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: G4PEEffectFluoModel.cc 73607 2013-09-02 10:04:03Z gcosmo $
27 //
28 // -------------------------------------------------------------------
29 //
30 // GEANT4 Class file
31 //
32 //
33 // File name: G4PEEffectFluoModel
34 //
35 // Author: Vladimir Ivanchenko on base of G4PEEffectModel
36 //
37 // Creation date: 13.06.2010
38 //
39 // Modifications:
40 //
41 // Class Description:
42 // Implementation of the photo-electric effect with deexcitation
43 //
44 // -------------------------------------------------------------------
45 //
46 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
47 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
48 
49 #include "G4PEEffectFluoModel.hh"
50 #include "G4PhysicalConstants.hh"
51 #include "G4SystemOfUnits.hh"
52 #include "G4Electron.hh"
53 #include "G4Gamma.hh"
54 #include "Randomize.hh"
55 #include "G4Material.hh"
56 #include "G4DataVector.hh"
58 #include "G4VAtomDeexcitation.hh"
59 #include "G4LossTableManager.hh"
61 
62 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
63 
64 using namespace std;
65 
67  : G4VEmModel(nam)
68 {
69  theGamma = G4Gamma::Gamma();
70  theElectron = G4Electron::Electron();
71  fminimalEnergy = 1.0*eV;
72  SetDeexcitationFlag(true);
73  fParticleChange = 0;
74  fAtomDeexcitation = 0;
75 
76  fSandiaCof.resize(4,0.0);
77 
78  // default generator
80 }
81 
82 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
83 
85 {}
86 
87 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
88 
90  const G4DataVector&)
91 {
92  fAtomDeexcitation = G4LossTableManager::Instance()->AtomDeexcitation();
93  if(!fParticleChange) { fParticleChange = GetParticleChangeForGamma(); }
94 }
95 
96 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
97 
98 G4double
103 {
104  // This method may be used only if G4MaterialCutsCouple pointer
105  // has been set properly
106 
108  ->GetSandiaTable()->GetSandiaCofPerAtom((G4int)Z, energy, fSandiaCof);
109 
110  G4double energy2 = energy*energy;
111  G4double energy3 = energy*energy2;
112  G4double energy4 = energy2*energy2;
113 
114  return fSandiaCof[0]/energy + fSandiaCof[1]/energy2 +
115  fSandiaCof[2]/energy3 + fSandiaCof[3]/energy4;
116 }
117 
118 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
119 
120 G4double
122  const G4ParticleDefinition*,
125 {
126  const G4double* SandiaCof =
127  material->GetSandiaTable()->GetSandiaCofForMaterial(energy);
128 
129  G4double energy2 = energy*energy;
130  G4double energy3 = energy*energy2;
131  G4double energy4 = energy2*energy2;
132 
133  return SandiaCof[0]/energy + SandiaCof[1]/energy2 +
134  SandiaCof[2]/energy3 + SandiaCof[3]/energy4;
135 }
136 
137 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
138 
139 void
140 G4PEEffectFluoModel::SampleSecondaries(std::vector<G4DynamicParticle*>* fvect,
141  const G4MaterialCutsCouple* couple,
142  const G4DynamicParticle* aDynamicPhoton,
143  G4double,
144  G4double)
145 {
146  SetCurrentCouple(couple);
147  const G4Material* aMaterial = couple->GetMaterial();
148 
149  G4double energy = aDynamicPhoton->GetKineticEnergy();
150 
151  // select randomly one element constituing the material.
152  const G4Element* anElement = SelectRandomAtom(aMaterial,theGamma,energy);
153 
154  //
155  // Photo electron
156  //
157 
158  // Select atomic shell
159  G4int nShells = anElement->GetNbOfAtomicShells();
160  G4int i = 0;
161  for(; i<nShells; ++i) {
162  /*
163  G4cout << "i= " << i << " E(eV)= " << energy/eV
164  << " Eb(eV)= " << anElement->GetAtomicShell(i)/eV
165  << " " << anElement->GetName()
166  << G4endl;
167  */
168  if(energy >= anElement->GetAtomicShell(i)) { break; }
169  }
170 
171  G4double edep = energy;
172 
173  // Normally one shell is available
174  if (i < nShells) {
175 
176  G4double bindingEnergy = anElement->GetAtomicShell(i);
177  edep = bindingEnergy;
178  G4double esec = 0.0;
179 
180  // sample deexcitation
181  //
182  if(fAtomDeexcitation) {
183  G4int index = couple->GetIndex();
184  if(fAtomDeexcitation->CheckDeexcitationActiveRegion(index)) {
185  G4int Z = G4lrint(anElement->GetZ());
187  const G4AtomicShell* shell = fAtomDeexcitation->GetAtomicShell(Z, as);
188  G4double eshell = shell->BindingEnergy();
189  if(eshell > bindingEnergy && eshell <= energy) {
190  bindingEnergy = eshell;
191  edep = eshell;
192  }
193  size_t nbefore = fvect->size();
194  fAtomDeexcitation->GenerateParticles(fvect, shell, Z, index);
195  size_t nafter = fvect->size();
196  if(nafter > nbefore) {
197  for (size_t j=nbefore; j<nafter; ++j) {
198  G4double e = ((*fvect)[j])->GetKineticEnergy();
199  if(esec + e > edep) {
200  /*
201  G4cout << "### G4PEffectFluoModel Edep(eV)= " << edep/eV
202  << " Esec(eV)= " << esec/eV
203  << " E["<< j << "](eV)= " << e/eV
204  << " N= " << nafter
205  << " Z= " << Z << " shell= " << i
206  << " Ebind(keV)= " << bindingEnergy/keV
207  << " Eshell(keV)= " << shell->BindingEnergy()/keV
208  << G4endl;
209  */
210  for (size_t jj=j; jj<nafter; ++jj) { delete (*fvect)[jj]; }
211  for (size_t jj=j; jj<nafter; ++jj) { fvect->pop_back(); }
212  break;
213  }
214  esec += e;
215  }
216  }
217  edep -= esec;
218  }
219  }
220  // create photo electron
221  //
222  G4double elecKineEnergy = energy - bindingEnergy;
223  if (elecKineEnergy > fminimalEnergy) {
224  G4DynamicParticle* aParticle = new G4DynamicParticle(theElectron,
225  GetAngularDistribution()->SampleDirection(aDynamicPhoton,
226  elecKineEnergy,
227  i, couple->GetMaterial()),
228  elecKineEnergy);
229  fvect->push_back(aParticle);
230  } else {
231  edep += elecKineEnergy;
232  elecKineEnergy = 0.0;
233  }
234  if(fabs(energy - elecKineEnergy - esec - edep) > eV) {
235  G4cout << "### G4PEffectFluoModel dE(eV)= "
236  << (energy - elecKineEnergy - esec - edep)/eV
237  << " shell= " << i
238  << " E(keV)= " << energy/keV
239  << " Ebind(keV)= " << bindingEnergy/keV
240  << " Ee(keV)= " << elecKineEnergy/keV
241  << " Esec(keV)= " << esec/keV
242  << " Edep(keV)= " << edep/keV
243  << G4endl;
244  }
245  }
246 
247  // kill primary photon
248  fParticleChange->SetProposedKineticEnergy(0.);
249  fParticleChange->ProposeTrackStatus(fStopAndKill);
250  if(edep > 0.0) {
251  fParticleChange->ProposeLocalEnergyDeposit(edep);
252  }
253 }
254 
255 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
G4bool CheckDeexcitationActiveRegion(G4int coupleIndex)
static G4LossTableManager * Instance()
G4int GetNbOfAtomicShells() const
Definition: G4Element.hh:146
G4double GetKineticEnergy() const
Int_t index
Definition: macro.C:9
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)
G4SandiaTable * GetSandiaTable() const
Definition: G4Material.hh:227
string material
Definition: eplot.py:19
virtual const G4AtomicShell * GetAtomicShell(G4int Z, G4AtomicShellEnumerator shell)=0
const G4MaterialCutsCouple * CurrentCouple() const
Definition: G4VEmModel.hh:426
G4GLOB_DLL std::ostream G4cout
double energy
Definition: plottest35.C:25
Float_t Z
Definition: plot.C:39
G4double GetSandiaCofForMaterial(G4int, G4int)
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A, G4double, G4double)
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
void GetSandiaCofPerAtom(G4int Z, G4double energy, std::vector< G4double > &coeff)
int G4lrint(double ad)
Definition: templates.hh:163
void SetCurrentCouple(const G4MaterialCutsCouple *)
Definition: G4VEmModel.hh:419
void SetAngularDistribution(G4VEmAngularDistribution *)
Definition: G4VEmModel.hh:585
virtual void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)
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)
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &)
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)
G4PEEffectFluoModel(const G4String &nam="PhotoElectric")
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