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

#include <G4PenelopeCrossSection.hh>

Public Member Functions

 G4PenelopeCrossSection (size_t nOfEnergyPoints, size_t nOfShells=0)
 
 ~G4PenelopeCrossSection ()
 
G4double GetTotalCrossSection (G4double energy) const
 Returns total cross section at the given energy. More...
 
G4double GetHardCrossSection (G4double energy) const
 Returns hard cross section at the given energy. More...
 
G4double GetSoftStoppingPower (G4double energy) const
 Returns the total stopping power due to soft collisions. More...
 
G4double GetShellCrossSection (size_t shellID, G4double energy) const
 Returns the hard cross section for the given shell (per molecule) More...
 
G4double GetNormalizedShellCrossSection (size_t shellID, G4double energy) const
 Returns the hard cross section for the given shell (normalized to 1) More...
 
size_t GetNumberOfShells () const
 
void AddCrossSectionPoint (size_t binNumber, G4double energy, G4double XH0, G4double XH1, G4double XH2, G4double XS0, G4double XS1, G4double XS2)
 
void AddShellCrossSectionPoint (size_t binNumber, size_t shellID, G4double energy, G4double xs)
 
void NormalizeShellCrossSections ()
 

Detailed Description

Definition at line 72 of file G4PenelopeCrossSection.hh.

Constructor & Destructor Documentation

G4PenelopeCrossSection::G4PenelopeCrossSection ( size_t  nOfEnergyPoints,
size_t  nOfShells = 0 
)

Definition at line 45 of file G4PenelopeCrossSection.cc.

45  :
46  numberOfEnergyPoints(nPointsE),numberOfShells(nShells),softCrossSections(0),
47  hardCrossSections(0),shellCrossSections(0),shellNormalizedCrossSections(0)
48 {
49  //check the number of points is not zero
50  if (!numberOfEnergyPoints)
51  {
53  ed << "G4PenelopeCrossSection: invalid number of energy points " << G4endl;
54  G4Exception("G4PenelopeCrossSection::G4PenelopeCrossSection()",
55  "em2017",FatalException,ed);
56  }
57 
58  isNormalized = false;
59 
60  // 1) soft XS table
61  softCrossSections = new G4PhysicsTable();
62  //the table contains 3 G4PhysicsFreeVectors,
63  //(softCrossSections)[0] --> log XS0 vs. log E
64  //(softCrossSections)[1] --> log XS1 vs. log E
65  //(softCrossSections)[2] --> log XS2 vs. log E
66 
67  //everything is log-log
68  for (size_t i=0;i<3;i++)
69  softCrossSections->push_back(new G4PhysicsFreeVector(numberOfEnergyPoints));
70 
71  //2) hard XS table
72  hardCrossSections = new G4PhysicsTable();
73  //the table contains 3 G4PhysicsFreeVectors,
74  //(hardCrossSections)[0] --> log XH0 vs. log E
75  //(hardCrossSections)[1] --> log XH1 vs. log E
76  //(hardCrossSections)[2] --> log XH2 vs. log E
77 
78  //everything is log-log
79  for (size_t i=0;i<3;i++)
80  hardCrossSections->push_back(new G4PhysicsFreeVector(numberOfEnergyPoints));
81 
82  //3) shell XS table, if it is the case
83  if (numberOfShells)
84  {
85  shellCrossSections = new G4PhysicsTable();
86  shellNormalizedCrossSections = new G4PhysicsTable();
87  //the table has to contain numberofShells G4PhysicsFreeVectors,
88  //(theTable)[ishell] --> cross section for shell #ishell
89  for (size_t i=0;i<numberOfShells;i++)
90  {
91  shellCrossSections->push_back(new G4PhysicsFreeVector(numberOfEnergyPoints));
92  shellNormalizedCrossSections->push_back(new G4PhysicsFreeVector(numberOfEnergyPoints));
93  }
94  }
95 }
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
void push_back(G4PhysicsVector *)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

G4PenelopeCrossSection::~G4PenelopeCrossSection ( )

Definition at line 98 of file G4PenelopeCrossSection.cc.

99 {
100  //clean up tables
101  if (shellCrossSections)
102  {
103  //shellCrossSections->clearAndDestroy();
104  delete shellCrossSections;
105  }
106  if (shellNormalizedCrossSections)
107  {
108  //shellNormalizedCrossSections->clearAndDestroy();
109  delete shellNormalizedCrossSections;
110  }
111  if (softCrossSections)
112  {
113  //softCrossSections->clearAndDestroy();
114  delete softCrossSections;
115  }
116  if (hardCrossSections)
117  {
118  //hardCrossSections->clearAndDestroy();
119  delete hardCrossSections;
120  }
121 }

Member Function Documentation

void G4PenelopeCrossSection::AddCrossSectionPoint ( size_t  binNumber,
G4double  energy,
G4double  XH0,
G4double  XH1,
G4double  XH2,
G4double  XS0,
G4double  XS1,
G4double  XS2 
)

Public interface for the master thread

Definition at line 124 of file G4PenelopeCrossSection.cc.

129 {
130  if (!softCrossSections || !hardCrossSections)
131  {
132  G4cout << "Something wrong in G4PenelopeCrossSection::AddCrossSectionPoint" <<
133  G4endl;
134  G4cout << "Trying to fill un-initialized tables" << G4endl;
135  return;
136  }
137 
138  //fill vectors
139  G4PhysicsFreeVector* theVector = (G4PhysicsFreeVector*) (*softCrossSections)[0];
140 
141  if (binNumber >= numberOfEnergyPoints)
142  {
143  G4cout << "Something wrong in G4PenelopeCrossSection::AddCrossSectionPoint" <<
144  G4endl;
145  G4cout << "Trying to register more points than originally declared" << G4endl;
146  return;
147  }
148  G4double logEne = std::log(energy);
149 
150  //XS0
151  G4double val = std::log(std::max(XS0,1e-42*cm2)); //avoid log(0)
152  theVector->PutValue(binNumber,logEne,val);
153 
154  //XS1
155  theVector = (G4PhysicsFreeVector*) (*softCrossSections)[1];
156  val = std::log(std::max(XS1,1e-42*eV*cm2)); //avoid log(0)
157  theVector->PutValue(binNumber,logEne,val);
158 
159  //XS2
160  theVector = (G4PhysicsFreeVector*) (*softCrossSections)[2];
161  val = std::log(std::max(XS2,1e-42*eV*eV*cm2)); //avoid log(0)
162  theVector->PutValue(binNumber,logEne,val);
163 
164  //XH0
165  theVector = (G4PhysicsFreeVector*) (*hardCrossSections)[0];
166  val = std::log(std::max(XH0,1e-42*cm2)); //avoid log(0)
167  theVector->PutValue(binNumber,logEne,val);
168 
169  //XH1
170  theVector = (G4PhysicsFreeVector*) (*hardCrossSections)[1];
171  val = std::log(std::max(XH1,1e-42*eV*cm2)); //avoid log(0)
172  theVector->PutValue(binNumber,logEne,val);
173 
174  //XH2
175  theVector = (G4PhysicsFreeVector*) (*hardCrossSections)[2];
176  val = std::log(std::max(XH2,1e-42*eV*eV*cm2)); //avoid log(0)
177  theVector->PutValue(binNumber,logEne,val);
178 
179  return;
180 }
static constexpr double cm2
Definition: G4SIunits.hh:120
void PutValue(size_t index, G4double energy, G4double dataValue)
G4GLOB_DLL std::ostream G4cout
static constexpr double eV
Definition: G4SIunits.hh:215
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4double energy(const ThreeVector &p, const G4double m)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

void G4PenelopeCrossSection::AddShellCrossSectionPoint ( size_t  binNumber,
size_t  shellID,
G4double  energy,
G4double  xs 
)

Definition at line 184 of file G4PenelopeCrossSection.cc.

188 {
189  if (!shellCrossSections)
190  {
191  G4cout << "Something wrong in G4PenelopeCrossSection::AddShellCrossSectionPoint" <<
192  G4endl;
193  G4cout << "Trying to fill un-initialized table" << G4endl;
194  return;
195  }
196 
197  if (shellID >= numberOfShells)
198  {
199  G4cout << "Something wrong in G4PenelopeCrossSection::AddShellCrossSectionPoint" <<
200  G4endl;
201  G4cout << "Trying to fill shell #" << shellID << " while the maximum is "
202  << numberOfShells-1 << G4endl;
203  return;
204  }
205 
206  //fill vector
207  G4PhysicsFreeVector* theVector = (G4PhysicsFreeVector*) (*shellCrossSections)[shellID];
208 
209  if (binNumber >= numberOfEnergyPoints)
210  {
211  G4cout << "Something wrong in G4PenelopeCrossSection::AddShellCrossSectionPoint" <<
212  G4endl;
213  G4cout << "Trying to register more points than originally declared" << G4endl;
214  return;
215  }
216  G4double logEne = std::log(energy);
217  G4double val = std::log(std::max(xs,1e-42*cm2)); //avoid log(0)
218  theVector->PutValue(binNumber,logEne,val);
219 
220  return;
221 }
static constexpr double cm2
Definition: G4SIunits.hh:120
void PutValue(size_t index, G4double energy, G4double dataValue)
G4GLOB_DLL std::ostream G4cout
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4double energy(const ThreeVector &p, const G4double m)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4PenelopeCrossSection::GetHardCrossSection ( G4double  energy) const

Returns hard cross section at the given energy.

Definition at line 269 of file G4PenelopeCrossSection.cc.

270 {
271  G4double result = 0;
272  //take here XH0
273  if (!hardCrossSections)
274  {
275  G4cout << "Something wrong in G4PenelopeCrossSection::GetHardCrossSection" <<
276  G4endl;
277  G4cout << "Trying to retrieve from un-initialized tables" << G4endl;
278  return result;
279  }
280 
281  G4PhysicsFreeVector* theVector = (G4PhysicsFreeVector*) (*hardCrossSections)[0];
282  if (theVector->GetVectorLength() < numberOfEnergyPoints)
283  {
284  G4cout << "Something wrong in G4PenelopeCrossSection::GetHardCrossSection" <<
285  G4endl;
286  G4cout << "Hard cross section table looks not filled" << G4endl;
287  return result;
288  }
289  G4double logene = std::log(energy);
290  G4double logXS = theVector->Value(logene);
291  result = G4Exp(logXS);
292 
293  return result;
294 }
G4double G4ParticleHPJENDLHEData::G4double result
size_t GetVectorLength() const
G4GLOB_DLL std::ostream G4cout
G4double Value(G4double theEnergy, size_t &lastidx) const
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
G4double energy(const ThreeVector &p, const G4double m)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4PenelopeCrossSection::GetNormalizedShellCrossSection ( size_t  shellID,
G4double  energy 
) const

Returns the hard cross section for the given shell (normalized to 1)

Definition at line 364 of file G4PenelopeCrossSection.cc.

365 {
366  G4double result = 0;
367  if (!shellNormalizedCrossSections)
368  {
369  G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" <<
370  G4endl;
371  G4cout << "Trying to retrieve from un-initialized tables" << G4endl;
372  return result;
373  }
374 
375  if (!isNormalized)
376  {
377  G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" << G4endl;
378  G4cout << "The table of normalized cross section is not initialized" << G4endl;
379  }
380 
381 
382  if (shellID >= numberOfShells)
383  {
384  G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" <<
385  G4endl;
386  G4cout << "Trying to retrieve shell #" << shellID << " while the maximum is "
387  << numberOfShells-1 << G4endl;
388  return result;
389  }
390 
391  const G4PhysicsFreeVector* theVector =
392  (G4PhysicsFreeVector*) (*shellNormalizedCrossSections)[shellID];
393 
394  if (theVector->GetVectorLength() < numberOfEnergyPoints)
395  {
396  G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" <<
397  G4endl;
398  G4cout << "Shell cross section table looks not filled" << G4endl;
399  return result;
400  }
401  G4double logene = std::log(energy);
402  G4double logXS = theVector->Value(logene);
403  result = G4Exp(logXS);
404 
405  return result;
406 }
G4double G4ParticleHPJENDLHEData::G4double result
size_t GetVectorLength() const
G4GLOB_DLL std::ostream G4cout
G4double Value(G4double theEnergy, size_t &lastidx) const
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
G4double energy(const ThreeVector &p, const G4double m)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

size_t G4PenelopeCrossSection::GetNumberOfShells ( ) const
inline

Definition at line 94 of file G4PenelopeCrossSection.hh.

94 {return numberOfShells;};
G4double G4PenelopeCrossSection::GetShellCrossSection ( size_t  shellID,
G4double  energy 
) const

Returns the hard cross section for the given shell (per molecule)

Definition at line 328 of file G4PenelopeCrossSection.cc.

329 {
330  G4double result = 0;
331  if (!shellCrossSections)
332  {
333  G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" <<
334  G4endl;
335  G4cout << "Trying to retrieve from un-initialized tables" << G4endl;
336  return result;
337  }
338  if (shellID >= numberOfShells)
339  {
340  G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" <<
341  G4endl;
342  G4cout << "Trying to retrieve shell #" << shellID << " while the maximum is "
343  << numberOfShells-1 << G4endl;
344  return result;
345  }
346 
347  G4PhysicsFreeVector* theVector = (G4PhysicsFreeVector*) (*shellCrossSections)[shellID];
348 
349  if (theVector->GetVectorLength() < numberOfEnergyPoints)
350  {
351  G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" <<
352  G4endl;
353  G4cout << "Shell cross section table looks not filled" << G4endl;
354  return result;
355  }
356  G4double logene = std::log(energy);
357  G4double logXS = theVector->Value(logene);
358  result = G4Exp(logXS);
359 
360  return result;
361 }
G4double G4ParticleHPJENDLHEData::G4double result
size_t GetVectorLength() const
G4GLOB_DLL std::ostream G4cout
G4double Value(G4double theEnergy, size_t &lastidx) const
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
G4double energy(const ThreeVector &p, const G4double m)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4PenelopeCrossSection::GetSoftStoppingPower ( G4double  energy) const

Returns the total stopping power due to soft collisions.

Definition at line 299 of file G4PenelopeCrossSection.cc.

300 {
301  G4double result = 0;
302  //take here XH0
303  if (!softCrossSections)
304  {
305  G4cout << "Something wrong in G4PenelopeCrossSection::GetSoftStoppingPower" <<
306  G4endl;
307  G4cout << "Trying to retrieve from un-initialized tables" << G4endl;
308  return result;
309  }
310 
311  G4PhysicsFreeVector* theVector = (G4PhysicsFreeVector*) (*softCrossSections)[1];
312  if (theVector->GetVectorLength() < numberOfEnergyPoints)
313  {
314  G4cout << "Something wrong in G4PenelopeCrossSection::GetSoftStoppingPower" <<
315  G4endl;
316  G4cout << "Soft cross section table looks not filled" << G4endl;
317  return result;
318  }
319  G4double logene = std::log(energy);
320  G4double logXS = theVector->Value(logene);
321  result = G4Exp(logXS);
322 
323  return result;
324 }
G4double G4ParticleHPJENDLHEData::G4double result
size_t GetVectorLength() const
G4GLOB_DLL std::ostream G4cout
G4double Value(G4double theEnergy, size_t &lastidx) const
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
G4double energy(const ThreeVector &p, const G4double m)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4PenelopeCrossSection::GetTotalCrossSection ( G4double  energy) const

Returns total cross section at the given energy.

Definition at line 225 of file G4PenelopeCrossSection.cc.

226 {
227  G4double result = 0;
228  //take here XS0 + XH0
229  if (!softCrossSections || !hardCrossSections)
230  {
231  G4cout << "Something wrong in G4PenelopeCrossSection::GetTotalCrossSection" <<
232  G4endl;
233  G4cout << "Trying to retrieve from un-initialized tables" << G4endl;
234  return result;
235  }
236 
237  // 1) soft part
238  G4PhysicsFreeVector* theVector = (G4PhysicsFreeVector*) (*softCrossSections)[0];
239  if (theVector->GetVectorLength() < numberOfEnergyPoints)
240  {
241  G4cout << "Something wrong in G4PenelopeCrossSection::GetTotalCrossSection" <<
242  G4endl;
243  G4cout << "Soft cross section table looks not filled" << G4endl;
244  return result;
245  }
246  G4double logene = std::log(energy);
247  G4double logXS = theVector->Value(logene);
248  G4double softXS = G4Exp(logXS);
249 
250  // 2) hard part
251  theVector = (G4PhysicsFreeVector*) (*hardCrossSections)[0];
252  if (theVector->GetVectorLength() < numberOfEnergyPoints)
253  {
254  G4cout << "Something wrong in G4PenelopeCrossSection::GetTotalCrossSection" <<
255  G4endl;
256  G4cout << "Hard cross section table looks not filled" << G4endl;
257  return result;
258  }
259  logXS = theVector->Value(logene);
260  G4double hardXS = G4Exp(logXS);
261 
262  result = hardXS + softXS;
263  return result;
264 
265 }
G4double G4ParticleHPJENDLHEData::G4double result
size_t GetVectorLength() const
G4GLOB_DLL std::ostream G4cout
G4double Value(G4double theEnergy, size_t &lastidx) const
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
G4double energy(const ThreeVector &p, const G4double m)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

void G4PenelopeCrossSection::NormalizeShellCrossSections ( )

Definition at line 412 of file G4PenelopeCrossSection.cc.

413 {
414  if (isNormalized) //already done!
415  {
416  G4cout << "G4PenelopeCrossSection::NormalizeShellCrossSections()" << G4endl;
417  G4cout << "already invoked. Ignore it" << G4endl;
418  return;
419  }
420 
421  if (!shellNormalizedCrossSections)
422  {
423  G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" <<
424  G4endl;
425  G4cout << "Trying to retrieve from un-initialized tables" << G4endl;
426  return;
427  }
428 
429  for (size_t i=0;i<numberOfEnergyPoints;i++) //loop on energy
430  {
431  //energy grid is the same for all shells
432 
433  //Recalculate manually the XS factor, to avoid problems with
434  //underflows
435  G4double normFactor = 0.;
436  for (size_t shellID=0;shellID<numberOfShells;shellID++)
437  {
438  G4PhysicsFreeVector* theVec =
439  (G4PhysicsFreeVector*) (*shellCrossSections)[shellID];
440 
441  normFactor += G4Exp((*theVec)[i]);
442  }
443  G4double logNormFactor = std::log(normFactor);
444  //Normalize
445  for (size_t shellID=0;shellID<numberOfShells;shellID++)
446  {
447  G4PhysicsFreeVector* theVec =
448  (G4PhysicsFreeVector*) (*shellNormalizedCrossSections)[shellID];
449  G4PhysicsFreeVector* theFullVec =
450  (G4PhysicsFreeVector*) (*shellCrossSections)[shellID];
451  G4double previousValue = (*theFullVec)[i]; //log(XS)
452  G4double logEnergy = theFullVec->GetLowEdgeEnergy(i);
453  //log(XS/normFactor) = log(XS) - log(normFactor)
454  theVec->PutValue(i,logEnergy,previousValue-logNormFactor);
455  }
456  }
457 
458  isNormalized = true;
459 
460 
461  /*
462  //TESTING
463  for (size_t shellID=0;shellID<numberOfShells;shellID++)
464  {
465  G4cout << "SHELL " << shellID << G4endl;
466  G4PhysicsFreeVector* theVec =
467  (G4PhysicsFreeVector*) (*shellCrossSections)[shellID];
468  for (size_t i=0;i<numberOfEnergyPoints;i++) //loop on energy
469  {
470  G4double logene = theVec->GetLowEdgeEnergy(i);
471  G4cout << G4Exp(logene)/MeV << " " << G4Exp((*theVec)[i]) << G4endl;
472  }
473  }
474  */
475 
476  return;
477 }
void PutValue(size_t index, G4double energy, G4double dataValue)
G4double GetLowEdgeEnergy(size_t binNumber) const
G4GLOB_DLL std::ostream G4cout
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
#define G4endl
Definition: G4ios.hh:61
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: