Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4RDVCrossSectionHandler Class Referenceabstract

#include <G4RDVCrossSectionHandler.hh>

Inheritance diagram for G4RDVCrossSectionHandler:

Public Member Functions

 G4RDVCrossSectionHandler ()
 
 G4RDVCrossSectionHandler (G4RDVDataSetAlgorithm *interpolation, G4double minE=250 *CLHEP::eV, G4double maxE=100 *CLHEP::GeV, G4int nBins=200, G4double unitE=CLHEP::MeV, G4double unitData=CLHEP::barn, G4int minZ=1, G4int maxZ=99)
 
virtual ~G4RDVCrossSectionHandler ()
 
void Initialise (G4RDVDataSetAlgorithm *interpolation=0, G4double minE=250 *CLHEP::eV, G4double maxE=100 *CLHEP::GeV, G4int numberOfBins=200, G4double unitE=CLHEP::MeV, G4double unitData=CLHEP::barn, G4int minZ=1, G4int maxZ=99)
 
G4int SelectRandomAtom (const G4MaterialCutsCouple *couple, G4double e) const
 
const G4ElementSelectRandomElement (const G4MaterialCutsCouple *material, G4double e) const
 
G4int SelectRandomShell (G4int Z, G4double e) const
 
G4RDVEMDataSetBuildMeanFreePathForMaterials (const G4DataVector *energyCuts=0)
 
G4double FindValue (G4int Z, G4double e) const
 
G4double FindValue (G4int Z, G4double e, G4int shellIndex) const
 
G4double ValueForMaterial (const G4Material *material, G4double e) const
 
void LoadData (const G4String &dataFile)
 
void LoadShellData (const G4String &dataFile)
 
void PrintData () const
 
void Clear ()
 

Protected Member Functions

G4int NumberOfComponents (G4int Z) const
 
void ActiveElements ()
 
virtual std::vector
< G4RDVEMDataSet * > * 
BuildCrossSectionsForMaterials (const G4DataVector &energyVector, const G4DataVector *energyCuts=0)=0
 
virtual G4RDVDataSetAlgorithmCreateInterpolation ()
 
const G4RDVDataSetAlgorithmGetInterpolation () const
 

Detailed Description

Definition at line 64 of file G4RDVCrossSectionHandler.hh.

Constructor & Destructor Documentation

G4RDVCrossSectionHandler::G4RDVCrossSectionHandler ( )

Definition at line 59 of file G4RDVCrossSectionHandler.cc.

60 {
61  crossSections = 0;
62  interpolation = 0;
63  Initialise();
65 }
void Initialise(G4RDVDataSetAlgorithm *interpolation=0, G4double minE=250 *CLHEP::eV, G4double maxE=100 *CLHEP::GeV, G4int numberOfBins=200, G4double unitE=CLHEP::MeV, G4double unitData=CLHEP::barn, G4int minZ=1, G4int maxZ=99)

Here is the call graph for this function:

G4RDVCrossSectionHandler::G4RDVCrossSectionHandler ( G4RDVDataSetAlgorithm interpolation,
G4double  minE = 250*CLHEP::eV,
G4double  maxE = 100*CLHEP::GeV,
G4int  nBins = 200,
G4double  unitE = CLHEP::MeV,
G4double  unitData = CLHEP::barn,
G4int  minZ = 1,
G4int  maxZ = 99 
)

Definition at line 68 of file G4RDVCrossSectionHandler.cc.

76  : interpolation(algorithm), eMin(minE), eMax(maxE), nBins(bins),
77  unit1(unitE), unit2(unitData), zMin(minZ), zMax(maxZ)
78 {
79  crossSections = 0;
81 }

Here is the call graph for this function:

G4RDVCrossSectionHandler::~G4RDVCrossSectionHandler ( )
virtual

Definition at line 83 of file G4RDVCrossSectionHandler.cc.

84 {
85  delete interpolation;
86  interpolation = 0;
87  std::map<G4int,G4RDVEMDataSet*,std::less<G4int> >::iterator pos;
88 
89  for (pos = dataMap.begin(); pos != dataMap.end(); ++pos)
90  {
91  // The following is a workaround for STL ObjectSpace implementation,
92  // which does not support the standard and does not accept
93  // the syntax pos->second
94  // G4RDVEMDataSet* dataSet = pos->second;
95  G4RDVEMDataSet* dataSet = (*pos).second;
96  delete dataSet;
97  }
98 
99  if (crossSections != 0)
100  {
101  size_t n = crossSections->size();
102  for (size_t i=0; i<n; i++)
103  {
104  delete (*crossSections)[i];
105  }
106  delete crossSections;
107  crossSections = 0;
108  }
109 }
const G4int n
static const G4double pos

Member Function Documentation

void G4RDVCrossSectionHandler::ActiveElements ( )
protected

Definition at line 624 of file G4RDVCrossSectionHandler.cc.

625 {
626  const G4MaterialTable* materialTable = G4Material::GetMaterialTable();
627  if (materialTable == 0)
628  G4Exception("G4RDVCrossSectionHandler::ActiveElements",
629  "InvalidSetup", FatalException, "No MaterialTable found!");
630 
632 
633  for (G4int m=0; m<nMaterials; m++)
634  {
635  const G4Material* material= (*materialTable)[m];
636  const G4ElementVector* elementVector = material->GetElementVector();
637  const G4int nElements = material->GetNumberOfElements();
638 
639  for (G4int iEl=0; iEl<nElements; iEl++)
640  {
641  G4Element* element = (*elementVector)[iEl];
642  G4double Z = element->GetZ();
643  if (!(activeZ.contains(Z)) && Z >= zMin && Z <= zMax)
644  {
645  activeZ.push_back(Z);
646  }
647  }
648  }
649 }
std::vector< G4Element * > G4ElementVector
G4double GetZ() const
Definition: G4Element.hh:131
static G4MaterialTable * GetMaterialTable()
Definition: G4Material.cc:587
std::vector< G4Material * > G4MaterialTable
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:190
int G4int
Definition: G4Types.hh:78
string material
Definition: eplot.py:19
static constexpr double m
Definition: G4SIunits.hh:129
static size_t GetNumberOfMaterials()
Definition: G4Material.cc:594
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
size_t GetNumberOfElements() const
Definition: G4Material.hh:186
double G4double
Definition: G4Types.hh:76
G4bool contains(const G4double &) const

Here is the call graph for this function:

Here is the caller graph for this function:

virtual std::vector<G4RDVEMDataSet*>* G4RDVCrossSectionHandler::BuildCrossSectionsForMaterials ( const G4DataVector energyVector,
const G4DataVector energyCuts = 0 
)
protectedpure virtual

Implemented in G4RDeIonisationCrossSectionHandler, G4RDBremsstrahlungCrossSectionHandler, and G4RDCrossSectionHandler.

Here is the caller graph for this function:

G4RDVEMDataSet * G4RDVCrossSectionHandler::BuildMeanFreePathForMaterials ( const G4DataVector energyCuts = 0)

Definition at line 416 of file G4RDVCrossSectionHandler.cc.

417 {
418  // Builds a CompositeDataSet containing the mean free path for each material
419  // in the material table
420 
421  G4DataVector energyVector;
422  G4double dBin = std::log10(eMax/eMin) / nBins;
423 
424  for (G4int i=0; i<nBins+1; i++)
425  {
426  energyVector.push_back(std::pow(10., std::log10(eMin)+i*dBin));
427  }
428 
429  // Factory method to build cross sections in derived classes,
430  // related to the type of physics process
431 
432  if (crossSections != 0)
433  { // Reset the list of cross sections
434  std::vector<G4RDVEMDataSet*>::iterator mat;
435  if (! crossSections->empty())
436  {
437  for (mat = crossSections->begin(); mat!= crossSections->end(); ++mat)
438  {
439  G4RDVEMDataSet* set = *mat;
440  delete set;
441  set = 0;
442  }
443  crossSections->clear();
444  delete crossSections;
445  crossSections = 0;
446  }
447  }
448 
449  crossSections = BuildCrossSectionsForMaterials(energyVector,energyCuts);
450 
451  if (crossSections == 0)
452  G4Exception("G4RDVCrossSectionHandler::BuildMeanFreePathForMaterials()",
453  "InvalidCondition", FatalException, "CrossSections = 0!");
454 
456  G4RDVEMDataSet* materialSet = new G4RDCompositeEMDataSet(algo);
457 
458  G4DataVector* energies;
460 
461  const G4ProductionCutsTable* theCoupleTable=
463  size_t numOfCouples = theCoupleTable->GetTableSize();
464 
465 
466  for (size_t m=0; m<numOfCouples; m++)
467  {
468  energies = new G4DataVector;
469  data = new G4DataVector;
470  for (G4int bin=0; bin<nBins; bin++)
471  {
472  G4double energy = energyVector[bin];
473  energies->push_back(energy);
474  G4RDVEMDataSet* matCrossSet = (*crossSections)[m];
475  G4double materialCrossSection = 0.0;
476  G4int nElm = matCrossSet->NumberOfComponents();
477  for(G4int j=0; j<nElm; j++) {
478  materialCrossSection += matCrossSet->GetComponent(j)->FindValue(energy);
479  }
480 
481  if (materialCrossSection > 0.)
482  {
483  data->push_back(1./materialCrossSection);
484  }
485  else
486  {
487  data->push_back(DBL_MAX);
488  }
489  }
491  G4RDVEMDataSet* dataSet = new G4RDEMDataSet(m,energies,data,algo,1.,1.);
492  materialSet->AddComponent(dataSet);
493  }
494 
495  return materialSet;
496 }
tuple bin
Definition: plottest35.py:22
virtual const G4RDVEMDataSet * GetComponent(G4int componentId) const =0
int G4int
Definition: G4Types.hh:78
virtual size_t NumberOfComponents(void) const =0
const XML_Char const XML_Char * data
Definition: expat.h:268
static constexpr double m
Definition: G4SIunits.hh:129
virtual G4RDVDataSetAlgorithm * CreateInterpolation()
virtual G4double FindValue(G4double x, G4int componentId=0) const =0
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
static G4ProductionCutsTable * GetProductionCutsTable()
virtual std::vector< G4RDVEMDataSet * > * BuildCrossSectionsForMaterials(const G4DataVector &energyVector, const G4DataVector *energyCuts=0)=0
G4double energy(const ThreeVector &p, const G4double m)
double G4double
Definition: G4Types.hh:76
virtual void AddComponent(G4RDVEMDataSet *dataSet)=0
#define DBL_MAX
Definition: templates.hh:83

Here is the call graph for this function:

Here is the caller graph for this function:

void G4RDVCrossSectionHandler::Clear ( )

Definition at line 306 of file G4RDVCrossSectionHandler.cc.

307 {
308  // Reset the map of data sets: remove the data sets from the map
309  std::map<G4int,G4RDVEMDataSet*,std::less<G4int> >::iterator pos;
310 
311  if(! dataMap.empty())
312  {
313  for (pos = dataMap.begin(); pos != dataMap.end(); ++pos)
314  {
315  // The following is a workaround for STL ObjectSpace implementation,
316  // which does not support the standard and does not accept
317  // the syntax pos->first or pos->second
318  // G4RDVEMDataSet* dataSet = pos->second;
319  G4RDVEMDataSet* dataSet = (*pos).second;
320  delete dataSet;
321  dataSet = 0;
322  G4int i = (*pos).first;
323  dataMap[i] = 0;
324  }
325  dataMap.clear();
326  }
327 
328  activeZ.clear();
329  ActiveElements();
330 }
int G4int
Definition: G4Types.hh:78
static const G4double pos

Here is the call graph for this function:

Here is the caller graph for this function:

G4RDVDataSetAlgorithm * G4RDVCrossSectionHandler::CreateInterpolation ( )
protectedvirtual

Definition at line 651 of file G4RDVCrossSectionHandler.cc.

Here is the caller graph for this function:

G4double G4RDVCrossSectionHandler::FindValue ( G4int  Z,
G4double  e 
) const

Definition at line 332 of file G4RDVCrossSectionHandler.cc.

333 {
334  G4double value = 0.;
335 
336  std::map<G4int,G4RDVEMDataSet*,std::less<G4int> >::const_iterator pos;
337  pos = dataMap.find(Z);
338  if (pos!= dataMap.end())
339  {
340  // The following is a workaround for STL ObjectSpace implementation,
341  // which does not support the standard and does not accept
342  // the syntax pos->first or pos->second
343  // G4RDVEMDataSet* dataSet = pos->second;
344  G4RDVEMDataSet* dataSet = (*pos).second;
345  value = dataSet->FindValue(energy);
346  }
347  else
348  {
349  G4cout << "WARNING: G4RDVCrossSectionHandler::FindValue did not find Z = "
350  << Z << G4endl;
351  }
352  return value;
353 }
G4GLOB_DLL std::ostream G4cout
const XML_Char int const XML_Char * value
Definition: expat.h:331
virtual G4double FindValue(G4double x, G4int componentId=0) const =0
G4double energy(const ThreeVector &p, const G4double m)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
static const G4double pos

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4RDVCrossSectionHandler::FindValue ( G4int  Z,
G4double  e,
G4int  shellIndex 
) const

Definition at line 355 of file G4RDVCrossSectionHandler.cc.

357 {
358  G4double value = 0.;
359 
360  std::map<G4int,G4RDVEMDataSet*,std::less<G4int> >::const_iterator pos;
361  pos = dataMap.find(Z);
362  if (pos!= dataMap.end())
363  {
364  // The following is a workaround for STL ObjectSpace implementation,
365  // which does not support the standard and does not accept
366  // the syntax pos->first or pos->second
367  // G4RDVEMDataSet* dataSet = pos->second;
368  G4RDVEMDataSet* dataSet = (*pos).second;
369  if (shellIndex >= 0)
370  {
371  G4int nComponents = dataSet->NumberOfComponents();
372  if(shellIndex < nComponents)
373  // - MGP - Why doesn't it use G4RDVEMDataSet::FindValue directly?
374  value = dataSet->GetComponent(shellIndex)->FindValue(energy);
375  else
376  {
377  G4cout << "WARNING: G4RDVCrossSectionHandler::FindValue did not find"
378  << " shellIndex= " << shellIndex
379  << " for Z= "
380  << Z << G4endl;
381  }
382  } else {
383  value = dataSet->FindValue(energy);
384  }
385  }
386  else
387  {
388  G4cout << "WARNING: G4RDVCrossSectionHandler::FindValue did not find Z = "
389  << Z << G4endl;
390  }
391  return value;
392 }
virtual const G4RDVEMDataSet * GetComponent(G4int componentId) const =0
int G4int
Definition: G4Types.hh:78
virtual size_t NumberOfComponents(void) const =0
G4GLOB_DLL std::ostream G4cout
const XML_Char int const XML_Char * value
Definition: expat.h:331
virtual G4double FindValue(G4double x, G4int componentId=0) const =0
G4double energy(const ThreeVector &p, const G4double m)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
static const G4double pos

Here is the call graph for this function:

const G4RDVDataSetAlgorithm* G4RDVCrossSectionHandler::GetInterpolation ( ) const
inlineprotected

Definition at line 120 of file G4RDVCrossSectionHandler.hh.

120 { return interpolation; }
void G4RDVCrossSectionHandler::Initialise ( G4RDVDataSetAlgorithm interpolation = 0,
G4double  minE = 250*CLHEP::eV,
G4double  maxE = 100*CLHEP::GeV,
G4int  numberOfBins = 200,
G4double  unitE = CLHEP::MeV,
G4double  unitData = CLHEP::barn,
G4int  minZ = 1,
G4int  maxZ = 99 
)

Definition at line 111 of file G4RDVCrossSectionHandler.cc.

116 {
117  if (algorithm != 0)
118  {
119  delete interpolation;
120  interpolation = algorithm;
121  }
122  else
123  {
124  interpolation = CreateInterpolation();
125  }
126 
127  eMin = minE;
128  eMax = maxE;
129  nBins = numberOfBins;
130  unit1 = unitE;
131  unit2 = unitData;
132  zMin = minZ;
133  zMax = maxZ;
134 }
virtual G4RDVDataSetAlgorithm * CreateInterpolation()

Here is the call graph for this function:

Here is the caller graph for this function:

void G4RDVCrossSectionHandler::LoadData ( const G4String dataFile)

Definition at line 157 of file G4RDVCrossSectionHandler.cc.

158 {
159  size_t nZ = activeZ.size();
160  for (size_t i=0; i<nZ; i++)
161  {
162  G4int Z = (G4int) activeZ[i];
163 
164  // Build the complete string identifying the file with the data set
165 
166  char* path = getenv("G4LEDATA");
167  if (!path)
168  {
169  G4String excep = "G4LEDATA environment variable not set!";
170  G4Exception("G4RDVCrossSectionHandler::LoadData()",
171  "InvalidSetup", FatalException, excep);
172  }
173 
174  std::ostringstream ost;
175  ost << path << '/' << fileName << Z << ".dat";
176  std::ifstream file(ost.str().c_str());
177  std::filebuf* lsdp = file.rdbuf();
178 
179  if (! (lsdp->is_open()) )
180  {
181  G4String excep = "Data file: " + ost.str() + " not found!";
182  G4Exception("G4RDVCrossSectionHandler::LoadData()",
183  "DataNotFound", FatalException, excep);
184  }
185  G4double a = 0;
186  G4int k = 1;
187  G4DataVector* energies = new G4DataVector;
189  do
190  {
191  file >> a;
192  G4int nColumns = 2;
193  // The file is organized into two columns:
194  // 1st column is the energy
195  // 2nd column is the corresponding value
196  // The file terminates with the pattern: -1 -1
197  // -2 -2
198  if (a == -1 || a == -2)
199  {
200  }
201  else
202  {
203  if (k%nColumns != 0)
204  {
205  G4double e = a * unit1;
206  energies->push_back(e);
207  k++;
208  }
209  else if (k%nColumns == 0)
210  {
211  G4double value = a * unit2;
212  data->push_back(value);
213  k = 1;
214  }
215  }
216  } while (a != -2); // end of file
217 
218  file.close();
219  G4RDVDataSetAlgorithm* algo = interpolation->Clone();
220  G4RDVEMDataSet* dataSet = new G4RDEMDataSet(Z,energies,data,algo);
221  dataMap[Z] = dataSet;
222  }
223 }
virtual G4RDVDataSetAlgorithm * Clone() const =0
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
int G4int
Definition: G4Types.hh:78
const XML_Char const XML_Char * data
Definition: expat.h:268
const XML_Char int const XML_Char * value
Definition: expat.h:331
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

void G4RDVCrossSectionHandler::LoadShellData ( const G4String dataFile)

Definition at line 225 of file G4RDVCrossSectionHandler.cc.

226 {
227  size_t nZ = activeZ.size();
228  for (size_t i=0; i<nZ; i++)
229  {
230  G4int Z = (G4int) activeZ[i];
231 
232  // Riccardo Capra <capra@ge.infn.it>: PLEASE CHECK THE FOLLOWING PIECE OF CODE
233  // "energies" AND "data" G4DataVector ARE ALLOCATED, FILLED IN AND NEVER USED OR
234  // DELETED. WHATSMORE LOADING FILE OPERATIONS WERE DONE BY G4RDShellEMDataSet
235  // EVEN BEFORE THE CHANGES I DID ON THIS FILE. SO THE FOLLOWING CODE IN MY
236  // OPINION SHOULD BE USELESS AND SHOULD PRODUCE A MEMORY LEAK.
237 
238  // Build the complete string identifying the file with the data set
239 
240  char* path = getenv("G4LEDATA");
241  if (!path)
242  {
243  G4String excep = "G4LEDATA environment variable not set!";
244  G4Exception("G4RDVCrossSectionHandler::LoadShellData()",
245  "InvalidSetup", FatalException, excep);
246  }
247 
248  std::ostringstream ost;
249 
250  ost << path << '/' << fileName << Z << ".dat";
251 
252  std::ifstream file(ost.str().c_str());
253  std::filebuf* lsdp = file.rdbuf();
254 
255  if (! (lsdp->is_open()) )
256  {
257  G4String excep = "Data file: " + ost.str() + " not found!";
258  G4Exception("G4RDVCrossSectionHandler::LoadShellData()",
259  "DataNotFound", FatalException, excep);
260  }
261  G4double a = 0;
262  G4int k = 1;
263  G4DataVector* energies = new G4DataVector;
265  do
266  {
267  file >> a;
268  G4int nColumns = 2;
269  // The file is organized into two columns:
270  // 1st column is the energy
271  // 2nd column is the corresponding value
272  // The file terminates with the pattern: -1 -1
273  // -2 -2
274  if (a == -1 || a == -2)
275  {
276  }
277  else
278  {
279  if (k%nColumns != 0)
280  {
281  G4double e = a * unit1;
282  energies->push_back(e);
283  k++;
284  }
285  else if (k%nColumns == 0)
286  {
287  G4double value = a * unit2;
288  data->push_back(value);
289  k = 1;
290  }
291  }
292  } while (a != -2); // end of file
293 
294  file.close();
295 
296  // Riccardo Capra <capra@ge.infn.it>: END OF CODE THAT IN MY OPINION SHOULD BE
297  // REMOVED.
298 
299  G4RDVDataSetAlgorithm* algo = interpolation->Clone();
300  G4RDVEMDataSet* dataSet = new G4RDShellEMDataSet(Z, algo);
301  dataSet->LoadData(fileName);
302  dataMap[Z] = dataSet;
303  }
304 }
virtual G4RDVDataSetAlgorithm * Clone() const =0
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
int G4int
Definition: G4Types.hh:78
const XML_Char const XML_Char * data
Definition: expat.h:268
const XML_Char int const XML_Char * value
Definition: expat.h:331
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
double G4double
Definition: G4Types.hh:76
virtual G4bool LoadData(const G4String &fileName)=0

Here is the call graph for this function:

Here is the caller graph for this function:

G4int G4RDVCrossSectionHandler::NumberOfComponents ( G4int  Z) const
protected

Definition at line 657 of file G4RDVCrossSectionHandler.cc.

658 {
659  G4int n = 0;
660 
661  std::map<G4int,G4RDVEMDataSet*,std::less<G4int> >::const_iterator pos;
662  pos = dataMap.find(Z);
663  if (pos!= dataMap.end())
664  {
665  G4RDVEMDataSet* dataSet = (*pos).second;
666  n = dataSet->NumberOfComponents();
667  }
668  else
669  {
670  G4cout << "WARNING: G4RDVCrossSectionHandler::NumberOfComponents did not "
671  << "find Z = "
672  << Z << G4endl;
673  }
674  return n;
675 }
int G4int
Definition: G4Types.hh:78
virtual size_t NumberOfComponents(void) const =0
G4GLOB_DLL std::ostream G4cout
const G4int n
#define G4endl
Definition: G4ios.hh:61
static const G4double pos

Here is the call graph for this function:

Here is the caller graph for this function:

void G4RDVCrossSectionHandler::PrintData ( void  ) const

Definition at line 136 of file G4RDVCrossSectionHandler.cc.

137 {
138  std::map<G4int,G4RDVEMDataSet*,std::less<G4int> >::const_iterator pos;
139 
140  for (pos = dataMap.begin(); pos != dataMap.end(); pos++)
141  {
142  // The following is a workaround for STL ObjectSpace implementation,
143  // which does not support the standard and does not accept
144  // the syntax pos->first or pos->second
145  // G4int z = pos->first;
146  // G4RDVEMDataSet* dataSet = pos->second;
147  G4int z = (*pos).first;
148  G4RDVEMDataSet* dataSet = (*pos).second;
149  G4cout << "---- Data set for Z = "
150  << z
151  << G4endl;
152  dataSet->PrintData();
153  G4cout << "--------------------------------------------------" << G4endl;
154  }
155 }
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
virtual void PrintData(void) const =0
tuple z
Definition: test.py:28
#define G4endl
Definition: G4ios.hh:61
static const G4double pos

Here is the call graph for this function:

Here is the caller graph for this function:

G4int G4RDVCrossSectionHandler::SelectRandomAtom ( const G4MaterialCutsCouple couple,
G4double  e 
) const

Definition at line 498 of file G4RDVCrossSectionHandler.cc.

500 {
501  // Select randomly an element within the material, according to the weight
502  // determined by the cross sections in the data set
503 
504  const G4Material* material = couple->GetMaterial();
505  G4int nElements = material->GetNumberOfElements();
506 
507  // Special case: the material consists of one element
508  if (nElements == 1)
509  {
510  G4int Z = (G4int) material->GetZ();
511  return Z;
512  }
513 
514  // Composite material
515 
516  const G4ElementVector* elementVector = material->GetElementVector();
517  size_t materialIndex = couple->GetIndex();
518 
519  G4RDVEMDataSet* materialSet = (*crossSections)[materialIndex];
520  G4double materialCrossSection0 = 0.0;
521  G4DataVector cross;
522  cross.clear();
523  for ( G4int i=0; i < nElements; i++ )
524  {
525  G4double cr = materialSet->GetComponent(i)->FindValue(e);
526  materialCrossSection0 += cr;
527  cross.push_back(materialCrossSection0);
528  }
529 
530  G4double random = G4UniformRand() * materialCrossSection0;
531 
532  for (G4int k=0 ; k < nElements ; k++ )
533  {
534  if (random <= cross[k]) return (G4int) (*elementVector)[k]->GetZ();
535  }
536  // It should never get here
537  return 0;
538 }
G4double GetZ() const
Definition: G4Material.cc:623
std::vector< G4Element * > G4ElementVector
virtual const G4RDVEMDataSet * GetComponent(G4int componentId) const =0
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:190
int G4int
Definition: G4Types.hh:78
string material
Definition: eplot.py:19
#define G4UniformRand()
Definition: Randomize.hh:97
virtual G4double FindValue(G4double x, G4int componentId=0) const =0
size_t GetNumberOfElements() const
Definition: G4Material.hh:186
double G4double
Definition: G4Types.hh:76
const G4Material * GetMaterial() const

Here is the call graph for this function:

Here is the caller graph for this function:

const G4Element * G4RDVCrossSectionHandler::SelectRandomElement ( const G4MaterialCutsCouple material,
G4double  e 
) const

Definition at line 540 of file G4RDVCrossSectionHandler.cc.

542 {
543  // Select randomly an element within the material, according to the weight determined
544  // by the cross sections in the data set
545 
546  const G4Material* material = couple->GetMaterial();
547  G4Element* nullElement = 0;
548  G4int nElements = material->GetNumberOfElements();
549  const G4ElementVector* elementVector = material->GetElementVector();
550 
551  // Special case: the material consists of one element
552  if (nElements == 1)
553  {
554  G4Element* element = (*elementVector)[0];
555  return element;
556  }
557  else
558  {
559  // Composite material
560 
561  size_t materialIndex = couple->GetIndex();
562 
563  G4RDVEMDataSet* materialSet = (*crossSections)[materialIndex];
564  G4double materialCrossSection0 = 0.0;
565  G4DataVector cross;
566  cross.clear();
567  for (G4int i=0; i<nElements; i++)
568  {
569  G4double cr = materialSet->GetComponent(i)->FindValue(e);
570  materialCrossSection0 += cr;
571  cross.push_back(materialCrossSection0);
572  }
573 
574  G4double random = G4UniformRand() * materialCrossSection0;
575 
576  for (G4int k=0 ; k < nElements ; k++ )
577  {
578  if (random <= cross[k]) return (*elementVector)[k];
579  }
580  // It should never end up here
581  G4cout << "G4RDVCrossSectionHandler::SelectRandomElement - no element found" << G4endl;
582  return nullElement;
583  }
584 }
std::vector< G4Element * > G4ElementVector
virtual const G4RDVEMDataSet * GetComponent(G4int componentId) const =0
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:190
int G4int
Definition: G4Types.hh:78
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
size_t GetIndex() const
Definition: G4Element.hh:182
virtual G4double FindValue(G4double x, G4int componentId=0) const =0
#define G4endl
Definition: G4ios.hh:61
size_t GetNumberOfElements() const
Definition: G4Material.hh:186
double G4double
Definition: G4Types.hh:76
const G4Material * GetMaterial() const

Here is the call graph for this function:

Here is the caller graph for this function:

G4int G4RDVCrossSectionHandler::SelectRandomShell ( G4int  Z,
G4double  e 
) const

Definition at line 586 of file G4RDVCrossSectionHandler.cc.

587 {
588  // Select randomly a shell, according to the weight determined by the cross sections
589  // in the data set
590 
591  // Note for later improvement: it would be useful to add a cache mechanism for already
592  // used shells to improve performance
593 
594  G4int shell = 0;
595 
596  G4double totCrossSection = FindValue(Z,e);
597  G4double random = G4UniformRand() * totCrossSection;
598  G4double partialSum = 0.;
599 
600  G4RDVEMDataSet* dataSet = 0;
601  std::map<G4int,G4RDVEMDataSet*,std::less<G4int> >::const_iterator pos;
602  pos = dataMap.find(Z);
603  // The following is a workaround for STL ObjectSpace implementation,
604  // which does not support the standard and does not accept
605  // the syntax pos->first or pos->second
606  // if (pos != dataMap.end()) dataSet = pos->second;
607  if (pos != dataMap.end()) dataSet = (*pos).second;
608 
609  size_t nShells = dataSet->NumberOfComponents();
610  for (size_t i=0; i<nShells; i++)
611  {
612  const G4RDVEMDataSet* shellDataSet = dataSet->GetComponent(i);
613  if (shellDataSet != 0)
614  {
615  G4double value = shellDataSet->FindValue(e);
616  partialSum += value;
617  if (random <= partialSum) return i;
618  }
619  }
620  // It should never get here
621  return shell;
622 }
virtual const G4RDVEMDataSet * GetComponent(G4int componentId) const =0
G4double FindValue(G4int Z, G4double e) const
int G4int
Definition: G4Types.hh:78
virtual size_t NumberOfComponents(void) const =0
#define G4UniformRand()
Definition: Randomize.hh:97
const XML_Char int const XML_Char * value
Definition: expat.h:331
virtual G4double FindValue(G4double x, G4int componentId=0) const =0
double G4double
Definition: G4Types.hh:76
static const G4double pos

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4RDVCrossSectionHandler::ValueForMaterial ( const G4Material material,
G4double  e 
) const

Definition at line 395 of file G4RDVCrossSectionHandler.cc.

397 {
398  G4double value = 0.;
399 
400  const G4ElementVector* elementVector = material->GetElementVector();
401  const G4double* nAtomsPerVolume = material->GetVecNbOfAtomsPerVolume();
402  G4int nElements = material->GetNumberOfElements();
403 
404  for (G4int i=0 ; i<nElements ; i++)
405  {
406  G4int Z = (G4int) (*elementVector)[i]->GetZ();
407  G4double elementValue = FindValue(Z,energy);
408  G4double nAtomsVol = nAtomsPerVolume[i];
409  value += nAtomsVol * elementValue;
410  }
411 
412  return value;
413 }
std::vector< G4Element * > G4ElementVector
G4double FindValue(G4int Z, G4double e) const
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:190
int G4int
Definition: G4Types.hh:78
const G4double * GetVecNbOfAtomsPerVolume() const
Definition: G4Material.hh:206
const XML_Char int const XML_Char * value
Definition: expat.h:331
G4double energy(const ThreeVector &p, const G4double m)
size_t GetNumberOfElements() const
Definition: G4Material.hh:186
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:


The documentation for this class was generated from the following files: