CMS 3D CMS Logo

List of all members | Classes | Public Types | Public Member Functions | Protected Member Functions | Protected Attributes | Static Protected Attributes
ReferenceTrajectoryBase Class Referenceabstract

#include <ReferenceTrajectoryBase.h>

Inheritance diagram for ReferenceTrajectoryBase:
ReferenceCounted DualReferenceTrajectory ReferenceTrajectory TwoBodyDecayTrajectory DualBzeroReferenceTrajectory BzeroReferenceTrajectory

Classes

struct  Config
 

Public Types

enum  MaterialEffects {
  none, multipleScattering, energyLoss, combined,
  breakPoints, brokenLinesCoarse, brokenLinesFine, localGBL,
  curvlinGBL
}
 
typedef ReferenceCountingPointer< ReferenceTrajectoryBaseReferenceTrajectoryPtr
 

Public Member Functions

virtual ReferenceTrajectoryBaseclone () const =0
 
const AlgebraicMatrixderivatives () const
 
const Eigen::MatrixXd & gblExtDerivatives () const
 
const Eigen::VectorXd & gblExtMeasurements () const
 
const Eigen::VectorXd & gblExtPrecisions () const
 
std::vector< std::pair< std::vector< gbl::GblPoint >, Eigen::MatrixXd > > & gblInput ()
 
bool isValid ()
 
const AlgebraicMatrixlocalToTrajectory () const
 
const AlgebraicSymMatrixmeasurementErrors () const
 
const AlgebraicVectormeasurements () const
 
int nominalField () const
 
unsigned int numberOfHitMeas () const
 
unsigned int numberOfHits () const
 
unsigned int numberOfPar () const
 
unsigned int numberOfVirtualMeas () const
 
unsigned int numberOfVirtualPar () const
 
const AlgebraicSymMatrixparameterErrors () const
 
bool parameterErrorsAvailable () const
 
const AlgebraicVectorparameters () const
 
const TransientTrackingRecHit::ConstRecHitContainerrecHits () const
 
void setParameterErrors (const AlgebraicSymMatrix &error)
 
const AlgebraicSymMatrixtrajectoryPositionErrors () const
 
const AlgebraicVectortrajectoryPositions () const
 
const std::vector< TrajectoryStateOnSurface > & trajectoryStates () const
 
const AlgebraicMatrixtrajectoryToCurv () const
 
 ~ReferenceTrajectoryBase () override
 

Protected Member Functions

unsigned int numberOfUsedRecHits (const TransientTrackingRecHit::ConstRecHitContainer &recHits) const
 
 ReferenceTrajectoryBase (unsigned int nPar, unsigned int nHits, unsigned int nVirtualPar, unsigned int nVirtualMeas)
 
bool useRecHit (const TransientTrackingRecHit::ConstRecHitPointer &hitPtr) const
 

Protected Attributes

AlgebraicMatrix theDerivatives
 
Eigen::MatrixXd theGblExtDerivatives
 
Eigen::VectorXd theGblExtMeasurements
 
Eigen::VectorXd theGblExtPrecisions
 
std::vector< std::pair< std::vector< gbl::GblPoint >, Eigen::MatrixXd > > theGblInput
 
AlgebraicMatrix theInnerLocalToTrajectory
 
AlgebraicMatrix theInnerTrajectoryToCurvilinear
 
AlgebraicVector theMeasurements
 
AlgebraicSymMatrix theMeasurementsCov
 
int theNomField
 
unsigned int theNumberOfHits
 
unsigned int theNumberOfPars
 
unsigned int theNumberOfVirtualMeas
 
unsigned int theNumberOfVirtualPars
 
bool theParamCovFlag
 
AlgebraicSymMatrix theParameterCov
 
AlgebraicVector theParameters
 
TransientTrackingRecHit::ConstRecHitContainer theRecHits
 
AlgebraicSymMatrix theTrajectoryPositionCov
 
AlgebraicVector theTrajectoryPositions
 
std::vector< TrajectoryStateOnSurfacetheTsosVec
 
bool theValidityFlag
 

Static Protected Attributes

static constexpr unsigned int nMeasPerHit {2}
 

Detailed Description

Author : Gero Flucke (based on code for ORCA by Edmund Widl) date : 2006/09/17 last update:

Date
2010/09/10 12:11:39

by :

Author
mussgill

Base class for reference 'trajectories' of single- or multiparticles stated. Inheriting classes have to calculate all relevant quantities accessed through member functions of this base class:

The local measured x/y coordinates on all crossed detectors as vector:

m = (x1, y1, x2, y2, ..., xN, yN) [transposed vector shown]

their covariance matrix (possibly containing correlations between hits due to material effects which are not taken into account by the ideal trajectory parametrisation, cf. below),

similarly the local x/y cordinates of the reference trajectory with covariance,

the parameters of the (ideal) 'trajectory' (e.g. 5 values for a single track or 9 for a two-track-state with vertex constraint),

the derivatives of the local coordinates of the reference trajectory with respect to the initial 'trajectory'-parameters, e.g. for n parameters 'p' and N hits with 'x/y' coordinates:

D = ( dx1/dp1, dx1/dp2, ..., dx1/dpn,

  dy1/dp1, dy1/dp2, ..., dy1/dpn,

  dx2/dp1, dx2/dp2, ..., dx2/dpn,

  dy2/dp1, dy2/dp2, ..., dy2/dpn,

     .        .             .

     .        .             .

  dxN/dp1, dxN/dp2, ..., dxN/dpn,

  dyN/dp1, dyN/dp2, ..., dyN/dpn )

and finally the TrajectoryStateOnSurface's of the reference trajectory.

Take care to check validity of the calculation (isValid()) before using the results.

.............................................................................

090730 C. Kleinwort: 'Break Points' introduced for better description of multiple scattering (correlations)

For each detector layer (in the hit collection) two ortogonal scattering angles are introduced as new local track parameters ("break points"). To constrain those to the expected mean (=0.) and variance (~1/p^2 X/X0) corresponding measurements are added. Invalid hits may be used to produce break points too (UseInvalidHits = True).

Break points have been implemented for: ReferenceTrajectory, BzeroReferenceTrajectory, DualReferenceTrajectory and DualBzeroReferenceTrajectory For 'Dual' trajectories they make no sense as only the correlations between the hits
in a track half are accounted for and not those between the halves!

Break Points are selected by TrajectoryFactory.MaterialEffects = "BreakPoints"

090909 C. Kleinwort: 'Broken Lines' introduced for description of multiple scattering (V. Blobel, Nuclear Instruments and Methods A, 566 (2006), pp. 14-17) Fine Broken Lines are selected by TrajectoryFactory.MaterialEffects = "BrokenLinesFine" (exact derivatives) Coarse Broken Lines are selected by TrajectoryFactory.MaterialEffects = "BrokenLines[Coarse]" (approximate derivatives, closeby hits (ds<1cm) combined)

091127 C. Kleinwort: 1) 'Broken Lines' extended to PCA, required for beamspot constraint and TwoBodyDecayTrajectory, selected with "BrokenLines[Coarse]Pca" or "BrokenLinesFinePca" 2) For coarse Broken Lines linear interpolation is used for combined hits 3) TwoBodyDecayTrajectory implemented for break points and Broken Lines

141103 C. Kleinwort: 'General Broken Lines' introduced for description of multiple scattering (C. Kleinwort, Nuclear Instruments and Methods A, 673 (2012), pp. 107-110) needs GBL version >= V01-13-00 (from svnsrv.desy.de) Selected by TrajectoryFactory.MaterialEffects = "LocalGBL" or = "CurvlinGBL" (for trajectory constructed in local or curvilinear system)

Definition at line 105 of file ReferenceTrajectoryBase.h.

Member Typedef Documentation

◆ ReferenceTrajectoryPtr

Definition at line 107 of file ReferenceTrajectoryBase.h.

Member Enumeration Documentation

◆ MaterialEffects

Constructor & Destructor Documentation

◆ ~ReferenceTrajectoryBase()

ReferenceTrajectoryBase::~ReferenceTrajectoryBase ( )
inlineoverride

Definition at line 140 of file ReferenceTrajectoryBase.h.

140 {}

◆ ReferenceTrajectoryBase()

ReferenceTrajectoryBase::ReferenceTrajectoryBase ( unsigned int  nPar,
unsigned int  nHits,
unsigned int  nVirtualPar,
unsigned int  nVirtualMeas 
)
explicitprotected

Definition at line 3 of file ReferenceTrajectoryBase.cc.

References nHits, theRecHits, and theTsosVec.

7  : theValidityFlag(false),
8  theParamCovFlag(false),
11  theNumberOfVirtualMeas(nVirtualMeas),
12  theNumberOfVirtualPars(nVirtualPar),
13  theTsosVec(),
14  theRecHits(),
15  theMeasurements(nMeasPerHit * nHits + nVirtualMeas),
16  theMeasurementsCov(nMeasPerHit * nHits + nVirtualMeas, 0),
21  theDerivatives(nMeasPerHit * nHits + nVirtualMeas, nPar + nVirtualPar, 0),
23  theInnerLocalToTrajectory(5, 5, 0) {
24  theTsosVec.reserve(nHits);
25  theRecHits.reserve(nHits);
26 }
AlgebraicMatrix theInnerTrajectoryToCurvilinear
AlgebraicMatrix theInnerLocalToTrajectory
AlgebraicSymMatrix theTrajectoryPositionCov
TransientTrackingRecHit::ConstRecHitContainer theRecHits
AlgebraicSymMatrix theMeasurementsCov
static constexpr unsigned int nMeasPerHit
AlgebraicSymMatrix theParameterCov
std::vector< TrajectoryStateOnSurface > theTsosVec
TupleMultiplicity< TrackerTraits > const *__restrict__ uint32_t nHits

Member Function Documentation

◆ clone()

virtual ReferenceTrajectoryBase* ReferenceTrajectoryBase::clone ( ) const
pure virtual

◆ derivatives()

const AlgebraicMatrix& ReferenceTrajectoryBase::derivatives ( void  ) const
inline

Returns the derivatives of the local coordinates of the reference trajectory (i.e. trajectoryPositions) w.r.t. the initial 'trajectory'-parameters.

Definition at line 164 of file ReferenceTrajectoryBase.h.

References theDerivatives.

Referenced by TwoBodyDecayTrajectory::construct(), and DualReferenceTrajectory::construct().

164 { return theDerivatives; }

◆ gblExtDerivatives()

const Eigen::MatrixXd& ReferenceTrajectoryBase::gblExtDerivatives ( ) const
inline

Returns the GBL external derivatives.

Definition at line 179 of file ReferenceTrajectoryBase.h.

References theGblExtDerivatives.

179 { return theGblExtDerivatives; }

◆ gblExtMeasurements()

const Eigen::VectorXd& ReferenceTrajectoryBase::gblExtMeasurements ( ) const
inline

Returns the GBL external derivatives.

Definition at line 183 of file ReferenceTrajectoryBase.h.

References theGblExtMeasurements.

183 { return theGblExtMeasurements; }

◆ gblExtPrecisions()

const Eigen::VectorXd& ReferenceTrajectoryBase::gblExtPrecisions ( ) const
inline

Returns the GBL external derivatives.

Definition at line 187 of file ReferenceTrajectoryBase.h.

References theGblExtPrecisions.

187 { return theGblExtPrecisions; }

◆ gblInput()

std::vector<std::pair<std::vector<gbl::GblPoint>, Eigen::MatrixXd> >& ReferenceTrajectoryBase::gblInput ( )
inline

Returns the GBL input

Definition at line 175 of file ReferenceTrajectoryBase.h.

References theGblInput.

Referenced by TwoBodyDecayTrajectory::construct().

175 { return theGblInput; }
std::vector< std::pair< std::vector< gbl::GblPoint >, Eigen::MatrixXd > > theGblInput

◆ isValid()

bool ReferenceTrajectoryBase::isValid ( void  )
inline

◆ localToTrajectory()

const AlgebraicMatrix& ReferenceTrajectoryBase::localToTrajectory ( ) const
inline

Returns the transformation of local to tracjectory parameters

Definition at line 171 of file ReferenceTrajectoryBase.h.

References theInnerLocalToTrajectory.

Referenced by TwoBodyDecayTrajectory::construct().

171 { return theInnerLocalToTrajectory; }
AlgebraicMatrix theInnerLocalToTrajectory

◆ measurementErrors()

const AlgebraicSymMatrix& ReferenceTrajectoryBase::measurementErrors ( ) const
inline

Returns the full covariance matrix of the measurements. ORCA-equivalent: covariance()

Definition at line 150 of file ReferenceTrajectoryBase.h.

References theMeasurementsCov.

Referenced by TwoBodyDecayTrajectory::construct(), and DualReferenceTrajectory::construct().

150 { return theMeasurementsCov; }
AlgebraicSymMatrix theMeasurementsCov

◆ measurements()

const AlgebraicVector& ReferenceTrajectoryBase::measurements ( ) const
inline

Returns the measurements in local coordinates.

Definition at line 146 of file ReferenceTrajectoryBase.h.

References theMeasurements.

Referenced by TwoBodyDecayTrajectory::construct(), and DualReferenceTrajectory::construct().

146 { return theMeasurements; }

◆ nominalField()

int ReferenceTrajectoryBase::nominalField ( ) const
inline

Definition at line 222 of file ReferenceTrajectoryBase.h.

References theNomField.

Referenced by TwoBodyDecayTrajectory::construct().

◆ numberOfHitMeas()

unsigned int ReferenceTrajectoryBase::numberOfHitMeas ( ) const
inline

◆ numberOfHits()

unsigned int ReferenceTrajectoryBase::numberOfHits ( ) const
inline

Definition at line 217 of file ReferenceTrajectoryBase.h.

References theNumberOfHits.

Referenced by TwoBodyDecayTrajectory::construct().

217 { return theNumberOfHits; }

◆ numberOfPar()

unsigned int ReferenceTrajectoryBase::numberOfPar ( ) const
inline

Definition at line 218 of file ReferenceTrajectoryBase.h.

References theNumberOfPars.

Referenced by TwoBodyDecayTrajectory::construct().

218 { return theNumberOfPars; }

◆ numberOfUsedRecHits()

unsigned int ReferenceTrajectoryBase::numberOfUsedRecHits ( const TransientTrackingRecHit::ConstRecHitContainer recHits) const
protected

Definition at line 28 of file ReferenceTrajectoryBase.cc.

References recHits(), and useRecHit().

29  {
30  unsigned int nUsedHits = 0;
31  TransientTrackingRecHit::ConstRecHitContainer::const_iterator itHit;
32  for (itHit = recHits.begin(); itHit != recHits.end(); ++itHit)
33  if (useRecHit(*itHit))
34  ++nUsedHits;
35  return nUsedHits;
36 }
const TransientTrackingRecHit::ConstRecHitContainer & recHits() const
bool useRecHit(const TransientTrackingRecHit::ConstRecHitPointer &hitPtr) const

◆ numberOfVirtualMeas()

unsigned int ReferenceTrajectoryBase::numberOfVirtualMeas ( ) const
inline

◆ numberOfVirtualPar()

unsigned int ReferenceTrajectoryBase::numberOfVirtualPar ( ) const
inline

Definition at line 220 of file ReferenceTrajectoryBase.h.

References theNumberOfVirtualPars.

Referenced by TwoBodyDecayTrajectory::construct().

220 { return theNumberOfVirtualPars; }

◆ parameterErrors()

const AlgebraicSymMatrix& ReferenceTrajectoryBase::parameterErrors ( ) const
inline

Returns the covariance matrix of the 'track'-parameters.

Definition at line 206 of file ReferenceTrajectoryBase.h.

References theParameterCov.

206 { return theParameterCov; }
AlgebraicSymMatrix theParameterCov

◆ parameterErrorsAvailable()

bool ReferenceTrajectoryBase::parameterErrorsAvailable ( ) const
inline

Returns true if the covariance matrix of the 'track'-parameters is set.

Definition at line 195 of file ReferenceTrajectoryBase.h.

References theParamCovFlag.

◆ parameters()

const AlgebraicVector& ReferenceTrajectoryBase::parameters ( void  ) const
inline

Returns the set of 'track'-parameters.

Definition at line 191 of file ReferenceTrajectoryBase.h.

References theParameters.

Referenced by ReferenceTrajectory::fillDerivatives().

191 { return theParameters; }

◆ recHits()

const TransientTrackingRecHit::ConstRecHitContainer& ReferenceTrajectoryBase::recHits ( ) const
inline

Returns the TransientTrackingRecHits (as ConstRecHitPointer's), order might depend on concrete implementation of inheriting class

Definition at line 215 of file ReferenceTrajectoryBase.h.

References theRecHits.

Referenced by BzeroReferenceTrajectory::BzeroReferenceTrajectory(), TwoBodyDecayTrajectory::construct(), DualBzeroReferenceTrajectory::construct(), DualReferenceTrajectory::construct(), ReferenceTrajectory::construct(), numberOfUsedRecHits(), ReferenceTrajectory::ReferenceTrajectory(), and TwoBodyDecayTrajectory::TwoBodyDecayTrajectory().

215 { return theRecHits; }
TransientTrackingRecHit::ConstRecHitContainer theRecHits

◆ setParameterErrors()

void ReferenceTrajectoryBase::setParameterErrors ( const AlgebraicSymMatrix error)
inline

Set the covariance matrix of the 'track'-parameters.

Definition at line 199 of file ReferenceTrajectoryBase.h.

References relativeConstraints::error, theParamCovFlag, and theParameterCov.

◆ trajectoryPositionErrors()

const AlgebraicSymMatrix& ReferenceTrajectoryBase::trajectoryPositionErrors ( ) const
inline

Returns the covariance matrix of the reference trajectory.

Definition at line 159 of file ReferenceTrajectoryBase.h.

References theTrajectoryPositionCov.

Referenced by DualReferenceTrajectory::construct().

159 { return theTrajectoryPositionCov; }
AlgebraicSymMatrix theTrajectoryPositionCov

◆ trajectoryPositions()

const AlgebraicVector& ReferenceTrajectoryBase::trajectoryPositions ( ) const
inline

Returns the local coordinates of the reference trajectory. ORCA-equivalent: referenceTrack()

Definition at line 155 of file ReferenceTrajectoryBase.h.

References theTrajectoryPositions.

Referenced by TwoBodyDecayTrajectory::construct(), and DualReferenceTrajectory::construct().

155 { return theTrajectoryPositions; }

◆ trajectoryStates()

const std::vector<TrajectoryStateOnSurface>& ReferenceTrajectoryBase::trajectoryStates ( ) const
inline

Returns the Tsos at the surfaces of the hits, parallel to recHits()

Definition at line 210 of file ReferenceTrajectoryBase.h.

References theTsosVec.

Referenced by TwoBodyDecayTrajectory::construct(), DualReferenceTrajectory::construct(), and TwoBodyDecayTrajectory::constructTsosVecWithErrors().

210 { return theTsosVec; }
std::vector< TrajectoryStateOnSurface > theTsosVec

◆ trajectoryToCurv()

const AlgebraicMatrix& ReferenceTrajectoryBase::trajectoryToCurv ( ) const
inline

Returns the transformation of tracjectory to curvilinear parameters

Definition at line 168 of file ReferenceTrajectoryBase.h.

References theInnerTrajectoryToCurvilinear.

AlgebraicMatrix theInnerTrajectoryToCurvilinear

◆ useRecHit()

bool ReferenceTrajectoryBase::useRecHit ( const TransientTrackingRecHit::ConstRecHitPointer hitPtr) const
protected

Definition at line 38 of file ReferenceTrajectoryBase.cc.

Referenced by ReferenceTrajectory::construct(), ReferenceTrajectory::getHitProjectionMatrix(), and numberOfUsedRecHits().

38  {
39  return hitPtr->isValid();
40 }

Member Data Documentation

◆ nMeasPerHit

constexpr unsigned int ReferenceTrajectoryBase::nMeasPerHit {2}
staticprotected

◆ theDerivatives

AlgebraicMatrix ReferenceTrajectoryBase::theDerivatives
protected

◆ theGblExtDerivatives

Eigen::MatrixXd ReferenceTrajectoryBase::theGblExtDerivatives
protected

◆ theGblExtMeasurements

Eigen::VectorXd ReferenceTrajectoryBase::theGblExtMeasurements
protected

◆ theGblExtPrecisions

Eigen::VectorXd ReferenceTrajectoryBase::theGblExtPrecisions
protected

◆ theGblInput

std::vector<std::pair<std::vector<gbl::GblPoint>, Eigen::MatrixXd> > ReferenceTrajectoryBase::theGblInput
protected

◆ theInnerLocalToTrajectory

AlgebraicMatrix ReferenceTrajectoryBase::theInnerLocalToTrajectory
protected

◆ theInnerTrajectoryToCurvilinear

AlgebraicMatrix ReferenceTrajectoryBase::theInnerTrajectoryToCurvilinear
protected

◆ theMeasurements

AlgebraicVector ReferenceTrajectoryBase::theMeasurements
protected

◆ theMeasurementsCov

AlgebraicSymMatrix ReferenceTrajectoryBase::theMeasurementsCov
protected

◆ theNomField

int ReferenceTrajectoryBase::theNomField
protected

◆ theNumberOfHits

unsigned int ReferenceTrajectoryBase::theNumberOfHits
protected

◆ theNumberOfPars

unsigned int ReferenceTrajectoryBase::theNumberOfPars
protected

◆ theNumberOfVirtualMeas

unsigned int ReferenceTrajectoryBase::theNumberOfVirtualMeas
protected

◆ theNumberOfVirtualPars

unsigned int ReferenceTrajectoryBase::theNumberOfVirtualPars
protected

◆ theParamCovFlag

bool ReferenceTrajectoryBase::theParamCovFlag
protected

Definition at line 236 of file ReferenceTrajectoryBase.h.

Referenced by parameterErrorsAvailable(), and setParameterErrors().

◆ theParameterCov

AlgebraicSymMatrix ReferenceTrajectoryBase::theParameterCov
protected

Definition at line 253 of file ReferenceTrajectoryBase.h.

Referenced by parameterErrors(), and setParameterErrors().

◆ theParameters

AlgebraicVector ReferenceTrajectoryBase::theParameters
protected

◆ theRecHits

TransientTrackingRecHit::ConstRecHitContainer ReferenceTrajectoryBase::theRecHits
protected

◆ theTrajectoryPositionCov

AlgebraicSymMatrix ReferenceTrajectoryBase::theTrajectoryPositionCov
protected

◆ theTrajectoryPositions

AlgebraicVector ReferenceTrajectoryBase::theTrajectoryPositions
protected

◆ theTsosVec

std::vector<TrajectoryStateOnSurface> ReferenceTrajectoryBase::theTsosVec
protected

◆ theValidityFlag

bool ReferenceTrajectoryBase::theValidityFlag
protected