TwoBodyDecayTrajectory Class Reference

#include <TwoBodyDecayTrajectory.h>

Inheritance diagram for TwoBodyDecayTrajectory:

Public Types
typedef std::pair < ConstRecHitContainer, ConstRecHitContainer >	ConstRecHitCollection
typedef TransientTrackingRecHit::ConstRecHitContainer	ConstRecHitContainer
Public Types inherited from ReferenceTrajectoryBase
enum	MaterialEffects {   none, multipleScattering, energyLoss, combined,   breakPoints, brokenLinesCoarse, brokenLinesFine, localGBL,   curvlinGBL }
typedef ReferenceCountingPointer < ReferenceTrajectoryBase >	ReferenceTrajectoryPtr

Public Member Functions
virtual TwoBodyDecayTrajectory *	clone (void) const
const std::pair< int, int >	numberOfRecHits (void)
	TwoBodyDecayTrajectory (const TwoBodyDecayTrajectoryState &trajectoryState, const ConstRecHitCollection &recHits, const MagneticField *magField, MaterialEffects materialEffects, PropagationDirection propDir, bool hitsAreReverse, const reco::BeamSpot &beamSpot, bool useRefittedState, bool constructTsosWithErrors)
	TwoBodyDecayTrajectory (void)
	~TwoBodyDecayTrajectory (void)
Public Member Functions inherited from ReferenceTrajectoryBase
const AlgebraicMatrix &	derivatives () const
const TMatrixD &	gblExtDerivatives () const
const TVectorD &	gblExtMeasurements () const
const TVectorD &	gblExtPrecisions () const
std::vector< std::pair < std::vector< gbl::GblPoint > , TMatrixD > > &	gblInput ()
bool	isValid ()
const AlgebraicMatrix &	localToTrajectory () const
const AlgebraicSymMatrix &	measurementErrors () const
const AlgebraicVector &	measurements () 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 AlgebraicSymMatrix &	parameterErrors () const
bool	parameterErrorsAvailable () const
const AlgebraicVector &	parameters () const
const TransientTrackingRecHit::ConstRecHitContainer &	recHits () const
void	setParameterErrors (const AlgebraicSymMatrix &error)
const AlgebraicSymMatrix &	trajectoryPositionErrors () const
const AlgebraicVector &	trajectoryPositions () const
const std::vector < TrajectoryStateOnSurface > &	trajectoryStates () const
const AlgebraicMatrix &	trajectoryToCurv () const
virtual	~ReferenceTrajectoryBase ()

Private Member Functions
bool	construct (const TwoBodyDecayTrajectoryState &state, const ConstRecHitCollection &recHits, const MagneticField *field, MaterialEffects materialEffects, PropagationDirection propDir, const reco::BeamSpot &beamSpot, bool useRefittedState, bool constructTsosWithErrors)
void	constructSingleTsosWithErrors (const TrajectoryStateOnSurface &tsos, int iTsos, const MagneticField *field)
void	constructTsosVecWithErrors (const ReferenceTrajectory &traj1, const ReferenceTrajectory &traj2, const MagneticField *field)

Private Attributes
std::pair< int, int >	theNumberOfRecHits

Additional Inherited Members
Protected Member Functions inherited from ReferenceTrajectoryBase
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 inherited from ReferenceTrajectoryBase
AlgebraicMatrix	theDerivatives
TMatrixD	theGblExtDerivatives
TVectorD	theGblExtMeasurements
TVectorD	theGblExtPrecisions
std::vector< std::pair < std::vector< gbl::GblPoint > , TMatrixD > >	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 < TrajectoryStateOnSurface >	theTsosVec
bool	theValidityFlag
Static Protected Attributes inherited from ReferenceTrajectoryBase
static const unsigned int	nMeasPerHit = 2

Detailed Description

Definition at line 14 of file TwoBodyDecayTrajectory.h.

Member Typedef Documentation

typedef std::pair< ConstRecHitContainer, ConstRecHitContainer > TwoBodyDecayTrajectory::ConstRecHitCollection

Definition at line 20 of file TwoBodyDecayTrajectory.h.

typedef TransientTrackingRecHit::ConstRecHitContainer TwoBodyDecayTrajectory::ConstRecHitContainer

Definition at line 19 of file TwoBodyDecayTrajectory.h.

Constructor & Destructor Documentation

TwoBodyDecayTrajectory::TwoBodyDecayTrajectory	(	const TwoBodyDecayTrajectoryState &	trajectoryState,
		const ConstRecHitCollection &	recHits,
		const MagneticField *	magField,
		MaterialEffects	materialEffects,
		PropagationDirection	propDir,
		bool	hitsAreReverse,
		const reco::BeamSpot &	beamSpot,
		bool	useRefittedState,
		bool	constructTsosWithErrors
	)

Definition at line 15 of file TwoBodyDecayTrajectory.cc.

References construct(), and ReferenceTrajectoryBase::theValidityFlag.

  : ReferenceTrajectoryBase( 
     TwoBodyDecayParameters::dimension, recHits.first.size() + recHits.second.size(),
  (materialEffects >= breakPoints) ? 2*(recHits.first.size() + recHits.second.size())-4 : 0,
  (materialEffects >= breakPoints) ? 2*(recHits.first.size() + recHits.second.size())-3 : 1 )
{
  if ( hitsAreReverse )
  {
    TransientTrackingRecHit::ConstRecHitContainer::const_reverse_iterator itRecHits;
    ConstRecHitCollection fwdRecHits;

    fwdRecHits.first.reserve( recHits.first.size() );
    for ( itRecHits = recHits.first.rbegin(); itRecHits != recHits.first.rend(); ++itRecHits )
    {
      fwdRecHits.first.push_back( *itRecHits );
    }

    fwdRecHits.second.reserve( recHits.second.size() );
    for ( itRecHits = recHits.second.rbegin(); itRecHits != recHits.second.rend(); ++itRecHits )
    {
      fwdRecHits.second.push_back( *itRecHits );
    }

    theValidityFlag = this->construct( trajectoryState, fwdRecHits, magField, materialEffects, propDir,
                       beamSpot, useRefittedState, constructTsosWithErrors );
  }
  else
  {
    theValidityFlag = this->construct( trajectoryState, recHits, magField, materialEffects, propDir,
                       beamSpot, useRefittedState, constructTsosWithErrors );
  }
}

TwoBodyDecayTrajectory::TwoBodyDecayTrajectory

(

void

)

Definition at line 58 of file TwoBodyDecayTrajectory.cc.

Referenced by clone().

  : ReferenceTrajectoryBase( 0, 0, 0, 0)
{}

TwoBodyDecayTrajectory::~TwoBodyDecayTrajectory ( void  )

inline

Definition at line 34 of file TwoBodyDecayTrajectory.h.

{}

Member Function Documentation

virtual TwoBodyDecayTrajectory* TwoBodyDecayTrajectory::clone ( void  ) const

inlinevirtual

Implements ReferenceTrajectoryBase.

Definition at line 36 of file TwoBodyDecayTrajectory.h.

References TwoBodyDecayTrajectory().

    { return new TwoBodyDecayTrajectory( *this ); }

bool TwoBodyDecayTrajectory::construct ( const TwoBodyDecayTrajectoryState &  state, const ConstRecHitCollection &  recHits, const MagneticField *  field, MaterialEffects  materialEffects, PropagationDirection  propDir, const reco::BeamSpot &  beamSpot, bool  useRefittedState, bool  constructTsosWithErrors  )

private

Definition at line 63 of file TwoBodyDecayTrajectory.cc.

References cuy::col, constructTsosVecWithErrors(), TwoBodyDecayParameters::covariance(), TwoBodyDecayTrajectoryState::decayParameters(), TwoBodyDecayTrajectoryState::derivatives(), ReferenceTrajectoryBase::derivatives(), ReferenceTrajectoryBase::gblInput(), ReferenceTrajectoryBase::isValid(), ReferenceTrajectoryBase::localGBL, ReferenceTrajectoryBase::localToTrajectory(), TwoBodyDecayParameters::mass, ReferenceTrajectoryBase::measurementErrors(), ReferenceTrajectoryBase::measurements(), ReferenceTrajectoryBase::nominalField(), ReferenceTrajectoryBase::numberOfHitMeas(), ReferenceTrajectoryBase::numberOfHits(), ReferenceTrajectoryBase::numberOfPar(), ReferenceTrajectoryBase::numberOfVirtualMeas(), ReferenceTrajectoryBase::numberOfVirtualPar(), TwoBodyDecayParameters::parameters(), TwoBodyDecayTrajectoryState::particleMass(), TwoBodyDecayTrajectoryState::primaryMass(), TwoBodyDecayTrajectoryState::primaryWidth(), ReferenceTrajectoryBase::theDerivatives, ReferenceTrajectoryBase::theGblExtDerivatives, ReferenceTrajectoryBase::theGblExtMeasurements, ReferenceTrajectoryBase::theGblExtPrecisions, ReferenceTrajectoryBase::theGblInput, ReferenceTrajectoryBase::theMeasurements, ReferenceTrajectoryBase::theMeasurementsCov, ReferenceTrajectoryBase::theNomField, ReferenceTrajectoryBase::theNumberOfHits, ReferenceTrajectoryBase::theNumberOfPars, theNumberOfRecHits, ReferenceTrajectoryBase::theNumberOfVirtualMeas, ReferenceTrajectoryBase::theNumberOfVirtualPars, ReferenceTrajectoryBase::theParameters, ReferenceTrajectoryBase::theRecHits, ReferenceTrajectoryBase::theTrajectoryPositionCov, ReferenceTrajectoryBase::theTrajectoryPositions, ReferenceTrajectoryBase::theTsosVec, ReferenceTrajectoryBase::trajectoryPositions(), TwoBodyDecayTrajectoryState::trajectoryStates(), and ReferenceTrajectoryBase::trajectoryStates().

Referenced by TwoBodyDecayTrajectory().

{  
  const TwoBodyDecayTrajectoryState::TsosContainer& tsos = state.trajectoryStates( useRefittedState );
  const TwoBodyDecayTrajectoryState::Derivatives& deriv = state.derivatives();
  double mass = state.particleMass();

  //
  // first track
  //

  // construct a trajectory (hits should be already in correct order)
  ReferenceTrajectory trajectory1( tsos.first, recHits.first, false, field, materialEffects,
                   propDir, mass, false, beamSpot);

  // check if construction of trajectory was successful
  if ( !trajectory1.isValid() ) return false;

  //
  // second track
  //

  ReferenceTrajectory trajectory2( tsos.second, recHits.second, false, field, materialEffects,
                   propDir, mass, false, beamSpot );

  if ( !trajectory2.isValid() ) return false;
  
  //
  // combine both tracks
  //
  unsigned int nLocal = deriv.first.num_row();
  unsigned int nTbd   = deriv.first.num_col();

  if (materialEffects >= localGBL) {
    // GBL trajectory inputs
    // convert to TMatrix
    TMatrixD tbdToLocal1(nLocal, nTbd);
    for (unsigned int row = 0; row < nLocal; ++row) {
      for (unsigned int col = 0; col < nTbd; ++col) {
        tbdToLocal1(row,col) = deriv.first[row][col];
      }
    }
    // add first body
    theGblInput.push_back(std::make_pair(trajectory1.gblInput().front().first, 
                                         trajectory1.gblInput().front().second*tbdToLocal1));
    // convert to TMatrix
    TMatrixD tbdToLocal2(nLocal, nTbd);
    for (unsigned int row = 0; row < nLocal; ++row) {
      for (unsigned int col = 0; col < nTbd; ++col) {
        tbdToLocal2(row,col) = deriv.second[row][col];
      }
    }
    // add second body
    theGblInput.push_back(std::make_pair(trajectory2.gblInput().front().first, 
                                         trajectory2.gblInput().front().second*tbdToLocal2));
    // add virtual mass measurement
    theGblExtDerivatives.ResizeTo(1,nTbd);
    theGblExtDerivatives(0,TwoBodyDecayParameters::mass) = 1.0;
    theGblExtMeasurements.ResizeTo(1);
    theGblExtMeasurements(0) = state.primaryMass() - state.decayParameters()[TwoBodyDecayParameters::mass];
    theGblExtPrecisions.ResizeTo(1);
    theGblExtPrecisions(0) = 1.0 / (state.primaryWidth() * state.primaryWidth());
    // nominal field
    theNomField = trajectory1.nominalField();
  } else {
    unsigned int nHitMeas1     = trajectory1.numberOfHitMeas();
    unsigned int nVirtualMeas1 = trajectory1.numberOfVirtualMeas(); 
    unsigned int nPar1         = trajectory1.numberOfPar();
    unsigned int nVirtualPar1  = trajectory1.numberOfVirtualPar(); 
     
    // derivatives of the trajectory w.r.t. to the decay parameters
    AlgebraicMatrix fullDeriv1 = trajectory1.derivatives().sub(1,nHitMeas1+nVirtualMeas1,1,nLocal) * trajectory1.localToTrajectory() * deriv.first;

    unsigned int nHitMeas2     = trajectory2.numberOfHitMeas();
    unsigned int nVirtualMeas2 = trajectory2.numberOfVirtualMeas();  
    unsigned int nPar2         = trajectory2.numberOfPar();
    unsigned int nVirtualPar2  = trajectory2.numberOfVirtualPar();

    AlgebraicMatrix fullDeriv2 = trajectory2.derivatives().sub(1,nHitMeas2+nVirtualMeas2,1,nLocal) * trajectory2.localToTrajectory() * deriv.second;

    theNumberOfRecHits.first = recHits.first.size();
    theNumberOfRecHits.second = recHits.second.size();

    theNumberOfHits = trajectory1.numberOfHits() + trajectory2.numberOfHits(); 
    theNumberOfPars = nPar1 + nPar2;
    theNumberOfVirtualPars = nVirtualPar1 + nVirtualPar2;
    theNumberOfVirtualMeas = nVirtualMeas1 + nVirtualMeas2 + 1; // add virtual mass measurement
  
    // hit measurements from trajectory 1
    int rowOffset = 1;
    int colOffset = 1; 
    theDerivatives.sub( rowOffset, colOffset,                fullDeriv1.sub(            1, nHitMeas1,                     1, nTbd ) );
    colOffset += nTbd;
    theDerivatives.sub( rowOffset, colOffset, trajectory1.derivatives().sub(            1, nHitMeas1,            nLocal + 1, nPar1 + nVirtualPar1 ) );
    // hit measurements from trajectory 2
    rowOffset += nHitMeas1;
    colOffset = 1; 
    theDerivatives.sub( rowOffset, colOffset,                fullDeriv2.sub(            1, nHitMeas2,                     1, nTbd ) );
    colOffset += (nPar1 + nVirtualPar1 + nTbd - nLocal);
    theDerivatives.sub( rowOffset, colOffset, trajectory2.derivatives().sub(            1, nHitMeas2,            nLocal + 1, nPar2 + nVirtualPar2 ) );  
    // MS measurements from trajectory 1
    rowOffset += nHitMeas2;  
    colOffset = 1; 
    theDerivatives.sub( rowOffset, colOffset,                fullDeriv1.sub(nHitMeas1 + 1, nHitMeas1 + nVirtualMeas1,          1, nTbd ) );  
    colOffset += nTbd;
    theDerivatives.sub( rowOffset, colOffset, trajectory1.derivatives().sub(nHitMeas1 + 1, nHitMeas1 + nVirtualMeas1, nLocal + 1, nPar1 + nVirtualPar1 ) ); 
    // MS measurements from trajectory 2
    rowOffset += nVirtualMeas1;  
    colOffset = 1; 
    theDerivatives.sub( rowOffset, colOffset,                fullDeriv2.sub(nHitMeas2 + 1, nHitMeas2 + nVirtualMeas2,          1, nTbd ) );
    colOffset += (nPar1 + nVirtualPar1 + nTbd - nLocal);  
    theDerivatives.sub( rowOffset, colOffset, trajectory2.derivatives().sub(nHitMeas2 + 1, nHitMeas2 + nVirtualMeas2, nLocal + 1, nPar2 + nVirtualPar2 ) ); 
        
    theMeasurements.sub(                                         1, trajectory1.measurements().sub(            1, nHitMeas1 ) );
    theMeasurements.sub( nHitMeas1                             + 1, trajectory2.measurements().sub(            1, nHitMeas2 ) );
    theMeasurements.sub( nHitMeas1 + nHitMeas2                 + 1, trajectory1.measurements().sub(nHitMeas1 + 1, nHitMeas1 + nVirtualMeas1 ) );
    theMeasurements.sub( nHitMeas1 + nHitMeas2 + nVirtualMeas1 + 1, trajectory2.measurements().sub(nHitMeas2 + 1, nHitMeas2 + nVirtualMeas2 ) );

    theMeasurementsCov.sub(                                         1, trajectory1.measurementErrors().sub(            1, nHitMeas1 ) );
    theMeasurementsCov.sub( nHitMeas1                             + 1, trajectory2.measurementErrors().sub(            1, nHitMeas2 ) );
    theMeasurementsCov.sub( nHitMeas1 + nHitMeas2                 + 1, trajectory1.measurementErrors().sub(nHitMeas1 + 1, nHitMeas1 + nVirtualMeas1 ) );
    theMeasurementsCov.sub( nHitMeas1 + nHitMeas2 + nVirtualMeas1 + 1, trajectory2.measurementErrors().sub(nHitMeas2 + 1, nHitMeas2 + nVirtualMeas2 ) );

    theTrajectoryPositions.sub(             1, trajectory1.trajectoryPositions() );
    theTrajectoryPositions.sub( nHitMeas1 + 1, trajectory2.trajectoryPositions() );

    theTrajectoryPositionCov = state.decayParameters().covariance().similarity( theDerivatives.sub(1, nHitMeas1 + nHitMeas2, 1, 9) );

    theParameters = state.decayParameters().parameters();

    // add virtual mass measurement
    rowOffset += nVirtualMeas2;
    int indMass = rowOffset-1;
    theMeasurements[indMass] = state.primaryMass() - state.decayParameters()[TwoBodyDecayParameters::mass];
    theMeasurementsCov[indMass][indMass] = state.primaryWidth() * state.primaryWidth();
    theDerivatives[indMass][TwoBodyDecayParameters::mass] = 1.0;
  }

  theRecHits.insert( theRecHits.end(), recHits.first.begin(), recHits.first.end() );
  theRecHits.insert( theRecHits.end(), recHits.second.begin(), recHits.second.end() );

  if ( constructTsosWithErrors )
  {
    constructTsosVecWithErrors( trajectory1, trajectory2, field );
  }
  else
  {
    theTsosVec.insert( theTsosVec.end(),
               trajectory1.trajectoryStates().begin(),
               trajectory1.trajectoryStates().end() );

    theTsosVec.insert( theTsosVec.end(),
               trajectory2.trajectoryStates().begin(),
               trajectory2.trajectoryStates().end() );
  }

  return true;
}

void TwoBodyDecayTrajectory::constructSingleTsosWithErrors ( const TrajectoryStateOnSurface &  tsos, int  iTsos, const MagneticField *  field  )

private

Definition at line 251 of file TwoBodyDecayTrajectory.cc.

References alignCSCRings::e, TrajectoryStateOnSurface::localParameters(), LocalTrajectoryParameters::momentum(), ReferenceTrajectoryBase::nMeasPerHit, AlCaHLTBitMon_ParallelJobs::p, LocalTrajectoryParameters::signedInverseMomentum(), TrajectoryStateOnSurface::surface(), TrajectoryStateOnSurface::surfaceSide(), ReferenceTrajectoryBase::theTrajectoryPositionCov, ReferenceTrajectoryBase::theTsosVec, PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by constructTsosVecWithErrors().

{
  AlgebraicSymMatrix55 localErrors;
  const double coeff = 1e-2;

  double invP = tsos.localParameters().signedInverseMomentum();
  LocalVector p = tsos.localParameters().momentum();

  // rough estimate for the errors of q/p, dx/dz and dy/dz, assuming that
  // sigma(px) = sigma(py) = sigma(pz) = coeff*p.
  float dpinv = coeff*( fabs(p.x()) + fabs(p.y()) + fabs(p.z()) )*invP*invP;
  float dxdir = coeff*( fabs(p.x()) + fabs(p.z()) )/p.z()/p.z();
  float dydir = coeff*( fabs(p.y()) + fabs(p.z()) )/p.z()/p.z();
  localErrors[0][0] = dpinv*dpinv;
  localErrors[1][1] = dxdir*dxdir;
  localErrors[2][2] = dydir*dydir;

  // exact values for the errors on local x and y
  localErrors[3][3] = theTrajectoryPositionCov[nMeasPerHit*iTsos][nMeasPerHit*iTsos];
  localErrors[3][4] = theTrajectoryPositionCov[nMeasPerHit*iTsos][nMeasPerHit*iTsos+1];
  localErrors[4][4] = theTrajectoryPositionCov[nMeasPerHit*iTsos+1][nMeasPerHit*iTsos+1];

  // construct tsos with local errors
  theTsosVec[iTsos] =  TrajectoryStateOnSurface( tsos.localParameters(),
                         LocalTrajectoryError( localErrors ),
                         tsos.surface(),
                         field,
                         tsos.surfaceSide() );
}

void TwoBodyDecayTrajectory::constructTsosVecWithErrors ( const ReferenceTrajectory &  traj1, const ReferenceTrajectory &  traj2, const MagneticField *  field  )

private

Definition at line 229 of file TwoBodyDecayTrajectory.cc.

References constructSingleTsosWithErrors(), and ReferenceTrajectoryBase::trajectoryStates().

Referenced by construct().

{
  int iTsos = 0;

  std::vector< TrajectoryStateOnSurface >::const_iterator itTsos;

  for ( itTsos = traj1.trajectoryStates().begin(); itTsos != traj1.trajectoryStates().end(); itTsos++ )
  {
    constructSingleTsosWithErrors( *itTsos, iTsos, field );
    iTsos++;
  }

  for ( itTsos = traj2.trajectoryStates().begin(); itTsos != traj2.trajectoryStates().end(); itTsos++ )
  {
    constructSingleTsosWithErrors( *itTsos, iTsos, field );
    iTsos++;
  }
}

const std::pair< int, int > TwoBodyDecayTrajectory::numberOfRecHits ( void  )

inline

Number of RecHits belonging to the first and second track.

Definition at line 41 of file TwoBodyDecayTrajectory.h.

References theNumberOfRecHits.

{ return theNumberOfRecHits; }

Member Data Documentation

std::pair< int, int > TwoBodyDecayTrajectory::theNumberOfRecHits

private

Definition at line 62 of file TwoBodyDecayTrajectory.h.

Referenced by construct(), and numberOfRecHits().