#include <G4CoupledTransportation.hh>

Inheritance diagram for G4CoupledTransportation:

Collaboration diagram for G4CoupledTransportation:

Public Member Functions
	G4CoupledTransportation (G4int verbosityLevel=0)

	~G4CoupledTransportation ()

G4double	AlongStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety, G4GPILSelection *selection)

G4VParticleChange *	AlongStepDoIt (const G4Track &track, const G4Step &stepData)

G4VParticleChange *	PostStepDoIt (const G4Track &track, const G4Step &stepData)

G4double	PostStepGetPhysicalInteractionLength (const G4Track &, G4double previousStepSize, G4ForceCondition *pForceCond)

G4PropagatorInField *	GetPropagatorInField ()

void	SetPropagatorInField (G4PropagatorInField *pFieldPropagator)

void	SetVerboseLevel (G4int verboseLevel)

G4int	GetVerboseLevel () const

G4double	GetThresholdWarningEnergy () const

G4double	GetThresholdImportantEnergy () const

G4int	GetThresholdTrials () const

void	SetThresholdWarningEnergy (G4double newEnWarn)

void	SetThresholdImportantEnergy (G4double newEnImp)

void	SetThresholdTrials (G4int newMaxTrials)

G4double	GetMaxEnergyKilled () const

G4double	GetSumEnergyKilled () const

void	ResetKilledStatistics (G4int report=1)

void	StartTracking (G4Track *aTrack)

void	EndTracking ()

G4double	AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)

G4VParticleChange *	AtRestDoIt (const G4Track &, const G4Step &)

Public Member Functions inherited from G4VProcess
	G4VProcess (const G4String &aName="NoName", G4ProcessType aType=fNotDefined)

	G4VProcess (const G4VProcess &right)

virtual	~G4VProcess ()

G4int	operator== (const G4VProcess &right) const

G4int	operator!= (const G4VProcess &right) const

G4double	GetCurrentInteractionLength () const

void	SetPILfactor (G4double value)

G4double	GetPILfactor () const

G4double	AlongStepGPIL (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)

G4double	AtRestGPIL (const G4Track &track, G4ForceCondition *condition)

G4double	PostStepGPIL (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)

virtual G4bool	IsApplicable (const G4ParticleDefinition &)

virtual void	BuildPhysicsTable (const G4ParticleDefinition &)

virtual void	PreparePhysicsTable (const G4ParticleDefinition &)

virtual G4bool	StorePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)

virtual G4bool	RetrievePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)

const G4String &	GetPhysicsTableFileName (const G4ParticleDefinition *, const G4String &directory, const G4String &tableName, G4bool ascii=false)

const G4String &	GetProcessName () const

G4ProcessType	GetProcessType () const

void	SetProcessType (G4ProcessType)

G4int	GetProcessSubType () const

void	SetProcessSubType (G4int)

virtual void	SetProcessManager (const G4ProcessManager *)

virtual const G4ProcessManager *	GetProcessManager ()

virtual void	ResetNumberOfInteractionLengthLeft ()

G4double	GetNumberOfInteractionLengthLeft () const

G4double	GetTotalNumberOfInteractionLengthTraversed () const

G4bool	isAtRestDoItIsEnabled () const

G4bool	isAlongStepDoItIsEnabled () const

G4bool	isPostStepDoItIsEnabled () const

virtual void	DumpInfo () const

void	SetVerboseLevel (G4int value)

G4int	GetVerboseLevel () const

virtual void	SetMasterProcess (G4VProcess *masterP)

const G4VProcess *	GetMasterProcess () const

virtual void	BuildWorkerPhysicsTable (const G4ParticleDefinition &part)

virtual void	PrepareWorkerPhysicsTable (const G4ParticleDefinition &)

Static Public Member Functions
static G4bool	EnableUseMagneticMoment (G4bool useMoment=true)

Static Public Member Functions inherited from G4VProcess
static const G4String &	GetProcessTypeName (G4ProcessType)

Protected Member Functions
G4bool	DoesGlobalFieldExist ()

void	ReportInexactEnergy (G4double startEnergy, G4double endEnergy)

void	ReportMove (G4ThreeVector OldVector, G4ThreeVector NewVector, const G4String &Quantity)

Protected Member Functions inherited from G4VProcess
void	SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)

void	ClearNumberOfInteractionLengthLeft ()

Private Attributes
G4Navigator *	fMassNavigator

G4PathFinder *	fPathFinder

G4int	fNavigatorId

G4PropagatorInField *	fFieldPropagator

G4bool	fGlobalFieldExists

G4ThreeVector	fTransportEndPosition

G4ThreeVector	fTransportEndMomentumDir

G4double	fTransportEndKineticEnergy

G4ThreeVector	fTransportEndSpin

G4bool	fMomentumChanged

G4bool	fEndGlobalTimeComputed

G4double	fCandidateEndGlobalTime

G4bool	fParticleIsLooping

G4ThreeVector	fPreviousSftOrigin

G4double	fPreviousMassSafety

G4double	fPreviousFullSafety

G4TouchableHandle	fCurrentTouchableHandle

G4bool	fMassGeometryLimitedStep

G4bool	fAnyGeometryLimitedStep

G4ParticleChangeForTransport	fParticleChange

G4double	fEndpointDistance

G4double	fThreshold_Warning_Energy

G4double	fThreshold_Important_Energy

G4int	fThresholdTrials

G4int	fNoLooperTrials

G4double	fSumEnergyKilled

G4double	fMaxEnergyKilled

G4SafetyHelper *	fpSafetyHelper

G4int	fVerboseLevel

Static Private Attributes
static G4bool	fUseMagneticMoment =false

Friends
class	G4Transportation

Additional Inherited Members
Protected Attributes inherited from G4VProcess
const G4ProcessManager *	aProcessManager

G4VParticleChange *	pParticleChange

G4ParticleChange	aParticleChange

G4double	theNumberOfInteractionLengthLeft

G4double	currentInteractionLength

G4double	theInitialNumberOfInteractionLength

G4String	theProcessName

G4String	thePhysicsTableFileName

G4ProcessType	theProcessType

G4int	theProcessSubType

G4double	thePILfactor

G4bool	enableAtRestDoIt

G4bool	enableAlongStepDoIt

G4bool	enablePostStepDoIt

G4int	verboseLevel

Detailed Description

Definition at line 64 of file G4CoupledTransportation.hh.

Constructor & Destructor Documentation

◆ G4CoupledTransportation()

G4CoupledTransportation::G4CoupledTransportation ( G4int verbosityLevel = 0 )

Definition at line 67 of file G4CoupledTransportation.cc.

   : G4VProcess( G4String("CoupledTransportation"), fTransportation ),
     fTransportEndPosition(0.0, 0.0, 0.0),
     fTransportEndMomentumDir(0.0, 0.0, 0.0),
     fTransportEndKineticEnergy(0.0), 
     fTransportEndSpin(0.0, 0.0, 0.0), // fTransportEndPolarization(0.0, 0.0, 0.0),
     fMomentumChanged(false), 
     fEndGlobalTimeComputed(false),
     fCandidateEndGlobalTime(0.0),
     fParticleIsLooping( false ),
     fPreviousSftOrigin( 0.,0.,0. ),
     fPreviousMassSafety( 0.0 ),
     fPreviousFullSafety( 0.0 ),
     fMassGeometryLimitedStep( false ), 
     fAnyGeometryLimitedStep( false ), 
     fEndpointDistance( -1.0 ), 
     fThreshold_Warning_Energy( 100 * MeV ),  
     fThreshold_Important_Energy( 250 * MeV ), 
     fThresholdTrials( 10 ), 
     fNoLooperTrials( 0 ),
     fSumEnergyKilled( 0.0 ), fMaxEnergyKilled( 0.0 ), 
     fVerboseLevel( verbosity )
 {
   // set Process Sub Type
   SetProcessSubType(static_cast<G4int>(COUPLED_TRANSPORTATION));
 
   G4TransportationManager* transportMgr ; 
 
   transportMgr = G4TransportationManager::GetTransportationManager() ; 
 
   fMassNavigator = transportMgr->GetNavigatorForTracking() ; 
   fFieldPropagator = transportMgr->GetPropagatorInField() ;
   // fGlobalFieldMgr = transportMgr->GetFieldManager() ;
   fNavigatorId= transportMgr->ActivateNavigator( fMassNavigator ); 
   if( fVerboseLevel > 0 )
   {
     G4cout << " G4CoupledTransportation constructor: ----- " << G4endl;
     G4cout << " Verbose level is " << fVerboseLevel << G4endl;
     G4cout << " Navigator Id obtained in G4CoupledTransportation constructor " 
            << fNavigatorId << G4endl;
   }
   fPathFinder=  G4PathFinder::GetInstance(); 
   fpSafetyHelper = transportMgr->GetSafetyHelper();  // New 
 
   // Following assignment is to fix small memory leak from simple use of 'new'
   static G4ThreadLocal G4TouchableHandle* pNullTouchableHandle = 0;
   if ( !pNullTouchableHandle)  { pNullTouchableHandle = new G4TouchableHandle; }
   fCurrentTouchableHandle = *pNullTouchableHandle;
     // Points to (G4VTouchable*) 0
 
   G4FieldManager  *globalFieldMgr= transportMgr->GetFieldManager();
   fGlobalFieldExists= globalFieldMgr ? globalFieldMgr->GetDetectorField() : 0 ; 
 }

Here is the call graph for this function:

◆ ~G4CoupledTransportation()

G4CoupledTransportation::~G4CoupledTransportation ( )

Definition at line 123 of file G4CoupledTransportation.cc.

 {
   // fCurrentTouchableHandle is a data member - no deletion required
 
   if( (fVerboseLevel > 0) || (fSumEnergyKilled > 0.0 ) )
   { 
     G4cout << " G4CoupledTransportation: Statistics for looping particles " << G4endl;
     G4cout << "   Sum of energy of loopers killed: " <<  fSumEnergyKilled << G4endl;
     G4cout << "   Max energy of loopers killed: " <<  fMaxEnergyKilled << G4endl;
   } 
 }

Here is the call graph for this function:

Member Function Documentation

◆ AlongStepDoIt()

G4VParticleChange * G4CoupledTransportation::AlongStepDoIt	(	const G4Track &	track,
		const G4Step &	stepData
	)

virtual

Implements G4VProcess.

Definition at line 557 of file G4CoupledTransportation.cc.

 {
   static G4ThreadLocal G4int noCalls=0;
   noCalls++;
 
   fParticleChange.Initialize(track) ;
       // sets all its members to the value of corresponding members in G4Track
 
   //  Code specific for Transport
   //
   fParticleChange.ProposePosition(fTransportEndPosition) ;
   // G4cout << " G4CoupledTransportation::AlongStepDoIt" 
   //     << " proposes position = " << fTransportEndPosition  
   //     << " and end momentum direction  = " << fTransportEndMomentumDir <<  G4endl;
   fParticleChange.ProposeMomentumDirection(fTransportEndMomentumDir) ;
   fParticleChange.ProposeEnergy(fTransportEndKineticEnergy) ;
   fParticleChange.SetMomentumChanged(fMomentumChanged) ;
 
   fParticleChange.ProposePolarization(fTransportEndSpin);
   
   G4double deltaTime = 0.0 ;
 
   // Calculate  Lab Time of Flight (ONLY if field Equations used it!)
      // G4double endTime   = fCandidateEndGlobalTime;
      // G4double delta_time = endTime - startTime;
 
   G4double startTime = track.GetGlobalTime() ;
   
   if (!fEndGlobalTimeComputed)
   {
      G4double finalInverseVel= DBL_MAX, initialInverseVel=DBL_MAX; 
 
      // The time was not integrated .. make the best estimate possible
      //
      G4double finalVelocity   = track.GetVelocity() ;
      if( finalVelocity > 0.0 ) { finalInverseVel= 1.0 / finalVelocity; }
      G4double initialVelocity = stepData.GetPreStepPoint()->GetVelocity() ;
      if( initialVelocity > 0.0 ) { initialInverseVel= 1.0 / initialVelocity; }
      G4double stepLength      = track.GetStepLength() ;
 
      if (finalVelocity > 0.0)
      {
         // deltaTime = stepLength/finalVelocity ;
         G4double meanInverseVelocity = 0.5 * ( initialInverseVel + finalInverseVel );
         deltaTime = stepLength * meanInverseVelocity ;
         // G4cout << " dt = s * mean(1/v) , with " << "  s = " << stepLength
         //     << "  mean(1/v)= "  << meanInverseVelocity << G4endl;
      }
      else
      {
         deltaTime = stepLength * initialInverseVel ;
         // G4cout << " dt = s / initV "  << "  s = "   << stepLength
         //        << " 1 / initV= " << initialInverseVel << G4endl; 
      }  //  Could do with better estimate for final step (finalVelocity = 0) ?
 
      fCandidateEndGlobalTime   = startTime + deltaTime ;
      fParticleChange.ProposeLocalTime(  track.GetLocalTime() + deltaTime) ;
 
      // G4cout << " Calculated global time from start = " << startTime << " and "
      //        << " delta time = " << deltaTime << G4endl;
   }
   else
   {
      deltaTime = fCandidateEndGlobalTime - startTime ;
      fParticleChange.ProposeGlobalTime( fCandidateEndGlobalTime ) ;
      // G4cout << " Calculated global time from candidate end time = "
      //    << fCandidateEndGlobalTime << " and start time = " << startTime << G4endl;
   }
 
   // G4cout << " G4CoupledTransportation::AlongStepDoIt  "
   // << " flag whether computed time = " << fEndGlobalTimeComputed  << " and " 
   // << " is proposes end time " << fCandidateEndGlobalTime << G4endl; 
 
   // Now Correct by Lorentz factor to get "proper" deltaTime
   
   G4double  restMass       = track.GetDynamicParticle()->GetMass() ;
   G4double deltaProperTime = deltaTime*( restMass/track.GetTotalEnergy() ) ;
 
   fParticleChange.ProposeProperTime(track.GetProperTime() + deltaProperTime) ;
   //fParticleChange. ProposeTrueStepLength( track.GetStepLength() ) ;
 
   // If the particle is caught looping or is stuck (in very difficult
   // boundaries) in a magnetic field (doing many steps) 
   //   THEN this kills it ...
   //
   if ( fParticleIsLooping )
   {
      G4double endEnergy= fTransportEndKineticEnergy;
 
      if( (endEnergy < fThreshold_Important_Energy) 
           || (fNoLooperTrials >= fThresholdTrials ) )
      {
         // Kill the looping particle 
         //
         fParticleChange.ProposeTrackStatus( fStopAndKill )  ;
 
         // 'Bare' statistics
         fSumEnergyKilled += endEnergy; 
         if( endEnergy > fMaxEnergyKilled) { fMaxEnergyKilled= endEnergy; }
 
 #ifdef G4VERBOSE
         if((fVerboseLevel > 1) && ( endEnergy > fThreshold_Warning_Energy ))
         { 
           G4cout << " G4CoupledTransportation is killing track that is looping or stuck " << G4endl
                  << "   This track has " << track.GetKineticEnergy() / MeV
                  << " MeV energy." << G4endl;
         }
         if( fVerboseLevel > 0 )
         { 
           G4cout << "   Steps by this track: " << track.GetCurrentStepNumber() << G4endl;
         }
 #endif
         fNoLooperTrials=0; 
       }
       else
       { 
         fNoLooperTrials ++; 
 #ifdef G4VERBOSE
         if( (fVerboseLevel > 2) )
         {
           G4cout << "  ** G4CoupledTransportation::AlongStepDoIt(): Particle looping -  " << G4endl
                  << "   Number of consecutive problem step (this track) = " << fNoLooperTrials << G4endl
                  << "   Steps by this track: " << track.GetCurrentStepNumber() << G4endl
                  << "   Total no of calls to this method (all tracks) = " << noCalls << G4endl;
         }
 #endif
       }
   }
   else
   { 
       fNoLooperTrials=0; 
   }
 
   // Another (sometimes better way) is to use a user-limit maximum Step size
   // to alleviate this problem .. 
 
   // Add smooth curved trajectories to particle-change
   //
   // fParticleChange.SetPointerToVectorOfAuxiliaryPoints
   //   (fFieldPropagator->GimmeTrajectoryVectorAndForgetIt() );
 
   return &fParticleChange ;
 }

Here is the call graph for this function:

◆ AlongStepGetPhysicalInteractionLength()

G4double G4CoupledTransportation::AlongStepGetPhysicalInteractionLength	(	const G4Track &	track,
		G4double	previousStepSize,
		G4double	currentMinimumStep,
		G4double &	currentSafety,
		G4GPILSelection *	selection
	)

virtual

Implements G4VProcess.

Definition at line 143 of file G4CoupledTransportation.cc.

 {
   G4double geometryStepLength; 
   G4double startMassSafety= 0.0;   //  estimated safety for start point (mass geometry)
   G4double startFullSafety= 0.0;   //  estimated safety for start point (all geometries)
   G4double safetyProposal= -1.0;   //  local copy of proposal 
 
   G4ThreeVector  EndUnitMomentum ;
   G4double       lengthAlongCurve=0.0 ;
  
   fParticleIsLooping = false ;
 
   // Initial actions moved to  StartTrack()   
   // --------------------------------------
   // Note: in case another process changes touchable handle
   //    it will be necessary to add here (for all steps)   
   // fCurrentTouchableHandle = aTrack->GetTouchableHandle();
 
   // GPILSelection is set to defaule value of CandidateForSelection
   // It is a return value
   //
   *selection = CandidateForSelection ;
 
   // Get initial Energy/Momentum of the track
   //
   const G4DynamicParticle*    pParticle  = track.GetDynamicParticle() ;
   const G4ParticleDefinition* pParticleDef   = pParticle->GetDefinition() ;
   G4ThreeVector startMomentumDir       = pParticle->GetMomentumDirection() ;
   G4ThreeVector startPosition          = track.GetPosition() ;
   G4VPhysicalVolume* currentVolume= track.GetVolume(); 
 
 #ifdef G4DEBUG_TRANSPORT
   if( fVerboseLevel > 1 )
   {
     G4cout << "G4CoupledTransportation::AlongStepGPIL> called in volume " 
            << currentVolume->GetName() << G4endl; 
   }
 #endif
   // G4double   theTime        = track.GetGlobalTime() ;
 
   // The Step Point safety can be limited by other geometries and/or the 
   // assumptions of any process - it's not always the geometrical safety.
   // We calculate the starting point's isotropic safety here.
   //
   G4ThreeVector OriginShift = startPosition - fPreviousSftOrigin ;
   G4double      MagSqShift  = OriginShift.mag2() ;
   startMassSafety = 0.0; 
   startFullSafety= 0.0; 
 
   //  Recall that FullSafety <= MassSafety 
   // Original: if( MagSqShift < sqr(fPreviousMassSafety) ) {
   if( MagSqShift < sqr(fPreviousFullSafety) )   // Revision proposed by Alex H, 2 Oct 07
   {
      G4double mag_shift= std::sqrt(MagSqShift); 
      startMassSafety = std::max( (fPreviousMassSafety - mag_shift), 0.0); 
      startFullSafety = std::max( (fPreviousFullSafety - mag_shift), 0.0);
        // Need to be consistent between full safety with Mass safety
        //   in order reproduce results in simple case  --> use same calculation method
 
      // Only compute full safety if massSafety > 0.  Else it remains 0
      //   startFullSafety = fPathFinder->ComputeSafety( startPosition ); 
   }
 
   // Is the particle charged or has it a magnetic moment?
   //
   G4double particleCharge = pParticle->GetCharge() ;
   G4double magneticMoment = pParticle->GetMagneticMoment() ;
   G4double       restMass = pParticle->GetMass() ; 
 
   fMassGeometryLimitedStep = false ; //  Set default - alex
   fAnyGeometryLimitedStep = false; 
 
   // fEndGlobalTimeComputed = false ;
 
   // There is no need to locate the current volume. It is Done elsewhere:
   //   On track construction 
   //   By the tracking, after all AlongStepDoIts, in "Relocation"
 
   // Check if the particle has a force, EM or gravitational, exerted on it
   //
   G4FieldManager* fieldMgr=0;
   G4bool          fieldExertsForce = false ;
 
   G4bool gravityOn = false;
   const G4Field* ptrField= 0;
 
   fieldMgr = fFieldPropagator->FindAndSetFieldManager( track.GetVolume() );
   if( fieldMgr != 0 )
   {
      // Message the field Manager, to configure it for this track
      fieldMgr->ConfigureForTrack( &track );
      // Here it can transition from a null field-ptr to a finite field 
 
      // If the field manager has no field ptr, the field is zero 
      //     by definition ( = there is no field ! )
      ptrField= fieldMgr->GetDetectorField();
  
      if( ptrField != 0)
      { 
         gravityOn= ptrField->IsGravityActive();
         if(  (particleCharge != 0.0) 
              || (fUseMagneticMoment && (magneticMoment != 0.0) )
              || (gravityOn && (restMass != 0.0))
           )
         {
            fieldExertsForce = true;
         }
      }
   }
   G4double momentumMagnitude = pParticle->GetTotalMomentum() ;
 
   if( fieldExertsForce )
   {
      G4EquationOfMotion*    equationOfMotion = 0; 
 
      // equationOfMotion = 
      //     (fFieldPropagator->GetChordFinder()->GetIntegrationDriver()->GetStepper())
      //      ->GetEquationOfMotion();
 
      // Consolidate into auxiliary method G4EquationOfMotion* GetEquationOfMotion() 
      //  equationOfMotion= fFieldPropagator->GetCurrentEquationOfMotion();
      G4MagIntegratorStepper*  pStepper= 0;
 
      G4ChordFinder*    pChordFinder= fFieldPropagator->GetChordFinder();
      if( pChordFinder ) 
      {
         G4MagInt_Driver*  pIntDriver= 0; 
 
         pIntDriver= pChordFinder->GetIntegrationDriver(); 
         if( pIntDriver )
         {
            pStepper= pIntDriver->GetStepper(); 
         }
         if( pStepper )
         {
            equationOfMotion= pStepper->GetEquationOfMotion();
         }
      }
      // End of proto GetEquationOfMotion()
 
      G4ChargeState chargeState(particleCharge,             // The charge can change (dynamic)
                                magneticMoment,
                                pParticleDef->GetPDGSpin() ); 
      // For insurance, could set it again
      // chargeState.SetPDGSpin( pParticleDef->GetPDGSpin() );   // Newly/provisionally in same object
 
      if( equationOfMotion )
      {
         equationOfMotion->SetChargeMomentumMass( chargeState,
                                                  momentumMagnitude,
                                                  restMass );
      }else{
         G4cerr << " ERROR > Cannot find valid Equation of motion - Unable to pass Charge, Momentum and Mass "  << G4endl;
      }
   }
 
   G4ThreeVector polarizationVec  = track.GetPolarization() ;
   G4FieldTrack  aFieldTrack = G4FieldTrack( startPosition, 
                                             track.GetGlobalTime(), // Lab.
                                             // track.GetProperTime(), // Particle rest frame
                                             track.GetMomentumDirection(),
                                             track.GetKineticEnergy(),
                                             restMass,
                                             particleCharge, 
                                             polarizationVec, 
                                             pParticleDef->GetPDGMagneticMoment(),
                                             0.0,                    // Length along track
                                             pParticleDef->GetPDGSpin()
      ) ;
   G4int stepNo= track.GetCurrentStepNumber(); 
 
   ELimited limitedStep; 
   G4FieldTrack endTrackState('a');  //  Default values
 
   fMassGeometryLimitedStep = false ;    //  default 
   fAnyGeometryLimitedStep  = false ;
   if( currentMinimumStep > 0 )
   {
       G4double newMassSafety= 0.0;     //  temp. for recalculation
 
       // Do the Transport in the field (non recti-linear)
       //
       lengthAlongCurve = fPathFinder->ComputeStep( aFieldTrack,
                                                    currentMinimumStep, 
                                                    fNavigatorId,
                                                    stepNo,
                                                    newMassSafety,
                                                    limitedStep,
                                                    endTrackState,
                                                    currentVolume ) ;
       // G4cout << " PathFinder ComputeStep returns " << lengthAlongCurve << G4endl; 
 
       G4double newFullSafety= fPathFinder->GetCurrentSafety();  
                // this was estimated already in step above
       // G4double newFullStep= fPathFinder->GetMinimumStep(); 
 
       if( limitedStep == kUnique || limitedStep == kSharedTransport )
       {
          fMassGeometryLimitedStep = true ;
       }
 
       fAnyGeometryLimitedStep = (fPathFinder->GetNumberGeometriesLimitingStep() != 0); 
 
 //#ifdef G4DEBUG_TRANSPORT
       if( fMassGeometryLimitedStep && !fAnyGeometryLimitedStep )
       {
          G4cerr << " Error in determining geometries limiting the step" << G4endl;
          G4cerr << "  Limiting:  mass=" << fMassGeometryLimitedStep
                 << " any= " << fAnyGeometryLimitedStep << G4endl;
          G4Exception("G4CoupledTransportation::AlongStepGetPhysicalInteractionLength()", 
                      "PathFinderConfused", FatalException, 
                      "Incompatible conditions - was limited by a geometry?");
       }
 //#endif
 
       // Other potential 
       // fAnyGeometryLimitedStep = newFullStep < currentMinimumStep; 
       //                                      ^^^ Not good enough; 
       //          Must compare with maximum requested step size
       //           (eg in case another process requested bigger, got this!)
 
       geometryStepLength = std::min( lengthAlongCurve, currentMinimumStep); 
 
       // Momentum:  Magnitude and direction can be changed too now ...
       //
       fMomentumChanged         = true ; 
       fTransportEndMomentumDir = endTrackState.GetMomentumDir() ;
 
       // Remember last safety origin & value.
       fPreviousSftOrigin  = startPosition ;
       fPreviousMassSafety = newMassSafety ;         
       fPreviousFullSafety = newFullSafety ; 
       // fpSafetyHelper->SetCurrentSafety( newFullSafety, startPosition);
 
 #ifdef G4DEBUG_TRANSPORT
       if( fVerboseLevel > 1 )
       {
         G4cout << "G4Transport:CompStep> " 
                << " called the pathfinder for a new step at " << startPosition
                << " and obtained step = " << lengthAlongCurve << G4endl;
         G4cout << "  New safety (preStep) = " << newMassSafety 
                << " versus precalculated = "  << startMassSafety << G4endl; 
       }
 #endif
 
       // Store as best estimate value
       startMassSafety    = newMassSafety ; 
       startFullSafety    = newFullSafety ; 
 
       // Get the End-Position and End-Momentum (Dir-ection)
       fTransportEndPosition = endTrackState.GetPosition() ;
       fTransportEndKineticEnergy  = endTrackState.GetKineticEnergy() ; 
   }
   else
   {
       geometryStepLength   = lengthAlongCurve= 0.0 ;
       fMomentumChanged         = false ; 
       // fMassGeometryLimitedStep = false ;   //  --- ???
       // fAnyGeometryLimitedStep = true;
       fTransportEndMomentumDir = track.GetMomentumDirection();
       fTransportEndKineticEnergy  = track.GetKineticEnergy();
 
       fTransportEndPosition = startPosition;
 
       endTrackState= aFieldTrack;  // Ensures that time is updated
 
       // If the step length requested is 0, and we are on a boundary
       //   then a boundary will also limit the step.
       if( startMassSafety == 0.0 )
       {
          fMassGeometryLimitedStep = true ;
          fAnyGeometryLimitedStep = true;
       }
       //   TODO:  Add explicit logical status for being at a boundary
   }
   // G4FieldTrack aTrackState(endTrackState);  
 
   if( !fieldExertsForce ) 
   { 
       fParticleIsLooping         = false ; 
       fMomentumChanged           = false ; 
       fEndGlobalTimeComputed     = false ; 
       // G4cout << " global time is false " << G4endl; 
   } 
   else 
   { 
   
 #ifdef G4DEBUG_TRANSPORT
       if( fVerboseLevel > 1 )
       {
         G4cout << " G4CT::CS End Position = "  << fTransportEndPosition << G4endl; 
         G4cout << " G4CT::CS End Direction = " << fTransportEndMomentumDir << G4endl; 
       }
 #endif
       if( fFieldPropagator->GetCurrentFieldManager()->DoesFieldChangeEnergy() )
       {
           // If the field can change energy, then the time must be integrated
           //    - so this should have been updated
           //
           fCandidateEndGlobalTime   = endTrackState.GetLabTimeOfFlight();
           fEndGlobalTimeComputed    = true;
   
           // was ( fCandidateEndGlobalTime != track.GetGlobalTime() );
           // a cleaner way is to have FieldTrack knowing whether time is updated.
       }
       else
       {
           // The energy should be unchanged by field transport,
           //    - so the time changed will be calculated elsewhere
           //
           fEndGlobalTimeComputed = false;
   
           // Check that the integration preserved the energy 
           //     -  and if not correct this!
           G4double  startEnergy= track.GetKineticEnergy();
           G4double  endEnergy= fTransportEndKineticEnergy; 
       
           static G4ThreadLocal G4int no_inexact_steps=0; // , no_large_ediff;
           G4double absEdiff = std::fabs(startEnergy- endEnergy);
           if( absEdiff > perMillion * endEnergy )
           {
             no_inexact_steps++;
             // Possible statistics keeping here ...
           }
 #ifdef G4VERBOSE
           if( (fVerboseLevel > 1) && ( absEdiff > perThousand * endEnergy) )
           {
             ReportInexactEnergy(startEnergy, endEnergy); 
           }  // end of if (fVerboseLevel)
 #endif
           // Correct the energy for fields that conserve it
           //  This - hides the integration error
           //       - but gives a better physical answer
           fTransportEndKineticEnergy= track.GetKineticEnergy(); 
       }
   }
 
   fEndpointDistance   = (fTransportEndPosition - startPosition).mag() ;
   fParticleIsLooping = fFieldPropagator->IsParticleLooping() ;
 
   fTransportEndSpin = endTrackState.GetSpin();
 
   // Calculate the safety 
   safetyProposal= startFullSafety;   // used to be startMassSafety
      // Changed to accomodate processes that cannot update the safety -- JA 22 Nov 06
 
   // Update safety for the end-point, if becomes negative at the end-point.
 
   if(   (startFullSafety < fEndpointDistance ) 
         && ( particleCharge != 0.0 ) )        //  Only needed to prepare for Mult Scat.
    //   && !fAnyGeometryLimitedStep )          // To-Try:  No safety update if at a boundary
   {
       G4double endFullSafety =
         fPathFinder->ComputeSafety( fTransportEndPosition); 
         // Expected mission -- only mass geometry's safety
         //   fMassNavigator->ComputeSafety( fTransportEndPosition) ;
         // Yet discrete processes only have poststep -- and this cannot 
         //   currently revise the safety  
         //   ==> so we use the all-geometry safety as a precaution
 
       fpSafetyHelper->SetCurrentSafety( endFullSafety, fTransportEndPosition);
         // Pushing safety to Helper avoids recalculation at this point
 
       G4ThreeVector centerPt= G4ThreeVector(0.0, 0.0, 0.0);  // Used for return value
       G4double endMassSafety= fPathFinder->ObtainSafety( fNavigatorId, centerPt); 
         //  Retrieves the mass value from PathFinder (it calculated it)
 
       fPreviousMassSafety = endMassSafety ; 
       fPreviousFullSafety = endFullSafety; 
       fPreviousSftOrigin = fTransportEndPosition ;
 
       // The convention (Stepping Manager's) is safety from the start point
       //
       safetyProposal = endFullSafety + fEndpointDistance;
           //  --> was endMassSafety
       // Changed to accomodate processes that cannot update the safety -- JA 22 Nov 06
 
       // #define G4DEBUG_TRANSPORT 1
 
 #ifdef G4DEBUG_TRANSPORT 
       G4int prec= G4cout.precision(12) ;
       G4cout << "***Transportation::AlongStepGPIL ** " << G4endl  ;
       G4cout << "  Revised Safety at endpoint "  << fTransportEndPosition
              << "   give safety values: Mass= " << endMassSafety 
              << "  All= " << endFullSafety << G4endl ; 
       G4cout << "  Adding endpoint distance " << fEndpointDistance 
              << "   to obtain pseudo-safety= " << safetyProposal << G4endl ; 
       G4cout.precision(prec); 
   }  
   else
   {
       G4int prec= G4cout.precision(12) ;
       G4cout << "***Transportation::AlongStepGPIL ** " << G4endl  ;
       G4cout << "  Quick Safety estimate at endpoint "  << fTransportEndPosition
              << "   gives safety endpoint value = " << startFullSafety - fEndpointDistance
              << "  using start-point value " << startFullSafety 
              << "  and endpointDistance " << fEndpointDistance << G4endl; 
       G4cout.precision(prec); 
 #endif
   }          
 
   proposedSafetyForStart= safetyProposal; 
   fParticleChange.ProposeTrueStepLength(geometryStepLength) ;
 
   return geometryStepLength ;
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ AtRestDoIt()

G4VParticleChange* G4CoupledTransportation::AtRestDoIt	(	const G4Track &	,
		const G4Step &
	)

inlinevirtual

Implements G4VProcess.

Definition at line 141 of file G4CoupledTransportation.hh.

144 {return 0;};

Here is the call graph for this function:

◆ AtRestGetPhysicalInteractionLength()

G4double G4CoupledTransportation::AtRestGetPhysicalInteractionLength	(	const G4Track &	,
		G4ForceCondition *
	)

inlinevirtual

Implements G4VProcess.

Definition at line 135 of file G4CoupledTransportation.hh.

138 { return -1.0; };

◆ DoesGlobalFieldExist()

G4bool G4CoupledTransportation::DoesGlobalFieldExist ( )

protected

Here is the caller graph for this function:

◆ EnableUseMagneticMoment()

G4bool G4CoupledTransportation::EnableUseMagneticMoment ( G4bool useMoment = true )

static

Definition at line 1006 of file G4CoupledTransportation.cc.

 {
   G4bool lastValue= fUseMagneticMoment;
   fUseMagneticMoment= useMoment;
   G4Transportation::fUseMagneticMoment= useMoment;
   return lastValue;
 }

◆ EndTracking()

void G4CoupledTransportation::EndTracking ( )

virtual

Reimplemented from G4VProcess.

Definition at line 959 of file G4CoupledTransportation.cc.

 {
   G4TransportationManager::GetTransportationManager()->InactivateAll();
   fPathFinder->EndTrack();   //  Resets TransportationManager to use ordinary Navigator
 }

Here is the call graph for this function:

◆ GetMaxEnergyKilled()

G4double G4CoupledTransportation::GetMaxEnergyKilled ( ) const

inline

◆ GetPropagatorInField()

G4PropagatorInField* G4CoupledTransportation::GetPropagatorInField ( )

◆ GetSumEnergyKilled()

G4double G4CoupledTransportation::GetSumEnergyKilled ( ) const

inline

◆ GetThresholdImportantEnergy()

G4double G4CoupledTransportation::GetThresholdImportantEnergy ( ) const

inline

◆ GetThresholdTrials()

G4int G4CoupledTransportation::GetThresholdTrials ( ) const

inline

◆ GetThresholdWarningEnergy()

G4double G4CoupledTransportation::GetThresholdWarningEnergy ( ) const

inline

◆ GetVerboseLevel()

G4int G4CoupledTransportation::GetVerboseLevel ( ) const

inline

◆ PostStepDoIt()

G4VParticleChange * G4CoupledTransportation::PostStepDoIt	(	const G4Track &	track,
		const G4Step &	stepData
	)

virtual

Implements G4VProcess.

Definition at line 735 of file G4CoupledTransportation.cc.

 {
   G4TouchableHandle retCurrentTouchable ;   // The one to return
 
   // Initialize ParticleChange  (by setting all its members equal
   //                             to corresponding members in G4Track)
   // fParticleChange.Initialize(track) ;  // To initialise TouchableChange
 
   fParticleChange.ProposeTrackStatus(track.GetTrackStatus()) ;
 
   // Check that the end position and direction are preserved 
   // since call to AlongStepDoIt
 
 #ifdef G4DEBUG_TRANSPORT
   if( ( fVerboseLevel > 0 )
      && ((fTransportEndPosition - track.GetPosition()).mag2() >= 1.0e-16) )
   {
      ReportMove( track.GetPosition(), fTransportEndPosition, "End of Step Position" ); 
      G4cerr << " Problem in G4CoupledTransportation::PostStepDoIt " << G4endl; 
   }
 
   // If the Step was determined by the volume boundary, relocate the particle
   // The pathFinder will know that the geometry limited the step (!?)
 
   if( fVerboseLevel > 0 )
   {
      G4cout << " Calling PathFinder::Locate() from " 
             << " G4CoupledTransportation::PostStepDoIt " << G4endl;
      G4cout << "  fAnyGeometryLimitedStep is " << fAnyGeometryLimitedStep << G4endl;
 
   }
 #endif
 
   if(fAnyGeometryLimitedStep)
   {  
     fPathFinder->Locate( track.GetPosition(), 
                        track.GetMomentumDirection(),
                        true); 
 
     // fCurrentTouchable will now become the previous touchable, 
     // and what was the previous will be freed.
     // (Needed because the preStepPoint can point to the previous touchable)
 
     fCurrentTouchableHandle= 
       fPathFinder->CreateTouchableHandle( fNavigatorId );
 
 #ifdef G4DEBUG_TRANSPORT
     if( fVerboseLevel > 0 )
     {
       G4cout << "G4CoupledTransportation::PostStepDoIt --- fNavigatorId = " 
              << fNavigatorId << G4endl;
     }
     if( fVerboseLevel > 1 )
     {
        G4VPhysicalVolume* vol= fCurrentTouchableHandle->GetVolume(); 
        G4cout << "CHECK !!!!!!!!!!! fCurrentTouchableHandle->GetVolume() = " << vol;
        if( vol ) { G4cout << "Name=" << vol->GetName(); }
        G4cout << G4endl;
     }
 #endif
 
     // Check whether the particle is out of the world volume 
     // If so it has exited and must be killed.
     //
     if( fCurrentTouchableHandle->GetVolume() == 0 )
     {
        fParticleChange.ProposeTrackStatus( fStopAndKill ) ;
     }
     retCurrentTouchable = fCurrentTouchableHandle ;
     // fParticleChange.SetTouchableHandle( fCurrentTouchableHandle ) ;
 
     // Notify particle change that this is last step in volume
     fParticleChange.ProposeLastStepInVolume(true);
   }
   else                 // fAnyGeometryLimitedStep  is false
   { 
 #ifdef G4DEBUG_TRANSPORT
     if( fVerboseLevel > 1 )
     {
        G4cout << "G4CoupledTransportation::PostStepDoIt -- "
               << " fAnyGeometryLimitedStep  = " << fAnyGeometryLimitedStep  
               << " must be false " << G4endl;
     }
 #endif
     // This serves only to move each of the Navigator's location
     //
     // fLinearNavigator->LocateGlobalPointWithinVolume( track.GetPosition() ) ;
 
     // G4cout << "G4CoupledTransportation calling PathFinder::ReLocate() " << G4endl;
     fPathFinder->ReLocate( track.GetPosition() );
                            // track.GetMomentumDirection() ); 
 
     // Keep the value of the track's current Touchable is retained,
     //  and use it to overwrite the (unset) one in particle change.
     // Expect this must be fCurrentTouchable too
     //   - could it be different, eg at the start of a step ?
     //
     retCurrentTouchable = track.GetTouchableHandle() ;
     // fParticleChange.SetTouchableHandle( track.GetTouchableHandle() ) ;
 
     // Have not reached a boundary
     fParticleChange.ProposeLastStepInVolume(false);
   }         // endif ( fAnyGeometryLimitedStep ) 
 
   const G4VPhysicalVolume* pNewVol = retCurrentTouchable->GetVolume() ;
   const G4Material* pNewMaterial   = 0 ;
   const G4VSensitiveDetector* pNewSensitiveDetector   = 0 ;
                                                                                        
   if( pNewVol != 0 )
   {
     pNewMaterial= pNewVol->GetLogicalVolume()->GetMaterial();
     pNewSensitiveDetector= pNewVol->GetLogicalVolume()->GetSensitiveDetector();
   }
 
   // ( const_cast<G4Material *> pNewMaterial ) ;
   // ( const_cast<G4VSensitiveDetetor *> pNewSensitiveDetector) ;
 
   fParticleChange.SetMaterialInTouchable( (G4Material *) pNewMaterial ) ;
   fParticleChange.SetSensitiveDetectorInTouchable( (G4VSensitiveDetector *) pNewSensitiveDetector ) ;
              // "temporarily" until Get/Set Material of ParticleChange, 
              // and StepPoint can be made const. 
 
   const G4MaterialCutsCouple* pNewMaterialCutsCouple = 0;
   if( pNewVol != 0 )
   {
     pNewMaterialCutsCouple=pNewVol->GetLogicalVolume()->GetMaterialCutsCouple();
     if( pNewMaterialCutsCouple!=0 
         && pNewMaterialCutsCouple->GetMaterial()!=pNewMaterial )
       {
         // for parametrized volume
         //
         pNewMaterialCutsCouple =
           G4ProductionCutsTable::GetProductionCutsTable()
                        ->GetMaterialCutsCouple(pNewMaterial,
                                                pNewMaterialCutsCouple->GetProductionCuts());
       }
   }
   fParticleChange.SetMaterialCutsCoupleInTouchable( pNewMaterialCutsCouple );
 
   // Must always set the touchable in ParticleChange, whether relocated or not
   fParticleChange.SetTouchableHandle(retCurrentTouchable) ;
 
   return &fParticleChange ;
 }

Here is the call graph for this function:

◆ PostStepGetPhysicalInteractionLength()

G4double G4CoupledTransportation::PostStepGetPhysicalInteractionLength	(	const G4Track &	,
		G4double	previousStepSize,
		G4ForceCondition *	pForceCond
	)

virtual

Implements G4VProcess.

Definition at line 709 of file G4CoupledTransportation.cc.

 { 
   // Must act as PostStep action -- to relocate particle
   *pForceCond = Forced ;    
   return DBL_MAX ;
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ ReportInexactEnergy()

void G4CoupledTransportation::ReportInexactEnergy	(	G4double	startEnergy,
		G4double	endEnergy
	)

protected

Definition at line 967 of file G4CoupledTransportation.cc.

 {
   static G4ThreadLocal G4int no_warnings= 0, warnModulo=1,  moduloFactor= 10, no_large_ediff= 0; 
 
   if( std::fabs(startEnergy- endEnergy) > perThousand * endEnergy )
   {
     no_large_ediff ++;
     if( (no_large_ediff% warnModulo) == 0 )
     {
       no_warnings++;
       G4cout << "WARNING - G4CoupledTransportation::AlongStepGetPIL() " 
              << "   Energy change in Step is above 1^-3 relative value. " << G4endl
              << "   Relative change in 'tracking' step = " 
              << std::setw(15) << (endEnergy-startEnergy)/startEnergy << G4endl
              << "     Starting E= " << std::setw(12) << startEnergy / MeV << " MeV " << G4endl
              << "     Ending   E= " << std::setw(12) << endEnergy   / MeV << " MeV " << G4endl;       
       G4cout << " Energy has been corrected -- however, review"
              << " field propagation parameters for accuracy."  << G4endl;
       if( (fVerboseLevel > 2 ) || (no_warnings<4) || (no_large_ediff == warnModulo * moduloFactor) )
       {
         G4cout << " These include EpsilonStepMax(/Min) in G4FieldManager "
                << " which determine fractional error per step for integrated quantities. " << G4endl
                << " Note also the influence of the permitted number of integration steps."
                << G4endl;
       }
       G4cerr << "ERROR - G4CoupledTransportation::AlongStepGetPIL()" << G4endl
              << "        Bad 'endpoint'. Energy change detected"
              << " and corrected. " 
              << " Has occurred already "
              << no_large_ediff << " times." << G4endl;
       if( no_large_ediff == warnModulo * moduloFactor )
       {
         warnModulo *= moduloFactor;
       }
     }
   }
 }

Here is the caller graph for this function:

◆ ReportMove()

void G4CoupledTransportation::ReportMove	(	G4ThreeVector	OldVector,
		G4ThreeVector	NewVector,
		const G4String &	Quantity
	)

protected

Definition at line 719 of file G4CoupledTransportation.cc.

 {
     G4ThreeVector moveVec = ( NewVector - OldVector );
 
     G4cerr << G4endl
            << "**************************************************************" << G4endl;
     G4cerr << "Endpoint has moved between value expected from TransportEndPosition "
            << " and value from Track in PostStepDoIt. " << G4endl
            << "Change of " << Quantity << " is " << moveVec.mag() / mm << " mm long, "
            << " and its vector is " << (1.0/mm) * moveVec << " mm " << G4endl
            << "Endpoint of ComputeStep was " << OldVector
            << " and current position to locate is " << NewVector << G4endl;
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ ResetKilledStatistics()

void G4CoupledTransportation::ResetKilledStatistics ( G4int report = 1 )

inline

◆ SetPropagatorInField()

void G4CoupledTransportation::SetPropagatorInField ( G4PropagatorInField * pFieldPropagator )

◆ SetThresholdImportantEnergy()

void G4CoupledTransportation::SetThresholdImportantEnergy ( G4double newEnImp )

inline

◆ SetThresholdTrials()

void G4CoupledTransportation::SetThresholdTrials ( G4int newMaxTrials )

inline

◆ SetThresholdWarningEnergy()

void G4CoupledTransportation::SetThresholdWarningEnergy ( G4double newEnWarn )

inline

◆ SetVerboseLevel()

void G4CoupledTransportation::SetVerboseLevel ( G4int verboseLevel )

inline

◆ StartTracking()

void G4CoupledTransportation::StartTracking ( G4Track * aTrack )

virtual

Reimplemented from G4VProcess.

Definition at line 887 of file G4CoupledTransportation.cc.

 {
 
   G4TransportationManager* transportMgr =
     G4TransportationManager::GetTransportationManager();
 
   // G4VProcess::StartTracking(aTrack);
 
   //  The 'initialising' actions
   //     once taken in AlongStepGPIL -- if ( track.GetCurrentStepNumber()==1 )
 
   // fStartedNewTrack= true; 
 
   fMassNavigator = transportMgr->GetNavigatorForTracking() ; 
   fNavigatorId= transportMgr->ActivateNavigator( fMassNavigator );  // Confirm it!
 
   // if( fVerboseLevel > 1 ){
   //  G4cout << " Navigator Id obtained in StartTracking " << fNavigatorId << G4endl;
   // }
   G4ThreeVector position = aTrack->GetPosition(); 
   G4ThreeVector direction = aTrack->GetMomentumDirection();
 
   // if( fVerboseLevel > 1 ){
   //   G4cout << " Calling PathFinder::PrepareNewTrack from    " 
   //   << " G4CoupledTransportation::StartTracking -- which calls Locate()" << G4endl;
   // }
   fPathFinder->PrepareNewTrack( position, direction); 
   // This implies a call to fPathFinder->Locate( position, direction ); 
 
   // Global field, if any, must exist before tracking is started
   fGlobalFieldExists= DoesGlobalFieldExist(); 
   // reset safety value and center
   //
   fPreviousMassSafety  = 0.0 ; 
   fPreviousFullSafety  = 0.0 ; 
   fPreviousSftOrigin = G4ThreeVector(0.,0.,0.) ;
   
   // reset looping counter -- for motion in field  
   fNoLooperTrials= 0; 
   // Must clear this state .. else it depends on last track's value
   //  --> a better solution would set this from state of suspended track TODO ? 
   // Was if( aTrack->GetCurrentStepNumber()==1 ) { .. }
 
   // ChordFinder reset internal state
   //
   if( fGlobalFieldExists )
   {
      fFieldPropagator->ClearPropagatorState();   
        // Resets safety values, in case of overlaps.  
 
      G4ChordFinder* chordF= fFieldPropagator->GetChordFinder();
      if( chordF )  { chordF->ResetStepEstimate(); }
   }
 
   // Clear the chord finders of all fields (ie managers) derived objects
   //
   G4FieldManagerStore* fieldMgrStore = G4FieldManagerStore::GetInstance();
   fieldMgrStore->ClearAllChordFindersState(); 
 
 #ifdef G4DEBUG_TRANSPORT
   if( fVerboseLevel > 1 )
   {
     G4cout << " Returning touchable handle " << fCurrentTouchableHandle << G4endl;
   }
 #endif
 
   // Update the current touchable handle  (from the track's)
   //
   fCurrentTouchableHandle = aTrack->GetTouchableHandle();  
 }