BLCoordinates Class Reference

#include <BLCoordinates.hh>

Inheritance diagram for BLCoordinates:

BLTrackInfo

List of all members.


Detailed Description

class BLCoordinates manages the centerline coordinates.

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 BLCoordinateTransformgetCurrentTransform ()
 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
BLCoordinateTransformreferenceTransformUsed

Static Private Attributes

static G4Allocator< BLCoordinatesallocator
static std::vector< CLSegmentsegmentCLVector
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< RCSegmentsegmentRCVector
static int currentRC
static G4double tolerance = 1.0e-4
static G4double mostRecentReferenceZ = 0

Classes

struct  CLSegment
struct  RCSegment
class  ReferenceCoordinates


Member Enumeration Documentation

anonymous enum [private]

Enumerator:
VERIFY 
00095 { VERIFY=0x13579246 };

enum BLCoordinates::State [private]

Enumerator:
UNUSED 
CONSTRUCTING 
TRACKING 


Constructor & Destructor Documentation

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 }

virtual BLCoordinates::~BLCoordinates (  )  [inline, virtual]

Destructor.

00115 { }


Member Function Documentation

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

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                 G4bool ret = findCLSegment(_global);
00119                 if(ret && validRC) updateReferenceCoordinates();
00120                 return ret;
00121         }
00122 
00123         segmentCLVector[segmentCL].transform.getLocal(centerline,_global);
00124 //if(steppingVerbose) 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);
00125         while(centerline[2] > segmentCLVector[segmentCL].maxZ ||
00126            (segmentCLVector[segmentCL].radiusCut > 0.0 &&
00127             (int)segmentCL != segmentCLRadiusFailure &&
00128             sqrt(centerline[0]*centerline[0]+centerline[1]*centerline[1]) >
00129                                         segmentCLVector[segmentCL].radiusCut)) {
00130                 if(++segmentCL >= segmentCLVector.size()) {
00131                         if(ring && findCLSegment(_global)) {
00132                                 segmentCLVector[segmentCL].transform.getLocal(centerline,
00133                                                                         _global);
00134                                 if(validRC) updateReferenceCoordinates();
00135                                 return true;
00136                         }
00137                         --segmentCL;
00138                         return false;
00139                 }
00140                 segmentCLVector[segmentCL].transform.getLocal(centerline,_global);
00141 //if(steppingVerbose) 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);
00142 /***
00143                 // note no check against minZ -- segmentCLs are ordered
00144                 if(centerline[2] > segmentCLVector[segmentCL].maxZ ||
00145                    (segmentCLVector[segmentCL].radiusCut > 0.0 &&
00146                     segmentCL != segmentCLRadiusFailure &&
00147                     sqrt(centerline[0]*centerline[0]+centerline[1]*centerline[1]) >
00148                                         segmentCLVector[segmentCL].radiusCut)) {
00149                         if(ring && findCLSegment(_global)) {
00150                                 segmentCLVector[segmentCL].transform.getLocal(centerline,
00151                                                                         _global);
00152                                 if(validRC) updateReferenceCoordinates();
00153                                 return true;
00154                         }
00155                         return false;
00156                 }
00157 ***/
00158         }
00159 
00160         if(validRC) updateReferenceCoordinates();
00161         return true;
00162 }

G4bool BLCoordinates::setGlobal ( const G4ThreeVector &  pos,
G4double  time 
)

References global, and setGlobal().

00165 {
00166         global[0]=pos[0],global[1]=pos[1],global[2]=pos[2];
00167         global[3]=time;
00168         return setGlobal(global);
00169 }

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(), 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().

00190 {
00191         prevType = type;
00192 
00193         switch(type) {
00194         case BLCOORD_GLOBAL:
00195                 local[0] = global[0]; local[1] = global[1];
00196                 local[2] = global[2]; local[3] = global[3];
00197                 return;
00198         case BLCOORD_CENTERLINE:
00199                 assert(segmentCL < segmentCLVector.size());
00200                 local[0] = centerline[0]; local[1] = centerline[1];
00201                 local[2] = centerline[2]; local[3] = centerline[3];
00202                 return;
00203         case BLCOORD_REFERENCE:
00204                 assert(validRC);
00205                 // interpolate
00206                 if(reference[2] < segmentRCVector[segmentRC].interpolatePrev &&
00207                                                            segmentRC > 0) {
00208                         referenceTransformUsed = 
00209                                 &segmentRCVector[segmentRC-1].interpolate;
00210                         referenceTransformUsed->getLocal(reference,global);
00211                 } else if(reference[2] > 
00212                                 segmentRCVector[segmentRC].interpolateThis) {
00213                         referenceTransformUsed = 
00214                                 &segmentRCVector[segmentRC].interpolate;
00215                         referenceTransformUsed->getLocal(reference,global);
00216                 }
00217 //G4ThreeVector axis = referenceTransformUsed->getRotation().getAxis();
00218 //G4double angle = referenceTransformUsed->getRotation().getDelta();
00219 //printf("segmentRC=%d minZ=%.4g maxZ=%.4g global=%.1f,%.1f,%.1f ref=%.1f,%.1f,%.3f 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);
00220                 local[0] = reference[0]; local[1] = reference[1];
00221                 local[2] = reference[2]; local[3] = reference[3];
00222                 return;
00223         }
00224         G4Exception("BLCoordinates","Invalid Coordinate Type",FatalException,
00225                                                                         "");
00226 }

void BLCoordinates::getCoords ( BLCoordinateType  type,
G4ThreeVector &  local 
)

References getCoords().

00229 {
00230         G4double tmp[4];
00231         getCoords(type,tmp);
00232         local[0] = tmp[0];
00233         local[1] = tmp[1];
00234         local[2] = tmp[2];
00235 }

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]

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.

00238 {
00239         switch(type) {
00240         case BLCOORD_GLOBAL:
00241                 return IdentityMatrix;
00242         case BLCOORD_CENTERLINE:
00243                 assert(segmentCL < segmentCLVector.size());
00244                 return segmentCLVector[segmentCL].transform.getRotation();
00245         case BLCOORD_REFERENCE:
00246                 assert(validRC && referenceTransformUsed != 0);
00247                 return referenceTransformUsed->getRotation();
00248         }
00249         G4Exception("BLCoordinates","Invalid Coordinate Type",FatalException,
00250                                                                         "");
00251         return IdentityMatrix;
00252 }

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().

00172 {
00173         // ensure there is a BLCoordiantes linked into the track
00174         BLCoordinates *c = (BLCoordinates *)track->GetUserInformation();
00175         if(!c || !c->isValid()) {
00176                 c = new BLCoordinates();
00177                 track->SetUserInformation(c);
00178         }
00179 
00180         // now update the object and apply radiusCut
00181         if(!c->setGlobal(track->GetPosition(),track->GetGlobalTime())) {
00182                 track->SetTrackStatus(fStopAndKill);
00183                 if(steppingVerbose)
00184                         printf("Track Centerline radius exceeds radiusCut\n");
00185                 c->segmentCLRadiusFailure = c->segmentCL;
00186         }
00187 }

virtual void BLCoordinates::Print (  )  const [inline, virtual]

Print from G4VUserTrackInformation.

00216 { }

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().

00255 {
00256         if(state == TRACKING) {
00257                 G4Exception("BLCoordinates","Invalid init during tracking",
00258                                                         FatalException, "");
00259         }
00260         if(state == UNUSED) {
00261                 state = CONSTRUCTING;
00262                 segmentCLVector.push_back(CLSegment(-DBL_MAX, +DBL_MAX, 0.0,
00263                                                 BLCoordinateTransform()));
00264                 currentCL = 0;
00265                 steppingVerbose = 
00266                         BLManager::getObject()->getSteppingVerbose() != 0;
00267         }
00268 }

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().

00273 {
00274         if(state != UNUSED || currentCL != 0) {
00275                 G4Exception("BLCoordinates","Invalid start",FatalException, "");
00276         }
00277         init();
00278         
00279         G4ThreeVector pos(0,0,-z);
00280         if(rotation)
00281                 pos = *rotation * pos;
00282         pos[0] += _global[0];
00283         pos[1] += _global[1];
00284         pos[2] += _global[2];
00285         segmentCLVector[0].transform = (rotation ?
00286                                         BLCoordinateTransform(rotation,pos) : 
00287                                         BLCoordinateTransform(pos));
00288         segmentCLVector[0].radiusCut = radiusCut;
00289         ring = _ring;
00290         steppingVerbose = BLManager::getObject()->getSteppingVerbose() != 0;
00291 }

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().

00295 {
00296         init();
00297         
00298         G4double local[4], glob[4];
00299         local[0]=local[1]=local[3]=0.0;
00300         local[2] = z;
00301         segmentCLVector[currentCL].transform.getGlobal(local,glob);
00302 
00303         segmentCLVector[currentCL].maxZ = z;
00304         G4RotationMatrix *p=0;
00305         if(rotation)
00306                 // if C=currentCL rot, and R=this rot, we want R referenced to
00307                 // the CURRENT centerline coordinates, so we need C R C^-1.
00308                 // But we compose this with C (on the right), and get C R.
00309                 p = new G4RotationMatrix(*getCurrentRotation() * *rotation);
00310         else if(getCurrentTransform().isRotated())
00311                 p = getCurrentRotation();
00312 
00313         G4ThreeVector pos(0,0,-z);
00314         if(p)
00315                 pos = *p * pos;
00316         pos[0] += glob[0];
00317         pos[1] += glob[1];
00318         pos[2] += glob[2];
00319         segmentCLVector.push_back(CLSegment(z, +DBL_MAX, radiusCut,
00320                                         BLCoordinateTransform(p,pos)));
00321         currentCL = segmentCLVector.size() - 1;
00322 }

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().

00325 {
00326         init();
00327         return segmentCLVector[currentCL].transform;
00328 }

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().

00331 {
00332         init();
00333         return &segmentCLVector[currentCL].transform.getRotationInverse();
00334 }

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().

00337 {
00338         init();
00339         segmentCLVector[currentCL].transform.getGlobal(local,_global);
00340 }

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.

00343 {
00344         init();
00345         segmentCLVector[currentCL].transform.getGlobal(local,_global);
00346 }

void BLCoordinates::getGlobalAnywhere ( const G4ThreeVector &  local,
G4ThreeVector &  _global 
) [static]

getGlobalAnywhere() returns global coordinates from centerline coordinates, using the lowest segmentCL containing local.

References segmentCLVector.

00350 {
00351         for(unsigned int seg=0; seg<segmentCLVector.size(); ++seg) {
00352                 // note no check against minZ -- segmentCLs are ordered
00353                 if(local[2] <= segmentCLVector[seg].maxZ &&
00354                    (segmentCLVector[seg].radiusCut <= 0.0 ||
00355                     sqrt(local[0]*local[0]+local[1]*local[1]) <=
00356                                         segmentCLVector[seg].radiusCut)) {
00357                         segmentCLVector[seg].transform.getGlobal(local,_global);
00358                         return;
00359                 }
00360         }
00361 }

G4RotationMatrix * BLCoordinates::getGlobalRotation ( const G4ThreeVector &  local  )  [static]

getGlobalRotation() returns the local->global rotation matrix coordinates, using the lowest segmentCL containing local.

References segmentCLVector.

00364 {
00365         for(unsigned int seg=0; seg<segmentCLVector.size(); ++seg) {
00366                 // note no check against minZ -- segmentCLs are ordered
00367                 if(local[2] <= segmentCLVector[seg].maxZ &&
00368                    (segmentCLVector[seg].radiusCut <= 0.0 ||
00369                     sqrt(local[0]*local[0]+local[1]*local[1]) <=
00370                                         segmentCLVector[seg].radiusCut)) {
00371                         return &segmentCLVector[seg].transform.getRotation();
00372                 }
00373         }
00374 
00375         return 0;
00376 }

G4double BLCoordinates::getRadiusCut ( const G4ThreeVector &  local  )  [static]

getRadiusCut() returns the radius cut for the lowest segmentCL containing local.

References segmentCLVector.

00379 {
00380         for(unsigned int seg=0; seg<segmentCLVector.size(); ++seg) {
00381                 // note no check against minZ -- segmentCLs are ordered
00382                 if(local[2] <= segmentCLVector[seg].maxZ &&
00383                    (segmentCLVector[seg].radiusCut <= 0.0 ||
00384                     sqrt(local[0]*local[0]+local[1]*local[1]) <=
00385                                         segmentCLVector[seg].radiusCut)) {
00386                         return segmentCLVector[seg].radiusCut;
00387                 }
00388         }
00389 
00390         return 0;
00391 }

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().

00394 {
00395         char c = tolower(s[0u]);
00396         switch(c) {
00397         case 'g':       return BLCOORD_GLOBAL;
00398         case 'c':       return BLCOORD_CENTERLINE;
00399         case 'r':       useReferenceCoordinates();
00400                         return BLCOORD_REFERENCE;
00401         default:        BLCommand::printError("Unknown coordinate type '%s',"
00402                                 " Centerline used.\n",s.c_str());
00403                         return BLCOORD_CENTERLINE;
00404         }
00405 }

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().

00418 {
00419         // if we are lost, use CL coord segment as a starting point
00420         if(segmentRC >= segmentRCVector.size()) {
00421                 segmentRC = segmentCLVector[segmentCL].firstRCSegment;
00422                 if(segmentRC > 0) --segmentRC; // try 1 previous
00423         }
00424         if(segmentRC >= segmentRCVector.size())
00425                 segmentRC = 0;
00426 
00427         for( ; segmentRC<segmentRCVector.size(); ++segmentRC) {
00428                 segmentRCVector[segmentRC].transform.getLocal(reference,global);
00429 //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]);
00430                 referenceTransformUsed = &segmentRCVector[segmentRC].transform;
00431                 reference[2] += (segmentRCVector[segmentRC].correction - 1.0) *
00432                                 (reference[2]-segmentRCVector[segmentRC].minZ);
00433                 if(reference[2] > segmentRCVector[segmentRC].maxZ)
00434                         continue;
00435                 return;
00436         }
00437         G4Exception("BLCoordinates","Cannot determine reference coordinate segment",
00438                                                 FatalException,"");
00439 }

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().

00408 {
00409 
00410         if(!useRC) {
00411                 useRC = true;
00412                 new ReferenceCoordinates();
00413         }
00414 }

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; } 


Member Data Documentation

int BLCoordinates::verify [private]

Referenced by BLCoordinates(), and isValid().

unsigned int BLCoordinates::segmentCL [private]

unsigned int BLCoordinates::segmentRC [private]

Referenced by getCoords(), and getRotation().

G4double BLCoordinates::global[4] [private]

G4double BLCoordinates::centerline[4] [private]

G4double BLCoordinates::reference[4] [private]

G4Allocator< BLCoordinates > BLCoordinates::allocator [static, private]

Referenced by operator delete(), and operator new().

std::vector< BLCoordinates::CLSegment > BLCoordinates::segmentCLVector [static, private]

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]

Referenced by init(), start(), and update().

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 BLCoordinates().


The documentation for this class was generated from the following files:
g4beamline