Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-03-17 11:18:32

0001 #include <memory>
0002 
0003 #include "FWCore/Framework/interface/Frameworkfwd.h"
0004 #include "FWCore/Framework/interface/global/EDProducer.h"
0005 
0006 #include "FWCore/Framework/interface/Event.h"
0007 
0008 #include "FWCore/Framework/interface/MakerMacros.h"
0009 
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011 
0012 #include "DataFormats/JetReco/interface/PileupJetIdentifier.h"
0013 #include "DataFormats/JetReco/interface/Jet.h"
0014 #include "DataFormats/JetReco/interface/PFJet.h"
0015 #include "DataFormats/JetReco/interface/PFJetCollection.h"
0016 #include "DataFormats/JetReco/interface/JetCollection.h"
0017 #include "FWCore/Utilities/interface/EDGetToken.h"
0018 
0019 #include "DataFormats/Common/interface/ValueMap.h"
0020 class PUFilter : public edm::global::EDProducer<> {
0021 public:
0022   explicit PUFilter(const edm::ParameterSet&);
0023 
0024   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0025 
0026 private:
0027   const edm::EDGetTokenT<edm::View<reco::PFJet>> jetsToken_;
0028   const edm::EDGetTokenT<edm::ValueMap<int>> jetPuIdToken_;
0029   void produce(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
0030 };
0031 
0032 PUFilter::PUFilter(const edm::ParameterSet& iConfig)
0033     : jetsToken_(consumes<edm::View<reco::PFJet>>(iConfig.getParameter<edm::InputTag>("Jets"))),
0034       jetPuIdToken_(consumes<edm::ValueMap<int>>(iConfig.getParameter<edm::InputTag>("JetPUID"))) {
0035   produces<std::vector<reco::PFJet>>();
0036 }
0037 
0038 void PUFilter::produce(edm::StreamID sid, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
0039   using namespace edm;
0040 
0041   Handle<edm::View<reco::PFJet>> jetsH;
0042   Handle<edm::ValueMap<int>> id_decisions;
0043 
0044   iEvent.getByToken(jetsToken_, jetsH);
0045   iEvent.getByToken(jetPuIdToken_, id_decisions);
0046 
0047   auto goodjets = std::make_unique<std::vector<reco::PFJet>>();
0048   for (size_t i = 0; i < jetsH->size(); ++i) {
0049     auto jet = jetsH->refAt(i);
0050     if ((*id_decisions)[jet])
0051       goodjets->push_back(*jet);
0052   }
0053   iEvent.put(std::move(goodjets));
0054 }
0055 
0056 void PUFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0057   edm::ParameterSetDescription desc;
0058   desc.add<edm::InputTag>("Jets", edm::InputTag("hltAK4PFJetsCorrected"));
0059   desc.add<edm::InputTag>("JetPUID", edm::InputTag("MVAJetPuIdProducer", "CATEv0Id"));
0060   descriptions.add("PUFilter", desc);
0061   desc.setUnknown();
0062 }
0063 DEFINE_FWK_MODULE(PUFilter);