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

#include <G4AntiNuclElastic.hh>

Inheritance diagram for G4AntiNuclElastic:
Inheritance graph
[legend]
Collaboration diagram for G4AntiNuclElastic:
Collaboration graph
[legend]

Public Member Functions

 G4AntiNuclElastic ()
 
virtual ~G4AntiNuclElastic ()
 
virtual G4double SampleInvariantT (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
G4double SampleThetaCMS (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
G4double SampleThetaLab (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
G4double CalculateParticleBeta (const G4ParticleDefinition *particle, G4double momentum)
 
G4double CalculateZommerfeld (G4double beta, G4double Z1, G4double Z2)
 
G4double CalculateAm (G4double momentum, G4double n, G4double Z)
 
G4double DampFactor (G4double z)
 
G4double BesselJzero (G4double z)
 
G4double BesselJone (G4double z)
 
G4double BesselOneByArg (G4double z)
 
G4double GetcosTeta1 (G4double plab, G4int A)
 
G4ComponentAntiNuclNuclearXS * GetComponentCrossSection ()
 
- Public Member Functions inherited from G4HadronElastic
 G4HadronElastic (const G4String &name="hElasticLHEP")
 
virtual ~G4HadronElastic ()
 
virtual G4HadFinalState * ApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
void SetLowestEnergyLimit (G4double value)
 
G4double LowestEnergyLimit () const
 
G4double ComputeMomentumCMS (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
virtual void ModelDescription (std::ostream &) const
 
- Public Member Functions inherited from G4HadronicInteraction
 G4HadronicInteraction (const G4String &modelName="HadronicModel")
 
virtual ~G4HadronicInteraction ()
 
virtual G4bool IsApplicable (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
G4double GetMinEnergy () const
 
G4double GetMinEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
void SetMinEnergy (G4double anEnergy)
 
void SetMinEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMinEnergy (G4double anEnergy, const G4Material *aMaterial)
 
G4double GetMaxEnergy () const
 
G4double GetMaxEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
void SetMaxEnergy (const G4double anEnergy)
 
void SetMaxEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMaxEnergy (G4double anEnergy, const G4Material *aMaterial)
 
G4int GetVerboseLevel () const
 
void SetVerboseLevel (G4int value)
 
const G4String & GetModelName () const
 
void DeActivateFor (const G4Material *aMaterial)
 
void ActivateFor (const G4Material *aMaterial)
 
void DeActivateFor (const G4Element *anElement)
 
void ActivateFor (const G4Element *anElement)
 
G4bool IsBlocked (const G4Material *aMaterial) const
 
G4bool IsBlocked (const G4Element *anElement) const
 
void SetRecoilEnergyThreshold (G4double val)
 
G4double GetRecoilEnergyThreshold () const
 
virtual const std::pair
< G4double, G4double > 		GetFatalEnergyCheckLevels () const
 
virtual std::pair< G4double,
G4double > 		GetEnergyMomentumCheckLevels () const
 
void SetEnergyMomentumCheckLevels (G4double relativeLevel, G4double absoluteLevel)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void InitialiseModel ()
 

Additional Inherited Members

- Protected Member Functions inherited from G4HadronicInteraction
void SetModelName (const G4String &nam)
 
G4bool IsBlocked () const
 
void Block ()
 
- Protected Attributes inherited from G4HadronicInteraction
G4HadFinalState theParticleChange
 
G4int verboseLevel
 
G4double theMinEnergy
 
G4double theMaxEnergy
 
G4bool isBlocked
 

Detailed Description

Definition at line 47 of file G4AntiNuclElastic.hh.

Constructor & Destructor Documentation

G4AntiNuclElastic::G4AntiNuclElastic ( )

Definition at line 57 of file G4AntiNuclElastic.cc.

58  : G4HadronElastic("AntiAElastic")
59 {
60  //V.Ivanchenko commented out
61  //SetMinEnergy( 0.1*GeV );
62  //SetMaxEnergy( 10.*TeV );
63 
64 
65  theAProton = G4AntiProton::AntiProton();
66  theANeutron = G4AntiNeutron::AntiNeutron();
67  theADeuteron = G4AntiDeuteron::AntiDeuteron();
68  theATriton = G4AntiTriton::AntiTriton();
69  theAAlpha = G4AntiAlpha::AntiAlpha();
70  theAHe3 = G4AntiHe3::AntiHe3();
71 
72  theProton = G4Proton::Proton();
73  theNeutron = G4Neutron::Neutron();
74  theDeuteron = G4Deuteron::Deuteron();
75  theAlpha = G4Alpha::Alpha();
76 
77 
78  cs = new G4ComponentAntiNuclNuclearXS();
79  fParticle = 0;
80  fWaveVector = 0.;
81  fBeta = 0.;
82  fZommerfeld = 0.;
83  fAm = 0.;
84  fTetaCMS = 0.;
85  fRa = 0.;
86  fRef = 0.;
87  fceff = 0.;
88  fptot = 0.;
89  fTmax = 0.;
90  fThetaLab = 0.;
91 }
G4AntiHe3::AntiHe3
static G4AntiHe3 * AntiHe3()
Definition: G4AntiHe3.cc:94
G4HadronElastic::G4HadronElastic
G4HadronElastic(const G4String &name="hElasticLHEP")
Definition: G4HadronElastic.cc:48
G4AntiDeuteron::AntiDeuteron
static G4AntiDeuteron * AntiDeuteron()
Definition: G4AntiDeuteron.cc:94
G4AntiAlpha::AntiAlpha
static G4AntiAlpha * AntiAlpha()
Definition: G4AntiAlpha.cc:89
G4AntiProton::AntiProton
static G4AntiProton * AntiProton()
Definition: G4AntiProton.cc:93
G4Proton::Proton
static G4Proton * Proton()
Definition: G4Proton.cc:93
G4Neutron::Neutron
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
G4Deuteron::Deuteron
static G4Deuteron * Deuteron()
Definition: G4Deuteron.cc:94
G4Alpha::Alpha
static G4Alpha * Alpha()
Definition: G4Alpha.cc:89
G4AntiTriton::AntiTriton
static G4AntiTriton * AntiTriton()
Definition: G4AntiTriton.cc:94
G4AntiNeutron::AntiNeutron
static G4AntiNeutron * AntiNeutron()
Definition: G4AntiNeutron.cc:103

Here is the call graph for this function:

G4AntiNuclElastic::~G4AntiNuclElastic ( )
virtual

Definition at line 94 of file G4AntiNuclElastic.cc.

95 {
96  delete cs;
97 }

Member Function Documentation

G4double G4AntiNuclElastic::BesselJone ( G4double  z)

Definition at line 592 of file G4AntiNuclElastic.cc.

593 {
594  G4double modvalue, value2, fact1, fact2, arg, shift, bessel;
595 
596  modvalue = std::fabs(value);
597 
598  if ( modvalue < 8.0 )
599  {
600  value2 = value*value;
601  fact1 = value*(72362614232.0 + value2*(-7895059235.0
602  + value2*( 242396853.1
603  + value2*(-2972611.439
604  + value2*( 15704.48260
605  + value2*(-30.16036606 ) ) ) ) ) );
606 
607  fact2 = 144725228442.0 + value2*(2300535178.0
608  + value2*(18583304.74
609  + value2*(99447.43394
610  + value2*(376.9991397
611  + value2*1.0 ) ) ) );
612  bessel = fact1/fact2;
613  }
614  else
615  {
616  arg = 8.0/modvalue;
617  value2 = arg*arg;
618 
619  shift = modvalue - 2.356194491;
620 
621  fact1 = 1.0 + value2*( 0.183105e-2
622  + value2*(-0.3516396496e-4
623  + value2*(0.2457520174e-5
624  + value2*(-0.240337019e-6 ) ) ) );
625 
626  fact2 = 0.04687499995 + value2*(-0.2002690873e-3
627  + value2*( 0.8449199096e-5
628  + value2*(-0.88228987e-6
629  + value2*0.105787412e-6 ) ) );
630 
631  bessel = std::sqrt( 0.636619772/modvalue)*(std::cos(shift)*fact1 - arg*std::sin(shift)*fact2);
632  if (value < 0.0) bessel = -bessel;
633  }
634  return bessel;
635 }
value
const XML_Char int const XML_Char * value
Definition: expat.h:331
G4double
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:

G4double G4AntiNuclElastic::BesselJzero ( G4double  z)

Definition at line 541 of file G4AntiNuclElastic.cc.

542 {
543  G4double modvalue, value2, fact1, fact2, arg, shift, bessel;
544 
545  modvalue = std::fabs(value);
546 
547  if ( value < 8.0 && value > -8.0 )
548  {
549  value2 = value*value;
550 
551  fact1 = 57568490574.0 + value2*(-13362590354.0
552  + value2*( 651619640.7
553  + value2*(-11214424.18
554  + value2*( 77392.33017
555  + value2*(-184.9052456 ) ) ) ) );
556 
557  fact2 = 57568490411.0 + value2*( 1029532985.0
558  + value2*( 9494680.718
559  + value2*(59272.64853
560  + value2*(267.8532712
561  + value2*1.0 ) ) ) );
562 
563  bessel = fact1/fact2;
564  }
565  else
566  {
567  arg = 8.0/modvalue;
568 
569  value2 = arg*arg;
570 
571  shift = modvalue-0.785398164;
572 
573  fact1 = 1.0 + value2*(-0.1098628627e-2
574  + value2*(0.2734510407e-4
575  + value2*(-0.2073370639e-5
576  + value2*0.2093887211e-6 ) ) );
577  fact2 = -0.1562499995e-1 + value2*(0.1430488765e-3
578  + value2*(-0.6911147651e-5
579  + value2*(0.7621095161e-6
580  - value2*0.934945152e-7 ) ) );
581 
582  bessel = std::sqrt(0.636619772/modvalue)*(std::cos(shift)*fact1 - arg*std::sin(shift)*fact2);
583  }
584  return bessel;
585 }
value
const XML_Char int const XML_Char * value
Definition: expat.h:331
G4double
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:

G4double G4AntiNuclElastic::BesselOneByArg ( G4double  z)

Definition at line 639 of file G4AntiNuclElastic.cc.

640 {
641  G4double x2, result;
642 
643  if( std::fabs(x) < 0.01 )
644  {
645  x *= 0.5;
646  x2 = x*x;
647  result = (2.- x2 + x2*x2/6.)/4.;
648  }
649  else
650  {
651  result = BesselJone(x)/x;
652  }
653  return result;
654 }
result
G4double G4ParticleHPJENDLHEData::G4double result
Definition: G4ParticleHPJENDLHEData.cc:257
G4AntiNuclElastic::BesselJone
G4double BesselJone(G4double z)
Definition: G4AntiNuclElastic.cc:592
G4double
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4AntiNuclElastic::CalculateAm ( G4double  momentum,
G4double  n,
G4double  Z 
)

Definition at line 525 of file G4AntiNuclElastic.cc.

526 {
527  G4double k = momentum/hbarc;
528  G4double ch = 1.13 + 3.76*n*n;
529  G4double zn = 1.77*k/G4Pow::GetInstance()->A13(Z)*Bohr_radius;
530  G4double zn2 = zn*zn;
531  fAm = ch/zn2;
532 
533  return fAm;
534 }
G4Pow::GetInstance
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
CLHEP::Bohr_radius
static constexpr double Bohr_radius
Definition: PhysicalConstants.h:99
CLHEP::hbarc
static constexpr double hbarc
Definition: PhysicalConstants.h:65
G4Pow::A13
G4double A13(G4double A) const
Definition: G4Pow.hh:132
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:

G4double G4AntiNuclElastic::CalculateParticleBeta ( const G4ParticleDefinition *  particle,
G4double  momentum 
)

Definition at line 502 of file G4AntiNuclElastic.cc.

504 {
505  G4double mass = particle->GetPDGMass();
506  G4double a = momentum/mass;
507  fBeta = a/std::sqrt(1+a*a);
508 
509  return fBeta;
510 }
G4ParticleDefinition::GetPDGMass
G4double GetPDGMass() const
Definition: G4ParticleDefinition.hh:122
G4double
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4AntiNuclElastic::CalculateZommerfeld ( G4double  beta,
G4double  Z1,
G4double  Z2 
)

Definition at line 516 of file G4AntiNuclElastic.cc.

517 {
518  fZommerfeld = fine_structure_const*Z1*Z2/beta;
519 
520  return fZommerfeld;
521 }
CLHEP::fine_structure_const
static constexpr double fine_structure_const
Definition: PhysicalConstants.h:96

Here is the caller graph for this function:

G4double G4AntiNuclElastic::DampFactor ( G4double  z)

Definition at line 482 of file G4AntiNuclElastic.cc.

483 {
484  G4double df;
485  G4double f3 = 6.; // first factorials
486 
487  if( std::fabs(x) < 0.01 )
488  {
489  df=1./(1.+x*x/f3);
490  }
491  else
492  {
493  df = x/std::sinh(x);
494  }
495  return df;
496 }
G4double
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:

G4ComponentAntiNuclNuclearXS * G4AntiNuclElastic::GetComponentCrossSection ( )
inline

Definition at line 123 of file G4AntiNuclElastic.hh.

124 {
125  return cs;
126 }

Here is the caller graph for this function:

G4double G4AntiNuclElastic::GetcosTeta1 ( G4double  plab,
G4int  A 
)

Definition at line 658 of file G4AntiNuclElastic.cc.

659 {
660 
661 // G4double p0 =G4LossTableManager::Instance()->FactorForAngleLimit()*CLHEP::hbarc/CLHEP::fermi;
662  G4double p0 = 1.*hbarc/fermi;
663 //G4double cteta1 = 1.0 - p0*p0/2.0 * pow(A,2./3.)/(plab*plab);
664  G4double cteta1 = 1.0 - p0*p0/2.0 * G4Pow::GetInstance()->Z23(A)/(plab*plab);
666  if(cteta1 < -1.) cteta1 = -1.0;
667  return cteta1;
668 }
G4Pow::GetInstance
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
CLHEP::hbarc
static constexpr double hbarc
Definition: PhysicalConstants.h:65
A
double A(double temperature)
Definition: G4DNAElectronHoleRecombination.cc:60
G4Pow::Z23
G4double Z23(G4int Z) const
Definition: G4Pow.hh:154
G4double
double G4double
Definition: G4Types.hh:76
fermi
static constexpr double fermi
Definition: G4SIunits.hh:103

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4AntiNuclElastic::SampleInvariantT ( const G4ParticleDefinition *  p,
G4double  plab,
G4int  Z,
G4int  A 
)
virtual

Reimplemented from G4HadronElastic.

Definition at line 101 of file G4AntiNuclElastic.cc.

103 {
104  G4double T;
105  G4double Mproj = particle->GetPDGMass();
106  G4LorentzVector Pproj(0.,0.,Plab,std::sqrt(Plab*Plab+Mproj*Mproj));
107  G4double ctet1 = GetcosTeta1(Plab, A);
108 
109  G4double energy=Pproj.e()-Mproj;
110 
111  const G4ParticleDefinition* theParticle = particle;
112 
113  G4ParticleDefinition * theDef = 0;
114 
115  if(Z == 1 && A == 1) theDef = theProton;
116  else if (Z == 1 && A == 2) theDef = theDeuteron;
117  else if (Z == 1 && A == 3) theDef = G4Triton::Triton();
118  else if (Z == 2 && A == 3) theDef = G4He3::He3();
119  else if (Z == 2 && A == 4) theDef = theAlpha;
120 
121 
122  G4double TargMass =G4NucleiProperties::GetNuclearMass(A,Z);
123 
124  //transform to CMS
125 
126  G4LorentzVector lv(0.0,0.0,0.0,TargMass);
127  lv += Pproj;
128  G4double S = lv.mag2()/GeV/GeV;
129 
130  G4ThreeVector bst = lv.boostVector();
131  Pproj.boost(-bst);
132 
133  G4ThreeVector p1 = Pproj.vect();
134  G4double ptot = p1.mag();
135 
136  fbst = bst;
137  fptot= ptot;
138  fTmax = 4.0*ptot*ptot;
139 
140  if(Plab/std::abs(particle->GetBaryonNumber()) < 100.*MeV) // Uzhi 24 Nov. 2011
141  {return fTmax*G4UniformRand();} // Uzhi 24 Nov. 2011
142 
143  G4double Z1 = particle->GetPDGCharge();
144  G4double Z2 = Z;
145 
146  G4double beta = CalculateParticleBeta(particle, ptot);
147  G4double n = CalculateZommerfeld( beta, Z1, Z2 );
148  G4double Am = CalculateAm( ptot, n, Z2 );
149  fWaveVector = ptot; // /hbarc;
150 
151  G4LorentzVector Fproj(0.,0.,0.,0.);
152  G4double XsCoulomb = sqr(n/fWaveVector)*pi*(1+ctet1)/(1.+Am)/(1.+2.*Am-ctet1);
153  XsCoulomb=XsCoulomb*0.38938e+6;
154 
155 
156  G4double XsElastHad =cs->GetElasticElementCrossSection(particle, energy, Z, (G4double)A);
157  G4double XstotalHad =cs->GetTotalElementCrossSection(particle, energy, Z, (G4double)A);
158 
159 
160  XsElastHad/=millibarn; XstotalHad/=millibarn;
161 
162 
163  G4double CoulombProb = XsCoulomb/(XsCoulomb+XsElastHad);
164 
165 // G4cout<<" XselastHadron " << XsElastHad << " XsCol "<< XsCoulomb <<G4endl;
166 // G4cout <<" XsTotal" << XstotalHad <<G4endl;
167 // G4cout<<"XsInel"<< XstotalHad-XsElastHad<<G4endl;
168 
169  if(G4UniformRand() < CoulombProb)
170  { // Simulation of Coulomb scattering
171 
172  G4double phi = twopi * G4UniformRand();
173  G4double Ksi = G4UniformRand();
174 
175  G4double par1 = 2.*(1.+Am)/(1.+ctet1);
176 
177 // ////sample ThetaCMS in Coulomb part
178 
179  G4double cosThetaCMS = (par1*ctet1- Ksi*(1.+2.*Am))/(par1-Ksi);
180 
181  G4double PtZ=ptot*cosThetaCMS;
182  Fproj.setPz(PtZ);
183  G4double PtProjCMS = ptot*std::sqrt(1.0 - cosThetaCMS*cosThetaCMS);
184  G4double PtX= PtProjCMS * std::cos(phi);
185  G4double PtY= PtProjCMS * std::sin(phi);
186  Fproj.setPx(PtX);
187  Fproj.setPy(PtY);
188  Fproj.setE(std::sqrt(PtX*PtX+PtY*PtY+PtZ*PtZ+Mproj*Mproj));
189  T = -(Pproj-Fproj).mag2();
190  } else
191 
192  {
194 
195 // G4double Qmax = 2.*ptot*197.33; // in fm^-1
196  G4double Qmax = 2.*3.0*197.33; // in fm^-1
197  G4double Amag = 70*70; // A1 in Magora funct:A1*exp(-q*A2)
198  G4double SlopeMag = 2.*3.0; // A2 in Magora funct:A1*exp(-q*A2)
199 
200  G4double sig_pbarp= cs->GetAntiHadronNucleonTotCrSc(particle,energy);
201 
202  fRa = 1.113*G4Pow::GetInstance()->Z13(A) -
203  0.227/G4Pow::GetInstance()->Z13(A);
204  if(A == 3) fRa=1.81;
205  if(A == 4) fRa=1.37;
206 
207  if((A>=12.) && (A<27) ) fRa=fRa*0.85;
208  if((A>=27.) && (A<48) ) fRa=fRa*0.90;
209  if((A>=48.) && (A<65) ) fRa=fRa*0.95;
210 
211  G4double Ref2 = 0;
212  G4double ceff2 =0;
213  G4double rho = 0;
214  if ((theParticle == theAProton) || (theParticle == theANeutron))
215  {
216  if(theDef == theProton)
217  {
218 // G4double Mp2=sqr(theDef->GetPDGMass()/GeV );
219 
220 // change 30 October
221 
222  if(Plab < 610.)
223  { rho = 1.3347-10.342*Plab/1000.+22.277*Plab/1000.*Plab/1000.-
224  13.634*Plab/1000.*Plab/1000.*Plab/1000. ;}
225  if((Plab < 5500.)&&(Plab >= 610.) )
226  { rho = 0.22; }
227  if((Plab >= 5500.)&&(Plab < 12300.) )
228  { rho = -0.32; }
229  if( Plab >= 12300.)
230  { rho = 0.135-2.26/(std::sqrt(S)) ;}
231 
232  Ref2 = 0.35 + 0.9/std::sqrt(std::sqrt(S-4.*0.88))+0.04*G4Log(S) ;
233  ceff2 = 0.375 - 2./S + 0.44/(sqr(S-4.)+1.5) ;
234 
235 /*
236  Ref2=0.8/std::sqrt(std::sqrt(S-4.*Mp2)) + 0.55;
237  if(S>1000.) Ref2=0.62+0.02*G4Log(S) ;
238  ceff2 = 0.035/(sqr(S-4.3)+0.4) + 0.085 * G4Log(S) ;
239  if(S>1000.) ceff2 = 0.005 * G4Log(S) + 0.29;
240 */
241 
242  Ref2=Ref2*Ref2;
243  ceff2 = ceff2*ceff2;
244 
245  SlopeMag = 0.5; // Uzhi
246  Amag= 1.; // Uzhi
247  }
248 
249  if(Z>2)
250  { Ref2 = fRa*fRa +2.48*0.01*sig_pbarp*fRa - 2.23e-6*sig_pbarp*sig_pbarp*fRa*fRa;
251  ceff2 = 0.16+3.3e-4*sig_pbarp+0.35*G4Exp(-0.03*sig_pbarp);
252  }
253  if( (Z==2)&&(A==4) )
254  { Ref2 = fRa*fRa -0.46 +0.03*sig_pbarp - 2.98e-6*sig_pbarp*sig_pbarp;
255  ceff2= 0.078 + 6.657e-4*sig_pbarp + 0.3359*G4Exp(-0.03*sig_pbarp);
256  }
257  if( (Z==1)&&(A==3) )
258  { Ref2 = fRa*fRa - 1.36 + 0.025 * sig_pbarp - 3.69e-7 * sig_pbarp*sig_pbarp;
259  ceff2 = 0.149 + 7.091e-04*sig_pbarp + 0.3743*G4Exp(-0.03*sig_pbarp);
260  }
261  if( (Z==2)&&(A==3) )
262  { Ref2 = fRa*fRa - 1.36 + 0.025 * sig_pbarp - 3.69e-7 * sig_pbarp*sig_pbarp;
263  ceff2 = 0.149 + 7.091e-04*sig_pbarp + 0.3743*G4Exp(-0.03*sig_pbarp);
264  }
265  if( (Z==1)&&(A==2) )
266  {
267  Ref2 = fRa*fRa - 0.28 + 0.019 * sig_pbarp + 2.06e-6 * sig_pbarp*sig_pbarp;
268  ceff2 = 0.297 + 7.853e-04*sig_pbarp + 0.2899*G4Exp(-0.03*sig_pbarp);
269  }
270  }
271 
272  if (theParticle == theADeuteron)
273  {
274  sig_pbarp= cs->GetAntiHadronNucleonTotCrSc(particle,energy/2.);
275  Ref2 = XstotalHad/10./2./pi ;
276  if(Z>2)
277  {
278  ceff2 = 0.38 + 2.0e-4 *sig_pbarp + 0.5 * G4Exp(-0.03*sig_pbarp);
279  }
280  if(theDef == theProton)
281  {
282  ceff2 = 0.297 + 7.853e-04*sig_pbarp + 0.2899*G4Exp(-0.03*sig_pbarp);
283  }
284  if(theDef == theDeuteron)
285  {
286  ceff2 = 0.65 + 3.0e-4*sig_pbarp + 0.55 * G4Exp(-0.03*sig_pbarp);
287  }
288  if( (theDef == G4Triton::Triton()) || (theDef == G4He3::He3() ) )
289  {
290  ceff2 = 0.57 + 2.5e-4*sig_pbarp + 0.65 * G4Exp(-0.02*sig_pbarp);
291  }
292  if(theDef == theAlpha)
293  {
294  ceff2 = 0.40 + 3.5e-4 *sig_pbarp + 0.45 * G4Exp(-0.02*sig_pbarp);
295  }
296  }
297 
298  if( (theParticle ==theAHe3) || (theParticle ==theATriton) )
299  {
300  sig_pbarp = cs->GetAntiHadronNucleonTotCrSc(particle,energy/3.);
301  Ref2 = XstotalHad/10./2./pi ;
302  if(Z>2)
303  {
304  ceff2 = 0.26 + 2.2e-4*sig_pbarp + 0.33*G4Exp(-0.03*sig_pbarp);
305  }
306  if(theDef == theProton)
307  {
308  ceff2 = 0.149 + 7.091e-04*sig_pbarp + 0.3743*G4Exp(-0.03*sig_pbarp);
309  }
310  if(theDef == theDeuteron)
311  {
312  ceff2 = 0.57 + 2.5e-4*sig_pbarp + 0.65 * G4Exp(-0.02*sig_pbarp);
313  }
314  if( (theDef == G4Triton::Triton()) || (theDef == G4He3::He3() ) )
315  {
316  ceff2 = 0.39 + 2.7e-4*sig_pbarp + 0.7 * G4Exp(-0.02*sig_pbarp);
317  }
318  if(theDef == theAlpha)
319  {
320  ceff2 = 0.24 + 3.5e-4*sig_pbarp + 0.75 * G4Exp(-0.03*sig_pbarp);
321  }
322  }
323 
324 
325  if (theParticle == theAAlpha)
326  {
327  sig_pbarp = cs->GetAntiHadronNucleonTotCrSc(particle,energy/3.);
328  Ref2 = XstotalHad/10./2./pi ;
329  if(Z>2)
330  {
331  ceff2 = 0.22 + 2.0e-4*sig_pbarp + 0.2 * G4Exp(-0.03*sig_pbarp);
332  }
333  if(theDef == theProton)
334  {
335  ceff2= 0.078 + 6.657e-4*sig_pbarp + 0.3359*G4Exp(-0.03*sig_pbarp);
336  }
337  if(theDef == theDeuteron)
338  {
339  ceff2 = 0.40 + 3.5e-4 *sig_pbarp + 0.45 * G4Exp(-0.02*sig_pbarp);
340  }
341  if( (theDef == G4Triton::Triton()) || (theDef == G4He3::He3() ) )
342  {
343  ceff2 = 0.24 + 3.5e-4*sig_pbarp + 0.75 * G4Exp(-0.03*sig_pbarp);
344  }
345  if(theDef == theAlpha)
346  {
347  ceff2 = 0.17 + 3.5e-4*sig_pbarp + 0.45 * G4Exp(-0.03*sig_pbarp);
348  }
349  }
350 
351  fRef=std::sqrt(Ref2);
352  fceff = std::sqrt(ceff2);
353 // G4cout<<" Ref "<<fRef<<" c_eff "<<fceff<< " rho "<< rho<<G4endl;
354 
355 
356  G4double Q = 0.0 ;
357  G4double BracFunct;
358  const G4int maxNumberOfLoops = 10000;
359  G4int loopCounter = 0;
360  do
361  {
362  Q = -G4Log(1.-(1.- G4Exp(-SlopeMag * Qmax))* G4UniformRand() )/SlopeMag;
363  G4double x = fRef * Q;
364  BracFunct = ( ( sqr(BesselOneByArg(x))+sqr(rho/2. * BesselJzero(x)) )
365 * sqr(DampFactor(pi*fceff*Q))) /(Amag*G4Exp(-SlopeMag*Q));
366 
367  BracFunct = BracFunct * Q * sqr(sqr(fRef));
368  }
369  while ( (G4UniformRand()>BracFunct) &&
370  ++loopCounter < maxNumberOfLoops ); /* Loop checking, 10.08.2015, A.Ribon */
371  if ( loopCounter >= maxNumberOfLoops ) {
372  fTetaCMS = 0.0;
373  return 0.0;
374  }
375 
376  T= sqr(Q);
377  T*=3.893913e+4; // fm -> MeV^2
378  }
379 
380  G4double cosTet=1.0-T/(2.*ptot*ptot);
381  if(cosTet > 1.0 ) cosTet= 1.; // Uzhi 30 Nov.
382  if(cosTet < -1.0 ) cosTet=-1.; // Uzhi 30 Nov.
383  fTetaCMS=std::acos(cosTet);
384 
385  return T;
386 }
G4Pow::GetInstance
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
G4NucleiProperties::GetNuclearMass
static G4double GetNuclearMass(const G4double A, const G4double Z)
Definition: G4NucleiProperties.cc:64
G4ComponentAntiNuclNuclearXS::GetElasticElementCrossSection
virtual G4double GetElasticElementCrossSection(const G4ParticleDefinition *aParticle, G4double kinEnergy, G4int Z, G4double A)
Definition: G4ComponentAntiNuclNuclearXS.cc:269
S
double S(double temp)
Definition: G4DNAElectronHoleRecombination.cc:77
Q
static double Q[]
Definition: nf_gammaFunctions.cc:91
G4AntiNuclElastic::BesselOneByArg
G4double BesselOneByArg(G4double z)
Definition: G4AntiNuclElastic.cc:639
G4int
int G4int
Definition: G4Types.hh:78
twopi
static constexpr double twopi
Definition: G4SIunits.hh:76
G4AntiNuclElastic::CalculateParticleBeta
G4double CalculateParticleBeta(const G4ParticleDefinition *particle, G4double momentum)
Definition: G4AntiNuclElastic.cc:502
G4ComponentAntiNuclNuclearXS::GetAntiHadronNucleonTotCrSc
G4double GetAntiHadronNucleonTotCrSc(const G4ParticleDefinition *aParticle, G4double kinEnergy)
Definition: G4ComponentAntiNuclNuclearXS.cc:292
G4UniformRand
#define G4UniformRand()
Definition: Randomize.hh:97
G4Pow::Z13
G4double Z13(G4int Z) const
Definition: G4Pow.hh:127
A
double A(double temperature)
Definition: G4DNAElectronHoleRecombination.cc:60
G4Triton::Triton
static G4Triton * Triton()
Definition: G4Triton.cc:95
G4ComponentAntiNuclNuclearXS::GetTotalElementCrossSection
virtual G4double GetTotalElementCrossSection(const G4ParticleDefinition *aParticle, G4double kinEnergy, G4int Z, G4double A)
Definition: G4ComponentAntiNuclNuclearXS.cc:89
G4Log
G4double G4Log(G4double x)
Definition: G4Log.hh:230
G4Exp
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
G4AntiNuclElastic::BesselJzero
G4double BesselJzero(G4double z)
Definition: G4AntiNuclElastic.cc:541
G4AntiNuclElastic::CalculateAm
G4double CalculateAm(G4double momentum, G4double n, G4double Z)
Definition: G4AntiNuclElastic.cc:525
G4AntiNuclElastic::DampFactor
G4double DampFactor(G4double z)
Definition: G4AntiNuclElastic.cc:482
G4INCL::KinematicsUtils::energy
G4double energy(const ThreeVector &p, const G4double m)
Definition: G4INCLKinematicsUtils.cc:157
GeV
static constexpr double GeV
Definition: G4SIunits.hh:217
MeV
static constexpr double MeV
Definition: G4SIunits.hh:214
pi
static constexpr double pi
Definition: G4SIunits.hh:75
sqr
T sqr(const T &x)
Definition: templates.hh:145
G4AntiNuclElastic::CalculateZommerfeld
G4double CalculateZommerfeld(G4double beta, G4double Z1, G4double Z2)
Definition: G4AntiNuclElastic.cc:516
G4double
double G4double
Definition: G4Types.hh:76
G4He3::He3
static G4He3 * He3()
Definition: G4He3.cc:94
CLHEP::Hep3Vector::mag
double mag() const
G4AntiNuclElastic::GetcosTeta1
G4double GetcosTeta1(G4double plab, G4int A)
Definition: G4AntiNuclElastic.cc:658
millibarn
static constexpr double millibarn
Definition: G4SIunits.hh:106
Z
G4int Z
Definition: G4ParticleHPJENDLHEData.cc:262

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4AntiNuclElastic::SampleThetaCMS ( const G4ParticleDefinition *  p,
G4double  plab,
G4int  Z,
G4int  A 
)

Definition at line 390 of file G4AntiNuclElastic.cc.

392 {
393  G4double T;
394  T = SampleInvariantT( p, plab, Z, A);
395 
396  // NaN finder
397  if(!(T < 0.0 || T >= 0.0))
398  {
399  if (verboseLevel > 0)
400  {
401  G4cout << "G4DiffuseElastic:WARNING: A = " << A
402  << " mom(GeV)= " << plab/GeV
403  << " S-wave will be sampled"
404  << G4endl;
405  }
406  T = G4UniformRand()*fTmax;
407 
408  }
409 
410  if(fptot > 0.) // Uzhi 24 Nov. 2011
411  {
412  G4double cosTet=1.0-T/(2.*fptot*fptot);
413  if(cosTet > 1.0 ) cosTet= 1.; // Uzhi 30 Nov.
414  if(cosTet < -1.0 ) cosTet=-1.; // Uzhi 30 Nov.
415  fTetaCMS=std::acos(cosTet);
416  return fTetaCMS;
417  } else // Uzhi 24 Nov. 2011
418  { // Uzhi 24 Nov. 2011
419  return 2.*G4UniformRand()-1.; // Uzhi 24 Nov. 2011
420  } // Uzhi 24 Nov. 2011
421 }
G4AntiNuclElastic::SampleInvariantT
virtual G4double SampleInvariantT(const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
Definition: G4AntiNuclElastic.cc:101
G4UniformRand
#define G4UniformRand()
Definition: Randomize.hh:97
G4cout
G4GLOB_DLL std::ostream G4cout
A
double A(double temperature)
Definition: G4DNAElectronHoleRecombination.cc:60
GeV
static constexpr double GeV
Definition: G4SIunits.hh:217
G4endl
#define G4endl
Definition: G4ios.hh:61
G4double
double G4double
Definition: G4Types.hh:76
Z
G4int Z
Definition: G4ParticleHPJENDLHEData.cc:262

Here is the call graph for this function:

G4double G4AntiNuclElastic::SampleThetaLab ( const G4ParticleDefinition *  p,
G4double  plab,
G4int  Z,
G4int  A 
)

Definition at line 426 of file G4AntiNuclElastic.cc.

428 {
429  G4double T;
430  T = SampleInvariantT( p, plab, Z, A);
431 
432  // NaN finder
433  if(!(T < 0.0 || T >= 0.0))
434  {
435  if (verboseLevel > 0)
436  {
437  G4cout << "G4DiffuseElastic:WARNING: A = " << A
438  << " mom(GeV)= " << plab/GeV
439  << " S-wave will be sampled"
440  << G4endl;
441  }
442  T = G4UniformRand()*fTmax;
443  }
444 
445  G4double phi = G4UniformRand()*twopi;
446 
447  G4double cost(1.);
448  if(fTmax > 0.) {cost = 1. - 2.0*T/fTmax;} // Uzhi 24 Nov. 2011
449 
450  G4double sint;
451  if( cost >= 1.0 )
452  {
453  cost = 1.0;
454  sint = 0.0;
455  }
456  else if( cost <= -1.0)
457  {
458  cost = -1.0;
459  sint = 0.0;
460  }
461  else
462  {
463  sint = std::sqrt((1.0-cost)*(1.0+cost));
464  }
465 
466  G4double m1 = p->GetPDGMass();
467  G4ThreeVector v(sint*std::cos(phi),sint*std::sin(phi),cost);
468  v *= fptot;
469  G4LorentzVector nlv(v.x(),v.y(),v.z(),std::sqrt(fptot*fptot + m1*m1));
470 
471  nlv.boost(fbst);
472 
473  G4ThreeVector np = nlv.vect();
474  G4double theta = np.theta();
475  fThetaLab = theta;
476 
477  return theta;
478 }
G4AntiNuclElastic::SampleInvariantT
virtual G4double SampleInvariantT(const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
Definition: G4AntiNuclElastic.cc:101
twopi
static constexpr double twopi
Definition: G4SIunits.hh:76
G4UniformRand
#define G4UniformRand()
Definition: Randomize.hh:97
G4cout
G4GLOB_DLL std::ostream G4cout
A
double A(double temperature)
Definition: G4DNAElectronHoleRecombination.cc:60
CLHEP::HepLorentzVector::boost
HepLorentzVector & boost(double, double, double)
Definition: LorentzVector.cc:59
CLHEP::Hep3Vector::theta
double theta() const
G4ParticleDefinition::GetPDGMass
G4double GetPDGMass() const
Definition: G4ParticleDefinition.hh:122
GeV
static constexpr double GeV
Definition: G4SIunits.hh:217
G4endl
#define G4endl
Definition: G4ios.hh:61
G4double
double G4double
Definition: G4Types.hh:76
Z
G4int Z
Definition: G4ParticleHPJENDLHEData.cc:262

Here is the call graph for this function:

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