Geant4  9.6.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4PolarizedMollerCrossSection.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 // $Id: G4PolarizedMollerCrossSection.cc 69847 2013-05-16 09:36:18Z gcosmo $
27 // -------------------------------------------------------------------
28 //
29 // GEANT4 Class file
30 //
31 //
32 // File name: G4PolarizedMollerCrossSection
33 //
34 // Author: Andreas Schaelicke
35 //
36 // Creation date: 12.01.2006
37 //
38 // Modifications:
39 // 16-01-06 included cross section as calculated by P.Starovoitov
40 //
41 // Class Description:
42 // * calculates the differential cross section
43 // incomming electron K1(along positive z direction) scatters at an electron K2 at rest
44 // * phi denotes the angle between the scattering plane (defined by the
45 // outgoing electron) and X-axis
46 // * all stokes vectors refer to spins in the Global System (X,Y,Z)
47 //
48 
49 #include "G4PolarizedMollerCrossSection.hh"
50 #include "G4PhysicalConstants.hh"
51 
52 G4PolarizedMollerCrossSection::G4PolarizedMollerCrossSection() :
53  phi0(0.)
54 {
55  SetXmax(.5);
56 }
57 G4PolarizedMollerCrossSection::~G4PolarizedMollerCrossSection() {}
58 void G4PolarizedMollerCrossSection::Initialize(
59  G4double e,
60  G4double gamma,
61  G4double /*phi*/,
62  const G4StokesVector & pol0,
63  const G4StokesVector & pol1,
64  G4int flag)
65 {
66  G4double re2 = classic_electr_radius * classic_electr_radius;
67  G4double gamma2=gamma*gamma;
68  G4double gmo = (gamma - 1.);
69  G4double gmo2 = (gamma - 1.)*(gamma - 1.);
70  G4double gpo = (gamma + 1.);
71  G4double pref = gamma2*re2/(gmo2*(gamma + 1.0));
72  G4double sqrttwo=std::sqrt(2.);
73  G4double f = (-1. + e);
74  G4double e2 = e*e;
75  G4double f2 = f*f;
76  // G4double w = e*(1. - e);
77 
78  G4bool polarized=(!pol0.IsZero())||(!pol1.IsZero());
79 
80  if (flag==0) polarized=false;
81  // Unpolarised part of XS
82  phi0 = 0.;
83  phi0+= gmo2/gamma2;
84  phi0+= ((1. - 2.*gamma)/gamma2)*(1./e + 1./(1.-e));
85  phi0+= 1./(e*e) + 1./((1. - e)*(1. - e));
86  phi0*=0.25;
87  // Initial state polarisarion dependence
88  if (polarized) {
89  G4double usephi=1.;
90  if (flag<=1) usephi=0.;
91  // G4cout<<"Polarized differential moller cross section"<<G4endl;
92  // G4cout<<"Initial state polarisation contributions"<<G4endl;
93  // G4cout<<"Diagonal Matrix Elements"<<G4endl;
94  G4double xx = (gamma - f*e*gmo*(3 + gamma))/(4*f*e*gamma2);
95  G4double yy = (-1 + f*e*gmo2 + 2*gamma)/(4*f*e*gamma2);
96  G4double zz = (-(e*gmo*(3 + gamma)) + e2*gmo*(3 + gamma) +
97  gamma*(-1 + 2*gamma))/(4*f*e*gamma2);
98 
99  phi0 += xx*pol0.x()*pol1.x() + yy*pol0.y()*pol1.y() + zz*pol0.z()*pol1.z();
100 
101  if (usephi==1.) {
102  // G4cout<<"Non-diagonal Matrix Elements"<<G4endl;
103  G4double xy = 0;
104  G4double xz = -((-1 + 2*e)*gmo)/(2*sqrttwo*gamma2*
105  std::sqrt(-((f*e)/gpo)));
106  G4double yx = 0;
107  G4double yz = 0;
108  G4double zx = -((-1 + 2*e)*gmo)/(2*sqrttwo*gamma2*
109  std::sqrt(-((f*e)/gpo)));
110  G4double zy = 0;
111  phi0+=yx*pol0.y()*pol1.x() + xy*pol0.x()*pol1.y();
112  phi0+=zx*pol0.z()*pol1.x() + xz*pol0.x()*pol1.z();
113  phi0+=zy*pol0.z()*pol1.y() + yz*pol0.y()*pol1.z();
114  }
115  }
116  // Final state polarisarion dependence
117  phi2=G4ThreeVector();
118  phi3=G4ThreeVector();
119 
120  if (flag>=1) {
121  //
122  // Final Electron P1
123  //
124 
125  // initial electron K1
126  if (!pol0.IsZero()) {
127  G4double xxP1K1 = (std::sqrt(gpo/(1 + e2*gmo + gamma - 2*e*gamma))*
128  (gamma - e*gpo))/(4*e2*gamma);
129  G4double xyP1K1 = 0;
130  G4double xzP1K1 = (-1 + 2*e*gamma)/(2*sqrttwo*f*gamma*
131  std::sqrt(e*e2*(1 + e + gamma - e*gamma)));
132  G4double yxP1K1 = 0;
133  G4double yyP1K1 = (-gamma2 + e*(-1 + gamma*(2 + gamma)))/(4*f*e2*gamma2);
134  G4double yzP1K1 = 0;
135  G4double zxP1K1 = (1 + 2*e2*gmo - 2*e*gamma)/(2*sqrttwo*f*e*gamma*
136  std::sqrt(e*(1 + e + gamma - e*gamma)));
137  G4double zyP1K1 = 0;
138  G4double zzP1K1 = (-gamma + e*(1 - 2*e*gmo + gamma))/(4*f*e2*gamma*
139  std::sqrt(1 - (2*e)/(f*gpo)));
140  phi2[0] += xxP1K1*pol0.x() + xyP1K1*pol0.y() + xzP1K1*pol0.z();
141  phi2[1] += yxP1K1*pol0.x() + yyP1K1*pol0.y() + yzP1K1*pol0.z();
142  phi2[2] += zxP1K1*pol0.x() + zyP1K1*pol0.y() + zzP1K1*pol0.z();
143  }
144  // initial electron K2
145  if (!pol1.IsZero()) {
146  G4double xxP1K2 = ((1 + e*(-3 + gamma))*std::sqrt(gpo/(1 + e2*gmo + gamma -
147  2*e*gamma)))/(4*f*e*gamma);
148  G4double xyP1K2 = 0;
149  G4double xzP1K2 = (-2 + 2*e + gamma)/(2*sqrttwo*f2*gamma*
150  std::sqrt(e*(1 + e + gamma - e*gamma)));
151  G4double yxP1K2 = 0;
152  G4double yyP1K2 = (1 - 2*gamma + e*(-1 + gamma*(2 + gamma)))/(4*f2*e*gamma2);
153  G4double yzP1K2 = 0;
154  G4double zxP1K2 = (2*e*(1 + e*gmo - 2*gamma) + gamma)/(2*sqrttwo*f2*gamma*
155  std::sqrt(e*(1 + e + gamma - e*gamma)));
156  G4double zyP1K2 = 0;
157  G4double zzP1K2 = (1 - 2*gamma + e*(-1 - 2*e*gmo + 3*gamma))/
158  (4*f2*e*gamma*std::sqrt(1 - (2*e)/(f*gpo)));
159  phi2[0] += xxP1K2*pol1.x() + xyP1K2*pol1.y() + xzP1K2*pol1.z();
160  phi2[1] += yxP1K2*pol1.x() + yyP1K2*pol1.y() + yzP1K2*pol1.z();
161  phi2[2] += zxP1K2*pol1.x() + zyP1K2*pol1.y() + zzP1K2*pol1.z();
162  }
163  //
164  // Final Electron P2
165  //
166 
167  // initial electron K1
168  if (!pol0.IsZero()) {
169 
170 
171  G4double xxP2K1 = (-1 + e + e*gamma)/(4*f2*gamma*
172  std::sqrt((e*(2 + e*gmo))/gpo));
173  G4double xyP2K1 = 0;
174  G4double xzP2K1 = -((1 + 2*f*gamma)*std::sqrt(f/(-2 + e - e*gamma)))/
175  (2*sqrttwo*f2*e*gamma);
176  G4double yxP2K1 = 0;
177  G4double yyP2K1 = (1 - 2*gamma + e*(-1 + gamma*(2 + gamma)))/(4*f2*e*gamma2);
178  G4double yzP2K1 = 0;
179  G4double zxP2K1 = (1 + 2*e*(-2 + e + gamma - e*gamma))/(2*sqrttwo*f*e*
180  std::sqrt(-(f*(2 + e*gmo)))*gamma);
181  G4double zyP2K1 = 0;
182  G4double zzP2K1 = (std::sqrt((e*gpo)/(2 + e*gmo))*
183  (-3 + e*(5 + 2*e*gmo - 3*gamma) + 2*gamma))/(4*f2*e*gamma);
184 
185  phi3[0] += xxP2K1*pol0.x() + xyP2K1*pol0.y() + xzP2K1*pol0.z();
186  phi3[1] += yxP2K1*pol0.x() + yyP2K1*pol0.y() + yzP2K1*pol0.z();
187  phi3[2] += zxP2K1*pol0.x() + zyP2K1*pol0.y() + zzP2K1*pol0.z();
188  }
189  // initial electron K2
190  if (!pol1.IsZero()) {
191 
192  G4double xxP2K2 = (-2 - e*(-3 + gamma) + gamma)/
193  (4*f*e*gamma* std::sqrt((e*(2 + e*gmo))/gpo));
194  G4double xyP2K2 = 0;
195  G4double xzP2K2 = ((-2*e + gamma)*std::sqrt(f/(-2 + e - e*gamma)))/
196  (2*sqrttwo*f*e2*gamma);
197  G4double yxP2K2 = 0;
198  G4double yyP2K2 = (-gamma2 + e*(-1 + gamma*(2 + gamma)))/(4*f*e2*gamma2);
199  G4double yzP2K2 = 0;
200  G4double zxP2K2 = (gamma + 2*e*(-1 + e - e*gamma))/
201  (2*sqrttwo*e2* std::sqrt(-(f*(2 + e*gmo)))*gamma);
202  G4double zyP2K2 = 0;
203  G4double zzP2K2 = (std::sqrt((e*gpo)/(2 + e*gmo))*
204  (-2 + e*(3 + 2*e*gmo - gamma) + gamma))/(4*f*e2*gamma);
205  phi3[0] += xxP2K2*pol1.x() + xyP2K2*pol1.y() + xzP2K2*pol1.z();
206  phi3[1] += yxP2K2*pol1.x() + yyP2K2*pol1.y() + yzP2K2*pol1.z();
207  phi3[2] += zxP2K2*pol1.x() + zyP2K2*pol1.y() + zzP2K2*pol1.z();
208  }
209  }
210  phi0 *= pref;
211  phi2 *= pref;
212  phi3 *= pref;
213 }
214 
215 G4double G4PolarizedMollerCrossSection::XSection(const G4StokesVector & pol2,
216  const G4StokesVector & pol3)
217 {
218  G4double xs=0.;
219  xs+=phi0;
220 
221  G4bool polarized=(!pol2.IsZero())||(!pol3.IsZero());
222  if (polarized) {
223  xs+=phi2*pol2 + phi3*pol3;
224  }
225  return xs;
226 }
227 
228 G4double G4PolarizedMollerCrossSection::TotalXSection(
229  G4double xmin, G4double xmax, G4double gamma,
230  const G4StokesVector & pol0,const G4StokesVector & pol1)
231 {
232  G4double xs=0.;
233 
234  G4double x=xmin;
235 
236  if (xmax != 1./2.) G4cout<<" warning xmax expected to be 1/2 but is "<<xmax<< G4endl;
237 
238  // re -> electron radius^2;
239  G4double re2 = classic_electr_radius * classic_electr_radius;
240  G4double gamma2=gamma*gamma;
241  G4double gmo2 = (gamma - 1.)*(gamma - 1.);
242  G4double logMEM = std::log(1./x - 1.);
243  G4double pref = twopi*gamma2*re2/(gmo2*(gamma + 1.0));
244  // unpolarise XS
245  G4double sigma0 = 0.;
246  sigma0 += (gmo2/gamma2)*(0.5 - x);
247  sigma0 += ((1. - 2.*gamma)/gamma2)*logMEM;
248  sigma0 += 1./x - 1./(1. - x);
249  // longitudinal part
250  G4double sigma2=0.;
251  sigma2 += ((gamma2 + 2.*gamma - 3.)/gamma2)*(0.5 - x);
252  sigma2 += (1./gamma - 2.)*logMEM;
253  // transverse part
254  G4double sigma3=0.;
255  sigma3 += (2.*(1. - gamma)/gamma2)*(0.5 - x);
256  sigma3 += (1. - 3.*gamma)/(2.*gamma2)*logMEM;
257  // total cross section
258  xs+=pref*(sigma0 + sigma2*pol0.z()*pol1.z() + sigma3*(pol0.x()*pol1.x()+pol0.y()*pol1.y()));
259 
260  return xs;
261 }
262 
263 
264 G4StokesVector G4PolarizedMollerCrossSection::GetPol2()
265 {
266  // Note, mean polarization can not contain correlation
267  // effects.
268  return 1./phi0 * phi2;
269 }
270 G4StokesVector G4PolarizedMollerCrossSection::GetPol3()
271 {
272  // Note, mean polarization can not contain correlation
273  // effects.
274  return 1./phi0 * phi3;
275 }
Generated on Sat May 25 2013 14:33:36 for Geant4 by   doxygen 1.8.4