JetMET/plugins/HLTScoutingPFProducer.cc

0001 // -*- C++ -*-
0002 //
0003 // Package:    HLTrigger/JetMET
0004 // Class:      HLTScoutingPFProducer
0005 //
0006 /**\class HLTScoutingPFProducer HLTScoutingPFProducer.cc HLTrigger/JetMET/plugins/HLTScoutingPFProducer.cc
0007
0008 Description: Producer for ScoutingPFJets from reco::PFJet objects, ScoutingVertexs from reco::Vertexs and ScoutingParticles from reco::PFCandidates
0009
0010 */
0011 //
0012 // Original Author:  Dustin James Anderson
0013 //         Created:  Fri, 12 Jun 2015 15:49:20 GMT
0014 //
0015 //
0016
0017 // system include files
0018 #include <memory>
0019
0020 // user include files
0021 #include "FWCore/Framework/interface/Frameworkfwd.h"
0022 #include "FWCore/Framework/interface/global/EDProducer.h"
0023 #include "FWCore/Framework/interface/Event.h"
0024 #include "FWCore/Framework/interface/MakerMacros.h"
0025 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0026
0027 #include "DataFormats/JetReco/interface/PFJet.h"
0028 #include "DataFormats/METReco/interface/PFMET.h"
0029 #include "DataFormats/METReco/interface/PFMETCollection.h"
0030 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidateFwd.h"
0031 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0032 #include "DataFormats/VertexReco/interface/Vertex.h"
0033 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0034 #include "DataFormats/BTauReco/interface/JetTag.h"
0035
0036 #include "DataFormats/Scouting/interface/Run3ScoutingPFJet.h"
0037 #include "DataFormats/Scouting/interface/Run3ScoutingParticle.h"
0038 #include "DataFormats/Scouting/interface/Run3ScoutingVertex.h"
0039
0040 #include "DataFormats/Math/interface/deltaR.h"
0041
0042 #include "DataFormats/Math/interface/libminifloat.h"
0043
0044 #include "RecoBTag/FeatureTools/interface/deep_helpers.h"
0045
0046 class HLTScoutingPFProducer : public edm::global::EDProducer<> {
0047 public:
0048   explicit HLTScoutingPFProducer(const edm::ParameterSet &);
0049   ~HLTScoutingPFProducer() override;
0050
0051   static void fillDescriptions(edm::ConfigurationDescriptions &descriptions);
0052
0053 private:
0054   void produce(edm::StreamID sid, edm::Event &iEvent, edm::EventSetup const &setup) const final;
0055
0056   const edm::EDGetTokenT<reco::PFJetCollection> pfJetCollection_;
0057   const edm::EDGetTokenT<reco::JetTagCollection> pfJetTagCollection_;
0058   const edm::EDGetTokenT<reco::PFCandidateCollection> pfCandidateCollection_;
0059   const edm::EDGetTokenT<reco::VertexCollection> vertexCollection_;
0060   const edm::EDGetTokenT<reco::PFMETCollection> metCollection_;
0061   const edm::EDGetTokenT<double> rho_;
0062
0063   const double pfJetPtCut_;
0064   const double pfJetEtaCut_;
0065   const double pfCandidatePtCut_;
0066   const double pfCandidateEtaCut_;
0067   const int mantissaPrecision_;
0068
0069   const bool doJetTags_;
0070   const bool doCandidates_;
0071   const bool doMet_;
0072   const bool doTrackVars_;
0073   const bool relativeTrackVars_;
0074   const bool doCandIndsForJets_;
0075 };
0076
0077 //
0078 // constructors and destructor
0079 //
0080 HLTScoutingPFProducer::HLTScoutingPFProducer(const edm::ParameterSet &iConfig)
0081     : pfJetCollection_(consumes(iConfig.getParameter<edm::InputTag>("pfJetCollection"))),
0082       pfJetTagCollection_(consumes(iConfig.getParameter<edm::InputTag>("pfJetTagCollection"))),
0083       pfCandidateCollection_(consumes(iConfig.getParameter<edm::InputTag>("pfCandidateCollection"))),
0084       vertexCollection_(consumes(iConfig.getParameter<edm::InputTag>("vertexCollection"))),
0085       metCollection_(consumes(iConfig.getParameter<edm::InputTag>("metCollection"))),
0086       rho_(consumes(iConfig.getParameter<edm::InputTag>("rho"))),
0087       pfJetPtCut_(iConfig.getParameter<double>("pfJetPtCut")),
0088       pfJetEtaCut_(iConfig.getParameter<double>("pfJetEtaCut")),
0089       pfCandidatePtCut_(iConfig.getParameter<double>("pfCandidatePtCut")),
0090       pfCandidateEtaCut_(iConfig.getParameter<double>("pfCandidateEtaCut")),
0091       mantissaPrecision_(iConfig.getParameter<int>("mantissaPrecision")),
0092       doJetTags_(iConfig.getParameter<bool>("doJetTags")),
0093       doCandidates_(iConfig.getParameter<bool>("doCandidates")),
0094       doMet_(iConfig.getParameter<bool>("doMet")),
0095       doTrackVars_(iConfig.getParameter<bool>("doTrackVars")),
0096       relativeTrackVars_(iConfig.getParameter<bool>("relativeTrackVars")),
0097       doCandIndsForJets_(iConfig.getParameter<bool>("doCandIndsForJets")) {
0098   //register products
0099   produces<Run3ScoutingPFJetCollection>();
0100   produces<Run3ScoutingParticleCollection>();
0101   produces<double>("rho");
0102   produces<double>("pfMetPt");
0103   produces<double>("pfMetPhi");
0104 }
0105
0106 HLTScoutingPFProducer::~HLTScoutingPFProducer() = default;
0107
0108 // ------------ method called to produce the data  ------------
0109 void HLTScoutingPFProducer::produce(edm::StreamID sid, edm::Event &iEvent, edm::EventSetup const &setup) const {
0110   using namespace edm;
0111
0112   //get vertices
0113   Handle<reco::VertexCollection> vertexCollection;
0114   auto outVertices = std::make_unique<Run3ScoutingVertexCollection>();
0115   if (iEvent.getByToken(vertexCollection_, vertexCollection)) {
0116     for (auto const &vtx : *vertexCollection) {
0117       outVertices->emplace_back(MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.x(), mantissaPrecision_),
0118                                 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.y(), mantissaPrecision_),
0119                                 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.z(), mantissaPrecision_),
0120                                 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.zError(), mantissaPrecision_),
0121                                 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.xError(), mantissaPrecision_),
0122                                 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.yError(), mantissaPrecision_),
0123                                 vtx.tracksSize(),
0124                                 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.chi2(), mantissaPrecision_),
0125                                 vtx.ndof(),
0126                                 vtx.isValid());
0127     }
0128   }
0129
0130   //get rho
0131   Handle<double> rho;
0132   auto outRho = std::make_unique<double>(-999);
0133   if (iEvent.getByToken(rho_, rho)) {
0134     *outRho = *rho;
0135   }
0136
0137   //get MET
0138   Handle<reco::PFMETCollection> metCollection;
0139   auto outMetPt = std::make_unique<double>(-999);
0140   auto outMetPhi = std::make_unique<double>(-999);
0141   if (doMet_ && iEvent.getByToken(metCollection_, metCollection)) {
0142     outMetPt = std::make_unique<double>(metCollection->front().pt());
0143     outMetPhi = std::make_unique<double>(metCollection->front().phi());
0144   }
0145
0146   //get PF candidates
0147   Handle<reco::PFCandidateCollection> pfCandidateCollection;
0148   auto outPFCandidates = std::make_unique<Run3ScoutingParticleCollection>();
0149   if (doCandidates_ && iEvent.getByToken(pfCandidateCollection_, pfCandidateCollection)) {
0150     for (auto const &cand : *pfCandidateCollection) {
0151       if (cand.pt() > pfCandidatePtCut_ && std::abs(cand.eta()) < pfCandidateEtaCut_) {
0152         int vertex_index = -1;
0153         int index_counter = 0;
0154         double dr2 = 0.0001;
0155         for (auto const &vtx : *outVertices) {
0156           double tmp_dr2 = pow(vtx.x() - cand.vx(), 2) + pow(vtx.y() - cand.vy(), 2) + pow(vtx.z() - cand.vz(), 2);
0157           if (tmp_dr2 < dr2) {
0158             dr2 = tmp_dr2;
0159             vertex_index = index_counter;
0160           }
0161           if (dr2 == 0.0)
0162             break;
0163           ++index_counter;
0164         }
0165         float normchi2{0}, dz{0}, dxy{0}, dzError{0}, dxyError{0}, trk_pt{0}, trk_eta{0}, trk_phi{0};
0166         uint8_t lostInnerHits{0}, quality{0};
0167         if (doTrackVars_) {
0168           const auto *trk = cand.bestTrack();
0169           if (trk != nullptr) {
0170             normchi2 = MiniFloatConverter::reduceMantissaToNbitsRounding(trk->normalizedChi2(), mantissaPrecision_);
0171             lostInnerHits = btagbtvdeep::lost_inner_hits_from_pfcand(cand);
0172             quality = btagbtvdeep::quality_from_pfcand(cand);
0173             if (relativeTrackVars_) {
0174               trk_pt = MiniFloatConverter::reduceMantissaToNbitsRounding(trk->pt() - cand.pt(), mantissaPrecision_);
0175               trk_eta = MiniFloatConverter::reduceMantissaToNbitsRounding(trk->eta() - cand.eta(), mantissaPrecision_);
0176               trk_phi = MiniFloatConverter::reduceMantissaToNbitsRounding(trk->phi() - cand.phi(), mantissaPrecision_);
0177             } else {
0178               trk_pt = MiniFloatConverter::reduceMantissaToNbitsRounding(trk->pt(), mantissaPrecision_);
0179               trk_eta = MiniFloatConverter::reduceMantissaToNbitsRounding(trk->eta(), mantissaPrecision_);
0180               trk_phi = MiniFloatConverter::reduceMantissaToNbitsRounding(trk->phi(), mantissaPrecision_);
0181             }
0182             if (not vertexCollection->empty()) {
0183               const reco::Vertex &pv = (*vertexCollection)[0];
0184
0185               dz = trk->dz(pv.position());
0186               dzError = MiniFloatConverter::reduceMantissaToNbitsRounding(dz / trk->dzError(), mantissaPrecision_);
0187               dz = MiniFloatConverter::reduceMantissaToNbitsRounding(dz, mantissaPrecision_);
0188
0189               dxy = trk->dxy(pv.position());
0190               dxyError = MiniFloatConverter::reduceMantissaToNbitsRounding(dxy / trk->dxyError(), mantissaPrecision_);
0191               dxy = MiniFloatConverter::reduceMantissaToNbitsRounding(dxy, mantissaPrecision_);
0192             }
0193           } else {
0194             normchi2 = MiniFloatConverter::reduceMantissaToNbitsRounding(999, mantissaPrecision_);
0195           }
0196         }
0197         outPFCandidates->emplace_back(MiniFloatConverter::reduceMantissaToNbitsRounding(cand.pt(), mantissaPrecision_),
0198                                       MiniFloatConverter::reduceMantissaToNbitsRounding(cand.eta(), mantissaPrecision_),
0199                                       MiniFloatConverter::reduceMantissaToNbitsRounding(cand.phi(), mantissaPrecision_),
0200                                       cand.pdgId(),
0201                                       vertex_index,
0202                                       normchi2,
0203                                       dz,
0204                                       dxy,
0205                                       dzError,
0206                                       dxyError,
0207                                       lostInnerHits,
0208                                       quality,
0209                                       trk_pt,
0210                                       trk_eta,
0211                                       trk_phi,
0212                                       relativeTrackVars_);
0213       }
0214     }
0215   }
0216
0217   //get PF jets
0218   Handle<reco::PFJetCollection> pfJetCollection;
0219   auto outPFJets = std::make_unique<Run3ScoutingPFJetCollection>();
0220   if (iEvent.getByToken(pfJetCollection_, pfJetCollection)) {
0221     //get PF jet tags
0222     Handle<reco::JetTagCollection> pfJetTagCollection;
0223     bool haveJetTags = false;
0224     if (doJetTags_ && iEvent.getByToken(pfJetTagCollection_, pfJetTagCollection)) {
0225       haveJetTags = true;
0226     }
0227
0228     for (auto const &jet : *pfJetCollection) {
0229       if (jet.pt() < pfJetPtCut_ || std::abs(jet.eta()) > pfJetEtaCut_)
0230         continue;
0231       //find the jet tag corresponding to the jet
0232       float tagValue = -20;
0233       float minDR2 = 0.01;
0234       if (haveJetTags) {
0235         for (auto const &tag : *pfJetTagCollection) {
0236           float dR2 = reco::deltaR2(jet, *(tag.first));
0237           if (dR2 < minDR2) {
0238             minDR2 = dR2;
0239             tagValue = tag.second;
0240           }
0241         }
0242       }
0243       //get the PF constituents of the jet
0244       std::vector<int> candIndices;
0245       if (doCandidates_ && doCandIndsForJets_) {
0246         for (auto const &cand : jet.getPFConstituents()) {
0247           if (not(cand.isNonnull() and cand.isAvailable())) {
0248             throw cms::Exception("HLTScoutingPFProducer")
0249                 << "invalid reference to reco::PFCandidate from reco::PFJet::getPFConstituents()";
0250           }
0251           if (cand->pt() > pfCandidatePtCut_ && std::abs(cand->eta()) < pfCandidateEtaCut_) {
0252             //search for the candidate in the collection
0253             float minDR2 = 0.0001;
0254             int matchIndex = -1;
0255             int outIndex = 0;
0256             for (auto &outCand : *outPFCandidates) {
0257               auto const dR2 = reco::deltaR2(*cand, outCand);
0258               if (dR2 < minDR2) {
0259                 minDR2 = dR2;
0260                 matchIndex = outIndex;
0261               }
0262               if (minDR2 == 0) {
0263                 break;
0264               }
0265               outIndex++;
0266             }
0267             candIndices.push_back(matchIndex);
0268           }
0269         }
0270       }
0271       outPFJets->emplace_back(
0272           MiniFloatConverter::reduceMantissaToNbitsRounding(jet.pt(), mantissaPrecision_),
0273           MiniFloatConverter::reduceMantissaToNbitsRounding(jet.eta(), mantissaPrecision_),
0274           MiniFloatConverter::reduceMantissaToNbitsRounding(jet.phi(), mantissaPrecision_),
0275           MiniFloatConverter::reduceMantissaToNbitsRounding(jet.mass(), mantissaPrecision_),
0276           MiniFloatConverter::reduceMantissaToNbitsRounding(jet.jetArea(), mantissaPrecision_),
0277           MiniFloatConverter::reduceMantissaToNbitsRounding(jet.chargedHadronEnergy(), mantissaPrecision_),
0278           MiniFloatConverter::reduceMantissaToNbitsRounding(jet.neutralHadronEnergy(), mantissaPrecision_),
0279           MiniFloatConverter::reduceMantissaToNbitsRounding(jet.photonEnergy(), mantissaPrecision_),
0280           MiniFloatConverter::reduceMantissaToNbitsRounding(jet.electronEnergy(), mantissaPrecision_),
0281           MiniFloatConverter::reduceMantissaToNbitsRounding(jet.muonEnergy(), mantissaPrecision_),
0282           MiniFloatConverter::reduceMantissaToNbitsRounding(jet.HFHadronEnergy(), mantissaPrecision_),
0283           MiniFloatConverter::reduceMantissaToNbitsRounding(jet.HFEMEnergy(), mantissaPrecision_),
0284           jet.chargedHadronMultiplicity(),
0285           jet.neutralHadronMultiplicity(),
0286           jet.photonMultiplicity(),
0287           jet.electronMultiplicity(),
0288           jet.muonMultiplicity(),
0289           jet.HFHadronMultiplicity(),
0290           jet.HFEMMultiplicity(),
0291           MiniFloatConverter::reduceMantissaToNbitsRounding(jet.hoEnergy(), mantissaPrecision_),
0292           MiniFloatConverter::reduceMantissaToNbitsRounding(tagValue, mantissaPrecision_),
0293           MiniFloatConverter::reduceMantissaToNbitsRounding(0.0, mantissaPrecision_),
0294           candIndices);
0295     }
0296   }
0297
0298   //put output
0299   iEvent.put(std::move(outPFCandidates));
0300   iEvent.put(std::move(outPFJets));
0301   iEvent.put(std::move(outRho), "rho");
0302   iEvent.put(std::move(outMetPt), "pfMetPt");
0303   iEvent.put(std::move(outMetPhi), "pfMetPhi");
0304 }
0305
0306 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0307 void HLTScoutingPFProducer::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
0308   edm::ParameterSetDescription desc;
0309   desc.add<edm::InputTag>("pfJetCollection", edm::InputTag("hltAK4PFJets"));
0310   desc.add<edm::InputTag>("pfJetTagCollection", edm::InputTag("hltDeepCombinedSecondaryVertexBJetTagsPF"));
0311   desc.add<edm::InputTag>("pfCandidateCollection", edm::InputTag("hltParticleFlow"));
0312   desc.add<edm::InputTag>("vertexCollection", edm::InputTag("hltPixelVertices"));
0313   desc.add<edm::InputTag>("metCollection", edm::InputTag("hltPFMETProducer"));
0314   desc.add<edm::InputTag>("rho", edm::InputTag("hltFixedGridRhoFastjetAll"));
0315   desc.add<double>("pfJetPtCut", 20.0);
0316   desc.add<double>("pfJetEtaCut", 3.0);
0317   desc.add<double>("pfCandidatePtCut", 0.6);
0318   desc.add<double>("pfCandidateEtaCut", 5.0);
0319   desc.add<int>("mantissaPrecision", 10)->setComment("default float16, change to 23 for float32");
0320   desc.add<bool>("doJetTags", true);
0321   desc.add<bool>("doCandidates", true);
0322   desc.add<bool>("doMet", true);
0323   desc.add<bool>("doTrackVars", true);
0324   desc.add<bool>("relativeTrackVars", true);
0325   desc.add<bool>("doCandIndsForJets", false);
0326   descriptions.addWithDefaultLabel(desc);
0327 }
0328
0329 // declare this class as a framework plugin
0330 DEFINE_FWK_MODULE(HLTScoutingPFProducer);