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.

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

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