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 #ifndef G4QElectronNuclearCrossSection_h
00040 #define G4QElectronNuclearCrossSection_h 1
00041
00042 #include <vector>
00043 #include "Randomize.hh"
00044 #include "G4VQCrossSection.hh"
00045
00046 class G4QElectronNuclearCrossSection : public G4VQCrossSection
00047 {
00048 protected:
00049
00050 G4QElectronNuclearCrossSection() {}
00051
00052 public:
00053
00054 ~G4QElectronNuclearCrossSection();
00055
00056 static G4VQCrossSection* GetPointer();
00057
00058 G4double ThresholdEnergy(G4int Z, G4int N, G4int PDG=11);
00059
00060
00061 virtual G4double GetCrossSection(G4bool fCS, G4double pMom, G4int tgZ, G4int tgN,
00062 G4int pPDG=0);
00063
00064 G4double CalculateCrossSection(G4bool CS, G4int F, G4int I, G4int PDG, G4int Z, G4int N,
00065 G4double Momentum);
00066
00067 G4int GetExchangePDGCode();
00068
00069 G4double GetExchangeEnergy();
00070
00071 G4double GetVirtualFactor(G4double nu, G4double Q2);
00072
00073 G4double GetExchangeQ2(G4double nu);
00074
00075 private:
00076 G4int GetFunctions(G4double a, G4double* x, G4double* y, G4double* z);
00077 G4double HighEnergyJ1(G4double lE);
00078 G4double HighEnergyJ2(G4double lE);
00079 G4double HighEnergyJ3(G4double lE);
00080 G4double SolveTheEquation(G4double f);
00081 G4double Fun(G4double x);
00082 G4double DFun(G4double x);
00083
00084
00085 private:
00086 static G4bool onlyCS;
00087 static G4double lastSig;
00088 static G4int lastL;
00089 static G4double lastE;
00090 static G4int lastF;
00091 static G4double lastG;
00092 static G4double lastH;
00093 static G4double* lastJ1;
00094 static G4double* lastJ2;
00095 static G4double* lastJ3;
00096 static G4int lastPDG;
00097 static G4int lastN;
00098 static G4int lastZ;
00099 static G4double lastP;
00100 static G4double lastTH;
00101 static G4double lastCS;
00102 static G4int lastI;
00103 static std::vector <G4double*>* J1;
00104 static std::vector <G4double*>* J2;
00105 static std::vector <G4double*>* J3;
00106 };
00107
00108 inline G4double G4QElectronNuclearCrossSection::DFun(G4double x)
00109 {
00110 static const G4double shd=1.0734;
00111 static const G4double poc=0.0375;
00112 static const G4double pos=16.5;
00113 static const G4double reg=.11;
00114 static const G4double mel=0.5109989;
00115 static const G4double lmel=std::log(mel);
00116 G4double y=std::exp(x-lastG-lmel);
00117 G4double flux=lastG*(2.-y*(2.-y))-1.;
00118 return (poc*(x-pos)+shd*std::exp(-reg*x))*flux;
00119 }
00120
00121 inline G4double G4QElectronNuclearCrossSection::Fun(G4double x)
00122 {
00123 G4double dlg1=lastG+lastG-1.;
00124 G4double lgoe=lastG/lastE;
00125 G4double HE2=HighEnergyJ2(x);
00126 return dlg1*HighEnergyJ1(x)-lgoe*(HE2+HE2-HighEnergyJ3(x)/lastE);
00127 }
00128
00129 inline G4double G4QElectronNuclearCrossSection::HighEnergyJ1(G4double lEn)
00130 {
00131 static const G4double le=std::log(50000.);
00132 static const G4double le2=le*le;
00133 static const G4double a=.0375;
00134 static const G4double ha=a*.5;
00135 static const G4double ab=a*16.5;
00136 static const G4double d=0.11;
00137 static const G4double cd=1.0734/d;
00138 static const G4double ele=std::exp(-d*le);
00139 return ha*(lEn*lEn-le2)-ab*(lEn-le)-cd*(std::exp(-d*lEn)-ele);
00140 }
00141
00142 inline G4double G4QElectronNuclearCrossSection::HighEnergyJ2(G4double lEn)
00143 {
00144 static const G4double e=50000.;
00145 static const G4double le=std::log(e);
00146 static const G4double le1=(le-1.)*e;
00147 static const G4double a=.0375;
00148 static const G4double ab=a*16.5;
00149 static const G4double d=1.-0.11;
00150 static const G4double cd=1.0734/d;
00151 static const G4double ele=std::exp(d*le);
00152 G4double En=std::exp(lEn);
00153 return a*((lEn-1.)*En-le1)-ab*(En-e)+cd*(std::exp(d*lEn)-ele);
00154 }
00155
00156 inline G4double G4QElectronNuclearCrossSection::HighEnergyJ3(G4double lEn)
00157 {
00158 static const G4double e=50000.;
00159 static const G4double le=std::log(e);
00160 static const G4double e2=e*e;
00161 static const G4double leh=(le-.5)*e2;
00162 static const G4double ha=.0375*.5;
00163 static const G4double hab=ha*16.5;
00164 static const G4double d=2.-.11;
00165 static const G4double cd=1.0734/d;
00166 static const G4double ele=std::exp(d*le);
00167 G4double lastE2=std::exp(lEn+lEn);
00168 return ha*((lEn-.5)*lastE2-leh)-hab*(lastE2-e2)+cd*(std::exp(d*lEn)-ele);
00169 }
00170
00171 #endif