Geant4  10.02.p03
G4ParallelWorldProcess Class Reference

#include <G4ParallelWorldProcess.hh>

Inheritance diagram for G4ParallelWorldProcess:
Collaboration diagram for G4ParallelWorldProcess:

Public Member Functions

 G4ParallelWorldProcess (const G4String &processName="ParaWorld", G4ProcessType theType=fParallel)
 
virtual ~G4ParallelWorldProcess ()
 
void SetParallelWorld (G4String parallelWorldName)
 
void SetParallelWorld (G4VPhysicalVolume *parallelWorld)
 
void StartTracking (G4Track *)
 
G4double AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)
 
G4VParticleChange * AtRestDoIt (const G4Track &, const G4Step &)
 
G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
 
G4VParticleChange * AlongStepDoIt (const G4Track &, const G4Step &)
 
G4double PostStepGetPhysicalInteractionLength (const G4Track &, G4double, G4ForceCondition *)
 
G4VParticleChange * PostStepDoIt (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
 
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 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
 
virtual void SetMasterProcess (G4VProcess *masterP)
 
const G4VProcessGetMasterProcess () const
 
virtual void BuildWorkerPhysicsTable (const G4ParticleDefinition &part)
 
virtual void PrepareWorkerPhysicsTable (const G4ParticleDefinition &)
 

Static Public Member Functions

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

Private Member Functions

void CopyStep (const G4Step &step)
 
void SwitchMaterial (G4StepPoint *)
 

Private Attributes

G4Step * fGhostStep
 
G4StepPoint * fGhostPreStepPoint
 
G4StepPoint * fGhostPostStepPoint
 
G4VParticleChange aDummyParticleChange
 
G4ParticleChange xParticleChange
 
G4TransportationManagerfTransportationManager
 
G4PathFinderfPathFinder
 
G4String fGhostWorldName
 
G4VPhysicalVolumefGhostWorld
 
G4NavigatorfGhostNavigator
 
G4int fNavigatorID
 
G4TouchableHandle fOldGhostTouchable
 
G4TouchableHandle fNewGhostTouchable
 
G4FieldTrack fFieldTrack
 
G4double fGhostSafety
 
G4bool fOnBoundary
 
G4bool layeredMaterialFlag
 
G4int iParallelWorld
 

Static Private Attributes

static G4ThreadLocal G4Step * fpHyperStep = 0
 
static G4ThreadLocal G4int nParallelWorlds = 0
 
static G4ThreadLocal G4int fNavIDHyp = 0
 

Additional Inherited Members

- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 
- Protected Attributes inherited from G4VProcess
const G4ProcessManageraProcessManager
 
G4VParticleChange * pParticleChange
 
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 = fParallel 
)

Definition at line 62 of file G4ParallelWorldProcess.cc.

63 :G4VProcess(processName,theType),fGhostWorld(nullptr),fGhostNavigator(nullptr),
65  layeredMaterialFlag(false)
66 {
67  SetProcessSubType(491);
68  if(!fpHyperStep) fpHyperStep = new G4Step();
70 
72 
73  fGhostStep = new G4Step();
74  fGhostPreStepPoint = fGhostStep->GetPreStepPoint();
75  fGhostPostStepPoint = fGhostStep->GetPostStepPoint();
76 
80 
81  fGhostWorldName = "** NotDefined **";
83 
84  if (verboseLevel>0)
85  {
86  G4cout << GetProcessName() << " is created " << G4endl;
87  }
88 }
static G4PathFinder * GetInstance()
Definition: G4PathFinder.cc:57
G4int verboseLevel
Definition: G4VProcess.hh:368
G4VProcess(const G4String &aName="NoName", G4ProcessType aType=fNotDefined)
Definition: G4VProcess.cc:52
G4VPhysicalVolume * fGhostWorld
static G4ParallelWorldProcessStore * GetInstance()
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
G4GLOB_DLL std::ostream G4cout
static G4ThreadLocal G4int nParallelWorlds
void SetParallelWorld(G4ParallelWorldProcess *proc, G4String parallelWorldName)
void SetProcessSubType(G4int)
Definition: G4VProcess.hh:432
G4Navigator * GetNavigatorForTracking() const
static G4TransportationManager * GetTransportationManager()
G4VParticleChange aDummyParticleChange
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
void SetPushVerbosity(G4bool mode)
#define G4endl
Definition: G4ios.hh:61
static G4ThreadLocal G4Step * fpHyperStep
G4TransportationManager * fTransportationManager
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ~G4ParallelWorldProcess()

G4ParallelWorldProcess::~G4ParallelWorldProcess ( )
virtual

Definition at line 90 of file G4ParallelWorldProcess.cc.

91 {
92  delete fGhostStep;
94  if(nParallelWorlds==0)
95  {
96  delete fpHyperStep;
97  fpHyperStep = 0;
98  }
99 }
static G4ThreadLocal G4int nParallelWorlds
static G4ThreadLocal G4Step * fpHyperStep
Here is the call graph for this function:

Member Function Documentation

◆ AlongStepDoIt()

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

Implements G4VProcess.

Definition at line 338 of file G4ParallelWorldProcess.cc.

340 {
341  pParticleChange->Initialize(track);
342  return pParticleChange;
343 }
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283

◆ AlongStepGetPhysicalInteractionLength()

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

Implements G4VProcess.

Definition at line 265 of file G4ParallelWorldProcess.cc.

268 {
269  static G4ThreadLocal G4FieldTrack *endTrack_G4MT_TLS_ = 0 ; if (!endTrack_G4MT_TLS_) endTrack_G4MT_TLS_ = new G4FieldTrack ('0') ; G4FieldTrack &endTrack = *endTrack_G4MT_TLS_;
270  //static ELimited eLimited;
271  ELimited eLimited;
272  ELimited eLim = kUndefLimited;
273 
274  *selection = NotCandidateForSelection;
275  G4double returnedStep = DBL_MAX;
276 
277  if (previousStepSize > 0.)
278  { fGhostSafety -= previousStepSize; }
279  if (fGhostSafety < 0.) fGhostSafety = 0.0;
280 
281  if (currentMinimumStep <= fGhostSafety && currentMinimumStep > 0.)
282  {
283  // I have no chance to limit
284  returnedStep = currentMinimumStep;
285  fOnBoundary = false;
286  proposedSafety = fGhostSafety - currentMinimumStep;
287  eLim = kDoNot;
288  }
289  else
290  {
292 
293 #ifdef G4DEBUG_PARALLEL_WORLD_PROCESS
294  if( verboseLevel > 0 ){
295  int localVerb = verboseLevel-1;
296 
297  if( localVerb == 1 ) {
298  G4cout << " Pll Wrl proc::AlongStepGPIL " << this->GetProcessName() << G4endl;
299  }else if( localVerb > 1 ) {
300  G4cout << "----------------------------------------------" << G4endl;
301  G4cout << " ParallelWorldProcess: field Track set to : " << G4endl;
302  G4cout << "----------------------------------------------" << G4endl;
303  G4cout << fFieldTrack << G4endl;
304  G4cout << "----------------------------------------------" << G4endl;
305  }
306  }
307 #endif
308 
309  returnedStep
310  = fPathFinder->ComputeStep(fFieldTrack,currentMinimumStep,fNavigatorID,
311  track.GetCurrentStepNumber(),fGhostSafety,eLimited,
312  endTrack,track.GetVolume());
313  if(eLimited == kDoNot)
314  {
315  fOnBoundary = false;
316  fGhostSafety = fGhostNavigator->ComputeSafety(endTrack.GetPosition());
317  }
318  else
319  {
320  fOnBoundary = true;
321  // fGhostSafetyEnd = 0.0; // At end-point of expected step only
322  }
323  proposedSafety = fGhostSafety;
324  if(eLimited == kUnique || eLimited == kSharedOther) {
325  *selection = CandidateForSelection;
326  }
327  else if (eLimited == kSharedTransport) {
328  returnedStep *= (1.0 + 1.0e-9);
329  }
330  eLim = eLimited;
331  }
332 
334  if(eLim == kUnique || eLim == kSharedOther) fNavIDHyp = fNavigatorID;
335  return returnedStep;
336 }
G4int verboseLevel
Definition: G4VProcess.hh:368
ELimited
#define G4ThreadLocal
Definition: tls.hh:89
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
G4GLOB_DLL std::ostream G4cout
static G4ThreadLocal G4int nParallelWorlds
static G4ThreadLocal G4int fNavIDHyp
G4double ComputeStep(const G4FieldTrack &pFieldTrack, G4double pCurrentProposedStepLength, G4int navigatorId, G4int stepNo, G4double &pNewSafety, ELimited &limitedStep, G4FieldTrack &EndState, G4VPhysicalVolume *currentVolume)
c1 Update()
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
virtual G4double ComputeSafety(const G4ThreeVector &globalpoint, const G4double pProposedMaxLength=DBL_MAX, const G4bool keepState=true)
#define DBL_MAX
Definition: templates.hh:83
Here is the call graph for this function:

◆ AtRestDoIt()

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

Implements G4VProcess.

Definition at line 171 of file G4ParallelWorldProcess.cc.

174 {
175 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
176 // At Rest must be registered ONLY for the particle which has other At Rest
177 // process(es).
178 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
179  fOldGhostTouchable = fGhostPostStepPoint->GetTouchableHandle();
180  G4VSensitiveDetector* aSD = 0;
183  fOnBoundary = false;
184  if(aSD)
185  {
186  CopyStep(step);
187  fGhostPreStepPoint->SetSensitiveDetector(aSD);
188 
190 
191  fGhostPreStepPoint->SetTouchableHandle(fOldGhostTouchable);
192  fGhostPostStepPoint->SetTouchableHandle(fNewGhostTouchable);
194  {
195  fGhostPostStepPoint->SetSensitiveDetector(
197  }
198  else
199  { fGhostPostStepPoint->SetSensitiveDetector(0); }
200 
201  aSD->Hit(fGhostStep);
202  }
203 
204  pParticleChange->Initialize(track);
205  return pParticleChange;
206 }
void CopyStep(const G4Step &step)
G4TouchableHandle fNewGhostTouchable
G4bool Hit(G4Step *aStep)
virtual G4VPhysicalVolume * GetVolume(G4int depth=0) const
Definition: G4VTouchable.cc:44
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
G4VSensitiveDetector * GetSensitiveDetector() const
G4TouchableHandle fOldGhostTouchable
G4LogicalVolume * GetLogicalVolume() const
Here is the call graph for this function:

◆ AtRestGetPhysicalInteractionLength()

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

Implements G4VProcess.

Definition at line 159 of file G4ParallelWorldProcess.cc.

162 {
163 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
164 // At Rest must be registered ONLY for the particle which has other At Rest
165 // process(es).
166 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
167  *condition = Forced;
168  return DBL_MAX;
169 }
G4double condition(const G4ErrorSymMatrix &m)
#define DBL_MAX
Definition: templates.hh:83

◆ CopyStep()

void G4ParallelWorldProcess::CopyStep ( const G4Step &  step)
private

Definition at line 345 of file G4ParallelWorldProcess.cc.

346 {
347  G4StepStatus prevStat = fGhostPostStepPoint->GetStepStatus();
348 
349  fGhostStep->SetTrack(step.GetTrack());
350  fGhostStep->SetStepLength(step.GetStepLength());
351  fGhostStep->SetTotalEnergyDeposit(step.GetTotalEnergyDeposit());
352  fGhostStep->SetNonIonizingEnergyDeposit(step.GetNonIonizingEnergyDeposit());
353  fGhostStep->SetControlFlag(step.GetControlFlag());
354 
355  *fGhostPreStepPoint = *(step.GetPreStepPoint());
356  *fGhostPostStepPoint = *(step.GetPostStepPoint());
357 
358  fGhostPreStepPoint->SetStepStatus(prevStat);
359  if(fOnBoundary)
360  { fGhostPostStepPoint->SetStepStatus(fGeomBoundary); }
361  else if(fGhostPostStepPoint->GetStepStatus()==fGeomBoundary)
362  { fGhostPostStepPoint->SetStepStatus(fPostStepDoItProc); }
363 
364  if(iParallelWorld==1)
365  {
366  G4StepStatus prevStatHyp = fpHyperStep->GetPostStepPoint()->GetStepStatus();
367 
368  fpHyperStep->SetTrack(step.GetTrack());
369  fpHyperStep->SetStepLength(step.GetStepLength());
370  fpHyperStep->SetTotalEnergyDeposit(step.GetTotalEnergyDeposit());
371  fpHyperStep->SetNonIonizingEnergyDeposit(step.GetNonIonizingEnergyDeposit());
372  fpHyperStep->SetControlFlag(step.GetControlFlag());
373 
374  *(fpHyperStep->GetPreStepPoint()) = *(fpHyperStep->GetPostStepPoint());
375  *(fpHyperStep->GetPostStepPoint()) = *(step.GetPostStepPoint());
376 
377  fpHyperStep->GetPreStepPoint()->SetStepStatus(prevStatHyp);
378  }
379 
380  if(fOnBoundary)
381  { fpHyperStep->GetPostStepPoint()->SetStepStatus(fGeomBoundary); }
382 }
static G4ThreadLocal G4Step * fpHyperStep
Here is the caller graph for this function:

◆ GetHyperStep()

const G4Step * G4ParallelWorldProcess::GetHyperStep ( )
static

Definition at line 56 of file G4ParallelWorldProcess.cc.

57 { return fpHyperStep; }
static G4ThreadLocal G4Step * fpHyperStep
Here is the caller graph for this function:

◆ GetHypNavigatorID()

G4int G4ParallelWorldProcess::GetHypNavigatorID ( )
static

Definition at line 58 of file G4ParallelWorldProcess.cc.

59 { return fNavIDHyp; }
static G4ThreadLocal G4int fNavIDHyp
Here is the call graph for this function:
Here is the caller graph for this function:

◆ GetLayeredMaterialFlag()

G4bool G4ParallelWorldProcess::GetLayeredMaterialFlag ( ) const
inline

Definition at line 127 of file G4ParallelWorldProcess.hh.

Here is the call graph for this function:

◆ IsAtRestRequired()

G4bool G4ParallelWorldProcess::IsAtRestRequired ( G4ParticleDefinition partDef)

Definition at line 424 of file G4ParallelWorldProcess.cc.

425 {
426  G4int pdgCode = partDef->GetPDGEncoding();
427  if(pdgCode==0)
428  {
429  G4String partName = partDef->GetParticleName();
430  if(partName=="opticalphoton") return false;
431  if(partName=="geantino") return false;
432  if(partName=="chargedgeantino") return false;
433  }
434  else
435  {
436  if(pdgCode==22) return false; // gamma
437  if(pdgCode==11) return false; // electron
438  if(pdgCode==2212) return false; // proton
439  if(pdgCode==-12) return false; // anti_nu_e
440  if(pdgCode==12) return false; // nu_e
441  if(pdgCode==-14) return false; // anti_nu_mu
442  if(pdgCode==14) return false; // nu_mu
443  if(pdgCode==-16) return false; // anti_nu_tau
444  if(pdgCode==16) return false; // nu_tau
445  }
446  return true;
447 }
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
Here is the call graph for this function:
Here is the caller graph for this function:

◆ PostStepDoIt()

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

Implements G4VProcess.

Definition at line 218 of file G4ParallelWorldProcess.cc.

221 {
222  fOldGhostTouchable = fGhostPostStepPoint->GetTouchableHandle();
223  G4VSensitiveDetector* aSD = 0;
226  CopyStep(step);
227  fGhostPreStepPoint->SetSensitiveDetector(aSD);
228 
229  if(fOnBoundary)
230  {
232  }
233  else
234  {
236  }
237 
238  fGhostPreStepPoint->SetTouchableHandle(fOldGhostTouchable);
239  fGhostPostStepPoint->SetTouchableHandle(fNewGhostTouchable);
240 
242  {
243  fGhostPostStepPoint->SetSensitiveDetector(
245  }
246  else
247  { fGhostPostStepPoint->SetSensitiveDetector(0); }
248 
249  G4VSensitiveDetector* sd = fGhostPreStepPoint->GetSensitiveDetector();
250  if(sd)
251  {
252  sd->Hit(fGhostStep);
253  }
254 
255  pParticleChange->Initialize(track);
257  {
258  G4StepPoint* realWorldPostStepPoint =
259  ((G4Step*)(track.GetStep()))->GetPostStepPoint();
260  SwitchMaterial(realWorldPostStepPoint);
261  }
262  return pParticleChange;
263 }
void CopyStep(const G4Step &step)
G4TouchableHandle CreateTouchableHandle(G4int navId) const
G4TouchableHandle fNewGhostTouchable
G4bool Hit(G4Step *aStep)
virtual G4VPhysicalVolume * GetVolume(G4int depth=0) const
Definition: G4VTouchable.cc:44
void SwitchMaterial(G4StepPoint *)
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
G4VSensitiveDetector * GetSensitiveDetector() const
G4TouchableHandle fOldGhostTouchable
G4LogicalVolume * GetLogicalVolume() const
Here is the call graph for this function:

◆ PostStepGetPhysicalInteractionLength()

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

Implements G4VProcess.

Definition at line 209 of file G4ParallelWorldProcess.cc.

213 {
214  *condition = StronglyForced;
215  return DBL_MAX;
216 }
G4double condition(const G4ErrorSymMatrix &m)
#define DBL_MAX
Definition: templates.hh:83

◆ SetLayeredMaterialFlag()

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

Definition at line 125 of file G4ParallelWorldProcess.hh.

Here is the caller graph for this function:

◆ SetParallelWorld() [1/2]

void G4ParallelWorldProcess::SetParallelWorld ( G4String  parallelWorldName)

Definition at line 102 of file G4ParallelWorldProcess.cc.

103 {
104  fGhostWorldName = parallelWorldName;
108 }
G4VPhysicalVolume * fGhostWorld
G4VPhysicalVolume * GetParallelWorld(const G4String &worldName)
G4Navigator * GetNavigator(const G4String &worldName)
void SetPushVerbosity(G4bool mode)
G4TransportationManager * fTransportationManager
Here is the call graph for this function:
Here is the caller graph for this function:

◆ SetParallelWorld() [2/2]

void G4ParallelWorldProcess::SetParallelWorld ( G4VPhysicalVolume parallelWorld)

Definition at line 111 of file G4ParallelWorldProcess.cc.

112 {
113  fGhostWorldName = parallelWorld->GetName();
114  fGhostWorld = parallelWorld;
117 }
G4VPhysicalVolume * fGhostWorld
const G4String & GetName() const
G4Navigator * GetNavigator(const G4String &worldName)
void SetPushVerbosity(G4bool mode)
G4TransportationManager * fTransportationManager
Here is the call graph for this function:

◆ StartTracking()

void G4ParallelWorldProcess::StartTracking ( G4Track *  trk)
virtual

Reimplemented from G4VProcess.

Definition at line 119 of file G4ParallelWorldProcess.cc.

120 {
121  if(fGhostNavigator)
123  else
124  {
125  G4Exception("G4ParallelWorldProcess::StartTracking",
126  "ProcParaWorld000",FatalException,
127  "G4ParallelWorldProcess is used for tracking without having a parallel world assigned");
128  }
129  fPathFinder->PrepareNewTrack(trk->GetPosition(),trk->GetMomentumDirection());
130 
132  fGhostPreStepPoint->SetTouchableHandle(fOldGhostTouchable);
134  fGhostPostStepPoint->SetTouchableHandle(fNewGhostTouchable);
135 
136  fGhostSafety = -1.;
137  fOnBoundary = false;
138  fGhostPreStepPoint->SetStepStatus(fUndefined);
139  fGhostPostStepPoint->SetStepStatus(fUndefined);
140 
141 // G4VPhysicalVolume* thePhys = fNewGhostTouchable->GetVolume();
142 // if(thePhys)
143 // {
144 // G4Material* ghostMaterial = thePhys->GetLogicalVolume()->GetMaterial();
145 // if(ghostMaterial)
146 // { G4cout << " --- Material : " << ghostMaterial->GetName() << G4endl; }
147 // }
148 
149  *(fpHyperStep->GetPostStepPoint()) = *(trk->GetStep()->GetPostStepPoint());
151  {
152  G4StepPoint* realWorldPostStepPoint = trk->GetStep()->GetPostStepPoint();
153  SwitchMaterial(realWorldPostStepPoint);
154  }
155  *(fpHyperStep->GetPreStepPoint()) = *(fpHyperStep->GetPostStepPoint());
156 }
void PrepareNewTrack(const G4ThreeVector &position, const G4ThreeVector &direction, G4VPhysicalVolume *massStartVol=0)
G4TouchableHandle CreateTouchableHandle(G4int navId) const
G4TouchableHandle fNewGhostTouchable
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
void SwitchMaterial(G4StepPoint *)
G4int ActivateNavigator(G4Navigator *aNavigator)
G4TouchableHandle fOldGhostTouchable
static G4ThreadLocal G4Step * fpHyperStep
G4TransportationManager * fTransportationManager
Here is the call graph for this function:

◆ SwitchMaterial()

void G4ParallelWorldProcess::SwitchMaterial ( G4StepPoint *  realWorldStepPoint)
private

Definition at line 384 of file G4ParallelWorldProcess.cc.

385 {
386  if(realWorldStepPoint->GetStepStatus()==fWorldBoundary) return;
388  if(thePhys)
389  {
390  G4Material* ghostMaterial = thePhys->GetLogicalVolume()->GetMaterial();
391  if(ghostMaterial)
392  {
393  G4Region* ghostRegion = thePhys->GetLogicalVolume()->GetRegion();
394  G4ProductionCuts* prodCuts =
395  realWorldStepPoint->GetMaterialCutsCouple()->GetProductionCuts();
396  if(ghostRegion)
397  {
398  G4ProductionCuts* ghostProdCuts = ghostRegion->GetProductionCuts();
399  if(ghostProdCuts) prodCuts = ghostProdCuts;
400  }
401  const G4MaterialCutsCouple* ghostMCCouple =
403  ->GetMaterialCutsCouple(ghostMaterial,prodCuts);
404  if(ghostMCCouple)
405  {
406  realWorldStepPoint->SetMaterial(ghostMaterial);
407  realWorldStepPoint->SetMaterialCutsCouple(ghostMCCouple);
408  *(fpHyperStep->GetPostStepPoint()) = *(fGhostPostStepPoint);
409  fpHyperStep->GetPostStepPoint()->SetMaterial(ghostMaterial);
410  fpHyperStep->GetPostStepPoint()->SetMaterialCutsCouple(ghostMCCouple);
411  }
412  else
413  {
414  G4cout << "!!! MaterialCutsCouple is not found for "
415  << ghostMaterial->GetName() << "." << G4endl
416  << " Material in real world ("
417  << realWorldStepPoint->GetMaterial()->GetName()
418  << ") is used." << G4endl;
419  }
420  }
421  }
422 }
G4Material * GetMaterial() const
G4TouchableHandle fNewGhostTouchable
G4GLOB_DLL std::ostream G4cout
virtual G4VPhysicalVolume * GetVolume(G4int depth=0) const
Definition: G4VTouchable.cc:44
static G4ProductionCutsTable * GetProductionCutsTable()
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
G4Region * GetRegion() const
#define G4endl
Definition: G4ios.hh:61
G4ProductionCuts * GetProductionCuts() const
G4LogicalVolume * GetLogicalVolume() const
const G4String & GetName() const
Definition: G4Material.hh:178
static G4ThreadLocal G4Step * fpHyperStep
const std::vector< G4double > & GetProductionCuts() const
void SetMaterial(const G4Material *)
Here is the call graph for this function:
Here is the caller graph for this function:

Member Data Documentation

◆ aDummyParticleChange

G4VParticleChange G4ParallelWorldProcess::aDummyParticleChange
private

Definition at line 145 of file G4ParallelWorldProcess.hh.

◆ fFieldTrack

G4FieldTrack G4ParallelWorldProcess::fFieldTrack
private

Definition at line 160 of file G4ParallelWorldProcess.hh.

◆ fGhostNavigator

G4Navigator* G4ParallelWorldProcess::fGhostNavigator
private

Definition at line 156 of file G4ParallelWorldProcess.hh.

◆ fGhostPostStepPoint

G4StepPoint* G4ParallelWorldProcess::fGhostPostStepPoint
private

Definition at line 143 of file G4ParallelWorldProcess.hh.

◆ fGhostPreStepPoint

G4StepPoint* G4ParallelWorldProcess::fGhostPreStepPoint
private

Definition at line 142 of file G4ParallelWorldProcess.hh.

◆ fGhostSafety

G4double G4ParallelWorldProcess::fGhostSafety
private

Definition at line 161 of file G4ParallelWorldProcess.hh.

◆ fGhostStep

G4Step* G4ParallelWorldProcess::fGhostStep
private

Definition at line 141 of file G4ParallelWorldProcess.hh.

◆ fGhostWorld

G4VPhysicalVolume* G4ParallelWorldProcess::fGhostWorld
private

Definition at line 155 of file G4ParallelWorldProcess.hh.

◆ fGhostWorldName

G4String G4ParallelWorldProcess::fGhostWorldName
private

Definition at line 154 of file G4ParallelWorldProcess.hh.

◆ fNavIDHyp

G4ThreadLocal G4int G4ParallelWorldProcess::fNavIDHyp = 0
staticprivate

Definition at line 178 of file G4ParallelWorldProcess.hh.

◆ fNavigatorID

G4int G4ParallelWorldProcess::fNavigatorID
private

Definition at line 157 of file G4ParallelWorldProcess.hh.

◆ fNewGhostTouchable

G4TouchableHandle G4ParallelWorldProcess::fNewGhostTouchable
private

Definition at line 159 of file G4ParallelWorldProcess.hh.

◆ fOldGhostTouchable

G4TouchableHandle G4ParallelWorldProcess::fOldGhostTouchable
private

Definition at line 158 of file G4ParallelWorldProcess.hh.

◆ fOnBoundary

G4bool G4ParallelWorldProcess::fOnBoundary
private

Definition at line 162 of file G4ParallelWorldProcess.hh.

◆ fPathFinder

G4PathFinder* G4ParallelWorldProcess::fPathFinder
private

Definition at line 149 of file G4ParallelWorldProcess.hh.

◆ fpHyperStep

G4ThreadLocal G4Step * G4ParallelWorldProcess::fpHyperStep = 0
staticprivate

Definition at line 176 of file G4ParallelWorldProcess.hh.

◆ fTransportationManager

G4TransportationManager* G4ParallelWorldProcess::fTransportationManager
private

Definition at line 148 of file G4ParallelWorldProcess.hh.

◆ iParallelWorld

G4int G4ParallelWorldProcess::iParallelWorld
private

Definition at line 179 of file G4ParallelWorldProcess.hh.

◆ layeredMaterialFlag

G4bool G4ParallelWorldProcess::layeredMaterialFlag
private

Definition at line 167 of file G4ParallelWorldProcess.hh.

◆ nParallelWorlds

G4ThreadLocal G4int G4ParallelWorldProcess::nParallelWorlds = 0
staticprivate

Definition at line 177 of file G4ParallelWorldProcess.hh.

◆ xParticleChange

G4ParticleChange G4ParallelWorldProcess::xParticleChange
private

Definition at line 146 of file G4ParallelWorldProcess.hh.


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