ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4PolarizedComptonCrossSection.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4PolarizedComptonCrossSection.cc
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 //
27 // GEANT4 Class file
28 //
29 //
30 // File name: G4PolarizedComptonCrossSection
31 //
32 // Author: Andreas Schaelicke
33 //
34 // Creation date: 15.05.2005
35 //
36 // Modifications:
37 //
38 // Class Description:
39 // determine the polarization of the final state
40 // in a Compton scattering process employing the differential
41 // cross section by F.W.Lipps & H.A.Tolhoek
42 // ( Physica 20 (1954) 395 )
43 // recalculated by P.Starovoitov
44 //
45 #include "G4PolarizedComptonCrossSection.hh"
46 #include "G4PhysicalConstants.hh"
47 #include "Randomize.hh"
48 
49 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
50 G4PolarizedComptonCrossSection::G4PolarizedComptonCrossSection()
51  : gammaPol2(false), electronPol3(false)
52 {
53  SetYmin(0.);
54 
55  // G4cout<<"G4PolarizedComptonCrossSection() init\n";
56 
57  re2 = classic_electr_radius * classic_electr_radius * sqr(4*pi/hbarc);
58  // G4double unit_conversion = hbarc_squared ;
59  // G4cout<<" (keV)^2* m^2 ="<<unit_conversion<<"\n";
60  phi0 = 0.; polXS = 0.; unpXS = 0.;
61  phi2 = G4ThreeVector(0., 0., 0.);
62  phi3 = G4ThreeVector(0., 0., 0.);
63  polxx = polyy = polzz = polxz = polzx = polyz = polzy = polxy = polyx = 0.;
64  diffXSFactor = 1.;
65  totalXSFactor = 1.;
66 }
67 
68 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
69 G4PolarizedComptonCrossSection::~G4PolarizedComptonCrossSection()
70 {}
71 
72 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
73 void G4PolarizedComptonCrossSection::Initialize(G4double eps, G4double X, G4double , // phi
74  const G4StokesVector & pol0,
75  const G4StokesVector & pol1,
76  G4int flag)
77 {
78  G4double cosT = 1. - (1./eps - 1.)/X;
79  if(cosT > 1.+1.e-8) cosT = 1.;
80  if(cosT < -1.-1.e-8) cosT = -1.;
81  G4double cosT2 = cosT*cosT;
82  G4double cosT3 = cosT2*cosT;
83  G4double sinT2 = 1. - cosT2;
84  if(sinT2 > 1. + 1.e-8) sinT2 = 1.;
85  if(sinT2 < 0.) sinT2 = 0.;
86  G4double sinT = std::sqrt(sinT2);
87  G4double cos2T = 2.*cosT2 - 1.;
88  G4double sin2T = 2.*sinT*cosT;
89  G4double eps2 = sqr(eps);
90  DefineCoefficients(pol0,pol1);
91  diffXSFactor = re2/(4.*X);
92 
93  // unpolarized Cross Section
94  unpXS = (eps2 + 1. - eps*sinT2)/(2.*eps);
95  // initial polarization dependence
96  polXS = -sinT2*pol0.x() + (1. - eps)*sinT*polzx + ((eps2 - 1.)/eps)*cosT*polzz;
97  polXS *= 0.5;
98 
99  phi0 = unpXS + polXS;
100 
101  if (flag == 2 ){
102  // polarization of outgoing photon
103  G4double PHI21 = -sinT2 + 0.5*(cos2T + 3.)*pol0.x() - ((1. - eps)/eps)*sinT*polzx;
104  PHI21 *= 0.5;
105  G4double PHI22 = cosT*pol0.y() + ((1. - eps)/(2.*eps))*sinT*polzy;
106  G4double PHI23 = ((eps2 + 1.)/eps)*cosT*pol0.z() - ((1. - eps)/eps)*(eps*cosT2 + 1.)*pol1.z();
107  PHI23 += 0.5*(1. - eps)*sin2T*pol1.x();
108  PHI23 += (eps - 1.)*(-sinT2*polxz + sinT*polyy - 0.5*sin2T*polxx);
109  PHI23 *= 0.5;
110  phi2 = G4ThreeVector(PHI21, PHI22, PHI23);
111 
112  // polarization of outgoing electron
113  G4double PHI32 = -sinT2*polxy + ((1. - eps)/eps)*sinT*polyz + 0.5*(cos2T + 3.)*pol1.y();
114  PHI32 *= 0.5;
115 
116  G4double PHI31 = 0., PHI31add = 0., PHI33 = 0., PHI33add = 0.;
117 
118  if ((1. - eps) > 1.e-12){
119  G4double helpVar = std::sqrt(eps2 - 2.*cosT*eps + 1.);
120 
121  PHI31 = (1. - eps)*(1. + cosT)*sinT*pol0.z();
122  PHI31 += (-eps*cosT3 + eps*cosT2 + (eps - 2.)*cosT + eps)*pol1.x();
123  PHI31 += -(eps*cosT2 - eps*cosT + cosT + 1.)*sinT*pol1.z();
124  PHI31 /= 2.*helpVar;
125 
126  PHI31add = -eps*sqr(1. - cosT)*(1. + cosT)*polxx;
127  PHI31add += (1. - eps)*sinT2*polyy;
128  PHI31add += -(-eps2 + cosT*(cosT*eps - eps + 1.)*eps + eps - 1.)*sinT*polxz/eps;
129  PHI31add /= 2.*helpVar;
130 
131  PHI33 = ((1. - eps)/eps)*(-eps*cosT2 + eps*(eps + 1.)*cosT - 1.)*pol0.z();
132  PHI33 += -(eps*cosT2 + (1. - eps)*eps*cosT + 1.)*sinT*pol1.x();
133  PHI33 += -(-eps2*cosT3 + eps*(eps2 - eps + 1.)*cosT2 - cosT + eps2)*pol1.z()/eps;
134  PHI33 /= -2.*helpVar;
135 
136  PHI33add = (eps*(eps - cosT - 1.)*cosT + 1.)*sinT*polxx;
137  PHI33add += -(-eps2 + cosT*eps + eps - 1.)*sinT2*polxz;
138  PHI33add += (eps - 1.)*(cosT - eps)*sinT*polyy;
139  PHI33add /= -2.*helpVar;
140  }else{
141  PHI31 = -pol1.z() - (X - 1.)*std::sqrt(1. - eps)*pol1.x()/std::sqrt(2.*X);
142  PHI31add = -(-X*X*pol1.z() - 2.*X*(2.*pol0.z() - pol1.z()) - (4.*pol0.x() + 5.)*pol1.z())*(1. - eps)/(4.*X);
143 
144  PHI33 = pol1.x() - (X - 1.)*std::sqrt(1. - eps)*pol1.z()/std::sqrt(2.*X);
145  PHI33add = -(X*X - 2.*X + 4.*pol0.x() + 5.)*(1. - eps)*pol1.x()/(4.*X);
146  }
147  phi3 = G4ThreeVector(PHI31 + PHI31add, PHI32, PHI33 + PHI33add);
148 
149  }
150  unpXS *= diffXSFactor;
151  polXS *= diffXSFactor;
152  phi0 *= diffXSFactor;
153  phi2 *= diffXSFactor;
154  phi3 *= diffXSFactor;
155 
156 }
157 
158 G4double G4PolarizedComptonCrossSection::XSection(const G4StokesVector & pol2,const G4StokesVector & pol3)
159 {
160  gammaPol2 = !(pol2==G4StokesVector::ZERO);
161  electronPol3 = !(pol3==G4StokesVector::ZERO);
162 
163  G4double phi = 0.;
164  // polarization independent part
165  phi += phi0;
166 
167 
168  if (gammaPol2) {
169  // part depending on the polarization of the final photon
170  phi += phi2*pol2;
171  }
172 
173  if (electronPol3) {
174  // part depending on the polarization of the final electron
175  phi += phi3*pol3;
176  }
177 
178  // return cross section.
179  return phi;
180 }
181 
182 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
183 
184 
185 G4double G4PolarizedComptonCrossSection::TotalXSection(G4double /*xmin*/, G4double /*xmax*/, G4double k0,
186  const G4StokesVector & pol0,
187  const G4StokesVector & pol1)
188 {
189 
190  // G4double k0 = gammaEnergy / electron_mass_c2 ;
191  G4double k1 = 1. + 2.*k0 ;
192 
193 // // pi*re^2
194 // G4double re=2.81794e-15; //m
195 // G4double barn=1.e-28; //m^2
196  G4double Z=theZ;
197 
198  G4double unit = Z*pi*classic_electr_radius
199  * classic_electr_radius ; // *1./barn;
200 
201  G4double pre = unit/(sqr(k0)*sqr(1.+2.*k0));
202 
203  G4double xs_0 = ((k0 - 2.)*k0 -2.)*sqr(k1)*std::log(k1) + 2.*k0*(k0*(k0 + 1.)*(k0 + 8.) + 2.);
204  G4double xs_pol = (k0 + 1.)*sqr(k1)*std::log(k1) - 2.*k0*(5.*sqr(k0) + 4.*k0 + 1.);
205 
206  return pre*(xs_0/k0 + pol0.p3()*pol1.z()*xs_pol);
207 }
208 
209 
210 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
211 
212 
213 G4StokesVector G4PolarizedComptonCrossSection::GetPol2()
214 {
215  // Note, mean polarization can not contain correlation
216  // effects.
217  return 1./phi0 * phi2;
218 }
219 
220 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
221 
222 G4StokesVector G4PolarizedComptonCrossSection::GetPol3()
223 {
224  // Note, mean polarization can not contain correlation
225  // effects.
226  return 1./phi0 * phi3;
227 }
228 
229 
230 
231 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
232 void G4PolarizedComptonCrossSection::DefineCoefficients(const G4StokesVector & pol0,
233  const G4StokesVector & pol1)
234 {
235  polxx=pol0.x()*pol1.x();
236  polyy=pol0.y()*pol1.y();
237  polzz=pol0.z()*pol1.z();
238 
239  polxz=pol0.x()*pol1.z();
240  polzx=pol0.z()*pol1.x();
241 
242  polyz=pol0.y()*pol1.z();
243  polzy=pol0.z()*pol1.y();
244 
245  polxy=pol0.x()*pol1.y();
246  polyx=pol0.y()*pol1.x();
247 }
248 
249 
250 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......