[9_4_X] Egamma Fall17 V2 IDs

Configuration/Eras/python/Modifier_run2_miniAOD_devel_cff.py

@@ -0,0 +1,3 @@
+import FWCore.ParameterSet.Config as cms
+
+run2_miniAOD_devel = cms.Modifier()

Configuration/StandardSequences/python/Eras.py

@@ -46,6 +46,7 @@ def __init__(self):
                            'trackingLowPU', 'trackingPhase1', 'trackingPhase1QuadProp', 'ctpps_2016', 'trackingPhase2PU140',
                            'tracker_apv_vfp30_2016', 'run2_miniAOD_80XLegacy', 'run2_miniAOD_94XFall17', 'run2_nanoAOD_92X',
                            'run2_nanoAOD_94XMiniAODv1', 'run2_nanoAOD_94XMiniAODv2', 'run2_nanoAOD_94X2016',
+                           'run2_miniAOD_devel',
                            'hcalHardcodeConditions', 'hcalSkipPacker']
         internalUseModChains = ['run2_2017_noTrackingModifier']
 

PhysicsTools/PatAlgos/python/slimming/miniAOD_tools.py

@@ -268,7 +268,8 @@ def miniAOD_customizeCommon(process):
     process.slimmedPhotons.modifierConfig.modifications   = egamma_modifications
 
     #VID Electron IDs
-    electron_ids = ['RecoEgamma.ElectronIdentification.Identification.cutBasedElectronID_PHYS14_PU20bx25_V2_cff',
+    electron_id_config = cms.PSet(electron_ids = cms.vstring([
+                    'RecoEgamma.ElectronIdentification.Identification.cutBasedElectronID_PHYS14_PU20bx25_V2_cff',
                     'RecoEgamma.ElectronIdentification.Identification.cutBasedElectronID_Spring15_25ns_V1_cff',
                     'RecoEgamma.ElectronIdentification.Identification.cutBasedElectronID_Spring15_50ns_V2_cff',
                     'RecoEgamma.ElectronIdentification.Identification.heepElectronID_HEEPV60_cff',
@@ -282,15 +283,23 @@ def miniAOD_customizeCommon(process):
                     'RecoEgamma.ElectronIdentification.Identification.cutBasedElectronID_Summer16_80X_V1_cff',
                     'RecoEgamma.ElectronIdentification.Identification.mvaElectronID_Spring16_GeneralPurpose_V1_cff',
                     'RecoEgamma.ElectronIdentification.Identification.mvaElectronID_Spring16_HZZ_V1_cff',
-                    ]
+                    ]))
+
+    from Configuration.Eras.Modifier_run2_miniAOD_devel_cff import run2_miniAOD_devel
+    run2_miniAOD_devel.toModify(electron_id_config, electron_ids = electron_id_config.electron_ids + [\
+                    'RecoEgamma.ElectronIdentification.Identification.cutBasedElectronID_Fall17_94X_V2_cff',
+                    'RecoEgamma.ElectronIdentification.Identification.mvaElectronID_Fall17_noIso_V2_cff', 
+                    'RecoEgamma.ElectronIdentification.Identification.mvaElectronID_Fall17_iso_V2_cff',
+                    ])
+
     switchOnVIDElectronIdProducer(process,DataFormat.MiniAOD, task)
     process.egmGsfElectronIDs.physicsObjectSrc = \
         cms.InputTag("reducedEgamma","reducedGedGsfElectrons")
     process.electronMVAValueMapProducer.src = \
         cms.InputTag('reducedEgamma','reducedGedGsfElectrons')
     process.electronRegressionValueMapProducer.src = \
         cms.InputTag('reducedEgamma','reducedGedGsfElectrons')
-    for idmod in electron_ids:
+    for idmod in electron_id_config.electron_ids.value():
         setupAllVIDIdsInModule(process,idmod,setupVIDElectronSelection,None,False,task)
 
     #heepIDVarValueMaps only exists if HEEP V6.1 or HEEP 7.0 ID has already been loaded
@@ -306,15 +315,23 @@ def miniAOD_customizeCommon(process):
                 break
 
     #VID Photon IDs
-    photon_ids = ['RecoEgamma.PhotonIdentification.Identification.cutBasedPhotonID_Spring15_25ns_V1_cff',
+    photon_id_config = cms.PSet(photon_ids = cms.vstring([
+                  'RecoEgamma.PhotonIdentification.Identification.cutBasedPhotonID_Spring15_25ns_V1_cff',
                   'RecoEgamma.PhotonIdentification.Identification.cutBasedPhotonID_Spring15_50ns_V1_cff',
                   'RecoEgamma.PhotonIdentification.Identification.mvaPhotonID_Spring15_25ns_nonTrig_V2p1_cff',
                   'RecoEgamma.PhotonIdentification.Identification.mvaPhotonID_Spring15_50ns_nonTrig_V2p1_cff',
                   'RecoEgamma.PhotonIdentification.Identification.cutBasedPhotonID_Fall17_94X_V1_TrueVtx_cff',
-                  'RecoEgamma.PhotonIdentification.Identification.mvaPhotonID_Fall17_94X_V1_cff',
+                  'RecoEgamma.PhotonIdentification.Identification.mvaPhotonID_Fall17_94X_V1_cff', 
                   'RecoEgamma.PhotonIdentification.Identification.mvaPhotonID_Fall17_94X_V1p1_cff', 
                   'RecoEgamma.PhotonIdentification.Identification.cutBasedPhotonID_Spring16_V2p2_cff',
-                  'RecoEgamma.PhotonIdentification.Identification.mvaPhotonID_Spring16_nonTrig_V1_cff']
+                  'RecoEgamma.PhotonIdentification.Identification.mvaPhotonID_Spring16_nonTrig_V1_cff',
+                  ]))
+
+    from Configuration.Eras.Modifier_run2_miniAOD_devel_cff import run2_miniAOD_devel
+    run2_miniAOD_devel.toModify(photon_id_config, photon_ids = photon_id_config.photon_ids + \
+            ["RecoEgamma.PhotonIdentification.Identification.mvaPhotonID_Fall17_94X_V2_cff",
+             "RecoEgamma.PhotonIdentification.Identification.cutBasedPhotonID_Fall17_94X_V2_cff"])
+
     switchOnVIDPhotonIdProducer(process,DataFormat.AOD, task) 
     process.egmPhotonIsolation.srcToIsolate = \
         cms.InputTag("reducedEgamma","reducedGedPhotons")  
@@ -331,7 +348,7 @@ def miniAOD_customizeCommon(process):
         cms.InputTag("reducedEgamma","reducedPhotonPfCandMap")
     process.photonMVAValueMapProducer.src = \
         cms.InputTag('reducedEgamma','reducedGedPhotons')
-    for idmod in photon_ids:
+    for idmod in photon_id_config.photon_ids.value():
         setupAllVIDIdsInModule(process,idmod,setupVIDPhotonSelection,None,False,task)
 
     #add the cut base IDs bitmaps of which cuts passed

RecoEgamma/EgammaTools/interface/MVAValueMapProducer.h

@@ -56,6 +56,7 @@ class MVAValueMapProducer : public edm::stream::EDProducer< edm::GlobalCache<ega
 
   // Value map names
   std::vector <std::string> mvaValueMapNames_;
+  std::vector <std::string> mvaRawValueMapNames_;
   std::vector <std::string> mvaCategoriesMapNames_;
 
 };
@@ -84,12 +85,15 @@ MVAValueMapProducer<ParticleType>::MVAValueMapProducer(const edm::ParameterSet&
     const std::string full_name = ( currentEstimator->getName() + 
                                     currentEstimator->getTag()    );
     std::string thisValueMapName = full_name + "Values";
+    std::string thisRawValueMapName = full_name + "RawValues";
     std::string thisCategoriesMapName = full_name + "Categories";    
     mvaValueMapNames_.push_back( thisValueMapName );
+    mvaRawValueMapNames_.push_back( thisRawValueMapName );
     mvaCategoriesMapNames_.push_back( thisCategoriesMapName );
 
     // Declare the maps to the framework
     produces<edm::ValueMap<float> >(thisValueMapName);  
+    produces<edm::ValueMap<float> >(thisRawValueMapName);
     produces<edm::ValueMap<int> >(thisCategoriesMapName);
   }
 
@@ -130,16 +134,20 @@ void MVAValueMapProducer<ParticleType>::produce(edm::Event& iEvent, const edm::E
     const auto& thisEstimator = mva_itr->second;
 
     std::vector<float> mvaValues;
+    std::vector<float> mvaRawValues;
     std::vector<int> mvaCategories;
 
     // Loop over particles
     for (size_t i = 0; i < src->size(); ++i){
       auto iCand = src->ptrAt(i);      
-      mvaValues.push_back( thisEstimator->mvaValue( iCand, iEvent ) );
+      const float response = thisEstimator->mvaValue( iCand, iEvent );
+      mvaRawValues.push_back( response ); // The MVA score
+      mvaValues.push_back( 2.0/(1.0+exp(-2.0*response))-1 ); // MVA output between -1 and 1
       mvaCategories.push_back( thisEstimator->findCategory( iCand ) );
     } // end loop over particles
 
     writeValueMap(iEvent, src, mvaValues, mvaValueMapNames_[iEstimator] );  
+    writeValueMap(iEvent, src, mvaRawValues, mvaRawValueMapNames_[iEstimator] );
     writeValueMap(iEvent, src, mvaCategories, mvaCategoriesMapNames_[iEstimator] );
   } // end loop over estimators
 

RecoEgamma/ElectronIdentification/data/Fall17/effAreaElectrons_cone03_pfNeuHadronsAndPhotons_94X.txt

@@ -0,0 +1,18 @@
+# This file contains Effective Area constants for 
+# computing pile-up corrections for the neutral hadron and photon
+# isolation for an electron object.
+#   Documentation:
+# 
+#   https://indico.cern.ch/event/697576/contributions/2940576/attachments/1620927/2578913/eleIdTuning.pdf
+#   (slides 3 to 5)
+# 
+#  The effective areas are based on 90% efficient contours
+#
+# |eta| min   |eta| max   effective area
+0.000        1.000       0.1440
+1.000        1.479       0.1562
+1.479        2.000       0.1032
+2.000        2.200       0.0859
+2.200        2.300       0.1116
+2.300        2.400       0.1321
+2.400        2.500       0.1654

RecoEgamma/ElectronIdentification/interface/ElectronMVAEstimatorRun2Fall17.h

@@ -128,6 +128,8 @@ class ElectronMVAEstimatorRun2Fall17 : public AnyMVAEstimatorRun2Base{
   bool withIso_;
 
   bool debug_;
+
+  const bool absTrackClusterMatching_;
 };
 
 #endif

RecoEgamma/ElectronIdentification/plugins/ElectronMVAEstimatorRun2Phys14NonTrig.cc

@@ -43,7 +43,7 @@ mvaValue( const edm::Ptr<reco::Candidate>& particle, const edm::Event& evt) cons
 
   const int iCategory = findCategory( particle );
   const std::vector<float> vars = fillMVAVariables( particle, evt );
-  const float result = _gbrForests.at(iCategory)->GetClassifier(vars.data());
+  const float result = _gbrForests.at(iCategory)->GetResponse(vars.data()); // The BDT score
 
   constexpr bool debug = false;
   if(debug) {

RecoEgamma/ElectronIdentification/plugins/ElectronMVAEstimatorRun2Spring15NonTrig.cc

@@ -63,7 +63,7 @@ mvaValue( const edm::Ptr<reco::Candidate>& particle, const edm::Event& iEvent) c
 
   const int iCategory = findCategory( particle );
   const std::vector<float> vars = fillMVAVariables( particle, iEvent );
-  const float result = _gbrForests.at(iCategory)->GetClassifier(vars.data());
+  const float result = _gbrForests.at(iCategory)->GetResponse(vars.data()); // The BDT score
 
   const bool debug = false;
   if(debug) {

RecoEgamma/ElectronIdentification/plugins/ElectronMVAEstimatorRun2Spring15Trig.cc

@@ -63,7 +63,7 @@ mvaValue( const edm::Ptr<reco::Candidate>& particle, const edm::Event& iEvent) c
 
   const int iCategory = findCategory( particle );
   const std::vector<float> vars = fillMVAVariables( particle, iEvent );
-  const float result = _gbrForests.at(iCategory)->GetClassifier(vars.data());
+  const float result = _gbrForests.at(iCategory)->GetResponse(vars.data()); // The BDT score
 
   const bool debug = false;
   if(debug) {

RecoEgamma/ElectronIdentification/plugins/cuts/GsfEleRelPFIsoScaledCut.cc

@@ -0,0 +1,77 @@
+#include "PhysicsTools/SelectorUtils/interface/CutApplicatorWithEventContentBase.h"
+#include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
+#include "RecoEgamma/EgammaTools/interface/EffectiveAreas.h"
+
+
+class GsfEleRelPFIsoScaledCut : public CutApplicatorWithEventContentBase {
+  public:
+    GsfEleRelPFIsoScaledCut(const edm::ParameterSet& c);
+
+    result_type operator()(const reco::GsfElectronPtr&) const final;
+
+    void setConsumes(edm::ConsumesCollector&) final;
+    void getEventContent(const edm::EventBase&) final;
+
+    double value(const reco::CandidatePtr& cand) const final;
+
+    CandidateType candidateType() const final { 
+      return ELECTRON; 
+    }
+
+  private:
+    const float barrelC0_, endcapC0_, barrelCpt_, endcapCpt_, barrelCutOff_;
+    EffectiveAreas effectiveAreas_;
+    edm::Handle<double> rhoHandle_;
+};
+
+DEFINE_EDM_PLUGIN(CutApplicatorFactory, GsfEleRelPFIsoScaledCut, "GsfEleRelPFIsoScaledCut");
+
+GsfEleRelPFIsoScaledCut::GsfEleRelPFIsoScaledCut(const edm::ParameterSet& c) :
+  CutApplicatorWithEventContentBase(c),
+  barrelC0_(c.getParameter<double>("barrelC0")),
+  endcapC0_(c.getParameter<double>("endcapC0")),
+  barrelCpt_(c.getParameter<double>("barrelCpt")),
+  endcapCpt_(c.getParameter<double>("endcapCpt")),
+  barrelCutOff_(c.getParameter<double>("barrelCutOff")),
+  effectiveAreas_((c.getParameter<edm::FileInPath>("effAreasConfigFile")).fullPath())
+{
+  edm::InputTag rhoTag = c.getParameter<edm::InputTag>("rho");    
+  contentTags_.emplace("rho",rhoTag);  
+}
+
+void GsfEleRelPFIsoScaledCut::setConsumes(edm::ConsumesCollector& cc){
+  auto rho = cc.consumes<double>(contentTags_["rho"]);
+  contentTokens_.emplace("rho", rho);
+}
+
+void GsfEleRelPFIsoScaledCut::getEventContent(const edm::EventBase& ev){
+  ev.getByLabel(contentTags_["rho"], rhoHandle_);
+}
+
+CutApplicatorBase::result_type GsfEleRelPFIsoScaledCut::operator()(const reco::GsfElectronPtr& cand) const {  
+  // Establish the cut value
+  double absEta = std::abs(cand->superCluster()->eta());
+
+  const float C0     = (absEta < barrelCutOff_ ? barrelC0_  : endcapC0_);
+  const float Cpt    = (absEta < barrelCutOff_ ? barrelCpt_ : endcapCpt_);
+  const float isoCut = C0+Cpt/cand->pt();
+
+  return value(cand) < isoCut;
+}
+
+double GsfEleRelPFIsoScaledCut::value(const reco::CandidatePtr& cand) const {
+
+  // Establish the cut value
+  reco::GsfElectronPtr ele(cand);
+  double absEta = std::abs(ele->superCluster()->eta());
+
+  // Compute the combined isolation with effective area correction
+  auto pfIso = ele->pfIsolationVariables();
+  const float chad = pfIso.sumChargedHadronPt;
+  const float nhad = pfIso.sumNeutralHadronEt;
+  const float pho  = pfIso.sumPhotonEt;
+  const float  eA  = effectiveAreas_.getEffectiveArea(absEta);
+  const float rho  = rhoHandle_.isValid() ? (float)(*rhoHandle_) : 0; // std::max likes float arguments
+  const float iso  = chad + std::max(0.0f, nhad + pho - rho*eA);
+  return iso/cand->pt();
+}

RecoEgamma/ElectronIdentification/python/ElectronMVAValueMapProducer_cfi.py

@@ -31,6 +31,14 @@
     import mvaEleID_Fall17_iso_V1_producer_config
 mvaConfigsForEleProducer.append( mvaEleID_Fall17_iso_V1_producer_config )
 
+from RecoEgamma.ElectronIdentification.Identification.mvaElectronID_Fall17_noIso_V2_cff \
+    import mvaEleID_Fall17_noIso_V2_producer_config
+mvaConfigsForEleProducer.append( mvaEleID_Fall17_noIso_V2_producer_config )
+
+from RecoEgamma.ElectronIdentification.Identification.mvaElectronID_Fall17_iso_V2_cff \
+    import mvaEleID_Fall17_iso_V2_producer_config
+mvaConfigsForEleProducer.append( mvaEleID_Fall17_iso_V2_producer_config )
+
 
 electronMVAValueMapProducer = cms.EDProducer('ElectronMVAValueMapProducer',
                                              # The module automatically detects AOD vs miniAOD, so we configure both