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.

00380 { 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.

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

void TrackerPlaneInstance::setTrueMode (  )  [inline]

References BLTRACKER_TRUE, mode, and setNotHit().

00385 { 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.

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

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

bool TrackerPlaneInstance::isHit (  )  [inline]

References haveBeenHit.

Referenced by UserSteppingAction().

00390 { 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