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
00050
00051
00052
00053
00054 #ifndef G4EmCorrections_h
00055 #define G4EmCorrections_h 1
00056
00057 #include <CLHEP/Units/PhysicalConstants.h>
00058
00059 #include "globals.hh"
00060 #include "G4ionEffectiveCharge.hh"
00061 #include "G4Material.hh"
00062 #include "G4ParticleDefinition.hh"
00063 #include "G4NistManager.hh"
00064
00065 class G4VEmModel;
00066 class G4PhysicsVector;
00067 class G4IonTable;
00068 class G4MaterialCutsCouple;
00069 class G4LPhysicsFreeVector;
00070
00071 class G4EmCorrections
00072 {
00073
00074 public:
00075
00076 G4EmCorrections();
00077
00078 virtual ~G4EmCorrections();
00079
00080 G4double HighOrderCorrections(const G4ParticleDefinition*,
00081 const G4Material*,
00082 G4double kineticEnergy,
00083 G4double cutEnergy);
00084
00085 G4double IonHighOrderCorrections(const G4ParticleDefinition*,
00086 const G4MaterialCutsCouple*,
00087 G4double kineticEnergy);
00088
00089 G4double ComputeIonCorrections(const G4ParticleDefinition*,
00090 const G4Material*,
00091 G4double kineticEnergy);
00092
00093 G4double IonBarkasCorrection(const G4ParticleDefinition*,
00094 const G4Material*,
00095 G4double kineticEnergy);
00096
00097 G4double Bethe(const G4ParticleDefinition*,
00098 const G4Material*,
00099 G4double kineticEnergy);
00100
00101 G4double SpinCorrection(const G4ParticleDefinition*,
00102 const G4Material*,
00103 G4double kineticEnergy);
00104
00105 G4double KShellCorrection(const G4ParticleDefinition*,
00106 const G4Material*,
00107 G4double kineticEnergy);
00108
00109 G4double LShellCorrection(const G4ParticleDefinition*,
00110 const G4Material*,
00111 G4double kineticEnergy);
00112
00113 G4double ShellCorrection(const G4ParticleDefinition*,
00114 const G4Material*,
00115 G4double kineticEnergy);
00116
00117 G4double ShellCorrectionSTD(const G4ParticleDefinition*,
00118 const G4Material*,
00119 G4double kineticEnergy);
00120
00121 G4double DensityCorrection(const G4ParticleDefinition*,
00122 const G4Material*,
00123 G4double kineticEnergy);
00124
00125 G4double BarkasCorrection(const G4ParticleDefinition*,
00126 const G4Material*,
00127 G4double kineticEnergy);
00128
00129 G4double BlochCorrection(const G4ParticleDefinition*,
00130 const G4Material*,
00131 G4double kineticEnergy);
00132
00133 G4double MottCorrection(const G4ParticleDefinition*,
00134 const G4Material*,
00135 G4double kineticEnergy);
00136
00137 G4double NuclearDEDX(const G4ParticleDefinition*,
00138 const G4Material*,
00139 G4double kineticEnergy,
00140 G4bool fluct = true);
00141
00142 void AddStoppingData(G4int Z, G4int A, const G4String& materialName,
00143 G4PhysicsVector* dVector);
00144
00145 void InitialiseForNewRun();
00146
00147
00148 G4double EffectiveChargeCorrection(const G4ParticleDefinition*,
00149 const G4Material*,
00150 G4double kineticEnergy);
00151
00152
00153 inline G4double GetParticleCharge(const G4ParticleDefinition*,
00154 const G4Material*,
00155 G4double kineticEnergy);
00156
00157 inline
00158 G4double EffectiveChargeSquareRatio(const G4ParticleDefinition*,
00159 const G4Material*,
00160 G4double kineticEnergy);
00161
00162
00163 inline void SetIonisationModels(G4VEmModel* m1 = 0, G4VEmModel* m2 = 0);
00164
00165 inline G4int GetNumberOfStoppingVectors();
00166
00167 private:
00168
00169 void Initialise();
00170
00171 void BuildCorrectionVector();
00172
00173 void SetupKinematics(const G4ParticleDefinition*,
00174 const G4Material*,
00175 G4double kineticEnergy);
00176
00177 G4double KShell(G4double theta, G4double eta);
00178
00179 G4double LShell(G4double theta, G4double eta);
00180
00181 G4int Index(G4double x, G4double* y, G4int n);
00182
00183 G4double Value(G4double xv, G4double x1, G4double x2, G4double y1, G4double y2);
00184
00185 G4double Value2(G4double xv, G4double yv, G4double x1, G4double x2,
00186 G4double y1, G4double y2,
00187 G4double z11, G4double z21, G4double z12, G4double z22);
00188
00189 G4double NuclearStoppingPower(G4double e, G4double z1, G4double z2,
00190 G4double m1, G4double m2);
00191
00192
00193 G4EmCorrections & operator=(const G4EmCorrections &right);
00194 G4EmCorrections(const G4EmCorrections&);
00195
00196 G4double ed[104];
00197 G4double a[104];
00198 G4double theZieglerFactor;
00199 G4double alpha2;
00200 G4bool lossFlucFlag;
00201
00202 G4int verbose;
00203
00204 G4int nK;
00205 G4int nL;
00206 G4int nEtaK;
00207 G4int nEtaL;
00208
00209 G4double COSEB[14];
00210 G4double COSXI[14];
00211 G4double ZD[11];
00212
00213 G4double TheK[20];
00214 G4double SK[20];
00215 G4double TK[20];
00216 G4double UK[20];
00217 G4double VK[20];
00218 G4double ZK[20];
00219
00220 G4double TheL[26];
00221 G4double SL[26];
00222 G4double TL[26];
00223 G4double UL[26];
00224 G4double VL[26];
00225
00226 G4double Eta[29];
00227 G4double CK[20][29];
00228 G4double CL[26][28];
00229 G4double HM[53];
00230 G4double HN[31];
00231 G4double Z23[100];
00232
00233 G4LPhysicsFreeVector* BarkasCorr;
00234 G4LPhysicsFreeVector* ThetaK;
00235 G4LPhysicsFreeVector* ThetaL;
00236
00237 std::vector<const G4Material*> currmat;
00238 std::map< G4int, std::vector<G4double> > thcorr;
00239 size_t ncouples;
00240
00241 const G4ParticleDefinition* particle;
00242 const G4ParticleDefinition* curParticle;
00243 const G4Material* material;
00244 const G4Material* curMaterial;
00245 const G4ElementVector* theElementVector;
00246 const G4double* atomDensity;
00247
00248 G4int numberOfElements;
00249 G4double kinEnergy;
00250 G4double mass;
00251 G4double massFactor;
00252 G4double formfact;
00253 G4double eth;
00254 G4double tau;
00255 G4double gamma;
00256 G4double bg2;
00257 G4double beta2;
00258 G4double beta;
00259 G4double ba2;
00260 G4double tmax;
00261 G4double charge;
00262 G4double q2;
00263 G4double eCorrMin;
00264 G4double eCorrMax;
00265 G4int nbinCorr;
00266
00267 G4ionEffectiveCharge effCharge;
00268
00269 G4NistManager* nist;
00270 const G4IonTable* ionTable;
00271 G4VEmModel* ionLEModel;
00272 G4VEmModel* ionHEModel;
00273
00274
00275 G4int nIons;
00276 G4int idx;
00277 G4int currentZ;
00278 std::vector<G4int> Zion;
00279 std::vector<G4int> Aion;
00280 std::vector<G4String> materialName;
00281
00282 std::vector<const G4ParticleDefinition*> ionList;
00283
00284 std::vector<const G4Material*> materialList;
00285 std::vector<G4PhysicsVector*> stopData;
00286 G4PhysicsVector* curVector;
00287 };
00288
00289 inline G4int G4EmCorrections::Index(G4double x, G4double* y, G4int n)
00290 {
00291 G4int iddd = n-1;
00292 do {iddd--;} while (iddd>0 && x<y[iddd]);
00293 return iddd;
00294 }
00295
00296 inline G4double G4EmCorrections::Value(G4double xv, G4double x1, G4double x2,
00297 G4double y1, G4double y2)
00298 {
00299 return y1 + (y2 - y1)*(xv - x1)/(x2 - x1);
00300 }
00301
00302 inline G4double G4EmCorrections::Value2(G4double xv, G4double yv,
00303 G4double x1, G4double x2,
00304 G4double y1, G4double y2,
00305 G4double z11, G4double z21,
00306 G4double z12, G4double z22)
00307 {
00308 return (z11*(x2-xv)*(y2-yv) + z22*(xv-x1)*(yv-y1) +
00309 0.5*(z12*((x2-xv)*(yv-y1)+(xv-x1)*(y2-yv))+
00310 z21*((xv-x1)*(y2-yv)+(yv-y1)*(x2-xv))))
00311 / ((x2-x1)*(y2-y1));
00312 }
00313
00314 inline void
00315 G4EmCorrections::SetIonisationModels(G4VEmModel* mod1, G4VEmModel* mod2)
00316 {
00317 if(mod1) { ionLEModel = mod1; }
00318 if(mod2) { ionHEModel = mod2; }
00319 }
00320
00321 inline G4int G4EmCorrections::GetNumberOfStoppingVectors()
00322 {
00323 return nIons;
00324 }
00325
00326 inline G4double
00327 G4EmCorrections::GetParticleCharge(const G4ParticleDefinition* p,
00328 const G4Material* mat,
00329 G4double kineticEnergy)
00330 {
00331 return effCharge.EffectiveCharge(p,mat,kineticEnergy);
00332 }
00333
00334 inline G4double
00335 G4EmCorrections::EffectiveChargeSquareRatio(const G4ParticleDefinition* p,
00336 const G4Material* mat,
00337 G4double kineticEnergy)
00338 {
00339 return effCharge.EffectiveChargeSquareRatio(p,mat,kineticEnergy);
00340 }
00341
00342 inline void G4EmCorrections::SetupKinematics(const G4ParticleDefinition* p,
00343 const G4Material* mat,
00344 G4double kineticEnergy)
00345 {
00346 if(kineticEnergy != kinEnergy || p != particle) {
00347 particle = p;
00348 kinEnergy = kineticEnergy;
00349 mass = p->GetPDGMass();
00350 tau = kineticEnergy / mass;
00351 gamma = 1.0 + tau;
00352 bg2 = tau * (tau+2.0);
00353 beta2 = bg2/(gamma*gamma);
00354 beta = std::sqrt(beta2);
00355 ba2 = beta2/alpha2;
00356 G4double ratio = CLHEP::electron_mass_c2/mass;
00357 tmax = 2.0*CLHEP::electron_mass_c2*bg2 /(1. + 2.0*gamma*ratio + ratio*ratio);
00358 charge = p->GetPDGCharge()/CLHEP::eplus;
00359
00360
00361
00362
00363
00364 if(charge > 1.5) { charge = effCharge.EffectiveCharge(p,mat,kinEnergy); }
00365 q2 = charge*charge;
00366 }
00367 if(mat != material) {
00368 material = mat;
00369 theElementVector = material->GetElementVector();
00370 atomDensity = material->GetAtomicNumDensityVector();
00371 numberOfElements = material->GetNumberOfElements();
00372 }
00373 }
00374
00375
00376
00377 #endif