Skip to content

Commit

Permalink
Merge pull request cms-sw#186 from gpetruc/CMGTools-from-CMSSW_7_2_3
Browse files Browse the repository at this point in the history
Upstream updates from cms-sw#185
  • Loading branch information
gpetruc committed Dec 19, 2014
2 parents 9e0742b + e55fbe1 commit 2440eb3
Show file tree
Hide file tree
Showing 19 changed files with 630 additions and 599 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -20,8 +20,6 @@ def __init__(self, cfg_ana, cfg_comp, looperName):
if not getattr(self.cfg_ana, 'saveTLorentzVectors', False):
fourVectorType.removeVariable("p4")

## Declare how we store floats by default
self.tree.setDefaultFloatType("F"); # otherwise it's "D"

self.collections = {}
self.globalObjects = {}
Expand All @@ -33,6 +31,11 @@ def __init__(self, cfg_ana, cfg_comp, looperName):
if hasattr(cfg_ana,"globalVariables"):
self.globalVariables=cfg_ana.globalVariables

def beginLoop(self, setup) :
super(AutoFillTreeProducer, self).beginLoop(setup)
## Declare how we store floats by default
self.tree.setDefaultFloatType("F"); # otherwise it's "D"

def declareHandles(self):
super(AutoFillTreeProducer, self).declareHandles()
# self.handles['TriggerResults'] = AutoHandle( ('TriggerResults','','HLT'), 'edm::TriggerResults' )
Expand Down
20 changes: 13 additions & 7 deletions PhysicsTools/Heppy/python/analyzers/core/TreeAnalyzerNumpy.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,17 +11,22 @@ class TreeAnalyzerNumpy( Analyzer ):

def __init__(self, cfg_ana, cfg_comp, looperName):
super(TreeAnalyzerNumpy,self).__init__(cfg_ana, cfg_comp, looperName)
fileName = '/'.join([self.dirName,
'tree.root'])

isCompressed = self.cfg_ana.isCompressed if hasattr(cfg_ana,'isCompressed') else 1
print 'Compression', isCompressed

self.file = TFile( fileName, 'recreate', '', isCompressed )
self.tree = Tree('tree', self.name)

def beginLoop(self, setup) :
super(TreeAnalyzerNumpy, self).beginLoop(setup)
print setup.services
if "outputfile" in setup.services:
print "Using outputfile given in setup.outputfile"
self.file = setup.services["outputfile"].file
else :
fileName = '/'.join([self.dirName,
'tree.root'])
isCompressed = self.cfg_ana.isCompressed if hasattr(self.cfg_ana,'isCompressed') else 1
print 'Compression', isCompressed
self.file = TFile( fileName, 'recreate', '', isCompressed )
self.tree = Tree('tree', self.name)
self.declareVariables(setup)

def declareVariables(self,setup):
Expand All @@ -30,5 +35,6 @@ def declareVariables(self,setup):

def write(self, setup):
super(TreeAnalyzerNumpy, self).write(setup)
self.file.Write()
if "outputfile" not in setup.services:
self.file.Write()

55 changes: 54 additions & 1 deletion PhysicsTools/Heppy/python/analyzers/core/TriggerBitAnalyzer.py
Original file line number Diff line number Diff line change
Expand Up @@ -3,12 +3,15 @@
from PhysicsTools.Heppy.analyzers.core.Analyzer import Analyzer
from PhysicsTools.Heppy.analyzers.core.AutoHandle import AutoHandle
from PhysicsTools.Heppy.analyzers.core.AutoFillTreeProducer import NTupleVariable
import PhysicsTools.HeppyCore.framework.config as cfg

class TriggerBitAnalyzer( Analyzer ):
def __init__(self, cfg_ana, cfg_comp, looperName ):
super(TriggerBitAnalyzer,self).__init__(cfg_ana,cfg_comp,looperName)
self.processName = getattr(self.cfg_ana,"processName","HLT")
self.outprefix = getattr(self.cfg_ana,"outprefix", self.processName)
self.unrollbits = ( hasattr(self.cfg_ana,"unrollbits") and self.cfg_ana.unrollbits )


def declareHandles(self):
super(TriggerBitAnalyzer, self).declareHandles()
Expand All @@ -17,22 +20,72 @@ def declareHandles(self):
def beginLoop(self, setup):
super(TriggerBitAnalyzer,self).beginLoop(setup)
self.triggerBitCheckers = []
if self.unrollbits :
self.allPaths = set()
self.triggerBitCheckersSingleBits = []

for T, TL in self.cfg_ana.triggerBits.iteritems():
trigVec = ROOT.vector(ROOT.string)()
for TP in TL:
trigVec.push_back(TP)
outname="%s_%s"%(self.outprefix,T)
if self.unrollbits :
if TP not in self.allPaths :
self.allPaths.update(TP)
trigVecBit = ROOT.vector(ROOT.string)()
trigVecBit.push_back(TP)
outname="%s_BIT_%s"%(self.outprefix,TP)
if not hasattr(setup ,"globalVariables") :
setup.globalVariables = []
setup.globalVariables.append( NTupleVariable(outname, eval("lambda ev: ev.%s" % outname), help="Trigger bit %s"%TP) )
self.triggerBitCheckersSingleBits.append( (TP, ROOT.heppy.TriggerBitChecker(trigVecBit)) )

outname="%s_%s"%(self.outprefix,T)
if not hasattr(setup ,"globalVariables") :
setup.globalVariables = []
setup.globalVariables.append( NTupleVariable(outname, eval("lambda ev: ev.%s" % outname), help="OR of %s"%TL) )
self.triggerBitCheckers.append( (T, ROOT.heppy.TriggerBitChecker(trigVec)) )


def process(self, event):
self.readCollections( event.input )
triggerResults = self.handles['TriggerResults'].product()
for T,TC in self.triggerBitCheckers:
outname="%s_%s"%(self.outprefix,T)
setattr(event,outname, TC.check(event.input.object(), triggerResults))
if self.unrollbits :
for TP,TC in self.triggerBitCheckersSingleBits:
outname="%s_BIT_%s"%(self.outprefix,TP)
setattr(event,outname, TC.check(event.input.object(), triggerResults))


return True


setattr(TriggerBitAnalyzer,"defaultConfig",cfg.Analyzer(
TriggerBitAnalyzer, name="TriggerFlags",
processName = 'HLT',
triggerBits = {
# "<name>" : [ 'HLT_<Something>_v*', 'HLT_<SomethingElse>_v*' ]
}
)
)
setattr(TriggerBitAnalyzer,"defaultEventFlagsConfig",cfg.Analyzer(
TriggerBitAnalyzer, name="EventFlags",
processName = 'PAT',
outprefix = 'Flag',
triggerBits = {
"HBHENoiseFilter" : [ "Flag_HBHENoiseFilter" ],
"CSCTightHaloFilter" : [ "Flag_CSCTightHaloFilter" ],
"hcalLaserEventFilter" : [ "Flag_hcalLaserEventFilter" ],
"EcalDeadCellTriggerPrimitiveFilter" : [ "Flag_EcalDeadCellTriggerPrimitiveFilter" ],
"goodVertices" : [ "Flag_goodVertices" ],
"trackingFailureFilter" : [ "Flag_trackingFailureFilter" ],
"eeBadScFilter" : [ "Flag_eeBadScFilter" ],
"ecalLaserCorrFilter" : [ "Flag_ecalLaserCorrFilter" ],
"trkPOGFilters" : [ "Flag_trkPOGFilters" ],
"trkPOG_manystripclus53X" : [ "Flag_trkPOG_manystripclus53X" ],
"trkPOG_toomanystripclus53X" : [ "Flag_trkPOG_toomanystripclus53X" ],
"trkPOG_logErrorTooManyClusters" : [ "Flag_trkPOG_logErrorTooManyClusters" ],
}
)
)
2 changes: 1 addition & 1 deletion PhysicsTools/Heppy/python/analyzers/core/autovars.py
Original file line number Diff line number Diff line change
Expand Up @@ -33,7 +33,7 @@ def __init__(self,name,baseObjectTypes=[],mcOnly=[],variables=[]):
self.variables = variables
def ownVars(self,isMC):
"""Return only my vars, not including the ones from the bases"""
return [ v for v in self.variables if (isMC or not var.mcOnly) ]
return [ v for v in self.variables if (isMC or not v.mcOnly) ]
def allVars(self,isMC):
"""Return all vars, including the base ones. Duplicate bases are not added twice"""
ret = []; names = {}
Expand Down
1 change: 1 addition & 0 deletions PhysicsTools/Heppy/python/analyzers/objects/JetAnalyzer.py
Original file line number Diff line number Diff line change
Expand Up @@ -321,6 +321,7 @@ def smearJets(self, event, jets):
minLepPt = 10,
relaxJetId = False,
doPuId = False, # Not commissioned in 7.0.X
doQG = False,
recalibrateJets = False,
shiftJEC = 0, # set to +1 or -1 to get +/-1 sigma shifts
smearJets = True,
Expand Down
200 changes: 195 additions & 5 deletions PhysicsTools/Heppy/python/physicsobjects/Electron.py
Original file line number Diff line number Diff line change
@@ -1,6 +1,196 @@
from classes.Electron import Electron
#comment for old
from PhysicsTools.Heppy.physicsobjects.Lepton import Lepton
from PhysicsTools.Heppy.physicsutils.ElectronMVAID import ElectronMVAID_Trig, ElectronMVAID_NonTrig, ElectronMVAID_TrigNoIP
import ROOT
import FastObjects
FastObjects.decorate(ROOT.pat.Electron,Electron)
Electron=FastObjects.AddPhysObjAndCallInit

class Electron( Lepton ):

def __init__(self, *args, **kwargs):
'''Initializing tightIdResult to None. The user is responsible
for setting this attribute externally if he wants to use the tightId
function.'''
super(Electron, self).__init__(*args, **kwargs)
self._physObjInit()

def _physObjInit(self):
self.tightIdResult = None
self.associatedVertex = None
self.rho = None
self._mvaNonTrigV0 = {True:None, False:None}
self._mvaTrigV0 = {True:None, False:None}
self._mvaTrigNoIPV0 = {True:None, False:None}

def electronID( self, id, vertex=None, rho=None ):
if id is None or id == "": return True
if vertex == None and hasattr(self,'associatedVertex') and self.associatedVertex != None: vertex = self.associatedVertex
if rho == None and hasattr(self,'rho') and self.rho != None: rho = self.rho
if id == "POG_MVA_ID_NonTrig": return self.mvaIDLoose()
elif id == "POG_MVA_ID_Trig": return self.mvaIDTight()
elif id == "POG_MVA_ID_NonTrig_full5x5": return self.mvaIDLoose(full5x5=True)
elif id == "POG_MVA_ID_Trig_full5x5": return self.mvaIDTight(full5x5=True)
elif id.startswith("POG_Cuts_ID_"):
return self.cutBasedId(id.replace("POG_Cuts_ID_","POG_"))
for ID in self.electronIDs():
if ID.first == id:
return ID.second
raise RuntimeError, "Electron id '%s' not yet implemented in Electron.py" % id

def cutBasedId(self, wp, showerShapes="auto"):
if "_full5x5" in wp:
showerShapes = "full5x5"
wp = wp.replace("_full5x5","")
elif showerShapes == "auto":
if "POG_CSA14_25ns_v1" in wp or "POG_CSA14_50ns_v1" in wp:
showerShapes = "full5x5"
vars = {
'dEtaIn' : abs(self.physObj.deltaEtaSuperClusterTrackAtVtx()),
'dPhiIn' : abs(self.physObj.deltaPhiSuperClusterTrackAtVtx()),
'sigmaIEtaIEta' : self.physObj.full5x5_sigmaIetaIeta() if showerShapes == "full5x5" else self.physObj.sigmaIetaIeta(),
'H/E' : self.physObj.hadronicOverEm(),
#'1/E-1/p' : abs(1.0/self.physObj.ecalEnergy() - self.physObj.eSuperClusterOverP()/self.physObj.ecalEnergy()),
'1/E-1/p' : abs(1.0/self.physObj.ecalEnergy() - self.physObj.eSuperClusterOverP()/self.physObj.ecalEnergy()) if self.physObj.ecalEnergy()>0. else 9e9,
}
WP = {
## ------- https://twiki.cern.ch/twiki/bin/viewauth/CMS/EgammaCutBasedIdentification?rev=31
'POG_2012_Veto' : [('dEtaIn', [0.007, 0.01]), ('dPhiIn', [0.8, 0.7 ]), ('sigmaIEtaIEta', [0.01, 0.03]), ('H/E', [0.15, 9e9]), ('1/E-1/p', [9e9, 9e9])],
'POG_2012_Loose' : [('dEtaIn', [0.007, 0.009]), ('dPhiIn', [0.15, 0.1 ]), ('sigmaIEtaIEta', [0.01, 0.03]), ('H/E', [0.12, 0.1]), ('1/E-1/p', [0.05, 0.05])],
'POG_2012_Medium' : [('dEtaIn', [0.004, 0.007]), ('dPhiIn', [0.06, 0.03]), ('sigmaIEtaIEta', [0.01, 0.03]), ('H/E', [0.12, 0.1]), ('1/E-1/p', [0.05, 0.05])],
'POG_2012_Tight' : [('dEtaIn', [0.004, 0.005]), ('dPhiIn', [0.03, 0.02]), ('sigmaIEtaIEta', [0.01, 0.03]), ('H/E', [0.12, 0.1]), ('1/E-1/p', [0.05, 0.05])],
## ------- https://indico.cern.cH/Event/298242/contribution/1/material/slides/5.pdf (slide 5)
'POG_CSA14_25ns_v1_Veto' : [('dEtaIn', [0.012, 0.015]), ('dPhiIn', [0.8, 0.7 ]), ('sigmaIEtaIEta', [0.01 , 0.033]), ('H/E', [0.15, 9e9 ]), ('1/E-1/p', [9e9, 9e9])],
'POG_CSA14_25ns_v1_Loose' : [('dEtaIn', [0.012, 0.014]), ('dPhiIn', [0.15, 0.1 ]), ('sigmaIEtaIEta', [0.01 , 0.033]), ('H/E', [0.12, 0.12]), ('1/E-1/p', [0.05, 0.05])],
'POG_CSA14_25ns_v1_Medium' : [('dEtaIn', [0.009, 0.012]), ('dPhiIn', [0.06, 0.03]), ('sigmaIEtaIEta', [0.01 , 0.031]), ('H/E', [0.12, 0.12]), ('1/E-1/p', [0.05, 0.05])],
'POG_CSA14_25ns_v1_Tight' : [('dEtaIn', [0.009, 0.010]), ('dPhiIn', [0.03, 0.02]), ('sigmaIEtaIEta', [0.01 , 0.031]), ('H/E', [0.12, 0.12]), ('1/E-1/p', [0.05, 0.05])],
'POG_CSA14_50ns_v1_Veto' : [('dEtaIn', [0.012, 0.022]), ('dPhiIn', [0.8, 0.7 ]), ('sigmaIEtaIEta', [0.012, 0.033]), ('H/E', [0.15, 9e9 ]), ('1/E-1/p', [9e9, 9e9])],
'POG_CSA14_50ns_v1_Loose' : [('dEtaIn', [0.012, 0.021]), ('dPhiIn', [0.15, 0.1 ]), ('sigmaIEtaIEta', [0.012, 0.033]), ('H/E', [0.12, 0.12]), ('1/E-1/p', [0.05, 0.05])],
'POG_CSA14_50ns_v1_Medium' : [('dEtaIn', [0.009, 0.019]), ('dPhiIn', [0.06, 0.03]), ('sigmaIEtaIEta', [0.01 , 0.031]), ('H/E', [0.12, 0.12]), ('1/E-1/p', [0.05, 0.05])],
'POG_CSA14_50ns_v1_Tight' : [('dEtaIn', [0.009, 0.017]), ('dPhiIn', [0.03, 0.02]), ('sigmaIEtaIEta', [0.01 , 0.031]), ('H/E', [0.12, 0.12]), ('1/E-1/p', [0.05, 0.05])],
}
if wp not in WP:
raise RuntimeError, "Working point '%s' not yet implemented in Electron.py" % wp
for (cut_name,(cut_eb,cut_ee)) in WP[wp]:
if vars[cut_name] >= (cut_eb if self.physObj.isEB() else cut_ee):
return False
return True

def absEffAreaIso(self,rho,effectiveAreas):
'''MIKE, missing doc.
Should have the same name as the function in the mother class.
Can call the mother class function with super.
'''
return self.absIsoFromEA(rho,self.superCluster().eta(),effectiveAreas.eGamma)

def mvaId( self ):
return self.mvaNonTrigV0()

def tightId( self ):
return self.tightIdResult

def mvaNonTrigV0( self, full5x5=False, debug = False ):
if self._mvaNonTrigV0[full5x5] == None:
if self.associatedVertex == None: raise RuntimeError, "You need to set electron.associatedVertex before calling any MVA"
if self.rho == None: raise RuntimeError, "You need to set electron.rho before calling any MVA"
self._mvaNonTrigV0[full5x5] = ElectronMVAID_NonTrig(self.physObj, self.associatedVertex, self.rho, full5x5, debug)
return self._mvaNonTrigV0[full5x5]

def mvaTrigV0( self, full5x5=False, debug = False ):
if self._mvaTrigV0[full5x5] == None:
if self.associatedVertex == None: raise RuntimeError, "You need to set electron.associatedVertex before calling any MVA"
if self.rho == None: raise RuntimeError, "You need to set electron.rho before calling any MVA"
self._mvaTrigV0[full5x5] = ElectronMVAID_Trig(self.physObj, self.associatedVertex, self.rho, full5x5, debug)
return self._mvaTrigV0[full5x5]

def mvaTrigNoIPV0( self, full5x5=False, debug = False ):
if self._mvaTrigNoIPV0[full5x5] == None:
if self.associatedVertex == None: raise RuntimeError, "You need to set electron.associatedVertex before calling any MVA"
if self.rho == None: raise RuntimeError, "You need to set electron.rho before calling any MVA"
self._mvaTrigNoIPV0[full5x5] = ElectronMVAID_TrigNoIP(self.physObj, self.associatedVertex, self.rho, full5x5, debug)
return self._mvaTrigNoIPV0[full5x5]


def mvaIDTight(self, full5x5=False):
eta = abs(self.superCluster().eta())
if self.pt() < 20:
if (eta < 0.8) : return self.mvaTrigV0(full5x5) > +0.00;
elif (eta < 1.479): return self.mvaTrigV0(full5x5) > +0.10;
else : return self.mvaTrigV0(full5x5) > +0.62;
else:
if (eta < 0.8) : return self.mvaTrigV0(full5x5) > +0.94;
elif (eta < 1.479): return self.mvaTrigV0(full5x5) > +0.85;
else : return self.mvaTrigV0(full5x5) > +0.92;

def mvaIDLoose(self, full5x5=False):
eta = abs(self.superCluster().eta())
if self.pt() < 10:
if (eta < 0.8) : return self.mvaNonTrigV0(full5x5) > +0.47;
elif (eta < 1.479): return self.mvaNonTrigV0(full5x5) > +0.004;
else : return self.mvaNonTrigV0(full5x5) > +0.295;
else:
if (eta < 0.8) : return self.mvaNonTrigV0(full5x5) > -0.34;
elif (eta < 1.479): return self.mvaNonTrigV0(full5x5) > -0.65;
else : return self.mvaNonTrigV0(full5x5) > +0.60;

def mvaIDZZ(self):
return self.mvaIDLoose() and (self.gsfTrack().trackerExpectedHitsInner().numberOfLostHits()<=1)

def chargedHadronIsoR(self,R=0.4):
if R == 0.3: return self.physObj.pfIsolationVariables().sumChargedHadronPt
elif R == 0.4: return self.physObj.chargedHadronIso()
raise RuntimeError, "Electron chargedHadronIso missing for R=%s" % R

def neutralHadronIsoR(self,R=0.4):
if R == 0.3: return self.physObj.pfIsolationVariables().sumNeutralHadronEt
elif R == 0.4: return self.physObj.neutralHadronIso()
raise RuntimeError, "Electron neutralHadronIso missing for R=%s" % R

def photonIsoR(self,R=0.4):
if R == 0.3: return self.physObj.pfIsolationVariables().sumPhotonEt
elif R == 0.4: return self.physObj.photonIso()
raise RuntimeError, "Electron photonIso missing for R=%s" % R

def chargedAllIso(self,R=0.4):
if R == 0.3: return self.physObj.pfIsolationVariables().sumChargedParticlePt
raise RuntimeError, "Electron chargedAllIso missing for R=%s" % R

def puChargedHadronIsoR(self,R=0.4):
if R == 0.3: return self.physObj.pfIsolationVariables().sumPUPt
elif R == 0.4: return self.physObj.puChargedHadronIso()
raise RuntimeError, "Electron chargedHadronIso missing for R=%s" % R




def chargedAllIso(self):
'''This function is used in the isolation, see Lepton class.
Here, we replace the all charged isolation by the all charged isolation with cone veto'''
return self.chargedAllIsoWithConeVeto()


def dxy(self, vertex=None):
'''Returns dxy.
Computed using vertex (or self.associatedVertex if vertex not specified),
and the gsf track.
'''
if vertex is None:
vertex = self.associatedVertex
return self.gsfTrack().dxy( vertex.position() )
def p4(self):
return ROOT.reco.Candidate.p4(self.physObj)

# def p4(self):
# return self.physObj.p4(self.physObj.candidateP4Kind()) # if kind == None else kind)

def dz(self, vertex=None):
'''Returns dz.
Computed using vertex (or self.associatedVertex if vertex not specified),
and the gsf track.
'''
if vertex is None:
vertex = self.associatedVertex
return self.gsfTrack().dz( vertex.position() )

def lostInner(self) :
if hasattr(self.gsfTrack(),"trackerExpectedHitsInner") :
return self.gsfTrack().trackerExpectedHitsInner().numberOfLostHits()
else :
return self.gsfTrack().hitPattern().numberOfHits(ROOT.reco.HitPattern.MISSING_INNER_HITS)

Loading

0 comments on commit 2440eb3

Please sign in to comment.