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

#include <G4RDeBremsstrahlungSpectrum.hh>

Inheritance diagram for G4RDeBremsstrahlungSpectrum:
Collaboration diagram for G4RDeBremsstrahlungSpectrum:

Public Member Functions

 G4RDeBremsstrahlungSpectrum (const G4DataVector &bins, const G4String &name)
 
 ~G4RDeBremsstrahlungSpectrum ()
 
G4double Probability (G4int Z, G4double tMin, G4double tMax, G4double kineticEnergy, G4int shell=0, const G4ParticleDefinition *pd=0) const
 
G4double AverageEnergy (G4int Z, G4double tMin, G4double tMax, G4double kineticEnergy, G4int shell=0, const G4ParticleDefinition *pd=0) const
 
G4double SampleEnergy (G4int Z, G4double tMin, G4double tMax, G4double kineticEnergy, G4int shell=0, const G4ParticleDefinition *pd=0) const
 
G4double MaxEnergyOfSecondaries (G4double kineticEnergy, G4int Z=0, const G4ParticleDefinition *pd=0) const
 
G4double Excitation (G4int Z, G4double kineticEnergy) const
 
void PrintData () const
 
- Public Member Functions inherited from G4RDVEnergySpectrum
 G4RDVEnergySpectrum ()
 
virtual ~G4RDVEnergySpectrum ()
 

Detailed Description

Definition at line 66 of file G4RDeBremsstrahlungSpectrum.hh.

Constructor & Destructor Documentation

G4RDeBremsstrahlungSpectrum::G4RDeBremsstrahlungSpectrum ( const G4DataVector bins,
const G4String name 
)

Definition at line 56 of file G4RDeBremsstrahlungSpectrum.cc.

58  lowestE(0.1*eV),
59  xp(bins)
60 {
61  length = xp.size();
62  theBRparam = new G4RDBremsstrahlungParameters(name,length+1);
63 
64  verbose = 0;
65 }
static constexpr double eV
Definition: G4SIunits.hh:215
G4RDeBremsstrahlungSpectrum::~G4RDeBremsstrahlungSpectrum ( )

Definition at line 68 of file G4RDeBremsstrahlungSpectrum.cc.

69 {
70  delete theBRparam;
71 }

Member Function Documentation

G4double G4RDeBremsstrahlungSpectrum::AverageEnergy ( G4int  Z,
G4double  tMin,
G4double  tMax,
G4double  kineticEnergy,
G4int  shell = 0,
const G4ParticleDefinition pd = 0 
) const
virtual

Implements G4RDVEnergySpectrum.

Definition at line 121 of file G4RDeBremsstrahlungSpectrum.cc.

127 {
128  G4double tm = std::min(tmax, e);
129  G4double t0 = std::max(tmin, lowestE);
130  if(t0 >= tm) return 0.0;
131 
132  t0 /= e;
133  tm /= e;
134 
135  G4double z0 = lowestE/e;
136 
137  G4DataVector p;
138 
139  // Access parameters
140  for (size_t i=0; i<=length; i++) {
141  p.push_back(theBRparam->Parameter(i, Z, e));
142  }
143 
144  G4double x = AverageValue(t0, tm, p);
145  G4double y = IntSpectrum(z0, 1.0, p);
146 
147  // Add integrant over lowest energies
148  G4double zmin = tmin/e;
149  if(zmin < t0) {
150  G4double c = std::sqrt(theBRparam->ParameterC(Z));
151  x += p[0]*(t0 - zmin - c*(std::atan(t0/c) - std::atan(zmin/c)));
152  }
153  x *= e;
154 
155  if(1 < verbose) {
156  G4cout << "tcut(MeV)= " << tmin/MeV
157  << "; tMax(MeV)= " << tmax/MeV
158  << "; e(MeV)= " << e/MeV
159  << "; t0= " << t0
160  << "; tm= " << tm
161  << "; y= " << y
162  << "; x= " << x
163  << G4endl;
164  }
165  p.clear();
166 
167  if(y > 0.0) x /= y;
168  else x = 0.0;
169  // if(x < 0.0) x = 0.0;
170 
171  return x;
172 }
const char * p
Definition: xmltok.h:285
tuple x
Definition: test.py:50
G4GLOB_DLL std::ostream G4cout
G4double ParameterC(G4int index) const
G4double Parameter(G4int parameterIndex, G4int Z, G4double energy) const
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
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
double G4double
Definition: G4Types.hh:76
tuple c
Definition: test.py:13

Here is the call graph for this function:

G4double G4RDeBremsstrahlungSpectrum::Excitation ( G4int  Z,
G4double  kineticEnergy 
) const
virtual

Implements G4RDVEnergySpectrum.

Definition at line 303 of file G4RDeBremsstrahlungSpectrum.cc.

304 {
305  return 0.0;
306 }
G4double G4RDeBremsstrahlungSpectrum::MaxEnergyOfSecondaries ( G4double  kineticEnergy,
G4int  Z = 0,
const G4ParticleDefinition pd = 0 
) const
virtual

Implements G4RDVEnergySpectrum.

Definition at line 308 of file G4RDeBremsstrahlungSpectrum.cc.

311 {
312  return kineticEnergy;
313 }
void G4RDeBremsstrahlungSpectrum::PrintData ( void  ) const
virtual

Implements G4RDVEnergySpectrum.

Definition at line 300 of file G4RDeBremsstrahlungSpectrum.cc.

301 { theBRparam->PrintData(); }

Here is the call graph for this function:

G4double G4RDeBremsstrahlungSpectrum::Probability ( G4int  Z,
G4double  tMin,
G4double  tMax,
G4double  kineticEnergy,
G4int  shell = 0,
const G4ParticleDefinition pd = 0 
) const
virtual

Implements G4RDVEnergySpectrum.

Definition at line 74 of file G4RDeBremsstrahlungSpectrum.cc.

80 {
81  G4double tm = std::min(tmax, e);
82  G4double t0 = std::max(tmin, lowestE);
83  if(t0 >= tm) return 0.0;
84 
85  t0 /= e;
86  tm /= e;
87 
88  G4double z0 = lowestE/e;
90 
91  // Access parameters
92  for (size_t i=0; i<=length; i++) {
93  p.push_back(theBRparam->Parameter(i, Z, e));
94  }
95 
96  G4double x = IntSpectrum(t0, tm, p);
97  G4double y = IntSpectrum(z0, 1.0, p);
98 
99 
100  if(1 < verbose) {
101  G4cout << "tcut(MeV)= " << tmin/MeV
102  << "; tMax(MeV)= " << tmax/MeV
103  << "; t0= " << t0
104  << "; tm= " << tm
105  << "; xp[0]= " << xp[0]
106  << "; z= " << z0
107  << "; val= " << x
108  << "; nor= " << y
109  << G4endl;
110  }
111  p.clear();
112 
113  if(y > 0.0) x /= y;
114  else x = 0.0;
115  // if(x < 0.0) x = 0.0;
116 
117  return x;
118 }
const char * p
Definition: xmltok.h:285
tuple x
Definition: test.py:50
G4GLOB_DLL std::ostream G4cout
G4double Parameter(G4int parameterIndex, G4int Z, G4double energy) const
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
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4double G4RDeBremsstrahlungSpectrum::SampleEnergy ( G4int  Z,
G4double  tMin,
G4double  tMax,
G4double  kineticEnergy,
G4int  shell = 0,
const G4ParticleDefinition pd = 0 
) const
virtual

Implements G4RDVEnergySpectrum.

Definition at line 175 of file G4RDeBremsstrahlungSpectrum.cc.

181 {
182  G4double tm = std::min(tmax, e);
183  G4double t0 = std::max(tmin, lowestE);
184  if(t0 >= tm) return 0.0;
185 
186  t0 /= e;
187  tm /= e;
188 
189  G4DataVector p;
190 
191  for (size_t i=0; i<=length; i++) {
192  p.push_back(theBRparam->Parameter(i, Z, e));
193  }
194  G4double amaj = std::max(p[length], 1. - (p[1] - p[0])*xp[0]/(xp[1] - xp[0]) );
195 
196  G4double amax = std::log(tm);
197  G4double amin = std::log(t0);
198  G4double tgam, q, fun;
199 
200  do {
201  G4double x = amin + G4UniformRand()*(amax - amin);
202  tgam = std::exp(x);
203  fun = Function(tgam, p);
204 
205  if(fun > amaj) {
206  G4cout << "WARNING in G4RDeBremsstrahlungSpectrum::SampleEnergy:"
207  << " Majoranta " << amaj
208  << " < " << fun
209  << G4endl;
210  }
211 
212  q = amaj * G4UniformRand();
213  } while (q > fun);
214 
215  tgam *= e;
216 
217  p.clear();
218 
219  return tgam;
220 }
const char * p
Definition: xmltok.h:285
tuple x
Definition: test.py:50
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
G4double Parameter(G4int parameterIndex, G4int Z, G4double energy) const
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
#define G4endl
Definition: G4ios.hh:61
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: