Geant4 11.1.1
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4OpWLS2 Class Reference

#include <G4OpWLS2.hh>

+ Inheritance diagram for G4OpWLS2:

Public Member Functions

 G4OpWLS2 (const G4String &processName="OpWLS2", G4ProcessType type=fOptical)
 
virtual ~G4OpWLS2 ()
 
virtual G4bool IsApplicable (const G4ParticleDefinition &aParticleType) override
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &aParticleType) override
 
virtual G4double GetMeanFreePath (const G4Track &aTrack, G4double, G4ForceCondition *) override
 
virtual G4VParticleChangePostStepDoIt (const G4Track &aTrack, const G4Step &aStep) override
 
virtual G4PhysicsTableGetIntegralTable () const
 
virtual void DumpPhysicsTable () const
 
virtual void UseTimeProfile (const G4String name)
 
virtual void PreparePhysicsTable (const G4ParticleDefinition &) override
 
virtual void Initialise ()
 
void SetVerboseLevel (G4int)
 
- Public Member Functions inherited from G4VDiscreteProcess
 G4VDiscreteProcess (const G4String &aName, G4ProcessType aType=fNotDefined)
 
 G4VDiscreteProcess (G4VDiscreteProcess &)
 
virtual ~G4VDiscreteProcess ()
 
G4VDiscreteProcessoperator= (const G4VDiscreteProcess &)=delete
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4VParticleChangePostStepDoIt (const G4Track &, const G4Step &)
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &, const G4Step &)
 
virtual G4double GetCrossSection (const G4double, const G4MaterialCutsCouple *)
 
virtual G4double MinPrimaryEnergy (const G4ParticleDefinition *, const G4Material *)
 
- Public Member Functions inherited from G4VProcess
 G4VProcess (const G4String &aName="NoName", G4ProcessType aType=fNotDefined)
 
 G4VProcess (const G4VProcess &right)
 
virtual ~G4VProcess ()
 
G4VProcessoperator= (const G4VProcess &)=delete
 
G4bool operator== (const G4VProcess &right) const
 
G4bool operator!= (const G4VProcess &right) const
 
virtual G4VParticleChangePostStepDoIt (const G4Track &track, const G4Step &stepData)=0
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &track, const G4Step &stepData)=0
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &track, const G4Step &stepData)=0
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)=0
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &track, G4ForceCondition *condition)=0
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)=0
 
G4double GetCurrentInteractionLength () const
 
void SetPILfactor (G4double value)
 
G4double GetPILfactor () const
 
G4double AlongStepGPIL (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)
 
G4double AtRestGPIL (const G4Track &track, G4ForceCondition *condition)
 
G4double PostStepGPIL (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4bool IsApplicable (const G4ParticleDefinition &)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void PreparePhysicsTable (const G4ParticleDefinition &)
 
virtual G4bool StorePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
virtual G4bool RetrievePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
const G4StringGetPhysicsTableFileName (const G4ParticleDefinition *, const G4String &directory, const G4String &tableName, G4bool ascii=false)
 
const G4StringGetProcessName () const
 
G4ProcessType GetProcessType () const
 
void SetProcessType (G4ProcessType)
 
G4int GetProcessSubType () const
 
void SetProcessSubType (G4int)
 
virtual const G4VProcessGetCreatorProcess () const
 
virtual void StartTracking (G4Track *)
 
virtual void EndTracking ()
 
virtual void SetProcessManager (const G4ProcessManager *)
 
virtual const G4ProcessManagerGetProcessManager ()
 
virtual void ResetNumberOfInteractionLengthLeft ()
 
G4double GetNumberOfInteractionLengthLeft () const
 
G4double GetTotalNumberOfInteractionLengthTraversed () const
 
G4bool isAtRestDoItIsEnabled () const
 
G4bool isAlongStepDoItIsEnabled () const
 
G4bool isPostStepDoItIsEnabled () const
 
virtual void DumpInfo () const
 
virtual void ProcessDescription (std::ostream &outfile) const
 
void SetVerboseLevel (G4int value)
 
G4int GetVerboseLevel () const
 
virtual void SetMasterProcess (G4VProcess *masterP)
 
const G4VProcessGetMasterProcess () const
 
virtual void BuildWorkerPhysicsTable (const G4ParticleDefinition &part)
 
virtual void PrepareWorkerPhysicsTable (const G4ParticleDefinition &)
 

Protected Attributes

G4VWLSTimeGeneratorProfileWLSTimeGeneratorProfile
 
G4PhysicsTabletheIntegralTable
 
- Protected Attributes inherited from G4VProcess
const G4ProcessManageraProcessManager = nullptr
 
G4VParticleChangepParticleChange = nullptr
 
G4ParticleChange aParticleChange
 
G4double theNumberOfInteractionLengthLeft = -1.0
 
G4double currentInteractionLength = -1.0
 
G4double theInitialNumberOfInteractionLength = -1.0
 
G4String theProcessName
 
G4String thePhysicsTableFileName
 
G4ProcessType theProcessType = fNotDefined
 
G4int theProcessSubType = -1
 
G4double thePILfactor = 1.0
 
G4int verboseLevel = 0
 
G4bool enableAtRestDoIt = true
 
G4bool enableAlongStepDoIt = true
 
G4bool enablePostStepDoIt = true
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 
virtual G4double GetMeanFreePath (const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)=0
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double prevStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 

Detailed Description

Definition at line 51 of file G4OpWLS2.hh.

Constructor & Destructor Documentation

◆ G4OpWLS2()

G4OpWLS2::G4OpWLS2 ( const G4String processName = "OpWLS2",
G4ProcessType  type = fOptical 
)
explicit

Definition at line 54 of file G4OpWLS2.cc.

55 : G4VDiscreteProcess(processName, type)
56{
58 Initialise();
60 theIntegralTable = nullptr;
61
62 if(verboseLevel > 0)
63 G4cout << GetProcessName() << " is created " << G4endl;
64}
@ fOpWLS
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
virtual void Initialise()
Definition: G4OpWLS2.cc:84
G4PhysicsTable * theIntegralTable
Definition: G4OpWLS2.hh:91
G4VWLSTimeGeneratorProfile * WLSTimeGeneratorProfile
Definition: G4OpWLS2.hh:90
G4int verboseLevel
Definition: G4VProcess.hh:360
void SetProcessSubType(G4int)
Definition: G4VProcess.hh:410
const G4String & GetProcessName() const
Definition: G4VProcess.hh:386

◆ ~G4OpWLS2()

G4OpWLS2::~G4OpWLS2 ( )
virtual

Definition at line 67 of file G4OpWLS2.cc.

68{
70 {
72 delete theIntegralTable;
73 }
75}
void clearAndDestroy()

Member Function Documentation

◆ BuildPhysicsTable()

void G4OpWLS2::BuildPhysicsTable ( const G4ParticleDefinition aParticleType)
overridevirtual

Reimplemented from G4VProcess.

Definition at line 233 of file G4OpWLS2.cc.

234{
236 {
238 delete theIntegralTable;
239 theIntegralTable = nullptr;
240 }
241
242 const G4MaterialTable* materialTable = G4Material::GetMaterialTable();
243 std::size_t numOfMaterials = G4Material::GetNumberOfMaterials();
244 theIntegralTable = new G4PhysicsTable(numOfMaterials);
245
246 // loop for materials
247 for(std::size_t i = 0; i < numOfMaterials; ++i)
248 {
249 auto physVector = new G4PhysicsFreeVector();
250
251 // Retrieve vector of WLS2 wavelength intensity for
252 // the material from the material's optical properties table.
254 (*materialTable)[i]->GetMaterialPropertiesTable();
255 if(MPT)
256 {
258 if(wlsVector)
259 {
260 // Retrieve the first intensity point in vector
261 // of (photon energy, intensity) pairs
262 G4double currentIN = (*wlsVector)[0];
263 if(currentIN >= 0.0)
264 {
265 // Create first (photon energy)
266 G4double currentPM = wlsVector->Energy(0);
267 G4double currentCII = 0.0;
268 physVector->InsertValues(currentPM, currentCII);
269
270 // Set previous values to current ones prior to loop
271 G4double prevPM = currentPM;
272 G4double prevCII = currentCII;
273 G4double prevIN = currentIN;
274
275 // loop over all (photon energy, intensity)
276 // pairs stored for this material
277 for(std::size_t j = 1; j < wlsVector->GetVectorLength(); ++j)
278 {
279 currentPM = wlsVector->Energy(j);
280 currentIN = (*wlsVector)[j];
281 currentCII =
282 prevCII + 0.5 * (currentPM - prevPM) * (prevIN + currentIN);
283
284 physVector->InsertValues(currentPM, currentCII);
285
286 prevPM = currentPM;
287 prevCII = currentCII;
288 prevIN = currentIN;
289 }
290 }
291 }
292 }
293 theIntegralTable->insertAt(i, physVector);
294 }
295}
std::vector< G4Material * > G4MaterialTable
double G4double
Definition: G4Types.hh:83
G4MaterialPropertyVector * GetProperty(const char *key) const
static size_t GetNumberOfMaterials()
Definition: G4Material.cc:684
static G4MaterialTable * GetMaterialTable()
Definition: G4Material.cc:677
void insertAt(std::size_t, G4PhysicsVector *)
G4double Energy(const std::size_t index) const
std::size_t GetVectorLength() const

◆ DumpPhysicsTable()

void G4OpWLS2::DumpPhysicsTable ( ) const
inlinevirtual

Definition at line 114 of file G4OpWLS2.hh.

115{
116 std::size_t PhysicsTableSize = theIntegralTable->entries();
118
119 for(std::size_t i = 0; i < PhysicsTableSize; ++i)
120 {
122 v->DumpValues();
123 }
124}
std::size_t entries() const
void DumpValues(G4double unitE=1.0, G4double unitV=1.0) const

◆ GetIntegralTable()

G4PhysicsTable * G4OpWLS2::GetIntegralTable ( ) const
inlinevirtual

Definition at line 109 of file G4OpWLS2.hh.

110{
111 return theIntegralTable;
112}

◆ GetMeanFreePath()

G4double G4OpWLS2::GetMeanFreePath ( const G4Track aTrack,
G4double  ,
G4ForceCondition  
)
overridevirtual

Implements G4VDiscreteProcess.

Definition at line 298 of file G4OpWLS2.cc.

300{
301 G4double thePhotonEnergy = aTrack.GetDynamicParticle()->GetTotalEnergy();
302 G4double attLength = DBL_MAX;
305
306 if(MPT)
307 {
309 if(attVector)
310 {
311 attLength = attVector->Value(thePhotonEnergy, idx_wls2);
312 }
313 }
314 return attLength;
315}
G4double GetTotalEnergy() const
G4MaterialPropertiesTable * GetMaterialPropertiesTable() const
Definition: G4Material.hh:251
G4double Value(const G4double energy, std::size_t &lastidx) const
G4Material * GetMaterial() const
const G4DynamicParticle * GetDynamicParticle() const
#define DBL_MAX
Definition: templates.hh:62

◆ Initialise()

void G4OpWLS2::Initialise ( )
virtual

Definition at line 84 of file G4OpWLS2.cc.

85{
89}
void SetVerboseLevel(G4int)
Definition: G4OpWLS2.cc:343
virtual void UseTimeProfile(const G4String name)
Definition: G4OpWLS2.cc:318
G4int GetWLS2VerboseLevel() const
G4String GetWLS2TimeProfile() const
static G4OpticalParameters * Instance()

Referenced by G4OpWLS2(), and PreparePhysicsTable().

◆ IsApplicable()

G4bool G4OpWLS2::IsApplicable ( const G4ParticleDefinition aParticleType)
inlineoverridevirtual

Reimplemented from G4VProcess.

Definition at line 104 of file G4OpWLS2.hh.

105{
106 return (&aParticleType == G4OpticalPhoton::OpticalPhoton());
107}
static G4OpticalPhoton * OpticalPhoton()

◆ PostStepDoIt()

G4VParticleChange * G4OpWLS2::PostStepDoIt ( const G4Track aTrack,
const G4Step aStep 
)
overridevirtual

Reimplemented from G4VDiscreteProcess.

Definition at line 92 of file G4OpWLS2.cc.

94{
95 std::vector<G4Track*> proposedSecondaries;
98
99 if(verboseLevel > 1)
100 {
101 G4cout << "\n** G4OpWLS2: Photon absorbed! **" << G4endl;
102 }
103
104 G4StepPoint* pPostStepPoint = aStep.GetPostStepPoint();
107 if(!MPT)
108 {
109 return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
110 }
111 if(!MPT->GetProperty(kWLSCOMPONENT2))
112 {
113 return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
114 }
115
116 G4int NumPhotons = 1;
118 {
119 G4double MeanNumberOfPhotons =
121 NumPhotons = G4int(G4Poisson(MeanNumberOfPhotons));
122 if(NumPhotons <= 0)
123 {
124 // return unchanged particle and no secondaries
126 return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
127 }
128 }
129
130 // Retrieve the WLS Integral for this material
131 // new G4PhysicsFreeVector allocated to hold CII's
132 G4double primaryEnergy = aTrack.GetDynamicParticle()->GetKineticEnergy();
133 G4double WLSTime = 0.;
134 G4PhysicsFreeVector* WLSIntegral = nullptr;
135
136 WLSTime = MPT->GetConstProperty(kWLSTIMECONSTANT2);
137 WLSIntegral = (G4PhysicsFreeVector*) ((*theIntegralTable)(
138 aTrack.GetMaterial()->GetIndex()));
139
140 // Max WLS Integral
141 G4double CIImax = WLSIntegral->GetMaxValue();
142 G4int NumberOfPhotons = NumPhotons;
143
144 for(G4int i = 0; i < NumPhotons; ++i)
145 {
146 G4double sampledEnergy;
147 // Make sure the energy of the secondary is less than that of the primary
148 for(G4int j = 1; j <= 100; ++j)
149 {
150 // Determine photon energy
151 G4double CIIvalue = G4UniformRand() * CIImax;
152 sampledEnergy = WLSIntegral->GetEnergy(CIIvalue);
153 if(sampledEnergy <= primaryEnergy)
154 break;
155 }
156 // If no such energy can be sampled, return one less secondary, or none
157 if(sampledEnergy > primaryEnergy)
158 {
159 if(verboseLevel > 1)
160 {
161 G4cout << " *** G4OpWLS2: One less WLS2 photon will be returned ***"
162 << G4endl;
163 }
164 NumberOfPhotons--;
165 if(NumberOfPhotons == 0)
166 {
167 if(verboseLevel > 1)
168 {
169 G4cout << " *** G4OpWLS2: No WLS2 photon can be sampled for this "
170 "primary ***"
171 << G4endl;
172 }
173 // return unchanged particle and no secondaries
175 return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
176 }
177 continue;
178 }
179 else if(verboseLevel > 1)
180 {
181 G4cout << "G4OpWLS2: Created photon with energy: " << sampledEnergy
182 << G4endl;
183 }
184
185 // Generate random photon direction
186 G4double cost = 1. - 2. * G4UniformRand();
187 G4double sint = std::sqrt((1. - cost) * (1. + cost));
188 G4double phi = twopi * G4UniformRand();
189 G4double sinp = std::sin(phi);
190 G4double cosp = std::cos(phi);
191 G4ParticleMomentum photonMomentum(sint * cosp, sint * sinp, cost);
192
193 G4ThreeVector photonPolarization(cost * cosp, cost * sinp, -sint);
194 G4ThreeVector perp = photonMomentum.cross(photonPolarization);
195
196 phi = twopi * G4UniformRand();
197 sinp = std::sin(phi);
198 cosp = std::cos(phi);
199 photonPolarization = (cosp * photonPolarization + sinp * perp).unit();
200
201 // Generate a new photon:
202 auto sec_dp =
204 sec_dp->SetPolarization(photonPolarization);
205 sec_dp->SetKineticEnergy(sampledEnergy);
206
207 G4double secTime = pPostStepPoint->GetGlobalTime() +
209 G4ThreeVector secPos = pPostStepPoint->GetPosition();
210 G4Track* secTrack = new G4Track(sec_dp, secTime, secPos);
211
212 secTrack->SetTouchableHandle(aTrack.GetTouchableHandle());
213 secTrack->SetParentID(aTrack.GetTrackID());
214
215 proposedSecondaries.push_back(secTrack);
216 }
217
218 aParticleChange.SetNumberOfSecondaries((G4int)proposedSecondaries.size());
219 for(auto sec : proposedSecondaries)
220 {
222 }
223 if(verboseLevel > 1)
224 {
225 G4cout << "\n Exiting from G4OpWLS2::DoIt -- NumberOfSecondaries = "
227 }
228
229 return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
230}
@ kWLSMEANNUMBERPHOTONS2
G4long G4Poisson(G4double mean)
Definition: G4Poisson.hh:50
@ fStopAndKill
int G4int
Definition: G4Types.hh:85
#define G4UniformRand()
Definition: Randomize.hh:52
Hep3Vector cross(const Hep3Vector &) const
G4double GetKineticEnergy() const
G4bool ConstPropertyExists(const G4String &key) const
G4double GetConstProperty(const G4String &key) const
size_t GetIndex() const
Definition: G4Material.hh:255
void AddSecondary(G4Track *aSecondary)
void Initialize(const G4Track &) override
G4double GetEnergy(const G4double value) const
G4double GetMaxValue() const
G4double GetGlobalTime() const
const G4ThreeVector & GetPosition() const
G4StepPoint * GetPostStepPoint() const
G4int GetTrackID() const
void SetTouchableHandle(const G4TouchableHandle &apValue)
const G4TouchableHandle & GetTouchableHandle() const
void SetParentID(const G4int aValue)
virtual G4VParticleChange * PostStepDoIt(const G4Track &, const G4Step &)
void ProposeTrackStatus(G4TrackStatus status)
G4int GetNumberOfSecondaries() const
void SetNumberOfSecondaries(G4int totSecondaries)
G4ParticleChange aParticleChange
Definition: G4VProcess.hh:331
virtual G4double GenerateTime(const G4double time_constant)=0

◆ PreparePhysicsTable()

void G4OpWLS2::PreparePhysicsTable ( const G4ParticleDefinition )
overridevirtual

Reimplemented from G4VProcess.

Definition at line 78 of file G4OpWLS2.cc.

79{
80 Initialise();
81}

◆ SetVerboseLevel()

void G4OpWLS2::SetVerboseLevel ( G4int  verbose)

Definition at line 343 of file G4OpWLS2.cc.

Referenced by Initialise().

◆ UseTimeProfile()

void G4OpWLS2::UseTimeProfile ( const G4String  name)
virtual

Definition at line 318 of file G4OpWLS2.cc.

319{
321 {
323 WLSTimeGeneratorProfile = nullptr;
324 }
325 if(name == "delta")
326 {
328 }
329 else if(name == "exponential")
330 {
332 new G4WLSTimeGeneratorProfileExponential("exponential");
333 }
334 else
335 {
336 G4Exception("G4OpWLS::UseTimeProfile", "em0202", FatalException,
337 "generator does not exist");
338 }
340}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:59
void SetWLS2TimeProfile(const G4String &)

Referenced by Initialise().

Member Data Documentation

◆ theIntegralTable

G4PhysicsTable* G4OpWLS2::theIntegralTable
protected

◆ WLSTimeGeneratorProfile

G4VWLSTimeGeneratorProfile* G4OpWLS2::WLSTimeGeneratorProfile
protected

Definition at line 90 of file G4OpWLS2.hh.

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


The documentation for this class was generated from the following files: