TwoBodyDecayLinearizationPointFinder Class Reference

#include <TwoBodyDecayLinearizationPointFinder.h>

Public Types
typedef PerigeeLinearizedTrackState::RefCountedLinearizedTrackState	RefCountedLinearizedTrackState

Public Member Functions
virtual TwoBodyDecayLinearizationPointFinder *	clone (void) const
virtual const TwoBodyDecayParameters	getLinearizationPoint (const std::vector< RefCountedLinearizedTrackState > &tracks, const double primaryMass, const double secondaryMass) const
	TwoBodyDecayLinearizationPointFinder (const edm::ParameterSet &config)
virtual	~TwoBodyDecayLinearizationPointFinder (void)

Detailed Description

Class TwoBodyDecayLinearizationPointFinder computes a rough estimate of the parameters of a decay. This serves as linearization point for TwoBodyDecayEstimator.

/author Edmund Widl

Definition at line 16 of file TwoBodyDecayLinearizationPointFinder.h.

Member Typedef Documentation

typedef PerigeeLinearizedTrackState::RefCountedLinearizedTrackState TwoBodyDecayLinearizationPointFinder::RefCountedLinearizedTrackState

Definition at line 21 of file TwoBodyDecayLinearizationPointFinder.h.

Constructor & Destructor Documentation

TwoBodyDecayLinearizationPointFinder::TwoBodyDecayLinearizationPointFinder ( const edm::ParameterSet &  config )

inline

Definition at line 23 of file TwoBodyDecayLinearizationPointFinder.h.

Referenced by clone().

{}

virtual TwoBodyDecayLinearizationPointFinder::~TwoBodyDecayLinearizationPointFinder ( void  )

inlinevirtual

Definition at line 25 of file TwoBodyDecayLinearizationPointFinder.h.

{}

Member Function Documentation

virtual TwoBodyDecayLinearizationPointFinder* TwoBodyDecayLinearizationPointFinder::clone ( void  ) const

inlinevirtual

Definition at line 32 of file TwoBodyDecayLinearizationPointFinder.h.

References TwoBodyDecayLinearizationPointFinder().

{ return new TwoBodyDecayLinearizationPointFinder( *this ); }

const TwoBodyDecayParameters TwoBodyDecayLinearizationPointFinder::getLinearizationPoint ( const std::vector< RefCountedLinearizedTrackState > &  tracks, const double  primaryMass, const double  secondaryMass  ) const

virtual

Definition at line 7 of file TwoBodyDecayLinearizationPointFinder.cc.

References TwoBodyDecayParameters::dimension, TwoBodyDecayParameters::mass, TrajectoryStateClosestToPoint::momentum(), AlCaHLTBitMon_ParallelJobs::p, p1, p2, perp2(), TwoBodyDecayParameters::phi, pi, PerigeeLinearizedTrackState::predictedState(), TwoBodyDecayParameters::px, TwoBodyDecayParameters::py, TwoBodyDecayParameters::pz, TwoBodyDecayModel::rotationMatrix(), mathSSE::sqrt(), TwoBodyDecayParameters::theta, TwoBodyDecayParameters::x, PV3DBase< T, PVType, FrameType >::x(), TwoBodyDecayParameters::y, PV3DBase< T, PVType, FrameType >::y(), TwoBodyDecayParameters::z, and PV3DBase< T, PVType, FrameType >::z().

Referenced by TwoBodyDecayFitter::estimate().

{
  GlobalPoint linPoint = tracks[0]->linearizationPoint();
  PerigeeLinearizedTrackState* linTrack1 = dynamic_cast<PerigeeLinearizedTrackState*>( tracks[0].get() );
  PerigeeLinearizedTrackState* linTrack2 = dynamic_cast<PerigeeLinearizedTrackState*>( tracks[1].get() );
  if (!linTrack1 || !linTrack2) return TwoBodyDecayParameters();

  GlobalVector firstMomentum = linTrack1->predictedState().momentum();
  GlobalVector secondMomentum = linTrack2->predictedState().momentum();

  AlgebraicVector secondaryMomentum1( 3 );
  secondaryMomentum1[0] = firstMomentum.x();
  secondaryMomentum1[1] = firstMomentum.y();
  secondaryMomentum1[2] = firstMomentum.z();

  AlgebraicVector secondaryMomentum2( 3 );
  secondaryMomentum2[0] = secondMomentum.x();
  secondaryMomentum2[1] = secondMomentum.y();
  secondaryMomentum2[2] = secondMomentum.z();

  AlgebraicVector primaryMomentum = secondaryMomentum1 + secondaryMomentum2;

  TwoBodyDecayModel decayModel( primaryMass, secondaryMass );
  AlgebraicMatrix rotMat = decayModel.rotationMatrix( primaryMomentum[0], primaryMomentum[1], primaryMomentum[2] );
  AlgebraicMatrix invRotMat = rotMat.T();

  double p = primaryMomentum.norm();
  double pSquared = p*p;
  double gamma = sqrt( pSquared + primaryMass*primaryMass )/primaryMass;
  double betaGamma = p/primaryMass;
  AlgebraicSymMatrix lorentzTransformation( 4, 1 );
  lorentzTransformation[0][0] = gamma;
  lorentzTransformation[3][3] = gamma;
  lorentzTransformation[0][3] = -betaGamma;

  double p1 = secondaryMomentum1.norm();
  AlgebraicVector boostedLorentzMomentum1( 4 );
  boostedLorentzMomentum1[0] = sqrt( p1*p1 + secondaryMass*secondaryMass );
  boostedLorentzMomentum1.sub( 2, invRotMat*secondaryMomentum1 );

  AlgebraicVector restFrameLorentzMomentum1 = lorentzTransformation*boostedLorentzMomentum1;
  AlgebraicVector restFrameMomentum1 = restFrameLorentzMomentum1.sub( 2, 4 );
  double perp1 = sqrt( restFrameMomentum1[0]*restFrameMomentum1[0] + restFrameMomentum1[1]*restFrameMomentum1[1] );
  double theta1 = atan2( perp1, restFrameMomentum1[2] );
  double phi1 = atan2( restFrameMomentum1[1], restFrameMomentum1[0] );

  double p2 = secondaryMomentum2.norm();
  AlgebraicVector boostedLorentzMomentum2( 4 );
  boostedLorentzMomentum2[0] = sqrt( p2*p2 + secondaryMass*secondaryMass );
  boostedLorentzMomentum2.sub( 2, invRotMat*secondaryMomentum2 );

  AlgebraicVector restFrameLorentzMomentum2 = lorentzTransformation*boostedLorentzMomentum2;
  AlgebraicVector restFrameMomentum2 = restFrameLorentzMomentum2.sub( 2, 4 );
  double perp2 = sqrt( restFrameMomentum2[0]*restFrameMomentum2[0] + restFrameMomentum2[1]*restFrameMomentum2[1] );
  double theta2 = atan2( perp2, restFrameMomentum2[2] );
  double phi2 = atan2( restFrameMomentum2[1], restFrameMomentum2[0] );

  double pi = 3.141592654;
  double relSign = -1.;

  if ( phi1 < 0 ) phi1 += 2*pi;
  if ( phi2 < 0 ) phi2 += 2*pi;
  if ( phi1 > phi2 ) relSign = 1.;

  double momentumSquared1 = secondaryMomentum1.normsq();
  double energy1 = sqrt( secondaryMass*secondaryMass + momentumSquared1 );
  double momentumSquared2 = secondaryMomentum2.normsq();
  double energy2 = sqrt( secondaryMass*secondaryMass + momentumSquared2 );
  double sumMomentaSquared = ( secondaryMomentum1 + secondaryMomentum2 ).normsq();
  double sumEnergiesSquared = ( energy1 + energy2 )*( energy1 + energy2 );
  double estimatedPrimaryMass = sqrt( sumEnergiesSquared - sumMomentaSquared );

  AlgebraicVector linParam( TwoBodyDecayParameters::dimension, 0 );
  linParam[TwoBodyDecayParameters::x] = linPoint.x();
  linParam[TwoBodyDecayParameters::y] = linPoint.y();
  linParam[TwoBodyDecayParameters::z] = linPoint.z();
  linParam[TwoBodyDecayParameters::px] = primaryMomentum[0];
  linParam[TwoBodyDecayParameters::py] = primaryMomentum[1];
  linParam[TwoBodyDecayParameters::pz] = primaryMomentum[2];
  linParam[TwoBodyDecayParameters::theta] = 0.5*( theta1 - theta2 + pi ) ;
  linParam[TwoBodyDecayParameters::phi] = 0.5*( phi1 + phi2 + relSign*pi );
  linParam[TwoBodyDecayParameters::mass] = estimatedPrimaryMass;

  return TwoBodyDecayParameters( linParam );
}