File indexing completed on 2023-03-17 11:16:04
0001
0002
0003 #include <memory>
0004
0005
0006 #include "FWCore/Framework/interface/Frameworkfwd.h"
0007 #include "FWCore/Framework/interface/stream/EDProducer.h"
0008
0009 #include "FWCore/Framework/interface/Event.h"
0010 #include "FWCore/Framework/interface/MakerMacros.h"
0011
0012 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0013
0014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0015 #include "FWCore/Utilities/interface/StreamID.h"
0016
0017 #include "DataFormats/JetReco/interface/GenJetCollection.h"
0018 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
0019
0020 #include "DataFormats/Common/interface/ValueMap.h"
0021 #include "CommonTools/Utils/interface/StringCutObjectSelector.h"
0022 #include "CommonTools/Utils/interface/StringObjectFunction.h"
0023
0024
0025
0026
0027
0028 class GenJetGenPartMerger : public edm::stream::EDProducer<> {
0029 public:
0030 explicit GenJetGenPartMerger(const edm::ParameterSet&);
0031 ~GenJetGenPartMerger() override;
0032 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0033
0034 private:
0035 void beginStream(edm::StreamID) override;
0036 void produce(edm::Event&, const edm::EventSetup&) override;
0037 void endStream() override;
0038
0039 const edm::EDGetTokenT<reco::GenJetCollection> jetToken_;
0040 const edm::EDGetTokenT<reco::GenParticleCollection> partToken_;
0041 const StringCutObjectSelector<reco::Candidate> cut_;
0042 const edm::EDGetTokenT<edm::ValueMap<bool>> tauAncToken_;
0043 };
0044
0045
0046
0047
0048
0049
0050
0051
0052
0053
0054
0055
0056 GenJetGenPartMerger::GenJetGenPartMerger(const edm::ParameterSet& iConfig)
0057 : jetToken_(consumes<reco::GenJetCollection>(iConfig.getParameter<edm::InputTag>("srcJet"))),
0058 partToken_(consumes<reco::GenParticleCollection>(iConfig.getParameter<edm::InputTag>("srcPart"))),
0059 cut_(iConfig.getParameter<std::string>("cut")),
0060 tauAncToken_(consumes<edm::ValueMap<bool>>(iConfig.getParameter<edm::InputTag>("hasTauAnc"))) {
0061 produces<reco::GenJetCollection>("merged");
0062 produces<edm::ValueMap<bool>>("hasTauAnc");
0063 }
0064
0065 GenJetGenPartMerger::~GenJetGenPartMerger() {}
0066
0067 void GenJetGenPartMerger::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0068 edm::ParameterSetDescription desc;
0069 desc.add<edm::InputTag>("srcJet")->setComment("reco::GenJetCollection input collection");
0070 desc.add<edm::InputTag>("srcPart")->setComment("reco::GenParticleCollection input collection");
0071 desc.add<std::string>("cut")->setComment("a selection to apply to GenJet");
0072 desc.add<edm::InputTag>("hasTauAnc")->setComment("value map defining GenJet tau origin");
0073 descriptions.addWithDefaultLabel(desc);
0074 }
0075
0076
0077
0078
0079
0080 void GenJetGenPartMerger::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0081 using namespace edm;
0082 auto merged = std::make_unique<reco::GenJetCollection>();
0083
0084 std::vector<bool> hasTauAncValues;
0085
0086 auto jetHandle = iEvent.getHandle(jetToken_);
0087 const auto& partProd = iEvent.get(partToken_);
0088 const auto& tauAncProd = iEvent.get(tauAncToken_);
0089
0090 for (unsigned int ijet = 0; ijet < jetHandle->size(); ++ijet) {
0091 auto jet = jetHandle->at(ijet);
0092 if (cut_(jet)) {
0093 merged->push_back(reco::GenJet(jet));
0094 reco::GenJetRef jetRef(jetHandle, ijet);
0095 hasTauAncValues.push_back(tauAncProd[jetRef]);
0096 }
0097 }
0098
0099 for (const auto& part : partProd) {
0100 reco::GenJet jet;
0101 jet.setP4(part.p4());
0102 jet.setPdgId(part.pdgId());
0103 jet.setCharge(part.charge());
0104 merged->push_back(jet);
0105 hasTauAncValues.push_back(false);
0106 }
0107
0108 auto newmerged = iEvent.put(std::move(merged), "merged");
0109
0110 auto out = std::make_unique<edm::ValueMap<bool>>();
0111 edm::ValueMap<bool>::Filler filler(*out);
0112 filler.insert(newmerged, hasTauAncValues.begin(), hasTauAncValues.end());
0113 filler.fill();
0114
0115 iEvent.put(std::move(out), "hasTauAnc");
0116 }
0117
0118
0119 void GenJetGenPartMerger::beginStream(edm::StreamID) {}
0120
0121
0122 void GenJetGenPartMerger::endStream() {}
0123
0124
0125 DEFINE_FWK_MODULE(GenJetGenPartMerger);