File indexing completed on 2023-03-17 11:25:44
0001 #include "DataFormats/Math/interface/deltaPhi.h"
0002 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0003 #include "SimTracker/TrackAssociation/interface/trackAssociationChi2.h"
0004 #include "TrackingTools/PatternTools/interface/trackingParametersAtClosestApproachToBeamSpot.h"
0005
0006 namespace track_associator {
0007 constexpr double invalidChi2 = 10000000000.;
0008
0009 double trackAssociationChi2(const reco::TrackBase::ParameterVector &rParameters,
0010 const reco::TrackBase::CovarianceMatrix &recoTrackCovMatrix,
0011 const reco::TrackBase::ParameterVector &sParameters) {
0012 double chi2 = invalidChi2;
0013
0014 reco::TrackBase::ParameterVector diffParameters = rParameters - sParameters;
0015 diffParameters[2] = reco::deltaPhi(diffParameters[2], 0.f);
0016 chi2 = ROOT::Math::Dot(diffParameters * recoTrackCovMatrix, diffParameters);
0017 chi2 /= 5;
0018
0019 LogDebug("TrackAssociator") << "====NEW RECO TRACK WITH PT=" << sin(rParameters[1]) / rParameters[0] << "====\n"
0020 << "qoverp sim: " << sParameters[0] << "\n"
0021 << "lambda sim: " << sParameters[1] << "\n"
0022 << "phi sim: " << sParameters[2] << "\n"
0023 << "dxy sim: " << sParameters[3] << "\n"
0024 << "dsz sim: " << sParameters[4] << "\n"
0025 << ": " << "\n"
0026 << "qoverp rec: " << rParameters[0] << "\n"
0027 << "lambda rec: " << rParameters[1] << "\n"
0028 << "phi rec: " << rParameters[2] << "\n"
0029 << "dxy rec: " << rParameters[3] << "\n"
0030 << "dsz rec: " << rParameters[4] << "\n"
0031 << ": " << "\n"
0032 << "chi2: " << chi2 << "\n";
0033 return chi2;
0034 }
0035
0036 double trackAssociationChi2(const reco::TrackBase::ParameterVector &rParameters,
0037 const reco::TrackBase::CovarianceMatrix &recoTrackCovMatrix,
0038 const Basic3DVector<double> &momAtVtx,
0039 const Basic3DVector<double> &vert,
0040 int charge,
0041 const MagneticField &magfield,
0042 const reco::BeamSpot &bs) {
0043 double chi2 = invalidChi2;
0044
0045 std::pair<bool, reco::TrackBase::ParameterVector> params =
0046 reco::trackingParametersAtClosestApproachToBeamSpot(vert, momAtVtx, charge, magfield, bs);
0047 if (params.first) {
0048 return trackAssociationChi2(rParameters, recoTrackCovMatrix, params.second);
0049 }
0050 return chi2;
0051 }
0052
0053 double trackAssociationChi2(const reco::TrackBase::ParameterVector &rParameters,
0054 const reco::TrackBase::CovarianceMatrix &recoTrackCovMatrix,
0055 const TrackingParticle &trackingParticle,
0056 const MagneticField &magfield,
0057 const reco::BeamSpot &bs) {
0058 const int charge = trackingParticle.charge();
0059 if (charge == 0)
0060 return invalidChi2;
0061
0062 const auto tpMom = trackingParticle.momentum();
0063 Basic3DVector<double> momAtVtx(tpMom.x(), tpMom.y(), tpMom.z());
0064 Basic3DVector<double> vert(trackingParticle.vertex());
0065
0066 return trackAssociationChi2(rParameters, recoTrackCovMatrix, momAtVtx, vert, charge, magfield, bs);
0067 }
0068
0069 double trackAssociationChi2(const reco::TrackBase &track,
0070 const TrackingParticle &trackingParticle,
0071 const MagneticField &magfield,
0072 const reco::BeamSpot &bs) {
0073 auto cov = track.covariance();
0074 cov.Invert();
0075
0076 return trackAssociationChi2(track.parameters(), cov, trackingParticle, magfield, bs);
0077 }
0078 }