Geant4  9.6.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 
46 
48 : G4VComponentCrossSection("AntiAGlauber"), fUpperLimit(10000*GeV),
49  fLowerLimit(10*MeV), fRadiusEff(0.0), fRadiusNN2(0.0),
50  fTotalXsc(0.0), fElasticXsc(0.0), fInelasticXsc(0.0),
51  fAntiHadronNucleonTotXsc(0.0), fAntiHadronNucleonElXsc(0.0),
52  Elab(0.0), S(0.0), SqrtS(0)
53 {
54  theAProton = G4AntiProton::AntiProton();
55  theANeutron = G4AntiNeutron::AntiNeutron();
56  theADeuteron = G4AntiDeuteron::AntiDeuteron();
57  theATriton = G4AntiTriton::AntiTriton();
58  theAAlpha = G4AntiAlpha::AntiAlpha();
59  theAHe3 = G4AntiHe3::AntiHe3();
60 
61  Mn = 0.93827231; // GeV
62  b0 = 11.92; // GeV^(-2)
63  b2 = 0.3036; // GeV^(-2)
64  SqrtS0 = 20.74; // GeV
65  S0 = 33.0625; // GeV^2
66 
67 }
68 
70 //
71 //
72 
74 {
75 }
76 
78 
79 
80 
82 //
83 // Calculation of total CrossSection of Anti-Nucleus - Nucleus
84 
85 
87 (const G4ParticleDefinition* aParticle, G4double kinEnergy, G4int Z, G4double A)
88 {
89  G4double xsection, sigmaTotal, sigmaElastic;
90 
91  const G4ParticleDefinition* theParticle = aParticle;
92 
93  sigmaTotal = GetAntiHadronNucleonTotCrSc(theParticle,kinEnergy);
94  sigmaElastic = GetAntiHadronNucleonElCrSc(theParticle,kinEnergy);
95 
96 // calculation of squared radius of NN-collision
97  fRadiusNN2=sigmaTotal*sigmaTotal*0.1/(8.*sigmaElastic*pi) ; //fm^2
98 
99 // calculation of effective nuclear radius for Pbar and Nbar interactions (can be changed)
100 
101  //A.R. 29-Jan-2013 : use antiprotons/antineutrons as the default case,
102  // to be used for instance, as first approximation
103  // without validation, for anti-hyperons.
104  //if ( (theParticle == theAProton) || (theParticle == theANeutron) ) {
105  if(A==1)
106  { fTotalXsc = sigmaTotal * millibarn;
107  return fTotalXsc; }
108 
109  fRadiusEff = 1.34*std::pow(A,0.23)+1.35/std::pow(A,1./3.); //fm
110 
111  if( (Z==1) && (A==2) ) fRadiusEff = 3.800; //fm
112  if( (Z==1) && (A==3) ) fRadiusEff = 3.300;
113  if( (Z==2) && (A==3) ) fRadiusEff = 3.300;
114  if( (Z==2) && (A==4) ) fRadiusEff = 2.376;
115  //}
116 
117 //calculation of effective nuclear radius for AntiDeuteron interaction (can be changed)
118  if (theParticle == theADeuteron)
119  { fRadiusEff = 1.46 * std::pow(A,0.21) + 1.45 / std::pow(A,1./3.);
120 
121  if( (Z==1) && (A==2) ) fRadiusEff = 3.238; //fm
122  if( (Z==1) && (A==3) ) fRadiusEff = 3.144;
123  if( (Z==2) && (A==3) ) fRadiusEff = 3.144;
124  if( (Z==2) && (A==4) ) fRadiusEff = 2.544;
125  }
126 // calculation of effective nuclear radius for AntiHe3 interaction (can be changed)
127 
128  if( (theParticle ==theAHe3) || (theParticle ==theATriton) )
129  { fRadiusEff = 1.40* std::pow(A,0.21)+1.63/std::pow(A,1./3.);
130 
131  if( (Z==1) && (A==2) ) fRadiusEff = 3.144; //fm
132  if( (Z==1) && (A==3) ) fRadiusEff = 3.075;
133  if( (Z==2) && (A==3) ) fRadiusEff = 3.075;
134  if( (Z==2) && (A==4) ) fRadiusEff = 2.589;
135  }
136 
137 //calculation of effective nuclear radius for AntiAlpha interaction (can be changed)
138 
139  if (theParticle == theAAlpha)
140  {
141  fRadiusEff = 1.35* std::pow(A,0.21)+1.1/std::pow(A,1./3.);
142 
143  if( (Z==1) && (A==2) ) fRadiusEff = 2.544; //fm
144  if( (Z==1) && (A==3) ) fRadiusEff = 2.589;
145  if( (Z==2) && (A==3) ) fRadiusEff = 2.589;
146  if( (Z==2) && (A==4) ) fRadiusEff = 2.241;
147 
148  }
149 
150  G4double R2 = fRadiusEff*fRadiusEff;
151  G4double REf2 = R2+fRadiusNN2;
152  G4double ApAt = std::abs(theParticle->GetBaryonNumber()) * A;
153 
154  xsection = 2*pi*REf2*10.*std::log(1+(ApAt*sigmaTotal/(2*pi*REf2*10.))); //mb
155  xsection =xsection *millibarn;
156  fTotalXsc = xsection;
157 
158  return fTotalXsc;
159 }
160 
161 
163 //
164 // Calculation of total CrossSection of Anti-Nucleus - Nucleus
167 (const G4ParticleDefinition* aParticle, G4double kinEnergy, G4int Z, G4int A )
168 { return GetTotalElementCrossSection(aParticle, kinEnergy, Z, (G4double) A); }
169 
171 // Calculation of inelastic CrossSection of Anti-Nucleus - Nucleus
173 
175 (const G4ParticleDefinition* aParticle, G4double kinEnergy, G4int Z, G4double A)
176 {
177  G4double inelxsection, sigmaTotal, sigmaElastic;
178 
179  const G4ParticleDefinition* theParticle = aParticle;
180 
181  sigmaTotal = GetAntiHadronNucleonTotCrSc(theParticle,kinEnergy);
182  sigmaElastic = GetAntiHadronNucleonElCrSc(theParticle,kinEnergy);
183 
184 // calculation of sqr of radius NN-collision
185  fRadiusNN2=sigmaTotal*sigmaTotal*0.1/(8.*sigmaElastic*pi); // fm^2
186 
187 
188 // calculation of effective nuclear radius for Pbar and Nbar interaction (can be changed)
189 
190  //A.R. 29-Jan-2013 : use antiprotons/antineutrons as the default case,
191  // to be used for instance, as first approximation
192  // without validation, for anti-hyperons.
193  //if ( (theParticle == theAProton) || (theParticle == theANeutron) ) {
194  if (A==1)
195  { fInelasticXsc = (sigmaTotal - sigmaElastic) * millibarn;
196  return fInelasticXsc;
197  }
198  fRadiusEff = 1.31*std::pow(A, 0.22)+0.9/std::pow(A, 1./3.); //fm
199 
200  if( (Z==1) && (A==2) ) fRadiusEff = 3.582; //fm
201  if( (Z==1) && (A==3) ) fRadiusEff = 3.105;
202  if( (Z==2) && (A==3) ) fRadiusEff = 3.105;
203  if( (Z==2) && (A==4) ) fRadiusEff = 2.209;
204  //}
205 
206 //calculation of effective nuclear radius for AntiDeuteron interaction (can be changed)
207 
208  if (theParticle ==theADeuteron)
209 {
210  fRadiusEff = 1.38*std::pow(A, 0.21)+1.55/std::pow(A, 1./3.);
211 
212  if( (Z==1) && (A==2) ) fRadiusEff = 3.169; //fm
213  if( (Z==1) && (A==3) ) fRadiusEff = 3.066;
214  if( (Z==2) && (A==3) ) fRadiusEff = 3.066;
215  if( (Z==2) && (A==4) ) fRadiusEff = 2.498;
216  }
217 
218 //calculation of effective nuclear radius for AntiHe3 interaction (can be changed)
219 
220  if( (theParticle ==theAHe3) || (theParticle ==theATriton) )
221  {
222  fRadiusEff = 1.34 * std::pow(A, 0.21)+1.51/std::pow(A, 1./3.);
223 
224  if( (Z==1) && (A==2) ) fRadiusEff = 3.066; //fm
225  if( (Z==1) && (A==3) ) fRadiusEff = 2.973;
226  if( (Z==2) && (A==3) ) fRadiusEff = 2.973;
227  if( (Z==2) && (A==4) ) fRadiusEff = 2.508;
228 
229  }
230 
231 //calculation of effective nuclear radius for AntiAlpha interaction (can be changed)
232 
233  if (theParticle == theAAlpha)
234  {
235  fRadiusEff = 1.3*std::pow(A, 0.21)+1.05/std::pow(A, 1./3.);
236 
237  if( (Z==1) && (A==2) ) fRadiusEff = 2.498; //fm
238  if( (Z==1) && (A==3) ) fRadiusEff = 2.508;
239  if( (Z==2) && (A==3) ) fRadiusEff = 2.508;
240  if( (Z==2) && (A==4) ) fRadiusEff = 2.158;
241  }
242  G4double R2 = fRadiusEff*fRadiusEff;
243  G4double REf2 = R2+fRadiusNN2;
244  G4double ApAt= std::abs(theParticle->GetBaryonNumber()) * A;
245 
246  inelxsection = pi*REf2 *10* std::log(1+(ApAt*sigmaTotal/(pi*REf2*10.))); //mb
247  inelxsection = inelxsection * millibarn;
248  fInelasticXsc = inelxsection;
249  return fInelasticXsc;
250 }
251 
253 //
254 // Calculates Inelastic Anti-nucleus-Nucleus cross-section
255 //
257 (const G4ParticleDefinition* aParticle, G4double kinEnergy, G4int Z, G4int A)
258 {return GetInelasticElementCrossSection(aParticle, kinEnergy, Z, (G4double) A); }
259 
260 
261 
263 //
264 // Calculates elastic Anti-nucleus-Nucleus cross-section as Total - Inelastic
265 //
267 (const G4ParticleDefinition* aParticle, G4double kinEnergy, G4int Z, G4double A)
268 {
269  fElasticXsc = GetTotalElementCrossSection(aParticle, kinEnergy, Z, A)-
270  GetInelasticElementCrossSection(aParticle, kinEnergy, Z, A);
271 
272  if (fElasticXsc < 0.) fElasticXsc = 0.;
273 
274  return fElasticXsc;
275 }
276 
278 //
279 // Calculates elastic Anti-nucleus-Nucleus cross-section
280 //
282 (const G4ParticleDefinition* aParticle, G4double kinEnergy, G4int Z, G4int A)
283 { return GetElasticElementCrossSection(aParticle, kinEnergy, Z, (G4double) A); }
284 
285 
287 // Calculation of Antihadron - hadron Total Cross-section
288 
290 (const G4ParticleDefinition* aParticle, G4double kinEnergy)
291 {
292  G4double xsection, Pmass, Energy, momentum;
293  const G4ParticleDefinition* theParticle = aParticle;
294  Pmass=theParticle->GetPDGMass();
295  Energy=Pmass+kinEnergy;
296  momentum=std::sqrt(Energy*Energy-Pmass*Pmass)/std::abs(theParticle->GetBaryonNumber());
297  G4double Plab = momentum / GeV;
298 
299  G4double B, SigAss;
300  G4double C, d1, d2, d3 ;
301 
302  Elab = std::sqrt(Mn*Mn + Plab*Plab); // GeV
303  S = 2.*Mn*Mn + 2. *Mn*Elab; // GeV^2
304  SqrtS = std::sqrt(S); // GeV
305 
306  B = b0+b2*std::log(SqrtS/SqrtS0)*std::log(SqrtS/SqrtS0); //GeV^(-2)
307  SigAss = 36.04 +0.304*std::log(S/S0)*std::log(S/S0); //mb
308  R0 = std::sqrt(0.40874044*SigAss - B); //GeV^(-2)
309 
310  C = 13.55;
311  d1 = -4.47;
312  d2 = 12.38;
313  d3 = -12.43;
314  xsection = SigAss*(1 + 1./(std::sqrt(S-4.*Mn*Mn)) / (std::pow(R0, 3.))
315  *C* (1+d1/SqrtS+d2/(std::pow(SqrtS,2.))+d3/(std::pow(SqrtS,3.)) ));
316 
317 // xsection *= millibarn;
318 
319  fAntiHadronNucleonTotXsc = xsection;
320  return fAntiHadronNucleonTotXsc;
321 }
322 
323 
324 //
325 // /////////////////////////////////////////////////////////////////////////////////
326 // Calculation of Antihadron - hadron Elastic Cross-section
327 
330 {
331  G4double xsection;
332 
333  G4double SigAss;
334  G4double C, d1, d2, d3 ;
335 
336  GetAntiHadronNucleonTotCrSc(aParticle,kinEnergy);
337 
338  SigAss = 4.5 + 0.101*std::log(S/S0)*std::log(S/S0); //mb
339 
340  C = 59.27;
341  d1 = -6.95;
342  d2 = 23.54;
343  d3 = -25.34;
344 
345  xsection = SigAss* (1 + 1. / (std::sqrt(S-4.*Mn*Mn)) / (std::pow(R0, 3.))
346  *C* ( 1+d1/SqrtS+d2/(std::pow(SqrtS,2.))+d3/(std::pow(SqrtS,3.)) ));
347 
348 // xsection *= millibarn;
349 
350  fAntiHadronNucleonElXsc = xsection;
351  return fAntiHadronNucleonElXsc;
352 }
353 
355 {
356  outFile << "The G4ComponentAntiNuclNuclearXS calculates total,\n"
357  << "inelastic, elastic cross sections of anti-nucleons and light \n"
358  << "anti-nucleus interactions with nuclei using Glauber's approach.\n"
359  << "It uses parametrizations of antiproton-proton total and elastic \n"
360  << "cross sections and Wood-Saxon distribution of nuclear density.\n"
361  << "The lower limit is 10 MeV, the upper limit is 10 TeV. \n"
362  << "See details in Phys.Lett. B705 (2011) 235. \n";
363 }