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
00056
00057
00058
00059
00060
00061
00062
00063
00064
00065
00066
00067
00068
00069
00070
00071 #ifndef G4VMultipleScattering_h
00072 #define G4VMultipleScattering_h 1
00073
00074 #include "G4VContinuousDiscreteProcess.hh"
00075 #include "globals.hh"
00076 #include "G4Material.hh"
00077 #include "G4ParticleChangeForMSC.hh"
00078 #include "G4Track.hh"
00079 #include "G4Step.hh"
00080 #include "G4EmModelManager.hh"
00081 #include "G4VMscModel.hh"
00082 #include "G4MscStepLimitType.hh"
00083
00084 class G4ParticleDefinition;
00085 class G4VEnergyLossProcess;
00086 class G4LossTableManager;
00087 class G4SafetyHelper;
00088
00089
00090
00091 class G4VMultipleScattering : public G4VContinuousDiscreteProcess
00092 {
00093 public:
00094
00095 G4VMultipleScattering(const G4String& name = "msc",
00096 G4ProcessType type = fElectromagnetic);
00097
00098 virtual ~G4VMultipleScattering();
00099
00100
00101
00102
00103
00104 virtual G4bool IsApplicable(const G4ParticleDefinition& p) = 0;
00105
00106 virtual void PrintInfo() = 0;
00107
00108 protected:
00109
00110 virtual void InitialiseProcess(const G4ParticleDefinition*) = 0;
00111
00112 public:
00113
00114
00115
00116
00117
00118
00119 void PreparePhysicsTable(const G4ParticleDefinition&);
00120
00121
00122 void BuildPhysicsTable(const G4ParticleDefinition&);
00123
00124
00125 void PrintInfoDefinition();
00126
00127
00128
00129 G4bool StorePhysicsTable(const G4ParticleDefinition*,
00130 const G4String& directory,
00131 G4bool ascii = false);
00132
00133
00134
00135
00136
00137
00138 G4bool RetrievePhysicsTable(const G4ParticleDefinition*,
00139 const G4String& directory,
00140 G4bool ascii);
00141
00142
00143 void StartTracking(G4Track*);
00144
00145
00146
00147
00148 G4double AlongStepGetPhysicalInteractionLength(
00149 const G4Track&,
00150 G4double previousStepSize,
00151 G4double currentMinimalStep,
00152 G4double& currentSafety,
00153 G4GPILSelection* selection);
00154
00155
00156
00157 G4double PostStepGetPhysicalInteractionLength(
00158 const G4Track&,
00159 G4double previousStepSize,
00160 G4ForceCondition* condition);
00161
00162
00163 G4VParticleChange* AlongStepDoIt(const G4Track&, const G4Step&);
00164
00165
00166 G4VParticleChange* PostStepDoIt(const G4Track&, const G4Step&);
00167
00168
00169 G4double ContinuousStepLimit(const G4Track& track,
00170 G4double previousStepSize,
00171 G4double currentMinimalStep,
00172 G4double& currentSafety);
00173
00174
00175
00176
00177
00178
00179 inline G4VEmModel* SelectModel(G4double kinEnergy, size_t idx);
00180
00181 public:
00182
00183
00184
00185 void AddEmModel(G4int order, G4VEmModel*, const G4Region* region = 0);
00186
00187
00188 void SetModel(G4VMscModel*, G4int index = 1);
00189
00190
00191 G4VMscModel* Model(G4int index = 1);
00192
00193
00194 void SetEmModel(G4VMscModel*, G4int index = 1);
00195
00196
00197 G4VMscModel* EmModel(G4int index = 1);
00198
00199
00200 G4VEmModel* GetModelByIndex(G4int idx = 0, G4bool ver = false) const;
00201
00202
00203
00204
00205
00206 void SetIonisation(G4VEnergyLossProcess*);
00207
00208 inline G4bool LateralDisplasmentFlag() const;
00209 inline void SetLateralDisplasmentFlag(G4bool val);
00210
00211 inline G4double Skin() const;
00212 inline void SetSkin(G4double val);
00213
00214 inline G4double RangeFactor() const;
00215 inline void SetRangeFactor(G4double val);
00216
00217 inline G4double GeomFactor() const;
00218 inline void SetGeomFactor(G4double val);
00219
00220 inline G4double PolarAngleLimit() const;
00221 inline void SetPolarAngleLimit(G4double val);
00222
00223 inline G4MscStepLimitType StepLimitType() const;
00224 inline void SetStepLimitType(G4MscStepLimitType val);
00225
00226 inline const G4ParticleDefinition* FirstParticle() const;
00227
00228
00229
00230
00231
00232 protected:
00233
00234
00235 G4double GetMeanFreePath(const G4Track& track,
00236 G4double,
00237 G4ForceCondition* condition);
00238
00239
00240 G4double GetContinuousStepLimit(const G4Track& track,
00241 G4double previousStepSize,
00242 G4double currentMinimalStep,
00243 G4double& currentSafety);
00244
00245 private:
00246
00247
00248 G4VMultipleScattering(G4VMultipleScattering &);
00249 G4VMultipleScattering & operator=(const G4VMultipleScattering &right);
00250
00251
00252
00253 G4EmModelManager* modelManager;
00254 G4LossTableManager* emManager;
00255 G4double geomMin;
00256
00257
00258
00259 G4SafetyHelper* safetyHelper;
00260
00261 std::vector<G4VMscModel*> mscModels;
00262 G4int numberOfModels;
00263
00264 const G4ParticleDefinition* firstParticle;
00265 const G4ParticleDefinition* currParticle;
00266
00267 G4MscStepLimitType stepLimit;
00268
00269 G4double skin;
00270 G4double facrange;
00271 G4double facgeom;
00272 G4double polarAngleLimit;
00273 G4double lowestKinEnergy;
00274
00275 G4bool latDisplasment;
00276 G4bool isIon;
00277
00278
00279
00280 protected:
00281
00282 G4GPILSelection valueGPILSelectionMSC;
00283 G4ParticleChangeForMSC fParticleChange;
00284
00285 private:
00286
00287
00288 G4VMscModel* currentModel;
00289 G4VEnergyLossProcess* fIonisation;
00290
00291 G4double physStepLimit;
00292 G4double tPathLength;
00293 G4double gPathLength;
00294
00295 G4ThreeVector fNewPosition;
00296 G4bool fPositionChanged;
00297 G4bool isActive;
00298
00299 G4int warn;
00300 };
00301
00302
00303
00304
00305
00306 inline G4VEmModel*
00307 G4VMultipleScattering::SelectModel(G4double kinEnergy, size_t coupleIndex)
00308 {
00309 return modelManager->SelectModel(kinEnergy, coupleIndex);
00310 }
00311
00312
00313
00314 inline G4bool G4VMultipleScattering::LateralDisplasmentFlag() const
00315 {
00316 return latDisplasment;
00317 }
00318
00319
00320
00321 inline void G4VMultipleScattering::SetLateralDisplasmentFlag(G4bool val)
00322 {
00323 latDisplasment = val;
00324 }
00325
00326
00327
00328 inline G4double G4VMultipleScattering::Skin() const
00329 {
00330 return skin;
00331 }
00332
00333
00334
00335 inline void G4VMultipleScattering::SetSkin(G4double val)
00336 {
00337 if(val < 1.0) { skin = 0.0; }
00338 else { skin = val; }
00339 }
00340
00341
00342
00343 inline G4double G4VMultipleScattering::RangeFactor() const
00344 {
00345 return facrange;
00346 }
00347
00348
00349
00350 inline void G4VMultipleScattering::SetRangeFactor(G4double val)
00351 {
00352 if(val > 0.0) facrange = val;
00353 }
00354
00355
00356
00357 inline G4double G4VMultipleScattering::GeomFactor() const
00358 {
00359 return facgeom;
00360 }
00361
00362
00363
00364 inline void G4VMultipleScattering::SetGeomFactor(G4double val)
00365 {
00366 if(val > 0.0) facgeom = val;
00367 }
00368
00369
00370
00371 inline G4double G4VMultipleScattering::PolarAngleLimit() const
00372 {
00373 return polarAngleLimit;
00374 }
00375
00376
00377
00378 inline void G4VMultipleScattering::SetPolarAngleLimit(G4double val)
00379 {
00380 if(val < 0.0) { polarAngleLimit = 0.0; }
00381 else if(val > CLHEP::pi) { polarAngleLimit = CLHEP::pi; }
00382 else { polarAngleLimit = val; }
00383 }
00384
00385
00386
00387 inline G4MscStepLimitType G4VMultipleScattering::StepLimitType() const
00388 {
00389 return stepLimit;
00390 }
00391
00392
00393
00394 inline void G4VMultipleScattering::SetStepLimitType(G4MscStepLimitType val)
00395 {
00396 stepLimit = val;
00397 if(val == fMinimal) { facrange = 0.2; }
00398 }
00399
00400
00401
00402 inline const G4ParticleDefinition* G4VMultipleScattering::FirstParticle() const
00403 {
00404 return firstParticle;
00405 }
00406
00407
00408
00409 #endif