Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Static Public Member Functions | List of all members
G4KalbachCrossSection Class Reference

#include <G4KalbachCrossSection.hh>

Static Public Member Functions

static G4double ComputePowerParameter (G4int resA, G4int idx)
 
static G4double ComputeCrossSection (G4double K, G4double cb, G4double resA13, G4double amu1, G4int idx, G4int Z, G4int A, G4int resA)
 

Detailed Description

Definition at line 74 of file G4KalbachCrossSection.hh.

Member Function Documentation

static G4double G4KalbachCrossSection::ComputeCrossSection ( G4double  K,
G4double  cb,
G4double  resA13,
G4double  amu1,
G4int  idx,
G4int  Z,
G4int  A,
G4int  resA 
)
inlinestatic

Definition at line 83 of file G4KalbachCrossSection.hh.

87  {
88  G4double sig = 0.0;
89  G4double signor = 1.0;
90  G4double lambda, mu, nu;
91  G4double ec = 0.5;
92  if(0 < Z) { ec = cb; }
93  //JMQ 13.02.2009 tuning for improving cluster emission ddxs
94  // (spallation benchmark)
95  /*
96  G4double xx = 1.7;
97  if(1 == A) { xx = 1.5; }
98  ec = 1.44 * Z * resZ / (xx*resA13 + paramK[idx][10]);
99  }
100  */
101  G4double ecsq = ec*ec;
102  G4double elab = K * (A + resA) / G4double(resA);
103 
104  if(idx == 0) { // parameterization for neutron
105 
106  if(resA < 40) { signor =0.7 + resA*0.0075; }
107  else if(resA > 210) { signor = 1. + (resA-210)*0.004; }
108  lambda = paramK[idx][3]/resA13 + paramK[idx][4];
109  mu = (paramK[idx][5] + paramK[idx][6]*resA13)*resA13;
110  // JMQ 20.11.2008 very low energy behaviour corrected
111  // (problem for A (apprx.)>60) fix for avoiding
112  // neutron xs going to zero at very low energies
113  nu = std::abs((paramK[idx][7]*resA + paramK[idx][8]*resA13)*resA13
114  + paramK[idx][9]);
115 
116  } else { // parameterization for charged
117  // proton correction
118  if(idx == 1) {
119  if (resA <= 60) { signor = 0.92; }
120  else if (resA < 100) { signor = 0.8 + resA*0.002; }
121  }
122  lambda = paramK[idx][3]*resA + paramK[idx][4];
123  mu = paramK[idx][5]*amu1;
124  nu = amu1* (paramK[idx][7] + paramK[idx][8]*ec + paramK[idx][9]*ecsq);
125  }
126  /*
127  G4cout << "## idx= " << idx << " K= " << K << " elab= " << elab << " ec= " << ec
128  << " lambda= " << lambda << " mu= " << mu << " nu= " << nu << G4endl;
129  */
130  // threashold cross section
131  if(elab < ec) {
132  G4double p = paramK[idx][0];
133  if(0 < Z) { p += paramK[idx][1]/ec + paramK[idx][2]/ecsq; }
134  G4double a = -2*p*ec + lambda - nu/ecsq;
135  G4double b = p*ecsq + mu + 2*nu/ec;
136  G4double ecut;
137  G4double det = a*a - 4*p*b;
138  if (det > 0.0) { ecut = (std::sqrt(det) - a)/(2*p); }
139  else { ecut = -a/(2*p); }
140 
141  //G4cout << " elab= " << elab << " ecut= " << ecut << " sig= " << sig
142  // << " sig1= " << (p*elab*elab + a*elab + b)*signor << G4endl;
143  // If ecut>0, sig=0 at elab=ecut
144  if(0 == idx) {
145  sig = (lambda*ec + mu + nu/ec)*signor*std::sqrt(elab/ec);
146  } else if(elab >= ecut) {
147  sig = (p*elab*elab + a*elab + b)*signor;
148 
149  // extra proton correction
150  if(1 == idx) {
151  // c and w are for global correction factor for
152  // they are scaled down for light targets where ec is low.
153  G4double cc = std::min(3.15, ec*0.5);
154  G4double signor2 = (ec - elab - cc) *3.15/ (0.7*cc);
155  sig /= (1. + G4Exp(signor2));
156  }
157  }
158  //G4cout << " ecut= " << ecut << " a= " << a << " b= " << b
159  // << " signor= " << signor << " sig= " << sig << G4endl;
160 
161  // high energy cross section
162  } else {
163  // etest is the energy above which the rxn cross section is
164  // compared with the geometrical limit and the max taken.
165 
166  // neutron parameters
167  G4double etest = 32.;
168  G4double xnulam = 1.0;
169 
170  // parameters for charged
171  static const G4double flow = 1.e-18;
172  static const G4double spill= 1.e+18;
173  if(0 < Z) {
174  etest = 0.0;
175  xnulam = nu / lambda;
176  xnulam = std::min(xnulam, spill);
177  if (xnulam >= flow) {
178  if(1 == idx) { etest = std::sqrt(xnulam) + 7.; }
179  else { etest = 1.2 *std::sqrt(xnulam); }
180  }
181  }
182  // ** For xnulam.gt.0, sig reaches a maximum at sqrt(xnulam).
183  sig = (lambda*elab + mu + nu/elab)*signor;
184  if (xnulam >= flow && elab >= etest) {
185  G4double geom = std::sqrt(A*K);
186  geom = 1.23*resA13 + paramK[idx][10] + 4.573/geom;
187  geom = 31.416 * geom * geom;
188  sig = std::max(sig, geom);
189  }
190  }
191  sig = std::max(sig, 0.0);
192  //G4cout << " ---- sig= " << sig << G4endl;
193  return sig;
194  }
a
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
p
const char * p
Definition: xmltok.h:285
test.b
tuple b
Definition: test.py:12
A
double A(double temperature)
Definition: G4DNAElectronHoleRecombination.cc:60
G4Exp
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
G4INCL::Math::max
T max(const T t1, const T t2)
brief Return the largest of the two arguments
Definition: G4INCLGlobals.hh:112
G4INCL::Math::min
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
Definition: G4INCLGlobals.hh:117
G4double
double G4double
Definition: G4Types.hh:76
paramK
static const G4double paramK[6][11]
Definition: G4KalbachCrossSection.hh:58
Z
G4int Z
Definition: G4ParticleHPJENDLHEData.cc:262

Here is the call graph for this function:

Here is the caller graph for this function:

static G4double G4KalbachCrossSection::ComputePowerParameter ( G4int  resA,
G4int  idx 
)
inlinestatic

Definition at line 78 of file G4KalbachCrossSection.hh.

79  {
80  return G4Pow::GetInstance()->powZ(resA, paramK[idx][6]);
81  }
G4Pow::GetInstance
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
G4Pow::powZ
G4double powZ(G4int Z, G4double y) const
Definition: G4Pow.hh:254
paramK
static const G4double paramK[6][11]
Definition: G4KalbachCrossSection.hh:58

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 file: