00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049 #include "G4eeCrossSections.hh"
00050 #include "G4PhysicalConstants.hh"
00051 #include "G4SystemOfUnits.hh"
00052 #include "G4PionPlus.hh"
00053 #include "G4PionMinus.hh"
00054 #include "G4PionZero.hh"
00055 #include "G4Eta.hh"
00056 #include "G4KaonPlus.hh"
00057 #include "G4KaonMinus.hh"
00058 #include "G4KaonZeroLong.hh"
00059 #include "G4PhysicsLinearVector.hh"
00060
00061 #include <iostream>
00062 #include <fstream>
00063
00064
00065
00066 using namespace std;
00067
00068 G4eeCrossSections::G4eeCrossSections()
00069 {
00070 Initialise();
00071 }
00072
00073
00074
00075 G4eeCrossSections::~G4eeCrossSections()
00076 {}
00077
00078
00079
00080 void G4eeCrossSections::Initialise()
00081 {
00082 MsPi = G4PionPlus::PionPlus()->GetPDGMass();
00083 MsPi0= G4PionZero::PionZero()->GetPDGMass();
00084 MsEta= G4Eta::Eta()->GetPDGMass();
00085 MsEtap=957.78*MeV;
00086 MsKs = G4KaonZeroLong::KaonZeroLong()->GetPDGMass();
00087 MsKc = G4KaonPlus::KaonPlus()->GetPDGMass();
00088 MsRho= 775.5*MeV;
00089 MsOm = 782.62*MeV;
00090 MsF0 = 980.0*MeV;
00091 MsA0 = 984.7*MeV;
00092 MsPhi= 1019.46*MeV;
00093 MsK892 = 891.66*MeV;
00094 MsK0892 = 896.0*MeV;
00095 GRho = 149.4*MeV;
00096 GOm = 8.49*MeV;
00097 GPhi = 4.26*MeV;
00098 GK892 = 50.8*MeV;
00099 GK0892 = 50.3*MeV;
00100 PhRho = 0.0;
00101 PhOm = 0.0;
00102 PhPhi = 155.0*degree;
00103 PhRhoPi = 186.0*degree;
00104
00105 BrRhoPiG = 4.5e-4;
00106 BrRhoPi0G= 6.8e-4;
00107 BrRhoEtaG= 2.95e-4;
00108 BrRhoEe = 4.7e-5;
00109 BrOm3Pi = 0.891;
00110 BrOmPi0G= 0.089;
00111 BrOmEtaG= 4.9e-4;
00112 BrOm2Pi = 0.017;
00113 PhOm2Pi = 90.0;
00114 BrOmEe = 7.18e-5;
00115 BrPhi2Kc = 0.492;
00116 BrPhiKsKl= 0.34;
00117 BrPhi3Pi = 0.153;
00118 BrPhiPi0G= 1.25e-3;
00119 BrPhiEtaG= 1.301e-2;
00120 BrPhi2Pi = 7.3e-5;
00121 PhPhi2Pi = -20.0*degree;
00122 BrPhiEe = 2.97e-4;
00123
00124 MsRho3 = MsRho*MsRho*MsRho;
00125 MsOm3 = MsOm*MsOm*MsOm;
00126 MsPhi3 = MsPhi*MsPhi*MsPhi;
00127
00128 MeVnb = 3.8938e+11*nanobarn;
00129 Alpha = fine_structure_const;
00130
00131 AOmRho = 3.0;
00132 ARhoPRho = 0.72;
00133 cterm=0.;
00134 mssig = 600.*MeV;
00135 gsig = 500.*MeV;
00136 brsigpipi = 1.;
00137
00138 msrho1450 = 1459.*MeV;
00139 msrho1700 = 1688.8*MeV;
00140 grho1450 = 171.*MeV;
00141 grho1700 = 161.*MeV;
00142 arhoompi0 = 1.;
00143 arho1450ompi0 = 1.;
00144 arho1700ompi0 = 1.;
00145 phrhoompi0 = 0.;
00146 phrho1450ompi0 = pi;
00147 phrho1700ompi0 = 0.;
00148 aomrhopi0 = 1.;
00149 phomrhopi0 = 0.;
00150 arhopi0pi0g = 0.;
00151 aompi0pi0g = 0.;
00152 phrhopi0pi0g = 0.;
00153 phompi0pi0g = 0.;
00154 brrho1450ompi0 = 0.02;
00155 brrho1450pipi = 0.50;
00156 brrho1700ompi0 = 1.0;
00157 brrho1700pipi = 0.02;
00158 aphirhopi0 = 1.;
00159 phphirhopi0 = pi;
00160 arhosigg = 0.;
00161 phrhosigg = 0.;
00162 aomsigg = 0.;
00163 phomsigg = 0.;
00164
00165 G4String w0, w1, w2;
00166 ph3p = 0;
00167
00168
00169
00170
00171
00172
00173
00174
00175
00176
00177
00178
00179
00180
00181
00182
00183
00184
00185
00186
00187
00188
00189 }
00190
00191
00192
00193 G4double G4eeCrossSections::CrossSection2pi(G4double e)
00194 {
00195 complex<G4double> xr(cos(PhRho),sin(PhRho));
00196 complex<G4double> xo(cos(PhOm2Pi),sin(PhOm2Pi));
00197 complex<G4double> xf(cos(PhPhi2Pi),sin(PhPhi2Pi));
00198
00199 G4double s_inv = e*e;
00200 complex<G4double> drho = DpRho(e);
00201 complex<G4double> dom = DpOm(e);
00202 complex<G4double> dphi = DpPhi(e);
00203
00204 complex<G4double> amp =
00205 sqrt(Width2p(s_inv,MsRho,GRho,1.0,MsPi)*MsRho3*BrRhoEe*GRho)*xr/drho
00206 + sqrt(Width2p(s_inv,MsOm,GOm,BrOm2Pi,MsPi)*MsOm3*BrOmEe*GOm)*xo/dom
00207 + sqrt(Width2p(s_inv,MsPhi,GPhi,BrPhi2Pi,MsPi)*MsPhi3*BrPhiEe*GPhi)*xf/dphi;
00208
00209 G4double cross = 12.0*pi*MeVnb*norm(amp)/(e*s_inv);
00210
00211 return cross;
00212 }
00213
00214
00215
00216 G4double G4eeCrossSections::CrossSection3pi(G4double e)
00217 {
00218 complex<G4double> xf(cos(PhPhi2Pi),sin(PhPhi));
00219
00220 G4double s_inv = e*e;
00221 complex<G4double> dom = DpOm(e);
00222 complex<G4double> dphi = DpPhi(e);
00223
00224 complex<G4double> amp =
00225 sqrt(Width3p(s_inv,MsOm,GOm,BrOm3Pi)*MsOm3*BrOmEe*GOm)/dom
00226 + sqrt(Width3p(s_inv,MsPhi,GPhi,BrPhi3Pi)*MsPhi3*BrPhiEe*GPhi)*xf/dphi;
00227
00228 G4double cross = 12.0*pi*MeVnb*norm(amp)/(e*s_inv);
00229
00230 return cross;
00231 }
00232
00233
00234
00235 G4double G4eeCrossSections::CrossSectionPi0G(G4double e)
00236 {
00237 complex<G4double> xf(cos(PhPhi),sin(PhPhi));
00238
00239 G4double s_inv = e*e;
00240 complex<G4double> drho = DpRho(e);
00241 complex<G4double> dom = DpOm(e);
00242 complex<G4double> dphi = DpPhi(e);
00243
00244 complex<G4double> amp =
00245 sqrt(WidthPg(s_inv,MsRho,GRho,BrRhoPi0G,MsPi0)*MsRho3*BrRhoEe*GRho)/drho
00246 + sqrt(WidthPg(s_inv,MsOm,GOm,BrOmPi0G,MsPi0)*MsOm3*BrOmEe*GOm)/dom
00247 + sqrt(WidthPg(s_inv,MsPhi,GPhi,BrPhiPi0G,MsPi0)*MsPhi3*BrPhiEe*GPhi)*xf/dphi;
00248
00249 G4double cross = 12.0*pi*MeVnb*norm(amp)/(e*s_inv);
00250
00251 return cross;
00252 }
00253
00254
00255
00256 G4double G4eeCrossSections::CrossSectionEtaG(G4double e)
00257 {
00258 complex<G4double> xf(cos(PhPhi),sin(PhPhi));
00259
00260 G4double s_inv = e*e;
00261 complex<G4double> drho = DpRho(e);
00262 complex<G4double> dom = DpOm(e);
00263 complex<G4double> dphi = DpPhi(e);
00264
00265 complex<G4double> amp =
00266 sqrt(WidthPg(s_inv,MsRho,GRho,BrRhoEtaG,MsEta)*MsRho3*BrRhoEe*GRho)/drho
00267 + sqrt(WidthPg(s_inv,MsOm,GOm,BrOmEtaG,MsEta)*MsOm3*BrOmEe*GOm)/dom
00268 + sqrt(WidthPg(s_inv,MsPhi,GPhi,BrPhiEtaG,MsEta)*MsPhi3*BrPhiEe*GPhi)*xf/dphi;
00269
00270 G4double cross = 12.0*pi*MeVnb*norm(amp)/(e*s_inv);
00271
00272 return cross;
00273 }
00274
00275
00276
00277 G4double G4eeCrossSections::CrossSection2Kcharged(G4double e)
00278 {
00279 G4double s_inv = e*e;
00280 complex<G4double> dphi = DpPhi(e);
00281
00282 complex<G4double> amp =
00283 sqrt(Width2p(s_inv,MsPhi,GPhi,BrPhi2Kc,MsKc)*MsPhi3*BrPhiEe*GPhi)/dphi;
00284
00285 G4double cross = 12.0*pi*MeVnb*norm(amp)/(e*s_inv);
00286
00287 return cross;
00288 }
00289
00290
00291
00292 G4double G4eeCrossSections::CrossSection2Kneutral(G4double e)
00293 {
00294 G4double s_inv = e*e;
00295 complex<G4double> dphi = DpPhi(e);
00296
00297 complex<G4double> amp =
00298 sqrt(Width2p(s_inv,MsPhi,GPhi,BrPhiKsKl,MsKs)*MsPhi3*BrPhiEe*GPhi)/dphi;
00299
00300 G4double cross = 12.0*pi*MeVnb*norm(amp)/(e*s_inv);
00301
00302 return cross;
00303 }
00304
00305
00306
00307 G4double G4eeCrossSections::Width2p(G4double s_inv, G4double mres,
00308 G4double gconst, G4double br, G4double mp)
00309 {
00310 G4double mp2 = 4.0*mp*mp;
00311 G4double s0 = mres*mres;
00312 G4double f = (s_inv - mp2)/(s0 - mp2);
00313 if(f < 0.0) f = 0.0;
00314 return gconst*br*sqrt(f)*f*s0/s_inv;
00315 }
00316
00317
00318
00319 G4double G4eeCrossSections::Width3p(G4double s_inv, G4double mres,
00320 G4double gconst, G4double br)
00321 {
00322 G4double w = PhaseSpace3p(sqrt(s_inv));
00323 G4double w0= PhaseSpace3p(mres);
00324 G4double x = gconst*br*w/w0;
00325 return x;
00326 }
00327
00328
00329
00330 G4double G4eeCrossSections::PhaseSpace3p(G4double e)
00331 {
00332
00333
00334
00335
00336
00337
00338 G4double x = 1.0;
00339 G4double emev = e/MeV;
00340 G4double y = 414.12/emev;
00341 x *= pow(e/MsOm, 5.0) * pow(emev*0.1, 3.0)*(1.0 - y*y);
00342 return x;
00343 }
00344
00345
00346
00347 G4double G4eeCrossSections::WidthPg(G4double s_inv, G4double mres,
00348 G4double gconst, G4double br, G4double mp)
00349 {
00350 G4double mp2 = mp*mp;
00351 G4double s0 = mres*mres;
00352 G4double f = (s_inv - mp2)*mres/((s0 - mp2)*sqrt(s_inv));
00353 if(f < 0.0) f = 0.0;
00354 return gconst*br*f*f*f;
00355 }
00356
00357
00358
00359 G4double G4eeCrossSections::WidthRho(G4double e)
00360 {
00361 G4double w = Width2p(e*e, MsRho, GRho, 1.0, MsPi);
00362 return w;
00363 }
00364
00365
00366
00367 G4double G4eeCrossSections::WidthOm(G4double e)
00368 {
00369 G4double s_inv = e*e;
00370 G4double w = (Width3p(s_inv, MsOm, GOm, BrOm3Pi) +
00371 WidthPg(s_inv, MsOm, GOm, BrOmPi0G, MsPi0) +
00372 WidthPg(s_inv, MsOm, GOm, BrOmEtaG, MsEta) +
00373 Width2p(s_inv, MsOm, GOm, BrOm2Pi, MsPi)) /
00374 (BrOm3Pi+BrOmPi0G+BrOmEtaG+BrOm2Pi);
00375 return w;
00376 }
00377
00378
00379
00380 G4double G4eeCrossSections::WidthPhi(G4double e)
00381 {
00382 G4double s_inv = e*e;
00383 G4double w = (Width3p(s_inv, MsPhi, GPhi, BrPhi3Pi) +
00384 WidthPg(s_inv, MsPhi, GPhi, BrPhiPi0G, MsPi0) +
00385 WidthPg(s_inv, MsPhi, GPhi, BrPhiEtaG, MsEta) +
00386 Width2p(s_inv, MsPhi, GPhi, BrPhi2Kc, MsKc) +
00387 Width2p(s_inv, MsPhi, GPhi, BrPhiKsKl, MsKs)) /
00388 (BrPhi3Pi+BrPhiPi0G+BrPhiEtaG+BrPhi2Kc+BrPhiKsKl);
00389 return w;
00390 }
00391
00392
00393
00394 complex<G4double> G4eeCrossSections::DpRho(G4double e)
00395 {
00396 complex<G4double> d(MsRho*MsRho - e*e, -e*WidthRho(e));
00397 return d;
00398 }
00399
00400
00401
00402 complex<G4double> G4eeCrossSections::DpOm(G4double e)
00403 {
00404 complex<G4double> d(MsOm*MsOm - e*e, -e*WidthOm(e));
00405 return d;
00406 }
00407
00408
00409
00410 complex<G4double> G4eeCrossSections::DpPhi(G4double e)
00411 {
00412 complex<G4double> d(MsPhi*MsPhi - e*e, -e*WidthPhi(e));
00413 return d;
00414 }
00415
00416