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

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

BLCMDtrackerplane::~BLCMDtrackerplane (  )  [virtual]

Destructor.

References solid.

01311 {
01312         if(solid) delete solid;
01313 }

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.

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


Member Function Documentation

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

clone()

Implements BLElement.

References BLCMDtrackerplane().

00309 { return new BLCMDtrackerplane(*this); }

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

commandName() returns "trackerplane".

Implements BLCommand.

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

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

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.

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

void BLCMDtrackerplane::argChanged (  )  [virtual]

argChanged() does internal computations after some arg changed

Reimplemented from BLCommand.

References height, radius, and width.

01398 {
01399         if(radius > 0.0) width = height = radius*2.0;
01400 }

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.

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

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

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

Implements BLElement.

References length.

00336 { return length; }

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

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

Implements BLElement.

References width.

00339 { return width; }

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

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

Implements BLElement.

References height.

00342 { return height; }

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

isOK() returns true.

Implements BLElement.

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

01470 {
01471         if(radius > 0.0)
01472                 generateTubs(npoints, innerRadius, radius, 0.0, 360.0*deg,
01473                                                                 length, v);
01474         else
01475                 generateBox(npoints,width,height,length,v);
01476 }

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

isOutside() from BLElement

Implements BLElement.

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

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


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