Merge pull request #33817 from bainbrid/LowPtElectrons_nanoAOD_integr…

…ation

LowPtElectrons: NanoAOD integration

DataFormats/EgammaCandidates/interface/GsfElectron.h

@@ -787,6 +787,7 @@ namespace reco {
     // setters
     void setCorrectedEcalEnergyError(float newEnergyError);
     void setCorrectedEcalEnergy(float newEnergy);
+    void setCorrectedEcalEnergy(float newEnergy, bool rescaleDependentValues);
     void setTrackMomentumError(float trackMomentumError);
     void setP4(P4Kind kind, const LorentzVector &p4, float p4Error, bool setCandidate);
     using RecoCandidate::setP4;
@@ -812,8 +813,8 @@ namespace reco {
     //bool isMomentumCorrected() const { return corrections_.isMomentumCorrected ; }
     float caloEnergy() const { return correctedEcalEnergy(); }
     bool isEnergyScaleCorrected() const { return isEcalEnergyCorrected(); }
-    void correctEcalEnergy(float newEnergy, float newEnergyError) {
-      setCorrectedEcalEnergy(newEnergy);
+    void correctEcalEnergy(float newEnergy, float newEnergyError, bool corrEovP = true) {
+      setCorrectedEcalEnergy(newEnergy, corrEovP);
       setEcalEnergyError(newEnergyError);
     }
     void correctMomentum(const LorentzVector &p4, float trackMomentumError, float p4Error) {

DataFormats/EgammaCandidates/src/GsfElectron.cc

@@ -171,14 +171,18 @@ void GsfElectron::setCorrectedEcalEnergyError(float energyError) {
   corrections_.correctedEcalEnergyError = energyError;
 }
 
-void GsfElectron::setCorrectedEcalEnergy(float newEnergy) {
+void GsfElectron::setCorrectedEcalEnergy(float newEnergy) { setCorrectedEcalEnergy(newEnergy, true); }
+
+void GsfElectron::setCorrectedEcalEnergy(float newEnergy, bool rescaleDependentValues) {
   math::XYZTLorentzVectorD momentum = p4();
   momentum *= newEnergy / momentum.e();
   setP4(momentum);
-  showerShape_.hcalDepth1OverEcal *= corrections_.correctedEcalEnergy / newEnergy;
-  showerShape_.hcalDepth2OverEcal *= corrections_.correctedEcalEnergy / newEnergy;
-  trackClusterMatching_.eSuperClusterOverP *= newEnergy / corrections_.correctedEcalEnergy;
-  corrections_.correctedEcalEnergyError *= newEnergy / corrections_.correctedEcalEnergy;
+  if (corrections_.correctedEcalEnergy > 0. && rescaleDependentValues) {
+    showerShape_.hcalDepth1OverEcal *= corrections_.correctedEcalEnergy / newEnergy;
+    showerShape_.hcalDepth2OverEcal *= corrections_.correctedEcalEnergy / newEnergy;
+    trackClusterMatching_.eSuperClusterOverP *= newEnergy / corrections_.correctedEcalEnergy;
+    corrections_.correctedEcalEnergyError *= newEnergy / corrections_.correctedEcalEnergy;
+  }
   corrections_.correctedEcalEnergy = newEnergy;
   corrections_.isEcalEnergyCorrected = true;
 }

PhysicsTools/NanoAOD/python/lowPtElectrons_cff.py

@@ -0,0 +1,203 @@
+import FWCore.ParameterSet.Config as cms
+from PhysicsTools.NanoAOD.nano_eras_cff import *
+from PhysicsTools.NanoAOD.common_cff import *
+
+################################################################################
+# Modules
+################################################################################
+
+from RecoEgamma.EgammaTools.lowPtElectronModifier_cfi import lowPtElectronModifier
+from RecoEgamma.EgammaElectronProducers.lowPtGsfElectrons_cfi import lowPtRegressionModifier
+modifiedLowPtElectrons = cms.EDProducer(
+    "ModifiedElectronProducer",
+    src = cms.InputTag("slimmedLowPtElectrons"),
+    modifierConfig = cms.PSet(
+        modifications = cms.VPSet(lowPtElectronModifier,lowPtRegressionModifier)
+    )
+)
+
+from RecoEgamma.EgammaElectronProducers.lowPtGsfElectronID_cfi import lowPtGsfElectronID
+lowPtPATElectronID = lowPtGsfElectronID.clone(
+    usePAT = True,
+    electrons = "modifiedLowPtElectrons",
+    unbiased = "",
+    ModelWeights = [
+        'RecoEgamma/ElectronIdentification/data/LowPtElectrons/LowPtElectrons_ID_2020Nov28.root',
+    ],
+)
+
+isoForLowPtEle = cms.EDProducer(
+    "EleIsoValueMapProducer",
+    src = cms.InputTag("modifiedLowPtElectrons"),
+    relative = cms.bool(False),
+    rho_MiniIso = cms.InputTag("fixedGridRhoFastjetAll"),
+    rho_PFIso = cms.InputTag("fixedGridRhoFastjetAll"),
+    EAFile_MiniIso = cms.FileInPath("RecoEgamma/ElectronIdentification/data/Fall17/effAreaElectrons_cone03_pfNeuHadronsAndPhotons_94X.txt"),
+    EAFile_PFIso = cms.FileInPath("RecoEgamma/ElectronIdentification/data/Fall17/effAreaElectrons_cone03_pfNeuHadronsAndPhotons_94X.txt"),
+)
+
+updatedLowPtElectronsWithUserData = cms.EDProducer(
+    "PATElectronUserDataEmbedder",
+    src = cms.InputTag("modifiedLowPtElectrons"),
+    userFloats = cms.PSet(
+        ID = cms.InputTag("lowPtPATElectronID"),
+        miniIsoChg = cms.InputTag("isoForLowPtEle:miniIsoChg"),
+        miniIsoAll = cms.InputTag("isoForLowPtEle:miniIsoAll"),
+    ),
+    userIntFromBools = cms.PSet(),
+    userInts = cms.PSet(),
+    userCands = cms.PSet(),
+)
+
+finalLowPtElectrons = cms.EDFilter(
+    "PATElectronRefSelector",
+    src = cms.InputTag("updatedLowPtElectronsWithUserData"),
+    cut = cms.string("pt > 1. && userFloat('ID') > -0.25"),
+)
+
+################################################################################
+# electronTable 
+################################################################################
+
+lowPtElectronTable = cms.EDProducer(
+    "SimpleCandidateFlatTableProducer",
+    src = cms.InputTag("finalLowPtElectrons"),
+    cut = cms.string(""),
+    name= cms.string("LowPtElectron"),
+    doc = cms.string("slimmedLowPtElectrons after basic selection (" + finalLowPtElectrons.cut.value()+")"),
+    singleton = cms.bool(False), # the number of entries is variable
+    extension = cms.bool(False), # this is the main table for the electrons
+    variables = cms.PSet(
+        # Basic variables
+        CandVars,
+        # BDT scores and WPs
+        embeddedID = Var("electronID('ID')",float,doc="ID, BDT (raw) score"),
+        ID = Var("userFloat('ID')",float,doc="New ID, BDT (raw) score"),
+        unbiased = Var("electronID('unbiased')",float,doc="ElectronSeed, pT- and dxy- agnostic BDT (raw) score"),
+        ptbiased = Var("electronID('ptbiased')",float,doc="ElectronSeed, pT- and dxy- dependent BDT (raw) score"),
+        # Isolation
+        miniPFRelIso_chg = Var("userFloat('miniIsoChg')/pt",float,
+                               doc="mini PF relative isolation, charged component"),
+        miniPFRelIso_all = Var("userFloat('miniIsoAll')/pt",float,
+                               doc="mini PF relative isolation, total (with scaled rho*EA PU corrections)"),
+        # Conversions
+        convVeto = Var("passConversionVeto()",bool,doc="pass conversion veto"),
+        convWP = Var("userInt('convOpen')*1 + userInt('convLoose')*2 + userInt('convTight')*4",
+                     int,doc="conversion flag bit map: 1=Veto, 2=Loose, 3=Tight"),
+        convVtxRadius = Var("userFloat('convVtxRadius')",float,doc="conversion vertex radius (cm)",precision=7),
+        # Tracking
+        lostHits = Var("gsfTrack.hitPattern.numberOfLostHits('MISSING_INNER_HITS')","uint8",doc="number of missing inner hits"),
+        # Cluster-related
+        energyErr = Var("p4Error('P4_COMBINATION')",float,doc="energy error of the cluster-track combination",precision=6),
+        deltaEtaSC = Var("superCluster().eta()-eta()",float,doc="delta eta (SC,ele) with sign",precision=10),
+        r9 = Var("full5x5_r9()",float,doc="R9 of the SC, calculated with full 5x5 region",precision=10),
+        sieie = Var("full5x5_sigmaIetaIeta()",float,doc="sigma_IetaIeta of the SC, calculated with full 5x5 region",precision=10),
+        eInvMinusPInv = Var("(1-eSuperClusterOverP())/ecalEnergy()",float,doc="1/E_SC - 1/p_trk",precision=10),
+        scEtOverPt = Var("(superCluster().energy()/(pt*cosh(superCluster().eta())))-1",float,doc="(SC energy)/pt-1",precision=8),
+        hoe = Var("hadronicOverEm()",float,doc="H over E",precision=8),
+        # Displacement
+        dxy = Var("dB('PV2D')",float,doc="dxy (with sign) wrt first PV, in cm",precision=10),
+        dxyErr = Var("edB('PV2D')",float,doc="dxy uncertainty, in cm",precision=6),
+        dz = Var("dB('PVDZ')",float,doc="dz (with sign) wrt first PV, in cm",precision=10),
+        dzErr = Var("abs(edB('PVDZ'))",float,doc="dz uncertainty, in cm",precision=6),
+        ip3d = Var("abs(dB('PV3D'))",float,doc="3D impact parameter wrt first PV, in cm",precision=10),
+        sip3d = Var("abs(dB('PV3D')/edB('PV3D'))",float,doc="3D impact parameter significance wrt first PV, in cm",precision=10),
+        # Cross-referencing
+        #jetIdx
+        #photonIdx
+    ),
+)
+
+################################################################################
+# electronTable (MC)
+################################################################################
+
+from PhysicsTools.NanoAOD.particlelevel_cff import particleLevel
+particleLevelForMatchingLowPt = particleLevel.clone(
+    lepMinPt = cms.double(1.),
+    phoMinPt = cms.double(1),
+)
+
+tautaggerForMatchingLowPt = cms.EDProducer(
+    "GenJetTauTaggerProducer",
+    src = cms.InputTag('particleLevelForMatchingLowPt:leptons')
+)
+
+matchingLowPtElecPhoton = cms.EDProducer(
+    "GenJetGenPartMerger",
+    srcJet =cms.InputTag("particleLevelForMatchingLowPt:leptons"),
+    srcPart=cms.InputTag("particleLevelForMatchingLowPt:photons"),
+    hasTauAnc=cms.InputTag("tautaggerForMatchingLowPt"),
+)
+
+lowPtElectronsMCMatchForTableAlt = cms.EDProducer(
+    "GenJetMatcherDRPtByDR",                # cut on deltaR, deltaPt/Pt; pick best by deltaR
+    src         = lowPtElectronTable.src,   # final reco collection
+    matched     = cms.InputTag("matchingLowPtElecPhoton:merged"), # final mc-truth particle collection
+    mcPdgId     = cms.vint32(11,22),        # one or more PDG ID (11 = el, 22 = pho); absolute values (see below)
+    checkCharge = cms.bool(False),          # True = require RECO and MC objects to have the same charge
+    mcStatus    = cms.vint32(),
+    maxDeltaR   = cms.double(0.3),          # Minimum deltaR for the match
+    maxDPtRel   = cms.double(0.5),          # Minimum deltaPt/Pt for the match
+    resolveAmbiguities    = cms.bool(True), # Forbid two RECO objects to match to the same GEN object
+    resolveByMatchQuality = cms.bool(True), # False = just match input in order; True = pick lowest deltaR pair first
+) 
+
+lowPtElectronsMCMatchForTable = cms.EDProducer(
+    "MCMatcher",                            # cut on deltaR, deltaPt/Pt; pick best by deltaR
+    src         = lowPtElectronTable.src,   # final reco collection
+    matched     = cms.InputTag("finalGenParticles"), # final mc-truth particle collection
+    mcPdgId     = cms.vint32(11),           # one or more PDG ID (11 = ele); absolute values (see below)
+    checkCharge = cms.bool(False),          # True = require RECO and MC objects to have the same charge
+    mcStatus    = cms.vint32(1),            # PYTHIA status code (1 = stable, 2 = shower, 3 = hard scattering)
+    maxDeltaR   = cms.double(0.3),          # Minimum deltaR for the match
+    maxDPtRel   = cms.double(0.5),          # Minimum deltaPt/Pt for the match
+    resolveAmbiguities    = cms.bool(True), # Forbid two RECO objects to match to the same GEN object
+    resolveByMatchQuality = cms.bool(True), # False = just match input in order; True = pick lowest deltaR pair first
+)
+
+from PhysicsTools.NanoAOD.electrons_cff import electronMCTable
+lowPtElectronMCTable = cms.EDProducer(
+    "CandMCMatchTableProducer",
+    src = lowPtElectronTable.src,
+    mcMapDressedLep = cms.InputTag("lowPtElectronsMCMatchForTableAlt"),
+    mcMap = cms.InputTag("lowPtElectronsMCMatchForTable"),
+    mapTauAnc = cms.InputTag("matchingLowPtElecPhoton:hasTauAnc"),
+    objName = lowPtElectronTable.name,
+    objType = electronMCTable.objType,
+    branchName = cms.string("genPart"),
+    docString = cms.string("MC matching to status==1 electrons or photons"),
+    genparticles = cms.InputTag("finalGenParticles"), 
+)
+
+################################################################################
+# Sequences
+################################################################################
+
+lowPtElectronSequence = cms.Sequence(modifiedLowPtElectrons
+                                     +lowPtPATElectronID
+                                     +isoForLowPtEle
+                                     +updatedLowPtElectronsWithUserData
+                                     +finalLowPtElectrons)
+lowPtElectronTables = cms.Sequence(lowPtElectronTable)
+lowPtElectronMC = cms.Sequence(
+    particleLevelForMatchingLowPt
+    +tautaggerForMatchingLowPt
+    +matchingLowPtElecPhoton
+    +lowPtElectronsMCMatchForTable
+    +lowPtElectronsMCMatchForTableAlt
+    +lowPtElectronMCTable)
+
+################################################################################
+# Modifiers
+################################################################################
+
+_modifiers = ( run2_miniAOD_80XLegacy |
+               run2_nanoAOD_94XMiniAODv1 |
+               run2_nanoAOD_94XMiniAODv2 |
+               run2_nanoAOD_94X2016 |
+               run2_nanoAOD_102Xv1 |
+               run2_nanoAOD_106Xv1 )
+(_modifiers).toReplaceWith(lowPtElectronSequence,cms.Sequence())
+(_modifiers).toReplaceWith(lowPtElectronTables,cms.Sequence())
+(_modifiers).toReplaceWith(lowPtElectronMC,cms.Sequence())

PhysicsTools/NanoAOD/python/nanoDQM_cff.py

@@ -108,6 +108,7 @@
 ## MC
 nanoDQMMC = nanoDQM.clone()
 nanoDQMMC.vplots.Electron.sels.Prompt = cms.string("genPartFlav == 1")
+nanoDQMMC.vplots.LowPtElectron.sels.Prompt = cms.string("genPartFlav == 1")
 nanoDQMMC.vplots.Muon.sels.Prompt = cms.string("genPartFlav == 1")
 nanoDQMMC.vplots.Photon.sels.Prompt = cms.string("genPartFlav == 1")
 nanoDQMMC.vplots.Tau.sels.Prompt = cms.string("genPartFlav == 5")

PhysicsTools/NanoAOD/python/nanoDQM_cfi.py

@@ -112,6 +112,53 @@
                 Plot1D('dEsigmaDown', 'dEsigmaDown', 100, -0.1, 0.1, '#Delta E sigmaDown'),
             )
         ),
+
+        LowPtElectron = cms.PSet(
+            sels = cms.PSet(
+                Good = cms.string('pt > 1. && ID > 5.')
+            ),
+            plots = cms.VPSet(
+                #
+                Count1D('_size', 8, -0.5, 7.5, 'slimmedLowPtElectrons after basic selection'),
+                # CandVars
+                Plot1D('charge', 'charge', 3, -1.5, 1.5, 'electric charge'),
+                Plot1D('eta', 'eta', 20, -3., 3., 'eta'),
+                NoPlot('mass'),
+                Plot1D('pdgId', 'pdgId', 101, -50.5, 50.5, 'PDG code assigned by the event reconstruction (not by MC truth)'),
+                Plot1D('phi', 'phi', 20, -3.14159, 3.14159, 'phi'),
+                Plot1D('pt', 'pt', 40, 0., 20., 'pt (corrected)'),
+                # BDT scores and WPs
+                Plot1D('embeddedID', 'embeddedID', 40, -10., 10., 'Embedded ID, BDT (raw) score'),
+                Plot1D('ID', 'ID', 40, -10., 10., 'ID, BDT (raw) score'),
+                Plot1D('unbiased', 'unbiased', 40, -10., 10., 'ElectronSeed, pT- and dxy- agnostic BDT (raw) score'),
+                Plot1D('ptbiased', 'ptbiased', 40, -10., 10., 'ElectronSeed, pT- and dxy- dependent BDT (raw) score'),
+                # Isolation
+                Plot1D('miniPFRelIso_chg', 'miniPFRelIso_chg', 20, 0, 2, 'mini PF relative isolation, charged component'),
+                Plot1D('miniPFRelIso_all', 'miniPFRelIso_all', 20, 0, 2, 'mini PF relative isolation, total (with scaled rho*EA PU corrections)'),
+                # Conversions
+                Plot1D('convVeto', 'convVeto', 2, -0.5, 1.5, 'pass conversion veto'),
+                Plot1D('convWP', 'convWP', 8, -0.5, 7.5, 'conversion flag bit map: 1=Veto, 2=Loose, 3=Tight'),
+                Plot1D('convVtxRadius', 'convVtxRadius', 40, 0., 20.0, 'conversion vertex radius (cm)'),
+                # Tracking
+                Plot1D('lostHits', 'lostHits', 4, -0.5, 3.5, 'number of missing inner hits'),
+                # Cluster-related
+                Plot1D('energyErr', 'energyErr', 40, 0., 20., 'energy error of the cluster from regression'),
+                Plot1D('deltaEtaSC', 'deltaEtaSC', 20, -0.2, 0.2, 'delta eta (SC,ele) with sign'),
+                Plot1D('r9', 'r9', 20, 0, 1.1, 'R9 of the supercluster, calculated with full 5x5 region'),
+                Plot1D('sieie', 'sieie', 20, 0, 0.05, 'sigma_IetaIeta of the supercluster, calculated with full 5x5 region'),
+                Plot1D('eInvMinusPInv', 'eInvMinusPInv', 20, -0.1, 0.1, '1/E_SC - 1/p_trk'),
+                Plot1D('scEtOverPt', 'scEtOverPt', 20, -0.5, 0.5, '(supercluster transverse energy)/pt - 1'),
+                Plot1D('hoe', 'hoe', 20, 0, 0.6, 'H over E'),
+                # Displacement
+                Plot1D('dxy', 'dxy', 20, -0.1, 0.1, 'dxy (with sign) wrt first PV, in cm'),
+                Plot1D('dz', 'dz', 20, -0.3, 0.3, 'dz (with sign) wrt first PV, in cm'),
+                Plot1D('dxyErr', 'dxyErr', 20, 0., 0.2, 'dxy uncertainty, in cm'),
+                Plot1D('dzErr', 'dzErr', 20, 0., 0.2, 'dz uncertainty, in cm'),
+                Plot1D('ip3d', 'ip3d', 20, 0., 0.2, '3D impact parameter wrt first PV, in cm'),
+                Plot1D('sip3d', 'sip3d', 20, 0., 20., '3D impact parameter significance wrt first PV, in cm'),
+            ),
+        ),
+
         FatJet = cms.PSet(
             sels = cms.PSet(),
             plots = cms.VPSet(

PhysicsTools/NanoAOD/python/nano_cff.py

@@ -7,6 +7,7 @@
 from PhysicsTools.NanoAOD.taus_cff import *
 from PhysicsTools.NanoAOD.boostedTaus_cff import *
 from PhysicsTools.NanoAOD.electrons_cff import *
+from PhysicsTools.NanoAOD.lowPtElectrons_cff import *
 from PhysicsTools.NanoAOD.photons_cff import *
 from PhysicsTools.NanoAOD.globals_cff import *
 from PhysicsTools.NanoAOD.extraflags_cff import *
@@ -106,10 +107,10 @@
 (run2_miniAOD_80XLegacy | run2_nanoAOD_94X2016 | run2_nanoAOD_94XMiniAODv1 | run2_nanoAOD_94XMiniAODv2 | run2_nanoAOD_102Xv1).toModify(l1bits, storeUnprefireableBit=False)
 
 nanoSequenceCommon = cms.Sequence(
-        nanoMetadata + jetSequence + muonSequence + tauSequence + boostedTauSequence + electronSequence+photonSequence+vertexSequence+
+        nanoMetadata + jetSequence + muonSequence + tauSequence + boostedTauSequence + electronSequence + lowPtElectronSequence + photonSequence+vertexSequence+
         isoTrackSequence + jetLepSequence + # must be after all the leptons
         linkedObjects  +
-        jetTables + muonTables + tauTables + boostedTauTables + electronTables + photonTables +  globalTables +vertexTables+ metTables+simpleCleanerTable + isoTrackTables
+        jetTables + muonTables + tauTables + boostedTauTables + electronTables + lowPtElectronTables + photonTables +  globalTables +vertexTables+ metTables+simpleCleanerTable + isoTrackTables
         )
 #remove boosted tau from previous eras
 (run2_miniAOD_80XLegacy | run2_nanoAOD_92X | run2_nanoAOD_94XMiniAODv1 | run2_nanoAOD_94X2016 | run2_nanoAOD_94XMiniAODv2 | run2_nanoAOD_102Xv1 | run2_nanoAOD_106Xv1).toReplaceWith(nanoSequenceCommon, nanoSequenceCommon.copyAndExclude([boostedTauSequence, boostedTauTables]))
@@ -119,7 +120,7 @@
 
 nanoSequence = cms.Sequence(nanoSequenceCommon + nanoSequenceOnlyData + nanoSequenceOnlyFullSim)
 
-nanoSequenceFS = cms.Sequence(genParticleSequence + genVertexTables + particleLevelSequence + nanoSequenceCommon + jetMC + muonMC + electronMC + photonMC + tauMC + boostedTauMC + metMC + ttbarCatMCProducers +  globalTablesMC + btagWeightTable + genWeightsTable + genVertexTable + genParticleTables + particleLevelTables + lheInfoTable  + ttbarCategoryTable )
+nanoSequenceFS = cms.Sequence(genParticleSequence + genVertexTables + particleLevelSequence + nanoSequenceCommon + jetMC + muonMC + electronMC + lowPtElectronMC + photonMC + tauMC + boostedTauMC + metMC + ttbarCatMCProducers +  globalTablesMC + btagWeightTable + genWeightsTable + genVertexTable + genParticleTables + particleLevelTables + lheInfoTable  + ttbarCategoryTable )
 
 (run2_nanoAOD_92X | run2_miniAOD_80XLegacy | run2_nanoAOD_94X2016 | run2_nanoAOD_94X2016 | \
     run2_nanoAOD_94XMiniAODv1 | run2_nanoAOD_94XMiniAODv2 | \

PhysicsTools/PatAlgos/python/slimming/slimmedLowPtElectrons_cfi.py

@@ -1,5 +1,7 @@
 import FWCore.ParameterSet.Config as cms
 
+from RecoEgamma.EgammaTools.lowPtElectronModifier_cfi import lowPtElectronModifier
+
 slimmedLowPtElectrons = cms.EDProducer("PATElectronSlimmer",
    src = cms.InputTag("selectedPatLowPtElectrons"),                                  
    dropSuperCluster = cms.string("0"), # you can put a cut to slim selectively, e.g. pt < 10
@@ -39,6 +41,7 @@
                 modifierName = cms.string('EGExtraInfoModifierFromPackedCandPtrValueMaps'),
                 photon_config = cms.PSet()
             ),
+            lowPtElectronModifier,
         )
    )
 )