CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | Private Types | Private Member Functions | Private Attributes
AnalyticalPropagator Class Reference

#include <AnalyticalPropagator.h>

Inheritance diagram for AnalyticalPropagator:
Propagator

Public Member Functions

 AnalyticalPropagator (const MagneticField *field, PropagationDirection dir=alongMomentum, float maxDPhi=1.6, bool isOld=true)
 
virtual AnalyticalPropagatorclone () const
 
TrajectoryStateOnSurface propagate (const FreeTrajectoryState &fts, const Plane &plane) const
 propagation to plane More...
 
TrajectoryStateOnSurface propagate (const FreeTrajectoryState &fts, const Cylinder &cylinder) const
 propagation to cylinder More...
 
std::pair
< TrajectoryStateOnSurface,
double > 
propagateWithPath (const FreeTrajectoryState &fts, const Plane &plane) const
 propagation to plane with path length More...
 
std::pair
< TrajectoryStateOnSurface,
double > 
propagateWithPath (const FreeTrajectoryState &fts, const Cylinder &cylinder) const
 propagation to cylinder with path length More...
 
virtual bool setMaxDirectionChange (float phiMax)
 
void setMaxRelativeChangeInBz (const float maxDBz)
 
 ~AnalyticalPropagator ()
 
- Public Member Functions inherited from Propagator
virtual TrajectoryStateOnSurface propagate (const FreeTrajectoryState &, const Surface &) const
 
virtual TrajectoryStateOnSurface propagate (const TrajectoryStateOnSurface &, const Surface &) const
 
virtual TrajectoryStateOnSurface propagate (const TrajectoryStateOnSurface &, const Plane &) const
 
virtual TrajectoryStateOnSurface propagate (const TrajectoryStateOnSurface &, const Cylinder &) const
 
virtual FreeTrajectoryState propagate (const FreeTrajectoryState &, const reco::BeamSpot &) const
 
virtual std::pair
< TrajectoryStateOnSurface,
double > 
propagateWithPath (const FreeTrajectoryState &, const Surface &) const
 
virtual std::pair
< TrajectoryStateOnSurface,
double > 
propagateWithPath (const TrajectoryStateOnSurface &, const Surface &) const
 
virtual std::pair
< TrajectoryStateOnSurface,
double > 
propagateWithPath (const TrajectoryStateOnSurface &, const Plane &) const
 
virtual std::pair
< TrajectoryStateOnSurface,
double > 
propagateWithPath (const TrajectoryStateOnSurface &, const Cylinder &) const
 
virtual std::pair
< FreeTrajectoryState, double > 
propagateWithPath (const FreeTrajectoryState &, const GlobalPoint &, const GlobalPoint &) const
 
virtual PropagationDirection propagationDirection () const
 
 Propagator (PropagationDirection dir=alongMomentum)
 
virtual void setPropagationDirection (PropagationDirection dir) const
 
virtual ~Propagator ()
 

Private Types

typedef std::pair
< TrajectoryStateOnSurface,
double > 
TsosWP
 

Private Member Functions

virtual const MagneticFieldmagneticField () const
 
std::pair
< TrajectoryStateOnSurface,
double > 
propagatedStateWithPath (const FreeTrajectoryState &fts, const Surface &surface, const GlobalTrajectoryParameters &gtp, const double &s) const dso_internal
 propagation of errors (if needed) and generation of a new TSOS More...
 
bool propagateParametersOnCylinder (const FreeTrajectoryState &fts, const Cylinder &cylinder, GlobalPoint &x, GlobalVector &p, double &s) const dso_internal
 parameter propagation to cylinder (returns position, momentum and path length) More...
 
bool propagateParametersOnPlane (const FreeTrajectoryState &fts, const Plane &plane, GlobalPoint &x, GlobalVector &p, double &s) const dso_internal
 parameter propagation to plane (returns position, momentum and path length) More...
 
bool propagateWithHelixCrossing (HelixPlaneCrossing &, const Plane &, const float, GlobalPoint &, GlobalVector &, double &s) const dso_internal
 helix parameter propagation to a plane using HelixPlaneCrossing More...
 
bool propagateWithLineCrossing (const GlobalPoint &, const GlobalVector &, const Plane &, GlobalPoint &, double &) const dso_internal
 straight line parameter propagation to a plane More...
 
bool propagateWithLineCrossing (const GlobalPoint &, const GlobalVector &, const Cylinder &, GlobalPoint &, double &) const dso_internal
 straight line parameter propagation to a cylinder More...
 

Private Attributes

bool isOldPropagationType
 
const MagneticFieldtheField
 
float theMaxDBzRatio
 
float theMaxDPhi2
 

Detailed Description

(Mostly) analytical helix propagation to cylindrical or planar surfaces. Based on GtfGeometricalPropagator with successive replacement of components (currently: propagation to arbitrary plane).

Definition at line 23 of file AnalyticalPropagator.h.

Member Typedef Documentation

typedef std::pair<TrajectoryStateOnSurface,double> AnalyticalPropagator::TsosWP
private

Definition at line 127 of file AnalyticalPropagator.h.

Constructor & Destructor Documentation

AnalyticalPropagator::AnalyticalPropagator ( const MagneticField field,
PropagationDirection  dir = alongMomentum,
float  maxDPhi = 1.6,
bool  isOld = true 
)
inline

Definition at line 27 of file AnalyticalPropagator.h.

Referenced by clone().

29  :
30  Propagator(dir),
31  theMaxDPhi2(maxDPhi*maxDPhi),
32  theMaxDBzRatio(0.5),
33  theField(field),
34  isOldPropagationType(isOld) {}
Propagator(PropagationDirection dir=alongMomentum)
Definition: Propagator.h:41
const MagneticField * theField
dbl *** dir
Definition: mlp_gen.cc:35
AnalyticalPropagator::~AnalyticalPropagator ( )
inline

Definition at line 36 of file AnalyticalPropagator.h.

36 {}

Member Function Documentation

virtual AnalyticalPropagator* AnalyticalPropagator::clone ( void  ) const
inlinevirtual

Implements Propagator.

Definition at line 75 of file AnalyticalPropagator.h.

References AnalyticalPropagator().

79  {
80  return new AnalyticalPropagator(*this);
81  }
AnalyticalPropagator(const MagneticField *field, PropagationDirection dir=alongMomentum, float maxDPhi=1.6, bool isOld=true)
virtual const MagneticField* AnalyticalPropagator::magneticField ( ) const
inlineprivatevirtual

Implements Propagator.

Definition at line 124 of file AnalyticalPropagator.h.

References theField.

124 {return theField;}
const MagneticField * theField
TrajectoryStateOnSurface AnalyticalPropagator::propagate ( const FreeTrajectoryState fts,
const Plane plane 
) const
inlinevirtual
TrajectoryStateOnSurface AnalyticalPropagator::propagate ( const FreeTrajectoryState fts,
const Cylinder cylinder 
) const
inlinevirtual

propagation to cylinder

Implements Propagator.

Definition at line 58 of file AnalyticalPropagator.h.

References propagateWithPath().

59  {
60  return propagateWithPath(fts,cylinder).first;
61  }
std::pair< TrajectoryStateOnSurface, double > propagateWithPath(const FreeTrajectoryState &fts, const Plane &plane) const
propagation to plane with path length
std::pair< TrajectoryStateOnSurface, double > AnalyticalPropagator::propagatedStateWithPath ( const FreeTrajectoryState fts,
const Surface surface,
const GlobalTrajectoryParameters gtp,
const double &  s 
) const
private

propagation of errors (if needed) and generation of a new TSOS

Definition at line 115 of file AnalyticalPropagator.cc.

References SurfaceSideDefinition::afterSurface, alongMomentum, SurfaceSideDefinition::beforeSurface, FreeTrajectoryState::curvilinearError(), FreeTrajectoryState::hasError(), AnalyticalCurvilinearJacobian::jacobian(), CurvilinearTrajectoryError::matrix(), GlobalTrajectoryParameters::momentum(), FreeTrajectoryState::parameters(), GlobalTrajectoryParameters::position(), and alignCSCRings::s.

119 {
120  //
121  // for forward propagation: state is before surface,
122  // for backward propagation: state is after surface
123  //
126  //
127  //
128  // error propagation (if needed) and conversion to a TrajectoryStateOnSurface
129  //
130  if (fts.hasError()) {
131  //
132  // compute jacobian
133  //
134  AnalyticalCurvilinearJacobian analyticalJacobian(fts.parameters(), gtp.position(), gtp.momentum(), s);
135  const AlgebraicMatrix55 &jacobian = analyticalJacobian.jacobian();
136  CurvilinearTrajectoryError cte(ROOT::Math::Similarity(jacobian, fts.curvilinearError().matrix()));
137  return TsosWP(TrajectoryStateOnSurface(gtp,cte,surface,side),s);
138  }
139  else {
140  //
141  // return state without errors
142  //
143  return TsosWP(TrajectoryStateOnSurface(gtp,surface,side),s);
144  }
145 }
virtual PropagationDirection propagationDirection() const
Definition: Propagator.h:143
const GlobalTrajectoryParameters & parameters() const
const AlgebraicMatrix55 & jacobian() const
const CurvilinearTrajectoryError & curvilinearError() const
const AlgebraicSymMatrix55 & matrix() const
std::pair< TrajectoryStateOnSurface, double > TsosWP
ROOT::Math::SMatrix< double, 5, 5, ROOT::Math::MatRepStd< double, 5, 5 > > AlgebraicMatrix55
bool AnalyticalPropagator::propagateParametersOnCylinder ( const FreeTrajectoryState fts,
const Cylinder cylinder,
GlobalPoint x,
GlobalVector p,
double &  s 
) const
private

parameter propagation to cylinder (returns position, momentum and path length)

Definition at line 147 of file AnalyticalPropagator.cc.

References alignCSCRings::e, PV3DBase< T, PVType, FrameType >::mag(), FreeTrajectoryState::momentum(), AlCaHLTBitMon_ParallelJobs::p, HelixBarrelCylinderCrossing::pathLength(), GloballyPositioned< T >::position(), FreeTrajectoryState::position(), rho, alignCSCRings::s, FreeTrajectoryState::transverseCurvature(), Vector3DBase< T, FrameTag >::unit(), x, PV3DBase< T, PVType, FrameType >::x(), and PV3DBase< T, PVType, FrameType >::y().

150 {
151 
152  GlobalPoint const & sp = cylinder.position();
153  if (sp.x()!=0. || sp.y()!=0.) {
154  throw PropagationException("Cannot propagate to an arbitrary cylinder");
155  }
156  // preset output
157  x = fts.position();
158  p = fts.momentum();
159  s = 0;
160  // (transverse) curvature
161  double rho = fts.transverseCurvature();
162  //
163  // Straight line approximation? |rho|<1.e-10 equivalent to ~ 1um
164  // difference in transversal position at 10m.
165  //
166  if( fabs(rho)<1.e-10 )
167  return propagateWithLineCrossing(fts.position(),p,cylinder,x,s);
168  //
169  // Helix case
170  //
171  // check for possible intersection
172  const double tolerance = 1.e-4; // 1 micron distance
173  double rdiff = x.perp() - cylinder.radius();
174  if ( fabs(rdiff) < tolerance ) return true;
175  //
176  // Instantiate HelixBarrelCylinderCrossing and get solutions
177  //
178  HelixBarrelCylinderCrossing cylinderCrossing(fts.position(),fts.momentum(),rho,
179  propagationDirection(),cylinder);
180  if ( !cylinderCrossing.hasSolution() ) return false;
181  // path length
182  s = cylinderCrossing.pathLength();
183  // point
184  x = cylinderCrossing.position();
185  // direction (renormalised)
186  p = cylinderCrossing.direction().unit()*fts.momentum().mag();
187  return true;
188 }
virtual PropagationDirection propagationDirection() const
Definition: Propagator.h:143
Definition: DDAxes.h:10
T y() const
Definition: PV3DBase.h:62
T mag() const
Definition: PV3DBase.h:66
Common base class.
GlobalVector momentum() const
Vector3DBase unit() const
Definition: Vector3DBase.h:57
GlobalPoint position() const
double transverseCurvature() const
Definition: DDAxes.h:10
T x() const
Definition: PV3DBase.h:61
const PositionType & position() const
bool propagateWithLineCrossing(const GlobalPoint &, const GlobalVector &, const Plane &, GlobalPoint &, double &) const dso_internal
straight line parameter propagation to a plane
bool AnalyticalPropagator::propagateParametersOnPlane ( const FreeTrajectoryState fts,
const Plane plane,
GlobalPoint x,
GlobalVector p,
double &  s 
) const
private

parameter propagation to plane (returns position, momentum and path length)

Definition at line 191 of file AnalyticalPropagator.cc.

References HelixArbitraryPlaneCrossing::direction(), alignCSCRings::e, LogDebug, PV3DBase< T, PVType, FrameType >::mag(), Basic3DVector< T >::mag(), FreeTrajectoryState::momentum(), oppositeToMomentum, AlCaHLTBitMon_ParallelJobs::p, perp(), PV3DBase< T, PVType, FrameType >::phi(), HelixArbitraryPlaneCrossing::position(), FreeTrajectoryState::position(), rho, alignCSCRings::s, indexGen::s2, FreeTrajectoryState::transverseCurvature(), and x.

196 {
197  // initialisation of position, momentum and path length
198  x = fts.position();
199  p = fts.momentum();
200  s = 0;
201  // (transverse) curvature
202  double rho = fts.transverseCurvature();
203  //
204  // Straight line approximation? |rho|<1.e-10 equivalent to ~ 1um
205  // difference in transversal position at 10m.
206  //
207  if( fabs(rho)<1.e-10 )
208  return propagateWithLineCrossing(fts.position(),p,plane,x,s);
209  //
210  // Helix case
211  //
212  GlobalVector u = plane.normalVector();
213  const double small = 1.e-6; // for orientation of planes
214  //
215  // Frame-independant point and vector are created explicitely to
216  // avoid confusing gcc (refuses to compile with temporary objects
217  // in the constructor).
218  //
221  if (isOldPropagationType && fabs(u.z()) < small) {
222  // barrel plane:
223  // instantiate HelixBarrelPlaneCrossing, get vector of solutions and check for existance
224  HelixBarrelPlaneCrossingByCircle planeCrossing(helixPos,helixDir,rho,propagationDirection());
225  return propagateWithHelixCrossing(planeCrossing,plane,fts.momentum().mag(),x,p,s);
226  }
227  if (isOldPropagationType && fabs(u.x()) < small && fabs(u.y()) < small) {
228  // forward plane:
229  // instantiate HelixForwardPlaneCrossing, get vector of solutions and check for existance
230  HelixForwardPlaneCrossing planeCrossing(helixPos,helixDir,rho,propagationDirection());
231  return propagateWithHelixCrossing(planeCrossing,plane,fts.momentum().mag(),x,p,s);
232  }
233  else {
234  // arbitrary plane:
235  // instantiate HelixArbitraryPlaneCrossing, get vector of solutions and check for existance
237  HelixArbitraryPlaneCrossing planeCrossing(helixPos,helixDir,rho,propagationDirection());
238  return propagateWithHelixCrossing(planeCrossing,plane,fts.momentum().mag(),x,p,s);
239  }else{
240  //--- Alternative implementation to be used for the propagation of the parameters of looping
241  // particles that cross twice the (infinite) surface of the plane. It is not trivial to determine
242  // which of the two intersections has to be returned.
243 
244  //---- FIXME: WHAT FOLLOWS HAS TO BE REWRITTEN IN A CLEANER (AND CPU-OPTIMIZED) WAY ---------
245  LogDebug("AnalyticalPropagator") << "In AnaliticalProp, calling HAPC " << "\n"
246  << "plane is centered in xyz: "
247  << plane.position().x() << " , "
248  << plane.position().y() << " , "
249  << plane.position().z() << "\n";
250 
251 
252  GlobalPoint gp1 = fts.position();
253  GlobalVector gm1 = fts.momentum();
254  double s1 = 0;
255  double rho1 = fts.transverseCurvature();
256  HelixPlaneCrossing::PositionType helixPos1(gp1);
257  HelixPlaneCrossing::DirectionType helixDir1(gm1);
258  LogDebug("AnalyticalPropagator") << "gp1 before calling planeCrossing1: " << gp1 << "\n";
259  HelixArbitraryPlaneCrossing planeCrossing1(helixPos1,helixDir1,rho1,propagationDirection());
260 
263 
264  double tolerance(0.0050);
266  tolerance *=-1;
267 
268  bool check1 = propagateWithHelixCrossing(planeCrossing1,plane,fts.momentum().mag(),gp1,gm1,s1);
269  double dphi1 = fabs(fts.momentum().phi()-gm1.phi());
270  LogDebug("AnalyticalPropagator") << "check1, s1, dphi, gp1: "
271  << check1 << " , "
272  << s1 << " , "
273  << dphi1 << " , "
274  << gp1 << "\n";
275 
276  //move forward a bit to avoid that the propagator doesn't propagate because the state is already on surface.
277  //we want to go to the other point of intersection between the helix and the plane
278  xGen = planeCrossing1.position(s1+tolerance);
279  pGen = planeCrossing1.direction(s1+tolerance);
280 
281  /*
282  if(!check1 || s1>170 ){
283  //PropagationDirection newDir = (propagationDirection() == alongMomentum) ? oppositeToMomentum : alongMomentum;
284  PropagationDirection newDir = anyDirection;
285  HelixArbitraryPlaneCrossing planeCrossing1B(helixPos1,helixDir1,rho1,newDir);
286  check1 = propagateWithHelixCrossing(planeCrossing1B,plane,fts.momentum().mag(),gp1,gm1,s1);
287  LogDebug("AnalyticalPropagator") << "after second attempt, check1, s1,gp1: "
288  << check1 << " , "
289  << s1 << " , " << gp1 << "\n";
290 
291  xGen = planeCrossing1B.position(s1+tolerance);
292  pGen = planeCrossing1B.direction(s1+tolerance);
293  }
294  */
295 
296  if(!check1){
297  LogDebug("AnalyticalPropagator") << "failed also second attempt. No idea what to do, then bailout" << "\n";
298  }
299 
300 
301  pGen *= gm1.mag()/pGen.mag();
302  GlobalPoint gp2(xGen);
303  GlobalVector gm2(pGen);
304  double s2 = 0;
305  double rho2 = rho1;
306  HelixPlaneCrossing::PositionType helixPos2(gp2);
307  HelixPlaneCrossing::DirectionType helixDir2(gm2);
308  HelixArbitraryPlaneCrossing planeCrossing2(helixPos2,helixDir2,rho2,propagationDirection());
309 
310  bool check2 = propagateWithHelixCrossing(planeCrossing2,plane,gm2.mag(),gp2,gm2,s2);
311 
312  if(!check2){
313  x = gp1;
314  p = gm1;
315  s = s1;
316  return check1;
317  }
318 
319  if(!check1){
320  edm::LogError("AnalyticalPropagator") << "LOGIC ERROR: I should not have entered here!" << "\n";
321  return false;
322  }
323 
324 
325  LogDebug("AnalyticalPropagator") << "check2, s2, gp2: "
326  << check2 << " , "
327  << s2 << " , " << gp2 << "\n";
328 
329 
330  double dist1 = (plane.position()-gp1).perp();
331  double dist2 = (plane.position()-gp2).perp();
332 
333 
334  LogDebug("AnalyticalPropagator") << "propDir, dist1, dist2: "
335  << propagationDirection() << " , "
336  << dist1 << " , "
337  << dist2 << "\n";
338 
339  //If there are two solutions, the one which is the closest to the module's center is chosen
340  if(dist1<2*dist2){
341  x = gp1;
342  p = gm1;
343  s = s1;
344  return check1;
345  }else if(dist2<2*dist1){
346  x = gp2;
347  p = gm2;
348  s = s1+s2+tolerance;
349  return check2;
350  }else{
351  if(fabs(s1)<fabs(s2)){
352  x = gp1;
353  p = gm1;
354  s = s1;
355  return check1;
356  }else{
357  x = gp2;
358  p = gm2;
359  s = s1+s2+tolerance;
360  return check2;
361  }
362  }
363 
364  //-------- END of ugly piece of code ---------------
365  }
366 
367  }
368 }
#define LogDebug(id)
virtual PropagationDirection propagationDirection() const
Definition: Propagator.h:143
T mag() const
The vector magnitude. Equivalent to sqrt(vec.mag2())
Definition: DDAxes.h:10
Geom::Phi< T > phi() const
Definition: PV3DBase.h:68
tuple s2
Definition: indexGen.py:106
T mag() const
Definition: PV3DBase.h:66
GlobalVector momentum() const
GlobalPoint position() const
double transverseCurvature() const
T perp() const
Magnitude of transverse component.
bool propagateWithHelixCrossing(HelixPlaneCrossing &, const Plane &, const float, GlobalPoint &, GlobalVector &, double &s) const dso_internal
helix parameter propagation to a plane using HelixPlaneCrossing
Definition: DDAxes.h:10
bool propagateWithLineCrossing(const GlobalPoint &, const GlobalVector &, const Plane &, GlobalPoint &, double &) const dso_internal
straight line parameter propagation to a plane
bool AnalyticalPropagator::propagateWithHelixCrossing ( HelixPlaneCrossing planeCrossing,
const Plane plane,
const float  pmag,
GlobalPoint x,
GlobalVector p,
double &  s 
) const
private

helix parameter propagation to a plane using HelixPlaneCrossing

Definition at line 421 of file AnalyticalPropagator.cc.

References HelixPlaneCrossing::direction(), Basic3DVector< T >::mag(), HelixPlaneCrossing::pathLength(), and HelixPlaneCrossing::position().

426  {
427  // get solution
428  std::pair<bool,double> propResult = planeCrossing.pathLength(plane);
429  if ( !propResult.first ) return false;
430  s = propResult.second;
431  // point (reconverted to GlobalPoint)
432  HelixPlaneCrossing::PositionType xGen = planeCrossing.position(s);
433  x = GlobalPoint(xGen);
434  // direction (reconverted to GlobalVector, renormalised)
435  HelixPlaneCrossing::DirectionType pGen = planeCrossing.direction(s);
436  pGen *= pmag/pGen.mag();
437  p = GlobalVector(pGen);
438  //
439  return true;
440 }
T mag() const
The vector magnitude. Equivalent to sqrt(vec.mag2())
virtual DirectionType direction(double s) const =0
Global3DPoint GlobalPoint
Definition: GlobalPoint.h:10
virtual std::pair< bool, double > pathLength(const Plane &)=0
virtual PositionType position(double s) const =0
Global3DVector GlobalVector
Definition: GlobalVector.h:10
bool AnalyticalPropagator::propagateWithLineCrossing ( const GlobalPoint x0,
const GlobalVector p0,
const Plane plane,
GlobalPoint x,
double &  s 
) const
private

straight line parameter propagation to a plane

Definition at line 371 of file AnalyticalPropagator.cc.

References dir, StraightLinePlaneCrossing::pathLength(), pos, and StraightLinePlaneCrossing::position().

374  {
375  //
376  // Instantiate auxiliary object for finding intersection.
377  // Frame-independant point and vector are created explicitely to
378  // avoid confusing gcc (refuses to compile with temporary objects
379  // in the constructor).
380  //
384  //
385  // get solution
386  //
387  std::pair<bool,double> propResult = planeCrossing.pathLength(plane);
388  if ( !propResult.first ) return false;
389  s = propResult.second;
390  // point (reconverted to GlobalPoint)
391  x = GlobalPoint(planeCrossing.position(s));
392  //
393  return true;
394 }
virtual PropagationDirection propagationDirection() const
Definition: Propagator.h:143
Global3DPoint GlobalPoint
Definition: GlobalPoint.h:10
dbl *** dir
Definition: mlp_gen.cc:35
bool AnalyticalPropagator::propagateWithLineCrossing ( const GlobalPoint x0,
const GlobalVector p0,
const Cylinder cylinder,
GlobalPoint x,
double &  s 
) const
private

straight line parameter propagation to a cylinder

Definition at line 397 of file AnalyticalPropagator.cc.

References StraightLineBarrelCylinderCrossing::pathLength(), and StraightLineBarrelCylinderCrossing::position().

400  {
401  //
402  // Instantiate auxiliary object for finding intersection.
403  // Frame-independant point and vector are created explicitely to
404  // avoid confusing gcc (refuses to compile with temporary objects
405  // in the constructor).
406  //
408  //
409  // get solution
410  //
411  std::pair<bool,double> propResult = cylCrossing.pathLength(cylinder);
412  if ( !propResult.first ) return false;
413  s = propResult.second;
414  // point (reconverted to GlobalPoint)
415  x = cylCrossing.position(s);
416  //
417  return true;
418 }
virtual PropagationDirection propagationDirection() const
Definition: Propagator.h:143
std::pair< TrajectoryStateOnSurface, double > AnalyticalPropagator::propagateWithPath ( const FreeTrajectoryState fts,
const Plane plane 
) const
virtual

propagation to plane with path length

Implements Propagator.

Definition at line 28 of file AnalyticalPropagator.cc.

References FreeTrajectoryState::charge(), alignCSCRings::e, Plane::localZ(), LogDebug, PV3DBase< T, PVType, FrameType >::mag(), PV3DBase< T, PVType, FrameType >::mag2(), FreeTrajectoryState::momentum(), AlCaHLTBitMon_ParallelJobs::p, PV3DBase< T, PVType, FrameType >::perp2(), GloballyPositioned< T >::position(), FreeTrajectoryState::position(), rho, alignCSCRings::s, FreeTrajectoryState::transverseCurvature(), and x.

Referenced by ReferenceTrajectory::construct(), propagate(), DualTrajectoryFactory::propagateExternal(), DualBzeroTrajectoryFactory::propagateExternal(), and TwoBodyDecayTrajectoryState::propagateSingleState().

30 {
31  // check curvature
32  double rho = fts.transverseCurvature();
33 
34  // propagate parameters
35  GlobalPoint x;
37  double s;
38 
39  // check if already on plane
40  const float maxDistToPlane(0.1e-4);
41  const float numericalPrecision(5.e-7);
42  float maxDz = numericalPrecision*plane.position().mag();
43  if ( fabs(plane.localZ(fts.position()))>(maxDistToPlane>maxDz?maxDistToPlane:maxDz) ) {
44  // propagate
45  bool parametersOK = this->propagateParametersOnPlane(fts, plane, x, p, s);
46  // check status and deltaPhi limit
47  float dphi2 = s*rho;
48  dphi2 = dphi2*dphi2*fts.momentum().perp2()/fts.momentum().mag2();
49  if ( !parametersOK || dphi2>theMaxDPhi2 ) return TsosWP(TrajectoryStateOnSurface(),0.);
50  }
51  else {
52  LogDebug("AnalyticalPropagator")<<"not going anywhere. Already on surface.\n"
53  <<"plane.localZ(fts.position()): "<<plane.localZ(fts.position())<<"\n"
54  <<"maxDistToPlane: "<<maxDistToPlane<<"\n"
55  <<"maxDz: "<<maxDz<<"\n"
56  <<"plane.position().mag(): "<<plane.position().mag();
57  x = fts.position();
58  p = fts.momentum();
59  s = 0.;
60  }
61  //
62  // Compute propagated state and check change in curvature
63  //
65  if ( fabs(rho)>1.e-10 && fabs((gtp.transverseCurvature()-rho)/rho)>theMaxDBzRatio )
66  return TsosWP(TrajectoryStateOnSurface(),0.);
67  //
68  // construct TrajectoryStateOnSurface
69  //
70  return propagatedStateWithPath(fts,plane,gtp,s);
71 }
#define LogDebug(id)
T mag2() const
Definition: PV3DBase.h:65
bool propagateParametersOnPlane(const FreeTrajectoryState &fts, const Plane &plane, GlobalPoint &x, GlobalVector &p, double &s) const dso_internal
parameter propagation to plane (returns position, momentum and path length)
std::pair< TrajectoryStateOnSurface, double > propagatedStateWithPath(const FreeTrajectoryState &fts, const Surface &surface, const GlobalTrajectoryParameters &gtp, const double &s) const dso_internal
propagation of errors (if needed) and generation of a new TSOS
Definition: DDAxes.h:10
float localZ(const GlobalPoint &gp) const
Fast access to distance from plane for a point.
Definition: Plane.h:52
T perp2() const
Definition: PV3DBase.h:70
TrackCharge charge() const
T mag() const
Definition: PV3DBase.h:66
GlobalVector momentum() const
GlobalPoint position() const
double transverseCurvature() const
std::pair< TrajectoryStateOnSurface, double > TsosWP
const MagneticField * theField
Definition: DDAxes.h:10
const PositionType & position() const
std::pair< TrajectoryStateOnSurface, double > AnalyticalPropagator::propagateWithPath ( const FreeTrajectoryState fts,
const Cylinder cylinder 
) const
virtual

propagation to cylinder with path length

Implements Propagator.

Definition at line 75 of file AnalyticalPropagator.cc.

References FreeTrajectoryState::charge(), PV3DBase< T, PVType, FrameType >::mag2(), FreeTrajectoryState::momentum(), AlCaHLTBitMon_ParallelJobs::p, PV3DBase< T, PVType, FrameType >::perp2(), rho, alignCSCRings::s, Cylinder::tangentPlane(), FreeTrajectoryState::transverseCurvature(), and x.

77 {
78  // check curvature
79  double rho = fts.transverseCurvature();
80 
81  // propagate parameters
82  GlobalPoint x;
84  double s = 0;
85 
86  bool parametersOK = this->propagateParametersOnCylinder(fts, cylinder, x, p, s);
87  // check status and deltaPhi limit
88  float dphi2 = s*rho;
89  dphi2 = dphi2*dphi2*fts.momentum().perp2()/fts.momentum().mag2();
90  if ( !parametersOK || dphi2>theMaxDPhi2 ) return TsosWP(TrajectoryStateOnSurface(),0.);
91  //
92  // Compute propagated state and check change in curvature
93  //
95  if ( fabs(rho)>1.e-10 && fabs((gtp.transverseCurvature()-rho)/rho)>theMaxDBzRatio )
96  return TsosWP(TrajectoryStateOnSurface(),0.);
97  //
98  // create result TSOS on TangentPlane (local parameters & errors are better defined)
99  //
100 
101  //try {
102  ReferenceCountingPointer<TangentPlane> plane(cylinder.tangentPlane(x)); // need to be here until tsos is created!
103  return propagatedStateWithPath(fts,*plane,gtp,s);
104  /*
105  } catch(...) {
106  std::cout << "wrong tangent to cylinder " << x
107  << " pos, rad " << cylinder.position() << " " << cylinder.radius()
108  << std::endl;
109  return TsosWP(TrajectoryStateOnSurface(),0.);
110  }
111  */
112 }
T mag2() const
Definition: PV3DBase.h:65
virtual ReferenceCountingPointer< TangentPlane > tangentPlane(const GlobalPoint &) const
tangent plane to surface from global point
Definition: Cylinder.cc:23
std::pair< TrajectoryStateOnSurface, double > propagatedStateWithPath(const FreeTrajectoryState &fts, const Surface &surface, const GlobalTrajectoryParameters &gtp, const double &s) const dso_internal
propagation of errors (if needed) and generation of a new TSOS
Definition: DDAxes.h:10
T perp2() const
Definition: PV3DBase.h:70
TrackCharge charge() const
bool propagateParametersOnCylinder(const FreeTrajectoryState &fts, const Cylinder &cylinder, GlobalPoint &x, GlobalVector &p, double &s) const dso_internal
parameter propagation to cylinder (returns position, momentum and path length)
GlobalVector momentum() const
double transverseCurvature() const
std::pair< TrajectoryStateOnSurface, double > TsosWP
const MagneticField * theField
Definition: DDAxes.h:10
virtual bool AnalyticalPropagator::setMaxDirectionChange ( float  phiMax)
inlinevirtual

limitation of change in transverse direction (to avoid loops).

Reimplemented from Propagator.

Definition at line 69 of file AnalyticalPropagator.h.

References jptDQMConfig_cff::phiMax, and theMaxDPhi2.

69  {
71  return true;
72  }
void AnalyticalPropagator::setMaxRelativeChangeInBz ( const float  maxDBz)
inline

Set the maximum relative change in Bz (Bz_at_end-Bz_at_start)/Bz_at_start for a single propagation. The default is no limit. NB: this propagator assumes constant, non-zero magnetic field parallel to the z-axis!

Definition at line 87 of file AnalyticalPropagator.h.

References theMaxDBzRatio.

87  {
88  theMaxDBzRatio = maxDBz;
89  }

Member Data Documentation

bool AnalyticalPropagator::isOldPropagationType
private

Definition at line 131 of file AnalyticalPropagator.h.

const MagneticField* AnalyticalPropagator::theField
private

Definition at line 130 of file AnalyticalPropagator.h.

Referenced by magneticField().

float AnalyticalPropagator::theMaxDBzRatio
private

Definition at line 129 of file AnalyticalPropagator.h.

Referenced by setMaxRelativeChangeInBz().

float AnalyticalPropagator::theMaxDPhi2
private

Definition at line 128 of file AnalyticalPropagator.h.

Referenced by setMaxDirectionChange().