File indexing completed on 2024-04-06 12:09:56
0001 #include "DQMOffline/Trigger/interface/HLTTauRefProducer.h"
0002 #include "TLorentzVector.h"
0003
0004 #include "DataFormats/TauReco/interface/PFTau.h"
0005 #include "DataFormats/TauReco/interface/PFTauDiscriminator.h"
0006 #include "DataFormats/TauReco/interface/TauDiscriminatorContainer.h"
0007
0008 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
0009 #include "DataFormats/EgammaCandidates/interface/GsfElectronFwd.h"
0010 #include "DataFormats/EgammaCandidates/interface/Electron.h"
0011 #include "DataFormats/EgammaCandidates/interface/ElectronFwd.h"
0012 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0013 #include "DataFormats/TrackReco/interface/Track.h"
0014
0015 #include "DataFormats/MuonReco/interface/Muon.h"
0016 #include "DataFormats/JetReco/interface/CaloJet.h"
0017 #include "TLorentzVector.h"
0018 #include "DataFormats/Math/interface/deltaR.h"
0019
0020 #include "DataFormats/CaloTowers/interface/CaloTower.h"
0021 #include "DataFormats/CaloTowers/interface/CaloTowerCollection.h"
0022 #include "DataFormats/CaloTowers/interface/CaloTowerDefs.h"
0023 #include "Math/GenVector/VectorUtil.h"
0024
0025 using namespace edm;
0026 using namespace reco;
0027 using namespace std;
0028
0029 HLTTauRefProducer::HLTTauRefProducer(const edm::ParameterSet& iConfig) {
0030
0031 {
0032 auto const& pfTau = iConfig.getUntrackedParameter<edm::ParameterSet>("PFTaus");
0033 PFTaus_ = consumes<reco::PFTauCollection>(pfTau.getUntrackedParameter<InputTag>("PFTauProducer"));
0034 auto discs = pfTau.getUntrackedParameter<vector<InputTag>>("PFTauDiscriminators");
0035 auto discConts = pfTau.getUntrackedParameter<vector<InputTag>>("PFTauDiscriminatorContainers");
0036 PFTauDisContWPs_ = pfTau.getUntrackedParameter<vector<std::string>>("PFTauDiscriminatorContainerWPs");
0037 if (discConts.size() != PFTauDisContWPs_.size())
0038 throw cms::Exception("Configuration") << "HLTTauRefProducer: Input parameters PFTauDiscriminatorContainers and "
0039 "PFTauDiscriminatorContainerWPs must have the same number of entries!\n";
0040 for (auto const& tag : discs) {
0041 PFTauDis_.push_back(consumes<reco::PFTauDiscriminator>(tag));
0042 }
0043 for (auto const& tag : discConts) {
0044 PFTauDisCont_.push_back(consumes<reco::TauDiscriminatorContainer>(tag));
0045 }
0046 doPFTaus_ = pfTau.getUntrackedParameter<bool>("doPFTaus", false);
0047 ptMinPFTau_ = pfTau.getUntrackedParameter<double>("ptMin", 15.);
0048 etaMinPFTau_ = pfTau.getUntrackedParameter<double>("etaMin", -2.5);
0049 etaMaxPFTau_ = pfTau.getUntrackedParameter<double>("etaMax", 2.5);
0050 phiMinPFTau_ = pfTau.getUntrackedParameter<double>("phiMin", -3.15);
0051 phiMaxPFTau_ = pfTau.getUntrackedParameter<double>("phiMax", 3.15);
0052 }
0053
0054 {
0055 auto const& electrons = iConfig.getUntrackedParameter<edm::ParameterSet>("Electrons");
0056 Electrons_ = consumes<reco::GsfElectronCollection>(electrons.getUntrackedParameter<InputTag>("ElectronCollection"));
0057 doElectrons_ = electrons.getUntrackedParameter<bool>("doElectrons", false);
0058 e_ctfTrackCollectionSrc_ = electrons.getUntrackedParameter<InputTag>("TrackCollection");
0059 e_ctfTrackCollection_ = consumes<reco::TrackCollection>(e_ctfTrackCollectionSrc_);
0060 ptMinElectron_ = electrons.getUntrackedParameter<double>("ptMin", 15.);
0061 e_doTrackIso_ = electrons.getUntrackedParameter<bool>("doTrackIso", false);
0062 e_trackMinPt_ = electrons.getUntrackedParameter<double>("ptMinTrack", 1.5);
0063 e_lipCut_ = electrons.getUntrackedParameter<double>("lipMinTrack", 1.5);
0064 e_minIsoDR_ = electrons.getUntrackedParameter<double>("InnerConeDR", 0.02);
0065 e_maxIsoDR_ = electrons.getUntrackedParameter<double>("OuterConeDR", 0.6);
0066 e_isoMaxSumPt_ = electrons.getUntrackedParameter<double>("MaxIsoVar", 0.02);
0067 }
0068
0069 {
0070 auto const& muons = iConfig.getUntrackedParameter<edm::ParameterSet>("Muons");
0071 Muons_ = consumes<reco::MuonCollection>(muons.getUntrackedParameter<InputTag>("MuonCollection"));
0072 doMuons_ = muons.getUntrackedParameter<bool>("doMuons", false);
0073 ptMinMuon_ = muons.getUntrackedParameter<double>("ptMin", 15.);
0074 }
0075
0076 {
0077 auto const& jets = iConfig.getUntrackedParameter<edm::ParameterSet>("Jets");
0078 Jets_ = consumes<reco::CaloJetCollection>(jets.getUntrackedParameter<InputTag>("JetCollection"));
0079 doJets_ = jets.getUntrackedParameter<bool>("doJets");
0080 ptMinJet_ = jets.getUntrackedParameter<double>("etMin");
0081 }
0082
0083 {
0084 auto const& towers = iConfig.getUntrackedParameter<edm::ParameterSet>("Towers");
0085 Towers_ = consumes<CaloTowerCollection>(towers.getUntrackedParameter<InputTag>("TowerCollection"));
0086 doTowers_ = towers.getUntrackedParameter<bool>("doTowers");
0087 ptMinTower_ = towers.getUntrackedParameter<double>("etMin");
0088 towerIsol_ = towers.getUntrackedParameter<double>("towerIsolation");
0089 }
0090
0091 {
0092 auto const& photons = iConfig.getUntrackedParameter<edm::ParameterSet>("Photons");
0093 Photons_ = consumes<reco::PhotonCollection>(photons.getUntrackedParameter<InputTag>("PhotonCollection"));
0094 doPhotons_ = photons.getUntrackedParameter<bool>("doPhotons");
0095 ptMinPhoton_ = photons.getUntrackedParameter<double>("etMin");
0096 photonEcalIso_ = photons.getUntrackedParameter<double>("ECALIso");
0097 }
0098
0099 {
0100 auto const& met = iConfig.getUntrackedParameter<edm::ParameterSet>("MET");
0101 MET_ = consumes<reco::CaloMETCollection>(met.getUntrackedParameter<InputTag>("METCollection"));
0102 doMET_ = met.getUntrackedParameter<bool>("doMET", false);
0103 ptMinMET_ = met.getUntrackedParameter<double>("ptMin", 15.);
0104 }
0105
0106 etaMin_ = iConfig.getUntrackedParameter<double>("EtaMin", -2.5);
0107 etaMax_ = iConfig.getUntrackedParameter<double>("EtaMax", 2.5);
0108 phiMin_ = iConfig.getUntrackedParameter<double>("PhiMin", -3.15);
0109 phiMax_ = iConfig.getUntrackedParameter<double>("PhiMax", 3.15);
0110
0111
0112 produces<LorentzVectorCollection>("PFTaus");
0113 produces<LorentzVectorCollection>("Electrons");
0114 produces<LorentzVectorCollection>("Muons");
0115 produces<LorentzVectorCollection>("Jets");
0116 produces<LorentzVectorCollection>("Photons");
0117 produces<LorentzVectorCollection>("Towers");
0118 produces<LorentzVectorCollection>("MET");
0119 }
0120
0121 void HLTTauRefProducer::produce(edm::StreamID iID, edm::Event& iEvent, edm::EventSetup const&) const {
0122 if (doPFTaus_)
0123 doPFTaus(iID, iEvent);
0124 if (doElectrons_)
0125 doElectrons(iEvent);
0126 if (doMuons_)
0127 doMuons(iEvent);
0128 if (doJets_)
0129 doJets(iEvent);
0130 if (doPhotons_)
0131 doPhotons(iEvent);
0132 if (doTowers_)
0133 doTowers(iEvent);
0134 if (doMET_)
0135 doMET(iEvent);
0136 }
0137
0138 void HLTTauRefProducer::doPFTaus(edm::StreamID iID, edm::Event& iEvent) const {
0139 auto product_PFTaus = make_unique<LorentzVectorCollection>();
0140
0141 edm::Handle<PFTauCollection> pftaus;
0142 if (iEvent.getByToken(PFTaus_, pftaus)) {
0143
0144 if (streamCache(iID)->first != iEvent.processHistoryID()) {
0145 streamCache(iID)->first = iEvent.processHistoryID();
0146 streamCache(iID)->second.resize(PFTauDisContWPs_.size());
0147 for (size_t i = 0; i < PFTauDisCont_.size(); ++i) {
0148 auto const aHandle = iEvent.getHandle(PFTauDisCont_[i]);
0149 auto const aProv = aHandle.provenance();
0150 if (aProv == nullptr)
0151 aHandle.whyFailed()->raise();
0152 const auto& psetsFromProvenance = edm::parameterSet(aProv->stable(), iEvent.processHistory());
0153 if (psetsFromProvenance.exists("workingPoints")) {
0154 auto const idlist = psetsFromProvenance.getParameter<std::vector<std::string>>("workingPoints");
0155 bool found = false;
0156 for (size_t j = 0; j < idlist.size(); ++j) {
0157 if (PFTauDisContWPs_[i] == idlist[j]) {
0158 found = true;
0159 streamCache(iID)->second[i] = j;
0160 }
0161 }
0162 if (!found)
0163 throw cms::Exception("Configuration")
0164 << "HLTTauRefProducer: Requested working point '" << PFTauDisContWPs_[i] << "' not found!\n";
0165 } else if (psetsFromProvenance.exists("IDWPdefinitions")) {
0166 auto const idlist = psetsFromProvenance.getParameter<std::vector<edm::ParameterSet>>("IDWPdefinitions");
0167 bool found = false;
0168 for (size_t j = 0; j < idlist.size(); ++j) {
0169 if (PFTauDisContWPs_[i] == idlist[j].getParameter<std::string>("IDname")) {
0170 found = true;
0171 streamCache(iID)->second[i] = j;
0172 }
0173 }
0174 if (!found)
0175 throw cms::Exception("Configuration")
0176 << "HLTTauRefProducer: Requested working point '" << PFTauDisContWPs_[i] << "' not found!\n";
0177 } else
0178 throw cms::Exception("Configuration")
0179 << "HLTTauRefProducer: No suitable ID list found in provenace config!\n";
0180 }
0181 }
0182 for (unsigned int i = 0; i < pftaus->size(); ++i) {
0183 auto const& pftau = (*pftaus)[i];
0184 if (pftau.pt() > ptMinPFTau_ && pftau.eta() > etaMinPFTau_ && pftau.eta() < etaMaxPFTau_ &&
0185 pftau.phi() > phiMinPFTau_ && pftau.phi() < phiMaxPFTau_) {
0186 reco::PFTauRef thePFTau{pftaus, i};
0187 bool passAll{true};
0188
0189 for (auto const& token : PFTauDis_) {
0190 edm::Handle<reco::PFTauDiscriminator> pftaudis;
0191 if (iEvent.getByToken(token, pftaudis)) {
0192 if ((*pftaudis)[thePFTau] < 0.5) {
0193 passAll = false;
0194 break;
0195 }
0196 } else {
0197 passAll = false;
0198 break;
0199 }
0200 }
0201
0202 int idx = 0;
0203 for (auto const& token : PFTauDisCont_) {
0204 edm::Handle<reco::TauDiscriminatorContainer> pftaudis;
0205 if (iEvent.getByToken(token, pftaudis)) {
0206
0207 if ((*pftaudis)[thePFTau].workingPoints.empty() ||
0208 !(*pftaudis)[thePFTau].workingPoints.at(streamCache(iID)->second[idx])) {
0209 passAll = false;
0210 break;
0211 }
0212 } else {
0213 passAll = false;
0214 break;
0215 }
0216 idx++;
0217 }
0218 if (passAll) {
0219 product_PFTaus->emplace_back(pftau.px(), pftau.py(), pftau.pz(), pftau.energy());
0220 }
0221 }
0222 }
0223 }
0224 iEvent.put(std::move(product_PFTaus), "PFTaus");
0225 }
0226
0227 void HLTTauRefProducer::doElectrons(edm::Event& iEvent) const {
0228 auto product_Electrons = make_unique<LorentzVectorCollection>();
0229
0230 edm::Handle<reco::TrackCollection> pCtfTracks;
0231 if (!iEvent.getByToken(e_ctfTrackCollection_, pCtfTracks)) {
0232 edm::LogInfo("") << "Error! Can't get " << e_ctfTrackCollectionSrc_.label() << " by label. ";
0233 iEvent.put(std::move(product_Electrons), "Electrons");
0234 return;
0235 }
0236
0237 edm::Handle<GsfElectronCollection> electrons;
0238 if (iEvent.getByToken(Electrons_, electrons)) {
0239 for (size_t i = 0; i < electrons->size(); ++i) {
0240 edm::Ref<reco::GsfElectronCollection> electronRef(electrons, i);
0241 auto const& electron = (*electrons)[i];
0242 if (electron.pt() > ptMinElectron_ && fabs(electron.eta()) < etaMax_) {
0243 if (e_doTrackIso_) {
0244 double sum_of_pt_ele{};
0245 for (auto const& tr : *pCtfTracks) {
0246 double const lip{electron.gsfTrack()->dz() - tr.dz()};
0247 if (tr.pt() > e_trackMinPt_ && fabs(lip) < e_lipCut_) {
0248 double dphi{fabs(tr.phi() - electron.trackMomentumAtVtx().phi())};
0249 if (dphi > acos(-1.)) {
0250 dphi = 2 * acos(-1.) - dphi;
0251 }
0252 double const deta{fabs(tr.eta() - electron.trackMomentumAtVtx().eta())};
0253 double const dr_ctf_ele{sqrt(deta * deta + dphi * dphi)};
0254 if ((dr_ctf_ele > e_minIsoDR_) && (dr_ctf_ele < e_maxIsoDR_)) {
0255 double const cft_pt_2{tr.pt() * tr.pt()};
0256 sum_of_pt_ele += cft_pt_2;
0257 }
0258 }
0259 }
0260 double const isolation_value_ele{sum_of_pt_ele /
0261 (electron.trackMomentumAtVtx().Rho() * electron.trackMomentumAtVtx().Rho())};
0262 if (isolation_value_ele < e_isoMaxSumPt_) {
0263 product_Electrons->emplace_back(electron.px(), electron.py(), electron.pz(), electron.energy());
0264 }
0265 } else {
0266 product_Electrons->emplace_back(electron.px(), electron.py(), electron.pz(), electron.energy());
0267 }
0268 }
0269 }
0270 }
0271 iEvent.put(std::move(product_Electrons), "Electrons");
0272 }
0273
0274 void HLTTauRefProducer::doMuons(edm::Event& iEvent) const {
0275 auto product_Muons = make_unique<LorentzVectorCollection>();
0276
0277 edm::Handle<MuonCollection> muons;
0278 if (iEvent.getByToken(Muons_, muons)) {
0279 for (auto const& muon : *muons) {
0280 if (muon.pt() > ptMinMuon_ && muon.eta() > etaMin_ && muon.eta() < etaMax_ && muon.phi() > phiMin_ &&
0281 muon.phi() < phiMax_) {
0282 product_Muons->emplace_back(muon.px(), muon.py(), muon.pz(), muon.energy());
0283 }
0284 }
0285 }
0286 iEvent.put(std::move(product_Muons), "Muons");
0287 }
0288
0289 void HLTTauRefProducer::doJets(edm::Event& iEvent) const {
0290 auto product_Jets = make_unique<LorentzVectorCollection>();
0291
0292 edm::Handle<CaloJetCollection> jets;
0293 if (iEvent.getByToken(Jets_, jets)) {
0294 for (auto const& jet : *jets) {
0295 if (jet.et() > ptMinJet_ && jet.eta() > etaMin_ && jet.eta() < etaMax_ && jet.phi() > phiMin_ &&
0296 jet.phi() < phiMax_) {
0297 product_Jets->emplace_back(jet.px(), jet.py(), jet.pz(), jet.energy());
0298 }
0299 }
0300 }
0301 iEvent.put(std::move(product_Jets), "Jets");
0302 }
0303
0304 void HLTTauRefProducer::doTowers(edm::Event& iEvent) const {
0305 auto product_Towers = make_unique<LorentzVectorCollection>();
0306
0307 edm::Handle<CaloTowerCollection> towers;
0308 if (iEvent.getByToken(Towers_, towers)) {
0309 for (auto const& tower1 : *towers) {
0310 if (tower1.pt() > ptMinTower_ && tower1.eta() > etaMin_ && tower1.eta() < etaMax_ && tower1.phi() > phiMin_ &&
0311 tower1.phi() < phiMax_) {
0312
0313 double isolET{};
0314 for (auto const& tower2 : *towers) {
0315 if (ROOT::Math::VectorUtil::DeltaR(tower1.p4(), tower2.p4()) < 0.5) {
0316 isolET += tower2.pt();
0317 }
0318 isolET -= tower1.pt();
0319 }
0320 if (isolET < towerIsol_) {
0321 product_Towers->emplace_back(tower1.px(), tower1.py(), tower1.pz(), tower1.energy());
0322 }
0323 }
0324 }
0325 }
0326 iEvent.put(std::move(product_Towers), "Towers");
0327 }
0328
0329 void HLTTauRefProducer::doPhotons(edm::Event& iEvent) const {
0330 auto product_Gammas = make_unique<LorentzVectorCollection>();
0331
0332 edm::Handle<PhotonCollection> photons;
0333 if (iEvent.getByToken(Photons_, photons)) {
0334 for (auto const& photon : *photons) {
0335 if (photon.ecalRecHitSumEtConeDR04() < photonEcalIso_ && photon.et() > ptMinPhoton_ && photon.eta() > etaMin_ &&
0336 photon.eta() < etaMax_ && photon.phi() > phiMin_ && photon.phi() < phiMax_) {
0337 product_Gammas->emplace_back(photon.px(), photon.py(), photon.pz(), photon.energy());
0338 }
0339 }
0340 }
0341 iEvent.put(std::move(product_Gammas), "Photons");
0342 }
0343
0344 void HLTTauRefProducer::doMET(edm::Event& iEvent) const {
0345 auto product_MET = make_unique<LorentzVectorCollection>();
0346
0347 edm::Handle<reco::CaloMETCollection> met;
0348 if (iEvent.getByToken(MET_, met) && !met->empty()) {
0349 auto const& metMom = met->front().p4();
0350 product_MET->emplace_back(metMom.Px(), metMom.Py(), 0, metMom.Pt());
0351 }
0352 iEvent.put(std::move(product_MET), "MET");
0353 }
0354
0355 #include "FWCore/Framework/interface/MakerMacros.h"
0356 DEFINE_FWK_MODULE(HLTTauRefProducer);
