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
00031
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049
00050
00051
00052
00054
00055 #ifndef G4Scintillation_h
00056 #define G4Scintillation_h 1
00057
00059
00061
00062 #include "globals.hh"
00063 #include "templates.hh"
00064 #include "Randomize.hh"
00065 #include "G4Poisson.hh"
00066 #include "G4ThreeVector.hh"
00067 #include "G4ParticleMomentum.hh"
00068 #include "G4Step.hh"
00069 #include "G4VRestDiscreteProcess.hh"
00070 #include "G4OpticalPhoton.hh"
00071 #include "G4DynamicParticle.hh"
00072 #include "G4Material.hh"
00073 #include "G4PhysicsTable.hh"
00074 #include "G4MaterialPropertiesTable.hh"
00075 #include "G4PhysicsOrderedFreeVector.hh"
00076
00077 #include "G4EmSaturation.hh"
00078
00079
00080
00081
00082
00083
00085
00087
00088 class G4Scintillation : public G4VRestDiscreteProcess
00089 {
00090
00091 public:
00092
00094
00096
00097 G4Scintillation(const G4String& processName = "Scintillation",
00098 G4ProcessType type = fElectromagnetic);
00099 ~G4Scintillation();
00100
00101 private:
00102
00103 G4Scintillation(const G4Scintillation &right);
00104
00106
00108
00109 G4Scintillation& operator=(const G4Scintillation &right);
00110
00111 public:
00112
00114
00116
00117
00118
00119
00120
00121 G4bool IsApplicable(const G4ParticleDefinition& aParticleType);
00122
00123
00124
00125 G4double GetMeanFreePath(const G4Track& aTrack,
00126 G4double ,
00127 G4ForceCondition* );
00128
00129
00130
00131
00132 G4double GetMeanLifeTime(const G4Track& aTrack,
00133 G4ForceCondition* );
00134
00135
00136
00137
00138 G4VParticleChange* PostStepDoIt(const G4Track& aTrack,
00139 const G4Step& aStep);
00140 G4VParticleChange* AtRestDoIt (const G4Track& aTrack,
00141 const G4Step& aStep);
00142
00143
00144
00145 void SetTrackSecondariesFirst(const G4bool state);
00146
00147
00148
00149
00150 void SetFiniteRiseTime(const G4bool state);
00151
00152
00153
00154 G4bool GetTrackSecondariesFirst() const;
00155
00156
00157 G4bool GetFiniteRiseTime() const;
00158
00159
00160 void SetScintillationYieldFactor(const G4double yieldfactor);
00161
00162
00163
00164
00165 G4double GetScintillationYieldFactor() const;
00166
00167
00168 void SetScintillationExcitationRatio(const G4double excitationratio);
00169
00170
00171
00172
00173
00174 G4double GetScintillationExcitationRatio() const;
00175
00176
00177 G4PhysicsTable* GetFastIntegralTable() const;
00178
00179
00180 G4PhysicsTable* GetSlowIntegralTable() const;
00181
00182
00183 void AddSaturation(G4EmSaturation* sat) { emSaturation = sat; }
00184
00185
00186 void RemoveSaturation() { emSaturation = NULL; }
00187
00188
00189 G4EmSaturation* GetSaturation() const { return emSaturation; }
00190
00191
00192 void SetScintillationByParticleType(const G4bool );
00193
00194
00195
00196 G4bool GetScintillationByParticleType() const
00197 { return scintillationByParticleType; }
00198
00199
00200
00201 void DumpPhysicsTable() const;
00202
00203
00204 protected:
00205
00206 void BuildThePhysicsTable();
00207
00208
00209
00211
00213
00214 G4PhysicsTable* theSlowIntegralTable;
00215 G4PhysicsTable* theFastIntegralTable;
00216
00217 G4bool fTrackSecondariesFirst;
00218 G4bool fFiniteRiseTime;
00219
00220 G4double YieldFactor;
00221
00222 G4double ExcitationRatio;
00223
00224 G4bool scintillationByParticleType;
00225
00226 private:
00227
00228 G4double single_exp(G4double t, G4double tau2);
00229 G4double bi_exp(G4double t, G4double tau1, G4double tau2);
00230
00231
00232 G4double sample_time(G4double tau1, G4double tau2);
00233
00234 G4EmSaturation* emSaturation;
00235
00236 };
00237
00239
00241
00242 inline
00243 G4bool G4Scintillation::IsApplicable(const G4ParticleDefinition& aParticleType)
00244 {
00245 if (aParticleType.GetParticleName() == "opticalphoton") return false;
00246 if (aParticleType.IsShortLived()) return false;
00247
00248 return true;
00249 }
00250
00251 inline
00252 void G4Scintillation::SetTrackSecondariesFirst(const G4bool state)
00253 {
00254 fTrackSecondariesFirst = state;
00255 }
00256
00257 inline
00258 void G4Scintillation::SetFiniteRiseTime(const G4bool state)
00259 {
00260 fFiniteRiseTime = state;
00261 }
00262
00263 inline
00264 G4bool G4Scintillation::GetTrackSecondariesFirst() const
00265 {
00266 return fTrackSecondariesFirst;
00267 }
00268
00269 inline
00270 G4bool G4Scintillation::GetFiniteRiseTime() const
00271 {
00272 return fFiniteRiseTime;
00273 }
00274
00275 inline
00276 void G4Scintillation::SetScintillationYieldFactor(const G4double yieldfactor)
00277 {
00278 YieldFactor = yieldfactor;
00279 }
00280
00281 inline
00282 G4double G4Scintillation::GetScintillationYieldFactor() const
00283 {
00284 return YieldFactor;
00285 }
00286
00287 inline
00288 void G4Scintillation::SetScintillationExcitationRatio(const G4double excitationratio)
00289 {
00290 ExcitationRatio = excitationratio;
00291 }
00292
00293 inline
00294 G4double G4Scintillation::GetScintillationExcitationRatio() const
00295 {
00296 return ExcitationRatio;
00297 }
00298
00299 inline
00300 G4PhysicsTable* G4Scintillation::GetSlowIntegralTable() const
00301 {
00302 return theSlowIntegralTable;
00303 }
00304
00305 inline
00306 G4PhysicsTable* G4Scintillation::GetFastIntegralTable() const
00307 {
00308 return theFastIntegralTable;
00309 }
00310
00311 inline
00312 void G4Scintillation::DumpPhysicsTable() const
00313 {
00314 if (theFastIntegralTable) {
00315 G4int PhysicsTableSize = theFastIntegralTable->entries();
00316 G4PhysicsOrderedFreeVector *v;
00317
00318 for (G4int i = 0 ; i < PhysicsTableSize ; i++ )
00319 {
00320 v = (G4PhysicsOrderedFreeVector*)(*theFastIntegralTable)[i];
00321 v->DumpValues();
00322 }
00323 }
00324
00325 if (theSlowIntegralTable) {
00326 G4int PhysicsTableSize = theSlowIntegralTable->entries();
00327 G4PhysicsOrderedFreeVector *v;
00328
00329 for (G4int i = 0 ; i < PhysicsTableSize ; i++ )
00330 {
00331 v = (G4PhysicsOrderedFreeVector*)(*theSlowIntegralTable)[i];
00332 v->DumpValues();
00333 }
00334 }
00335 }
00336
00337 inline
00338 G4double G4Scintillation::single_exp(G4double t, G4double tau2)
00339 {
00340 return std::exp(-1.0*t/tau2)/tau2;
00341 }
00342
00343 inline
00344 G4double G4Scintillation::bi_exp(G4double t, G4double tau1, G4double tau2)
00345 {
00346 return std::exp(-1.0*t/tau2)*(1-std::exp(-1.0*t/tau1))/tau2/tau2*(tau1+tau2);
00347 }
00348
00349 #endif