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