HLTrigger/src/HLTGenValObjectMgr.cc

0001 #include "Validation/HLTrigger/interface/HLTGenValObjectMgr.h"
0002 #include "PhysicsTools/JetMCUtils/interface/JetMCTag.h"
0003
0004 HLTGenValObjectMgr::HLTGenValObjectMgr(const edm::ParameterSet& iConfig, edm::ConsumesCollector cc)
0005     : genParticleToken_(cc.consumes<reco::GenParticleCollection>(iConfig.getParameter<edm::InputTag>("genParticles"))),
0006       genMETToken_(cc.consumes<reco::GenMETCollection>(iConfig.getParameter<edm::InputTag>("genMET"))),
0007       ak4GenJetToken_(cc.consumes<reco::GenJetCollection>(iConfig.getParameter<edm::InputTag>("ak4GenJets"))),
0008       ak8GenJetToken_(cc.consumes<reco::GenJetCollection>(iConfig.getParameter<edm::InputTag>("ak8GenJets"))),
0009       tauGenJetToken_(cc.consumes<reco::GenJetCollection>(iConfig.getParameter<edm::InputTag>("tauGenJets"))),
0010       maxPromptGenJetFrac_(iConfig.getParameter<double>("maxPromptGenJetFrac")),
0011       minPtForGenHT_(iConfig.getParameter<double>("minPtForGenHT")),
0012       maxAbsEtaForGenHT_(iConfig.getParameter<double>("maxAbsEtaForGenHT")) {}
0013
0014 edm::ParameterSetDescription HLTGenValObjectMgr::makePSetDescription() {
0015   edm::ParameterSetDescription desc;
0016   desc.add<edm::InputTag>("genParticles", edm::InputTag("genParticles"));
0017   desc.add<edm::InputTag>("genMET", edm::InputTag("genMetTrue"));
0018   desc.add<edm::InputTag>("ak4GenJets", edm::InputTag("ak4GenJetsNoNu"));
0019   desc.add<edm::InputTag>("ak8GenJets", edm::InputTag("ak8GenJetsNoNu"));
0020   desc.add<edm::InputTag>("tauGenJets", edm::InputTag("tauGenJets"));
0021   desc.add<double>("maxPromptGenJetFrac", 0.1);
0022   desc.add<double>("minPtForGenHT", 30);
0023   desc.add<double>("maxAbsEtaForGenHT", 2.5);
0024   return desc;
0025 }
0026
0027 // this method handles the different object types and collections that can be used for efficiency calculation
0028 std::vector<HLTGenValObject> HLTGenValObjectMgr::getGenValObjects(const edm::Event& iEvent,
0029                                                                   const std::string& objType) {
0030   std::vector<HLTGenValObject> objects;  // the vector of objects to be filled
0031
0032   // handle object type
0033   std::vector<std::string> implementedGenParticles = {"ele", "pho", "mu", "tau"};
0034   if (std::find(implementedGenParticles.begin(), implementedGenParticles.end(), objType) !=
0035       implementedGenParticles.end()) {
0036     objects = getGenParticles(iEvent, objType);
0037   } else if (objType == "AK4jet") {  // ak4 jets, using the ak4GenJets collection
0038     const auto& genJets = iEvent.getHandle(ak4GenJetToken_);
0039     for (size_t i = 0; i < genJets->size(); i++) {
0040       const reco::GenJet p = (*genJets)[i];
0041       if (passGenJetID(p)) {
0042         objects.emplace_back(p);
0043       }
0044     }
0045   } else if (objType == "AK8jet") {  // ak8 jets, using the ak8GenJets collection
0046     const auto& genJets = iEvent.getHandle(ak8GenJetToken_);
0047     for (size_t i = 0; i < genJets->size(); i++) {
0048       const reco::GenJet p = (*genJets)[i];
0049       if (passGenJetID(p)) {
0050         objects.emplace_back(p);
0051       }
0052     }
0053   } else if (objType == "AK4HT") {  // ak4-based HT, using the ak4GenJets collection
0054     const auto& genJets = iEvent.getHandle(ak4GenJetToken_);
0055     if (!genJets->empty()) {
0056       double HTsum = 0.;
0057       for (const auto& genJet : *genJets) {
0058         if (genJet.pt() > minPtForGenHT_ && std::abs(genJet.eta()) < maxAbsEtaForGenHT_ && passGenJetID(genJet)) {
0059           HTsum += genJet.pt();
0060         }
0061       }
0062       if (HTsum > 0)
0063         objects.emplace_back(reco::Candidate::PolarLorentzVector(HTsum, 0, 0, 0));
0064     }
0065   } else if (objType == "AK8HT") {  // ak8-based HT, using the ak8GenJets collection
0066     const auto& genJets = iEvent.getHandle(ak8GenJetToken_);
0067     if (!genJets->empty()) {
0068       double HTsum = 0.;
0069       for (const auto& genJet : *genJets) {
0070         if (genJet.pt() > minPtForGenHT_ && std::abs(genJet.eta()) < maxAbsEtaForGenHT_ && passGenJetID(genJet)) {
0071           HTsum += genJet.pt();
0072         }
0073       }
0074       if (HTsum > 0)
0075         objects.emplace_back(reco::Candidate::PolarLorentzVector(HTsum, 0, 0, 0));
0076     }
0077   } else if (objType == "MET") {  // MET, using genMET
0078     const auto& genMET = iEvent.getHandle(genMETToken_);
0079     if (!genMET->empty()) {
0080       auto genMETpt = (*genMET)[0].pt();
0081       objects.emplace_back(reco::Candidate::PolarLorentzVector(genMETpt, 0, 0, 0));
0082     }
0083   } else if (objType == "tauHAD") {
0084     const auto& tauJets = iEvent.getHandle(tauGenJetToken_);
0085     for (const auto& tauJet : *tauJets) {
0086       const std::string& decayMode = JetMCTagUtils::genTauDecayMode(tauJet);
0087       if (decayMode != "electron" && decayMode != "muon") {
0088         objects.emplace_back(tauJet);
0089       }
0090     }
0091   } else
0092     throw cms::Exception("InputError") << "Generator-level validation is not available for type " << objType << ".\n"
0093                                        << "Please check for a potential spelling error.\n";
0094
0095   return objects;
0096 }
0097
0098 // in case of GenParticles, a subset of the entire collection needs to be chosen
0099 std::vector<HLTGenValObject> HLTGenValObjectMgr::getGenParticles(const edm::Event& iEvent, const std::string& objType) {
0100   std::vector<HLTGenValObject> objects;  // vector to be filled
0101
0102   const auto& genParticles = iEvent.getHandle(genParticleToken_);  // getting all GenParticles
0103
0104   // we need to ge the ID corresponding to the desired GenParticle type
0105   int pdgID = -1;  // setting to -1 should not be needed, but prevents the compiler warning :)
0106   if (objType == "ele")
0107     pdgID = 11;
0108   else if (objType == "pho")
0109     pdgID = 22;
0110   else if (objType == "mu")
0111     pdgID = 13;
0112   else if (objType == "tau")
0113     pdgID = 15;
0114
0115   // main loop over GenParticles
0116   for (size_t i = 0; i < genParticles->size(); ++i) {
0117     const reco::GenParticle p = (*genParticles)[i];
0118
0119     // only GenParticles with correct ID
0120     if (std::abs(p.pdgId()) != pdgID)
0121       continue;
0122
0123     // checking if particle comes from "hard process"
0124     if (p.isHardProcess()) {
0125       // depending on the particle type, last particle before or after FSR is chosen
0126       if ((objType == "ele") || (objType == "pho"))
0127         objects.emplace_back(getLastCopyPreFSR(p));
0128       else if ((objType == "mu") || (objType == "tau"))
0129         objects.emplace_back(getLastCopy(p));
0130     }
0131   }
0132
0133   return objects;
0134 }
0135
0136 // function returning the last GenParticle in a decay chain before FSR
0137 reco::GenParticle HLTGenValObjectMgr::getLastCopyPreFSR(reco::GenParticle part) {
0138   const auto& daughters = part.daughterRefVector();
0139   if (daughters.size() == 1 && daughters.at(0)->pdgId() == part.pdgId())
0140     return getLastCopyPreFSR(*daughters.at(0).get());  // recursion, whooo
0141   else
0142     return part;
0143 }
0144
0145 // function returning the last GenParticle in a decay chain
0146 reco::GenParticle HLTGenValObjectMgr::getLastCopy(reco::GenParticle part) {
0147   for (const auto& daughter : part.daughterRefVector()) {
0148     if (daughter->pdgId() == part.pdgId())
0149       return getLastCopy(*daughter.get());
0150   }
0151   return part;
0152 }
0153
0154 bool HLTGenValObjectMgr::passGenJetID(const reco::GenJet& jet) {
0155   float promptPt = 0;
0156   for (const auto& genPart : jet.getGenConstituents()) {
0157     if (genPart->fromHardProcessFinalState()) {
0158       promptPt += genPart->pt();
0159     }
0160   }
0161   return promptPt < jet.pt() * maxPromptGenJetFrac_;
0162 }