File indexing completed on 2024-04-06 12:18:32
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018 #include <memory>
0019
0020
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
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
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
0109 void HLTScoutingPFProducer::produce(edm::StreamID sid, edm::Event &iEvent, edm::EventSetup const &setup) const {
0110 using namespace edm;
0111
0112
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(
0118 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.x(), mantissaPrecision_),
0119 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.y(), mantissaPrecision_),
0120 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.z(), mantissaPrecision_),
0121 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.zError(), mantissaPrecision_),
0122 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.xError(), mantissaPrecision_),
0123 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.yError(), mantissaPrecision_),
0124 vtx.tracksSize(),
0125 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.chi2(), mantissaPrecision_),
0126 vtx.ndof(),
0127 vtx.isValid(),
0128 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.covariance(0, 1), mantissaPrecision_),
0129 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.covariance(0, 2), mantissaPrecision_),
0130 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.covariance(1, 2), mantissaPrecision_));
0131 }
0132 }
0133
0134
0135 Handle<double> rho;
0136 auto outRho = std::make_unique<double>(-999);
0137 if (iEvent.getByToken(rho_, rho)) {
0138 *outRho = *rho;
0139 }
0140
0141
0142 Handle<reco::PFMETCollection> metCollection;
0143 auto outMetPt = std::make_unique<double>(-999);
0144 auto outMetPhi = std::make_unique<double>(-999);
0145 if (doMet_ && iEvent.getByToken(metCollection_, metCollection)) {
0146 outMetPt = std::make_unique<double>(metCollection->front().pt());
0147 outMetPhi = std::make_unique<double>(metCollection->front().phi());
0148 }
0149
0150
0151 Handle<reco::PFCandidateCollection> pfCandidateCollection;
0152 auto outPFCandidates = std::make_unique<Run3ScoutingParticleCollection>();
0153 if (doCandidates_ && iEvent.getByToken(pfCandidateCollection_, pfCandidateCollection)) {
0154 for (auto const &cand : *pfCandidateCollection) {
0155 if (cand.pt() > pfCandidatePtCut_ && std::abs(cand.eta()) < pfCandidateEtaCut_) {
0156 int vertex_index = -1;
0157 int index_counter = 0;
0158 double dr2 = 0.0001;
0159 for (auto const &vtx : *outVertices) {
0160 double tmp_dr2 = pow(vtx.x() - cand.vx(), 2) + pow(vtx.y() - cand.vy(), 2) + pow(vtx.z() - cand.vz(), 2);
0161 if (tmp_dr2 < dr2) {
0162 dr2 = tmp_dr2;
0163 vertex_index = index_counter;
0164 }
0165 if (dr2 == 0.0)
0166 break;
0167 ++index_counter;
0168 }
0169 float normchi2{0}, dz{0}, dxy{0}, dzError{0}, dxyError{0}, trk_pt{0}, trk_eta{0}, trk_phi{0};
0170 uint8_t lostInnerHits{0}, quality{0};
0171 if (doTrackVars_) {
0172 const auto *trk = cand.bestTrack();
0173 if (trk != nullptr) {
0174 normchi2 = MiniFloatConverter::reduceMantissaToNbitsRounding(trk->normalizedChi2(), mantissaPrecision_);
0175 lostInnerHits = btagbtvdeep::lost_inner_hits_from_pfcand(cand);
0176 quality = btagbtvdeep::quality_from_pfcand(cand);
0177 if (relativeTrackVars_) {
0178 trk_pt = MiniFloatConverter::reduceMantissaToNbitsRounding(trk->pt() - cand.pt(), mantissaPrecision_);
0179 trk_eta = MiniFloatConverter::reduceMantissaToNbitsRounding(trk->eta() - cand.eta(), mantissaPrecision_);
0180 trk_phi = MiniFloatConverter::reduceMantissaToNbitsRounding(trk->phi() - cand.phi(), mantissaPrecision_);
0181 } else {
0182 trk_pt = MiniFloatConverter::reduceMantissaToNbitsRounding(trk->pt(), mantissaPrecision_);
0183 trk_eta = MiniFloatConverter::reduceMantissaToNbitsRounding(trk->eta(), mantissaPrecision_);
0184 trk_phi = MiniFloatConverter::reduceMantissaToNbitsRounding(trk->phi(), mantissaPrecision_);
0185 }
0186 if (not vertexCollection->empty()) {
0187 const reco::Vertex &pv = (*vertexCollection)[0];
0188
0189 dz = trk->dz(pv.position());
0190 dzError = MiniFloatConverter::reduceMantissaToNbitsRounding(dz / trk->dzError(), mantissaPrecision_);
0191 dz = MiniFloatConverter::reduceMantissaToNbitsRounding(dz, mantissaPrecision_);
0192
0193 dxy = trk->dxy(pv.position());
0194 dxyError = MiniFloatConverter::reduceMantissaToNbitsRounding(dxy / trk->dxyError(), mantissaPrecision_);
0195 dxy = MiniFloatConverter::reduceMantissaToNbitsRounding(dxy, mantissaPrecision_);
0196 }
0197 } else {
0198 normchi2 = MiniFloatConverter::reduceMantissaToNbitsRounding(999, mantissaPrecision_);
0199 }
0200 }
0201 outPFCandidates->emplace_back(MiniFloatConverter::reduceMantissaToNbitsRounding(cand.pt(), mantissaPrecision_),
0202 MiniFloatConverter::reduceMantissaToNbitsRounding(cand.eta(), mantissaPrecision_),
0203 MiniFloatConverter::reduceMantissaToNbitsRounding(cand.phi(), mantissaPrecision_),
0204 cand.pdgId(),
0205 vertex_index,
0206 normchi2,
0207 dz,
0208 dxy,
0209 dzError,
0210 dxyError,
0211 lostInnerHits,
0212 quality,
0213 trk_pt,
0214 trk_eta,
0215 trk_phi,
0216 relativeTrackVars_);
0217 }
0218 }
0219 }
0220
0221
0222 Handle<reco::PFJetCollection> pfJetCollection;
0223 auto outPFJets = std::make_unique<Run3ScoutingPFJetCollection>();
0224 if (iEvent.getByToken(pfJetCollection_, pfJetCollection)) {
0225
0226 Handle<reco::JetTagCollection> pfJetTagCollection;
0227 bool haveJetTags = false;
0228 if (doJetTags_ && iEvent.getByToken(pfJetTagCollection_, pfJetTagCollection)) {
0229 haveJetTags = true;
0230 }
0231
0232 for (auto const &jet : *pfJetCollection) {
0233 if (jet.pt() < pfJetPtCut_ || std::abs(jet.eta()) > pfJetEtaCut_)
0234 continue;
0235
0236 float tagValue = -20;
0237 float minDR2 = 0.01;
0238 if (haveJetTags) {
0239 for (auto const &tag : *pfJetTagCollection) {
0240 float dR2 = reco::deltaR2(jet, *(tag.first));
0241 if (dR2 < minDR2) {
0242 minDR2 = dR2;
0243 tagValue = tag.second;
0244 }
0245 }
0246 }
0247
0248 std::vector<int> candIndices;
0249 if (doCandidates_ && doCandIndsForJets_) {
0250 for (auto const &cand : jet.getPFConstituents()) {
0251 if (not(cand.isNonnull() and cand.isAvailable())) {
0252 throw cms::Exception("HLTScoutingPFProducer")
0253 << "invalid reference to reco::PFCandidate from reco::PFJet::getPFConstituents()";
0254 }
0255 if (cand->pt() > pfCandidatePtCut_ && std::abs(cand->eta()) < pfCandidateEtaCut_) {
0256
0257 float minDR2 = 0.0001;
0258 int matchIndex = -1;
0259 int outIndex = 0;
0260 for (auto &outCand : *outPFCandidates) {
0261 auto const dR2 = reco::deltaR2(*cand, outCand);
0262 if (dR2 < minDR2) {
0263 minDR2 = dR2;
0264 matchIndex = outIndex;
0265 }
0266 if (minDR2 == 0) {
0267 break;
0268 }
0269 outIndex++;
0270 }
0271 candIndices.push_back(matchIndex);
0272 }
0273 }
0274 }
0275 outPFJets->emplace_back(
0276 MiniFloatConverter::reduceMantissaToNbitsRounding(jet.pt(), mantissaPrecision_),
0277 MiniFloatConverter::reduceMantissaToNbitsRounding(jet.eta(), mantissaPrecision_),
0278 MiniFloatConverter::reduceMantissaToNbitsRounding(jet.phi(), mantissaPrecision_),
0279 MiniFloatConverter::reduceMantissaToNbitsRounding(jet.mass(), mantissaPrecision_),
0280 MiniFloatConverter::reduceMantissaToNbitsRounding(jet.jetArea(), mantissaPrecision_),
0281 MiniFloatConverter::reduceMantissaToNbitsRounding(jet.chargedHadronEnergy(), mantissaPrecision_),
0282 MiniFloatConverter::reduceMantissaToNbitsRounding(jet.neutralHadronEnergy(), mantissaPrecision_),
0283 MiniFloatConverter::reduceMantissaToNbitsRounding(jet.photonEnergy(), mantissaPrecision_),
0284 MiniFloatConverter::reduceMantissaToNbitsRounding(jet.electronEnergy(), mantissaPrecision_),
0285 MiniFloatConverter::reduceMantissaToNbitsRounding(jet.muonEnergy(), mantissaPrecision_),
0286 MiniFloatConverter::reduceMantissaToNbitsRounding(jet.HFHadronEnergy(), mantissaPrecision_),
0287 MiniFloatConverter::reduceMantissaToNbitsRounding(jet.HFEMEnergy(), mantissaPrecision_),
0288 jet.chargedHadronMultiplicity(),
0289 jet.neutralHadronMultiplicity(),
0290 jet.photonMultiplicity(),
0291 jet.electronMultiplicity(),
0292 jet.muonMultiplicity(),
0293 jet.HFHadronMultiplicity(),
0294 jet.HFEMMultiplicity(),
0295 MiniFloatConverter::reduceMantissaToNbitsRounding(jet.hoEnergy(), mantissaPrecision_),
0296 MiniFloatConverter::reduceMantissaToNbitsRounding(tagValue, mantissaPrecision_),
0297 MiniFloatConverter::reduceMantissaToNbitsRounding(0.0, mantissaPrecision_),
0298 candIndices);
0299 }
0300 }
0301
0302
0303 iEvent.put(std::move(outPFCandidates));
0304 iEvent.put(std::move(outPFJets));
0305 iEvent.put(std::move(outRho), "rho");
0306 iEvent.put(std::move(outMetPt), "pfMetPt");
0307 iEvent.put(std::move(outMetPhi), "pfMetPhi");
0308 }
0309
0310
0311 void HLTScoutingPFProducer::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
0312 edm::ParameterSetDescription desc;
0313 desc.add<edm::InputTag>("pfJetCollection", edm::InputTag("hltAK4PFJets"));
0314 desc.add<edm::InputTag>("pfJetTagCollection", edm::InputTag("hltDeepCombinedSecondaryVertexBJetTagsPF"));
0315 desc.add<edm::InputTag>("pfCandidateCollection", edm::InputTag("hltParticleFlow"));
0316 desc.add<edm::InputTag>("vertexCollection", edm::InputTag("hltPixelVertices"));
0317 desc.add<edm::InputTag>("metCollection", edm::InputTag("hltPFMETProducer"));
0318 desc.add<edm::InputTag>("rho", edm::InputTag("hltFixedGridRhoFastjetAll"));
0319 desc.add<double>("pfJetPtCut", 20.0);
0320 desc.add<double>("pfJetEtaCut", 3.0);
0321 desc.add<double>("pfCandidatePtCut", 0.6);
0322 desc.add<double>("pfCandidateEtaCut", 5.0);
0323 desc.add<int>("mantissaPrecision", 10)->setComment("default float16, change to 23 for float32");
0324 desc.add<bool>("doJetTags", true);
0325 desc.add<bool>("doCandidates", true);
0326 desc.add<bool>("doMet", true);
0327 desc.add<bool>("doTrackVars", true);
0328 desc.add<bool>("relativeTrackVars", true);
0329 desc.add<bool>("doCandIndsForJets", false);
0330 descriptions.addWithDefaultLabel(desc);
0331 }
0332
0333
0334 DEFINE_FWK_MODULE(HLTScoutingPFProducer);