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
00055 #ifndef G4EmCalculator_h
00056 #define G4EmCalculator_h 1
00057
00058 #include <vector>
00059 #include "globals.hh"
00060 #include "G4DataVector.hh"
00061 #include "G4DynamicParticle.hh"
00062 #include "G4VAtomDeexcitation.hh"
00063
00064 class G4LossTableManager;
00065 class G4Material;
00066 class G4MaterialCutsCouple;
00067 class G4ParticleDefinition;
00068 class G4PhysicsTable;
00069 class G4VEmModel;
00070 class G4VEnergyLossProcess;
00071 class G4VEmProcess;
00072 class G4VMultipleScattering;
00073 class G4VProcess;
00074 class G4ionEffectiveCharge;
00075 class G4Region;
00076 class G4Element;
00077 class G4EmCorrections;
00078
00079 class G4EmCalculator
00080 {
00081
00082 public:
00083
00084 G4EmCalculator();
00085
00086 ~G4EmCalculator();
00087
00088
00089
00090
00091
00092
00093 G4double GetDEDX(G4double kinEnergy, const G4ParticleDefinition*,
00094 const G4Material*,
00095 const G4Region* r = 0);
00096 G4double GetDEDX(G4double kinEnergy, const G4String& part,
00097 const G4String& mat,
00098 const G4String& s = "world");
00099
00100 G4double GetRangeFromRestricteDEDX(G4double kinEnergy,
00101 const G4ParticleDefinition*,
00102 const G4Material*,
00103 const G4Region* r = 0);
00104 G4double GetRangeFromRestricteDEDX(G4double kinEnergy, const G4String& part,
00105 const G4String& mat,
00106 const G4String& s = "world");
00107
00108 G4double GetCSDARange(G4double kinEnergy, const G4ParticleDefinition*,
00109 const G4Material*,
00110 const G4Region* r = 0);
00111 G4double GetCSDARange(G4double kinEnergy, const G4String& part,
00112 const G4String& mat,
00113 const G4String& s = "world");
00114
00115 G4double GetRange(G4double kinEnergy, const G4ParticleDefinition*,
00116 const G4Material*,
00117 const G4Region* r = 0);
00118 G4double GetRange(G4double kinEnergy, const G4String& part,
00119 const G4String& mat,
00120 const G4String& s = "world");
00121
00122 G4double GetKinEnergy(G4double range, const G4ParticleDefinition*,
00123 const G4Material*,
00124 const G4Region* r = 0);
00125 G4double GetKinEnergy(G4double range, const G4String& part,
00126 const G4String& mat,
00127 const G4String& s = "world");
00128
00129 G4double GetCrossSectionPerVolume(
00130 G4double kinEnergy, const G4ParticleDefinition*,
00131 const G4String& processName, const G4Material*,
00132 const G4Region* r = 0);
00133 G4double GetCrossSectionPerVolume(
00134 G4double kinEnergy, const G4String& part, const G4String& proc,
00135 const G4String& mat, const G4String& s = "world");
00136
00137 G4double GetShellIonisationCrossSectionPerAtom(
00138 const G4String& part, G4int Z,
00139 G4AtomicShellEnumerator shell,
00140 G4double kinEnergy);
00141
00142 G4double GetMeanFreePath(G4double kinEnergy, const G4ParticleDefinition*,
00143 const G4String& processName, const G4Material*,
00144 const G4Region* r = 0);
00145 G4double GetMeanFreePath(G4double kinEnergy, const G4String& part,
00146 const G4String& proc,
00147 const G4String& mat, const G4String& s = "world");
00148
00149 void PrintDEDXTable(const G4ParticleDefinition*);
00150
00151 void PrintRangeTable(const G4ParticleDefinition*);
00152
00153 void PrintInverseRangeTable(const G4ParticleDefinition*);
00154
00155
00156
00157
00158
00159
00160 G4double ComputeDEDX(G4double kinEnergy, const G4ParticleDefinition*,
00161 const G4String& processName, const G4Material*,
00162 G4double cut = DBL_MAX);
00163 G4double ComputeDEDX(G4double kinEnergy, const G4String& part,
00164 const G4String& proc,
00165 const G4String& mat, G4double cut = DBL_MAX);
00166
00167 G4double ComputeElectronicDEDX(G4double kinEnergy,
00168 const G4ParticleDefinition*,
00169 const G4Material* mat, G4double cut = DBL_MAX);
00170 G4double ComputeElectronicDEDX(G4double kinEnergy, const G4String& part,
00171 const G4String& mat, G4double cut = DBL_MAX);
00172
00173 G4double ComputeNuclearDEDX(G4double kinEnergy, const G4ParticleDefinition*,
00174 const G4Material*);
00175 G4double ComputeNuclearDEDX(G4double kinEnergy, const G4String& part,
00176 const G4String& mat);
00177
00178 G4double ComputeTotalDEDX(G4double kinEnergy, const G4ParticleDefinition*,
00179 const G4Material*, G4double cut = DBL_MAX);
00180 G4double ComputeTotalDEDX(G4double kinEnergy, const G4String& part,
00181 const G4String& mat, G4double cut = DBL_MAX);
00182
00183 G4double ComputeCrossSectionPerVolume(
00184 G4double kinEnergy, const G4ParticleDefinition*,
00185 const G4String& processName, const G4Material*,
00186 G4double cut = 0.0);
00187 G4double ComputeCrossSectionPerVolume(
00188 G4double kinEnergy, const G4String& part,
00189 const G4String& proc,
00190 const G4String& mat, G4double cut = 0.0);
00191
00192 G4double ComputeCrossSectionPerAtom(
00193 G4double kinEnergy, const G4ParticleDefinition*,
00194 const G4String& processName, G4double Z, G4double A,
00195 G4double cut = 0.0);
00196 G4double ComputeCrossSectionPerAtom(
00197 G4double kinEnergy, const G4String& part,
00198 const G4String& processName, const G4Element*,
00199 G4double cut = 0.0);
00200
00201 G4double ComputeGammaAttenuationLength(G4double kinEnergy,
00202 const G4Material*);
00203
00204 G4double ComputeShellIonisationCrossSectionPerAtom(
00205 const G4String& part, G4int Z,
00206 G4AtomicShellEnumerator shell,
00207 G4double kinEnergy,
00208 const G4Material* mat = 0);
00209
00210 G4double ComputeMeanFreePath(
00211 G4double kinEnergy, const G4ParticleDefinition*,
00212 const G4String& processName, const G4Material*,
00213 G4double cut = 0.0);
00214 G4double ComputeMeanFreePath(
00215 G4double kinEnergy, const G4String&, const G4String&,
00216 const G4String& processName, G4double cut = 0.0);
00217
00218 G4double ComputeEnergyCutFromRangeCut(
00219 G4double range, const G4ParticleDefinition*,
00220 const G4Material*);
00221 G4double ComputeEnergyCutFromRangeCut(
00222 G4double range, const G4String&,
00223 const G4String&);
00224
00225
00226
00227
00228
00229 const G4ParticleDefinition* FindParticle(const G4String&);
00230
00231 const G4ParticleDefinition* FindIon(G4int Z, G4int A);
00232
00233 const G4Material* FindMaterial(const G4String&);
00234
00235 const G4Region* FindRegion(const G4String&);
00236
00237 const G4MaterialCutsCouple* FindCouple(const G4Material*,
00238 const G4Region* r = 0);
00239
00240 void SetVerbose(G4int val);
00241
00242
00243
00244
00245
00246 private:
00247
00248 G4bool UpdateParticle(const G4ParticleDefinition*, G4double kinEnergy);
00249
00250 G4bool UpdateCouple(const G4Material*, G4double cut);
00251
00252 void FindLambdaTable(const G4ParticleDefinition*,
00253 const G4String& processName,
00254 G4double kinEnergy);
00255
00256 G4bool FindEmModel(const G4ParticleDefinition*,
00257 const G4String& processName,
00258 G4double kinEnergy);
00259
00260 G4VEnergyLossProcess* FindEnergyLossProcess(const G4ParticleDefinition*);
00261
00262 G4VEnergyLossProcess* FindEnLossProcess(const G4ParticleDefinition*,
00263 const G4String& processName);
00264
00265 G4VEmProcess* FindDiscreteProcess(const G4ParticleDefinition*,
00266 const G4String& processName);
00267
00268 G4VMultipleScattering* FindMscProcess(const G4ParticleDefinition*,
00269 const G4String& processName);
00270
00271 G4bool ActiveForParticle(const G4ParticleDefinition* part,
00272 G4VProcess* proc);
00273
00274 G4EmCalculator & operator=(const G4EmCalculator &right);
00275 G4EmCalculator(const G4EmCalculator&);
00276
00277 std::vector<const G4Material*> localMaterials;
00278 std::vector<const G4MaterialCutsCouple*> localCouples;
00279
00280 G4LossTableManager* manager;
00281 G4EmCorrections* corr;
00282 G4DataVector localCuts;
00283 G4int nLocalMaterials;
00284
00285 G4int verbose;
00286
00287
00288 G4int currentCoupleIndex;
00289 const G4MaterialCutsCouple* currentCouple;
00290 const G4Material* currentMaterial;
00291 const G4ParticleDefinition* currentParticle;
00292 const G4ParticleDefinition* lambdaParticle;
00293 const G4ParticleDefinition* baseParticle;
00294 const G4PhysicsTable* currentLambda;
00295 G4VEmModel* currentModel;
00296 G4VEmModel* loweModel;
00297 G4VEnergyLossProcess* currentProcess;
00298
00299 const G4ParticleDefinition* theGenericIon;
00300 G4ionEffectiveCharge* ionEffCharge;
00301 G4DynamicParticle dynParticle;
00302
00303 G4String currentName;
00304 G4String lambdaName;
00305 G4double currentCut;
00306 G4double chargeSquare;
00307 G4double massRatio;
00308 G4double mass;
00309 G4bool isIon;
00310 G4bool isApplicable;
00311
00312 G4String currentParticleName;
00313 G4String currentMaterialName;
00314 G4String currentProcessName;
00315 };
00316
00317
00318
00319 #endif