Geant4-11
Public Member Functions | Static Public Member Functions | Protected Member Functions | Protected Attributes | Static Protected Attributes | Private Member Functions | Static Private Member Functions
G4RToEConvForGamma Class Reference

#include <G4RToEConvForGamma.hh>

Inheritance diagram for G4RToEConvForGamma:
G4VRangeToEnergyConverter

Public Member Functions

virtual G4double Convert (const G4double rangeCut, const G4Material *material)
 
 G4RToEConvForGamma ()
 
const G4ParticleDefinitionGetParticleType () const
 
G4int GetVerboseLevel () const
 
G4bool operator!= (const G4VRangeToEnergyConverter &r) const =delete
 
G4bool operator== (const G4VRangeToEnergyConverter &r) const =delete
 
void SetVerboseLevel (G4int value)
 
virtual ~G4RToEConvForGamma ()
 

Static Public Member Functions

static G4double GetHighEdgeEnergy ()
 
static G4double GetLowEdgeEnergy ()
 
static G4double GetMaxEnergyCut ()
 
static void SetEnergyRange (const G4double lowedge, const G4double highedge)
 
static void SetMaxEnergyCut (const G4double value)
 

Protected Member Functions

G4double ComputeValue (const G4int Z, const G4double kinEnergy) final
 

Protected Attributes

G4int fPDG = 0
 
G4bool isFirstInstance = false
 
const G4ParticleDefinitiontheParticle = nullptr
 
G4int verboseLevel = 1
 

Static Protected Attributes

static G4double Emax = 0.0
 
static G4double Emin = 0.0
 
static std::vector< G4double > * Energy = nullptr
 
static G4int Nbin = 350
 
static G4int NbinPerDecade = 50
 

Private Member Functions

G4double ConvertForElectron (const G4double rangeCut, const G4Material *material)
 
G4double ConvertForGamma (const G4double rangeCut, const G4Material *material)
 
G4double LiniearInterpolation (const G4double e1, const G4double e2, const G4double r1, const G4double r2, const G4double r)
 

Static Private Member Functions

static void FillEnergyVector (const G4double emin, const G4double emax)
 

Detailed Description

Definition at line 40 of file G4RToEConvForGamma.hh.

Constructor & Destructor Documentation

◆ G4RToEConvForGamma()

G4RToEConvForGamma::G4RToEConvForGamma ( )
explicit

Definition at line 40 of file G4RToEConvForGamma.cc.

42{
44 if (theParticle == nullptr)
45 {
46#ifdef G4VERBOSE
47 if (GetVerboseLevel()>0)
48 {
49 G4cout << " G4RToEConvForGamma::G4RToEConvForGamma() - ";
50 G4cout << "Gamma is not defined !!" << G4endl;
51 }
52#endif
53 }
54 else
55 {
57 }
58}
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
static G4ParticleTable * GetParticleTable()
const G4ParticleDefinition * theParticle

References G4ParticleTable::FindParticle(), G4VRangeToEnergyConverter::fPDG, G4cout, G4endl, G4ParticleTable::GetParticleTable(), G4ParticleDefinition::GetPDGEncoding(), G4VRangeToEnergyConverter::GetVerboseLevel(), and G4VRangeToEnergyConverter::theParticle.

◆ ~G4RToEConvForGamma()

G4RToEConvForGamma::~G4RToEConvForGamma ( )
virtual

Definition at line 61 of file G4RToEConvForGamma.cc.

62{}

Member Function Documentation

◆ ComputeValue()

G4double G4RToEConvForGamma::ComputeValue ( const G4int  Z,
const G4double  kinEnergy 
)
finalprotectedvirtual

Implements G4VRangeToEnergyConverter.

Definition at line 65 of file G4RToEConvForGamma.cc.

67{
68 // Compute the "absorption" cross-section of the photon "absorption".
69 // Cross-section means here the sum of the cross-sections of the
70 // pair production, Compton scattering and photoelectric processes
71
72 const G4double t1keV = 1.*CLHEP::keV;
73 const G4double t200keV = 200.*CLHEP::keV;
74 const G4double t100MeV = 100.*CLHEP::MeV;
75
76 G4double Zsquare = Z*Z;
78 G4double Zlogsquare = Zlog*Zlog;
79
80 G4double tmin = (0.552+218.5/Z+557.17/Zsquare)*CLHEP::MeV;
81 G4double tlow = 0.2*G4Exp(-7.355/std::sqrt(Z))*CLHEP::MeV;
82
83 G4double smin = (0.01239+0.005585*Zlog-0.000923*Zlogsquare)*G4Exp(1.5*Zlog);
84 G4double s200keV = (0.2651-0.1501*Zlog+0.02283*Zlogsquare)*Zsquare;
85
86 G4double cminlog = G4Log(tmin/t200keV);
87 G4double cmin = G4Log(s200keV/smin)/(cminlog*cminlog);
88
89 G4double slowlog = G4Log(t200keV/tlow);
90 G4double slow = s200keV * G4Exp(0.042*Z*slowlog*slowlog);
91 G4double logtlow = G4Log(tlow/t1keV);
92 G4double clow = G4Log(300.*Zsquare/slow)/logtlow;
93 G4double chigh = (7.55e-5 - 0.0542e-5*Z)*Zsquare*Z/G4Log(t100MeV/tmin);
94
95 // Calculate the cross-section (using an approximate empirical formula)
96 G4double xs;
97 if ( energy < tlow )
98 {
99 xs = (energy < t1keV) ? slow*G4Exp(clow*logtlow) :
100 slow*G4Exp(clow*G4Log(tlow/energy));
101 }
102 else if ( energy < t200keV )
103 {
104 G4double x = G4Log(t200keV/energy);
105 xs = s200keV * G4Exp(0.042*Z*x*x);
106 }
107 else if( energy<tmin )
108 {
109 const G4double x = G4Log(tmin/energy);
110 xs = smin * G4Exp(cmin*x*x);
111 }
112 else
113 {
114 xs = smin + chigh*G4Log(energy/tmin);
115 }
116 return xs * CLHEP::barn;
117}
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:179
G4double G4Log(G4double x)
Definition: G4Log.hh:226
double G4double
Definition: G4Types.hh:83
const G4int Z[17]
static const G4double tlow
static G4Pow * GetInstance()
Definition: G4Pow.cc:41
G4double logZ(G4int Z) const
Definition: G4Pow.hh:137
static constexpr double barn
Definition: SystemOfUnits.h:86
static constexpr double keV
static constexpr double MeV
G4double energy(const ThreeVector &p, const G4double m)

References CLHEP::barn, G4INCL::KinematicsUtils::energy(), G4Exp(), G4Log(), G4Pow::GetInstance(), CLHEP::keV, G4Pow::logZ(), CLHEP::MeV, tlow, and Z.

◆ Convert()

G4double G4VRangeToEnergyConverter::Convert ( const G4double  rangeCut,
const G4Material material 
)
virtualinherited

Reimplemented in G4RToEConvForProton.

Definition at line 87 of file G4VRangeToEnergyConverter.cc.

89{
90#ifdef G4VERBOSE
91 if (GetVerboseLevel()>3)
92 {
93 G4cout << "G4VRangeToEnergyConverter::Convert() - ";
94 G4cout << "Convert for " << material->GetName()
95 << " with Range Cut " << rangeCut/mm << "[mm]" << G4endl;
96 }
97#endif
98
99 G4double cut = 0.0;
100 if(fPDG == 22)
101 {
102 cut = ConvertForGamma(rangeCut, material);
103 }
104 else
105 {
106 cut = ConvertForElectron(rangeCut, material);
107
108 const G4double tune = 0.025*CLHEP::mm*CLHEP::g/CLHEP::cm3;
109 const G4double lowen = 30.*CLHEP::keV;
110 if(cut < lowen)
111 {
112 // corr. should be switched on smoothly
113 cut /= (1.+(1.-cut/lowen)*tune/(rangeCut*material->GetDensity()));
114 }
115 }
116
117 cut = std::max(Emin, std::min(cut, Emax));
118 return cut;
119}
static constexpr double mm
Definition: G4SIunits.hh:95
G4double ConvertForGamma(const G4double rangeCut, const G4Material *material)
G4double ConvertForElectron(const G4double rangeCut, const G4Material *material)
static constexpr double mm
Definition: SystemOfUnits.h:96
static constexpr double cm3
static constexpr double g
T max(const T t1, const T t2)
brief Return the largest of the two arguments
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
string material
Definition: eplot.py:19

References CLHEP::cm3, G4VRangeToEnergyConverter::ConvertForElectron(), G4VRangeToEnergyConverter::ConvertForGamma(), G4VRangeToEnergyConverter::Emax, G4VRangeToEnergyConverter::Emin, G4VRangeToEnergyConverter::fPDG, CLHEP::g, G4cout, G4endl, G4VRangeToEnergyConverter::GetVerboseLevel(), CLHEP::keV, eplot::material, G4INCL::Math::max(), G4INCL::Math::min(), CLHEP::mm, and mm.

Referenced by G4ProductionCutsTable::ConvertRangeToEnergy(), and G4ProductionCutsTable::UpdateCoupleTable().

◆ ConvertForElectron()

G4double G4VRangeToEnergyConverter::ConvertForElectron ( const G4double  rangeCut,
const G4Material material 
)
privateinherited

Definition at line 225 of file G4VRangeToEnergyConverter.cc.

227{
228 const G4ElementVector* elm = material->GetElementVector();
229 const G4double* dens = material->GetAtomicNumDensityVector();
230
231 // fill absorption length vector
232 G4int nelm = material->GetNumberOfElements();
233 G4double dedx1 = 0.0;
234 G4double dedx2 = 0.0;
235 G4double range1 = 0.0;
236 G4double range2 = 0.0;
237 G4double e1 = 0.0;
238 G4double e2 = 0.0;
239 G4double range = 0.;
240 for (G4int i=0; i<Nbin; ++i)
241 {
242 e2 = (*Energy)[i];
243 dedx2 = 0.0;
244 for (G4int j=0; j<nelm; ++j)
245 {
246 dedx2 += dens[j]*ComputeValue((*elm)[j]->GetZasInt(), e2);
247 }
248 range += (dedx1 + dedx2 > 0.0) ? 2*(e2 - e1)/(dedx1 + dedx2) : 0.0;
249 range2 = range;
250 if(range2 < rangeCut)
251 {
252 e1 = e2;
253 dedx1 = dedx2;
254 range1 = range2;
255 }
256 else
257 {
258 break;
259 }
260 }
261 return LiniearInterpolation(e1, e2, range1, range2, rangeCut);
262}
static const G4double e1[44]
static const G4double e2[44]
std::vector< const G4Element * > G4ElementVector
int G4int
Definition: G4Types.hh:85
virtual G4double ComputeValue(const G4int Z, const G4double kinEnergy)=0
G4double LiniearInterpolation(const G4double e1, const G4double e2, const G4double r1, const G4double r2, const G4double r)

References G4VRangeToEnergyConverter::ComputeValue(), e1, e2, G4VRangeToEnergyConverter::LiniearInterpolation(), eplot::material, and G4VRangeToEnergyConverter::Nbin.

Referenced by G4VRangeToEnergyConverter::Convert().

◆ ConvertForGamma()

G4double G4VRangeToEnergyConverter::ConvertForGamma ( const G4double  rangeCut,
const G4Material material 
)
privateinherited

Definition at line 188 of file G4VRangeToEnergyConverter.cc.

190{
191 const G4ElementVector* elm = material->GetElementVector();
192 const G4double* dens = material->GetAtomicNumDensityVector();
193
194 // fill absorption length vector
195 G4int nelm = material->GetNumberOfElements();
196 G4double range1 = 0.0;
197 G4double range2 = 0.0;
198 G4double e1 = 0.0;
199 G4double e2 = 0.0;
200 for (G4int i=0; i<Nbin; ++i)
201 {
202 e2 = (*Energy)[i];
203 G4double sig = 0.;
204
205 for (G4int j=0; j<nelm; ++j)
206 {
207 sig += dens[j]*ComputeValue((*elm)[j]->GetZasInt(), e2);
208 }
209 range2 = (sig > 0.0) ? 5./sig : DBL_MAX;
210 if(i == 0 || range2 < rangeCut)
211 {
212 e1 = e2;
213 range1 = range2;
214 }
215 else
216 {
217 break;
218 }
219 }
220 return LiniearInterpolation(e1, e2, range1, range2, rangeCut);
221}
#define DBL_MAX
Definition: templates.hh:62

References G4VRangeToEnergyConverter::ComputeValue(), DBL_MAX, e1, e2, G4VRangeToEnergyConverter::LiniearInterpolation(), eplot::material, and G4VRangeToEnergyConverter::Nbin.

Referenced by G4VRangeToEnergyConverter::Convert().

◆ FillEnergyVector()

void G4VRangeToEnergyConverter::FillEnergyVector ( const G4double  emin,
const G4double  emax 
)
staticprivateinherited

Definition at line 161 of file G4VRangeToEnergyConverter.cc.

163{
164 if(emin != Emin || emax != Emax)
165 {
166#ifdef G4MULTITHREADED
167 G4MUTEXLOCK(&theMutex);
168 if(emin != Emin || emax != Emax)
169 {
170#endif
171 Emin = emin;
172 Emax = emax;
173 Nbin = NbinPerDecade*static_cast<G4int>(std::log10(emax/emin));
174 Energy->resize(Nbin + 1);
175 (*Energy)[0] = emin;
176 (*Energy)[Nbin] = emax;
177 G4double fact = G4Log(emax/emin)/Nbin;
178 for(G4int i=1; i<Nbin; ++i) { (*Energy)[i] = emin*G4Exp(i * fact); }
179#ifdef G4MULTITHREADED
180 }
181 G4MUTEXUNLOCK(&theMutex);
182#endif
183 }
184}
static const G4double emax
#define G4MUTEXLOCK(mutex)
Definition: G4Threading.hh:251
#define G4MUTEXUNLOCK(mutex)
Definition: G4Threading.hh:254
static std::vector< G4double > * Energy

References G4VRangeToEnergyConverter::Emax, emax, G4VRangeToEnergyConverter::Emin, G4VRangeToEnergyConverter::Energy, G4Exp(), G4Log(), G4MUTEXLOCK, G4MUTEXUNLOCK, G4VRangeToEnergyConverter::Nbin, and G4VRangeToEnergyConverter::NbinPerDecade.

Referenced by G4VRangeToEnergyConverter::G4VRangeToEnergyConverter(), G4VRangeToEnergyConverter::SetEnergyRange(), and G4VRangeToEnergyConverter::SetMaxEnergyCut().

◆ GetHighEdgeEnergy()

G4double G4VRangeToEnergyConverter::GetHighEdgeEnergy ( )
staticinherited

Definition at line 139 of file G4VRangeToEnergyConverter.cc.

140{
141 return Emax;
142}

References G4VRangeToEnergyConverter::Emax.

Referenced by G4ProductionCutsTable::GetHighEdgeEnergy().

◆ GetLowEdgeEnergy()

G4double G4VRangeToEnergyConverter::GetLowEdgeEnergy ( )
staticinherited

Definition at line 133 of file G4VRangeToEnergyConverter.cc.

134{
135 return Emin;
136}

References G4VRangeToEnergyConverter::Emin.

Referenced by G4ProductionCutsTable::GetLowEdgeEnergy().

◆ GetMaxEnergyCut()

G4double G4VRangeToEnergyConverter::GetMaxEnergyCut ( )
staticinherited

Definition at line 146 of file G4VRangeToEnergyConverter.cc.

147{
148 return Emax;
149}

References G4VRangeToEnergyConverter::Emax.

Referenced by G4ProductionCutsTable::GetMaxEnergyCut().

◆ GetParticleType()

const G4ParticleDefinition * G4VRangeToEnergyConverter::GetParticleType ( ) const
inlineinherited

Definition at line 140 of file G4VRangeToEnergyConverter.hh.

141{
142 return theParticle;
143}

References G4VRangeToEnergyConverter::theParticle.

◆ GetVerboseLevel()

G4int G4VRangeToEnergyConverter::GetVerboseLevel ( ) const
inlineinherited

◆ LiniearInterpolation()

G4double G4VRangeToEnergyConverter::LiniearInterpolation ( const G4double  e1,
const G4double  e2,
const G4double  r1,
const G4double  r2,
const G4double  r 
)
inlineprivateinherited

Definition at line 145 of file G4VRangeToEnergyConverter.hh.

148{
149 return (r1 == r2) ? e1 : e1 + (e2 - e1)*(r - r1)/(r2 - r1);
150}

References e1, and e2.

Referenced by G4VRangeToEnergyConverter::ConvertForElectron(), and G4VRangeToEnergyConverter::ConvertForGamma().

◆ operator!=()

G4bool G4VRangeToEnergyConverter::operator!= ( const G4VRangeToEnergyConverter r) const
deleteinherited

◆ operator==()

G4bool G4VRangeToEnergyConverter::operator== ( const G4VRangeToEnergyConverter r) const
deleteinherited

◆ SetEnergyRange()

void G4VRangeToEnergyConverter::SetEnergyRange ( const G4double  lowedge,
const G4double  highedge 
)
staticinherited

Definition at line 122 of file G4VRangeToEnergyConverter.cc.

124{
125 G4double ehigh = std::min(Emax, highedge);
126 if(ehigh > lowedge)
127 {
128 FillEnergyVector(lowedge, ehigh);
129 }
130}
static void FillEnergyVector(const G4double emin, const G4double emax)

References G4VRangeToEnergyConverter::Emax, G4VRangeToEnergyConverter::FillEnergyVector(), and G4INCL::Math::min().

Referenced by G4ProductionCutsTable::SetEnergyRange().

◆ SetMaxEnergyCut()

void G4VRangeToEnergyConverter::SetMaxEnergyCut ( const G4double  value)
staticinherited

Definition at line 152 of file G4VRangeToEnergyConverter.cc.

153{
154 if(value > Emin)
155 {
156 FillEnergyVector(Emin, value);
157 }
158}

References G4VRangeToEnergyConverter::Emin, and G4VRangeToEnergyConverter::FillEnergyVector().

Referenced by G4ProductionCutsTable::SetMaxEnergyCut().

◆ SetVerboseLevel()

void G4VRangeToEnergyConverter::SetVerboseLevel ( G4int  value)
inlineinherited

Field Documentation

◆ Emax

G4double G4VRangeToEnergyConverter::Emax = 0.0
staticprotectedinherited

◆ Emin

G4double G4VRangeToEnergyConverter::Emin = 0.0
staticprotectedinherited

◆ Energy

std::vector< G4double > * G4VRangeToEnergyConverter::Energy = nullptr
staticprotectedinherited

◆ fPDG

G4int G4VRangeToEnergyConverter::fPDG = 0
protectedinherited

◆ isFirstInstance

G4bool G4VRangeToEnergyConverter::isFirstInstance = false
protectedinherited

◆ Nbin

G4int G4VRangeToEnergyConverter::Nbin = 350
staticprotectedinherited

◆ NbinPerDecade

G4int G4VRangeToEnergyConverter::NbinPerDecade = 50
staticprotectedinherited

◆ theParticle

const G4ParticleDefinition* G4VRangeToEnergyConverter::theParticle = nullptr
protectedinherited

◆ verboseLevel

G4int G4VRangeToEnergyConverter::verboseLevel = 1
protectedinherited

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