CMS 3D CMS Logo

List of all members | Public Member Functions | Protected Member Functions | Protected Attributes | Private Attributes
BaseParticlePropagator Class Reference

#include <BaseParticlePropagator.h>

Inheritance diagram for BaseParticlePropagator:
ParticlePropagator

Public Member Functions

bool backPropagate ()
 
 BaseParticlePropagator ()
 Default c'tor. More...
 
 BaseParticlePropagator (const RawParticle &myPart, double r, double z, double B)
 
 BaseParticlePropagator (const RawParticle &myPart, double r, double z, double B, double t)
 
double getMagneticField () const
 Get the magnetic field. More...
 
int getSuccess () const
 Has propagation been performed and was barrel or endcap reached ? More...
 
bool hasDecayed () const
 Has the particle decayed while propagated ? More...
 
double helixCentreDistToAxis () const
 The distance between the cylinder and the helix axes. More...
 
double helixCentreDistToAxis (double xC, double yC) const
 
double helixCentrePhi () const
 The azimuth if the vector joining the cylinder and the helix axes. More...
 
double helixCentrePhi (double xC, double yC) const
 
double helixCentreX () const
 The x coordinate of the helix axis. More...
 
double helixCentreX (double radius, double phi) const
 
double helixCentreY () const
 The y coordinate of the helix axis. More...
 
double helixCentreY (double radius, double phi) const
 
double helixRadius () const
 The helix Radius. More...
 
double helixRadius (double pT) const
 
double helixStartPhi () const
 The azimuth of the momentum at the vertex. More...
 
void increaseRCyl (double delta)
 Just an internal trick. More...
 
void init ()
 Initialize internal switches and quantities. More...
 
bool inside () const
 Is the vertex inside the cylinder ? (stricly inside : true) More...
 
bool inside (double rPos2) const
 
bool onBarrel () const
 Is the vertex already on the cylinder barrel ? More...
 
bool onBarrel (double rPos2) const
 
bool onEndcap () const
 Is the vertex already on the cylinder endcap ? More...
 
bool onEndcap (double rPos2) const
 
bool onFiducial () const
 Is the vertex on some material ? More...
 
bool onSurface () const
 Is the vertex already on the cylinder surface ? More...
 
bool onSurface (double rPos2) const
 
RawParticle const & particle () const
 The particle being propagated. More...
 
RawParticleparticle ()
 
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;. More...
 
void setMagneticField (double b)
 Set the magnetic field. More...
 
void setParticle (RawParticle const &iParticle)
 
void setPropagationConditions (double r, double z, bool firstLoop=true)
 Set the propagation characteristics (rCyl, zCyl and first loop only) More...
 
void setProperDecayTime (double t)
 Set the proper decay time. More...
 
double xyImpactParameter (double x0=0., double y0=0.) const
 Transverse impact parameter. More...
 
double zImpactParameter (double x0=0, double y0=0.) const
 Longitudinal impact parameter. More...
 

Protected Member Functions

double c_light () const
 The speed of light in mm/ns (!) without clhep (yeaaahhh!) More...
 

Protected Attributes

bool fiducial
 The particle traverses some real material. More...
 
int success
 0:propagation still be done, 1:reached 'barrel', 2:reached 'endcaps' More...
 

Private Attributes

double bField
 Magnetic field in the cylinder, oriented along the Z axis. More...
 
bool debug
 The debug level. More...
 
bool decayed
 The particle decayed while propagated ! More...
 
bool firstLoop
 Do only the first half-loop. More...
 
RawParticle particle_
 
int propDir
 The propagation direction. More...
 
double properDecayTime
 The proper decay time of the particle. More...
 
double properTime
 The proper time of the particle. More...
 
double rCyl
 Simulated particle that is to be resp has been propagated. More...
 
double rCyl2
 
double zCyl
 Half-height of the cylinder (centred at 0,0,0) to which propagation is done. More...
 

Detailed Description

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:

The method

propagated()

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

backPropagate()

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.

Author
Patrick Janot $Date : 19-Aug-2002, with subsequent modification for FAMOS
Version
15-Dec-2003

Definition at line 82 of file BaseParticlePropagator.h.

Constructor & Destructor Documentation

BaseParticlePropagator::BaseParticlePropagator ( )

Default c'tor.

Definition at line 5 of file BaseParticlePropagator.cc.

5  :
6  rCyl(0.), rCyl2(0.), zCyl(0.), bField(0), properDecayTime(1E99)
7 {;}
double zCyl
Half-height of the cylinder (centred at 0,0,0) to which propagation is done.
double rCyl
Simulated particle that is to be resp has been propagated.
double properDecayTime
The proper decay time of the particle.
double bField
Magnetic field in the cylinder, oriented along the Z axis.
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 16 of file BaseParticlePropagator.cc.

References init().

17  :
18  particle_(myPart), rCyl(R), rCyl2(R*R), zCyl(Z), bField(B), properDecayTime(1E99)
19 {
20  init();
21 }
double zCyl
Half-height of the cylinder (centred at 0,0,0) to which propagation is done.
double rCyl
Simulated particle that is to be resp has been propagated.
void init()
Initialize internal switches and quantities.
double properDecayTime
The proper decay time of the particle.
static const std::string B
double bField
Magnetic field in the cylinder, oriented along the Z axis.
BaseParticlePropagator::BaseParticlePropagator ( const RawParticle myPart,
double  r,
double  z,
double  B,
double  t 
)

Definition at line 9 of file BaseParticlePropagator.cc.

References init().

10  :
11  particle_(myPart), rCyl(R), rCyl2(R*R), zCyl(Z), bField(B), properDecayTime(t)
12 {
13  init();
14 }
double zCyl
Half-height of the cylinder (centred at 0,0,0) to which propagation is done.
double rCyl
Simulated particle that is to be resp has been propagated.
void init()
Initialize internal switches and quantities.
double properDecayTime
The proper decay time of the particle.
static const std::string B
double bField
Magnetic field in the cylinder, oriented along the Z axis.

Member Function Documentation

bool BaseParticlePropagator::backPropagate ( )

Update the current instance, after the back-propagation of the particle to the surface of the cylinder

Definition at line 321 of file BaseParticlePropagator.cc.

References RawParticle::charge(), RawParticle::E(), particle_, propagate(), propDir, RawParticle::Px(), RawParticle::Py(), RawParticle::Pz(), RawParticle::setCharge(), and RawParticle::setMomentum().

Referenced by propagated(), and propagateToClosestApproach().

321  {
322 
323  // Backpropagate
326  propDir = -1;
327  bool done = propagate();
330  propDir = +1;
331 
332  return done;
333 }
void setMomentum(const XYZTLorentzVector &vtx)
set the momentum
Definition: RawParticle.h:347
void setCharge(float q)
set the MEASURED charge
Definition: RawParticle.cc:67
double Py() const
y of the momentum
Definition: RawParticle.h:319
double Pz() const
z of the momentum
Definition: RawParticle.h:322
double charge() const
get the MEASURED charge
Definition: RawParticle.h:313
int propDir
The propagation direction.
double Px() const
x of the momentum
Definition: RawParticle.h:316
double E() const
energy of the momentum
Definition: RawParticle.h:325
double BaseParticlePropagator::c_light ( ) const
inlineprotected

The speed of light in mm/ns (!) without clhep (yeaaahhh!)

Definition at line 154 of file BaseParticlePropagator.h.

Referenced by helixRadius(), and propagateToBeamCylinder().

154 { return 299.792458; }
double BaseParticlePropagator::getMagneticField ( ) const
inline

Get the magnetic field.

Definition at line 309 of file BaseParticlePropagator.h.

References bField.

Referenced by TrajectoryManager::createPSimHits(), and ConvBremSeedProducer::produce().

309 { return bField; }
double bField
Magnetic field in the cylinder, oriented along the Z axis.
int BaseParticlePropagator::getSuccess ( ) const
inline

Has propagation been performed and was barrel or endcap reached ?

Definition at line 303 of file BaseParticlePropagator.h.

References success.

Referenced by PFTrackTransformer::addPoints(), PFTrackTransformer::addPointsAndBrems(), AntiElectronIDMVA6::atECalEntrance(), reco::tau::atECALEntrance(), FBaseSimEvent::fill(), ConvBremSeedProducer::GoodCluster(), ConvBremSeedProducer::produce(), TrajectoryManager::propagateToCalorimeters(), and TrajectoryManager::propagateToLayer().

303 { return success; }
int success
0:propagation still be done, 1:reached &#39;barrel&#39;, 2:reached &#39;endcaps&#39;
bool BaseParticlePropagator::hasDecayed ( ) const
inline

Has the particle decayed while propagated ?

Definition at line 300 of file BaseParticlePropagator.h.

References decayed.

Referenced by ParticlePropagator::propagateToBoundSurface(), and TrajectoryManager::propagateToCalorimeters().

300 { return decayed; }
bool decayed
The particle decayed while propagated !
double BaseParticlePropagator::helixCentreDistToAxis ( ) const
inline

The distance between the cylinder and the helix axes.

Definition at line 235 of file BaseParticlePropagator.h.

References helixCentreX(), helixCentreY(), and mathSSE::sqrt().

Referenced by propagate(), propagateToBeamCylinder(), propagateToClosestApproach(), and xyImpactParameter().

235  {
236  // The distance between the cylinder and the helix axes
237  double xC = helixCentreX();
238  double yC = helixCentreY();
239  return std::sqrt( xC*xC + yC*yC );
240  }
double helixCentreX() const
The x coordinate of the helix axis.
T sqrt(T t)
Definition: SSEVec.h:18
double helixCentreY() const
The y coordinate of the helix axis.
double BaseParticlePropagator::helixCentreDistToAxis ( double  xC,
double  yC 
) const
inline

Definition at line 242 of file BaseParticlePropagator.h.

References mathSSE::sqrt().

242  {
243  // Faster version of helixCentreDistToAxis
244  return std::sqrt( xC*xC + yC*yC );
245  }
T sqrt(T t)
Definition: SSEVec.h:18
double BaseParticlePropagator::helixCentrePhi ( ) const
inline

The azimuth if the vector joining the cylinder and the helix axes.

Definition at line 248 of file BaseParticlePropagator.h.

References helixCentreX(), and helixCentreY().

Referenced by propagate().

248  {
249  // The azimuth if the vector joining the cylinder and the helix axes
250  double xC = helixCentreX();
251  double yC = helixCentreY();
252  return xC == 0.0 && yC == 0.0 ? 0.0 : std::atan2(yC,xC);
253  }
double helixCentreX() const
The x coordinate of the helix axis.
double helixCentreY() const
The y coordinate of the helix axis.
double BaseParticlePropagator::helixCentrePhi ( double  xC,
double  yC 
) const
inline

Definition at line 255 of file BaseParticlePropagator.h.

255  {
256  // Faster version of helixCentrePhi()
257  return xC == 0.0 && yC == 0.0 ? 0.0 : std::atan2(yC,xC);
258  }
double BaseParticlePropagator::helixCentreX ( ) const
inline

The x coordinate of the helix axis.

Definition at line 213 of file BaseParticlePropagator.h.

References helixRadius(), helixStartPhi(), funct::sin(), and RawParticle::X().

Referenced by helixCentreDistToAxis(), helixCentrePhi(), propagate(), propagateToBeamCylinder(), propagateToClosestApproach(), and xyImpactParameter().

213  {
214  // The x coordinate of the helix axis
215  return particle_.X() - helixRadius() * std::sin ( helixStartPhi() );
216  }
double helixRadius() const
The helix Radius.
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
double X() const
x of vertex
Definition: RawParticle.h:305
double helixStartPhi() const
The azimuth of the momentum at the vertex.
double BaseParticlePropagator::helixCentreX ( double  radius,
double  phi 
) const
inline

Definition at line 218 of file BaseParticlePropagator.h.

References funct::sin(), and RawParticle::X().

218  {
219  // Fast version of helixCentreX()
220  return particle_.X() - radius * std::sin (phi);
221  }
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
double X() const
x of vertex
Definition: RawParticle.h:305
double BaseParticlePropagator::helixCentreY ( ) const
inline

The y coordinate of the helix axis.

Definition at line 224 of file BaseParticlePropagator.h.

References funct::cos(), helixRadius(), helixStartPhi(), and RawParticle::Y().

Referenced by helixCentreDistToAxis(), helixCentrePhi(), propagate(), propagateToBeamCylinder(), propagateToClosestApproach(), and xyImpactParameter().

224  {
225  // The y coordinate of the helix axis
226  return particle_.Y() + helixRadius() * std::cos ( helixStartPhi() );
227 }
double helixRadius() const
The helix Radius.
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
double Y() const
y of vertex
Definition: RawParticle.h:306
double helixStartPhi() const
The azimuth of the momentum at the vertex.
double BaseParticlePropagator::helixCentreY ( double  radius,
double  phi 
) const
inline

Definition at line 229 of file BaseParticlePropagator.h.

References funct::cos(), and RawParticle::Y().

229  {
230  // Fast version of helixCentreX()
231  return particle_.Y() + radius * std::cos (phi);
232  }
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
double Y() const
y of vertex
Definition: RawParticle.h:306
double BaseParticlePropagator::helixRadius ( ) const
inline

The helix Radius.

Definition at line 189 of file BaseParticlePropagator.h.

References c_light(), RawParticle::charge(), MillePedeFileConverter_cfg::e, and RawParticle::Pt().

Referenced by helixCentreX(), helixCentreY(), propagate(), propagateToClosestApproach(), propagateToNominalVertex(), and xyImpactParameter().

189  {
190  // The helix Radius
191  //
192  // The helix' Radius sign accounts for the orientation of the magnetic field
193  // (+ = along z axis) and the sign of the electric charge of the particle.
194  // It signs the rotation of the (charged) particle around the z axis:
195  // Positive means anti-clockwise, negative means clockwise rotation.
196  //
197  // The radius is returned in cm to match the units in RawParticle.
198  return particle_.charge() == 0 ? 0.0 : - particle_.Pt() / ( c_light() * 1e-5 * bField * particle_.charge() );
199  }
double c_light() const
The speed of light in mm/ns (!) without clhep (yeaaahhh!)
double Pt() const
transverse momentum
Definition: RawParticle.h:327
double charge() const
get the MEASURED charge
Definition: RawParticle.h:313
double bField
Magnetic field in the cylinder, oriented along the Z axis.
double BaseParticlePropagator::helixRadius ( double  pT) const
inline

Definition at line 201 of file BaseParticlePropagator.h.

References c_light(), RawParticle::charge(), and MillePedeFileConverter_cfg::e.

201  {
202  // a faster version of helixRadius, once Perp() has been computed
203  return particle_.charge() == 0 ? 0.0 : - pT / ( c_light() * 1e-5 * bField * particle_.charge() );
204  }
double c_light() const
The speed of light in mm/ns (!) without clhep (yeaaahhh!)
double charge() const
get the MEASURED charge
Definition: RawParticle.h:313
double bField
Magnetic field in the cylinder, oriented along the Z axis.
double BaseParticlePropagator::helixStartPhi ( ) const
inline

The azimuth of the momentum at the vertex.

Definition at line 207 of file BaseParticlePropagator.h.

References RawParticle::Px(), and RawParticle::Py().

Referenced by helixCentreX(), helixCentreY(), propagate(), propagateToClosestApproach(), and xyImpactParameter().

207  {
208  // The azimuth of the momentum at the vertex
209  return particle_.Px() == 0.0 && particle_.Py() == 0.0 ? 0.0 : std::atan2(particle_.Py(),particle_.Px());
210  }
double Py() const
y of the momentum
Definition: RawParticle.h:319
double Px() const
x of the momentum
Definition: RawParticle.h:316
void BaseParticlePropagator::increaseRCyl ( double  delta)
inline

Just an internal trick.

Definition at line 177 of file BaseParticlePropagator.h.

References delta, rCyl, and xyImpactParameter().

177 {rCyl = rCyl + delta; rCyl2 = rCyl*rCyl; }
dbl * delta
Definition: mlp_gen.cc:36
double rCyl
Simulated particle that is to be resp has been propagated.
void BaseParticlePropagator::init ( void  )

Initialize internal switches and quantities.

Definition at line 24 of file BaseParticlePropagator.cc.

References debug, decayed, fiducial, firstLoop, propDir, properTime, and success.

Referenced by BaseParticlePropagator().

24  {
25 
26  // The status of the propagation
27  success = 0;
28  // Propagate only for the first half-loop
29  firstLoop = true;
30  //
31  fiducial = true;
32  // The particle has not yet decayed
33  decayed = false;
34  // The proper time is set to zero
35  properTime = 0.;
36  // The propagation direction is along the momentum
37  propDir = 1;
38  //
39  debug = false;
40 
41 }
bool firstLoop
Do only the first half-loop.
int success
0:propagation still be done, 1:reached &#39;barrel&#39;, 2:reached &#39;endcaps&#39;
double properTime
The proper time of the particle.
int propDir
The propagation direction.
bool decayed
The particle decayed while propagated !
bool debug
The debug level.
bool fiducial
The particle traverses some real material.
bool BaseParticlePropagator::inside ( ) const
inline

Is the vertex inside the cylinder ? (stricly inside : true)

Definition at line 261 of file BaseParticlePropagator.h.

References RawParticle::R2(), and RawParticle::Z().

Referenced by propagate().

261  {
262  return (particle_.R2()<rCyl2-0.00001*rCyl && fabs(particle_.Z())<zCyl-0.00001);}
double zCyl
Half-height of the cylinder (centred at 0,0,0) to which propagation is done.
double rCyl
Simulated particle that is to be resp has been propagated.
double R2() const
vertex radius**2
Definition: RawParticle.h:310
double Z() const
z of vertex
Definition: RawParticle.h:307
bool BaseParticlePropagator::inside ( double  rPos2) const
inline

Definition at line 264 of file BaseParticlePropagator.h.

References RawParticle::Z().

264  {
265  return (rPos2<rCyl2-0.00001*rCyl && fabs(particle_.Z())<zCyl-0.00001);}
double zCyl
Half-height of the cylinder (centred at 0,0,0) to which propagation is done.
double rCyl
Simulated particle that is to be resp has been propagated.
double Z() const
z of vertex
Definition: RawParticle.h:307
bool BaseParticlePropagator::onBarrel ( ) const
inline

Is the vertex already on the cylinder barrel ?

Definition at line 278 of file BaseParticlePropagator.h.

References RawParticle::R2(), and RawParticle::Z().

Referenced by onSurface(), and propagate().

278  {
279  double rPos2 = particle_.R2();
280  return ( fabs(rPos2-rCyl2) < 0.00001*rCyl && fabs(particle_.Z()) <= zCyl );
281  }
double zCyl
Half-height of the cylinder (centred at 0,0,0) to which propagation is done.
double rCyl
Simulated particle that is to be resp has been propagated.
double R2() const
vertex radius**2
Definition: RawParticle.h:310
double Z() const
z of vertex
Definition: RawParticle.h:307
bool BaseParticlePropagator::onBarrel ( double  rPos2) const
inline

Definition at line 283 of file BaseParticlePropagator.h.

References RawParticle::Z().

283  {
284  return ( fabs(rPos2-rCyl2) < 0.00001*rCyl && fabs(particle_.Z()) <= zCyl );
285  }
double zCyl
Half-height of the cylinder (centred at 0,0,0) to which propagation is done.
double rCyl
Simulated particle that is to be resp has been propagated.
double Z() const
z of vertex
Definition: RawParticle.h:307
bool BaseParticlePropagator::onEndcap ( ) const
inline

Is the vertex already on the cylinder endcap ?

Definition at line 288 of file BaseParticlePropagator.h.

References RawParticle::R2(), rCyl2, and RawParticle::Z().

Referenced by onSurface(), and propagate().

288  {
289  return ( fabs(fabs(particle_.Z())-zCyl) < 0.00001 && particle_.R2() <= rCyl2 );
290  }
double zCyl
Half-height of the cylinder (centred at 0,0,0) to which propagation is done.
double R2() const
vertex radius**2
Definition: RawParticle.h:310
double Z() const
z of vertex
Definition: RawParticle.h:307
bool BaseParticlePropagator::onEndcap ( double  rPos2) const
inline

Definition at line 292 of file BaseParticlePropagator.h.

References rCyl2, and RawParticle::Z().

292  {
293  return ( fabs(fabs(particle_.Z())-zCyl) < 0.00001 && rPos2 <= rCyl2 );
294  }
double zCyl
Half-height of the cylinder (centred at 0,0,0) to which propagation is done.
double Z() const
z of vertex
Definition: RawParticle.h:307
bool BaseParticlePropagator::onFiducial ( ) const
inline

Is the vertex on some material ?

Definition at line 297 of file BaseParticlePropagator.h.

References fiducial.

Referenced by TrajectoryManager::propagateToLayer().

297 { return fiducial; }
bool fiducial
The particle traverses some real material.
bool BaseParticlePropagator::onSurface ( ) const
inline

Is the vertex already on the cylinder surface ?

Definition at line 269 of file BaseParticlePropagator.h.

References onBarrel(), and onEndcap().

Referenced by propagate().

269  {
270  return ( onBarrel() || onEndcap() );
271  }
bool onEndcap() const
Is the vertex already on the cylinder endcap ?
bool onBarrel() const
Is the vertex already on the cylinder barrel ?
bool BaseParticlePropagator::onSurface ( double  rPos2) const
inline

Definition at line 273 of file BaseParticlePropagator.h.

References onBarrel(), and onEndcap().

273  {
274  return ( onBarrel(rPos2) || onEndcap(rPos2) );
275  }
bool onEndcap() const
Is the vertex already on the cylinder endcap ?
bool onBarrel() const
Is the vertex already on the cylinder barrel ?
RawParticle const& BaseParticlePropagator::particle ( ) const
inline
RawParticle& BaseParticlePropagator::particle ( )
inline

Definition at line 170 of file BaseParticlePropagator.h.

References particle_.

170 { return particle_;}
bool BaseParticlePropagator::propagate ( )

Update the current instance, after the propagation of the particle to the surface of the cylinder

Definition at line 44 of file BaseParticlePropagator.cc.

References bField, RawParticle::charge(), funct::cos(), decayed, delta, hiPixelPairStep_cff::deltaPhi, RawParticle::E(), firstLoop, helixCentreDistToAxis(), helixCentrePhi(), helixCentreX(), helixCentreY(), helixRadius(), helixStartPhi(), inside(), M_PI, RawParticle::mass(), min(), onBarrel(), onEndcap(), onSurface(), particle_, RawParticle::Perp2(), propDir, properDecayTime, properTime, PVValHelper::pT, RawParticle::Pt(), RawParticle::Px(), RawParticle::Py(), RawParticle::Pz(), RawParticle::R2(), TCMET_cfi::radius, rCyl, rCyl2, RawParticle::setMomentum(), RawParticle::setVertex(), funct::sin(), mathSSE::sqrt(), success, RawParticle::T(), RawParticle::X(), RawParticle::Y(), RawParticle::Z(), and zCyl.

Referenced by PFTrackTransformer::addPoints(), PFTrackTransformer::addPointsAndBrems(), backPropagate(), ConvBremSeedProducer::produce(), propagated(), ParticlePropagator::propagateToBoundSurface(), propagateToEcal(), propagateToEcalEntrance(), propagateToHcalEntrance(), propagateToHcalExit(), propagateToHOLayer(), propagateToPreshowerLayer1(), propagateToPreshowerLayer2(), and propagateToVFcalEntrance().

44  {
45 
46 
47  //
48  // Check that the particle is not already on the cylinder surface
49  //
50  double rPos2 = particle_.R2();
51 
52  if ( onBarrel(rPos2) ) {
53  success = 1;
54  return true;
55  }
56  //
57  if ( onEndcap(rPos2) ) {
58  success = 2;
59  return true;
60  }
61  //
62  // Treat neutral particles (or no magnetic field) first
63  //
64  if ( fabs(particle_.charge()) < 1E-12 || bField < 1E-12 ) {
65 
66  //
67  // First check that the particle crosses the cylinder
68  //
69  double pT2 = particle_.Perp2();
70  // double pT2 = pT*pT;
71  // double b2 = rCyl * pT;
72  double b2 = rCyl2 * pT2;
73  double ac = fabs ( particle_.X()*particle_.Py() - particle_.Y()*particle_.Px() );
74  double ac2 = ac*ac;
75  //
76  // The particle never crosses (or never crossed) the cylinder
77  //
78  // if ( ac > b2 ) return false;
79  if ( ac2 > b2 ) return false;
80 
81  // double delta = std::sqrt(b2*b2 - ac*ac);
82  double delta = std::sqrt(b2 - ac2);
83  double tplus = -( particle_.X()*particle_.Px()+particle_.Y()*particle_.Py() ) + delta;
84  double tminus = -( particle_.X()*particle_.Px()+particle_.Y()*particle_.Py() ) - delta;
85  //
86  // Find the first (time-wise) intersection point with the cylinder
87  //
88 
89  double solution = tminus < 0 ? tplus/pT2 : tminus/pT2;
90  if ( solution < 0. ) return false;
91  //
92  // Check that the particle does (did) not exit from the endcaps first.
93  //
94  double zProp = particle_.Z() + particle_.Pz() * solution;
95  if ( fabs(zProp) > zCyl ) {
96  tplus = ( zCyl - particle_.Z() ) / particle_.Pz();
97  tminus = (-zCyl - particle_.Z() ) / particle_.Pz();
98  solution = tminus < 0 ? tplus : tminus;
99  if ( solution < 0. ) return false;
100  zProp = particle_.Z() + particle_.Pz() * solution;
101  success = 2;
102  }
103  else {
104  success = 1;
105  }
106 
107  //
108  // Check the decay time
109  //
110  double delTime = propDir * particle_.mass() * solution;
111  double factor = 1.;
112  properTime += delTime;
113  if ( properTime > properDecayTime ) {
114  factor = 1.-(properTime-properDecayTime)/delTime;
116  decayed = true;
117  }
118 
119  //
120  // Compute the coordinates of the RawParticle after propagation
121  //
122  double xProp = particle_.X() + particle_.Px() * solution * factor;
123  double yProp = particle_.Y() + particle_.Py() * solution * factor;
124  zProp = particle_.Z() + particle_.Pz() * solution * factor;
125  double tProp = particle_.T() + particle_.E() * solution * factor;
126 
127  //
128  // Last check : Although propagated to the endcaps, R could still be
129  // larger than rCyl because the cylinder was too short...
130  //
131  if ( success == 2 && xProp*xProp+yProp*yProp > rCyl2 ) {
132  success = -1;
133  return true;
134  }
135 
136  //
137  // ... and update the particle with its propagated coordinates
138  //
139  particle_.setVertex( XYZTLorentzVector(xProp,yProp,zProp,tProp) );
140 
141  return true;
142  }
143  //
144  // Treat charged particles with magnetic field next.
145  //
146  else {
147  //
148  // First check that the particle can cross the cylinder
149  //
150  double pT = particle_.Pt();
151  double radius = helixRadius(pT);
152  double phi0 = helixStartPhi();
153  double xC = helixCentreX(radius,phi0);
154  double yC = helixCentreY(radius,phi0);
155  double dist = helixCentreDistToAxis(xC,yC);
156  //
157  // The helix either is outside or includes the cylinder -> no intersection
158  //
159  if ( fabs ( fabs(radius) - dist ) > rCyl ) return false;
160  //
161  // The particle is already away from the endcaps, and will never come back
162  // Could be back-propagated, so warn the user that it is so.
163  //
164  if ( particle_.Z() * particle_.Pz() > zCyl * fabs(particle_.Pz()) ) {
165  success = -2;
166  return true;
167  }
168 
169  double pZ = particle_.Pz();
170  double phiProp, zProp;
171  //
172  // Does the helix cross the cylinder barrel ? If not, do the first half-loop
173  //
174  double rProp = std::min(fabs(radius)+dist-0.000001, rCyl);
175 
176  // Check for rounding errors in the ArcSin.
177  double sinPhiProp =
178  (rProp*rProp-radius*radius-dist*dist)/( 2.*dist*radius);
179 
180  //
181  double deltaPhi = 1E99;
182 
183  // Taylor development up to third order for large momenta
184  if ( 1.-fabs(sinPhiProp) < 1E-9 ) {
185 
186  double cphi0 = std::cos(phi0);
187  double sphi0 = std::sin(phi0);
188  double r0 = (particle_.X()*cphi0 + particle_.Y()*sphi0)/radius;
189  double q0 = (particle_.X()*sphi0 - particle_.Y()*cphi0)/radius;
190  double rcyl2 = (rCyl2 - particle_.X()*particle_.X() - particle_.Y()*particle_.Y())/(radius*radius);
191  double delta = r0*r0 + rcyl2*(1.-q0);
192 
193  // This is a solution of a second order equation, assuming phi = phi0 + epsilon
194  // and epsilon is small.
195  deltaPhi = radius > 0 ?
196  ( -r0 + std::sqrt(delta ) ) / (1.-q0) :
197  ( -r0 - std::sqrt(delta ) ) / (1.-q0);
198 
199  }
200 
201  // Use complete calculation otherwise, or if the delta phi is large anyway.
202  if ( fabs(deltaPhi) > 1E-3 ) {
203 
204  // phiProp = fabs(sinPhiProp) < 0.99999999 ?
205  // std::asin(sinPhiProp) : M_PI/2. * sinPhiProp;
206  phiProp = std::asin(sinPhiProp);
207 
208  //
209  // Compute phi after propagation (two solutions, but the asin definition
210  // (between -pi/2 and +pi/2) takes care of the wrong one.
211  //
212  phiProp = helixCentrePhi(xC,yC) +
213  ( inside(rPos2) || onSurface(rPos2) ? phiProp : M_PI-phiProp );
214 
215  //
216  // Solve the obvious two-pi ambiguities: more than one turn!
217  //
218  if ( fabs(phiProp - phi0) > 2.*M_PI )
219  phiProp += ( phiProp > phi0 ? -2.*M_PI : +2.*M_PI );
220 
221  //
222  // Check that the phi difference has the right sign, according to Q*B
223  // (Another two pi ambuiguity, slightly more subtle)
224  //
225  if ( (phiProp-phi0)*radius < 0.0 )
226  radius > 0.0 ? phiProp += 2 * M_PI : phiProp -= 2 * M_PI;
227 
228  } else {
229 
230  // Use Taylor
231  phiProp = phi0 + deltaPhi;
232 
233  }
234 
235  //
236  // Check that the particle does not exit from the endcaps first.
237  //
238  zProp = particle_.Z() + ( phiProp - phi0 ) * pZ * radius / pT ;
239  if ( fabs(zProp) > zCyl ) {
240 
241  zProp = zCyl * fabs(pZ)/pZ;
242  phiProp = phi0 + (zProp - particle_.Z()) / radius * pT / pZ;
243  success = 2;
244 
245  } else {
246 
247  //
248  // If the particle does not cross the barrel, either process the
249  // first-half loop only or propagate all the way down to the endcaps
250  //
251  if ( rProp < rCyl ) {
252 
253  if ( firstLoop || fabs(pZ)/pT < 1E-10 ) {
254 
255  success = 0;
256 
257  } else {
258 
259  zProp = zCyl * fabs(pZ)/pZ;
260  phiProp = phi0 + (zProp - particle_.Z()) / radius * pT / pZ;
261  success = 2;
262 
263  }
264  //
265  // The particle crossed the barrel
266  //
267  } else {
268 
269  success = 1;
270 
271  }
272  }
273  //
274  // Compute the coordinates of the RawParticle after propagation
275  //
276  //
277  // Check the decay time
278  //
279  double delTime = propDir * (phiProp-phi0)*radius*particle_.mass() / pT;
280  double factor = 1.;
281  properTime += delTime;
282  if ( properTime > properDecayTime ) {
283  factor = 1.-(properTime-properDecayTime)/delTime;
285  decayed = true;
286  }
287 
288  zProp = particle_.Z() + (zProp-particle_.Z())*factor;
289 
290  phiProp = phi0 + (phiProp-phi0)*factor;
291 
292  double sProp = std::sin(phiProp);
293  double cProp = std::cos(phiProp);
294  double xProp = xC + radius * sProp;
295  double yProp = yC - radius * cProp;
296  double tProp = particle_.T() + (phiProp-phi0)*radius*particle_.E()/pT;
297  double pxProp = pT * cProp;
298  double pyProp = pT * sProp;
299 
300  //
301  // Last check : Although propagated to the endcaps, R could still be
302  // larger than rCyl because the cylinder was too short...
303  //
304  if ( success == 2 && xProp*xProp+yProp*yProp > rCyl2 ) {
305  success = -1;
306  return true;
307  }
308  //
309  // ... and update the particle vertex and momentum
310  //
311  particle_.setVertex( XYZTLorentzVector(xProp,yProp,zProp,tProp) );
312  particle_.setMomentum(pxProp,pyProp,pZ,particle_.E());
313  // particle_.SetPx(pxProp);
314  // particle_.SetPy(pyProp);
315  return true;
316 
317  }
318 }
dbl * delta
Definition: mlp_gen.cc:36
double zCyl
Half-height of the cylinder (centred at 0,0,0) to which propagation is done.
void setMomentum(const XYZTLorentzVector &vtx)
set the momentum
Definition: RawParticle.h:347
double helixRadius() const
The helix Radius.
double helixCentrePhi() const
The azimuth if the vector joining the cylinder and the helix axes.
double rCyl
Simulated particle that is to be resp has been propagated.
double Perp2() const
perpendicular momentum squared
Definition: RawParticle.h:330
bool firstLoop
Do only the first half-loop.
int success
0:propagation still be done, 1:reached &#39;barrel&#39;, 2:reached &#39;endcaps&#39;
double properTime
The proper time of the particle.
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
double mass() const
get the MEASURED mass
Definition: RawParticle.h:314
double R2() const
vertex radius**2
Definition: RawParticle.h:310
double helixCentreX() const
The x coordinate of the helix axis.
double helixCentreDistToAxis() const
The distance between the cylinder and the helix axes.
double properDecayTime
The proper decay time of the particle.
T sqrt(T t)
Definition: SSEVec.h:18
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
double Y() const
y of vertex
Definition: RawParticle.h:306
double Py() const
y of the momentum
Definition: RawParticle.h:319
double Pt() const
transverse momentum
Definition: RawParticle.h:327
double Z() const
z of vertex
Definition: RawParticle.h:307
T min(T a, T b)
Definition: MathUtil.h:58
double Pz() const
z of the momentum
Definition: RawParticle.h:322
double charge() const
get the MEASURED charge
Definition: RawParticle.h:313
int propDir
The propagation direction.
#define M_PI
bool onEndcap() const
Is the vertex already on the cylinder endcap ?
bool decayed
The particle decayed while propagated !
double helixCentreY() const
The y coordinate of the helix axis.
bool inside() const
Is the vertex inside the cylinder ? (stricly inside : true)
bool onSurface() const
Is the vertex already on the cylinder surface ?
double X() const
x of vertex
Definition: RawParticle.h:305
double helixStartPhi() const
The azimuth of the momentum at the vertex.
double Px() const
x of the momentum
Definition: RawParticle.h:316
double E() const
energy of the momentum
Definition: RawParticle.h:325
double bField
Magnetic field in the cylinder, oriented along the Z axis.
double T() const
vertex time
Definition: RawParticle.h:308
bool onBarrel() const
Is the vertex already on the cylinder barrel ?
void setVertex(const XYZTLorentzVector &vtx)
set the vertex
Definition: RawParticle.h:344
math::XYZTLorentzVector XYZTLorentzVector
Definition: RawParticle.h:27
BaseParticlePropagator BaseParticlePropagator::propagated ( ) const

Return a new instance, corresponding to the particle propagated to the surface of the cylinder

Definition at line 336 of file BaseParticlePropagator.cc.

References backPropagate(), firstLoop, propagate(), and success.

Referenced by ParticlePropagator::propagated().

336  {
337  // Copy the input particle in a new instance
338  BaseParticlePropagator myPropPart(*this);
339  // Allow for many loops
340  myPropPart.firstLoop=false;
341  // Propagate the particle !
342  myPropPart.propagate();
343  // Back propagate if the forward propagation was not a success
344  if (myPropPart.success < 0 )
345  myPropPart.backPropagate();
346  // Return the new instance
347  return myPropPart;
348 }
bool BaseParticlePropagator::propagateToBeamCylinder ( const XYZTLorentzVector v,
double  radius = 0. 
)

Definition at line 631 of file BaseParticlePropagator.cc.

References a, b, bField, EnergyCorrector::c, c_light(), RawParticle::charge(), funct::cos(), gather_cfg::cout, PVValHelper::dx, PVValHelper::dy, PVValHelper::dz, MillePedeFileConverter_cfg::e, helixCentreDistToAxis(), helixCentreX(), helixCentreY(), mps_fire::i, RawParticle::momentum(), particle_, propagateToClosestApproach(), PVValHelper::pT, RawParticle::Pt(), RawParticle::Px(), RawParticle::Py(), alignCSCRings::r, diffTwoXMLs::r2, TCMET_cfi::radius, RawParticle::setCharge(), RawParticle::SetE(), RawParticle::setMomentum(), funct::sin(), mathSSE::sqrt(), RawParticle::X(), RawParticle::Y(), and RawParticle::Z().

632 {
633 
634  // For neutral (or BField = 0, simply check that the track passes through the cylinder
635  if ( particle_.charge() == 0. || bField == 0.)
636  return fabs( particle_.Px() * particle_.Y() - particle_.Py() * particle_.X() ) / particle_.Pt() < radius;
637 
638  // Now go to the charged particles
639 
640  // A few needed intermediate variables (to make the code understandable)
641  // The inner cylinder
642  double r = radius;
643  double r2 = r*r;
644  double r4 = r2*r2;
645  // The two hits
646  double dx = particle_.X()-v.X();
647  double dy = particle_.Y()-v.Y();
648  double dz = particle_.Z()-v.Z();
649  double Sxy = particle_.X()*v.X() + particle_.Y()*v.Y();
650  double Dxy = particle_.Y()*v.X() - particle_.X()*v.Y();
651  double Dxy2 = Dxy*Dxy;
652  double Sxy2 = Sxy*Sxy;
653  double SDxy = dx*dx + dy*dy;
654  double SSDxy = std::sqrt(SDxy);
655  double ATxy = std::atan2(dy,dx);
656  // The coefficients of the second order equation for the trajectory radius
657  double a = r2 - Dxy2/SDxy;
658  double b = r * (r2 - Sxy);
659  double c = r4 - 2.*Sxy*r2 + Sxy2 + Dxy2;
660 
661  // Here are the four possible solutions for
662  // 1) The trajectory radius
663  std::array<double,4> helixRs = {{0.}};
664  helixRs[0] = (b - std::sqrt(b*b - a*c))/(2.*a);
665  helixRs[1] = (b + std::sqrt(b*b - a*c))/(2.*a);
666  helixRs[2] = -helixRs[0];
667  helixRs[3] = -helixRs[1];
668  // 2) The azimuthal direction at the second point
669  std::array<double,4> helixPhis = {{0.}};
670  helixPhis[0] = std::asin ( SSDxy/(2.*helixRs[0]) );
671  helixPhis[1] = std::asin ( SSDxy/(2.*helixRs[1]) );
672  helixPhis[2] = -helixPhis[0];
673  helixPhis[3] = -helixPhis[1];
674  // Only two solutions are valid though
675  double solution1=0.;
676  double solution2=0.;
677  double phi1=0.;
678  double phi2=0.;
679  // Loop over the four possibilities
680  for ( unsigned int i=0; i<4; ++i ) {
681  helixPhis[i] = ATxy + helixPhis[i];
682  // The centre of the helix
683  double xC = helixCentreX(helixRs[i],helixPhis[i]);
684  double yC = helixCentreY(helixRs[i],helixPhis[i]);
685  // The distance of that centre to the origin
686  double dist = helixCentreDistToAxis(xC, yC);
687  /*
688  std::cout << "Solution "<< i << " = " << helixRs[i] << " accepted ? "
689  << fabs(fabs(helixRs[i]) - dist ) << " - " << radius << " = "
690  << fabs(fabs(helixRs[i]) - dist ) - fabs(radius) << " "
691  << ( fabs( fabs(fabs(helixRs[i]) - dist) -fabs(radius) ) < 1e-6) << std::endl;
692  */
693  // Check that the propagation will lead to a trajectroy tangent to the inner cylinder
694  if ( fabs( fabs(fabs(helixRs[i]) - dist) -fabs(radius) ) < 1e-6 ) {
695  // Fill first solution
696  if ( solution1 == 0. ) {
697  solution1 = helixRs[i];
698  phi1 = helixPhis[i];
699  // Fill second solution (ordered)
700  } else if ( solution2 == 0. ) {
701  if ( helixRs[i] < solution1 ) {
702  solution2 = solution1;
703  solution1 = helixRs[i];
704  phi2 = phi1;
705  phi1 = helixPhis[i];
706  } else {
707  solution2 = helixRs[i];
708  phi2 = helixPhis[i];
709  }
710  // Must not happen!
711  } else {
712  std::cout << "warning !!! More than two solutions for this track !!! " << std::endl;
713  }
714  }
715  }
716 
717  // Find the largest possible pT compatible with coming from within the inne cylinder
718  double pT = 0.;
719  double helixPhi = 0.;
720  double helixR = 0.;
721  // This should not happen
722  if ( solution1*solution2 == 0. ) {
723  std::cout << "warning !!! Less than two solution for this track! "
724  << solution1 << " " << solution2 << std::endl;
725  return false;
726  // If the two solutions have opposite sign, it means that pT = infinity is a possibility.
727  } else if ( solution1*solution2 < 0. ) {
728  pT = 1000.;
729  double norm = pT/SSDxy;
730  particle_.setCharge(+1.);
731  particle_.setMomentum(dx*norm,dy*norm,dz*norm,0.);
732  particle_.SetE(std::sqrt(particle_.momentum().Vect().Mag2()));
733  // Otherwise take the solution that gives the largest transverse momentum
734  } else {
735  if (solution1<0.) {
736  helixR = solution1;
737  helixPhi = phi1;
738  particle_.setCharge(+1.);
739  } else {
740  helixR = solution2;
741  helixPhi = phi2;
742  particle_.setCharge(-1.);
743  }
744  pT = fabs(helixR) * 1e-5 * c_light() *bField;
745  double norm = pT/SSDxy;
746  particle_.setMomentum(pT*std::cos(helixPhi),pT*std::sin(helixPhi),dz*norm,0.);
747  particle_.SetE(std::sqrt(particle_.momentum().Vect().Mag2()));
748  }
749 
750  // Propagate to closest approach to get the Z value (a bit of an overkill)
752 
753 }
void setMomentum(const XYZTLorentzVector &vtx)
set the momentum
Definition: RawParticle.h:347
void setCharge(float q)
set the MEASURED charge
Definition: RawParticle.cc:67
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
void SetE(double)
Definition: RawParticle.h:353
double c_light() const
The speed of light in mm/ns (!) without clhep (yeaaahhh!)
bool propagateToClosestApproach(double x0=0., double y0=0, bool first=true)
double helixCentreX() const
The x coordinate of the helix axis.
const XYZTLorentzVector & momentum() const
the momentum fourvector
Definition: RawParticle.h:340
double helixCentreDistToAxis() const
The distance between the cylinder and the helix axes.
T sqrt(T t)
Definition: SSEVec.h:18
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
double Y() const
y of vertex
Definition: RawParticle.h:306
double Py() const
y of the momentum
Definition: RawParticle.h:319
double Pt() const
transverse momentum
Definition: RawParticle.h:327
double Z() const
z of vertex
Definition: RawParticle.h:307
double charge() const
get the MEASURED charge
Definition: RawParticle.h:313
double helixCentreY() const
The y coordinate of the helix axis.
double b
Definition: hdecay.h:120
double X() const
x of vertex
Definition: RawParticle.h:305
double a
Definition: hdecay.h:121
double Px() const
x of the momentum
Definition: RawParticle.h:316
double bField
Magnetic field in the cylinder, oriented along the Z axis.
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)

Definition at line 359 of file BaseParticlePropagator.cc.

References backPropagate(), bField, RawParticle::charge(), helixCentreDistToAxis(), helixCentreX(), helixCentreY(), helixRadius(), helixStartPhi(), RawParticle::momentum(), particle_, PVValHelper::pT, RawParticle::Pt(), RawParticle::Px(), RawParticle::Py(), TCMET_cfi::radius, RawParticle::setMomentum(), setPropagationConditions(), RawParticle::setVertex(), mathSSE::sqrt(), RawParticle::vertex(), RawParticle::X(), RawParticle::Y(), and z.

Referenced by propagateToBeamCylinder(), ParticlePropagator::propagateToClosestApproach(), and propagateToNominalVertex().

359  {
360 
361  // Pre-computed quantities
362  double pT = particle_.Pt();
363  double radius = helixRadius(pT);
364  double phi0 = helixStartPhi();
365 
366  // The Helix centre coordinates are computed wrt to the z axis
367  // Actually the z axis might be centred on 0,0: it's the beam spot position (x0,y0)!
368  double xC = helixCentreX(radius,phi0);
369  double yC = helixCentreY(radius,phi0);
370  double distz = helixCentreDistToAxis(xC-x0,yC-y0);
371  double dist0 = helixCentreDistToAxis(xC,yC);
372 
373  //
374  // Propagation to point of clostest approach to z axis
375  //
376  double rz,r0,z;
377  if ( particle_.charge() != 0.0 && bField != 0.0 ) {
378  rz = fabs ( fabs(radius) - distz ) + std::sqrt(x0*x0+y0*y0) + 0.0000001;
379  r0 = fabs ( fabs(radius) - dist0 ) + 0.0000001;
380  } else {
381  rz = fabs( particle_.Px() * (particle_.Y()-y0) - particle_.Py() * (particle_.X()-x0) ) / particle_.Pt();
382  r0 = fabs( particle_.Px() * particle_.Y() - particle_.Py() * particle_.X() ) / particle_.Pt();
383  }
384 
385  z = 999999.;
386 
387  // Propagate to the first interesection point
388  // with cylinder of radius sqrt(x0**2+y0**2)
390  bool done = backPropagate();
391 
392  // Unsuccessful propagation - should not happen!
393  if ( !done ) return done;
394 
395  // The z axis is (0,0) - no need to go further
396  if ( fabs(rz-r0) < 1E-10 ) return done;
397  double dist1 = (particle_.X()-x0)*(particle_.X()-x0) + (particle_.Y()-y0)*(particle_.Y()-y0);
398 
399  // We are already at closest approach - no need to go further
400  if ( dist1 < 1E-10 ) return done;
401 
402  // Keep for later if it happens to be the right solution
403  XYZTLorentzVector vertex1 = particle_.vertex();
404  XYZTLorentzVector momentum1 = particle_.momentum();
405 
406  // Propagate to the distance of closest approach to (0,0)
408  done = backPropagate();
409  if ( !done ) return done;
410 
411  // Propagate to the first interesection point
412  // with cylinder of radius sqrt(x0**2+y0**2)
414  done = backPropagate();
415  if ( !done ) return done;
416  double dist2 = (particle_.X()-x0)*(particle_.X()-x0) + (particle_.Y()-y0)*(particle_.Y()-y0);
417 
418  // Keep the good solution.
419  if ( dist2 > dist1 ) {
420  particle_.setVertex(vertex1);
421  particle_.setMomentum(momentum1.X(),momentum1.Y(),momentum1.Z(),momentum1.E());
422  dist2 = dist1;
423  }
424 
425  // Done
426  return done;
427 
428 }
void setMomentum(const XYZTLorentzVector &vtx)
set the momentum
Definition: RawParticle.h:347
double helixRadius() const
The helix Radius.
void setPropagationConditions(double r, double z, bool firstLoop=true)
Set the propagation characteristics (rCyl, zCyl and first loop only)
double helixCentreX() const
The x coordinate of the helix axis.
const XYZTLorentzVector & momentum() const
the momentum fourvector
Definition: RawParticle.h:340
double helixCentreDistToAxis() const
The distance between the cylinder and the helix axes.
T sqrt(T t)
Definition: SSEVec.h:18
double Y() const
y of vertex
Definition: RawParticle.h:306
double Py() const
y of the momentum
Definition: RawParticle.h:319
double Pt() const
transverse momentum
Definition: RawParticle.h:327
double charge() const
get the MEASURED charge
Definition: RawParticle.h:313
const XYZTLorentzVector & vertex() const
the vertex fourvector
Definition: RawParticle.h:339
double helixCentreY() const
The y coordinate of the helix axis.
double X() const
x of vertex
Definition: RawParticle.h:305
double helixStartPhi() const
The azimuth of the momentum at the vertex.
double Px() const
x of the momentum
Definition: RawParticle.h:316
double bField
Magnetic field in the cylinder, oriented along the Z axis.
void setVertex(const XYZTLorentzVector &vtx)
set the vertex
Definition: RawParticle.h:344
math::XYZTLorentzVector XYZTLorentzVector
Definition: RawParticle.h:27
bool BaseParticlePropagator::propagateToEcal ( bool  first = true)

Definition at line 431 of file BaseParticlePropagator.cc.

References propagate(), and setPropagationConditions().

431  {
432  //
433  // Propagation to Ecal (including preshower) after the
434  // last Tracker layer
435  // TODO: include proper geometry
436  // Geometry taken from CMS ECAL TDR
437  //
438  // First propagate to global barrel / endcap cylinder
439  // setPropagationConditions(1290. , 3045 , first);
440  // New position of the preshower as of CMSSW_1_7_X
441  setPropagationConditions(129.0 , 303.353 , first);
442  return propagate();
443 
444 }
void setPropagationConditions(double r, double z, bool firstLoop=true)
Set the propagation characteristics (rCyl, zCyl and first loop only)
bool BaseParticlePropagator::propagateToEcalEntrance ( bool  first = true)

Definition at line 488 of file BaseParticlePropagator.cc.

References RawParticle::cos2ThetaV(), particle_, propagate(), setPropagationConditions(), and success.

Referenced by PFTrackTransformer::addPointsAndBrems(), AntiElectronIDMVA6::atECalEntrance(), reco::tau::atECALEntrance(), FBaseSimEvent::fill(), ConvBremSeedProducer::GoodCluster(), reco::tau::PFRecoTauChargedHadronFromGenericTrackPlugin< TrackClass >::operator()(), GoodSeedProducer::produce(), TrajectoryManager::propagateToCalorimeters(), and LowPtGsfElectronSeedProducer::propagateTrackToCalo().

488  {
489  //
490  // Propagation to ECAL entrance
491  // TODO: include proper geometry
492  // Geometry taken from CMS ECAL TDR
493  //
494  // First propagate to global barrel / endcap cylinder
495  setPropagationConditions(129.0 , 320.9,first);
496  bool done = propagate();
497 
498  // Go to endcap cylinder in the "barrel cut corner"
499  // eta = 1.479 -> cos^2(theta) = 0.81230
500  // if ( done && eta > 1.479 && success == 1 ) {
501  if ( done && particle_.cos2ThetaV() > 0.81230 && success == 1 ) {
502  setPropagationConditions(152.6 , 320.9, first);
503  done = propagate();
504  }
505 
506  // We are not in the ECAL acceptance
507  // eta = 3.0 -> cos^2(theta) = 0.99013
508  if ( particle_.cos2ThetaV() > 0.99014 ) success = 0;
509 
510  return done;
511 }
int success
0:propagation still be done, 1:reached &#39;barrel&#39;, 2:reached &#39;endcaps&#39;
void setPropagationConditions(double r, double z, bool firstLoop=true)
Set the propagation characteristics (rCyl, zCyl and first loop only)
double cos2ThetaV() const
Definition: RawParticle.h:300
bool BaseParticlePropagator::propagateToHcalEntrance ( bool  first = true)

Definition at line 514 of file BaseParticlePropagator.cc.

References RawParticle::cos2ThetaV(), particle_, propagate(), propDir, setPropagationConditions(), and success.

Referenced by PFTrackTransformer::addPointsAndBrems(), FBaseSimEvent::fill(), and TrajectoryManager::propagateToCalorimeters().

514  {
515  //
516  // Propagation to HCAL entrance
517  // TODO: include proper geometry
518  // Geometry taken from CMSSW_3_1_X xml (Sunanda)
519  //
520 
521  // First propagate to global barrel / endcap cylinder
522  setPropagationConditions(177.5 , 335.0, first);
523  propDir = 0;
524  bool done = propagate();
525  propDir = 1;
526 
527  // If went through the bottom of HB cylinder -> re-propagate to HE surface
528  if (done && success == 2) {
529  setPropagationConditions(300.0, 400.458, first);
530  propDir = 0;
531  done = propagate();
532  propDir = 1;
533  }
534 
535 
536  // out of the HB/HE acceptance
537  // eta = 3.0 -> cos^2(theta) = 0.99014
538  if ( done && particle_.cos2ThetaV() > 0.99014 ) success = 0;
539 
540  return done;
541 }
int success
0:propagation still be done, 1:reached &#39;barrel&#39;, 2:reached &#39;endcaps&#39;
void setPropagationConditions(double r, double z, bool firstLoop=true)
Set the propagation characteristics (rCyl, zCyl and first loop only)
double cos2ThetaV() const
Definition: RawParticle.h:300
int propDir
The propagation direction.
bool BaseParticlePropagator::propagateToHcalExit ( bool  first = true)

Definition at line 567 of file BaseParticlePropagator.cc.

References propagate(), propDir, and setPropagationConditions().

Referenced by PFTrackTransformer::addPointsAndBrems(), and FBaseSimEvent::fill().

567  {
568  //
569  // Propagation to HCAL exit
570  // TODO: include proper geometry
571  // Geometry taken from CMSSW_3_1_X xml (Sunanda)
572  //
573 
574  // Approximate it to a single cylinder as it is not that crucial.
575  setPropagationConditions(263.9 , 554.1, first);
576  // this->rawPart().particle_.setCharge(0.0); ?? Shower Propagation ??
577  propDir = 0;
578  bool done = propagate();
579  propDir = 1;
580 
581  return done;
582 }
void setPropagationConditions(double r, double z, bool firstLoop=true)
Set the propagation characteristics (rCyl, zCyl and first loop only)
int propDir
The propagation direction.
bool BaseParticlePropagator::propagateToHOLayer ( bool  first = true)

Definition at line 585 of file BaseParticlePropagator.cc.

References funct::abs(), particle_, propagate(), propDir, setPropagationConditions(), success, and RawParticle::Z().

Referenced by PFTrackTransformer::addPointsAndBrems(), and FBaseSimEvent::fill().

585  {
586  //
587  // Propagation to Layer0 of HO (Ring-0)
588  // TODO: include proper geometry
589  // Geometry taken from CMSSW_3_1_X xml (Sunanda)
590  //
591 
592  // Approximate it to a single cylinder as it is not that crucial.
593  setPropagationConditions(387.6, 1000.0, first); //Dia is the middle of HO \sqrt{384.8**2+46.2**2}
594  // setPropagationConditions(410.2, 1000.0, first); //sqrt{406.6**2+54.2**2} //for outer layer
595  // this->rawPart().setCharge(0.0); ?? Shower Propagation ??
596  propDir = 0;
597  bool done = propagate();
598  propDir = 1;
599 
600  if ( done && std::abs(particle_.Z()) > 700.25) success = 0;
601 
602  return done;
603 }
int success
0:propagation still be done, 1:reached &#39;barrel&#39;, 2:reached &#39;endcaps&#39;
void setPropagationConditions(double r, double z, bool firstLoop=true)
Set the propagation characteristics (rCyl, zCyl and first loop only)
double Z() const
z of vertex
Definition: RawParticle.h:307
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
int propDir
The propagation direction.
bool BaseParticlePropagator::propagateToNominalVertex ( const XYZTLorentzVector hit2 = XYZTLorentzVector(0.,0.,0.,0.))

Definition at line 607 of file BaseParticlePropagator.cc.

References bField, RawParticle::charge(), funct::cos(), PVValHelper::dx, PVValHelper::dy, helixRadius(), particle_, phi, propagateToClosestApproach(), PVValHelper::pT, RawParticle::pt(), RawParticle::SetPx(), RawParticle::SetPy(), funct::sin(), mathSSE::sqrt(), RawParticle::X(), and RawParticle::Y().

Referenced by ParticlePropagator::propagateToNominalVertex().

608 {
609 
610 
611  // Not implemented for neutrals (used for electrons only)
612  if ( particle_.charge() == 0. || bField == 0.) return false;
613 
614  // Define the proper pT direction to meet the point (vx,vy) in (x,y)
615  double dx = particle_.X()-v.X();
616  double dy = particle_.Y()-v.Y();
617  double phi = std::atan2(dy,dx) + std::asin ( std::sqrt(dx*dx+dy*dy)/(2.*helixRadius()) );
618 
619  // The absolute pT (kept to the original value)
620  double pT = particle_.pt();
621 
622  // Set the new pT
623  particle_.SetPx(pT*std::cos(phi));
624  particle_.SetPy(pT*std::sin(phi));
625 
626  return propagateToClosestApproach(v.X(),v.Y());
627 
628 }
double helixRadius() const
The helix Radius.
double pt() const
transverse momentum
Definition: RawParticle.h:328
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
bool propagateToClosestApproach(double x0=0., double y0=0, bool first=true)
T sqrt(T t)
Definition: SSEVec.h:18
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
double Y() const
y of vertex
Definition: RawParticle.h:306
double charge() const
get the MEASURED charge
Definition: RawParticle.h:313
void SetPx(double)
Definition: RawParticle.h:350
void SetPy(double)
Definition: RawParticle.h:351
double X() const
x of vertex
Definition: RawParticle.h:305
double bField
Magnetic field in the cylinder, oriented along the Z axis.
bool BaseParticlePropagator::propagateToPreshowerLayer1 ( bool  first = true)

Definition at line 447 of file BaseParticlePropagator.cc.

References particle_, propagate(), RawParticle::R2(), setPropagationConditions(), and success.

Referenced by PFTrackTransformer::addPointsAndBrems(), FBaseSimEvent::fill(), and TrajectoryManager::propagateToCalorimeters().

447  {
448  //
449  // Propagation to Preshower Layer 1
450  // TODO: include proper geometry
451  // Geometry taken from CMS ECAL TDR
452  //
453  // First propagate to global barrel / endcap cylinder
454  // setPropagationConditions(1290., 3045 , first);
455  // New position of the preshower as of CMSSW_1_7_X
456  setPropagationConditions(129.0, 303.353 , first);
457  bool done = propagate();
458 
459  // Check that were are on the Layer 1
460  if ( done && (particle_.R2() > 125.0*125.0 || particle_.R2() < 45.0*45.0) )
461  success = 0;
462 
463  return done;
464 }
int success
0:propagation still be done, 1:reached &#39;barrel&#39;, 2:reached &#39;endcaps&#39;
void setPropagationConditions(double r, double z, bool firstLoop=true)
Set the propagation characteristics (rCyl, zCyl and first loop only)
double R2() const
vertex radius**2
Definition: RawParticle.h:310
bool BaseParticlePropagator::propagateToPreshowerLayer2 ( bool  first = true)

Definition at line 467 of file BaseParticlePropagator.cc.

References particle_, propagate(), RawParticle::R2(), setPropagationConditions(), and success.

Referenced by PFTrackTransformer::addPointsAndBrems(), FBaseSimEvent::fill(), and TrajectoryManager::propagateToCalorimeters().

467  {
468  //
469  // Propagation to Preshower Layer 2
470  // TODO: include proper geometry
471  // Geometry taken from CMS ECAL TDR
472  //
473  // First propagate to global barrel / endcap cylinder
474  // setPropagationConditions(1290. , 3090 , first);
475  // New position of the preshower as of CMSSW_1_7_X
476  setPropagationConditions(129.0 , 307.838 , first);
477  bool done = propagate();
478 
479  // Check that we are on Layer 2
480  if ( done && (particle_.R2() > 125.0*125.0 || particle_.R2() < 45.0*45.0 ) )
481  success = 0;
482 
483  return done;
484 
485 }
int success
0:propagation still be done, 1:reached &#39;barrel&#39;, 2:reached &#39;endcaps&#39;
void setPropagationConditions(double r, double z, bool firstLoop=true)
Set the propagation characteristics (rCyl, zCyl and first loop only)
double R2() const
vertex radius**2
Definition: RawParticle.h:310
bool BaseParticlePropagator::propagateToVFcalEntrance ( bool  first = true)

Definition at line 544 of file BaseParticlePropagator.cc.

References RawParticle::cos2ThetaV(), particle_, propagate(), propDir, setPropagationConditions(), and success.

Referenced by FBaseSimEvent::fill(), and TrajectoryManager::propagateToCalorimeters().

544  {
545  //
546  // Propagation to VFCAL entrance
547  // TODO: include proper geometry
548  // Geometry taken from DAQ TDR Chapter 13
549 
550  setPropagationConditions(400.0 , 1110.0, first);
551  propDir = 0;
552  bool done = propagate();
553  propDir = 1;
554 
555  if (!done) success = 0;
556 
557  // We are not in the VFCAL acceptance
558  // eta = 3.0 -> cos^2(theta) = 0.99014
559  // eta = 5.2 -> cos^2(theta) = 0.9998755
560  double c2teta = particle_.cos2ThetaV();
561  if ( done && ( c2teta < 0.99014 || c2teta > 0.9998755 ) ) success = 0;
562 
563  return done;
564 }
int success
0:propagation still be done, 1:reached &#39;barrel&#39;, 2:reached &#39;endcaps&#39;
void setPropagationConditions(double r, double z, bool firstLoop=true)
Set the propagation characteristics (rCyl, zCyl and first loop only)
double cos2ThetaV() const
Definition: RawParticle.h:300
int propDir
The propagation direction.
void BaseParticlePropagator::resetDebug ( )
inline

Definition at line 313 of file BaseParticlePropagator.h.

313 { debug = false; }
bool debug
The debug level.
void BaseParticlePropagator::setDebug ( )
inline

Set the debug leve;.

Definition at line 312 of file BaseParticlePropagator.h.

312 { debug = true; }
bool debug
The debug level.
void BaseParticlePropagator::setMagneticField ( double  b)
inline
void BaseParticlePropagator::setParticle ( RawParticle const &  iParticle)
inline

Definition at line 171 of file BaseParticlePropagator.h.

171 { particle_=iParticle;}
void BaseParticlePropagator::setPropagationConditions ( double  r,
double  z,
bool  firstLoop = true 
)

Set the propagation characteristics (rCyl, zCyl and first loop only)

Definition at line 351 of file BaseParticlePropagator.cc.

References f, firstLoop, dttmaxenums::R, rCyl, rCyl2, DOFs::Z, and zCyl.

Referenced by PFTrackTransformer::addPoints(), PFTrackTransformer::addPointsAndBrems(), propagateToClosestApproach(), propagateToEcal(), propagateToEcalEntrance(), propagateToHcalEntrance(), propagateToHcalExit(), propagateToHOLayer(), propagateToPreshowerLayer1(), propagateToPreshowerLayer2(), propagateToVFcalEntrance(), and ParticlePropagator::setPropagationConditions().

351  {
352  rCyl = R;
353  rCyl2 = R*R;
354  zCyl = Z;
355  firstLoop = f;
356 }
double zCyl
Half-height of the cylinder (centred at 0,0,0) to which propagation is done.
double rCyl
Simulated particle that is to be resp has been propagated.
bool firstLoop
Do only the first half-loop.
double f[11][100]
void BaseParticlePropagator::setProperDecayTime ( double  t)
inline

Set the proper decay time.

Definition at line 174 of file BaseParticlePropagator.h.

References lumiQTWidget::t.

Referenced by ParticlePropagator::initProperDecayTime(), and ParticlePropagator::ParticlePropagator().

174 { properDecayTime = t; }
double properDecayTime
The proper decay time of the particle.
double BaseParticlePropagator::xyImpactParameter ( double  x0 = 0.,
double  y0 = 0. 
) const

Transverse impact parameter.

Definition at line 756 of file BaseParticlePropagator.cc.

References bField, RawParticle::charge(), helixCentreDistToAxis(), helixCentreX(), helixCentreY(), helixRadius(), helixStartPhi(), particle_, PVValHelper::pT, RawParticle::Pt(), RawParticle::Px(), RawParticle::Py(), TCMET_cfi::radius, RawParticle::X(), and RawParticle::Y().

Referenced by increaseRCyl().

756  {
757 
758  double ip=0.;
759  double pT = particle_.Pt();
760 
761  if ( particle_.charge() != 0.0 && bField != 0.0 ) {
762  double radius = helixRadius(pT);
763  double phi0 = helixStartPhi();
764 
765  // The Helix centre coordinates are computed wrt to the z axis
766  double xC = helixCentreX(radius,phi0);
767  double yC = helixCentreY(radius,phi0);
768  double distz = helixCentreDistToAxis(xC-x0,yC-y0);
769  ip = distz - fabs(radius);
770  } else {
771  ip = fabs( particle_.Px() * (particle_.Y()-y0) - particle_.Py() * (particle_.X()-x0) ) / pT;
772  }
773 
774  return ip;
775 
776 }
double helixRadius() const
The helix Radius.
double helixCentreX() const
The x coordinate of the helix axis.
double helixCentreDistToAxis() const
The distance between the cylinder and the helix axes.
double Y() const
y of vertex
Definition: RawParticle.h:306
double Py() const
y of the momentum
Definition: RawParticle.h:319
double Pt() const
transverse momentum
Definition: RawParticle.h:327
double charge() const
get the MEASURED charge
Definition: RawParticle.h:313
double helixCentreY() const
The y coordinate of the helix axis.
double X() const
x of vertex
Definition: RawParticle.h:305
double helixStartPhi() const
The azimuth of the momentum at the vertex.
double Px() const
x of the momentum
Definition: RawParticle.h:316
double bField
Magnetic field in the cylinder, oriented along the Z axis.
double BaseParticlePropagator::zImpactParameter ( double  x0 = 0,
double  y0 = 0. 
) const
inline

Longitudinal impact parameter.

Definition at line 183 of file BaseParticlePropagator.h.

References RawParticle::Perp2(), RawParticle::Pz(), mathSSE::sqrt(), RawParticle::X(), RawParticle::Y(), and RawParticle::Z().

183  {
184  // Longitudinal impact parameter
185  return particle_.Z() - particle_.Pz() * std::sqrt( ((particle_.X()-x0)*(particle_.X()-x0) + (particle_.Y()-y0)*(particle_.Y()-y0) ) / particle_.Perp2());
186  }
double Perp2() const
perpendicular momentum squared
Definition: RawParticle.h:330
T sqrt(T t)
Definition: SSEVec.h:18
double Y() const
y of vertex
Definition: RawParticle.h:306
double Z() const
z of vertex
Definition: RawParticle.h:307
double Pz() const
z of the momentum
Definition: RawParticle.h:322
double X() const
x of vertex
Definition: RawParticle.h:305

Member Data Documentation

double BaseParticlePropagator::bField
private

Magnetic field in the cylinder, oriented along the Z axis.

Definition at line 141 of file BaseParticlePropagator.h.

Referenced by getMagneticField(), propagate(), propagateToBeamCylinder(), propagateToClosestApproach(), propagateToNominalVertex(), and xyImpactParameter().

bool BaseParticlePropagator::debug
private
bool BaseParticlePropagator::decayed
private

The particle decayed while propagated !

Definition at line 160 of file BaseParticlePropagator.h.

Referenced by hasDecayed(), init(), and propagate().

bool BaseParticlePropagator::fiducial
protected

The particle traverses some real material.

Definition at line 151 of file BaseParticlePropagator.h.

Referenced by init(), onFiducial(), and ParticlePropagator::propagateToBoundSurface().

bool BaseParticlePropagator::firstLoop
private

Do only the first half-loop.

Definition at line 158 of file BaseParticlePropagator.h.

Referenced by init(), propagate(), propagated(), ParticlePropagator::setPropagationConditions(), and setPropagationConditions().

RawParticle BaseParticlePropagator::particle_
private
int BaseParticlePropagator::propDir
private
double BaseParticlePropagator::properDecayTime
private

The proper decay time of the particle.

Definition at line 143 of file BaseParticlePropagator.h.

Referenced by ParticlePropagator::initProperDecayTime(), and propagate().

double BaseParticlePropagator::properTime
private

The proper time of the particle.

Definition at line 162 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 136 of file BaseParticlePropagator.h.

Referenced by increaseRCyl(), propagate(), and setPropagationConditions().

double BaseParticlePropagator::rCyl2
private

Definition at line 137 of file BaseParticlePropagator.h.

Referenced by onEndcap(), propagate(), and setPropagationConditions().

int BaseParticlePropagator::success
protected
double BaseParticlePropagator::zCyl
private

Half-height of the cylinder (centred at 0,0,0) to which propagation is done.

Definition at line 139 of file BaseParticlePropagator.h.

Referenced by propagate(), and setPropagationConditions().