Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4ErrorPropagator.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: G4ErrorPropagator.cc 91809 2015-08-06 13:00:09Z gcosmo $
27 //
28 // ------------------------------------------------------------
29 // GEANT 4 class implementation file
30 // ------------------------------------------------------------
31 //
32 
33 #include "G4ErrorPropagator.hh"
34 #include "G4ErrorPropagatorData.hh"
35 #include "G4ErrorFreeTrajState.hh"
38 #include "G4ErrorSurfaceTarget.hh"
39 
40 #include "G4SystemOfUnits.hh"
41 #include "G4DynamicParticle.hh"
42 #include "G4Track.hh"
43 #include "G4SteppingManager.hh"
44 #include "G4EventManager.hh"
45 #include "G4TrackingManager.hh"
46 #include "G4ParticleTable.hh"
47 #include "G4StateManager.hh"
48 
49 #include "G4VPhysicalVolume.hh"
50 #include "G4PhysicalVolumeStore.hh"
51 #include "G4UnitsTable.hh"
52 
53 #include <vector>
54 
55 
56 //---------------------------------------------------------------------------
58  : theStepLength(0.), theInitialTrajState(0), theStepN(0), theG4Track(0)
59 {
61 #ifdef G4EVERBOSE
62  if(verbose >= 5) { G4cout << "G4ErrorPropagator " << this << G4endl; }
63 #endif
64 
65  fpSteppingManager = G4EventManager::GetEventManager()
67  thePropIsInitialized = false;
68 }
69 
70 
71 //-----------------------------------------------------------------------
73  const G4ErrorTarget* target, G4ErrorMode mode )
74 {
75  // to start ierror is set to 1 (= OK)
76  //
77  G4int ierr = 1;
78 
79  G4ErrorPropagatorData* g4edata =
81 
82  //----- Do not propagate zero or too low energy particles
83  //
84  if( currentTS->GetMomentum().mag() < 1.E-9*MeV )
85  {
86  std::ostringstream message;
87  message << "Energy too low to be propagated: "
88  << G4BestUnit(currentTS->GetMomentum().mag(),"Energy");
89  G4Exception("G4ErrorPropagator::Propagate()", "GEANT4e-Notification",
90  JustWarning, message);
91  return -3;
92  }
93 
94  g4edata->SetMode( mode );
95 
96 #ifdef G4EVERBOSE
97  if( verbose >= 1 )
98  {
99  G4cout << " =====> starting GEANT4E tracking for "
100  << currentTS->GetParticleType()
101  << " Forwards= " << g4edata->GetMode() << G4endl;
102  }
103  if(verbose >= 1 )
104  {
105  G4cout << G4endl << "@@@@@@@@@@@@@@@@@@@@@@@@@ NEW STEP " << G4endl;
106  }
107 
108  if( verbose >= 3 )
109  {
110  G4cout << " G4ErrorPropagator::Propagate initialTS ";
111  G4cout << *currentTS << G4endl;
112  target->Dump(G4String(" to target "));
113  }
114 #endif
115 
116  g4edata->SetTarget( target );
117 
118  //----- Create a track
119  //
120  if( theG4Track != 0 ) { delete theG4Track; }
121  theG4Track = InitG4Track( *currentTS );
122 
123  //----- Create a G4ErrorFreeTrajState
124  //
125  G4ErrorFreeTrajState* currentTS_FREE = InitFreeTrajState( currentTS );
126 
127  //----- Track the particle
128  //
129  ierr = MakeSteps( currentTS_FREE );
130 
131  //------ Tracking ended, check if target has been reached
132  // if target not found
133  //
134  if( g4edata->GetState() != G4ErrorState_StoppedAtTarget )
135  {
136  if( theG4Track->GetKineticEnergy() > 0. )
137  {
138  ierr = -ierr - 10;
139  }
140  else
141  {
142  ierr = -ierr - 20;
143  }
144  *currentTS = *currentTS_FREE;
145  if(verbose >= 0 )
146  {
147  std::ostringstream message;
148  message << "Particle does not reach target: " << *currentTS;
149  G4Exception("G4ErrorPropagator::Propagate()", "GEANT4e-Notification",
150  JustWarning, message);
151  }
152  }
153  else
154  {
155  GetFinalTrajState( currentTS, currentTS_FREE, target );
156  }
157 
158 #ifdef G4EVERBOSE
159  if( verbose >= 1 )
160  {
161  G4cout << " G4ErrorPropagator: propagation ended " << G4endl;
162  }
163  if( verbose >= 2 )
164  {
165  G4cout << " Current TrajState " << currentTS << G4endl;
166  }
167 #endif
168 
169  // Inform end of tracking to physics processes
170  //
171  theG4Track->GetDefinition()->GetProcessManager()->EndTracking();
172 
173  InvokePostUserTrackingAction( theG4Track );
174 
175  // delete currentTS_FREE;
176 
177  return ierr;
178 }
179 
180 
181 //-----------------------------------------------------------------------
183 {
184  G4ErrorPropagatorData* g4edata =
186 
187  if ( (g4edata->GetState() == G4ErrorState_PreInit)
189  != G4State_GeomClosed) )
190  {
191  std::ostringstream message;
192  message << "Called before initialization is done for this track!";
193  G4Exception("G4ErrorPropagator::PropagateOneStep()",
194  "InvalidCall", FatalException, message,
195  "Please call G4ErrorPropagatorManager::InitGeant4e().");
196  }
197 
198  // to start ierror is set to 0 (= OK)
199  //
200  G4int ierr = 0;
201 
202  //--- Do not propagate zero or too low energy particles
203  //
204  if( currentTS->GetMomentum().mag() < 1.E-9*MeV )
205  {
206  std::ostringstream message;
207  message << "Energy too low to be propagated: "
208  << G4BestUnit(currentTS->GetMomentum().mag(),"Energy");
209  G4Exception("G4ErrorPropagator::PropagateOneStep()",
210  "GEANT4e-Notification", JustWarning, message);
211  return -3;
212  }
213 
214 #ifdef G4EVERBOSE
215  if( verbose >= 1 )
216  {
217  G4cout << " =====> starting GEANT4E tracking for "
218  << currentTS->GetParticleType()
219  << " Forwards= " << g4edata->GetMode() << G4endl;
220  }
221  if( verbose >= 3 )
222  {
223  G4cout << " G4ErrorPropagator::Propagate initialTS ";
224  G4cout << *currentTS << G4endl;
225  }
226 #endif
227 
228  //----- If it is the first step, create a track
229  //
230  if( theStepN == 0 )
231  {
232  if( theG4Track != 0 ) { delete theG4Track; }
233  theG4Track = InitG4Track( *currentTS );
234  }
235  // set to 0 by the initialization in G4ErrorPropagatorManager
236  theStepN++;
237 
238  //----- Create a G4ErrorFreeTrajState
239  //
240  G4ErrorFreeTrajState* currentTS_FREE = InitFreeTrajState( currentTS );
241 
242  //----- Track the particle one step
243  //
244  ierr = MakeOneStep( currentTS_FREE );
245 
246  //----- Get the state on target
247  //
248  GetFinalTrajState( currentTS, currentTS_FREE, g4edata->GetTarget() );
249 
250  return ierr;
251 }
252 
253 
254 //---------------------------------------------------------------------------
256 {
257  if( verbose >= 5 ) { G4cout << "InitG4Track " << G4endl; }
258 
259  //----- Create Particle
260  //
261  const G4String partType = initialTS.GetParticleType();
263  G4ParticleDefinition* particle = particleTable->FindParticle(partType);
264  if( particle == 0)
265  {
266  std::ostringstream message;
267  message << "Particle type not defined: " << partType;
268  G4Exception("G4ErrorPropagator::InitG4Track()", "InvalidSetup",
269  FatalException, message);
270  }
271 
272  G4DynamicParticle* DP =
273  new G4DynamicParticle(particle,initialTS.GetMomentum());
274 
275  DP->SetPolarization(0.,0.,0.);
276 
277  // Set Charge
278  //
279  if( particle->GetPDGCharge() < 0 )
280  {
281  DP->SetCharge(-eplus);
282  }
283  else
284  {
285  DP->SetCharge(eplus);
286  }
287 
288  //----- Create Track
289  //
290  theG4Track = new G4Track(DP, 0., initialTS.GetPosition() );
291  theG4Track->SetParentID(0);
292 
293 #ifdef G4EVERBOSE
294  if(verbose >= 3)
295  {
296  G4cout << " G4ErrorPropagator new track of energy: "
297  << theG4Track->GetKineticEnergy() << G4endl;
298  }
299 #endif
300 
301  //---- Reproduce G4TrackingManager::ProcessOneTrack initialization
302  InvokePreUserTrackingAction( theG4Track );
303 
304  if( fpSteppingManager == 0 )
305  {
306  G4Exception("G4ErrorPropagator::InitG4Track()", "InvalidSetup",
307  FatalException, "G4SteppingManager not initialized yet!");
308  }
309  else
310  {
311  fpSteppingManager->SetInitialStep(theG4Track);
312  }
313 
314  // Give SteppingManger the maximum number of processes
315  //
316  fpSteppingManager->GetProcessNumber();
317 
318  // Give track the pointer to the Step
319  //
320  theG4Track->SetStep(fpSteppingManager->GetStep());
321 
322  // Inform beginning of tracking to physics processes
323  //
324  theG4Track->GetDefinition()->GetProcessManager()->StartTracking(theG4Track);
325 
326  initialTS.SetG4Track( theG4Track );
327 
328  return theG4Track;
329 }
330 
331 
332 //-----------------------------------------------------------------------
333 G4int G4ErrorPropagator::MakeSteps( G4ErrorFreeTrajState* currentTS_FREE )
334 {
335  G4int ierr = 0;
336 
337  //----- Track the particle Step-by-Step while it is alive
338  //
339  theStepLength = 0.;
340 
341  while( (theG4Track->GetTrackStatus() == fAlive) ||
342  (theG4Track->GetTrackStatus() == fStopButAlive) )
343  {
344  ierr = MakeOneStep( currentTS_FREE );
345  if( ierr != 0 ) { break; }
346 
347  //----- Check if last step for error propagation
348  //
349  if( CheckIfLastStep( theG4Track ) )
350  {
351  break;
352  }
353  } // Loop checking, 06.08.2015, G.Cosmo
354  return ierr;
355 }
356 
357 
358 //-----------------------------------------------------------------------
360 {
361  G4ErrorPropagatorData* g4edata =
363  G4int ierr = 0;
364 
365  //---------- Track one step
366 #ifdef G4EVERBOSE
367  if(verbose >= 2 )
368  {
369  G4cout << G4endl
370  << "@@@@@@@@@@@@@@@@@@@@@@@@@ NEW STEP " << G4endl;
371  }
372 #endif
373 
374  theG4Track->IncrementCurrentStepNumber();
375 
376  fpSteppingManager->Stepping();
377 
378  //---------- Check if Target has been reached (and then set G4ErrorState)
379 
380  // G4ErrorPropagationNavigator limits the step if target is closer than
381  // boundary (but the winner process is always "Transportation": then
382  // error propagator will stop the track
383 
384  if( theG4Track->GetStep()->GetPostStepPoint()
385  ->GetProcessDefinedStep()->GetProcessName() == "Transportation" )
386  {
387  if( g4edata->GetState()
389  { // target or step length reached
390 
391 #ifdef G4EVERBOSE
392  if(verbose >= 5 )
393  {
394  G4cout << " transportation determined by geant4e " << G4endl;
395  }
396 #endif
398  }
399  else if( g4edata->GetTarget()->GetType() == G4ErrorTarget_GeomVolume )
400  {
402  (G4ErrorGeomVolumeTarget*)(g4edata->GetTarget());
403  if( static_cast<G4ErrorGeomVolumeTarget*>( target )
404  ->TargetReached( theG4Track->GetStep() ) )
405  {
407  }
408  }
409  }
410  else if( theG4Track->GetStep()->GetPostStepPoint()->GetProcessDefinedStep()
411  ->GetProcessName() == "G4ErrorTrackLengthTarget" )
412  {
414  }
415 
416  //---------- Propagate error
417 
418 #ifdef G4EVERBOSE
419  if(verbose >= 2 )
420  {
421  G4cout << " propagating error " << *currentTS_FREE << G4endl;
422  }
423 #endif
424  const G4Track* cTrack = const_cast<G4Track*>(theG4Track);
425  ierr = currentTS_FREE->PropagateError( cTrack );
426 
427 #ifdef G4EVERBOSE
428  if(verbose >= 3 )
429  {
430  G4cout << " PropagateError returns " << ierr << G4endl;
431  }
432 #endif
433 
434  currentTS_FREE->Update( cTrack );
435 
436  theStepLength += theG4Track->GetStepLength();
437 
438  if(ierr != 0 )
439  {
440  std::ostringstream message;
441  message << "Error returned: " << ierr;
442  G4Exception("G4ErrorPropagator::MakeOneStep()",
443  "GEANT4e-Notification", JustWarning, message,
444  "Geant4 tracking will be stopped !");
445  }
446 
447  return ierr;
448 }
449 
450 
451 //-----------------------------------------------------------------------
454 {
455  G4ErrorFreeTrajState* currentTS_FREE = 0;
456 
457  //----- Transform the TrajState to Free coordinates if it is OnSurface
458  //
459  if( currentTS->GetTSType() == G4eTS_OS )
460  {
462  static_cast<G4ErrorSurfaceTrajState*>(currentTS);
463  currentTS_FREE = new G4ErrorFreeTrajState( *tssd );
464  }
465  else if( currentTS->GetTSType() == G4eTS_FREE )
466  {
467  currentTS_FREE = static_cast<G4ErrorFreeTrajState*>(currentTS);
468  }
469  else
470  {
471  std::ostringstream message;
472  message << "Wrong trajectory state: " << currentTS->GetTSType();
473  G4Exception("G4ErrorPropagator::InitFreeTrajState()", "InvalidState",
474  FatalException, message);
475  }
476  return currentTS_FREE;
477 }
478 
479 
480 //-----------------------------------------------------------------------
482  G4ErrorFreeTrajState* currentTS_FREE,
483  const G4ErrorTarget* target )
484 {
485  G4ErrorPropagatorData* g4edata =
487 
488 #ifdef G4EVERBOSE
489  if(verbose >= 1 )
490  {
491  G4cout << " G4ErrorPropagator::Propagate: final state "
492  << G4int(g4edata->GetState()) << " TSType "
493  << currentTS->GetTSType() << G4endl;
494  }
495 #endif
496 
497  if( (currentTS->GetTSType() == G4eTS_FREE) ||
498  (g4edata->GetState() != G4ErrorState_StoppedAtTarget) )
499  {
500  currentTS = currentTS_FREE;
501  }
502  else if( currentTS->GetTSType() == G4eTS_OS )
503  {
504  if( target->GetType() == G4ErrorTarget_TrkL )
505  {
506  G4Exception("G4ErrorPropagator:GetFinalTrajState()",
507  "InvalidSetup", FatalException,
508  "Using a G4ErrorSurfaceTrajState with wrong target");
509  }
510  const G4ErrorTanPlaneTarget* targetWTP =
511  static_cast<const G4ErrorTanPlaneTarget*>(target);
512  *currentTS = G4ErrorSurfaceTrajState(
513  *(static_cast<G4ErrorFreeTrajState*>(currentTS_FREE)),
514  targetWTP->GetTangentPlane( currentTS_FREE->GetPosition() ) );
515 #ifdef G4EVERBOSE
516  if(verbose >= 1 )
517  {
518  G4cout << currentTS << " returning tssd " << *currentTS << G4endl;
519  }
520 #endif
521  delete currentTS_FREE;
522  }
523 }
524 
525 
526 //-------------------------------------------------------------------------
528 {
529  G4bool exception = 0;
530  G4bool lastG4eStep = false;
531  G4ErrorPropagatorData* g4edata =
533 
534 #ifdef G4EVERBOSE
535  if( verbose >= 4 )
536  {
537  G4cout << " G4ErrorPropagator::CheckIfLastStep G4ErrorState= "
538  << G4int(g4edata->GetState()) << G4endl;
539  }
540 #endif
541 
542  //----- Check if this is the last step: track has reached the target
543  // or the end of world
544  //
546  {
547  lastG4eStep = true;
548 #ifdef G4EVERBOSE
549  if(verbose >= 5 )
550  {
551  G4cout << " G4ErrorPropagator::CheckIfLastStep " << lastG4eStep
552  << " " << G4int(g4edata->GetState()) << G4endl;
553  }
554 #endif
555  }
556  else if( aTrack->GetNextVolume() == 0 )
557  {
558  //----- If particle is out of world, without finding the G4ErrorTarget,
559  // give a n error/warning
560  //
561  lastG4eStep = true;
562  if( exception )
563  {
564  std::ostringstream message;
565  message << "Track extrapolated until end of World" << G4endl
566  << "without finding the defined target!";
567  G4Exception("G4ErrorPropagator::CheckIfLastStep()",
568  "InvalidSetup", FatalException, message);
569  }
570  else
571  {
572  if( verbose >= 1 )
573  {
574  std::ostringstream message;
575  message << "Track extrapolated until end of World" << G4endl
576  << "without finding the defined target.";
577  G4Exception("G4ErrorPropagator::CheckIfLastStep()",
578  "GEANT4e-Notification", JustWarning, message);
579  }
580  }
581  } //----- not last step from G4e, but track is stopped (energy exhausted)
582  else if( aTrack->GetTrackStatus() == fStopAndKill )
583  {
584  if( exception )
585  {
586  std::ostringstream message;
587  message << "Track extrapolated until energy is exhausted" << G4endl
588  << "without finding the defined target!";
589  G4Exception("G4ErrorPropagator::CheckIfLastStep()",
590  "InvalidSetup", FatalException, message);
591  }
592  else
593  {
594  if( verbose >= 1 )
595  {
596  std::ostringstream message;
597  message << "Track extrapolated until energy is exhausted" << G4endl
598  << "without finding the defined target.";
599  G4Exception("G4ErrorPropagator::CheckIfLastStep()",
600  "GEANT4e-Notification", JustWarning, message);
601  }
602  lastG4eStep = 1;
603  }
604  }
605 
606 #ifdef G4EVERBOSE
607  if( verbose >= 5 )
608  {
609  G4cout << " return CheckIfLastStep " << lastG4eStep << G4endl;
610  }
611 #endif
612 
613  return lastG4eStep;
614 }
615 
616 
617 //---------------------------------------------------------------------------
619 {
620  const G4UserTrackingAction* fpUserTrackingAction =
622  if( fpUserTrackingAction != 0 )
623  {
624  const_cast<G4UserTrackingAction*>(fpUserTrackingAction)
625  ->PreUserTrackingAction((fpTrack) );
626  }
627 }
628 
629 
630 //---------------------------------------------------------------------------
632 {
633  const G4UserTrackingAction* fpUserTrackingAction =
635  if( fpUserTrackingAction != 0 )
636  {
637  const_cast<G4UserTrackingAction*>(fpUserTrackingAction)
638  ->PostUserTrackingAction((fpTrack) );
639  }
640 }
G4ParticleDefinition * GetDefinition() const
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
const XML_Char * target
Definition: expat.h:268
G4SteppingManager * GetSteppingManager() const
G4int MakeOneStep(G4ErrorFreeTrajState *currentTS_FREE)
G4int Propagate(G4ErrorTrajState *currentTS, const G4ErrorTarget *target, G4ErrorMode mode=G4ErrorMode_PropForwards)
void InvokePostUserTrackingAction(G4Track *fpTrack)
void SetG4Track(G4Track *trk)
G4int PropagateOneStep(G4ErrorTrajState *currentTS)
G4TrackStatus GetTrackStatus() const
G4Point3D GetPosition() const
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
G4ErrorTargetType GetType() const
const G4Step * GetStep() const
int G4int
Definition: G4Types.hh:78
G4VPhysicalVolume * GetNextVolume() const
virtual void Dump(const G4String &msg) const =0
void SetState(G4ErrorState sta)
static G4StateManager * GetStateManager()
G4double GetKineticEnergy() const
G4TrackingManager * GetTrackingManager() const
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
void SetMode(G4ErrorMode mode)
G4Vector3D GetMomentum() const
static constexpr double eplus
Definition: G4SIunits.hh:199
void GetFinalTrajState(G4ErrorTrajState *currentTS, G4ErrorFreeTrajState *currentTS_FREE, const G4ErrorTarget *target)
G4ErrorMode GetMode() const
void InvokePreUserTrackingAction(G4Track *fpTrack)
G4ApplicationState GetCurrentState() const
void SetPolarization(G4double polX, G4double polY, G4double polZ)
G4StepStatus Stepping()
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
G4Track * InitG4Track(G4ErrorTrajState &initialTS)
void SetInitialStep(G4Track *valueTrack)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
virtual G4Plane3D GetTangentPlane(const G4ThreeVector &point) const =0
void SetTarget(const G4ErrorTarget *target)
const G4VProcess * GetProcessDefinedStep() const
void IncrementCurrentStepNumber()
static G4ParticleTable * GetParticleTable()
void StartTracking(G4Track *aTrack=0)
const G4String & GetParticleType() const
G4ProcessManager * GetProcessManager() const
G4bool CheckIfLastStep(G4Track *aTrack)
G4UserTrackingAction * GetUserTrackingAction()
void SetParentID(const G4int aValue)
G4StepPoint * GetPostStepPoint() const
virtual G4eTSType GetTSType() const
static G4EventManager * GetEventManager()
const G4ErrorTarget * GetTarget(G4bool mustExist=0) const
virtual G4int Update(const G4Track *aTrack)
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
static G4ErrorPropagatorData * GetErrorPropagatorData()
virtual G4int PropagateError(const G4Track *aTrack)
void SetCharge(G4double charge)
G4double GetPDGCharge() const
G4ErrorFreeTrajState * InitFreeTrajState(G4ErrorTrajState *currentTS)
void SetStep(const G4Step *aValue)
G4ErrorState GetState() const
G4double GetStepLength() const
G4Step * GetStep() const