From cbc38e0cc8bfcbb932d3b9917ea558148292d7ce Mon Sep 17 00:00:00 2001 From: Ilya Kravchenko Date: Mon, 14 Sep 2015 00:18:16 +0200 Subject: [PATCH 1/6] New version of this ID, WP Tight EB is revised only, to be tighter than HLT --- ...cutBasedElectronID_Spring15_50ns_V2_cff.py | 189 ++++++++++++++++++ 1 file changed, 189 insertions(+) create mode 100644 RecoEgamma/ElectronIdentification/python/Identification/cutBasedElectronID_Spring15_50ns_V2_cff.py diff --git a/RecoEgamma/ElectronIdentification/python/Identification/cutBasedElectronID_Spring15_50ns_V2_cff.py b/RecoEgamma/ElectronIdentification/python/Identification/cutBasedElectronID_Spring15_50ns_V2_cff.py new file mode 100644 index 0000000000000..e92987e7c851e --- /dev/null +++ b/RecoEgamma/ElectronIdentification/python/Identification/cutBasedElectronID_Spring15_50ns_V2_cff.py @@ -0,0 +1,189 @@ +from PhysicsTools.SelectorUtils.centralIDRegistry import central_id_registry + +# Common functions and classes for ID definition are imported here: +from RecoEgamma.ElectronIdentification.Identification.cutBasedElectronID_tools import * + +# +# This is the first round of Spring15 50ns cuts, optimized on Spring15 50ns samples. +# +# The ID cuts below are optimized IDs for Spring15 Scenario with 50ns bunch spacing +# The cut values are taken from the twiki: +# https://twiki.cern.ch/twiki/bin/view/CMS/CutBasedElectronIdentificationRun2 +# (where they may not stay, if a newer version of cuts becomes available for these +# conditions) +# See also the presentation explaining these working points (this will not change): +# https://indico.cern.ch/event/369245/contribution/2/attachments/1153005/1655984/Rami_eleCB_ID_50ns_V2.pdf +# NOTE: this V2 version is different from the V1 version for these conditions only +# by the cuts of the WP Tight for the barrel. All other WP and endcap are the same as V1. +# The changes was needed to make the WP Tight EB tighter than the HLT. For reference, the full V1 talk is here: +# https://indico.cern.ch/event/369239/contribution/6/attachments/1134836/1623383/Rami_eleCB_ID_50ns.pdf +# +# First, define cut values +# + +# Veto working point Barrel and Endcap +idName = "cutBasedElectronID-Spring15-50ns-V2-standalone-veto" +WP_Veto_EB = EleWorkingPoint_V2( + idName , # idName + 0.0126 , # dEtaInCut + 0.107 , # dPhiInCut + 0.012 , # full5x5_sigmaIEtaIEtaCut + 0.186 , # hOverECut + 0.0621 , # dxyCut + 0.613 , # dzCut + 0.239 , # absEInverseMinusPInverseCut + 0.161 , # relCombIsolationWithEALowPtCut + 0.161 , # relCombIsolationWithEAHighPtCut + # conversion veto cut needs no parameters, so not mentioned + 2 # missingHitsCut + ) + +WP_Veto_EE = EleWorkingPoint_V2( + idName , # idName + 0.0109 , # dEtaInCut + 0.219 , # dPhiInCut + 0.0339 , # full5x5_sigmaIEtaIEtaCut + 0.0962 , # hOverECut + 0.279 , # dxyCut + 0.947 , # dzCut + 0.141 , # absEInverseMinusPInverseCut + 0.193 , # relCombIsolationWithEALowPtCut + 0.193 , # relCombIsolationWithEAHighPtCut + # conversion veto cut needs no parameters, so not mentioned + 3 # missingHitsCut + ) + +# Loose working point Barrel and Endcap +idName = "cutBasedElectronID-Spring15-50ns-V2-standalone-loose" +WP_Loose_EB = EleWorkingPoint_V2( + idName , # idName + 0.0098 , # dEtaInCut + 0.0929 , # dPhiInCut + 0.0105 , # full5x5_sigmaIEtaIEtaCut + 0.0765 , # hOverECut + 0.0227 , # dxyCut + 0.379 , # dzCut + 0.184 , # absEInverseMinusPInverseCut + 0.118 , # relCombIsolationWithEALowPtCut + 0.118 , # relCombIsolationWithEAHighPtCut + # conversion veto cut needs no parameters, so not mentioned + 2 # missingHitsCut + ) + +WP_Loose_EE = EleWorkingPoint_V2( + idName , # idName + 0.00950 , # dEtaInCut + 0.181 , # dPhiInCut + 0.0318 , # full5x5_sigmaIEtaIEtaCut + 0.0824 , # hOverECut + 0.242 , # dxyCut + 0.921 , # dzCut + 0.125 , # absEInverseMinusPInverseCut + 0.118 , # relCombIsolationWithEALowPtCut + 0.118 , # relCombIsolationWithEAHighPtCut + # conversion veto cut needs no parameters, so not mentioned + 1 # missingHitsCut + ) + +# Medium working point Barrel and Endcap +idName = "cutBasedElectronID-Spring15-50ns-V2-standalone-medium" +WP_Medium_EB = EleWorkingPoint_V2( + idName , # idName + 0.00945 , # dEtaInCut + 0.0296 , # dPhiInCut + 0.0101 , # full5x5_sigmaIEtaIEtaCut + 0.0372 , # hOverECut + 0.0151 , # dxyCut + 0.238 , # dzCut + 0.118 , # absEInverseMinusPInverseCut + 0.0987 , # relCombIsolationWithEALowPtCut + 0.0987 , # relCombIsolationWithEAHighPtCut + # conversion veto cut needs no parameters, so not mentioned + 2 # missingHitsCut + ) + +WP_Medium_EE = EleWorkingPoint_V2( + idName , # idName + 0.00773 , # dEtaInCut + 0.148 , # dPhiInCut + 0.0287 , # full5x5_sigmaIEtaIEtaCut + 0.0546 , # hOverECut + 0.0535 , # dxyCut + 0.572 , # dzCut + 0.104 , # absEInverseMinusPInverseCut + 0.0902 , # relCombIsolationWithEALowPtCut + 0.0902 , # relCombIsolationWithEAHighPtCut + # conversion veto cut needs no parameters, so not mentioned + 1 # missingHitsCut + ) + +# Tight working point Barrel and Endcap +idName = "cutBasedElectronID-Spring15-50ns-V2-standalone-tight" +WP_Tight_EB = EleWorkingPoint_V2( + idName , # idName + 0.00864 , # dEtaInCut + 0.0286 , # dPhiInCut + 0.0101 , # full5x5_sigmaIEtaIEtaCut + 0.0342 , # hOverECut + 0.0103 , # dxyCut + 0.170 , # dzCut + 0.0116 , # absEInverseMinusPInverseCut + 0.0591 , # relCombIsolationWithEALowPtCut + 0.0591 , # relCombIsolationWithEAHighPtCut + # conversion veto cut needs no parameters, so not mentioned + 2 # missingHitsCut + ) + +WP_Tight_EE = EleWorkingPoint_V2( + idName , # idName + 0.00762 , # dEtaInCut + 0.0439 , # dPhiInCut + 0.0287 , # full5x5_sigmaIEtaIEtaCut + 0.0544 , # hOverECut + 0.0377 , # dxyCut + 0.571 , # dzCut + 0.0100 , # absEInverseMinusPInverseCut + 0.0759 , # relCombIsolationWithEALowPtCut + 0.0759 , # relCombIsolationWithEAHighPtCut + # conversion veto cut needs no parameters, so not mentioned + 1 # missingHitsCut + ) + +# Second, define what effective areas to use for pile-up correction +isoInputs = IsolationCutInputs_V2( + # phoIsolationEffAreas + "RecoEgamma/ElectronIdentification/data/Spring15/effAreaElectrons_cone03_pfNeuHadronsAndPhotons_50ns.txt" +) + + +# +# Set up VID configuration for all cuts and working points +# + +cutBasedElectronID_Spring15_50ns_V2_standalone_veto = configureVIDCutBasedEleID_V2(WP_Veto_EB, WP_Veto_EE, isoInputs) +cutBasedElectronID_Spring15_50ns_V2_standalone_loose = configureVIDCutBasedEleID_V2(WP_Loose_EB, WP_Loose_EE, isoInputs) +cutBasedElectronID_Spring15_50ns_V2_standalone_medium = configureVIDCutBasedEleID_V2(WP_Medium_EB, WP_Medium_EE, isoInputs) +cutBasedElectronID_Spring15_50ns_V2_standalone_tight = configureVIDCutBasedEleID_V2(WP_Tight_EB, WP_Tight_EE, isoInputs) + +# The MD5 sum numbers below reflect the exact set of cut variables +# and values above. If anything changes, one has to +# 1) comment out the lines below about the registry, +# 2) run "calculateMD5 +# 3) update the MD5 sum strings below and uncomment the lines again. +# + +central_id_registry.register(cutBasedElectronID_Spring15_50ns_V2_standalone_veto.idName, + '3edb92d2aee3dcf1e366a927a5660155') +central_id_registry.register(cutBasedElectronID_Spring15_50ns_V2_standalone_loose.idName, + '527c1b1ddcf9061dc4093dc95590e2bb') +central_id_registry.register(cutBasedElectronID_Spring15_50ns_V2_standalone_medium.idName, + '6837f1ac82974f19d2a15041a2e52ebb') +central_id_registry.register(cutBasedElectronID_Spring15_50ns_V2_standalone_tight.idName, + 'ab050ff2bd30d832881bfac03c9d3a8a') + + +# for now until we have a database... +cutBasedElectronID_Spring15_50ns_V2_standalone_veto.isPOGApproved = cms.untracked.bool(True) +cutBasedElectronID_Spring15_50ns_V2_standalone_loose.isPOGApproved = cms.untracked.bool(True) +cutBasedElectronID_Spring15_50ns_V2_standalone_medium.isPOGApproved = cms.untracked.bool(True) +cutBasedElectronID_Spring15_50ns_V2_standalone_tight.isPOGApproved = cms.untracked.bool(True) From 240c6edaad899466dfdac3f441e00f8e4892654e Mon Sep 17 00:00:00 2001 From: Ilya Kravchenko Date: Fri, 18 Sep 2015 15:30:07 -0500 Subject: [PATCH 2/6] First implementation of the triggering electron MVA for Spring15 in C++ --- .../ElectronMVAEstimatorRun2Spring15Trig.h | 147 +++++++ .../ElectronMVAEstimatorRun2Spring15Trig.cc | 414 ++++++++++++++++++ 2 files changed, 561 insertions(+) create mode 100644 RecoEgamma/ElectronIdentification/interface/ElectronMVAEstimatorRun2Spring15Trig.h create mode 100644 RecoEgamma/ElectronIdentification/plugins/ElectronMVAEstimatorRun2Spring15Trig.cc diff --git a/RecoEgamma/ElectronIdentification/interface/ElectronMVAEstimatorRun2Spring15Trig.h b/RecoEgamma/ElectronIdentification/interface/ElectronMVAEstimatorRun2Spring15Trig.h new file mode 100644 index 0000000000000..cdc5f9cda90f4 --- /dev/null +++ b/RecoEgamma/ElectronIdentification/interface/ElectronMVAEstimatorRun2Spring15Trig.h @@ -0,0 +1,147 @@ +#ifndef RecoEgamma_ElectronIdentification_ElectronMVAEstimatorRun2Spring15Trig_H +#define RecoEgamma_ElectronIdentification_ElectronMVAEstimatorRun2Spring15Trig_H + +#include "RecoEgamma/EgammaTools/interface/AnyMVAEstimatorRun2Base.h" + +#include "DataFormats/EgammaCandidates/interface/GsfElectron.h" + +#include "DataFormats/BeamSpot/interface/BeamSpot.h" + +#include "DataFormats/EgammaCandidates/interface/ConversionFwd.h" +#include "DataFormats/EgammaCandidates/interface/Conversion.h" +#include "RecoEgamma/EgammaTools/interface/ConversionTools.h" + +#include "CondFormats/EgammaObjects/interface/GBRForest.h" + +#include +#include +#include +#include "TMVA/Factory.h" +#include "TMVA/Tools.h" +#include "TMVA/Reader.h" + +class ElectronMVAEstimatorRun2Spring15Trig : public AnyMVAEstimatorRun2Base{ + + public: + + // Define here the number and the meaning of the categories + // for this specific MVA + const int nCategories = 3; + enum mvaCategories { + UNDEFINED = -1, + CAT_EB1 = 0, + CAT_EB2 = 1, + CAT_EE = 2 + }; + + // Define the struct that contains all necessary for MVA variables + // Note: all variables have to be floats for TMVA Reader, even if + // the training was done with ints. + struct AllVariables { + // Pure ECAL -> shower shapes + float see; // 0 + float spp; // 1 + float OneMinusE1x5E5x5; // 2 + float R9; // 3 + float etawidth; // 4 + float phiwidth; // 5 + float HoE; // 6 + + // Endcap only variables + float PreShowerOverRaw; // 7 + + //Pure tracking variables + float kfhits; // 8 + float kfchi2; // 9 + float gsfchi2; // 10 + // Energy matching + float fbrem; // 11 + + float gsfhits; // 12 + float expectedMissingInnerHits; // 13 + float convVtxFitProbability; // 14 + + float EoP; // 15 + float eleEoPout; // 16 + float IoEmIoP; // 17 + // Geometrical matchings + float deta; // 18 + float dphi; // 19 + float detacalo; // 20 + + // Spectator variables + // ... none in this version ... + + // Other variables. These ones are not needed for this version + // of MVA, but kept in this data structure for convenience. + float pt; // 21 + float SCeta; //22 + + + }; + + // Constructor and destructor + ElectronMVAEstimatorRun2Spring15Trig(const edm::ParameterSet& conf); + ~ElectronMVAEstimatorRun2Spring15Trig(); + + // Calculation of the MVA value + float mvaValue( const edm::Ptr& particle, const edm::Event&) const override; + + // Utility functions + std::unique_ptr createSingleReader(const int iCategory, + const edm::FileInPath &weightFile); + + virtual int getNCategories() const override { return nCategories; } + bool isEndcapCategory( int category ) const; + virtual const std::string& getName() const override final { return _name; } + virtual const std::string& getTag() const override final { return _tag; } + + // Functions that should work on both pat and reco electrons + // (use the fact that pat::Electron inherits from reco::GsfElectron) + std::vector fillMVAVariables(const edm::Ptr& particle, const edm::Event&) const override; + int findCategory( const edm::Ptr& particle) const override; + // The function below ensures that the variables passed to MVA are + // within reasonable bounds + void constrainMVAVariables(AllVariables&) const; + + // Call this function once after the constructor to declare + // the needed event content pieces to the framework + void setConsumes(edm::ConsumesCollector&&) const override final; + // Call this function once per event to retrieve all needed + // event content pices + + private: + + // MVA name. This is a unique name for this MVA implementation. + // It will be used as part of ValueMap names. + // For simplicity, keep it set to the class name. + const std::string _name = "ElectronMVAEstimatorRun2Spring15Trig"; + // MVA tag. This is an additional string variable to distinguish + // instances of the estimator of this class configured with different + // weight files. + const std::string _tag; + + // Data members + std::vector< std::unique_ptr > _gbrForests; + + // All variables needed by this MVA + const std::string _MethodName; + AllVariables _allMVAVars; + + // + // Declare all tokens that will be needed to retrieve misc + // data from the event content required by this MVA + // + const edm::InputTag _beamSpotLabel; + // Conversions in AOD and miniAOD have different names + const edm::InputTag _conversionsLabelAOD; + const edm::InputTag _conversionsLabelMiniAOD; + + +}; + +DEFINE_EDM_PLUGIN(AnyMVAEstimatorRun2Factory, + ElectronMVAEstimatorRun2Spring15Trig, + "ElectronMVAEstimatorRun2Spring15Trig"); + +#endif diff --git a/RecoEgamma/ElectronIdentification/plugins/ElectronMVAEstimatorRun2Spring15Trig.cc b/RecoEgamma/ElectronIdentification/plugins/ElectronMVAEstimatorRun2Spring15Trig.cc new file mode 100644 index 0000000000000..f37079d8ab714 --- /dev/null +++ b/RecoEgamma/ElectronIdentification/plugins/ElectronMVAEstimatorRun2Spring15Trig.cc @@ -0,0 +1,414 @@ +#include "RecoEgamma/ElectronIdentification/interface/ElectronMVAEstimatorRun2Spring15Trig.h" + +#include "DataFormats/TrackReco/interface/Track.h" +#include "DataFormats/GsfTrackReco/interface/GsfTrack.h" + +#include "DataFormats/PatCandidates/interface/Electron.h" + +#include "FWCore/ParameterSet/interface/FileInPath.h" + +#include "TMath.h" +#include "TMVA/MethodBDT.h" + +ElectronMVAEstimatorRun2Spring15Trig::ElectronMVAEstimatorRun2Spring15Trig(const edm::ParameterSet& conf): + AnyMVAEstimatorRun2Base(conf), + _tag(conf.getParameter("mvaTag")), + _MethodName("BDTG method"), + _beamSpotLabel(conf.getParameter("beamSpot")), + _conversionsLabelAOD(conf.getParameter("conversionsAOD")), + _conversionsLabelMiniAOD(conf.getParameter("conversionsMiniAOD")) { + + const std::vector weightFileNames + = conf.getParameter >("weightFileNames"); + + if( (int)(weightFileNames.size()) != nCategories ) + throw cms::Exception("MVA config failure: ") + << "wrong number of weightfiles" << std::endl; + + _gbrForests.clear(); + // Create a TMVA reader object for each category + for(int i=0; i(_beamSpotLabel); + + // Conversions collection (different names in AOD and miniAOD) + cc.mayConsume(_conversionsLabelAOD); + cc.mayConsume(_conversionsLabelMiniAOD); + + +} + +float ElectronMVAEstimatorRun2Spring15Trig:: +mvaValue( const edm::Ptr& particle, const edm::Event& iEvent) const { + + const int iCategory = findCategory( particle ); + const std::vector vars = std::move( fillMVAVariables( particle, iEvent ) ); + const float result = _gbrForests.at(iCategory)->GetClassifier(vars.data()); + + const bool debug = false; + if(debug) { + std::cout << " *** Inside the class _MethodName " << _MethodName << std::endl; + std::cout << " bin " << iCategory + << " fbrem " << vars[11] + << " kfchi2 " << vars[9] + << " mykfhits " << vars[8] + << " gsfchi2 " << vars[10] + << " deta " << vars[18] + << " dphi " << vars[19] + << " detacalo " << vars[20] + << " see " << vars[0] + << " spp " << vars[1] + << " etawidth " << vars[4] + << " phiwidth " << vars[5] + << " OneMinusE1x5E5x5 " << vars[2] + << " R9 " << vars[3] + << " HoE " << vars[6] + << " EoP " << vars[15] + << " IoEmIoP " << vars[17] + << " eleEoPout " << vars[16] + << " eta " << vars[22] + << " pt " << vars[21] << std::endl; + std::cout << " ### MVA " << result << std::endl; + } + + return result; +} + +int ElectronMVAEstimatorRun2Spring15Trig::findCategory( const edm::Ptr& particle) const { + + // Try to cast the particle into a reco particle. + // This should work for both reco and pat. + const edm::Ptr eleRecoPtr = ( edm::Ptr )particle; + if( eleRecoPtr.get() == NULL ) + throw cms::Exception("MVA failure: ") + << " given particle is expected to be reco::GsfElectron or pat::Electron," << std::endl + << " but appears to be neither" << std::endl; + + float eta = eleRecoPtr->superCluster()->eta(); + + // + // Determine the category + // + int iCategory = UNDEFINED; + const float ebSplit = 0.800;// barrel is split into two regions + const float ebeeSplit = 1.479; // division between barrel and endcap + + if (std::abs(eta) < ebSplit) + iCategory = CAT_EB1; + + if (std::abs(eta) >= ebSplit && std::abs(eta) < ebeeSplit) + iCategory = CAT_EB2; + + if (std::abs(eta) >= ebeeSplit) + iCategory = CAT_EE; + + return iCategory; +} + +bool ElectronMVAEstimatorRun2Spring15Trig:: +isEndcapCategory(int category ) const { + + bool isEndcap = false; + if( category == CAT_EE ) + isEndcap = true; + + return isEndcap; +} + + +std::unique_ptr ElectronMVAEstimatorRun2Spring15Trig:: +createSingleReader(const int iCategory, const edm::FileInPath &weightFile){ + + // + // Create the reader + // + TMVA::Reader tmpTMVAReader( "!Color:Silent:!Error" ); + + // + // Configure all variables and spectators. Note: the order and names + // must match what is found in the xml weights file! + // + // Pure ECAL -> shower shapes + tmpTMVAReader.AddVariable("ele_oldsigmaietaieta", &_allMVAVars.see); + tmpTMVAReader.AddVariable("ele_oldsigmaiphiiphi", &_allMVAVars.spp); + tmpTMVAReader.AddVariable("ele_oldcircularity", &_allMVAVars.OneMinusE1x5E5x5); + tmpTMVAReader.AddVariable("ele_oldr9", &_allMVAVars.R9); + tmpTMVAReader.AddVariable("ele_scletawidth", &_allMVAVars.etawidth); + tmpTMVAReader.AddVariable("ele_sclphiwidth", &_allMVAVars.phiwidth); + tmpTMVAReader.AddVariable("ele_he", &_allMVAVars.HoE); + // Endcap only variables + if( isEndcapCategory(iCategory) ) + tmpTMVAReader.AddVariable("ele_psEoverEraw", &_allMVAVars.PreShowerOverRaw); + + //Pure tracking variables + tmpTMVAReader.AddVariable("ele_kfhits", &_allMVAVars.kfhits); + tmpTMVAReader.AddVariable("ele_kfchi2", &_allMVAVars.kfchi2); + tmpTMVAReader.AddVariable("ele_gsfchi2", &_allMVAVars.gsfchi2); + + // Energy matching + tmpTMVAReader.AddVariable("ele_fbrem", &_allMVAVars.fbrem); + + tmpTMVAReader.AddVariable("ele_gsfhits", &_allMVAVars.gsfhits); + tmpTMVAReader.AddVariable("ele_expected_inner_hits", &_allMVAVars.expectedMissingInnerHits); + tmpTMVAReader.AddVariable("ele_conversionVertexFitProbability", &_allMVAVars.convVtxFitProbability); + + tmpTMVAReader.AddVariable("ele_ep", &_allMVAVars.EoP); + tmpTMVAReader.AddVariable("ele_eelepout", &_allMVAVars.eleEoPout); + tmpTMVAReader.AddVariable("ele_IoEmIop", &_allMVAVars.IoEmIoP); + + // Geometrical matchings + tmpTMVAReader.AddVariable("ele_deltaetain", &_allMVAVars.deta); + tmpTMVAReader.AddVariable("ele_deltaphiin", &_allMVAVars.dphi); + tmpTMVAReader.AddVariable("ele_deltaetaseed", &_allMVAVars.detacalo); + + // Spectator variables + // .... none ... + + // + // Book the method and set up the weights file + // + std::unique_ptr temp( tmpTMVAReader.BookMVA(_MethodName , weightFile.fullPath() ) ); + + return std::unique_ptr ( new GBRForest( dynamic_cast( tmpTMVAReader.FindMVA(_MethodName) ) ) ); +} + +// A function that should work on both pat and reco objects +std::vector ElectronMVAEstimatorRun2Spring15Trig:: +fillMVAVariables(const edm::Ptr& particle, + const edm::Event& iEvent ) const { + + // + // Declare all value maps corresponding to the products we defined earlier + // + edm::Handle theBeamSpot; + edm::Handle conversions; + + // Get data needed for conversion rejection + iEvent.getByLabel(_beamSpotLabel, theBeamSpot); + + // Conversions in miniAOD and AOD have different names, + // but the same type, so we use the same handle with different tokens. + iEvent.getByLabel(_conversionsLabelAOD, conversions); + if( !conversions.isValid() ) + iEvent.getByLabel(_conversionsLabelMiniAOD, conversions); + + // Make sure everything is retrieved successfully + if(! (theBeamSpot.isValid() + && conversions.isValid() ) + ) + throw cms::Exception("MVA failure: ") + << "Failed to retrieve event content needed for this MVA" + << std::endl + << "Check python MVA configuration file." + << std::endl; + + // Try to cast the particle into a reco particle. + // This should work for both reco and pat. + const edm::Ptr eleRecoPtr = ( edm::Ptr )particle; + if( eleRecoPtr.get() == NULL ) + throw cms::Exception("MVA failure: ") + << " given particle is expected to be reco::GsfElectron or pat::Electron," << std::endl + << " but appears to be neither" << std::endl; + + // Both pat and reco particles have exactly the same accessors, so we use a reco ptr + // throughout the code, with a single exception as of this writing, handled separately below. + auto superCluster = eleRecoPtr->superCluster(); + + AllVariables allMVAVars; + + // Pure ECAL -> shower shapes + allMVAVars.see = eleRecoPtr->full5x5_sigmaIetaIeta(); + allMVAVars.spp = eleRecoPtr->full5x5_sigmaIphiIphi(); + allMVAVars.OneMinusE1x5E5x5 = 1. - eleRecoPtr->full5x5_e1x5() / eleRecoPtr->full5x5_e5x5(); + allMVAVars.R9 = eleRecoPtr->full5x5_r9(); + allMVAVars.etawidth = superCluster->etaWidth(); + allMVAVars.phiwidth = superCluster->phiWidth(); + allMVAVars.HoE = eleRecoPtr->hadronicOverEm(); + // Endcap only variables + allMVAVars.PreShowerOverRaw = superCluster->preshowerEnergy() / superCluster->rawEnergy(); + + // To get to CTF track information in pat::Electron, we have to have the pointer + // to pat::Electron, it is not accessible from the pointer to reco::GsfElectron. + // This behavior is reported and is expected to change in the future (post-7.4.5 some time). + bool validKF= false; + reco::TrackRef myTrackRef = eleRecoPtr->closestCtfTrackRef(); + const edm::Ptr elePatPtr(eleRecoPtr); + // Check if this is really a pat::Electron, and if yes, get the track ref from this new + // pointer instead + if( elePatPtr.get() != NULL ) + myTrackRef = elePatPtr->closestCtfTrackRef(); + validKF = (myTrackRef.isAvailable() && (myTrackRef.isNonnull()) ); + + //Pure tracking variables + allMVAVars.kfhits = (validKF) ? myTrackRef->hitPattern().trackerLayersWithMeasurement() : -1. ; + allMVAVars.kfchi2 = (validKF) ? myTrackRef->normalizedChi2() : 0; + allMVAVars.gsfchi2 = eleRecoPtr->gsfTrack()->normalizedChi2(); + + // Energy matching + allMVAVars.fbrem = eleRecoPtr->fbrem(); + + allMVAVars.gsfhits = eleRecoPtr->gsfTrack()->hitPattern().trackerLayersWithMeasurement(); + allMVAVars.expectedMissingInnerHits = eleRecoPtr->gsfTrack() + ->hitPattern().numberOfHits(reco::HitPattern::MISSING_INNER_HITS); + + reco::ConversionRef conv_ref = ConversionTools::matchedConversion(*eleRecoPtr, + conversions, + theBeamSpot->position()); + double vertexFitProbability = -1.; + if(!conv_ref.isNull()) { + const reco::Vertex &vtx = conv_ref.get()->conversionVertex(); if (vtx.isValid()) { + vertexFitProbability = TMath::Prob( vtx.chi2(), vtx.ndof()); + } + } + allMVAVars.convVtxFitProbability = vertexFitProbability; + + allMVAVars.EoP = eleRecoPtr->eSuperClusterOverP(); + allMVAVars.eleEoPout = eleRecoPtr->eEleClusterOverPout(); + float pAtVertex = eleRecoPtr->trackMomentumAtVtx().R(); + allMVAVars.IoEmIoP = (1.0/eleRecoPtr->ecalEnergy()) - (1.0 / pAtVertex ); + + // Geometrical matchings + allMVAVars.deta = eleRecoPtr->deltaEtaSuperClusterTrackAtVtx(); + allMVAVars.dphi = eleRecoPtr->deltaPhiSuperClusterTrackAtVtx(); + allMVAVars.detacalo = eleRecoPtr->deltaEtaSeedClusterTrackAtCalo(); + + // Spectator variables + allMVAVars.pt = eleRecoPtr->pt(); + float scEta = superCluster->eta(); + allMVAVars.SCeta = scEta; + + constrainMVAVariables(allMVAVars); + + std::vector vars; + + if( isEndcapCategory( findCategory( particle ) ) ) { + vars = std::move( packMVAVariables(allMVAVars.see, + allMVAVars.spp, + allMVAVars.OneMinusE1x5E5x5, + allMVAVars.R9, + allMVAVars.etawidth, + allMVAVars.phiwidth, + allMVAVars.HoE, + // Endcap only variables + allMVAVars.PreShowerOverRaw, + //Pure tracking variables + allMVAVars.kfhits, + allMVAVars.kfchi2, + allMVAVars.gsfchi2, + // Energy matching + allMVAVars.fbrem, + allMVAVars.gsfhits, + allMVAVars.expectedMissingInnerHits, + allMVAVars.convVtxFitProbability, + allMVAVars.EoP, + allMVAVars.eleEoPout, + allMVAVars.IoEmIoP, + // Geometrical matchings + allMVAVars.deta, + allMVAVars.dphi, + allMVAVars.detacalo, + // Spectator variables + allMVAVars.pt, + allMVAVars.SCeta) + ); + } else { + vars = std::move( packMVAVariables(allMVAVars.see, + allMVAVars.spp, + allMVAVars.OneMinusE1x5E5x5, + allMVAVars.R9, + allMVAVars.etawidth, + allMVAVars.phiwidth, + allMVAVars.HoE, + //Pure tracking variables + allMVAVars.kfhits, + allMVAVars.kfchi2, + allMVAVars.gsfchi2, + // Energy matching + allMVAVars.fbrem, + allMVAVars.gsfhits, + allMVAVars.expectedMissingInnerHits, + allMVAVars.convVtxFitProbability, + allMVAVars.EoP, + allMVAVars.eleEoPout, + allMVAVars.IoEmIoP, + // Geometrical matchings + allMVAVars.deta, + allMVAVars.dphi, + allMVAVars.detacalo, + // Spectator variables + allMVAVars.pt, + allMVAVars.SCeta) + ); + } + return vars; +} + +void ElectronMVAEstimatorRun2Spring15Trig::constrainMVAVariables(AllVariables& allMVAVars) const { + + // Check that variables do not have crazy values + + if(allMVAVars.fbrem < -1.) + allMVAVars.fbrem = -1.; + + allMVAVars.deta = fabs(allMVAVars.deta); + if(allMVAVars.deta > 0.06) + allMVAVars.deta = 0.06; + + + allMVAVars.dphi = fabs(allMVAVars.dphi); + if(allMVAVars.dphi > 0.6) + allMVAVars.dphi = 0.6; + + + if(allMVAVars.EoP > 20.) + allMVAVars.EoP = 20.; + + if(allMVAVars.eleEoPout > 20.) + allMVAVars.eleEoPout = 20.; + + + allMVAVars.detacalo = fabs(allMVAVars.detacalo); + if(allMVAVars.detacalo > 0.2) + allMVAVars.detacalo = 0.2; + + if(allMVAVars.OneMinusE1x5E5x5 < -1.) + allMVAVars.OneMinusE1x5E5x5 = -1; + + if(allMVAVars.OneMinusE1x5E5x5 > 2.) + allMVAVars.OneMinusE1x5E5x5 = 2.; + + + + if(allMVAVars.R9 > 5) + allMVAVars.R9 = 5; + + if(allMVAVars.gsfchi2 > 200.) + allMVAVars.gsfchi2 = 200; + + + if(allMVAVars.kfchi2 > 10.) + allMVAVars.kfchi2 = 10.; + + +} + From f3fadc697054abbeda528eb9bb22921eaa180219 Mon Sep 17 00:00:00 2001 From: Ilya Kravchenko Date: Fri, 18 Sep 2015 15:31:10 -0500 Subject: [PATCH 3/6] Python configuration and addition to VID of the triggering electron MVA V1 for both 50ns and 25ns of Spring15 --- .../python/ElectronMVAValueMapProducer_cfi.py | 6 ++ ...mvaElectronID_Spring15_25ns_Trig_V1_cff.py | 100 ++++++++++++++++++ ...mvaElectronID_Spring15_50ns_Trig_V1_cff.py | 100 ++++++++++++++++++ .../Identification/mvaElectronID_tools.py | 24 +++++ 4 files changed, 230 insertions(+) create mode 100644 RecoEgamma/ElectronIdentification/python/Identification/mvaElectronID_Spring15_25ns_Trig_V1_cff.py create mode 100644 RecoEgamma/ElectronIdentification/python/Identification/mvaElectronID_Spring15_50ns_Trig_V1_cff.py diff --git a/RecoEgamma/ElectronIdentification/python/ElectronMVAValueMapProducer_cfi.py b/RecoEgamma/ElectronIdentification/python/ElectronMVAValueMapProducer_cfi.py index 9023e9dd60616..2899d50b17b8a 100644 --- a/RecoEgamma/ElectronIdentification/python/ElectronMVAValueMapProducer_cfi.py +++ b/RecoEgamma/ElectronIdentification/python/ElectronMVAValueMapProducer_cfi.py @@ -9,6 +9,12 @@ from RecoEgamma.ElectronIdentification.Identification.mvaElectronID_Spring15_25ns_nonTrig_V1_cff import * mvaConfigsForEleProducer.append( mvaEleID_Spring15_25ns_nonTrig_V1_producer_config ) +from RecoEgamma.ElectronIdentification.Identification.mvaElectronID_Spring15_50ns_Trig_V1_cff import * +mvaConfigsForEleProducer.append( mvaEleID_Spring15_50ns_Trig_V1_producer_config ) + +from RecoEgamma.ElectronIdentification.Identification.mvaElectronID_Spring15_25ns_Trig_V1_cff import * +mvaConfigsForEleProducer.append( mvaEleID_Spring15_25ns_Trig_V1_producer_config ) + electronMVAValueMapProducer = cms.EDProducer('ElectronMVAValueMapProducer', # The module automatically detects AOD vs miniAOD, so we configure both # diff --git a/RecoEgamma/ElectronIdentification/python/Identification/mvaElectronID_Spring15_25ns_Trig_V1_cff.py b/RecoEgamma/ElectronIdentification/python/Identification/mvaElectronID_Spring15_25ns_Trig_V1_cff.py new file mode 100644 index 0000000000000..c00b0bdfc0a0d --- /dev/null +++ b/RecoEgamma/ElectronIdentification/python/Identification/mvaElectronID_Spring15_25ns_Trig_V1_cff.py @@ -0,0 +1,100 @@ +from PhysicsTools.SelectorUtils.centralIDRegistry import central_id_registry + +import FWCore.ParameterSet.Config as cms + +# +# In this file we define the locations of the MVA weights, cuts on the MVA values +# for specific working points, and configure those cuts in VID +# + +# +# The following MVA is derived for 25ns Spring15 MC samples for triggering electrons. +# See more documentation in this presentation (P.Pigard): +# https://indico.cern.ch/event/369245/contribution/3/attachments/1153011/1655996/20150910_EID_POG_vAsPresented.pdf +# + +# This MVA implementation class name +mvaSpring15TrigClassName = "ElectronMVAEstimatorRun2Spring15Trig" +# The tag is an extra string attached to the names of the products +# such as ValueMaps that needs to distinguish cases when the same MVA estimator +# class is used with different tuning/weights +mvaTag = "25nsV1" + +# There are 3 categories in this MVA. They have to be configured in this strict order +# (cuts and weight files order): +# 0 EB1 (eta<0.8) +# 1 EB2 (eta>=0.8) +# 2 EE + +mvaSpring15TrigWeightFiles_V1 = cms.vstring( + "RecoEgamma/ElectronIdentification/data/Spring15/EIDmva_EB1_10_oldTrigSpring15_25ns_data_1_VarD_TMVA412_Sig6BkgAll_MG_noSpec_BDT.weights.xml", + "RecoEgamma/ElectronIdentification/data/Spring15/EIDmva_EB2_10_oldTrigSpring15_25ns_data_1_VarD_TMVA412_Sig6BkgAll_MG_noSpec_BDT.weights.xml", + "RecoEgamma/ElectronIdentification/data/Spring15/EIDmva_EE_10_oldTrigSpring15_25ns_data_1_VarD_TMVA412_Sig6BkgAll_MG_noSpec_BDT.weights.xml" + ) + +# Load some common definitions for MVA machinery +from RecoEgamma.ElectronIdentification.Identification.mvaElectronID_tools import * + +# The locatoins of value maps with the actual MVA values and categories +# for all particles. +# The names for the maps are ":Values" +# and ":Categories" +mvaProducerModuleLabel = "electronMVAValueMapProducer" +mvaValueMapName = mvaProducerModuleLabel + ":" + mvaSpring15TrigClassName + mvaTag + "Values" +mvaCategoriesMapName = mvaProducerModuleLabel + ":" + mvaSpring15TrigClassName + mvaTag + "Categories" + +# The working point for this MVA that is expected to have about 90% signal +# efficiency in each category +idName90 = "mvaEleID-Spring15-25ns-Trig-V1-wp90" +MVA_WP90 = EleMVA_3Categories_WP( + idName = idName90, + mvaValueMapName = mvaValueMapName, # map with MVA values for all particles + mvaCategoriesMapName = mvaCategoriesMapName, # map with category index for all particles + cutCategory0 = 0.972153, # EB1 + cutCategory1 = 0.922126, # EB2 + cutCategory2 = 0.610764 # EE + ) + +idName80 = "mvaEleID-Spring15-25ns-Trig-V1-wp80" +MVA_WP80 = EleMVA_3Categories_WP( + idName = idName80, + mvaValueMapName = mvaValueMapName, # map with MVA values for all particles + mvaCategoriesMapName = mvaCategoriesMapName, # map with category index for all particles + cutCategory0 = 0.988153, # EB1 + cutCategory1 = 0.96791 , # EB2 + cutCategory2 = 0.841729 # EE + ) + +# +# Finally, set up VID configuration for all cuts +# + +# Create the PSet that will be fed to the MVA value map producer +mvaEleID_Spring15_25ns_Trig_V1_producer_config = cms.PSet( + mvaName = cms.string(mvaSpring15TrigClassName), + mvaTag = cms.string(mvaTag), + # This MVA uses conversion info, so configure several data items on that + beamSpot = cms.InputTag('offlineBeamSpot'), + conversionsAOD = cms.InputTag('allConversions'), + conversionsMiniAOD = cms.InputTag('reducedEgamma:reducedConversions'), + # + weightFileNames = mvaSpring15TrigWeightFiles_V1 + ) +# Create the VPset's for VID cuts +mvaEleID_Spring15_25ns_Trig_V1_wp90 = configureVIDMVAEleID_V1( MVA_WP90 ) +mvaEleID_Spring15_25ns_Trig_V1_wp80 = configureVIDMVAEleID_V1( MVA_WP80 ) + +# The MD5 sum numbers below reflect the exact set of cut variables +# and values above. If anything changes, one has to +# 1) comment out the lines below about the registry, +# 2) run "calculateMD5 +# 3) update the MD5 sum strings below and uncomment the lines again. +# + +central_id_registry.register(mvaEleID_Spring15_25ns_Trig_V1_wp90.idName, + 'bb430b638bf3a4d970627021b0da63ae') +central_id_registry.register(mvaEleID_Spring15_25ns_Trig_V1_wp80.idName, + '81046ab478185af337be1be9b30948ae') + +mvaEleID_Spring15_25ns_Trig_V1_wp90.isPOGApproved = cms.untracked.bool(True) +mvaEleID_Spring15_25ns_Trig_V1_wp80.isPOGApproved = cms.untracked.bool(True) diff --git a/RecoEgamma/ElectronIdentification/python/Identification/mvaElectronID_Spring15_50ns_Trig_V1_cff.py b/RecoEgamma/ElectronIdentification/python/Identification/mvaElectronID_Spring15_50ns_Trig_V1_cff.py new file mode 100644 index 0000000000000..11dac559af53c --- /dev/null +++ b/RecoEgamma/ElectronIdentification/python/Identification/mvaElectronID_Spring15_50ns_Trig_V1_cff.py @@ -0,0 +1,100 @@ +from PhysicsTools.SelectorUtils.centralIDRegistry import central_id_registry + +import FWCore.ParameterSet.Config as cms + +# +# In this file we define the locations of the MVA weights, cuts on the MVA values +# for specific working points, and configure those cuts in VID +# + +# +# The following MVA is derived for 50ns Spring15 MC samples for triggering electrons. +# See more documentation in this presentation (P.Pigard): +# https://indico.cern.ch/event/369245/contribution/3/attachments/1153011/1655996/20150910_EID_POG_vAsPresented.pdf +# + +# This MVA implementation class name +mvaSpring15TrigClassName = "ElectronMVAEstimatorRun2Spring15Trig" +# The tag is an extra string attached to the names of the products +# such as ValueMaps that needs to distinguish cases when the same MVA estimator +# class is used with different tuning/weights +mvaTag = "50nsV1" + +# There are 3 categories in this MVA. They have to be configured in this strict order +# (cuts and weight files order): +# 0 EB1 (eta<0.8) +# 1 EB2 (eta>=0.8) +# 2 EE + +mvaSpring15TrigWeightFiles_V1 = cms.vstring( + "RecoEgamma/ElectronIdentification/data/Spring15/EIDmva_EB1_10_oldTrigSpring15_50ns_data_1_VarD_TMVA412_Sig6BkgAll_MG_noSpec_BDT.weights.xml", + "RecoEgamma/ElectronIdentification/data/Spring15/EIDmva_EB2_10_oldTrigSpring15_50ns_data_1_VarD_TMVA412_Sig6BkgAll_MG_noSpec_BDT.weights.xml", + "RecoEgamma/ElectronIdentification/data/Spring15/EIDmva_EE_10_oldTrigSpring15_50ns_data_1_VarD_TMVA412_Sig6BkgAll_MG_noSpec_BDT.weights.xml" + ) + +# Load some common definitions for MVA machinery +from RecoEgamma.ElectronIdentification.Identification.mvaElectronID_tools import * + +# The locatoins of value maps with the actual MVA values and categories +# for all particles. +# The names for the maps are ":Values" +# and ":Categories" +mvaProducerModuleLabel = "electronMVAValueMapProducer" +mvaValueMapName = mvaProducerModuleLabel + ":" + mvaSpring15TrigClassName + mvaTag + "Values" +mvaCategoriesMapName = mvaProducerModuleLabel + ":" + mvaSpring15TrigClassName + mvaTag + "Categories" + +# The working point for this MVA that is expected to have about 90% signal +# efficiency in each category +idName90 = "mvaEleID-Spring15-50ns-Trig-V1-wp90" +MVA_WP90 = EleMVA_3Categories_WP( + idName = idName90, + mvaValueMapName = mvaValueMapName, # map with MVA values for all particles + mvaCategoriesMapName = mvaCategoriesMapName, # map with category index for all particles + cutCategory0 = 0.923654, # EB1 + cutCategory1 = 0.777484, # EB2 + cutCategory2 = 0.924466 # EE + ) + +idName80 = "mvaEleID-Spring15-50ns-Trig-V1-wp80" +MVA_WP80 = EleMVA_3Categories_WP( + idName = idName80, + mvaValueMapName = mvaValueMapName, # map with MVA values for all particles + mvaCategoriesMapName = mvaCategoriesMapName, # map with category index for all particles + cutCategory0 = 0.971645, # EB1 + cutCategory1 = 0.913198, # EB2 + cutCategory2 = 0.966658 # EE + ) + +# +# Finally, set up VID configuration for all cuts +# + +# Create the PSet that will be fed to the MVA value map producer +mvaEleID_Spring15_50ns_Trig_V1_producer_config = cms.PSet( + mvaName = cms.string(mvaSpring15TrigClassName), + mvaTag = cms.string(mvaTag), + # This MVA uses conversion info, so configure several data items on that + beamSpot = cms.InputTag('offlineBeamSpot'), + conversionsAOD = cms.InputTag('allConversions'), + conversionsMiniAOD = cms.InputTag('reducedEgamma:reducedConversions'), + # + weightFileNames = mvaSpring15TrigWeightFiles_V1 + ) +# Create the VPset's for VID cuts +mvaEleID_Spring15_50ns_Trig_V1_wp90 = configureVIDMVAEleID_V1( MVA_WP90 ) +mvaEleID_Spring15_50ns_Trig_V1_wp80 = configureVIDMVAEleID_V1( MVA_WP80 ) + +# The MD5 sum numbers below reflect the exact set of cut variables +# and values above. If anything changes, one has to +# 1) comment out the lines below about the registry, +# 2) run "calculateMD5 +# 3) update the MD5 sum strings below and uncomment the lines again. +# + +central_id_registry.register(mvaEleID_Spring15_50ns_Trig_V1_wp90.idName, + '8a781409e526a974b260dfcd9c2035c2') +central_id_registry.register(mvaEleID_Spring15_50ns_Trig_V1_wp80.idName, + 'efbd46fb1259325f5d2d49df5c4cd323') + +mvaEleID_Spring15_50ns_Trig_V1_wp90.isPOGApproved = cms.untracked.bool(True) +mvaEleID_Spring15_50ns_Trig_V1_wp80.isPOGApproved = cms.untracked.bool(True) diff --git a/RecoEgamma/ElectronIdentification/python/Identification/mvaElectronID_tools.py b/RecoEgamma/ElectronIdentification/python/Identification/mvaElectronID_tools.py index 1fb4d555010eb..1b0584dc3b96d 100644 --- a/RecoEgamma/ElectronIdentification/python/Identification/mvaElectronID_tools.py +++ b/RecoEgamma/ElectronIdentification/python/Identification/mvaElectronID_tools.py @@ -4,6 +4,30 @@ # Define simple containers for MVA cut values and related # ======================================================= +class EleMVA_3Categories_WP: + """ + This is a container class to hold MVA cut values for a 3-category MVA + as well as the names of the value maps that contain the MVA values computed + for all particles in a producer upstream. + """ + def __init__(self, + idName, + mvaValueMapName, + mvaCategoriesMapName, + cutCategory0, + cutCategory1, + cutCategory2 + ): + self.idName = idName + self.mvaValueMapName = mvaValueMapName + self.mvaCategoriesMapName = mvaCategoriesMapName + self.cutCategory0 = cutCategory0 + self.cutCategory1 = cutCategory1 + self.cutCategory2 = cutCategory2 + + def getCutValues(self): + return [self.cutCategory0, self.cutCategory1, self.cutCategory2] + class EleMVA_6Categories_WP: """ This is a container class to hold MVA cut values for a 6-category MVA From 4a86cfa434a27a0731a458bb89f9022eb3d8ac32 Mon Sep 17 00:00:00 2001 From: Ilya Kravchenko Date: Mon, 21 Sep 2015 14:20:16 -0500 Subject: [PATCH 4/6] Fixed the wp90 and wp80 example cut values for the 50ns electron triggering MVA --- .../mvaElectronID_Spring15_50ns_Trig_V1_cff.py | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/RecoEgamma/ElectronIdentification/python/Identification/mvaElectronID_Spring15_50ns_Trig_V1_cff.py b/RecoEgamma/ElectronIdentification/python/Identification/mvaElectronID_Spring15_50ns_Trig_V1_cff.py index 11dac559af53c..8372ed5926796 100644 --- a/RecoEgamma/ElectronIdentification/python/Identification/mvaElectronID_Spring15_50ns_Trig_V1_cff.py +++ b/RecoEgamma/ElectronIdentification/python/Identification/mvaElectronID_Spring15_50ns_Trig_V1_cff.py @@ -50,9 +50,9 @@ idName = idName90, mvaValueMapName = mvaValueMapName, # map with MVA values for all particles mvaCategoriesMapName = mvaCategoriesMapName, # map with category index for all particles - cutCategory0 = 0.923654, # EB1 - cutCategory1 = 0.777484, # EB2 - cutCategory2 = 0.924466 # EE + cutCategory0 = 0.953843, # EB1 + cutCategory1 = 0.849994, # EB2 + cutCategory2 = 0.514118 # EE ) idName80 = "mvaEleID-Spring15-50ns-Trig-V1-wp80" @@ -60,9 +60,9 @@ idName = idName80, mvaValueMapName = mvaValueMapName, # map with MVA values for all particles mvaCategoriesMapName = mvaCategoriesMapName, # map with category index for all particles - cutCategory0 = 0.971645, # EB1 - cutCategory1 = 0.913198, # EB2 - cutCategory2 = 0.966658 # EE + cutCategory0 = 0.981841, # EB1 + cutCategory1 = 0.946762, # EB2 + cutCategory2 = 0.79704 # EE ) # @@ -92,9 +92,9 @@ # central_id_registry.register(mvaEleID_Spring15_50ns_Trig_V1_wp90.idName, - '8a781409e526a974b260dfcd9c2035c2') + 'b2c90992b008ddca2daf74fe578e7e1b') central_id_registry.register(mvaEleID_Spring15_50ns_Trig_V1_wp80.idName, - 'efbd46fb1259325f5d2d49df5c4cd323') + 'e9127ba0788d88aa64abe68b13749519') mvaEleID_Spring15_50ns_Trig_V1_wp90.isPOGApproved = cms.untracked.bool(True) mvaEleID_Spring15_50ns_Trig_V1_wp80.isPOGApproved = cms.untracked.bool(True) From 7ae911ea0aa6e36fc55863a3c239fca4d2ee8536 Mon Sep 17 00:00:00 2001 From: Ilya Kravchenko Date: Fri, 9 Oct 2015 19:02:19 +0200 Subject: [PATCH 5/6] A fix of 50ns cut-based photon ID: by mistake instead of the exponential, a linear scaling was used for neutral hadron isolation cut for the endcap. Fixed now. --- .../PhoAnyPFIsoWithEAAndExpoScalingCut.cc | 175 ++++++++++++++++++ .../cutBasedPhotonID_Spring15_50ns_V1_cff.py | 14 +- .../Identification/cutBasedPhotonID_tools.py | 55 ++++++ 3 files changed, 237 insertions(+), 7 deletions(-) create mode 100644 RecoEgamma/PhotonIdentification/plugins/cuts/PhoAnyPFIsoWithEAAndExpoScalingCut.cc diff --git a/RecoEgamma/PhotonIdentification/plugins/cuts/PhoAnyPFIsoWithEAAndExpoScalingCut.cc b/RecoEgamma/PhotonIdentification/plugins/cuts/PhoAnyPFIsoWithEAAndExpoScalingCut.cc new file mode 100644 index 0000000000000..0366a6a6e33e4 --- /dev/null +++ b/RecoEgamma/PhotonIdentification/plugins/cuts/PhoAnyPFIsoWithEAAndExpoScalingCut.cc @@ -0,0 +1,175 @@ +#include "PhysicsTools/SelectorUtils/interface/CutApplicatorWithEventContentBase.h" +#include "DataFormats/EgammaCandidates/interface/Photon.h" +#include "RecoEgamma/EgammaTools/interface/EffectiveAreas.h" + + +class PhoAnyPFIsoWithEAAndExpoScalingCut : public CutApplicatorWithEventContentBase { +public: + PhoAnyPFIsoWithEAAndExpoScalingCut(const edm::ParameterSet& c); + + result_type operator()(const reco::PhotonPtr&) const override final; + + void setConsumes(edm::ConsumesCollector&) override final; + void getEventContent(const edm::EventBase&) override final; + + double value(const reco::CandidatePtr& cand) const override final; + + CandidateType candidateType() const override final { + return PHOTON; + } + +private: + // Cut values + float _C1_EB; + float _C2_EB; + float _C3_EB; + float _C1_EE; + float _C2_EE; + float _C3_EE; + // Configuration + float _barrelCutOff; + bool _useRelativeIso; + // Effective area constants + EffectiveAreas _effectiveAreas; + // The isolations computed upstream + edm::Handle > _anyPFIsoMap; + // The rho + edm::Handle< double > _rhoHandle; + + constexpr static char anyPFIsoWithEA_[] = "anyPFIsoWithEA"; + constexpr static char rhoString_ [] = "rho"; +}; + +constexpr char PhoAnyPFIsoWithEAAndExpoScalingCut::anyPFIsoWithEA_[]; +constexpr char PhoAnyPFIsoWithEAAndExpoScalingCut::rhoString_[]; + +DEFINE_EDM_PLUGIN(CutApplicatorFactory, + PhoAnyPFIsoWithEAAndExpoScalingCut, + "PhoAnyPFIsoWithEAAndExpoScalingCut"); + +PhoAnyPFIsoWithEAAndExpoScalingCut::PhoAnyPFIsoWithEAAndExpoScalingCut(const edm::ParameterSet& c) : + CutApplicatorWithEventContentBase(c), + _C1_EB(c.getParameter("C1_EB")), + _C2_EB(c.getParameter("C2_EB")), + _C3_EB(c.getParameter("C3_EB")), + _C1_EE(c.getParameter("C1_EE")), + _C2_EE(c.getParameter("C2_EE")), + _C3_EE(c.getParameter("C3_EE")), + _barrelCutOff(c.getParameter("barrelCutOff")), + _useRelativeIso(c.getParameter("useRelativeIso")), + _effectiveAreas( (c.getParameter("effAreasConfigFile")).fullPath()) +{ + + edm::InputTag maptag = c.getParameter("anyPFIsoMap"); + contentTags_.emplace(anyPFIsoWithEA_,maptag); + + edm::InputTag rhoTag = c.getParameter("rho"); + contentTags_.emplace(rhoString_,rhoTag); + +} + +void PhoAnyPFIsoWithEAAndExpoScalingCut::setConsumes(edm::ConsumesCollector& cc) { + auto anyPFIsoWithEA = + cc.consumes >(contentTags_[anyPFIsoWithEA_]); + contentTokens_.emplace(anyPFIsoWithEA_,anyPFIsoWithEA); + + auto rho = cc.consumes(contentTags_[rhoString_]); + contentTokens_.emplace(rhoString_, rho); +} + +void PhoAnyPFIsoWithEAAndExpoScalingCut::getEventContent(const edm::EventBase& ev) { + ev.getByLabel(contentTags_[anyPFIsoWithEA_],_anyPFIsoMap); + ev.getByLabel(contentTags_[rhoString_],_rhoHandle); +} + +CutApplicatorBase::result_type +PhoAnyPFIsoWithEAAndExpoScalingCut:: +operator()(const reco::PhotonPtr& cand) const{ + + // in case we are by-value + const std::string& inst_name = contentTags_.find(anyPFIsoWithEA_)->second.instance(); + edm::Ptr pat(cand); + float anyisoval = -1.0; + if( _anyPFIsoMap.isValid() && _anyPFIsoMap->contains( cand.id() ) ) { + anyisoval = (*_anyPFIsoMap)[cand]; + } else if ( _anyPFIsoMap.isValid() && _anyPFIsoMap->idSize() == 1 && + cand.id() == edm::ProductID() ) { + // in case we have spoofed a ptr + //note this must be a 1:1 valuemap (only one product input) + anyisoval = _anyPFIsoMap->begin()[cand.key()]; + } else if ( _anyPFIsoMap.isValid() ){ // throw an exception + anyisoval = (*_anyPFIsoMap)[cand]; + } + + // Figure out the cut value + // The value is generally pt-dependent: C1 + pt * C2 + const float pt = cand->pt(); + + // In this version of the isolation cut we apply + // exponential pt scaling to the barrel isolation cut, + // and linear pt scaling to the endcap isolation cut. + double absEta = std::abs(cand->superCluster()->eta()); + const float isolationCutValue = + ( absEta < _barrelCutOff ? + _C1_EB + exp( pt*_C2_EB + _C3_EB) + : _C1_EE + exp( pt*_C2_EE + _C3_EE) ); + + // Retrieve the variable value for this particle + float anyPFIso = _anyPFIsoMap.isValid() ? anyisoval : pat->userFloat(inst_name); + + // Apply pile-up correction + double eA = _effectiveAreas.getEffectiveArea( absEta ); + double rho = *_rhoHandle; + float anyPFIsoWithEA = std::max(0.0, anyPFIso - rho * eA); + + // Divide by pT if the relative isolation is requested + if( _useRelativeIso ) + anyPFIsoWithEA /= pt; + + // Apply the cut and return the result + return anyPFIsoWithEA < isolationCutValue; +} + +double PhoAnyPFIsoWithEAAndExpoScalingCut:: +value(const reco::CandidatePtr& cand) const { + reco::PhotonPtr pho(cand); + + // in case we are by-value + const std::string& inst_name = contentTags_.find(anyPFIsoWithEA_)->second.instance(); + edm::Ptr pat(cand); + float anyisoval = -1.0; + if( _anyPFIsoMap.isValid() && _anyPFIsoMap->contains( cand.id() ) ) { + anyisoval = (*_anyPFIsoMap)[cand]; + } else if ( _anyPFIsoMap.isValid() && _anyPFIsoMap->idSize() == 1 && + cand.id() == edm::ProductID() ) { + // in case we have spoofed a ptr + //note this must be a 1:1 valuemap (only one product input) + anyisoval = _anyPFIsoMap->begin()[cand.key()]; + } else if ( _anyPFIsoMap.isValid() ){ // throw an exception + anyisoval = (*_anyPFIsoMap)[cand]; + } + + // Figure out the cut value + // The value is generally pt-dependent: C1 + pt * C2 + const float pt = pho->pt(); + + // In this version of the isolation cut we apply + // exponential pt scaling to the barrel isolation cut, + // and linear pt scaling to the endcap isolation cut. + double absEta = std::abs(pho->superCluster()->eta()); + + // Retrieve the variable value for this particle + float anyPFIso = _anyPFIsoMap.isValid() ? anyisoval : pat->userFloat(inst_name); + + // Apply pile-up correction + double eA = _effectiveAreas.getEffectiveArea( absEta ); + double rho = *_rhoHandle; + float anyPFIsoWithEA = std::max(0.0, anyPFIso - rho * eA); + + // Divide by pT if the relative isolation is requested + if( _useRelativeIso ) + anyPFIsoWithEA /= pt; + + // Apply the cut and return the result + return anyPFIsoWithEA; +} diff --git a/RecoEgamma/PhotonIdentification/python/Identification/cutBasedPhotonID_Spring15_50ns_V1_cff.py b/RecoEgamma/PhotonIdentification/python/Identification/cutBasedPhotonID_Spring15_50ns_V1_cff.py index a70d75f428c28..1b5517009ea44 100644 --- a/RecoEgamma/PhotonIdentification/python/Identification/cutBasedPhotonID_Spring15_50ns_V1_cff.py +++ b/RecoEgamma/PhotonIdentification/python/Identification/cutBasedPhotonID_Spring15_50ns_V1_cff.py @@ -43,7 +43,7 @@ 0.00 , # absPFChaHadIsoWithEACut_C2 4.00 , # absPFNeuHadIsoWithEACut_C1 0.0040 , # absPFNeuHadIsoWithEACut_C2 - 0.9402 , # absPFNeuHadIsoWithEACut_C2 + 0.9402 , # absPFNeuHadIsoWithEACut_C3 2.15 , # absPFPhoIsoWithEACut_C1 0.0041 # absPFPhoIsoWithEACut_C2 ) @@ -129,9 +129,9 @@ # # Finally, set up VID configuration for all cuts # -cutBasedPhotonID_Spring15_50ns_V1_standalone_loose = configureVIDCutBasedPhoID_V3 ( WP_Loose_EB, WP_Loose_EE, isoInputs) -cutBasedPhotonID_Spring15_50ns_V1_standalone_medium = configureVIDCutBasedPhoID_V3 ( WP_Medium_EB, WP_Medium_EE, isoInputs) -cutBasedPhotonID_Spring15_50ns_V1_standalone_tight = configureVIDCutBasedPhoID_V3 ( WP_Tight_EB, WP_Tight_EE, isoInputs) +cutBasedPhotonID_Spring15_50ns_V1_standalone_loose = configureVIDCutBasedPhoID_V4 ( WP_Loose_EB, WP_Loose_EE, isoInputs) +cutBasedPhotonID_Spring15_50ns_V1_standalone_medium = configureVIDCutBasedPhoID_V4 ( WP_Medium_EB, WP_Medium_EE, isoInputs) +cutBasedPhotonID_Spring15_50ns_V1_standalone_tight = configureVIDCutBasedPhoID_V4 ( WP_Tight_EB, WP_Tight_EE, isoInputs) ## The MD5 sum numbers below reflect the exact set of cut variables # and values above. If anything changes, one has to @@ -141,11 +141,11 @@ # central_id_registry.register(cutBasedPhotonID_Spring15_50ns_V1_standalone_loose.idName, - '3d50a36a9fe1a807fefffe0e6712210a') + '0ecdbf2d9e48f89051f6dd7f9e6547e8') central_id_registry.register(cutBasedPhotonID_Spring15_50ns_V1_standalone_medium.idName, - '63a4ab695fabdae62764db5c55f57b10') + '9efd6037fab2ff3bb65b7592fac5efde') central_id_registry.register(cutBasedPhotonID_Spring15_50ns_V1_standalone_tight.idName, - 'cb046b1400392c9f5db251b5316a87cb') + 'e6a9a77984738862cda3a2fd966fc05e') cutBasedPhotonID_Spring15_50ns_V1_standalone_loose.isPOGApproved = cms.untracked.bool(True) cutBasedPhotonID_Spring15_50ns_V1_standalone_medium.isPOGApproved = cms.untracked.bool(True) diff --git a/RecoEgamma/PhotonIdentification/python/Identification/cutBasedPhotonID_tools.py b/RecoEgamma/PhotonIdentification/python/Identification/cutBasedPhotonID_tools.py index afd09f620dfc7..60e590c4666fe 100644 --- a/RecoEgamma/PhotonIdentification/python/Identification/cutBasedPhotonID_tools.py +++ b/RecoEgamma/PhotonIdentification/python/Identification/cutBasedPhotonID_tools.py @@ -234,6 +234,32 @@ def psetNeuHadIsoWithEAExpoScalingEBCut(wpEB, wpEE, isoInputs): effAreasConfigFile = cms.FileInPath( isoInputs.neuHadIsolationEffAreas ) ) +# Configure the cut on the neutral hadron isolation that uses +# the exponential pt scaling for both barrel and endcap +def psetNeuHadIsoWithEAExpoScalingCut(wpEB, wpEE, isoInputs): + """ + Arguments: two containers of working point cut values of the type WorkingPoint_* + The third argument contains data for isolation calculation. + """ + return cms.PSet( + cutName = cms.string('PhoAnyPFIsoWithEAAndExpoScalingCut'), # Neutral hadrons isolation block + # Barrel: cut = c1 + expo(pt*c2+c3) + C1_EB = cms.double( wpEB.absPFNeuHadIsoWithEACut_C1 ), + C2_EB = cms.double( wpEB.absPFNeuHadIsoWithEACut_C2 ), + C3_EB = cms.double( wpEB.absPFNeuHadIsoWithEACut_C3 ), + # Endcap: cut = cut = c1 + expo(pt*c2+c3) + C1_EE = cms.double( wpEE.absPFNeuHadIsoWithEACut_C1 ), + C2_EE = cms.double( wpEE.absPFNeuHadIsoWithEACut_C2 ), + C3_EE = cms.double( wpEE.absPFNeuHadIsoWithEACut_C3 ), + anyPFIsoMap = cms.InputTag( isoInputs.neuHadIsolationMapName ), + barrelCutOff = cms.double(ebCutOff), + useRelativeIso = cms.bool(False), + needsAdditionalProducts = cms.bool(True), + isIgnored = cms.bool(False), + rho = cms.InputTag("fixedGridRhoFastjetAll"), + effAreasConfigFile = cms.FileInPath( isoInputs.neuHadIsolationEffAreas ) + ) + # Configure the cut on the photon isolation that uses # the linear pt scaling for barrel and endcap def psetPhoIsoWithEALinScalingCut(wpEB, wpEE, isoInputs): @@ -343,3 +369,32 @@ def configureVIDCutBasedPhoID_V3( wpEB, wpEE, isoInputs ): # return parameterSet +def configureVIDCutBasedPhoID_V4( wpEB, wpEE, isoInputs ): + """ + This function configures the full cms.PSet for a VID ID and returns it. + The inputs: first object is of the type WorkingPoint_V2, second object + is of the type WorkingPoint_V2 as well, first containing the cuts for the + Barrel (EB) and the other one for the Endcap (EE). + The third argument contains data for isolation calculation. + + The V4 with respect to V3 has one change: both barrel and endcap + use the exponential scaling for the neutral hadron isolation cut + (in V3 it was only done for the barrel) + """ + # print "VID: Configuring cut set %s" % wpEB.idName + parameterSet = cms.PSet( + # + idName = cms.string( wpEB.idName ), # same name stored in the _EB and _EE objects + cutFlow = cms.VPSet( + psetMinPtCut(), # pt cut + psetPhoSCEtaMultiRangeCut(), # eta cut + psetPhoSingleTowerHadOverEmCut(wpEB,wpEE), # H/E cut + psetPhoFull5x5SigmaIEtaIEtaCut(wpEB,wpEE), # full 5x5 sigmaIEtaIEta cut + psetChHadIsoWithEALinScalingCut(wpEB,wpEE,isoInputs), # charged hadron isolation cut + psetNeuHadIsoWithEAExpoScalingCut(wpEB,wpEE,isoInputs), # neutral hadron isolation cut + psetPhoIsoWithEALinScalingCut(wpEB,wpEE,isoInputs) # photon isolation cut + ) + ) + # + return parameterSet + From 2ba0e1f8cd060ac1d23646e4bb9cdbdd7694651e Mon Sep 17 00:00:00 2001 From: Lindsey Gray Date: Wed, 14 Oct 2015 09:44:40 +0200 Subject: [PATCH 6/6] update electron unit tests --- RecoEgamma/ElectronIdentification/test/runtests.sh | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/RecoEgamma/ElectronIdentification/test/runtests.sh b/RecoEgamma/ElectronIdentification/test/runtests.sh index d678321462770..b338f0294970f 100755 --- a/RecoEgamma/ElectronIdentification/test/runtests.sh +++ b/RecoEgamma/ElectronIdentification/test/runtests.sh @@ -3,9 +3,11 @@ function die { echo $1: status $2 ; exit $2; } ids_to_test=( "RecoEgamma.ElectronIdentification.Identification.cutBasedElectronID_PHYS14_PU20bx25_V2_cff" "RecoEgamma.ElectronIdentification.Identification.cutBasedElectronID_Spring15_25ns_V1_cff" -"RecoEgamma.ElectronIdentification.Identification.cutBasedElectronID_Spring15_50ns_V1_cff" +"RecoEgamma.ElectronIdentification.Identification.cutBasedElectronID_Spring15_50ns_V2_cff" "RecoEgamma.ElectronIdentification.Identification.heepElectronID_HEEPV60_cff" "RecoEgamma.ElectronIdentification.Identification.mvaElectronID_Spring15_25ns_nonTrig_V1_cff" +"RecoEgamma.ElectronIdentification.Identification.mvaElectronID_Spring15_25ns_Trig_V1_cff" +"RecoEgamma.ElectronIdentification.Identification.mvaElectronID_Spring15_50ns_Trig_V1_cff" ) for id_set in "${ids_to_test[@]}"; do