Geant4  10.03.p01
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
G4VCrossSectionHandler Class Referenceabstract

#include <G4VCrossSectionHandler.hh>

Inheritance diagram for G4VCrossSectionHandler:

Public Member Functions

 G4VCrossSectionHandler ()
 
 G4VCrossSectionHandler (G4VDataSetAlgorithm *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 ~G4VCrossSectionHandler ()
 
void Initialise (G4VDataSetAlgorithm *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
 
G4VEMDataSetBuildMeanFreePathForMaterials (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 LoadNonLogData (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
< G4VEMDataSet * > * 
BuildCrossSectionsForMaterials (const G4DataVector &energyVector, const G4DataVector *energyCuts=0)=0
 
virtual G4VDataSetAlgorithmCreateInterpolation ()
 
const G4VDataSetAlgorithmGetInterpolation () const
 

Detailed Description

Definition at line 64 of file G4VCrossSectionHandler.hh.

Constructor & Destructor Documentation

G4VCrossSectionHandler::G4VCrossSectionHandler ( )

Definition at line 88 of file G4VCrossSectionHandler.cc.

89 {
90  crossSections = 0;
91  interpolation = 0;
92  Initialise();
94 }
void Initialise(G4VDataSetAlgorithm *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:

G4VCrossSectionHandler::G4VCrossSectionHandler ( G4VDataSetAlgorithm 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 97 of file G4VCrossSectionHandler.cc.

105  : interpolation(algorithm), eMin(minE), eMax(maxE), nBins(bins),
106  unit1(unitE), unit2(unitData), zMin(minZ), zMax(maxZ)
107 {
108  crossSections = 0;
109  ActiveElements();
110 }

Here is the call graph for this function:

G4VCrossSectionHandler::~G4VCrossSectionHandler ( )
virtual

Definition at line 112 of file G4VCrossSectionHandler.cc.

113 {
114  delete interpolation;
115  interpolation = 0;
116  std::map<G4int,G4VEMDataSet*,std::less<G4int> >::iterator pos;
117 
118  for (pos = dataMap.begin(); pos != dataMap.end(); ++pos)
119  {
120  // The following is a workaround for STL ObjectSpace implementation,
121  // which does not support the standard and does not accept
122  // the syntax pos->second
123  // G4VEMDataSet* dataSet = pos->second;
124  G4VEMDataSet* dataSet = (*pos).second;
125  delete dataSet;
126  }
127 
128  if (crossSections != 0)
129  {
130  size_t n = crossSections->size();
131  for (size_t i=0; i<n; i++)
132  {
133  delete (*crossSections)[i];
134  }
135  delete crossSections;
136  crossSections = 0;
137  }
138 }
static const G4double pos

Member Function Documentation

void G4VCrossSectionHandler::ActiveElements ( )
protected

Definition at line 695 of file G4VCrossSectionHandler.cc.

696 {
697  const G4MaterialTable* materialTable = G4Material::GetMaterialTable();
698  if (materialTable == 0)
699  G4Exception("G4VCrossSectionHandler::ActiveElements",
700  "em1001",FatalException,"no MaterialTable found");
701 
703 
704  for (G4int mLocal2=0; mLocal2<nMaterials; mLocal2++)
705  {
706  const G4Material* material= (*materialTable)[mLocal2];
707  const G4ElementVector* elementVector = material->GetElementVector();
708  const G4int nElements = material->GetNumberOfElements();
709 
710  for (G4int iEl=0; iEl<nElements; iEl++)
711  {
712  G4Element* element = (*elementVector)[iEl];
713  G4double Z = element->GetZ();
714  if (!(activeZ.contains(Z)) && Z >= zMin && Z <= zMax)
715  {
716  activeZ.push_back(Z);
717  }
718  }
719  }
720 }
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
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<G4VEMDataSet*>* G4VCrossSectionHandler::BuildCrossSectionsForMaterials ( const G4DataVector energyVector,
const G4DataVector energyCuts = 0 
)
protectedpure virtual
G4VEMDataSet * G4VCrossSectionHandler::BuildMeanFreePathForMaterials ( const G4DataVector energyCuts = 0)

Definition at line 463 of file G4VCrossSectionHandler.cc.

464 {
465  // Builds a CompositeDataSet containing the mean free path for each material
466  // in the material table
467 
468  G4DataVector energyVector;
469  G4double dBin = std::log10(eMax/eMin) / nBins;
470 
471  for (G4int i=0; i<nBins+1; i++)
472  {
473  energyVector.push_back(std::pow(10., std::log10(eMin)+i*dBin));
474  }
475 
476  // Factory method to build cross sections in derived classes,
477  // related to the type of physics process
478 
479  if (crossSections != 0)
480  { // Reset the list of cross sections
481  std::vector<G4VEMDataSet*>::iterator mat;
482  if (! crossSections->empty())
483  {
484  for (mat = crossSections->begin(); mat!= crossSections->end(); ++mat)
485  {
486  G4VEMDataSet* set = *mat;
487  delete set;
488  set = 0;
489  }
490  crossSections->clear();
491  delete crossSections;
492  crossSections = 0;
493  }
494  }
495 
496  crossSections = BuildCrossSectionsForMaterials(energyVector,energyCuts);
497 
498  if (crossSections == 0)
499  {
500  G4Exception("G4VCrossSectionHandler::BuildMeanFreePathForMaterials",
501  "em1010",FatalException,"crossSections = 0");
502  return 0;
503  }
504 
506  G4VEMDataSet* materialSet = new G4CompositeEMDataSet(algo);
507  //G4cout << "G4VCrossSectionHandler new dataset " << materialSet << G4endl;
508 
509  G4DataVector* energies;
511  G4DataVector* log_energies;
512  G4DataVector* log_data;
513 
514 
515  const G4ProductionCutsTable* theCoupleTable=
517  size_t numOfCouples = theCoupleTable->GetTableSize();
518 
519 
520  for (size_t mLocal=0; mLocal<numOfCouples; mLocal++)
521  {
522  energies = new G4DataVector;
523  data = new G4DataVector;
524  log_energies = new G4DataVector;
525  log_data = new G4DataVector;
526  for (G4int bin=0; bin<nBins; bin++)
527  {
528  G4double energy = energyVector[bin];
529  energies->push_back(energy);
530  log_energies->push_back(std::log10(energy));
531  G4VEMDataSet* matCrossSet = (*crossSections)[mLocal];
532  G4double materialCrossSection = 0.0;
533  G4int nElm = matCrossSet->NumberOfComponents();
534  for(G4int j=0; j<nElm; j++) {
535  materialCrossSection += matCrossSet->GetComponent(j)->FindValue(energy);
536  }
537 
538  if (materialCrossSection > 0.)
539  {
540  data->push_back(1./materialCrossSection);
541  log_data->push_back(std::log10(1./materialCrossSection));
542  }
543  else
544  {
545  data->push_back(DBL_MAX);
546  log_data->push_back(std::log10(DBL_MAX));
547  }
548  }
550 
551  //G4VEMDataSet* dataSet = new G4EMDataSet(m,energies,data,algo,1.,1.);
552 
553  G4VEMDataSet* dataSet = new G4EMDataSet(mLocal,energies,data,log_energies,log_data,algoLocal,1.,1.);
554 
555  materialSet->AddComponent(dataSet);
556  }
557 
558  return materialSet;
559 }
virtual G4double FindValue(G4double x, G4int componentId=0) const =0
virtual std::vector< G4VEMDataSet * > * BuildCrossSectionsForMaterials(const G4DataVector &energyVector, const G4DataVector *energyCuts=0)=0
int G4int
Definition: G4Types.hh:78
const XML_Char const XML_Char * data
Definition: expat.h:268
virtual const G4VEMDataSet * GetComponent(G4int componentId) const =0
virtual size_t NumberOfComponents(void) const =0
virtual void AddComponent(G4VEMDataSet *dataSet)=0
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
static G4ProductionCutsTable * GetProductionCutsTable()
G4double energy(const ThreeVector &p, const G4double m)
virtual G4VDataSetAlgorithm * CreateInterpolation()
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VCrossSectionHandler::Clear ( )

Definition at line 353 of file G4VCrossSectionHandler.cc.

354 {
355  // Reset the map of data sets: remove the data sets from the map
356  std::map<G4int,G4VEMDataSet*,std::less<G4int> >::iterator pos;
357 
358  if(! dataMap.empty())
359  {
360  for (pos = dataMap.begin(); pos != dataMap.end(); ++pos)
361  {
362  // The following is a workaround for STL ObjectSpace implementation,
363  // which does not support the standard and does not accept
364  // the syntax pos->first or pos->second
365  // G4VEMDataSet* dataSet = pos->second;
366  G4VEMDataSet* dataSet = (*pos).second;
367  delete dataSet;
368  dataSet = 0;
369  G4int i = (*pos).first;
370  dataMap[i] = 0;
371  }
372  dataMap.clear();
373  }
374 
375  activeZ.clear();
376  ActiveElements();
377 }
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:

G4VDataSetAlgorithm * G4VCrossSectionHandler::CreateInterpolation ( )
protectedvirtual

Definition at line 722 of file G4VCrossSectionHandler.cc.

723 {
725  return algorithm;
726 }

Here is the caller graph for this function:

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

Definition at line 379 of file G4VCrossSectionHandler.cc.

380 {
381  G4double value = 0.;
382 
383  std::map<G4int,G4VEMDataSet*,std::less<G4int> >::const_iterator pos;
384  pos = dataMap.find(Z);
385  if (pos!= dataMap.end())
386  {
387  // The following is a workaround for STL ObjectSpace implementation,
388  // which does not support the standard and does not accept
389  // the syntax pos->first or pos->second
390  // G4VEMDataSet* dataSet = pos->second;
391  G4VEMDataSet* dataSet = (*pos).second;
392  value = dataSet->FindValue(energy);
393  }
394  else
395  {
396  G4cout << "WARNING: G4VCrossSectionHandler::FindValue did not find Z = "
397  << Z << G4endl;
398  }
399  return value;
400 }
virtual G4double FindValue(G4double x, G4int componentId=0) const =0
G4GLOB_DLL std::ostream G4cout
const XML_Char int const XML_Char * value
Definition: expat.h:331
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 G4VCrossSectionHandler::FindValue ( G4int  Z,
G4double  e,
G4int  shellIndex 
) const

Definition at line 402 of file G4VCrossSectionHandler.cc.

404 {
405  G4double value = 0.;
406 
407  std::map<G4int,G4VEMDataSet*,std::less<G4int> >::const_iterator pos;
408  pos = dataMap.find(Z);
409  if (pos!= dataMap.end())
410  {
411  // The following is a workaround for STL ObjectSpace implementation,
412  // which does not support the standard and does not accept
413  // the syntax pos->first or pos->second
414  // G4VEMDataSet* dataSet = pos->second;
415  G4VEMDataSet* dataSet = (*pos).second;
416  if (shellIndex >= 0)
417  {
418  G4int nComponents = dataSet->NumberOfComponents();
419  if(shellIndex < nComponents)
420  // - MGP - Why doesn't it use G4VEMDataSet::FindValue directly?
421  value = dataSet->GetComponent(shellIndex)->FindValue(energy);
422  else
423  {
424  G4cout << "WARNING: G4VCrossSectionHandler::FindValue did not find"
425  << " shellIndex= " << shellIndex
426  << " for Z= "
427  << Z << G4endl;
428  }
429  } else {
430  value = dataSet->FindValue(energy);
431  }
432  }
433  else
434  {
435  G4cout << "WARNING: G4VCrossSectionHandler::FindValue did not find Z = "
436  << Z << G4endl;
437  }
438  return value;
439 }
virtual G4double FindValue(G4double x, G4int componentId=0) const =0
int G4int
Definition: G4Types.hh:78
virtual const G4VEMDataSet * GetComponent(G4int componentId) const =0
G4GLOB_DLL std::ostream G4cout
const XML_Char int const XML_Char * value
Definition: expat.h:331
virtual size_t NumberOfComponents(void) 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 G4VDataSetAlgorithm* G4VCrossSectionHandler::GetInterpolation ( ) const
inlineprotected

Definition at line 126 of file G4VCrossSectionHandler.hh.

126 { return interpolation; }
void G4VCrossSectionHandler::Initialise ( G4VDataSetAlgorithm 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 140 of file G4VCrossSectionHandler.cc.

145 {
146  if (algorithm != 0)
147  {
148  delete interpolation;
149  interpolation = algorithm;
150  }
151  else
152  {
153  delete interpolation;
154  interpolation = CreateInterpolation();
155  }
156 
157  eMin = minE;
158  eMax = maxE;
159  nBins = numberOfBins;
160  unit1 = unitE;
161  unit2 = unitData;
162  zMin = minZ;
163  zMax = maxZ;
164 }
virtual G4VDataSetAlgorithm * CreateInterpolation()

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VCrossSectionHandler::LoadData ( const G4String dataFile)

Definition at line 187 of file G4VCrossSectionHandler.cc.

188 {
189  size_t nZ = activeZ.size();
190  for (size_t i=0; i<nZ; i++)
191  {
192  G4int Z = (G4int) activeZ[i];
193 
194  // Build the complete string identifying the file with the data set
195 
196  char* path = getenv("G4LEDATA");
197  if (!path)
198  {
199  G4Exception("G4VCrossSectionHandler::LoadData",
200  "em0006",FatalException,"G4LEDATA environment variable not set");
201  return;
202  }
203 
204  std::ostringstream ost;
205  ost << path << '/' << fileName << Z << ".dat";
206  std::ifstream file(ost.str().c_str());
207  std::filebuf* lsdp = file.rdbuf();
208 
209  if (! (lsdp->is_open()) )
210  {
211  G4String excep = "data file: " + ost.str() + " not found";
212  G4Exception("G4VCrossSectionHandler::LoadData",
213  "em0003",FatalException,excep);
214  }
215  G4double a = 0;
216  G4int k = 0;
217  G4int nColumns = 2;
218 
219  G4DataVector* orig_reg_energies = new G4DataVector;
220  G4DataVector* orig_reg_data = new G4DataVector;
221  G4DataVector* log_reg_energies = new G4DataVector;
222  G4DataVector* log_reg_data = new G4DataVector;
223 
224  do
225  {
226  file >> a;
227 
228  if (a==0.) a=1e-300;
229 
230  // The file is organized into four columns:
231  // 1st column contains the values of energy
232  // 2nd column contains the corresponding data value
233  // The file terminates with the pattern: -1 -1
234  // -2 -2
235  //
236  if (a != -1 && a != -2)
237  {
238  if (k%nColumns == 0)
239  {
240  orig_reg_energies->push_back(a*unit1);
241  log_reg_energies->push_back(std::log10(a)+std::log10(unit1));
242  }
243  else if (k%nColumns == 1)
244  {
245  orig_reg_data->push_back(a*unit2);
246  log_reg_data->push_back(std::log10(a)+std::log10(unit2));
247  }
248  k++;
249  }
250  }
251  while (a != -2); // End of File
252 
253  file.close();
254  G4VDataSetAlgorithm* algo = interpolation->Clone();
255 
256  G4VEMDataSet* dataSet = new G4EMDataSet(Z,orig_reg_energies,orig_reg_data,log_reg_energies,log_reg_data,algo);
257 
258  dataMap[Z] = dataSet;
259 
260  }
261 }
int G4int
Definition: G4Types.hh:78
virtual G4VDataSetAlgorithm * Clone() const =0
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 G4VCrossSectionHandler::LoadNonLogData ( const G4String dataFile)

Definition at line 264 of file G4VCrossSectionHandler.cc.

265 {
266  size_t nZ = activeZ.size();
267  for (size_t i=0; i<nZ; i++)
268  {
269  G4int Z = (G4int) activeZ[i];
270 
271  // Build the complete string identifying the file with the data set
272 
273  char* path = getenv("G4LEDATA");
274  if (!path)
275  {
276  G4Exception("G4VCrossSectionHandler::LoadNonLogData",
277  "em0006",FatalException,"G4LEDATA environment variable not set");
278  return;
279  }
280 
281  std::ostringstream ost;
282  ost << path << '/' << fileName << Z << ".dat";
283  std::ifstream file(ost.str().c_str());
284  std::filebuf* lsdp = file.rdbuf();
285 
286  if (! (lsdp->is_open()) )
287  {
288  G4String excep = "data file: " + ost.str() + " not found";
289  G4Exception("G4VCrossSectionHandler::LoadNonLogData",
290  "em0003",FatalException,excep);
291  }
292  G4double a = 0;
293  G4int k = 0;
294  G4int nColumns = 2;
295 
296  G4DataVector* orig_reg_energies = new G4DataVector;
297  G4DataVector* orig_reg_data = new G4DataVector;
298 
299  do
300  {
301  file >> a;
302 
303  // The file is organized into four columns:
304  // 1st column contains the values of energy
305  // 2nd column contains the corresponding data value
306  // The file terminates with the pattern: -1 -1
307  // -2 -2
308  //
309  if (a != -1 && a != -2)
310  {
311  if (k%nColumns == 0)
312  {
313  orig_reg_energies->push_back(a*unit1);
314  }
315  else if (k%nColumns == 1)
316  {
317  orig_reg_data->push_back(a*unit2);
318  }
319  k++;
320  }
321  }
322  while (a != -2); // End of File
323 
324  file.close();
325  G4VDataSetAlgorithm* algo = interpolation->Clone();
326 
327  G4VEMDataSet* dataSet = new G4EMDataSet(Z,orig_reg_energies,orig_reg_data,algo);
328 
329  dataMap[Z] = dataSet;
330 
331  }
332 }
int G4int
Definition: G4Types.hh:78
virtual G4VDataSetAlgorithm * Clone() const =0
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:

void G4VCrossSectionHandler::LoadShellData ( const G4String dataFile)

Definition at line 334 of file G4VCrossSectionHandler.cc.

335 {
336  size_t nZ = activeZ.size();
337  for (size_t i=0; i<nZ; i++)
338  {
339  G4int Z = (G4int) activeZ[i];
340 
341  G4VDataSetAlgorithm* algo = interpolation->Clone();
342  G4VEMDataSet* dataSet = new G4ShellEMDataSet(Z, algo);
343 
344  dataSet->LoadData(fileName);
345 
346  dataMap[Z] = dataSet;
347  }
348 }
int G4int
Definition: G4Types.hh:78
virtual G4VDataSetAlgorithm * Clone() const =0
virtual G4bool LoadData(const G4String &fileName)=0

Here is the call graph for this function:

Here is the caller graph for this function:

G4int G4VCrossSectionHandler::NumberOfComponents ( G4int  Z) const
protected

Definition at line 728 of file G4VCrossSectionHandler.cc.

729 {
730  G4int n = 0;
731 
732  std::map<G4int,G4VEMDataSet*,std::less<G4int> >::const_iterator pos;
733  pos = dataMap.find(Z);
734  if (pos!= dataMap.end())
735  {
736  G4VEMDataSet* dataSet = (*pos).second;
737  n = dataSet->NumberOfComponents();
738  }
739  else
740  {
741  G4cout << "WARNING: G4VCrossSectionHandler::NumberOfComponents did not "
742  << "find Z = "
743  << Z << G4endl;
744  }
745  return n;
746 }
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
virtual size_t NumberOfComponents(void) const =0
#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 G4VCrossSectionHandler::PrintData ( void  ) const

Definition at line 166 of file G4VCrossSectionHandler.cc.

167 {
168  std::map<G4int,G4VEMDataSet*,std::less<G4int> >::const_iterator pos;
169 
170  for (pos = dataMap.begin(); pos != dataMap.end(); pos++)
171  {
172  // The following is a workaround for STL ObjectSpace implementation,
173  // which does not support the standard and does not accept
174  // the syntax pos->first or pos->second
175  // G4int z = pos->first;
176  // G4VEMDataSet* dataSet = pos->second;
177  G4int z = (*pos).first;
178  G4VEMDataSet* dataSet = (*pos).second;
179  G4cout << "---- Data set for Z = "
180  << z
181  << G4endl;
182  dataSet->PrintData();
183  G4cout << "--------------------------------------------------" << G4endl;
184  }
185 }
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
virtual void PrintData(void) const =0
#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 G4VCrossSectionHandler::SelectRandomAtom ( const G4MaterialCutsCouple couple,
G4double  e 
) const

Definition at line 562 of file G4VCrossSectionHandler.cc.

564 {
565  // Select randomly an element within the material, according to the weight
566  // determined by the cross sections in the data set
567 
568  const G4Material* material = couple->GetMaterial();
569  G4int nElements = material->GetNumberOfElements();
570 
571  // Special case: the material consists of one element
572  if (nElements == 1)
573  {
574  G4int Z = (G4int) material->GetZ();
575  return Z;
576  }
577 
578  // Composite material
579 
580  const G4ElementVector* elementVector = material->GetElementVector();
581  size_t materialIndex = couple->GetIndex();
582 
583  G4VEMDataSet* materialSet = (*crossSections)[materialIndex];
584  G4double materialCrossSection0 = 0.0;
585  G4DataVector cross;
586  cross.clear();
587  for ( G4int i=0; i < nElements; i++ )
588  {
589  G4double cr = materialSet->GetComponent(i)->FindValue(e);
590  materialCrossSection0 += cr;
591  cross.push_back(materialCrossSection0);
592  }
593 
594  G4double random = G4UniformRand() * materialCrossSection0;
595 
596  for (G4int k=0 ; k < nElements ; k++ )
597  {
598  if (random <= cross[k]) return (G4int) (*elementVector)[k]->GetZ();
599  }
600  // It should never get here
601  return 0;
602 }
virtual G4double FindValue(G4double x, G4int componentId=0) const =0
G4double GetZ() const
Definition: G4Material.cc:623
std::vector< G4Element * > G4ElementVector
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:190
int G4int
Definition: G4Types.hh:78
virtual const G4VEMDataSet * GetComponent(G4int componentId) const =0
#define G4UniformRand()
Definition: Randomize.hh:97
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 * G4VCrossSectionHandler::SelectRandomElement ( const G4MaterialCutsCouple material,
G4double  e 
) const

Definition at line 604 of file G4VCrossSectionHandler.cc.

606 {
607  // Select randomly an element within the material, according to the weight determined
608  // by the cross sections in the data set
609 
610  const G4Material* material = couple->GetMaterial();
611  G4Element* nullElement = 0;
612  G4int nElements = material->GetNumberOfElements();
613  const G4ElementVector* elementVector = material->GetElementVector();
614 
615  // Special case: the material consists of one element
616  if (nElements == 1)
617  {
618  G4Element* element = (*elementVector)[0];
619  return element;
620  }
621  else
622  {
623  // Composite material
624 
625  size_t materialIndex = couple->GetIndex();
626 
627  G4VEMDataSet* materialSet = (*crossSections)[materialIndex];
628  G4double materialCrossSection0 = 0.0;
629  G4DataVector cross;
630  cross.clear();
631  for (G4int i=0; i<nElements; i++)
632  {
633  G4double cr = materialSet->GetComponent(i)->FindValue(e);
634  materialCrossSection0 += cr;
635  cross.push_back(materialCrossSection0);
636  }
637 
638  G4double random = G4UniformRand() * materialCrossSection0;
639 
640  for (G4int k=0 ; k < nElements ; k++ )
641  {
642  if (random <= cross[k]) return (*elementVector)[k];
643  }
644  // It should never end up here
645  G4cout << "G4VCrossSectionHandler::SelectRandomElement - no element found" << G4endl;
646  return nullElement;
647  }
648 }
virtual G4double FindValue(G4double x, G4int componentId=0) const =0
std::vector< G4Element * > G4ElementVector
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:190
int G4int
Definition: G4Types.hh:78
virtual const G4VEMDataSet * GetComponent(G4int componentId) const =0
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
size_t GetIndex() const
Definition: G4Element.hh:182
#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:

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

Definition at line 650 of file G4VCrossSectionHandler.cc.

651 {
652  // Select randomly a shell, according to the weight determined by the cross sections
653  // in the data set
654 
655  // Note for later improvement: it would be useful to add a cache mechanism for already
656  // used shells to improve performance
657 
658  G4int shell = 0;
659 
660  G4double totCrossSection = FindValue(Z,e);
661  G4double random = G4UniformRand() * totCrossSection;
662  G4double partialSum = 0.;
663 
664  G4VEMDataSet* dataSet = 0;
665  std::map<G4int,G4VEMDataSet*,std::less<G4int> >::const_iterator pos;
666  pos = dataMap.find(Z);
667  // The following is a workaround for STL ObjectSpace implementation,
668  // which does not support the standard and does not accept
669  // the syntax pos->first or pos->second
670  // if (pos != dataMap.end()) dataSet = pos->second;
671  if (pos != dataMap.end())
672  dataSet = (*pos).second;
673  else
674  {
675  G4Exception("G4VCrossSectionHandler::SelectRandomShell",
676  "em1011",FatalException,"unable to load the dataSet");
677  return 0;
678  }
679 
680  size_t nShells = dataSet->NumberOfComponents();
681  for (size_t i=0; i<nShells; i++)
682  {
683  const G4VEMDataSet* shellDataSet = dataSet->GetComponent(i);
684  if (shellDataSet != 0)
685  {
686  G4double value = shellDataSet->FindValue(e);
687  partialSum += value;
688  if (random <= partialSum) return i;
689  }
690  }
691  // It should never get here
692  return shell;
693 }
virtual G4double FindValue(G4double x, G4int componentId=0) const =0
int G4int
Definition: G4Types.hh:78
virtual const G4VEMDataSet * GetComponent(G4int componentId) const =0
G4double FindValue(G4int Z, G4double e) const
#define G4UniformRand()
Definition: Randomize.hh:97
const XML_Char int const XML_Char * value
Definition: expat.h:331
virtual size_t NumberOfComponents(void) const =0
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
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 G4VCrossSectionHandler::ValueForMaterial ( const G4Material material,
G4double  e 
) const

Definition at line 442 of file G4VCrossSectionHandler.cc.

444 {
445  G4double value = 0.;
446 
447  const G4ElementVector* elementVector = material->GetElementVector();
448  const G4double* nAtomsPerVolume = material->GetVecNbOfAtomsPerVolume();
449  G4int nElements = material->GetNumberOfElements();
450 
451  for (G4int i=0 ; i<nElements ; i++)
452  {
453  G4int Z = (G4int) (*elementVector)[i]->GetZ();
454  G4double elementValue = FindValue(Z,energy);
455  G4double nAtomsVol = nAtomsPerVolume[i];
456  value += nAtomsVol * elementValue;
457  }
458 
459  return value;
460 }
std::vector< G4Element * > G4ElementVector
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:190
int G4int
Definition: G4Types.hh:78
G4double FindValue(G4int Z, G4double e) const
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:


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