#include <BaseParticlePropagator.h>
Public Member Functions | |
bool | backPropagate () |
BaseParticlePropagator (const RawParticle &myPart, double r, double z, double B) | |
BaseParticlePropagator () | |
Default c'tor. | |
BaseParticlePropagator (const RawParticle &myPart, double r, double z, double B, double t) | |
double | getMagneticField () const |
Get the magnetic field. | |
int | getSuccess () const |
Has propagation been performed and was barrel or endcap reached ? | |
bool | hasDecayed () const |
Has the particle decayed while propagated ? | |
double | helixCentreDistToAxis (double xC, double yC) const |
double | helixCentreDistToAxis () const |
The distance between the cylinder and the helix axes. | |
double | helixCentrePhi () const |
The azimuth if the vector joining the cylinder and the helix axes. | |
double | helixCentrePhi (double xC, double yC) const |
double | helixCentreX () const |
The x coordinate of the helix axis. | |
double | helixCentreX (double radius, double phi) const |
double | helixCentreY () const |
The y coordinate of the helix axis. | |
double | helixCentreY (double radius, double phi) const |
double | helixRadius () const |
The helix Radius. | |
double | helixRadius (double pT) const |
double | helixStartPhi () const |
The azimuth of the momentum at the vertex. | |
void | increaseRCyl (double delta) |
Just an internal trick. | |
void | init () |
Initialize internal switches and quantities. | |
bool | inside () const |
Is the vertex inside the cylinder ? (stricly inside : true) | |
bool | inside (double rPos2) const |
bool | onBarrel () const |
Is the vertex already on the cylinder barrel ? | |
bool | onBarrel (double rPos2) const |
bool | onEndcap () const |
Is the vertex already on the cylinder endcap ? | |
bool | onEndcap (double rPos2) const |
bool | onFiducial () const |
Is the vertex on some material ? | |
bool | onSurface () const |
Is the vertex already on the cylinder surface ? | |
bool | onSurface (double rPos2) const |
bool | propagate () |
BaseParticlePropagator | propagated () const |
bool | propagateToBeamCylinder (const XYZTLorentzVector &v, double radius=0.) |
bool | propagateToClosestApproach (double x0=0., double y0=0, bool first=true) |
bool | propagateToEcal (bool first=true) |
bool | propagateToEcalEntrance (bool first=true) |
bool | propagateToHcalEntrance (bool first=true) |
bool | propagateToHcalExit (bool first=true) |
bool | propagateToHOLayer (bool first=true) |
bool | propagateToNominalVertex (const XYZTLorentzVector &hit2=XYZTLorentzVector(0., 0., 0., 0.)) |
bool | propagateToPreshowerLayer1 (bool first=true) |
bool | propagateToPreshowerLayer2 (bool first=true) |
bool | propagateToVFcalEntrance (bool first=true) |
void | resetDebug () |
void | setDebug () |
Set the debug leve;. | |
void | setMagneticField (double b) |
Set the magnetic field. | |
void | setPropagationConditions (double r, double z, bool firstLoop=true) |
Set the propagation characteristics (rCyl, zCyl and first loop only) | |
void | setProperDecayTime (double t) |
Set the proper decay time. | |
double | xyImpactParameter (double x0=0., double y0=0.) const |
Transverse impact parameter. | |
double | zImpactParameter (double x0=0, double y0=0.) const |
Longitudinal impact parameter. | |
Protected Member Functions | |
double | c_light () const |
The speed of light in mm/ns (!) without clhep (yeaaahhh!) | |
Protected Attributes | |
bool | fiducial |
The particle traverses some real material. | |
int | success |
0:propagation still be done, 1:reached 'barrel', 2:reached 'endcaps' | |
Private Attributes | |
double | bField |
Magnetic field in the cylinder, oriented along the Z axis. | |
bool | debug |
The debug level. | |
bool | decayed |
The particle decayed while propagated ! | |
bool | firstLoop |
Do only the first half-loop. | |
int | propDir |
The propagation direction. | |
double | properDecayTime |
The proper decay time of the particle. | |
double | properTime |
The proper time of the particle. | |
double | rCyl |
Simulated particle that is to be resp has been propagated. | |
double | rCyl2 |
double | zCyl |
Half-height of the cylinder (centred at 0,0,0) to which propagation is done. |
This class is aimed at propagating charged and neutral particles (yet under the form of a RawParticle) from a given vertex to a cylinder, defined by a radius rCyl and a length 2*zCyl, centered in (0,0,0) and whose axis is parallel to the B field (B is oriented along z, by definition of the z axis).
The main method
bool propagate()
returns true if an intersection between the propagated RawParticle and the cylinder is found. The location of the intersection is given by the value of success:
* success = 2 : an intersection is found forward in the cylinder endcaps. The RawParticle vertex and momentum are overwritten by the values at this intersection.
* success = 1 : an intersection is found forward in the cylinder barrel. The RawParticle vertex and momentum are overwritten by the values at this intersection.
* success =-2 : an intersection is found backward in the cylinder endcaps. The RawParticle vertex and momentum are NOT overwritten by the values at this intersection.
* success =-1 : an intersection is found backward in the cylinder barrel. The RawParticle vertex and momentum are NOT overwritten by the values at this intersection.
* success = 0 : No intersection is found after half a loop (only possible if firstLoop=true). The vertex and momentum are overwritten by the values after this half loop.
The method
returns a new RawParticle with the propagated coordinates, if overwriting is not considered an advantage by the user.
Particles can be propagated backward with the method
Member functions
o propagateToPreshowerLayer1(), o propagateToPreshowerLayer2(), o propagateToEcalEntrance(), o propagateToHcalEntrance(), o propagateToHcalExit(), o propagateToClosestApproach(),
only need a momentum, a vertex and an electric charge to operate. Radii, half-lengths and default B field (4T) are set therein by default.
As of today, no average loss of energy (dE/dx, Brems) is considered in the propagation. No uncertainty (multiple scattering, dE/dx, Brems) is yet implemented.
Definition at line 82 of file BaseParticlePropagator.h.
BaseParticlePropagator::BaseParticlePropagator | ( | ) |
Default c'tor.
Definition at line 4 of file BaseParticlePropagator.cc.
: RawParticle(), rCyl(0.), rCyl2(0.), zCyl(0.), bField(0), properDecayTime(1E99) {;}
BaseParticlePropagator::BaseParticlePropagator | ( | const RawParticle & | myPart, |
double | r, | ||
double | z, | ||
double | B | ||
) |
Constructors taking as arguments a RawParticle, as well as the radius, half-height and magnetic field defining the cylinder for which propagation is to be performed, and optionally, the proper decay time
Definition at line 15 of file BaseParticlePropagator.cc.
References init().
BaseParticlePropagator::BaseParticlePropagator | ( | const RawParticle & | myPart, |
double | r, | ||
double | z, | ||
double | B, | ||
double | t | ||
) |
bool BaseParticlePropagator::backPropagate | ( | ) |
Update the current instance, after the back-propagation of the particle to the surface of the cylinder
Definition at line 320 of file BaseParticlePropagator.cc.
References RawParticle::charge(), run_regression::done, propagate(), propDir, and RawParticle::setCharge().
Referenced by propagated(), and propagateToClosestApproach().
double BaseParticlePropagator::c_light | ( | ) | const [inline, protected] |
The speed of light in mm/ns (!) without clhep (yeaaahhh!)
Definition at line 153 of file BaseParticlePropagator.h.
Referenced by helixRadius(), and propagateToBeamCylinder().
{ return 299.792458; }
double BaseParticlePropagator::getMagneticField | ( | ) | const [inline] |
Get the magnetic field.
Definition at line 303 of file BaseParticlePropagator.h.
References bField.
Referenced by TrajectoryManager::createPSimHits(), and ConvBremSeedProducer::produce().
{ return bField; }
int BaseParticlePropagator::getSuccess | ( | ) | const [inline] |
Has propagation been performed and was barrel or endcap reached ?
Definition at line 297 of file BaseParticlePropagator.h.
References success.
Referenced by PFTrackTransformer::addPointsAndBrems(), FBaseSimEvent::fill(), ConvBremSeedProducer::GoodCluster(), GoodSeedProducer::produce(), ConvBremSeedProducer::produce(), TrajectoryManager::propagateToCalorimeters(), and TrajectoryManager::propagateToLayer().
{ return success; }
bool BaseParticlePropagator::hasDecayed | ( | ) | const [inline] |
Has the particle decayed while propagated ?
Definition at line 294 of file BaseParticlePropagator.h.
References decayed.
Referenced by ParticlePropagator::propagateToBoundSurface(), and TrajectoryManager::propagateToCalorimeters().
{ return decayed; }
double BaseParticlePropagator::helixCentreDistToAxis | ( | ) | const [inline] |
The distance between the cylinder and the helix axes.
Definition at line 229 of file BaseParticlePropagator.h.
References helixCentreX(), helixCentreY(), and mathSSE::sqrt().
Referenced by propagate(), propagateToBeamCylinder(), propagateToClosestApproach(), and xyImpactParameter().
{ // The distance between the cylinder and the helix axes double xC = helixCentreX(); double yC = helixCentreY(); return std::sqrt( xC*xC + yC*yC ); }
double BaseParticlePropagator::helixCentreDistToAxis | ( | double | xC, |
double | yC | ||
) | const [inline] |
Definition at line 236 of file BaseParticlePropagator.h.
References mathSSE::sqrt().
{ // Faster version of helixCentreDistToAxis return std::sqrt( xC*xC + yC*yC ); }
double BaseParticlePropagator::helixCentrePhi | ( | double | xC, |
double | yC | ||
) | const [inline] |
Definition at line 249 of file BaseParticlePropagator.h.
{ // Faster version of helixCentrePhi() return xC == 0.0 && yC == 0.0 ? 0.0 : std::atan2(yC,xC); }
double BaseParticlePropagator::helixCentrePhi | ( | ) | const [inline] |
The azimuth if the vector joining the cylinder and the helix axes.
Definition at line 242 of file BaseParticlePropagator.h.
References helixCentreX(), and helixCentreY().
Referenced by propagate().
{ // The azimuth if the vector joining the cylinder and the helix axes double xC = helixCentreX(); double yC = helixCentreY(); return xC == 0.0 && yC == 0.0 ? 0.0 : std::atan2(yC,xC); }
double BaseParticlePropagator::helixCentreX | ( | double | radius, |
double | phi | ||
) | const [inline] |
Definition at line 212 of file BaseParticlePropagator.h.
References funct::sin(), and RawParticle::X().
double BaseParticlePropagator::helixCentreX | ( | ) | const [inline] |
The x coordinate of the helix axis.
Definition at line 207 of file BaseParticlePropagator.h.
References helixRadius(), helixStartPhi(), funct::sin(), and RawParticle::X().
Referenced by helixCentreDistToAxis(), helixCentrePhi(), propagate(), propagateToBeamCylinder(), propagateToClosestApproach(), and xyImpactParameter().
{ // The x coordinate of the helix axis return X() - helixRadius() * std::sin ( helixStartPhi() ); }
double BaseParticlePropagator::helixCentreY | ( | ) | const [inline] |
The y coordinate of the helix axis.
Definition at line 218 of file BaseParticlePropagator.h.
References funct::cos(), helixRadius(), helixStartPhi(), and RawParticle::Y().
Referenced by helixCentreDistToAxis(), helixCentrePhi(), propagate(), propagateToBeamCylinder(), propagateToClosestApproach(), and xyImpactParameter().
{ // The y coordinate of the helix axis return Y() + helixRadius() * std::cos ( helixStartPhi() ); }
double BaseParticlePropagator::helixCentreY | ( | double | radius, |
double | phi | ||
) | const [inline] |
Definition at line 223 of file BaseParticlePropagator.h.
References funct::cos(), and RawParticle::Y().
double BaseParticlePropagator::helixRadius | ( | ) | const [inline] |
The helix Radius.
Definition at line 183 of file BaseParticlePropagator.h.
References bField, c_light(), RawParticle::charge(), alignCSCRings::e, and reco::tau::disc::Pt().
Referenced by helixCentreX(), helixCentreY(), propagate(), propagateToClosestApproach(), propagateToNominalVertex(), and xyImpactParameter().
{ // The helix Radius // // The helix' Radius sign accounts for the orientation of the magnetic field // (+ = along z axis) and the sign of the electric charge of the particle. // It signs the rotation of the (charged) particle around the z axis: // Positive means anti-clockwise, negative means clockwise rotation. // // The radius is returned in cm to match the units in RawParticle. return charge() == 0 ? 0.0 : - Pt() / ( c_light() * 1e-5 * bField * charge() ); }
double BaseParticlePropagator::helixRadius | ( | double | pT | ) | const [inline] |
Definition at line 195 of file BaseParticlePropagator.h.
References bField, c_light(), RawParticle::charge(), and alignCSCRings::e.
double BaseParticlePropagator::helixStartPhi | ( | ) | const [inline] |
The azimuth of the momentum at the vertex.
Definition at line 201 of file BaseParticlePropagator.h.
Referenced by helixCentreX(), helixCentreY(), propagate(), propagateToClosestApproach(), and xyImpactParameter().
{ // The azimuth of the momentum at the vertex return Px() == 0.0 && Py() == 0.0 ? 0.0 : std::atan2(Py(),Px()); }
void BaseParticlePropagator::increaseRCyl | ( | double | delta | ) | [inline] |
void BaseParticlePropagator::init | ( | void | ) |
Initialize internal switches and quantities.
Reimplemented from RawParticle.
Definition at line 23 of file BaseParticlePropagator.cc.
References debug, decayed, fiducial, firstLoop, propDir, properTime, and success.
Referenced by BaseParticlePropagator().
{ // The status of the propagation success = 0; // Propagate only for the first half-loop firstLoop = true; // fiducial = true; // The particle has not yet decayed decayed = false; // The proper time is set to zero properTime = 0.; // The propagation direction is along the momentum propDir = 1; // debug = false; }
bool BaseParticlePropagator::inside | ( | ) | const [inline] |
Is the vertex inside the cylinder ? (stricly inside : true)
Definition at line 255 of file BaseParticlePropagator.h.
References RawParticle::R2(), rCyl, rCyl2, RawParticle::Z(), and zCyl.
Referenced by propagate().
bool BaseParticlePropagator::inside | ( | double | rPos2 | ) | const [inline] |
Definition at line 258 of file BaseParticlePropagator.h.
References rCyl, rCyl2, RawParticle::Z(), and zCyl.
bool BaseParticlePropagator::onBarrel | ( | ) | const [inline] |
Is the vertex already on the cylinder barrel ?
Definition at line 272 of file BaseParticlePropagator.h.
References RawParticle::R2(), rCyl, rCyl2, RawParticle::Z(), and zCyl.
Referenced by onSurface(), and propagate().
bool BaseParticlePropagator::onBarrel | ( | double | rPos2 | ) | const [inline] |
Definition at line 277 of file BaseParticlePropagator.h.
References rCyl, rCyl2, RawParticle::Z(), and zCyl.
bool BaseParticlePropagator::onEndcap | ( | double | rPos2 | ) | const [inline] |
Definition at line 286 of file BaseParticlePropagator.h.
References rCyl2, RawParticle::Z(), and zCyl.
bool BaseParticlePropagator::onEndcap | ( | ) | const [inline] |
Is the vertex already on the cylinder endcap ?
Definition at line 282 of file BaseParticlePropagator.h.
References RawParticle::R2(), rCyl2, RawParticle::Z(), and zCyl.
Referenced by onSurface(), and propagate().
bool BaseParticlePropagator::onFiducial | ( | ) | const [inline] |
Is the vertex on some material ?
Definition at line 291 of file BaseParticlePropagator.h.
References fiducial.
Referenced by TrajectoryManager::propagateToLayer().
{ return fiducial; }
bool BaseParticlePropagator::onSurface | ( | double | rPos2 | ) | const [inline] |
Definition at line 267 of file BaseParticlePropagator.h.
References onBarrel(), and onEndcap().
bool BaseParticlePropagator::onSurface | ( | ) | const [inline] |
Is the vertex already on the cylinder surface ?
Definition at line 263 of file BaseParticlePropagator.h.
References onBarrel(), and onEndcap().
Referenced by propagate().
bool BaseParticlePropagator::propagate | ( | ) |
Update the current instance, after the propagation of the particle to the surface of the cylinder
Definition at line 43 of file BaseParticlePropagator.cc.
References bField, RawParticle::charge(), funct::cos(), decayed, delta, SiPixelRawToDigiRegional_cfi::deltaPhi, firstLoop, helixCentreDistToAxis(), helixCentrePhi(), helixCentreX(), helixCentreY(), helixRadius(), helixStartPhi(), inside(), M_PI, RawParticle::mass(), min, onBarrel(), onEndcap(), onSurface(), propDir, properDecayTime, properTime, reco::tau::disc::Pt(), RawParticle::R2(), CosmicsPD_Skims::radius, rCyl, rCyl2, RawParticle::setVertex(), funct::sin(), mathSSE::sqrt(), success, RawParticle::T(), RawParticle::X(), RawParticle::Y(), RawParticle::Z(), and zCyl.
Referenced by PFTrackTransformer::addPointsAndBrems(), backPropagate(), ConvBremSeedProducer::produce(), propagated(), ParticlePropagator::propagateToBoundSurface(), propagateToEcal(), propagateToEcalEntrance(), propagateToHcalEntrance(), propagateToHcalExit(), propagateToHOLayer(), FastPixelHitMatcher::propagateToLayer(), propagateToPreshowerLayer1(), propagateToPreshowerLayer2(), and propagateToVFcalEntrance().
{ // // Check that the particle is not already on the cylinder surface // double rPos2 = R2(); if ( onBarrel(rPos2) ) { success = 1; return true; } // if ( onEndcap(rPos2) ) { success = 2; return true; } // // Treat neutral particles (or no magnetic field) first // if ( fabs(charge()) < 1E-12 || bField < 1E-12 ) { // // First check that the particle crosses the cylinder // double pT2 = Perp2(); // double pT2 = pT*pT; // double b2 = rCyl * pT; double b2 = rCyl2 * pT2; double ac = fabs ( X()*Py() - Y()*Px() ); double ac2 = ac*ac; // // The particle never crosses (or never crossed) the cylinder // // if ( ac > b2 ) return false; if ( ac2 > b2 ) return false; // double delta = std::sqrt(b2*b2 - ac*ac); double delta = std::sqrt(b2 - ac2); double tplus = -( X()*Px()+Y()*Py() ) + delta; double tminus = -( X()*Px()+Y()*Py() ) - delta; // // Find the first (time-wise) intersection point with the cylinder // double solution = tminus < 0 ? tplus/pT2 : tminus/pT2; if ( solution < 0. ) return false; // // Check that the particle does (did) not exit from the endcaps first. // double zProp = Z() + Pz() * solution; if ( fabs(zProp) > zCyl ) { tplus = ( zCyl - Z() ) / Pz(); tminus = (-zCyl - Z() ) / Pz(); solution = tminus < 0 ? tplus : tminus; if ( solution < 0. ) return false; zProp = Z() + Pz() * solution; success = 2; } else { success = 1; } // // Check the decay time // double delTime = propDir * mass() * solution; double factor = 1.; properTime += delTime; if ( properTime > properDecayTime ) { factor = 1.-(properTime-properDecayTime)/delTime; properTime = properDecayTime; decayed = true; } // // Compute the coordinates of the RawParticle after propagation // double xProp = X() + Px() * solution * factor; double yProp = Y() + Py() * solution * factor; zProp = Z() + Pz() * solution * factor; double tProp = T() + E() * solution * factor; // // Last check : Although propagated to the endcaps, R could still be // larger than rCyl because the cylinder was too short... // if ( success == 2 && xProp*xProp+yProp*yProp > rCyl2 ) { success = -1; return true; } // // ... and update the particle with its propagated coordinates // setVertex( XYZTLorentzVector(xProp,yProp,zProp,tProp) ); return true; } // // Treat charged particles with magnetic field next. // else { // // First check that the particle can cross the cylinder // double pT = Pt(); double radius = helixRadius(pT); double phi0 = helixStartPhi(); double xC = helixCentreX(radius,phi0); double yC = helixCentreY(radius,phi0); double dist = helixCentreDistToAxis(xC,yC); // // The helix either is outside or includes the cylinder -> no intersection // if ( fabs ( fabs(radius) - dist ) > rCyl ) return false; // // The particle is already away from the endcaps, and will never come back // Could be back-propagated, so warn the user that it is so. // if ( Z() * Pz() > zCyl * fabs(Pz()) ) { success = -2; return true; } double pZ = Pz(); double phiProp, zProp; // // Does the helix cross the cylinder barrel ? If not, do the first half-loop // double rProp = std::min(fabs(radius)+dist-0.000001, rCyl); // Check for rounding errors in the ArcSin. double sinPhiProp = (rProp*rProp-radius*radius-dist*dist)/( 2.*dist*radius); // double deltaPhi = 1E99; // Taylor development up to third order for large momenta if ( 1.-fabs(sinPhiProp) < 1E-9 ) { double cphi0 = std::cos(phi0); double sphi0 = std::sin(phi0); double r0 = (X()*cphi0 + Y()*sphi0)/radius; double q0 = (X()*sphi0 - Y()*cphi0)/radius; double rcyl2 = (rCyl2 - X()*X() - Y()*Y())/(radius*radius); double delta = r0*r0 + rcyl2*(1.-q0); // This is a solution of a second order equation, assuming phi = phi0 + epsilon // and epsilon is small. deltaPhi = radius > 0 ? ( -r0 + std::sqrt(delta ) ) / (1.-q0) : ( -r0 - std::sqrt(delta ) ) / (1.-q0); } // Use complete calculation otherwise, or if the delta phi is large anyway. if ( fabs(deltaPhi) > 1E-3 ) { // phiProp = fabs(sinPhiProp) < 0.99999999 ? // std::asin(sinPhiProp) : M_PI/2. * sinPhiProp; phiProp = std::asin(sinPhiProp); // // Compute phi after propagation (two solutions, but the asin definition // (between -pi/2 and +pi/2) takes care of the wrong one. // phiProp = helixCentrePhi(xC,yC) + ( inside(rPos2) || onSurface(rPos2) ? phiProp : M_PI-phiProp ); // // Solve the obvious two-pi ambiguities: more than one turn! // if ( fabs(phiProp - phi0) > 2.*M_PI ) phiProp += ( phiProp > phi0 ? -2.*M_PI : +2.*M_PI ); // // Check that the phi difference has the right sign, according to Q*B // (Another two pi ambuiguity, slightly more subtle) // if ( (phiProp-phi0)*radius < 0.0 ) radius > 0.0 ? phiProp += 2 * M_PI : phiProp -= 2 * M_PI; } else { // Use Taylor phiProp = phi0 + deltaPhi; } // // Check that the particle does not exit from the endcaps first. // zProp = Z() + ( phiProp - phi0 ) * pZ * radius / pT ; if ( fabs(zProp) > zCyl ) { zProp = zCyl * fabs(pZ)/pZ; phiProp = phi0 + (zProp - Z()) / radius * pT / pZ; success = 2; } else { // // If the particle does not cross the barrel, either process the // first-half loop only or propagate all the way down to the endcaps // if ( rProp < rCyl ) { if ( firstLoop || fabs(pZ)/pT < 1E-10 ) { success = 0; } else { zProp = zCyl * fabs(pZ)/pZ; phiProp = phi0 + (zProp - Z()) / radius * pT / pZ; success = 2; } // // The particle crossed the barrel // } else { success = 1; } } // // Compute the coordinates of the RawParticle after propagation // // // Check the decay time // double delTime = propDir * (phiProp-phi0)*radius*mass() / pT; double factor = 1.; properTime += delTime; if ( properTime > properDecayTime ) { factor = 1.-(properTime-properDecayTime)/delTime; properTime = properDecayTime; decayed = true; } zProp = Z() + (zProp-Z())*factor; phiProp = phi0 + (phiProp-phi0)*factor; double sProp = std::sin(phiProp); double cProp = std::cos(phiProp); double xProp = xC + radius * sProp; double yProp = yC - radius * cProp; double tProp = T() + (phiProp-phi0)*radius*E()/pT; double pxProp = pT * cProp; double pyProp = pT * sProp; // // Last check : Although propagated to the endcaps, R could still be // larger than rCyl because the cylinder was too short... // if ( success == 2 && xProp*xProp+yProp*yProp > rCyl2 ) { success = -1; return true; } // // ... and update the particle vertex and momentum // setVertex( XYZTLorentzVector(xProp,yProp,zProp,tProp) ); SetXYZT(pxProp,pyProp,pZ,E()); // SetPx(pxProp); // SetPy(pyProp); return true; } }
BaseParticlePropagator BaseParticlePropagator::propagated | ( | ) | const |
Return a new instance, corresponding to the particle propagated to the surface of the cylinder
Reimplemented in ParticlePropagator.
Definition at line 335 of file BaseParticlePropagator.cc.
References backPropagate(), firstLoop, propagate(), and success.
Referenced by ParticlePropagator::propagated().
{ // Copy the input particle in a new instance BaseParticlePropagator myPropPart(*this); // Allow for many loops myPropPart.firstLoop=false; // Propagate the particle ! myPropPart.propagate(); // Back propagate if the forward propagation was not a success if (myPropPart.success < 0 ) myPropPart.backPropagate(); // Return the new instance return myPropPart; }
bool BaseParticlePropagator::propagateToBeamCylinder | ( | const XYZTLorentzVector & | v, |
double | radius = 0. |
||
) |
Definition at line 630 of file BaseParticlePropagator.cc.
References a, b, bField, trackerHits::c, c_light(), RawParticle::charge(), funct::cos(), gather_cfg::cout, alignCSCRings::e, helixCentreDistToAxis(), helixCentreX(), helixCentreY(), i, propagateToClosestApproach(), reco::tau::disc::Pt(), RawParticle::r(), RawParticle::r2(), CosmicsPD_Skims::radius, RawParticle::setCharge(), funct::sin(), mathSSE::sqrt(), RawParticle::X(), RawParticle::Y(), and RawParticle::Z().
Referenced by TrajectorySeedProducer::compatibleWithBeamAxis().
{ // For neutral (or BField = 0, simply check that the track passes through the cylinder if ( charge() == 0. || bField == 0.) return fabs( Px() * Y() - Py() * X() ) / Pt() < radius; // Now go to the charged particles // A few needed intermediate variables (to make the code understandable) // The inner cylinder double r = radius; double r2 = r*r; double r4 = r2*r2; // The two hits double dx = X()-v.X(); double dy = Y()-v.Y(); double dz = Z()-v.Z(); double Sxy = X()*v.X() + Y()*v.Y(); double Dxy = Y()*v.X() - X()*v.Y(); double Dxy2 = Dxy*Dxy; double Sxy2 = Sxy*Sxy; double SDxy = dx*dx + dy*dy; double SSDxy = std::sqrt(SDxy); double ATxy = std::atan2(dy,dx); // The coefficients of the second order equation for the trajectory radius double a = r2 - Dxy2/SDxy; double b = r * (r2 - Sxy); double c = r4 - 2.*Sxy*r2 + Sxy2 + Dxy2; // Here are the four possible solutions for // 1) The trajectory radius std::vector<double> helixRs(4,static_cast<double>(0.)); helixRs[0] = (b - std::sqrt(b*b - a*c))/(2.*a); helixRs[1] = (b + std::sqrt(b*b - a*c))/(2.*a); helixRs[2] = -helixRs[0]; helixRs[3] = -helixRs[1]; // 2) The azimuthal direction at the second point std::vector<double> helixPhis(4,static_cast<double>(0.)); helixPhis[0] = std::asin ( SSDxy/(2.*helixRs[0]) ); helixPhis[1] = std::asin ( SSDxy/(2.*helixRs[1]) ); helixPhis[2] = -helixPhis[0]; helixPhis[3] = -helixPhis[1]; // Only two solutions are valid though double solution1=0.; double solution2=0.; double phi1=0.; double phi2=0.; // Loop over the four possibilities for ( unsigned int i=0; i<4; ++i ) { helixPhis[i] = ATxy + helixPhis[i]; // The centre of the helix double xC = helixCentreX(helixRs[i],helixPhis[i]); double yC = helixCentreY(helixRs[i],helixPhis[i]); // The distance of that centre to the origin double dist = helixCentreDistToAxis(xC, yC); /* std::cout << "Solution "<< i << " = " << helixRs[i] << " accepted ? " << fabs(fabs(helixRs[i]) - dist ) << " - " << radius << " = " << fabs(fabs(helixRs[i]) - dist ) - fabs(radius) << " " << ( fabs( fabs(fabs(helixRs[i]) - dist) -fabs(radius) ) < 1e-6) << std::endl; */ // Check that the propagation will lead to a trajectroy tangent to the inner cylinder if ( fabs( fabs(fabs(helixRs[i]) - dist) -fabs(radius) ) < 1e-6 ) { // Fill first solution if ( solution1 == 0. ) { solution1 = helixRs[i]; phi1 = helixPhis[i]; // Fill second solution (ordered) } else if ( solution2 == 0. ) { if ( helixRs[i] < solution1 ) { solution2 = solution1; solution1 = helixRs[i]; phi2 = phi1; phi1 = helixPhis[i]; } else { solution2 = helixRs[i]; phi2 = helixPhis[i]; } // Must not happen! } else { std::cout << "warning !!! More than two solutions for this track !!! " << std::endl; } } } // Find the largest possible pT compatible with coming from within the inne cylinder double pT = 0.; double helixPhi = 0.; double helixR = 0.; // This should not happen if ( solution1*solution2 == 0. ) { std::cout << "warning !!! Less than two solution for this track! " << solution1 << " " << solution2 << std::endl; return false; // If the two solutions have opposite sign, it means that pT = infinity is a possibility. } else if ( solution1*solution2 < 0. ) { pT = 1000.; double norm = pT/SSDxy; setCharge(+1.); SetXYZT(dx*norm,dy*norm,dz*norm,0.); SetE(std::sqrt(Vect().Mag2())); // Otherwise take the solution that gives the largest transverse momentum } else { if (solution1<0.) { helixR = solution1; helixPhi = phi1; setCharge(+1.); } else { helixR = solution2; helixPhi = phi2; setCharge(-1.); } pT = fabs(helixR) * 1e-5 * c_light() *bField; double norm = pT/SSDxy; SetXYZT(pT*std::cos(helixPhi),pT*std::sin(helixPhi),dz*norm,0.); SetE(std::sqrt(Vect().Mag2())); } // Propagate to closest approach to get the Z value (a bit of an overkill) return propagateToClosestApproach(); }
bool BaseParticlePropagator::propagateToClosestApproach | ( | double | x0 = 0. , |
double | y0 = 0 , |
||
bool | first = true |
||
) |
Update the particle after propagation to the closest approach from Z axis, to the preshower layer 1 & 2, to the ECAL entrance, to the HCAL entrance, the HCAL 2nd and 3rd layer (not coded yet), the VFCAL entrance, or any BoundSurface(disk or cylinder)
Reimplemented in ParticlePropagator.
Definition at line 358 of file BaseParticlePropagator.cc.
References backPropagate(), bField, RawParticle::charge(), run_regression::done, helixCentreDistToAxis(), helixCentreX(), helixCentreY(), helixRadius(), helixStartPhi(), RawParticle::momentum(), reco::tau::disc::Pt(), CosmicsPD_Skims::radius, setPropagationConditions(), RawParticle::setVertex(), mathSSE::sqrt(), RawParticle::vertex(), RawParticle::X(), RawParticle::Y(), and RawParticle::z().
Referenced by propagateToBeamCylinder(), and propagateToNominalVertex().
{ // Pre-computed quantities double pT = Pt(); double radius = helixRadius(pT); double phi0 = helixStartPhi(); // The Helix centre coordinates are computed wrt to the z axis // Actually the z axis might be centred on 0,0: it's the beam spot position (x0,y0)! double xC = helixCentreX(radius,phi0); double yC = helixCentreY(radius,phi0); double distz = helixCentreDistToAxis(xC-x0,yC-y0); double dist0 = helixCentreDistToAxis(xC,yC); // // Propagation to point of clostest approach to z axis // double rz,r0,z; if ( charge() != 0.0 && bField != 0.0 ) { rz = fabs ( fabs(radius) - distz ) + std::sqrt(x0*x0+y0*y0) + 0.0000001; r0 = fabs ( fabs(radius) - dist0 ) + 0.0000001; } else { rz = fabs( Px() * (Y()-y0) - Py() * (X()-x0) ) / Pt(); r0 = fabs( Px() * Y() - Py() * X() ) / Pt(); } z = 999999.; // Propagate to the first interesection point // with cylinder of radius sqrt(x0**2+y0**2) setPropagationConditions(rz , z, first); bool done = backPropagate(); // Unsuccessful propagation - should not happen! if ( !done ) return done; // The z axis is (0,0) - no need to go further if ( fabs(rz-r0) < 1E-10 ) return done; double dist1 = (X()-x0)*(X()-x0) + (Y()-y0)*(Y()-y0); // We are already at closest approach - no need to go further if ( dist1 < 1E-10 ) return done; // Keep for later if it happens to be the right solution XYZTLorentzVector vertex1 = vertex(); XYZTLorentzVector momentum1 = momentum(); // Propagate to the distance of closest approach to (0,0) setPropagationConditions(r0 , z, first); done = backPropagate(); if ( !done ) return done; // Propagate to the first interesection point // with cylinder of radius sqrt(x0**2+y0**2) setPropagationConditions(rz , z, first); done = backPropagate(); if ( !done ) return done; double dist2 = (X()-x0)*(X()-x0) + (Y()-y0)*(Y()-y0); // Keep the good solution. if ( dist2 > dist1 ) { setVertex(vertex1); SetXYZT(momentum1.X(),momentum1.Y(),momentum1.Z(),momentum1.E()); dist2 = dist1; } // Done return done; }
bool BaseParticlePropagator::propagateToEcal | ( | bool | first = true | ) |
Definition at line 430 of file BaseParticlePropagator.cc.
References propagate(), and setPropagationConditions().
{ // // Propagation to Ecal (including preshower) after the // last Tracker layer // TODO: include proper geometry // Geometry taken from CMS ECAL TDR // // First propagate to global barrel / endcap cylinder // setPropagationConditions(1290. , 3045 , first); // New position of the preshower as of CMSSW_1_7_X setPropagationConditions(129.0 , 303.353 , first); return propagate(); }
bool BaseParticlePropagator::propagateToEcalEntrance | ( | bool | first = true | ) |
Definition at line 487 of file BaseParticlePropagator.cc.
References RawParticle::cos2ThetaV(), run_regression::done, propagate(), setPropagationConditions(), and success.
Referenced by PFTrackTransformer::addPointsAndBrems(), FBaseSimEvent::fill(), ConvBremSeedProducer::GoodCluster(), GoodSeedProducer::produce(), and TrajectoryManager::propagateToCalorimeters().
{ // // Propagation to ECAL entrance // TODO: include proper geometry // Geometry taken from CMS ECAL TDR // // First propagate to global barrel / endcap cylinder setPropagationConditions(129.0 , 320.9,first); bool done = propagate(); // Go to endcap cylinder in the "barrel cut corner" // eta = 1.479 -> cos^2(theta) = 0.81230 // if ( done && eta > 1.479 && success == 1 ) { if ( done && cos2ThetaV() > 0.81230 && success == 1 ) { setPropagationConditions(152.6 , 320.9, first); done = propagate(); } // We are not in the ECAL acceptance // eta = 3.0 -> cos^2(theta) = 0.99013 if ( cos2ThetaV() > 0.99014 ) success = 0; return done; }
bool BaseParticlePropagator::propagateToHcalEntrance | ( | bool | first = true | ) |
Definition at line 513 of file BaseParticlePropagator.cc.
References RawParticle::cos2ThetaV(), run_regression::done, propagate(), propDir, setPropagationConditions(), and success.
Referenced by PFTrackTransformer::addPointsAndBrems(), FBaseSimEvent::fill(), and TrajectoryManager::propagateToCalorimeters().
{ // // Propagation to HCAL entrance // TODO: include proper geometry // Geometry taken from CMSSW_3_1_X xml (Sunanda) // // First propagate to global barrel / endcap cylinder setPropagationConditions(177.5 , 335.0, first); propDir = 0; bool done = propagate(); propDir = 1; // If went through the bottom of HB cylinder -> re-propagate to HE surface if (done && success == 2) { setPropagationConditions(300.0, 400.458, first); propDir = 0; done = propagate(); propDir = 1; } // out of the HB/HE acceptance // eta = 3.0 -> cos^2(theta) = 0.99014 if ( done && cos2ThetaV() > 0.99014 ) success = 0; return done; }
bool BaseParticlePropagator::propagateToHcalExit | ( | bool | first = true | ) |
Definition at line 566 of file BaseParticlePropagator.cc.
References run_regression::done, propagate(), propDir, and setPropagationConditions().
Referenced by PFTrackTransformer::addPointsAndBrems(), and FBaseSimEvent::fill().
{ // // Propagation to HCAL exit // TODO: include proper geometry // Geometry taken from CMSSW_3_1_X xml (Sunanda) // // Approximate it to a single cylinder as it is not that crucial. setPropagationConditions(263.9 , 554.1, first); // this->rawPart().setCharge(0.0); ?? Shower Propagation ?? propDir = 0; bool done = propagate(); propDir = 1; return done; }
bool BaseParticlePropagator::propagateToHOLayer | ( | bool | first = true | ) |
Definition at line 584 of file BaseParticlePropagator.cc.
References abs, run_regression::done, propagate(), propDir, setPropagationConditions(), success, and RawParticle::Z().
Referenced by PFTrackTransformer::addPointsAndBrems(), and FBaseSimEvent::fill().
{ // // Propagation to Layer0 of HO (Ring-0) // TODO: include proper geometry // Geometry taken from CMSSW_3_1_X xml (Sunanda) // // Approximate it to a single cylinder as it is not that crucial. setPropagationConditions(387.6, 1000.0, first); //Dia is the middle of HO \sqrt{384.8**2+46.2**2} // setPropagationConditions(410.2, 1000.0, first); //sqrt{406.6**2+54.2**2} //for outer layer // this->rawPart().setCharge(0.0); ?? Shower Propagation ?? propDir = 0; bool done = propagate(); propDir = 1; if ( done && abs(Z()) > 700.25) success = 0; return done; }
bool BaseParticlePropagator::propagateToNominalVertex | ( | const XYZTLorentzVector & | hit2 = XYZTLorentzVector(0.,0.,0.,0.) | ) |
Reimplemented in ParticlePropagator.
Definition at line 606 of file BaseParticlePropagator.cc.
References bField, RawParticle::charge(), funct::cos(), helixRadius(), phi, propagateToClosestApproach(), funct::sin(), mathSSE::sqrt(), RawParticle::X(), and RawParticle::Y().
{ // Not implemented for neutrals (used for electrons only) if ( charge() == 0. || bField == 0.) return false; // Define the proper pT direction to meet the point (vx,vy) in (x,y) double dx = X()-v.X(); double dy = Y()-v.Y(); double phi = std::atan2(dy,dx) + std::asin ( std::sqrt(dx*dx+dy*dy)/(2.*helixRadius()) ); // The absolute pT (kept to the original value) double pT = pt(); // Set the new pT SetPx(pT*std::cos(phi)); SetPy(pT*std::sin(phi)); return propagateToClosestApproach(v.X(),v.Y()); }
bool BaseParticlePropagator::propagateToPreshowerLayer1 | ( | bool | first = true | ) |
Definition at line 446 of file BaseParticlePropagator.cc.
References run_regression::done, propagate(), RawParticle::R2(), setPropagationConditions(), and success.
Referenced by PFTrackTransformer::addPointsAndBrems(), FBaseSimEvent::fill(), TrajectoryManager::propagateToCalorimeters(), and GoodSeedProducer::PSforTMVA().
{ // // Propagation to Preshower Layer 1 // TODO: include proper geometry // Geometry taken from CMS ECAL TDR // // First propagate to global barrel / endcap cylinder // setPropagationConditions(1290., 3045 , first); // New position of the preshower as of CMSSW_1_7_X setPropagationConditions(129.0, 303.353 , first); bool done = propagate(); // Check that were are on the Layer 1 if ( done && (R2() > 125.0*125.0 || R2() < 45.0*45.0) ) success = 0; return done; }
bool BaseParticlePropagator::propagateToPreshowerLayer2 | ( | bool | first = true | ) |
Definition at line 466 of file BaseParticlePropagator.cc.
References run_regression::done, propagate(), RawParticle::R2(), setPropagationConditions(), and success.
Referenced by PFTrackTransformer::addPointsAndBrems(), FBaseSimEvent::fill(), and TrajectoryManager::propagateToCalorimeters().
{ // // Propagation to Preshower Layer 2 // TODO: include proper geometry // Geometry taken from CMS ECAL TDR // // First propagate to global barrel / endcap cylinder // setPropagationConditions(1290. , 3090 , first); // New position of the preshower as of CMSSW_1_7_X setPropagationConditions(129.0 , 307.838 , first); bool done = propagate(); // Check that we are on Layer 2 if ( done && (R2() > 125.0*125.0 || R2() < 45.0*45.0 ) ) success = 0; return done; }
bool BaseParticlePropagator::propagateToVFcalEntrance | ( | bool | first = true | ) |
Definition at line 543 of file BaseParticlePropagator.cc.
References RawParticle::cos2ThetaV(), run_regression::done, propagate(), propDir, setPropagationConditions(), and success.
Referenced by FBaseSimEvent::fill(), and TrajectoryManager::propagateToCalorimeters().
{ // // Propagation to VFCAL entrance // TODO: include proper geometry // Geometry taken from DAQ TDR Chapter 13 setPropagationConditions(400.0 , 1110.0, first); propDir = 0; bool done = propagate(); propDir = 1; if (!done) success = 0; // We are not in the VFCAL acceptance // eta = 3.0 -> cos^2(theta) = 0.99014 // eta = 5.2 -> cos^2(theta) = 0.9998755 double c2teta = cos2ThetaV(); if ( done && ( c2teta < 0.99014 || c2teta > 0.9998755 ) ) success = 0; return done; }
void BaseParticlePropagator::resetDebug | ( | ) | [inline] |
void BaseParticlePropagator::setDebug | ( | ) | [inline] |
Set the debug leve;.
Definition at line 306 of file BaseParticlePropagator.h.
References debug.
{ debug = true; }
void BaseParticlePropagator::setMagneticField | ( | double | b | ) | [inline] |
Set the magnetic field.
Definition at line 300 of file BaseParticlePropagator.h.
Referenced by PFTrackTransformer::addPointsAndBrems(), FBaseSimEvent::fill(), ParticlePropagator::ParticlePropagator(), ParticlePropagator::propagateToBoundSurface(), ParticlePropagator::propagateToClosestApproach(), and ParticlePropagator::propagateToNominalVertex().
void BaseParticlePropagator::setPropagationConditions | ( | double | r, |
double | z, | ||
bool | firstLoop = true |
||
) |
Set the propagation characteristics (rCyl, zCyl and first loop only)
Definition at line 350 of file BaseParticlePropagator.cc.
References f, firstLoop, RawParticle::R(), rCyl, rCyl2, RawParticle::Z(), and zCyl.
Referenced by PFTrackTransformer::addPointsAndBrems(), propagateToClosestApproach(), propagateToEcal(), propagateToEcalEntrance(), propagateToHcalEntrance(), propagateToHcalExit(), propagateToHOLayer(), FastPixelHitMatcher::propagateToLayer(), propagateToPreshowerLayer1(), propagateToPreshowerLayer2(), and propagateToVFcalEntrance().
void BaseParticlePropagator::setProperDecayTime | ( | double | t | ) | [inline] |
Set the proper decay time.
Definition at line 168 of file BaseParticlePropagator.h.
References properDecayTime, and RawParticle::t().
Referenced by ParticlePropagator::initProperDecayTime(), and ParticlePropagator::ParticlePropagator().
{ properDecayTime = t; }
double BaseParticlePropagator::xyImpactParameter | ( | double | x0 = 0. , |
double | y0 = 0. |
||
) | const |
Transverse impact parameter.
Definition at line 755 of file BaseParticlePropagator.cc.
References bField, RawParticle::charge(), helixCentreDistToAxis(), helixCentreX(), helixCentreY(), helixRadius(), helixStartPhi(), reco::tau::disc::Pt(), CosmicsPD_Skims::radius, RawParticle::X(), and RawParticle::Y().
Referenced by TrajectorySeedProducer::produce().
{ double ip=0.; double pT = Pt(); if ( charge() != 0.0 && bField != 0.0 ) { double radius = helixRadius(pT); double phi0 = helixStartPhi(); // The Helix centre coordinates are computed wrt to the z axis double xC = helixCentreX(radius,phi0); double yC = helixCentreY(radius,phi0); double distz = helixCentreDistToAxis(xC-x0,yC-y0); ip = distz - fabs(radius); } else { ip = fabs( Px() * (Y()-y0) - Py() * (X()-x0) ) / pT; } return ip; }
double BaseParticlePropagator::zImpactParameter | ( | double | x0 = 0 , |
double | y0 = 0. |
||
) | const [inline] |
Longitudinal impact parameter.
Definition at line 177 of file BaseParticlePropagator.h.
References mathSSE::sqrt(), RawParticle::X(), RawParticle::Y(), and RawParticle::Z().
Referenced by TrajectorySeedProducer::produce().
double BaseParticlePropagator::bField [private] |
Magnetic field in the cylinder, oriented along the Z axis.
Definition at line 140 of file BaseParticlePropagator.h.
Referenced by getMagneticField(), helixRadius(), propagate(), propagateToBeamCylinder(), propagateToClosestApproach(), propagateToNominalVertex(), setMagneticField(), and xyImpactParameter().
bool BaseParticlePropagator::debug [private] |
The debug level.
Definition at line 144 of file BaseParticlePropagator.h.
Referenced by init(), resetDebug(), and setDebug().
bool BaseParticlePropagator::decayed [private] |
The particle decayed while propagated !
Definition at line 159 of file BaseParticlePropagator.h.
Referenced by hasDecayed(), init(), and propagate().
bool BaseParticlePropagator::fiducial [protected] |
The particle traverses some real material.
Definition at line 150 of file BaseParticlePropagator.h.
Referenced by init(), onFiducial(), and ParticlePropagator::propagateToBoundSurface().
bool BaseParticlePropagator::firstLoop [private] |
Do only the first half-loop.
Definition at line 157 of file BaseParticlePropagator.h.
Referenced by init(), propagate(), propagated(), ParticlePropagator::setPropagationConditions(), and setPropagationConditions().
int BaseParticlePropagator::propDir [private] |
The propagation direction.
Definition at line 163 of file BaseParticlePropagator.h.
Referenced by backPropagate(), init(), propagate(), propagateToHcalEntrance(), propagateToHcalExit(), propagateToHOLayer(), and propagateToVFcalEntrance().
double BaseParticlePropagator::properDecayTime [private] |
The proper decay time of the particle.
Definition at line 142 of file BaseParticlePropagator.h.
Referenced by ParticlePropagator::initProperDecayTime(), propagate(), and setProperDecayTime().
double BaseParticlePropagator::properTime [private] |
The proper time of the particle.
Definition at line 161 of file BaseParticlePropagator.h.
Referenced by init(), and propagate().
double BaseParticlePropagator::rCyl [private] |
Simulated particle that is to be resp has been propagated.
Radius of the cylinder (centred at 0,0,0) to which propagation is done
Definition at line 135 of file BaseParticlePropagator.h.
Referenced by increaseRCyl(), inside(), onBarrel(), propagate(), and setPropagationConditions().
double BaseParticlePropagator::rCyl2 [private] |
Definition at line 136 of file BaseParticlePropagator.h.
Referenced by increaseRCyl(), inside(), onBarrel(), onEndcap(), propagate(), and setPropagationConditions().
int BaseParticlePropagator::success [protected] |
0:propagation still be done, 1:reached 'barrel', 2:reached 'endcaps'
Definition at line 148 of file BaseParticlePropagator.h.
Referenced by getSuccess(), init(), propagate(), propagated(), ParticlePropagator::propagateToBoundSurface(), propagateToEcalEntrance(), propagateToHcalEntrance(), propagateToHOLayer(), propagateToPreshowerLayer1(), propagateToPreshowerLayer2(), and propagateToVFcalEntrance().
double BaseParticlePropagator::zCyl [private] |
Half-height of the cylinder (centred at 0,0,0) to which propagation is done.
Definition at line 138 of file BaseParticlePropagator.h.
Referenced by inside(), onBarrel(), onEndcap(), propagate(), and setPropagationConditions().