CMS 3D CMS Logo

PFJetIDSelectionFunctor.h
Go to the documentation of this file.
1 #ifndef PhysicsTools_SelectorUtils_interface_PFJetIDSelectionFunctor_h
2 #define PhysicsTools_SelectorUtils_interface_PFJetIDSelectionFunctor_h
3 
15 #ifndef __GCCXML__
17 #endif
20 
22 
23 #include <TMath.h>
24 class PFJetIDSelectionFunctor : public Selector<pat::Jet> {
25 
26  public: // interface
27 
30 
32 
33 #ifndef __GCCXML__
36  {}
37 #endif
38 
40  {
41  std::string versionStr = params.getParameter<std::string>("version");
42  std::string qualityStr = params.getParameter<std::string>("quality");
43 
44  if ( versionStr == "FIRSTDATA" )
46  else if( versionStr == "RUNIISTARTUP")
48  // WINTER16 implements most recent (as of Feb 2017) JetID criteria
49  // See: https://twiki.cern.ch/twiki/bin/view/CMS/JetID13TeVRun2016
50  else if( versionStr == "WINTER16")
51  version_ = WINTER16;
52  else version_ = WINTER16;//set WINTER16 as default
53 
54 
55  if ( qualityStr == "LOOSE") quality_ = LOOSE;
56  else if ( qualityStr == "TIGHT") quality_ = TIGHT;
57  else quality_ = LOOSE;
58 
59  push_back("CHF" );
60  push_back("NHF" );
61  push_back("CEF" );
62  push_back("NEF" );
63  push_back("NCH" );
64  push_back("nConstituents");
65  if(version_ == RUNIISTARTUP ){
66  push_back("NEF_FW");
67  push_back("nNeutrals_FW");
68  }
69  if(version_ == WINTER16 ){
70  push_back("NHF_EC");
71  push_back("NEF_EC");
72  push_back("nNeutrals_EC");
73  push_back("NEF_FW");
74  push_back("nNeutrals_FW");
75  }
76 
77 
78 
79  // Set some default cuts for LOOSE, TIGHT
80  if ( quality_ == LOOSE ) {
81  set("CHF", 0.0);
82  set("NHF", 0.99);
83  set("CEF", 0.99);
84  set("NEF", 0.99);
85  set("NCH", 0);
86  set("nConstituents", 1);
87  if(version_ == RUNIISTARTUP){
88  set("NEF_FW",0.90);
89  set("nNeutrals_FW",10);
90  }
91  if(version_ == WINTER16){
92  set("NHF_EC",0.98);
93  set("NEF_EC",0.01);
94  set("nNeutrals_EC",2);
95  set("NEF_FW",0.90);
96  set("nNeutrals_FW",10);
97  }
98 
99 
100  } else if ( quality_ == TIGHT ) {
101  set("CHF", 0.0);
102  set("NHF", 0.9);
103  set("CEF", 0.99);
104  set("NEF", 0.9);
105  set("NCH", 0);
106  set("nConstituents", 1);
107  if(version_ == RUNIISTARTUP){
108  set("NEF_FW",0.90);
109  set("nNeutrals_FW",10);
110  }
111  if(version_ == WINTER16){
112  set("NHF_EC",0.98);
113  set("NEF_EC",0.01);
114  set("nNeutrals_EC",2);
115  set("NEF_FW",0.90);
116  set("nNeutrals_FW",10);
117  }
118 
119  }
120 
121 
122  // Now check the configuration to see if the user changed anything
123  if ( params.exists("CHF") ) set("CHF", params.getParameter<double>("CHF") );
124  if ( params.exists("NHF") ) set("NHF", params.getParameter<double>("NHF") );
125  if ( params.exists("CEF") ) set("CEF", params.getParameter<double>("CEF") );
126  if ( params.exists("NEF") ) set("NEF", params.getParameter<double>("NEF") );
127  if ( params.exists("NCH") ) set("NCH", params.getParameter<int> ("NCH") );
128  if ( params.exists("nConstituents") ) set("nConstituents", params.getParameter<int> ("nConstituents") );
129  if(version_ == RUNIISTARTUP){
130  if ( params.exists("NEF_FW") ) set("NEF_FW", params.getParameter<double> ("NEF_FW") );
131  if ( params.exists("nNeutrals_FW") ) set("nNeutrals_FW", params.getParameter<int> ("nNeutrals_FW") );
132  }
133  if(version_ == WINTER16){
134  if ( params.exists("NHF_EC") ) set("NHF_EC", params.getParameter<int> ("NHF_EC") );
135  if ( params.exists("NEF_EC") ) set("NEF_EC", params.getParameter<int> ("NEF_EC") );
136  if ( params.exists("nNeutrals_EC") ) set("nNeutrals_EC", params.getParameter<int> ("nNeutrals_EC") );
137  if ( params.exists("NEF_FW") ) set("NEF_FW", params.getParameter<double> ("NEF_FW") );
138  if ( params.exists("nNeutrals_FW") ) set("nNeutrals_FW", params.getParameter<int> ("nNeutrals_FW") );
139  }
140 
141 
142  if ( params.exists("cutsToIgnore") )
143  setIgnoredCuts( params.getParameter<std::vector<std::string> >("cutsToIgnore") );
144 
145 
146  indexNConstituents_ = index_type (&bits_, "nConstituents");
147  indexNEF_ = index_type (&bits_, "NEF");
148  indexNHF_ = index_type (&bits_, "NHF");
149  indexCEF_ = index_type (&bits_, "CEF");
150  indexCHF_ = index_type (&bits_, "CHF");
151  indexNCH_ = index_type (&bits_, "NCH");
152  if(version_ == RUNIISTARTUP){
153  indexNEF_FW_ = index_type (&bits_, "NEF_FW");
154  indexNNeutrals_FW_ = index_type (&bits_, "nNeutrals_FW");
155  }
156  if(version_ == WINTER16){
157  indexNHF_EC_ = index_type (&bits_, "NHF_EC");
158  indexNEF_EC_ = index_type (&bits_, "NEF_EC");
159  indexNNeutrals_EC_ = index_type (&bits_, "nNeutrals_EC");
160  indexNEF_FW_ = index_type (&bits_, "NEF_FW");
161  indexNNeutrals_FW_ = index_type (&bits_, "nNeutrals_FW");
162  }
163 
165 
166  }
167 
168 
170  Quality_t quality ) :
171  version_(version), quality_(quality)
172  {
173 
174  push_back("CHF" );
175  push_back("NHF" );
176  push_back("CEF" );
177  push_back("NEF" );
178  push_back("NCH" );
179  push_back("nConstituents");
180  if(version_ == RUNIISTARTUP){
181  push_back("NEF_FW");
182  push_back("nNeutrals_FW");
183  }
184  if(version_ == WINTER16){
185  push_back("NHF_EC");
186  push_back("NEF_EC");
187  push_back("nNeutrals_EC");
188  push_back("NEF_FW");
189  push_back("nNeutrals_FW");
190  }
191 
192 
193  // Set some default cuts for LOOSE, TIGHT
194  if ( quality_ == LOOSE ) {
195  set("CHF", 0.0);
196  set("NHF", 0.99);
197  set("CEF", 0.99);
198  set("NEF", 0.99);
199  set("NCH", 0);
200  set("nConstituents", 1);
201  if(version_ == RUNIISTARTUP){
202  set("NEF_FW",0.90);
203  set("nNeutrals_FW",10);
204  }
205  if(version_ == WINTER16){
206  set("NHF_EC",0.98);
207  set("NEF_EC",0.01);
208  set("nNeutrals_EC",2);
209  set("NEF_FW",0.90);
210  set("nNeutrals_FW",10);
211  }
212 
213  } else if ( quality_ == TIGHT ) {
214  set("CHF", 0.0);
215  set("NHF", 0.9);
216  set("CEF", 0.99);
217  set("NEF", 0.9);
218  set("NCH", 0);
219  set("nConstituents", 1);
220  if(version_ == RUNIISTARTUP){
221  set("NEF_FW",0.90);
222  set("nNeutrals_FW",10);
223  }
224  if(version_ == WINTER16){
225  set("NHF_EC",0.98);
226  set("NEF_EC",0.01);
227  set("nNeutrals_EC",2);
228  set("NEF_FW",0.90);
229  set("nNeutrals_FW",10);
230  }
231  }
232 
233 
234  indexNConstituents_ = index_type (&bits_, "nConstituents");
235  indexNEF_ = index_type (&bits_, "NEF");
236  indexNHF_ = index_type (&bits_, "NHF");
237  indexCEF_ = index_type (&bits_, "CEF");
238  indexCHF_ = index_type (&bits_, "CHF");
239  indexNCH_ = index_type (&bits_, "NCH");
240  if(version_ == RUNIISTARTUP){
241  indexNEF_FW_ = index_type (&bits_, "NEF_FW");
242  indexNNeutrals_FW_ = index_type (&bits_, "nNeutrals_FW");
243  }
244  if(version_ == WINTER16){
245  indexNHF_EC_ = index_type (&bits_, "NHF_EC");
246  indexNEF_EC_ = index_type (&bits_, "NEF_EC");
247  indexNNeutrals_EC_ = index_type (&bits_, "nNeutrals_EC");
248  indexNEF_FW_ = index_type (&bits_, "NEF_FW");
249  indexNNeutrals_FW_ = index_type (&bits_, "nNeutrals_FW");
250  }
251 
252 
254  }
255 
256 
257  //
258  // Accessor from PAT jets
259  //
260  bool operator()( const pat::Jet & jet, pat::strbitset & ret )
261  {
263  if ( jet.currentJECLevel() == "Uncorrected" || !jet.jecSetsAvailable() )
264  return firstDataCuts( jet, ret, version_);
265  else
266  return firstDataCuts( jet.correctedJet("Uncorrected"), ret, version_ );
267  }
268  else {
269  return false;
270  }
271  }
273 
274  //
275  // Accessor from *CORRECTED* 4-vector, EMF, and Jet ID.
276  // This can be used with reco quantities.
277  //
278  bool operator()( const reco::PFJet & jet, pat::strbitset & ret )
279  {
280  if ( version_ == FIRSTDATA || version_ == RUNIISTARTUP || version_ == WINTER16 ){ return firstDataCuts( jet, ret, version_);
281  }
282  else {
283  return false;
284  }
285  }
286 
287  bool operator()( const reco::PFJet & jet )
288  {
289  retInternal_.set(false);
290  operator()(jet, retInternal_);
292  return (bool)retInternal_;
293  }
294 
295  //
296  // cuts based on craft 08 analysis.
297  //
298  bool firstDataCuts( reco::Jet const & jet,
300  {
301  ret.set(false);
302 
303  // cache some variables
304  double chf = 0.0;
305  double nhf = 0.0;
306  double cef = 0.0;
307  double nef = 0.0;
308  int nch = 0;
309  int nconstituents = 0;
310  int nneutrals = 0;
311 
312  // Have to do this because pat::Jet inherits from reco::Jet but not reco::PFJet
313  reco::PFJet const * pfJet = dynamic_cast<reco::PFJet const *>(&jet);
314  pat::Jet const * patJet = dynamic_cast<pat::Jet const *>(&jet);
315  reco::BasicJet const * basicJet = dynamic_cast<reco::BasicJet const *>(&jet);
316 
317  if ( patJet != 0 ) {
318  if ( patJet->isPFJet() ) {
319  chf = patJet->chargedHadronEnergyFraction();
320  nhf = patJet->neutralHadronEnergyFraction();
321  cef = patJet->chargedEmEnergyFraction();
322  nef = patJet->neutralEmEnergyFraction();
323  nch = patJet->chargedMultiplicity();
324  nconstituents = patJet->numberOfDaughters();
325  nneutrals = patJet->neutralMultiplicity();
326  }
327  // Handle the special case where this is a composed jet for
328  // subjet analyses
329  else if ( patJet->isBasicJet() ) {
330  double e_chf = 0.0;
331  double e_nhf = 0.0;
332  double e_cef = 0.0;
333  double e_nef = 0.0;
334  nch = 0;
335  nconstituents = 0;
336  nneutrals = 0;
337 
338  for ( reco::Jet::const_iterator ibegin = patJet->begin(),
339  iend = patJet->end(), isub = ibegin;
340  isub != iend; ++isub ) {
341  reco::PFJet const * pfsub = dynamic_cast<reco::PFJet const *>( &*isub );
342  e_chf += pfsub->chargedHadronEnergy();
343  e_nhf += pfsub->neutralHadronEnergy();
344  e_cef += pfsub->chargedEmEnergy();
345  e_nef += pfsub->neutralEmEnergy();
346  nch += pfsub->chargedMultiplicity();
347  nconstituents += pfsub->numberOfDaughters();
348  nneutrals += pfsub->neutralMultiplicity();
349  }
350  double e = patJet->energy();
351  if ( e > 0.000001 ) {
352  chf = e_chf / e;
353  nhf = e_nhf / e;
354  cef = e_cef / e;
355  nef = e_nef / e;
356  } else {
357  chf = nhf = cef = nef = 0.0;
358  }
359  }
360  } // end if pat jet
361  else if ( pfJet != 0 ) {
362  // CV: need to compute energy fractions in a way that works for corrected as well as for uncorrected PFJets
363  double jetEnergyUncorrected =
364  pfJet->chargedHadronEnergy()
365  + pfJet->neutralHadronEnergy()
366  + pfJet->photonEnergy()
367  + pfJet->electronEnergy()
368  + pfJet->muonEnergy()
369  + pfJet->HFEMEnergy();
370  if ( jetEnergyUncorrected > 0. ) {
371  chf = pfJet->chargedHadronEnergy() / jetEnergyUncorrected;
372  nhf = pfJet->neutralHadronEnergy() / jetEnergyUncorrected;
373  cef = pfJet->chargedEmEnergy() / jetEnergyUncorrected;
374  nef = pfJet->neutralEmEnergy() / jetEnergyUncorrected;
375  }
376  nch = pfJet->chargedMultiplicity();
377  nconstituents = pfJet->numberOfDaughters();
378  nneutrals = pfJet->neutralMultiplicity();
379  } // end if PF jet
380  // Handle the special case where this is a composed jet for
381  // subjet analyses
382  else if ( basicJet != 0 ) {
383  double e_chf = 0.0;
384  double e_nhf = 0.0;
385  double e_cef = 0.0;
386  double e_nef = 0.0;
387  nch = 0;
388  nconstituents = 0;
389  for ( reco::Jet::const_iterator ibegin = basicJet->begin(),
390  iend = patJet->end(), isub = ibegin;
391  isub != iend; ++isub ) {
392  reco::PFJet const * pfsub = dynamic_cast<reco::PFJet const *>( &*isub );
393  e_chf += pfsub->chargedHadronEnergy();
394  e_nhf += pfsub->neutralHadronEnergy();
395  e_cef += pfsub->chargedEmEnergy();
396  e_nef += pfsub->neutralEmEnergy();
397  nch += pfsub->chargedMultiplicity();
398  nconstituents += pfsub->numberOfDaughters();
399  nneutrals += pfsub->neutralMultiplicity();
400  }
401  double e = basicJet->energy();
402  if ( e > 0.000001 ) {
403  chf = e_chf / e;
404  nhf = e_nhf / e;
405  cef = e_cef / e;
406  nef = e_nef / e;
407  }
408  } // end if basic jet
409 
410 
411 
412  // Cuts for |eta| < 2.4 for FIRSTDATA, RUNIISTARTUP and WINTER16
413  if ( ignoreCut(indexCEF_) || ( cef < cut(indexCEF_, double()) || std::abs(jet.eta()) > 2.4 ) ) passCut( ret, indexCEF_);
414  if ( ignoreCut(indexCHF_) || ( chf > cut(indexCHF_, double()) || std::abs(jet.eta()) > 2.4 ) ) passCut( ret, indexCHF_);
415  if ( ignoreCut(indexNCH_) || ( nch > cut(indexNCH_, int()) || std::abs(jet.eta()) > 2.4 ) ) passCut( ret, indexNCH_);
416 
417  if(version_ == FIRSTDATA){// Cuts for all eta for FIRSTDATA
418  if ( ignoreCut(indexNConstituents_) || ( nconstituents > cut(indexNConstituents_, int()) ) ) passCut( ret, indexNConstituents_);
419  if ( ignoreCut(indexNEF_) || ( nef < cut(indexNEF_, double()) ) ) passCut( ret, indexNEF_);
420  if ( ignoreCut(indexNHF_) || ( nhf < cut(indexNHF_, double()) ) ) passCut( ret, indexNHF_);
421  }else if(version_ == RUNIISTARTUP){
422  // Cuts for |eta| <= 3.0 for RUNIISTARTUP scenario
423  if ( ignoreCut(indexNConstituents_) || ( nconstituents > cut(indexNConstituents_, int()) || std::abs(jet.eta()) > 3.0 ) ) passCut( ret, indexNConstituents_);
424  if ( ignoreCut(indexNEF_) || ( nef < cut(indexNEF_, double()) || std::abs(jet.eta()) > 3.0 ) ) passCut( ret, indexNEF_);
425  if ( ignoreCut(indexNHF_) || ( nhf < cut(indexNHF_, double()) || std::abs(jet.eta()) > 3.0 ) ) passCut( ret, indexNHF_);
426  // Cuts for |eta| > 3.0 for RUNIISTARTUP scenario
427  if ( ignoreCut(indexNEF_FW_) || ( nef < cut(indexNEF_FW_, double()) || std::abs(jet.eta()) <= 3.0 ) ) passCut( ret, indexNEF_FW_);
428  if ( ignoreCut(indexNNeutrals_FW_) || ( nneutrals > cut(indexNNeutrals_FW_, int()) || std::abs(jet.eta()) <= 3.0 ) ) passCut( ret, indexNNeutrals_FW_);
429  }
430  else if(version_ == WINTER16){
431  // Cuts for |eta| <= 2.7 for WINTER16 scenario
432  if ( ignoreCut(indexNConstituents_) || ( nconstituents > cut(indexNConstituents_, int()) || std::abs(jet.eta()) > 2.7 ) ) passCut( ret, indexNConstituents_);
433  if ( ignoreCut(indexNEF_) || ( nef < cut(indexNEF_, double()) || std::abs(jet.eta()) > 2.7 ) ) passCut( ret, indexNEF_);
434  if ( ignoreCut(indexNHF_) || ( nhf < cut(indexNHF_, double()) || std::abs(jet.eta()) > 2.7 ) ) passCut( ret, indexNHF_);
435 
436  // Cuts for 2.7 < |eta| <= 3.0 for WINTER16 scenario
437  if ( ignoreCut(indexNHF_EC_) || ( nhf < cut(indexNHF_EC_, double()) || std::abs(jet.eta()) <= 2.7 || std::abs(jet.eta()) > 3.0) ) passCut( ret, indexNHF_EC_);
438  if ( ignoreCut(indexNEF_EC_) || ( nef > cut(indexNEF_EC_, double()) || std::abs(jet.eta()) <= 2.7 || std::abs(jet.eta()) > 3.0) ) passCut( ret, indexNEF_EC_);
439  if ( ignoreCut(indexNNeutrals_EC_) || ( nneutrals > cut(indexNNeutrals_EC_, int()) || std::abs(jet.eta()) <= 2.7 || std::abs(jet.eta()) > 3.0) ) passCut( ret, indexNNeutrals_EC_);
440 
441  // Cuts for |eta| > 3.0 for WINTER16 scenario
442  if ( ignoreCut(indexNEF_FW_) || ( nef < cut(indexNEF_FW_, double()) || std::abs(jet.eta()) <= 3.0 ) ) passCut( ret, indexNEF_FW_);
443  if ( ignoreCut(indexNNeutrals_FW_) || ( nneutrals > cut(indexNNeutrals_FW_, int()) || std::abs(jet.eta()) <= 3.0 ) ) passCut( ret, indexNNeutrals_FW_);
444  }
445 
446 
447  //std::cout << "<PFJetIDSelectionFunctor::firstDataCuts>:" << std::endl;
448  //std::cout << " jet: Pt = " << jet.pt() << ", eta = " << jet.eta() << ", phi = " << jet.phi() << std::endl;
449  //ret.print(std::cout);
450 
451  setIgnored( ret );
452  return (bool)ret;
453  }
454 
455  private: // member variables
456 
459 
466 
469 
473 
474 
475 };
476 
477 #endif
T getParameter(std::string const &) const
float photonEnergy() const
photonEnergy
Definition: PFJet.h:106
float neutralHadronEnergyFraction() const
neutralHadronEnergyFraction (relative to uncorrected jet energy)
Definition: Jet.h:370
bool jecSetsAvailable() const
Definition: Jet.h:129
float muonEnergy() const
muonEnergy
Definition: PFJet.h:114
float chargedEmEnergy() const
chargedEmEnergy
Definition: PFJet.h:142
virtual double eta() const final
momentum pseudorapidity
float chargedHadronEnergyFraction() const
chargedHadronEnergyFraction (relative to uncorrected jet energy)
Definition: Jet.h:368
PFJetIDSelectionFunctor(edm::ParameterSet const &params, edm::ConsumesCollector &iC)
Base class for all types of Jets.
Definition: Jet.h:20
float chargedEmEnergyFraction() const
chargedEmEnergyFraction (relative to uncorrected jet energy)
Definition: Jet.h:372
bool exists(std::string const &parameterName) const
checks if a parameter exists
pat::strbitset::index_type index_type
Definition: Selector.h:29
void setIgnored(pat::strbitset &ret)
set ignored bits
Definition: Selector.h:224
int neutralMultiplicity() const
neutralMultiplicity
Definition: Jet.h:418
pat::strbitset retInternal_
internal ret if users don&#39;t care about return bits
Definition: Selector.h:287
Jets made from CaloTowers.
Definition: BasicJet.h:20
int chargedMultiplicity() const
chargedMultiplicity
Definition: PFJet.h:155
Jets made from PFObjects.
Definition: PFJet.h:21
PFJetIDSelectionFunctor(Version_t version, Quality_t quality)
float neutralEmEnergy() const
neutralEmEnergy
Definition: PFJet.h:150
bool operator()(const reco::PFJet &jet, pat::strbitset &ret)
float electronEnergy() const
electronEnergy
Definition: PFJet.h:110
PFJetIDSelectionFunctor(edm::ParameterSet const &params)
pat::strbitset bits_
the bitset indexed by strings
Definition: Selector.h:286
virtual size_t numberOfDaughters() const
number of daughters
void passCut(pat::strbitset &ret, std::string const &s)
Passing cuts.
Definition: Selector.h:176
float HFEMEnergy() const
HFEMEnergy.
Definition: PFJet.h:122
bool ignoreCut(std::string const &s) const
ignore the cut at index "s"
Definition: Selector.h:159
const_iterator end() const
last daughter const_iterator
Definition: Candidate.h:146
bool operator()(const pat::Jet &jet, pat::strbitset &ret)
virtual double energy() const final
energy
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
virtual void push_back(std::string const &s)
This is the registration of an individual cut string.
Definition: Selector.h:46
bool operator()(const reco::PFJet &jet)
Functor that operates on <T>
Definition: Selector.h:24
PF Jet selector for pat::Jets.
bool isPFJet() const
check to see if the jet is a reco::PFJet
Definition: Jet.h:253
int neutralMultiplicity() const
neutralMultiplicity
Definition: PFJet.h:157
virtual size_t numberOfDaughters() const
strbitset & set(bool val=true)
set method of all bits
Definition: strbitset.h:144
std::string currentJECLevel() const
return the name of the current step of jet energy corrections
Definition: Jet.h:139
Analysis-level calorimeter jet class.
Definition: Jet.h:78
pat::strbitset getBitTemplate() const
Get an empty bitset with the proper names.
Definition: Selector.h:212
bool isBasicJet() const
check to see if the jet is no more than a reco::BasicJet
Definition: Jet.h:255
const_iterator begin() const
first daughter const_iterator
Definition: Candidate.h:144
float neutralHadronEnergy() const
neutralHadronEnergy
Definition: PFJet.h:102
void setIgnoredCuts(std::vector< std::string > const &bitsToIgnore)
set the bits to ignore from a vector
Definition: Selector.h:167
float neutralEmEnergyFraction() const
neutralEmEnergyFraction (relative to uncorrected jet energy)
Definition: Jet.h:374
Jet correctedJet(const std::string &level, const std::string &flavor="none", const std::string &set="") const
int chargedMultiplicity() const
chargedMultiplicity
Definition: Jet.h:659
bool firstDataCuts(reco::Jet const &jet, pat::strbitset &ret, Version_t version_)
float chargedHadronEnergy() const
chargedHadronEnergy
Definition: PFJet.h:98
int cut(index_type const &i, int val) const
Access the int cut values at index "s".
Definition: Selector.h:194