CMS 3D CMS Logo

trackAssociationChi2.cc
Go to the documentation of this file.
5 
6 namespace track_associator {
7  constexpr double invalidChi2 = 10000000000.;
8 
10  const reco::TrackBase::CovarianceMatrix &recoTrackCovMatrix,
11  const Basic3DVector<double> &momAtVtx,
12  const Basic3DVector<double> &vert,
13  int charge,
14  const MagneticField &magfield,
15  const reco::BeamSpot &bs) {
16  double chi2 = invalidChi2;
17 
18  std::pair<bool, reco::TrackBase::ParameterVector> params =
19  reco::trackingParametersAtClosestApproachToBeamSpot(vert, momAtVtx, charge, magfield, bs);
20  if (params.first) {
21  reco::TrackBase::ParameterVector sParameters = params.second;
22 
23  reco::TrackBase::ParameterVector diffParameters = rParameters - sParameters;
24  diffParameters[2] = reco::deltaPhi(diffParameters[2], 0.f);
25  chi2 = ROOT::Math::Dot(diffParameters * recoTrackCovMatrix, diffParameters);
26  chi2 /= 5;
27 
28  LogDebug("TrackAssociator") << "====NEW RECO TRACK WITH PT="
29  << sin(rParameters[1]) * float(charge) / rParameters[0] << "====\n"
30  << "qoverp sim: " << sParameters[0] << "\n"
31  << "lambda sim: " << sParameters[1] << "\n"
32  << "phi sim: " << sParameters[2] << "\n"
33  << "dxy sim: " << sParameters[3] << "\n"
34  << "dsz sim: " << sParameters[4] << "\n"
35  << ": " /*<< */ << "\n"
36  << "qoverp rec: " << rParameters[0] << "\n"
37  << "lambda rec: " << rParameters[1] << "\n"
38  << "phi rec: " << rParameters[2] << "\n"
39  << "dxy rec: " << rParameters[3] << "\n"
40  << "dsz rec: " << rParameters[4] << "\n"
41  << ": " /*<< */ << "\n"
42  << "chi2: " << chi2 << "\n";
43  }
44  return chi2;
45  }
46 
48  const reco::TrackBase::CovarianceMatrix &recoTrackCovMatrix,
49  const TrackingParticle &trackingParticle,
50  const MagneticField &magfield,
51  const reco::BeamSpot &bs) {
52  const int charge = trackingParticle.charge();
53  if (charge == 0)
54  return invalidChi2;
55 
56  const auto tpMom = trackingParticle.momentum();
57  Basic3DVector<double> momAtVtx(tpMom.x(), tpMom.y(), tpMom.z());
58  Basic3DVector<double> vert(trackingParticle.vertex());
59 
60  return trackAssociationChi2(rParameters, recoTrackCovMatrix, momAtVtx, vert, charge, magfield, bs);
61  }
62 
64  const TrackingParticle &trackingParticle,
65  const MagneticField &magfield,
66  const reco::BeamSpot &bs) {
67  auto cov = track.covariance();
68  cov.Invert();
69 
70  return trackAssociationChi2(track.parameters(), cov, trackingParticle, magfield, bs);
71  }
72 } // namespace track_associator
#define LogDebug(id)
constexpr double deltaPhi(double phi1, double phi2)
Definition: deltaPhi.h:22
Vector momentum() const
spatial momentum vector
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
math::Vector< dimension >::type ParameterVector
parameter vector
Definition: TrackBase.h:74
float charge() const
Electric charge. Note this is taken from the first SimTrack only.
CovarianceMatrix covariance() const
return track covariance matrix
Definition: TrackBase.h:782
double f[11][100]
std::pair< bool, reco::TrackBase::ParameterVector > trackingParametersAtClosestApproachToBeamSpot(const Basic3DVector< double > &vertex, const Basic3DVector< double > &momAtVtx, float charge, const MagneticField &magField, const BeamSpot &bs)
Point vertex() const
Parent vertex position.
ParameterVector parameters() const
Track parameters with one-to-one correspondence to the covariance matrix.
Definition: TrackBase.h:776
Monte Carlo truth information used for tracking validation.
double trackAssociationChi2(const reco::TrackBase::ParameterVector &rParameters, const reco::TrackBase::CovarianceMatrix &recoTrackCovMatrix, const Basic3DVector< double > &momAtVtx, const Basic3DVector< double > &vert, int charge, const MagneticField &magfield, const reco::BeamSpot &bs)
basic method where chi2 is computed
#define constexpr
math::Error< dimension >::type CovarianceMatrix
5 parameter covariance matrix
Definition: TrackBase.h:77