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 #include "G4PolarizedPairProductionCrossSection.hh"
00041 #include "G4PhysicalConstants.hh"
00042
00043 G4bool G4PolarizedPairProductionCrossSection::scrnInitialized=false;
00044 G4double G4PolarizedPairProductionCrossSection::SCRN [3][20];
00045
00046
00047 void G4PolarizedPairProductionCrossSection::InitializeMe()
00048 {
00049 if (!scrnInitialized) {
00050 SCRN [1][1]= 0.5 ; SCRN [2][1] = 0.0145;
00051 SCRN [1][2]= 1.0 ; SCRN [2][2] = 0.0490;
00052 SCRN [1][3]= 2.0 ; SCRN [2][3] = 0.1400;
00053 SCRN [1][4]= 4.0 ; SCRN [2][4] = 0.3312;
00054 SCRN [1][5]= 8.0 ; SCRN [2][5] = 0.6758;
00055 SCRN [1][6]= 15.0 ; SCRN [2][6] = 1.126;
00056 SCRN [1][7]= 20.0 ; SCRN [2][7] = 1.367;
00057 SCRN [1][8]= 25.0 ; SCRN [2][8] = 1.564;
00058 SCRN [1][9]= 30.0 ; SCRN [2][9] = 1.731;
00059 SCRN [1][10]= 35.0 ; SCRN [2][10]= 1.875;
00060 SCRN [1][11]= 40.0 ; SCRN [2][11]= 2.001;
00061 SCRN [1][12]= 45.0 ; SCRN [2][12]= 2.114;
00062 SCRN [1][13]= 50.0 ; SCRN [2][13]= 2.216;
00063 SCRN [1][14]= 60.0 ; SCRN [2][14]= 2.393;
00064 SCRN [1][15]= 70.0 ; SCRN [2][15]= 2.545;
00065 SCRN [1][16]= 80.0 ; SCRN [2][16]= 2.676;
00066 SCRN [1][17]= 90.0 ; SCRN [2][17]= 2.793;
00067 SCRN [1][18]= 100.0 ; SCRN [2][18]= 2.897;
00068 SCRN [1][19]= 120.0 ; SCRN [2][19]= 3.078;
00069
00070 scrnInitialized=true;
00071 }
00072 }
00073
00074 G4PolarizedPairProductionCrossSection::G4PolarizedPairProductionCrossSection()
00075 {
00076 InitializeMe();
00077 }
00078
00079
00080 void G4PolarizedPairProductionCrossSection::Initialize(G4double aGammaE,
00081 G4double aLept0E,
00082 G4double sintheta,
00083 const G4StokesVector & beamPol,
00084 const G4StokesVector & ,
00085 G4int )
00086 {
00087 G4double aLept1E = aGammaE - aLept0E;
00088
00089 G4double Stokes_P3 = beamPol.z() ;
00090
00091
00092 G4double m0_c2 = electron_mass_c2;
00093 G4double Lept0E = aLept0E/m0_c2+1., Lept0E2 = Lept0E * Lept0E ;
00094 G4double GammaE = aGammaE/m0_c2;
00095 G4double Lept1E = aLept1E/m0_c2-1., Lept1E2 = Lept1E * Lept1E ;
00096
00097
00098
00099
00100
00101
00102 G4double TMom = std::sqrt(Lept0E2 -1.)* sintheta, u = TMom , u2 =u * u ;
00103 G4double Xsi = 1./(1.+u2) , Xsi2 = Xsi * Xsi ;
00104
00105
00106
00107
00108 G4double delta = 12. * std::pow(theZ, 1./3.) * Lept0E * Lept1E * Xsi / (121. * GammaE);
00109 G4double GG=0.;
00110
00111 if(delta < 0.5) {
00112 GG = std::log(2.* Lept0E * Lept1E / GammaE) - 2. - fCoul;
00113 }
00114 else if ( delta < 120) {
00115 for (G4int j=2; j<=19; j++) {
00116 if(SCRN[1][j] >= delta) {
00117 GG =std::log(2 * Lept0E * Lept1E / GammaE) - 2 - fCoul
00118 -(SCRN[2][j-1]+(delta-SCRN[1][j-1])*(SCRN[2][j]-SCRN[2][j-1])/(SCRN[1][j]-SCRN[1][j-1]));
00119 break;
00120 }
00121 }
00122 }
00123 else {
00124 G4double alpha_sc = (111 * std::pow(theZ, -1./3.)) / Xsi;
00125 GG = std::log(alpha_sc)- 2 - fCoul;
00126 }
00127
00128 if(GG<-1) GG=-1;
00129
00130
00131 G4double I_Lepton = (Lept0E2 + Lept1E2)*(3+2*GG) + 2 * Lept0E * Lept1E * (1 + 4 * u2 * Xsi2 * GG);
00132
00133
00134
00135 G4double L_Lepton1 = GammaE * ((Lept0E - Lept1E) * (3. + 2. * GG)+2 * Lept1E * (1 + 4 * u2 * Xsi2 * GG))/I_Lepton;
00136
00137 G4double T_Lepton1 = 4 * GammaE * Lept1E * Xsi * u * (1. - 2. * Xsi) * GG / I_Lepton ;
00138
00139
00140 G4double Stokes_S1 = (Stokes_P3 * T_Lepton1) ;
00141 G4double Stokes_S2 = 0;
00142 G4double Stokes_S3 = (Stokes_P3 * L_Lepton1) ;
00143
00144
00145 theFinalElectronPolarization.setX(Stokes_S1);
00146 theFinalElectronPolarization.setY(Stokes_S2);
00147 theFinalElectronPolarization.setZ(Stokes_S3);
00148
00149 if(theFinalElectronPolarization.mag2()>1) {
00150 G4cout<<" WARNING in pol-conv theFinalElectronPolarization \n";
00151 G4cout
00152 <<"\t"<<theFinalElectronPolarization
00153 <<"\t GG\t"<<GG
00154 <<"\t delta\t"<<delta
00155 <<G4endl;
00156 theFinalElectronPolarization.setX(0);
00157 theFinalElectronPolarization.setY(0);
00158 theFinalElectronPolarization.setZ(Stokes_S3);
00159 if(Stokes_S3>1) theFinalElectronPolarization.setZ(1);
00160 }
00161
00162
00163 G4double L_Lepton2 = GammaE * ((Lept1E - Lept0E) * (3. + 2. * GG)+2 * Lept0E * (1 + 4 * u2 * Xsi2 * GG))/I_Lepton;
00164
00165 G4double T_Lepton2 = 4 * GammaE * Lept0E * Xsi * u * (1. - 2. * Xsi) * GG / I_Lepton ;
00166
00167 G4double Stokes_SS1 = (Stokes_P3 * T_Lepton2) ;
00168 G4double Stokes_SS2 = 0;
00169 G4double Stokes_SS3 = (Stokes_P3 * L_Lepton2) ;
00170
00171 theFinalPositronPolarization.SetPhoton();
00172
00173 theFinalPositronPolarization.setX(Stokes_SS1);
00174 theFinalPositronPolarization.setY(Stokes_SS2);
00175 theFinalPositronPolarization.setZ(Stokes_SS3);
00176
00177 if(theFinalPositronPolarization.mag2()>1) {
00178 G4cout<<" WARNING in pol-conv theFinalPositronPolarization \n";
00179 G4cout
00180 <<"\t"<<theFinalPositronPolarization
00181 <<"\t GG\t"<<GG
00182 <<"\t delta\t"<<delta
00183 <<G4endl;
00184 }
00185 }
00186
00187 G4double G4PolarizedPairProductionCrossSection::XSection(const G4StokesVector & ,
00188 const G4StokesVector & )
00189 {
00190 G4cout<<"ERROR dummy routine G4PolarizedPairProductionCrossSection::XSection called \n";
00191 return 0;
00192 }
00193
00194
00195 G4StokesVector G4PolarizedPairProductionCrossSection::GetPol2()
00196 {
00197
00198 return theFinalElectronPolarization;
00199 }
00200 G4StokesVector G4PolarizedPairProductionCrossSection::GetPol3()
00201 {
00202
00203 return theFinalPositronPolarization;;
00204 }
00205
00206