TrackerPlaneInstance Class Reference

Inheritance diagram for TrackerPlaneInstance:

BLManager::SteppingAction

List of all members.


Detailed Description

class TrackerPlaneInstance implements one placement of a tracker plane.

Public Member Functions

 TrackerPlaneInstance (G4String _name, G4VPhysicalVolume *pv, BLCMDtrackerplane *plane)
 constructor.
G4String getName ()
void UserSteppingAction (const G4Step *step)
 UserSteppingAction() from BLManager::SteppingAction.
void setTrueMode ()
void setIgnoreMode ()
void setFitMode (int w, double t)
void setNotHit ()
bool isHit ()

Private Attributes

G4String name
G4VPhysicalVolume * thisVol
BLTrackerMode mode
BLCMDtrackertracker
int planeID
int hitWire
double hitTime
bool haveBeenHit
double sinTheta
double cosTheta
double wireSpacing
double wireOffset
double sigmaWire
double sigmaTime
double sinThetaErr
double cosThetaErr
double wireSpacingErr
double wireOffsetErr

Constructor & Destructor Documentation

TrackerPlaneInstance::TrackerPlaneInstance ( G4String  _name,
G4VPhysicalVolume *  pv,
BLCMDtrackerplane plane 
)


Member Function Documentation

G4String TrackerPlaneInstance::getName (  )  [inline]

References name.

00384 { return name; }

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

UserSteppingAction() from BLManager::SteppingAction.

Implements BLManager::SteppingAction.

References BEAM, BLCOORD_CENTERLINE, BLTRACKER_FIT, BLTRACKER_IGNORE, BLTRACKER_TRUE, cosTheta, cosThetaErr, BLCMDtracker::delta(), BLCoordinates::getCoords(), BLManager::getObject(), BLCoordinates::getRotation(), BLManager::getState(), haveBeenHit, hitTime, hitWire, isHit(), BLCoordinates::isValid(), mode, name, NO_HIT, BLCMDtracker::planeHit(), planeID, sigmaTime, sigmaWire, sinTheta, sinThetaErr, thisVol, tracker, BLCMDtracker::verbose, wireOffset, wireOffsetErr, wireSpacing, and wireSpacingErr.

01527 {
01528         BLManagerState state = BLManager::getObject()->getState();
01529         if(state != BEAM || mode == BLTRACKER_IGNORE) return;
01530 
01531         if(tracker->verbose >= 3)
01532                 printf("trackerplane %s UserSteppingAction\n", name.c_str());
01533         if(isHit()) {
01534                 if(tracker->verbose >= 3)
01535                         printf("trackerplane %s already hit\n", name.c_str());
01536                 return;
01537         }
01538 
01539         // get basic physical-volume info
01540         G4StepPoint *prePoint = step->GetPreStepPoint();
01541         if(!prePoint) return;
01542         G4VPhysicalVolume *preVol = prePoint->GetPhysicalVolume();
01543         if(!preVol) return;
01544         G4StepPoint *postPoint = step->GetPostStepPoint();
01545         if(!postPoint) return;
01546         G4VPhysicalVolume *postVol = postPoint->GetPhysicalVolume();
01547         if(!postVol) return;
01548         
01549         // return if not entering thisVol
01550         if(preVol == postVol || postVol != thisVol) return;
01551 
01552         G4Track *track = step->GetTrack();
01553 
01554         G4ThreeVector position = track->GetPosition();
01555         G4double time = track->GetGlobalTime();
01556         G4ThreeVector momentum = track->GetMomentum();
01557 
01558         // transform to centerline coordinates, if available
01559         BLCoordinates *c = (BLCoordinates *)track->GetUserInformation();
01560         if(c && c->isValid()) {
01561                 c->getCoords(BLCOORD_CENTERLINE,position);
01562                 momentum = c->getRotation() * momentum;
01563         }
01564 
01565         // Remember that isHit() is known to be false.
01566         haveBeenHit = true;
01567 
01568         if(mode == BLTRACKER_TRUE) {
01569                 hitWire = 0; // fake wire hit if no wires
01570                 if(wireSpacing > 0.0) {
01571                         double x = cosTheta*position[0] + sinTheta*position[1];
01572                         hitWire = (int)floor((x-wireOffset)/wireSpacing+0.5);
01573                 }
01574                 hitTime = 0.0; // fake time if no timing
01575                 if(sigmaTime > 0.0)
01576                         hitTime = time + CLHEP::RandGauss::shoot(0.0,sigmaTime);
01577                 tracker->planeHit(planeID,hitWire,hitTime);
01578                 if(tracker->verbose >= 3)
01579                         printf("trackerplane %s hitWire=%d hitTIme=%.3f\n",
01580                                                 name.c_str(),hitWire,hitTime);
01581         } else if(mode == BLTRACKER_FIT) {
01582                 if(hitWire == NO_HIT) {
01583                         if(tracker->verbose >= 3)
01584                                 printf("trackerplane %s no true hit\n",
01585                                                                 name.c_str());
01586                         return;
01587                 }
01588                 double dw=0.0, dt=0.0;
01589                 bool validWire=false, validTime=false;
01590                 if(wireSpacing > 0.0) {
01591                         double x = cosThetaErr*position[0] +
01592                                                         sinThetaErr*position[1];
01593                         dw = (x-wireSpacingErr*hitWire-wireOffsetErr)/sigmaWire;
01594                         validWire = true;
01595                 }
01596                 if(sigmaTime > 0.0) {
01597                         dt = (time-hitTime)/sigmaTime;
01598                         validTime = true;
01599                 }
01600                 if(tracker->verbose >= 3)
01601                         printf("trackerplane %s hit dw=%.3f dt=%.3f\n",
01602                                                         name.c_str(),dw,dt);
01603                 if(!tracker->delta(planeID,dw,dt,validWire,validTime)) {
01604                         track->SetTrackStatus(fStopAndKill);
01605                         if(tracker->verbose >= 3) 
01606                                 printf("trackerplane %s kills fitting track\n",
01607                                                                 name.c_str());
01608                 }
01609         }
01610 }

void TrackerPlaneInstance::setTrueMode (  )  [inline]

References BLTRACKER_TRUE, mode, and setNotHit().

00389 { mode = BLTRACKER_TRUE; setNotHit(); }

void TrackerPlaneInstance::setIgnoreMode (  )  [inline]

References BLTRACKER_IGNORE, mode, and setNotHit().

void TrackerPlaneInstance::setFitMode ( int  w,
double  t 
) [inline]

References BLTRACKER_FIT, haveBeenHit, hitTime, hitWire, and mode.

00392           { mode = BLTRACKER_FIT; hitWire = w; hitTime = t; haveBeenHit=false; }

void TrackerPlaneInstance::setNotHit (  )  [inline]

References haveBeenHit, hitTime, hitWire, and NO_HIT.

Referenced by setIgnoreMode(), setTrueMode(), and TrackerPlaneInstance().

00393 { hitWire=NO_HIT; hitTime=0.0; haveBeenHit=false; }

bool TrackerPlaneInstance::isHit (  )  [inline]

References haveBeenHit.

Referenced by UserSteppingAction().

00394 { return haveBeenHit; }


Member Data Documentation

G4String TrackerPlaneInstance::name [private]

G4VPhysicalVolume* TrackerPlaneInstance::thisVol [private]


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