CMS 3D CMS Logo

LHEJetFilter.cc
Go to the documentation of this file.
1 #include <vector>
2 
7 
10 
12 
14 
15 #include "fastjet/PseudoJet.hh"
16 #include "fastjet/JetDefinition.hh"
17 #include "fastjet/ClusterSequence.hh"
18 #include "fastjet/Selector.hh"
19 
20 using namespace std;
21 using namespace fastjet;
22 
24 public:
25  explicit LHEJetFilter(const edm::ParameterSet&);
26  ~LHEJetFilter() override {}
27 
28 private:
29  bool filter(edm::StreamID strid, edm::Event& evt, const edm::EventSetup& params) const override;
30 
32  double jetPtMin_;
33  JetDefinition jetdef_;
34 };
35 
37  : tokenLHEEvent_(consumes<LHEEventProduct>(params.getParameter<edm::InputTag>("src"))),
38  jetPtMin_(params.getParameter<double>("jetPtMin")),
39  jetdef_(antikt_algorithm, params.getParameter<double>("jetR")) {}
40 
41 bool LHEJetFilter::filter(edm::StreamID strid, edm::Event& evt, const edm::EventSetup& params) const {
43  evt.getByToken(tokenLHEEvent_, lheinfo);
44 
45  if (!lheinfo.isValid()) {
46  return true;
47  }
48 
49  vector<PseudoJet> jetconsts;
50  jetconsts.reserve(10);
51  const lhef::HEPEUP& hepeup = lheinfo->hepeup();
52  for (size_t p = 0; p < hepeup.IDUP.size(); ++p) {
53  if (hepeup.ISTUP[p] == 1) {
54  const lhef::HEPEUP::FiveVector& mom = hepeup.PUP[p];
55  jetconsts.emplace_back(mom[0], mom[1], mom[2], mom[3]);
56  }
57  }
58 
59  ClusterSequence cs(jetconsts, jetdef_);
60  vector<PseudoJet> jets = cs.inclusive_jets(jetPtMin_);
61 
62  return !jets.empty();
63 }
64 
65 // Define module as a plug-in
const lhef::HEPEUP & hepeup() const
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:517
unique_ptr< ClusterSequence > cs
edm::EDGetTokenT< LHEEventProduct > tokenLHEEvent_
Definition: LHEJetFilter.cc:31
JetDefinition jetdef_
Definition: LHEJetFilter.cc:33
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
double jetPtMin_
Definition: LHEJetFilter.cc:32
LHEJetFilter(const edm::ParameterSet &)
Definition: LHEJetFilter.cc:36
std::vector< FiveVector > PUP
Definition: LesHouches.h:261
vector< PseudoJet > jets
std::vector< int > ISTUP
Definition: LesHouches.h:243
bool isValid() const
Definition: HandleBase.h:74
std::vector< int > IDUP
Definition: LesHouches.h:238
bool filter(edm::StreamID strid, edm::Event &evt, const edm::EventSetup &params) const override
Definition: LHEJetFilter.cc:41
HLT enums.
~LHEJetFilter() override
Definition: LHEJetFilter.cc:26