BLCoordinates::ReferenceCoordinates Class Reference

Inheritance diagram for BLCoordinates::ReferenceCoordinates:

BLManager::TrackingAction BLManager::SteppingAction

List of all members.

Public Member Functions

 ReferenceCoordinates ()
virtual ~ReferenceCoordinates ()
void PreUserTrackingAction (const G4Track *track)
void PostUserTrackingAction (const G4Track *track)
void UserSteppingAction (const G4Step *step)
void addStep (const G4ThreeVector &pos, G4double z, bool final=false)

Private Attributes

BLManagermanager
G4ThreeVector prevDir
G4double prevZ
G4ThreeVector prevPos
G4double stepZ
G4ThreeVector stepPos
G4ThreeVector prevStepPos

Static Private Attributes

static G4ThreeVector unitZ


Constructor & Destructor Documentation

BLCoordinates::ReferenceCoordinates::ReferenceCoordinates (  ) 

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

00039 { }


Member Function Documentation

void BLCoordinates::ReferenceCoordinates::PreUserTrackingAction ( const G4Track *  track  )  [virtual]

Implements BLManager::TrackingAction.

References BEAM, BLCOORD_CENTERLINE, BLCoordinates::getCoords(), BLManager::getObject(), BLCoordinates::isValid(), prevDir, prevPos, prevStepPos, prevZ, REFERENCE, BLCoordinates::segmentRC, BLCoordinates::setMostRecentReferenceZ(), stepPos, stepZ, and BLCoordinates::validRC.

00464 {
00465         BLCoordinates *coord = (BLCoordinates *)track->GetUserInformation();
00466         if(!coord || !coord->isValid()) {
00467                 G4Exception("BLCoordinates","Coordinate Object Missing in Reference Track",
00468                                                         FatalException, "");
00469         }
00470 
00471         coord->segmentRC = 99999999;
00472 
00473         if(BLManager::getObject()->getState() == BEAM) {
00474                 if(!validRC) {
00475                         G4Exception("BLCoordinates","Reference Coordinates not available",
00476                                                         FatalException, "");
00477                 }
00478                 return;
00479         }
00480 
00481         if(BLManager::getObject()->getState() != REFERENCE) return;
00482 
00483         // @@@ need to handle multiple reference particles !!!!
00484 
00485         prevPos = track->GetPosition();
00486         prevDir = track->GetMomentumDirection();
00487 
00488         // initial Z is same as Centerline Coord value.
00489         G4ThreeVector cl;
00490         coord->getCoords(BLCOORD_CENTERLINE,cl);
00491         prevZ = cl[2];
00492         coord->setMostRecentReferenceZ(prevZ);
00493 
00494         // fake initial step of 1 mm 
00495         // (ensures sanity if reference bends in first step)
00496         stepPos = prevPos - prevDir;
00497         prevStepPos = stepPos - prevDir;
00498         stepZ = prevZ - 1.0;
00499 }

void BLCoordinates::ReferenceCoordinates::PostUserTrackingAction ( const G4Track *  track  )  [virtual]

Implements BLManager::TrackingAction.

References addStep(), BLManager::getObject(), BLCoordinates::global, REFERENCE, BLCoordinates::segmentRCVector, and BLCoordinates::validRC.

00502 {
00503         if(BLManager::getObject()->getState() != REFERENCE) return;
00504 
00505         // add 1 mm to final step
00506         // (ensures sanity if reference bends in final step)
00507         G4ThreeVector pos = track->GetPosition()+track->GetMomentumDirection();
00508         G4double z = track->GetTrackLength() + 1.0;
00509         addStep(pos,z,true);
00510 
00511         validRC = true;
00512 
00513         // for debugging -- dump segmentRCVector
00514         printf("\n# segmentRCVector\n#N minZ maxZ globX globY globZ axisX axisY axisZ angleDeg corr\n");
00515         G4double local[4];
00516         local[0] = local[1] = local[2] = local[3] = 0.0;
00517         for(unsigned i=0; i<segmentRCVector.size(); ++i) {
00518                 G4double global[4];
00519                 local[2] = (i==0 ? 0.0 : segmentRCVector[i].minZ);
00520                 segmentRCVector[i].transform.getGlobal(local,global);
00521                 G4ThreeVector axis = 
00522                         segmentRCVector[i].transform.getRotation().getAxis();
00523                 printf("%d %.4g %.4g %.4g %.4g %.4g %.4f %.4f %.4f %.1f %.3f\n",i,
00524                         segmentRCVector[i].minZ,segmentRCVector[i].maxZ,
00525                         global[0],global[1],global[2],axis[0],axis[1],axis[2],
00526                         segmentRCVector[i].transform.getRotation().getDelta()/deg,
00527                         segmentRCVector[i].correction);
00528         }
00529         printf("\n");
00530         fflush(stdout);
00531 }

void BLCoordinates::ReferenceCoordinates::UserSteppingAction ( const G4Step *  step  )  [virtual]

Implements BLManager::SteppingAction.

References addStep(), BLManager::getObject(), prevDir, prevPos, prevZ, REFERENCE, BLCoordinates::setMostRecentReferenceZ(), and BLCoordinates::tolerance.

00534 {
00535         // handle all states
00536         if(BLManager::getObject()->getState() != REFERENCE) return;
00537 
00538         G4Track *track = step->GetTrack();
00539         G4ThreeVector dir = track->GetMomentumDirection();
00540         if((dir-prevDir).mag2() > tolerance*tolerance)
00541                 addStep(prevPos, prevZ);
00542 
00543         prevDir = dir;
00544         prevPos = track->GetPosition();
00545         prevZ = track->GetTrackLength();
00546         BLCoordinates::setMostRecentReferenceZ(prevZ);
00547 }

void BLCoordinates::ReferenceCoordinates::addStep ( const G4ThreeVector &  pos,
G4double  z,
bool  final = false 
)

References BLCoordinates::currentCL, M_PI, prevStepPos, BLCoordinates::segmentCLVector, BLCoordinates::segmentRCVector, stepPos, stepZ, BLCoordinates::tolerance, and unitZ.

Referenced by PostUserTrackingAction(), and UserSteppingAction().

00551 {
00552         G4ThreeVector delta = pos - stepPos;
00553         G4double magDelta = delta.mag();
00554 
00555         // protect against delta being too small (should never happen,
00556         // as momentum direction changed significantly or first/final step
00557         // is >= 1 mm).
00558         if(magDelta < tolerance) return;
00559 
00560         // logic for rings is in CL coordinates; when we are lost, use CL
00561         // segment as starting point
00562         if(segmentCLVector[currentCL].firstRCSegment > segmentRCVector.size())
00563                 segmentCLVector[currentCL].firstRCSegment = 
00564                                                         segmentRCVector.size();
00565 
00566         G4ThreeVector rot = unitZ.cross(delta.unit());
00567         G4double a = asin(rot.mag());
00568         rot = rot.unit();
00569         // handle all 4 quadrants
00570         G4double b = unitZ.dot(delta.unit());
00571         if(a >= 0.0) {
00572                 if(b >= 0.0)
00573                         a = a;
00574                 else
00575                         a = M_PI - a;
00576         } else {
00577                 if(b >= 0.0)
00578                         a = a;
00579                 else
00580                         a = -M_PI + a;
00581         }
00582 
00583         G4RotationMatrix *p = 0;
00584         if(a > tolerance)
00585                 p = new G4RotationMatrix(rot,a);
00586         G4ThreeVector tmp(0,0,-stepZ);
00587         if(p)
00588                 tmp = *p * tmp;
00589         tmp += stepPos;
00590 
00591         G4double minZ = (segmentRCVector.size()==0 ? -DBL_MAX : stepZ);
00592         G4double maxZ = (final ? DBL_MAX : z);
00593         G4double interpPrev = stepZ + (z-stepZ)/4.0;
00594         G4double interpThis = z - (z-stepZ)/4.0;
00595         G4double correction = (z-stepZ)/magDelta;
00596 //printf("addStep minZ=%.4g maxZ=%.4g start x,y,z=%.3f,%.3f,%.3f  axis=%.3f,%.3f,%.3f aDeg=%.1f  corr=%.6f\n",minZ,maxZ,stepPos.x(),stepPos.y(),stepPos.z(),rot[0],rot[1],rot[2],a/deg,correction);
00597         segmentRCVector.push_back(RCSegment(minZ,maxZ,
00598                 BLCoordinateTransform(p,tmp),correction,interpPrev,interpThis));
00599 
00600         // set previous interpolate
00601         if(segmentRCVector.size() >= 2) {
00602                 delta = (pos-prevStepPos);
00603                 rot = unitZ.cross(delta.unit());
00604                 a = asin(rot.mag());
00605                 rot = rot.unit();
00606                 // handle all 4 quadrants
00607                 b = unitZ.dot(delta.unit());
00608                 if(a >= 0.0) {
00609                         if(b >= 0.0)
00610                                 a = a;
00611                         else
00612                                 a = M_PI - a;
00613                 } else {
00614                         if(b >= 0.0)
00615                                 a = a;
00616                         else
00617                                 a = -M_PI + a;
00618                 }
00619                 p = 0;
00620                 if(a > tolerance)
00621                         p = new G4RotationMatrix(rot,a);
00622                 tmp = G4ThreeVector(0,0,-stepZ);
00623                 if(p)
00624                         tmp = *p * tmp;
00625                 tmp += stepPos;
00626 //printf("addStep prev interpolate x,y,z=%.3f,%.3f,%.3f  axis=%.3f,%.3f,%.3f aDeg=%.1f\n",stepPos.x(),stepPos.y(),stepPos.z(),rot[0],rot[1],rot[2],a/deg);
00627                 segmentRCVector[segmentRCVector.size()-2].interpolate
00628                                                 = BLCoordinateTransform(p,tmp);
00629         }
00630 
00631         // update after this segment
00632         prevStepPos = stepPos;
00633         stepPos = pos;
00634         stepZ = z;
00635 }


Member Data Documentation

Referenced by ReferenceCoordinates().

Referenced by addStep(), and PreUserTrackingAction().

Referenced by addStep(), and PreUserTrackingAction().

G4ThreeVector BLCoordinates::ReferenceCoordinates::unitZ [static, private]

Referenced by addStep().


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