File indexing completed on 2024-04-06 12:27:48
0001 #include "FWCore/Utilities/interface/Exception.h"
0002
0003 #include "DataFormats/MuonReco/interface/Muon.h"
0004 #include "DataFormats/MuonReco/interface/MuonFwd.h"
0005 #include "DataFormats/MuonReco/interface/MuonChamberMatch.h"
0006 #include "DataFormats/MuonReco/interface/MuonSelectors.h"
0007 #include "DataFormats/TauReco/interface/PFTau.h"
0008 #include "DataFormats/TrackReco/interface/HitPattern.h"
0009 #include "DataFormats/TrackReco/interface/Track.h"
0010 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0011 #include "DataFormats/Math/interface/deltaR.h"
0012
0013 #include "RecoTauTag/RecoTau/plugins/PFRecoTauDiscriminationAgainstMuon2Helper.h"
0014 #include "RecoTauTag/RecoTau/interface/RecoTauMuonTools.h"
0015 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0016
0017 #include <iostream>
0018
0019 using reco::tau::format_vint;
0020
0021 PFRecoTauDiscriminationAgainstMuon2Helper::PFRecoTauDiscriminationAgainstMuon2Helper(
0022 const bool& verbosity,
0023 const std::string& moduleLabel,
0024 const bool srcMuons_label_empty,
0025 const double& minPtMatchedMuon,
0026 const double& dRmuonMatch,
0027 const bool& dRmuonMatchLimitedToJetArea,
0028 std::atomic<unsigned int>& numWarnings,
0029 const unsigned int& maxWarnings,
0030 const std::vector<int>& maskMatchesDT,
0031 const std::vector<int>& maskMatchesCSC,
0032 const std::vector<int>& maskMatchesRPC,
0033 const std::vector<int>& maskHitsDT,
0034 const std::vector<int>& maskHitsCSC,
0035 const std::vector<int>& maskHitsRPC,
0036 const edm::Handle<reco::MuonCollection>& muons,
0037 const reco::PFTauRef& pfTau,
0038 const reco::PFCandidatePtr& pfCand)
0039 : pfLeadChargedHadron_{pfCand} {
0040 if (verbosity) {
0041 edm::LogPrint("PFTauAgainstMuon2") << "<PFRecoTauDiscriminationAgainstMuon2Container::discriminate>:";
0042 edm::LogPrint("PFTauAgainstMuon2") << " moduleLabel = " << moduleLabel;
0043 edm::LogPrint("PFTauAgainstMuon2") << "tau #" << pfTau.key() << ": Pt = " << pfTau->pt()
0044 << ", eta = " << pfTau->eta() << ", phi = " << pfTau->phi();
0045 }
0046
0047 std::vector<int> numMatchesDT(4);
0048 std::vector<int> numMatchesCSC(4);
0049 std::vector<int> numMatchesRPC(4);
0050 std::vector<int> numHitsDT(4);
0051 std::vector<int> numHitsCSC(4);
0052 std::vector<int> numHitsRPC(4);
0053 for (int iStation = 0; iStation < 4; ++iStation) {
0054 numMatchesDT[iStation] = 0;
0055 numMatchesCSC[iStation] = 0;
0056 numMatchesRPC[iStation] = 0;
0057 numHitsDT[iStation] = 0;
0058 numHitsCSC[iStation] = 0;
0059 numHitsRPC[iStation] = 0;
0060 }
0061
0062
0063 if (pfLeadChargedHadron_.isNonnull()) {
0064 reco::MuonRef muonRef = pfLeadChargedHadron_->muonRef();
0065 if (muonRef.isNonnull()) {
0066 if (verbosity)
0067 edm::LogPrint("PFTauAgainstMuon2") << " has muonRef.";
0068 reco::tau::countMatches(*muonRef, numMatchesDT, numMatchesCSC, numMatchesRPC);
0069 reco::tau::countHits(*muonRef, numHitsDT, numHitsCSC, numHitsRPC);
0070 }
0071 }
0072
0073 if (!srcMuons_label_empty) {
0074 size_t numMuons = muons->size();
0075 for (size_t idxMuon = 0; idxMuon < numMuons; ++idxMuon) {
0076 reco::MuonRef muon(muons, idxMuon);
0077 if (verbosity)
0078 edm::LogPrint("PFTauAgainstMuon2") << "muon #" << muon.key() << ": Pt = " << muon->pt()
0079 << ", eta = " << muon->eta() << ", phi = " << muon->phi();
0080 if (!(muon->pt() > minPtMatchedMuon)) {
0081 if (verbosity) {
0082 edm::LogPrint("PFTauAgainstMuon2") << " fails Pt cut --> skipping it.";
0083 }
0084 continue;
0085 }
0086 if (pfLeadChargedHadron_.isNonnull()) {
0087 reco::MuonRef muonRef = pfLeadChargedHadron_->muonRef();
0088 if (muonRef.isNonnull() && muon == pfLeadChargedHadron_->muonRef()) {
0089 if (verbosity) {
0090 edm::LogPrint("PFTauAgainstMuon2") << " matches muonRef of tau --> skipping it.";
0091 }
0092 continue;
0093 }
0094 }
0095 double dR = deltaR(muon->p4(), pfTau->p4());
0096 double dRmatch = dRmuonMatch;
0097 if (dRmuonMatchLimitedToJetArea) {
0098 double jetArea = 0.;
0099 if (pfTau->jetRef().isNonnull())
0100 jetArea = pfTau->jetRef()->jetArea();
0101 if (jetArea > 0.) {
0102 dRmatch = std::min(dRmatch, std::sqrt(jetArea / M_PI));
0103 } else {
0104 if (numWarnings < maxWarnings) {
0105 edm::LogInfo("PFRecoTauDiscriminationAgainstMuon2Container::discriminate")
0106 << "Jet associated to Tau: Pt = " << pfTau->pt() << ", eta = " << pfTau->eta()
0107 << ", phi = " << pfTau->phi() << " has area = " << jetArea << " !!";
0108 ++numWarnings;
0109 }
0110 dRmatch = 0.1;
0111 }
0112 }
0113 if (dR < dRmatch) {
0114 if (verbosity)
0115 edm::LogPrint("PFTauAgainstMuon2") << " overlaps with tau, dR = " << dR;
0116 reco::tau::countMatches(*muon, numMatchesDT, numMatchesCSC, numMatchesRPC);
0117 reco::tau::countHits(*muon, numHitsDT, numHitsCSC, numHitsRPC);
0118 }
0119 }
0120 }
0121
0122 for (int iStation = 0; iStation < 4; ++iStation) {
0123 if (numMatchesDT[iStation] > 0 && !maskMatchesDT[iStation])
0124 ++numStationsWithMatches_;
0125 if (numMatchesCSC[iStation] > 0 && !maskMatchesCSC[iStation])
0126 ++numStationsWithMatches_;
0127 if (numMatchesRPC[iStation] > 0 && !maskMatchesRPC[iStation])
0128 ++numStationsWithMatches_;
0129 }
0130
0131 for (int iStation = 2; iStation < 4; ++iStation) {
0132 if (numHitsDT[iStation] > 0 && !maskHitsDT[iStation])
0133 ++numLast2StationsWithHits_;
0134 if (numHitsCSC[iStation] > 0 && !maskHitsCSC[iStation])
0135 ++numLast2StationsWithHits_;
0136 if (numHitsRPC[iStation] > 0 && !maskHitsRPC[iStation])
0137 ++numLast2StationsWithHits_;
0138 }
0139
0140 if (verbosity) {
0141 edm::LogPrint("PFTauAgainstMuon2") << "numMatchesDT = " << format_vint(numMatchesDT);
0142 edm::LogPrint("PFTauAgainstMuon2") << "numMatchesCSC = " << format_vint(numMatchesCSC);
0143 edm::LogPrint("PFTauAgainstMuon2") << "numMatchesRPC = " << format_vint(numMatchesRPC);
0144 edm::LogPrint("PFTauAgainstMuon2") << " --> numStationsWithMatches_ = " << numStationsWithMatches_;
0145 edm::LogPrint("PFTauAgainstMuon2") << "numHitsDT = " << format_vint(numHitsDT);
0146 edm::LogPrint("PFTauAgainstMuon2") << "numHitsCSC = " << format_vint(numHitsCSC);
0147 edm::LogPrint("PFTauAgainstMuon2") << "numHitsRPC = " << format_vint(numHitsRPC);
0148 edm::LogPrint("PFTauAgainstMuon2") << " --> numLast2StationsWithHits_ = " << numLast2StationsWithHits_;
0149 }
0150
0151 if (pfLeadChargedHadron_.isNonnull()) {
0152 energyECALplusHCAL_ = pfLeadChargedHadron_->ecalEnergy() + pfLeadChargedHadron_->hcalEnergy();
0153 if (verbosity) {
0154 if (pfLeadChargedHadron_->trackRef().isNonnull()) {
0155 edm::LogPrint("PFTauAgainstMuon2")
0156 << "decayMode = " << pfTau->decayMode() << ", energy(ECAL+HCAL) = " << energyECALplusHCAL_
0157 << ", leadPFChargedHadronP = " << pfLeadChargedHadron_->trackRef()->p();
0158 } else if (pfLeadChargedHadron_->gsfTrackRef().isNonnull()) {
0159 edm::LogPrint("PFTauAgainstMuon2")
0160 << "decayMode = " << pfTau->decayMode() << ", energy(ECAL+HCAL) = " << energyECALplusHCAL_
0161 << ", leadPFChargedHadronP = " << pfLeadChargedHadron_->gsfTrackRef()->p();
0162 }
0163 }
0164 if (pfLeadChargedHadron_->trackRef().isNonnull())
0165 leadTrack_ = pfLeadChargedHadron_->trackRef().get();
0166 else if (pfLeadChargedHadron_->gsfTrackRef().isNonnull())
0167 leadTrack_ = pfLeadChargedHadron_->gsfTrackRef().get();
0168 }
0169 }
0170
0171 bool PFRecoTauDiscriminationAgainstMuon2Helper::eval(const PFRecoTauDiscriminationAgainstMuonConfigSet& config,
0172 const reco::PFTauRef& pfTau) const {
0173 bool passesCaloMuonVeto = true;
0174 if (pfLeadChargedHadron_.isNonnull()) {
0175 if (pfTau->decayMode() == 0 && leadTrack_ && energyECALplusHCAL_ < (config.hop * leadTrack_->p()))
0176 passesCaloMuonVeto = false;
0177 }
0178
0179 bool discriminatorValue = false;
0180 if (config.discriminatorOption == PFRecoTauDiscriminationAgainstMuonConfigSet::kLoose &&
0181 numStationsWithMatches_ <= config.maxNumberOfMatches)
0182 discriminatorValue = true;
0183 else if (config.discriminatorOption == PFRecoTauDiscriminationAgainstMuonConfigSet::kMedium &&
0184 numStationsWithMatches_ <= config.maxNumberOfMatches &&
0185 numLast2StationsWithHits_ <= config.maxNumberOfHitsLast2Stations)
0186 discriminatorValue = true;
0187 else if (config.discriminatorOption == PFRecoTauDiscriminationAgainstMuonConfigSet::kTight &&
0188 numStationsWithMatches_ <= config.maxNumberOfMatches &&
0189 numLast2StationsWithHits_ <= config.maxNumberOfHitsLast2Stations && passesCaloMuonVeto)
0190 discriminatorValue = true;
0191 else if (config.discriminatorOption == PFRecoTauDiscriminationAgainstMuonConfigSet::kCustom) {
0192 discriminatorValue = true;
0193 if (config.maxNumberOfMatches >= 0 && numStationsWithMatches_ > config.maxNumberOfMatches)
0194 discriminatorValue = false;
0195 if (config.maxNumberOfHitsLast2Stations >= 0 && numLast2StationsWithHits_ > config.maxNumberOfHitsLast2Stations)
0196 discriminatorValue = false;
0197 if (config.doCaloMuonVeto && !passesCaloMuonVeto)
0198 discriminatorValue = false;
0199 }
0200 return discriminatorValue;
0201 }
