Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2026-05-02 07:42:28

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 /// \file G4MonopoleTransportation.cc
0027 /// \brief Implementation of the G4MonopoleTransportation class
0028 
0029 // This class is a process responsible for the transportation of
0030 // magnetic monopoles, ie the geometrical propagation that encounters the
0031 // geometrical sub-volumes of the detectors.
0032 //
0033 // For monopoles, uses a different equation of motion and ignores energy
0034 // conservation.
0035 //
0036 
0037 // =======================================================================
0038 // Created:  3 May 2010, J. Apostolakis, B. Bozsogi
0039 // =======================================================================
0040 
0041 #include "G4MonopoleTransportation.hh"
0042 
0043 #include "DetectorConstruction.hh"
0044 
0045 #include "G4ChordFinder.hh"
0046 #include "G4FieldManagerStore.hh"
0047 #include "G4Monopole.hh"
0048 #include "G4ParticleTable.hh"
0049 #include "G4ProductionCutsTable.hh"
0050 #include "G4RunManager.hh"
0051 #include "G4SafetyHelper.hh"
0052 #include "G4SystemOfUnits.hh"
0053 #include "G4TransportationProcessType.hh"
0054 
0055 class G4VSensitiveDetector;
0056 
0057 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0058 
0059 G4MonopoleTransportation::G4MonopoleTransportation(const G4Monopole* mpl, G4int verb)
0060   : G4VProcess(G4String("MonopoleTransportation"), fTransportation),
0061     fParticleDef(mpl),
0062     fMagSetup(0),
0063     fLinearNavigator(0),
0064     fFieldPropagator(0),
0065     fParticleIsLooping(false),
0066     fPreviousSftOrigin(0., 0., 0.),
0067     fPreviousSafety(0.0),
0068     fThreshold_Warning_Energy(100 * MeV),
0069     fThreshold_Important_Energy(250 * MeV),
0070     fThresholdTrials(10),
0071     // fUnimportant_Energy( 1 * MeV ),
0072     fNoLooperTrials(0),
0073     fSumEnergyKilled(0.0),
0074     fMaxEnergyKilled(0.0),
0075     fShortStepOptimisation(false),  // Old default: true (=fast short steps)
0076     fpSafetyHelper(0),
0077     noCalls(0)
0078 {
0079   verboseLevel = verb;
0080 
0081   // set Process Sub Type
0082   SetProcessSubType(TRANSPORTATION);
0083 
0084   // Do not finalize the G4MonopoleTransportation class
0085   if (G4Threading::IsMasterThread() && G4Threading::IsMultithreadedApplication()) {
0086     return;
0087   }
0088 
0089   const DetectorConstruction* detector = static_cast<const DetectorConstruction*>(
0090     G4RunManager::GetRunManager()->GetUserDetectorConstruction());
0091 
0092   fMagSetup = detector->GetMonopoleFieldSetup();
0093 
0094   G4TransportationManager* transportMgr = G4TransportationManager::GetTransportationManager();
0095 
0096   fLinearNavigator = transportMgr->GetNavigatorForTracking();
0097 
0098   fFieldPropagator = transportMgr->GetPropagatorInField();
0099   fpSafetyHelper = transportMgr->GetSafetyHelper();
0100 
0101   // New
0102 
0103   // Cannot determine whether a field exists here,
0104   //  because it would only work if the field manager has informed
0105   //  about the detector's field before this transportation process
0106   //  is constructed.
0107   // Instead later the method DoesGlobalFieldExist() is called
0108   fCurrentTouchableHandle = nullptr;
0109 
0110   fEndGlobalTimeComputed = false;
0111   fCandidateEndGlobalTime = 0;
0112 }
0113 
0114 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0115 
0116 G4MonopoleTransportation::~G4MonopoleTransportation()
0117 {
0118   if ((verboseLevel > 0) && (fSumEnergyKilled > 0.0)) {
0119     G4cout << " G4MonopoleTransportation: Statistics for looping particles " << G4endl;
0120     G4cout << "   Sum of energy of loopers killed: " << fSumEnergyKilled << G4endl;
0121     G4cout << "   Max energy of loopers killed: " << fMaxEnergyKilled << G4endl;
0122   }
0123 }
0124 
0125 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0126 //
0127 // Responsibilities:
0128 //    Find whether the geometry limits the Step, and to what length
0129 //    Calculate the new value of the safety and return it.
0130 //    Store the final time, position and momentum.
0131 
0132 G4double G4MonopoleTransportation::AlongStepGetPhysicalInteractionLength(
0133   const G4Track& track,
0134   G4double,  //  previousStepSize
0135   G4double currentMinimumStep, G4double& currentSafety, G4GPILSelection* selection)
0136 {
0137   fMagSetup->SetStepperAndChordFinder(1);
0138   // change to monopole equation
0139 
0140   G4double geometryStepLength, newSafety;
0141   fParticleIsLooping = false;
0142 
0143   // Initial actions moved to  StartTrack()
0144   // --------------------------------------
0145   // Note: in case another process changes touchable handle
0146   //    it will be necessary to add here (for all steps)
0147   // fCurrentTouchableHandle = aTrack->GetTouchableHandle();
0148 
0149   // GPILSelection is set to defaule value of CandidateForSelection
0150   // It is a return value
0151   //
0152   *selection = CandidateForSelection;
0153 
0154   // Get initial Energy/Momentum of the track
0155   //
0156   const G4DynamicParticle* pParticle = track.GetDynamicParticle();
0157   G4ThreeVector startMomentumDir = pParticle->GetMomentumDirection();
0158   G4ThreeVector startPosition = track.GetPosition();
0159 
0160   // G4double   theTime        = track.GetGlobalTime() ;
0161 
0162   // The Step Point safety can be limited by other geometries and/or the
0163   // assumptions of any process - it's not always the geometrical safety.
0164   // We calculate the starting point's isotropic safety here.
0165   //
0166   G4ThreeVector OriginShift = startPosition - fPreviousSftOrigin;
0167   G4double MagSqShift = OriginShift.mag2();
0168   if (MagSqShift >= sqr(fPreviousSafety)) {
0169     currentSafety = 0.0;
0170   }
0171   else {
0172     currentSafety = fPreviousSafety - std::sqrt(MagSqShift);
0173   }
0174 
0175   // Is the monopole charged ?
0176   //
0177   G4double particleMagneticCharge = fParticleDef->MagneticCharge();
0178   G4double particleElectricCharge = pParticle->GetCharge();
0179 
0180   fGeometryLimitedStep = false;
0181   // fEndGlobalTimeComputed = false ;
0182 
0183   // There is no need to locate the current volume. It is Done elsewhere:
0184   //   On track construction
0185   //   By the tracking, after all AlongStepDoIts, in "Relocation"
0186 
0187   // Check whether the particle have an (EM) field force exerting upon it
0188   //
0189   G4FieldManager* fieldMgr = 0;
0190   G4bool fieldExertsForce = false;
0191 
0192   if ((particleMagneticCharge != 0.0)) {
0193     fieldMgr = fFieldPropagator->FindAndSetFieldManager(track.GetVolume());
0194     if (fieldMgr != 0) {
0195       // Message the field Manager, to configure it for this track
0196       fieldMgr->ConfigureForTrack(&track);
0197       // Moved here, in order to allow a transition
0198       //   from a zero-field  status (with fieldMgr->(field)0
0199       //   to a finite field  status
0200 
0201       // If the field manager has no field, there is no field !
0202       fieldExertsForce = (fieldMgr->GetDetectorField() != 0);
0203     }
0204   }
0205 
0206   // G4cout << " G4Transport:  field exerts force= " << fieldExertsForce
0207   //          << "  fieldMgr= " << fieldMgr << G4endl;
0208 
0209   // Choose the calculation of the transportation: Field or not
0210   //
0211   if (!fieldExertsForce) {
0212     G4double linearStepLength;
0213     if (fShortStepOptimisation && (currentMinimumStep <= currentSafety)) {
0214       // The Step is guaranteed to be taken
0215       //
0216       geometryStepLength = currentMinimumStep;
0217       fGeometryLimitedStep = false;
0218     }
0219     else {
0220       //  Find whether the straight path intersects a volume
0221       //
0222       linearStepLength = fLinearNavigator->ComputeStep(startPosition, startMomentumDir,
0223                                                        currentMinimumStep, newSafety);
0224       // Remember last safety origin & value.
0225       //
0226       fPreviousSftOrigin = startPosition;
0227       fPreviousSafety = newSafety;
0228       // fpSafetyHelper->SetCurrentSafety( newSafety, startPosition);
0229 
0230       // The safety at the initial point has been re-calculated:
0231       //
0232       currentSafety = newSafety;
0233 
0234       fGeometryLimitedStep = (linearStepLength <= currentMinimumStep);
0235       if (fGeometryLimitedStep) {
0236         // The geometry limits the Step size (an intersection was found.)
0237         geometryStepLength = linearStepLength;
0238       }
0239       else {
0240         // The full Step is taken.
0241         geometryStepLength = currentMinimumStep;
0242       }
0243     }
0244     endpointDistance = geometryStepLength;
0245 
0246     // Calculate final position
0247     //
0248     fTransportEndPosition = startPosition + geometryStepLength * startMomentumDir;
0249 
0250     // Momentum direction, energy and polarisation are unchanged by transport
0251     //
0252     fTransportEndMomentumDir = startMomentumDir;
0253     fTransportEndKineticEnergy = track.GetKineticEnergy();
0254     fTransportEndSpin = track.GetPolarization();
0255     fParticleIsLooping = false;
0256     fMomentumChanged = false;
0257     fEndGlobalTimeComputed = false;
0258   }
0259   else  //  A field exerts force
0260   {
0261     G4double momentumMagnitude = pParticle->GetTotalMomentum();
0262     G4ThreeVector EndUnitMomentum;
0263     G4double lengthAlongCurve;
0264     G4double restMass = fParticleDef->GetPDGMass();
0265 
0266     G4ChargeState chargeState(particleElectricCharge,  // The charge can change
0267                               fParticleDef->GetPDGSpin(),
0268                               0,  //  Magnetic moment:
0269                               0,  //  Electric Dipole moment
0270                               particleMagneticCharge);  // in Mev/c
0271 
0272     G4EquationOfMotion* equationOfMotion =
0273       fFieldPropagator->GetChordFinder()->GetIntegrationDriver()->GetEquationOfMotion();
0274 
0275     equationOfMotion->SetChargeMomentumMass(chargeState, momentumMagnitude, restMass);
0276     // SetChargeMomentumMass now passes both the electric and magnetic
0277     // charge in chargeState
0278 
0279     G4ThreeVector spin = track.GetPolarization();
0280     G4FieldTrack aFieldTrack = G4FieldTrack(startPosition, track.GetMomentumDirection(), 0.0,
0281                                             track.GetKineticEnergy(), restMass, track.GetVelocity(),
0282                                             track.GetGlobalTime(),  // Lab.
0283                                             track.GetProperTime(),  // Part.
0284                                             &spin);
0285     if (currentMinimumStep > 0) {
0286       // Do the Transport in the field (non recti-linear)
0287       //
0288       lengthAlongCurve = fFieldPropagator->ComputeStep(aFieldTrack, currentMinimumStep,
0289                                                        currentSafety, track.GetVolume());
0290       fGeometryLimitedStep = lengthAlongCurve < currentMinimumStep;
0291       if (fGeometryLimitedStep) {
0292         geometryStepLength = lengthAlongCurve;
0293       }
0294       else {
0295         geometryStepLength = currentMinimumStep;
0296       }
0297     }
0298     else {
0299       geometryStepLength = lengthAlongCurve = 0.0;
0300       fGeometryLimitedStep = false;
0301     }
0302 
0303     // Remember last safety origin & value.
0304     //
0305     fPreviousSftOrigin = startPosition;
0306     fPreviousSafety = currentSafety;
0307     // fpSafetyHelper->SetCurrentSafety( newSafety, startPosition);
0308 
0309     // Get the End-Position and End-Momentum (Dir-ection)
0310     //
0311     fTransportEndPosition = aFieldTrack.GetPosition();
0312 
0313     // Momentum:  Magnitude and direction can be changed too now ...
0314     //
0315     fMomentumChanged = true;
0316     fTransportEndMomentumDir = aFieldTrack.GetMomentumDir();
0317 
0318     fTransportEndKineticEnergy = aFieldTrack.GetKineticEnergy();
0319 
0320     fCandidateEndGlobalTime = aFieldTrack.GetLabTimeOfFlight();
0321     fEndGlobalTimeComputed = true;
0322 
0323     fTransportEndSpin = aFieldTrack.GetSpin();
0324     fParticleIsLooping = fFieldPropagator->IsParticleLooping();
0325     endpointDistance = (fTransportEndPosition - startPosition).mag();
0326   }
0327 
0328   // If we are asked to go a step length of 0, and we are on a boundary
0329   // then a boundary will also limit the step -> we must flag this.
0330   //
0331   if (currentMinimumStep == 0.0) {
0332     if (currentSafety == 0.0) fGeometryLimitedStep = true;
0333   }
0334 
0335   // Update the safety starting from the end-point,
0336   // if it will become negative at the end-point.
0337   //
0338   if (currentSafety < endpointDistance) {
0339     // if( particleMagneticCharge == 0.0 )
0340     // G4cout  << "  Avoiding call to ComputeSafety : charge = 0.0 " << G4endl;
0341 
0342     if (particleMagneticCharge != 0.0) {
0343       G4double endSafety = fLinearNavigator->ComputeSafety(fTransportEndPosition);
0344       currentSafety = endSafety;
0345       fPreviousSftOrigin = fTransportEndPosition;
0346       fPreviousSafety = currentSafety;
0347       fpSafetyHelper->SetCurrentSafety(currentSafety, fTransportEndPosition);
0348 
0349       // Because the Stepping Manager assumes it is from the start point,
0350       //  add the StepLength
0351       //
0352       currentSafety += endpointDistance;
0353 
0354 #ifdef G4DEBUG_TRANSPORT
0355       G4cout.precision(12);
0356       G4cout << "***G4MonopoleTransportation::AlongStepGPIL ** " << G4endl;
0357       G4cout << "  Called Navigator->ComputeSafety at " << fTransportEndPosition
0358              << "    and it returned safety= " << endSafety << G4endl;
0359       G4cout << "  Adding endpoint distance " << endpointDistance
0360              << "   to obtain pseudo-safety= " << currentSafety << G4endl;
0361 #endif
0362     }
0363   }
0364 
0365   fParticleChange.ProposeTrueStepLength(geometryStepLength);
0366 
0367   fMagSetup->SetStepperAndChordFinder(0);
0368   // change back to usual equation
0369 
0370   return geometryStepLength;
0371 }
0372 
0373 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0374 //
0375 //   Initialize ParticleChange  (by setting all its members equal
0376 //                               to corresponding members in G4Track)
0377 
0378 G4VParticleChange* G4MonopoleTransportation::AlongStepDoIt(const G4Track& track,
0379                                                            const G4Step& stepData)
0380 {
0381   ++noCalls;
0382 
0383   if (fGeometryLimitedStep) {
0384     stepData.GetPostStepPoint()->SetStepStatus(fGeomBoundary);
0385   }
0386 
0387   fParticleChange.Initialize(track);
0388 
0389   //  Code for specific process
0390   //
0391   fParticleChange.ProposePosition(fTransportEndPosition);
0392   fParticleChange.ProposeMomentumDirection(fTransportEndMomentumDir);
0393   fParticleChange.ProposeEnergy(fTransportEndKineticEnergy);
0394   fParticleChange.SetMomentumChanged(fMomentumChanged);
0395 
0396   fParticleChange.ProposePolarization(fTransportEndSpin);
0397 
0398   G4double deltaTime = 0.0;
0399 
0400   // Calculate  Lab Time of Flight (ONLY if field Equations used it!)
0401 
0402   G4double startTime = track.GetGlobalTime();
0403 
0404   if (!fEndGlobalTimeComputed) {
0405     // The time was not integrated .. make the best estimate possible
0406     //
0407     G4double finalVelocity = track.GetVelocity();
0408     G4double initialVelocity = stepData.GetPreStepPoint()->GetVelocity();
0409     G4double stepLength = track.GetStepLength();
0410 
0411     deltaTime = 0.0;  // in case initialVelocity = 0
0412     if (finalVelocity > 0.0) {
0413       G4double meanInverseVelocity;
0414       // deltaTime = stepLength/finalVelocity ;
0415       meanInverseVelocity = 0.5 * (1.0 / initialVelocity + 1.0 / finalVelocity);
0416       deltaTime = stepLength * meanInverseVelocity;
0417     }
0418     else if (initialVelocity > 0.0) {
0419       deltaTime = stepLength / initialVelocity;
0420     }
0421     fCandidateEndGlobalTime = startTime + deltaTime;
0422   }
0423   else {
0424     deltaTime = fCandidateEndGlobalTime - startTime;
0425   }
0426 
0427   fParticleChange.ProposeGlobalTime(fCandidateEndGlobalTime);
0428 
0429   // Now Correct by Lorentz factor to get "proper" deltaTime
0430 
0431   G4double restMass = track.GetDynamicParticle()->GetMass();
0432   G4double deltaProperTime = deltaTime * (restMass / track.GetTotalEnergy());
0433 
0434   fParticleChange.ProposeProperTime(track.GetProperTime() + deltaProperTime);
0435 
0436   // If the particle is caught looping or is stuck (in very difficult
0437   // boundaries) in a magnetic field (doing many steps)
0438   //   THEN this kills it ...
0439   //
0440   if (fParticleIsLooping) {
0441     G4double endEnergy = fTransportEndKineticEnergy;
0442 
0443     if ((endEnergy < fThreshold_Important_Energy) || (fNoLooperTrials >= fThresholdTrials)) {
0444       // Kill the looping particle
0445       //
0446       fParticleChange.ProposeTrackStatus(fStopAndKill);
0447 
0448       // 'Bare' statistics
0449       fSumEnergyKilled += endEnergy;
0450       if (endEnergy > fMaxEnergyKilled) {
0451         fMaxEnergyKilled = endEnergy;
0452       }
0453 
0454 #ifdef G4VERBOSE
0455       if ((verboseLevel > 1) || (endEnergy > fThreshold_Warning_Energy)) {
0456         G4cout << " G4MonopoleTransportation is killing track "
0457                << "that is looping or stuck " << G4endl << "   This track has "
0458                << track.GetKineticEnergy() / MeV << " MeV energy." << G4endl;
0459         G4cout << "   Number of trials = " << fNoLooperTrials
0460                << "   No of calls to AlongStepDoIt = " << noCalls << G4endl;
0461       }
0462 #endif
0463       fNoLooperTrials = 0;
0464     }
0465     else {
0466       fNoLooperTrials++;
0467 #ifdef G4VERBOSE
0468       if ((verboseLevel > 2)) {
0469         G4cout << "   G4MonopoleTransportation::AlongStepDoIt(): "
0470                << "Particle looping - "
0471                << "   Number of trials = " << fNoLooperTrials << "   No of calls to  = " << noCalls
0472                << G4endl;
0473       }
0474 #endif
0475     }
0476   }
0477   else {
0478     fNoLooperTrials = 0;
0479   }
0480 
0481   // Another (sometimes better way) is to use a user-limit maximum Step size
0482   // to alleviate this problem ..
0483 
0484   // Introduce smooth curved trajectories to particle-change
0485   //
0486   fParticleChange.SetPointerToVectorOfAuxiliaryPoints(
0487     fFieldPropagator->GimmeTrajectoryVectorAndForgetIt());
0488 
0489   return &fParticleChange;
0490 }
0491 
0492 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0493 //
0494 //  This ensures that the PostStep action is always called,
0495 //  so that it can do the relocation if it is needed.
0496 //
0497 
0498 G4double
0499 G4MonopoleTransportation::PostStepGetPhysicalInteractionLength(const G4Track&,
0500                                                                G4double,  // previousStepSize
0501                                                                G4ForceCondition* pForceCond)
0502 {
0503   *pForceCond = Forced;
0504   return DBL_MAX;  // was kInfinity ; but convention now is DBL_MAX
0505 }
0506 
0507 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0508 
0509 G4VParticleChange* G4MonopoleTransportation::PostStepDoIt(const G4Track& track, const G4Step&)
0510 {
0511   G4TouchableHandle retCurrentTouchable;  // The one to return
0512 
0513   // Initialize ParticleChange  (by setting all its members equal
0514   //                             to corresponding members in G4Track)
0515   // fParticleChange.Initialize(track) ;  // To initialise TouchableChange
0516 
0517   fParticleChange.ProposeTrackStatus(track.GetTrackStatus());
0518 
0519   // If the Step was determined by the volume boundary,
0520   // logically relocate the particle
0521 
0522   if (fGeometryLimitedStep) {
0523     // fCurrentTouchable will now become the previous touchable,
0524     // and what was the previous will be freed.
0525     // (Needed because the preStepPoint can point to the previous touchable)
0526 
0527     fLinearNavigator->SetGeometricallyLimitedStep();
0528     fLinearNavigator->LocateGlobalPointAndUpdateTouchableHandle(
0529       track.GetPosition(), track.GetMomentumDirection(), fCurrentTouchableHandle, true);
0530     // Check whether the particle is out of the world volume
0531     // If so it has exited and must be killed.
0532     //
0533     if (fCurrentTouchableHandle->GetVolume() == 0) {
0534       fParticleChange.ProposeTrackStatus(fStopAndKill);
0535     }
0536     retCurrentTouchable = fCurrentTouchableHandle;
0537     fParticleChange.SetTouchableHandle(fCurrentTouchableHandle);
0538   }
0539   else  // fGeometryLimitedStep  is false
0540   {
0541     // This serves only to move the Navigator's location
0542     //
0543     fLinearNavigator->LocateGlobalPointWithinVolume(track.GetPosition());
0544 
0545     // The value of the track's current Touchable is retained.
0546     // (and it must be correct because we must use it below to
0547     // overwrite the (unset) one in particle change)
0548     //  It must be fCurrentTouchable too ??
0549     //
0550     fParticleChange.SetTouchableHandle(track.GetTouchableHandle());
0551     retCurrentTouchable = track.GetTouchableHandle();
0552   }  // endif ( fGeometryLimitedStep )
0553 
0554   const G4VPhysicalVolume* pNewVol = retCurrentTouchable->GetVolume();
0555   const G4Material* pNewMaterial = 0;
0556   G4VSensitiveDetector* pNewSensitiveDetector = 0;
0557   if (pNewVol != 0) {
0558     pNewMaterial = pNewVol->GetLogicalVolume()->GetMaterial();
0559     pNewSensitiveDetector = pNewVol->GetLogicalVolume()->GetSensitiveDetector();
0560   }
0561 
0562   // ( <const_cast> pNewMaterial ) ;
0563 
0564   fParticleChange.SetMaterialInTouchable((G4Material*)pNewMaterial);
0565   fParticleChange.SetSensitiveDetectorInTouchable(pNewSensitiveDetector);
0566 
0567   const G4MaterialCutsCouple* pNewMaterialCutsCouple = 0;
0568   if (pNewVol != 0) {
0569     pNewMaterialCutsCouple = pNewVol->GetLogicalVolume()->GetMaterialCutsCouple();
0570   }
0571 
0572   if (pNewVol != 0 && pNewMaterialCutsCouple != 0
0573       && pNewMaterialCutsCouple->GetMaterial() != pNewMaterial)
0574   {
0575     // for parametrized volume
0576     //
0577     pNewMaterialCutsCouple = G4ProductionCutsTable::GetProductionCutsTable()->GetMaterialCutsCouple(
0578       pNewMaterial, pNewMaterialCutsCouple->GetProductionCuts());
0579   }
0580   fParticleChange.SetMaterialCutsCoupleInTouchable(pNewMaterialCutsCouple);
0581 
0582   // temporarily until Get/Set Material of ParticleChange,
0583   // and StepPoint can be made const.
0584   // Set the touchable in ParticleChange
0585   // this must always be done because the particle change always
0586   // uses this value to overwrite the current touchable pointer.
0587   //
0588   fParticleChange.SetTouchableHandle(retCurrentTouchable);
0589 
0590   return &fParticleChange;
0591 }
0592 
0593 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0594 
0595 // New method takes over the responsibility to reset the state
0596 // of G4MonopoleTransportation object at the start of a new track
0597 // or the resumption of a suspended track.
0598 
0599 void G4MonopoleTransportation::StartTracking(G4Track* aTrack)
0600 {
0601   G4VProcess::StartTracking(aTrack);
0602 
0603   // The actions here are those that were taken in AlongStepGPIL
0604   //   when track.GetCurrentStepNumber()==1
0605 
0606   // reset safety value and center
0607   //
0608   fPreviousSafety = 0.0;
0609   fPreviousSftOrigin = G4ThreeVector(0., 0., 0.);
0610 
0611   // reset looping counter -- for motion in field
0612   fNoLooperTrials = 0;
0613   // Must clear this state .. else it depends on last track's value
0614   //  --> a better solution would set this from state of suspended track TODO ?
0615   // Was if( aTrack->GetCurrentStepNumber()==1 ) { .. }
0616 
0617   // ChordFinder reset internal state
0618   //
0619   if (DoesGlobalFieldExist()) {
0620     fFieldPropagator->ClearPropagatorState();
0621     // Resets all state of field propagator class (ONLY)
0622     // including safety values
0623     // in case of overlaps and to wipe for first track).
0624 
0625     // G4ChordFinder* chordF= fFieldPropagator->GetChordFinder();
0626     // if( chordF ) chordF->ResetStepEstimate();
0627   }
0628 
0629   // Make sure to clear the chord finders of all fields (ie managers)
0630   G4FieldManagerStore* fieldMgrStore = G4FieldManagerStore::GetInstance();
0631   fieldMgrStore->ClearAllChordFindersState();
0632 
0633   // Update the current touchable handle  (from the track's)
0634   //
0635   fCurrentTouchableHandle = aTrack->GetTouchableHandle();
0636 }
0637 
0638 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......