BLCMDtrackerplane Class Reference

Inheritance diagram for BLCMDtrackerplane:

BLElement BLCommand

List of all members.


Detailed Description

class BLCMDtrackerplane - implements one plane of a tracker.

A tracker plane is a detector consisting of either: A) a large number of parallel wires arranged in a plane B) a counter that gives timing information C) both (A) and (B) For a given track hitting the plane, one and only one wire is hit. Arguments to the command specify the wire spacing, orientation, and offset; the number of wires is implicit from their spacing and orientation and the size of the plane. if wireSpacing=0, no wires are present.

If sigmaTime>0 this plane measures the time of the track. A Gaussian random number with that sigma is added to the true value.

This class works with the BLCMDtracker class to permit the fitting of tracks to the wire hits. modes: true reports the hit wire # to the BLCMDtracker instance (this is G4beamline tracking "true" tracks) fit computes the plane's contribution to Chisq for the track and reports it to the BLCMDtracker instance ignore track is ignored

Public Member Functions

 BLCMDtrackerplane ()
 Default constructor.
virtual ~BLCMDtrackerplane ()
 Destructor.
BLElementclone ()
 clone()
 BLCMDtrackerplane (const BLCMDtrackerplane &r)
 Copy constructor.
G4String commandName ()
 commandName() returns "trackerplane".
int command (BLArgumentVector &argv, BLArgumentMap &namedArgs)
 command() implements the trackerplane command.
void defineNamedArgs ()
 defineNamedArgs() defines the named arguments for the command.
void argChanged ()
 argChanged() does internal computations after some arg changed
virtual void construct (G4RotationMatrix *relativeRotation, G4ThreeVector relativePosition, G4LogicalVolume *parent, G4String parentName, G4RotationMatrix *parentRotation, G4ThreeVector parentPosition)
 construct() will construct the trackerplane. Used for normal placements of a tracker plane object.
virtual G4double getLength ()
 getLength() returns this element's Length along the Z axis.
virtual G4double getWidth ()
 getWidth() returns this element's Width along the X axis.
virtual G4double getHeight ()
 getHeight() returns this element's height along the Y axis.
virtual G4bool isOK ()
 isOK() returns true.
void generatePoints (int npoints, std::vector< G4ThreeVector > &v)
 generatePoints() from BLElement
G4bool isOutside (G4ThreeVector &local, G4double tolerance)
 isOutside() from BLElement

Private Attributes

G4String tracker
G4double radius
G4double innerRadius
G4double height
G4double width
G4double length
G4String material
G4String color
G4VSolid * solid
G4double maxStep
G4double theta
G4double wireSpacing
G4double wireOffset
G4String errType
G4double errTheta
G4double errSpacing
G4double errOffset
G4double sigmaTime

Friends

class TrackerPlaneInstance


Constructor & Destructor Documentation

BLCMDtrackerplane::BLCMDtrackerplane (  ) 

Default constructor.

References BLCMDTYPE_DATA, color, errOffset, errSpacing, errTheta, errType, height, innerRadius, length, material, maxStep, radius, BLCommand::registerCommand(), BLCommand::setDescription(), BLCommand::setSynopsis(), sigmaTime, solid, theta, tracker, width, wireOffset, and wireSpacing.

Referenced by clone(), and command().

01251                                      : BLElement()
01252 {
01253         registerCommand(BLCMDTYPE_DATA);
01254         setSynopsis("Construct a tracker plane.");
01255         setDescription("A trackerplane belongs to a specific tracker, and "
01256                 "represents one measuring element of the tracker. "
01257                 "While the term 'wire' is used, a trackerplane can model "
01258                 "any planar measuring device that measures one dimension "
01259                 "using equally spaced detectors that are either on or off "
01260                 "(hit or not hit). A trackerplane can be "
01261                 "circular (specify radius and possibly innerRadius) or "
01262                 "rectangular (specify height and width).\n\n"
01263                 "The wires are at angle theta from the vertical, so theta=0 "
01264                 "means vertical wires that measure x; theta=90 means "
01265                 "horizontal wires that measure y; theta=180 also "
01266                 "measures x, but increasing x means decreasing wire #.\n\n"
01267                 "Setting wireSpacing=0 means there are no wires, which is "
01268                 "useful for a trackerplane that models a scintillator used "
01269                 "for track timing (set sigmaTime >= 0 to indicate that).\n\n"
01270                 "sigmaTime>0 means this plane can measure the time of the "
01271                 "track with that resolution. A Gaussian random number is added "
01272                 "to the true track's time at this plane when reading the "
01273                 "TrackHit NTuple.\n\n"
01274                 "Survey errors can be modeled using the err-arguments -- "
01275                 "the values they specify are applied during trackfitting "
01276                 "(but not for true tracks). "
01277                 "errType: 'fixed' means error values given are the actual values, "
01278                 "'gaussian' means error values are the sigma of a Gaussian random "
01279                 "number, 'rect' means error values are the half-width of a uniform "
01280                 "random number. The random number is picked before the run "
01281                 "begins; the random-number seed is set from the clock so "
01282                 "every run will have different random errors.\n\n"
01283                 "trackerplane has 3 modes:\n"
01284                 "  true    The trackerplane reports the hit wire and time to\n"
01285                 "          the tracker;\n"
01286                 "  fit     The trackerplane reports the Chisq contribution of\n"
01287                 "          the fit track distance to the true track's hit wire\n"
01288                 "          center, plus survey errors (if any). The track time\n"
01289                 "          also contributes to the chisq.\n"
01290                 "  ignore  Any track is ignored."
01291         );
01292         
01293         // initial field values
01294         tracker = "";
01295         radius = 0.0;
01296         innerRadius = 0.0;
01297         height = 0.0;
01298         width = 0.0;
01299         length = 1.0*mm;
01300         material = "";
01301         color = "1,1,1";
01302         solid = 0;
01303         maxStep = -99.0;
01304         theta = 0.0;
01305         wireSpacing = -1.0;
01306         wireOffset = 0.0;
01307         errType = "fixed";
01308         errTheta = 0.0;
01309         errSpacing = 0.0;
01310         errOffset = 0.0;
01311         sigmaTime = -1.0;
01312 }

BLCMDtrackerplane::~BLCMDtrackerplane (  )  [virtual]

Destructor.

References solid.

01315 {
01316         if(solid) delete solid;
01317 }

BLCMDtrackerplane::BLCMDtrackerplane ( const BLCMDtrackerplane r  ) 

Copy constructor.

References color, errOffset, errSpacing, errTheta, errType, height, innerRadius, length, material, maxStep, radius, sigmaTime, solid, theta, tracker, width, wireOffset, and wireSpacing.

01319                                                                : BLElement(r)
01320 {
01321         tracker = r.tracker;
01322         radius = r.radius;
01323         innerRadius = r.innerRadius;
01324         height = r.height;
01325         width = r.width;
01326         length = r.length;
01327         material = r.material;
01328         color = r.color;
01329         solid = r.solid;
01330         maxStep = r.maxStep;
01331         theta = r.theta;
01332         wireSpacing = r.wireSpacing;
01333         wireOffset = r.wireOffset;
01334         errType = r.errType;
01335         errTheta = r.errTheta;
01336         errSpacing = r.errSpacing;
01337         errOffset = r.errOffset;
01338         sigmaTime = r.sigmaTime;
01339 }


Member Function Documentation

BLElement* BLCMDtrackerplane::clone (  )  [inline, virtual]

clone()

Implements BLElement.

References BLCMDtrackerplane().

00313 { return new BLCMDtrackerplane(*this); }

G4String BLCMDtrackerplane::commandName (  )  [inline, virtual]

commandName() returns "trackerplane".

Implements BLCommand.

00319 { return "trackerplane"; }

int BLCMDtrackerplane::command ( BLArgumentVector argv,
BLArgumentMap namedArgs 
) [virtual]

command() implements the trackerplane command.

Implements BLCommand.

References BLCMDtrackerplane(), errOffset, errSpacing, errTheta, errType, BLCMDtracker::findTracker(), BLParam::getDouble(), BLCommand::getMaterial(), BLCommand::handleNamedArgs(), material, maxStep, Param, BLCommand::print(), BLCommand::printError(), BLElement::setName(), and tracker.

01342 {
01343         if(argv.size() != 1) {
01344                 printError("trackerplane: Invalid command, must have name");
01345                 return -1;
01346         }
01347 
01348         if(argv[0] == "default") {
01349                 return defaultTrackerPlane.handleNamedArgs(namedArgs);
01350         }
01351 
01352         BLCMDtrackerplane *t = new BLCMDtrackerplane(defaultTrackerPlane);
01353         t->setName(argv[0]);
01354         int retval = t->handleNamedArgs(namedArgs);
01355 
01356         if(t->maxStep < 0.0) t->maxStep = Param.getDouble("maxStep");
01357 
01358         if(BLCMDtracker::findTracker(t->tracker) == 0)
01359                 printError("trackerplane: invalid tracker");
01360 
01361         // handle errType
01362         if(t->errType.at(0) == 'g' || t->errType.at(0) == 'G') {
01363                 t->errTheta = t->errTheta*CLHEP::RandGauss::shoot();
01364                 t->errSpacing = t->errSpacing*CLHEP::RandGauss::shoot();
01365                 t->errOffset = t->errOffset*CLHEP::RandGauss::shoot();
01366         } else if(t->errType.at(0) == 'r' || t->errType.at(0) == 'R') {
01367                 t->errTheta = t->errTheta - t->errTheta*2.0*G4UniformRand();
01368                 t->errSpacing = t->errSpacing - t->errSpacing*2.0*G4UniformRand();
01369                 t->errOffset = t->errOffset - t->errOffset*2.0*G4UniformRand();
01370         }
01371 
01372         // check material exists
01373         if(t->material.size() > 0) getMaterial(t->material,false);
01374 
01375         t->print(argv[0]);
01376 
01377         return retval;
01378 }

void BLCMDtrackerplane::defineNamedArgs (  )  [virtual]

defineNamedArgs() defines the named arguments for the command.

Reimplemented from BLCommand.

References BLCommand::argDouble(), BLCommand::argString(), color, errOffset, errSpacing, errTheta, errType, height, innerRadius, length, material, maxStep, radius, sigmaTime, theta, tracker, width, wireOffset, and wireSpacing.

01381 {
01382         argString(tracker,"tracker","The tracker to which this plane belongs; REQUIRED.");
01383         argDouble(radius,"radius","The radius of the circular tracker plane (mm).");
01384         argDouble(innerRadius,"innerRadius","The inner radius of the circular tracker plane (0 mm).");
01385         argDouble(height,"height","The height of the rectangular tracker plane (mm).");
01386         argDouble(width,"width","The width of the rectangular tracker plane (mm).");
01387         argDouble(length,"length","The length of the tracker plane (mm).");
01388         argDouble(theta,"theta","Wire angle in X-Y plane (deg). 0=>x, 90=>y...",deg);
01389         argDouble(wireSpacing,"wireSpacing","Wire spacing (mm).");
01390         argDouble(wireOffset,"wireOffset","Wire # 0 offset (mm).");
01391         argString(errType,"errType","Error type: fixed, gaussian, rect.");
01392         argDouble(errTheta,"errTheta","Error in wire angle (deg).",deg);
01393         argDouble(errSpacing,"errSpacing","Error in wire spacing (mm).");
01394         argDouble(errOffset,"errOffset","Error in wire offset (mm).");
01395         argDouble(sigmaTime,"sigmaTime","Sigma for timing plane (ns), default=-1.");
01396         argDouble(maxStep,"maxStep","The maximum stepsize in the element (mm).");
01397         argString(material,"material","The material of the tracker plane.");
01398         argString(color,"color","The color of the tracker plane (''=invisible).");
01399 }

void BLCMDtrackerplane::argChanged (  )  [virtual]

argChanged() does internal computations after some arg changed

Reimplemented from BLCommand.

References height, radius, and width.

01402 {
01403         if(radius > 0.0) width = height = radius*2.0;
01404 }

void BLCMDtrackerplane::construct ( G4RotationMatrix *  relativeRotation,
G4ThreeVector  relativePosition,
G4LogicalVolume *  parent,
G4String  parentName,
G4RotationMatrix *  parentRotation,
G4ThreeVector  parentPosition 
) [virtual]

construct() will construct the trackerplane. Used for normal placements of a tracker plane object.

Implements BLElement.

References color, BLParam::getDouble(), BLCommand::getMaterial(), BLElement::getName(), BLManager::getObject(), BLCommand::getVisAttrib(), height, innerRadius, length, material, maxStep, Param, BLCommand::printError(), radius, BLManager::registerBeamStep(), BLManager::registerReferenceParticleStep(), solid, TrackerPlaneInstance, and width.

01412 {
01413         G4Material *mat;
01414         if(material != "")
01415                 mat = getMaterial(material);
01416         else
01417                 mat = parent->GetMaterial();
01418 
01419         G4String thisname = parentName+getName();
01420 
01421         if(!solid) {
01422                 if(radius > 0.0) {
01423                         solid = new G4Tubs(thisname+"Tubs", innerRadius, radius,
01424                                         length/2.0, 0.0, 2.0*pi);
01425                 } else if(height > 0.0 && width > 0.0) {
01426                         solid = new G4Box(thisname+"Box",width/2.0,
01427                                         height/2.0,length/2.0);
01428                 } else {
01429                         printError("trackerplane::construct %s INVALID - no "
01430                                 "radius or height&width",thisname.c_str());
01431                         return;
01432                 }
01433         }
01434         G4LogicalVolume *lv = new G4LogicalVolume(solid,mat, thisname+"LogVol");
01435         lv->SetVisAttributes(getVisAttrib(color));
01436         if(maxStep < 0.0) maxStep = Param.getDouble("maxStep");
01437         lv->SetUserLimits(new G4UserLimits(maxStep));
01438 
01439         // geant4 rotation convention is backwards from g4beamline
01440         G4RotationMatrix *g4rot = 0;
01441         if(relativeRotation)
01442                 g4rot = new G4RotationMatrix(relativeRotation->inverse());
01443 
01444         G4VPhysicalVolume *pv = new G4PVPlacement(g4rot, relativePosition,
01445                                         lv,thisname,parent,false,0);
01446 
01447         // get globalRotation and globalPosition
01448         G4RotationMatrix *globalRotation = 0;
01449         if(relativeRotation && parentRotation) {
01450                 globalRotation = 
01451                     new G4RotationMatrix(*parentRotation * *relativeRotation);
01452         } else if(relativeRotation) {
01453                 globalRotation = relativeRotation;
01454         } else if(parentRotation) {
01455                 globalRotation = parentRotation;
01456         }
01457         G4ThreeVector globalPosition(relativePosition + parentPosition);
01458         if(globalRotation)
01459                 globalPosition = *globalRotation * relativePosition +
01460                                 parentPosition;
01461 
01462         TrackerPlaneInstance *tpi = new TrackerPlaneInstance(thisname,pv,this);
01463 
01464         BLManager::getObject()->registerBeamStep(pv,tpi);
01465         // I'm not sure if this needs reference step any more...
01466         BLManager::getObject()->registerReferenceParticleStep(pv,tpi);
01467 
01468         printf("BLCMDtrackerplane::Construct %s parent=%s relZ=%.1f globZ=%.1f\n",
01469                 thisname.c_str(),parentName.c_str(),relativePosition[2],
01470                 globalPosition[2]);
01471 }

virtual G4double BLCMDtrackerplane::getLength (  )  [inline, virtual]

getLength() returns this element's Length along the Z axis.

Implements BLElement.

References length.

00340 { return length; }

virtual G4double BLCMDtrackerplane::getWidth (  )  [inline, virtual]

getWidth() returns this element's Width along the X axis.

Implements BLElement.

References width.

00343 { return width; }

virtual G4double BLCMDtrackerplane::getHeight (  )  [inline, virtual]

getHeight() returns this element's height along the Y axis.

Implements BLElement.

References height.

00346 { return height; }

virtual G4bool BLCMDtrackerplane::isOK (  )  [inline, virtual]

isOK() returns true.

Implements BLElement.

00349 { return true; }

void BLCMDtrackerplane::generatePoints ( int  npoints,
std::vector< G4ThreeVector > &  v 
) [virtual]

generatePoints() from BLElement

Implements BLElement.

References BLElement::generateBox(), BLElement::generateTubs(), height, innerRadius, length, radius, and width.

01474 {
01475         if(radius > 0.0)
01476                 generateTubs(npoints, innerRadius, radius, 0.0, 360.0*deg,
01477                                                                 length, v);
01478         else
01479                 generateBox(npoints,width,height,length,v);
01480 }

G4bool BLCMDtrackerplane::isOutside ( G4ThreeVector &  local,
G4double  tolerance 
) [virtual]

isOutside() from BLElement

Implements BLElement.

References height, innerRadius, length, radius, and width.

01483 {
01484         if(radius > 0.0) {
01485                 G4double r = sqrt(local[0]*local[0]+local[1]*local[1]);
01486                 return r > radius-tolerance || r < innerRadius+tolerance ||
01487                         fabs(local[2]) > length/2.0-tolerance;
01488         } else {
01489                 return fabs(local[0]) > width/2.0-tolerance ||
01490                         fabs(local[1]) > height/2.0-tolerance ||
01491                         fabs(local[2]) > length/2.0-tolerance;
01492         }
01493 }


Friends And Related Function Documentation

friend class TrackerPlaneInstance [friend]

Referenced by construct().


Member Data Documentation

G4String BLCMDtrackerplane::tracker [private]

G4double BLCMDtrackerplane::radius [private]

G4double BLCMDtrackerplane::innerRadius [private]

G4double BLCMDtrackerplane::height [private]

G4double BLCMDtrackerplane::width [private]

G4double BLCMDtrackerplane::length [private]

G4String BLCMDtrackerplane::material [private]

G4String BLCMDtrackerplane::color [private]

G4VSolid* BLCMDtrackerplane::solid [private]

G4double BLCMDtrackerplane::maxStep [private]

G4double BLCMDtrackerplane::theta [private]

G4double BLCMDtrackerplane::wireSpacing [private]

G4double BLCMDtrackerplane::wireOffset [private]

G4String BLCMDtrackerplane::errType [private]

G4double BLCMDtrackerplane::errTheta [private]

G4double BLCMDtrackerplane::errSpacing [private]

G4double BLCMDtrackerplane::errOffset [private]

G4double BLCMDtrackerplane::sigmaTime [private]


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