SimpleCutBasedElectronIDSelectionFunctor.h

Go to the documentation of this file.

#ifndef PhysicsTools_PatUtils_interface_SimpleCutBasedElectronIDSelectionFunctor_h
#define PhysicsTools_PatUtils_interface_SimpleCutBasedElectronIDSelectionFunctor_h

#include "DataFormats/PatCandidates/interface/Electron.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"

#include "PhysicsTools/SelectorUtils/interface/Selector.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"

#include "FWCore/Utilities/interface/InputTag.h"

/*
___________________________________________________________________________________

Description:
^^^^^^^^^^^^
    This is a class that implements the Simple Cut Based Electron 
    Identification cuts.  A more detailed description of the cuts
    and the tuning method can be found on this twiki:
    
    https://twiki.cern.ch/twiki/bin/view/CMS/SimpleCutBasedEleID

    For more information on how to calculate the magnetic field
    look here:

    https://twiki.cern.ch/twiki/bin/viewauth/CMS/ConversionBackgroundRejection
___________________________________________________________________________________

How to use:
^^^^^^^^^^^
    From CMSSW39 onwards you can simply define an instance of this class:

      SimpleCutBasedElectronIDSelectionFunctor patSele95
      (SimpleCutBasedElectronIDSelectionFunctor::relIso95);

    and get the decision with the following method:
      pat::Electron *myElec = .....;
      bool pass = patSele90(*myElec);

    The various options are listed in the enumeration Version_t. There
    is also the option to enter as a constructor argument a PSet
    with your favorite cuts.
___________________________________________________________________________________

    Contacts: Nikolaos Rompotis and Chris Seez
    Nikolaos dot Rompotis at Cern dot ch
    Chris    dot Seez     at Cern dot ch

    Author:    Nikolaos Rompotis
               many thanks to Sal Rappoccio
    Imperial College London
    7 June 2010, first commit for CMSSW_3_6_1_patchX
    11July 2010, implementing the ICHEP Egamma recommendation for 
                 removing the Delta Eta cut in the endcaps
    30Sept 2010, simplification of conversion rejection in CMSSW39X
___________________________________________________________________________________

*/


class SimpleCutBasedElectronIDSelectionFunctor : public Selector<pat::Electron>  {

 public: // interface  
  
  enum Version_t { relIso95=0, cIso95,  relIso90, cIso90, relIso85, cIso85, 
           relIso80, cIso80,  relIso70, cIso70, relIso60, cIso60, NONE };
  
  SimpleCutBasedElectronIDSelectionFunctor() {}
  
  // initialize it by inserting directly the cut values in a parameter set
  SimpleCutBasedElectronIDSelectionFunctor(edm::ParameterSet const & parameters)
    {
      // get the cuts from the PS
      initialize( parameters.getParameter<Double_t>("trackIso_EB"), 
          parameters.getParameter<Double_t>("ecalIso_EB"), 
          parameters.getParameter<Double_t>("hcalIso_EB"), 
          parameters.getParameter<Double_t>("sihih_EB"), 
          parameters.getParameter<Double_t>("dphi_EB"), 
          parameters.getParameter<Double_t>("deta_EB"), 
          parameters.getParameter<Double_t>("hoe_EB"), 
          parameters.getParameter<Double_t>("cIso_EB"), 
          parameters.getParameter<Double_t>("trackIso_EE"), 
          parameters.getParameter<Double_t>("ecalIso_EE"), 
          parameters.getParameter<Double_t>("hcalIso_EE"), 
          parameters.getParameter<Double_t>("sihih_EE"), 
          parameters.getParameter<Double_t>("dphi_EE"), 
          parameters.getParameter<Double_t>("deta_EE"), 
          parameters.getParameter<Double_t>("hoe_EE"), 
          parameters.getParameter<Double_t>("cIso_EE"), 
          parameters.getParameter<Int_t>("conversionRejection"), 
          parameters.getParameter<Int_t>("maxNumberOfExpectedMissingHits"));
      retInternal_ = getBitTemplate();
    }
  // initialize it by using only the version name
  SimpleCutBasedElectronIDSelectionFunctor(Version_t  version)
    {
      if (version == NONE) {
    std::cout << "SimpleCutBasedElectronIDSelectionFunctor: If you want to use version NONE "
          << "then you have also to provide the selection cuts by yourself " << std::endl;
    std::cout << "SimpleCutBasedElectronIDSelectionFunctor: ID Version is changed to 80cIso "
          << std::endl;
    version = cIso80;
      }
      initialize(version);
      retInternal_ = getBitTemplate();
    }

  void initialize( Version_t version ) 
  {
    version_ = version;
    // push back the variables
    push_back("trackIso_EB");
    push_back("ecalIso_EB" );
    push_back("hcalIso_EB" );
    push_back("sihih_EB"   );
    push_back("dphi_EB"    );
    push_back("deta_EB"    );
    push_back("hoe_EB"     );
    push_back("cIso_EB"    );
    
    push_back("trackIso_EE");
    push_back("ecalIso_EE" );
    push_back("hcalIso_EE" );
    push_back("sihih_EE"   );
    push_back("dphi_EE"    );
    push_back("deta_EE"    );
    push_back("hoe_EE"     );
    push_back("cIso_EE"    );
    
    push_back("conversionRejection"            );
    push_back("maxNumberOfExpectedMissingHits" );
    
    
    
    
    if (version_ == relIso95) {
      set("trackIso_EB", 1.5e-01);
      set("ecalIso_EB",  2.0e+00);
      set("hcalIso_EB",  1.2e-01);
      set("sihih_EB",    1.0e-02);
      set("dphi_EB",     8.0e-01);
      set("deta_EB",     7.0e-03);
      set("hoe_EB",      1.5e-01);
      set("cIso_EB",     10000. );
      
      set("trackIso_EE", 8.0e-02);
      set("ecalIso_EE",  6.0e-02);
      set("hcalIso_EE",  5.0e-02);
      set("sihih_EE",    3.0e-02);
      set("dphi_EE",     7.0e-01);
      set("deta_EE",     1.0e-02);
      set("hoe_EE",      7.0e-02);
      set("cIso_EE",     10000. );
      
      set("conversionRejection",            0);
      set("maxNumberOfExpectedMissingHits", 1);

    }
    else if (version_ == cIso95) {
      set("trackIso_EB", 100000.);
      set("ecalIso_EB",  100000.);
      set("hcalIso_EB",  100000.);
      set("sihih_EB",    1.0e-02);
      set("dphi_EB",     8.0e-01);
      set("deta_EB",     7.0e-03);
      set("hoe_EB",      1.5e-01);
      set("cIso_EB",     1.5e-01);
                                              
      set("trackIso_EE", 100000.);
      set("ecalIso_EE",  100000.);
      set("hcalIso_EE",  100000.);
      set("sihih_EE",    3.0e-02);
      set("dphi_EE",     7.0e-01);
      set("deta_EE",     1.0e-02);
      set("hoe_EE",      7.0e-02);
      set("cIso_EE",     1.0e-01);
      
      set("conversionRejection",            0);
      set("maxNumberOfExpectedMissingHits", 1);
      
    }
    else if (version_ == relIso90) {
      set("trackIso_EB", 1.2e-01);
      set("ecalIso_EB",  9.0e-02);
      set("hcalIso_EB",  1.0e-01);
      set("sihih_EB",    1.0e-02);
      set("dphi_EB",     8.0e-01);
      set("deta_EB",     7.0e-03);
      set("hoe_EB",      1.2e-01);
      set("cIso_EB",     10000. );

      set("trackIso_EE", 5.0e-02);
      set("ecalIso_EE",  6.0e-02);
      set("hcalIso_EE",  3.0e-02);
      set("sihih_EE",    3.0e-02);
      set("dphi_EE",     7.0e-01);
      set("deta_EE",     9.0e-03);
      set("hoe_EE",      5.0e-02);
      set("cIso_EE",     10000. );
      
      set("conversionRejection",            1);
      set("maxNumberOfExpectedMissingHits", 1);
    }
    else if (version_ == cIso90) {
      set("trackIso_EB", 100000.);
      set("ecalIso_EB",  100000.);
      set("hcalIso_EB",  100000.);
      set("sihih_EB",    1.0e-02);
      set("dphi_EB",     8.0e-01);
      set("deta_EB",     7.0e-03);
      set("hoe_EB",      1.2e-01);
      set("cIso_EB",     1.0e-01);

      set("trackIso_EE", 100000.);
      set("ecalIso_EE",  100000.);
      set("hcalIso_EE",  100000.);
      set("sihih_EE",    3.0e-02);
      set("dphi_EE",     7.0e-01);
      set("deta_EE",     9.0e-03);
      set("hoe_EE",      5.0e-02);
      set("cIso_EE",     7.0e-02);
      
      set("conversionRejection",            1);
      set("maxNumberOfExpectedMissingHits", 1);
    }
    else if (version_ == relIso85) {
      set("trackIso_EB", 9.0e-02);
      set("ecalIso_EB",  8.0e-02);
      set("hcalIso_EB",  1.0e-01);
      set("sihih_EB",    1.0e-02);
      set("dphi_EB",     6.0e-02);
      set("deta_EB",     6.0e-03);
      set("hoe_EB",      4.0e-02);
      set("cIso_EB",     10000. );

      set("trackIso_EE", 5.0e-02);
      set("ecalIso_EE",  5.0e-02);
      set("hcalIso_EE",  2.5e-02);
      set("sihih_EE",    3.0e-02);
      set("dphi_EE",     4.0e-02);
      set("deta_EE",     7.0e-03);
      set("hoe_EE",      2.5e-02);
      set("cIso_EE",     10000. );
      
      set("conversionRejection",            1);
      set("maxNumberOfExpectedMissingHits", 1);
    }
    else if (version_ == cIso85) {
      set("trackIso_EB", 100000.);
      set("ecalIso_EB",  100000.);
      set("hcalIso_EB",  100000.);
      set("sihih_EB",    1.0e-02);
      set("dphi_EB",     6.0e-02);
      set("deta_EB",     6.0e-03);
      set("hoe_EB",      4.0e-02);
      set("cIso_EB",     9.0e-02);

      set("trackIso_EE", 100000.);
      set("ecalIso_EE",  100000.);
      set("hcalIso_EE",  100000.);
      set("sihih_EE",    3.0e-02);
      set("dphi_EE",     4.0e-02);
      set("deta_EE",     7.0e-03);
      set("hoe_EE",      2.5e-02);
      set("cIso_EE",     6.0e-02);
      
      set("conversionRejection",            1);
      set("maxNumberOfExpectedMissingHits", 1);
    }
    else if (version_ == relIso80) {
      set("trackIso_EB", 9.0e-02);
      set("ecalIso_EB",  7.0e-02);
      set("hcalIso_EB",  1.0e-01);
      set("sihih_EB",    1.0e-02);
      set("dphi_EB",     6.0e-02);
      set("deta_EB",     4.0e-03);
      set("hoe_EB",      4.0e-02);
      set("cIso_EB",     100000.);

      set("trackIso_EE", 4.0e-02);
      set("ecalIso_EE",  5.0e-02);
      set("hcalIso_EE",  2.5e-02);
      set("sihih_EE",    3.0e-02);
      set("dphi_EE",     3.0e-02);
      set("deta_EE",     7.0e-03);
      set("hoe_EE",      2.5e-02);
      set("cIso_EE",     100000.);
      
      set("conversionRejection",            1);
      set("maxNumberOfExpectedMissingHits", 0);
    }
    else if (version_ == cIso80) {
      set("trackIso_EB", 100000.);
      set("ecalIso_EB",  100000.);
      set("hcalIso_EB",  100000.);
      set("sihih_EB",    1.0e-02);
      set("dphi_EB",     6.0e-02);
      set("deta_EB",     4.0e-03);
      set("hoe_EB",      4.0e-02);
      set("cIso_EB",     7.0e-02);

      set("trackIso_EE", 100000.);
      set("ecalIso_EE",  100000.);
      set("hcalIso_EE",  100000.);
      set("sihih_EE",    3.0e-02);
      set("dphi_EE",     3.0e-02);
      set("deta_EE",     7.0e-03);
      set("hoe_EE",      2.5e-02);
      set("cIso_EE",     6.0e-02);
      
      set("conversionRejection",            1);
      set("maxNumberOfExpectedMissingHits", 0);
    }
    else if (version_ == relIso70) {
      set("trackIso_EB", 5.0e-02);
      set("ecalIso_EB",  6.0e-02);
      set("hcalIso_EB",  3.0e-02);
      set("sihih_EB",    1.0e-02);
      set("dphi_EB",     3.0e-02);
      set("deta_EB",     4.0e-03);
      set("hoe_EB",      2.5e-02);
      set("cIso_EB",     100000.);

      set("trackIso_EE", 2.5e-02);
      set("ecalIso_EE",  2.5e-02);
      set("hcalIso_EE",  2.0e-02);
      set("sihih_EE",    3.0e-02);
      set("dphi_EE",     2.0e-02);
      set("deta_EE",     5.0e-03);
      set("hoe_EE",      2.5e-02);
      set("cIso_EE",     100000.);
      
      set("conversionRejection",            1);
      set("maxNumberOfExpectedMissingHits", 0);
    }
    else if (version_ == cIso70) {
      set("trackIso_EB", 100000.);
      set("ecalIso_EB",  100000.);
      set("hcalIso_EB",  100000.);
      set("sihih_EB",    1.0e-02);
      set("dphi_EB",     3.0e-02);
      set("deta_EB",     4.0e-03);
      set("hoe_EB",      2.5e-02);
      set("cIso_EB",     4.0e-02);

      set("trackIso_EE", 100000.);
      set("ecalIso_EE",  100000.);
      set("hcalIso_EE",  100000.);
      set("sihih_EE",    3.0e-02);
      set("dphi_EE",     2.0e-02);
      set("deta_EE",     5.0e-03);
      set("hoe_EE",      2.5e-02);
      set("cIso_EE",     3.0e-02);
      
      set("conversionRejection",            1);
      set("maxNumberOfExpectedMissingHits", 0);
    }
    else if (version_ == relIso60) {
      set("trackIso_EB", 4.0e-02);
      set("ecalIso_EB",  4.0e-02);
      set("hcalIso_EB",  3.0e-02);
      set("sihih_EB",    1.0e-02);
      set("dphi_EB",     2.5e-02);
      set("deta_EB",     4.0e-03);
      set("hoe_EB",      2.5e-02);
      set("cIso_EB",     100000.);

      set("trackIso_EE", 2.5e-02);
      set("ecalIso_EE",  2.0e-02);
      set("hcalIso_EE",  2.0e-02);
      set("sihih_EE",    3.0e-02);
      set("dphi_EE",     2.0e-02);
      set("deta_EE",     5.0e-03);
      set("hoe_EE",      2.5e-02);
      set("cIso_EE",     100000.);
      
      set("conversionRejection",            1);
      set("maxNumberOfExpectedMissingHits", 0);
    }
    else if (version_ == cIso60) {
      set("trackIso_EB", 100000.);
      set("ecalIso_EB",  100000.);
      set("hcalIso_EB",  100000.);
      set("sihih_EB",    1.0e-02);
      set("dphi_EB",     2.5e-02);
      set("deta_EB",     4.0e-03);
      set("hoe_EB",      2.5e-02);
      set("cIso_EB",     3.0e-02);

      set("trackIso_EE", 100000.);
      set("ecalIso_EE",  100000.);
      set("hcalIso_EE",  100000.);
      set("sihih_EE",    3.0e-02);
      set("dphi_EE",     2.0e-02);
      set("deta_EE",     5.0e-03);
      set("hoe_EE",      2.5e-02);
      set("cIso_EE",     2.0e-02);
      
      set("conversionRejection",            1);
      set("maxNumberOfExpectedMissingHits", 0);
    }
  }

  void initialize(Double_t trackIso_EB, Double_t ecalIso_EB, Double_t hcalIso_EB,
          Double_t sihih_EB, Double_t  dphi_EB, Double_t deta_EB, Double_t hoe_EB,
          Double_t cIso_EB,
          Double_t trackIso_EE, Double_t ecalIso_EE, Double_t hcalIso_EE,
          Double_t sihih_EE, Double_t  dphi_EE, Double_t deta_EE, Double_t hoe_EE,
          Double_t cIso_EE, Int_t conversionRejection, 
          Int_t maxNumberOfExpectedMissingHits)
  {
    version_ = NONE;
    push_back("trackIso_EB");
    push_back("ecalIso_EB" );
    push_back("hcalIso_EB" );
    push_back("sihih_EB"   );
    push_back("dphi_EB"    );
    push_back("deta_EB"    );
    push_back("hoe_EB"     );
    push_back("cIso_EB"    );
    
    push_back("trackIso_EE");
    push_back("ecalIso_EE" );
    push_back("hcalIso_EE" );
    push_back("sihih_EE"   );
    push_back("dphi_EE"    );
    push_back("deta_EE"    );
    push_back("hoe_EE"     );
    push_back("cIso_EE"    );
    
    push_back("conversionRejection"            );
    push_back("maxNumberOfExpectedMissingHits" );
    
   
    set("trackIso_EB", trackIso_EB);
    set("ecalIso_EB",  ecalIso_EB);
    set("hcalIso_EB",  hcalIso_EB);
    set("sihih_EB",    sihih_EB);
    set("dphi_EB",     dphi_EB);
    set("deta_EB",     deta_EB);
    set("hoe_EB",      hoe_EB);
    set("cIso_EB",     cIso_EB);
    
    set("trackIso_EE", trackIso_EE);
    set("ecalIso_EE",  ecalIso_EE);
    set("hcalIso_EE",  hcalIso_EE);
    set("sihih_EE",    sihih_EE);
    set("dphi_EE",     dphi_EE);
    set("deta_EE",     deta_EE);
    set("hoe_EE",      hoe_EE);
    set("cIso_EE",     cIso_EE);
    
    set("conversionRejection",            conversionRejection);
    set("maxNumberOfExpectedMissingHits", maxNumberOfExpectedMissingHits);
    
  }

  bool operator()( const pat::Electron & electron, pat::strbitset & ret ) 
  {
    // for the time being only Spring10 variable definition
    return spring10Variables(electron, ret);
  }
  using Selector<pat::Electron>::operator();
  // function with the Spring10 variable definitions
  bool spring10Variables( const pat::Electron & electron, pat::strbitset & ret) 
  {
    ret.set(false);
    //
    Double_t eleET = electron.p4().Pt();
    Double_t trackIso = electron.dr03TkSumPt()/eleET;
    Double_t ecalIso = electron.dr03EcalRecHitSumEt()/eleET;
    Double_t hcalIso = electron.dr03HcalTowerSumEt()/eleET;
    Double_t sihih   = electron.sigmaIetaIeta();
    Double_t Dphi    = electron.deltaPhiSuperClusterTrackAtVtx();
    Double_t Deta    = electron.deltaEtaSuperClusterTrackAtVtx();
    Double_t HoE     = electron.hadronicOverEm();
    Double_t cIso    = 0;
    if (electron.isEB()) { cIso = 
    ( electron.dr03TkSumPt() + std::max(0.,electron.dr03EcalRecHitSumEt() -1.) 
      + electron.dr03HcalTowerSumEt() ) / eleET;
    }
    else {
      cIso = ( electron.dr03TkSumPt()+electron.dr03EcalRecHitSumEt()+
           electron.dr03HcalTowerSumEt()  ) / eleET;
    }
    Int_t innerHits = electron.gsfTrack()->trackerExpectedHitsInner().numberOfHits();
    // in 39 conversion rejection variables are accessible from Gsf electron
    Double_t dist = electron.convDist(); // default value is -9999 if conversion partner not found
    Double_t dcot = electron.convDcot(); // default value is -9999 if conversion partner not found
    Bool_t isConv = fabs(dist) < 0.02 && fabs(dcot) < 0.02;
    // now apply the cuts
    if (electron.isEB()) { // BARREL case
      // check the EB cuts
      if ( trackIso   <  cut("trackIso_EB", double()) || ignoreCut("trackIso_EB")) passCut(ret, "trackIso_EB");
      if ( ecalIso    <  cut("ecalIso_EB",  double()) || ignoreCut("ecalIso_EB") ) passCut(ret, "ecalIso_EB");
      if ( hcalIso    <  cut("hcalIso_EB",  double()) || ignoreCut("hcalIso_EB") ) passCut(ret, "hcalIso_EB");
      if ( sihih      <  cut("sihih_EB",    double()) || ignoreCut("sihih_EB")   ) passCut(ret, "sihih_EB");
      if ( fabs(Dphi) <  cut("dphi_EB",     double()) || ignoreCut("dphi_EB")    ) passCut(ret, "dphi_EB");
      if ( fabs(Deta) <  cut("deta_EB",     double()) || ignoreCut("deta_EB")    ) passCut(ret, "deta_EB");
      if ( HoE        <  cut("hoe_EB",      double()) || ignoreCut("hoe_EB")     ) passCut(ret, "hoe_EB");
      if ( cIso       <  cut("cIso_EB",     double()) || ignoreCut("cIso_EB")    ) passCut(ret, "cIso_EB");
      // pass all the EE cuts
      passCut(ret, "trackIso_EE");  
      passCut(ret, "ecalIso_EE");   
      passCut(ret, "hcalIso_EE");   
      passCut(ret, "sihih_EE"); 
      passCut(ret, "dphi_EE");  
      passCut(ret, "deta_EE");  
      passCut(ret, "hoe_EE");   
      passCut(ret, "cIso_EE");     
    } else {  // ENDCAPS case
      // check the EE cuts
      if ( trackIso   <  cut("trackIso_EE", double()) || ignoreCut("trackIso_EE")) passCut(ret, "trackIso_EE");
      if ( ecalIso    <  cut("ecalIso_EE",  double()) || ignoreCut("ecalIso_EE") ) passCut(ret, "ecalIso_EE");
      if ( hcalIso    <  cut("hcalIso_EE",  double()) || ignoreCut("hcalIso_EE") ) passCut(ret, "hcalIso_EE");
      if ( sihih      <  cut("sihih_EE",    double()) || ignoreCut("sihih_EE")   ) passCut(ret, "sihih_EE");
      if ( fabs(Dphi) <  cut("dphi_EE",     double()) || ignoreCut("dphi_EE")    ) passCut(ret, "dphi_EE");
      if ( fabs(Deta) <  cut("deta_EE",     double()) || ignoreCut("deta_EE")    ) passCut(ret, "deta_EE");
      if ( HoE        <  cut("hoe_EE",      double()) || ignoreCut("hoe_EE")     ) passCut(ret, "hoe_EE");
      if ( cIso       <  cut("cIso_EE",     double()) || ignoreCut("cIso_EE")    ) passCut(ret, "cIso_EE");     
      // pass all the EB cuts
      passCut(ret, "trackIso_EB");  
      passCut(ret, "ecalIso_EB");   
      passCut(ret, "hcalIso_EB");   
      passCut(ret, "sihih_EB"); 
      passCut(ret, "dphi_EB");  
      passCut(ret, "deta_EB");  
      passCut(ret, "hoe_EB");   
      passCut(ret, "cIso_EB");     
    }

    // conversion rejection common for EB and EE
    if ( innerHits  <=  cut("maxNumberOfExpectedMissingHits", int())) 
      passCut(ret, "maxNumberOfExpectedMissingHits");    
    if ( 0==cut("conversionRejection", int()) || isConv==false)
      passCut(ret, "conversionRejection");
    setIgnored(ret);   
    return (bool)ret;
  }



 private: // member variables
  // version of the cuts  
  Version_t version_;
};


#endif