source: trunk/source/geometry/navigation/src/G4PathFinder.cc@ 834

Last change on this file since 834 was 831, checked in by garnier, 17 years ago

import all except CVS
File size: 46.8 KB
Line 
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// $Id: G4PathFinder.cc,v 1.58.2.1 2008/04/29 15:37:16 gcosmo Exp $
28// GEANT4 tag $ Name: $
29//
30// class G4PathFinder Implementation
31//
32// Original author: John Apostolakis, April 2006
33//
34// --------------------------------------------------------------------
35
36#include "G4PathFinder.hh"
37
38#include <iomanip>
39
40#include "G4GeometryTolerance.hh"
41#include "G4Navigator.hh"
42#include "G4PropagatorInField.hh"
43#include "G4TransportationManager.hh"
44#include "G4MultiNavigator.hh"
45#include "G4SafetyHelper.hh"
46
47// Initialise the static instance of the singleton
48//
49G4PathFinder* G4PathFinder::fpPathFinder=0;
50
51// ----------------------------------------------------------------------------
52// GetInstance()
53//
54// Retrieve the static instance of the singleton
55//
56G4PathFinder* G4PathFinder::GetInstance()
57{
58 static G4PathFinder theInstance;
59 if( ! fpPathFinder )
60 {
61 fpPathFinder = &theInstance;
62 }
63 return fpPathFinder;
64}
65
66// ----------------------------------------------------------------------------
67// Constructor
68//
69G4PathFinder::G4PathFinder()
70 : fEndState( G4ThreeVector(), G4ThreeVector(), 0., 0., 0., 0., 0.),
71 fRelocatedPoint(true),
72 fLastStepNo(-1),
73 fVerboseLevel(0)
74{
75 fpMultiNavigator= new G4MultiNavigator();
76
77 fpTransportManager= G4TransportationManager::GetTransportationManager();
78 fpFieldPropagator = fpTransportManager->GetPropagatorInField();
79
80 kCarTolerance = G4GeometryTolerance::GetInstance()->GetSurfaceTolerance();
81
82 fNoActiveNavigators= 0;
83 G4ThreeVector Big3Vector( kInfinity, kInfinity, kInfinity );
84 fLastLocatedPosition= Big3Vector;
85 fSafetyLocation= Big3Vector;
86 fPreSafetyLocation= Big3Vector;
87 fPreStepLocation= Big3Vector;
88
89 fPreSafetyMinValue= -1.0;
90 fMinSafety_PreStepPt= -1.0;
91 fMinSafety_atSafLocation= -1.0;
92 fMinStep= -1.0;
93 fNewTrack= false;
94 fNoGeometriesLimiting= 0;
95
96 for( register int num=0; num<= fMaxNav; ++num )
97 {
98 fpNavigator[num] = 0;
99 fLimitTruth[num] = false;
100 fLimitedStep[num] = kUndefLimited;
101 fCurrentStepSize[num] = -1.0;
102 fLocatedVolume[num] = 0;
103 fPreSafetyValues[num]= -1.0;
104 fCurrentPreStepSafety[num] = -1.0;
105 fNewSafetyComputed[num]= -1.0;
106 }
107}
108
109// ----------------------------------------------------------------------------
110// Destructor
111//
112G4PathFinder::~G4PathFinder()
113{
114 delete fpMultiNavigator;
115}
116
117// ----------------------------------------------------------------------------
118//
119void
120G4PathFinder::EnableParallelNavigation(G4bool enableChoice)
121{
122 G4Navigator *navigatorForPropagation=0, *massNavigator=0;
123
124 massNavigator= fpTransportManager->GetNavigatorForTracking();
125 if( enableChoice )
126 {
127 navigatorForPropagation= fpMultiNavigator;
128
129 // Enable SafetyHelper to use PF
130 //
131 fpTransportManager->GetSafetyHelper()->EnableParallelNavigation(true);
132 }
133 else
134 {
135 navigatorForPropagation= massNavigator;
136
137 // Disable SafetyHelper to use PF
138 //
139 fpTransportManager->GetSafetyHelper()->EnableParallelNavigation(false);
140 }
141 fpFieldPropagator->SetNavigatorForPropagating(navigatorForPropagation);
142}
143
144// ----------------------------------------------------------------------------
145//
146G4double
147G4PathFinder::ComputeStep( const G4FieldTrack &InitialFieldTrack,
148 G4double proposedStepLength,
149 G4int navigatorNo,
150 G4int stepNo, // find next step
151 G4double &pNewSafety, // for this geom
152 ELimited &limitedStep,
153 G4FieldTrack &EndState,
154 G4VPhysicalVolume* currentVolume)
155{
156 G4double possibleStep= -1.0;
157
158#ifdef G4DEBUG_PATHFINDER
159 if( fVerboseLevel > 2 )
160 {
161 G4cout << " -------------------------" << G4endl;
162 G4cout << " G4PathFinder::ComputeStep - entered " << G4endl;
163 G4cout << " - stepNo = " << std::setw(4) << stepNo << " "
164 << " navigatorId = " << std::setw(2) << navigatorNo << " "
165 << " proposed step len = " << proposedStepLength << " " << G4endl;
166 G4cout << " PF::ComputeStep requested step "
167 << " from " << InitialFieldTrack.GetPosition()
168 << " dir " << InitialFieldTrack.GetMomentumDirection() << G4endl;
169 }
170#endif
171#ifdef G4VERBOSE
172 if( navigatorNo >= fNoActiveNavigators )
173 {
174 G4cerr << "ERROR - G4PathFinder::ComputeStep()" << G4endl
175 << " Requested Navigator ID = " << navigatorNo << G4endl
176 << " Number of active navigators = " << fNoActiveNavigators
177 << G4endl;
178 G4Exception("G4PathFinder::ComputeStep()", "InvalidSetup",
179 FatalException, "Bad Navigator ID !");
180 }
181#endif
182
183 if( fNewTrack || (stepNo != fLastStepNo) )
184 {
185 // This is a new track or a new step, so we must make the step
186 // ( else we can simply retrieve its results for this Navigator Id )
187
188 G4FieldTrack currentState= InitialFieldTrack;
189
190 fCurrentStepNo = stepNo;
191
192 // Check whether a process shifted the position
193 // since the last step -- by physics processes
194 //
195 G4ThreeVector newPosition = InitialFieldTrack.GetPosition();
196 G4ThreeVector moveVector= newPosition - fLastLocatedPosition;
197 G4double moveLenSq= moveVector.mag2();
198 if( moveLenSq > kCarTolerance * kCarTolerance )
199 {
200 G4ThreeVector newDirection = InitialFieldTrack.GetMomentumDirection();
201#ifdef G4DEBUG_PATHFINDER
202 if( fVerboseLevel > 2 )
203 {
204 G4double moveLen= std::sqrt( moveLenSq );
205 G4cout << " G4PathFinder::ComputeStep : Point moved since last step "
206 << " -- at step # = " << stepNo << G4endl
207 << " by " << moveLen << " to " << newPosition << G4endl;
208 }
209#endif
210 MovePoint(); // Unintentional changed -- ????
211
212 // Relocate to cope with this move -- else could abort !?
213 //
214 Locate( newPosition, newDirection );
215 }
216
217 // Check whether the particle have an (EM) field force exerting upon it
218 //
219 G4double particleCharge= currentState.GetCharge();
220
221 G4FieldManager* fieldMgr=0;
222 G4bool fieldExertsForce = false ;
223 if( (particleCharge != 0.0) )
224 {
225 fieldMgr= fpFieldPropagator->FindAndSetFieldManager( currentVolume );
226
227 // Protect for case where field manager has no field (= field is zero)
228 //
229 fieldExertsForce = (fieldMgr != 0)
230 && (fieldMgr->GetDetectorField() != 0);
231 }
232 fFieldExertedForce = fieldExertsForce; // Store for use in later calls
233 // referring to this 'step'.
234
235 fNoGeometriesLimiting= -1; // At start of track, no process limited step
236 if( fieldExertsForce )
237 {
238 DoNextCurvedStep( currentState, proposedStepLength, currentVolume );
239 //--------------
240 }else{
241 DoNextLinearStep( currentState, proposedStepLength );
242 //--------------
243 }
244 fLastStepNo= stepNo;
245
246 if ( (fNoGeometriesLimiting < 0)
247 || (fNoGeometriesLimiting > fNoActiveNavigators) )
248 {
249 G4cout << "ERROR - G4PathFinder::ComputeStep()" << G4endl
250 << " Number of geometries limiting step = "
251 << fNoGeometriesLimiting << G4endl;
252 G4Exception("G4PathFinder::ComputeStep()",
253 "NumGeometriesOutOfRange", FatalException,
254 "Number of geometries limiting the step not set.");
255 }
256 }
257#ifdef G4DEBUG_PATHFINDER
258 else
259 {
260 if( proposedStepLength < fTrueMinStep ) // For 2nd+ geometry
261 {
262 G4cout << "ERROR - G4PathFinder::ComputeStep()" << G4endl
263 << " Problem in step size request." << G4endl
264 << " Being requested to make a step which is shorter"
265 << " than the minimum Step " << G4endl
266 << " already computed for any Navigator/geometry during"
267 << " this tracking-step: " << G4endl;
268 G4cout << " This can happen due to an error in process ordering."
269 << G4endl;
270 G4cout << " Check that all physics processes are registered"
271 << G4endl
272 << " before all processes with a navigator/geometry."
273 << G4endl;
274 G4cout << " If using pre-packaged physics list and/or"
275 << G4endl
276 << " functionality, please report this error."
277 << G4endl << G4endl;
278 G4cout << " Additional information for problem: " << G4endl
279 << " Steps request/proposed = " << proposedStepLength
280 << G4endl
281 << " MinimumStep (true) = " << fTrueMinStep
282 << G4endl
283 << " MinimumStep (navraw) = " << fMinStep
284 << G4endl
285 << " Navigator raw return value" << G4endl
286 << " Requested step now = " << proposedStepLength
287 << G4endl
288 << " Difference min-req = "
289 << fTrueMinStep-proposedStepLength << G4endl;
290 G4cout << " -- Step info> stepNo= " << stepNo
291 << " last= " << fLastStepNo
292 << " newTr= " << fNewTrack << G4endl;
293 G4cerr << "ERROR - G4PathFinder::ComputeStep()" << G4endl
294 << " Problem in step size request. " << G4endl
295 << " Error can be caused by incorrect process ordering."
296 << G4endl
297 << " Please see more information in standard output."
298 << G4endl;
299 G4Exception("G4PathFinder::ComputeStep()",
300 "ReductionOfRequestedStepSizeBelowMinimum",
301 FatalException,
302 "Not part of specification - not implemented.");
303 }
304 else
305 {
306 // This is neither a new track nor a new step -- just another
307 // client accessing information for the current track, step
308 // We will simply retrieve the results of the synchronous
309 // stepping for this Navigator Id below.
310 //
311 if( fVerboseLevel > 1 )
312 {
313 G4cout << " G4P::CS -> Not calling DoNextLinearStep: "
314 << " stepNo= " << stepNo << " last= " << fLastStepNo
315 << " new= " << fNewTrack << " Step already done" << G4endl;
316 }
317 }
318 }
319#endif
320
321 fNewTrack= false;
322
323 // Prepare the information to return
324
325 pNewSafety = fCurrentPreStepSafety[ navigatorNo ];
326 limitedStep = fLimitedStep[ navigatorNo ];
327 fRelocatedPoint= false;
328
329 possibleStep= std::min(proposedStepLength, fCurrentStepSize[ navigatorNo ]);
330 EndState = fEndState; // now corrected for smaller step, if needed
331
332#ifdef G4DEBUG_PATHFINDER
333 if( fVerboseLevel > 0 )
334 {
335 G4cout << " G4PathFinder::ComputeStep returns "
336 << fCurrentStepSize[ navigatorNo ]
337 << " for Navigator " << navigatorNo
338 << " Limited step = " << limitedStep
339 << " Safety(mm) = " << pNewSafety / mm
340 << G4endl;
341 }
342#endif
343
344 return possibleStep;
345}
346
347// ----------------------------------------------------------------------
348
349void
350G4PathFinder::PrepareNewTrack( const G4ThreeVector& position,
351 const G4ThreeVector& direction,
352 G4VPhysicalVolume* massStartVol)
353{
354 // Key purposes:
355 // - Check and cache set of active navigators
356 // - Reset state for new track
357
358 G4int num=0;
359
360 EnableParallelNavigation(true);
361 // Switch PropagatorInField to use MultiNavigator
362
363 fpTransportManager->GetSafetyHelper()->InitialiseHelper();
364 // Reinitialise state of safety helper -- avoid problems with overlaps
365
366 fNewTrack= true;
367 this->MovePoint(); // Signal further that the last status is wiped
368
369 // Message the G4NavigatorPanel / Dispatcher to find active navigators
370 //
371 std::vector<G4Navigator*>::iterator pNavigatorIter;
372
373 fNoActiveNavigators= fpTransportManager-> GetNoActiveNavigators();
374 if( fNoActiveNavigators > fMaxNav )
375 {
376 G4cerr << "ERROR - G4PathFinder::PrepareNewTrack()" << G4endl
377 << " Too many active Navigators. G4PathFinder fails."
378 << G4endl
379 << " Transportation Manager has "
380 << fNoActiveNavigators << " active navigators." << G4endl
381 << " This is more than the number allowed = "
382 << fMaxNav << " !" << G4endl;
383 G4Exception("G4PathFinder::PrepareNewTrack()", "TooManyNavigators",
384 FatalException, "Too many active Navigators / worlds");
385 }
386
387 fpMultiNavigator->PrepareNavigators();
388 //------------------------------------
389
390 pNavigatorIter= fpTransportManager->GetActiveNavigatorsIterator();
391 for( num=0; num< fNoActiveNavigators; ++pNavigatorIter,++num )
392 {
393 // Keep information in C-array ... for creating touchables - at least
394
395 fpNavigator[num] = *pNavigatorIter;
396 fLimitTruth[num] = false;
397 fLimitedStep[num] = kDoNot;
398 fCurrentStepSize[num] = 0.0;
399 fLocatedVolume[num] = 0;
400 }
401 fNoGeometriesLimiting= 0; // At start of track, no process limited step
402
403 // In case of one geometry, the tracking will have done the locating!!
404
405 if( fNoActiveNavigators > 1 )
406 {
407 Locate( position, direction, false );
408 }
409 else
410 {
411 // Update state -- depending on the tracking's call to Mass Navigator
412
413 fLastLocatedPosition= position;
414 fLocatedVolume[0]= massStartVol; // This information must be given
415 // by transportation
416 fLimitedStep[0] = kDoNot;
417 fCurrentStepSize[0] = 0.0;
418 }
419
420 // Reset Safety Information -- as in case of overlaps this can cause
421 // inconsistencies ...
422 //
423 fMinSafety_PreStepPt= fPreSafetyMinValue= fMinSafety_atSafLocation= 0.0;
424
425 for( num=0; num< fNoActiveNavigators; ++num )
426 {
427 fPreSafetyValues[num]= 0.0;
428 fNewSafetyComputed[num]= 0.0;
429 fCurrentPreStepSafety[num] = 0.0;
430 }
431
432 // The first location for each Navigator must be non-relative
433 // or else call ResetStackAndState() for each Navigator
434
435 fRelocatedPoint= false;
436}
437
438void G4PathFinder::ReportMove( const G4ThreeVector& OldVector,
439 const G4ThreeVector& NewVector,
440 const G4String& Quantity ) const
441{
442 G4ThreeVector moveVec = ( NewVector - OldVector );
443
444 G4int prc= G4cerr.precision(12);
445 G4cerr << G4endl
446 << "WARNING - G4PathFinder::ReportMove()" << G4endl
447 << " Endpoint moved between value returned by ComputeStep()"
448 << " and call to Locate(). " << G4endl
449 << " Change of " << Quantity << " is "
450 << moveVec.mag() / mm << " mm long" << G4endl
451 << " and its vector is "
452 << (1.0/mm) * moveVec << " mm " << G4endl
453 << " Endpoint of ComputeStep() was " << OldVector << G4endl
454 << " and current position to locate is " << NewVector
455 << G4endl;
456 G4cerr.precision(prc);
457}
458
459void
460G4PathFinder::Locate( const G4ThreeVector& position,
461 const G4ThreeVector& direction,
462 G4bool relative)
463{
464 // Locate the point in each geometry
465
466 std::vector<G4Navigator*>::iterator pNavIter=
467 fpTransportManager->GetActiveNavigatorsIterator();
468
469 G4ThreeVector lastEndPosition= fEndState.GetPosition();
470 G4ThreeVector moveVec = (position - lastEndPosition );
471 G4double moveLenSq= moveVec.mag2();
472 if( (!fNewTrack) && (!fRelocatedPoint)
473 && ( moveLenSq> kCarTolerance*kCarTolerance ) )
474 {
475 ReportMove( position, lastEndPosition, "Position" );
476 G4Exception("G4PathFinder::Locate", "201-LocateUnexpectedPoint",
477 JustWarning,
478 "Location is not where last ComputeStep() ended.");
479 }
480 fLastLocatedPosition= position;
481
482#ifdef G4DEBUG_PATHFINDER
483 if( fVerboseLevel > 2 )
484 {
485 G4cout << G4endl;
486 G4cout << " G4PathFinder::Locate : entered " << G4endl;
487 G4cout << " -------------------- -------" << G4endl;
488 G4cout << " Locating at position " << position
489 << " with direction " << direction
490 << " relative= " << relative << G4endl;
491 if ( (fVerboseLevel > 1) || ( moveLenSq > 0.0) )
492 {
493 G4cout << " lastEndPosition = " << lastEndPosition
494 << " moveVec = " << moveVec
495 << " newTr = " << fNewTrack
496 << " relocated = " << fRelocatedPoint << G4endl;
497 }
498
499 G4cout << " Located at " << position ;
500 if( fNoActiveNavigators > 1 ) { G4cout << G4endl; }
501 }
502#endif
503
504 for ( register G4int num=0; num< fNoActiveNavigators ; ++pNavIter,++num )
505 {
506 // ... who limited the step ....
507
508 if( fLimitTruth[num] ) { (*pNavIter)->SetGeometricallyLimitedStep(); }
509
510 G4VPhysicalVolume *pLocated=
511 (*pNavIter)->LocateGlobalPointAndSetup( position, &direction,
512 relative,
513 false);
514 // Set the state related to the location
515 //
516 fLocatedVolume[num] = pLocated;
517
518 // Clear state related to the step
519 //
520 fLimitedStep[num] = kDoNot;
521 fCurrentStepSize[num] = 0.0;
522
523#ifdef G4DEBUG_PATHFINDER
524 if( fVerboseLevel > 2 )
525 {
526 G4cout << " - In world " << num << " geomLimStep= " << fLimitTruth[num]
527 << " gives volume= " << pLocated ;
528 if( pLocated )
529 {
530 G4cout << " name = '" << pLocated->GetName() << "'";
531 G4cout << " - CopyNo= " << pLocated->GetCopyNo();
532 }
533 G4cout << G4endl;
534 }
535#endif
536 }
537
538 fRelocatedPoint= false;
539}
540
541void G4PathFinder::ReLocate( const G4ThreeVector& position )
542{
543 // Locate the point in each geometry
544
545 std::vector<G4Navigator*>::iterator pNavIter=
546 fpTransportManager->GetActiveNavigatorsIterator();
547
548 // Check that this relocation does not violate safety
549 // - at endpoint (computed from start point) AND
550 // - at last safety location (likely just called)
551
552 G4ThreeVector lastEndPosition= fEndState.GetPosition();
553
554 // Calculate end-point safety ...
555 //
556 G4double DistanceStartEnd= (lastEndPosition - fPreStepLocation).mag();
557 G4double endPointSafety_raw = fMinSafety_PreStepPt - DistanceStartEnd;
558 G4double endPointSafety_Est1 = std::max( 0.0, endPointSafety_raw );
559
560 // ... and check move from endpoint against this endpoint safety
561 //
562 G4ThreeVector moveVecEndPos = position - lastEndPosition;
563 G4double moveLenEndPosSq = moveVecEndPos.mag2();
564
565 // Check that move from endpoint of last step is within safety
566 // -- or check against last location or relocation ??
567 //
568 G4ThreeVector moveVecSafety= position - fSafetyLocation;
569 G4double moveLenSafSq= moveVecSafety.mag2();
570
571 G4double distCheckEnd_sq= ( moveLenEndPosSq - endPointSafety_Est1
572 *endPointSafety_Est1 );
573 G4double distCheckSaf_sq= ( moveLenSafSq - fMinSafety_atSafLocation
574 *fMinSafety_atSafLocation );
575
576 G4bool longMoveEnd = distCheckEnd_sq > 0.0;
577 G4bool longMoveSaf = distCheckSaf_sq > 0.0;
578
579 G4double revisedSafety= 0.0;
580
581 if( (!fNewTrack) && ( longMoveEnd && longMoveSaf ) )
582 {
583 G4ThreeVector LastSafetyLocation;
584 // Copy to keep last value - and restore
585
586 LastSafetyLocation= fSafetyLocation;
587
588 // Recompute ComputeSafety for end position
589 //
590 revisedSafety= ComputeSafety(lastEndPosition);
591
592 // Reset the state of last call to ComputeSafety
593 //
594 ComputeSafety( LastSafetyLocation );
595
596#ifdef G4DEBUG_PATHFINDER
597
598 const G4double kRadTolerance =
599 G4GeometryTolerance::GetInstance()->GetRadialTolerance();
600 const G4double cErrorTolerance=1e-12;
601 // Maximum relative error from roundoff of arithmetic
602
603 G4double distCheckRevisedEnd= moveLenEndPosSq-revisedSafety*revisedSafety;
604
605 G4bool longMoveRevisedEnd= ( distCheckRevisedEnd > 0. ) ;
606
607 G4double moveMinusSafety= 0.0;
608 G4double moveLenEndPosition= std::sqrt( moveLenEndPosSq );
609 moveMinusSafety = moveLenEndPosition - revisedSafety;
610
611 if ( longMoveRevisedEnd && (moveMinusSafety > 0.0 )
612 && ( revisedSafety > 0.0 ) )
613 {
614 // Take into account possibility of roundoff error causing
615 // this apparent move further than safety
616
617 if( fVerboseLevel > 0 )
618 {
619 G4cout << " G4PF:Relocate> Ratio to revised safety is "
620 << std::fabs(moveMinusSafety)/revisedSafety << G4endl;
621 }
622
623 G4double absMoveMinusSafety= std::fabs(moveMinusSafety);
624 G4bool smallRatio= absMoveMinusSafety < kRadTolerance * revisedSafety ;
625 G4double maxCoordPos = std::max(
626 std::max( std::fabs(position.x()),
627 std::fabs(position.y())),
628 std::fabs(position.z()) );
629 G4bool smallValue= absMoveMinusSafety < cErrorTolerance * maxCoordPos;
630 if( ! (smallRatio || smallValue) )
631 {
632 G4cout << " G4PF:Relocate> Ratio to revised safety is "
633 << std::fabs(moveMinusSafety)/revisedSafety << G4endl;
634 G4cout << " Difference of move and safety is not very small."
635 << G4endl;
636 }
637 else
638 {
639 moveMinusSafety = 0.0;
640 longMoveRevisedEnd = false; // Numerical issue -- not too long!
641
642 G4cout << " Difference of move & safety is very small in magnitude, "
643 << absMoveMinusSafety << G4endl;
644 if( smallRatio )
645 {
646 G4cout << " ratio to safety " << revisedSafety
647 << " is " << absMoveMinusSafety / revisedSafety
648 << "smaller than " << kRadTolerance << " of safety ";
649 }
650 else
651 {
652 G4cout << " as fraction " << absMoveMinusSafety / maxCoordPos
653 << " of position vector max-coord " << maxCoordPos
654 << " smaller than " << cErrorTolerance ;
655 }
656 G4cout << " -- reset moveMinusSafety to "
657 << moveMinusSafety << G4endl;
658 }
659 }
660
661 if ( longMoveEnd && longMoveSaf
662 && longMoveRevisedEnd && (moveMinusSafety>0.0) )
663 {
664 G4int oldPrec= G4cout.precision(9);
665 G4cout << " Problem in G4PathFinder::Relocate() " << G4endl;
666 G4cout << " Moved from last endpoint by " << moveLenEndPosition
667 << " compared to end safety (from preStep point) = "
668 << endPointSafety_Est1 << G4endl;
669
670 G4cout << " --> last PreSafety Location was " << fPreSafetyLocation
671 << G4endl;
672 G4cout << " safety value = " << fPreSafetyMinValue << G4endl;
673 G4cout << " --> last PreStep Location was " << fPreStepLocation
674 << G4endl;
675 G4cout << " safety value = " << fMinSafety_PreStepPt << G4endl;
676 G4cout << " --> last EndStep Location was " << lastEndPosition
677 << G4endl;
678 G4cout << " safety value = " << endPointSafety_Est1
679 << " raw-value = " << endPointSafety_raw << G4endl;
680 G4cout << " --> Calling again at this endpoint, we get "
681 << revisedSafety << " as safety value." << G4endl;
682 G4cout << " --> last position for safety " << fSafetyLocation
683 << G4endl;
684 G4cout << " its safety value = " << fMinSafety_atSafLocation
685 << G4endl;
686 G4cout << " move from safety location = "
687 << std::sqrt(moveLenSafSq) << G4endl
688 << " again= " << moveVecSafety.mag() << G4endl;
689 G4cout << " safety - Move-from-end= "
690 << revisedSafety - moveLenEndPosition
691 << " (negative is Bad.)" << G4endl;
692
693 G4cout << " Debug: distCheckRevisedEnd = "
694 << distCheckRevisedEnd << G4endl;
695
696 ReportMove( lastEndPosition, position, "Position" );
697 G4Exception( "G4PathFinder::ReLocate", "205-RelocatePointTooFar",
698 FatalException,
699 "ReLocation is further than end-safety value.");
700 G4cout.precision(oldPrec);
701 }
702
703#endif
704 }
705
706#ifdef G4DEBUG_PATHFINDER
707 if( fVerboseLevel > 2 )
708 {
709 G4cout << G4endl;
710 G4cout << " G4PathFinder::ReLocate : entered " << G4endl;
711 G4cout << " ---------------------- -------" << G4endl;
712 G4cout << " *Re*Locating at position " << position << G4endl;
713 // << " with direction " << direction
714 // << " relative= " << relative << G4endl;
715 if ( (fVerboseLevel > -1) || ( moveLenEndPosSq > 0.0) )
716 {
717 G4cout << " lastEndPosition = " << lastEndPosition
718 << " moveVec from step-end = " << moveVecEndPos
719 << " is new Track = " << fNewTrack
720 << " relocated = " << fRelocatedPoint << G4endl;
721 }
722 }
723#endif
724
725 for ( register G4int num=0; num< fNoActiveNavigators ; ++pNavIter,++num )
726 {
727 // ... none limited the step
728
729 (*pNavIter)->LocateGlobalPointWithinVolume( position );
730
731 // Clear state related to the step
732 //
733 fLimitedStep[num] = kDoNot;
734 fCurrentStepSize[num] = 0.0;
735 fLimitTruth[num] = false;
736 }
737
738 fLastLocatedPosition= position;
739 fRelocatedPoint= false;
740
741#ifdef G4DEBUG_PATHFINDER
742 if( fVerboseLevel > 2 )
743 {
744 G4cout << " G4PathFinder::ReLocate : exiting "
745 << " at position " << fLastLocatedPosition << G4endl << G4endl;
746 }
747#endif
748}
749
750// -----------------------------------------------------------------------------
751
752G4double G4PathFinder::ComputeSafety( const G4ThreeVector& position )
753{
754 // Recompute safety for the relevant point
755
756 G4double minSafety= kInfinity;
757
758 std::vector<G4Navigator*>::iterator pNavigatorIter;
759 pNavigatorIter= fpTransportManager-> GetActiveNavigatorsIterator();
760
761 for( register G4int num=0; num<fNoActiveNavigators; ++pNavigatorIter,++num )
762 {
763 G4double safety = (*pNavigatorIter)->ComputeSafety( position );
764 if( safety < minSafety ) { minSafety = safety; }
765 fNewSafetyComputed[num]= safety;
766 }
767
768 fSafetyLocation= position;
769 fMinSafety_atSafLocation = minSafety;
770
771#ifdef G4DEBUG_PATHFINDER
772 if( fVerboseLevel > 1 )
773 {
774 G4cout << " G4PathFinder::ComputeSafety - returns "
775 << minSafety << " at location " << position << G4endl;
776 }
777#endif
778 return minSafety;
779}
780
781
782// -----------------------------------------------------------------------------
783
784G4TouchableHandle
785G4PathFinder::CreateTouchableHandle( G4int navId ) const
786{
787#ifdef G4DEBUG_PATHFINDER
788 if( fVerboseLevel > 2 )
789 {
790 G4cout << "G4PathFinder::CreateTouchableHandle : navId = "
791 << navId << " -- " << GetNavigator(navId) << G4endl;
792 }
793#endif
794
795 G4TouchableHistory* touchHist;
796 touchHist= GetNavigator(navId) -> CreateTouchableHistory();
797
798 G4VPhysicalVolume* locatedVolume= fLocatedVolume[navId];
799 if( locatedVolume == 0 )
800 {
801 // Workaround to ensure that the touchable is fixed !! // TODO: fix
802
803 touchHist->UpdateYourself( locatedVolume, touchHist->GetHistory() );
804 }
805
806#ifdef G4DEBUG_PATHFINDER
807 if( fVerboseLevel > 2 )
808 {
809 G4String VolumeName("None");
810 if( locatedVolume ) { VolumeName= locatedVolume->GetName(); }
811 G4cout << " Touchable History created at address " << touchHist
812 << " volume = " << locatedVolume << " name= " << VolumeName
813 << G4endl;
814 }
815#endif
816
817 return G4TouchableHandle(touchHist);
818}
819
820G4double
821G4PathFinder::DoNextLinearStep( const G4FieldTrack &initialState,
822 G4double proposedStepLength )
823{
824 std::vector<G4Navigator*>::iterator pNavigatorIter;
825 G4double safety= 0.0, step=0.0;
826 G4double minSafety= kInfinity, minStep;
827
828 const G4int IdTransport= 0; // Id of Mass Navigator !!
829 register G4int num=0;
830
831#ifdef G4DEBUG_PATHFINDER
832 if( fVerboseLevel > 2 )
833 {
834 G4cout << " G4PathFinder::DoNextLinearStep : entered " << G4endl;
835 G4cout << " Input field track= " << initialState << G4endl;
836 G4cout << " Requested step= " << proposedStepLength << G4endl;
837 }
838#endif
839
840 G4ThreeVector initialPosition= initialState.GetPosition();
841 G4ThreeVector initialDirection= initialState.GetMomentumDirection();
842
843 G4ThreeVector OriginShift = initialPosition - fPreSafetyLocation;
844 G4double MagSqShift = OriginShift.mag2() ;
845 G4double MagShift; // Only given value if it larger than minimum safety
846
847 G4double fullSafety; // For all geometries, for prestep point
848
849 // Potential optimisation using Maximum Value of safety!
850 // if( MagSqShift >= sqr(fPreSafetyMaxValue ) ){
851 // MagShift= kInfinity; // Not a useful value -- all will not use/ignore
852 // else
853 // MagShift= std::sqrt(MagSqShift) ;
854
855 MagShift= std::sqrt(MagSqShift) ;
856 if( MagSqShift >= sqr(fPreSafetyMinValue ) )
857 {
858 fullSafety = 0.0 ;
859 }
860 else
861 {
862 fullSafety = fPreSafetyMinValue - MagShift;
863 }
864
865#ifdef G4PATHFINDER_OPTIMISATION
866
867 if( proposedStepLength < fullSafety )
868 {
869 // Move is smaller than all safeties
870 // -> so we do not have to move the safety center
871
872 fPreStepCenterRenewed= false;
873
874 for( num=0; num< fNoActiveNavigators; ++num )
875 {
876 fCurrentStepSize[num]= kInfinity;
877 safety = std::max( 0.0, fPreSafetyValues[num] - MagShift);
878 minSafety= std::min( safety, minSafety );
879 fCurrentPreStepSafety[num]= safety;
880 }
881 minStep= kInfinity;
882
883#ifdef G4DEBUG_PATHFINDER
884 if( fVerboseLevel > 2 )
885 {
886 G4cout << "G4PathFinder::DoNextLinearStep : Quick Stepping. " << G4endl
887 << " proposedStepLength " << proposedStepLength
888 << " < (full) safety = " << fullSafety
889 << " at " << initialPosition
890 << G4endl;
891 }
892#endif
893 }
894 else
895#endif // End of G4PATHFINDER_OPTIMISATION 1
896 {
897 // Move is larger than at least one of the safeties
898 // -> so we must move the safety center!
899
900 fPreStepCenterRenewed= true;
901 pNavigatorIter= fpTransportManager-> GetActiveNavigatorsIterator();
902
903 minStep= kInfinity; // Not proposedStepLength;
904
905 for( num=0; num< fNoActiveNavigators; ++pNavigatorIter,++num )
906 {
907 safety = std::max( 0.0, fPreSafetyValues[num] - MagShift);
908
909#ifdef G4PATHFINDER_OPTIMISATION
910 if( proposedStepLength <= safety ) // Should be just < safety ?
911 {
912 // The Step is guaranteed to be taken
913
914 step= kInfinity; // ComputeStep Would return this
915
916#ifdef G4DEBUG_PATHFINDER
917 G4cout.precision(8);
918 G4cout << "PathFinder::ComputeStep> small proposed step = "
919 << proposedStepLength
920 << " <= safety = " << safety << " for nav " << num
921 << " Step fully taken. " << G4endl;
922#endif
923 }
924 else
925#endif // End of G4PATHFINDER_OPTIMISATION 2
926 {
927#ifdef G4DEBUG_PATHFINDER
928 G4double previousSafety= safety;
929#endif
930 step= (*pNavigatorIter)->ComputeStep( initialPosition,
931 initialDirection,
932 proposedStepLength,
933 safety );
934 minStep = std::min( step, minStep);
935
936 // TODO: consider whether/how to reduce the proposed step
937 // to the latest minStep value - to reduce calculations
938
939#ifdef G4DEBUG_PATHFINDER
940 if( fVerboseLevel > 0)
941 {
942 G4cout.precision(8);
943 G4cout << "PathFinder::ComputeStep> long proposed step = "
944 << proposedStepLength
945 << " > safety = " << previousSafety
946 << " for nav " << num
947 << " . New safety = " << safety << " step= " << step
948 << G4endl;
949 }
950#endif
951 }
952 fCurrentStepSize[num] = step;
953
954 // Save safety value, must be done for all geometries "together"
955 // (even if not recomputed using call to ComputeStep)
956 // since they share the fPreSafetyLocation
957
958 fPreSafetyValues[num]= safety;
959 fCurrentPreStepSafety[num]= safety;
960
961 minSafety= std::min( safety, minSafety );
962
963#ifdef G4DEBUG_PATHFINDER
964 if( fVerboseLevel > 2 )
965 {
966 G4cout << "G4PathFinder::DoNextLinearStep : Navigator ["
967 << num << "] -- step size " << step << G4endl;
968 }
969#endif
970 }
971
972 // Only change these when safety is recalculated
973 // it is good/relevant only for safety calculations
974
975 fPreSafetyLocation= initialPosition;
976 fPreSafetyMinValue= minSafety;
977 } // end of else for if( proposedStepLength <= fullSafety)
978
979 // For use in Relocation, need PreStep point location, min-safety
980 //
981 fPreStepLocation= initialPosition;
982 fMinSafety_PreStepPt= minSafety;
983
984 fMinStep= minStep;
985
986 if( fMinStep == kInfinity )
987 {
988 minStep = proposedStepLength; // Use this below for endpoint !!
989 }
990 fTrueMinStep = minStep;
991
992 // Set the EndState
993
994 G4ThreeVector endPosition;
995
996 fEndState= initialState;
997 endPosition= initialPosition + minStep * initialDirection ;
998
999#ifdef G4DEBUG_PATHFINDER
1000 if( fVerboseLevel > 1 )
1001 {
1002 G4cout << "G4PathFinder::DoNextLinearStep : "
1003 << " initialPosition = " << initialPosition
1004 << " and endPosition = " << endPosition<< G4endl;
1005 }
1006#endif
1007
1008 fEndState.SetPosition( endPosition );
1009 fEndState.SetProperTimeOfFlight( -1.000 ); // Not defined YET
1010
1011 if( fNoActiveNavigators == 1 )
1012 {
1013 G4bool transportLimited = (fMinStep!= kInfinity);
1014 fLimitTruth[IdTransport] = transportLimited;
1015 fLimitedStep[IdTransport] = transportLimited ? kUnique : kDoNot;
1016
1017 // Set fNoGeometriesLimiting - as WhichLimited does
1018 fNoGeometriesLimiting = transportLimited ? 1 : 0;
1019 }
1020 else
1021 {
1022 WhichLimited();
1023 }
1024
1025#ifdef G4DEBUG_PATHFINDER
1026 if( fVerboseLevel > 2 )
1027 {
1028 G4cout << " G4PathFinder::DoNextLinearStep : exits returning "
1029 << minStep << G4endl;
1030 G4cout << " Endpoint values = " << fEndState << G4endl;
1031 G4cout << G4endl;
1032 }
1033#endif
1034
1035 return minStep;
1036}
1037
1038void G4PathFinder::WhichLimited()
1039{
1040 // Flag which processes limited the step
1041
1042 G4int num=-1, last=-1;
1043 G4int noLimited=0;
1044 ELimited shared= kSharedOther;
1045
1046 const G4int IdTransport= 0; // Id of Mass Navigator !!
1047
1048 // Assume that [IdTransport] is Mass / Transport
1049 //
1050 G4bool transportLimited = (fCurrentStepSize[IdTransport] == fMinStep)
1051 && ( fMinStep!= kInfinity) ;
1052
1053 if( transportLimited ) {
1054 shared= kSharedTransport;
1055 }
1056
1057 for ( num= 0; num < fNoActiveNavigators; num++ )
1058 {
1059 G4bool limitedStep;
1060
1061 G4double step= fCurrentStepSize[num];
1062
1063 limitedStep = ( step == fMinStep ) && ( step != kInfinity);
1064
1065 fLimitTruth[ num ] = limitedStep;
1066 if( limitedStep )
1067 {
1068 noLimited++;
1069 fLimitedStep[num] = shared;
1070 last= num;
1071 }
1072 else
1073 {
1074 fLimitedStep[num] = kDoNot;
1075 }
1076 }
1077 fNoGeometriesLimiting= noLimited; // Save # processes limiting step
1078
1079 if( (last > -1) && (noLimited == 1 ) )
1080 {
1081 fLimitedStep[ last ] = kUnique;
1082 }
1083
1084#ifdef G4DEBUG_PATHFINDER
1085 if( fVerboseLevel > 1 )
1086 {
1087 PrintLimited(); // --> for tracing
1088 if( fVerboseLevel > 4 ) {
1089 G4cout << " G4PathFinder::WhichLimited - exiting. " << G4endl;
1090 }
1091 }
1092#endif
1093}
1094
1095void G4PathFinder::PrintLimited()
1096{
1097 // Report results -- for checking
1098
1099 G4cout << "G4PathFinder::PrintLimited reports: " ;
1100 G4cout << " Minimum step (true)= " << fTrueMinStep
1101 << " reported min = " << fMinStep
1102 << G4endl;
1103 if( (fCurrentStepNo <= 2) || (fVerboseLevel>=2) )
1104 {
1105 G4cout << std::setw(5) << " Step#" << " "
1106 << std::setw(5) << " NavId" << " "
1107 << std::setw(12) << " step-size " << " "
1108 << std::setw(12) << " raw-size " << " "
1109 << std::setw(12) << " pre-safety " << " "
1110 << std::setw(15) << " Limited / flag" << " "
1111 << std::setw(15) << " World " << " "
1112 << G4endl;
1113 }
1114 G4int num;
1115 for ( num= 0; num < fNoActiveNavigators; num++ )
1116 {
1117 G4double rawStep = fCurrentStepSize[num];
1118 G4double stepLen = fCurrentStepSize[num];
1119 if( stepLen > fTrueMinStep )
1120 {
1121 stepLen = fTrueMinStep; // did not limit (went as far as asked)
1122 }
1123 G4int oldPrec= G4cout.precision(9);
1124
1125 G4cout << std::setw(5) << fCurrentStepNo << " "
1126 << std::setw(5) << num << " "
1127 << std::setw(12) << stepLen << " "
1128 << std::setw(12) << rawStep << " "
1129 << std::setw(12) << fCurrentPreStepSafety[num] << " "
1130 << std::setw(5) << (fLimitTruth[num] ? "YES" : " NO") << " ";
1131 G4String limitedStr= LimitedString(fLimitedStep[num]);
1132 G4cout << " " << std::setw(15) << limitedStr << " ";
1133 G4cout.precision(oldPrec);
1134
1135 G4Navigator *pNav= GetNavigator( num );
1136 G4String WorldName( "Not-Set" );
1137 if (pNav)
1138 {
1139 G4VPhysicalVolume *pWorld= pNav->GetWorldVolume();
1140 if( pWorld )
1141 {
1142 WorldName = pWorld->GetName();
1143 }
1144 }
1145 G4cout << " " << WorldName ;
1146 G4cout << G4endl;
1147 }
1148
1149 if( fVerboseLevel > 4 )
1150 {
1151 G4cout << " G4PathFinder::PrintLimited - exiting. " << G4endl;
1152 }
1153}
1154
1155G4double
1156G4PathFinder::DoNextCurvedStep( const G4FieldTrack &initialState,
1157 G4double proposedStepLength,
1158 G4VPhysicalVolume* pCurrentPhysicalVolume )
1159{
1160 const G4double toleratedRelativeError= 1.0e-10;
1161 G4double minStep= kInfinity, newSafety=0.0;
1162 G4int numNav;
1163 G4FieldTrack fieldTrack= initialState;
1164 G4ThreeVector startPoint= initialState.GetPosition();
1165
1166#ifdef G4DEBUG_PATHFINDER
1167 G4int prc= G4cout.precision(9);
1168 if( fVerboseLevel > 2 )
1169 {
1170 G4cout << " G4PathFinder::DoNextCurvedStep ****** " << G4endl;
1171 G4cout << " Initial value of field track is " << fieldTrack
1172 << " and proposed step= " << proposedStepLength << G4endl;
1173 }
1174#endif
1175
1176 fPreStepCenterRenewed= true; // Always update PreSafety with PreStep point
1177
1178 if( fNoActiveNavigators > 1 )
1179 {
1180 // Calculate the safety values before making the step
1181
1182 G4double minSafety= kInfinity, safety;
1183 for( numNav=0; numNav < fNoActiveNavigators; ++numNav )
1184 {
1185 safety= fpNavigator[numNav]->ComputeSafety( startPoint );
1186 fPreSafetyValues[numNav]= safety;
1187 fCurrentPreStepSafety[numNav]= safety;
1188 minSafety = std::min( safety, minSafety );
1189 }
1190
1191 // Save safety value, related position
1192
1193 fPreSafetyLocation= startPoint;
1194 fPreSafetyMinValue= minSafety;
1195 fPreStepLocation= startPoint;
1196 fMinSafety_PreStepPt= minSafety;
1197 }
1198
1199 // Allow Propagator In Field to do the hard work, calling G4MultiNavigator
1200 //
1201 minStep= fpFieldPropagator->ComputeStep( fieldTrack,
1202 proposedStepLength,
1203 newSafety,
1204 pCurrentPhysicalVolume );
1205
1206 // fieldTrack now contains the endpoint information
1207 //
1208 fEndState= fieldTrack;
1209 fMinStep= minStep;
1210 fTrueMinStep = std::min( minStep, proposedStepLength );
1211
1212 if( fNoActiveNavigators== 1 )
1213 {
1214 // Update the 'PreSafety' sphere - as any ComputeStep was called
1215 // (must be done anyway in field)
1216
1217 fPreSafetyValues[0]= newSafety;
1218 fPreSafetyLocation= startPoint;
1219 fPreSafetyMinValue= newSafety;
1220
1221 // Update the current 'PreStep' point's values - mandatory
1222 //
1223 fCurrentPreStepSafety[0]= newSafety;
1224 fPreStepLocation= startPoint;
1225 fMinSafety_PreStepPt= newSafety;
1226 }
1227
1228#ifdef G4DEBUG_PATHFINDER
1229 if( fVerboseLevel > 2 )
1230 {
1231 G4cout << "G4PathFinder::DoNextCurvedStep : " << G4endl
1232 << " initialState = " << initialState << G4endl
1233 << " and endState = " << fEndState << G4endl;
1234 G4cout << "G4PathFinder::DoNextCurvedStep : "
1235 << " minStep = " << minStep
1236 << " proposedStepLength " << proposedStepLength
1237 << " safety = " << newSafety << G4endl;
1238 }
1239#endif
1240 G4double currentStepSize; // = 0.0;
1241 if( minStep < proposedStepLength ) // if == , then a boundary found at end ??
1242 {
1243 // Recover the remaining information from MultiNavigator
1244 // especially regarding which Navigator limited the step
1245
1246 G4int noLimited= 0; // No geometries limiting step
1247 for( numNav=0; numNav < fNoActiveNavigators; ++numNav )
1248 {
1249 G4double finalStep, lastPreSafety=0.0, minStepLast;
1250 ELimited didLimit;
1251 G4bool limited;
1252
1253 finalStep= fpMultiNavigator->ObtainFinalStep( numNav, lastPreSafety,
1254 minStepLast, didLimit );
1255
1256 // Calculate the step for this geometry, using the
1257 // final step (the only one which can differ.)
1258
1259 currentStepSize = fTrueMinStep;
1260 G4double diffStep= 0.0;
1261 if( (minStepLast != kInfinity) )
1262 {
1263 diffStep = (finalStep-minStepLast);
1264 if ( std::abs(diffStep) <= toleratedRelativeError * finalStep )
1265 {
1266 diffStep = 0.0;
1267 }
1268 currentStepSize += diffStep;
1269 }
1270 fCurrentStepSize[numNav] = currentStepSize;
1271
1272 // TODO: could refine the way to obtain safeties for > 1 geometries
1273 // - for pre step safety
1274 // notify MultiNavigator about new set of sub-steps
1275 // allow it to return this value in ObtainFinalStep
1276 // instead of lastPreSafety (or as well?)
1277 // - for final step start (available)
1278 // get final Step start from MultiNavigator
1279 // and corresponding safety values
1280 // and/or ALSO calculate ComputeSafety at endpoint
1281 // endSafety= fpNavigator[numNav]->ComputeSafety( endPoint );
1282
1283 fLimitedStep[numNav] = didLimit;
1284 fLimitTruth[numNav] = limited = (didLimit != kDoNot );
1285 if( limited ) { noLimited++; }
1286
1287#ifdef G4DEBUG_PATHFINDER
1288 G4bool StepError= (currentStepSize < 0)
1289 || ( (minStepLast != kInfinity) && (diffStep < 0) ) ;
1290 if( StepError || (fVerboseLevel > 2) )
1291 {
1292 G4String limitedString= LimitedString( fLimitedStep[numNav] );
1293
1294 G4cout << " G4PathFinder::ComputeStep. Geometry " << numNav
1295 << " step= " << fCurrentStepSize[numNav]
1296 << " from final-step= " << finalStep
1297 << " fTrueMinStep= " << fTrueMinStep
1298 << " minStepLast= " << minStepLast
1299 << " limited = " << (fLimitTruth[numNav] ? "YES" : " NO")
1300 << " ";
1301 G4cout << " status = " << limitedString << " #= " << didLimit
1302 << G4endl;
1303
1304 if( StepError )
1305 {
1306 G4cerr << " currentStepSize = " << currentStepSize
1307 << " diffStep= " << diffStep << G4endl;
1308 G4cerr << "ERROR in computing step size for this navigator."
1309 << G4endl;
1310 G4Exception("G4PathFinder::DoNextCurvedStep",
1311 "207-StepGoingBackwards", FatalException,
1312 "Incorrect calculation of step size for one navigator");
1313 }
1314 }
1315#endif
1316 } // for num Navigators
1317
1318 fNoGeometriesLimiting= noLimited; // Save # processes limiting step
1319 }
1320 else if ( (minStep == proposedStepLength)
1321 || (minStep == kInfinity)
1322 || ( std::abs(minStep-proposedStepLength)
1323 < toleratedRelativeError * proposedStepLength ) )
1324 {
1325 // In case the step was not limited, use default responses
1326 // --> all Navigators
1327 // Also avoid problems in case of PathFinder using safety to optimise
1328 // - it is possible that the Navigators were not called
1329 // if the safety was already satisfactory.
1330 // (In that case calling ObtainFinalStep gives invalid results.)
1331
1332 currentStepSize= minStep;
1333 for( numNav=0; numNav < fNoActiveNavigators; ++numNav )
1334 {
1335 fCurrentStepSize[numNav] = minStep;
1336 // Safety for endpoint ?? // Can eventuall improve it -- see TODO above
1337 fLimitedStep[numNav] = kDoNot;
1338 fLimitTruth[numNav] = false;
1339 }
1340 fNoGeometriesLimiting= 0; // Save # processes limiting step
1341 }
1342 else // (minStep > proposedStepLength) and not (minStep == kInfinity)
1343 {
1344 G4cerr << G4endl;
1345 G4cerr << "ERROR - G4PathFinder::DoNextCurvedStep()" << G4endl
1346 << " currentStepSize = " << minStep << " is larger than "
1347 << " proposed StepSize = " << proposedStepLength << "." << G4endl;
1348 G4Exception("G4PathFinder::DoNextCurvedStep()",
1349 "208-StepLongerThanRequested", FatalException,
1350 "Incorrect calculation of step size for one navigator.");
1351 }
1352
1353#ifdef G4DEBUG_PATHFINDER
1354 if( fVerboseLevel > 2 )
1355 {
1356 G4cout << " Exiting G4PathFinder::DoNextCurvedStep " << G4endl;
1357 PrintLimited();
1358 }
1359 G4cout.precision(prc);
1360#endif
1361
1362 return minStep;
1363}
1364
1365G4String& G4PathFinder::LimitedString( ELimited lim )
1366{
1367 static G4String StrDoNot("DoNot"),
1368 StrUnique("Unique"),
1369 StrUndefined("Undefined"),
1370 StrSharedTransport("SharedTransport"),
1371 StrSharedOther("SharedOther");
1372
1373 G4String* limitedStr;
1374 switch ( lim )
1375 {
1376 case kDoNot: limitedStr= &StrDoNot; break;
1377 case kUnique: limitedStr = &StrUnique; break;
1378 case kSharedTransport: limitedStr= &StrSharedTransport; break;
1379 case kSharedOther: limitedStr = &StrSharedOther; break;
1380 default: limitedStr = &StrUndefined; break;
1381 }
1382 return *limitedStr;
1383}
1384
1385void G4PathFinder::PushPostSafetyToPreSafety()
1386{
1387 fPreSafetyLocation= fSafetyLocation;
1388 fPreSafetyMinValue= fMinSafety_atSafLocation;
1389 for( register G4int nav=0; nav < fNoActiveNavigators; ++nav )
1390 {
1391 fPreSafetyValues[nav]= fNewSafetyComputed[nav];
1392 }
1393}
Note: See TracBrowser for help on using the repository browser.