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 G4HIMPACTIONISATION
00056 #define G4HIMPACTIONISATION 1
00057
00058 #include <map>
00059 #include <CLHEP/Units/PhysicalConstants.h>
00060
00061 #include "globals.hh"
00062 #include "G4hRDEnergyLoss.hh"
00063 #include "G4DataVector.hh"
00064 #include "G4AtomicDeexcitation.hh"
00065 #include "G4PixeCrossSectionHandler.hh"
00066
00067 class G4VLowEnergyModel;
00068 class G4VParticleChange;
00069 class G4ParticleDefinition;
00070 class G4PhysicsTable;
00071 class G4MaterialCutsCouple;
00072 class G4Track;
00073 class G4Step;
00074
00075 class G4hImpactIonisation : public G4hRDEnergyLoss
00076 {
00077 public:
00078
00079 G4hImpactIonisation(const G4String& processName = "hImpactIoni");
00080
00081
00082
00083 ~G4hImpactIonisation();
00084
00085
00086 G4bool IsApplicable(const G4ParticleDefinition&);
00087
00088
00089 void BuildPhysicsTable(const G4ParticleDefinition& aParticleType) ;
00090
00091
00092 G4double GetMeanFreePath(const G4Track& track,
00093 G4double previousStepSize,
00094 enum G4ForceCondition* condition );
00095
00096
00097 void PrintInfoDefinition() const;
00098
00099
00100 void SetHighEnergyForProtonParametrisation(G4double energy) {protonHighEnergy = energy;} ;
00101
00102
00103
00104
00105 void SetLowEnergyForProtonParametrisation(G4double energy) {protonLowEnergy = energy;} ;
00106
00107
00108
00109
00110 void SetHighEnergyForAntiProtonParametrisation(G4double energy) {antiprotonHighEnergy = energy;} ;
00111
00112
00113
00114
00115 void SetLowEnergyForAntiProtonParametrisation(G4double energy) {antiprotonLowEnergy = energy;} ;
00116
00117
00118
00119
00120 G4double GetContinuousStepLimit(const G4Track& track,
00121 G4double previousStepSize,
00122 G4double currentMinimumStep,
00123 G4double& currentSafety);
00124
00125
00126 void SetElectronicStoppingPowerModel(const G4ParticleDefinition* aParticle,
00127 const G4String& dedxTable);
00128
00129
00130
00131 void SetNuclearStoppingPowerModel(const G4String& dedxTable)
00132 {theNuclearTable = dedxTable; SetNuclearStoppingOn();};
00133
00134
00135
00136
00137
00138
00139 void SetNuclearStoppingOn() {nStopping = true;};
00140
00141
00142 void SetNuclearStoppingOff() {nStopping = false;};
00143
00144
00145 void SetBarkasOn() {theBarkas = true;};
00146
00147
00148 void SetBarkasOff() {theBarkas = false;};
00149
00150
00151 void SetPixe(const G4bool ) {pixeIsActive = true;};
00152
00153
00154 G4VParticleChange* AlongStepDoIt(const G4Track& trackData ,
00155 const G4Step& stepData ) ;
00156
00157
00158 G4VParticleChange* PostStepDoIt(const G4Track& track,
00159 const G4Step& Step ) ;
00160
00161
00162 G4double ComputeDEDX(const G4ParticleDefinition* aParticle,
00163 const G4MaterialCutsCouple* couple,
00164 G4double kineticEnergy);
00165
00166
00167 void SetCutForSecondaryPhotons(G4double cut);
00168
00169
00170 void SetCutForAugerElectrons(G4double cut);
00171
00172
00173 void ActivateAugerElectronProduction(G4bool val);
00174
00175
00176
00177 void SetPixeCrossSectionK(const G4String& name) { modelK = name; }
00178 void SetPixeCrossSectionL(const G4String& name) { modelL = name; }
00179 void SetPixeCrossSectionM(const G4String& name) { modelM = name; }
00180 void SetPixeProjectileMinEnergy(G4double energy) { eMinPixe = energy; }
00181 void SetPixeProjectileMaxEnergy(G4double energy) { eMaxPixe = energy; }
00182
00183 protected:
00184
00185 private:
00186
00187 void InitializeMe();
00188 void InitializeParametrisation();
00189 void BuildLossTable(const G4ParticleDefinition& aParticleType);
00190
00191 void BuildLambdaTable(const G4ParticleDefinition& aParticleType);
00192 void SetProtonElectronicStoppingPowerModel(const G4String& dedxTable)
00193 {protonTable = dedxTable ;};
00194
00195
00196 void SetAntiProtonElectronicStoppingPowerModel(const G4String& dedxTable)
00197 {antiprotonTable = dedxTable;};
00198
00199 G4double MicroscopicCrossSection(const G4ParticleDefinition& aParticleType,
00200 G4double kineticEnergy,
00201 G4double atomicNumber,
00202 G4double deltaCutInEnergy) const;
00203
00204 G4double GetConstraints(const G4DynamicParticle* particle,
00205 const G4MaterialCutsCouple* couple);
00206
00207
00208 G4double ProtonParametrisedDEDX(const G4MaterialCutsCouple* couple,
00209 G4double kineticEnergy) const;
00210
00211 G4double AntiProtonParametrisedDEDX(const G4MaterialCutsCouple* couple,
00212 G4double kineticEnergy) const;
00213
00214 G4double DeltaRaysEnergy(const G4MaterialCutsCouple* couple,
00215 G4double kineticEnergy,
00216 G4double particleMass) const;
00217
00218
00219
00220 G4double BarkasTerm(const G4Material* material,
00221 G4double kineticEnergy) const;
00222
00223
00224 G4double BlochTerm(const G4Material* material,
00225 G4double kineticEnergy,
00226 G4double cSquare) const;
00227
00228
00229 G4double ElectronicLossFluctuation(const G4DynamicParticle* particle,
00230 const G4MaterialCutsCouple* material,
00231 G4double meanLoss,
00232 G4double step) const;
00233
00234
00235
00236 G4hImpactIonisation & operator=(const G4hImpactIonisation &right);
00237 G4hImpactIonisation(const G4hImpactIonisation&);
00238
00239 private:
00240
00241 G4VLowEnergyModel* betheBlochModel;
00242 G4VLowEnergyModel* protonModel;
00243 G4VLowEnergyModel* antiprotonModel;
00244 G4VLowEnergyModel* theIonEffChargeModel;
00245 G4VLowEnergyModel* theNuclearStoppingModel;
00246 G4VLowEnergyModel* theIonChuFluctuationModel;
00247 G4VLowEnergyModel* theIonYangFluctuationModel;
00248
00249
00250
00251
00252 G4String protonTable;
00253 G4String antiprotonTable;
00254 G4String theNuclearTable;
00255
00256
00257 G4double protonLowEnergy;
00258 G4double protonHighEnergy;
00259 G4double antiprotonLowEnergy;
00260 G4double antiprotonHighEnergy;
00261
00262
00263 G4bool nStopping;
00264 G4bool theBarkas;
00265
00266 G4DataVector cutForDelta;
00267 G4DataVector cutForGamma;
00268 G4double minGammaEnergy;
00269 G4double minElectronEnergy;
00270 G4PhysicsTable* theMeanFreePathTable;
00271
00272 const G4double paramStepLimit;
00273
00274 G4double fdEdx;
00275 G4double fRangeNow ;
00276 G4double charge;
00277 G4double chargeSquare;
00278 G4double initialMass;
00279 G4double fBarkas;
00280
00281 G4PixeCrossSectionHandler* pixeCrossSectionHandler;
00282 G4AtomicDeexcitation atomicDeexcitation;
00283 G4String modelK;
00284 G4String modelL;
00285 G4String modelM;
00286 G4double eMinPixe;
00287 G4double eMaxPixe;
00288
00289 G4bool pixeIsActive;
00290
00291 };
00292
00293
00294 inline G4double G4hImpactIonisation::GetContinuousStepLimit(const G4Track& track,
00295 G4double,
00296 G4double currentMinimumStep,
00297 G4double&)
00298 {
00299 G4double step = GetConstraints(track.GetDynamicParticle(),track.GetMaterialCutsCouple()) ;
00300
00301
00302
00303
00304
00305 if ((step > 0.) && (step < currentMinimumStep)) currentMinimumStep = step ;
00306
00307 return step ;
00308 }
00309
00310
00311 inline G4bool G4hImpactIonisation::IsApplicable(const G4ParticleDefinition& particle)
00312 {
00313
00314
00315
00316 return (particle.GetPDGCharge() != 0.0 && particle.GetPDGMass() > CLHEP::proton_mass_c2*0.1);
00317 }
00318
00319 #endif
00320
00321
00322
00323
00324
00325
00326
00327