CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
HLTPFJetIDProducer.cc
Go to the documentation of this file.
1 
10 
15 
16 // Constructor
18  : minPt_(iConfig.getParameter<double>("minPt")),
19  maxEta_(iConfig.getParameter<double>("maxEta")),
20  CHF_(iConfig.getParameter<double>("CHF")),
21  NHF_(iConfig.getParameter<double>("NHF")),
22  CEF_(iConfig.getParameter<double>("CEF")),
23  NEF_(iConfig.getParameter<double>("NEF")),
24  maxCF_(iConfig.getParameter<double>("maxCF")),
25  NCH_(iConfig.getParameter<int>("NCH")),
26  NTOT_(iConfig.getParameter<int>("NTOT")),
27  inputTag_(iConfig.getParameter<edm::InputTag>("jetsInput")) {
28  m_thePFJetToken = consumes<reco::PFJetCollection>(inputTag_);
29 
30  // Register the products
31  produces<reco::PFJetCollection>();
32 }
33 
34 // Destructor
36 
37 // Fill descriptions
40  desc.add<double>("minPt", 20.);
41  desc.add<double>("maxEta", 1e99);
42  desc.add<double>("CHF", -99.);
43  desc.add<double>("NHF", 99.);
44  desc.add<double>("CEF", 99.);
45  desc.add<double>("NEF", 99.);
46  desc.add<double>("maxCF", 99.);
47  desc.add<int>("NCH", 0);
48  desc.add<int>("NTOT", 0);
49  desc.add<edm::InputTag>("jetsInput", edm::InputTag("hltAntiKT4PFJets"));
50  descriptions.add("hltPFJetIDProducer", desc);
51 }
52 
53 // Produce the products
55  // Create a pointer to the products
56  std::unique_ptr<reco::PFJetCollection> result(new reco::PFJetCollection());
57 
59  iEvent.getByToken(m_thePFJetToken, pfjets);
60 
61  for (auto const& j : *pfjets) {
62  bool pass = false;
63  double pt = j.pt();
64  double eta = j.eta();
65  double abseta = std::abs(eta);
66 
67  if (!(pt > 0.))
68  continue; // skip jets with zero or negative pt
69 
70  if (pt < minPt_) {
71  pass = true;
72 
73  } else if (abseta >= maxEta_) {
74  pass = true;
75 
76  } else {
77  double chf = j.chargedHadronEnergyFraction();
78  //double nhf = j->neutralHadronEnergyFraction() + j->HFHadronEnergyFraction();
79  double nhf = j.neutralHadronEnergyFraction();
80  double cef = j.chargedEmEnergyFraction();
81  double nef = j.neutralEmEnergyFraction();
82  double cftot = chf + cef + j.chargedMuEnergyFraction();
83  int nch = j.chargedMultiplicity();
84  int ntot = j.numberOfDaughters();
85 
86  pass = true;
87  pass = pass && (ntot > NTOT_);
88  pass = pass && (nef < NEF_);
89  pass = pass && (nhf < NHF_ || abseta >= 2.4); //NHF-cut does not work in HF anymore with recent PF
90  pass = pass && (cef < CEF_ || abseta >= 2.4);
91  pass = pass && (chf > CHF_ || abseta >= 2.4);
92  pass = pass && (nch > NCH_ || abseta >= 2.4);
93  pass = pass && (cftot < maxCF_ || abseta >= 2.4);
94  }
95 
96  if (pass)
97  result->push_back(j);
98  }
99 
100  // Put the products into the Event
101  iEvent.put(std::move(result));
102 }
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
Definition: Event.h:133
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:539
HLTPFJetIDProducer(const edm::ParameterSet &iConfig)
~HLTPFJetIDProducer() override
tuple result
Definition: mps_fire.py:311
int iEvent
Definition: GenABIO.cc:224
edm::EDGetTokenT< reco::PFJetCollection > m_thePFJetToken
def move
Definition: eostools.py:511
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
ParameterDescriptionBase * add(U const &iLabel, T const &value)
double NEF_
neutral EM fraction
int NTOT_
number of constituents
double CHF_
charged hadron fraction
int NCH_
number of charged constituents
edm::InputTag inputTag_
input PFJet collection
void add(std::string const &label, ParameterSetDescription const &psetDescription)
std::vector< PFJet > PFJetCollection
collection of PFJet objects
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
void produce(edm::Event &iEvent, const edm::EventSetup &iSetup) override
__shared__ uint32_t ntot