Geant4_10
G4ITStepProcessor.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id: G4ITStepProcessor.cc 74551 2013-10-14 12:59:14Z gcosmo $
27 //
28 // Author: Mathieu Karamitros (kara (AT) cenbg . in2p3 . fr)
29 //
30 // History:
31 // -----------
32 // 10 Oct 2011 M.Karamitros created
33 //
34 // -------------------------------------------------------------------
35 
36 #include "G4ITStepProcessor.hh"
37 #include "G4UImanager.hh"
38 #include "G4ForceCondition.hh"
39 #include "G4GPILSelection.hh"
41 // #include "G4VSensitiveDetector.hh" // Include from 'hits/digi'
42 #include "G4GeometryTolerance.hh"
43 #include "G4ParticleTable.hh"
44 #include "G4ITTrackingManager.hh"
45 #include "G4TrackingInformation.hh"
46 #include "G4IT.hh"
47 #include "G4ITNavigator.hh" // Include from 'geometry'
48 
49 #include "G4VITProcess.hh"
50 #include "G4VProcess.hh"
51 #include "G4ITTransportation.hh"
52 
53 #include <iomanip> // Include from 'system'
54 #include <vector> // Include from 'system'
55 
56 using namespace std;
57 
58 static const size_t SizeOfSelectedDoItVector=100;
59 //static const size_t& gMaxNProcesses(G4VITProcess::GetMaxProcessIndex());
60 
61 //____________________________________________________________________________________
62 
64 {
65  verboseLevel = 0 ;
66  // fpUserSteppingAction = 0 ;
67  fStoreTrajectory = 0;
68  fpTrackingManager = 0;
69  fpNavigator = 0;
70  kCarTolerance = -1.;
71  fInitialized = false;
72  fPreviousTimeStep = DBL_MAX;
73  CleanProcessor();
74  ResetSecondaries();
75 }
76 
77 G4ITStepProcessor::G4ITStepProcessorState::G4ITStepProcessorState() :
79 // fSelectedAtRestDoItVector (gMaxNProcesses,0),
80 // fSelectedPostStepDoItVector (gMaxNProcesses,0)
81  fSelectedAtRestDoItVector (G4VITProcess::GetMaxProcessIndex(),0),
82  fSelectedPostStepDoItVector (G4VITProcess::GetMaxProcessIndex(),0)
83 {
84  fPhysicalStep = -1.;
85  fPreviousStepSize = -1.;
86 
87  fSafety = -1.;
88  proposedSafety = -1.;
89  endpointSafety = -1;
90 
91  fStepStatus = fUndefined;
92 
93  fTouchableHandle = 0;
94 }
95 
96 // should not be used
97 G4ITStepProcessor::G4ITStepProcessorState::G4ITStepProcessorState(const G4ITStepProcessorState& ) :
99 // fSelectedAtRestDoItVector (gMaxNProcesses,0),
100 // fSelectedPostStepDoItVector (gMaxNProcesses,0)
101  fSelectedAtRestDoItVector (G4VITProcess::GetMaxProcessIndex(),0),
102  fSelectedPostStepDoItVector (G4VITProcess::GetMaxProcessIndex(),0)
103 {
104  fPhysicalStep = -1.;
105  fPreviousStepSize = -1.;
106 
107  fSafety = -1.;
108  proposedSafety = -1.;
109  endpointSafety = -1;
110 
111  fStepStatus = fUndefined;
112 
113  fTouchableHandle = 0;
114 }
115 
116 // should not be used
117 G4ITStepProcessor::G4ITStepProcessorState& G4ITStepProcessor::G4ITStepProcessorState::operator=(const G4ITStepProcessorState& rhs)
118 {
119  if(this == &rhs) return *this;
120 
121  fSelectedAtRestDoItVector.clear();
122 // fSelectedAtRestDoItVector.resize(gMaxNProcesses,0);
123  fSelectedAtRestDoItVector.resize(G4VITProcess::GetMaxProcessIndex(),0);
124  fSelectedPostStepDoItVector.clear();
125 // fSelectedPostStepDoItVector.resize(gMaxNProcesses,0);
126  fSelectedPostStepDoItVector.resize(G4VITProcess::GetMaxProcessIndex(),0);
127 
128  fPhysicalStep = -1.;
129  fPreviousStepSize = -1.;
130 
131  fSafety = -1.;
132  proposedSafety = -1.;
133  endpointSafety = -1;
134 
135  fStepStatus = fUndefined;
136 
137  fTouchableHandle = 0;
138  return *this;
139 }
140 //____________________________________________________________________________________
141 
142 G4ITStepProcessor::G4ITStepProcessorState::~G4ITStepProcessorState()
143 {;}
144 //____________________________________________________________________________________
145 
147 {
148  std::map<const G4ParticleDefinition*, ProcessGeneralInfo*> ::iterator it;
149 
150  for(it = fProcessGeneralInfoMap.begin();it != fProcessGeneralInfoMap.end();it++)
151  {
152  if(it->second)
153  {
154  delete it->second;
155  it->second = 0;
156  }
157  }
158 
159  fProcessGeneralInfoMap.clear();
160 }
161 
162 //____________________________________________________________________________________
163 
165 {
166  fInitialized = false;
168  Initialize();
169 }
170 
171 //____________________________________________________________________________________
172 
174 {
175  CleanProcessor();
176  if(fInitialized) return;
177  // ActiveOnlyITProcess();
178 
180  ->GetNavigatorForTracking());
181 
182  fPhysIntLength = DBL_MAX;
183  kCarTolerance = 0.5*G4GeometryTolerance::GetInstance()->GetSurfaceTolerance();
184 
185  fInitialized = true;
186 }
187 //______________________________________________________________________________
188 
190 {
191  if(fpStep)
192  {
193  fpStep->DeleteSecondaryVector();
194  delete fpStep;
195  }
196 
197  if(fpSecondary) delete fpSecondary;
200 
201  // if(fpUserSteppingAction) delete fpUserSteppingAction;
202 }
203 //______________________________________________________________________________
204 // should not be used
206 {
207  verboseLevel = rhs.verboseLevel ;
208  fStoreTrajectory = rhs.fStoreTrajectory ;
209 
210  // fpUserSteppingAction = 0 ;
211  fpTrackingManager = 0;
212  fpNavigator = 0;
213  fInitialized = false;
214 
215  kCarTolerance = rhs.kCarTolerance;
216  fInitialized = false;
217  fPreviousTimeStep = DBL_MAX;
218 
219  CleanProcessor();
221 }
222 //______________________________________________________________________________
223 
225 {
226  if (this == &rhs) return *this; // handle self assignment
227  //assignment operator
228  return *this;
229 }
230 // ******************************************************************
231 
233 {
234  // Method not used for the time being
235 #ifdef debug
236  G4cout<<"G4ITStepProcessor::CloneProcesses: is called"<<G4endl;
237 #endif
238 
241 
242  theParticleIterator->reset();
243  // TODO : Ne faire la boucle que sur les IT **** !!!
244  while( (*theParticleIterator)() )
245  {
246  G4ParticleDefinition* particle = theParticleIterator->value();
247  G4ProcessManager* pm= particle->GetProcessManager();
248 
249  if(!pm)
250  {
251  G4cerr << "ERROR - G4ITStepProcessor::GetProcessNumber()" << G4endl
252  << " ProcessManager is NULL for particle = "
253  << particle->GetParticleName() << ", PDG_code = "
254  << particle->GetPDGEncoding() << G4endl;
255  G4Exception("G4ITStepProcessor::GetProcessNumber()", "ITStepProcessor0001",
256  FatalException, "Process Manager is not found.");
257  return;
258  }
259 
261  }
262 }
263 // ******************************************************************
264 
266 {
267  // Method not used for the time being
268  G4ProcessVector* processVector = processManager->GetProcessList();
269 
270  G4VITProcess* itProcess = 0 ;
271  for(int i = 0 ; i < processVector->size() ; i++)
272  {
273  G4VProcess* base_process = (*processVector)[i];
274  itProcess = dynamic_cast<G4VITProcess*>(base_process);
275 
276  if(!itProcess)
277  {
278  processManager->SetProcessActivation(base_process, false);
279  }
280  }
281 }
282 // ******************************************************************
284  G4ProcessManager* pm)
285 {
286 
287 #ifdef debug
288  G4cout<<"G4ITStepProcessor::GetProcessNumber: is called track"<<G4endl;
289 #endif
290  if(!pm)
291  {
292  G4cerr << "ERROR - G4SteppingManager::GetProcessNumber()" << G4endl
293  << " ProcessManager is NULL for particle = "
294  << particle->GetParticleName() << ", PDG_code = "
295  << particle->GetPDGEncoding() << G4endl;
296  G4Exception("G4SteppingManager::GetProcessNumber()", "ITStepProcessor0002",
297  FatalException, "Process Manager is not found.");
298  return;
299  }
300 
301  std::map<const G4ParticleDefinition*, ProcessGeneralInfo*>::iterator it = fProcessGeneralInfoMap.find(particle);
302  if(it != fProcessGeneralInfoMap.end())
303  {
304  G4Exception("G4SteppingManager::SetupGeneralProcessInfo()", "ITStepProcessor0003",
305  FatalException, "Process info already registered.");
306  return;
307  }
308 
309  // here used as temporary
310  fpProcessInfo = new ProcessGeneralInfo();
311 
312  // AtRestDoits
313  fpProcessInfo->MAXofAtRestLoops = pm->GetAtRestProcessVector()->entries();
314  fpProcessInfo->fpAtRestDoItVector = pm->GetAtRestProcessVector(typeDoIt);
315  fpProcessInfo->fpAtRestGetPhysIntVector = pm->GetAtRestProcessVector(typeGPIL);
316 #ifdef debug
317  G4cout << "G4ITStepProcessor::GetProcessNumber: #ofAtRest="
318  << fpProcessInfo->MAXofAtRestLoops << G4endl;
319 #endif
320 
321  // AlongStepDoits
322  fpProcessInfo->MAXofAlongStepLoops = pm->GetAlongStepProcessVector()->entries();
323  fpProcessInfo->fpAlongStepDoItVector = pm->GetAlongStepProcessVector(typeDoIt);
324  fpProcessInfo->fpAlongStepGetPhysIntVector = pm->GetAlongStepProcessVector(typeGPIL);
325 #ifdef debug
326  G4cout << "G4ITStepProcessor::GetProcessNumber:#ofAlongStp="
327  << fpProcessInfo->MAXofAlongStepLoops << G4endl;
328 #endif
329 
330  // PostStepDoits
331  fpProcessInfo->MAXofPostStepLoops = pm->GetPostStepProcessVector()->entries();
332  fpProcessInfo->fpPostStepDoItVector = pm->GetPostStepProcessVector(typeDoIt);
333  fpProcessInfo->fpPostStepGetPhysIntVector = pm->GetPostStepProcessVector(typeGPIL);
334 #ifdef debug
335  G4cout << "G4ITStepProcessor::GetProcessNumber: #ofPostStep="
336  << fpProcessInfo->MAXofPostStepLoops << G4endl;
337 #endif
338 
339  if (SizeOfSelectedDoItVector<fpProcessInfo->MAXofAtRestLoops ||
340  SizeOfSelectedDoItVector<fpProcessInfo->MAXofAlongStepLoops ||
341  SizeOfSelectedDoItVector<fpProcessInfo->MAXofPostStepLoops )
342  {
343  G4cerr << "ERROR - G4ITStepProcessor::GetProcessNumber()" << G4endl
344  << " SizeOfSelectedDoItVector= " << SizeOfSelectedDoItVector
345  << " ; is smaller then one of MAXofAtRestLoops= "
346  << fpProcessInfo->MAXofAtRestLoops << G4endl
347  << " or MAXofAlongStepLoops= " << fpProcessInfo->MAXofAlongStepLoops
348  << " or MAXofPostStepLoops= " << fpProcessInfo->MAXofPostStepLoops << G4endl;
349  G4Exception("G4ITStepProcessor::GetProcessNumber()",
350  "ITStepProcessor0004", FatalException,
351  "The array size is smaller than the actual No of processes.");
352  }
353 
354  if(!fpProcessInfo->fpAtRestDoItVector &&
355  !fpProcessInfo->fpAlongStepDoItVector &&
356  !fpProcessInfo->fpPostStepDoItVector)
357  {
358  G4ExceptionDescription exceptionDescription ;
359  exceptionDescription << "No DoIt process found " ;
360  G4Exception("G4ITStepProcessor::DoStepping","ITStepProcessor0005",
361  FatalErrorInArgument,exceptionDescription);
362  return ;
363  }
364 
365  if(fpProcessInfo->fpAlongStepGetPhysIntVector && fpProcessInfo->MAXofAlongStepLoops>0)
366  {
367  fpProcessInfo->fpTransportation = dynamic_cast<G4ITTransportation*>
368  ((*fpProcessInfo->fpAlongStepGetPhysIntVector)[fpProcessInfo->MAXofAlongStepLoops-1]);
369 
370  if(fpProcessInfo->fpTransportation == 0)
371  {
372  G4ExceptionDescription exceptionDescription ;
373  exceptionDescription << "No transportation process found " ;
374  G4Exception("G4ITStepProcessor::SetupGeneralProcessInfo","ITStepProcessor0006",
375  FatalErrorInArgument,exceptionDescription);
376  }
377  }
378  fProcessGeneralInfoMap[particle] = fpProcessInfo;
379  // fpProcessInfo = 0;
380 }
381 
382 // ******************************************************************
383 
385 {
386  fpTrack = track ;
387  if(fpTrack)
388  {
389  fpITrack = GetIT(fpTrack) ;
390  fpStep = const_cast<G4Step*>(fpTrack -> GetStep());
391 
392  if(fpITrack)
393  {
394  fpTrackingInfo = fpITrack->GetTrackingInfo() ;
395  }
396  else
397  {
398  fpTrackingInfo = 0;
399  G4cerr << "Track ID : " << fpTrack->GetTrackID() << G4endl;
400 
401  G4ExceptionDescription exceptionDescription ("No IT pointer was attached to the track you try to process.");
402  G4Exception("G4ITStepProcessor::SetTrack","ITStepProcessor0007",
403  FatalErrorInArgument,exceptionDescription);
404  }
405  }
406  else
407  {
408  fpITrack = 0;
409  fpStep = 0 ;
410  }
411 }
412 //______________________________________________________________________________
413 
415 {
416  G4ParticleDefinition* particle = fpTrack->GetDefinition();
417  std::map<const G4ParticleDefinition*, ProcessGeneralInfo*>::iterator it = fProcessGeneralInfoMap.find(particle);
418 
419  if(it == fProcessGeneralInfoMap.end())
420  {
422  if(fpProcessInfo == 0)
423  {
424  G4ExceptionDescription exceptionDescription ("...");
425  G4Exception("G4ITStepProcessor::GetProcessNumber","ITStepProcessor0008",
426  FatalErrorInArgument,exceptionDescription);
427  return;
428  }
429  }
430  else
431  {
432  fpProcessInfo = it->second;
433  }
434 }
435 //______________________________________________________________________________
436 
438 {
439  fpSecondary = fpStep->GetfSecondary();
440  fpPreStepPoint = fpStep->GetPreStepPoint();
441  fpPostStepPoint = fpStep->GetPostStepPoint();
442 
443  fpState = (G4ITStepProcessorState*) fpITrack->GetTrackingInfo()->GetStepProcessorState();
444 
445  GetProcessInfo();
447 }
448 //______________________________________________________________________________
449 
451 {
452  // Reset the secondary particles
453  fN2ndariesAtRestDoIt = 0;
454  fN2ndariesAlongStepDoIt = 0;
455  fN2ndariesPostStepDoIt = 0;
456 }
457 //______________________________________________________________________________
458 
460 {
461  // Select the rest process which has the shortest time before
462  // it is invoked. In rest processes, GPIL()
463  // returns the time before a process occurs.
464  G4double lifeTime (DBL_MAX), shortestLifeTime (DBL_MAX);
465 
466  fAtRestDoItProcTriggered = 0;
467  shortestLifeTime = DBL_MAX;
468 
469  unsigned int NofInactiveProc=0;
470 
471  for( size_t ri=0 ; ri < fpProcessInfo->MAXofAtRestLoops ; ri++ )
472  {
473  fpCurrentProcess = (G4VITProcess*) (*fpProcessInfo->fpAtRestGetPhysIntVector)[ri];
474  if (fpCurrentProcess== 0)
475  {
476  (fpState->fSelectedAtRestDoItVector)[ri] = InActivated;
477  NofInactiveProc++;
478  continue;
479  } // NULL means the process is inactivated by a user on fly.
480 
481  fCondition=NotForced;
482  fpCurrentProcess->SetProcessState(fpTrackingInfo->GetProcessState(fpCurrentProcess->GetProcessID()));
483  lifeTime = fpCurrentProcess->AtRestGPIL( *fpTrack, &fCondition );
484  fpCurrentProcess->SetProcessState(0);
485 
486  if(fCondition==Forced)
487  {
488  (fpState->fSelectedAtRestDoItVector)[ri] = Forced;
489  }
490  else
491  {
492  (fpState->fSelectedAtRestDoItVector)[ri] = InActivated;
493  if(lifeTime < shortestLifeTime )
494  {
495  shortestLifeTime = lifeTime;
496  fAtRestDoItProcTriggered = G4int(ri);
497  (fpState->fSelectedAtRestDoItVector)[fAtRestDoItProcTriggered] = NotForced;
498  }
499  }
500  }
501 
502  fTimeStep = shortestLifeTime ;
503 
504  // at least one process is necessary to destroy the particle
505  // exit with warning
506  if(NofInactiveProc==fpProcessInfo->MAXofAtRestLoops)
507  {
508  G4cerr << "ERROR - G4ITStepProcessor::InvokeAtRestDoItProcs()" << G4endl
509  << " No AtRestDoIt process is active!" << G4endl;
510  }
511 }
512 //___________________________________________________________________________
513 
515 {
516  SetTrack(track);
518 }
519 //______________________________________________________________________________
520 
521 
523 {
524  // DEBUG
525  // G4cout << "SetInitialStep for : " << fpITrack-> GetName() << G4endl;
526  //________________________________________________________
527  // Initialize geometry
528 
529 
530  if ( ! fpTrack->GetTouchableHandle())
531  {
532  G4ThreeVector direction= fpTrack->GetMomentumDirection();
533  fpNavigator->LocateGlobalPointAndSetup( fpTrack->GetPosition(),
534  &direction, false, false );
535  fpState->fTouchableHandle = fpNavigator->CreateTouchableHistory();
536 
537  fpTrack->SetTouchableHandle( fpState->fTouchableHandle );
538  fpTrack->SetNextTouchableHandle( fpState->fTouchableHandle );
539  }
540  else
541  {
542  fpState->fTouchableHandle = fpTrack->GetTouchableHandle();
543  fpTrack->SetNextTouchableHandle( fpState->fTouchableHandle );
544  G4VPhysicalVolume* oldTopVolume= fpTrack->GetTouchableHandle()->GetVolume();
545  G4VPhysicalVolume* newTopVolume=
546  fpNavigator->ResetHierarchyAndLocate( fpTrack->GetPosition(),
547  fpTrack->GetMomentumDirection(),
548  *((G4TouchableHistory*)fpTrack->GetTouchableHandle()()) );
549  if(newTopVolume != oldTopVolume || oldTopVolume->GetRegularStructureId() == 1 )
550  {
551  fpState->fTouchableHandle = fpNavigator->CreateTouchableHistory();
552  fpTrack->SetTouchableHandle( fpState->fTouchableHandle );
553  fpTrack->SetNextTouchableHandle( fpState->fTouchableHandle );
554  }
555  }
556 
557  fpCurrentVolume = fpState->fTouchableHandle->GetVolume();
558 
559  //________________________________________________________
560  // If the primary track has 'Suspend' or 'PostponeToNextEvent' state,
561  // set the track state to 'Alive'.
562  if( (fpTrack->GetTrackStatus()==fSuspend) ||
563  (fpTrack->GetTrackStatus()==fPostponeToNextEvent) )
564  {
565  fpTrack->SetTrackStatus(fAlive);
566  }
567 
568  // If the primary track has 'zero' kinetic energy, set the track
569  // state to 'StopButAlive'.
570  if(fpTrack->GetKineticEnergy() <= 0.0)
571  {
572  fpTrack->SetTrackStatus( fStopButAlive );
573  }
574  //________________________________________________________
575  // Set vertex information of G4Track at here
576  if ( fpTrack->GetCurrentStepNumber() == 0 )
577  {
578  fpTrack->SetVertexPosition( fpTrack->GetPosition() );
579  fpTrack->SetVertexMomentumDirection( fpTrack->GetMomentumDirection() );
580  fpTrack->SetVertexKineticEnergy( fpTrack->GetKineticEnergy() );
581  fpTrack->SetLogicalVolumeAtVertex( fpTrack->GetVolume()->GetLogicalVolume() );
582  }
583  //________________________________________________________
584  // If track is already outside the world boundary, kill it
585  if( fpCurrentVolume==0 )
586  {
587  // If the track is a primary, stop processing
588  if(fpTrack->GetParentID()==0)
589  {
590  G4cerr << "ERROR - G4ITStepProcessor::SetInitialStep()" << G4endl
591  << " Primary particle starting at - "
592  << fpTrack->GetPosition()
593  << " - is outside of the world volume." << G4endl;
594  G4Exception("G4ITStepProcessor::SetInitialStep()", "ITStepProcessor0011",
595  FatalException, "Primary vertex outside of the world!");
596  }
597 
598  fpTrack->SetTrackStatus( fStopAndKill );
599  G4cout << "WARNING - G4ITStepProcessor::SetInitialStep()" << G4endl
600  << " Initial track position is outside world! - "
601  << fpTrack->GetPosition() << G4endl;
602  }
603  else{
604  // Initial set up for attribues of 'Step'
605  fpStep->InitializeStep( fpTrack );
606  }
607 
608 
609  if( fpTrack->GetTrackStatus() == fStopAndKill ) return ;
610 
611  fpTrackingManager->StartTracking(fpTrack);
612 
613  fpState->fStepStatus = fUndefined;
614 }
615 //______________________________________________________________________________
616 
618 {
619 
620  if(!fpStep)
621  {
622 
623  // Create new Step and give it to the track
624  fpStep = new G4Step();
625  fpTrack->SetStep(fpStep);
626  fpSecondary = fpStep->NewSecondaryVector();
627 
628  // Create new state and set it in the trackingInfo
629  fpState = new G4ITStepProcessorState();
631 
632  SetupMembers();
633  fpNavigator->NewNavigatorState();
634 
635  SetInitialStep();
636  }
637  else
638  {
639  SetupMembers();
640 
641  fpState->fPreviousStepSize = fpTrack->GetStepLength();
642 /***
643  // Send G4Step information to Hit/Dig if the volume is sensitive
644  fpCurrentVolume = fpStep->GetPreStepPoint()->GetPhysicalVolume();
645  StepControlFlag = fpStep->GetControlFlag();
646  if( fpCurrentVolume != 0 && StepControlFlag != AvoidHitInvocation)
647  {
648  fpSensitive = fpStep->GetPreStepPoint()->
649  GetSensitiveDetector();
650 
651  // if( fSensitive != 0 ) {
652  // fSensitive->Hit(fStep);
653  // }
654  }
655 ***/
656  // Store last PostStepPoint to PreStepPoint, and swap current and next
657  // volume information of G4Track. Reset total energy deposit in one Step.
658  fpStep->CopyPostToPreStepPoint();
659  fpStep->ResetTotalEnergyDeposit();
660 
661  //JA Set the volume before it is used (in DefineStepLength() for User Limit)
662  fpCurrentVolume = fpStep->GetPreStepPoint()->GetPhysicalVolume();
663 /*
664  G4cout << G4endl;
665  G4cout << "!!!!!!!!!!!!!!!!!!!!!!!!!!!" << G4endl;
666  G4cout << "PreStepPoint Volume : " << fpCurrentVolume->GetName() << G4endl;
667  G4cout << "Track Touchable : " << fpTrack->GetTouchableHandle()->GetVolume()->GetName() << G4endl;
668  G4cout << "Track NextTouchable : " << fpTrack->GetNextTouchableHandle()->GetVolume()->GetName() << G4endl;
669 */
670  // Reset the step's auxiliary points vector pointer
672 
673  // Switch next touchable in track to current one
674  fpTrack->SetTouchableHandle(fpTrack->GetNextTouchableHandle());
675  fpState->fTouchableHandle = fpTrack->GetTouchableHandle();
676  fpTrack->SetNextTouchableHandle( fpState->fTouchableHandle );
677  G4VPhysicalVolume* oldTopVolume= fpTrack->GetTouchableHandle()->GetVolume();
678  fpNavigator->SetNavigatorState(fpITrack->GetTrackingInfo()->GetNavigatorState());
679 
680  G4VPhysicalVolume* newTopVolume=
681  fpNavigator->ResetHierarchyAndLocate( fpTrack->GetPosition(),
682  fpTrack->GetMomentumDirection(),
683  *((G4TouchableHistory*)fpTrack->GetTouchableHandle()()) );
684 
685  // G4cout << "New Top Volume : " << newTopVolume->GetName() << G4endl;
686 
687  if(newTopVolume != oldTopVolume || oldTopVolume->GetRegularStructureId() == 1 )
688  {
689  fpState->fTouchableHandle = fpNavigator->CreateTouchableHistory();
690  fpTrack->SetTouchableHandle( fpState->fTouchableHandle );
691  fpTrack->SetNextTouchableHandle( fpState->fTouchableHandle );
692  }
693 
694  fpNavigator->SetNavigatorState(fpITrack->GetTrackingInfo()->GetNavigatorState());
695  }
696 }
697 
698 //______________________________________________________________________________
699 
700 
701 // ************************************************************************
702 // Compute Interaction Length
703 // ************************************************************************
705 {
706 
707  InitDefineStep();
708 
709  G4TrackStatus trackStatus = fpTrack -> GetTrackStatus() ;
710 
711  if(trackStatus == fStopAndKill)
712  {
713  return ;
714  }
715 
716  if(trackStatus == fStopButAlive)
717  {
718  fpITrack->GetTrackingInfo()->SetNavigatorState(fpNavigator->GetNavigatorState());
719  fpNavigator->SetNavigatorState(0);
720  return GetAtRestIL() ;
721  }
722 
723 
724  // Find minimum Step length and corresponding time
725  // demanded by active disc./cont. processes
726 
727  // ReSet the counter etc.
728  fpState->fPhysicalStep = DBL_MAX; // Initialize by a huge number
729  fPhysIntLength = DBL_MAX; // Initialize by a huge number
730 
731  double proposedTimeStep = DBL_MAX;
732  G4VProcess* processWithPostStepGivenByTimeStep(0);
733 
734  // GPIL for PostStep
735  fPostStepDoItProcTriggered = fpProcessInfo->MAXofPostStepLoops;
736  fPostStepAtTimeDoItProcTriggered = fpProcessInfo->MAXofPostStepLoops;
737 
738  // G4cout << "fpProcessInfo->MAXofPostStepLoops : " << fpProcessInfo->MAXofPostStepLoops
739  // << " mol : " << fpITrack -> GetName() << " id : " << fpTrack->GetTrackID()
740  // << G4endl;
741 
742  for(size_t np=0; np < fpProcessInfo->MAXofPostStepLoops; np++)
743  {
744  fpCurrentProcess = (G4VITProcess*) (*fpProcessInfo->fpPostStepGetPhysIntVector)[np];
745  if (fpCurrentProcess== 0)
746  {
747  (fpState->fSelectedPostStepDoItVector)[np] = InActivated;
748  continue;
749  } // NULL means the process is inactivated by a user on fly.
750 
751  fCondition=NotForced;
752  fpCurrentProcess->SetProcessState(fpTrackingInfo->GetProcessState(fpCurrentProcess->GetProcessID()));
753 
754  // G4cout << "Is going to call : " << fpCurrentProcess -> GetProcessName() << G4endl;
755  fPhysIntLength = fpCurrentProcess->
756  PostStepGPIL( *fpTrack,
757  fpState->fPreviousStepSize,
758  &fCondition );
759  fpCurrentProcess->SetProcessState(0);
760 
761  switch (fCondition)
762  {
763  case ExclusivelyForced: // Will need special treatment
764  (fpState->fSelectedPostStepDoItVector)[np] = ExclusivelyForced;
765  fpState->fStepStatus = fExclusivelyForcedProc;
766  fpStep->GetPostStepPoint()
767  ->SetProcessDefinedStep(fpCurrentProcess);
768  break;
769 
770  case Conditionally:
771  // (fpState->fSelectedPostStepDoItVector)[np] = Conditionally;
772  G4Exception("G4ITStepProcessor::DefinePhysicalStepLength()", "ITStepProcessor0008",
773  FatalException, "This feature is no more supported");
774  break;
775 
776  case Forced:
777  (fpState->fSelectedPostStepDoItVector)[np] = Forced;
778  break;
779 
780  case StronglyForced:
781  (fpState->fSelectedPostStepDoItVector)[np] = StronglyForced;
782  break;
783 
784  default:
785  (fpState->fSelectedPostStepDoItVector)[np] = InActivated;
786  break;
787  }
788 
789  if (fCondition==ExclusivelyForced)
790  {
791  for(size_t nrest=np+1; nrest < fpProcessInfo->MAXofPostStepLoops; nrest++)
792  {
793  (fpState->fSelectedPostStepDoItVector)[nrest] = InActivated;
794  }
795  return; // Please note the 'return' at here !!!
796  }
797  else
798  {
799  if(fPhysIntLength < fpState->fPhysicalStep )
800  {
801  // To avoid checking whether the process is actually
802  // proposing a time step, the returned time steps are
803  // negative (just for tagging)
804  if(fpCurrentProcess->ProposesTimeStep())
805  {
806  fPhysIntLength *= -1;
807  if(fPhysIntLength < proposedTimeStep)
808  {
809  proposedTimeStep = fPhysIntLength;
810  fPostStepAtTimeDoItProcTriggered = np;
811  processWithPostStepGivenByTimeStep = fpCurrentProcess;
812  }
813  }
814  else
815  {
816  fpState->fPhysicalStep = fPhysIntLength;
817  fpState->fStepStatus = fPostStepDoItProc;
818  fPostStepDoItProcTriggered = G4int(np);
819  fpStep->GetPostStepPoint()
820  ->SetProcessDefinedStep(fpCurrentProcess);
821  }
822  }
823  }
824  }
825 
826  // GPIL for AlongStep
827  fpState->proposedSafety = DBL_MAX;
828  G4double safetyProposedToAndByProcess = fpState->proposedSafety;
829 
830  for(size_t kp=0; kp < fpProcessInfo->MAXofAlongStepLoops; kp++)
831  {
832  fpCurrentProcess = (G4VITProcess*) (*fpProcessInfo->fpAlongStepGetPhysIntVector)[kp];
833  if (fpCurrentProcess== 0) continue;
834  // NULL means the process is inactivated by a user on fly.
835 
836  fpCurrentProcess->SetProcessState(fpTrackingInfo->GetProcessState(fpCurrentProcess->GetProcessID()));
837  fPhysIntLength = fpCurrentProcess-> AlongStepGPIL( *fpTrack,
838  fpState->fPreviousStepSize,
839  fpState->fPhysicalStep,
840  safetyProposedToAndByProcess,
841  &fGPILSelection );
842 
843  if(fPhysIntLength < fpState->fPhysicalStep)
844  {
845  fpState->fPhysicalStep = fPhysIntLength;
846  // Should save PS and TS in IT
847 
848  // Check if the process wants to be the GPIL winner. For example,
849  // multi-scattering proposes Step limit, but won't be the winner.
850  if(fGPILSelection==CandidateForSelection)
851  {
852  fpState->fStepStatus = fAlongStepDoItProc;
853  fpStep->GetPostStepPoint()
854  ->SetProcessDefinedStep(fpCurrentProcess);
855  }
856 
857  // Transportation is assumed to be the last process in the vector
858  if(kp == fpProcessInfo->MAXofAlongStepLoops-1)
859  {
860  fpTransportation = dynamic_cast<G4ITTransportation*>(fpCurrentProcess);
861 
862  if(! fpTransportation)
863  {
864  G4ExceptionDescription exceptionDescription ;
865  exceptionDescription << "No transportation process found " ;
866  G4Exception("G4ITStepProcessor::DoDefinePhysicalStepLength","ITStepProcessor0009",
867  FatalErrorInArgument,exceptionDescription);
868  }
869 
870  fTimeStep = fpTransportation->GetInteractionTimeLeft();
871 
872 
873  if (fpTrack->GetNextVolume() != 0)
874  fpState->fStepStatus = fGeomBoundary;
875  else
876  fpState->fStepStatus = fWorldBoundary;
877  }
878  }
879  else
880  {
881  if(kp == fpProcessInfo->MAXofAlongStepLoops-1)
882  {
883  fpTransportation = dynamic_cast<G4ITTransportation*>(fpCurrentProcess);
884 
885  if(! fpTransportation)
886  {
887  G4ExceptionDescription exceptionDescription ;
888  exceptionDescription << "No transportation process found " ;
889  G4Exception("G4ITStepProcessor::DoDefinePhysicalStepLength","ITStepProcessor0010",
890  FatalErrorInArgument,exceptionDescription);
891  }
892 
893  fTimeStep = fpTransportation->GetInteractionTimeLeft();
894  }
895  }
896 
897  if(proposedTimeStep < fTimeStep)
898  {
899  if(fPostStepAtTimeDoItProcTriggered<fpProcessInfo->MAXofPostStepLoops)
900  {
901  if ((fpState->fSelectedPostStepDoItVector)[fPostStepAtTimeDoItProcTriggered] ==
902  InActivated)
903  {
904  (fpState->fSelectedPostStepDoItVector)[fPostStepAtTimeDoItProcTriggered] = NotForced;
905  // (fpState->fSelectedPostStepDoItVector)[fPostStepDoItProcTriggered] = InActivated;
906 
907  fpState->fStepStatus = fPostStepDoItProc;
908  fpStep->GetPostStepPoint()->SetProcessDefinedStep(processWithPostStepGivenByTimeStep);
909 
910  fTimeStep = proposedTimeStep;
911 
912  fpTransportation->ComputeStep(*fpTrack,*fpStep,fTimeStep,fpState->fPhysicalStep);
913  }
914  }
915  }
916  else
917  {
918  if (fPostStepDoItProcTriggered<fpProcessInfo->MAXofPostStepLoops)
919  {
920  if ((fpState->fSelectedPostStepDoItVector)[fPostStepDoItProcTriggered] ==
921  InActivated)
922  {
923  (fpState->fSelectedPostStepDoItVector)[fPostStepDoItProcTriggered] =
924  NotForced;
925  }
926  }
927  }
928 
929  fpCurrentProcess->SetProcessState(0);
930 
931  // Make sure to check the safety, even if Step is not limited
932  // by this process. J. Apostolakis, June 20, 1998
933  //
934  if (safetyProposedToAndByProcess < fpState->proposedSafety)
935  // proposedSafety keeps the smallest value:
936  fpState->proposedSafety = safetyProposedToAndByProcess;
937  else
938  // safetyProposedToAndByProcess always proposes a valid safety:
939  safetyProposedToAndByProcess = fpState->proposedSafety;
940 
941  }
942 
943  fpITrack->GetTrackingInfo()->SetNavigatorState(fpNavigator->GetNavigatorState());
944  fpNavigator->SetNavigatorState(0);
945 }
946 
947 //______________________________________________________________________________
void SetTrackStatus(const G4TrackStatus aTrackStatus)
G4ParticleDefinition * GetDefinition() const
G4int GetParentID() const
virtual G4VPhysicalVolume * ResetHierarchyAndLocate(const G4ThreeVector &point, const G4ThreeVector &direction, const G4TouchableHistory &h)
void SetVertexMomentumDirection(const G4ThreeVector &aValue)
void SetProcessDefinedStep(const G4VProcess *aValue)
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
static const size_t & GetMaxProcessIndex()
void DeleteSecondaryVector()
void SetStepProcessorState(G4ITStepProcessorState_Lock *)
void SetNavigator(G4ITNavigator *value)
const G4ThreeVector & GetPosition() const
G4double GetSurfaceTolerance() const
G4TrackStatus GetTrackStatus() const
void SetNextTouchableHandle(const G4TouchableHandle &apValue)
void SetTouchableHandle(const G4TouchableHandle &apValue)
virtual G4VPhysicalVolume * LocateGlobalPointAndSetup(const G4ThreeVector &point, const G4ThreeVector *direction=0, const G4bool pRelativeSearch=true, const G4bool ignoreDirection=true)
void SetPointerToVectorOfAuxiliaryPoints(std::vector< G4ThreeVector > *theNewVectorPointer)
Definition: G4Step.hh:237
G4ProcessManager * GetProcessManager() const
int G4int
Definition: G4Types.hh:78
virtual void StartTracking(G4Track *)
const G4String & GetParticleName() const
G4VProcess * SetProcessActivation(G4VProcess *aProcess, G4bool fActive)
G4VPhysicalVolume * GetNextVolume() const
void SetTrack(G4Track *)
void DefinePhysicalStepLength(G4Track *)
G4double AtRestGPIL(const G4Track &track, G4ForceCondition *condition)
Definition: G4VProcess.hh:490
G4StepPoint * GetPreStepPoint() const
G4ITNavigatorState_Lock * GetNavigatorState()
G4IT * GetIT(const G4Track *track)
Definition: G4IT.cc:48
G4ProcessVector * GetAtRestProcessVector(G4ProcessVectorTypeIndex typ=typeGPIL) const
G4double GetKineticEnergy() const
G4int entries() const
static G4ITTransportationManager * GetTransportationManager()
G4GLOB_DLL std::ostream G4cout
void reset(G4bool ifSkipIon=true)
G4int GetCurrentStepNumber() const
G4TrackVector * NewSecondaryVector()
G4VPhysicalVolume * GetPhysicalVolume() const
virtual void ComputeStep(const G4Track &, const G4Step &, const double timeStep, double &spaceStep)
virtual G4int GetRegularStructureId() const =0
void ResetTotalEnergyDeposit()
G4ITStepProcessorState_Lock * GetStepProcessorState()
void SetVertexKineticEnergy(const G4double aValue)
const G4TouchableHandle & GetNextTouchableHandle() const
G4ITNavigatorState_Lock * GetNavigatorState() const
Definition: G4Step.hh:76
G4int GetTrackID() const
G4double GetInteractionTimeLeft()
void SetProcessState(G4ProcessState_Lock *aProcInfo)
Definition: G4VITProcess.hh:90
const G4TouchableHandle & GetTouchableHandle() const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
virtual void Initialize()
G4int size() const
G4TrackingInformation * GetTrackingInfo()
Definition: G4IT.hh:134
const G4ThreeVector & GetMomentumDirection() const
G4LogicalVolume * GetLogicalVolume() const
G4TrackVector * GetfSecondary()
static G4ParticleTable * GetParticleTable()
G4ProcessState_Lock * GetProcessState(size_t index)
virtual G4VPhysicalVolume * GetVolume(G4int depth=0) const
Definition: G4VTouchable.cc:44
void SetVertexPosition(const G4ThreeVector &aValue)
G4StepPoint * GetPostStepPoint() const
G4ITStepProcessor & operator=(const G4ITStepProcessor &other)
void SetNavigatorState(G4ITNavigatorState_Lock *)
G4VPhysicalVolume * GetVolume() const
void SetNavigatorState(G4ITNavigatorState_Lock *)
#define G4endl
Definition: G4ios.hh:61
void InitializeStep(G4Track *aValue)
double G4double
Definition: G4Types.hh:76
G4TouchableHistory * CreateTouchableHistory() const
void SetupGeneralProcessInfo(G4ParticleDefinition *, G4ProcessManager *)
G4ProcessVector * GetAlongStepProcessVector(G4ProcessVectorTypeIndex typ=typeGPIL) const
G4bool ProposesTimeStep() const
void NewNavigatorState()
#define DBL_MAX
Definition: templates.hh:83
size_t GetProcessID() const
Definition: G4VITProcess.hh:80
G4TrackStatus
G4PTblDicIterator * GetIterator() const
void SetStep(const G4Step *aValue)
#define theParticleIterator
G4double GetStepLength() const
static G4GeometryTolerance * GetInstance()
G4ProcessVector * GetProcessList() const
G4GLOB_DLL std::ostream G4cerr
void SetLogicalVolumeAtVertex(const G4LogicalVolume *)
void CopyPostToPreStepPoint()
G4ProcessVector * GetPostStepProcessVector(G4ProcessVectorTypeIndex typ=typeGPIL) const