BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
|
Extra G4Navigator to get coordinate transforms. More...
#include <BDSAuxiliaryNavigator.hh>
Public Member Functions | |
G4VPhysicalVolume * | LocateGlobalPointAndSetup (const G4ThreeVector &point, const G4ThreeVector *direction=nullptr, const G4bool pRelativeSearch=true, const G4bool ignoreDirection=true, G4bool useCurvilinear=true) const |
A wrapper for the underlying static navigator instance located within this class. | |
G4VPhysicalVolume * | LocateGlobalPointAndSetup (G4Step const *const step, G4bool useCurvilinear=true) const |
BDSStep | ConvertToLocal (G4Step const *const step, G4bool useCurvilinear=true) const |
BDSStep | ConvertToLocal (const G4ThreeVector &globalPosition, const G4ThreeVector &globalDirection, const G4double stepLength=0, const G4bool useCurvilinear=true, const G4double marginLength=1) const |
BDSStep | ConvertToGlobalStep (const G4ThreeVector &localPosition, const G4ThreeVector &localDirection, const G4bool useCurvilinear=true) const |
G4ThreeVector | ConvertToLocal (const G4double globalPoint[3], const G4bool useCurvilinear=true) const |
G4ThreeVector | ConvertToLocal (const G4ThreeVector &globalPosition, const G4bool useCurvilinear=true) const |
Vector version - see notes above. | |
G4ThreeVector | ConvertToLocalNoSetup (const G4ThreeVector &globalPosition, const G4bool useCurvilinear=true) const |
Similar to above function but does NOT initialise the transforms. | |
G4ThreeVector | ConvertAxisToLocal (const G4ThreeVector &globalAxis, const G4bool useCurvilinear=true) const |
G4ThreeVector | ConvertAxisToLocal (const G4double globalPoint[3], const G4double globalAxis[3], const G4bool useCurvilinear=true) const |
G4ThreeVector | ConvertAxisToLocal (const G4ThreeVector &globalPoint, const G4ThreeVector &globalAxis, const G4bool useCurvilinear=true) const |
Vector version. | |
G4ThreeVector | ConvertAxisToGlobal (const G4ThreeVector &localAxis, const G4bool useCurvilinear=true) const |
std::pair< G4ThreeVector, G4ThreeVector > | ConvertAxisToGlobal (const std::pair< G4ThreeVector, G4ThreeVector > &localAxis, const G4bool useCurvilinear=true) const |
G4ThreeVector | ConvertToGlobal (const G4ThreeVector &localPosition, const G4bool useCurvilinear=true) const |
G4ThreeVector | ConvertAxisToGlobal (const G4ThreeVector &globalPosition, const G4ThreeVector &localAxis, const G4bool useCurvilinear=true) const |
G4ThreeVector | ConvertToGlobal (const G4ThreeVector &globalPosition, const G4ThreeVector &localPosition, const G4bool useCurvilinear=true) const |
BDSStep | GlobalToCurvilinear (const G4double fieldArcLength, const G4ThreeVector &unitField, const G4double angle, const G4ThreeVector &position, const G4ThreeVector &unitMomentum, const G4double h, const G4bool useCurvilinearWorld, const G4double FCof, const G4double tilt=0) |
BDSStep | GlobalToCurvilinear (const G4ThreeVector &position, const G4ThreeVector &unitMomentum, const G4double h, const G4bool useCurvilinearWorld) |
BDSStep | CurvilinearToGlobal (const G4ThreeVector &localPosition, const G4ThreeVector &localMomentum, const G4bool useCurvilinearWorld) |
BDSStep | CurvilinearToGlobal (const G4double fieldArcLength, const G4ThreeVector &unitField, const G4double angle, const G4ThreeVector &CLPosition, const G4ThreeVector &CLMomentum, const G4bool useCurvilinearWorld, const G4double FCof, const G4double tilt=0) |
Static Public Member Functions | |
static void | AttachWorldVolumeToNavigator (G4VPhysicalVolume *worldPVIn) |
Setup the navigator w.r.t. to a world volume - typically real world. | |
static void | AttachWorldVolumeToNavigatorCL (G4VPhysicalVolume *curvilinearWorldPVIn) |
static void | RegisterCurvilinearBridgeWorld (G4VPhysicalVolume *curvilinearBridgeWorldPVIn) |
static void | ResetNavigatorStates () |
Protected Attributes | |
G4AffineTransform | globalToLocal |
G4AffineTransform | localToGlobal |
G4AffineTransform | globalToLocalCL |
G4AffineTransform | localToGlobalCL |
G4bool | bridgeVolumeWasUsed |
Static Protected Attributes | |
static G4Navigator * | auxNavigator = new G4Navigator() |
static G4Navigator * | auxNavigatorCL = new G4Navigator() |
static G4Navigator * | auxNavigatorCLB = new G4Navigator() |
Private Member Functions | |
G4Navigator * | Navigator (G4bool curvilinear) const |
Utility function to select appropriate navigator. | |
void | InitialiseTransform (const G4bool massworld=true, const G4bool curvilinearWorld=true) const |
void | InitialiseTransform (const G4ThreeVector &globalPosition) const |
void | InitialiseTransform (const G4ThreeVector &globalPosition, const G4ThreeVector &globalMomentum, const G4double stepLength) |
const G4AffineTransform & | GlobalToLocal (G4bool curvilinear) const |
Utility function to select appropriate transform. | |
const G4AffineTransform & | LocalToGlobal (G4bool curvilinear) const |
Utility function to select appropriate transform. | |
Private Attributes | |
G4double | volumeMargin |
Static Private Attributes | |
static G4int | numberOfInstances = 0 |
static G4VPhysicalVolume * | worldPV = nullptr |
Cache of world PV to test if we're getting the wrong volume for the transform. | |
static G4VPhysicalVolume * | curvilinearWorldPV = nullptr |
Cache of world PV to test if we're getting the wrong volume for the transform. | |
static G4VPhysicalVolume * | curvilinearBridgeWorldPV = nullptr |
Cache of world PV to test if we're getting the wrong volume for the transform. | |
Extra G4Navigator to get coordinate transforms.
All BDSIM steppers and magnetic fields require the ability to convert from global to local coordinates. The prescribed method to do this is by using a G4Navigator instance. There is of course the main navigator for tracking, but requesting information on a global point changes the state of the navigator - ie the particle being tracked by the navigator is now that point. This is VERY dangerous.
This class provides a static auxiliary navigator that each derived class can use. Making the auxiliary navigator static is not done to reduce memory usage but because navigating from an unknown place to anywhere in the geometry is much more costly than a relative move in the geometry. If we only use one auxiliary navigator, it will always be relatively close in the geometry even if a different stepper has been used and therefore more efficient. This is important as Geant4 may use the steppers at least two or three times to estimate the best next step and the stepper itself may make three steps (full, and two half) to estimate the error in the integration.
This class in fact provides access two static navigators. The intended use is one for the real world and one for the read out geometry / world for curvilinear coordinates. All functions have an optional last argument to select which navigator is required - the default is the curvilinear one.
Definition at line 64 of file BDSAuxiliaryNavigator.hh.
BDSAuxiliaryNavigator::BDSAuxiliaryNavigator | ( | ) |
Definition at line 38 of file BDSAuxiliaryNavigator.cc.
BDSAuxiliaryNavigator::~BDSAuxiliaryNavigator | ( | ) |
Definition at line 49 of file BDSAuxiliaryNavigator.cc.
|
inlinestatic |
Setup the navigator w.r.t. to a world volume - typically real world.
Definition at line 71 of file BDSAuxiliaryNavigator.hh.
References auxNavigator, and worldPV.
Referenced by BDSDetectorConstruction::BuildWorld().
|
inlinestatic |
Setup the navigator w.r.t. to the read out world / geometry to provide curvilinear coordinates.
Definition at line 76 of file BDSAuxiliaryNavigator.hh.
References auxNavigatorCL, and curvilinearWorldPV.
Referenced by BDSParallelWorldCurvilinear::Construct().
G4ThreeVector BDSAuxiliaryNavigator::ConvertAxisToGlobal | ( | const G4ThreeVector & | globalPosition, |
const G4ThreeVector & | localAxis, | ||
const G4bool | useCurvilinear = true |
||
) | const |
Convert a vector (axis) from local to global coordinates. This uses a global position to ensure the transform is initialised.
Definition at line 278 of file BDSAuxiliaryNavigator.cc.
References LocalToGlobal().
G4ThreeVector BDSAuxiliaryNavigator::ConvertAxisToGlobal | ( | const G4ThreeVector & | localAxis, |
const G4bool | useCurvilinear = true |
||
) | const |
Convert a vector (axis) from local to global coordinates. NOTE this function must only be used once the instance of this class has been initialised, setting up the transforms. It is up to the developer to ensure this, otherwise you'll find a bad access.
Definition at line 261 of file BDSAuxiliaryNavigator.cc.
References LocalToGlobal().
Referenced by BDSIntegratorDipoleFringe::BaseStepper(), BDSFieldEGlobal::GetField(), BDSFieldEMGlobal::GetField(), BDSFieldMagGlobal::GetField(), and BDSIntegratorDipoleFringe::MultipoleStep().
std::pair< G4ThreeVector, G4ThreeVector > BDSAuxiliaryNavigator::ConvertAxisToGlobal | ( | const std::pair< G4ThreeVector, G4ThreeVector > & | localAxis, |
const G4bool | useCurvilinear = true |
||
) | const |
Convert a vector (axis) from local to global coordinates. Note this function must only be used once the instance of this class has been initialised, setting up the transforms. It is up to the developer to ensure this. This utility function operates on two threevectors in a pair.
Definition at line 265 of file BDSAuxiliaryNavigator.cc.
References LocalToGlobal().
G4ThreeVector BDSAuxiliaryNavigator::ConvertAxisToLocal | ( | const G4double | globalPoint[3], |
const G4double | globalAxis[3], | ||
const G4bool | useCurvilinear = true |
||
) | const |
Calculate the local version of a global vector (axis). This is done w.r.t. a point - ie, the point is used to initialise the transform if not done already.
Definition at line 244 of file BDSAuxiliaryNavigator.cc.
References ConvertAxisToLocal().
G4ThreeVector BDSAuxiliaryNavigator::ConvertAxisToLocal | ( | const G4ThreeVector & | globalAxis, |
const G4bool | useCurvilinear = true |
||
) | const |
Convert an axis to curvilinear coordinates using the existing cached transforms. Therefore, this should only be used if you have converted a point or axis already in the current volume. Provided for the situation where multiple axis conversions are required without relocating the volume in the hierarchy, which may introduce ambiguities.
Definition at line 238 of file BDSAuxiliaryNavigator.cc.
References GlobalToLocal().
Referenced by ConvertAxisToLocal(), and BDSSDCollimator::ProcessHitsOrdered().
G4ThreeVector BDSAuxiliaryNavigator::ConvertAxisToLocal | ( | const G4ThreeVector & | globalPoint, |
const G4ThreeVector & | globalAxis, | ||
const G4bool | useCurvilinear = true |
||
) | const |
Vector version.
Definition at line 253 of file BDSAuxiliaryNavigator.cc.
References GlobalToLocal().
G4ThreeVector BDSAuxiliaryNavigator::ConvertToGlobal | ( | const G4ThreeVector & | globalPosition, |
const G4ThreeVector & | localPosition, | ||
const G4bool | useCurvilinear = true |
||
) | const |
Convert a position in local coordinates to global coordinates. This uses a global position to ensure the transform is initialised.
Definition at line 296 of file BDSAuxiliaryNavigator.cc.
References LocalToGlobal().
G4ThreeVector BDSAuxiliaryNavigator::ConvertToGlobal | ( | const G4ThreeVector & | localPosition, |
const G4bool | useCurvilinear = true |
||
) | const |
Convert a position in local coordinates to global coordinates. NOTE a similar caution to ConvertAxisToGlobal applies.
Definition at line 274 of file BDSAuxiliaryNavigator.cc.
References LocalToGlobal().
BDSStep BDSAuxiliaryNavigator::ConvertToGlobalStep | ( | const G4ThreeVector & | localPosition, |
const G4ThreeVector & | localDirection, | ||
const G4bool | useCurvilinear = true |
||
) | const |
Convert back to global coordinates. This DOES NOT update the transforms and uses the existing transforms inside this class - ie this relies on ConvertToLocal being used beforehand to initialise the transforms. This is done as we'd need a look up point behind this point.
Definition at line 286 of file BDSAuxiliaryNavigator.cc.
References LocalToGlobal().
Referenced by BDSIntegratorMag::ConvertToGlobal().
G4ThreeVector BDSAuxiliaryNavigator::ConvertToLocal | ( | const G4double | globalPoint[3], |
const G4bool | useCurvilinear = true |
||
) | const |
Calculate the local coordinates of a global point. Initialises the transforms by looking up the global point. Note, this is done without looking up the volume with direction, so it could be inaccurate on boundaries.
Definition at line 218 of file BDSAuxiliaryNavigator.cc.
References ConvertToLocal().
G4ThreeVector BDSAuxiliaryNavigator::ConvertToLocal | ( | const G4ThreeVector & | globalPosition, |
const G4bool | useCurvilinear = true |
||
) | const |
Vector version - see notes above.
Definition at line 225 of file BDSAuxiliaryNavigator.cc.
References GlobalToLocal().
BDSStep BDSAuxiliaryNavigator::ConvertToLocal | ( | const G4ThreeVector & | globalPosition, |
const G4ThreeVector & | globalDirection, | ||
const G4double | stepLength = 0 , |
||
const G4bool | useCurvilinear = true , |
||
const G4double | marginLength = 1 |
||
) | const |
Calculate the local coordinates for a position and direction along a step length. This is similar to the same function but for a G4Step but split apart. The direction vector can be used as the momentum vector without being a unit vector. The 'post step' vector in the BDSStep instance will be the direction vector (of same magnitude) but rotated to the local frame. This uses LocateGlobalPointAndSetup. Can control default lookahead distance default (1 -> 1mm).
Definition at line 182 of file BDSAuxiliaryNavigator.cc.
References GlobalToLocal(), and LocateGlobalPointAndSetup().
BDSStep BDSAuxiliaryNavigator::ConvertToLocal | ( | G4Step const *const | step, |
G4bool | useCurvilinear = true |
||
) | const |
Calculate the local coordinates for both a pre and post step point. The mid point of the step is used for the volume (and therefore transform) lookup which should ensure the correct volume is found - avoiding potential boundary issues between parallel worlds. This uses LocateGlobalPointAndSetup.
Definition at line 166 of file BDSAuxiliaryNavigator.cc.
References GlobalToLocal(), and LocateGlobalPointAndSetup().
Referenced by BDSIntegratorDipoleRodrigues::AdvanceHelix(), BDSIntegratorDecapole::AdvanceHelix(), BDSIntegratorOctupole::AdvanceHelix(), BDSIntegratorSextupole::AdvanceHelix(), BDSTrajectoryPoint::BDSTrajectoryPoint(), ConvertToLocal(), BDSFieldEGlobal::GetField(), BDSFieldEMGlobal::GetField(), BDSFieldMagGlobal::GetField(), GlobalToCurvilinear(), BDSSDApertureImpacts::ProcessHits(), BDSSDEnergyDeposition::ProcessHits(), BDSSDCollimator::ProcessHitsOrdered(), BDSSDEnergyDeposition::ProcessHitsTrack(), BDSIntegratorKickerThin::Stepper(), and BDSIntegratorMultipoleThin::Stepper().
G4ThreeVector BDSAuxiliaryNavigator::ConvertToLocalNoSetup | ( | const G4ThreeVector & | globalPosition, |
const G4bool | useCurvilinear = true |
||
) | const |
Similar to above function but does NOT initialise the transforms.
Definition at line 232 of file BDSAuxiliaryNavigator.cc.
References GlobalToLocal().
BDSStep BDSAuxiliaryNavigator::CurvilinearToGlobal | ( | const G4double | fieldArcLength, |
const G4ThreeVector & | unitField, | ||
const G4double | angle, | ||
const G4ThreeVector & | CLPosition, | ||
const G4ThreeVector & | CLMomentum, | ||
const G4bool | useCurvilinearWorld, | ||
const G4double | FCof, | ||
const G4double | tilt = 0 |
||
) |
Definition at line 395 of file BDSAuxiliaryNavigator.cc.
BDSStep BDSAuxiliaryNavigator::CurvilinearToGlobal | ( | const G4ThreeVector & | localPosition, |
const G4ThreeVector & | localMomentum, | ||
const G4bool | useCurvilinearWorld | ||
) |
Definition at line 388 of file BDSAuxiliaryNavigator.cc.
BDSStep BDSAuxiliaryNavigator::GlobalToCurvilinear | ( | const G4double | fieldArcLength, |
const G4ThreeVector & | unitField, | ||
const G4double | angle, | ||
const G4ThreeVector & | position, | ||
const G4ThreeVector & | unitMomentum, | ||
const G4double | h, | ||
const G4bool | useCurvilinearWorld, | ||
const G4double | FCof, | ||
const G4double | tilt = 0 |
||
) |
Convert to curvilinear coordinates. Angle supplied separately in case of over/underpowered magnets.
Definition at line 312 of file BDSAuxiliaryNavigator.cc.
References ConvertToLocal(), BDS::IsFinite(), BDSStep::PostStepPoint(), BDSStep::PreStepPoint(), and BDSStep::rotateZ().
Referenced by BDSIntegratorDipoleFringe::BaseStepper(), BDSIntegratorDipoleFringe::MultipoleStep(), BDSIntegratorDipoleQuadrupole::Stepper(), BDSIntegratorQuadrupole::Stepper(), and BDSIntegratorSolenoid::Stepper().
BDSStep BDSAuxiliaryNavigator::GlobalToCurvilinear | ( | const G4ThreeVector & | position, |
const G4ThreeVector & | unitMomentum, | ||
const G4double | h, | ||
const G4bool | useCurvilinearWorld | ||
) |
Definition at line 304 of file BDSAuxiliaryNavigator.cc.
|
private |
Utility function to select appropriate transform.
Definition at line 460 of file BDSAuxiliaryNavigator.cc.
Referenced by ConvertAxisToLocal(), ConvertToLocal(), and ConvertToLocalNoSetup().
|
private |
Definition at line 470 of file BDSAuxiliaryNavigator.cc.
|
private |
Locate the supplied point the in the geometry and get and store the transform to that volume in the member variable. This function has to be const as it's called the first time in GetField which is a pure virtual const function from G4MagneticField that we have to implement and have to keep const. This function doesn't change the const pointer but does change the contents of what it points to.
Definition at line 493 of file BDSAuxiliaryNavigator.cc.
References auxNavigator, and auxNavigatorCL.
|
private |
This is used to forcibly initialise the transforms using a position, momentum vector and step length. The free drift of the particle is calculated and the average of the two points is used to locate and initialise the transforms (in global coordinates).
Definition at line 503 of file BDSAuxiliaryNavigator.cc.
|
private |
Utility function to select appropriate transform.
Definition at line 465 of file BDSAuxiliaryNavigator.cc.
Referenced by ConvertAxisToGlobal(), ConvertToGlobal(), and ConvertToGlobalStep().
G4VPhysicalVolume * BDSAuxiliaryNavigator::LocateGlobalPointAndSetup | ( | const G4ThreeVector & | point, |
const G4ThreeVector * | direction = nullptr , |
||
const G4bool | pRelativeSearch = true , |
||
const G4bool | ignoreDirection = true , |
||
G4bool | useCurvilinear = true |
||
) | const |
A wrapper for the underlying static navigator instance located within this class.
Definition at line 68 of file BDSAuxiliaryNavigator.cc.
References auxNavigatorCLB, curvilinearWorldPV, Navigator(), volumeMargin, and worldPV.
Referenced by ConvertToLocal(), and BDSOutputROOTEventTrajectory::FillIndividualTrajectory().
G4VPhysicalVolume * BDSAuxiliaryNavigator::LocateGlobalPointAndSetup | ( | G4Step const *const | step, |
G4bool | useCurvilinear = true |
||
) | const |
A safe way to locate and setup a point inside a volume. This uses the mid point of the step to locate the volume rather than the edges which may lie on a boundary and typically find the world or previous volume.
Definition at line 117 of file BDSAuxiliaryNavigator.cc.
References auxNavigatorCLB, curvilinearWorldPV, Navigator(), volumeMargin, and worldPV.
|
private |
Utility function to select appropriate navigator.
Definition at line 454 of file BDSAuxiliaryNavigator.cc.
References auxNavigator, and auxNavigatorCL.
Referenced by LocateGlobalPointAndSetup().
|
inlinestatic |
Definition at line 79 of file BDSAuxiliaryNavigator.hh.
|
static |
Definition at line 61 of file BDSAuxiliaryNavigator.cc.
|
staticprotected |
Navigator object for safe navigation in the real (mass) world without affecting tracking of the particle.
Definition at line 229 of file BDSAuxiliaryNavigator.hh.
Referenced by AttachWorldVolumeToNavigator(), InitialiseTransform(), Navigator(), and BDSIntegratorEulerOld::Stepper().
|
staticprotected |
Navigator object for curvilinear world that contains simple cylinders for each element whose local coordinates represent the curvilinear coordinate system.
Definition at line 234 of file BDSAuxiliaryNavigator.hh.
Referenced by AttachWorldVolumeToNavigatorCL(), InitialiseTransform(), and Navigator().
|
staticprotected |
Navigator object for bridge world. This contains bridging volumes for the gaps in the curivlinear world. It therefore acts as a fall back if we find the world volume when we know we really shouldn't.
Definition at line 239 of file BDSAuxiliaryNavigator.hh.
Referenced by LocateGlobalPointAndSetup().
|
mutableprotected |
Definition at line 225 of file BDSAuxiliaryNavigator.hh.
|
staticprivate |
Cache of world PV to test if we're getting the wrong volume for the transform.
Definition at line 276 of file BDSAuxiliaryNavigator.hh.
|
staticprivate |
Cache of world PV to test if we're getting the wrong volume for the transform.
Definition at line 275 of file BDSAuxiliaryNavigator.hh.
Referenced by AttachWorldVolumeToNavigatorCL(), and LocateGlobalPointAndSetup().
|
mutableprotected |
Definition at line 221 of file BDSAuxiliaryNavigator.hh.
|
mutableprotected |
Definition at line 223 of file BDSAuxiliaryNavigator.hh.
|
mutableprotected |
Definition at line 222 of file BDSAuxiliaryNavigator.hh.
|
mutableprotected |
Definition at line 224 of file BDSAuxiliaryNavigator.hh.
|
staticprivate |
Counter to keep track of when the last instance of the class is deleted and therefore when the navigators can be safely deleted without affecting
Definition at line 271 of file BDSAuxiliaryNavigator.hh.
|
private |
Margin by which to advance the point along the step direction if the world volume is found for transforms. This is in an attempt to find a real curvilinear volume. Therefore, this should be greater than lengthSafety or the geometrical tolerance.
Definition at line 283 of file BDSAuxiliaryNavigator.hh.
Referenced by LocateGlobalPointAndSetup().
|
staticprivate |
Cache of world PV to test if we're getting the wrong volume for the transform.
Definition at line 274 of file BDSAuxiliaryNavigator.hh.
Referenced by AttachWorldVolumeToNavigator(), and LocateGlobalPointAndSetup().