Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros
G4PropagatorInField.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//
27//
28//
29// This class implements an algorithm to track a particle in a
30// non-uniform magnetic field. It utilises an ODE solver (with
31// the Runge - Kutta method) to evolve the particle, and drives it
32// until the particle has traveled a set distance or it enters a new
33// volume.
34//
35// 14.10.96 John Apostolakis, design and implementation
36// 17.03.97 John Apostolakis, renaming new set functions being added
37//
38// $Id$
39// GEANT4 tag $ Name: $
40// ---------------------------------------------------------------------------
41
42#include <iomanip>
43
45#include "G4ios.hh"
46#include "G4SystemOfUnits.hh"
47#include "G4ThreeVector.hh"
48#include "G4VPhysicalVolume.hh"
49#include "G4Navigator.hh"
52#include "G4ChordFinder.hh"
54
55///////////////////////////////////////////////////////////////////////////
56//
57// Constructors and destructor
58
60 G4FieldManager *detectorFieldMgr,
61 G4VIntersectionLocator *vLocator )
62 :
63 fMax_loop_count(1000),
64 fUseSafetyForOptimisation(true), // (false) is less sensitive to incorrect safety
65 fZeroStepThreshold( 0.0 ), // length of what is recognised as 'zero' step
66 fDetectorFieldMgr(detectorFieldMgr),
67 fpTrajectoryFilter( 0 ),
68 fNavigator(theNavigator),
69 fCurrentFieldMgr(detectorFieldMgr),
70 fSetFieldMgr(false),
71 fCharge(0.0), fInitialMomentumModulus(0.0), fMass(0.0),
72 End_PointAndTangent(G4ThreeVector(0.,0.,0.),
73 G4ThreeVector(0.,0.,0.),0.0,0.0,0.0,0.0,0.0),
74 fParticleIsLooping(false),
75 fNoZeroStep(0),
76 fVerboseLevel(0)
77{
78 if(fDetectorFieldMgr) { fEpsilonStep = fDetectorFieldMgr->GetMaximumEpsilonStep();}
79 else { fEpsilonStep= 1.0e-5; }
80 fActionThreshold_NoZeroSteps = 2;
81 fSevereActionThreshold_NoZeroSteps = 10;
82 fAbandonThreshold_NoZeroSteps = 50;
83 fFull_CurveLen_of_LastAttempt = -1;
84 fLast_ProposedStepLength = -1;
85 fLargestAcceptableStep = 1000.0 * meter;
86
87 fPreviousSftOrigin= G4ThreeVector(0.,0.,0.);
88 fPreviousSafety= 0.0;
90 fZeroStepThreshold= std::max( 1.0e5 * kCarTolerance, 1.0e-1 * micrometer );
91
92#ifdef G4DEBUG_FIELD
93 G4cout << " PiF: Zero Step Threshold set to "
94 << fZeroStepThreshold / millimeter
95 << " mm." << G4endl;
96 G4cout << " PiF: Value of kCarTolerance = "
97 << kCarTolerance / millimeter
98 << " mm. " << G4endl;
99#endif
100
101 // Definding Intersection Locator and his parameters
102 if(vLocator==0){
103 fIntersectionLocator= new G4MultiLevelLocator(theNavigator);
104 fAllocatedLocator=true;
105 }else{
106 fIntersectionLocator=vLocator;
107 fAllocatedLocator=false;
108 }
109 RefreshIntersectionLocator(); // Copy all relevant parameters
110}
111
113{
114 if(fAllocatedLocator)delete fIntersectionLocator;
115}
116
117// Update the IntersectionLocator with current parameters
118void
120{
121 fIntersectionLocator->SetEpsilonStepFor(fEpsilonStep);
122 fIntersectionLocator->SetDeltaIntersectionFor(fCurrentFieldMgr->GetDeltaIntersection());
123 fIntersectionLocator->SetChordFinderFor(GetChordFinder());
124 fIntersectionLocator->SetSafetyParametersFor( fUseSafetyForOptimisation);
125}
126///////////////////////////////////////////////////////////////////////////
127//
128// Compute the next geometric Step
129
132 G4FieldTrack& pFieldTrack,
133 G4double CurrentProposedStepLength,
134 G4double& currentSafety, // IN/OUT
135 G4VPhysicalVolume* pPhysVol)
136{
137 // If CurrentProposedStepLength is too small for finding Chords
138 // then return with no action (for now - TODO: some action)
139 //
140 if(CurrentProposedStepLength<kCarTolerance)
141 {
142 return kInfinity;
143 }
144
145 // Introducing smooth trajectory display (jacek 01/11/2002)
146 //
147 if (fpTrajectoryFilter)
148 {
149 fpTrajectoryFilter->CreateNewTrajectorySegment();
150 }
151
152 // Parameters for adaptive Runge-Kutta integration
153
154 G4double h_TrialStepSize; // 1st Step Size
155 G4double TruePathLength = CurrentProposedStepLength;
156 G4double StepTaken = 0.0;
157 G4double s_length_taken, epsilon ;
158 G4bool intersects;
159 G4bool first_substep = true;
160
161 G4double NewSafety;
162 fParticleIsLooping = false;
163
164 // If not yet done,
165 // Set the field manager to the local one if the volume has one,
166 // or to the global one if not
167 //
168 if( !fSetFieldMgr ) fCurrentFieldMgr= FindAndSetFieldManager( pPhysVol );
169 // For the next call, the field manager must again be set
170 fSetFieldMgr= false;
171
172 GetChordFinder()->SetChargeMomentumMass(fCharge,fInitialMomentumModulus,fMass);
173
174 // Values for Intersection Locator has to be updated on each call for the
175 // case that CurrentFieldManager has changed from the one of previous step
177
178 G4FieldTrack CurrentState(pFieldTrack);
179 G4FieldTrack OriginalState = CurrentState;
180
181 // If the Step length is "infinite", then an approximate-maximum Step
182 // length (used to calculate the relative accuracy) must be guessed.
183 //
184 if( CurrentProposedStepLength >= fLargestAcceptableStep )
185 {
186 G4ThreeVector StartPointA, VelocityUnit;
187 StartPointA = pFieldTrack.GetPosition();
188 VelocityUnit = pFieldTrack.GetMomentumDir();
189
190 G4double trialProposedStep = 1.e2 * ( 10.0 * cm +
191 fNavigator->GetWorldVolume()->GetLogicalVolume()->
192 GetSolid()->DistanceToOut(StartPointA, VelocityUnit) );
193 CurrentProposedStepLength= std::min( trialProposedStep,
194 fLargestAcceptableStep );
195 }
196 epsilon = fCurrentFieldMgr->GetDeltaOneStep() / CurrentProposedStepLength;
197 // G4double raw_epsilon= epsilon;
198 G4double epsilonMin= fCurrentFieldMgr->GetMinimumEpsilonStep();
199 G4double epsilonMax= fCurrentFieldMgr->GetMaximumEpsilonStep();;
200 if( epsilon < epsilonMin ) epsilon = epsilonMin;
201 if( epsilon > epsilonMax ) epsilon = epsilonMax;
202 SetEpsilonStep( epsilon );
203
204 // G4cout << "G4PiF: Epsilon of current step - raw= " << raw_epsilon
205 // << " final= " << epsilon << G4endl;
206
207 // Shorten the proposed step in case of earlier problems (zero steps)
208 //
209 if( fNoZeroStep > fActionThreshold_NoZeroSteps )
210 {
211 G4double stepTrial;
212
213 stepTrial= fFull_CurveLen_of_LastAttempt;
214 if( (stepTrial <= 0.0) && (fLast_ProposedStepLength > 0.0) )
215 stepTrial= fLast_ProposedStepLength;
216
217 G4double decreaseFactor = 0.9; // Unused default
218 if( (fNoZeroStep < fSevereActionThreshold_NoZeroSteps)
219 && (stepTrial > 100.0*fZeroStepThreshold) )
220 {
221 // Attempt quick convergence
222 //
223 decreaseFactor= 0.25;
224 }
225 else
226 {
227 // We are in significant difficulties, probably at a boundary that
228 // is either geometrically sharp or between very different materials.
229 // Careful decreases to cope with tolerance are required.
230 //
231 if( stepTrial > 100.0*fZeroStepThreshold )
232 decreaseFactor = 0.35; // Try decreasing slower
233 else if( stepTrial > 30.0*fZeroStepThreshold )
234 decreaseFactor= 0.5; // Try yet slower decrease
235 else if( stepTrial > 10.0*fZeroStepThreshold )
236 decreaseFactor= 0.75; // Try even slower decreases
237 else
238 decreaseFactor= 0.9; // Try very slow decreases
239 }
240 stepTrial *= decreaseFactor;
241
242#ifdef G4DEBUG_FIELD
243 G4cout << " G4PropagatorInField::ComputeStep(): " << G4endl
244 << " Decreasing step - "
245 << " decreaseFactor= " << std::setw(8) << decreaseFactor
246 << " stepTrial = " << std::setw(18) << stepTrial << " "
247 << " fZeroStepThreshold = " << fZeroStepThreshold << G4endl;
248 PrintStepLengthDiagnostic(CurrentProposedStepLength, decreaseFactor,
249 stepTrial, pFieldTrack);
250#endif
251 if( stepTrial == 0.0 ) // Change to make it < 0.1 * kCarTolerance ??
252 {
253 std::ostringstream message;
254 message << "Particle abandoned due to lack of progress in field."
255 << G4endl
256 << " Properties : " << pFieldTrack << G4endl
257 << " Attempting a zero step = " << stepTrial << G4endl
258 << " while attempting to progress after " << fNoZeroStep
259 << " trial steps. Will abandon step.";
260 G4Exception("G4PropagatorInField::ComputeStep()", "GeomNav1002",
261 JustWarning, message);
262 fParticleIsLooping= true;
263 return 0; // = stepTrial;
264 }
265 if( stepTrial < CurrentProposedStepLength )
266 CurrentProposedStepLength = stepTrial;
267 }
268 fLast_ProposedStepLength = CurrentProposedStepLength;
269
270 G4int do_loop_count = 0;
271 do
272 {
273 G4FieldTrack SubStepStartState = CurrentState;
274 G4ThreeVector SubStartPoint = CurrentState.GetPosition();
275
276 if( !first_substep) {
277 fNavigator->LocateGlobalPointWithinVolume( SubStartPoint );
278 }
279
280 // How far to attempt to move the particle !
281 //
282 h_TrialStepSize = CurrentProposedStepLength - StepTaken;
283
284 // Integrate as far as "chord miss" rule allows.
285 //
286 s_length_taken = GetChordFinder()->AdvanceChordLimited(
287 CurrentState, // Position & velocity
288 h_TrialStepSize,
289 fEpsilonStep,
290 fPreviousSftOrigin,
291 fPreviousSafety
292 );
293 // CurrentState is now updated with the final position and velocity.
294
295 fFull_CurveLen_of_LastAttempt = s_length_taken;
296
297 G4ThreeVector EndPointB = CurrentState.GetPosition();
298 G4ThreeVector InterSectionPointE;
299 G4double LinearStepLength;
300
301 // Intersect chord AB with geometry
302 intersects= IntersectChord( SubStartPoint, EndPointB,
303 NewSafety, LinearStepLength,
304 InterSectionPointE );
305 // E <- Intersection Point of chord AB and either volume A's surface
306 // or a daughter volume's surface ..
307
308 if( first_substep ) {
309 currentSafety = NewSafety;
310 } // Updating safety in other steps is potential future extention
311
312 if( intersects )
313 {
314 G4FieldTrack IntersectPointVelct_G(CurrentState); // FT-Def-Construct
315
316 // Find the intersection point of AB true path with the surface
317 // of vol(A), if it exists. Start with point E as first "estimate".
318 G4bool recalculatedEndPt= false;
319
320 G4bool found_intersection = fIntersectionLocator->
321 EstimateIntersectionPoint( SubStepStartState, CurrentState,
322 InterSectionPointE, IntersectPointVelct_G,
323 recalculatedEndPt,fPreviousSafety,fPreviousSftOrigin);
324 intersects = intersects && found_intersection;
325 if( found_intersection ) {
326 End_PointAndTangent= IntersectPointVelct_G; // G is our EndPoint ...
327 StepTaken = TruePathLength = IntersectPointVelct_G.GetCurveLength()
328 - OriginalState.GetCurveLength();
329 } else {
330 // intersects= false; // "Minor" chords do not intersect
331 if( recalculatedEndPt ){
332 CurrentState= IntersectPointVelct_G;
333 }
334 }
335 }
336 if( !intersects )
337 {
338 StepTaken += s_length_taken;
339 // For smooth trajectory display (jacek 01/11/2002)
340 if (fpTrajectoryFilter) {
341 fpTrajectoryFilter->TakeIntermediatePoint(CurrentState.GetPosition());
342 }
343 }
344 first_substep = false;
345
346#ifdef G4DEBUG_FIELD
347 if( fNoZeroStep > fActionThreshold_NoZeroSteps ) {
348 printStatus( SubStepStartState, // or OriginalState,
349 CurrentState, CurrentProposedStepLength,
350 NewSafety, do_loop_count, pPhysVol );
351 }
352 if( (fVerboseLevel > 1) && (do_loop_count > fMax_loop_count-10 )) {
353 if( do_loop_count == fMax_loop_count-9 ){
354 G4cout << " G4PropagatorInField::ComputeStep(): " << G4endl
355 << " Difficult track - taking many sub steps." << G4endl;
356 }
357 printStatus( SubStepStartState, CurrentState, CurrentProposedStepLength,
358 NewSafety, do_loop_count, pPhysVol );
359 }
360#endif
361
362 do_loop_count++;
363
364 } while( (!intersects )
365 && (StepTaken + kCarTolerance < CurrentProposedStepLength)
366 && ( do_loop_count < fMax_loop_count ) );
367
368 if( do_loop_count >= fMax_loop_count )
369 {
370 fParticleIsLooping = true;
371
372 if ( fVerboseLevel > 0 )
373 {
374 G4cout << " G4PropagateInField::ComputeStep(): " << G4endl
375 << " Killing looping particle "
376 // << " of " << energy << " energy "
377 << " after " << do_loop_count << " field substeps "
378 << " totaling " << StepTaken / mm << " mm " ;
379 if( pPhysVol )
380 G4cout << " in volume " << pPhysVol->GetName() ;
381 else
382 G4cout << " in unknown or null volume. " ;
383 G4cout << G4endl;
384 }
385 }
386
387 if( !intersects )
388 {
389 // Chord AB or "minor chords" do not intersect
390 // B is the endpoint Step of the current Step.
391 //
392 End_PointAndTangent = CurrentState;
393 TruePathLength = StepTaken;
394 }
395
396 // Set pFieldTrack to the return value
397 //
398 pFieldTrack = End_PointAndTangent;
399
400#ifdef G4VERBOSE
401 // Check that "s" is correct
402 //
403 if( std::fabs(OriginalState.GetCurveLength() + TruePathLength
404 - End_PointAndTangent.GetCurveLength()) > 3.e-4 * TruePathLength )
405 {
406 std::ostringstream message;
407 message << "Curve length mis-match between original state "
408 << "and proposed endpoint of propagation." << G4endl
409 << " The curve length of the endpoint should be: "
410 << OriginalState.GetCurveLength() + TruePathLength << G4endl
411 << " and it is instead: "
412 << End_PointAndTangent.GetCurveLength() << "." << G4endl
413 << " A difference of: "
414 << OriginalState.GetCurveLength() + TruePathLength
415 - End_PointAndTangent.GetCurveLength() << G4endl
416 << " Original state = " << OriginalState << G4endl
417 << " Proposed state = " << End_PointAndTangent;
418 G4Exception("G4PropagatorInField::ComputeStep()",
419 "GeomNav0003", FatalException, message);
420 }
421#endif
422
423 // In particular anomalous cases, we can get repeated zero steps
424 // In order to correct this efficiently, we identify these cases
425 // and only take corrective action when they occur.
426 //
427 if( ( (TruePathLength < fZeroStepThreshold)
428 && ( TruePathLength+kCarTolerance < CurrentProposedStepLength )
429 )
430 || ( TruePathLength < 0.5*kCarTolerance )
431 )
432 {
433 fNoZeroStep++;
434 }
435 else{
436 fNoZeroStep = 0;
437 }
438
439 if( fNoZeroStep > fAbandonThreshold_NoZeroSteps )
440 {
441 fParticleIsLooping = true;
442 std::ostringstream message;
443 message << "Particle is stuck; it will be killed." << G4endl
444 << " Zero progress for " << fNoZeroStep << " attempted steps."
445 << G4endl
446 << " Proposed Step is " << CurrentProposedStepLength
447 << " but Step Taken is "<< fFull_CurveLen_of_LastAttempt << G4endl
448 << " For Particle with Charge = " << fCharge
449 << " Momentum = "<< fInitialMomentumModulus
450 << " Mass = " << fMass << G4endl;
451 if( pPhysVol )
452 message << " in volume " << pPhysVol->GetName() ;
453 else
454 message << " in unknown or null volume. " ;
455 G4Exception("G4PropagatorInField::ComputeStep()",
456 "GeomNav1002", JustWarning, message);
457 fNoZeroStep = 0;
458 }
459
460 return TruePathLength;
461}
462
463///////////////////////////////////////////////////////////////////////////
464//
465// Dumps status of propagator.
466
467void
469 const G4FieldTrack& CurrentFT,
470 G4double requestStep,
471 G4double safety,
472 G4int stepNo,
473 G4VPhysicalVolume* startVolume)
474{
475 const G4int verboseLevel=fVerboseLevel;
476 const G4ThreeVector StartPosition = StartFT.GetPosition();
477 const G4ThreeVector StartUnitVelocity = StartFT.GetMomentumDir();
478 const G4ThreeVector CurrentPosition = CurrentFT.GetPosition();
479 const G4ThreeVector CurrentUnitVelocity = CurrentFT.GetMomentumDir();
480
481 G4double step_len = CurrentFT.GetCurveLength() - StartFT.GetCurveLength();
482
483 G4int oldprec; // cout/cerr precision settings
484
485 if( ((stepNo == 0) && (verboseLevel <3)) || (verboseLevel >= 3) )
486 {
487 oldprec = G4cout.precision(4);
488 G4cout << std::setw( 6) << " "
489 << std::setw( 25) << " Current Position and Direction" << " "
490 << G4endl;
491 G4cout << std::setw( 5) << "Step#"
492 << std::setw(10) << " s " << " "
493 << std::setw(10) << "X(mm)" << " "
494 << std::setw(10) << "Y(mm)" << " "
495 << std::setw(10) << "Z(mm)" << " "
496 << std::setw( 7) << " N_x " << " "
497 << std::setw( 7) << " N_y " << " "
498 << std::setw( 7) << " N_z " << " " ;
499 G4cout << std::setw( 7) << " Delta|N|" << " "
500 << std::setw( 9) << "StepLen" << " "
501 << std::setw(12) << "StartSafety" << " "
502 << std::setw( 9) << "PhsStep" << " ";
503 if( startVolume )
504 { G4cout << std::setw(18) << "NextVolume" << " "; }
505 G4cout.precision(oldprec);
506 G4cout << G4endl;
507 }
508 if((stepNo == 0) && (verboseLevel <=3))
509 {
510 // Recurse to print the start values
511 //
512 printStatus( StartFT, StartFT, -1.0, safety, -1, startVolume);
513 }
514 if( verboseLevel <= 3 )
515 {
516 if( stepNo >= 0)
517 { G4cout << std::setw( 4) << stepNo << " "; }
518 else
519 { G4cout << std::setw( 5) << "Start" ; }
520 oldprec = G4cout.precision(8);
521 G4cout << std::setw(10) << CurrentFT.GetCurveLength() << " ";
522 G4cout.precision(8);
523 G4cout << std::setw(10) << CurrentPosition.x() << " "
524 << std::setw(10) << CurrentPosition.y() << " "
525 << std::setw(10) << CurrentPosition.z() << " ";
526 G4cout.precision(4);
527 G4cout << std::setw( 7) << CurrentUnitVelocity.x() << " "
528 << std::setw( 7) << CurrentUnitVelocity.y() << " "
529 << std::setw( 7) << CurrentUnitVelocity.z() << " ";
530 G4cout.precision(3);
531 G4cout << std::setw( 7)
532 << CurrentFT.GetMomentum().mag()-StartFT.GetMomentum().mag() << " ";
533 G4cout << std::setw( 9) << step_len << " ";
534 G4cout << std::setw(12) << safety << " ";
535 if( requestStep != -1.0 )
536 { G4cout << std::setw( 9) << requestStep << " "; }
537 else
538 { G4cout << std::setw( 9) << "Init/NotKnown" << " "; }
539 if( startVolume != 0)
540 { G4cout << std::setw(12) << startVolume->GetName() << " "; }
541 G4cout.precision(oldprec);
542 G4cout << G4endl;
543 }
544 else // if( verboseLevel > 3 )
545 {
546 // Multi-line output
547
548 G4cout << "Step taken was " << step_len
549 << " out of PhysicalStep = " << requestStep << G4endl;
550 G4cout << "Final safety is: " << safety << G4endl;
551 G4cout << "Chord length = " << (CurrentPosition-StartPosition).mag()
552 << G4endl;
553 G4cout << G4endl;
554 }
555}
556
557///////////////////////////////////////////////////////////////////////////
558//
559// Prints Step diagnostics
560
561void
563 G4double CurrentProposedStepLength,
564 G4double decreaseFactor,
565 G4double stepTrial,
566 const G4FieldTrack& )
567{
568 G4int iprec= G4cout.precision(8);
569 G4cout << " " << std::setw(12) << " PiF: NoZeroStep "
570 << " " << std::setw(20) << " CurrentProposed len "
571 << " " << std::setw(18) << " Full_curvelen_last"
572 << " " << std::setw(18) << " last proposed len "
573 << " " << std::setw(18) << " decrease factor "
574 << " " << std::setw(15) << " step trial "
575 << G4endl;
576
577 G4cout << " " << std::setw(10) << fNoZeroStep << " "
578 << " " << std::setw(20) << CurrentProposedStepLength
579 << " " << std::setw(18) << fFull_CurveLen_of_LastAttempt
580 << " " << std::setw(18) << fLast_ProposedStepLength
581 << " " << std::setw(18) << decreaseFactor
582 << " " << std::setw(15) << stepTrial
583 << G4endl;
584 G4cout.precision( iprec );
585
586}
587
588// Access the points which have passed through the filter. The
589// points are stored as ThreeVectors for the initial impelmentation
590// only (jacek 30/10/2002)
591// Responsibility for deleting the points lies with
592// SmoothTrajectoryPoint, which is the points' final
593// destination. The points pointer is set to NULL, to ensure that
594// the points are not re-used in subsequent steps, therefore THIS
595// METHOD MUST BE CALLED EXACTLY ONCE PER STEP. (jacek 08/11/2002)
596
597std::vector<G4ThreeVector>*
599{
600 // NB, GimmeThePointsAndForgetThem really forgets them, so it can
601 // only be called (exactly) once for each step.
602
603 if (fpTrajectoryFilter)
604 {
605 return fpTrajectoryFilter->GimmeThePointsAndForgetThem();
606 }
607 else
608 {
609 return 0;
610 }
611}
612
613void
615{
616 fpTrajectoryFilter = filter;
617}
618
620{
621 // Goal: Clear all memory of previous steps, cached information
622
623 fParticleIsLooping= false;
624 fNoZeroStep= 0;
625
626 End_PointAndTangent= G4FieldTrack( G4ThreeVector(0.,0.,0.),
627 G4ThreeVector(0.,0.,0.),
628 0.0,0.0,0.0,0.0,0.0);
629 fFull_CurveLen_of_LastAttempt = -1;
630 fLast_ProposedStepLength = -1;
631
632 fPreviousSftOrigin= G4ThreeVector(0.,0.,0.);
633 fPreviousSafety= 0.0;
634}
635
637FindAndSetFieldManager( G4VPhysicalVolume* pCurrentPhysicalVolume)
638{
639 G4FieldManager* currentFieldMgr;
640
641 currentFieldMgr = fDetectorFieldMgr;
642 if( pCurrentPhysicalVolume)
643 {
644 G4FieldManager *pRegionFieldMgr= 0, *localFieldMgr = 0;
645 G4LogicalVolume* pLogicalVol= pCurrentPhysicalVolume->GetLogicalVolume();
646
647 if( pLogicalVol ) {
648 // Value for Region, if any, Overrides
649 G4Region* pRegion= pLogicalVol->GetRegion();
650 if( pRegion ) {
651 pRegionFieldMgr= pRegion->GetFieldManager();
652 if( pRegionFieldMgr )
653 currentFieldMgr= pRegionFieldMgr;
654 }
655
656 // 'Local' Value from logical volume, if any, Overrides
657 localFieldMgr= pLogicalVol->GetFieldManager();
658 if ( localFieldMgr )
659 currentFieldMgr = localFieldMgr;
660 }
661 }
662 fCurrentFieldMgr= currentFieldMgr;
663
664 // Flag that field manager has been set.
665 fSetFieldMgr= true;
666
667 return currentFieldMgr;
668}
669
671{
672 G4int oldval= fVerboseLevel;
673 fVerboseLevel= level;
674
675 // Forward the verbose level 'reduced' to ChordFinder,
676 // MagIntegratorDriver ... ?
677 //
679 integrDriver->SetVerboseLevel( fVerboseLevel - 2 );
680 G4cout << "Set Driver verbosity to " << fVerboseLevel - 2 << G4endl;
681
682 return oldval;
683}
@ JustWarning
@ FatalException
CLHEP::Hep3Vector G4ThreeVector
double G4double
Definition: G4Types.hh:64
int G4int
Definition: G4Types.hh:66
bool G4bool
Definition: G4Types.hh:67
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cout
double z() const
double x() const
double y() const
double mag() const
void SetChargeMomentumMass(G4double pCharge, G4double pMomentum, G4double pMass)
G4double AdvanceChordLimited(G4FieldTrack &yCurrent, G4double stepInitial, G4double epsStep_Relative, const G4ThreeVector latestSafetyOrigin, G4double lasestSafetyRadius)
G4MagInt_Driver * GetIntegrationDriver()
G4double GetMinimumEpsilonStep() const
G4double GetMaximumEpsilonStep() const
G4double GetDeltaOneStep() const
G4double GetDeltaIntersection() const
const G4ThreeVector & GetMomentumDir() const
G4double GetCurveLength() const
G4ThreeVector GetPosition() const
G4ThreeVector GetMomentum() const
G4double GetSurfaceTolerance() const
static G4GeometryTolerance * GetInstance()
G4Region * GetRegion() const
G4FieldManager * GetFieldManager() const
void SetVerboseLevel(G4int newLevel)
virtual void LocateGlobalPointWithinVolume(const G4ThreeVector &position)
Definition: G4Navigator.cc:548
G4VPhysicalVolume * GetWorldVolume() const
void printStatus(const G4FieldTrack &startFT, const G4FieldTrack &currentFT, G4double requestStep, G4double safety, G4int step, G4VPhysicalVolume *startVolume)
G4double ComputeStep(G4FieldTrack &pFieldTrack, G4double pCurrentProposedStepLength, G4double &pNewSafety, G4VPhysicalVolume *pPhysVol=0)
G4FieldManager * FindAndSetFieldManager(G4VPhysicalVolume *pCurrentPhysVol)
G4int SetVerboseLevel(G4int verbose)
void PrintStepLengthDiagnostic(G4double currentProposedStepLength, G4double decreaseFactor, G4double stepTrial, const G4FieldTrack &aFieldTrack)
G4ChordFinder * GetChordFinder()
G4bool IntersectChord(const G4ThreeVector &StartPointA, const G4ThreeVector &EndPointB, G4double &NewSafety, G4double &LinearStepLength, G4ThreeVector &IntersectionPoint)
G4PropagatorInField(G4Navigator *theNavigator, G4FieldManager *detectorFieldMgr, G4VIntersectionLocator *vLocator=0)
std::vector< G4ThreeVector > * GimmeTrajectoryVectorAndForgetIt() const
void SetTrajectoryFilter(G4VCurvedTrajectoryFilter *filter)
void SetEpsilonStep(G4double newEps)
G4FieldManager * GetFieldManager() const
std::vector< G4ThreeVector > * GimmeThePointsAndForgetThem()
virtual void TakeIntermediatePoint(G4ThreeVector newPoint)=0
void SetDeltaIntersectionFor(G4double deltaIntersection)
void SetSafetyParametersFor(G4bool UseSafety)
void SetChordFinderFor(G4ChordFinder *fCFinder)
void SetEpsilonStepFor(G4double EpsilonStep)
G4LogicalVolume * GetLogicalVolume() const
const G4String & GetName() const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41