File indexing completed on 2024-06-06 04:27:07
0001 #include "PhysicsTools/PatAlgos/interface/SoftMuonMvaRun3Estimator.h"
0002 #include "DataFormats/PatCandidates/interface/Muon.h"
0003 #include "PhysicsTools/XGBoost/interface/XGBooster.h"
0004
0005 typedef std::pair<const reco::MuonChamberMatch*, const reco::MuonSegmentMatch*> MatchPair;
0006
0007 const MatchPair& getBetterMatch(const MatchPair& match1, const MatchPair& match2) {
0008
0009
0010
0011 if (match2.first->detector() == MuonSubdetId::DT and match1.first->detector() != MuonSubdetId::DT)
0012 return match2;
0013
0014
0015
0016
0017 if (abs(match1.first->x - match1.second->x) > abs(match2.first->x - match2.second->x))
0018 return match2;
0019
0020 return match1;
0021 }
0022
0023 float dX(const MatchPair& match) {
0024 if (match.first and match.second->hasPhi())
0025 return (match.first->x - match.second->x);
0026 else
0027 return 9999.;
0028 }
0029
0030 float pullX(const MatchPair& match) {
0031 if (match.first and match.second->hasPhi())
0032 return dX(match) / sqrt(pow(match.first->xErr, 2) + pow(match.second->xErr, 2));
0033 else
0034 return 9999.;
0035 }
0036
0037 float pullDxDz(const MatchPair& match) {
0038 if (match.first and match.second->hasPhi())
0039 return (match.first->dXdZ - match.second->dXdZ) /
0040 sqrt(pow(match.first->dXdZErr, 2) + pow(match.second->dXdZErr, 2));
0041 else
0042 return 9999.;
0043 }
0044
0045 float dY(const MatchPair& match) {
0046 if (match.first and match.second->hasZed())
0047 return (match.first->y - match.second->y);
0048 else
0049 return 9999.;
0050 }
0051
0052 float pullY(const MatchPair& match) {
0053 if (match.first and match.second->hasZed())
0054 return dY(match) / sqrt(pow(match.first->yErr, 2) + pow(match.second->yErr, 2));
0055 else
0056 return 9999.;
0057 }
0058
0059 float pullDyDz(const MatchPair& match) {
0060 if (match.first and match.second->hasZed())
0061 return (match.first->dYdZ - match.second->dYdZ) /
0062 sqrt(pow(match.first->dYdZErr, 2) + pow(match.second->dYdZErr, 2));
0063 else
0064 return 9999.;
0065 }
0066
0067 void fillMatchInfoForStation(std::string prefix, pat::XGBooster& booster, const MatchPair& match) {
0068 booster.set(prefix + "_dX", dX(match));
0069 booster.set(prefix + "_pullX", pullX(match));
0070 booster.set(prefix + "_pullDxDz", pullDxDz(match));
0071 booster.set(prefix + "_dY", dY(match));
0072 booster.set(prefix + "_pullY", pullY(match));
0073 booster.set(prefix + "_pullDyDz", pullDyDz(match));
0074 }
0075
0076 void fillMatchInfo(pat::XGBooster& booster, const pat::Muon& muon) {
0077
0078 const int n_stations = 2;
0079 std::vector<MatchPair> matches;
0080 for (unsigned int i = 0; i < n_stations; ++i)
0081 matches.push_back(std::pair(nullptr, nullptr));
0082
0083
0084 for (auto& chamberMatch : muon.matches()) {
0085 unsigned int station = chamberMatch.station() - 1;
0086 if (station >= n_stations)
0087 continue;
0088
0089
0090
0091
0092 for (auto& segmentMatch : chamberMatch.segmentMatches) {
0093 if (not segmentMatch.isMask(reco::MuonSegmentMatch::BestInStationByDR) ||
0094 not segmentMatch.isMask(reco::MuonSegmentMatch::BelongsToTrackByDR))
0095 continue;
0096
0097
0098
0099
0100 auto match_pair = MatchPair(&chamberMatch, &segmentMatch);
0101
0102 if (matches[station].first)
0103 matches[station] = getBetterMatch(matches[station], match_pair);
0104 else
0105 matches[station] = match_pair;
0106 }
0107 }
0108
0109
0110 fillMatchInfoForStation("match1", booster, matches[0]);
0111 fillMatchInfoForStation("match2", booster, matches[1]);
0112 }
0113
0114 float pat::computeSoftMvaRun3(pat::XGBooster& booster, const pat::Muon& muon) {
0115 if (!muon.isTrackerMuon() && !muon.isGlobalMuon())
0116 return 0;
0117
0118 fillMatchInfo(booster, muon);
0119
0120 booster.set("pt", muon.pt());
0121 booster.set("eta", muon.eta());
0122 booster.set("trkValidFrac", muon.innerTrack()->validFraction());
0123 booster.set("glbTrackProbability", muon.combinedQuality().glbTrackProbability);
0124 booster.set("nLostHitsInner",
0125 muon.innerTrack()->hitPattern().numberOfLostTrackerHits(reco::HitPattern::MISSING_INNER_HITS));
0126 booster.set("nLostHitsOuter",
0127 muon.innerTrack()->hitPattern().numberOfLostTrackerHits(reco::HitPattern::MISSING_OUTER_HITS));
0128 booster.set("trkKink", muon.combinedQuality().trkKink);
0129 booster.set("chi2LocalPosition", muon.combinedQuality().chi2LocalPosition);
0130 booster.set("nPixels", muon.innerTrack()->hitPattern().numberOfValidPixelHits());
0131 booster.set("nValidHits", muon.innerTrack()->hitPattern().numberOfValidTrackerHits());
0132 booster.set("nLostHitsOn", muon.innerTrack()->hitPattern().numberOfLostTrackerHits(reco::HitPattern::TRACK_HITS));
0133 booster.set("glbNormChi2", muon.isGlobalMuon() ? muon.globalTrack()->normalizedChi2() : 9999.);
0134 booster.set("trkLayers", muon.innerTrack()->hitPattern().trackerLayersWithMeasurement());
0135 booster.set("highPurity", muon.innerTrack()->quality(reco::Track::highPurity));
0136
0137 return booster.predict();
0138 }