Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
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: G4LivermorePhotoElectricModel.cc 90578 2015-06-04 09:55:07Z gcosmo $
27 //
28 //
29 // Author: Sebastien Incerti
30 // 30 October 2008
31 // on base of G4LowEnergyPhotoElectric developed by A.Forti and M.G.Pia
32 //
33 // 22 Oct 2012 A & V Ivanchenko Migration data structure to G4PhysicsVector
34 //
35 
37 #include "G4SystemOfUnits.hh"
38 #include "G4PhysicalConstants.hh"
39 #include "G4LossTableManager.hh"
40 #include "G4Electron.hh"
41 #include "G4Gamma.hh"
43 #include "G4CrossSectionHandler.hh"
44 #include "G4LPhysicsFreeVector.hh"
45 #include "G4VAtomDeexcitation.hh"
48 #include "G4AtomicShell.hh"
49 
50 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
51 
52 G4LPhysicsFreeVector* G4LivermorePolarizedPhotoElectricModel::fCrossSection[] = {nullptr};
53 G4LPhysicsFreeVector* G4LivermorePolarizedPhotoElectricModel::fCrossSectionLE[] = {nullptr};
54 std::vector<G4double>* G4LivermorePolarizedPhotoElectricModel::fParam[] = {nullptr};
55 G4int G4LivermorePolarizedPhotoElectricModel::fNShells[] = {0};
56 G4int G4LivermorePolarizedPhotoElectricModel::fNShellsUsed[] = {0};
57 G4ElementData* G4LivermorePolarizedPhotoElectricModel::fShellCrossSection = nullptr;
58 G4Material* G4LivermorePolarizedPhotoElectricModel::fWater = nullptr;
59 G4double G4LivermorePolarizedPhotoElectricModel::fWaterEnergyLimit = 0.0;
60 
61 using namespace std;
62 
63 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
64 
66  : G4VEmModel(nam),fParticleChange(nullptr),maxZ(99),
67  nShellLimit(100),fDeexcitationActive(false),isInitialised(false),
68  fAtomDeexcitation(nullptr)
69 {
70  verboseLevel= 0;
71  // Verbosity scale:
72  // 0 = nothing
73  // 1 = warning for energy non-conservation
74  // 2 = details of energy budget
75  // 3 = calculation of cross sections, file openings, sampling of atoms
76  // 4 = entering in methods
77 
78  theGamma = G4Gamma::Gamma();
79  theElectron = G4Electron::Electron();
80 
81  // default generator
82  // SetAngularDistribution(new G4SauterGavrilaAngularDistribution());
84 
85  if(verboseLevel>0) {
86  G4cout << "Livermore PhotoElectric is constructed "
87  << " nShellLimit= " << nShellLimit << G4endl;
88  }
89 
90  //Mark this model as "applicable" for atomic deexcitation
91  SetDeexcitationFlag(true);
92  fSandiaCof.resize(4,0.0);
93  fCurrSection = 0.0;
94 
95 }
96 
97 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
98 
100 {
101  if(IsMaster()) {
102  delete fShellCrossSection;
103  for(G4int i=0; i<maxZ; ++i) {
104  delete fParam[i];
105  fParam[i] = 0;
106  delete fCrossSection[i];
107  fCrossSection[i] = 0;
108  delete fCrossSectionLE[i];
109  fCrossSectionLE[i] = 0;
110  }
111  }
112 }
113 
114 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
115 
116 void
118  const G4DataVector&)
119 {
120  if (verboseLevel > 2) {
121  G4cout << "Calling G4LivermorePolarizedPhotoElectricModel::Initialise()" << G4endl;
122  }
123 
124  if(IsMaster()) {
125 
126  if(!fWater) {
127  fWater = G4Material::GetMaterial("G4_WATER", false);
128  if(fWater) { fWaterEnergyLimit = 13.6*eV; }
129  }
130 
131  if(!fShellCrossSection) { fShellCrossSection = new G4ElementData(); }
132 
133  char* path = getenv("G4LEDATA");
134 
135  G4ProductionCutsTable* theCoupleTable =
137  G4int numOfCouples = theCoupleTable->GetTableSize();
138 
139  for(G4int i=0; i<numOfCouples; ++i) {
140  const G4MaterialCutsCouple* couple =
141  theCoupleTable->GetMaterialCutsCouple(i);
142  const G4Material* material = couple->GetMaterial();
143  const G4ElementVector* theElementVector = material->GetElementVector();
144  G4int nelm = material->GetNumberOfElements();
145 
146  for (G4int j=0; j<nelm; ++j) {
147  G4int Z = (G4int)(*theElementVector)[j]->GetZ();
148  if(Z < 1) { Z = 1; }
149  else if(Z > maxZ) { Z = maxZ; }
150  if(!fCrossSection[Z]) { ReadData(Z, path); }
151  }
152  }
153  }
154  //
155  if (verboseLevel > 2) {
156  G4cout << "Loaded cross section files for LivermorePolarizedPhotoElectric model"
157  << G4endl;
158  }
159  if(!isInitialised) {
160  isInitialised = true;
162 
163  fAtomDeexcitation = G4LossTableManager::Instance()->AtomDeexcitation();
164  }
165  fDeexcitationActive = false;
166  if(fAtomDeexcitation) {
167  fDeexcitationActive = fAtomDeexcitation->IsFluoActive();
168  }
169 
170  if (verboseLevel > 0) {
171  G4cout << "LivermorePolarizedPhotoElectric model is initialized " << G4endl
172  << G4endl;
173  }
174 }
175 
176 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
177 
179  const G4Material* material,
180  const G4ParticleDefinition* p,
183 {
184  fCurrSection = 0.0;
185  if(fWater && (material == fWater ||
186  material->GetBaseMaterial() == fWater)) {
187  if(energy <= fWaterEnergyLimit) {
188  fWater->GetSandiaTable()->GetSandiaCofWater(energy, fSandiaCof);
189 
190  G4double energy2 = energy*energy;
191  G4double energy3 = energy*energy2;
192  G4double energy4 = energy2*energy2;
193 
194  fCurrSection = material->GetDensity()*
195  (fSandiaCof[0]/energy + fSandiaCof[1]/energy2 +
196  fSandiaCof[2]/energy3 + fSandiaCof[3]/energy4);
197  }
198  }
199  if(0.0 == fCurrSection) {
200  fCurrSection = G4VEmModel::CrossSectionPerVolume(material, p, energy);
201  }
202  return fCurrSection;
203 }
204 
205 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
206 
208  const G4ParticleDefinition*,
210  G4double ZZ, G4double,
212 {
213  if (verboseLevel > 3) {
214  G4cout << "G4LivermorePolarizedPhotoElectricModel::ComputeCrossSectionPerAtom():"
215  << " Z= " << ZZ << " R(keV)= " << energy/keV << G4endl;
216  }
217  G4double cs = 0.0;
218  G4int Z = G4lrint(ZZ);
219  if(Z < 1 || Z >= maxZ) { return cs; }
220 
221  // if element was not initialised
222  // do initialisation safely for MT mode
223  if(!fCrossSection[Z]) {
224  InitialiseForElement(0, Z);
225  if(!fCrossSection[Z]) { return cs; }
226  }
227 
228  G4int idx = fNShells[Z]*6 - 4;
229  if (energy < (*(fParam[Z]))[idx-1]) { energy = (*(fParam[Z]))[idx-1]; }
230 
231  G4double x1 = 1.0/energy;
232  G4double x2 = x1*x1;
233  G4double x3 = x2*x1;
234 
235  // parameterisation
236  if(energy >= (*(fParam[Z]))[0]) {
237  G4double x4 = x2*x2;
238  cs = x1*((*(fParam[Z]))[idx] + x1*(*(fParam[Z]))[idx+1]
239  + x2*(*(fParam[Z]))[idx+2] + x3*(*(fParam[Z]))[idx+3]
240  + x4*(*(fParam[Z]))[idx+4]);
241  // high energy part
242  } else if(energy >= (*(fParam[Z]))[1]) {
243  cs = x3*(fCrossSection[Z])->Value(energy);
244 
245  // low energy part
246  } else {
247  cs = x3*(fCrossSectionLE[Z])->Value(energy);
248  }
249  if (verboseLevel > 1) {
250  G4cout << "LivermorePolarizedPhotoElectricModel: E(keV)= " << energy/keV
251  << " Z= " << Z << " cross(barn)= " << cs/barn << G4endl;
252  }
253  return cs;
254 }
255 
256 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
257 
258 void
260  std::vector<G4DynamicParticle*>* fvect,
261  const G4MaterialCutsCouple* couple,
262  const G4DynamicParticle* aDynamicGamma,
263  G4double,
264  G4double)
265 {
266  G4double gammaEnergy = aDynamicGamma->GetKineticEnergy();
267  if (verboseLevel > 3) {
268  G4cout << "G4LivermorePolarizedPhotoElectricModel::SampleSecondaries() Egamma(keV)= "
269  << gammaEnergy/keV << G4endl;
270  }
271 
272  // kill incident photon
275 
276  // low-energy photo-effect in water - full absorption
277  const G4Material* material = couple->GetMaterial();
278  if(fWater && (material == fWater ||
279  material->GetBaseMaterial() == fWater)) {
280  if(gammaEnergy <= fWaterEnergyLimit) {
282  return;
283  }
284  }
285 
286  // Returns the normalized direction of the momentum
287  G4ThreeVector photonDirection = aDynamicGamma->GetMomentumDirection();
288 
289  // Select randomly one element in the current material
290  //G4cout << "Select random atom Egamma(keV)= " << gammaEnergy/keV << G4endl;
291  const G4Element* elm = SelectRandomAtom(material, theGamma, gammaEnergy);
292  G4int Z = G4lrint(elm->GetZ());
293 
294  // Select the ionised shell in the current atom according to shell
295  // cross sections
296  // G4cout << "Select random shell Z= " << Z << G4endl;
297 
298  if(Z >= maxZ) { Z = maxZ-1; }
299 
300  // element was not initialised gamma should be absorbed
301  if(!fCrossSection[Z]) {
303  return;
304  }
305 
306  // shell index
307  size_t shellIdx = 0;
308  size_t nn = fNShellsUsed[Z];
309 
310  if(nn > 1) {
311  if(gammaEnergy >= (*(fParam[Z]))[0]) {
312  G4double x1 = 1.0/gammaEnergy;
313  G4double x2 = x1*x1;
314  G4double x3 = x2*x1;
315  G4double x4 = x3*x1;
316  G4int idx = nn*6 - 4;
317  // when do sampling common factors are not taken into account
318  // so cross section is not real
319  G4double cs0 = G4UniformRand()*((*(fParam[Z]))[idx]
320  + x1*(*(fParam[Z]))[idx+1]
321  + x2*(*(fParam[Z]))[idx+2]
322  + x3*(*(fParam[Z]))[idx+3]
323  + x4*(*(fParam[Z]))[idx+4]);
324  for(shellIdx=0; shellIdx<nn; ++shellIdx) {
325  idx = shellIdx*6 + 2;
326  if(gammaEnergy > (*(fParam[Z]))[idx-1]) {
327  G4double cs = (*(fParam[Z]))[idx] + x1*(*(fParam[Z]))[idx+1]
328  + x2*(*(fParam[Z]))[idx+2] + x3*(*(fParam[Z]))[idx+3]
329  + x4*(*(fParam[Z]))[idx+4];
330  if(cs >= cs0) { break; }
331  }
332  }
333  if(shellIdx >= nn) { shellIdx = nn-1; }
334 
335  } else {
336 
337  // when do sampling common factors are not taken into account
338  // so cross section is not real
339  G4double cs = G4UniformRand();
340 
341  if(gammaEnergy >= (*(fParam[Z]))[1]) {
342  cs *= (fCrossSection[Z])->Value(gammaEnergy);
343  } else {
344  cs *= (fCrossSectionLE[Z])->Value(gammaEnergy);
345  }
346 
347  for(size_t j=0; j<nn; ++j) {
348  shellIdx = (size_t)fShellCrossSection->GetComponentID(Z, j);
349  if(gammaEnergy > (*(fParam[Z]))[6*shellIdx+1]) {
350  cs -= fShellCrossSection->GetValueForComponent(Z, j, gammaEnergy);
351  }
352  if(cs <= 0.0 || j+1 == nn) { break; }
353  }
354  }
355  }
356 
357  G4double bindingEnergy = (*(fParam[Z]))[shellIdx*6 + 1];
358  //G4cout << "Z= " << Z << " shellIdx= " << shellIdx
359  // << " nShells= " << fNShells[Z]
360  // << " Ebind(keV)= " << bindingEnergy/keV
361  // << " Egamma(keV)= " << gammaEnergy/keV << G4endl;
362 
363  const G4AtomicShell* shell = 0;
364 
365  // no de-excitation from the last shell
366  if(fDeexcitationActive && shellIdx + 1 < nn) {
368  shell = fAtomDeexcitation->GetAtomicShell(Z, as);
369  }
370 
371  // If binding energy of the selected shell is larger than photon energy
372  // do not generate secondaries
373  if(gammaEnergy < bindingEnergy) {
375  return;
376  }
377 
378  // Primary outcoming electron
379  G4double eKineticEnergy = gammaEnergy - bindingEnergy;
380  G4double edep = bindingEnergy;
381 
382  // Calculate direction of the photoelectron
383  G4ThreeVector electronDirection =
384  GetAngularDistribution()->SampleDirection(aDynamicGamma,
385  eKineticEnergy,
386  shellIdx,
387  couple->GetMaterial());
388 
389  // The electron is created
390  G4DynamicParticle* electron = new G4DynamicParticle (theElectron,
391  electronDirection,
392  eKineticEnergy);
393  fvect->push_back(electron);
394 
395  // Sample deexcitation
396  if(shell) {
397  G4int index = couple->GetIndex();
398  if(fAtomDeexcitation->CheckDeexcitationActiveRegion(index)) {
399  G4int nbefore = fvect->size();
400 
401  fAtomDeexcitation->GenerateParticles(fvect, shell, Z, index);
402  G4int nafter = fvect->size();
403  if(nafter > nbefore) {
404  G4double esec = 0.0;
405  for (G4int j=nbefore; j<nafter; ++j) {
406 
407  G4double e = ((*fvect)[j])->GetKineticEnergy();
408  if(esec + e > edep) {
409  // correct energy in order to have energy balance
410  e = edep - esec;
411  ((*fvect)[j])->SetKineticEnergy(e);
412  esec += e;
413  // delete the rest of secondaries (should not happens)
414  for (G4int jj=nafter-1; jj>j; --jj) {
415  delete (*fvect)[jj];
416  fvect->pop_back();
417  }
418  break;
419  }
420  esec += e;
421  }
422  edep -= esec;
423  }
424  }
425  }
426  // energy balance - excitation energy left
427  if(edep > 0.0) {
429  }
430 }
431 
432 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
433 
434 void
435 G4LivermorePolarizedPhotoElectricModel::ReadData(G4int Z, const char* path)
436 {
437  if (verboseLevel > 1)
438  {
439  G4cout << "Calling ReadData() of G4LivermoreGammaConversionModel"
440  << G4endl;
441  }
442 
443  if(fCrossSection[Z]) { return; }
444 
445  const char* datadir = path;
446 
447  if(!datadir)
448  {
449  datadir = getenv("G4LEDATA");
450  if(!datadir)
451  {
452  G4Exception("G4LivermorePolarizedPhotoElectricModel::ReadData()",
453  "em0006",FatalException,
454  "Environment variable G4LEDATA not defined");
455  return;
456  }
457  }
458 
459  // spline for photoeffect total x-section above K-shell
460  fCrossSection[Z] = new G4LPhysicsFreeVector();
461  fCrossSection[Z]->SetSpline(true);
462 
463  std::ostringstream ost;
464  ost << datadir << "/livermore/phot/pe-cs-" << Z <<".dat";
465  std::ifstream fin(ost.str().c_str());
466  if( !fin.is_open()) {
468  ed << "G4LivermorePolarizedPhotoElectricModel data file <" << ost.str().c_str()
469  << "> is not opened!" << G4endl;
470  G4Exception("G4LivermorePolarizedPhotoElectricModel::ReadData()",
471  "em0003",FatalException,
472  ed,"G4LEDATA version should be G4EMLOW6.32 or later.");
473  return;
474  } else {
475  if(verboseLevel > 3) { G4cout << "File " << ost.str().c_str()
476  << " is opened by G4LivermorePolarizedPhotoElectricModel" << G4endl;}
477  fCrossSection[Z]->Retrieve(fin, true);
478  fCrossSection[Z]->ScaleVector(MeV, barn);
479  fin.close();
480  }
481 
482  fParam[Z] = new std::vector<G4double>;
483 
484  // read fit parameters
485  G4int n1 = 0;
486  G4int n2 = 0;
487  G4double x;
488  std::ostringstream ost1;
489  ost1 << datadir << "/livermore/phot/pe-" << Z <<".dat";
490  std::ifstream fin1(ost1.str().c_str());
491  if( !fin1.is_open()) {
493  ed << "G4LivermorePolarizedPhotoElectricModel data file <" << ost1.str().c_str()
494  << "> is not opened!" << G4endl;
495  G4Exception("G4LivermorePolarizedPhotoElectricModel::ReadData()",
496  "em0003",FatalException,
497  ed,"G4LEDATA version should be G4EMLOW6.32 or later.");
498  return;
499  } else {
500  if(verboseLevel > 3) {
501  G4cout << "File " << ost1.str().c_str()
502  << " is opened by G4LivermorePolarizedPhotoElectricModel" << G4endl;
503  }
504  fin1 >> n1;
505  if(fin1.fail()) { return; }
506  if(0 > n1 || n1 >= INT_MAX) { n1 = 0; }
507 
508  fin1 >> n2;
509  if(fin1.fail()) { return; }
510  if(0 > n2 || n2 >= INT_MAX) { n2 = 0; }
511 
512  fin1 >> x;
513  if(fin1.fail()) { return; }
514 
515  fNShells[Z] = n1;
516  fParam[Z]->reserve(6*n1+1);
517  fParam[Z]->push_back(x*MeV);
518  for(G4int i=0; i<n1; ++i) {
519  for(G4int j=0; j<6; ++j) {
520  fin1 >> x;
521  if(0 == j) { x *= MeV; }
522  else { x *= barn; }
523  fParam[Z]->push_back(x);
524  }
525  }
526  fin1.close();
527  }
528  // there is a possibility to used only main shells
529  if(nShellLimit < n2) { n2 = nShellLimit; }
530  fShellCrossSection->InitialiseForComponent(Z, n2);
531  fNShellsUsed[Z] = n2;
532 
533  if(1 < n2) {
534  std::ostringstream ost2;
535  ost2 << datadir << "/livermore/phot/pe-ss-cs-" << Z <<".dat";
536  std::ifstream fin2(ost2.str().c_str());
537  if( !fin2.is_open()) {
539  ed << "G4LivermorePolarizedPhotoElectricModel data file <" << ost2.str().c_str()
540  << "> is not opened!" << G4endl;
541  G4Exception("G4LivermorePolarizedPhotoElectricModel::ReadData()",
542  "em0003",FatalException,
543  ed,"G4LEDATA version should be G4EMLOW6.32 or later.");
544  return;
545  } else {
546  if(verboseLevel > 3) {
547  G4cout << "File " << ost2.str().c_str()
548  << " is opened by G4LivermorePolarizedPhotoElectricModel" << G4endl;
549  }
550 
551  G4int n3, n4;
552  G4double y;
553  for(G4int i=0; i<n2; ++i) {
554  fin2 >> x >> y >> n3 >> n4;
556  for(G4int j=0; j<n3; ++j) {
557  fin2 >> x >> y;
558  v->PutValues(j, x*MeV, y*barn);
559  }
560  fShellCrossSection->AddComponent(Z, n4, v);
561  }
562  fin2.close();
563  }
564  }
565 
566  // no spline for photoeffect total x-section below K-shell
567  if(1 < fNShells[Z]) {
568  fCrossSectionLE[Z] = new G4LPhysicsFreeVector();
569 
570  std::ostringstream ost3;
571  ost3 << datadir << "/livermore/phot/pe-le-cs-" << Z <<".dat";
572  std::ifstream fin3(ost3.str().c_str());
573  if( !fin3.is_open()) {
575  ed << "G4LivermorePolarizedPhotoElectricModel data file <" << ost3.str().c_str()
576  << "> is not opened!" << G4endl;
577  G4Exception("G4LivermorePolarizedPhotoElectricModel::ReadData()",
578  "em0003",FatalException,
579  ed,"G4LEDATA version should be G4EMLOW6.32 or later.");
580  return;
581  } else {
582  if(verboseLevel > 3) {
583  G4cout << "File " << ost3.str().c_str()
584  << " is opened by G4LivermorePolarizedPhotoElectricModel" << G4endl;
585  }
586  fCrossSectionLE[Z]->Retrieve(fin3, true);
587  fCrossSectionLE[Z]->ScaleVector(MeV, barn);
588  fin3.close();
589  }
590  }
591 }
592 
593 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
594 
595 #include "G4AutoLock.hh"
596 namespace { G4Mutex LivermorePolarizedPhotoElectricModelMutex = G4MUTEX_INITIALIZER; }
597 
599  const G4ParticleDefinition*, G4int Z)
600 {
601  G4AutoLock l(&LivermorePolarizedPhotoElectricModelMutex);
602  // G4cout << "G4LivermorePolarizedPhotoElectricModel::InitialiseForElement Z= "
603  // << Z << G4endl;
604  if(!fCrossSection[Z]) { ReadData(Z); }
605  l.unlock();
606 }
607 
608 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:257
G4bool CheckDeexcitationActiveRegion(G4int coupleIndex)
static G4LossTableManager * Instance()
G4bool IsFluoActive() const
void AddComponent(G4int Z, G4int id, G4PhysicsVector *v)
std::vector< G4Element * > G4ElementVector
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4double GetKineticEnergy() const
virtual void InitialiseForElement(const G4ParticleDefinition *, G4int Z)
static G4Material * GetMaterial(const G4String &name, G4bool warning=true)
Definition: G4Material.cc:602
const char * p
Definition: xmltok.h:285
G4double GetZ() const
Definition: G4Element.hh:131
G4VEmAngularDistribution * GetAngularDistribution()
Definition: G4VEmModel.hh:619
G4double GetDensity() const
Definition: G4Material.hh:180
G4bool IsMaster() const
Definition: G4VEmModel.hh:725
tuple x
Definition: test.py:50
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:190
int G4int
Definition: G4Types.hh:78
#define G4MUTEX_INITIALIZER
Definition: G4Threading.hh:175
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double energy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
virtual void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)
void SetSpline(G4bool)
G4SandiaTable * GetSandiaTable() const
Definition: G4Material.hh:229
string material
Definition: eplot.py:19
virtual const G4AtomicShell * GetAtomicShell(G4int Z, G4AtomicShellEnumerator shell)=0
G4int GetComponentID(G4int Z, size_t idx)
G4LivermorePolarizedPhotoElectricModel(const G4String &nam="LivermorePolarizedPhotoElectric")
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
virtual G4ThreeVector & SampleDirection(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, const G4Material *)=0
const G4ThreeVector & GetMomentumDirection() const
static constexpr double eV
Definition: G4SIunits.hh:215
virtual void ScaleVector(G4double factorE, G4double factorV)
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
tuple v
Definition: test.py:18
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
#define INT_MAX
Definition: templates.hh:111
static G4ProductionCutsTable * GetProductionCutsTable()
G4int G4Mutex
Definition: G4Threading.hh:173
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
int G4lrint(double ad)
Definition: templates.hh:163
G4double energy(const ThreeVector &p, const G4double m)
void SetAngularDistribution(G4VEmAngularDistribution *)
Definition: G4VEmModel.hh:626
virtual G4bool Retrieve(std::ifstream &fIn, G4bool ascii=false)
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &)
G4double GetValueForComponent(G4int Z, size_t idx, G4double kinEnergy)
const G4Material * GetBaseMaterial() const
Definition: G4Material.hh:233
static G4Electron * Electron()
Definition: G4Electron.cc:94
void SetProposedKineticEnergy(G4double proposedKinEnergy)
void PutValues(size_t index, G4double e, G4double dataValue)
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
size_t GetNumberOfElements() const
Definition: G4Material.hh:186
G4VAtomDeexcitation * AtomDeexcitation()
void GetSandiaCofWater(G4double energy, std::vector< G4double > &coeff) const
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double energy, G4double Z, G4double A=0, G4double cut=0, G4double emax=DBL_MAX)
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:788
static constexpr double barn
Definition: G4SIunits.hh:105
G4double bindingEnergy(G4int A, G4int Z)
void InitialiseForComponent(G4int Z, G4int nComponents=0)
static constexpr double keV
Definition: G4SIunits.hh:216
G4AtomicShellEnumerator
const G4Element * SelectRandomAtom(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.hh:544
G4ParticleChangeForGamma * GetParticleChangeForGamma()
Definition: G4VEmModel.cc:132
virtual G4double Value(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy)
Definition: G4VEmModel.cc:377
const G4Material * GetMaterial() const