Geant4  10.03.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4ComponentAntiNuclNuclearXS.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 // Calculation of the total, elastic and inelastic cross-sections
27 // of anti-nucleon and anti-nucleus interactions with nuclei
28 // based on Glauber approach and V. Grishine formulaes for
29 // interpolations (ref. V.M.Grichine, Eur.Phys.J., C62(2009) 399;
30 // NIM, B267 (2009) 2460) and our parametrization of hadron-nucleon
31 // cross-sections
32 //
33 //
34 // Created by A.Galoyan and V. Uzhinsky, 18.11.2010
35 
36 
38 
39 #include "G4PhysicalConstants.hh"
40 #include "G4SystemOfUnits.hh"
41 #include "G4ParticleTable.hh"
42 #include "G4IonTable.hh"
43 #include "G4ParticleDefinition.hh"
44 #include "G4Pow.hh"
45 
47 
49 : G4VComponentCrossSection("AntiAGlauber"),
50 // fUpperLimit(10000*GeV), fLowerLimit(10*MeV),
51  fRadiusEff(0.0), fRadiusNN2(0.0),
52  fTotalXsc(0.0), fElasticXsc(0.0), fInelasticXsc(0.0),
53  fAntiHadronNucleonTotXsc(0.0), fAntiHadronNucleonElXsc(0.0),
54  Elab(0.0), S(0.0), SqrtS(0)
55 {
56  theAProton = G4AntiProton::AntiProton();
57  theANeutron = G4AntiNeutron::AntiNeutron();
58  theADeuteron = G4AntiDeuteron::AntiDeuteron();
59  theATriton = G4AntiTriton::AntiTriton();
60  theAAlpha = G4AntiAlpha::AntiAlpha();
61  theAHe3 = G4AntiHe3::AntiHe3();
62 
63  Mn = 0.93827231; // GeV
64  b0 = 11.92; // GeV^(-2)
65  b2 = 0.3036; // GeV^(-2)
66  SqrtS0 = 20.74; // GeV
67  S0 = 33.0625; // GeV^2
68  R0 = 1.0; // default value (V.Ivanchenko)
69 }
70 
72 //
73 //
74 
76 {
77 }
78 
80 
81 
82 
84 //
85 // Calculation of total CrossSection of Anti-Nucleus - Nucleus
86 
87 
89 (const G4ParticleDefinition* aParticle, G4double kinEnergy, G4int Z, G4double A)
90 {
91  G4double xsection, sigmaTotal, sigmaElastic;
92 
93  const G4ParticleDefinition* theParticle = aParticle;
94 
95  sigmaTotal = GetAntiHadronNucleonTotCrSc(theParticle,kinEnergy);
96  sigmaElastic = GetAntiHadronNucleonElCrSc(theParticle,kinEnergy);
97 
98 // calculation of squared radius of NN-collision
99  fRadiusNN2=sigmaTotal*sigmaTotal*0.1/(8.*sigmaElastic*pi) ; //fm^2
100 
101 // calculation of effective nuclear radius for Pbar and Nbar interactions (can be changed)
102 
103  //A.R. 29-Jan-2013 : use antiprotons/antineutrons as the default case,
104  // to be used for instance, as first approximation
105  // without validation, for anti-hyperons.
106  //if ( (theParticle == theAProton) || (theParticle == theANeutron) ) {
107  if(A==1)
108  { fTotalXsc = sigmaTotal * millibarn;
109  return fTotalXsc; }
110 
111  fRadiusEff = 1.34*G4Pow::GetInstance()->powA(A,0.23)+1.35/G4Pow::GetInstance()->powA(A,1./3.); //fm
112 
113  if( (Z==1) && (A==2) ) fRadiusEff = 3.800; //fm
114  if( (Z==1) && (A==3) ) fRadiusEff = 3.300;
115  if( (Z==2) && (A==3) ) fRadiusEff = 3.300;
116  if( (Z==2) && (A==4) ) fRadiusEff = 2.376;
117  //}
118 
119 //calculation of effective nuclear radius for AntiDeuteron interaction (can be changed)
120  if (theParticle == theADeuteron)
121  { fRadiusEff = 1.46 * G4Pow::GetInstance()->powA(A,0.21) + 1.45 / G4Pow::GetInstance()->powA(A,1./3.);
122 
123  if( (Z==1) && (A==2) ) fRadiusEff = 3.238; //fm
124  if( (Z==1) && (A==3) ) fRadiusEff = 3.144;
125  if( (Z==2) && (A==3) ) fRadiusEff = 3.144;
126  if( (Z==2) && (A==4) ) fRadiusEff = 2.544;
127  }
128 // calculation of effective nuclear radius for AntiHe3 interaction (can be changed)
129 
130  if( (theParticle ==theAHe3) || (theParticle ==theATriton) )
131  { fRadiusEff = 1.40* G4Pow::GetInstance()->powA(A,0.21)+1.63/G4Pow::GetInstance()->powA(A,1./3.);
132 
133  if( (Z==1) && (A==2) ) fRadiusEff = 3.144; //fm
134  if( (Z==1) && (A==3) ) fRadiusEff = 3.075;
135  if( (Z==2) && (A==3) ) fRadiusEff = 3.075;
136  if( (Z==2) && (A==4) ) fRadiusEff = 2.589;
137  }
138 
139 //calculation of effective nuclear radius for AntiAlpha interaction (can be changed)
140 
141  if (theParticle == theAAlpha)
142  {
143  fRadiusEff = 1.35* G4Pow::GetInstance()->powA(A,0.21)+1.1/G4Pow::GetInstance()->powA(A,1./3.);
144 
145  if( (Z==1) && (A==2) ) fRadiusEff = 2.544; //fm
146  if( (Z==1) && (A==3) ) fRadiusEff = 2.589;
147  if( (Z==2) && (A==3) ) fRadiusEff = 2.589;
148  if( (Z==2) && (A==4) ) fRadiusEff = 2.241;
149 
150  }
151 
152  G4double R2 = fRadiusEff*fRadiusEff;
153  G4double REf2 = R2+fRadiusNN2;
154  G4double ApAt = std::abs(theParticle->GetBaryonNumber()) * A;
155 
156  xsection = 2*pi*REf2*10.*G4Log(1+(ApAt*sigmaTotal/(2*pi*REf2*10.))); //mb
157  xsection =xsection *millibarn;
158  fTotalXsc = xsection;
159 
160  return fTotalXsc;
161 }
162 
163 
165 //
166 // Calculation of total CrossSection of Anti-Nucleus - Nucleus
169 (const G4ParticleDefinition* aParticle, G4double kinEnergy, G4int Z, G4int A )
170 { return GetTotalElementCrossSection(aParticle, kinEnergy, Z, (G4double) A); }
171 
173 // Calculation of inelastic CrossSection of Anti-Nucleus - Nucleus
175 
177 (const G4ParticleDefinition* aParticle, G4double kinEnergy, G4int Z, G4double A)
178 {
179  G4double inelxsection, sigmaTotal, sigmaElastic;
180 
181  const G4ParticleDefinition* theParticle = aParticle;
182 
183  sigmaTotal = GetAntiHadronNucleonTotCrSc(theParticle,kinEnergy);
184  sigmaElastic = GetAntiHadronNucleonElCrSc(theParticle,kinEnergy);
185 
186 // calculation of sqr of radius NN-collision
187  fRadiusNN2=sigmaTotal*sigmaTotal*0.1/(8.*sigmaElastic*pi); // fm^2
188 
189 
190 // calculation of effective nuclear radius for Pbar and Nbar interaction (can be changed)
191 
192  //A.R. 29-Jan-2013 : use antiprotons/antineutrons as the default case,
193  // to be used for instance, as first approximation
194  // without validation, for anti-hyperons.
195  //if ( (theParticle == theAProton) || (theParticle == theANeutron) ) {
196  if (A==1)
197  { fInelasticXsc = (sigmaTotal - sigmaElastic) * millibarn;
198  return fInelasticXsc;
199  }
200  fRadiusEff = 1.31*G4Pow::GetInstance()->powA(A, 0.22)+0.9/G4Pow::GetInstance()->powA(A, 1./3.); //fm
201 
202  if( (Z==1) && (A==2) ) fRadiusEff = 3.582; //fm
203  if( (Z==1) && (A==3) ) fRadiusEff = 3.105;
204  if( (Z==2) && (A==3) ) fRadiusEff = 3.105;
205  if( (Z==2) && (A==4) ) fRadiusEff = 2.209;
206  //}
207 
208 //calculation of effective nuclear radius for AntiDeuteron interaction (can be changed)
209 
210  if (theParticle ==theADeuteron)
211 {
212  fRadiusEff = 1.38*G4Pow::GetInstance()->powA(A, 0.21)+1.55/G4Pow::GetInstance()->powA(A, 1./3.);
213 
214  if( (Z==1) && (A==2) ) fRadiusEff = 3.169; //fm
215  if( (Z==1) && (A==3) ) fRadiusEff = 3.066;
216  if( (Z==2) && (A==3) ) fRadiusEff = 3.066;
217  if( (Z==2) && (A==4) ) fRadiusEff = 2.498;
218  }
219 
220 //calculation of effective nuclear radius for AntiHe3 interaction (can be changed)
221 
222  if( (theParticle ==theAHe3) || (theParticle ==theATriton) )
223  {
224  fRadiusEff = 1.34 * G4Pow::GetInstance()->powA(A, 0.21)+1.51/G4Pow::GetInstance()->powA(A, 1./3.);
225 
226  if( (Z==1) && (A==2) ) fRadiusEff = 3.066; //fm
227  if( (Z==1) && (A==3) ) fRadiusEff = 2.973;
228  if( (Z==2) && (A==3) ) fRadiusEff = 2.973;
229  if( (Z==2) && (A==4) ) fRadiusEff = 2.508;
230 
231  }
232 
233 //calculation of effective nuclear radius for AntiAlpha interaction (can be changed)
234 
235  if (theParticle == theAAlpha)
236  {
237  fRadiusEff = 1.3*G4Pow::GetInstance()->powA(A, 0.21)+1.05/G4Pow::GetInstance()->powA(A, 1./3.);
238 
239  if( (Z==1) && (A==2) ) fRadiusEff = 2.498; //fm
240  if( (Z==1) && (A==3) ) fRadiusEff = 2.508;
241  if( (Z==2) && (A==3) ) fRadiusEff = 2.508;
242  if( (Z==2) && (A==4) ) fRadiusEff = 2.158;
243  }
244  G4double R2 = fRadiusEff*fRadiusEff;
245  G4double REf2 = R2+fRadiusNN2;
246  G4double ApAt= std::abs(theParticle->GetBaryonNumber()) * A;
247 
248  inelxsection = pi*REf2 *10* G4Log(1+(ApAt*sigmaTotal/(pi*REf2*10.))); //mb
249  inelxsection = inelxsection * millibarn;
250  fInelasticXsc = inelxsection;
251  return fInelasticXsc;
252 }
253 
255 //
256 // Calculates Inelastic Anti-nucleus-Nucleus cross-section
257 //
259 (const G4ParticleDefinition* aParticle, G4double kinEnergy, G4int Z, G4int A)
260 {return GetInelasticElementCrossSection(aParticle, kinEnergy, Z, (G4double) A); }
261 
262 
263 
265 //
266 // Calculates elastic Anti-nucleus-Nucleus cross-section as Total - Inelastic
267 //
269 (const G4ParticleDefinition* aParticle, G4double kinEnergy, G4int Z, G4double A)
270 {
271  fElasticXsc = GetTotalElementCrossSection(aParticle, kinEnergy, Z, A)-
272  GetInelasticElementCrossSection(aParticle, kinEnergy, Z, A);
273 
274  if (fElasticXsc < 0.) fElasticXsc = 0.;
275 
276  return fElasticXsc;
277 }
278 
280 //
281 // Calculates elastic Anti-nucleus-Nucleus cross-section
282 //
284 (const G4ParticleDefinition* aParticle, G4double kinEnergy, G4int Z, G4int A)
285 { return GetElasticElementCrossSection(aParticle, kinEnergy, Z, (G4double) A); }
286 
287 
289 // Calculation of Antihadron - hadron Total Cross-section
290 
292 (const G4ParticleDefinition* aParticle, G4double kinEnergy)
293 {
294  G4double xsection, Pmass, Energy, momentum;
295  const G4ParticleDefinition* theParticle = aParticle;
296  Pmass=theParticle->GetPDGMass();
297  Energy=Pmass+kinEnergy;
298  momentum=std::sqrt(Energy*Energy-Pmass*Pmass)/std::abs(theParticle->GetBaryonNumber());
299  G4double Plab = momentum / GeV;
300 
301  G4double B, SigAss;
302  G4double C, d1, d2, d3 ;
303 
304  Elab = std::sqrt(Mn*Mn + Plab*Plab); // GeV
305  S = 2.*Mn*Mn + 2. *Mn*Elab; // GeV^2
306  SqrtS = std::sqrt(S); // GeV
307 
308  B = b0+b2*G4Log(SqrtS/SqrtS0)*G4Log(SqrtS/SqrtS0); //GeV^(-2)
309  SigAss = 36.04 +0.304*G4Log(S/S0)*G4Log(S/S0); //mb
310  R0 = std::sqrt(0.40874044*SigAss - B); //GeV^(-2)
311 
312  C = 13.55;
313  d1 = -4.47;
314  d2 = 12.38;
315  d3 = -12.43;
316  xsection = SigAss*(1 + 1./(std::sqrt(S-4.*Mn*Mn)) / (G4Pow::GetInstance()->powA(R0, 3.))
317  *C* (1+d1/SqrtS+d2/(G4Pow::GetInstance()->powA(SqrtS,2.))+d3/(G4Pow::GetInstance()->powA(SqrtS,3.)) ));
318 
319 // xsection *= millibarn;
320 
321  fAntiHadronNucleonTotXsc = xsection;
322  return fAntiHadronNucleonTotXsc;
323 }
324 
325 
326 //
327 // /////////////////////////////////////////////////////////////////////////////////
328 // Calculation of Antihadron - hadron Elastic Cross-section
329 
332 {
333  G4double xsection;
334 
335  G4double SigAss;
336  G4double C, d1, d2, d3 ;
337 
338  GetAntiHadronNucleonTotCrSc(aParticle,kinEnergy);
339 
340  SigAss = 4.5 + 0.101*G4Log(S/S0)*G4Log(S/S0); //mb
341 
342  C = 59.27;
343  d1 = -6.95;
344  d2 = 23.54;
345  d3 = -25.34;
346 
347  xsection = SigAss* (1 + 1. / (std::sqrt(S-4.*Mn*Mn)) / (G4Pow::GetInstance()->powA(R0, 3.))
348  *C* ( 1+d1/SqrtS+d2/(G4Pow::GetInstance()->powA(SqrtS,2.))+d3/(G4Pow::GetInstance()->powA(SqrtS,3.)) ));
349 
350 // xsection *= millibarn;
351 
352  fAntiHadronNucleonElXsc = xsection;
353  return fAntiHadronNucleonElXsc;
354 }
355 
357 {
358  outFile << "The G4ComponentAntiNuclNuclearXS calculates total,\n"
359  << "inelastic, elastic cross sections of anti-nucleons and light \n"
360  << "anti-nucleus interactions with nuclei using Glauber's approach.\n"
361  << "It uses parametrizations of antiproton-proton total and elastic \n"
362  << "cross sections and Wood-Saxon distribution of nuclear density.\n"
363  << "The lower limit is 10 MeV, the upper limit is 10 TeV. \n"
364  << "See details in Phys.Lett. B705 (2011) 235. \n";
365 }
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
static G4AntiHe3 * AntiHe3()
Definition: G4AntiHe3.cc:94
G4double powA(G4double A, G4double y) const
Definition: G4Pow.hh:259
virtual G4double GetElasticElementCrossSection(const G4ParticleDefinition *aParticle, G4double kinEnergy, G4int Z, G4double A)
double S(double temp)
static G4AntiDeuteron * AntiDeuteron()
static const G4double d2
virtual G4double GetInelasticElementCrossSection(const G4ParticleDefinition *aParticle, G4double kinEnergy, G4int Z, G4double A)
G4double GetAntiHadronNucleonElCrSc(const G4ParticleDefinition *aParticle, G4double kinEnergy)
double C(double temp)
double B(double temperature)
static G4AntiAlpha * AntiAlpha()
Definition: G4AntiAlpha.cc:89
int G4int
Definition: G4Types.hh:78
G4double GetAntiHadronNucleonTotCrSc(const G4ParticleDefinition *aParticle, G4double kinEnergy)
virtual void CrossSectionDescription(std::ostream &) const
double A(double temperature)
static G4AntiProton * AntiProton()
Definition: G4AntiProton.cc:93
static const G4double d1
virtual G4double GetTotalElementCrossSection(const G4ParticleDefinition *aParticle, G4double kinEnergy, G4int Z, G4double A)
G4double G4Log(G4double x)
Definition: G4Log.hh:230
G4double GetPDGMass() const
virtual G4double GetTotalIsotopeCrossSection(const G4ParticleDefinition *aParticle, G4double kinEnergy, G4int Z, G4int A)
virtual G4double GetInelasticIsotopeCrossSection(const G4ParticleDefinition *aParticle, G4double kinEnergy, G4int Z, G4int A)
static constexpr double GeV
Definition: G4SIunits.hh:217
static constexpr double pi
Definition: G4SIunits.hh:75
double G4double
Definition: G4Types.hh:76
static G4AntiTriton * AntiTriton()
Definition: G4AntiTriton.cc:94
static G4AntiNeutron * AntiNeutron()
static constexpr double millibarn
Definition: G4SIunits.hh:106
virtual G4double GetElasticIsotopeCrossSection(const G4ParticleDefinition *aParticle, G4double kinEnergy, G4int Z, G4int A)