GenFilters/plugins/LHEJetFilter.cc

0001 #include <vector>
0002
0003 #include "FWCore/Framework/interface/global/EDFilter.h"
0004 #include "FWCore/Framework/interface/Frameworkfwd.h"
0005 #include "FWCore/Framework/interface/Event.h"
0006 #include "FWCore/Framework/interface/MakerMacros.h"
0007
0008 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0009 #include "FWCore/Utilities/interface/InputTag.h"
0010
0011 #include "DataFormats/Common/interface/Handle.h"
0012
0013 #include "SimDataFormats/GeneratorProducts/interface/LHEEventProduct.h"
0014
0015 #include "fastjet/PseudoJet.hh"
0016 #include "fastjet/JetDefinition.hh"
0017 #include "fastjet/ClusterSequence.hh"
0018 #include "fastjet/Selector.hh"
0019
0020 using namespace std;
0021 using namespace fastjet;
0022
0023 class LHEJetFilter : public edm::global::EDFilter<> {
0024 public:
0025   explicit LHEJetFilter(const edm::ParameterSet&);
0026   ~LHEJetFilter() override {}
0027
0028 private:
0029   bool filter(edm::StreamID strid, edm::Event& evt, const edm::EventSetup& params) const override;
0030
0031   edm::EDGetTokenT<LHEEventProduct> tokenLHEEvent_;
0032   double jetPtMin_;
0033   JetDefinition jetdef_;
0034 };
0035
0036 LHEJetFilter::LHEJetFilter(const edm::ParameterSet& params)
0037     : tokenLHEEvent_(consumes<LHEEventProduct>(params.getParameter<edm::InputTag>("src"))),
0038       jetPtMin_(params.getParameter<double>("jetPtMin")),
0039       jetdef_(antikt_algorithm, params.getParameter<double>("jetR")) {}
0040
0041 bool LHEJetFilter::filter(edm::StreamID strid, edm::Event& evt, const edm::EventSetup& params) const {
0042   edm::Handle<LHEEventProduct> lheinfo;
0043   evt.getByToken(tokenLHEEvent_, lheinfo);
0044
0045   if (!lheinfo.isValid()) {
0046     return true;
0047   }
0048
0049   vector<PseudoJet> jetconsts;
0050   jetconsts.reserve(10);
0051   const lhef::HEPEUP& hepeup = lheinfo->hepeup();
0052   for (size_t p = 0; p < hepeup.IDUP.size(); ++p) {
0053     if (hepeup.ISTUP[p] == 1) {
0054       const lhef::HEPEUP::FiveVector& mom = hepeup.PUP[p];
0055       jetconsts.emplace_back(mom[0], mom[1], mom[2], mom[3]);
0056     }
0057   }
0058
0059   ClusterSequence cs(jetconsts, jetdef_);
0060   vector<PseudoJet> jets = cs.inclusive_jets(jetPtMin_);
0061
0062   return !jets.empty();
0063 }
0064
0065 // Define module as a plug-in
0066 DEFINE_FWK_MODULE(LHEJetFilter);