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

#include <G4ParallelWorldProcess.hh>

+ Inheritance diagram for G4ParallelWorldProcess:

Public Member Functions

 G4ParallelWorldProcess (const G4String &processName="ParaWorld", G4ProcessType theType=fParameterisation)
 
virtual ~G4ParallelWorldProcess ()
 
void SetParallelWorld (G4String parallelWorldName)
 
void SetParallelWorld (G4VPhysicalVolume *parallelWorld)
 
void StartTracking (G4Track *)
 
G4double AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)
 
G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
 
G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
 
G4VParticleChangeAlongStepDoIt (const G4Track &, const G4Step &)
 
G4double PostStepGetPhysicalInteractionLength (const G4Track &, G4double, G4ForceCondition *)
 
G4VParticleChangePostStepDoIt (const G4Track &, const G4Step &)
 
void SetLayeredMaterialFlag (G4bool flg=true)
 
G4bool GetLayeredMaterialFlag () const
 
G4bool IsAtRestRequired (G4ParticleDefinition *)
 
- Public Member Functions inherited from G4VProcess
 G4VProcess (const G4String &aName="NoName", G4ProcessType aType=fNotDefined)
 
 G4VProcess (const G4VProcess &right)
 
virtual ~G4VProcess ()
 
G4int operator== (const G4VProcess &right) const
 
G4int 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 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
 
void SetVerboseLevel (G4int value)
 
G4int GetVerboseLevel () const
 

Static Public Member Functions

static const G4StepGetHyperStep ()
 
static G4int GetHypNavigatorID ()
 
- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 

Additional Inherited Members

- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 
- Protected Attributes inherited from G4VProcess
const G4ProcessManageraProcessManager
 
G4VParticleChangepParticleChange
 
G4ParticleChange aParticleChange
 
G4double theNumberOfInteractionLengthLeft
 
G4double currentInteractionLength
 
G4double theInitialNumberOfInteractionLength
 
G4String theProcessName
 
G4String thePhysicsTableFileName
 
G4ProcessType theProcessType
 
G4int theProcessSubType
 
G4double thePILfactor
 
G4bool enableAtRestDoIt
 
G4bool enableAlongStepDoIt
 
G4bool enablePostStepDoIt
 
G4int verboseLevel
 

Detailed Description

Definition at line 73 of file G4ParallelWorldProcess.hh.

Constructor & Destructor Documentation

◆ G4ParallelWorldProcess()

G4ParallelWorldProcess::G4ParallelWorldProcess ( const G4String processName = "ParaWorld",
G4ProcessType  theType = fParameterisation 
)

Definition at line 60 of file G4ParallelWorldProcess.cc.

62:G4VProcess(processName,theType), fGhostNavigator(0), fNavigatorID(-1),
63 fFieldTrack('0'),layeredMaterialFlag(false)
64{
65 if(!fpHyperStep) fpHyperStep = new G4Step();
66 iParallelWorld = ++nParallelWorlds;
67
68 pParticleChange = &aDummyParticleChange;
69
70 fGhostStep = new G4Step();
71 fGhostPreStepPoint = fGhostStep->GetPreStepPoint();
72 fGhostPostStepPoint = fGhostStep->GetPostStepPoint();
73
75 fPathFinder = G4PathFinder::GetInstance();
76
77 if (verboseLevel>0)
78 {
79 G4cout << GetProcessName() << " is created " << G4endl;
80 }
81}
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cout
static G4PathFinder * GetInstance()
Definition: G4PathFinder.cc:57
Definition: G4Step.hh:78
G4StepPoint * GetPreStepPoint() const
G4StepPoint * GetPostStepPoint() const
static G4TransportationManager * GetTransportationManager()
G4int verboseLevel
Definition: G4VProcess.hh:368
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
const G4String & GetProcessName() const
Definition: G4VProcess.hh:379

◆ ~G4ParallelWorldProcess()

G4ParallelWorldProcess::~G4ParallelWorldProcess ( )
virtual

Definition at line 83 of file G4ParallelWorldProcess.cc.

84{
85 delete fGhostStep;
86 nParallelWorlds--;
87 if(nParallelWorlds==0)
88 {
89 delete fpHyperStep;
90 fpHyperStep = 0;
91 }
92}

Member Function Documentation

◆ AlongStepDoIt()

G4VParticleChange * G4ParallelWorldProcess::AlongStepDoIt ( const G4Track track,
const G4Step  
)
virtual

Implements G4VProcess.

Definition at line 313 of file G4ParallelWorldProcess.cc.

315{
317 return pParticleChange;
318}
virtual void Initialize(const G4Track &)

◆ AlongStepGetPhysicalInteractionLength()

G4double G4ParallelWorldProcess::AlongStepGetPhysicalInteractionLength ( const G4Track track,
G4double  previousStepSize,
G4double  currentMinimumStep,
G4double proposedSafety,
G4GPILSelection selection 
)
virtual

Implements G4VProcess.

Definition at line 258 of file G4ParallelWorldProcess.cc.

261{
262 static G4FieldTrack endTrack('0');
263 //static ELimited eLimited;
264 ELimited eLimited;
265 ELimited eLim = kUndefLimited;
266
267 *selection = NotCandidateForSelection;
268 G4double returnedStep = DBL_MAX;
269
270 if (previousStepSize > 0.)
271 { fGhostSafety -= previousStepSize; }
272 if (fGhostSafety < 0.) fGhostSafety = 0.0;
273
274 if (currentMinimumStep <= fGhostSafety && currentMinimumStep > 0.)
275 {
276 // I have no chance to limit
277 returnedStep = currentMinimumStep;
278 fOnBoundary = false;
279 proposedSafety = fGhostSafety - currentMinimumStep;
280 eLim = kDoNot;
281 }
282 else
283 {
284 G4FieldTrackUpdator::Update(&fFieldTrack,&track);
285 returnedStep
286 = fPathFinder->ComputeStep(fFieldTrack,currentMinimumStep,fNavigatorID,
287 track.GetCurrentStepNumber(),fGhostSafety,eLimited,
288 endTrack,track.GetVolume());
289 if(eLimited == kDoNot)
290 {
291 fOnBoundary = false;
292 fGhostSafety = fGhostNavigator->ComputeSafety(endTrack.GetPosition());
293 }
294 else
295 {
296 fOnBoundary = true;
297 }
298 proposedSafety = fGhostSafety;
299 if(eLimited == kUnique || eLimited == kSharedOther) {
300 *selection = CandidateForSelection;
301 }
302 else if (eLimited == kSharedTransport) {
303 returnedStep *= (1.0 + 1.0e-9);
304 }
305 eLim = eLimited;
306 }
307
308 if(iParallelWorld==nParallelWorlds) fNavIDHyp = 0;
309 if(eLim == kUnique || eLim == kSharedOther) fNavIDHyp = fNavigatorID;
310 return returnedStep;
311}
@ CandidateForSelection
@ NotCandidateForSelection
ELimited
@ kDoNot
@ kUndefLimited
@ kUnique
@ kSharedOther
@ kSharedTransport
double G4double
Definition: G4Types.hh:64
static void Update(G4FieldTrack *, const G4Track *)
virtual G4double ComputeSafety(const G4ThreeVector &globalpoint, const G4double pProposedMaxLength=DBL_MAX, const G4bool keepState=true)
G4double ComputeStep(const G4FieldTrack &pFieldTrack, G4double pCurrentProposedStepLength, G4int navigatorId, G4int stepNo, G4double &pNewSafety, ELimited &limitedStep, G4FieldTrack &EndState, G4VPhysicalVolume *currentVolume)
G4VPhysicalVolume * GetVolume() const
G4int GetCurrentStepNumber() const
#define DBL_MAX
Definition: templates.hh:83

◆ AtRestDoIt()

G4VParticleChange * G4ParallelWorldProcess::AtRestDoIt ( const G4Track track,
const G4Step step 
)
virtual

Implements G4VProcess.

Definition at line 164 of file G4ParallelWorldProcess.cc.

167{
168//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
169// At Rest must be registered ONLY for the particle which has other At Rest
170// process(es).
171//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
172 fOldGhostTouchable = fGhostPostStepPoint->GetTouchableHandle();
173 G4VSensitiveDetector* aSD = 0;
174 if(fOldGhostTouchable->GetVolume())
175 { aSD = fOldGhostTouchable->GetVolume()->GetLogicalVolume()->GetSensitiveDetector(); }
176 fOnBoundary = false;
177 if(aSD)
178 {
179 CopyStep(step);
180 fGhostPreStepPoint->SetSensitiveDetector(aSD);
181
182 fNewGhostTouchable = fOldGhostTouchable;
183
184 fGhostPreStepPoint->SetTouchableHandle(fOldGhostTouchable);
185 fGhostPostStepPoint->SetTouchableHandle(fNewGhostTouchable);
186 if(fNewGhostTouchable->GetVolume())
187 {
188 fGhostPostStepPoint->SetSensitiveDetector(
189 fNewGhostTouchable->GetVolume()->GetLogicalVolume()->GetSensitiveDetector());
190 }
191 else
192 { fGhostPostStepPoint->SetSensitiveDetector(0); }
193
194 aSD->Hit(fGhostStep);
195 }
196
198 return pParticleChange;
199}
G4VSensitiveDetector * GetSensitiveDetector() const
void SetSensitiveDetector(G4VSensitiveDetector *)
void SetTouchableHandle(const G4TouchableHandle &apValue)
const G4TouchableHandle & GetTouchableHandle() const
G4LogicalVolume * GetLogicalVolume() const
G4bool Hit(G4Step *aStep)
virtual G4VPhysicalVolume * GetVolume(G4int depth=0) const
Definition: G4VTouchable.cc:44

◆ AtRestGetPhysicalInteractionLength()

G4double G4ParallelWorldProcess::AtRestGetPhysicalInteractionLength ( const G4Track ,
G4ForceCondition condition 
)
virtual

Implements G4VProcess.

Definition at line 152 of file G4ParallelWorldProcess.cc.

155{
156//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
157// At Rest must be registered ONLY for the particle which has other At Rest
158// process(es).
159//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
160 *condition = Forced;
161 return DBL_MAX;
162}
G4double condition(const G4ErrorSymMatrix &m)
@ Forced

◆ GetHyperStep()

const G4Step * G4ParallelWorldProcess::GetHyperStep ( )
static

Definition at line 55 of file G4ParallelWorldProcess.cc.

56{ return fpHyperStep; }

◆ GetHypNavigatorID()

G4int G4ParallelWorldProcess::GetHypNavigatorID ( )
static

Definition at line 57 of file G4ParallelWorldProcess.cc.

58{ return fNavIDHyp; }

◆ GetLayeredMaterialFlag()

G4bool G4ParallelWorldProcess::GetLayeredMaterialFlag ( ) const
inline

Definition at line 127 of file G4ParallelWorldProcess.hh.

128 { return layeredMaterialFlag; }

◆ IsAtRestRequired()

G4bool G4ParallelWorldProcess::IsAtRestRequired ( G4ParticleDefinition partDef)

Definition at line 399 of file G4ParallelWorldProcess.cc.

400{
401 G4int pdgCode = partDef->GetPDGEncoding();
402 if(pdgCode==0)
403 {
404 G4String partName = partDef->GetParticleName();
405 if(partName=="opticalphoton") return false;
406 if(partName=="geantino") return false;
407 if(partName=="chargedgeantino") return false;
408 }
409 else
410 {
411 if(pdgCode==22) return false; // gamma
412 if(pdgCode==11) return false; // electron
413 if(pdgCode==2212) return false; // proton
414 if(pdgCode==-12) return false; // anti_nu_e
415 if(pdgCode==12) return false; // nu_e
416 if(pdgCode==-14) return false; // anti_nu_mu
417 if(pdgCode==14) return false; // nu_mu
418 if(pdgCode==-16) return false; // anti_nu_tau
419 if(pdgCode==16) return false; // nu_tau
420 }
421 return true;
422}
int G4int
Definition: G4Types.hh:66
const G4String & GetParticleName() const

◆ PostStepDoIt()

G4VParticleChange * G4ParallelWorldProcess::PostStepDoIt ( const G4Track track,
const G4Step step 
)
virtual

Implements G4VProcess.

Definition at line 211 of file G4ParallelWorldProcess.cc.

214{
215 fOldGhostTouchable = fGhostPostStepPoint->GetTouchableHandle();
216 G4VSensitiveDetector* aSD = 0;
217 if(fOldGhostTouchable->GetVolume())
218 { aSD = fOldGhostTouchable->GetVolume()->GetLogicalVolume()->GetSensitiveDetector(); }
219 CopyStep(step);
220 fGhostPreStepPoint->SetSensitiveDetector(aSD);
221
222 if(fOnBoundary)
223 {
224 fNewGhostTouchable = fPathFinder->CreateTouchableHandle(fNavigatorID);
225 }
226 else
227 {
228 fNewGhostTouchable = fOldGhostTouchable;
229 }
230
231 fGhostPreStepPoint->SetTouchableHandle(fOldGhostTouchable);
232 fGhostPostStepPoint->SetTouchableHandle(fNewGhostTouchable);
233
234 if(fNewGhostTouchable->GetVolume())
235 {
236 fGhostPostStepPoint->SetSensitiveDetector(
237 fNewGhostTouchable->GetVolume()->GetLogicalVolume()->GetSensitiveDetector());
238 }
239 else
240 { fGhostPostStepPoint->SetSensitiveDetector(0); }
241
242 G4VSensitiveDetector* sd = fGhostPreStepPoint->GetSensitiveDetector();
243 if(sd)
244 {
245 sd->Hit(fGhostStep);
246 }
247
249 if(layeredMaterialFlag)
250 {
251 G4StepPoint* realWorldPostStepPoint =
252 ((G4Step*)(track.GetStep()))->GetPostStepPoint();
253 SwitchMaterial(realWorldPostStepPoint);
254 }
255 return pParticleChange;
256}
G4TouchableHandle CreateTouchableHandle(G4int navId) const
G4VSensitiveDetector * GetSensitiveDetector() const
const G4Step * GetStep() const

◆ PostStepGetPhysicalInteractionLength()

G4double G4ParallelWorldProcess::PostStepGetPhysicalInteractionLength ( const G4Track ,
G4double  ,
G4ForceCondition condition 
)
virtual

Implements G4VProcess.

Definition at line 202 of file G4ParallelWorldProcess.cc.

206{
208 return DBL_MAX;
209}
@ StronglyForced

◆ SetLayeredMaterialFlag()

void G4ParallelWorldProcess::SetLayeredMaterialFlag ( G4bool  flg = true)
inline

Definition at line 125 of file G4ParallelWorldProcess.hh.

126 { layeredMaterialFlag = flg; }

◆ SetParallelWorld() [1/2]

void G4ParallelWorldProcess::SetParallelWorld ( G4String  parallelWorldName)

Definition at line 94 of file G4ParallelWorldProcess.cc.

96{
97 fGhostWorldName = parallelWorldName;
98 fGhostWorld = fTransportationManager->GetParallelWorld(fGhostWorldName);
99 fGhostNavigator = fTransportationManager->GetNavigator(fGhostWorld);
100}
G4VPhysicalVolume * GetParallelWorld(const G4String &worldName)
G4Navigator * GetNavigator(const G4String &worldName)

◆ SetParallelWorld() [2/2]

void G4ParallelWorldProcess::SetParallelWorld ( G4VPhysicalVolume parallelWorld)

Definition at line 102 of file G4ParallelWorldProcess.cc.

104{
105 fGhostWorldName = parallelWorld->GetName();
106 fGhostWorld = parallelWorld;
107 fGhostNavigator = fTransportationManager->GetNavigator(fGhostWorld);
108}
const G4String & GetName() const

◆ StartTracking()

void G4ParallelWorldProcess::StartTracking ( G4Track trk)
virtual

Reimplemented from G4VProcess.

Definition at line 110 of file G4ParallelWorldProcess.cc.

111{
112 if(fGhostNavigator)
113 { fNavigatorID = fTransportationManager->ActivateNavigator(fGhostNavigator); }
114 else
115 {
116 G4Exception("G4ParallelWorldProcess::StartTracking",
117 "ProcParaWorld000",FatalException,
118 "G4ParallelWorldProcess is used for tracking without having a parallel world assigned");
119 }
120 fPathFinder->PrepareNewTrack(trk->GetPosition(),trk->GetMomentumDirection());
121
122 fOldGhostTouchable = fPathFinder->CreateTouchableHandle(fNavigatorID);
123 fGhostPreStepPoint->SetTouchableHandle(fOldGhostTouchable);
124 fNewGhostTouchable = fOldGhostTouchable;
125 fGhostPostStepPoint->SetTouchableHandle(fNewGhostTouchable);
126
127 fGhostSafety = -1.;
128 fOnBoundary = false;
129 fGhostPreStepPoint->SetStepStatus(fUndefined);
130 fGhostPostStepPoint->SetStepStatus(fUndefined);
131
132// G4VPhysicalVolume* thePhys = fNewGhostTouchable->GetVolume();
133// G4cout << "======= G4ParallelWorldProcess::StartTracking() =======" << G4endl;
134// if(thePhys)
135// {
136// G4cout << " --- Parallel world volume : " << thePhys->GetName() << G4endl;
137// G4Material* ghostMaterial = thePhys->GetLogicalVolume()->GetMaterial();
138// if(ghostMaterial)
139// { G4cout << " --- Material : " << ghostMaterial->GetName() << G4endl; }
140// }
141
142 *(fpHyperStep->GetPostStepPoint()) = *(trk->GetStep()->GetPostStepPoint());
143 if(layeredMaterialFlag)
144 {
145 G4StepPoint* realWorldPostStepPoint = trk->GetStep()->GetPostStepPoint();
146 SwitchMaterial(realWorldPostStepPoint);
147 }
148 *(fpHyperStep->GetPreStepPoint()) = *(fpHyperStep->GetPostStepPoint());
149}
@ FatalException
@ fUndefined
Definition: G4StepStatus.hh:66
void PrepareNewTrack(const G4ThreeVector &position, const G4ThreeVector &direction, G4VPhysicalVolume *massStartVol=0)
void SetStepStatus(const G4StepStatus aValue)
const G4ThreeVector & GetPosition() const
const G4ThreeVector & GetMomentumDirection() const
G4int ActivateNavigator(G4Navigator *aNavigator)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41

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