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 #ifndef G4VProcess_h
00057 #define G4VProcess_h 1
00058
00059 #include "globals.hh"
00060 #include <cmath>
00061 #include "G4ios.hh"
00062
00063 class G4ParticleDefinition;
00064 class G4DynamicParticle;
00065 class G4Track;
00066 class G4Step;
00067
00068 #include "G4PhysicsTable.hh"
00069 #include "G4VParticleChange.hh"
00070 #include "G4ForceCondition.hh"
00071 #include "G4GPILSelection.hh"
00072 #include "G4ParticleChange.hh"
00073 #include "G4ProcessType.hh"
00074
00075 class G4VProcess
00076 {
00077
00078
00079
00080
00081 private:
00082
00083
00084
00085 G4VProcess & operator=(const G4VProcess &right);
00086
00087 public:
00088
00089 G4VProcess(const G4String& aName = "NoName",
00090 G4ProcessType aType = fNotDefined );
00091
00092
00093
00094 G4VProcess(const G4VProcess &right);
00095
00096 public:
00097
00098 virtual ~G4VProcess();
00099
00100
00101 G4int operator==(const G4VProcess &right) const;
00102 G4int operator!=(const G4VProcess &right) const;
00103
00104 public:
00106
00108 virtual G4VParticleChange* PostStepDoIt(
00109 const G4Track& track,
00110 const G4Step& stepData
00111 ) = 0;
00112
00113 virtual G4VParticleChange* AlongStepDoIt(
00114 const G4Track& track,
00115 const G4Step& stepData
00116 ) = 0;
00117 virtual G4VParticleChange* AtRestDoIt(
00118 const G4Track& track,
00119 const G4Step& stepData
00120 ) = 0;
00121
00122
00123
00124
00125
00126
00127
00128
00129
00130
00132
00134 virtual G4double AlongStepGetPhysicalInteractionLength(
00135 const G4Track& track,
00136 G4double previousStepSize,
00137 G4double currentMinimumStep,
00138 G4double& proposedSafety,
00139 G4GPILSelection* selection) = 0;
00140
00141 virtual G4double AtRestGetPhysicalInteractionLength(
00142 const G4Track& track,
00143 G4ForceCondition* condition
00144 ) = 0;
00145
00146 virtual G4double PostStepGetPhysicalInteractionLength(
00147 const G4Track& track,
00148 G4double previousStepSize,
00149 G4ForceCondition* condition
00150 ) = 0;
00151
00152
00153
00154
00155
00156
00157
00158
00159
00160
00161
00162
00163
00164
00165
00166
00167
00168
00169
00170
00171
00172
00173
00174
00175
00176
00177
00178 G4double GetCurrentInteractionLength() const;
00179
00180
00182 void SetPILfactor(G4double value);
00183 G4double GetPILfactor() const;
00184
00185
00186
00187
00188
00189
00190
00191 G4double AlongStepGPIL( const G4Track& track,
00192 G4double previousStepSize,
00193 G4double currentMinimumStep,
00194 G4double& proposedSafety,
00195 G4GPILSelection* selection );
00196
00197 G4double AtRestGPIL( const G4Track& track,
00198 G4ForceCondition* condition );
00199
00200 G4double PostStepGPIL( const G4Track& track,
00201 G4double previousStepSize,
00202 G4ForceCondition* condition );
00203
00205 virtual G4bool IsApplicable(const G4ParticleDefinition&){return true;}
00206
00207
00208
00209
00210 virtual void BuildPhysicsTable(const G4ParticleDefinition&){}
00211
00212
00213
00214
00215
00216
00217 virtual void PreparePhysicsTable(const G4ParticleDefinition&){}
00218
00219
00220
00221
00222
00223
00224
00225
00226
00227
00228
00229
00230
00231 virtual G4bool StorePhysicsTable(const G4ParticleDefinition* ,
00232 const G4String&, G4bool){return true;}
00233
00234
00235
00236 virtual G4bool RetrievePhysicsTable( const G4ParticleDefinition* ,
00237 const G4String&, G4bool){return false;}
00238
00239
00240
00241
00242
00243 const G4String& GetPhysicsTableFileName(const G4ParticleDefinition* ,
00244 const G4String& directory,
00245 const G4String& tableName,
00246 G4bool ascii =false);
00247
00248
00250 const G4String& GetProcessName() const;
00251
00252
00253 G4ProcessType GetProcessType() const;
00254
00255
00256 void SetProcessType(G4ProcessType );
00257
00258
00259 G4int GetProcessSubType() const;
00260
00261
00262 void SetProcessSubType(G4int );
00263
00264
00265 static const G4String& GetProcessTypeName(G4ProcessType );
00266
00267
00268 virtual void StartTracking(G4Track*);
00269 virtual void EndTracking();
00270
00271
00272 public:
00273 virtual void SetProcessManager(const G4ProcessManager*);
00274
00275
00276 virtual const G4ProcessManager* GetProcessManager();
00277
00278
00279 protected:
00280 const G4ProcessManager* aProcessManager;
00281
00282 protected:
00283 G4VParticleChange* pParticleChange;
00284
00285
00286
00287
00288
00289 G4ParticleChange aParticleChange;
00290
00291
00292
00293 G4double theNumberOfInteractionLengthLeft;
00294
00295
00296
00297 G4double currentInteractionLength;
00298
00299
00300 G4double theInitialNumberOfInteractionLength;
00301
00302
00303 public:
00304 virtual void ResetNumberOfInteractionLengthLeft();
00305
00306
00307 G4double GetNumberOfInteractionLengthLeft() const;
00308
00309
00310 G4double GetTotalNumberOfInteractionLengthTraversed() const;
00311
00312
00313
00314 protected:
00315 void SubtractNumberOfInteractionLengthLeft(
00316 G4double previousStepSize
00317 );
00318
00319
00320
00321 void ClearNumberOfInteractionLengthLeft();
00322
00323
00324
00325
00326 public:
00327
00328
00329
00330 G4bool isAtRestDoItIsEnabled() const;
00331 G4bool isAlongStepDoItIsEnabled() const;
00332 G4bool isPostStepDoItIsEnabled() const;
00333
00334 protected:
00335 G4String theProcessName;
00336
00337
00338 G4String thePhysicsTableFileName;
00339
00340 G4ProcessType theProcessType;
00341
00342
00343 G4int theProcessSubType;
00344
00345
00346 G4double thePILfactor;
00347
00348
00349
00350 G4bool enableAtRestDoIt;
00351 G4bool enableAlongStepDoIt;
00352 G4bool enablePostStepDoIt;
00353
00354 public:
00355 virtual void DumpInfo() const;
00356
00357
00358 public:
00359 void SetVerboseLevel(G4int value);
00360 G4int GetVerboseLevel() const;
00361
00362
00363
00364
00365
00366
00367 protected:
00368 G4int verboseLevel;
00369
00370
00371 };
00372
00373
00374
00375
00376 #include "Randomize.hh"
00377
00378 inline
00379 const G4String& G4VProcess::GetProcessName() const
00380 {
00381 return theProcessName;
00382 }
00383
00384 inline
00385 G4ProcessType G4VProcess::GetProcessType() const
00386 {
00387 return theProcessType;
00388 }
00389
00390 inline
00391 void G4VProcess::SetProcessType(G4ProcessType aType)
00392 {
00393 theProcessType = aType;
00394 }
00395
00396 inline
00397 G4int G4VProcess::GetProcessSubType() const
00398 {
00399 return theProcessSubType;
00400 }
00401
00402 inline
00403 void G4VProcess::SetProcessSubType(G4int value)
00404 {
00405 theProcessSubType = value;
00406 }
00407
00408 inline void G4VProcess::SetVerboseLevel(G4int value)
00409 {
00410 verboseLevel = value;
00411 }
00412
00413 inline G4int G4VProcess::GetVerboseLevel() const
00414 {
00415 return verboseLevel;
00416 }
00417
00418 inline void G4VProcess::ClearNumberOfInteractionLengthLeft()
00419 {
00420 theInitialNumberOfInteractionLength = -1.0;
00421 theNumberOfInteractionLengthLeft = -1.0;
00422 }
00423
00424 inline G4double G4VProcess::GetNumberOfInteractionLengthLeft() const
00425 {
00426 return theNumberOfInteractionLengthLeft;
00427 }
00428
00429 inline G4double G4VProcess::GetTotalNumberOfInteractionLengthTraversed() const
00430 {
00431 return theInitialNumberOfInteractionLength - theNumberOfInteractionLengthLeft;}
00432
00433 inline G4double G4VProcess::GetCurrentInteractionLength() const
00434 {
00435 return currentInteractionLength;
00436 }
00437
00438 inline void G4VProcess::SetPILfactor(G4double value)
00439 {
00440 if (value>0.) {
00441 thePILfactor = value;
00442 }
00443 }
00444
00445 inline G4double G4VProcess::GetPILfactor() const
00446 {
00447 return thePILfactor;
00448 }
00449
00450 inline G4double G4VProcess::AlongStepGPIL( const G4Track& track,
00451 G4double previousStepSize,
00452 G4double currentMinimumStep,
00453 G4double& proposedSafety,
00454 G4GPILSelection* selection )
00455 {
00456 G4double value
00457 =AlongStepGetPhysicalInteractionLength(track, previousStepSize, currentMinimumStep, proposedSafety, selection);
00458 return value;
00459 }
00460
00461 inline G4double G4VProcess::AtRestGPIL( const G4Track& track,
00462 G4ForceCondition* condition )
00463 {
00464 G4double value
00465 =AtRestGetPhysicalInteractionLength(track, condition);
00466 return thePILfactor*value;
00467 }
00468
00469 inline G4double G4VProcess::PostStepGPIL( const G4Track& track,
00470 G4double previousStepSize,
00471 G4ForceCondition* condition )
00472 {
00473 G4double value
00474 =PostStepGetPhysicalInteractionLength(track, previousStepSize, condition);
00475 return thePILfactor*value;
00476 }
00477
00478 inline
00479 void G4VProcess::SetProcessManager(const G4ProcessManager* procMan)
00480 {
00481 aProcessManager = procMan;
00482 }
00483
00484 inline
00485 const G4ProcessManager* G4VProcess::GetProcessManager()
00486 {
00487 return aProcessManager;
00488 }
00489
00490 inline
00491 G4bool G4VProcess::isAtRestDoItIsEnabled() const
00492 {
00493 return enableAtRestDoIt;
00494 }
00495
00496 inline
00497 G4bool G4VProcess::isAlongStepDoItIsEnabled() const
00498 {
00499 return enableAlongStepDoIt;
00500 }
00501
00502 inline
00503 G4bool G4VProcess::isPostStepDoItIsEnabled() const
00504 {
00505 return enablePostStepDoIt;
00506 }
00507
00508 #endif