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

#include <G4HadronNucleonXsc.hh>

Public Member Functions

 G4HadronNucleonXsc ()
 
virtual ~G4HadronNucleonXsc ()
 
virtual G4bool IsApplicable (const G4DynamicParticle *aDP, const G4Element *)
 
virtual G4bool IsIsoApplicable (const G4DynamicParticle *aDP, G4int Z, G4int A)
 
virtual void DumpPhysicsTable (const G4ParticleDefinition &)
 
void CrossSectionDescription (std::ostream &) const
 
G4double GetHadronNucleonXscEL (const G4DynamicParticle *, const G4ParticleDefinition *)
 
G4double GetHadronNucleonXscPDG (const G4DynamicParticle *, const G4ParticleDefinition *)
 
G4double GetHadronNucleonXscNS (const G4DynamicParticle *, const G4ParticleDefinition *)
 
G4double GetKaonNucleonXscGG (const G4DynamicParticle *, const G4ParticleDefinition *)
 
G4double GetHadronNucleonXscVU (const G4DynamicParticle *, const G4ParticleDefinition *)
 
G4double CalculateEcmValue (const G4double, const G4double, const G4double)
 
G4double CalcMandelstamS (const G4double, const G4double, const G4double)
 
G4double GetCoulombBarrier (const G4DynamicParticle *aParticle, const G4ParticleDefinition *nucleon)
 
G4double GetTotalHadronNucleonXsc ()
 
G4double GetElasticHadronNucleonXsc ()
 
G4double GetInelasticHadronNucleonXsc ()
 

Detailed Description

Definition at line 51 of file G4HadronNucleonXsc.hh.

Constructor & Destructor Documentation

G4HadronNucleonXsc::G4HadronNucleonXsc ( )

Definition at line 41 of file G4HadronNucleonXsc.cc.

42 :
43 // fUpperLimit( 10000 * GeV ),
44  fLowerLimit( 0.03 * MeV ),
45  fTotalXsc(0.0), fElasticXsc(0.0), fInelasticXsc(0.0)
46 // , fHadronNucleonXsc(0.0)
47 {
48  theGamma = G4Gamma::Gamma();
49  theProton = G4Proton::Proton();
50  theNeutron = G4Neutron::Neutron();
51  theAProton = G4AntiProton::AntiProton();
52  theANeutron = G4AntiNeutron::AntiNeutron();
53  thePiPlus = G4PionPlus::PionPlus();
54  thePiMinus = G4PionMinus::PionMinus();
55  thePiZero = G4PionZero::PionZero();
56  theKPlus = G4KaonPlus::KaonPlus();
57  theKMinus = G4KaonMinus::KaonMinus();
58  theK0S = G4KaonZeroShort::KaonZeroShort();
59  theK0L = G4KaonZeroLong::KaonZeroLong();
60  theL = G4Lambda::Lambda();
61  theAntiL = G4AntiLambda::AntiLambda();
62  theSPlus = G4SigmaPlus::SigmaPlus();
63  theASPlus = G4AntiSigmaPlus::AntiSigmaPlus();
64  theSMinus = G4SigmaMinus::SigmaMinus();
65  theASMinus = G4AntiSigmaMinus::AntiSigmaMinus();
66  theS0 = G4SigmaZero::SigmaZero();
67  theAS0 = G4AntiSigmaZero::AntiSigmaZero();
68  theXiMinus = G4XiMinus::XiMinus();
69  theXi0 = G4XiZero::XiZero();
70  theAXiMinus = G4AntiXiMinus::AntiXiMinus();
71  theAXi0 = G4AntiXiZero::AntiXiZero();
72  theOmega = G4OmegaMinus::OmegaMinus();
73  theAOmega = G4AntiOmegaMinus::AntiOmegaMinus();
74  theD = G4Deuteron::Deuteron();
75  theT = G4Triton::Triton();
76  theA = G4Alpha::Alpha();
77  theHe3 = G4He3::He3();
78 
79  // InitialiseKaonNucleonTotXsc();
80 }
G4AntiOmegaMinus::AntiOmegaMinus
static G4AntiOmegaMinus * AntiOmegaMinus()
Definition: G4AntiOmegaMinus.cc:110
G4OmegaMinus::OmegaMinus
static G4OmegaMinus * OmegaMinus()
Definition: G4OmegaMinus.cc:110
G4KaonZeroLong::KaonZeroLong
static G4KaonZeroLong * KaonZeroLong()
Definition: G4KaonZeroLong.cc:115
G4AntiSigmaPlus::AntiSigmaPlus
static G4AntiSigmaPlus * AntiSigmaPlus()
Definition: G4AntiSigmaPlus.cc:108
G4SigmaZero::SigmaZero
static G4SigmaZero * SigmaZero()
Definition: G4SigmaZero.cc:102
G4KaonMinus::KaonMinus
static G4KaonMinus * KaonMinus()
Definition: G4KaonMinus.cc:113
G4AntiSigmaMinus::AntiSigmaMinus
static G4AntiSigmaMinus * AntiSigmaMinus()
Definition: G4AntiSigmaMinus.cc:106
G4XiZero::XiZero
static G4XiZero * XiZero()
Definition: G4XiZero.cc:106
G4KaonZeroShort::KaonZeroShort
static G4KaonZeroShort * KaonZeroShort()
Definition: G4KaonZeroShort.cc:104
G4AntiProton::AntiProton
static G4AntiProton * AntiProton()
Definition: G4AntiProton.cc:93
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
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:108
G4Deuteron::Deuteron
static G4Deuteron * Deuteron()
Definition: G4Deuteron.cc:94
G4SigmaMinus::SigmaMinus
static G4SigmaMinus * SigmaMinus()
Definition: G4SigmaMinus.cc:106
G4AntiLambda::AntiLambda
static G4AntiLambda * AntiLambda()
Definition: G4AntiLambda.cc:108
G4PionMinus::PionMinus
static G4PionMinus * PionMinus()
Definition: G4PionMinus.cc:98
G4AntiSigmaZero::AntiSigmaZero
static G4AntiSigmaZero * AntiSigmaZero()
Definition: G4AntiSigmaZero.cc:103
G4AntiXiZero::AntiXiZero
static G4AntiXiZero * AntiXiZero()
Definition: G4AntiXiZero.cc:106
MeV
static constexpr double MeV
Definition: G4SIunits.hh:214
G4Alpha::Alpha
static G4Alpha * Alpha()
Definition: G4Alpha.cc:89
G4SigmaPlus::SigmaPlus
static G4SigmaPlus * SigmaPlus()
Definition: G4SigmaPlus.cc:108
G4Lambda::Lambda
static G4Lambda * Lambda()
Definition: G4Lambda.cc:108
G4KaonPlus::KaonPlus
static G4KaonPlus * KaonPlus()
Definition: G4KaonPlus.cc:113
G4He3::He3
static G4He3 * He3()
Definition: G4He3.cc:94
G4AntiNeutron::AntiNeutron
static G4AntiNeutron * AntiNeutron()
Definition: G4AntiNeutron.cc:103

Here is the call graph for this function:

G4HadronNucleonXsc::~G4HadronNucleonXsc ( )
virtual

Definition at line 83 of file G4HadronNucleonXsc.cc.

84 {}

Member Function Documentation

G4double G4HadronNucleonXsc::CalcMandelstamS ( const G4double  mp,
const G4double  mt,
const G4double  Plab 
)

Definition at line 1418 of file G4HadronNucleonXsc.cc.

1421 {
1422  G4double Elab = std::sqrt ( mp * mp + Plab * Plab );
1423  G4double sMand = mp*mp + mt*mt + 2*Elab*mt ;
1424 
1425  return sMand;
1426 }
G4double
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:

G4double G4HadronNucleonXsc::CalculateEcmValue ( const G4double  mp,
const G4double  mt,
const G4double  Plab 
)

Definition at line 1401 of file G4HadronNucleonXsc.cc.

1404 {
1405  G4double Elab = std::sqrt ( mp * mp + Plab * Plab );
1406  G4double Ecm = std::sqrt ( mp * mp + mt * mt + 2 * Elab * mt );
1407  // G4double Pcm = Plab * mt / Ecm;
1408  // G4double KEcm = std::sqrt ( Pcm * Pcm + mp * mp ) - mp;
1409 
1410  return Ecm ; // KEcm;
1411 }
G4double
double G4double
Definition: G4Types.hh:76
void G4HadronNucleonXsc::CrossSectionDescription ( std::ostream &  outFile) const

Definition at line 86 of file G4HadronNucleonXsc.cc.

87 {
88  outFile << "G4HadronNucleonXsc calculates the total, inelastic and elastic\n"
89  << "hadron-nucleon cross sections using several different\n"
90  << "parameterizations within the Glauber-Gribov approach. It is\n"
91  << "valid for all incident gammas and long-lived hadrons at\n"
92  << "energies above 30 keV. This is a cross section component which\n"
93  << "is to be used to build a cross section data set.\n";
94 }
virtual void G4HadronNucleonXsc::DumpPhysicsTable ( const G4ParticleDefinition &  )
inlinevirtual

Definition at line 65 of file G4HadronNucleonXsc.hh.

66  {G4cout << "G4HadronNucleonXsc: uses parametrisation"<<G4endl;}
G4cout
G4GLOB_DLL std::ostream G4cout
G4endl
#define G4endl
Definition: G4ios.hh:61
G4double G4HadronNucleonXsc::GetCoulombBarrier ( const G4DynamicParticle *  aParticle,
const G4ParticleDefinition *  nucleon 
)

Definition at line 1433 of file G4HadronNucleonXsc.cc.

1435 {
1436  G4double ratio;
1437 
1438  G4double tR = 0.895*fermi, pR;
1439 
1440  if ( aParticle->GetDefinition() == theProton ) pR = 0.895*fermi;
1441  else if( aParticle->GetDefinition() == thePiPlus ) pR = 0.663*fermi;
1442  else if( aParticle->GetDefinition() == theKPlus ) pR = 0.340*fermi;
1443  else pR = 0.500*fermi;
1444 
1445  G4double pZ = aParticle->GetDefinition()->GetPDGCharge();
1446  G4double tZ = nucleon->GetPDGCharge();
1447 
1448  G4double pTkin = aParticle->GetKineticEnergy();
1449 
1450  G4double pM = aParticle->GetDefinition()->GetPDGMass();
1451  G4double tM = nucleon->GetPDGMass();
1452 
1453  G4double pElab = pTkin + pM;
1454 
1455  G4double totEcm = std::sqrt(pM*pM + tM*tM + 2.*pElab*tM);
1456 
1457  G4double totTcm = totEcm - pM -tM;
1458 
1459  G4double bC = fine_structure_const*hbarc*pZ*tZ;
1460  bC /= pR + tR;
1461  bC /= 2.; // 4., 2. parametrisation cof ??? vmg
1462 
1463  // G4cout<<"pTkin = "<<pTkin/GeV<<"; pPlab = "
1464  // <<pPlab/GeV<<"; bC = "<<bC/GeV<<"; pTcm = "<<pTcm/GeV<<G4endl;
1465 
1466  if( totTcm <= bC ) ratio = 0.;
1467  else ratio = 1. - bC/totTcm;
1468 
1469  // if(ratio < DBL_MIN) ratio = DBL_MIN;
1470  if( ratio < 0.) ratio = 0.;
1471 
1472  // G4cout <<"ratio = "<<ratio<<G4endl;
1473  return ratio;
1474 }
G4DynamicParticle::GetKineticEnergy
G4double GetKineticEnergy() const
CLHEP::hbarc
static constexpr double hbarc
Definition: PhysicalConstants.h:65
G4DynamicParticle::GetDefinition
G4ParticleDefinition * GetDefinition() const
G4ParticleDefinition::GetPDGMass
G4double GetPDGMass() const
Definition: G4ParticleDefinition.hh:122
G4double
double G4double
Definition: G4Types.hh:76
CLHEP::fine_structure_const
static constexpr double fine_structure_const
Definition: PhysicalConstants.h:96
fermi
static constexpr double fermi
Definition: G4SIunits.hh:103
G4ParticleDefinition::GetPDGCharge
G4double GetPDGCharge() const
Definition: G4ParticleDefinition.hh:124

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4HadronNucleonXsc::GetElasticHadronNucleonXsc ( )
inline

Definition at line 91 of file G4HadronNucleonXsc.hh.

91 { return fElasticXsc; };

Here is the caller graph for this function:

G4double G4HadronNucleonXsc::GetHadronNucleonXscEL ( const G4DynamicParticle *  aParticle,
const G4ParticleDefinition *  nucleon 
)

Definition at line 145 of file G4HadronNucleonXsc.cc.

147 {
148  G4double xsection;
149 
150 
151  G4double targ_mass = 0.939*GeV; // ~mean neutron and proton ???
152 
153  G4double proj_mass = aParticle->GetMass();
154  G4double proj_momentum = aParticle->GetMomentum().mag();
155  G4double sMand = CalcMandelstamS ( proj_mass , targ_mass , proj_momentum );
156 
157  sMand /= GeV*GeV; // in GeV for parametrisation
158  proj_momentum /= GeV;
159 
160  const G4ParticleDefinition* theParticle = aParticle->GetDefinition();
161 
162  G4bool pORn = (nucleon == theProton || nucleon == theNeutron );
163 
164 
165  if(theParticle == theGamma && pORn )
166  {
167  xsection = (0.0677*std::pow(sMand,0.0808) + 0.129*std::pow(sMand,-0.4525));
168  }
169  else if(theParticle == theNeutron && pORn ) // as proton ???
170  {
171  xsection = (21.70*std::pow(sMand,0.0808) + 56.08*std::pow(sMand,-0.4525));
172  }
173  else if(theParticle == theProton && pORn )
174  {
175  xsection = (21.70*std::pow(sMand,0.0808) + 56.08*std::pow(sMand,-0.4525));
176 
177  // xsection = At*( 49.51*std::pow(sMand,-0.097) + 0.314*G4Log(sMand)*G4Log(sMand) );
178  // xsection = At*( 38.4 + 0.85*std::abs(std::pow(log(sMand),1.47)) );
179  }
180  else if(theParticle == theAProton && pORn )
181  {
182  xsection = ( 21.70*std::pow(sMand,0.0808) + 98.39*std::pow(sMand,-0.4525));
183  }
184  else if(theParticle == thePiPlus && pORn )
185  {
186  xsection = (13.63*std::pow(sMand,0.0808) + 27.56*std::pow(sMand,-0.4525));
187  }
188  else if(theParticle == thePiMinus && pORn )
189  {
190  // xsection = At*( 55.2*std::pow(sMand,-0.255) + 0.346*G4Log(sMand)*G4Log(sMand) );
191  xsection = (13.63*std::pow(sMand,0.0808) + 36.02*std::pow(sMand,-0.4525));
192  }
193  else if(theParticle == theKPlus && pORn )
194  {
195  xsection = (11.82*std::pow(sMand,0.0808) + 8.15*std::pow(sMand,-0.4525));
196  }
197  else if(theParticle == theKMinus && pORn )
198  {
199  xsection = (11.82*std::pow(sMand,0.0808) + 26.36*std::pow(sMand,-0.4525));
200  }
201  else // as proton ???
202  {
203  xsection = (21.70*std::pow(sMand,0.0808) + 56.08*std::pow(sMand,-0.4525));
204  }
205  xsection *= millibarn;
206 
207  fTotalXsc = xsection;
208  fInelasticXsc = 0.83*xsection;
209  fElasticXsc = fTotalXsc - fInelasticXsc;
210  if (fElasticXsc < 0.)fElasticXsc = 0.;
211 
212  return xsection;
213 }
G4DynamicParticle::GetDefinition
G4ParticleDefinition * GetDefinition() const
G4DynamicParticle::GetMass
G4double GetMass() const
G4bool
bool G4bool
Definition: G4Types.hh:79
GeV
static constexpr double GeV
Definition: G4SIunits.hh:217
G4double
double G4double
Definition: G4Types.hh:76
CLHEP::Hep3Vector::mag
double mag() const
millibarn
static constexpr double millibarn
Definition: G4SIunits.hh:106
G4DynamicParticle::GetMomentum
G4ThreeVector GetMomentum() const
G4HadronNucleonXsc::CalcMandelstamS
G4double CalcMandelstamS(const G4double, const G4double, const G4double)
Definition: G4HadronNucleonXsc.cc:1418

Here is the call graph for this function:

G4double G4HadronNucleonXsc::GetHadronNucleonXscNS ( const G4DynamicParticle *  aParticle,
const G4ParticleDefinition *  nucleon 
)

pi- ////////////////////////////////////////////

Definition at line 367 of file G4HadronNucleonXsc.cc.

369 {
370  G4double xsection(0);
371 
372  G4double A0, B0;
373  G4double hpXsc(0);
374  G4double hnXsc(0);
375 
376 
377  G4double tM = 0.939*GeV; // ~mean neutron and proton ???
378 
379  G4double pM = aParticle->GetMass();
380  G4double pE = aParticle->GetTotalEnergy();
381  G4double pLab = aParticle->GetMomentum().mag();
382 
383  G4double sMand = CalcMandelstamS ( pM , tM , pLab );
384 
385  sMand /= GeV*GeV; // in GeV for parametrisation
386  pLab /= GeV;
387  pE /= GeV;
388  pM /= GeV;
389 
390  G4double logP = G4Log(pLab);
391 
392 
393  // General PDG fit constants
394 
395  G4double s0 = 5.38*5.38; // in Gev^2
396  G4double eta1 = 0.458;
397  G4double eta2 = 0.458;
398  G4double B = 0.308;
399  G4double minLogP = 3.5; // min of (lnP-minLogP)^2
400  G4double cofLogE = .0557; // elastic (lnP-minLogP)^2
401  G4double cofLogT = .3; // total (lnP-minLogP)^2
402  G4double pMin = .1; // fast LE calculation
403  G4double pMax = 1000.; // fast HE calculation
404 
405 
406  const G4ParticleDefinition* theParticle = aParticle->GetDefinition();
407 
408  G4bool pORn = (nucleon == theProton || nucleon == theNeutron );
409  G4bool proton = (nucleon == theProton);
410  G4bool neutron = (nucleon == theNeutron);
411 
412  if( theParticle == theNeutron && pORn )
413  {
414  if( pLab >= 373.)
415  {
416  xsection = GetHadronNucleonXscPDG(aParticle, nucleon)/millibarn;
417 
418  fElasticXsc = 6.5 + 0.308*std::pow(G4Log(sMand/400.),1.65) + 9.19*std::pow(sMand,-0.458);
419 
420  fTotalXsc = xsection;
421 
422  }
423  else if( pLab >= 100.)
424  {
425  B0 = 7.5;
426  A0 = 100. - B0*G4Log(3.0e7);
427 
428  xsection = A0 + B0*G4Log(pE) - 11
429  // + 103*std::pow(2*0.93827*pE + pM*pM+0.93827*0.93827,-0.165); // mb
430  + 103*std::pow(sMand,-0.165); // mb
431 
432  fElasticXsc = 5.53 + 0.308*std::pow(G4Log(sMand/28.9),1.1) + 9.19*std::pow(sMand,-0.458);
433 
434  fTotalXsc = xsection;
435  }
436  else if( pLab >= 10.)
437  {
438  B0 = 7.5;
439  A0 = 100. - B0*G4Log(3.0e7);
440 
441  xsection = A0 + B0*G4Log(pE) - 11
442  + 103*std::pow(2*0.93827*pE + pM*pM+
443  0.93827*0.93827,-0.165); // mb
444  fTotalXsc = xsection;
445  fElasticXsc = 6 + 20/( (logP-0.182)*(logP-0.182) + 1.0 );
446  }
447  else // pLab < 10 GeV/c
448  {
449  if( neutron ) // nn to be pp
450  {
451  if( pLab < 0.4 )
452  {
453  hnXsc = 23 + 50*( std::pow( G4Log(0.73/pLab), 3.5 ) );
454  fElasticXsc = hnXsc;
455  }
456  else if( pLab < 0.73 )
457  {
458  hnXsc = 23 + 50*( std::pow( G4Log(0.73/pLab), 3.5 ) );
459  fElasticXsc = hnXsc;
460  }
461  else if( pLab < 1.05 )
462  {
463  hnXsc = 23 + 40*(G4Log(pLab/0.73))*
464  (G4Log(pLab/0.73));
465  fElasticXsc = 23 + 20*(G4Log(pLab/0.73))*
466  (G4Log(pLab/0.73));
467  }
468  else // 1.05 - 10 GeV/c
469  {
470  hnXsc = 39.0+75*(pLab - 1.2)/(std::pow(pLab,3.0) + 0.15);
471 
472  fElasticXsc = 6 + 20/( (logP-0.182)*(logP-0.182) + 1.0 );
473  }
474  fTotalXsc = hnXsc;
475  }
476  if( proton ) // pn to be np
477  {
478  if( pLab < 0.02 )
479  {
480  hpXsc = 4100+30*std::pow(G4Log(1.3/pLab),3.6); // was as pLab < 0.8
481  fElasticXsc = hpXsc;
482  }
483  else if( pLab < 0.8 )
484  {
485  hpXsc = 33+30*std::pow(G4Log(pLab/1.3),4.0);
486  fElasticXsc = hpXsc;
487  }
488  else if( pLab < 1.05 )
489  {
490  hpXsc = 33+30*std::pow(G4Log(pLab/0.95),2.0);
491  fElasticXsc = 6 + 52/( G4Log(0.511/pLab)*G4Log(0.511/pLab) + 1.6 );
492  }
493  else if( pLab < 1.4 )
494  {
495  hpXsc = 33+30*std::pow(G4Log(pLab/0.95),2.0);
496  fElasticXsc = 6 + 52/( G4Log(0.511/pLab)*G4Log(0.511/pLab) + 1.6 );
497  }
498  else // 1.4 < pLab < 10. )
499  {
500  hpXsc = 33.3 + 20.8*(std::pow(pLab,2.0) - 1.35)/(std::pow(pLab,2.50) + 0.95);
501 
502  fElasticXsc = 6 + 20/( (logP-0.182)*(logP-0.182) + 1.0 );
503  }
504  fTotalXsc = hpXsc;
505  }
506  }
507  }
508  else if( theParticle == theProton && pORn )
509  {
510  if( pLab >= 373.) // pdg due to TOTEM data
511  {
512  xsection = GetHadronNucleonXscPDG(aParticle, nucleon)/millibarn;
513 
514  fElasticXsc = 6.5 + 0.308*std::pow(G4Log(sMand/400.),1.65) + 9.19*std::pow(sMand,-0.458);
515 
516  fTotalXsc = xsection;
517  }
518  else if( pLab >= 100.)
519  {
520  B0 = 7.5;
521  A0 = 100. - B0*G4Log(3.0e7);
522 
523  xsection = A0 + B0*G4Log(pE) - 11 + 103*std::pow(sMand,-0.165); // mb
524 
525  fElasticXsc = 5.53 + 0.308*std::pow(G4Log(sMand/28.9),1.1) + 9.19*std::pow(sMand,-0.458);
526 
527  fTotalXsc = xsection;
528  }
529  else if( pLab >= 10.)
530  {
531  B0 = 7.5;
532  A0 = 100. - B0*G4Log(3.0e7);
533 
534  xsection = A0 + B0*G4Log(pE) - 11 + 103*std::pow(sMand,-0.165); // mb
535 
536  fElasticXsc = 6 + 20/( (logP-0.182)*(logP-0.182) + 1.0 );
537 
538  fTotalXsc = xsection;
539  }
540  else
541  {
542  // pp
543 
544  if( proton )
545  {
546  if( pLab < 0.4 )
547  {
548  hpXsc = 23 + 50*( std::pow( G4Log(0.73/pLab), 3.5 ) );
549  fElasticXsc = hpXsc;
550  }
551  else if( pLab < 0.73 )
552  {
553  hpXsc = 23 + 50*( std::pow( G4Log(0.73/pLab), 3.5 ) );
554  fElasticXsc = hpXsc;
555  }
556  else if( pLab < 1.05 )
557  {
558  hpXsc = 23 + 40*(G4Log(pLab/0.73))*
559  (G4Log(pLab/0.73));
560  fElasticXsc = 23 + 20*(G4Log(pLab/0.73))*
561  (G4Log(pLab/0.73));
562  }
563  else // 1.05 - 10 GeV/c
564  {
565  hpXsc = 39.0+75*(pLab - 1.2)/(std::pow(pLab,3.0) + 0.15);
566 
567  fElasticXsc = 6 + 20/( (logP-0.182)*(logP-0.182) + 1.0 );
568  }
569  fTotalXsc = hpXsc;
570  }
571  if( neutron ) // pn to be np
572  {
573  if( pLab < 0.02 )
574  {
575  hnXsc = 4100+30*std::pow(G4Log(1.3/pLab),3.6); // was as pLab < 0.8
576  fElasticXsc = hnXsc;
577  }
578  else if( pLab < 0.8 )
579  {
580  hnXsc = 33+30*std::pow(G4Log(pLab/1.3),4.0);
581  fElasticXsc = hnXsc;
582  }
583  else if( pLab < 1.05 )
584  {
585  hnXsc = 33+30*std::pow(G4Log(pLab/0.95),2.0);
586  fElasticXsc = 6 + 52/( G4Log(0.511/pLab)*G4Log(0.511/pLab) + 1.6 );
587  }
588  else if( pLab < 1.4 )
589  {
590  hnXsc = 33+30*std::pow(G4Log(pLab/0.95),2.0);
591  fElasticXsc = 6 + 52/( G4Log(0.511/pLab)*G4Log(0.511/pLab) + 1.6 );
592  }
593  else // 1.4 < pLab < 10. )
594  {
595  hnXsc = 33.3 + 20.8*(std::pow(pLab,2.0) - 1.35)/(std::pow(pLab,2.50) + 0.95);
596 
597  fElasticXsc = 6 + 20/( (logP-0.182)*(logP-0.182) + 1.0 );
598  }
599  fTotalXsc = hnXsc;
600  }
601  }
602  }
603  else if( theParticle == theAProton && pORn )
604  {
605  if( proton )
606  {
607  xsection = 35.45 + B*std::pow(G4Log(sMand/s0),2.)
608  + 42.53*std::pow(sMand,-eta1) + 33.34*std::pow(sMand,-eta2);
609  }
610  if( neutron ) // ???
611  {
612  xsection = 35.80 + B*std::pow(G4Log(sMand/s0),2.)
613  + 40.15*std::pow(sMand,-eta1) + 30.*std::pow(sMand,-eta2);
614  }
615  fTotalXsc = xsection;
616  }
617  else if( theParticle == thePiPlus && pORn ) // pi+ /////////////////////////////////////////////
618  {
619  if( proton ) // pi+ p
620  {
621  if( pLab < 0.28 )
622  {
623  hpXsc = 10./((logP + 1.273)*(logP + 1.273) + 0.05);
624  fElasticXsc = hpXsc;
625  }
626  else if( pLab < 0.4 )
627  {
628  hpXsc = 14./( (logP + 1.273)*(logP + 1.273) + 0.07);
629  fElasticXsc = hpXsc;
630  }
631  else if( pLab < 0.68 )
632  {
633  hpXsc = 14./( (logP + 1.273)*(logP + 1.273) + 0.07);
634  fElasticXsc = hpXsc;
635  }
636  else if( pLab < 0.85 )
637  {
638  G4double Ex4 = 88*(G4Log(pLab/0.77))*(G4Log(pLab/0.77));
639  hpXsc = Ex4 + 14.9;
640  fElasticXsc = hpXsc*G4Exp(-3.*(pLab - 0.68));
641  }
642  else if( pLab < 1.15 )
643  {
644  G4double Ex4 = 88*(G4Log(pLab/0.77))*(G4Log(pLab/0.77));
645  hpXsc = Ex4 + 14.9;
646 
647  fElasticXsc = 6.0 + 1.4/(( pLab - 1.4)*( pLab - 1.4) + 0.1);
648  }
649  else if( pLab < 1.4) // ns original
650  {
651  G4double Ex1 = 3.2*G4Exp(-(pLab-2.55)*(pLab-2.55)/0.55/0.55);
652  G4double Ex2 = 12*G4Exp(-(pLab-1.47)*(pLab-1.47)/0.225/0.225);
653  hpXsc = Ex1 + Ex2 + 27.5;
654  fElasticXsc = 6.0 + 1.4/(( pLab - 1.4)*( pLab - 1.4) + 0.1);
655  }
656  else if( pLab < 2.0 ) // ns original
657  {
658  G4double Ex1 = 3.2*G4Exp(-(pLab-2.55)*(pLab-2.55)/0.55/0.55);
659  G4double Ex2 = 12*G4Exp(-(pLab-1.47)*(pLab-1.47)/0.225/0.225);
660  hpXsc = Ex1 + Ex2 + 27.5;
661  fElasticXsc = 3.0 + 1.36/( (logP - 0.336)*(logP - 0.336) + 0.08);
662  }
663  else if( pLab < 3.5 ) // ns original
664  {
665  G4double Ex1 = 3.2*G4Exp(-(pLab-2.55)*(pLab-2.55)/0.55/0.55);
666  G4double Ex2 = 12*G4Exp(-(pLab-1.47)*(pLab-1.47)/0.225/0.225);
667  hpXsc = Ex1 + Ex2 + 27.5;
668  fElasticXsc = 3.0 + 6.20/( (logP - 0.336)*(logP - 0.336) + 0.8);
669  }
670  else if( pLab < 200. ) // my
671  {
672  hpXsc = 10.6 + 2.*G4Log(pE) + 25*std::pow(pE, -0.43 ); // ns original
673  // hpXsc = GetHadronNucleonXscPDG(aParticle, nucleon )/millibarn;
674  fElasticXsc = 3.0 + 6.20/( (logP - 0.336)*(logP - 0.336) + 0.8);
675  }
676  else // pLab > 100 // my
677  {
678  hpXsc = GetHadronNucleonXscPDG(aParticle, nucleon )/millibarn;
679  fElasticXsc = 3.0 + 6.20/( (logP - 0.336)*(logP - 0.336) + 0.8);
680  }
681  fTotalXsc = hpXsc;
682  }
683  if( neutron ) // pi+ n = pi- p??
684  {
685  if( pLab < 0.28 )
686  {
687  hnXsc = 0.288/((pLab - 0.28)*(pLab - 0.28) + 0.004);
688  fElasticXsc = 1.8/((logP + 1.273)*(logP + 1.273) + 0.07);
689  }
690  else if( pLab < 0.395676 ) // first peak
691  {
692  hnXsc = 0.648/((pLab - 0.28)*(pLab - 0.28) + 0.009);
693  fElasticXsc = 0.257/((pLab - 0.28)*(pLab - 0.28) + 0.01);
694  }
695  else if( pLab < 0.5 )
696  {
697  hnXsc = 26 + 110*(G4Log(pLab/0.48))*(G4Log(pLab/0.48));
698  fElasticXsc = 0.37*hnXsc;
699  }
700  else if( pLab < 0.65 )
701  {
702  hnXsc = 26 + 110*(G4Log(pLab/0.48))*(G4Log(pLab/0.48));
703  fElasticXsc = 0.95/((pLab - 0.72)*(pLab - 0.72) + 0.049);
704  }
705  else if( pLab < 0.72 )
706  {
707  hnXsc = 36.1 + 10*G4Exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
708  24*G4Exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
709  fElasticXsc = 0.95/((pLab - 0.72)*(pLab - 0.72) + 0.049);
710  }
711  else if( pLab < 0.88 )
712  {
713  hnXsc = 36.1 + 10*G4Exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
714  24*G4Exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
715  fElasticXsc = 0.95/((pLab - 0.72)*(pLab - 0.72) + 0.049);
716  }
717  else if( pLab < 1.03 )
718  {
719  hnXsc = 36.1 + 10*G4Exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
720  24*G4Exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
721  fElasticXsc = 2.0 + 0.4/((pLab - 1.03)*(pLab - 1.03) + 0.016);
722  }
723  else if( pLab < 1.15 )
724  {
725  hnXsc = 36.1 + 10*G4Exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
726  24*G4Exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
727  fElasticXsc = 2.0 + 0.4/((pLab - 1.03)*(pLab - 1.03) + 0.016);
728  }
729  else if( pLab < 1.3 )
730  {
731  hnXsc = 36.1 + 10*G4Exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
732  24*G4Exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
733  fElasticXsc = 3. + 13./pLab;
734  }
735  else if( pLab < 2.6 ) // < 3.0) // ns original
736  {
737  hnXsc = 36.1 + 0.079-4.313*G4Log(pLab)+
738  3*G4Exp(-(pLab-2.1)*(pLab-2.1)/0.4/0.4)+
739  1.5*G4Exp(-(pLab-1.4)*(pLab-1.4)/0.12/0.12);
740  fElasticXsc = 3. + 13./pLab;
741  }
742  else if( pLab < 20. ) // < 3.0) // ns original
743  {
744  hnXsc = 36.1 + 0.079 - 4.313*G4Log(pLab)+
745  3*G4Exp(-(pLab-2.1)*(pLab-2.1)/0.4/0.4)+
746  1.5*G4Exp(-(pLab-1.4)*(pLab-1.4)/0.12/0.12);
747  fElasticXsc = 3. + 13./pLab;
748  }
749  else // mb
750  {
751  hnXsc = GetHadronNucleonXscPDG(aParticle, nucleon )/millibarn;
752  fElasticXsc = 3. + 13./pLab;
753  }
754  fTotalXsc = hnXsc;
755  }
756  }
757  else if( theParticle == thePiMinus && pORn )
758  {
759  if( neutron ) // pi- n = pi+ p??
760  {
761  if( pLab < 0.28 )
762  {
763  hnXsc = 10./((logP + 1.273)*(logP + 1.273) + 0.05);
764  fElasticXsc = hnXsc;
765  }
766  else if( pLab < 0.4 )
767  {
768  hnXsc = 14./( (logP + 1.273)*(logP + 1.273) + 0.07);
769  fElasticXsc = hnXsc;
770  }
771  else if( pLab < 0.68 )
772  {
773  hnXsc = 14./( (logP + 1.273)*(logP + 1.273) + 0.07);
774  fElasticXsc = hnXsc;
775  }
776  else if( pLab < 0.85 )
777  {
778  G4double Ex4 = 88*(G4Log(pLab/0.77))*(G4Log(pLab/0.77));
779  hnXsc = Ex4 + 14.9;
780  fElasticXsc = hnXsc*G4Exp(-3.*(pLab - 0.68));
781  }
782  else if( pLab < 1.15 )
783  {
784  G4double Ex4 = 88*(G4Log(pLab/0.77))*(G4Log(pLab/0.77));
785  hnXsc = Ex4 + 14.9;
786 
787  fElasticXsc = 6.0 + 1.4/(( pLab - 1.4)*( pLab - 1.4) + 0.1);
788  }
789  else if( pLab < 1.4) // ns original
790  {
791  G4double Ex1 = 3.2*G4Exp(-(pLab-2.55)*(pLab-2.55)/0.55/0.55);
792  G4double Ex2 = 12*G4Exp(-(pLab-1.47)*(pLab-1.47)/0.225/0.225);
793  hnXsc = Ex1 + Ex2 + 27.5;
794  fElasticXsc = 6.0 + 1.4/(( pLab - 1.4)*( pLab - 1.4) + 0.1);
795  }
796  else if( pLab < 2.0 ) // ns original
797  {
798  G4double Ex1 = 3.2*G4Exp(-(pLab-2.55)*(pLab-2.55)/0.55/0.55);
799  G4double Ex2 = 12*G4Exp(-(pLab-1.47)*(pLab-1.47)/0.225/0.225);
800  hnXsc = Ex1 + Ex2 + 27.5;
801  fElasticXsc = 3.0 + 1.36/( (logP - 0.336)*(logP - 0.336) + 0.08);
802  }
803  else if( pLab < 3.5 ) // ns original
804  {
805  G4double Ex1 = 3.2*G4Exp(-(pLab-2.55)*(pLab-2.55)/0.55/0.55);
806  G4double Ex2 = 12*G4Exp(-(pLab-1.47)*(pLab-1.47)/0.225/0.225);
807  hnXsc = Ex1 + Ex2 + 27.5;
808  fElasticXsc = 3.0 + 6.20/( (logP - 0.336)*(logP - 0.336) + 0.8);
809  }
810  else if( pLab < 200. ) // my
811  {
812  hnXsc = 10.6 + 2.*G4Log(pE) + 25*std::pow(pE, -0.43 ); // ns original
813  fElasticXsc = 3.0 + 6.20/( (logP - 0.336)*(logP - 0.336) + 0.8);
814  }
815  else // pLab > 100 // my
816  {
817  hnXsc = GetHadronNucleonXscPDG(aParticle, nucleon )/millibarn;
818  fElasticXsc = 3.0 + 6.20/( (logP - 0.336)*(logP - 0.336) + 0.8);
819  }
820  fTotalXsc = hnXsc;
821  }
822  if( proton ) // pi- p
823  {
824  if( pLab < 0.28 )
825  {
826  hpXsc = 0.288/((pLab - 0.28)*(pLab - 0.28) + 0.004);
827  fElasticXsc = 1.8/((logP + 1.273)*(logP + 1.273) + 0.07);
828  }
829  else if( pLab < 0.395676 ) // first peak
830  {
831  hpXsc = 0.648/((pLab - 0.28)*(pLab - 0.28) + 0.009);
832  fElasticXsc = 0.257/((pLab - 0.28)*(pLab - 0.28) + 0.01);
833  }
834  else if( pLab < 0.5 )
835  {
836  hpXsc = 26 + 110*(G4Log(pLab/0.48))*(G4Log(pLab/0.48));
837  fElasticXsc = 0.37*hpXsc;
838  }
839  else if( pLab < 0.65 )
840  {
841  hpXsc = 26 + 110*(G4Log(pLab/0.48))*(G4Log(pLab/0.48));
842  fElasticXsc = 0.95/((pLab - 0.72)*(pLab - 0.72) + 0.049);
843  }
844  else if( pLab < 0.72 )
845  {
846  hpXsc = 36.1+
847  10*G4Exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
848  24*G4Exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
849  fElasticXsc = 0.95/((pLab - 0.72)*(pLab - 0.72) + 0.049);
850  }
851  else if( pLab < 0.88 )
852  {
853  hpXsc = 36.1+
854  10*G4Exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
855  24*G4Exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
856  fElasticXsc = 0.95/((pLab - 0.72)*(pLab - 0.72) + 0.049);
857  }
858  else if( pLab < 1.03 )
859  {
860  hpXsc = 36.1+
861  10*G4Exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
862  24*G4Exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
863  fElasticXsc = 2.0 + 0.4/((pLab - 1.03)*(pLab - 1.03) + 0.016);
864  }
865  else if( pLab < 1.15 )
866  {
867  hpXsc = 36.1+
868  10*G4Exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
869  24*G4Exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
870  fElasticXsc = 2.0 + 0.4/((pLab - 1.03)*(pLab - 1.03) + 0.016);
871  }
872  else if( pLab < 1.3 )
873  {
874  hpXsc = 36.1+
875  10*G4Exp(-(pLab-0.72)*(pLab-0.72)/0.06/0.06)+
876  24*G4Exp(-(pLab-1.015)*(pLab-1.015)/0.075/0.075);
877  fElasticXsc = 3. + 13./pLab;
878  }
879  else if( pLab < 2.6 ) // < 3.0) // ns original
880  {
881  hpXsc = 36.1+0.079-4.313*G4Log(pLab)+
882  3*G4Exp(-(pLab-2.1)*(pLab-2.1)/0.4/0.4)+
883  1.5*G4Exp(-(pLab-1.4)*(pLab-1.4)/0.12/0.12);
884  fElasticXsc = 3. +13./pLab; // *G4Log(pLab*6.79);
885  }
886  else // mb
887  {
888  hpXsc = GetHadronNucleonXscPDG(aParticle, nucleon )/millibarn;
889  fElasticXsc = 3. + 13./pLab;
890  }
891  fTotalXsc = hpXsc;
892  }
893  }
894  else if( (theParticle == theKMinus || theParticle == theK0S) && proton ) // Kmp/K0p /////////////////////////////////
895  {
896  if( pLab < pMin)
897  {
898  G4double psp = pLab*std::sqrt(pLab);
899  fElasticXsc = 5.2/psp;
900  fTotalXsc = 14./psp;
901  }
902  else if( pLab > pMax )
903  {
904  G4double ld = G4Log(pLab) - minLogP;
905  G4double ld2 = ld*ld;
906  fElasticXsc = cofLogE*ld2 + 2.23;
907  fTotalXsc = 1.1*cofLogT*ld2 + 19.7;
908  }
909  else
910  {
911  G4double ld = G4Log(pLab) - minLogP;
912  G4double ld2 = ld*ld;
913  G4double sp = std::sqrt(pLab);
914  G4double psp = pLab*sp;
915  G4double p2 = pLab*pLab;
916  G4double p4 = p2*p2;
917  G4double lm = pLab - .39;
918  G4double md = lm*lm + .000356;
919 
920  G4double lh1 = pLab - 0.78;
921  G4double hd1 = lh1*lh1 + .00166;
922 
923  G4double lh = pLab - 1.01;
924  G4double hd = lh*lh + .011;
925 
926  G4double lh2 = pLab - 1.63;
927  G4double hd2 = lh2*lh2 + .007;
928 
929  fElasticXsc = 5.2/psp + (1.1*cofLogE*ld2 + 2.23)/(1. - .7/sp + .075/p4)
930  + .004/md + 0.005/hd1+ 0.01/hd2 +.15/hd; // small peaks were added
931 
932  fTotalXsc = 14./psp + (1.1*cofLogT*ld2 + 19.5)/(1. - .21/sp + .52/p4)
933  + .006/md + 0.01/hd1+ 0.02/hd2 + .20/hd ;
934  }
935  }
936  else if( (theParticle == theKMinus || theParticle == theK0S) && neutron ) // Kmn/K0n //////////////////////////////
937  {
938  if( pLab > pMax )
939  {
940  G4double ld = G4Log(pLab) - minLogP;
941  G4double ld2 = ld*ld;
942  fElasticXsc = cofLogE*ld2 + 2.23;
943  fTotalXsc = 1.1*cofLogT*ld2 + 19.7;
944  }
945  else
946  {
947 
948  G4double lh = pLab - 0.98;
949  G4double hd = lh*lh + .021;
950 
951  G4double LogPlab = G4Log( pLab );
952  G4double sqrLogPlab = LogPlab * LogPlab;
953 
954  fElasticXsc = // 5.2/psp + (cofLogE*ld2 + 2.23)/(1. - .7/sp + .075/p4) + .004/md
955  5.0 + 8.1*std::pow(pLab,-1.8 ) + 0.16*sqrLogPlab - 1.3*LogPlab + .15/hd;
956  fTotalXsc = // 14./psp +
957  // (1.1*cofLogT*ld2 + 19.5)/(1. - .21/sp + .52/p4)
958  // WP 25.2 + 0. *std::pow(pLab, 0. ) + 0.38*sqrLogPlab - 2.9*LogPlab
959  25.2 + 0.38*sqrLogPlab - 2.9*LogPlab
960  // + .006/md + 0.01/hd1+ 0.02/hd2
961  + 0.30/hd ;
962  }
963  }
964  else if( (theParticle == theKPlus || theParticle == theK0L) && proton ) // Kpp/aKp ////////////////////////
965  {
966  if( pLab < pMin )
967  {
968  G4double lr = pLab - .38;
969  G4double lm = pLab - 1.;
970  G4double md = lm*lm + .392;
971  fElasticXsc = .7/(lr*lr + .076) + 2./md;
972  fTotalXsc = .7/(lr*lr + .076) + 2.6/md;
973  }
974  else if( pLab > pMax )
975  {
976  G4double ld = G4Log(pLab) - minLogP;
977  G4double ld2 = ld*ld;
978  fElasticXsc = cofLogE*ld2 + 2.23;
979  fTotalXsc = cofLogT*ld2 + 19.2;
980  }
981  else
982  {
983  G4double ld = G4Log(pLab) - minLogP;
984  G4double ld2 = ld*ld;
985  G4double lr = pLab - .38;
986  G4double LE = .7/(lr*lr + .076);
987  G4double sp = std::sqrt(pLab);
988  G4double p2 = pLab*pLab;
989  G4double p4 = p2*p2;
990  G4double lm = pLab - 1.;
991  G4double md = lm*lm + .392;
992  fElasticXsc = LE + (cofLogE*ld2 + 2.23)/(1. - .7/sp + .1/p4) + 2./md;
993  fTotalXsc = LE + (cofLogT*ld2 + 19.5)/(1. + .46/sp + 1.6/p4) + 2.6/md;
994  }
995  }
996  else if( (theParticle == theKPlus || theParticle == theK0L) && neutron ) // Kpn/aKn ///////////////////////
997  {
998  if( pLab < pMin )
999  {
1000  G4double lm = pLab - 0.94;
1001  G4double md = lm*lm + .392;
1002  fElasticXsc = 2./md;
1003  fTotalXsc = 4.6/md;
1004  }
1005  else if( pLab > pMax )
1006  {
1007  G4double ld = G4Log(pLab) - minLogP;
1008  G4double ld2 = ld*ld;
1009  fElasticXsc = cofLogE*ld2 + 2.23;
1010  fTotalXsc = cofLogT*ld2 + 19.2;
1011  }
1012  else
1013  {
1014  G4double ld = G4Log(pLab) - minLogP;
1015  G4double ld2 = ld*ld;
1016  G4double sp = std::sqrt(pLab);
1017  G4double p2 = pLab*pLab;
1018  G4double p4 = p2*p2;
1019  G4double lm = pLab - 0.94;
1020  G4double md = lm*lm + .392;
1021  fElasticXsc = (cofLogE*ld2 + 2.23)/(1. - .7/sp + .1/p4) + 2./md;
1022  fTotalXsc = (cofLogT*ld2 + 19.5)/(1. + .46/sp + 1.6/p4) + 4.6/md;
1023  }
1024  }
1025  else if( theParticle == theSMinus && pORn )
1026  {
1027  xsection = 35.20 + B*std::pow(G4Log(sMand/s0),2.)
1028  - 199.*std::pow(sMand,-eta1) + 264.*std::pow(sMand,-eta2);
1029  }
1030  else if( theParticle == theGamma && pORn ) // modify later on
1031  {
1032  xsection = 0.0 + B*std::pow(G4Log(sMand/s0),2.)
1033  + 0.032*std::pow(sMand,-eta1); // WP - 0.0*std::pow(sMand,-eta2);
1034  fTotalXsc = xsection;
1035  }
1036  else // other then p,n,pi+,pi-,K+,K- as proton ???
1037  {
1038  if( proton )
1039  {
1040  xsection = 35.45 + B*std::pow(G4Log(sMand/s0),2.)
1041  + 42.53*std::pow(sMand,-eta1) - 33.34*std::pow(sMand,-eta2);
1042  }
1043  if( neutron )
1044  {
1045  xsection += 35.80 + B*std::pow(G4Log(sMand/s0),2.)
1046  + 40.15*std::pow(sMand,-eta1) - 30.*std::pow(sMand,-eta2);
1047  }
1048  fTotalXsc = xsection;
1049  }
1050  fTotalXsc *= millibarn; // parametrised in mb
1051  fElasticXsc *= millibarn; // parametrised in mb
1052 
1053  if( proton && aParticle->GetDefinition()->GetPDGCharge() > 0. )
1054  {
1055  G4double cB = GetCoulombBarrier(aParticle, nucleon);
1056  fTotalXsc *= cB;
1057  fElasticXsc *= cB;
1058  }
1059  fInelasticXsc = fTotalXsc - fElasticXsc;
1060  if( fInelasticXsc < 0. ) fInelasticXsc = 0.;
1061 
1062  // G4cout<<fTotalXsc/millibarn<<"; "<<fElasticXsc/millibarn<<"; "<<fInelasticXsc/millibarn<<G4endl;
1063 
1064  return fTotalXsc;
1065 }
LE
Definition: Evaluator.cc:66
G4HadronNucleonXsc::GetCoulombBarrier
G4double GetCoulombBarrier(const G4DynamicParticle *aParticle, const G4ParticleDefinition *nucleon)
Definition: G4HadronNucleonXsc.cc:1433
G4DynamicParticle::GetTotalEnergy
G4double GetTotalEnergy() const
neutron
Definition: G4DataQuestionaire.hh:35
G4DynamicParticle::GetDefinition
G4ParticleDefinition * GetDefinition() const
B
double B(double temperature)
Definition: G4DNAElectronHoleRecombination.cc:69
G4DynamicParticle::GetMass
G4double GetMass() const
G4bool
bool G4bool
Definition: G4Types.hh:79
G4Log
G4double G4Log(G4double x)
Definition: G4Log.hh:230
G4Exp
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
GeV
static constexpr double GeV
Definition: G4SIunits.hh:217
G4HadronNucleonXsc::GetHadronNucleonXscPDG
G4double GetHadronNucleonXscPDG(const G4DynamicParticle *, const G4ParticleDefinition *)
Definition: G4HadronNucleonXsc.cc:223
G4double
double G4double
Definition: G4Types.hh:76
G4ParticleDefinition::GetPDGCharge
G4double GetPDGCharge() const
Definition: G4ParticleDefinition.hh:124
CLHEP::Hep3Vector::mag
double mag() const
millibarn
static constexpr double millibarn
Definition: G4SIunits.hh:106
G4DynamicParticle::GetMomentum
G4ThreeVector GetMomentum() const
G4HadronNucleonXsc::CalcMandelstamS
G4double CalcMandelstamS(const G4double, const G4double, const G4double)
Definition: G4HadronNucleonXsc.cc:1418

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4HadronNucleonXsc::GetHadronNucleonXscPDG ( const G4DynamicParticle *  aParticle,
const G4ParticleDefinition *  nucleon 
)

Definition at line 223 of file G4HadronNucleonXsc.cc.

225 {
226  G4double xsection(0);
227  G4int Zt=1, Nt=1, At=1;
228 
229  G4double targ_mass = 0.939*GeV; // ~mean neutron and proton ???
230 
231  G4double proj_mass = aParticle->GetMass();
232  G4double proj_momentum = aParticle->GetMomentum().mag();
233 
234  G4double sMand = CalcMandelstamS ( proj_mass , targ_mass , proj_momentum );
235 
236  sMand /= GeV*GeV; // in GeV for parametrisation
237 
238  // General PDG fit constants
239 
240  G4double s0 = 5.38*5.38; // in Gev^2
241  G4double eta1 = 0.458;
242  G4double eta2 = 0.458;
243  G4double B = 0.308;
244 
245  const G4ParticleDefinition* theParticle = aParticle->GetDefinition();
246 
247  G4bool pORn = (nucleon == theProton || nucleon == theNeutron );
248  G4bool proton = (nucleon == theProton);
249  G4bool neutron = (nucleon == theNeutron);
250 
251  if(theParticle == theNeutron) // proton-neutron fit
252  {
253  if ( proton )
254  {
255  xsection = Zt*( 35.80 + B*std::pow(G4Log(sMand/s0),2.)
256  + 40.15*std::pow(sMand,-eta1) - 30.*std::pow(sMand,-eta2));// on p
257  }
258  if ( neutron )
259  {
260  xsection = Nt*( 35.45 + B*std::pow(G4Log(sMand/s0),2.)
261  + 42.53*std::pow(sMand,-eta1) - 33.34*std::pow(sMand,-eta2)); // on n pp for nn
262  }
263  }
264  else if(theParticle == theProton)
265  {
266  if ( proton )
267  {
268  xsection = Zt*( 35.45 + B*std::pow(G4Log(sMand/s0),2.)
269  + 42.53*std::pow(sMand,-eta1) - 33.34*std::pow(sMand,-eta2));
270  }
271  if ( neutron )
272  {
273  xsection = Nt*( 35.80 + B*std::pow(G4Log(sMand/s0),2.)
274  + 40.15*std::pow(sMand,-eta1) - 30.*std::pow(sMand,-eta2));
275  }
276  }
277  else if(theParticle == theAProton)
278  {
279  if ( proton )
280  {
281  xsection = Zt*( 35.45 + B*std::pow(G4Log(sMand/s0),2.)
282  + 42.53*std::pow(sMand,-eta1) + 33.34*std::pow(sMand,-eta2));
283  }
284  if ( neutron )
285  {
286  xsection = Nt*( 35.80 + B*std::pow(G4Log(sMand/s0),2.)
287  + 40.15*std::pow(sMand,-eta1) + 30.*std::pow(sMand,-eta2));
288  }
289  }
290  else if(theParticle == thePiPlus && pORn )
291  {
292  xsection = At*( 20.86 + B*std::pow(G4Log(sMand/s0),2.)
293  + 19.24*std::pow(sMand,-eta1) - 6.03*std::pow(sMand,-eta2));
294  }
295  else if(theParticle == thePiMinus && pORn )
296  {
297  xsection = At*( 20.86 + B*std::pow(G4Log(sMand/s0),2.)
298  + 19.24*std::pow(sMand,-eta1) + 6.03*std::pow(sMand,-eta2));
299  }
300  else if(theParticle == theKPlus)
301  {
302  if ( proton )
303  {
304  xsection = Zt*( 17.91 + B*std::pow(G4Log(sMand/s0),2.)
305  + 7.14*std::pow(sMand,-eta1) - 13.45*std::pow(sMand,-eta2));
306  }
307  if ( neutron )
308  {
309  xsection = Nt*( 17.87 + B*std::pow(G4Log(sMand/s0),2.)
310  + 5.17*std::pow(sMand,-eta1) - 7.23*std::pow(sMand,-eta2));
311  }
312  }
313  else if(theParticle == theKMinus)
314  {
315  if ( proton )
316  {
317  xsection = Zt*( 17.91 + B*std::pow(G4Log(sMand/s0),2.)
318  + 7.14*std::pow(sMand,-eta1) + 13.45*std::pow(sMand,-eta2));
319  }
320  if ( neutron )
321  {
322  xsection = Nt*( 17.87 + B*std::pow(G4Log(sMand/s0),2.)
323  + 5.17*std::pow(sMand,-eta1) + 7.23*std::pow(sMand,-eta2) );
324  }
325  }
326  else if(theParticle == theSMinus && pORn )
327  {
328  xsection = At*( 35.20 + B*std::pow(G4Log(sMand/s0),2.)
329  - 199.*std::pow(sMand,-eta1) + 264.*std::pow(sMand,-eta2) );
330  }
331  else if(theParticle == theGamma && pORn ) // modify later on
332  {
333  xsection = At*( 0.0 + B*std::pow(G4Log(sMand/s0),2.)
334  + 0.032*std::pow(sMand,-eta1) - 0.0*std::pow(sMand,-eta2) );
335 
336  }
337  else // as proton ???
338  {
339  if ( proton )
340  {
341  xsection = Zt*( 35.45 + B*std::pow(G4Log(sMand/s0),2.)
342  + 42.53*std::pow(sMand,-eta1) - 33.34*std::pow(sMand,-eta2) );
343  }
344  if ( neutron )
345  {
346  xsection = Nt*( 35.80 + B*std::pow(G4Log(sMand/s0),2.)
347  + 40.15*std::pow(sMand,-eta1) - 30.*std::pow(sMand,-eta2));
348  }
349  }
350  xsection *= millibarn; // parametrised in mb
351 
352  fTotalXsc = xsection;
353  fInelasticXsc = 0.75*xsection;
354  fElasticXsc = fTotalXsc - fInelasticXsc;
355  if (fElasticXsc < 0.) fElasticXsc = 0.;
356 
357  return xsection;
358 }
neutron
Definition: G4DataQuestionaire.hh:35
G4DynamicParticle::GetDefinition
G4ParticleDefinition * GetDefinition() const
B
double B(double temperature)
Definition: G4DNAElectronHoleRecombination.cc:69
G4int
int G4int
Definition: G4Types.hh:78
G4DynamicParticle::GetMass
G4double GetMass() const
G4bool
bool G4bool
Definition: G4Types.hh:79
G4Log
G4double G4Log(G4double x)
Definition: G4Log.hh:230
GeV
static constexpr double GeV
Definition: G4SIunits.hh:217
G4double
double G4double
Definition: G4Types.hh:76
CLHEP::Hep3Vector::mag
double mag() const
millibarn
static constexpr double millibarn
Definition: G4SIunits.hh:106
G4DynamicParticle::GetMomentum
G4ThreeVector GetMomentum() const
G4HadronNucleonXsc::CalcMandelstamS
G4double CalcMandelstamS(const G4double, const G4double, const G4double)
Definition: G4HadronNucleonXsc.cc:1418

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4HadronNucleonXsc::GetHadronNucleonXscVU ( const G4DynamicParticle *  aParticle,
const G4ParticleDefinition *  nucleon 
)

Definition at line 1238 of file G4HadronNucleonXsc.cc.

1240 {
1241  G4int PDGcode = aParticle->GetDefinition()->GetPDGEncoding();
1242  G4int absPDGcode = std::abs(PDGcode);
1243  G4double Elab = aParticle->GetTotalEnergy();
1244  // (s - 2*0.88*GeV*GeV)/(2*0.939*GeV)/GeV;
1245  G4double Plab = aParticle->GetMomentum().mag();
1246  // std::sqrt(Elab * Elab - 0.88);
1247 
1248  Elab /= GeV;
1249  Plab /= GeV;
1250 
1251  G4double LogPlab = G4Log( Plab );
1252  G4double sqrLogPlab = LogPlab * LogPlab;
1253 
1254  G4bool pORn = (nucleon == theProton || nucleon == theNeutron );
1255  G4bool proton = (nucleon == theProton);
1256  G4bool neutron = (nucleon == theNeutron);
1257 
1258 
1259  if( absPDGcode > 1000 && pORn ) //------Projectile is baryon -
1260  {
1261  if(proton)
1262  {
1263  // WP fTotalXsc = 48.0 + 0. *std::pow(Plab, 0. ) + 0.522*sqrLogPlab - 4.51*LogPlab;
1264  fTotalXsc = 48.0 + 0.522*sqrLogPlab - 4.51*LogPlab;
1265  fElasticXsc = 11.9 + 26.9*std::pow(Plab,-1.21) + 0.169*sqrLogPlab - 1.85*LogPlab;
1266  }
1267  if(neutron)
1268  {
1269  // WP fTotalXsc = 47.3 + 0. *std::pow(Plab, 0. ) + 0.513*sqrLogPlab - 4.27*LogPlab;
1270  fTotalXsc = 47.3 + 0.513*sqrLogPlab - 4.27*LogPlab;
1271  fElasticXsc = 11.9 + 26.9*std::pow(Plab,-1.21) + 0.169*sqrLogPlab - 1.85*LogPlab;
1272  }
1273  }
1274  else if( PDGcode == 211 && pORn ) //------Projectile is PionPlus ----
1275  {
1276  if(proton)
1277  {
1278  fTotalXsc = 16.4 + 19.3 *std::pow(Plab,-0.42) + 0.19 *sqrLogPlab - 0.0 *LogPlab;
1279  fElasticXsc = 0.0 + 11.4*std::pow(Plab,-0.40) + 0.079*sqrLogPlab - 0.0 *LogPlab;
1280  }
1281  if(neutron)
1282  {
1283  fTotalXsc = 33.0 + 14.0 *std::pow(Plab,-1.36) + 0.456*sqrLogPlab - 4.03*LogPlab;
1284  fElasticXsc = 1.76 + 11.2*std::pow(Plab,-0.64) + 0.043*sqrLogPlab - 0.0 *LogPlab;
1285  }
1286  }
1287  else if( PDGcode == -211 && pORn ) //------Projectile is PionMinus ----
1288  {
1289  if(proton)
1290  {
1291  fTotalXsc = 33.0 + 14.0 *std::pow(Plab,-1.36) + 0.456*sqrLogPlab - 4.03*LogPlab;
1292  fElasticXsc = 1.76 + 11.2*std::pow(Plab,-0.64) + 0.043*sqrLogPlab - 0.0 *LogPlab;
1293  }
1294  if(neutron)
1295  {
1296  fTotalXsc = 16.4 + 19.3 *std::pow(Plab,-0.42) + 0.19 *sqrLogPlab - 0.0 *LogPlab;
1297  fElasticXsc = 0.0 + 11.4*std::pow(Plab,-0.40) + 0.079*sqrLogPlab - 0.0 *LogPlab;
1298  }
1299  }
1300  else if( PDGcode == 111 && pORn ) //------Projectile is PionZero --
1301  {
1302  if(proton)
1303  {
1304  fTotalXsc = (16.4 + 19.3 *std::pow(Plab,-0.42) + 0.19 *sqrLogPlab - 0.0 *LogPlab + //Pi+
1305  33.0 + 14.0 *std::pow(Plab,-1.36) + 0.456*sqrLogPlab - 4.03*LogPlab)/2; //Pi-
1306 
1307  fElasticXsc = ( 0.0 + 11.4*std::pow(Plab,-0.40) + 0.079*sqrLogPlab - 0.0 *LogPlab + //Pi+
1308  1.76 + 11.2*std::pow(Plab,-0.64) + 0.043*sqrLogPlab - 0.0 *LogPlab)/2; //Pi-
1309 
1310  }
1311  if(neutron)
1312  {
1313  fTotalXsc = (33.0 + 14.0 *std::pow(Plab,-1.36) + 0.456*sqrLogPlab - 4.03*LogPlab + //Pi+
1314  16.4 + 19.3 *std::pow(Plab,-0.42) + 0.19 *sqrLogPlab - 0.0 *LogPlab)/2; //Pi-
1315  fElasticXsc = ( 1.76 + 11.2*std::pow(Plab,-0.64) + 0.043*sqrLogPlab - 0.0 *LogPlab + //Pi+
1316  0.0 + 11.4*std::pow(Plab,-0.40) + 0.079*sqrLogPlab - 0.0 *LogPlab)/2; //Pi-
1317  }
1318  }
1319  else if( PDGcode == 321 && pORn ) //------Projectile is KaonPlus --
1320  {
1321  if(proton)
1322  {
1323  // WP fTotalXsc = 18.1 + 0. *std::pow(Plab, 0. ) + 0.26 *sqrLogPlab - 1.0 *LogPlab;
1324  fTotalXsc = 18.1 + 0.26 *sqrLogPlab - 1.0 *LogPlab;
1325  fElasticXsc = 5.0 + 8.1*std::pow(Plab,-1.8 ) + 0.16 *sqrLogPlab - 1.3 *LogPlab;
1326  }
1327  if(neutron)
1328  {
1329  // WP fTotalXsc = 18.7 + 0. *std::pow(Plab, 0. ) + 0.21 *sqrLogPlab - 0.89*LogPlab;
1330  // WP fElasticXsc = 7.3 + 0. *std::pow(Plab,-0. ) + 0.29 *sqrLogPlab - 2.4 *LogPlab;
1331  fTotalXsc = 18.7 + 0.21 *sqrLogPlab - 0.89*LogPlab;
1332  fElasticXsc = 7.3 + 0.29 *sqrLogPlab - 2.4 *LogPlab;
1333  }
1334  }
1335  else if( PDGcode ==-321 && pORn ) //------Projectile is KaonMinus ----
1336  {
1337  if(proton)
1338  {
1339  // WP fTotalXsc = 32.1 + 0. *std::pow(Plab, 0. ) + 0.66*sqrLogPlab - 5.6*LogPlab;
1340  // WP fElasticXsc = 7.3 + 0. *std::pow(Plab,-0. ) + 0.29*sqrLogPlab - 2.4*LogPlab;
1341  fTotalXsc = 32.1 + 0.66*sqrLogPlab - 5.6*LogPlab;
1342  fElasticXsc = 7.3 + 0.29*sqrLogPlab - 2.4*LogPlab;
1343  }
1344  if(neutron)
1345  {
1346  // WP fTotalXsc = 25.2 + 0. *std::pow(Plab, 0. ) + 0.38*sqrLogPlab - 2.9*LogPlab;
1347  fTotalXsc = 25.2 + 0.38*sqrLogPlab - 2.9*LogPlab;
1348  fElasticXsc = 5.0 + 8.1*std::pow(Plab,-1.8 ) + 0.16*sqrLogPlab - 1.3*LogPlab;
1349  }
1350  }
1351  else if( PDGcode == 311 && pORn ) //------Projectile is KaonZero -----
1352  {
1353  if(proton)
1354  {
1355  // WP fTotalXsc = ( 18.1 + 0. *std::pow(Plab, 0. ) + 0.26 *sqrLogPlab - 1.0 *LogPlab + //K+
1356  // WP 32.1 + 0. *std::pow(Plab, 0. ) + 0.66 *sqrLogPlab - 5.6 *LogPlab)/2; //K-
1357  fTotalXsc = ( 18.1 + 0.26 *sqrLogPlab - 1.0 *LogPlab + //K+
1358  32.1 + 0.66 *sqrLogPlab - 5.6 *LogPlab)/2; //K-
1359  fElasticXsc = ( 5.0 + 8.1*std::pow(Plab,-1.8 ) + 0.16 *sqrLogPlab - 1.3 *LogPlab + //K+
1360  7.3 + 0.29 *sqrLogPlab - 2.4 *LogPlab)/2; //K-
1361  // WP 7.3 + 0. *std::pow(Plab,-0. ) + 0.29 *sqrLogPlab - 2.4 *LogPlab)/2; //K-
1362  }
1363  if(neutron)
1364  {
1365  // WP fTotalXsc = ( 18.7 + 0. *std::pow(Plab, 0. ) + 0.21 *sqrLogPlab - 0.89*LogPlab + //K+
1366  // WP 25.2 + 0. *std::pow(Plab, 0. ) + 0.38 *sqrLogPlab - 2.9 *LogPlab)/2; //K-
1367  fTotalXsc = ( 18.7 + 0.21 *sqrLogPlab - 0.89*LogPlab + //K+
1368  25.2 + 0.38 *sqrLogPlab - 2.9 *LogPlab)/2; //K-
1369  // WP fElasticXsc = ( 7.3 + 0. *std::pow(Plab,-0. ) + 0.29 *sqrLogPlab - 2.4 *LogPlab + //K+
1370  fElasticXsc = ( 7.3 + 0.29 *sqrLogPlab - 2.4 *LogPlab + //K+
1371  5.0 + 8.1*std::pow(Plab,-1.8 ) + 0.16 *sqrLogPlab - 1.3 *LogPlab)/2; //K-
1372  }
1373  }
1374  else //------Projectile is undefined, Nucleon assumed
1375  {
1376  if(proton)
1377  {
1378  // WP fTotalXsc = 48.0 + 0. *std::pow(Plab, 0. ) + 0.522*sqrLogPlab - 4.51*LogPlab;
1379  fTotalXsc = 48.0 + 0.522*sqrLogPlab - 4.51*LogPlab;
1380  fElasticXsc = 11.9 + 26.9*std::pow(Plab,-1.21) + 0.169*sqrLogPlab - 1.85*LogPlab;
1381  }
1382  if(neutron)
1383  {
1384  // WP fTotalXsc = 47.3 + 0. *std::pow(Plab, 0. ) + 0.513*sqrLogPlab - 4.27*LogPlab;
1385  fTotalXsc = 47.3 + 0.513*sqrLogPlab - 4.27*LogPlab;
1386  fElasticXsc = 11.9 + 26.9*std::pow(Plab,-1.21) + 0.169*sqrLogPlab - 1.85*LogPlab;
1387  }
1388  }
1389  fTotalXsc *= millibarn;
1390  fElasticXsc *= millibarn;
1391  fInelasticXsc = fTotalXsc - fElasticXsc;
1392  if (fInelasticXsc < 0.) fInelasticXsc = 0.;
1393 
1394  return fTotalXsc;
1395 }
G4DynamicParticle::GetTotalEnergy
G4double GetTotalEnergy() const
neutron
Definition: G4DataQuestionaire.hh:35
G4DynamicParticle::GetDefinition
G4ParticleDefinition * GetDefinition() const
G4int
int G4int
Definition: G4Types.hh:78
G4bool
bool G4bool
Definition: G4Types.hh:79
G4Log
G4double G4Log(G4double x)
Definition: G4Log.hh:230
GeV
static constexpr double GeV
Definition: G4SIunits.hh:217
G4double
double G4double
Definition: G4Types.hh:76
CLHEP::Hep3Vector::mag
double mag() const
millibarn
static constexpr double millibarn
Definition: G4SIunits.hh:106
G4DynamicParticle::GetMomentum
G4ThreeVector GetMomentum() const

Here is the call graph for this function:

G4double G4HadronNucleonXsc::GetInelasticHadronNucleonXsc ( )
inline

Definition at line 92 of file G4HadronNucleonXsc.hh.

92 { return fInelasticXsc; };

Here is the caller graph for this function:

G4double G4HadronNucleonXsc::GetKaonNucleonXscGG ( const G4DynamicParticle *  aParticle,
const G4ParticleDefinition *  nucleon 
)

Definition at line 1072 of file G4HadronNucleonXsc.cc.

1074 {
1075  G4double pLab = aParticle->GetMomentum().mag();
1076 
1077  pLab /= GeV;
1078  G4double LogPlab = G4Log( pLab );
1079  G4double sqrLogPlab = LogPlab * LogPlab;
1080 
1081  G4double minLogP = 3.5; // min of (lnP-minLogP)^2
1082  G4double cofLogE = .0557; // elastic (lnP-minLogP)^2
1083  G4double cofLogT = .3; // total (lnP-minLogP)^2
1084  G4double pMin = .1; // fast LE calculation
1085  G4double pMax = 1000.; // fast HE calculation
1086 
1087  const G4ParticleDefinition* theParticle = aParticle->GetDefinition();
1088 
1089  G4bool proton = (nucleon == theProton);
1090  G4bool neutron = (nucleon == theNeutron);
1091 
1092  if( (theParticle == theKMinus || theParticle == theK0S) && proton ) // (K-,K0)on p ////////////////////////////
1093  {
1094 
1095  if( pLab < pMin)
1096  {
1097  G4double psp = pLab*std::sqrt(pLab);
1098  fElasticXsc = 5.2/psp;
1099  fTotalXsc = 14./psp;
1100  }
1101  else if( pLab > pMax )
1102  {
1103  G4double ld = G4Log(pLab) - minLogP;
1104  G4double ld2 = ld*ld;
1105  fElasticXsc = cofLogE*ld2 + 2.23;
1106  fTotalXsc = 1.1*cofLogT*ld2 + 19.7;
1107  }
1108  else
1109  {
1110  G4double ld = G4Log(pLab) - minLogP;
1111  G4double ld2 = ld*ld;
1112  G4double sp = std::sqrt(pLab);
1113  G4double psp = pLab*sp;
1114  G4double p2 = pLab*pLab;
1115  G4double p4 = p2*p2;
1116 
1117  G4double lh = pLab - 0.98;
1118  G4double hd = lh*lh + .045;
1119 
1120 
1121  fElasticXsc = 5.2/psp + (cofLogE*ld2 + 2.23)/(1. - .7/sp + .075/p4) // + .004/md
1122  + .15/hd;
1123  fTotalXsc = 14./psp + (1.1*cofLogT*ld2 + 19.5)/(1. - .21/sp + .52/p4)
1124  // + .006/md + 0.01/hd1 + 0.02/hd2
1125  + .60/hd;
1126  }
1127  }
1128  else if( (theParticle == theKMinus || theParticle == theK0S) && neutron ) // Kmn/K0n /////////////////////////////
1129  {
1130  if( pLab > pMax )
1131  {
1132  G4double ld = G4Log(pLab) - minLogP;
1133  G4double ld2 = ld*ld;
1134  fElasticXsc = cofLogE*ld2 + 2.23;
1135  fTotalXsc = 1.1*cofLogT*ld2 + 19.7;
1136  }
1137  else
1138  {
1139 
1140  G4double lh = pLab - 0.98;
1141  G4double hd = lh*lh + .045;
1142 
1143  fElasticXsc = // 5.2/psp + (cofLogE*ld2 + 2.23)/(1. - .7/sp + .075/p4) + .004/md
1144  5.0 + 8.1*std::pow(pLab,-1.8 ) + 0.16*sqrLogPlab - 1.3*LogPlab + .15/hd;
1145  fTotalXsc = // 14./psp +
1146  // (1.1*cofLogT*ld2 + 19.5)/(1. - .21/sp + .52/p4)
1147  // WP 25.2 + 0. *std::pow(pLab, 0. ) + 0.38*sqrLogPlab - 2.9*LogPlab
1148  25.2 + 0.38*sqrLogPlab - 2.9*LogPlab
1149  // + .006/md + 0.01/hd1+ 0.02/hd2
1150  + 0.60/hd ;
1151  }
1152  }
1153  else if( (theParticle == theKPlus || theParticle == theK0L) && proton ) // Kpp/aKp //////////////////////
1154  {
1155  if( pLab < pMin )
1156  {
1157  G4double lr = pLab - .38;
1158  G4double lm = pLab - 1.;
1159  G4double md = lm*lm + .392;
1160  fElasticXsc = .7/(lr*lr + .076) + 2./md;
1161  fTotalXsc = // .7/(lr*lr + .076) +
1162  2.6/md;
1163  }
1164  else if( pLab > pMax )
1165  {
1166  G4double ld = G4Log(pLab) - minLogP;
1167  G4double ld2 = ld*ld;
1168  fElasticXsc = cofLogE*ld2 + 2.23;
1169  fTotalXsc = cofLogT*ld2 + 19.2;
1170  }
1171  else
1172  {
1173  G4double ld = G4Log(pLab) - minLogP;
1174  G4double ld2 = ld*ld;
1175  G4double lr = pLab - .38;
1176  G4double LE = .7/(lr*lr + .076);
1177  G4double sp = std::sqrt(pLab);
1178  G4double p2 = pLab*pLab;
1179  G4double p4 = p2*p2;
1180  G4double lm = pLab - 0.8;
1181  G4double md = lm*lm + .652;
1182  fElasticXsc = LE + (cofLogE*ld2 + 2.23)/(1. - .7/sp + .1/p4) + 2./md;
1183  fTotalXsc = (cofLogT*ld2 + 19.5)/(1. + .46/sp + 1.6/p4) + 7.6/md; // + LE;
1184  }
1185  }
1186  else if( (theParticle == theKPlus || theParticle == theK0L) && neutron ) // Kpn/aKn //////////////////////////////////
1187  {
1188  if( pLab < pMin )
1189  {
1190  G4double lm = pLab - 0.94;
1191  G4double md = lm*lm + .392;
1192  fElasticXsc = 2./md;
1193  fTotalXsc = 4.6/md;
1194  }
1195  else if( pLab > pMax )
1196  {
1197  G4double ld = G4Log(pLab) - minLogP;
1198  G4double ld2 = ld*ld;
1199  fElasticXsc = cofLogE*ld2 + 2.23;
1200  fTotalXsc = cofLogT*ld2 + 19.2;
1201  }
1202  else
1203  {
1204  G4double ld = G4Log(pLab) - minLogP;
1205  G4double ld2 = ld*ld;
1206  G4double sp = std::sqrt(pLab);
1207  G4double p2 = pLab*pLab;
1208  G4double p4 = p2*p2;
1209  G4double lm = pLab - 0.8;
1210  G4double md = lm*lm + .652;
1211  fElasticXsc = (cofLogE*ld2 + 2.23)/(1. - .7/sp + .1/p4) + 2./md;
1212  fTotalXsc = (cofLogT*ld2 + 19.5)/(1. + .46/sp + 1.6/p4) + 7.6/md;
1213  }
1214  }
1215  fTotalXsc *= millibarn; // parametrised in mb
1216  fElasticXsc *= millibarn; // parametrised in mb
1217 
1218  if( proton && aParticle->GetDefinition()->GetPDGCharge() > 0. )
1219  {
1220  G4double cB = GetCoulombBarrier(aParticle, nucleon);
1221  fTotalXsc *= cB;
1222  fElasticXsc *= cB;
1223  }
1224  fInelasticXsc = fTotalXsc - fElasticXsc;
1225  if( fInelasticXsc < 0. ) fInelasticXsc = 0.;
1226 
1227  // G4cout<<fTotalXsc/millibarn<<"; "<<fElasticXsc/millibarn<<"; "<<fInelasticXsc/millibarn<<G4endl;
1228 
1229  return fTotalXsc;
1230 }
LE
Definition: Evaluator.cc:66
G4HadronNucleonXsc::GetCoulombBarrier
G4double GetCoulombBarrier(const G4DynamicParticle *aParticle, const G4ParticleDefinition *nucleon)
Definition: G4HadronNucleonXsc.cc:1433
neutron
Definition: G4DataQuestionaire.hh:35
G4DynamicParticle::GetDefinition
G4ParticleDefinition * GetDefinition() const
G4bool
bool G4bool
Definition: G4Types.hh:79
G4Log
G4double G4Log(G4double x)
Definition: G4Log.hh:230
GeV
static constexpr double GeV
Definition: G4SIunits.hh:217
G4double
double G4double
Definition: G4Types.hh:76
G4ParticleDefinition::GetPDGCharge
G4double GetPDGCharge() const
Definition: G4ParticleDefinition.hh:124
CLHEP::Hep3Vector::mag
double mag() const
millibarn
static constexpr double millibarn
Definition: G4SIunits.hh:106
G4DynamicParticle::GetMomentum
G4ThreeVector GetMomentum() const

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4HadronNucleonXsc::GetTotalHadronNucleonXsc ( )
inline

Definition at line 90 of file G4HadronNucleonXsc.hh.

90 { return fTotalXsc; };
G4bool G4HadronNucleonXsc::IsApplicable ( const G4DynamicParticle *  aDP,
const G4Element *  anElement 
)
virtual

Definition at line 97 of file G4HadronNucleonXsc.cc.

99 {
100  G4int Z = G4lrint(anElement->GetZ());
101  G4int A = G4lrint(anElement->GetN());
102  return IsIsoApplicable(aDP, Z, A);
103 }
G4Element::GetN
G4double GetN() const
Definition: G4Element.hh:135
G4Element::GetZ
G4double GetZ() const
Definition: G4Element.hh:131
G4int
int G4int
Definition: G4Types.hh:78
A
double A(double temperature)
Definition: G4DNAElectronHoleRecombination.cc:60
G4HadronNucleonXsc::IsIsoApplicable
virtual G4bool IsIsoApplicable(const G4DynamicParticle *aDP, G4int Z, G4int A)
Definition: G4HadronNucleonXsc.cc:109
G4lrint
int G4lrint(double ad)
Definition: templates.hh:163
Z
G4int Z
Definition: G4ParticleHPJENDLHEData.cc:262

Here is the call graph for this function:

G4bool G4HadronNucleonXsc::IsIsoApplicable ( const G4DynamicParticle *  aDP,
G4int  Z,
G4int  A 
)
virtual

Definition at line 109 of file G4HadronNucleonXsc.cc.

111 {
112  G4bool applicable = false;
113  // G4int baryonNumber = aDP->GetDefinition()->GetBaryonNumber();
114  G4double kineticEnergy = aDP->GetKineticEnergy();
115 
116  const G4ParticleDefinition* theParticle = aDP->GetDefinition();
117 
118  if ( ( kineticEnergy >= fLowerLimit &&
119  Z > 1 && // >= He
120  ( theParticle == theAProton ||
121  theParticle == theGamma ||
122  theParticle == theKPlus ||
123  theParticle == theKMinus ||
124  theParticle == theSMinus) ) ||
125 
126  ( kineticEnergy >= 0.1*fLowerLimit &&
127  Z > 1 && // >= He
128  ( theParticle == theProton ||
129  theParticle == theNeutron ||
130  theParticle == thePiPlus ||
131  theParticle == thePiMinus ) ) ) applicable = true;
132 
133  return applicable;
134 }
G4DynamicParticle::GetKineticEnergy
G4double GetKineticEnergy() const
G4DynamicParticle::GetDefinition
G4ParticleDefinition * GetDefinition() const
G4bool
bool G4bool
Definition: G4Types.hh:79
G4double
double G4double
Definition: G4Types.hh:76
Z
G4int Z
Definition: G4ParticleHPJENDLHEData.cc:262

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