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 #include "G4OpticalPhysics.hh"
00040
00041 #include "G4LossTableManager.hh"
00042 #include "G4EmSaturation.hh"
00043
00044
00045
00046 #include "G4PhysicsConstructorFactory.hh"
00047
00048 G4_DECLARE_PHYSCONSTR_FACTORY(G4OpticalPhysics);
00049
00050
00051
00052
00053 G4OpticalPhysics::G4OpticalPhysics(G4int verbose, const G4String& name)
00054 : G4VPhysicsConstructor(name),
00055
00056 wasActivated(false),
00057
00058 fScintillationProcess(NULL),
00059 fCerenkovProcess(NULL),
00060 fOpWLSProcess(NULL),
00061 fOpAbsorptionProcess(NULL),
00062 fOpRayleighScatteringProcess(NULL),
00063 fOpMieHGScatteringProcess(NULL),
00064 fOpBoundaryProcess(NULL),
00065 fMaxNumPhotons(100),
00066 fMaxBetaChange(10.0),
00067 fYieldFactor(1.),
00068 fExcitationRatio(0.0),
00069 fProfile("delta"),
00070 fFiniteRiseTime(false),
00071 fScintillationByParticleType(false)
00072 {
00073 verboseLevel = verbose;
00074 fMessenger = new G4OpticalPhysicsMessenger(this);
00075
00076 for ( G4int i=0; i<kNoProcess; i++ ) {
00077 fProcesses.push_back(NULL);
00078 fProcessUse.push_back(true);
00079 fProcessVerbose.push_back(verbose);
00080 fProcessTrackSecondariesFirst.push_back(true);
00081 }
00082 }
00083
00084
00085
00086 G4OpticalPhysics::~G4OpticalPhysics()
00087 {
00088 delete fMessenger;
00089
00090 if (wasActivated) {
00091
00092 if (fScintillationProcess) delete fScintillationProcess;
00093 if (fCerenkovProcess) delete fCerenkovProcess;
00094 if (fOpWLSProcess) delete fOpWLSProcess;
00095
00096 if (fOpAbsorptionProcess) delete fOpAbsorptionProcess;
00097 if (fOpRayleighScatteringProcess) delete fOpRayleighScatteringProcess;
00098 if (fOpMieHGScatteringProcess) delete fOpMieHGScatteringProcess;
00099 if (fOpBoundaryProcess) delete fOpBoundaryProcess;
00100
00101 }
00102 }
00103
00104
00105
00106 void G4OpticalPhysics::PrintStatistics() const
00107 {
00108
00109
00110 for ( G4int i=0; i<kNoProcess; i++ ) {
00111 G4cout << " " << G4OpticalProcessName(i) << " process: ";
00112 if ( ! fProcessUse[i] ) {
00113 G4cout << "not used" << G4endl;
00114 }
00115 else {
00116 G4cout << "used" << G4endl;
00117 if ( i == kCerenkov ) {
00118 G4cout << " Max number of photons per step: " << fMaxNumPhotons << G4endl;
00119 G4cout << " Max beta change per step: " << fMaxBetaChange << G4endl;
00120 if ( fProcessTrackSecondariesFirst[kCerenkov] ) G4cout << " Track secondaries first: activated" << G4endl;
00121 }
00122 if ( i == kScintillation ) {
00123 if (fScintillationByParticleType)
00124 G4cout << " Scintillation by Particle Type: activated " << G4endl;
00125 G4cout << " Yield factor: " << fYieldFactor << G4endl;
00126 G4cout << " ExcitationRatio: " << fExcitationRatio << G4endl;
00127 if ( fProcessTrackSecondariesFirst[kScintillation] ) G4cout << " Track secondaries first: activated" << G4endl;
00128 }
00129 if ( i == kWLS ) {
00130 G4cout << " WLS process time profile: " << fProfile << G4endl;
00131 }
00132 }
00133 }
00134 }
00135
00136
00137
00138 #include "G4OpticalPhoton.hh"
00139
00140 void G4OpticalPhysics::ConstructParticle()
00141 {
00143
00144
00145 G4OpticalPhoton::OpticalPhotonDefinition();
00146 }
00147
00148
00149
00150 #include "G4ParticleDefinition.hh"
00151 #include "G4ProcessManager.hh"
00152
00153 void G4OpticalPhysics::ConstructProcess()
00154 {
00155
00156
00157 if (wasActivated) return;
00158
00159 if(verboseLevel>0)
00160 G4cout <<"G4OpticalPhysics:: Add Optical Physics Processes"<< G4endl;
00161
00162
00163
00164 fProcesses[kAbsorption] = fOpAbsorptionProcess = new G4OpAbsorption();
00165 fProcesses[kRayleigh] = fOpRayleighScatteringProcess = new G4OpRayleigh();
00166 fProcesses[kMieHG] = fOpMieHGScatteringProcess = new G4OpMieHG();
00167
00168 fProcesses[kBoundary] = fOpBoundaryProcess = new G4OpBoundaryProcess();
00169
00170 fProcesses[kWLS] = fOpWLSProcess = new G4OpWLS();
00171 fOpWLSProcess->UseTimeProfile(fProfile);
00172
00173 G4ProcessManager * pManager = 0;
00174 pManager = G4OpticalPhoton::OpticalPhoton()->GetProcessManager();
00175
00176 if (!pManager) {
00177 std::ostringstream o;
00178 o << "Optical Photon without a Process Manager";
00179 G4Exception("G4OpticalPhysics::ConstructProcess()","",
00180 FatalException,o.str().c_str());
00181 return;
00182 }
00183
00184 for ( G4int i=kAbsorption; i<=kWLS; i++ ) {
00185 if ( fProcessUse[i] ) {
00186 pManager->AddDiscreteProcess(fProcesses[i]);
00187 }
00188 }
00189
00190 fProcesses[kScintillation] = fScintillationProcess = new G4Scintillation();
00191 fScintillationProcess->SetScintillationYieldFactor(fYieldFactor);
00192 fScintillationProcess->SetScintillationExcitationRatio(fExcitationRatio);
00193 fScintillationProcess->SetFiniteRiseTime(fFiniteRiseTime);
00194 fScintillationProcess->SetScintillationByParticleType(fScintillationByParticleType);
00195 fScintillationProcess->
00196 SetTrackSecondariesFirst(fProcessTrackSecondariesFirst[kScintillation]);
00197
00198
00199
00200 G4EmSaturation* emSaturation = G4LossTableManager::Instance()->EmSaturation();
00201 fScintillationProcess->AddSaturation(emSaturation);
00202
00203 fProcesses[kCerenkov] = fCerenkovProcess = new G4Cerenkov();
00204 fCerenkovProcess->SetMaxNumPhotonsPerStep(fMaxNumPhotons);
00205 fCerenkovProcess->SetMaxBetaChangePerStep(fMaxBetaChange);
00206 fCerenkovProcess->
00207 SetTrackSecondariesFirst(fProcessTrackSecondariesFirst[kCerenkov]);
00208
00209 theParticleIterator->reset();
00210
00211 while( (*theParticleIterator)() ){
00212
00213 G4ParticleDefinition* particle = theParticleIterator->value();
00214 G4String particleName = particle->GetParticleName();
00215
00216 pManager = particle->GetProcessManager();
00217 if (!pManager) {
00218 std::ostringstream o;
00219 o << "Particle " << particleName << "without a Process Manager";
00220 G4Exception("G4OpticalPhysics::ConstructProcess()","",
00221 FatalException,o.str().c_str());
00222 return;
00223 }
00224
00225 if( fCerenkovProcess->IsApplicable(*particle) &&
00226 fProcessUse[kCerenkov] ) {
00227 pManager->AddProcess(fCerenkovProcess);
00228 pManager->SetProcessOrdering(fCerenkovProcess,idxPostStep);
00229 }
00230 if( fScintillationProcess->IsApplicable(*particle) &&
00231 fProcessUse[kScintillation] ){
00232 pManager->AddProcess(fScintillationProcess);
00233 pManager->SetProcessOrderingToLast(fScintillationProcess,idxAtRest);
00234 pManager->SetProcessOrderingToLast(fScintillationProcess,idxPostStep);
00235 }
00236
00237 }
00238
00239
00240 for ( G4int i=0; i<kNoProcess; i++ ) {
00241 fProcesses[i]->SetVerboseLevel(fProcessVerbose[i]);
00242 }
00243
00244 if (verboseLevel > 1) PrintStatistics();
00245 if (verboseLevel > 0)
00246 G4cout << "### " << namePhysics << " physics constructed." << G4endl;
00247
00248 wasActivated = true;
00249
00250 }
00251
00252 void G4OpticalPhysics::SetScintillationYieldFactor(G4double yieldFactor)
00253 {
00255
00256 fYieldFactor = yieldFactor;
00257
00258 if(fScintillationProcess)
00259 fScintillationProcess->SetScintillationYieldFactor(yieldFactor);
00260 }
00261
00262 void G4OpticalPhysics::SetScintillationExcitationRatio(G4double excitationRatio)
00263 {
00265
00266 fExcitationRatio = excitationRatio;
00267
00268 if(fScintillationProcess)
00269 fScintillationProcess->SetScintillationExcitationRatio(excitationRatio);
00270 }
00271
00272 void G4OpticalPhysics::SetMaxNumPhotonsPerStep(G4int maxNumPhotons)
00273 {
00275
00276 fMaxNumPhotons = maxNumPhotons;
00277
00278 if(fCerenkovProcess)
00279 fCerenkovProcess->SetMaxNumPhotonsPerStep(maxNumPhotons);
00280 }
00281
00282 void G4OpticalPhysics::SetMaxBetaChangePerStep(G4double maxBetaChange)
00283 {
00285
00286 fMaxBetaChange = maxBetaChange;
00287
00288 if(fCerenkovProcess)
00289 fCerenkovProcess->SetMaxBetaChangePerStep(maxBetaChange);
00290 }
00291
00292 void G4OpticalPhysics::SetWLSTimeProfile(G4String profile)
00293 {
00295
00296 fProfile = profile;
00297
00298 if(fOpWLSProcess)
00299 fOpWLSProcess->UseTimeProfile(fProfile);
00300 }
00301
00302 void G4OpticalPhysics::AddScintillationSaturation(G4EmSaturation* saturation)
00303 {
00305
00306 if(fScintillationProcess)
00307 fScintillationProcess->AddSaturation(saturation);
00308 }
00309
00310 void G4OpticalPhysics::SetScintillationByParticleType(G4bool scintillationByParticleType)
00311 {
00312 fScintillationByParticleType = scintillationByParticleType;
00313
00314 if (fScintillationProcess)
00315 fScintillationProcess->SetScintillationByParticleType(scintillationByParticleType);
00316 }
00317
00318 void G4OpticalPhysics::SetTrackSecondariesFirst(G4OpticalProcessIndex index,
00319 G4bool trackSecondariesFirst)
00320 {
00321 if ( index >= kNoProcess ) return;
00322 if ( fProcessTrackSecondariesFirst[index] == trackSecondariesFirst ) return;
00323 fProcessTrackSecondariesFirst[index] = trackSecondariesFirst;
00324
00325 if(fCerenkovProcess && index == kCerenkov )
00326 fCerenkovProcess->SetTrackSecondariesFirst(trackSecondariesFirst);
00327
00328 if(fScintillationProcess && index == kScintillation)
00329 fScintillationProcess->SetTrackSecondariesFirst(trackSecondariesFirst);
00330 }
00331
00332 void G4OpticalPhysics::SetFiniteRiseTime(G4bool finiteRiseTime)
00333 {
00334 fFiniteRiseTime = finiteRiseTime;
00335 if(fScintillationProcess)
00336 fScintillationProcess->SetFiniteRiseTime(finiteRiseTime);
00337 }
00338
00339 void G4OpticalPhysics::Configure(G4OpticalProcessIndex index, G4bool isUse)
00340 {
00341
00342
00343
00344
00345
00346
00347 if ( index >= kNoProcess ) return;
00348 if ( fProcessUse[index] == isUse ) return;
00349 fProcessUse[index] = isUse;
00350 }
00351
00352 void G4OpticalPhysics::SetProcessVerbose(G4int index,
00353 G4int inputVerboseLevel)
00354 {
00355
00356
00357 if ( index >= kNoProcess ) return;
00358 if ( fProcessVerbose[index] == inputVerboseLevel ) return;
00359
00360 fProcessVerbose[index] = inputVerboseLevel;
00361
00362 if ( fProcesses[index] ) fProcesses[index]->SetVerboseLevel(inputVerboseLevel);
00363 }
00364
00365