Geant4-11
Public Member Functions | Private Attributes
G4PenelopeCrossSection Class Reference

#include <G4PenelopeCrossSection.hh>

Public Member Functions

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)
 
 G4PenelopeCrossSection (const G4PenelopeCrossSection &)=delete
 
 G4PenelopeCrossSection (size_t nOfEnergyPoints, size_t nOfShells=0)
 
G4double GetHardCrossSection (G4double energy) const
 Returns hard cross section at the given energy. 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
 
G4double GetShellCrossSection (size_t shellID, G4double energy) const
 Returns the hard cross section for the given shell (per molecule) More...
 
G4double GetSoftStoppingPower (G4double energy) const
 Returns the total stopping power due to soft collisions. More...
 
G4double GetTotalCrossSection (G4double energy) const
 Returns total cross section at the given energy. More...
 
void NormalizeShellCrossSections ()
 
G4PenelopeCrossSectionoperator= (const G4PenelopeCrossSection &right)=delete
 
 ~G4PenelopeCrossSection ()
 

Private Attributes

G4PhysicsTablefHardCrossSections
 
G4bool fIsNormalized
 
size_t fNumberOfEnergyPoints
 
size_t fNumberOfShells
 
G4PhysicsTablefShellCrossSections
 
G4PhysicsTablefShellNormalizedCrossSections
 
G4PhysicsTablefSoftCrossSections
 

Detailed Description

Definition at line 71 of file G4PenelopeCrossSection.hh.

Constructor & Destructor Documentation

◆ G4PenelopeCrossSection() [1/2]

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

Definition at line 45 of file G4PenelopeCrossSection.cc.

45 :
46 fSoftCrossSections(nullptr),
49 fNumberOfEnergyPoints(nPointsE),fNumberOfShells(nShells)
50{
51 //check the number of points is not zero
53 {
55 ed << "G4PenelopeCrossSection: invalid number of energy points " << G4endl;
56 G4Exception("G4PenelopeCrossSection::G4PenelopeCrossSection()",
57 "em2017",FatalException,ed);
58 }
59
60 fIsNormalized = false;
61
62 // 1) soft XS table
64 //the table contains 3 G4PhysicsFreeVectors,
65 //(fSoftCrossSections)[0] --> log XS0 vs. log E
66 //(fSoftCrossSections)[1] --> log XS1 vs. log E
67 //(fSoftCrossSections)[2] --> log XS2 vs. log E
68
69 //everything is log-log
70 for (size_t i=0;i<3;i++)
72
73 //2) hard XS table
75 //the table contains 3 G4PhysicsFreeVectors,
76 //(fHardCrossSections)[0] --> log XH0 vs. log E
77 //(fHardCrossSections)[1] --> log XH1 vs. log E
78 //(fHardCrossSections)[2] --> log XH2 vs. log E
79
80 //everything is log-log
81 for (size_t i=0;i<3;i++)
83
84 //3) shell XS table, if it is the case
86 {
89 //the table has to contain numberofShells G4PhysicsFreeVectors,
90 //(theTable)[ishell] --> cross section for shell #ishell
91 for (size_t i=0;i<fNumberOfShells;i++)
92 {
95 }
96 }
97}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:40
#define G4endl
Definition: G4ios.hh:57
G4PhysicsTable * fShellNormalizedCrossSections
G4PhysicsTable * fShellCrossSections
void push_back(G4PhysicsVector *)

References FatalException, fHardCrossSections, fIsNormalized, fNumberOfEnergyPoints, fNumberOfShells, fShellCrossSections, fShellNormalizedCrossSections, fSoftCrossSections, G4endl, G4Exception(), and G4PhysicsTable::push_back().

◆ ~G4PenelopeCrossSection()

G4PenelopeCrossSection::~G4PenelopeCrossSection ( )

◆ G4PenelopeCrossSection() [2/2]

G4PenelopeCrossSection::G4PenelopeCrossSection ( const G4PenelopeCrossSection )
delete

Member Function Documentation

◆ AddCrossSectionPoint()

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 126 of file G4PenelopeCrossSection.cc.

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

References cm2, G4INCL::KinematicsUtils::energy(), eV, fHardCrossSections, fNumberOfEnergyPoints, fSoftCrossSections, G4cout, G4endl, G4Log(), G4INCL::Math::max(), and G4PhysicsFreeVector::PutValues().

Referenced by G4PenelopeIonisationXSHandler::BuildXSTable(), and G4PenelopeBremsstrahlungModel::BuildXSTable().

◆ AddShellCrossSectionPoint()

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

Definition at line 186 of file G4PenelopeCrossSection.cc.

190{
192 {
193 G4cout << "Something wrong in G4PenelopeCrossSection::AddShellCrossSectionPoint" <<
194 G4endl;
195 G4cout << "Trying to fill un-initialized table" << G4endl;
196 return;
197 }
198
199 if (shellID >= fNumberOfShells)
200 {
201 G4cout << "Something wrong in G4PenelopeCrossSection::AddShellCrossSectionPoint" <<
202 G4endl;
203 G4cout << "Trying to fill shell #" << shellID << " while the maximum is "
204 << fNumberOfShells-1 << G4endl;
205 return;
206 }
207
208 //fill vector
210
211 if (binNumber >= fNumberOfEnergyPoints)
212 {
213 G4cout << "Something wrong in G4PenelopeCrossSection::AddShellCrossSectionPoint" <<
214 G4endl;
215 G4cout << "Trying to register more points than originally declared" << G4endl;
216 return;
217 }
218 G4double logEne = G4Log(energy);
219 G4double val = G4Log(std::max(xs,1e-42*cm2)); //avoid log(0)
220 theVector->PutValues(binNumber,logEne,val);
221
222 return;
223}

References cm2, G4INCL::KinematicsUtils::energy(), fNumberOfEnergyPoints, fNumberOfShells, fShellCrossSections, G4cout, G4endl, G4Log(), G4INCL::Math::max(), and G4PhysicsFreeVector::PutValues().

Referenced by G4PenelopeIonisationXSHandler::BuildXSTable().

◆ GetHardCrossSection()

G4double G4PenelopeCrossSection::GetHardCrossSection ( G4double  energy) const

Returns hard cross section at the given energy.

Definition at line 270 of file G4PenelopeCrossSection.cc.

271{
272 G4double result = 0;
273 //take here XH0
275 {
276 G4cout << "Something wrong in G4PenelopeCrossSection::GetHardCrossSection" <<
277 G4endl;
278 G4cout << "Trying to retrieve from un-initialized tables" << G4endl;
279 return result;
280 }
281
283 if (theVector->GetVectorLength() < fNumberOfEnergyPoints)
284 {
285 G4cout << "Something wrong in G4PenelopeCrossSection::GetHardCrossSection" <<
286 G4endl;
287 G4cout << "Hard cross section table looks not filled" << G4endl;
288 return result;
289 }
290 G4double logene = G4Log(energy);
291 G4double logXS = theVector->Value(logene);
292 result = G4Exp(logXS);
293
294 return result;
295}
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:179
G4double Value(const G4double energy, std::size_t &lastidx) const
std::size_t GetVectorLength() const

References G4INCL::KinematicsUtils::energy(), fHardCrossSections, fNumberOfEnergyPoints, G4cout, G4endl, G4Exp(), G4Log(), G4PhysicsVector::GetVectorLength(), and G4PhysicsVector::Value().

Referenced by G4PenelopeBremsstrahlungModel::CrossSectionPerVolume(), and G4PenelopeIonisationModel::CrossSectionPerVolume().

◆ GetNormalizedShellCrossSection()

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;
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 (!fIsNormalized)
376 {
377 G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" << G4endl;
378 G4cout << "The table of normalized cross section is not initialized" << G4endl;
379 }
380
381 if (shellID >= fNumberOfShells)
382 {
383 G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" <<
384 G4endl;
385 G4cout << "Trying to retrieve shell #" << shellID << " while the maximum is "
386 << fNumberOfShells-1 << G4endl;
387 return result;
388 }
389
390 const G4PhysicsFreeVector* theVector =
392
393 if (theVector->GetVectorLength() < fNumberOfEnergyPoints)
394 {
395 G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" <<
396 G4endl;
397 G4cout << "Shell cross section table looks not filled" << G4endl;
398 return result;
399 }
400 G4double logene = G4Log(energy);
401 G4double logXS = theVector->Value(logene);
402 result = G4Exp(logXS);
403
404 return result;
405}

References G4INCL::KinematicsUtils::energy(), fIsNormalized, fNumberOfEnergyPoints, fNumberOfShells, fShellNormalizedCrossSections, G4cout, G4endl, G4Exp(), G4Log(), G4PhysicsVector::GetVectorLength(), and G4PhysicsVector::Value().

Referenced by G4PenelopeIonisationModel::SampleFinalStateElectron(), and G4PenelopeIonisationModel::SampleFinalStatePositron().

◆ GetNumberOfShells()

size_t G4PenelopeCrossSection::GetNumberOfShells ( ) const
inline

Definition at line 93 of file G4PenelopeCrossSection.hh.

93{return fNumberOfShells;};

References fNumberOfShells.

◆ GetShellCrossSection()

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;
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 >= fNumberOfShells)
339 {
340 G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" <<
341 G4endl;
342 G4cout << "Trying to retrieve shell #" << shellID << " while the maximum is "
343 << fNumberOfShells-1 << G4endl;
344 return result;
345 }
346
348
349 if (theVector->GetVectorLength() < fNumberOfEnergyPoints)
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 = G4Log(energy);
357 G4double logXS = theVector->Value(logene);
358 result = G4Exp(logXS);
359
360 return result;
361}

References G4INCL::KinematicsUtils::energy(), fNumberOfEnergyPoints, fNumberOfShells, fShellCrossSections, G4cout, G4endl, G4Exp(), G4Log(), G4PhysicsVector::GetVectorLength(), and G4PhysicsVector::Value().

Referenced by G4PenelopeIonisationCrossSection::CrossSection().

◆ GetSoftStoppingPower()

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
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
312 if (theVector->GetVectorLength() < fNumberOfEnergyPoints)
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 = G4Log(energy);
320 G4double logXS = theVector->Value(logene);
321 result = G4Exp(logXS);
322
323 return result;
324}

References G4INCL::KinematicsUtils::energy(), fNumberOfEnergyPoints, fSoftCrossSections, G4cout, G4endl, G4Exp(), G4Log(), G4PhysicsVector::GetVectorLength(), and G4PhysicsVector::Value().

Referenced by G4PenelopeBremsstrahlungModel::ComputeDEDXPerVolume(), and G4PenelopeIonisationModel::ComputeDEDXPerVolume().

◆ GetTotalCrossSection()

G4double G4PenelopeCrossSection::GetTotalCrossSection ( G4double  energy) const

Returns total cross section at the given energy.

Definition at line 227 of file G4PenelopeCrossSection.cc.

228{
229 G4double result = 0;
230 //take here XS0 + XH0
232 {
233 G4cout << "Something wrong in G4PenelopeCrossSection::GetTotalCrossSection" <<
234 G4endl;
235 G4cout << "Trying to retrieve from un-initialized tables" << G4endl;
236 return result;
237 }
238
239 // 1) soft part
241 if (theVector->GetVectorLength() < fNumberOfEnergyPoints)
242 {
243 G4cout << "Something wrong in G4PenelopeCrossSection::GetTotalCrossSection" <<
244 G4endl;
245 G4cout << "Soft cross section table looks not filled" << G4endl;
246 return result;
247 }
248 G4double logene = G4Log(energy);
249 G4double logXS = theVector->Value(logene);
250 G4double softXS = G4Exp(logXS);
251
252 // 2) hard part
253 theVector = (G4PhysicsFreeVector*) (*fHardCrossSections)[0];
254 if (theVector->GetVectorLength() < fNumberOfEnergyPoints)
255 {
256 G4cout << "Something wrong in G4PenelopeCrossSection::GetTotalCrossSection" <<
257 G4endl;
258 G4cout << "Hard cross section table looks not filled" << G4endl;
259 return result;
260 }
261 logXS = theVector->Value(logene);
262 G4double hardXS = G4Exp(logXS);
263
264 result = hardXS + softXS;
265 return result;
266}

References G4INCL::KinematicsUtils::energy(), fHardCrossSections, fNumberOfEnergyPoints, fSoftCrossSections, G4cout, G4endl, G4Exp(), G4Log(), G4PhysicsVector::GetVectorLength(), and G4PhysicsVector::Value().

Referenced by G4PenelopeIonisationModel::CrossSectionPerVolume().

◆ NormalizeShellCrossSections()

void G4PenelopeCrossSection::NormalizeShellCrossSections ( )

Definition at line 409 of file G4PenelopeCrossSection.cc.

410{
411 if (fIsNormalized) //already done!
412 {
413 G4cout << "G4PenelopeCrossSection::NormalizeShellCrossSections()" << G4endl;
414 G4cout << "already invoked. Ignore it" << G4endl;
415 return;
416 }
417
419 {
420 G4cout << "Something wrong in G4PenelopeCrossSection::GetShellCrossSection" <<
421 G4endl;
422 G4cout << "Trying to retrieve from un-initialized tables" << G4endl;
423 return;
424 }
425
426 for (size_t i=0;i<fNumberOfEnergyPoints;i++) //loop on energy
427 {
428 //energy grid is the same for all shells
429
430 //Recalculate manually the XS factor, to avoid problems with
431 //underflows
432 G4double normFactor = 0.;
433 for (size_t shellID=0;shellID<fNumberOfShells;shellID++)
434 {
435 G4PhysicsFreeVector* theVec =
437
438 normFactor += G4Exp((*theVec)[i]);
439 }
440 G4double logNormFactor = G4Log(normFactor);
441 //Normalize
442 for (size_t shellID=0;shellID<fNumberOfShells;shellID++)
443 {
444 G4PhysicsFreeVector* theVec =
446 G4PhysicsFreeVector* theFullVec =
448 G4double previousValue = (*theFullVec)[i]; //log(XS)
449 G4double logEnergy = theFullVec->GetLowEdgeEnergy(i);
450 //log(XS/normFactor) = log(XS) - log(normFactor)
451 theVec->PutValues(i,logEnergy,previousValue-logNormFactor);
452 }
453 }
454 fIsNormalized = true;
455
456 return;
457}
G4double GetLowEdgeEnergy(const std::size_t index) const

References fIsNormalized, fNumberOfEnergyPoints, fNumberOfShells, fShellCrossSections, fShellNormalizedCrossSections, G4cout, G4endl, G4Exp(), G4Log(), G4PhysicsVector::GetLowEdgeEnergy(), and G4PhysicsFreeVector::PutValues().

Referenced by G4PenelopeIonisationXSHandler::BuildXSTable().

◆ operator=()

G4PenelopeCrossSection & G4PenelopeCrossSection::operator= ( const G4PenelopeCrossSection right)
delete

Field Documentation

◆ fHardCrossSections

G4PhysicsTable* G4PenelopeCrossSection::fHardCrossSections
private

◆ fIsNormalized

G4bool G4PenelopeCrossSection::fIsNormalized
private

◆ fNumberOfEnergyPoints

size_t G4PenelopeCrossSection::fNumberOfEnergyPoints
private

◆ fNumberOfShells

size_t G4PenelopeCrossSection::fNumberOfShells
private

◆ fShellCrossSections

G4PhysicsTable* G4PenelopeCrossSection::fShellCrossSections
private

◆ fShellNormalizedCrossSections

G4PhysicsTable* G4PenelopeCrossSection::fShellNormalizedCrossSections
private

◆ fSoftCrossSections

G4PhysicsTable* G4PenelopeCrossSection::fSoftCrossSections
private

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