File indexing completed on 2024-09-07 04:37:19
0001
0002 #include <memory>
0003
0004
0005 #include "FWCore/Framework/interface/Frameworkfwd.h"
0006 #include "FWCore/Framework/interface/stream/EDProducer.h"
0007
0008 #include "FWCore/Framework/interface/Event.h"
0009
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0012 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0013
0014 #include "DataFormats/JetReco/interface/GenJet.h"
0015 #include "DataFormats/JetMatching/interface/JetFlavourInfo.h"
0016 #include "DataFormats/JetMatching/interface/JetFlavourInfoMatching.h"
0017 #include "DataFormats/Math/interface/deltaR.h"
0018
0019 #include "DataFormats/NanoAOD/interface/FlatTable.h"
0020
0021 #include "CommonTools/Utils/interface/StringCutObjectSelector.h"
0022 #include "CommonTools/Utils/interface/StringObjectFunction.h"
0023
0024 class GenJetFlavourTableProducer : public edm::stream::EDProducer<> {
0025 public:
0026 explicit GenJetFlavourTableProducer(const edm::ParameterSet& iConfig)
0027 : name_(iConfig.getParameter<std::string>("name")),
0028 src_(consumes<std::vector<reco::GenJet> >(iConfig.getParameter<edm::InputTag>("src"))),
0029 cut_(iConfig.getParameter<std::string>("cut"), true),
0030 deltaR_(iConfig.getParameter<double>("deltaR")),
0031 jetFlavourInfosToken_(
0032 consumes<reco::JetFlavourInfoMatchingCollection>(iConfig.getParameter<edm::InputTag>("jetFlavourInfos"))) {
0033 produces<nanoaod::FlatTable>();
0034 }
0035
0036 ~GenJetFlavourTableProducer() override {}
0037
0038 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0039 edm::ParameterSetDescription desc;
0040 desc.add<edm::InputTag>("src")->setComment("input genJet collection");
0041 desc.add<edm::InputTag>("jetFlavourInfos")->setComment("input flavour info collection");
0042 desc.add<std::string>("name")->setComment("name of the genJet FlatTable we are extending with flavour information");
0043 desc.add<std::string>("cut")->setComment("cut on input genJet collection");
0044 desc.add<double>("deltaR")->setComment("deltaR to match genjets");
0045 descriptions.add("genJetFlavourTable", desc);
0046 }
0047
0048 private:
0049 void produce(edm::Event&, edm::EventSetup const&) override;
0050
0051 std::string name_;
0052 edm::EDGetTokenT<std::vector<reco::GenJet> > src_;
0053 const StringCutObjectSelector<reco::GenJet> cut_;
0054 const double deltaR_;
0055 edm::EDGetTokenT<reco::JetFlavourInfoMatchingCollection> jetFlavourInfosToken_;
0056 };
0057
0058
0059 void GenJetFlavourTableProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0060 const auto& jetsProd = iEvent.get(src_);
0061 const auto& jetFlavourInfosProd = iEvent.get(jetFlavourInfosToken_);
0062
0063 unsigned int ncand = 0;
0064 std::vector<int16_t> partonFlavour;
0065 std::vector<uint8_t> hadronFlavour;
0066 std::vector<uint8_t> nBHadrons;
0067 std::vector<uint8_t> nCHadrons;
0068
0069 for (const reco::GenJet& jet : jetsProd) {
0070 if (!cut_(jet))
0071 continue;
0072 ++ncand;
0073 bool matched = false;
0074 for (const reco::JetFlavourInfoMatching& jetFlavourInfoMatching : jetFlavourInfosProd) {
0075 if (deltaR(jet.p4(), jetFlavourInfoMatching.first->p4()) < deltaR_) {
0076 partonFlavour.push_back(jetFlavourInfoMatching.second.getPartonFlavour());
0077 hadronFlavour.push_back(jetFlavourInfoMatching.second.getHadronFlavour());
0078 nBHadrons.push_back(jetFlavourInfoMatching.second.getbHadrons().size());
0079 nCHadrons.push_back(jetFlavourInfoMatching.second.getcHadrons().size());
0080 matched = true;
0081 break;
0082 }
0083 }
0084 if (!matched) {
0085 partonFlavour.push_back(0);
0086 hadronFlavour.push_back(0);
0087 nBHadrons.push_back(0);
0088 nCHadrons.push_back(0);
0089 }
0090 }
0091
0092 auto tab = std::make_unique<nanoaod::FlatTable>(ncand, name_, false, true);
0093 tab->addColumn<int16_t>("partonFlavour", partonFlavour, "flavour from parton matching");
0094 tab->addColumn<uint8_t>("hadronFlavour", hadronFlavour, "flavour from hadron ghost clustering");
0095 tab->addColumn<uint8_t>("nBHadrons", nBHadrons, "number of b-hadrons");
0096 tab->addColumn<uint8_t>("nCHadrons", nCHadrons, "number of c-hadrons");
0097
0098 iEvent.put(std::move(tab));
0099 }
0100
0101 #include "FWCore/Framework/interface/MakerMacros.h"
0102
0103 DEFINE_FWK_MODULE(GenJetFlavourTableProducer);