Geant4.10
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4HadronNucleonXsc.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // 14.03.07 V. Grichine - first implementation
27 //
28 
29 #include "G4HadronNucleonXsc.hh"
30 
31 #include "G4PhysicalConstants.hh"
32 #include "G4SystemOfUnits.hh"
33 #include "G4ParticleTable.hh"
34 #include "G4IonTable.hh"
35 #include "G4ParticleDefinition.hh"
36 #include "G4HadTmpUtil.hh"
37 
38 
39 G4HadronNucleonXsc::G4HadronNucleonXsc()
40 :
41 // fUpperLimit( 10000 * GeV ),
42  fLowerLimit( 0.03 * MeV ),
43  fTotalXsc(0.0), fElasticXsc(0.0), fInelasticXsc(0.0)
44 // , fHadronNucleonXsc(0.0)
45 {
46  theGamma = G4Gamma::Gamma();
47  theProton = G4Proton::Proton();
48  theNeutron = G4Neutron::Neutron();
49  theAProton = G4AntiProton::AntiProton();
50  theANeutron = G4AntiNeutron::AntiNeutron();
51  thePiPlus = G4PionPlus::PionPlus();
52  thePiMinus = G4PionMinus::PionMinus();
53  thePiZero = G4PionZero::PionZero();
54  theKPlus = G4KaonPlus::KaonPlus();
55  theKMinus = G4KaonMinus::KaonMinus();
56  theK0S = G4KaonZeroShort::KaonZeroShort();
57  theK0L = G4KaonZeroLong::KaonZeroLong();
58  theL = G4Lambda::Lambda();
59  theAntiL = G4AntiLambda::AntiLambda();
60  theSPlus = G4SigmaPlus::SigmaPlus();
61  theASPlus = G4AntiSigmaPlus::AntiSigmaPlus();
62  theSMinus = G4SigmaMinus::SigmaMinus();
63  theASMinus = G4AntiSigmaMinus::AntiSigmaMinus();
64  theS0 = G4SigmaZero::SigmaZero();
65  theAS0 = G4AntiSigmaZero::AntiSigmaZero();
66  theXiMinus = G4XiMinus::XiMinus();
67  theXi0 = G4XiZero::XiZero();
68  theAXiMinus = G4AntiXiMinus::AntiXiMinus();
69  theAXi0 = G4AntiXiZero::AntiXiZero();
70  theOmega = G4OmegaMinus::OmegaMinus();
71  theAOmega = G4AntiOmegaMinus::AntiOmegaMinus();
72  theD = G4Deuteron::Deuteron();
73  theT = G4Triton::Triton();
74  theA = G4Alpha::Alpha();
75  theHe3 = G4He3::He3();
76 
77  InitialiseKaonNucleonTotXsc();
78 }
79 
80 
81 G4HadronNucleonXsc::~G4HadronNucleonXsc()
82 {}
83 
84 void G4HadronNucleonXsc::CrossSectionDescription(std::ostream& outFile) const
85 {
86  outFile << "G4HadronNucleonXsc calculates the total, inelastic and elastic\n"
87  << "hadron-nucleon cross sections using several different\n"
88  << "parameterizations within the Glauber-Gribov approach. It is\n"
89  << "valid for all incident gammas and long-lived hadrons at\n"
90  << "energies above 30 keV. This is a cross section component which\n"
91  << "is to be used to build a cross section data set.\n";
92 }
93 
94 G4bool
95 G4HadronNucleonXsc::IsApplicable(const G4DynamicParticle* aDP,
96  const G4Element* anElement)
97 {
98  G4int Z = G4lrint(anElement->GetZ());
99  G4int A = G4lrint(anElement->GetN());
100  return IsIsoApplicable(aDP, Z, A);
101 }
102 
103 ////////////////////////////////////////////////////////////////////////////////////////
104 //
105 
106 G4bool
107 G4HadronNucleonXsc::IsIsoApplicable(const G4DynamicParticle* aDP,
108  G4int Z, G4int)
109 {
110  G4bool applicable = false;
111  // G4int baryonNumber = aDP->GetDefinition()->GetBaryonNumber();
112  G4double kineticEnergy = aDP->GetKineticEnergy();
113 
114  const G4ParticleDefinition* theParticle = aDP->GetDefinition();
115 
116  if ( ( kineticEnergy >= fLowerLimit &&
117  Z > 1 && // >= He
118  ( theParticle == theAProton ||
119  theParticle == theGamma ||
120  theParticle == theKPlus ||
121  theParticle == theKMinus ||
122  theParticle == theSMinus) ) ||
123 
124  ( kineticEnergy >= 0.1*fLowerLimit &&
125  Z > 1 && // >= He
126  ( theParticle == theProton ||
127  theParticle == theNeutron ||
128  theParticle == thePiPlus ||
129  theParticle == thePiMinus ) ) ) applicable = true;
130 
131  return applicable;
132 }
133 
134 
135 /////////////////////////////////////////////////////////////////////////////////////
136 //
137 // Returns hadron-nucleon Xsc according to differnt parametrisations:
138 // [2] E. Levin, hep-ph/9710546
139 // [3] U. Dersch, et al, hep-ex/9910052
140 // [4] M.J. Longo, et al, Phys.Rev.Lett. 33 (1974) 725
141 
142 G4double
143 G4HadronNucleonXsc::GetHadronNucleonXscEL(const G4DynamicParticle* aParticle,
144  const G4ParticleDefinition* nucleon )
145 {
146  G4double xsection;
147 
148 
149  G4double targ_mass = 0.939*GeV; // ~mean neutron and proton ???
150 
151  G4double proj_mass = aParticle->GetMass();
152  G4double proj_momentum = aParticle->GetMomentum().mag();
153  G4double sMand = CalcMandelstamS ( proj_mass , targ_mass , proj_momentum );
154 
155  sMand /= GeV*GeV; // in GeV for parametrisation
156  proj_momentum /= GeV;
157 
158  const G4ParticleDefinition* theParticle = aParticle->GetDefinition();
159 
160  G4bool pORn = (nucleon == theProton || nucleon == theNeutron );
161 
162 
163  if(theParticle == theGamma && pORn )
164  {
165  xsection = (0.0677*std::pow(sMand,0.0808) + 0.129*std::pow(sMand,-0.4525));
166  }
167  else if(theParticle == theNeutron && pORn ) // as proton ???
168  {
169  xsection = (21.70*std::pow(sMand,0.0808) + 56.08*std::pow(sMand,-0.4525));
170  }
171  else if(theParticle == theProton && pORn )
172  {
173  xsection = (21.70*std::pow(sMand,0.0808) + 56.08*std::pow(sMand,-0.4525));
174 
175  // xsection = At*( 49.51*std::pow(sMand,-0.097) + 0.314*std::log(sMand)*std::log(sMand) );
176  // xsection = At*( 38.4 + 0.85*std::abs(std::pow(log(sMand),1.47)) );
177  }
178  else if(theParticle == theAProton && pORn )
179  {
180  xsection = ( 21.70*std::pow(sMand,0.0808) + 98.39*std::pow(sMand,-0.4525));
181  }
182  else if(theParticle == thePiPlus && pORn )
183  {
184  xsection = (13.63*std::pow(sMand,0.0808) + 27.56*std::pow(sMand,-0.4525));
185  }
186  else if(theParticle == thePiMinus && pORn )
187  {
188  // xsection = At*( 55.2*std::pow(sMand,-0.255) + 0.346*std::log(sMand)*std::log(sMand) );
189  xsection = (13.63*std::pow(sMand,0.0808) + 36.02*std::pow(sMand,-0.4525));
190  }
191  else if(theParticle == theKPlus && pORn )
192  {
193  xsection = (11.82*std::pow(sMand,0.0808) + 8.15*std::pow(sMand,-0.4525));
194  }
195  else if(theParticle == theKMinus && pORn )
196  {
197  xsection = (11.82*std::pow(sMand,0.0808) + 26.36*std::pow(sMand,-0.4525));
198  }
199  else // as proton ???
200  {
201  xsection = (21.70*std::pow(sMand,0.0808) + 56.08*std::pow(sMand,-0.4525));
202  }
203  xsection *= millibarn;
204 
205  fTotalXsc = xsection;
206  fInelasticXsc = 0.83*xsection;
207  fElasticXsc = fTotalXsc - fInelasticXsc;
208  if (fElasticXsc < 0.)fElasticXsc = 0.;
209 
210  return xsection;
211 }
212 
213 
214 /////////////////////////////////////////////////////////////////////////////////////
215 //
216 // Returns hadron-nucleon Xsc according to PDG parametrisation (2005):
217 // http://pdg.lbl.gov/2006/reviews/hadronicrpp.pdf
218 // At = number of nucleons, Zt = number of protons
219 
220 G4double
221 G4HadronNucleonXsc::GetHadronNucleonXscPDG(const G4DynamicParticle* aParticle,
222  const G4ParticleDefinition* nucleon )
223 {
224  G4double xsection(0);
225  G4int Zt=1, Nt=1, At=1;
226 
227  G4double targ_mass = 0.939*GeV; // ~mean neutron and proton ???
228 
229  G4double proj_mass = aParticle->GetMass();
230  G4double proj_momentum = aParticle->GetMomentum().mag();
231 
232  G4double sMand = CalcMandelstamS ( proj_mass , targ_mass , proj_momentum );
233 
234  sMand /= GeV*GeV; // in GeV for parametrisation
235 
236  // General PDG fit constants
237 
238  G4double s0 = 5.38*5.38; // in Gev^2
239  G4double eta1 = 0.458;
240  G4double eta2 = 0.458;
241  G4double B = 0.308;
242 
243  const G4ParticleDefinition* theParticle = aParticle->GetDefinition();
244 
245  G4bool pORn = (nucleon == theProton || nucleon == theNeutron );
246  G4bool proton = (nucleon == theProton);
247  G4bool neutron = (nucleon == theNeutron);
248 
249  if(theParticle == theNeutron) // proton-neutron fit
250  {
251  if ( proton )
252  {
253  xsection = Zt*( 35.80 + B*std::pow(std::log(sMand/s0),2.)
254  + 40.15*std::pow(sMand,-eta1) - 30.*std::pow(sMand,-eta2));// on p
255  }
256  if ( neutron )
257  {
258  xsection = Nt*( 35.45 + B*std::pow(std::log(sMand/s0),2.)
259  + 42.53*std::pow(sMand,-eta1) - 33.34*std::pow(sMand,-eta2)); // on n pp for nn
260  }
261  }
262  else if(theParticle == theProton)
263  {
264  if ( proton )
265  {
266  xsection = Zt*( 35.45 + B*std::pow(std::log(sMand/s0),2.)
267  + 42.53*std::pow(sMand,-eta1) - 33.34*std::pow(sMand,-eta2));
268  }
269  if ( neutron )
270  {
271  xsection = Nt*( 35.80 + B*std::pow(std::log(sMand/s0),2.)
272  + 40.15*std::pow(sMand,-eta1) - 30.*std::pow(sMand,-eta2));
273  }
274  }
275  else if(theParticle == theAProton)
276  {
277  if ( proton )
278  {
279  xsection = Zt*( 35.45 + B*std::pow(std::log(sMand/s0),2.)
280  + 42.53*std::pow(sMand,-eta1) + 33.34*std::pow(sMand,-eta2));
281  }
282  if ( neutron )
283  {
284  xsection = Nt*( 35.80 + B*std::pow(std::log(sMand/s0),2.)
285  + 40.15*std::pow(sMand,-eta1) + 30.*std::pow(sMand,-eta2));
286  }
287  }
288  else if(theParticle == thePiPlus && pORn )
289  {
290  xsection = At*( 20.86 + B*std::pow(std::log(sMand/s0),2.)
291  + 19.24*std::pow(sMand,-eta1) - 6.03*std::pow(sMand,-eta2));
292  }
293  else if(theParticle == thePiMinus && pORn )
294  {
295  xsection = At*( 20.86 + B*std::pow(std::log(sMand/s0),2.)
296  + 19.24*std::pow(sMand,-eta1) + 6.03*std::pow(sMand,-eta2));
297  }
298  else if(theParticle == theKPlus)
299  {
300  if ( proton )
301  {
302  xsection = Zt*( 17.91 + B*std::pow(std::log(sMand/s0),2.)
303  + 7.14*std::pow(sMand,-eta1) - 13.45*std::pow(sMand,-eta2));
304  }
305  if ( neutron )
306  {
307  xsection = Nt*( 17.87 + B*std::pow(std::log(sMand/s0),2.)
308  + 5.17*std::pow(sMand,-eta1) - 7.23*std::pow(sMand,-eta2));
309  }
310  }
311  else if(theParticle == theKMinus)
312  {
313  if ( proton )
314  {
315  xsection = Zt*( 17.91 + B*std::pow(std::log(sMand/s0),2.)
316  + 7.14*std::pow(sMand,-eta1) + 13.45*std::pow(sMand,-eta2));
317  }
318  if ( neutron )
319  {
320  xsection = Nt*( 17.87 + B*std::pow(std::log(sMand/s0),2.)
321  + 5.17*std::pow(sMand,-eta1) + 7.23*std::pow(sMand,-eta2) );
322  }
323  }
324  else if(theParticle == theSMinus && pORn )
325  {
326  xsection = At*( 35.20 + B*std::pow(std::log(sMand/s0),2.)
327  - 199.*std::pow(sMand,-eta1) + 264.*std::pow(sMand,-eta2) );
328  }
329  else if(theParticle == theGamma && pORn ) // modify later on
330  {
331  xsection = At*( 0.0 + B*std::pow(std::log(sMand/s0),2.)
332  + 0.032*std::pow(sMand,-eta1) - 0.0*std::pow(sMand,-eta2) );
333 
334  }
335  else // as proton ???
336  {
337  if ( proton )
338  {
339  xsection = Zt*( 35.45 + B*std::pow(std::log(sMand/s0),2.)
340  + 42.53*std::pow(sMand,-eta1) - 33.34*std::pow(sMand,-eta2) );
341  }
342  if ( neutron )
343  {
344  xsection = Nt*( 35.80 + B*std::pow(std::log(sMand/s0),2.)
345  + 40.15*std::pow(sMand,-eta1) - 30.*std::pow(sMand,-eta2));
346  }
347  }
348  xsection *= millibarn; // parametrised in mb
349 
350  fTotalXsc = xsection;
351  fInelasticXsc = 0.75*xsection;
352  fElasticXsc = fTotalXsc - fInelasticXsc;
353  if (fElasticXsc < 0.) fElasticXsc = 0.;
354 
355  return xsection;
356 }
357 
358 
359 /////////////////////////////////////////////////////////////////////////////////////
360 //
361 // Returns hadron-nucleon cross-section based on N. Starkov parametrisation of
362 // data from mainly http://wwwppds.ihep.su:8001/c5-6A.html database
363 
364 G4double
365 G4HadronNucleonXsc::GetHadronNucleonXscNS(const G4DynamicParticle* aParticle,
366  const G4ParticleDefinition* nucleon )
367 {
368  G4double xsection(0);
369 
370  G4double A0, B0;
371  G4double hpXsc(0);
372  G4double hnXsc(0);
373 
374 
375  G4double tM = 0.939*GeV; // ~mean neutron and proton ???
376 
377  G4double pM = aParticle->GetMass();
378  G4double pE = aParticle->GetTotalEnergy();
379  G4double pLab = aParticle->GetMomentum().mag();
380 
381  G4double sMand = CalcMandelstamS ( pM , tM , pLab );
382 
383  sMand /= GeV*GeV; // in GeV for parametrisation
384  pLab /= GeV;
385  pE /= GeV;
386  pM /= GeV;
387 
388  G4double logP = std::log(pLab);
389 
390 
391  // General PDG fit constants
392 
393  G4double s0 = 5.38*5.38; // in Gev^2
394  G4double eta1 = 0.458;
395  G4double eta2 = 0.458;
396  G4double B = 0.308;
397  G4double minLogP = 3.5; // min of (lnP-minLogP)^2
398  G4double cofLogE = .0557; // elastic (lnP-minLogP)^2
399  G4double cofLogT = .3; // total (lnP-minLogP)^2
400  G4double pMin = .1; // fast LE calculation
401  G4double pMax = 1000.; // fast HE calculation
402 
403 
404  const G4ParticleDefinition* theParticle = aParticle->GetDefinition();
405 
406  G4bool pORn = (nucleon == theProton || nucleon == theNeutron );
407  G4bool proton = (nucleon == theProton);
408  G4bool neutron = (nucleon == theNeutron);
409 
410  if( theParticle == theNeutron && pORn )
411  {
412  if( pLab >= 373.)
413  {
414  xsection = GetHadronNucleonXscPDG(aParticle, nucleon)/millibarn;
415 
416  fElasticXsc = 6.5 + 0.308*std::pow(std::log(sMand/400.),1.65) + 9.19*std::pow(sMand,-0.458);
417 
418  fTotalXsc = xsection;
419 
420  }
421  else if( pLab >= 100.)
422  {
423  B0 = 7.5;
424  A0 = 100. - B0*std::log(3.0e7);
425 
426  xsection = A0 + B0*std::log(pE) - 11
427  // + 103*std::pow(2*0.93827*pE + pM*pM+0.93827*0.93827,-0.165); // mb
428  + 103*std::pow(sMand,-0.165); // mb
429 
430  fElasticXsc = 5.53 + 0.308*std::pow(std::log(sMand/28.9),1.1) + 9.19*std::pow(sMand,-0.458);
431 
432  fTotalXsc = xsection;
433  }
434  else if( pLab >= 10.)
435  {
436  B0 = 7.5;
437  A0 = 100. - B0*std::log(3.0e7);
438 
439  xsection = A0 + B0*std::log(pE) - 11
440  + 103*std::pow(2*0.93827*pE + pM*pM+
441  0.93827*0.93827,-0.165); // mb
442  fTotalXsc = xsection;
443  fElasticXsc = 6 + 20/( (logP-0.182)*(logP-0.182) + 1.0 );
444  }
445  else // pLab < 10 GeV/c
446  {
447  if( neutron ) // nn to be pp
448  {
449  if( pLab < 0.4 )
450  {
451  hnXsc = 23 + 50*( std::pow( std::log(0.73/pLab), 3.5 ) );
452  fElasticXsc = hnXsc;
453  }
454  else if( pLab < 0.73 )
455  {
456  hnXsc = 23 + 50*( std::pow( std::log(0.73/pLab), 3.5 ) );
457  fElasticXsc = hnXsc;
458  }
459  else if( pLab < 1.05 )
460  {
461  hnXsc = 23 + 40*(std::log(pLab/0.73))*
462  (std::log(pLab/0.73));
463  fElasticXsc = 23 + 20*(std::log(pLab/0.73))*
464  (std::log(pLab/0.73));
465  }
466  else // 1.05 - 10 GeV/c
467  {
468  hnXsc = 39.0+75*(pLab - 1.2)/(std::pow(pLab,3.0) + 0.15);
469 
470  fElasticXsc = 6 + 20/( (logP-0.182)*(logP-0.182) + 1.0 );
471  }
472  fTotalXsc = hnXsc;
473  }
474  if( proton ) // pn to be np
475  {
476  if( pLab < 0.02 )
477  {
478  hpXsc = 4100+30*std::pow(std::log(1.3/pLab),3.6); // was as pLab < 0.8
479  fElasticXsc = hpXsc;
480  }
481  else if( pLab < 0.8 )
482  {
483  hpXsc = 33+30*std::pow(std::log(pLab/1.3),4.0);
484  fElasticXsc = hpXsc;
485  }
486  else if( pLab < 1.05 )
487  {
488  hpXsc = 33+30*std::pow(std::log(pLab/0.95),2.0);
489  fElasticXsc = 6 + 52/( std::log(0.511/pLab)*std::log(0.511/pLab) + 1.6 );
490  }
491  else if( pLab < 1.4 )
492  {
493  hpXsc = 33+30*std::pow(std::log(pLab/0.95),2.0);
494  fElasticXsc = 6 + 52/( std::log(0.511/pLab)*std::log(0.511/pLab) + 1.6 );
495  }
496  else // 1.4 < pLab < 10. )
497  {
498  hpXsc = 33.3 + 20.8*(std::pow(pLab,2.0) - 1.35)/(std::pow(pLab,2.50) + 0.95);
499 
500  fElasticXsc = 6 + 20/( (logP-0.182)*(logP-0.182) + 1.0 );
501  }
502  fTotalXsc = hpXsc;
503  }
504  }
505  }
506  else if( theParticle == theProton && pORn ) ////// proton //////////////////////////////////////////////
507  {
508  if( pLab >= 373.) // pdg due to TOTEM data
509  {
510  xsection = GetHadronNucleonXscPDG(aParticle, nucleon)/millibarn;
511 
512  fElasticXsc = 6.5 + 0.308*std::pow(std::log(sMand/400.),1.65) + 9.19*std::pow(sMand,-0.458);
513 
514  fTotalXsc = xsection;
515  }
516  else if( pLab >= 100.)
517  {
518  B0 = 7.5;
519  A0 = 100. - B0*std::log(3.0e7);
520 
521  xsection = A0 + B0*std::log(pE) - 11 + 103*std::pow(sMand,-0.165); // mb
522 
523  fElasticXsc = 5.53 + 0.308*std::pow(std::log(sMand/28.9),1.1) + 9.19*std::pow(sMand,-0.458);
524 
525  fTotalXsc = xsection;
526  }
527  else if( pLab >= 10.)
528  {
529  B0 = 7.5;
530  A0 = 100. - B0*std::log(3.0e7);
531 
532  xsection = A0 + B0*std::log(pE) - 11 + 103*std::pow(sMand,-0.165); // mb
533 
534  fElasticXsc = 6 + 20/( (logP-0.182)*(logP-0.182) + 1.0 );
535 
536  fTotalXsc = xsection;
537  }
538  else
539  {
540  // pp
541 
542  if( proton )
543  {
544  if( pLab < 0.4 )
545  {
546  hpXsc = 23 + 50*( std::pow( std::log(0.73/pLab), 3.5 ) );
547  fElasticXsc = hpXsc;
548  }
549  else if( pLab < 0.73 )
550  {
551  hpXsc = 23 + 50*( std::pow( std::log(0.73/pLab), 3.5 ) );
552  fElasticXsc = hpXsc;
553  }
554  else if( pLab < 1.05 )
555  {
556  hpXsc = 23 + 40*(std::log(pLab/0.73))*
557  (std::log(pLab/0.73));
558  fElasticXsc = 23 + 20*(std::log(pLab/0.73))*
559  (std::log(pLab/0.73));
560  }
561  else // 1.05 - 10 GeV/c
562  {
563  hpXsc = 39.0+75*(pLab - 1.2)/(std::pow(pLab,3.0) + 0.15);
564 
565  fElasticXsc = 6 + 20/( (logP-0.182)*(logP-0.182) + 1.0 );
566  }
567  fTotalXsc = hpXsc;
568  }
569  if( neutron ) // pn to be np
570  {
571  if( pLab < 0.02 )
572  {
573  hnXsc = 4100+30*std::pow(std::log(1.3/pLab),3.6); // was as pLab < 0.8
574  fElasticXsc = hnXsc;
575  }
576  else if( pLab < 0.8 )
577  {
578  hnXsc = 33+30*std::pow(std::log(pLab/1.3),4.0);
579  fElasticXsc = hnXsc;
580  }
581  else if( pLab < 1.05 )
582  {
583  hnXsc = 33+30*std::pow(std::log(pLab/0.95),2.0);
584  fElasticXsc = 6 + 52/( std::log(0.511/pLab)*std::log(0.511/pLab) + 1.6 );
585  }
586  else if( pLab < 1.4 )
587  {
588  hnXsc = 33+30*std::pow(std::log(pLab/0.95),2.0);
589  fElasticXsc = 6 + 52/( std::log(0.511/pLab)*std::log(0.511/pLab) + 1.6 );
590  }
591  else // 1.4 < pLab < 10. )
592  {
593  hnXsc = 33.3 + 20.8*(std::pow(pLab,2.0) - 1.35)/(std::pow(pLab,2.50) + 0.95);
594 
595  fElasticXsc = 6 + 20/( (logP-0.182)*(logP-0.182) + 1.0 );
596  }
597  fTotalXsc = hnXsc;
598  }
599  }
600  }
601  else if( theParticle == theAProton && pORn ) /////////////////// p_bar ///////////////////////////
602  {
603  if( proton )
604  {
605  xsection = 35.45 + B*std::pow(std::log(sMand/s0),2.)
606  + 42.53*std::pow(sMand,-eta1) + 33.34*std::pow(sMand,-eta2);
607  }
608  if( neutron ) // ???
609  {
610  xsection = 35.80 + B*std::pow(std::log(sMand/s0),2.)
611  + 40.15*std::pow(sMand,-eta1) + 30.*std::pow(sMand,-eta2);
612  }
613  fTotalXsc = xsection;
614  }
615  else if( theParticle == thePiPlus && pORn ) // pi+ /////////////////////////////////////////////
616  {
617  if( proton ) // pi+ p
618  {
619  if( pLab < 0.28 )
620  {
621  hpXsc = 10./((logP + 1.273)*(logP + 1.273) + 0.05);
622  fElasticXsc = hpXsc;
623  }
624  else if( pLab < 0.4 )
625  {
626  hpXsc = 14./( (logP + 1.273)*(logP + 1.273) + 0.07);
627  fElasticXsc = hpXsc;
628  }
629  else if( pLab < 0.68 )
630  {
631  hpXsc = 14./( (logP + 1.273)*(logP + 1.273) + 0.07);
632  fElasticXsc = hpXsc;
633  }
634  else if( pLab < 0.85 )
635  {
636  G4double Ex4 = 88*(std::log(pLab/0.77))*(std::log(pLab/0.77));
637  hpXsc = Ex4 + 14.9;
638  fElasticXsc = hpXsc*std::exp(-3.*(pLab - 0.68));
639  }
640  else if( pLab < 1.15 )
641  {
642  G4double Ex4 = 88*(std::log(pLab/0.77))*(std::log(pLab/0.77));
643  hpXsc = Ex4 + 14.9;
644 
645  fElasticXsc = 6.0 + 1.4/(( pLab - 1.4)*( pLab - 1.4) + 0.1);
646  }
647  else if( pLab < 1.4) // ns original
648  {
649  G4double Ex1 = 3.2*std::exp(-(pLab-2.55)*(pLab-2.55)/0.55/0.55);
650  G4double Ex2 = 12*std::exp(-(pLab-1.47)*(pLab-1.47)/0.225/0.225);
651  hpXsc = Ex1 + Ex2 + 27.5;
652  fElasticXsc = 6.0 + 1.4/(( pLab - 1.4)*( pLab - 1.4) + 0.1);
653  }
654  else if( pLab < 2.0 ) // ns original
655  {
656  G4double Ex1 = 3.2*std::exp(-(pLab-2.55)*(pLab-2.55)/0.55/0.55);
657  G4double Ex2 = 12*std::exp(-(pLab-1.47)*(pLab-1.47)/0.225/0.225);
658  hpXsc = Ex1 + Ex2 + 27.5;
659  fElasticXsc = 3.0 + 1.36/( (logP - 0.336)*(logP - 0.336) + 0.08);
660  }
661  else if( pLab < 3.5 ) // ns original
662  {
663  G4double Ex1 = 3.2*std::exp(-(pLab-2.55)*(pLab-2.55)/0.55/0.55);
664  G4double Ex2 = 12*std::exp(-(pLab-1.47)*(pLab-1.47)/0.225/0.225);
665  hpXsc = Ex1 + Ex2 + 27.5;
666  fElasticXsc = 3.0 + 6.20/( (logP - 0.336)*(logP - 0.336) + 0.8);
667  }
668  else if( pLab < 200. ) // my
669  {
670  hpXsc = 10.6 + 2.*std::log(pE) + 25*std::pow(pE, -0.43 ); // ns original
671  // hpXsc = GetHadronNucleonXscPDG(aParticle, nucleon )/millibarn;
672  fElasticXsc = 3.0 + 6.20/( (logP - 0.336)*(logP - 0.336) + 0.8);
673  }
674  else // pLab > 100 // my
675  {
676  hpXsc = GetHadronNucleonXscPDG(aParticle, nucleon )/millibarn;
677  fElasticXsc = 3.0 + 6.20/( (logP - 0.336)*(logP - 0.336) + 0.8);
678  }
679  fTotalXsc = hpXsc;
680  }
681  if( neutron ) // pi+ n = pi- p??
682  {
683  if( pLab < 0.28 )
684  {
685  hnXsc = 0.288/((pLab - 0.28)*(pLab - 0.28) + 0.004);
686  fElasticXsc = 1.8/((logP + 1.273)*(logP + 1.273) + 0.07);
687  }
688  else if( pLab < 0.395676 ) // first peak
689  {
690  hnXsc = 0.648/((pLab - 0.28)*(pLab - 0.28) + 0.009);
691  fElasticXsc = 0.257/((pLab - 0.28)*(pLab - 0.28) + 0.01);
692  }
693  else if( pLab < 0.5 )
694  {
695  hnXsc = 26 + 110*(std::log(pLab/0.48))*(std::log(pLab/0.48));
696  fElasticXsc = 0.37*hnXsc;
697  }
698  else if( pLab < 0.65 )
699  {
700  hnXsc = 26 + 110*(std::log(pLab/0.48))*(std::log(pLab/0.48));
701  fElasticXsc = 0.95/((pLab - 0.72)*(pLab - 0.72) + 0.049);
702  }
703  else if( pLab < 0.72 )
704  {
705  hnXsc = 36.1 + 10*std::exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
706  24*std::exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
707  fElasticXsc = 0.95/((pLab - 0.72)*(pLab - 0.72) + 0.049);
708  }
709  else if( pLab < 0.88 )
710  {
711  hnXsc = 36.1 + 10*std::exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
712  24*std::exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
713  fElasticXsc = 0.95/((pLab - 0.72)*(pLab - 0.72) + 0.049);
714  }
715  else if( pLab < 1.03 )
716  {
717  hnXsc = 36.1 + 10*std::exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
718  24*std::exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
719  fElasticXsc = 2.0 + 0.4/((pLab - 1.03)*(pLab - 1.03) + 0.016);
720  }
721  else if( pLab < 1.15 )
722  {
723  hnXsc = 36.1 + 10*std::exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
724  24*std::exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
725  fElasticXsc = 2.0 + 0.4/((pLab - 1.03)*(pLab - 1.03) + 0.016);
726  }
727  else if( pLab < 1.3 )
728  {
729  hnXsc = 36.1 + 10*std::exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
730  24*std::exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
731  fElasticXsc = 3. + 13./pLab;
732  }
733  else if( pLab < 2.6 ) // < 3.0) // ns original
734  {
735  hnXsc = 36.1 + 0.079-4.313*std::log(pLab)+
736  3*std::exp(-(pLab-2.1)*(pLab-2.1)/0.4/0.4)+
737  1.5*std::exp(-(pLab-1.4)*(pLab-1.4)/0.12/0.12);
738  fElasticXsc = 3. + 13./pLab;
739  }
740  else if( pLab < 20. ) // < 3.0) // ns original
741  {
742  hnXsc = 36.1 + 0.079 - 4.313*std::log(pLab)+
743  3*std::exp(-(pLab-2.1)*(pLab-2.1)/0.4/0.4)+
744  1.5*std::exp(-(pLab-1.4)*(pLab-1.4)/0.12/0.12);
745  fElasticXsc = 3. + 13./pLab;
746  }
747  else // mb
748  {
749  hnXsc = GetHadronNucleonXscPDG(aParticle, nucleon )/millibarn;
750  fElasticXsc = 3. + 13./pLab;
751  }
752  fTotalXsc = hnXsc;
753  }
754  }
755  else if( theParticle == thePiMinus && pORn ) /// pi- ////////////////////////////////////////////
756  {
757  if( neutron ) // pi- n = pi+ p??
758  {
759  if( pLab < 0.28 )
760  {
761  hnXsc = 10./((logP + 1.273)*(logP + 1.273) + 0.05);
762  fElasticXsc = hnXsc;
763  }
764  else if( pLab < 0.4 )
765  {
766  hnXsc = 14./( (logP + 1.273)*(logP + 1.273) + 0.07);
767  fElasticXsc = hnXsc;
768  }
769  else if( pLab < 0.68 )
770  {
771  hnXsc = 14./( (logP + 1.273)*(logP + 1.273) + 0.07);
772  fElasticXsc = hnXsc;
773  }
774  else if( pLab < 0.85 )
775  {
776  G4double Ex4 = 88*(std::log(pLab/0.77))*(std::log(pLab/0.77));
777  hnXsc = Ex4 + 14.9;
778  fElasticXsc = hnXsc*std::exp(-3.*(pLab - 0.68));
779  }
780  else if( pLab < 1.15 )
781  {
782  G4double Ex4 = 88*(std::log(pLab/0.77))*(std::log(pLab/0.77));
783  hnXsc = Ex4 + 14.9;
784 
785  fElasticXsc = 6.0 + 1.4/(( pLab - 1.4)*( pLab - 1.4) + 0.1);
786  }
787  else if( pLab < 1.4) // ns original
788  {
789  G4double Ex1 = 3.2*std::exp(-(pLab-2.55)*(pLab-2.55)/0.55/0.55);
790  G4double Ex2 = 12*std::exp(-(pLab-1.47)*(pLab-1.47)/0.225/0.225);
791  hnXsc = Ex1 + Ex2 + 27.5;
792  fElasticXsc = 6.0 + 1.4/(( pLab - 1.4)*( pLab - 1.4) + 0.1);
793  }
794  else if( pLab < 2.0 ) // ns original
795  {
796  G4double Ex1 = 3.2*std::exp(-(pLab-2.55)*(pLab-2.55)/0.55/0.55);
797  G4double Ex2 = 12*std::exp(-(pLab-1.47)*(pLab-1.47)/0.225/0.225);
798  hnXsc = Ex1 + Ex2 + 27.5;
799  fElasticXsc = 3.0 + 1.36/( (logP - 0.336)*(logP - 0.336) + 0.08);
800  }
801  else if( pLab < 3.5 ) // ns original
802  {
803  G4double Ex1 = 3.2*std::exp(-(pLab-2.55)*(pLab-2.55)/0.55/0.55);
804  G4double Ex2 = 12*std::exp(-(pLab-1.47)*(pLab-1.47)/0.225/0.225);
805  hnXsc = Ex1 + Ex2 + 27.5;
806  fElasticXsc = 3.0 + 6.20/( (logP - 0.336)*(logP - 0.336) + 0.8);
807  }
808  else if( pLab < 200. ) // my
809  {
810  hnXsc = 10.6 + 2.*std::log(pE) + 25*std::pow(pE, -0.43 ); // ns original
811  fElasticXsc = 3.0 + 6.20/( (logP - 0.336)*(logP - 0.336) + 0.8);
812  }
813  else // pLab > 100 // my
814  {
815  hnXsc = GetHadronNucleonXscPDG(aParticle, nucleon )/millibarn;
816  fElasticXsc = 3.0 + 6.20/( (logP - 0.336)*(logP - 0.336) + 0.8);
817  }
818  fTotalXsc = hnXsc;
819  }
820  if( proton ) // pi- p
821  {
822  if( pLab < 0.28 )
823  {
824  hpXsc = 0.288/((pLab - 0.28)*(pLab - 0.28) + 0.004);
825  fElasticXsc = 1.8/((logP + 1.273)*(logP + 1.273) + 0.07);
826  }
827  else if( pLab < 0.395676 ) // first peak
828  {
829  hpXsc = 0.648/((pLab - 0.28)*(pLab - 0.28) + 0.009);
830  fElasticXsc = 0.257/((pLab - 0.28)*(pLab - 0.28) + 0.01);
831  }
832  else if( pLab < 0.5 )
833  {
834  hpXsc = 26 + 110*(std::log(pLab/0.48))*(std::log(pLab/0.48));
835  fElasticXsc = 0.37*hpXsc;
836  }
837  else if( pLab < 0.65 )
838  {
839  hpXsc = 26 + 110*(std::log(pLab/0.48))*(std::log(pLab/0.48));
840  fElasticXsc = 0.95/((pLab - 0.72)*(pLab - 0.72) + 0.049);
841  }
842  else if( pLab < 0.72 )
843  {
844  hpXsc = 36.1+
845  10*std::exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
846  24*std::exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
847  fElasticXsc = 0.95/((pLab - 0.72)*(pLab - 0.72) + 0.049);
848  }
849  else if( pLab < 0.88 )
850  {
851  hpXsc = 36.1+
852  10*std::exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
853  24*std::exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
854  fElasticXsc = 0.95/((pLab - 0.72)*(pLab - 0.72) + 0.049);
855  }
856  else if( pLab < 1.03 )
857  {
858  hpXsc = 36.1+
859  10*std::exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
860  24*std::exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
861  fElasticXsc = 2.0 + 0.4/((pLab - 1.03)*(pLab - 1.03) + 0.016);
862  }
863  else if( pLab < 1.15 )
864  {
865  hpXsc = 36.1+
866  10*std::exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
867  24*std::exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
868  fElasticXsc = 2.0 + 0.4/((pLab - 1.03)*(pLab - 1.03) + 0.016);
869  }
870  else if( pLab < 1.3 )
871  {
872  hpXsc = 36.1+
873  10*std::exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
874  24*std::exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
875  fElasticXsc = 3. + 13./pLab;
876  }
877  else if( pLab < 2.6 ) // < 3.0) // ns original
878  {
879  hpXsc = 36.1+0.079-4.313*std::log(pLab)+
880  3*std::exp(-(pLab-2.1)*(pLab-2.1)/0.4/0.4)+
881  1.5*std::exp(-(pLab-1.4)*(pLab-1.4)/0.12/0.12);
882  fElasticXsc = 3. +13./pLab; // *std::log(pLab*6.79);
883  }
884  else // mb
885  {
886  hpXsc = GetHadronNucleonXscPDG(aParticle, nucleon )/millibarn;
887  fElasticXsc = 3. + 13./pLab;
888  }
889  fTotalXsc = hpXsc;
890  }
891  }
892  else if( (theParticle == theKMinus || theParticle == theK0S) && proton ) // Kmp/K0p /////////////////////////////////
893  {
894  if( pLab < pMin)
895  {
896  G4double psp = pLab*std::sqrt(pLab);
897  fElasticXsc = 5.2/psp;
898  fTotalXsc = 14./psp;
899  }
900  else if( pLab > pMax )
901  {
902  G4double ld = std::log(pLab) - minLogP;
903  G4double ld2 = ld*ld;
904  fElasticXsc = cofLogE*ld2 + 2.23;
905  fTotalXsc = 1.1*cofLogT*ld2 + 19.7;
906  }
907  else
908  {
909  G4double ld = std::log(pLab) - minLogP;
910  G4double ld2 = ld*ld;
911  G4double sp = std::sqrt(pLab);
912  G4double psp = pLab*sp;
913  G4double p2 = pLab*pLab;
914  G4double p4 = p2*p2;
915  G4double lm = pLab - .39;
916  G4double md = lm*lm + .000356;
917 
918  G4double lh1 = pLab - 0.78;
919  G4double hd1 = lh1*lh1 + .00166;
920 
921  G4double lh = pLab - 1.01;
922  G4double hd = lh*lh + .011;
923 
924  G4double lh2 = pLab - 1.63;
925  G4double hd2 = lh2*lh2 + .007;
926 
927  fElasticXsc = 5.2/psp + (1.1*cofLogE*ld2 + 2.23)/(1. - .7/sp + .075/p4)
928  + .004/md + 0.005/hd1+ 0.01/hd2 +.15/hd; // small peaks were added
929 
930  fTotalXsc = 14./psp + (1.1*cofLogT*ld2 + 19.5)/(1. - .21/sp + .52/p4)
931  + .006/md + 0.01/hd1+ 0.02/hd2 + .20/hd ;
932  }
933  }
934  else if( (theParticle == theKMinus || theParticle == theK0S) && neutron ) // Kmn/K0n //////////////////////////////
935  {
936  if( pLab > pMax )
937  {
938  G4double ld = std::log(pLab) - minLogP;
939  G4double ld2 = ld*ld;
940  fElasticXsc = cofLogE*ld2 + 2.23;
941  fTotalXsc = 1.1*cofLogT*ld2 + 19.7;
942  }
943  else
944  {
945 
946  G4double lh = pLab - 0.98;
947  G4double hd = lh*lh + .021;
948 
949  G4double LogPlab = std::log( pLab );
950  G4double sqrLogPlab = LogPlab * LogPlab;
951 
952  fElasticXsc = // 5.2/psp + (cofLogE*ld2 + 2.23)/(1. - .7/sp + .075/p4) + .004/md
953  5.0 + 8.1*std::pow(pLab,-1.8 ) + 0.16*sqrLogPlab - 1.3*LogPlab + .15/hd;
954  fTotalXsc = // 14./psp +
955  // (1.1*cofLogT*ld2 + 19.5)/(1. - .21/sp + .52/p4)
956  25.2 + 0. *std::pow(pLab, 0. ) + 0.38*sqrLogPlab - 2.9*LogPlab
957  // + .006/md + 0.01/hd1+ 0.02/hd2
958  + 0.30/hd ;
959  }
960  }
961  else if( (theParticle == theKPlus || theParticle == theK0L) && proton ) // Kpp/aKp ////////////////////////
962  {
963  if( pLab < pMin )
964  {
965  G4double lr = pLab - .38;
966  G4double lm = pLab - 1.;
967  G4double md = lm*lm + .392;
968  fElasticXsc = .7/(lr*lr + .076) + 2./md;
969  fTotalXsc = .7/(lr*lr + .076) + 2.6/md;
970  }
971  else if( pLab > pMax )
972  {
973  G4double ld = std::log(pLab) - minLogP;
974  G4double ld2 = ld*ld;
975  fElasticXsc = cofLogE*ld2 + 2.23;
976  fTotalXsc = cofLogT*ld2 + 19.2;
977  }
978  else
979  {
980  G4double ld = std::log(pLab) - minLogP;
981  G4double ld2 = ld*ld;
982  G4double lr = pLab - .38;
983  G4double LE = .7/(lr*lr + .076);
984  G4double sp = std::sqrt(pLab);
985  G4double p2 = pLab*pLab;
986  G4double p4 = p2*p2;
987  G4double lm = pLab - 1.;
988  G4double md = lm*lm + .392;
989  fElasticXsc = LE + (cofLogE*ld2 + 2.23)/(1. - .7/sp + .1/p4) + 2./md;
990  fTotalXsc = LE + (cofLogT*ld2 + 19.5)/(1. + .46/sp + 1.6/p4) + 2.6/md;
991  }
992  }
993  else if( (theParticle == theKPlus || theParticle == theK0L) && neutron ) // Kpn/aKn ///////////////////////
994  {
995  if( pLab < pMin )
996  {
997  G4double lm = pLab - 0.94;
998  G4double md = lm*lm + .392;
999  fElasticXsc = 2./md;
1000  fTotalXsc = 4.6/md;
1001  }
1002  else if( pLab > pMax )
1003  {
1004  G4double ld = std::log(pLab) - minLogP;
1005  G4double ld2 = ld*ld;
1006  fElasticXsc = cofLogE*ld2 + 2.23;
1007  fTotalXsc = cofLogT*ld2 + 19.2;
1008  }
1009  else
1010  {
1011  G4double ld = std::log(pLab) - minLogP;
1012  G4double ld2 = ld*ld;
1013  G4double sp = std::sqrt(pLab);
1014  G4double p2 = pLab*pLab;
1015  G4double p4 = p2*p2;
1016  G4double lm = pLab - 0.94;
1017  G4double md = lm*lm + .392;
1018  fElasticXsc = (cofLogE*ld2 + 2.23)/(1. - .7/sp + .1/p4) + 2./md;
1019  fTotalXsc = (cofLogT*ld2 + 19.5)/(1. + .46/sp + 1.6/p4) + 4.6/md;
1020  }
1021  }
1022  else if( theParticle == theSMinus && pORn )
1023  {
1024  xsection = 35.20 + B*std::pow(std::log(sMand/s0),2.)
1025  - 199.*std::pow(sMand,-eta1) + 264.*std::pow(sMand,-eta2);
1026  }
1027  else if( theParticle == theGamma && pORn ) // modify later on
1028  {
1029  xsection = 0.0 + B*std::pow(std::log(sMand/s0),2.)
1030  + 0.032*std::pow(sMand,-eta1) - 0.0*std::pow(sMand,-eta2);
1031  fTotalXsc = xsection;
1032  }
1033  else // other then p,n,pi+,pi-,K+,K- as proton ???
1034  {
1035  if( proton )
1036  {
1037  xsection = 35.45 + B*std::pow(std::log(sMand/s0),2.)
1038  + 42.53*std::pow(sMand,-eta1) - 33.34*std::pow(sMand,-eta2);
1039  }
1040  if( neutron )
1041  {
1042  xsection += 35.80 + B*std::pow(std::log(sMand/s0),2.)
1043  + 40.15*std::pow(sMand,-eta1) - 30.*std::pow(sMand,-eta2);
1044  }
1045  fTotalXsc = xsection;
1046  }
1047  fTotalXsc *= millibarn; // parametrised in mb
1048  fElasticXsc *= millibarn; // parametrised in mb
1049 
1050  if( proton && aParticle->GetDefinition()->GetPDGCharge() > 0. )
1051  {
1052  G4double cB = GetCoulombBarrier(aParticle, nucleon);
1053  fTotalXsc *= cB;
1054  fElasticXsc *= cB;
1055  }
1056  fInelasticXsc = fTotalXsc - fElasticXsc;
1057  if( fInelasticXsc < 0. ) fInelasticXsc = 0.;
1058 
1059  // G4cout<<fTotalXsc/millibarn<<"; "<<fElasticXsc/millibarn<<"; "<<fInelasticXsc/millibarn<<G4endl;
1060 
1061  return fTotalXsc;
1062 }
1063 
1064 /////////////////////////////////////////////////////////////////////////////////////
1065 //
1066 // Returns kaon-nucleon cross-section based on smoothed NS for GG model
1067 
1068 G4double
1069 G4HadronNucleonXsc::GetKaonNucleonXscGG(const G4DynamicParticle* aParticle,
1070  const G4ParticleDefinition* nucleon )
1071 {
1072  G4double pLab = aParticle->GetMomentum().mag();
1073 
1074  pLab /= GeV;
1075  G4double LogPlab = std::log( pLab );
1076  G4double sqrLogPlab = LogPlab * LogPlab;
1077 
1078  G4double minLogP = 3.5; // min of (lnP-minLogP)^2
1079  G4double cofLogE = .0557; // elastic (lnP-minLogP)^2
1080  G4double cofLogT = .3; // total (lnP-minLogP)^2
1081  G4double pMin = .1; // fast LE calculation
1082  G4double pMax = 1000.; // fast HE calculation
1083 
1084  const G4ParticleDefinition* theParticle = aParticle->GetDefinition();
1085 
1086  G4bool proton = (nucleon == theProton);
1087  G4bool neutron = (nucleon == theNeutron);
1088 
1089  if( (theParticle == theKMinus || theParticle == theK0S) && proton ) // (K-,K0)on p ////////////////////////////
1090  {
1091 
1092  if( pLab < pMin)
1093  {
1094  G4double psp = pLab*std::sqrt(pLab);
1095  fElasticXsc = 5.2/psp;
1096  fTotalXsc = 14./psp;
1097  }
1098  else if( pLab > pMax )
1099  {
1100  G4double ld = std::log(pLab) - minLogP;
1101  G4double ld2 = ld*ld;
1102  fElasticXsc = cofLogE*ld2 + 2.23;
1103  fTotalXsc = 1.1*cofLogT*ld2 + 19.7;
1104  }
1105  else
1106  {
1107  G4double ld = std::log(pLab) - minLogP;
1108  G4double ld2 = ld*ld;
1109  G4double sp = std::sqrt(pLab);
1110  G4double psp = pLab*sp;
1111  G4double p2 = pLab*pLab;
1112  G4double p4 = p2*p2;
1113 
1114  G4double lh = pLab - 0.98;
1115  G4double hd = lh*lh + .045;
1116 
1117 
1118  fElasticXsc = 5.2/psp + (cofLogE*ld2 + 2.23)/(1. - .7/sp + .075/p4) // + .004/md
1119  + .15/hd;
1120  fTotalXsc = 14./psp + (1.1*cofLogT*ld2 + 19.5)/(1. - .21/sp + .52/p4)
1121  // + .006/md + 0.01/hd1 + 0.02/hd2
1122  + .60/hd;
1123  }
1124  }
1125  else if( (theParticle == theKMinus || theParticle == theK0S) && neutron ) // Kmn/K0n /////////////////////////////
1126  {
1127  if( pLab > pMax )
1128  {
1129  G4double ld = std::log(pLab) - minLogP;
1130  G4double ld2 = ld*ld;
1131  fElasticXsc = cofLogE*ld2 + 2.23;
1132  fTotalXsc = 1.1*cofLogT*ld2 + 19.7;
1133  }
1134  else
1135  {
1136 
1137  G4double lh = pLab - 0.98;
1138  G4double hd = lh*lh + .045;
1139 
1140  fElasticXsc = // 5.2/psp + (cofLogE*ld2 + 2.23)/(1. - .7/sp + .075/p4) + .004/md
1141  5.0 + 8.1*std::pow(pLab,-1.8 ) + 0.16*sqrLogPlab - 1.3*LogPlab + .15/hd;
1142  fTotalXsc = // 14./psp +
1143  // (1.1*cofLogT*ld2 + 19.5)/(1. - .21/sp + .52/p4)
1144  25.2 + 0. *std::pow(pLab, 0. ) + 0.38*sqrLogPlab - 2.9*LogPlab
1145  // + .006/md + 0.01/hd1+ 0.02/hd2
1146  + 0.60/hd ;
1147  }
1148  }
1149  else if( (theParticle == theKPlus || theParticle == theK0L) && proton ) // Kpp/aKp //////////////////////
1150  {
1151  if( pLab < pMin )
1152  {
1153  G4double lr = pLab - .38;
1154  G4double lm = pLab - 1.;
1155  G4double md = lm*lm + .392;
1156  fElasticXsc = .7/(lr*lr + .076) + 2./md;
1157  fTotalXsc = // .7/(lr*lr + .076) +
1158  2.6/md;
1159  }
1160  else if( pLab > pMax )
1161  {
1162  G4double ld = std::log(pLab) - minLogP;
1163  G4double ld2 = ld*ld;
1164  fElasticXsc = cofLogE*ld2 + 2.23;
1165  fTotalXsc = cofLogT*ld2 + 19.2;
1166  }
1167  else
1168  {
1169  G4double ld = std::log(pLab) - minLogP;
1170  G4double ld2 = ld*ld;
1171  G4double lr = pLab - .38;
1172  G4double LE = .7/(lr*lr + .076);
1173  G4double sp = std::sqrt(pLab);
1174  G4double p2 = pLab*pLab;
1175  G4double p4 = p2*p2;
1176  G4double lm = pLab - 0.8;
1177  G4double md = lm*lm + .652;
1178  fElasticXsc = LE + (cofLogE*ld2 + 2.23)/(1. - .7/sp + .1/p4) + 2./md;
1179  fTotalXsc = (cofLogT*ld2 + 19.5)/(1. + .46/sp + 1.6/p4) + 7.6/md; // + LE;
1180  }
1181  }
1182  else if( (theParticle == theKPlus || theParticle == theK0L) && neutron ) // Kpn/aKn //////////////////////////////////
1183  {
1184  if( pLab < pMin )
1185  {
1186  G4double lm = pLab - 0.94;
1187  G4double md = lm*lm + .392;
1188  fElasticXsc = 2./md;
1189  fTotalXsc = 4.6/md;
1190  }
1191  else if( pLab > pMax )
1192  {
1193  G4double ld = std::log(pLab) - minLogP;
1194  G4double ld2 = ld*ld;
1195  fElasticXsc = cofLogE*ld2 + 2.23;
1196  fTotalXsc = cofLogT*ld2 + 19.2;
1197  }
1198  else
1199  {
1200  G4double ld = std::log(pLab) - minLogP;
1201  G4double ld2 = ld*ld;
1202  G4double sp = std::sqrt(pLab);
1203  G4double p2 = pLab*pLab;
1204  G4double p4 = p2*p2;
1205  G4double lm = pLab - 0.8;
1206  G4double md = lm*lm + .652;
1207  fElasticXsc = (cofLogE*ld2 + 2.23)/(1. - .7/sp + .1/p4) + 2./md;
1208  fTotalXsc = (cofLogT*ld2 + 19.5)/(1. + .46/sp + 1.6/p4) + 7.6/md;
1209  }
1210  }
1211  fTotalXsc *= millibarn; // parametrised in mb
1212  fElasticXsc *= millibarn; // parametrised in mb
1213 
1214  if( proton && aParticle->GetDefinition()->GetPDGCharge() > 0. )
1215  {
1216  G4double cB = GetCoulombBarrier(aParticle, nucleon);
1217  fTotalXsc *= cB;
1218  fElasticXsc *= cB;
1219  }
1220  fInelasticXsc = fTotalXsc - fElasticXsc;
1221  if( fInelasticXsc < 0. ) fInelasticXsc = 0.;
1222 
1223  // G4cout<<fTotalXsc/millibarn<<"; "<<fElasticXsc/millibarn<<"; "<<fInelasticXsc/millibarn<<G4endl;
1224 
1225  return fTotalXsc;
1226 }
1227 
1228 /////////////////////////////////////////////////////////////////////////////////////
1229 //
1230 // Returns hadron-nucleon cross-section based on V. Uzjinsky parametrisation of
1231 // data from G4FTFCrossSection class
1232 
1233 G4double
1234 G4HadronNucleonXsc::GetHadronNucleonXscVU(const G4DynamicParticle* aParticle,
1235  const G4ParticleDefinition* nucleon )
1236 {
1237  G4int PDGcode = aParticle->GetDefinition()->GetPDGEncoding();
1238  G4int absPDGcode = std::abs(PDGcode);
1239  G4double Elab = aParticle->GetTotalEnergy();
1240  // (s - 2*0.88*GeV*GeV)/(2*0.939*GeV)/GeV;
1241  G4double Plab = aParticle->GetMomentum().mag();
1242  // std::sqrt(Elab * Elab - 0.88);
1243 
1244  Elab /= GeV;
1245  Plab /= GeV;
1246 
1247  G4double LogPlab = std::log( Plab );
1248  G4double sqrLogPlab = LogPlab * LogPlab;
1249 
1250  G4bool pORn = (nucleon == theProton || nucleon == theNeutron );
1251  G4bool proton = (nucleon == theProton);
1252  G4bool neutron = (nucleon == theNeutron);
1253 
1254 
1255  if( absPDGcode > 1000 && pORn ) //------Projectile is baryon -
1256  {
1257  if(proton)
1258  {
1259  fTotalXsc = 48.0 + 0. *std::pow(Plab, 0. ) + 0.522*sqrLogPlab - 4.51*LogPlab;
1260  fElasticXsc = 11.9 + 26.9*std::pow(Plab,-1.21) + 0.169*sqrLogPlab - 1.85*LogPlab;
1261  }
1262  if(neutron)
1263  {
1264  fTotalXsc = 47.3 + 0. *std::pow(Plab, 0. ) + 0.513*sqrLogPlab - 4.27*LogPlab;
1265  fElasticXsc = 11.9 + 26.9*std::pow(Plab,-1.21) + 0.169*sqrLogPlab - 1.85*LogPlab;
1266  }
1267  }
1268  else if( PDGcode == 211 && pORn ) //------Projectile is PionPlus ----
1269  {
1270  if(proton)
1271  {
1272  fTotalXsc = 16.4 + 19.3 *std::pow(Plab,-0.42) + 0.19 *sqrLogPlab - 0.0 *LogPlab;
1273  fElasticXsc = 0.0 + 11.4*std::pow(Plab,-0.40) + 0.079*sqrLogPlab - 0.0 *LogPlab;
1274  }
1275  if(neutron)
1276  {
1277  fTotalXsc = 33.0 + 14.0 *std::pow(Plab,-1.36) + 0.456*sqrLogPlab - 4.03*LogPlab;
1278  fElasticXsc = 1.76 + 11.2*std::pow(Plab,-0.64) + 0.043*sqrLogPlab - 0.0 *LogPlab;
1279  }
1280  }
1281  else if( PDGcode == -211 && pORn ) //------Projectile is PionMinus ----
1282  {
1283  if(proton)
1284  {
1285  fTotalXsc = 33.0 + 14.0 *std::pow(Plab,-1.36) + 0.456*sqrLogPlab - 4.03*LogPlab;
1286  fElasticXsc = 1.76 + 11.2*std::pow(Plab,-0.64) + 0.043*sqrLogPlab - 0.0 *LogPlab;
1287  }
1288  if(neutron)
1289  {
1290  fTotalXsc = 16.4 + 19.3 *std::pow(Plab,-0.42) + 0.19 *sqrLogPlab - 0.0 *LogPlab;
1291  fElasticXsc = 0.0 + 11.4*std::pow(Plab,-0.40) + 0.079*sqrLogPlab - 0.0 *LogPlab;
1292  }
1293  }
1294  else if( PDGcode == 111 && pORn ) //------Projectile is PionZero --
1295  {
1296  if(proton)
1297  {
1298  fTotalXsc = (16.4 + 19.3 *std::pow(Plab,-0.42) + 0.19 *sqrLogPlab - 0.0 *LogPlab + //Pi+
1299  33.0 + 14.0 *std::pow(Plab,-1.36) + 0.456*sqrLogPlab - 4.03*LogPlab)/2; //Pi-
1300 
1301  fElasticXsc = ( 0.0 + 11.4*std::pow(Plab,-0.40) + 0.079*sqrLogPlab - 0.0 *LogPlab + //Pi+
1302  1.76 + 11.2*std::pow(Plab,-0.64) + 0.043*sqrLogPlab - 0.0 *LogPlab)/2; //Pi-
1303 
1304  }
1305  if(neutron)
1306  {
1307  fTotalXsc = (33.0 + 14.0 *std::pow(Plab,-1.36) + 0.456*sqrLogPlab - 4.03*LogPlab + //Pi+
1308  16.4 + 19.3 *std::pow(Plab,-0.42) + 0.19 *sqrLogPlab - 0.0 *LogPlab)/2; //Pi-
1309  fElasticXsc = ( 1.76 + 11.2*std::pow(Plab,-0.64) + 0.043*sqrLogPlab - 0.0 *LogPlab + //Pi+
1310  0.0 + 11.4*std::pow(Plab,-0.40) + 0.079*sqrLogPlab - 0.0 *LogPlab)/2; //Pi-
1311  }
1312  }
1313  else if( PDGcode == 321 && pORn ) //------Projectile is KaonPlus --
1314  {
1315  if(proton)
1316  {
1317  fTotalXsc = 18.1 + 0. *std::pow(Plab, 0. ) + 0.26 *sqrLogPlab - 1.0 *LogPlab;
1318  fElasticXsc = 5.0 + 8.1*std::pow(Plab,-1.8 ) + 0.16 *sqrLogPlab - 1.3 *LogPlab;
1319  }
1320  if(neutron)
1321  {
1322  fTotalXsc = 18.7 + 0. *std::pow(Plab, 0. ) + 0.21 *sqrLogPlab - 0.89*LogPlab;
1323  fElasticXsc = 7.3 + 0. *std::pow(Plab,-0. ) + 0.29 *sqrLogPlab - 2.4 *LogPlab;
1324  }
1325  }
1326  else if( PDGcode ==-321 && pORn ) //------Projectile is KaonMinus ----
1327  {
1328  if(proton)
1329  {
1330  fTotalXsc = 32.1 + 0. *std::pow(Plab, 0. ) + 0.66*sqrLogPlab - 5.6*LogPlab;
1331  fElasticXsc = 7.3 + 0. *std::pow(Plab,-0. ) + 0.29*sqrLogPlab - 2.4*LogPlab;
1332  }
1333  if(neutron)
1334  {
1335  fTotalXsc = 25.2 + 0. *std::pow(Plab, 0. ) + 0.38*sqrLogPlab - 2.9*LogPlab;
1336  fElasticXsc = 5.0 + 8.1*std::pow(Plab,-1.8 ) + 0.16*sqrLogPlab - 1.3*LogPlab;
1337  }
1338  }
1339  else if( PDGcode == 311 && pORn ) //------Projectile is KaonZero -----
1340  {
1341  if(proton)
1342  {
1343  fTotalXsc = ( 18.1 + 0. *std::pow(Plab, 0. ) + 0.26 *sqrLogPlab - 1.0 *LogPlab + //K+
1344  32.1 + 0. *std::pow(Plab, 0. ) + 0.66 *sqrLogPlab - 5.6 *LogPlab)/2; //K-
1345  fElasticXsc = ( 5.0 + 8.1*std::pow(Plab,-1.8 ) + 0.16 *sqrLogPlab - 1.3 *LogPlab + //K+
1346  7.3 + 0. *std::pow(Plab,-0. ) + 0.29 *sqrLogPlab - 2.4 *LogPlab)/2; //K-
1347  }
1348  if(neutron)
1349  {
1350  fTotalXsc = ( 18.7 + 0. *std::pow(Plab, 0. ) + 0.21 *sqrLogPlab - 0.89*LogPlab + //K+
1351  25.2 + 0. *std::pow(Plab, 0. ) + 0.38 *sqrLogPlab - 2.9 *LogPlab)/2; //K-
1352  fElasticXsc = ( 7.3 + 0. *std::pow(Plab,-0. ) + 0.29 *sqrLogPlab - 2.4 *LogPlab + //K+
1353  5.0 + 8.1*std::pow(Plab,-1.8 ) + 0.16 *sqrLogPlab - 1.3 *LogPlab)/2; //K-
1354  }
1355  }
1356  else //------Projectile is undefined, Nucleon assumed
1357  {
1358  if(proton)
1359  {
1360  fTotalXsc = 48.0 + 0. *std::pow(Plab, 0. ) + 0.522*sqrLogPlab - 4.51*LogPlab;
1361  fElasticXsc = 11.9 + 26.9*std::pow(Plab,-1.21) + 0.169*sqrLogPlab - 1.85*LogPlab;
1362  }
1363  if(neutron)
1364  {
1365  fTotalXsc = 47.3 + 0. *std::pow(Plab, 0. ) + 0.513*sqrLogPlab - 4.27*LogPlab;
1366  fElasticXsc = 11.9 + 26.9*std::pow(Plab,-1.21) + 0.169*sqrLogPlab - 1.85*LogPlab;
1367  }
1368  }
1369  fTotalXsc *= millibarn;
1370  fElasticXsc *= millibarn;
1371  fInelasticXsc = fTotalXsc - fElasticXsc;
1372  if (fInelasticXsc < 0.) fInelasticXsc = 0.;
1373 
1374  return fTotalXsc;
1375 }
1376 
1377 ////////////////////////////////////////////////////////////////////////////////////
1378 //
1379 //
1380 
1381 G4double G4HadronNucleonXsc::CalculateEcmValue( const G4double mp ,
1382  const G4double mt ,
1383  const G4double Plab )
1384 {
1385  G4double Elab = std::sqrt ( mp * mp + Plab * Plab );
1386  G4double Ecm = std::sqrt ( mp * mp + mt * mt + 2 * Elab * mt );
1387  // G4double Pcm = Plab * mt / Ecm;
1388  // G4double KEcm = std::sqrt ( Pcm * Pcm + mp * mp ) - mp;
1389 
1390  return Ecm ; // KEcm;
1391 }
1392 
1393 
1394 ////////////////////////////////////////////////////////////////////////////////////
1395 //
1396 //
1397 
1398 G4double G4HadronNucleonXsc::CalcMandelstamS( const G4double mp ,
1399  const G4double mt ,
1400  const G4double Plab )
1401 {
1402  G4double Elab = std::sqrt ( mp * mp + Plab * Plab );
1403  G4double sMand = mp*mp + mt*mt + 2*Elab*mt ;
1404 
1405  return sMand;
1406 }
1407 
1408 
1409 ///////////////////////////////////////////////////////////////////////////////
1410 //
1411 //
1412 
1413 G4double G4HadronNucleonXsc::GetCoulombBarrier(const G4DynamicParticle* aParticle,
1414  const G4ParticleDefinition* nucleon )
1415 {
1416  G4double ratio;
1417 
1418  G4double tR = 0.895*fermi, pR;
1419 
1420  if ( aParticle->GetDefinition() == theProton ) pR = 0.895*fermi;
1421  else if( aParticle->GetDefinition() == thePiPlus ) pR = 0.663*fermi;
1422  else if( aParticle->GetDefinition() == theKPlus ) pR = 0.340*fermi;
1423  else pR = 0.500*fermi;
1424 
1425  G4double pZ = aParticle->GetDefinition()->GetPDGCharge();
1426  G4double tZ = nucleon->GetPDGCharge();
1427 
1428  G4double pTkin = aParticle->GetKineticEnergy();
1429 
1430  G4double pM = aParticle->GetDefinition()->GetPDGMass();
1431  G4double tM = nucleon->GetPDGMass();
1432 
1433  G4double pElab = pTkin + pM;
1434 
1435  G4double totEcm = std::sqrt(pM*pM + tM*tM + 2.*pElab*tM);
1436 
1437  G4double totTcm = totEcm - pM -tM;
1438 
1439  G4double bC = fine_structure_const*hbarc*pZ*tZ;
1440  bC /= pR + tR;
1441  bC /= 2.; // 4., 2. parametrisation cof ??? vmg
1442 
1443  // G4cout<<"pTkin = "<<pTkin/GeV<<"; pPlab = "
1444  // <<pPlab/GeV<<"; bC = "<<bC/GeV<<"; pTcm = "<<pTcm/GeV<<G4endl;
1445 
1446  if( totTcm <= bC ) ratio = 0.;
1447  else ratio = 1. - bC/totTcm;
1448 
1449  // if(ratio < DBL_MIN) ratio = DBL_MIN;
1450  if( ratio < 0.) ratio = 0.;
1451 
1452  // G4cout <<"ratio = "<<ratio<<G4endl;
1453  return ratio;
1454 }
1455 
1456 
1457 
1458 
1459 
1460 ////////////////////////////////////////////////////////////////////////////////////
1461 //
1462 // Initialaise K(p,m)-(p,n) total cross section vectors
1463 
1464 
1465 void G4HadronNucleonXsc::InitialiseKaonNucleonTotXsc()
1466 {
1467  G4int i = 0, iMax;
1468  G4double tmpxsc[106];
1469 
1470  // Kp-proton tot xsc
1471 
1472  iMax = 66;
1473  fKpProtonTotXscVector = G4LPhysicsFreeVector(iMax, fKpProtonTotTkin[0], fKpProtonTotTkin[iMax-1]);
1474  fKpProtonTotXscVector.SetSpline(true);
1475 
1476  for( i = 0; i < iMax; i++)
1477  {
1478  tmpxsc[i] = 0.;
1479 
1480  if( i == 0 || i == iMax-1 ) tmpxsc[i] = fKpProtonTotXsc[i];
1481  else tmpxsc[i] = 0.5*(fKpProtonTotXsc[i-1]+fKpProtonTotXsc[i+1]);
1482 
1483  fKpProtonTotXscVector.PutValues(size_t(i), fKpProtonTotTkin[i], tmpxsc[i]*millibarn);
1484  }
1485 
1486  // Kp-neutron tot xsc
1487 
1488  iMax = 75;
1489  fKpNeutronTotXscVector = G4LPhysicsFreeVector(iMax, fKpNeutronTotTkin[0], fKpNeutronTotTkin[iMax-1]);
1490  fKpNeutronTotXscVector.SetSpline(true);
1491 
1492  for( i = 0; i < iMax; i++)
1493  {
1494  tmpxsc[i] = 0.;
1495  if( i == 0 || i == iMax-1 ) tmpxsc[i] = fKpNeutronTotXsc[i];
1496  else tmpxsc[i] = 0.5*(fKpNeutronTotXsc[i-1]+fKpNeutronTotXsc[i+1]);
1497 
1498  fKpNeutronTotXscVector.PutValues(size_t(i), fKpNeutronTotTkin[i], tmpxsc[i]*millibarn);
1499  }
1500 
1501  // Km-proton tot xsc
1502 
1503  iMax = 106;
1504  fKmProtonTotXscVector = G4LPhysicsFreeVector(iMax, fKmProtonTotTkin[0], fKmProtonTotTkin[iMax-1]);
1505  fKmProtonTotXscVector.SetSpline(true);
1506 
1507  for( i = 0; i < iMax; i++)
1508  {
1509  tmpxsc[i] = 0.;
1510 
1511  if( i == 0 || i == iMax-1 ) tmpxsc[i] = fKmProtonTotXsc[i];
1512  else tmpxsc[i] = 0.5*(fKmProtonTotXsc[i-1]+fKmProtonTotXsc[i+1]);
1513 
1514  fKmProtonTotXscVector.PutValues(size_t(i), fKmProtonTotTkin[i], tmpxsc[i]*millibarn);
1515  }
1516 
1517  // Km-neutron tot xsc
1518 
1519  iMax = 68;
1520  fKmNeutronTotXscVector = G4LPhysicsFreeVector(iMax, fKmNeutronTotTkin[0], fKmNeutronTotTkin[iMax-1]);
1521  fKmNeutronTotXscVector.SetSpline(true);
1522 
1523  for( i = 0; i < iMax; i++)
1524  {
1525  tmpxsc[i] = 0.;
1526  if( i == 0 || i == iMax-1 ) tmpxsc[i] = fKmNeutronTotXsc[i];
1527  else tmpxsc[i] = 0.5*(fKmNeutronTotXsc[i-1]+fKmNeutronTotXsc[i+1]);
1528 
1529  fKmNeutronTotXscVector.PutValues(size_t(i), fKmNeutronTotTkin[i], tmpxsc[i]*millibarn);
1530  }
1531 }
1532 
1533 ///////////////////////////////////////////////////////
1534 //
1535 // K-nucleon tot xsc (mb) fit data, std::log(Tkin(MeV))
1536 
1537 const G4double G4HadronNucleonXsc::fKpProtonTotXsc[66] = {
1538 0.000000e+00, 1.592400e-01, 3.184700e-01, 7.961800e-01, 1.433120e+00, 2.070060e+00,
1539 2.866240e+00, 3.582800e+00, 4.378980e+00, 5.015920e+00, 5.573250e+00, 6.449040e+00,
1540 7.404460e+00, 8.200640e+00, 8.837580e+00, 9.713380e+00, 1.027070e+01, 1.090764e+01,
1541 1.130573e+01, 1.170382e+01, 1.242038e+01, 1.281847e+01, 1.321656e+01, 1.337580e+01,
1542 1.345541e+01, 1.329618e+01, 1.265924e+01, 1.242038e+01, 1.250000e+01, 1.305732e+01,
1543 1.369427e+01, 1.425159e+01, 1.544586e+01, 1.648089e+01, 1.751592e+01, 1.791401e+01,
1544 1.791401e+01, 1.775478e+01, 1.751592e+01, 1.735669e+01, 1.719745e+01, 1.711783e+01,
1545 1.703822e+01, 1.695860e+01, 1.695860e+01, 1.695860e+01, 1.695860e+01, 1.687898e+01,
1546 1.687898e+01, 1.703822e+01, 1.719745e+01, 1.735669e+01, 1.751592e+01, 1.767516e+01,
1547 1.783439e+01, 1.799363e+01, 1.815287e+01, 1.839172e+01, 1.855095e+01, 1.871019e+01,
1548 1.886943e+01, 1.918790e+01, 1.942675e+01, 1.966561e+01, 2.006369e+01, 2.054140e+01
1549 }; // 66
1550 
1551 
1552 const G4double G4HadronNucleonXsc::fKpProtonTotTkin[66] = {
1553 2.100000e+00, 2.180770e+00, 2.261540e+00, 2.396150e+00, 2.476920e+00, 2.557690e+00,
1554 2.557690e+00, 2.584620e+00, 2.638460e+00, 2.665380e+00, 2.719230e+00, 2.746150e+00,
1555 2.800000e+00, 2.853850e+00, 2.934620e+00, 3.042310e+00, 3.150000e+00, 3.311540e+00,
1556 3.446150e+00, 3.607690e+00, 3.930770e+00, 4.226920e+00, 4.361540e+00, 4.846150e+00,
1557 4.980770e+00, 5.088460e+00, 5.465380e+00, 5.653850e+00, 5.950000e+00, 6.084620e+00,
1558 6.246150e+00, 6.300000e+00, 6.380770e+00, 6.515380e+00, 6.730770e+00, 6.838460e+00,
1559 7.000000e+00, 7.161540e+00, 7.323080e+00, 7.457690e+00, 7.619230e+00, 7.780770e+00,
1560 7.915380e+00, 8.130770e+00, 8.265380e+00, 8.453850e+00, 8.642310e+00, 8.803850e+00,
1561 9.019230e+00, 9.234620e+00, 9.530770e+00, 9.773080e+00, 1.001538e+01, 1.017692e+01,
1562 1.033846e+01, 1.058077e+01, 1.082308e+01, 1.098462e+01, 1.114615e+01, 1.138846e+01,
1563 1.160385e+01, 1.173846e+01, 1.192692e+01, 1.216923e+01, 1.238461e+01, 1.257308e+01
1564 }; // 66
1565 
1566 const G4double G4HadronNucleonXsc::fKpNeutronTotXsc[75] = {
1567 3.980900e-01, 3.184700e-01, 3.184700e-01, 3.980900e-01, 3.980900e-01, 3.980900e-01,
1568 3.980900e-01, 3.980900e-01, 3.980900e-01, 4.777100e-01, 3.980900e-01, 3.980900e-01,
1569 4.777100e-01, 5.573200e-01, 1.035030e+00, 1.512740e+00, 2.149680e+00, 2.786620e+00,
1570 3.503180e+00, 4.219750e+00, 5.015920e+00, 5.652870e+00, 6.289810e+00, 7.245220e+00,
1571 8.121020e+00, 8.837580e+00, 9.633760e+00, 1.042994e+01, 1.114650e+01, 1.194268e+01,
1572 1.265924e+01, 1.329618e+01, 1.393312e+01, 1.449045e+01, 1.496815e+01, 1.552548e+01,
1573 1.592357e+01, 1.664013e+01, 1.727707e+01, 1.783439e+01, 1.831210e+01, 1.902866e+01,
1574 1.902866e+01, 1.878981e+01, 1.847134e+01, 1.831210e+01, 1.807325e+01, 1.791401e+01,
1575 1.783439e+01, 1.767516e+01, 1.759554e+01, 1.743631e+01, 1.743631e+01, 1.751592e+01,
1576 1.743631e+01, 1.735669e+01, 1.751592e+01, 1.759554e+01, 1.767516e+01, 1.783439e+01,
1577 1.783439e+01, 1.791401e+01, 1.815287e+01, 1.823248e+01, 1.847134e+01, 1.878981e+01,
1578 1.894905e+01, 1.902866e+01, 1.934713e+01, 1.966561e+01, 1.990446e+01, 2.014331e+01,
1579 2.030255e+01, 2.046178e+01, 2.085987e+01
1580 }; // 75
1581 
1582 const G4double G4HadronNucleonXsc::fKpNeutronTotTkin[75] = {
1583 2.692000e-02, 1.615400e-01, 2.961500e-01, 4.576900e-01, 6.461500e-01, 7.538500e-01,
1584 8.884600e-01, 1.103850e+00, 1.211540e+00, 1.400000e+00, 1.561540e+00, 1.776920e+00,
1585 1.992310e+00, 2.126920e+00, 2.342310e+00, 2.423080e+00, 2.557690e+00, 2.692310e+00,
1586 2.800000e+00, 2.988460e+00, 3.203850e+00, 3.365380e+00, 3.500000e+00, 3.688460e+00,
1587 3.850000e+00, 4.011540e+00, 4.173080e+00, 4.415380e+00, 4.630770e+00, 4.873080e+00,
1588 5.061540e+00, 5.276920e+00, 5.492310e+00, 5.707690e+00, 5.896150e+00, 6.030770e+00,
1589 6.138460e+00, 6.219230e+00, 6.273080e+00, 6.326920e+00, 6.407690e+00, 6.650000e+00,
1590 6.784620e+00, 7.026920e+00, 7.242310e+00, 7.350000e+00, 7.484620e+00, 7.619230e+00,
1591 7.807690e+00, 7.915380e+00, 8.050000e+00, 8.211540e+00, 8.453850e+00, 8.588460e+00,
1592 8.830770e+00, 9.073080e+00, 9.288460e+00, 9.476920e+00, 9.665380e+00, 9.826920e+00,
1593 1.004231e+01, 1.031154e+01, 1.052692e+01, 1.071538e+01, 1.095769e+01, 1.120000e+01,
1594 1.138846e+01, 1.155000e+01, 1.176538e+01, 1.190000e+01, 1.214231e+01, 1.222308e+01,
1595 1.238461e+01, 1.246538e+01, 1.265385e+01
1596 }; // 75
1597 
1598 const G4double G4HadronNucleonXsc::fKmProtonTotXsc[106] = {
1599 1.136585e+02, 9.749129e+01, 9.275262e+01, 8.885017e+01, 8.334146e+01, 7.955401e+01,
1600 7.504530e+01, 7.153658e+01, 6.858537e+01, 6.674913e+01, 6.525784e+01, 6.448781e+01,
1601 6.360279e+01, 6.255401e+01, 6.127526e+01, 6.032404e+01, 5.997910e+01, 5.443554e+01,
1602 5.376307e+01, 5.236934e+01, 5.113937e+01, 5.090941e+01, 4.967944e+01, 4.844948e+01,
1603 4.705575e+01, 4.638327e+01, 4.571080e+01, 4.475958e+01, 4.397213e+01, 4.257840e+01,
1604 4.102090e+01, 4.090592e+01, 3.906969e+01, 3.839721e+01, 3.756097e+01, 3.644599e+01,
1605 3.560976e+01, 3.533101e+01, 3.533101e+01, 3.644599e+01, 3.811847e+01, 3.839721e+01,
1606 3.979094e+01, 4.090592e+01, 4.257840e+01, 4.341463e+01, 4.425087e+01, 4.564460e+01,
1607 4.759582e+01, 4.703833e+01, 4.843206e+01, 4.787457e+01, 4.452962e+01, 4.202090e+01,
1608 4.034843e+01, 3.839721e+01, 3.616725e+01, 3.365854e+01, 3.170732e+01, 3.087108e+01,
1609 3.170732e+01, 3.254355e+01, 3.310104e+01, 3.254355e+01, 3.142857e+01, 3.059233e+01,
1610 2.947735e+01, 2.891986e+01, 2.836237e+01, 2.752613e+01, 2.696864e+01, 2.641115e+01,
1611 2.501742e+01, 2.473868e+01, 2.418118e+01, 2.362369e+01, 2.334495e+01, 2.278746e+01,
1612 2.250871e+01, 2.222997e+01, 2.167247e+01, 2.139373e+01, 2.139373e+01, 2.139373e+01,
1613 2.111498e+01, 2.083624e+01, 2.055749e+01, 2.083624e+01, 2.055749e+01, 2.083624e+01,
1614 2.083624e+01, 2.055749e+01, 2.055749e+01, 2.055749e+01, 2.027875e+01, 2.000000e+01,
1615 2.055749e+01, 2.027875e+01, 2.083624e+01, 2.083624e+01, 2.055749e+01, 2.083624e+01,
1616 2.083624e+01, 2.083624e+01, 2.139373e+01, 2.139373e+01
1617 }; // 106
1618 
1619 const G4double G4HadronNucleonXsc::fKmProtonTotTkin[106] = {
1620 4.017980e+00, 4.125840e+00, 4.179780e+00, 4.251690e+00, 4.287640e+00, 4.341570e+00,
1621 4.395510e+00, 4.467420e+00, 4.503370e+00, 4.575280e+00, 4.683150e+00, 4.737080e+00,
1622 4.773030e+00, 4.826970e+00, 4.880900e+00, 4.916850e+00, 4.952810e+00, 4.988760e+00,
1623 4.988760e+00, 5.006740e+00, 5.006740e+00, 5.042700e+00, 5.078650e+00, 5.114610e+00,
1624 5.132580e+00, 5.150560e+00, 5.186520e+00, 5.204490e+00, 5.276400e+00, 5.348310e+00,
1625 5.366290e+00, 5.384270e+00, 5.456180e+00, 5.564040e+00, 5.600000e+00, 5.671910e+00,
1626 5.743820e+00, 5.833710e+00, 5.905620e+00, 5.977530e+00, 6.085390e+00, 6.085390e+00,
1627 6.157300e+00, 6.175280e+00, 6.211240e+00, 6.229210e+00, 6.247190e+00, 6.337080e+00,
1628 6.391010e+00, 6.516850e+00, 6.462920e+00, 6.498880e+00, 6.570790e+00, 6.606740e+00,
1629 6.660670e+00, 6.678650e+00, 6.696630e+00, 6.732580e+00, 6.804490e+00, 6.876400e+00,
1630 6.948310e+00, 7.020220e+00, 7.074160e+00, 7.182020e+00, 7.235960e+00, 7.289890e+00,
1631 7.397750e+00, 7.523600e+00, 7.631460e+00, 7.757300e+00, 7.901120e+00, 8.062920e+00,
1632 8.260670e+00, 8.386520e+00, 8.530340e+00, 8.674160e+00, 8.817980e+00, 8.943820e+00,
1633 9.087640e+00, 9.267420e+00, 9.429210e+00, 9.573030e+00, 9.698880e+00, 9.896630e+00,
1634 1.002247e+01, 1.016629e+01, 1.031011e+01, 1.048989e+01, 1.063371e+01, 1.077753e+01,
1635 1.095730e+01, 1.108315e+01, 1.120899e+01, 1.135281e+01, 1.149663e+01, 1.162247e+01,
1636 1.174831e+01, 1.187416e+01, 1.200000e+01, 1.212584e+01, 1.221573e+01, 1.234157e+01,
1637 1.239551e+01, 1.250337e+01, 1.261124e+01, 1.273708e+01
1638 }; // 106
1639 
1640 const G4double G4HadronNucleonXsc::fKmNeutronTotXsc[68] = {
1641 2.621810e+01, 2.741123e+01, 2.868413e+01, 2.963889e+01, 3.067343e+01, 3.178759e+01,
1642 3.282148e+01, 3.417466e+01, 3.536778e+01, 3.552620e+01, 3.544576e+01, 3.496756e+01,
1643 3.433030e+01, 3.401166e+01, 3.313537e+01, 3.257772e+01, 3.178105e+01, 3.138264e+01,
1644 3.074553e+01, 2.970952e+01, 2.891301e+01, 2.827542e+01, 2.787700e+01, 2.715978e+01,
1645 2.660181e+01, 2.612394e+01, 2.564574e+01, 2.516721e+01, 2.421098e+01, 2.365235e+01,
1646 2.317366e+01, 2.261437e+01, 2.237389e+01, 2.205427e+01, 2.181395e+01, 2.165357e+01,
1647 2.149335e+01, 2.133297e+01, 2.109232e+01, 2.093128e+01, 2.069030e+01, 2.052992e+01,
1648 2.028927e+01, 2.012824e+01, 1.996737e+01, 1.996590e+01, 1.988530e+01, 1.964432e+01,
1649 1.948361e+01, 1.940236e+01, 1.940040e+01, 1.931882e+01, 1.947593e+01, 1.947429e+01,
1650 1.939320e+01, 1.939157e+01, 1.946922e+01, 1.962715e+01, 1.970481e+01, 1.970301e+01,
1651 1.993958e+01, 2.009669e+01, 2.025380e+01, 2.033178e+01, 2.049003e+01, 2.064747e+01,
1652 2.080540e+01, 2.096333e+01
1653 }; // 68
1654 
1655 const G4double G4HadronNucleonXsc::fKmNeutronTotTkin[68] = {
1656 5.708500e+00, 5.809560e+00, 5.896270e+00, 5.954120e+00, 5.997630e+00, 6.041160e+00,
1657 6.142160e+00, 6.171410e+00, 6.272470e+00, 6.344390e+00, 6.416230e+00, 6.459180e+00,
1658 6.487690e+00, 6.501940e+00, 6.544740e+00, 6.573280e+00, 6.616110e+00, 6.644710e+00,
1659 6.658840e+00, 6.744700e+00, 6.773150e+00, 6.830410e+00, 6.859010e+00, 6.916240e+00,
1660 6.973530e+00, 6.987730e+00, 7.030670e+00, 7.102360e+00, 7.173880e+00, 7.288660e+00,
1661 7.374720e+00, 7.547000e+00, 7.690650e+00, 7.791150e+00, 7.920420e+00, 8.020980e+00,
1662 8.107160e+00, 8.207720e+00, 8.365740e+00, 8.523790e+00, 8.710560e+00, 8.811110e+00,
1663 8.969140e+00, 9.127190e+00, 9.270860e+00, 9.400230e+00, 9.486440e+00, 9.673210e+00,
1664 9.802510e+00, 9.946220e+00, 1.011870e+01, 1.029116e+01, 1.047808e+01, 1.062181e+01,
1665 1.075114e+01, 1.089488e+01, 1.106739e+01, 1.118244e+01, 1.135496e+01, 1.151307e+01,
1666 1.171439e+01, 1.190130e+01, 1.208822e+01, 1.223199e+01, 1.231829e+01, 1.247646e+01,
1667 1.259150e+01, 1.270655e+01
1668 }; // 68
1669 
1670 
1671 
1672 
1673 
1674 //
1675 //
1676 ///////////////////////////////////////////////////////////////////////////////////////
G4LPhysicsFreeVector::PutValues
void PutValues(size_t binNumber, G4double binValue, G4double dataValue)
G4HadronNucleonXsc::GetCoulombBarrier
G4double GetCoulombBarrier(const G4DynamicParticle *aParticle, const G4ParticleDefinition *nucleon)
Definition: G4HadronNucleonXsc.cc:1413
G4AntiOmegaMinus::AntiOmegaMinus
static G4AntiOmegaMinus * AntiOmegaMinus()
Definition: G4AntiOmegaMinus.cc:110
G4DynamicParticle::GetKineticEnergy
G4double GetKineticEnergy() const
G4DynamicParticle::GetTotalEnergy
G4double GetTotalEnergy() const
outFile
std::ofstream outFile
Definition: GammaRayTel.cc:68
G4Element::GetN
G4double GetN() const
Definition: G4Element.hh:134
G4OmegaMinus::OmegaMinus
static G4OmegaMinus * OmegaMinus()
Definition: G4OmegaMinus.cc:110
G4Element::GetZ
G4double GetZ() const
Definition: G4Element.hh:131
G4KaonZeroLong::KaonZeroLong
static G4KaonZeroLong * KaonZeroLong()
Definition: G4KaonZeroLong.cc:115
neutron
Definition: G4DataQuestionaire.hh:35
G4DynamicParticle::GetDefinition
G4ParticleDefinition * GetDefinition() const
G4AntiSigmaPlus::AntiSigmaPlus
static G4AntiSigmaPlus * AntiSigmaPlus()
Definition: G4AntiSigmaPlus.cc:108
G4int
int G4int
Definition: G4Types.hh:78
python.hepunit.millibarn
int millibarn
Definition: hepunit.py:40
G4IonTable.hh
G4SigmaZero::SigmaZero
static G4SigmaZero * SigmaZero()
Definition: G4SigmaZero.cc:99
G4KaonMinus::KaonMinus
static G4KaonMinus * KaonMinus()
Definition: G4KaonMinus.cc:113
G4HadronNucleonXsc::GetHadronNucleonXscNS
G4double GetHadronNucleonXscNS(const G4DynamicParticle *, const G4ParticleDefinition *)
Definition: G4HadronNucleonXsc.cc:365
G4PhysicsVector::SetSpline
void SetSpline(G4bool)
G4AntiSigmaMinus::AntiSigmaMinus
static G4AntiSigmaMinus * AntiSigmaMinus()
Definition: G4AntiSigmaMinus.cc:106
G4XiZero::XiZero
static G4XiZero * XiZero()
Definition: G4XiZero.cc:106
G4ParticleTable.hh
G4HadronNucleonXsc::GetKaonNucleonXscGG
G4double GetKaonNucleonXscGG(const G4DynamicParticle *, const G4ParticleDefinition *)
Definition: G4HadronNucleonXsc.cc:1069
G4KaonZeroShort::KaonZeroShort
static G4KaonZeroShort * KaonZeroShort()
Definition: G4KaonZeroShort.cc:104
G4InuclParticleNames::nucleon
G4bool nucleon(G4int ityp)
Definition: G4InuclParticleNames.hh:95
G4AntiProton::AntiProton
static G4AntiProton * AntiProton()
Definition: G4AntiProton.cc:93
G4DynamicParticle::GetMass
G4double GetMass() const
G4bool
bool G4bool
Definition: G4Types.hh:79
G4HadronNucleonXsc::CalculateEcmValue
G4double CalculateEcmValue(const G4double, const G4double, const G4double)
Definition: G4HadronNucleonXsc.cc:1381
G4XiMinus::XiMinus
static G4XiMinus * XiMinus()
Definition: G4XiMinus.cc:106
G4AntiXiMinus::AntiXiMinus
static G4AntiXiMinus * AntiXiMinus()
Definition: G4AntiXiMinus.cc:105
G4Triton::Triton
static G4Triton * Triton()
Definition: G4Triton.cc:95
G4Proton::Proton
static G4Proton * Proton()
Definition: G4Proton.cc:93
G4PionPlus::PionPlus
static G4PionPlus * PionPlus()
Definition: G4PionPlus.cc:98
G4HadronNucleonXsc::IsIsoApplicable
virtual G4bool IsIsoApplicable(const G4DynamicParticle *aDP, G4int Z, G4int A)
Definition: G4HadronNucleonXsc.cc:107
G4ParticleDefinition.hh
G4HadronNucleonXsc::CrossSectionDescription
void CrossSectionDescription(std::ostream &) const
Definition: G4HadronNucleonXsc.cc:84
G4Neutron::Neutron
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
G4Gamma::Gamma
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
G4PionZero::PionZero
static G4PionZero * PionZero()
Definition: G4PionZero.cc:104
G4Deuteron::Deuteron
static G4Deuteron * Deuteron()
Definition: G4Deuteron.cc:94
G4PhysicalConstants.hh
G4HadronNucleonXsc::GetHadronNucleonXscVU
G4double GetHadronNucleonXscVU(const G4DynamicParticle *, const G4ParticleDefinition *)
Definition: G4HadronNucleonXsc.cc:1234
G4SigmaMinus::SigmaMinus
static G4SigmaMinus * SigmaMinus()
Definition: G4SigmaMinus.cc:106
G4HadronNucleonXsc::IsApplicable
virtual G4bool IsApplicable(const G4DynamicParticle *aDP, const G4Element *)
Definition: G4HadronNucleonXsc.cc:95
G4ParticleDefinition::GetPDGMass
G4double GetPDGMass() const
Definition: G4ParticleDefinition.hh:161
G4AntiLambda::AntiLambda
static G4AntiLambda * AntiLambda()
Definition: G4AntiLambda.cc:108
G4lrint
int G4lrint(double ad)
Definition: templates.hh:163
G4PionMinus::PionMinus
static G4PionMinus * PionMinus()
Definition: G4PionMinus.cc:98
G4AntiSigmaZero::AntiSigmaZero
static G4AntiSigmaZero * AntiSigmaZero()
Definition: G4AntiSigmaZero.cc:99
G4HadronNucleonXsc.hh
G4AntiXiZero::AntiXiZero
static G4AntiXiZero * AntiXiZero()
Definition: G4AntiXiZero.cc:106
G4HadronNucleonXsc::GetHadronNucleonXscPDG
G4double GetHadronNucleonXscPDG(const G4DynamicParticle *, const G4ParticleDefinition *)
Definition: G4HadronNucleonXsc.cc:221
G4Alpha::Alpha
static G4Alpha * Alpha()
Definition: G4Alpha.cc:89
LE
Definition: Evaluator.cc:66
G4SigmaPlus::SigmaPlus
static G4SigmaPlus * SigmaPlus()
Definition: G4SigmaPlus.cc:108
G4Lambda::Lambda
static G4Lambda * Lambda()
Definition: G4Lambda.cc:108
G4double
double G4double
Definition: G4Types.hh:76
G4SystemOfUnits.hh
G4KaonPlus::KaonPlus
static G4KaonPlus * KaonPlus()
Definition: G4KaonPlus.cc:113
G4ParticleDefinition::GetPDGCharge
G4double GetPDGCharge() const
Definition: G4ParticleDefinition.hh:163
G4He3::He3
static G4He3 * He3()
Definition: G4He3.cc:94
CLHEP::Hep3Vector::mag
double mag() const
G4HadronNucleonXsc::GetHadronNucleonXscEL
G4double GetHadronNucleonXscEL(const G4DynamicParticle *, const G4ParticleDefinition *)
Definition: G4HadronNucleonXsc.cc:143
G4AntiNeutron::AntiNeutron
static G4AntiNeutron * AntiNeutron()
Definition: G4AntiNeutron.cc:103
G4DynamicParticle::GetMomentum
G4ThreeVector GetMomentum() const
G4HadTmpUtil.hh
G4HadronNucleonXsc::CalcMandelstamS
G4double CalcMandelstamS(const G4double, const G4double, const G4double)
Definition: G4HadronNucleonXsc.cc:1398
Generated on Wed Apr 30 2014 15:54:48 for Geant4.10 by   doxygen 1.8.7