Geant4  10.01.p02
G4VIntersectionLocator.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: G4VIntersectionLocator.cc 90836 2015-06-10 09:31:06Z gcosmo $
27 //
28 // Class G4VIntersectionLocator implementation
29 //
30 // 27.10.08 - John Apostolakis, Tatiana Nikitina.
31 // ---------------------------------------------------------------------------
32 
33 #include <iomanip>
34 #include <sstream>
35 
36 #include "globals.hh"
37 #include "G4ios.hh"
38 #include "G4SystemOfUnits.hh"
40 #include "G4GeometryTolerance.hh"
41 
43 //
44 // Constructor
45 //
47  fUseNormalCorrection(false),
48  fiNavigator( theNavigator ),
49  fiChordFinder( 0 ), // Not set - overridden at each step
50  fiEpsilonStep( -1.0 ), // Out of range - overridden at each step
51  fiDeltaIntersection( -1.0 ), // Out of range - overridden at each step
52  fiUseSafety(false), // Default - overridden at each step
53  fpTouchable(0)
54 {
56  fVerboseLevel = 0;
58 }
59 
61 //
62 // Destructor.
63 //
65 {
66  delete fHelpingNavigator;
67  delete fpTouchable;
68 }
69 
71 //
72 // Dump status of propagator to cout (old method)
73 //
74 void
76  const G4FieldTrack& CurrentFT,
77  G4double requestStep,
78  G4double safety,
79  G4int stepNo)
80 {
81  std::ostringstream os;
82  printStatus( StartFT, CurrentFT, requestStep, safety, stepNo, os);
83  G4cout << os.str();
84 }
85 
87 //
88 // Dumps status of propagator.
89 //
90 void
92  const G4FieldTrack& CurrentFT,
93  G4double requestStep,
94  G4double safety,
95  G4int stepNo,
96  std::ostringstream& os )
97 {
98  const G4int verboseLevel= fVerboseLevel;
99  const G4ThreeVector StartPosition = StartFT.GetPosition();
100  const G4ThreeVector StartUnitVelocity = StartFT.GetMomentumDir();
101  const G4ThreeVector CurrentPosition = CurrentFT.GetPosition();
102  const G4ThreeVector CurrentUnitVelocity = CurrentFT.GetMomentumDir();
103 
104  G4double step_len = CurrentFT.GetCurveLength() - StartFT.GetCurveLength();
105  G4int oldprc; // cout/cerr precision settings
106 
107  if( ((stepNo == 0) && (verboseLevel <3)) || (verboseLevel >= 3) )
108  {
109  oldprc = os.precision(4);
110  os << std::setw( 6) << " "
111  << std::setw( 25) << " Current Position and Direction" << " "
112  << G4endl;
113  os << std::setw( 5) << "Step#"
114  << std::setw(10) << " s " << " "
115  << std::setw(10) << "X(mm)" << " "
116  << std::setw(10) << "Y(mm)" << " "
117  << std::setw(10) << "Z(mm)" << " "
118  << std::setw( 7) << " N_x " << " "
119  << std::setw( 7) << " N_y " << " "
120  << std::setw( 7) << " N_z " << " " ;
121  os << std::setw( 7) << " Delta|N|" << " "
122  << std::setw( 9) << "StepLen" << " "
123  << std::setw(12) << "StartSafety" << " "
124  << std::setw( 9) << "PhsStep" << " ";
125  os << G4endl;
126  os.precision(oldprc);
127  }
128  if((stepNo == 0) && (verboseLevel <=3))
129  {
130  // Recurse to print the start values
131  //
132  printStatus( StartFT, StartFT, -1.0, safety, -1);
133  }
134  if( verboseLevel <= 3 )
135  {
136  if( stepNo >= 0)
137  {
138  os << std::setw( 4) << stepNo << " ";
139  }
140  else
141  {
142  os << std::setw( 5) << "Start" ;
143  }
144  oldprc = os.precision(8);
145  os << std::setw(10) << CurrentFT.GetCurveLength() << " ";
146  os << std::setw(10) << CurrentPosition.x() << " "
147  << std::setw(10) << CurrentPosition.y() << " "
148  << std::setw(10) << CurrentPosition.z() << " ";
149  os.precision(4);
150  os << std::setw( 7) << CurrentUnitVelocity.x() << " "
151  << std::setw( 7) << CurrentUnitVelocity.y() << " "
152  << std::setw( 7) << CurrentUnitVelocity.z() << " ";
153  os.precision(3);
154  os << std::setw( 7)
155  << CurrentFT.GetMomentum().mag()- StartFT.GetMomentum().mag()
156  << " ";
157  os << std::setw( 9) << step_len << " ";
158  os << std::setw(12) << safety << " ";
159  if( requestStep != -1.0 )
160  {
161  os << std::setw( 9) << requestStep << " ";
162  }
163  else
164  {
165  os << std::setw( 9) << "Init/NotKnown" << " ";
166  }
167  os << G4endl;
168  os.precision(oldprc);
169  }
170  else // if( verboseLevel > 3 )
171  {
172  // Multi-line output
173 
174  os << "Step taken was " << step_len
175  << " out of PhysicalStep= " << requestStep << G4endl;
176  os << "Final safety is: " << safety << G4endl;
177  os << "Chord length = " << (CurrentPosition-StartPosition).mag()
178  << G4endl;
179  os << G4endl;
180  }
181 }
182 
184 //
185 // ReEstimateEndPoint.
186 //
188 ReEstimateEndpoint( const G4FieldTrack& CurrentStateA,
189  const G4FieldTrack& EstimatedEndStateB,
190  G4double linearDistSq,
191  G4double curveDist )
192 {
193  G4FieldTrack newEndPoint( CurrentStateA );
195 
196  G4FieldTrack retEndPoint( CurrentStateA );
197  G4bool goodAdvance;
198  G4int itrial=0;
199  const G4int no_trials= 20;
200 
201  G4double endCurveLen= EstimatedEndStateB.GetCurveLength();
202  do
203  {
204  G4double currentCurveLen= newEndPoint.GetCurveLength();
205  G4double advanceLength= endCurveLen - currentCurveLen ;
206  if (std::abs(advanceLength)<kCarTolerance)
207  {
208  goodAdvance=true;
209  }
210  else{
211  goodAdvance=
212  integrDriver->AccurateAdvance(newEndPoint, advanceLength,
214  }
215  }
216  while( !goodAdvance && (++itrial < no_trials) );
217 
218  if( goodAdvance )
219  {
220  retEndPoint= newEndPoint;
221  }
222  else
223  {
224  retEndPoint= EstimatedEndStateB; // Could not improve without major work !!
225  }
226 
227  // All the work is done
228  // below are some diagnostics only -- before the return!
229  //
230  static const G4String MethodName("G4VIntersectionLocator::ReEstimateEndpoint");
231 
232 #ifdef G4VERBOSE
233  G4int latest_good_trials=0;
234  if( itrial > 1)
235  {
236  if( fVerboseLevel > 0 )
237  {
238  G4cout << MethodName << " called - goodAdv= " << goodAdvance
239  << " trials = " << itrial
240  << " previous good= " << latest_good_trials
241  << G4endl;
242  }
243  latest_good_trials=0;
244  }
245  else
246  {
247  latest_good_trials++;
248  }
249 #endif
250 
251 #ifdef G4DEBUG_FIELD
252  G4double lengthDone = newEndPoint.GetCurveLength()
253  - CurrentStateA.GetCurveLength();
254  if( !goodAdvance )
255  {
256  if( fVerboseLevel >= 3 )
257  {
258  G4cout << MethodName << "> AccurateAdvance failed " ;
259  G4cout << " in " << itrial << " integration trials/steps. " << G4endl;
260  G4cout << " It went only " << lengthDone << " instead of " << curveDist
261  << " -- a difference of " << curveDist - lengthDone << G4endl;
262  G4cout << " ReEstimateEndpoint> Reset endPoint to original value!"
263  << G4endl;
264  }
265  }
266 
267  static G4int noInaccuracyWarnings = 0;
268  G4int maxNoWarnings = 10;
269  if ( (noInaccuracyWarnings < maxNoWarnings )
270  || (fVerboseLevel > 1) )
271  {
272  G4ThreeVector move = newEndPoint.GetPosition()
273  - EstimatedEndStateB.GetPosition();
274  std::ostringstream message;
275  message.precision(12);
276  message << " Integration inaccuracy requires"
277  << " an adjustment in the step's endpoint." << G4endl
278  << " Two mid-points are further apart than their"
279  << " curve length difference" << G4endl
280  << " Dist = " << std::sqrt(linearDistSq)
281  << " curve length = " << curveDist << G4endl;
282  message << " Correction applied is " << move.mag() << G4endl
283  << " Old Estimated B position= "
284  << EstimatedEndStateB.GetPosition() << G4endl
285  << " Recalculated Position= "
286  << newEndPoint.GetPosition() << G4endl
287  << " Change ( new - old ) = " << move;
288  G4Exception("G4VIntersectionLocator::ReEstimateEndpoint()",
289  "GeomNav1002", JustWarning, message);
290  }
291 #else
292  // Statistics on the RMS value of the corrections
293 
294  static G4ThreadLocal G4int noCorrections=0;
295  static G4ThreadLocal G4double sumCorrectionsSq = 0;
296  noCorrections++;
297  if( goodAdvance )
298  {
299  sumCorrectionsSq += (EstimatedEndStateB.GetPosition() -
300  newEndPoint.GetPosition()).mag2();
301  }
302  linearDistSq -= curveDist; // To use linearDistSq ... !
303 #endif
304 
305  return retEndPoint;
306 }
307 
309 //
310 // Method for finding SurfaceNormal of Intersecting Solid
311 //
313 GetLocalSurfaceNormal(const G4ThreeVector& CurrentE_Point, G4bool& validNormal)
314 {
315  G4ThreeVector Normal(G4ThreeVector(0.0,0.0,0.0));
316  G4VPhysicalVolume* located;
317 
318  validNormal = false;
319  fHelpingNavigator->SetWorldVolume(GetNavigatorFor()->GetWorldVolume());
320  located = fHelpingNavigator->LocateGlobalPointAndSetup( CurrentE_Point );
321 
322  delete fpTouchable;
324 
325  // To check if we can use GetGlobalExitNormal()
326  //
327  G4ThreeVector localPosition = fpTouchable->GetHistory()
328  ->GetTopTransform().TransformPoint(CurrentE_Point);
329 
330  // Issue: in the case of coincident surfaces, this version does not recognise
331  // which side you are located onto (can return vector with wrong sign.)
332  // TO-DO: use direction (of chord) to identify volume we will be "entering"
333 
334  if( located != 0)
335  {
336  G4LogicalVolume* pLogical= located->GetLogicalVolume();
337  G4VSolid* pSolid;
338 
339  if( (pLogical != 0) && ( (pSolid=pLogical->GetSolid()) !=0 ) )
340  {
341  // G4bool goodPoint, nearbyPoint;
342  // G4int numGoodPoints, numNearbyPoints; // --> use for stats
343  if ( ( pSolid->Inside(localPosition)==kSurface )
344  || ( pSolid->DistanceToOut(localPosition) < 1000.0 * kCarTolerance )
345  )
346  {
347  Normal = pSolid->SurfaceNormal(localPosition);
348  validNormal = true;
349 
350 #ifdef G4DEBUG_FIELD
351  if( std::fabs(Normal.mag2() - 1.0 ) > CLHEP::perThousand)
352  {
353  G4cerr << "PROBLEM in G4VIntersectionLocator::GetLocalSurfaceNormal."
354  << G4endl;
355  G4cerr << " Normal is not unit - mag=" << Normal.mag() << G4endl;
356  G4cerr << " at trial local point " << CurrentE_Point << G4endl;
357  G4cerr << " Solid is " << *pSolid << G4endl;
358  }
359 #endif
360  }
361  }
362  }
363 
364  return Normal;
365 }
366 
368 //
369 // Adjustment of Found Intersection
370 //
373  const G4ThreeVector& CurrentE_Point,
374  const G4ThreeVector& CurrentF_Point,
375  const G4ThreeVector& MomentumDir,
376  const G4bool IntersectAF,
377  G4ThreeVector& IntersectionPoint, // I/O
378  G4double& NewSafety, // I/O
379  G4double& fPreviousSafety, // I/O
381 {
382  G4double dist,lambda;
383  G4ThreeVector Normal, NewPoint, Point_G;
384  G4bool goodAdjust=false, Intersects_FP=false, validNormal=false;
385 
386  // Get SurfaceNormal of Intersecting Solid
387  //
388  Normal = GetGlobalSurfaceNormal(CurrentE_Point,validNormal);
389  if(!validNormal) { return false; }
390 
391  // Intersection between Line and Plane
392  //
393  G4double n_d_m = Normal.dot(MomentumDir);
394  if ( std::abs(n_d_m)>kCarTolerance )
395  {
396 #ifdef G4VERBOSE
397  if ( fVerboseLevel>1 )
398  {
399  G4cerr << "WARNING - "
400  << "G4VIntersectionLocator::AdjustementOfFoundIntersection()"
401  << G4endl
402  << " No intersection. Parallels lines!" << G4endl;
403  }
404 #endif
405  lambda =- Normal.dot(CurrentF_Point-CurrentE_Point)/n_d_m;
406 
407  // New candidate for Intersection
408  //
409  NewPoint = CurrentF_Point+lambda*MomentumDir;
410 
411  // Distance from CurrentF to Calculated Intersection
412  //
413  dist = std::abs(lambda);
414 
415  if ( dist<kCarTolerance*0.001 ) { return false; }
416 
417  // Calculation of new intersection point on the path.
418  //
419  if ( IntersectAF ) // First part intersects
420  {
421  G4double stepLengthFP;
422  G4ThreeVector Point_P = CurrentA_Point;
424  Intersects_FP = IntersectChord( Point_P, NewPoint, NewSafety,
425  fPreviousSafety, fPreviousSftOrigin,
426  stepLengthFP, Point_G );
427 
428  }
429  else // Second part intersects
430  {
431  G4double stepLengthFP;
433  Intersects_FP = IntersectChord( CurrentF_Point, NewPoint, NewSafety,
434  fPreviousSafety, fPreviousSftOrigin,
435  stepLengthFP, Point_G );
436  }
437  if ( Intersects_FP )
438  {
439  goodAdjust = true;
440  IntersectionPoint = Point_G;
441  }
442  }
443 
444  return goodAdjust;
445 }
446 
449  G4bool& validNormal) // const
450 {
451  G4ThreeVector NormalAtEntry; // ( -10. , -10., -10. );
452 
453  G4ThreeVector NormalAtEntryLast, NormalAtEntryGlobal, diffNormals;
454  G4bool validNormalLast;
455 
456  // Relies on a call to Navigator::ComputeStep in IntersectChord before
457  // this call
458  //
459  NormalAtEntryLast = GetLastSurfaceNormal( CurrentInt_Point, validNormalLast );
460  // May return valid=false in cases, including
461  // - if the candidate volume was not found (eg exiting world), or
462  // - a replica was involved -- determined the step size.
463  // (This list is not complete.)
464 
465 #ifdef G4DEBUG_FIELD
466  if ( validNormalLast
467  && ( std::fabs(NormalAtEntryLast.mag2() - 1.0) > perThousand ) )
468  {
469  std::ostringstream message;
470  message << "PROBLEM: Normal is not unit - magnitude = "
471  << NormalAtEntryLast.mag() << G4endl;
472  message << " at trial intersection point " << CurrentInt_Point << G4endl;
473  message << " Obtained from Get *Last* Surface Normal.";
474  G4Exception("G4VIntersectionLocator::GetSurfaceNormal()",
475  "GeomNav1002", JustWarning, message);
476  }
477 #endif
478 
479  if( validNormalLast )
480  {
481  NormalAtEntry=NormalAtEntryLast;
482  }
483  validNormal = validNormalLast;
484 
485  return NormalAtEntry;
486 }
487 
489 GetGlobalSurfaceNormal(const G4ThreeVector& CurrentE_Point,
490  G4bool& validNormal)
491 {
492  G4ThreeVector localNormal=
493  GetLocalSurfaceNormal( CurrentE_Point, validNormal );
494  G4AffineTransform localToGlobal= // Must use the same Navigator !!
496  G4ThreeVector globalNormal =
497  localToGlobal.TransformAxis( localNormal );
498 
499 #ifdef G4DEBUG_FIELD
500  if( validNormal && ( std::fabs(globalNormal.mag2() - 1.0) > perThousand ) )
501  {
502  std::ostringstream message;
503  message << "**************************************************************"
504  << G4endl;
505  message << " Bad Normal in G4VIntersectionLocator::GetGlobalSurfaceNormal "
506  << G4endl;
507  message << " * Constituents: " << G4endl;
508  message << " Local Normal= " << localNormal << G4endl;
509  message << " Transform: " << G4endl
510  << " Net Translation= " << localToGlobal.NetTranslation()
511  << G4endl
512  << " Net Rotation = " << localToGlobal.NetRotation()
513  << G4endl;
514  message << " * Result: " << G4endl;
515  message << " Global Normal= " << localNormal << G4endl;
516  message << "**************************************************************"
517  << G4endl;
518  G4Exception("G4VIntersectionLocator::GetGlobalSurfaceNormal()",
519  "GeomNav1002", JustWarning, message);
520  }
521 #endif
522 
523  return globalNormal;
524 }
525 
528  G4bool& normalIsValid) const
529 {
530  G4ThreeVector normalVec;
531  G4bool validNorm;
532  normalVec = fiNavigator->GetGlobalExitNormal( intersectPoint, &validNorm );
533  normalIsValid= validNorm;
534 
535  return normalVec;
536 }
537 
539  const G4ThreeVector& ChordAB_v,
540  const G4ThreeVector& ChordEF_v,
541  const G4ThreeVector& NewMomentumDir,
542  const G4ThreeVector& NormalAtEntry,
543  G4bool validNormal )
544 {
545  G4double ABchord_length = ChordAB_v.mag();
546  G4double MomDir_dot_Norm = NewMomentumDir.dot( NormalAtEntry ) ;
547  G4double MomDir_dot_ABchord;
548  MomDir_dot_ABchord= (1.0 / ABchord_length) * NewMomentumDir.dot( ChordAB_v );
549 
550  std::ostringstream outStream;
551  outStream // G4cout
552  << std::setw(6) << " Step# "
553  << std::setw(17) << " |ChordEF|(mag)" << " "
554  << std::setw(18) << " uMomentum.Normal" << " "
555  << std::setw(18) << " uMomentum.ABdir " << " "
556  << std::setw(16) << " AB-dist " << " "
557  << " Chord Vector (EF) "
558  << G4endl;
559  outStream.precision(7);
560  outStream // G4cout
561  << " " << std::setw(5) << step_no
562  << " " << std::setw(18) << ChordEF_v.mag()
563  << " " << std::setw(18) << MomDir_dot_Norm
564  << " " << std::setw(18) << MomDir_dot_ABchord
565  << " " << std::setw(12) << ABchord_length
566  << " " << ChordEF_v
567  << G4endl;
568  outStream // G4cout
569  << " MomentumDir= " << " " << NewMomentumDir
570  << " Normal at Entry E= " << NormalAtEntry
571  << " AB chord = " << ChordAB_v
572  << G4endl;
573  G4cout << outStream.str(); // ostr_verbose;
574 
575  if( ( std::fabs(NormalAtEntry.mag2() - 1.0) > perThousand ) )
576  {
577  G4cerr << " PROBLEM in G4VIntersectionLocator::ReportTrialStep " << G4endl
578  << " Normal is not unit - mag=" << NormalAtEntry.mag()
579  << " ValidNormalAtE = " << validNormal
580  << G4endl;
581  }
582  return;
583 }
584 
585 // Locate point using navigator: updates state of Navigator
586 // By default, it assumes that the point is inside the current volume,
587 // and returns true.
588 // In check mode, checks whether the point is *inside* the volume.
589 // If it is inside, it returns true
590 // If not, issues a warning and returns false.
591 //
592 G4bool
595 {
596  G4bool good= true;
598  const G4String MethodName("G4VIntersectionLocator::LocateGlobalPointWithinVolumeAndCheck()");
599 
600  if( fCheckMode )
601  {
602  G4bool navCheck= nav->IsCheckModeActive(); // Recover original value
603  nav->CheckMode(true);
604 
605  // Identify the current volume
606 
608  G4VPhysicalVolume* motherPhys= startTH->GetVolume();
609  G4VSolid* motherSolid= startTH->GetSolid();
610  G4AffineTransform transform = nav->GetGlobalToLocalTransform();
611  // GetLocalToGlobalTransform();
612  G4int motherCopyNo= motherPhys->GetCopyNo();
613 
614  // Let's check that the point is inside the current solid
615  G4ThreeVector localPosition = transform.TransformPoint(position);
616  EInside inMother= motherSolid->Inside( localPosition );
617  if( inMother != kInside )
618  {
619  G4cerr << " ERROR in " << MethodName << " Position located "
620  << ( inMother == kSurface ? " on Surface " : " outside " )
621  << "expected volume" << G4endl;
622  G4double safetyFromOut= motherSolid->DistanceToIn(localPosition);
623  G4cerr << " Safety (from Outside) = " << safetyFromOut << G4endl;
624  }
625 
626  // 1. Simple next step - quick relocation and check result.
627  // nav->LocateGlobalPointWithinVolume( position );
628 
629  // 2. Full relocation - to cross-check answer !
630  G4VPhysicalVolume* nextPhysical= nav->LocateGlobalPointAndSetup(position);
631  if( (nextPhysical != motherPhys)
632  || (nextPhysical->GetCopyNo() != motherCopyNo )
633  )
634  {
635  G4cerr << " ERROR in " << MethodName
636  << " Position located outside expected volume" << G4endl;
637  }
638  nav->CheckMode(navCheck); // Recover original value
639  }
640  else
641  {
642  nav->LocateGlobalPointWithinVolume( position );
643  }
644  return good;
645 }
646 
647 void
650  const G4String& CodeLocationInfo,
651  G4int CheckMode )
652 {
653  // Save value of Check mode first
654  G4bool oldCheck= GetCheckMode();
655 
657  if( !ok )
658  {
659  G4cerr << " ERROR occured in Intersection Locator" << G4endl;
660  G4cerr << " Code Location info: " << CodeLocationInfo << G4endl;
661  if( CheckMode > 1 ) {
662  // Additional information
663 
664  }
665  }
666 
667  SetCheckMode( oldCheck );
668 }
G4FieldTrack ReEstimateEndpoint(const G4FieldTrack &CurrentStateA, const G4FieldTrack &EstimtdEndStateB, G4double linearDistSq, G4double curveDist)
G4double GetCurveLength() const
CLHEP::Hep3Vector G4ThreeVector
G4ThreeVector GetGlobalSurfaceNormal(const G4ThreeVector &CurrentE_Point, G4bool &validNormal)
G4bool IsCheckModeActive() const
G4ThreeVector GetSurfaceNormal(const G4ThreeVector &CurrentInt_Point, G4bool &validNormal)
virtual G4TouchableHistoryHandle CreateTouchableHistoryHandle() const
G4bool IntersectChord(const G4ThreeVector &StartPointA, const G4ThreeVector &EndPointB, G4double &NewSafety, G4double &PreviousSafety, G4ThreeVector &PreviousSftOrigin, G4double &LinearStepLength, G4ThreeVector &IntersectionPoint, G4bool *calledNavigator=0)
G4double GetSurfaceTolerance() const
const G4ThreeVector & GetMomentumDir() const
G4ThreeVector GetLastSurfaceNormal(const G4ThreeVector &intersectPoint, G4bool &validNormal) const
G4ThreeVector NetTranslation() const
static const double perThousand
Definition: G4SIunits.hh:297
const G4AffineTransform GetLocalToGlobalTransform() const
#define G4ThreadLocal
Definition: tls.hh:89
int G4int
Definition: G4Types.hh:78
G4double GetEpsilonStepFor()
#define position
Definition: xmlparse.cc:605
G4ThreeVector GetPosition() const
G4GLOB_DLL std::ostream G4cout
virtual EInside Inside(const G4ThreeVector &p) const =0
bool G4bool
Definition: G4Types.hh:79
virtual G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const =0
void ReportTrialStep(G4int step_no, const G4ThreeVector &ChordAB_v, const G4ThreeVector &ChordEF_v, const G4ThreeVector &NewMomentumDir, const G4ThreeVector &NormalAtEntry, G4bool validNormal)
void LocateGlobalPointWithinVolumeCheckAndReport(const G4ThreeVector &pos, const G4String &CodeLocationInfo, G4int CheckMode)
virtual G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const =0
G4Navigator * GetNavigatorFor()
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4bool LocateGlobalPointWithinVolumeAndCheck(const G4ThreeVector &pos)
G4RotationMatrix NetRotation() const
G4ThreeVector TransformPoint(const G4ThreeVector &vec) const
G4TouchableHistory * CreateTouchableHistory() const
G4LogicalVolume * GetLogicalVolume() const
#define fPreviousSftOrigin
virtual G4ThreeVector GetGlobalExitNormal(const G4ThreeVector &point, G4bool *valid)
G4VIntersectionLocator(G4Navigator *theNavigator)
EInside
Definition: geomdefs.hh:58
#define fPreviousSafety
G4ThreeVector TransformAxis(const G4ThreeVector &axis) const
virtual G4int GetCopyNo() const =0
void printStatus(const G4FieldTrack &startFT, const G4FieldTrack &currentFT, G4double requestStep, G4double safety, G4int step)
void SetWorldVolume(G4VPhysicalVolume *pWorld)
const G4AffineTransform & GetGlobalToLocalTransform() const
virtual G4VPhysicalVolume * LocateGlobalPointAndSetup(const G4ThreeVector &point, const G4ThreeVector *direction=0, const G4bool pRelativeSearch=true, const G4bool ignoreDirection=true)
Definition: G4Navigator.cc:122
G4ThreeVector GetLocalSurfaceNormal(const G4ThreeVector &CurrentE_Point, G4bool &validNormal)
#define G4endl
Definition: G4ios.hh:61
const G4AffineTransform & GetTopTransform() const
void CheckMode(G4bool mode)
const G4NavigationHistory * GetHistory() const
G4ThreeVector GetMomentum() const
double G4double
Definition: G4Types.hh:76
G4TouchableHistory * fpTouchable
virtual G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=0, G4ThreeVector *n=0) const =0
G4MagInt_Driver * GetIntegrationDriver()
G4ChordFinder * GetChordFinderFor()
G4bool AccurateAdvance(G4FieldTrack &y_current, G4double hstep, G4double eps, G4double hinitial=0.0)
static G4GeometryTolerance * GetInstance()
virtual void LocateGlobalPointWithinVolume(const G4ThreeVector &position)
Definition: G4Navigator.cc:579
G4VSolid * GetSolid() const
G4bool AdjustmentOfFoundIntersection(const G4ThreeVector &A, const G4ThreeVector &CurrentE_Point, const G4ThreeVector &CurrentF_Point, const G4ThreeVector &MomentumDir, const G4bool IntersectAF, G4ThreeVector &IntersectionPoint, G4double &NewSafety, G4double &fPrevSafety, G4ThreeVector &fPrevSftOrigin)
G4GLOB_DLL std::ostream G4cerr