G4OpWLS Class Reference

#include <G4OpWLS.hh>

Inheritance diagram for G4OpWLS:

G4VDiscreteProcess G4VProcess

Public Member Functions

 G4OpWLS (const G4String &processName="OpWLS", G4ProcessType type=fOptical)
 ~G4OpWLS ()
G4bool IsApplicable (const G4ParticleDefinition &aParticleType)
G4double GetMeanFreePath (const G4Track &aTrack, G4double, G4ForceCondition *)
G4VParticleChangePostStepDoIt (const G4Track &aTrack, const G4Step &aStep)
G4PhysicsTableGetIntegralTable () const
void DumpPhysicsTable () const
void UseTimeProfile (const G4String name)

Protected Attributes

G4VWLSTimeGeneratorProfileWLSTimeGeneratorProfile
G4PhysicsTabletheIntegralTable

Detailed Description

Definition at line 80 of file G4OpWLS.hh.


Constructor & Destructor Documentation

G4OpWLS::G4OpWLS ( const G4String processName = "OpWLS",
G4ProcessType  type = fOptical 
)

Definition at line 64 of file G4OpWLS.cc.

References fOpWLS, G4cout, G4endl, G4VProcess::GetProcessName(), G4VProcess::SetProcessSubType(), theIntegralTable, G4VProcess::verboseLevel, and WLSTimeGeneratorProfile.

00065   : G4VDiscreteProcess(processName, type)
00066 {
00067   SetProcessSubType(fOpWLS);
00068 
00069   theIntegralTable = 0;
00070  
00071   if (verboseLevel>0) {
00072     G4cout << GetProcessName() << " is created " << G4endl;
00073   }
00074 
00075   WLSTimeGeneratorProfile = 
00076        new G4WLSTimeGeneratorProfileDelta("WLSTimeGeneratorProfileDelta");
00077 
00078   BuildThePhysicsTable();
00079 }

G4OpWLS::~G4OpWLS (  ) 

Definition at line 85 of file G4OpWLS.cc.

References G4PhysicsTable::clearAndDestroy(), theIntegralTable, and WLSTimeGeneratorProfile.

00086 {
00087   if (theIntegralTable != 0) {
00088     theIntegralTable->clearAndDestroy();
00089     delete theIntegralTable;
00090   }
00091   delete WLSTimeGeneratorProfile;
00092 }


Member Function Documentation

void G4OpWLS::DumpPhysicsTable (  )  const [inline]

Definition at line 161 of file G4OpWLS.hh.

References G4PhysicsOrderedFreeVector::DumpValues(), G4PhysicsTable::entries(), and theIntegralTable.

00162 {
00163   G4int PhysicsTableSize = theIntegralTable->entries();
00164   G4PhysicsOrderedFreeVector *v;
00165  
00166   for (G4int i = 0 ; i < PhysicsTableSize ; i++ )
00167     {
00168       v = (G4PhysicsOrderedFreeVector*)(*theIntegralTable)[i];
00169       v->DumpValues();
00170     }
00171 }

G4PhysicsTable * G4OpWLS::GetIntegralTable (  )  const [inline]

Definition at line 155 of file G4OpWLS.hh.

References theIntegralTable.

00156 {
00157   return theIntegralTable;
00158 }

G4double G4OpWLS::GetMeanFreePath ( const G4Track aTrack,
G4double  ,
G4ForceCondition  
) [virtual]

Implements G4VDiscreteProcess.

Definition at line 351 of file G4OpWLS.cc.

References DBL_MAX, G4Track::GetDynamicParticle(), G4Track::GetMaterial(), G4Material::GetMaterialPropertiesTable(), and G4DynamicParticle::GetTotalEnergy().

00354 {
00355   const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle();
00356   const G4Material* aMaterial = aTrack.GetMaterial();
00357 
00358   G4double thePhotonEnergy = aParticle->GetTotalEnergy();
00359 
00360   G4MaterialPropertiesTable* aMaterialPropertyTable;
00361   G4MaterialPropertyVector* AttenuationLengthVector;
00362         
00363   G4double AttenuationLength = DBL_MAX;
00364 
00365   aMaterialPropertyTable = aMaterial->GetMaterialPropertiesTable();
00366 
00367   if ( aMaterialPropertyTable ) {
00368     AttenuationLengthVector = aMaterialPropertyTable->
00369       GetProperty("WLSABSLENGTH");
00370     if ( AttenuationLengthVector ){
00371       AttenuationLength = AttenuationLengthVector->
00372         Value(thePhotonEnergy);
00373     }
00374     else {
00375       //             G4cout << "No WLS absorption length specified" << G4endl;
00376     }
00377   }
00378   else {
00379     //           G4cout << "No WLS absortion length specified" << G4endl;
00380   }
00381   
00382   return AttenuationLength;
00383 }

G4bool G4OpWLS::IsApplicable ( const G4ParticleDefinition aParticleType  )  [inline, virtual]

Reimplemented from G4VProcess.

Definition at line 149 of file G4OpWLS.hh.

References G4OpticalPhoton::OpticalPhoton().

00150 {
00151    return ( &aParticleType == G4OpticalPhoton::OpticalPhoton() );
00152 }

G4VParticleChange * G4OpWLS::PostStepDoIt ( const G4Track aTrack,
const G4Step aStep 
) [virtual]

Reimplemented from G4VDiscreteProcess.

Definition at line 102 of file G4OpWLS.cc.

References G4ParticleChange::AddSecondary(), G4VProcess::aParticleChange, G4MaterialPropertiesTable::ConstPropertyExists(), fStopAndKill, G4cout, G4endl, G4Poisson(), G4UniformRand, G4VWLSTimeGeneratorProfile::GenerateTime(), G4PhysicsOrderedFreeVector::GetEnergy(), G4StepPoint::GetGlobalTime(), G4Material::GetIndex(), G4Track::GetMaterial(), G4Material::GetMaterialPropertiesTable(), G4PhysicsOrderedFreeVector::GetMaxValue(), G4VParticleChange::GetNumberOfSecondaries(), G4StepPoint::GetPosition(), G4Step::GetPostStepPoint(), G4MaterialPropertiesTable::GetProperty(), G4ParticleChange::Initialize(), ns, G4OpticalPhoton::OpticalPhoton(), G4VDiscreteProcess::PostStepDoIt(), G4VParticleChange::ProposeTrackStatus(), G4DynamicParticle::SetKineticEnergy(), G4VParticleChange::SetNumberOfSecondaries(), G4Track::SetParentID(), G4DynamicParticle::SetPolarization(), G4Track::SetTouchableHandle(), theIntegralTable, G4VProcess::verboseLevel, and WLSTimeGeneratorProfile.

00103 {
00104   aParticleChange.Initialize(aTrack);
00105   
00106   aParticleChange.ProposeTrackStatus(fStopAndKill);
00107 
00108   if (verboseLevel>0) {
00109     G4cout << "\n** Photon absorbed! **" << G4endl;
00110   }
00111   
00112   const G4Material* aMaterial = aTrack.GetMaterial();
00113 
00114   G4StepPoint* pPostStepPoint = aStep.GetPostStepPoint();
00115     
00116   G4MaterialPropertiesTable* aMaterialPropertiesTable =
00117     aMaterial->GetMaterialPropertiesTable();
00118   if (!aMaterialPropertiesTable)
00119     return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
00120 
00121   const G4MaterialPropertyVector* WLS_Intensity = 
00122     aMaterialPropertiesTable->GetProperty("WLSCOMPONENT"); 
00123 
00124   if (!WLS_Intensity)
00125     return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
00126 
00127   G4int NumPhotons = 1;
00128 
00129   if (aMaterialPropertiesTable->ConstPropertyExists("WLSMEANNUMBERPHOTONS")) {
00130 
00131      G4double MeanNumberOfPhotons = aMaterialPropertiesTable->
00132                                     GetConstProperty("WLSMEANNUMBERPHOTONS");
00133 
00134      NumPhotons = G4int(G4Poisson(MeanNumberOfPhotons));
00135 
00136      if (NumPhotons <= 0) {
00137 
00138         // return unchanged particle and no secondaries
00139 
00140         aParticleChange.SetNumberOfSecondaries(0);
00141 
00142         return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
00143 
00144      }
00145 
00146   }
00147 
00148   aParticleChange.SetNumberOfSecondaries(NumPhotons);
00149 
00150   G4int materialIndex = aMaterial->GetIndex();
00151 
00152   // Retrieve the WLS Integral for this material
00153   // new G4PhysicsOrderedFreeVector allocated to hold CII's
00154 
00155   G4double WLSTime = 0.*ns;
00156   G4PhysicsOrderedFreeVector* WLSIntegral = 0;
00157 
00158   WLSTime   = aMaterialPropertiesTable->
00159     GetConstProperty("WLSTIMECONSTANT");
00160   WLSIntegral =
00161     (G4PhysicsOrderedFreeVector*)((*theIntegralTable)(materialIndex));
00162    
00163   // Max WLS Integral
00164   
00165   G4double CIImax = WLSIntegral->GetMaxValue();
00166   
00167   for (G4int i = 0; i < NumPhotons; i++) {
00168     
00169     // Determine photon energy
00170     
00171     G4double CIIvalue = G4UniformRand()*CIImax;
00172     G4double sampledEnergy = 
00173       WLSIntegral->GetEnergy(CIIvalue);
00174     
00175     if (verboseLevel>1) {
00176       G4cout << "sampledEnergy = " << sampledEnergy << G4endl;
00177       G4cout << "CIIvalue =        " << CIIvalue << G4endl;
00178     }
00179     
00180     // Generate random photon direction
00181     
00182     G4double cost = 1. - 2.*G4UniformRand();
00183     G4double sint = std::sqrt((1.-cost)*(1.+cost));
00184 
00185     G4double phi = twopi*G4UniformRand();
00186     G4double sinp = std::sin(phi);
00187     G4double cosp = std::cos(phi);
00188     
00189     G4double px = sint*cosp;
00190     G4double py = sint*sinp;
00191     G4double pz = cost;
00192     
00193     // Create photon momentum direction vector
00194     
00195     G4ParticleMomentum photonMomentum(px, py, pz);
00196     
00197     // Determine polarization of new photon
00198     
00199     G4double sx = cost*cosp;
00200     G4double sy = cost*sinp;
00201     G4double sz = -sint;
00202     
00203     G4ThreeVector photonPolarization(sx, sy, sz);
00204     
00205     G4ThreeVector perp = photonMomentum.cross(photonPolarization);
00206     
00207     phi = twopi*G4UniformRand();
00208     sinp = std::sin(phi);
00209     cosp = std::cos(phi);
00210     
00211     photonPolarization = cosp * photonPolarization + sinp * perp;
00212     
00213     photonPolarization = photonPolarization.unit();
00214     
00215     // Generate a new photon:
00216     
00217     G4DynamicParticle* aWLSPhoton =
00218       new G4DynamicParticle(G4OpticalPhoton::OpticalPhoton(),
00219                             photonMomentum);
00220     aWLSPhoton->SetPolarization
00221       (photonPolarization.x(),
00222        photonPolarization.y(),
00223        photonPolarization.z());
00224     
00225     aWLSPhoton->SetKineticEnergy(sampledEnergy);
00226     
00227     // Generate new G4Track object:
00228     
00229     // Must give position of WLS optical photon
00230 
00231     G4double TimeDelay = WLSTimeGeneratorProfile->GenerateTime(WLSTime);
00232     G4double aSecondaryTime = (pPostStepPoint->GetGlobalTime()) + TimeDelay;
00233 
00234     G4ThreeVector aSecondaryPosition = pPostStepPoint->GetPosition();
00235 
00236     G4Track* aSecondaryTrack = 
00237       new G4Track(aWLSPhoton,aSecondaryTime,aSecondaryPosition);
00238    
00239     aSecondaryTrack->SetTouchableHandle(aTrack.GetTouchableHandle()); 
00240     // aSecondaryTrack->SetTouchableHandle((G4VTouchable*)0);
00241     
00242     aSecondaryTrack->SetParentID(aTrack.GetTrackID());
00243     
00244     aParticleChange.AddSecondary(aSecondaryTrack);
00245   }
00246 
00247   if (verboseLevel>0) {
00248     G4cout << "\n Exiting from G4OpWLS::DoIt -- NumberOfSecondaries = " 
00249            << aParticleChange.GetNumberOfSecondaries() << G4endl;  
00250   }
00251   
00252   return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
00253 }

void G4OpWLS::UseTimeProfile ( const G4String  name  ) 

Definition at line 385 of file G4OpWLS.cc.

References FatalException, G4Exception(), and WLSTimeGeneratorProfile.

Referenced by G4OpticalPhysics::ConstructProcess(), and G4OpticalPhysics::SetWLSTimeProfile().

00386 {
00387   if (name == "delta")
00388     {
00389       delete WLSTimeGeneratorProfile;
00390       WLSTimeGeneratorProfile = 
00391              new G4WLSTimeGeneratorProfileDelta("delta");
00392     }
00393   else if (name == "exponential")
00394     {
00395       delete WLSTimeGeneratorProfile;
00396       WLSTimeGeneratorProfile =
00397              new G4WLSTimeGeneratorProfileExponential("exponential");
00398     }
00399   else
00400     {
00401       G4Exception("G4OpWLS::UseTimeProfile", "em0202",
00402                   FatalException,
00403                   "generator does not exist");
00404     }
00405 }


Field Documentation

G4PhysicsTable* G4OpWLS::theIntegralTable [protected]

Definition at line 140 of file G4OpWLS.hh.

Referenced by DumpPhysicsTable(), G4OpWLS(), GetIntegralTable(), PostStepDoIt(), and ~G4OpWLS().

G4VWLSTimeGeneratorProfile* G4OpWLS::WLSTimeGeneratorProfile [protected]

Definition at line 139 of file G4OpWLS.hh.

Referenced by G4OpWLS(), PostStepDoIt(), UseTimeProfile(), and ~G4OpWLS().


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:52:47 2013 for Geant4 by  doxygen 1.4.7