Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4ImportanceProcess Class Reference

#include <G4ImportanceProcess.hh>

Inheritance diagram for G4ImportanceProcess:
Collaboration diagram for G4ImportanceProcess:

Public Member Functions

 G4ImportanceProcess (const G4VImportanceAlgorithm &aImportanceAlgorithm, const G4VIStore &aIstore, const G4VTrackTerminator *TrackTerminator, const G4String &aName="ImportanceProcess", G4bool para=false)
 
virtual ~G4ImportanceProcess ()
 
void SetParallelWorld (const G4String &parallelWorldName)
 
void StartTracking (G4Track *)
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4VParticleChangePostStepDoIt (const G4Track &, const G4Step &)
 
virtual void KillTrack () const
 
virtual const G4StringGetName () const
 
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 &)
 
- 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 &)
 
- Public Member Functions inherited from G4VTrackTerminator
 G4VTrackTerminator ()
 
virtual ~G4VTrackTerminator ()
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 
- 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
 
- Protected Attributes inherited from G4VTrackTerminator
G4double kCarTolerance
 

Detailed Description

Definition at line 61 of file G4ImportanceProcess.hh.

Constructor & Destructor Documentation

G4ImportanceProcess::G4ImportanceProcess ( const G4VImportanceAlgorithm aImportanceAlgorithm,
const G4VIStore aIstore,
const G4VTrackTerminator TrackTerminator,
const G4String aName = "ImportanceProcess",
G4bool  para = false 
)

Definition at line 55 of file G4ImportanceProcess.cc.

59  : G4VProcess(aName),
60  fParticleChange(new G4ParticleChange),
61  fImportanceAlgorithm(aImportanceAlgorithm),
62  fIStore(aIstore),
63  fPostStepAction(0),
64  fGhostWorldName("NoParallelWorld"),fGhostWorld(0),
65  fGhostNavigator(0), fNavigatorID(-1), fFieldTrack('0'),
66  fParaflag(para), fEndTrack('0'), feLimited(kDoNot)
67 {
68  G4cout << G4endl << G4endl << G4endl;
69  G4cout << "G4ImportanceProcess:: Creating " << G4endl;
70  if (TrackTerminator)
71  {
72  fPostStepAction = new G4SamplingPostStepAction(*TrackTerminator);
73  }
74  else
75  {
76  fPostStepAction = new G4SamplingPostStepAction(*this);
77  }
78  if (!fParticleChange)
79  {
80  G4Exception("G4ImportanceProcess::G4ImportanceProcess()",
81  "FatalError", FatalException,
82  "Failed allocation of G4ParticleChange !");
83  }
84  G4VProcess::pParticleChange = fParticleChange;
85 
86 
87  fGhostStep = new G4Step();
88  fGhostPreStepPoint = fGhostStep->GetPreStepPoint();
89  fGhostPostStepPoint = fGhostStep->GetPostStepPoint();
90 
91  fTransportationManager = G4TransportationManager::GetTransportationManager();
92  fPathFinder = G4PathFinder::GetInstance();
93 
94  if (verboseLevel>0)
95  {
96  G4cout << GetProcessName() << " is created " << G4endl;
97  }
98 
99  G4cout << "G4ImportanceProcess:: importance process paraflag is: " << fParaflag << G4endl;
100 
101 }
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
G4StepPoint * GetPreStepPoint() const
G4GLOB_DLL std::ostream G4cout
Definition: G4Step.hh:76
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
static G4TransportationManager * GetTransportationManager()
G4StepPoint * GetPostStepPoint() const
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

G4ImportanceProcess::~G4ImportanceProcess ( )
virtual

Definition at line 103 of file G4ImportanceProcess.cc.

104 {
105 
106  delete fPostStepAction;
107  delete fParticleChange;
108  // delete fGhostStep;
109  // delete fGhostWorld;
110  // delete fGhostNavigator;
111 
112 }

Member Function Documentation

G4VParticleChange * G4ImportanceProcess::AlongStepDoIt ( const G4Track aTrack,
const G4Step  
)
virtual

Implements G4VProcess.

Definition at line 446 of file G4ImportanceProcess.cc.

447 {
448  // Dummy ParticleChange ie: does nothing
449  // Expecting G4Transportation to move the track
450  //AH G4cout << " along step do it " << G4endl;
451  pParticleChange->Initialize(aTrack);
452  return pParticleChange;
453  // return 0;
454 }
virtual void Initialize(const G4Track &)
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283

Here is the call graph for this function:

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

Implements G4VProcess.

Definition at line 357 of file G4ImportanceProcess.cc.

360 {
361  if(fParaflag) {
362  *selection = NotCandidateForSelection;
363  G4double returnedStep = DBL_MAX;
364 
365  if (previousStepSize > 0.)
366  { fGhostSafety -= previousStepSize; }
367  // else
368  // { fGhostSafety = -1.; }
369  if (fGhostSafety < 0.) fGhostSafety = 0.0;
370 
371  // ------------------------------------------
372  // Determination of the proposed STEP LENGTH:
373  // ------------------------------------------
374  if (currentMinimumStep <= fGhostSafety && currentMinimumStep > 0.)
375  {
376  // I have no chance to limit
377  returnedStep = currentMinimumStep;
378  fOnBoundary = false;
379  proposedSafety = fGhostSafety - currentMinimumStep;
380  //AH G4cout << " step not limited, why? " << G4endl;
381  }
382  else // (currentMinimumStep > fGhostSafety: I may limit the Step)
383  {
384  G4FieldTrackUpdator::Update(&fFieldTrack,&track);
385  //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
386  // ComputeStep
387  //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
388  returnedStep
389  = fPathFinder->ComputeStep(fFieldTrack,currentMinimumStep,fNavigatorID,
390  track.GetCurrentStepNumber(),fGhostSafety,feLimited,
391  fEndTrack,track.GetVolume());
392  //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
393  if(feLimited == kDoNot)
394  {
395  //AH G4cout << " computestep came back with not-boundary " << G4endl;
396  // Track is not on the boundary
397  fOnBoundary = false;
398  fGhostSafety = fGhostNavigator->ComputeSafety(fEndTrack.GetPosition());
399  }
400  else
401  {
402  // Track is on the boundary
403  //AH G4cout << " FOUND A BOUNDARY ! " << track.GetCurrentStepNumber() << G4endl;
404  fOnBoundary = true;
405  // proposedSafety = fGhostSafety;
406  }
407  proposedSafety = fGhostSafety;
408  if(feLimited == kUnique || feLimited == kSharedOther) {
409  *selection = CandidateForSelection;
410  }else if (feLimited == kSharedTransport) {
411  returnedStep *= (1.0 + 1.0e-9);
412  // Expand to disable its selection in Step Manager comparison
413  }
414 
415  }
416 
417  // ----------------------------------------------
418  // Returns the fGhostSafety as the proposedSafety
419  // The SteppingManager will take care of keeping
420  // the smallest one.
421  // ----------------------------------------------
422  return returnedStep;
423 
424  } else {
425 
426  return DBL_MAX;
427 
428  }
429 
430 }
G4int GetCurrentStepNumber() const
G4double ComputeStep(const G4FieldTrack &pFieldTrack, G4double pCurrentProposedStepLength, G4int navigatorId, G4int stepNo, G4double &pNewSafety, ELimited &limitedStep, G4FieldTrack &EndState, G4VPhysicalVolume *currentVolume)
G4VPhysicalVolume * GetVolume() const
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
static void Update(G4FieldTrack *, const G4Track *)

Here is the call graph for this function:

G4VParticleChange * G4ImportanceProcess::AtRestDoIt ( const G4Track ,
const G4Step  
)
virtual

Implements G4VProcess.

Definition at line 440 of file G4ImportanceProcess.cc.

441 {
442  return 0;
443 }
G4double G4ImportanceProcess::AtRestGetPhysicalInteractionLength ( const G4Track ,
G4ForceCondition  
)
virtual

Implements G4VProcess.

Definition at line 433 of file G4ImportanceProcess.cc.

435 {
436  return -1.0;
437 }
const G4String & G4ImportanceProcess::GetName ( void  ) const
virtual

Implements G4VTrackTerminator.

Definition at line 351 of file G4ImportanceProcess.cc.

352 {
353  return theProcessName;
354 }
G4String theProcessName
Definition: G4VProcess.hh:335
void G4ImportanceProcess::KillTrack ( ) const
virtual

Implements G4VTrackTerminator.

Definition at line 346 of file G4ImportanceProcess.cc.

347 {
348  fParticleChange->ProposeTrackStatus(fStopAndKill);
349 }
void ProposeTrackStatus(G4TrackStatus status)

Here is the call graph for this function:

G4VParticleChange * G4ImportanceProcess::PostStepDoIt ( const G4Track aTrack,
const G4Step aStep 
)
virtual

Implements G4VProcess.

Definition at line 211 of file G4ImportanceProcess.cc.

213 {
214  fParticleChange->Initialize(aTrack);
215 
216  if(fParaflag) {
217  fOldGhostTouchable = fGhostPostStepPoint->GetTouchableHandle();
218  //xbug? fOnBoundary = false;
219  CopyStep(aStep);
220 
221  if(fOnBoundary)
222  {
223 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
224 // Locate the point and get new touchable
225 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
226  //?? fPathFinder->Locate(step.GetPostStepPoint()->GetPosition(),
227  //?? step.GetPostStepPoint()->GetMomentumDirection());
228  fNewGhostTouchable = fPathFinder->CreateTouchableHandle(fNavigatorID);
229  //AH G4cout << " on boundary " << G4endl;
230 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
231  }
232  else
233  {
234  // Do I need this ??????????????????????????????????????????????????????????
235  // fGhostNavigator->LocateGlobalPointWithinVolume(track.GetPosition());
236  // ?????????????????????????????????????????????????????????????????????????
237 
238  // fPathFinder->ReLocate(track.GetPosition());
239 
240  // reuse the touchable
241  fNewGhostTouchable = fOldGhostTouchable;
242  //AH G4cout << " NOT on boundary " << G4endl;
243  }
244 
245  fGhostPreStepPoint->SetTouchableHandle(fOldGhostTouchable);
246  fGhostPostStepPoint->SetTouchableHandle(fNewGhostTouchable);
247 
248 // if ( (aStep.GetPostStepPoint()->GetStepStatus() == fGeomBoundary)
249 // && (aStep.GetStepLength() > kCarTolerance) )
250 // {
251 // if (aTrack.GetTrackStatus()==fStopAndKill)
252 // {
253 // G4cout << "WARNING - G4ImportanceProcess::PostStepDoIt()"
254 // << " StopAndKill track." << G4endl;
255 // }
256 
257 // G4StepPoint *prepoint = aStep.GetPreStepPoint();
258 // G4StepPoint *postpoint = aStep.GetPostStepPoint();
259 // G4GeometryCell prekey(*(prepoint->GetPhysicalVolume()),
260 // prepoint->GetTouchable()->GetReplicaNumber());
261 // G4GeometryCell postkey(*(postpoint->GetPhysicalVolume()),
262 // postpoint->GetTouchable()->GetReplicaNumber());
263 
264 
265  if ( (fGhostPostStepPoint->GetStepStatus() == fGeomBoundary)
266  && (aStep.GetStepLength() > kCarTolerance) )
267  {
268  if (aTrack.GetTrackStatus()==fStopAndKill)
269  {
270  G4cout << "WARNING - G4ImportanceProcess::PostStepDoIt()"
271  << " StopAndKill track. on boundary" << G4endl;
272  }
273 
274  G4GeometryCell prekey(*(fGhostPreStepPoint->GetPhysicalVolume()),
275  fGhostPreStepPoint->GetTouchable()->GetReplicaNumber());
276  G4GeometryCell postkey(*(fGhostPostStepPoint->GetPhysicalVolume()),
277  fGhostPostStepPoint->GetTouchable()->GetReplicaNumber());
278 
279  //AH
280  /*
281  G4cout << G4endl;
282  G4cout << G4endl;
283  G4cout << " inside parallel importance process " << aTrack.GetCurrentStepNumber() << G4endl;
284  G4cout << G4endl;
285  G4cout << G4endl;
286  G4cout << " prekey: " << fGhostPreStepPoint->GetPhysicalVolume()->GetName() << " replica:"
287  << fGhostPreStepPoint->GetTouchable()->GetReplicaNumber() << G4endl;
288  G4cout << " prekey ISTORE: " << fIStore.GetImportance(prekey) << G4endl;
289  G4cout << " postkey: " << G4endl;
290  G4cout << " postkey ISTORE: " << fIStore.GetImportance(postkey) << G4endl;
291  */
292  //AH
293  G4Nsplit_Weight nw = fImportanceAlgorithm.
294  Calculate(fIStore.GetImportance(prekey),
295  fIStore.GetImportance(postkey),
296  aTrack.GetWeight());
297  //AH
298  /*
299  G4cout << " prekey weight: " << fIStore.GetImportance(prekey)
300  << " postkey weight: " << fIStore.GetImportance(postkey)
301  << " split weight: " << nw << G4endl;
302  G4cout << " before poststepaction " << G4endl;
303  */
304  //AH
305  fPostStepAction->DoIt(aTrack, fParticleChange, nw);
306  //AH G4cout << " after post step do it " << G4endl;
307  }
308  } else {
309  if ( (aStep.GetPostStepPoint()->GetStepStatus() == fGeomBoundary)
310  && (aStep.GetStepLength() > kCarTolerance) )
311  {
312  //AH G4cout << " inside non-parallel importance process " << G4endl;
313  if (aTrack.GetTrackStatus()==fStopAndKill)
314  {
315  G4cout << "WARNING - G4ImportanceProcess::PostStepDoIt()"
316  << " StopAndKill track. on boundary non-parallel" << G4endl;
317  }
318 
319  G4StepPoint *prepoint = aStep.GetPreStepPoint();
320  G4StepPoint *postpoint = aStep.GetPostStepPoint();
321 
322  G4GeometryCell prekey(*(prepoint->GetPhysicalVolume()),
323  prepoint->GetTouchable()->GetReplicaNumber());
324  G4GeometryCell postkey(*(postpoint->GetPhysicalVolume()),
325  postpoint->GetTouchable()->GetReplicaNumber());
326 
327  G4Nsplit_Weight nw = fImportanceAlgorithm.
328  Calculate(fIStore.GetImportance(prekey),
329  fIStore.GetImportance(postkey),
330  aTrack.GetWeight());
331  //AH
332  /*
333  G4cout << " prekey weight: " << fIStore.GetImportance(prekey)
334  << " postkey weight: " << fIStore.GetImportance(postkey)
335  << " split weight: " << nw << G4endl;
336  G4cout << " before poststepaction 2 " << G4endl;
337  */
338  //AH
339  fPostStepAction->DoIt(aTrack, fParticleChange, nw);
340  //AH G4cout << " after poststepaction 2 " << G4endl;
341  }
342  }
343  return fParticleChange;
344 }
G4double GetStepLength() const
G4StepStatus GetStepStatus() const
G4TrackStatus GetTrackStatus() const
G4TouchableHandle CreateTouchableHandle(G4int navId) const
const G4VTouchable * GetTouchable() const
G4StepPoint * GetPreStepPoint() const
G4GLOB_DLL std::ostream G4cout
G4VPhysicalVolume * GetPhysicalVolume() const
virtual G4double GetImportance(const G4GeometryCell &gCell) const =0
virtual void Initialize(const G4Track &)
void DoIt(const G4Track &aTrack, G4ParticleChange *aParticleChange, const G4Nsplit_Weight &nw)
G4StepPoint * GetPostStepPoint() const
virtual G4int GetReplicaNumber(G4int depth=0) const
Definition: G4VTouchable.cc:58
G4double GetWeight() const
#define G4endl
Definition: G4ios.hh:61
void SetTouchableHandle(const G4TouchableHandle &apValue)
const G4TouchableHandle & GetTouchableHandle() const

Here is the call graph for this function:

G4double G4ImportanceProcess::PostStepGetPhysicalInteractionLength ( const G4Track aTrack,
G4double  previousStepSize,
G4ForceCondition condition 
)
virtual

Implements G4VProcess.

Definition at line 198 of file G4ImportanceProcess.cc.

201 {
202 // *condition = Forced;
203 // return kInfinity;
204 
205 // *condition = StronglyForced;
206  *condition = Forced;
207  return DBL_MAX;
208 }
G4double condition(const G4ErrorSymMatrix &m)
#define DBL_MAX
Definition: templates.hh:83
void G4ImportanceProcess::SetParallelWorld ( const G4String parallelWorldName)

Definition at line 122 of file G4ImportanceProcess.cc.

123 {
124  G4cout << G4endl << G4endl << G4endl;
125  G4cout << "G4ImportanceProcess:: SetParallelWorld name = " << parallelWorldName << G4endl;
126 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
127 // Get pointers of the parallel world and its navigator
128 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
129  fGhostWorldName = parallelWorldName;
130  fGhostWorld = fTransportationManager->GetParallelWorld(fGhostWorldName);
131  fGhostNavigator = fTransportationManager->GetNavigator(fGhostWorld);
132 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
133 }
G4GLOB_DLL std::ostream G4cout
G4VPhysicalVolume * GetParallelWorld(const G4String &worldName)
G4Navigator * GetNavigator(const G4String &worldName)
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

Here is the caller graph for this function:

void G4ImportanceProcess::StartTracking ( G4Track trk)
virtual

Reimplemented from G4VProcess.

Definition at line 156 of file G4ImportanceProcess.cc.

157 {
158 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
159 // Activate navigator and get the navigator ID
160 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
161 // G4cout << " G4ParallelWorldScoringProcess::StartTracking" << G4endl;
162 
163  if(fParaflag) {
164  if(fGhostNavigator)
165  { fNavigatorID = fTransportationManager->ActivateNavigator(fGhostNavigator); }
166  else
167  {
168  G4Exception("G4ImportanceProcess::StartTracking",
169  "ProcParaWorld000",FatalException,
170  "G4ImportanceProcess is used for tracking without having a parallel world assigned");
171  }
172 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
173 
174 // G4cout << "G4ParallelWorldScoringProcess::StartTracking <<<<<<<<<<<<<<<<<< " << G4endl;
175 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
176 // Let PathFinder initialize
177 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
178  fPathFinder->PrepareNewTrack(trk->GetPosition(),trk->GetMomentumDirection());
179 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
180 
181 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
182 // Setup initial touchables for the first step
183 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
184  fOldGhostTouchable = fPathFinder->CreateTouchableHandle(fNavigatorID);
185  fGhostPreStepPoint->SetTouchableHandle(fOldGhostTouchable);
186  fNewGhostTouchable = fOldGhostTouchable;
187  fGhostPostStepPoint->SetTouchableHandle(fNewGhostTouchable);
188 
189  // Initialize
190  fGhostSafety = -1.;
191  fOnBoundary = false;
192  }
193 
194 }
void PrepareNewTrack(const G4ThreeVector &position, const G4ThreeVector &direction, G4VPhysicalVolume *massStartVol=0)
const G4ThreeVector & GetPosition() const
G4TouchableHandle CreateTouchableHandle(G4int navId) const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4int ActivateNavigator(G4Navigator *aNavigator)
const G4ThreeVector & GetMomentumDirection() const
void SetTouchableHandle(const G4TouchableHandle &apValue)

Here is the call graph for this function:


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