Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-09-18 02:24:57

0001 /**
0002  *  @file     L1TTauOffline.cc
0003  *  @authors  Olivier Davignon (University of Bristol), C├ęcile Caillol (University of Wisconsin - Madison)
0004  *  @date     24/05/2017
0005  *  @version  1.1
0006  *
0007  */
0008 
0009 #include "DQMOffline/L1Trigger/interface/L1TTauOffline.h"
0010 
0011 #include "DataFormats/MuonReco/interface/MuonSelectors.h"
0012 #include "DataFormats/MuonReco/interface/MuonPFIsolation.h"
0013 #include "DataFormats/Math/interface/deltaR.h"
0014 
0015 #include "TLorentzVector.h"
0016 
0017 #include <iomanip>
0018 #include <cstdio>
0019 #include <string>
0020 #include <sstream>
0021 #include <cmath>
0022 #include <algorithm>
0023 
0024 using namespace reco;
0025 using namespace trigger;
0026 using namespace edm;
0027 using namespace std;
0028 
0029 TauL1TPair::TauL1TPair(const TauL1TPair& tauL1tPair) {
0030   m_tau = tauL1tPair.m_tau;
0031   m_regTau = tauL1tPair.m_regTau;
0032 
0033   m_eta = tauL1tPair.m_eta;
0034   m_phi_bar = tauL1tPair.m_phi_bar;
0035   m_phi_end = tauL1tPair.m_phi_end;
0036 }
0037 
0038 double TauL1TPair::dR() { return deltaR(m_regTau->eta(), m_regTau->phi(), eta(), phi()); }
0039 const std::map<std::string, unsigned int> L1TTauOffline::PlotConfigNames = {
0040     {"nVertex", PlotConfig::nVertex}, {"ETvsET", PlotConfig::ETvsET}, {"PHIvsPHI", PlotConfig::PHIvsPHI}};
0041 
0042 //
0043 // -------------------------------------- Constructor --------------------------------------------
0044 //
0045 L1TTauOffline::L1TTauOffline(const edm::ParameterSet& ps)
0046     : theTauCollection_(consumes<reco::PFTauCollection>(ps.getUntrackedParameter<edm::InputTag>("tauInputTag"))),
0047       AntiMuInputTag_(
0048           consumes<reco::TauDiscriminatorContainer>(ps.getUntrackedParameter<edm::InputTag>("antiMuInputTag"))),
0049       AntiMuWP_(ps.getUntrackedParameter<std::string>("antiMuWP")),
0050       AntiEleInputTag_(
0051           consumes<reco::TauDiscriminatorContainer>(ps.getUntrackedParameter<edm::InputTag>("antiEleInputTag"))),
0052       AntiEleWP_(ps.getUntrackedParameter<std::string>("antiEleWP")),
0053       DecayModeFindingInputTag_(
0054           consumes<reco::PFTauDiscriminator>(ps.getUntrackedParameter<edm::InputTag>("decayModeFindingInputTag"))),
0055       comb3TInputTag_(
0056           consumes<reco::TauDiscriminatorContainer>(ps.getUntrackedParameter<edm::InputTag>("comb3TInputTag"))),
0057       comb3TWP_(ps.getUntrackedParameter<std::string>("comb3TWP")),
0058       MuonInputTag_(consumes<reco::MuonCollection>(ps.getUntrackedParameter<edm::InputTag>("muonInputTag"))),
0059       MetInputTag_(consumes<reco::PFMETCollection>(ps.getUntrackedParameter<edm::InputTag>("metInputTag"))),
0060       VtxInputTag_(consumes<reco::VertexCollection>(ps.getUntrackedParameter<edm::InputTag>("vtxInputTag"))),
0061       BsInputTag_(consumes<reco::BeamSpot>(ps.getUntrackedParameter<edm::InputTag>("bsInputTag"))),
0062       triggerEvent_(consumes<trigger::TriggerEvent>(ps.getUntrackedParameter<edm::InputTag>("trigInputTag"))),
0063       trigProcess_(ps.getUntrackedParameter<string>("trigProcess")),
0064       triggerResults_(consumes<edm::TriggerResults>(ps.getUntrackedParameter<edm::InputTag>("trigProcess_token"))),
0065       triggerPath_(ps.getUntrackedParameter<vector<std::string>>("triggerNames")),
0066       histFolder_(ps.getParameter<std::string>("histFolder")),
0067       efficiencyFolder_(histFolder_ + "/efficiency_raw"),
0068       stage2CaloLayer2TauToken_(consumes<l1t::TauBxCollection>(ps.getUntrackedParameter<edm::InputTag>("l1tInputTag"))),
0069       tauEfficiencyThresholds_(ps.getParameter<std::vector<int>>("tauEfficiencyThresholds")),
0070       tauEfficiencyBins_(ps.getParameter<std::vector<double>>("tauEfficiencyBins")),
0071       histDefinitions_(dqmoffline::l1t::readHistDefinitions(ps.getParameterSet("histDefinitions"), PlotConfigNames)),
0072       m_TightMuons(),
0073       m_ProbeTaus(),
0074       m_TauL1tPairs(),
0075       m_RecoTaus(),
0076       m_L1tTaus(),
0077       m_RecoRecoTaus(),
0078       m_L1tL1tTaus(),
0079       m_L1tPtCuts(),
0080       m_MaxTauEta(99999),
0081       m_MaxL1tTauDR(99999),
0082       m_MaxHltTauDR(99999),
0083       m_trigIndices(),
0084       h_nVertex_(),
0085       h_tagAndProbeMass_(),
0086       h_L1TauETvsTauET_EB_(),
0087       h_L1TauETvsTauET_EE_(),
0088       h_L1TauETvsTauET_EB_EE_(),
0089       h_L1TauPhivsTauPhi_EB_(),
0090       h_L1TauPhivsTauPhi_EE_(),
0091       h_L1TauPhivsTauPhi_EB_EE_(),
0092       h_L1TauEtavsTauEta_(),
0093       h_resolutionTauET_EB_(),
0094       h_resolutionTauET_EE_(),
0095       h_resolutionTauET_EB_EE_(),
0096       h_resolutionTauPhi_EB_(),
0097       h_resolutionTauPhi_EE_(),
0098       h_resolutionTauPhi_EB_EE_(),
0099       h_resolutionTauEta_(),
0100       h_efficiencyIsoTauET_EB_pass_(),
0101       h_efficiencyIsoTauET_EE_pass_(),
0102       h_efficiencyIsoTauET_EB_EE_pass_(),
0103       h_efficiencyNonIsoTauET_EB_pass_(),
0104       h_efficiencyNonIsoTauET_EE_pass_(),
0105       h_efficiencyNonIsoTauET_EB_EE_pass_(),
0106       h_efficiencyIsoTauET_EB_total_(),
0107       h_efficiencyIsoTauET_EE_total_(),
0108       h_efficiencyIsoTauET_EB_EE_total_(),
0109       h_efficiencyNonIsoTauET_EB_total_(),
0110       h_efficiencyNonIsoTauET_EE_total_(),
0111       h_efficiencyNonIsoTauET_EB_EE_total_() {
0112   edm::LogInfo("L1TTauOffline") << "Constructor "
0113                                 << "L1TTauOffline::L1TTauOffline " << std::endl;
0114   mFieldToken_ = esConsumes();
0115 }
0116 
0117 //
0118 // -- Destructor
0119 //
0120 L1TTauOffline::~L1TTauOffline() {
0121   edm::LogInfo("L1TTauOffline") << "Destructor L1TTauOffline::~L1TTauOffline " << std::endl;
0122 }
0123 
0124 //
0125 // -------------------------------------- beginRun --------------------------------------------
0126 //
0127 void L1TTauOffline::dqmBeginRun(const edm::Run& run, const edm::EventSetup& iSetup)
0128 // void L1TTauOffline::dqmBeginRun(edm::Run const &, edm::EventSetup const &)
0129 {
0130   bool changed = true;
0131   m_hltConfig.init(run, iSetup, trigProcess_, changed);
0132 
0133   edm::LogInfo("L1TTauOffline") << "L1TTauOffline::beginRun" << std::endl;
0134 }
0135 
0136 //
0137 // -------------------------------------- bookHistos --------------------------------------------
0138 //
0139 void L1TTauOffline::bookHistograms(DQMStore::IBooker& ibooker, edm::Run const&, edm::EventSetup const&) {
0140   edm::LogInfo("L1TTauOffline") << "L1TTauOffline::bookHistograms" << std::endl;
0141 
0142   // book at beginRun
0143   bookTauHistos(ibooker);
0144 
0145   for (auto trigNamesIt = triggerPath_.begin(); trigNamesIt != triggerPath_.end(); trigNamesIt++) {
0146     std::string tNameTmp = (*trigNamesIt);
0147     std::string tNamePattern = "";
0148     std::size_t found0 = tNameTmp.find('*');
0149     if (found0 != std::string::npos)
0150       tNamePattern = tNameTmp.substr(0, tNameTmp.size() - 1);
0151     else
0152       tNamePattern = tNameTmp;
0153 
0154     int tIndex = -1;
0155 
0156     for (unsigned ipath = 0; ipath < m_hltConfig.size(); ++ipath) {
0157       std::string tmpName = m_hltConfig.triggerName(ipath);
0158 
0159       std::size_t found = tmpName.find(tNamePattern);
0160       if (found != std::string::npos) {
0161         tIndex = int(ipath);
0162         m_trigIndices.push_back(tIndex);
0163       }
0164     }
0165   }
0166 }
0167 
0168 //
0169 // -------------------------------------- Analyze --------------------------------------------
0170 //
0171 void L1TTauOffline::analyze(edm::Event const& e, edm::EventSetup const& eSetup) {
0172   m_MaxTauEta = 2.1;
0173   m_MaxL1tTauDR = 0.5;
0174   m_MaxHltTauDR = 0.5;
0175 
0176   edm::Handle<reco::PFTauCollection> taus;
0177   e.getByToken(theTauCollection_, taus);
0178 
0179   if (!taus.isValid()) {
0180     edm::LogWarning("L1TTauOffline") << "invalid collection: reco::PFTauCollection " << std::endl;
0181     return;
0182   }
0183 
0184   edm::Handle<reco::MuonCollection> muons;
0185   e.getByToken(MuonInputTag_, muons);
0186 
0187   if (!muons.isValid()) {
0188     edm::LogWarning("L1TTauOffline") << "invalid collection: reco::MuonCollection " << std::endl;
0189     return;
0190   }
0191 
0192   edm::Handle<reco::BeamSpot> beamSpot;
0193   e.getByToken(BsInputTag_, beamSpot);
0194 
0195   if (!beamSpot.isValid()) {
0196     edm::LogWarning("L1TTauOffline") << "invalid collection: reco::BeamSpot " << std::endl;
0197     return;
0198   }
0199 
0200   edm::Handle<reco::VertexCollection> vertex;
0201   e.getByToken(VtxInputTag_, vertex);
0202 
0203   if (!vertex.isValid()) {
0204     edm::LogWarning("L1TTauOffline") << "invalid collection: reco::VertexCollection " << std::endl;
0205     return;
0206   }
0207 
0208   edm::Handle<l1t::TauBxCollection> l1tCands;
0209   e.getByToken(stage2CaloLayer2TauToken_, l1tCands);
0210 
0211   if (!l1tCands.isValid()) {
0212     edm::LogWarning("L1TTauOffline") << "invalid collection: l1t::TauBxCollection " << std::endl;
0213     return;
0214   }
0215 
0216   edm::Handle<edm::TriggerResults> trigResults;
0217   e.getByToken(triggerResults_, trigResults);
0218 
0219   if (!trigResults.isValid()) {
0220     edm::LogWarning("L1TTauOffline") << "invalid collection: edm::TriggerResults " << std::endl;
0221     return;
0222   }
0223 
0224   edm::Handle<trigger::TriggerEvent> trigEvent;
0225   e.getByToken(triggerEvent_, trigEvent);
0226 
0227   if (!trigEvent.isValid()) {
0228     edm::LogWarning("L1TTauOffline") << "invalid collection: trigger::TriggerEvent " << std::endl;
0229     return;
0230   }
0231 
0232   edm::Handle<reco::PFMETCollection> mets;
0233   e.getByToken(MetInputTag_, mets);
0234 
0235   if (!mets.isValid()) {
0236     edm::LogWarning("L1TTauOffline") << "invalid collection: reco::PFMETCollection " << std::endl;
0237     return;
0238   }
0239 
0240   m_BField = eSetup.getHandle(mFieldToken_);
0241 
0242   const reco::Vertex primaryVertex = getPrimaryVertex(vertex, beamSpot);
0243 
0244   getTightMuons(muons, mets, primaryVertex, trigEvent);
0245   getProbeTaus(e, taus, muons, primaryVertex);
0246   getTauL1tPairs(l1tCands);
0247 
0248   vector<l1t::Tau> l1tContainer;
0249   l1tContainer.reserve(l1tCands->size() + 1);
0250 
0251   for (auto tau = l1tCands->begin(0); tau != l1tCands->end(0); ++tau) {
0252     l1tContainer.push_back(*tau);
0253   }
0254 
0255   for (auto tauL1tPairsIt = m_TauL1tPairs.begin(); tauL1tPairsIt != m_TauL1tPairs.end(); ++tauL1tPairsIt) {
0256     float eta = tauL1tPairsIt->eta();
0257     float phi = tauL1tPairsIt->phi();
0258     float pt = tauL1tPairsIt->pt();
0259 
0260     // unmatched gmt cands have l1tPt = -1.
0261     float l1tPt = tauL1tPairsIt->l1tPt();
0262 
0263     int counter = 0;
0264 
0265     for (auto threshold : tauEfficiencyThresholds_) {
0266       std::string str_threshold = std::to_string(threshold);
0267 
0268       int l1tPtCut = threshold;
0269       bool l1tAboveCut = (l1tPt >= l1tPtCut);
0270 
0271       stringstream ptCutToTag;
0272       ptCutToTag << l1tPtCut;
0273       string ptTag = ptCutToTag.str();
0274 
0275       if (fabs(eta) < m_MaxTauEta) {
0276         if (counter == 0) {
0277           if (fabs(eta) < 1.5) {
0278             h_L1TauETvsTauET_EB_->Fill(pt, l1tPt);
0279             h_L1TauPhivsTauPhi_EB_->Fill(phi, tauL1tPairsIt->l1tPhi());
0280             h_resolutionTauET_EB_->Fill((l1tPt - pt) / pt);
0281             h_resolutionTauPhi_EB_->Fill(tauL1tPairsIt->l1tPhi() - phi);
0282           } else {
0283             h_L1TauETvsTauET_EE_->Fill(pt, l1tPt);
0284             h_L1TauPhivsTauPhi_EE_->Fill(phi, tauL1tPairsIt->l1tPhi());
0285             h_resolutionTauET_EE_->Fill((l1tPt - pt) / pt);
0286             h_resolutionTauPhi_EE_->Fill(tauL1tPairsIt->l1tPhi() - phi);
0287           }
0288           h_L1TauETvsTauET_EB_EE_->Fill(pt, l1tPt);
0289           h_L1TauPhivsTauPhi_EB_EE_->Fill(phi, tauL1tPairsIt->l1tPhi());
0290           h_L1TauEtavsTauEta_->Fill(eta, tauL1tPairsIt->l1tEta());
0291           h_resolutionTauET_EB_EE_->Fill((l1tPt - pt) / pt);
0292           h_resolutionTauPhi_EB_EE_->Fill(tauL1tPairsIt->l1tPhi() - phi);
0293           h_resolutionTauEta_->Fill(tauL1tPairsIt->l1tEta() - eta);
0294 
0295           ++counter;
0296         }
0297 
0298         if (fabs(eta) < 1.5) {
0299           h_efficiencyNonIsoTauET_EB_total_[threshold]->Fill(pt);
0300           h_efficiencyIsoTauET_EB_total_[threshold]->Fill(pt);
0301         } else {
0302           h_efficiencyNonIsoTauET_EE_total_[threshold]->Fill(pt);
0303           h_efficiencyIsoTauET_EE_total_[threshold]->Fill(pt);
0304         }
0305         h_efficiencyNonIsoTauET_EB_EE_total_[threshold]->Fill(pt);
0306         h_efficiencyIsoTauET_EB_EE_total_[threshold]->Fill(pt);
0307 
0308         if (l1tAboveCut) {
0309           if (fabs(eta) < 1.5)
0310             h_efficiencyNonIsoTauET_EB_pass_[threshold]->Fill(pt);
0311           else
0312             h_efficiencyNonIsoTauET_EE_pass_[threshold]->Fill(pt);
0313           h_efficiencyNonIsoTauET_EB_EE_pass_[threshold]->Fill(pt);
0314 
0315           if (tauL1tPairsIt->l1tIso() > 0.5) {
0316             if (fabs(eta) < 1.5)
0317               h_efficiencyIsoTauET_EB_pass_[threshold]->Fill(pt);
0318             else
0319               h_efficiencyIsoTauET_EE_pass_[threshold]->Fill(pt);
0320             h_efficiencyIsoTauET_EB_EE_pass_[threshold]->Fill(pt);
0321           }
0322         }
0323       }
0324     }
0325   }  // loop over tau-L1 pairs
0326 }
0327 
0328 //
0329 // -------------------------------------- endRun --------------------------------------------
0330 //
0331 //
0332 // -------------------------------------- book histograms --------------------------------------------
0333 //
0334 void L1TTauOffline::bookTauHistos(DQMStore::IBooker& ibooker) {
0335   ibooker.cd();
0336   ibooker.setCurrentFolder(histFolder_);
0337   dqmoffline::l1t::HistDefinition nVertexDef = histDefinitions_[PlotConfig::nVertex];
0338   h_nVertex_ = ibooker.book1D(nVertexDef.name, nVertexDef.title, nVertexDef.nbinsX, nVertexDef.xmin, nVertexDef.xmax);
0339   h_tagAndProbeMass_ = ibooker.book1D("tagAndProbeMass", "Invariant mass of tag & probe pair", 100, 40, 140);
0340 
0341   dqmoffline::l1t::HistDefinition templateETvsET = histDefinitions_[PlotConfig::ETvsET];
0342   h_L1TauETvsTauET_EB_ = ibooker.book2D("L1TauETvsTauET_EB",
0343                                         "L1 Tau E_{T} vs PFTau E_{T} (EB); PFTau E_{T} (GeV); L1 Tau E_{T} (GeV)",
0344                                         templateETvsET.nbinsX,
0345                                         &templateETvsET.binsX[0],
0346                                         templateETvsET.nbinsY,
0347                                         &templateETvsET.binsY[0]);
0348   h_L1TauETvsTauET_EE_ = ibooker.book2D("L1TauETvsTauET_EE",
0349                                         "L1 Tau E_{T} vs PFTau E_{T} (EE); PFTau E_{T} (GeV); L1 Tau E_{T} (GeV)",
0350                                         templateETvsET.nbinsX,
0351                                         &templateETvsET.binsX[0],
0352                                         templateETvsET.nbinsY,
0353                                         &templateETvsET.binsY[0]);
0354   h_L1TauETvsTauET_EB_EE_ = ibooker.book2D("L1TauETvsTauET_EB_EE",
0355                                            "L1 Tau E_{T} vs PFTau E_{T} (EB+EE); PFTau E_{T} (GeV); L1 Tau E_{T} (GeV)",
0356                                            templateETvsET.nbinsX,
0357                                            &templateETvsET.binsX[0],
0358                                            templateETvsET.nbinsY,
0359                                            &templateETvsET.binsY[0]);
0360 
0361   dqmoffline::l1t::HistDefinition templatePHIvsPHI = histDefinitions_[PlotConfig::PHIvsPHI];
0362   h_L1TauPhivsTauPhi_EB_ =
0363       ibooker.book2D("L1TauPhivsTauPhi_EB",
0364                      "#phi_{tau}^{L1} vs #phi_{tau}^{offline} (EB); #phi_{tau}^{offline}; #phi_{tau}^{L1}",
0365                      templatePHIvsPHI.nbinsX,
0366                      templatePHIvsPHI.xmin,
0367                      templatePHIvsPHI.xmax,
0368                      templatePHIvsPHI.nbinsY,
0369                      templatePHIvsPHI.ymin,
0370                      templatePHIvsPHI.ymax);
0371   h_L1TauPhivsTauPhi_EE_ =
0372       ibooker.book2D("L1TauPhivsTauPhi_EE",
0373                      "#phi_{tau}^{L1} vs #phi_{tau}^{offline} (EE); #phi_{tau}^{offline}; #phi_{tau}^{L1}",
0374                      templatePHIvsPHI.nbinsX,
0375                      templatePHIvsPHI.xmin,
0376                      templatePHIvsPHI.xmax,
0377                      templatePHIvsPHI.nbinsY,
0378                      templatePHIvsPHI.ymin,
0379                      templatePHIvsPHI.ymax);
0380   h_L1TauPhivsTauPhi_EB_EE_ =
0381       ibooker.book2D("L1TauPhivsTauPhi_EB_EE",
0382                      "#phi_{tau}^{L1} vs #phi_{tau}^{offline} (EB+EE); #phi_{tau}^{offline}; #phi_{tau}^{L1}",
0383                      templatePHIvsPHI.nbinsX,
0384                      templatePHIvsPHI.xmin,
0385                      templatePHIvsPHI.xmax,
0386                      templatePHIvsPHI.nbinsY,
0387                      templatePHIvsPHI.ymin,
0388                      templatePHIvsPHI.ymax);
0389 
0390   h_L1TauEtavsTauEta_ =
0391       ibooker.book2D("L1TauEtavsTauEta", "L1 Tau #eta vs PFTau #eta; PFTau #eta; L1 Tau #eta", 100, -3, 3, 100, -3, 3);
0392 
0393   // tau resolutions
0394   h_resolutionTauET_EB_ = ibooker.book1D(
0395       "resolutionTauET_EB", "tau ET resolution (EB); (L1 Tau E_{T} - PFTau E_{T})/PFTau E_{T}; events", 50, -1, 1.5);
0396   h_resolutionTauET_EE_ = ibooker.book1D(
0397       "resolutionTauET_EE", "tau ET resolution (EE); (L1 Tau E_{T} - PFTau E_{T})/PFTau E_{T}; events", 50, -1, 1.5);
0398   h_resolutionTauET_EB_EE_ =
0399       ibooker.book1D("resolutionTauET_EB_EE",
0400                      "tau ET resolution (EB+EE); (L1 Tau E_{T} - PFTau E_{T})/PFTau E_{T}; events",
0401                      50,
0402                      -1,
0403                      1.5);
0404 
0405   h_resolutionTauPhi_EB_ = ibooker.book1D("resolutionTauPhi_EB",
0406                                           "#phi_{tau} resolution (EB); #phi_{tau}^{L1} - #phi_{tau}^{offline}; events",
0407                                           120,
0408                                           -0.3,
0409                                           0.3);
0410   h_resolutionTauPhi_EE_ = ibooker.book1D(
0411       "resolutionTauPhi_EE", "tau #phi resolution (EE); #phi_{tau}^{L1} - #phi_{tau}^{offline}; events", 120, -0.3, 0.3);
0412   h_resolutionTauPhi_EB_EE_ =
0413       ibooker.book1D("resolutionTauPhi_EB_EE",
0414                      "tau #phi resolution (EB+EE); #phi_{tau}^{L1} - #phi_{tau}^{offline}; events",
0415                      120,
0416                      -0.3,
0417                      0.3);
0418 
0419   h_resolutionTauEta_ =
0420       ibooker.book1D("resolutionTauEta", "tau #eta resolution  (EB); L1 Tau #eta - PFTau #eta; events", 120, -0.3, 0.3);
0421 
0422   // tau turn-ons
0423   ibooker.setCurrentFolder(efficiencyFolder_);
0424   std::vector<float> tauBins(tauEfficiencyBins_.begin(), tauEfficiencyBins_.end());
0425   int nBins = tauBins.size() - 1;
0426   float* tauBinArray = &(tauBins[0]);
0427 
0428   for (auto threshold : tauEfficiencyThresholds_) {
0429     std::string str_threshold = std::to_string(int(threshold));
0430     h_efficiencyIsoTauET_EB_pass_[threshold] =
0431         ibooker.book1D("efficiencyIsoTauET_EB_threshold_" + str_threshold + "_Num",
0432                        "iso tau efficiency (EB); PFTau E_{T} (GeV); events",
0433                        nBins,
0434                        tauBinArray);
0435     h_efficiencyIsoTauET_EE_pass_[threshold] =
0436         ibooker.book1D("efficiencyIsoTauET_EE_threshold_" + str_threshold + "_Num",
0437                        "iso tau efficiency (EE); PFTau E_{T} (GeV); events",
0438                        nBins,
0439                        tauBinArray);
0440     h_efficiencyIsoTauET_EB_EE_pass_[threshold] =
0441         ibooker.book1D("efficiencyIsoTauET_EB_EE_threshold_" + str_threshold + "_Num",
0442                        "iso tau efficiency (EB+EE); PFTau E_{T} (GeV); events",
0443                        nBins,
0444                        tauBinArray);
0445 
0446     h_efficiencyIsoTauET_EB_total_[threshold] =
0447         ibooker.book1D("efficiencyIsoTauET_EB_threshold_" + str_threshold + "_Den",
0448                        "iso tau efficiency (EB); PFTau E_{T} (GeV); events",
0449                        nBins,
0450                        tauBinArray);
0451     h_efficiencyIsoTauET_EE_total_[threshold] =
0452         ibooker.book1D("efficiencyIsoTauET_EE_threshold_" + str_threshold + "_Den",
0453                        "iso tau efficiency (EE); PFTau E_{T} (GeV); events",
0454                        nBins,
0455                        tauBinArray);
0456     h_efficiencyIsoTauET_EB_EE_total_[threshold] =
0457         ibooker.book1D("efficiencyIsoTauET_EB_EE_threshold_" + str_threshold + "_Den",
0458                        "iso tau efficiency (EB+EE); PFTau E_{T} (GeV); events",
0459                        nBins,
0460                        tauBinArray);
0461 
0462     // non iso
0463     h_efficiencyNonIsoTauET_EB_pass_[threshold] =
0464         ibooker.book1D("efficiencyNonIsoTauET_EB_threshold_" + str_threshold + "_Num",
0465                        "inclusive tau efficiency (EB); PFTau E_{T} (GeV); events",
0466                        nBins,
0467                        tauBinArray);
0468     h_efficiencyNonIsoTauET_EE_pass_[threshold] =
0469         ibooker.book1D("efficiencyNonIsoTauET_EE_threshold_" + str_threshold + "_Num",
0470                        "inclusive tau efficiency (EE); PFTau E_{T} (GeV); events",
0471                        nBins,
0472                        tauBinArray);
0473     h_efficiencyNonIsoTauET_EB_EE_pass_[threshold] =
0474         ibooker.book1D("efficiencyNonIsoTauET_EB_EE_threshold_" + str_threshold + "_Num",
0475                        "inclusive tau efficiency (EB+EE); PFTau E_{T} (GeV); events",
0476                        nBins,
0477                        tauBinArray);
0478 
0479     h_efficiencyNonIsoTauET_EB_total_[threshold] =
0480         ibooker.book1D("efficiencyNonIsoTauET_EB_threshold_" + str_threshold + "_Den",
0481                        "inclusive tau efficiency (EB); PFTau E_{T} (GeV); events",
0482                        nBins,
0483                        tauBinArray);
0484     h_efficiencyNonIsoTauET_EE_total_[threshold] =
0485         ibooker.book1D("efficiencyNonIsoTauET_EE_threshold_" + str_threshold + "_Den",
0486                        "inclusive tau efficiency (EE); PFTau E_{T} (GeV); events",
0487                        nBins,
0488                        tauBinArray);
0489     h_efficiencyNonIsoTauET_EB_EE_total_[threshold] =
0490         ibooker.book1D("efficiencyNonIsoTauET_EB_EE_threshold_" + str_threshold + "_Den",
0491                        "inclusive tau efficiency (EB+EE); PFTau E_{T} (GeV); events",
0492                        nBins,
0493                        tauBinArray);
0494   }
0495 
0496   ibooker.cd();
0497 
0498   return;
0499 }
0500 
0501 const reco::Vertex L1TTauOffline::getPrimaryVertex(edm::Handle<reco::VertexCollection> const& vertex,
0502                                                    edm::Handle<reco::BeamSpot> const& beamSpot) {
0503   reco::Vertex::Point posVtx;
0504   reco::Vertex::Error errVtx;
0505 
0506   bool hasPrimaryVertex = false;
0507 
0508   if (vertex.isValid()) {
0509     for (auto vertexIt = vertex->begin(); vertexIt != vertex->end(); ++vertexIt) {
0510       if (vertexIt->isValid() && !vertexIt->isFake()) {
0511         posVtx = vertexIt->position();
0512         errVtx = vertexIt->error();
0513         hasPrimaryVertex = true;
0514         break;
0515       }
0516     }
0517   }
0518 
0519   if (!hasPrimaryVertex) {
0520     posVtx = beamSpot->position();
0521     errVtx(0, 0) = beamSpot->BeamWidthX();
0522     errVtx(1, 1) = beamSpot->BeamWidthY();
0523     errVtx(2, 2) = beamSpot->sigmaZ();
0524   }
0525 
0526   const reco::Vertex primaryVertex(posVtx, errVtx);
0527 
0528   return primaryVertex;
0529 }
0530 
0531 bool L1TTauOffline::matchHlt(edm::Handle<trigger::TriggerEvent> const& triggerEvent, const reco::Muon* muon) {
0532   double matchDeltaR = 9999;
0533 
0534   trigger::TriggerObjectCollection trigObjs = triggerEvent->getObjects();
0535 
0536   for (auto trigIndexIt = m_trigIndices.begin(); trigIndexIt != m_trigIndices.end(); ++trigIndexIt) {
0537     const vector<string> moduleLabels(m_hltConfig.moduleLabels(*trigIndexIt));
0538     const unsigned moduleIndex = m_hltConfig.size((*trigIndexIt)) - 2;
0539 
0540     const unsigned hltFilterIndex = triggerEvent->filterIndex(InputTag(moduleLabels[moduleIndex], "", trigProcess_));
0541 
0542     if (hltFilterIndex < triggerEvent->sizeFilters()) {
0543       const Keys triggerKeys(triggerEvent->filterKeys(hltFilterIndex));
0544       const Vids triggerVids(triggerEvent->filterIds(hltFilterIndex));
0545 
0546       const unsigned nTriggers = triggerVids.size();
0547       for (size_t iTrig = 0; iTrig < nTriggers; ++iTrig) {
0548         const TriggerObject trigObject = trigObjs[triggerKeys[iTrig]];
0549 
0550         double dRtmp = deltaR((*muon), trigObject);
0551         if (dRtmp < matchDeltaR)
0552           matchDeltaR = dRtmp;
0553       }
0554     }
0555   }
0556 
0557   return (matchDeltaR < m_MaxHltTauDR);
0558 }
0559 
0560 void L1TTauOffline::getTauL1tPairs(edm::Handle<l1t::TauBxCollection> const& l1tCands) {
0561   m_TauL1tPairs.clear();
0562 
0563   vector<l1t::Tau> l1tContainer;
0564   l1tContainer.reserve(l1tCands->size() + 1);
0565 
0566   for (auto tau = l1tCands->begin(0); tau != l1tCands->end(0); ++tau) {
0567     l1tContainer.push_back(*tau);
0568   }
0569 
0570   for (auto probeTauIt = m_ProbeTaus.begin(); probeTauIt != m_ProbeTaus.end(); ++probeTauIt) {
0571     TauL1TPair pairBestCand((*probeTauIt), nullptr);
0572 
0573     for (auto l1tIt = l1tContainer.begin(); l1tIt != l1tContainer.end(); ++l1tIt) {
0574       TauL1TPair pairTmpCand((*probeTauIt), &(*l1tIt));
0575 
0576       if (pairTmpCand.dR() < m_MaxL1tTauDR && pairTmpCand.l1tPt() > pairBestCand.l1tPt())
0577         pairBestCand = pairTmpCand;
0578     }
0579 
0580     m_TauL1tPairs.push_back(pairBestCand);
0581   }
0582 }
0583 
0584 void L1TTauOffline::getTightMuons(edm::Handle<reco::MuonCollection> const& muons,
0585                                   edm::Handle<reco::PFMETCollection> const& mets,
0586                                   const reco::Vertex& vertex,
0587                                   edm::Handle<trigger::TriggerEvent> const& trigEvent) {
0588   m_TightMuons.clear();
0589 
0590   const reco::PFMET* pfmet = nullptr;
0591   pfmet = &(mets->front());
0592 
0593   int nb_mu = 0;
0594 
0595   for (auto muonIt2 = muons->begin(); muonIt2 != muons->end(); ++muonIt2) {
0596     if (fabs(muonIt2->eta()) < 2.4 && muonIt2->pt() > 10 && muon::isLooseMuon((*muonIt2)) &&
0597         (muonIt2->pfIsolationR04().sumChargedHadronPt +
0598          max(muonIt2->pfIsolationR04().sumNeutralHadronEt + muonIt2->pfIsolationR04().sumPhotonEt -
0599                  0.5 * muonIt2->pfIsolationR04().sumPUPt,
0600              0.0)) /
0601                 muonIt2->pt() <
0602             0.3) {
0603       ++nb_mu;
0604     }
0605   }
0606   bool foundTightMu = false;
0607   for (auto muonIt = muons->begin(); muonIt != muons->end(); ++muonIt) {
0608     if (!matchHlt(trigEvent, &(*muonIt)))
0609       continue;
0610     float muiso = (muonIt->pfIsolationR04().sumChargedHadronPt +
0611                    max(muonIt->pfIsolationR04().sumNeutralHadronEt + muonIt->pfIsolationR04().sumPhotonEt -
0612                            0.5 * muonIt->pfIsolationR04().sumPUPt,
0613                        0.0)) /
0614                   muonIt->pt();
0615 
0616     if (muiso < 0.1 && nb_mu < 2 && !foundTightMu && fabs(muonIt->eta()) < 2.1 && muonIt->pt() > 24 &&
0617         muon::isLooseMuon((*muonIt))) {
0618       float mt = sqrt(pow(muonIt->pt() + pfmet->pt(), 2) - pow(muonIt->px() + pfmet->px(), 2) -
0619                       pow(muonIt->py() + pfmet->py(), 2));
0620       if (mt < 30) {
0621         m_TightMuons.push_back(&(*muonIt));
0622         foundTightMu = true;
0623       }
0624     }
0625   }
0626 }
0627 
0628 void L1TTauOffline::getProbeTaus(const edm::Event& iEvent,
0629                                  edm::Handle<reco::PFTauCollection> const& taus,
0630                                  edm::Handle<reco::MuonCollection> const& muons,
0631                                  const reco::Vertex& vertex) {
0632   m_ProbeTaus.clear();
0633 
0634   edm::Handle<reco::TauDiscriminatorContainer> antimu;
0635   iEvent.getByToken(AntiMuInputTag_, antimu);
0636   if (!antimu.isValid()) {
0637     edm::LogWarning("L1TTauOffline") << "invalid collection: reco::PFTauDiscriminator " << std::endl;
0638     return;
0639   }
0640 
0641   edm::Handle<reco::PFTauDiscriminator> dmf;
0642   iEvent.getByToken(DecayModeFindingInputTag_, dmf);
0643   if (!dmf.isValid()) {
0644     edm::LogWarning("L1TTauOffline") << "invalid collection: reco::PFTauDiscriminator " << std::endl;
0645     return;
0646   }
0647 
0648   edm::Handle<reco::TauDiscriminatorContainer> antiele;
0649   iEvent.getByToken(AntiEleInputTag_, antiele);
0650   if (!antiele.isValid()) {
0651     edm::LogWarning("L1TTauOffline") << "invalid collection: reco::PFTauDiscriminator " << std::endl;
0652     return;
0653   }
0654 
0655   edm::Handle<reco::TauDiscriminatorContainer> comb3T;
0656   iEvent.getByToken(comb3TInputTag_, comb3T);
0657   if (!comb3T.isValid()) {
0658     edm::LogWarning("L1TTauOffline") << "invalid collection: reco::PFTauDiscriminator " << std::endl;
0659     return;
0660   }
0661 
0662   if (!m_TightMuons.empty()) {
0663     TLorentzVector mymu;
0664     mymu.SetPtEtaPhiE(m_TightMuons[0]->pt(), m_TightMuons[0]->eta(), m_TightMuons[0]->phi(), m_TightMuons[0]->energy());
0665     int iTau = 0;
0666 
0667     // load indices from input provenance config if process history changed, in particular for the first event
0668     if (iEvent.processHistoryID() != phID_) {
0669       phID_ = iEvent.processHistoryID();
0670       {
0671         const edm::Provenance* prov = antimu.provenance();
0672         const std::vector<edm::ParameterSet> psetsFromProvenance =
0673             edm::parameterSet(prov->stable(), iEvent.processHistory())
0674                 .getParameter<std::vector<edm::ParameterSet>>("IDWPdefinitions");
0675         for (uint i = 0; i < psetsFromProvenance.size(); i++) {
0676           if (psetsFromProvenance[i].getParameter<std::string>("IDname") == AntiMuWP_)
0677             AntiMuWPIndex_ = i;
0678         }
0679       }
0680       {
0681         const edm::Provenance* prov = antiele.provenance();
0682         const std::vector<std::string> psetsFromProvenance =
0683             edm::parameterSet(prov->stable(), iEvent.processHistory())
0684                 .getParameter<std::vector<std::string>>("workingPoints");
0685         for (uint i = 0; i < psetsFromProvenance.size(); i++) {
0686           if (psetsFromProvenance[i] == AntiEleWP_)
0687             AntiEleWPIndex_ = i;
0688         }
0689       }
0690       {
0691         const edm::Provenance* prov = comb3T.provenance();
0692         const std::vector<edm::ParameterSet> psetsFromProvenance =
0693             edm::parameterSet(prov->stable(), iEvent.processHistory())
0694                 .getParameter<std::vector<edm::ParameterSet>>("IDWPdefinitions");
0695         for (uint i = 0; i < psetsFromProvenance.size(); i++) {
0696           if (psetsFromProvenance[i].getParameter<std::string>("IDname") == comb3TWP_)
0697             comb3TWPIndex_ = i;
0698         }
0699       }
0700     }
0701 
0702     for (auto tauIt = taus->begin(); tauIt != taus->end(); ++tauIt, ++iTau) {
0703       reco::PFTauRef tauCandidate(taus, iTau);
0704       TLorentzVector mytau;
0705       mytau.SetPtEtaPhiE(tauIt->pt(), tauIt->eta(), tauIt->phi(), tauIt->energy());
0706 
0707       if ((*antimu)[tauCandidate].workingPoints.empty()) {
0708         edm::LogWarning("L1TTauOffline") << "This offline tau has no antimu discriminator, skipping" << std::endl;
0709         continue;
0710       }
0711       if ((*antiele)[tauCandidate].workingPoints.empty()) {
0712         edm::LogWarning("L1TTauOffline") << "This offline tau has no antiele discriminator, skipping" << std::endl;
0713         continue;
0714       }
0715       if ((*comb3T)[tauCandidate].workingPoints.empty()) {
0716         edm::LogWarning("L1TTauOffline") << "This offline tau has no comb3T discriminator, skipping" << std::endl;
0717         continue;
0718       }
0719 
0720       if (fabs(tauIt->charge()) == 1 && fabs(tauIt->eta()) < 2.1 && tauIt->pt() > 20 &&
0721           (*antimu)[tauCandidate].workingPoints[AntiMuWPIndex_] &&
0722           (*antiele)[tauCandidate].workingPoints[AntiEleWPIndex_] && (*dmf)[tauCandidate] > 0.5 &&
0723           (*comb3T)[tauCandidate].workingPoints[comb3TWPIndex_]) {
0724         if (mymu.DeltaR(mytau) > 0.5 && (mymu + mytau).M() > 40 && (mymu + mytau).M() < 80 &&
0725             m_TightMuons[0]->charge() * tauIt->charge() < 0) {
0726           m_ProbeTaus.push_back(&(*tauIt));
0727         }
0728       }
0729     }
0730   }
0731 }
0732 
0733 void L1TTauOffline::normalise2DHistogramsToBinArea() {
0734   std::vector<MonitorElement*> monElementstoNormalize = {h_L1TauETvsTauET_EB_,
0735                                                          h_L1TauETvsTauET_EE_,
0736                                                          h_L1TauETvsTauET_EB_EE_,
0737                                                          h_L1TauPhivsTauPhi_EB_,
0738                                                          h_L1TauPhivsTauPhi_EE_,
0739                                                          h_L1TauPhivsTauPhi_EB_EE_,
0740                                                          h_L1TauEtavsTauEta_};
0741 
0742   for (auto mon : monElementstoNormalize) {
0743     if (mon != nullptr) {
0744       auto h = mon->getTH2F();
0745       if (h != nullptr) {
0746         h->Scale(1, "width");
0747       }
0748     }
0749   }
0750 }
0751 // define this as a plug-in
0752 DEFINE_FWK_MODULE(L1TTauOffline);