#include <BLCoordinates.hh>
Whenever an element is placed into the world group, the current centerline coordinates are normally used, rather than global coordinates (the local coordinates of the world volume/group).
Centerline coordinates are X,Y,Z, with Z along the nominal beam direction, X is beam left, and Y is beam up. The centerline itself is always X=0,Y=0, and consists of a sequence of segmentCLs each of which is straight, and the centerline is continuous at each corner (bend in the centerline). Except within bending magnets the reference particle is supposed to travel directly down the centerline. Z is continuous along the centerline, including at corners; the first segmentCL extends to Z=-infinity, and the last segmentCL extends to Z=+infinity.
Note that the coordinate transform global->centerline is NOT unique -- near corners both segmentCLs may refer to a given point. The transform centerline->global is unique for a given segmentCL. The centerline itself is continuous at corners, but tracks will in general have a discontinuity in their centerline coordinates when transitioning from one segmentCL to the next. A track will transition from one segmentCL to the next when either its Z coordinate exceeds that of the initial segmentCL's end (corner), or when its radius exceeds the radiusCut for the initial segmentCL.
Initially the centerline coordinates are identical to global coordinates. The start() function will change the starting location and orientation of the first segmentCL, and the corner() function ends the current segmentCL and begins a new one. start() must not be called after any element has been placed in the first segmentCL, or after corner() has been called.
Instances of this class represent the current centerline coordinates of a track, and are linked into the track via G4Track::SetUserInformation() (an instance is registered with BLManager to manage the instance linked into the track, and they can also be used for individual points, if necessary). Static members and functions are used to access and create the centerline segmentCLs.
Because of the need to transition tracks from one segmentCL to the next, corners with acute angles are prohibited.
NOTE: BLManager expects G4Track::GetUserInformation() to return a pointer to a BLCoordinates object, or to an object of a class derived from it. This means that any G4VUserTrackInformation class used in a simulation using BLManager should be derived from this class, and not directly from G4VUserTrackInformation (but verify is used to detect if it is not a BLCoordinates object).
Note the BLManager::UserSteppingAction() manages the updates to the current track's BLCoordinates before doing anything else.
NOTE: "CL" in a vairable means it is for centerline coordinates; "RC" means it is for reference coordinates.
Public Member Functions | |
BLCoordinates () | |
Default constructor -- for UserSteppingAction or without point. | |
BLCoordinates (const G4double _global[4]) | |
Constructor -- finds the earliest segmentCL containing the point (global coordinates). | |
virtual | ~BLCoordinates () |
Destructor. | |
void * | operator new (size_t) |
operator new - using G4Allocator. NOTE: any derived class MUST re-implement this. | |
void | operator delete (void *p) |
operator delete - using G4Allocator. NOTE: any derived class MUST re-implement this. | |
G4bool | isValid () const |
isValid() verifies that this is a valid BLCoordinates instance. NOTE: this is not virtual, and cannot be made virtual! | |
G4bool | setGlobal (const G4double _global[4]) |
setGlobal() updates the centerline coordinates for the new point[] of the track; uses the current segmentCL #, and will update it if the track crosses into a new segmentCL. returns false if the track fails radiusCut in its current and next segmentCL (or all segmentCL if ring is true). Also updates reference coordinates, if validRC is true and the centerline radiusCut succeeded. | |
G4bool | setGlobal (const G4ThreeVector &pos, G4double time) |
void | getCenterlineCoords (G4double local[4]) const |
getCenterlineCoords() returns the centerline coordinates. DEPRECATED: use getCoords() instead. | |
void | getCenterlineCoords (G4ThreeVector &local) const |
void | getCoords (BLCoordinateType type, G4double local[4]) |
getCoords() returns the appropriate coordinate values. | |
void | getCoords (BLCoordinateType type, G4ThreeVector &local) |
G4double | getCLZ () |
getCLZ() returns the centerline z coordinate. | |
G4RotationMatrix & | getRotation () |
getRotation() returns the appropriate coordinate rotation for the previous call to getCoords(). Given a G4ThreeVector B in global coordinates, this transforms it to centerline components: coord->getCoords(BLCOORD_CENTERLINE,local); B = coord->getRotation() * B; | |
G4RotationMatrix & | getRotation (BLCoordinateType type) |
getRotation() returns the rotation for the specified coordinates. Given a G4ThreeVector B in global coordinates, this transforms it to centerline components: B = coord->getRotation(BLCOORD_CENTERLINE) * B; | |
G4RotationMatrix | getCenterlineRotation () |
getCenterlineRotation() returns the rotation matrix used for the last setGlobal() of this instance. This is the rotation of the centerline wrt global coordinate axes. DEPRECATED -- use getRotation(). | |
G4RotationMatrix | getTransformRotation () |
getTransformRotation() returns the rotation matrix of the global -> centerline coordinate transformation; it is the inverse of getCenterlineRotation(). Given a G4ThreeVector B in global coordinates, this transforms it to centerline components: B = getTransformRotation() * B; DEPRECATED -- use getRotation(). | |
int | getSegmentCL () |
getSegmentCL() re turns the CL segment # (for debugging). | |
virtual void | Print () const |
Print from G4VUserTrackInformation. | |
Static Public Member Functions | |
static double | getCurrentRadiusCut () |
getCurrentRadiusCut() returns the radiusCut for the current CL segment. | |
static void | update (G4Track *track) |
update() This function ensures that a UserTrackInformation object is linked into the track; if none is, it creates a new BLCoordinates instance and links it in. It also calls setGlobal() for the BLCoordinates linked into the track. | |
static void | start (const G4double point[4], G4double z, G4RotationMatrix *rotation=0, G4double radiusCut=0.0, G4bool _ring=false) |
start() sets the starting point and direction for the first segmentCL. Must not be called after corner() is called, or after any element has been placed using centerline coordinates. point[] is in global coordinates, z is the starting value of z in CENTERLINE coordinates, _ring indicates whether a ring is present (so jumps in segmentCL # are permitted), and rotation should be the rotation of the first segmentCL's centerline wrt the global coordinate axes (i.e. it is an active rotation of the centerline as an object wrt the global coordinate axes). | |
static void | corner (G4double z, G4RotationMatrix *rotation, G4double radiusCut) |
corner() creates a corner in the centerline by ending the current segmentCL and beginning a new one. z is the position in CENTERLINE coordinates of the current segmentCL, and rotation is relative to the current CENTERLINE coordinates (i.e. as an active rotation of the centerline as an object wrt the CENTERLINE coordinate axes of the current segmentCL). z must be greater than the centerline z of the previous corner or start. rotation can be an identity matrix (e.g. to merely set a new radiusCut). | |
static BLCoordinateTransform & | getCurrentTransform () |
getCurrentTransform() returns the current BLCoordinateTransform for the centerline coordinates (i.e. the coordinate transform global->centerline of the current segmentCL). | |
static G4RotationMatrix * | getCurrentRotation () |
getCurrentRotation() returns the current G4RotationMatrix for the centerline coordinates (i.e. the active rotation of the current segmentCL, as an object rotated wrt global coordinate axes). Given a G4ThreeVector B in global coordinates, this transforms it to centerline components: B = *getCurrentRotation() * B; | |
static void | getCurrentGlobal (const G4double local[4], G4double _global[4]) |
getCurrentGlobal() returns global coordinates from centerline coordinates, using the current segmentCL. | |
static void | getCurrentGlobal (const G4ThreeVector &local, G4ThreeVector &_global) |
getCurrentGlobal() returns global coordinates from centerline coordinates, using the current segmentCL. | |
static void | getGlobalAnywhere (const G4ThreeVector &local, G4ThreeVector &_global) |
getGlobalAnywhere() returns global coordinates from centerline coordinates, using the lowest segmentCL containing local. | |
static G4RotationMatrix * | getGlobalRotation (const G4ThreeVector &local) |
getGlobalRotation() returns the local->global rotation matrix coordinates, using the lowest segmentCL containing local. | |
static G4double | getRadiusCut (const G4ThreeVector &local) |
getRadiusCut() returns the radius cut for the lowest segmentCL containing local. | |
static BLCoordinateType | getCoordinateType (G4String s) |
getCoordinateType() converts a G4String name to a BLCoordinateType. Names are: global, centerline, reference; only the first letter is used, case insensitive. Also enables the generation of REference Coordinates, if they are selected by s. | |
static void | useReferenceCoordinates () |
useReferenceCoordinates() will arrange to fill segmentRCVector[] with reference coordinate segments as the reference particle is tracked. Must be called before the reference particle is tracked. Unless this is called, reference coordinates will not be available. Multiple calls are the same as one. | |
static void | setMostRecentReferenceZ (G4double z) |
setMostRecentReferenceZ() sets the most recent z value of the reference coords. VALID ONLY WHILE TRACKING REFERENCE PARTICLE! | |
static G4double | getMostRecentReferenceZ () |
getMostRecentReferenceZ() returns the most recent z value of the reference coords. VALID ONLY WHILE TRACKING REFERENCE PARTICLE! | |
Private Types | |
enum | { VERIFY = 0x13579246 } |
enum | State { UNUSED, CONSTRUCTING, TRACKING } |
Private Member Functions | |
G4bool | findCLSegment (const G4double _global[4]) |
void | updateReferenceCoordinates () |
updateReferenceCoordinates() will update the reference coordinates, given the centerline coordintes are already updated and valid. Called by setGlobal(). | |
Static Private Member Functions | |
static void | init () |
Private Attributes | |
int | verify |
unsigned int | segmentCL |
unsigned int | segmentRC |
BLCoordinateType | prevType |
G4double | global [4] |
G4double | centerline [4] |
G4double | reference [4] |
int | segmentCLRadiusFailure |
int | referenceCLsegment |
BLCoordinateTransform * | referenceTransformUsed |
Static Private Attributes | |
static G4Allocator< BLCoordinates > | allocator |
static std::vector< CLSegment > | segmentCLVector |
static unsigned | currentCL = 0 |
static G4bool | ring = false |
static State | state = UNUSED |
static G4bool | steppingVerbose = false |
static G4bool | useRC = false |
static G4bool | validRC = false |
static std::vector< RCSegment > | segmentRCVector |
static int | currentRC |
static G4double | tolerance = 1.0e-4 |
static G4double | mostRecentReferenceZ = 0 |
Classes | |
struct | CLSegment |
struct | RCSegment |
class | ReferenceCoordinates |
anonymous enum [private] |
enum BLCoordinates::State [private] |
BLCoordinates::BLCoordinates | ( | ) |
Default constructor -- for UserSteppingAction or without point.
References centerline, global, reference, referenceCLsegment, referenceTransformUsed, segmentCL, segmentCLRadiusFailure, segmentRC, VERIFY, and verify.
Referenced by update().
00064 : G4VUserTrackInformation() 00065 { 00066 // no state change, yet 00067 verify=VERIFY; 00068 segmentCL=99999999; 00069 segmentRC=99999999; 00070 referenceCLsegment = -1; 00071 referenceTransformUsed = 0; 00072 global[0]=global[1]=global[2]=global[3]=0.0; 00073 centerline[0]=centerline[1]=centerline[2]=centerline[3]=0.0; 00074 reference[0]=reference[1]=reference[2]=reference[3]=0.0; 00075 segmentCLRadiusFailure = -1; 00076 }
BLCoordinates::BLCoordinates | ( | const G4double | _global[4] | ) |
Constructor -- finds the earliest segmentCL containing the point (global coordinates).
References findCLSegment(), state, TRACKING, VERIFY, and verify.
00078 : 00079 G4VUserTrackInformation() 00080 { 00081 state = TRACKING; 00082 verify = VERIFY; 00083 findCLSegment(_global); 00084 }
G4bool BLCoordinates::findCLSegment | ( | const G4double | _global[4] | ) | [private] |
References centerline, segmentCL, segmentCLRadiusFailure, and segmentCLVector.
Referenced by BLCoordinates(), and setGlobal().
00093 { 00094 for(segmentCL=0; segmentCL<segmentCLVector.size(); ++segmentCL) { 00095 segmentCLVector[segmentCL].transform.getLocal(centerline,_global); 00096 //if(steppingVerbose) printf("FindCLSegment: seg=%d minZ=%.3f maxZ=%.3f Zcl=%.3f Rcl=%.3f rCut=%.3f\n",segmentCL,segmentCLVector[segmentCL].minZ,segmentCLVector[segmentCL].maxZ,centerline[2],sqrt(centerline[0]*centerline[0]+centerline[1]*centerline[1]),segmentCLVector[segmentCL].radiusCut); 00097 // note no check against minZ -- segmentCLs are ordered 00098 if(centerline[2] <= segmentCLVector[segmentCL].maxZ && 00099 (segmentCLVector[segmentCL].radiusCut <= 0.0 || 00100 (int)segmentCL == segmentCLRadiusFailure || 00101 sqrt(centerline[0]*centerline[0]+centerline[1]*centerline[1]) <= 00102 segmentCLVector[segmentCL].radiusCut)) { 00103 return true; 00104 } 00105 } 00106 --segmentCL; 00107 return false; 00108 }
void* BLCoordinates::operator new | ( | size_t | ) | [inline] |
operator new - using G4Allocator. NOTE: any derived class MUST re-implement this.
Reimplemented in BLTrackInfo.
References allocator.
00119 { 00120 return (void *)allocator.MallocSingle(); 00121 }
void BLCoordinates::operator delete | ( | void * | p | ) | [inline] |
operator delete - using G4Allocator. NOTE: any derived class MUST re-implement this.
Reimplemented in BLTrackInfo.
References allocator.
00125 { 00126 allocator.FreeSingle((BLCoordinates *)p); 00127 }
G4bool BLCoordinates::isValid | ( | ) | const |
isValid() verifies that this is a valid BLCoordinates instance. NOTE: this is not virtual, and cannot be made virtual!
References VERIFY, and verify.
Referenced by TraceNTuple::appendTrace(), BLCMDmovie::beamStep(), BLManager::getExternalParentID(), BLManager::getExternalTrackID(), BLCMDbeamlossntuple::PostUserTrackingAction(), BLManager::PreUserTrackingAction(), BLCoordinates::ReferenceCoordinates::PreUserTrackingAction(), BLCMDnewparticlentuple::PreUserTrackingAction(), BLCMDmovie::referenceStep(), BLManager::setExternalTrackID(), BLEvaluator::setTrackVariables(), BLManager::steppingVerbosePrint(), update(), BLManager::UserSteppingAction(), BLVirtualDetectorNTuple::UserSteppingAction(), TrackerPlaneInstance::UserSteppingAction(), TimeNTuple::UserSteppingAction(), BLCMDzntuple::Entry::UserZSteppingAction(), BLCMDtracker::UserZSteppingAction(), and BLCMDprofile::Entry::UserZSteppingAction().
G4bool BLCoordinates::setGlobal | ( | const G4double | _global[4] | ) |
setGlobal() updates the centerline coordinates for the new point[] of the track; uses the current segmentCL #, and will update it if the track crosses into a new segmentCL. returns false if the track fails radiusCut in its current and next segmentCL (or all segmentCL if ring is true). Also updates reference coordinates, if validRC is true and the centerline radiusCut succeeded.
References centerline, findCLSegment(), global, ring, segmentCL, segmentCLRadiusFailure, segmentCLVector, updateReferenceCoordinates(), and validRC.
Referenced by BLCMDmovie::dumpGeometry(), BLManager::PreUserTrackingAction(), setGlobal(), and update().
00111 { 00112 global[0] = _global[0]; 00113 global[1] = _global[1]; 00114 global[2] = _global[2]; 00115 global[3] = _global[3]; 00116 00117 if(segmentCL >= segmentCLVector.size()) { 00118 //printf("setGlobal-0: calling findCLSegment\n"); 00119 G4bool ret = findCLSegment(_global); 00120 if(ret && validRC) updateReferenceCoordinates(); 00121 //if(!ret) printf("setGlobal-0: findCLSegment FAILED\n"); 00122 return ret; 00123 } 00124 00125 segmentCLVector[segmentCL].transform.getLocal(centerline,_global); 00126 //printf("setGlobal-1A: global=%.1f,%.1f,%.1f cl=%.1f,%.1f,%.1f\n",_global[0],_global[1],_global[2],centerline[0],centerline[1],centerline[2]); 00127 //printf("setGlobal-1: seg=%d minZ=%.3f maxZ=%.3f Zcl=%.3f Rcl=%.3f rCut=%.3f\n",segmentCL,segmentCLVector[segmentCL].minZ,segmentCLVector[segmentCL].maxZ,centerline[2],sqrt(centerline[0]*centerline[0]+centerline[1]*centerline[1]),segmentCLVector[segmentCL].radiusCut); 00128 while(centerline[2] > segmentCLVector[segmentCL].maxZ || 00129 (segmentCLVector[segmentCL].radiusCut > 0.0 && 00130 (int)segmentCL != segmentCLRadiusFailure && 00131 sqrt(centerline[0]*centerline[0]+centerline[1]*centerline[1]) > 00132 segmentCLVector[segmentCL].radiusCut)) { 00133 if(++segmentCL >= segmentCLVector.size()) { 00134 if(ring && findCLSegment(_global)) { 00135 segmentCLVector[segmentCL].transform.getLocal(centerline, 00136 _global); 00137 if(validRC) updateReferenceCoordinates(); 00138 return true; 00139 } 00140 --segmentCL; 00141 //printf("setGlobal-1a: FAILED\n"); 00142 return false; 00143 } 00144 segmentCLVector[segmentCL].transform.getLocal(centerline,_global); 00145 //printf("setGlobal-2: seg=%d minZ=%.3f maxZ=%.3f Zcl=%.3f Rcl=%.3f rCut=%.3f\n",segmentCL,segmentCLVector[segmentCL].minZ,segmentCLVector[segmentCL].maxZ,centerline[2],sqrt(centerline[0]*centerline[0]+centerline[1]*centerline[1]),segmentCLVector[segmentCL].radiusCut); 00146 /*** 00147 // note no check against minZ -- segmentCLs are ordered 00148 if(centerline[2] > segmentCLVector[segmentCL].maxZ || 00149 (segmentCLVector[segmentCL].radiusCut > 0.0 && 00150 segmentCL != segmentCLRadiusFailure && 00151 sqrt(centerline[0]*centerline[0]+centerline[1]*centerline[1]) > 00152 segmentCLVector[segmentCL].radiusCut)) { 00153 if(ring && findCLSegment(_global)) { 00154 segmentCLVector[segmentCL].transform.getLocal(centerline, 00155 _global); 00156 if(validRC) updateReferenceCoordinates(); 00157 return true; 00158 } 00159 //printf("setGlobal-2a: FAILED\n"); 00160 return false; 00161 } 00162 ***/ 00163 } 00164 00165 if(validRC) updateReferenceCoordinates(); 00166 return true; 00167 }
G4bool BLCoordinates::setGlobal | ( | const G4ThreeVector & | pos, | |
G4double | time | |||
) |
References global, and setGlobal().
00170 { 00171 global[0]=pos[0],global[1]=pos[1],global[2]=pos[2]; 00172 global[3]=time; 00173 return setGlobal(global); 00174 }
void BLCoordinates::getCenterlineCoords | ( | G4double | local[4] | ) | const [inline] |
getCenterlineCoords() returns the centerline coordinates. DEPRECATED: use getCoords() instead.
References centerline.
00146 { local[0]=centerline[0]; local[1]=centerline[1]; 00147 local[2]=centerline[2]; local[3]=centerline[3]; }
void BLCoordinates::getCenterlineCoords | ( | G4ThreeVector & | local | ) | const [inline] |
References centerline.
00149 { local[0]=centerline[0]; local[1]=centerline[1]; 00150 local[2]=centerline[2]; }
void BLCoordinates::getCoords | ( | BLCoordinateType | type, | |
G4double | local[4] | |||
) |
getCoords() returns the appropriate coordinate values.
References BLCOORD_CENTERLINE, BLCOORD_GLOBAL, BLCOORD_REFERENCE, centerline, BLCoordinateTransform::getLocal(), BLCoordinateTransform::getRotation(), global, prevType, reference, referenceTransformUsed, segmentCL, segmentCLVector, segmentRC, segmentRCVector, and validRC.
Referenced by TraceNTuple::appendTrace(), BLCMDmovie::beamStep(), BLCMDmovie::dumpGeometry(), getCoords(), BLCMDbeamlossntuple::PostUserTrackingAction(), BLCoordinates::ReferenceCoordinates::PreUserTrackingAction(), BLCMDnewparticlentuple::PreUserTrackingAction(), BLCMDmovie::referenceStep(), BLEvaluator::setTrackVariables(), BLManager::steppingVerbosePrint(), BLVirtualDetectorNTuple::UserSteppingAction(), TrackerPlaneInstance::UserSteppingAction(), TimeNTuple::UserSteppingAction(), BLCMDzntuple::Entry::UserZSteppingAction(), BLCMDtune::UserZSteppingAction(), BLCMDtracker::UserZSteppingAction(), and BLCMDprofile::Entry::UserZSteppingAction().
00195 { 00196 prevType = type; 00197 00198 switch(type) { 00199 case BLCOORD_GLOBAL: 00200 local[0] = global[0]; local[1] = global[1]; 00201 local[2] = global[2]; local[3] = global[3]; 00202 //printf("getCoords(GLOBAL): %.1f,%.1f,%.1f\n",local[0],local[1],local[2]); 00203 return; 00204 case BLCOORD_CENTERLINE: 00205 assert(segmentCL < segmentCLVector.size()); 00206 local[0] = centerline[0]; local[1] = centerline[1]; 00207 local[2] = centerline[2]; local[3] = centerline[3]; 00208 //printf("getCoords(CL): %.1f,%.1f,%.1f\n",local[0],local[1],local[2]); 00209 return; 00210 case BLCOORD_REFERENCE: 00211 assert(validRC); 00212 // interpolate 00213 if(reference[2] < segmentRCVector[segmentRC].interpolatePrev && 00214 segmentRC > 0) { 00215 referenceTransformUsed = 00216 &segmentRCVector[segmentRC-1].interpolate; 00217 referenceTransformUsed->getLocal(reference,global); 00218 } else if(reference[2] > 00219 segmentRCVector[segmentRC].interpolateThis) { 00220 referenceTransformUsed = 00221 &segmentRCVector[segmentRC].interpolate; 00222 referenceTransformUsed->getLocal(reference,global); 00223 } 00224 G4ThreeVector axis = referenceTransformUsed->getRotation().getAxis(); 00225 G4double angle = referenceTransformUsed->getRotation().getDelta(); 00226 //printf("segmentRC=%d minZ=%.4g maxZ=%.4g global=%.1f,%.1f,%.1f ref=%.1f,%.1f,%.1f axis=%.1f,%.1f,%.1f angle(deg)=%.1f\n",segmentRC,segmentRCVector[segmentRC].minZ,segmentRCVector[segmentRC].maxZ,global[0],global[1],global[2],reference[0],reference[1],reference[2],axis[0],axis[1],axis[2],angle/deg); 00227 local[0] = reference[0]; local[1] = reference[1]; 00228 local[2] = reference[2]; local[3] = reference[3]; 00229 //printf("getCoords(REF): %.1f,%.1f,%.1f\n",local[0],local[1],local[2]); 00230 return; 00231 } 00232 G4Exception("BLCoordinates","Invalid Coordinate Type",FatalException, 00233 ""); 00234 }
void BLCoordinates::getCoords | ( | BLCoordinateType | type, | |
G4ThreeVector & | local | |||
) |
References getCoords().
00237 { 00238 G4double tmp[4]; 00239 getCoords(type,tmp); 00240 local[0] = tmp[0]; 00241 local[1] = tmp[1]; 00242 local[2] = tmp[2]; 00243 //printf("getCoords(): %.1f,%.1f,%.1f\n",local[0],local[1],local[2]); 00244 }
G4double BLCoordinates::getCLZ | ( | ) | [inline] |
getCLZ() returns the centerline z coordinate.
References centerline.
Referenced by BLManager::PreUserTrackingAction(), and BLManager::UserSteppingAction().
00157 { return centerline[2]; }
G4RotationMatrix& BLCoordinates::getRotation | ( | ) | [inline] |
getRotation() returns the appropriate coordinate rotation for the previous call to getCoords(). Given a G4ThreeVector B in global coordinates, this transforms it to centerline components: coord->getCoords(BLCOORD_CENTERLINE,local); B = coord->getRotation() * B;
References prevType.
Referenced by TraceNTuple::appendTrace(), BLCMDmovie::beamStep(), BLCMDbeamlossntuple::PostUserTrackingAction(), BLCMDnewparticlentuple::PreUserTrackingAction(), BLEvaluator::setTrackVariables(), BLManager::steppingVerbosePrint(), BLVirtualDetectorNTuple::UserSteppingAction(), TrackerPlaneInstance::UserSteppingAction(), TimeNTuple::UserSteppingAction(), BLCMDzntuple::Entry::UserZSteppingAction(), BLCMDtracker::UserZSteppingAction(), and BLCMDprofile::Entry::UserZSteppingAction().
00165 { return getRotation(prevType); }
G4RotationMatrix & BLCoordinates::getRotation | ( | BLCoordinateType | type | ) |
getRotation() returns the rotation for the specified coordinates. Given a G4ThreeVector B in global coordinates, this transforms it to centerline components: B = coord->getRotation(BLCOORD_CENTERLINE) * B;
References BLCOORD_CENTERLINE, BLCOORD_GLOBAL, BLCOORD_REFERENCE, BLCoordinateTransform::getRotation(), IdentityMatrix, referenceTransformUsed, segmentCL, segmentCLVector, and validRC.
00247 { 00248 switch(type) { 00249 case BLCOORD_GLOBAL: 00250 return IdentityMatrix; 00251 case BLCOORD_CENTERLINE: 00252 assert(segmentCL < segmentCLVector.size()); 00253 return segmentCLVector[segmentCL].transform.getRotation(); 00254 case BLCOORD_REFERENCE: 00255 assert(validRC && referenceTransformUsed != 0); 00256 return referenceTransformUsed->getRotation(); 00257 } 00258 G4Exception("BLCoordinates","Invalid Coordinate Type",FatalException, 00259 ""); 00260 return IdentityMatrix; 00261 }
G4RotationMatrix BLCoordinates::getCenterlineRotation | ( | ) | [inline] |
getCenterlineRotation() returns the rotation matrix used for the last setGlobal() of this instance. This is the rotation of the centerline wrt global coordinate axes. DEPRECATED -- use getRotation().
References segmentCL, and segmentCLVector.
00177 { 00178 if(segmentCL >= segmentCLVector.size()) 00179 return G4RotationMatrix(); 00180 return segmentCLVector[segmentCL].transform.getRotationInverse(); 00181 }
static double BLCoordinates::getCurrentRadiusCut | ( | ) | [inline, static] |
getCurrentRadiusCut() returns the radiusCut for the current CL segment.
References currentCL, and segmentCLVector.
Referenced by BLCMDcorner::command().
00185 { 00186 if(currentCL >= segmentCLVector.size()) 00187 return DBL_MAX; 00188 else 00189 return segmentCLVector[currentCL].radiusCut; 00190 }
G4RotationMatrix BLCoordinates::getTransformRotation | ( | ) | [inline] |
getTransformRotation() returns the rotation matrix of the global -> centerline coordinate transformation; it is the inverse of getCenterlineRotation(). Given a G4ThreeVector B in global coordinates, this transforms it to centerline components: B = getTransformRotation() * B; DEPRECATED -- use getRotation().
References segmentCL, and segmentCLVector.
00199 { 00200 if(segmentCL >= segmentCLVector.size()) 00201 return G4RotationMatrix(); 00202 return segmentCLVector[segmentCL].transform.getRotation(); 00203 }
int BLCoordinates::getSegmentCL | ( | ) | [inline] |
getSegmentCL() re turns the CL segment # (for debugging).
References segmentCL.
Referenced by BLManager::steppingVerbosePrint().
00206 { return segmentCL; }
void BLCoordinates::update | ( | G4Track * | track | ) | [static] |
update() This function ensures that a UserTrackInformation object is linked into the track; if none is, it creates a new BLCoordinates instance and links it in. It also calls setGlobal() for the BLCoordinates linked into the track.
References BLCoordinates(), isValid(), segmentCL, segmentCLRadiusFailure, setGlobal(), and steppingVerbose.
Referenced by BLManager::UserSteppingAction().
00177 { 00178 // ensure there is a BLCoordiantes linked into the track 00179 BLCoordinates *c = (BLCoordinates *)track->GetUserInformation(); 00180 if(!c || !c->isValid()) { 00181 c = new BLCoordinates(); 00182 track->SetUserInformation(c); 00183 } 00184 00185 // now update the object and apply radiusCut 00186 if(!c->setGlobal(track->GetPosition(),track->GetGlobalTime())) { 00187 track->SetTrackStatus(fStopAndKill); 00188 if(steppingVerbose) 00189 printf("Track Centerline radius exceeds radiusCut\n"); 00190 c->segmentCLRadiusFailure = c->segmentCL; 00191 } 00192 }
virtual void BLCoordinates::Print | ( | ) | const [inline, virtual] |
void BLCoordinates::init | ( | ) | [static, private] |
References CONSTRUCTING, currentCL, BLManager::getObject(), BLManager::getSteppingVerbose(), segmentCLVector, state, steppingVerbose, TRACKING, and UNUSED.
Referenced by corner(), getCurrentGlobal(), getCurrentRotation(), getCurrentTransform(), and start().
00264 { 00265 if(state == TRACKING) { 00266 G4Exception("BLCoordinates","Invalid init during tracking", 00267 FatalException, ""); 00268 } 00269 if(state == UNUSED) { 00270 state = CONSTRUCTING; 00271 segmentCLVector.push_back(CLSegment(-DBL_MAX, +DBL_MAX, 0.0, 00272 BLCoordinateTransform())); 00273 currentCL = 0; 00274 steppingVerbose = 00275 BLManager::getObject()->getSteppingVerbose() != 0; 00276 } 00277 }
void BLCoordinates::start | ( | const G4double | point[4], | |
G4double | z, | |||
G4RotationMatrix * | rotation = 0 , |
|||
G4double | radiusCut = 0.0 , |
|||
G4bool | _ring = false | |||
) | [static] |
start() sets the starting point and direction for the first segmentCL. Must not be called after corner() is called, or after any element has been placed using centerline coordinates. point[] is in global coordinates, z is the starting value of z in CENTERLINE coordinates, _ring indicates whether a ring is present (so jumps in segmentCL # are permitted), and rotation should be the rotation of the first segmentCL's centerline wrt the global coordinate axes (i.e. it is an active rotation of the centerline as an object wrt the global coordinate axes).
References currentCL, BLManager::getObject(), BLManager::getSteppingVerbose(), init(), ring, segmentCLVector, state, steppingVerbose, and UNUSED.
Referenced by BLCMDstart::command().
00282 { 00283 if(state != UNUSED || currentCL != 0) { 00284 G4Exception("BLCoordinates","Invalid start",FatalException, ""); 00285 } 00286 init(); 00287 00288 G4ThreeVector pos(0,0,-z); 00289 if(rotation) 00290 pos = *rotation * pos; 00291 pos[0] += _global[0]; 00292 pos[1] += _global[1]; 00293 pos[2] += _global[2]; 00294 segmentCLVector[0].transform = (rotation ? 00295 BLCoordinateTransform(rotation,pos) : 00296 BLCoordinateTransform(pos)); 00297 segmentCLVector[0].radiusCut = radiusCut; 00298 ring = _ring; 00299 steppingVerbose = BLManager::getObject()->getSteppingVerbose() != 0; 00300 }
void BLCoordinates::corner | ( | G4double | z, | |
G4RotationMatrix * | rotation, | |||
G4double | radiusCut | |||
) | [static] |
corner() creates a corner in the centerline by ending the current segmentCL and beginning a new one. z is the position in CENTERLINE coordinates of the current segmentCL, and rotation is relative to the current CENTERLINE coordinates (i.e. as an active rotation of the centerline as an object wrt the CENTERLINE coordinate axes of the current segmentCL). z must be greater than the centerline z of the previous corner or start. rotation can be an identity matrix (e.g. to merely set a new radiusCut).
References currentCL, getCurrentRotation(), getCurrentTransform(), init(), and segmentCLVector.
Referenced by BLCMDcornerarc::command(), and BLCMDcorner::command().
00304 { 00305 init(); 00306 00307 G4double local[4], glob[4]; 00308 local[0]=local[1]=local[3]=0.0; 00309 local[2] = z; 00310 segmentCLVector[currentCL].transform.getGlobal(local,glob); 00311 00312 segmentCLVector[currentCL].maxZ = z; 00313 G4RotationMatrix *p=0; 00314 if(rotation) 00315 // if C=currentCL rot, and R=this rot, we want R referenced to 00316 // the CURRENT centerline coordinates, so we need C R C^-1. 00317 // But we compose this with C (on the right), and get C R. 00318 p = new G4RotationMatrix(*getCurrentRotation() * *rotation); 00319 else if(getCurrentTransform().isRotated()) 00320 p = getCurrentRotation(); 00321 00322 G4ThreeVector pos(0,0,-z); 00323 if(p) 00324 pos = *p * pos; 00325 pos[0] += glob[0]; 00326 pos[1] += glob[1]; 00327 pos[2] += glob[2]; 00328 segmentCLVector.push_back(CLSegment(z, +DBL_MAX, radiusCut, 00329 BLCoordinateTransform(p,pos))); 00330 currentCL = segmentCLVector.size() - 1; 00331 }
BLCoordinateTransform & BLCoordinates::getCurrentTransform | ( | ) | [static] |
getCurrentTransform() returns the current BLCoordinateTransform for the centerline coordinates (i.e. the coordinate transform global->centerline of the current segmentCL).
References currentCL, init(), and segmentCLVector.
Referenced by corner().
00334 { 00335 init(); 00336 return segmentCLVector[currentCL].transform; 00337 }
G4RotationMatrix * BLCoordinates::getCurrentRotation | ( | ) | [static] |
getCurrentRotation() returns the current G4RotationMatrix for the centerline coordinates (i.e. the active rotation of the current segmentCL, as an object rotated wrt global coordinate axes). Given a G4ThreeVector B in global coordinates, this transforms it to centerline components: B = *getCurrentRotation() * B;
References currentCL, init(), and segmentCLVector.
Referenced by BLCMDtracker::command(), BLCMDreference::command(), BLCMDcorner::command(), BLCMDbeam::command(), corner(), and BLGroup::placeElement().
00340 { 00341 init(); 00342 return &segmentCLVector[currentCL].transform.getRotationInverse(); 00343 }
void BLCoordinates::getCurrentGlobal | ( | const G4double | local[4], | |
G4double | _global[4] | |||
) | [static] |
getCurrentGlobal() returns global coordinates from centerline coordinates, using the current segmentCL.
References currentCL, init(), and segmentCLVector.
Referenced by BLCMDtracker::command(), BLCMDreference::command(), BLCMDbeam::command(), and BLGroup::placeElement().
00346 { 00347 init(); 00348 segmentCLVector[currentCL].transform.getGlobal(local,_global); 00349 }
void BLCoordinates::getCurrentGlobal | ( | const G4ThreeVector & | local, | |
G4ThreeVector & | _global | |||
) | [static] |
getCurrentGlobal() returns global coordinates from centerline coordinates, using the current segmentCL.
References currentCL, init(), and segmentCLVector.
00352 { 00353 init(); 00354 segmentCLVector[currentCL].transform.getGlobal(local,_global); 00355 }
void BLCoordinates::getGlobalAnywhere | ( | const G4ThreeVector & | local, | |
G4ThreeVector & | _global | |||
) | [static] |
getGlobalAnywhere() returns global coordinates from centerline coordinates, using the lowest segmentCL containing local.
References segmentCLVector.
00359 { 00360 for(unsigned int seg=0; seg<segmentCLVector.size(); ++seg) { 00361 // note no check against minZ -- segmentCLs are ordered 00362 if(local[2] <= segmentCLVector[seg].maxZ && 00363 (segmentCLVector[seg].radiusCut <= 0.0 || 00364 sqrt(local[0]*local[0]+local[1]*local[1]) <= 00365 segmentCLVector[seg].radiusCut)) { 00366 segmentCLVector[seg].transform.getGlobal(local,_global); 00367 return; 00368 } 00369 } 00370 }
G4RotationMatrix * BLCoordinates::getGlobalRotation | ( | const G4ThreeVector & | local | ) | [static] |
getGlobalRotation() returns the local->global rotation matrix coordinates, using the lowest segmentCL containing local.
References segmentCLVector.
00373 { 00374 for(unsigned int seg=0; seg<segmentCLVector.size(); ++seg) { 00375 // note no check against minZ -- segmentCLs are ordered 00376 if(local[2] <= segmentCLVector[seg].maxZ && 00377 (segmentCLVector[seg].radiusCut <= 0.0 || 00378 sqrt(local[0]*local[0]+local[1]*local[1]) <= 00379 segmentCLVector[seg].radiusCut)) { 00380 return &segmentCLVector[seg].transform.getRotation(); 00381 } 00382 } 00383 00384 return 0; 00385 }
G4double BLCoordinates::getRadiusCut | ( | const G4ThreeVector & | local | ) | [static] |
getRadiusCut() returns the radius cut for the lowest segmentCL containing local.
References segmentCLVector.
00388 { 00389 for(unsigned int seg=0; seg<segmentCLVector.size(); ++seg) { 00390 // note no check against minZ -- segmentCLs are ordered 00391 if(local[2] <= segmentCLVector[seg].maxZ && 00392 (segmentCLVector[seg].radiusCut <= 0.0 || 00393 sqrt(local[0]*local[0]+local[1]*local[1]) <= 00394 segmentCLVector[seg].radiusCut)) { 00395 return segmentCLVector[seg].radiusCut; 00396 } 00397 } 00398 00399 return 0; 00400 }
BLCoordinateType BLCoordinates::getCoordinateType | ( | G4String | s | ) | [static] |
getCoordinateType() converts a G4String name to a BLCoordinateType. Names are: global, centerline, reference; only the first letter is used, case insensitive. Also enables the generation of REference Coordinates, if they are selected by s.
References BLCOORD_CENTERLINE, BLCOORD_GLOBAL, BLCOORD_REFERENCE, BLCommand::printError(), and useReferenceCoordinates().
Referenced by BLCMDzntuple::command(), BLCMDvirtualdetector::command(), BLCMDtrace::command(), BLCMDtimentuple::command(), BLCMDprofile::command(), BLCMDprintf::command(), BLCMDplace::command(), BLCMDnewparticlentuple::command(), BLCMDmovie::command(), and BLCMDbeamlossntuple::command().
00403 { 00404 char c = tolower(s[0u]); 00405 switch(c) { 00406 case 'g': return BLCOORD_GLOBAL; 00407 case 'c': return BLCOORD_CENTERLINE; 00408 case 'r': useReferenceCoordinates(); 00409 return BLCOORD_REFERENCE; 00410 default: BLCommand::printError("Unknown coordinate type '%s'," 00411 " Centerline used.\n",s.c_str()); 00412 return BLCOORD_CENTERLINE; 00413 } 00414 }
void BLCoordinates::updateReferenceCoordinates | ( | ) | [private] |
updateReferenceCoordinates() will update the reference coordinates, given the centerline coordintes are already updated and valid. Called by setGlobal().
References global, reference, referenceTransformUsed, segmentCL, segmentCLVector, segmentRC, and segmentRCVector.
Referenced by setGlobal().
00427 { 00428 // if we are lost, use CL coord segment as a starting point 00429 if(segmentRC >= segmentRCVector.size()) { 00430 segmentRC = segmentCLVector[segmentCL].firstRCSegment; 00431 if(segmentRC > 0) --segmentRC; // try 1 previous 00432 } 00433 if(segmentRC >= segmentRCVector.size()) 00434 segmentRC = 0; 00435 00436 for( ; segmentRC<segmentRCVector.size(); ++segmentRC) { 00437 segmentRCVector[segmentRC].transform.getLocal(reference,global); 00438 //printf("segmentRC=%d minZ=%.4g maxZ=%.4g global=%.1f,%.1f,%.1f ref=%.1f,%.1f,%.3f\n",segmentRC,segmentRCVector[segmentRC].minZ,segmentRCVector[segmentRC].maxZ,global[0],global[1],global[2],reference[0],reference[1],reference[2]); 00439 referenceTransformUsed = &segmentRCVector[segmentRC].transform; 00440 if(segmentRC != 0) // no correction (minZ = -DBL_MAX) 00441 reference[2] += (segmentRCVector[segmentRC].correction - 1.0) 00442 * (reference[2]-segmentRCVector[segmentRC].minZ); 00443 //printf(" corrected reference[2]=%.1f\n",reference[2]); 00444 if(reference[2] > segmentRCVector[segmentRC].maxZ) 00445 continue; 00446 return; 00447 } 00448 G4Exception("BLCoordinates","Cannot determine reference coordinate segment", 00449 FatalException,""); 00450 }
void BLCoordinates::useReferenceCoordinates | ( | ) | [static] |
useReferenceCoordinates() will arrange to fill segmentRCVector[] with reference coordinate segments as the reference particle is tracked. Must be called before the reference particle is tracked. Unless this is called, reference coordinates will not be available. Multiple calls are the same as one.
References useRC.
Referenced by BLCMDmovie::command(), and getCoordinateType().
00417 { 00418 00419 if(!useRC) { 00420 useRC = true; 00421 new ReferenceCoordinates(); 00422 } 00423 }
static void BLCoordinates::setMostRecentReferenceZ | ( | G4double | z | ) | [inline, static] |
setMostRecentReferenceZ() sets the most recent z value of the reference coords. VALID ONLY WHILE TRACKING REFERENCE PARTICLE!
References mostRecentReferenceZ.
Referenced by BLCoordinates::ReferenceCoordinates::PreUserTrackingAction(), and BLCoordinates::ReferenceCoordinates::UserSteppingAction().
00385 { mostRecentReferenceZ = z; }
static G4double BLCoordinates::getMostRecentReferenceZ | ( | ) | [inline, static] |
getMostRecentReferenceZ() returns the most recent z value of the reference coords. VALID ONLY WHILE TRACKING REFERENCE PARTICLE!
References mostRecentReferenceZ.
Referenced by BLCMDmovie::referenceStep().
00390 { return mostRecentReferenceZ; }
int BLCoordinates::verify [private] |
Referenced by BLCoordinates(), and isValid().
unsigned int BLCoordinates::segmentCL [private] |
unsigned int BLCoordinates::segmentRC [private] |
BLCoordinateType BLCoordinates::prevType [private] |
Referenced by getCoords(), and getRotation().
G4double BLCoordinates::global[4] [private] |
G4double BLCoordinates::centerline[4] [private] |
Referenced by BLCoordinates(), findCLSegment(), getCenterlineCoords(), getCLZ(), getCoords(), and setGlobal().
G4double BLCoordinates::reference[4] [private] |
Referenced by BLCoordinates(), getCoords(), and updateReferenceCoordinates().
int BLCoordinates::segmentCLRadiusFailure [private] |
Referenced by BLCoordinates(), findCLSegment(), setGlobal(), and update().
G4Allocator< BLCoordinates > BLCoordinates::allocator [static, private] |
Referenced by operator delete(), and operator new().
std::vector< BLCoordinates::CLSegment > BLCoordinates::segmentCLVector [static, private] |
Referenced by BLCoordinates::ReferenceCoordinates::addStep(), corner(), findCLSegment(), getCenterlineRotation(), getCoords(), getCurrentGlobal(), getCurrentRadiusCut(), getCurrentRotation(), getCurrentTransform(), getGlobalAnywhere(), getGlobalRotation(), getRadiusCut(), getRotation(), getTransformRotation(), init(), setGlobal(), start(), and updateReferenceCoordinates().
unsigned BLCoordinates::currentCL = 0 [static, private] |
G4bool BLCoordinates::ring = false [static, private] |
Referenced by setGlobal(), and start().
BLCoordinates::State BLCoordinates::state = UNUSED [static, private] |
Referenced by BLCoordinates(), init(), and start().
G4bool BLCoordinates::steppingVerbose = false [static, private] |
G4bool BLCoordinates::useRC = false [static, private] |
Referenced by useReferenceCoordinates().
G4bool BLCoordinates::validRC = false [static, private] |
std::vector< BLCoordinates::RCSegment > BLCoordinates::segmentRCVector [static, private] |
int BLCoordinates::currentRC [static, private] |
G4double BLCoordinates::tolerance = 1.0e-4 [static, private] |
G4double BLCoordinates::mostRecentReferenceZ = 0 [static, private] |
Referenced by getMostRecentReferenceZ(), and setMostRecentReferenceZ().
int BLCoordinates::referenceCLsegment [private] |
Referenced by BLCoordinates().
Referenced by BLCoordinates(), getCoords(), getRotation(), and updateReferenceCoordinates().