From 665c478f61865185743745da1528fc1c0a635499 Mon Sep 17 00:00:00 2001 From: Andrea Date: Wed, 17 Dec 2014 11:26:16 +0100 Subject: [PATCH 01/10] Heppy: use tfileservice in NumpyTreeProducer when available --- .../analyzers/core/AutoFillTreeProducer.py | 7 +++++-- .../analyzers/core/TreeAnalyzerNumpy.py | 20 ++++++++++++------- 2 files changed, 18 insertions(+), 9 deletions(-) diff --git a/PhysicsTools/Heppy/python/analyzers/core/AutoFillTreeProducer.py b/PhysicsTools/Heppy/python/analyzers/core/AutoFillTreeProducer.py index a1d256ac0b776..46d5e75e841a9 100644 --- a/PhysicsTools/Heppy/python/analyzers/core/AutoFillTreeProducer.py +++ b/PhysicsTools/Heppy/python/analyzers/core/AutoFillTreeProducer.py @@ -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 = {} @@ -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' ) diff --git a/PhysicsTools/Heppy/python/analyzers/core/TreeAnalyzerNumpy.py b/PhysicsTools/Heppy/python/analyzers/core/TreeAnalyzerNumpy.py index d8f4c9580710a..427c0f735e5c1 100644 --- a/PhysicsTools/Heppy/python/analyzers/core/TreeAnalyzerNumpy.py +++ b/PhysicsTools/Heppy/python/analyzers/core/TreeAnalyzerNumpy.py @@ -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): @@ -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() From 14785f6018964f88f1dabbaa2daf255c7b436f23 Mon Sep 17 00:00:00 2001 From: Andrea Date: Wed, 17 Dec 2014 12:25:55 +0100 Subject: [PATCH 02/10] add missing parameter in default config --- PhysicsTools/Heppy/python/analyzers/objects/JetAnalyzer.py | 1 + 1 file changed, 1 insertion(+) diff --git a/PhysicsTools/Heppy/python/analyzers/objects/JetAnalyzer.py b/PhysicsTools/Heppy/python/analyzers/objects/JetAnalyzer.py index 31cfbe77d65b7..55309728c36d7 100644 --- a/PhysicsTools/Heppy/python/analyzers/objects/JetAnalyzer.py +++ b/PhysicsTools/Heppy/python/analyzers/objects/JetAnalyzer.py @@ -317,6 +317,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, From 0254758e007bdb798af8f8442fe695b7781993d9 Mon Sep 17 00:00:00 2001 From: Giovanni Date: Thu, 18 Dec 2014 17:16:00 +0100 Subject: [PATCH 03/10] Heppy: fix application of jet energy correction, and allow specifying to recalibrate only data or only MC --- .../python/analyzers/objects/JetAnalyzer.py | 8 ++++-- .../python/physicsutils/JetReCalibrator.py | 28 +++++++++++-------- 2 files changed, 23 insertions(+), 13 deletions(-) diff --git a/PhysicsTools/Heppy/python/analyzers/objects/JetAnalyzer.py b/PhysicsTools/Heppy/python/analyzers/objects/JetAnalyzer.py index 31cfbe77d65b7..ce0d542b59f42 100644 --- a/PhysicsTools/Heppy/python/analyzers/objects/JetAnalyzer.py +++ b/PhysicsTools/Heppy/python/analyzers/objects/JetAnalyzer.py @@ -25,10 +25,14 @@ class JetAnalyzer( Analyzer ): """Taken from RootTools.JetAnalyzer, simplified, modified, added corrections """ def __init__(self, cfg_ana, cfg_comp, looperName): super(JetAnalyzer,self).__init__(cfg_ana, cfg_comp, looperName) - mcGT = cfg_ana.mcGT if hasattr(cfg_ana,'mcGT') else "PHYS14_25_V1" + mcGT = cfg_ana.mcGT if hasattr(cfg_ana,'mcGT') else "PHYS14_25_V2" dataGT = cfg_ana.dataGT if hasattr(cfg_ana,'dataGT') else "GR_70_V2_AN1" self.shiftJEC = self.cfg_ana.shiftJEC if hasattr(self.cfg_ana, 'shiftJEC') else 0 - self.doJEC = self.cfg_ana.recalibrateJets or (self.shiftJEC != 0) + self.recalibrateJets = self.cfg_ana.recalibrateJets + if self.recalibrateJets == "MC" : self.recalibrateJets = self.cfg_comp.isMC + elif self.recalibrateJets == "Data": self.recalibrateJets = not self.cfg_comp.isMC + elif self.recalibrateJets not in [True,False]: raise RuntimeError, "recalibrateJets must be any of { True, False, 'MC', 'Data' }, while it is %r " % self.recalibrateJets + self.doJEC = self.recalibrateJets or (self.shiftJEC != 0) if self.doJEC: if self.cfg_comp.isMC: self.jetReCalibrator = JetReCalibrator(mcGT,"AK4PFchs", False,cfg_ana.jecPath) diff --git a/PhysicsTools/Heppy/python/physicsutils/JetReCalibrator.py b/PhysicsTools/Heppy/python/physicsutils/JetReCalibrator.py index f0ed8269d782b..7f7893b15bd1d 100644 --- a/PhysicsTools/Heppy/python/physicsutils/JetReCalibrator.py +++ b/PhysicsTools/Heppy/python/physicsutils/JetReCalibrator.py @@ -23,7 +23,11 @@ def __init__(self,globalTag,jetFlavour,doResidualJECs,jecPath,upToLevel=3): self.vPar.push_back(self.ResJetPar); #Step3 (Construct a FactorizedJetCorrector object) self.JetCorrector = ROOT.FactorizedJetCorrector(self.vPar) - self.JetUncertainty = ROOT.JetCorrectionUncertainty("%s/%s_Uncertainty_%s.txt" % (path,globalTag,jetFlavour)); + if os.path.exists("%s/%s_Uncertainty_%s.txt" % (path,globalTag,jetFlavour)): + self.JetUncertainty = ROOT.JetCorrectionUncertainty("%s/%s_Uncertainty_%s.txt" % (path,globalTag,jetFlavour)); + else: + print 'Missing JEC uncertainty file "%s/%s_Uncertainty_%s.txt", so jet energy uncertainties will not be available' % (path,globalTag,jetFlavour) + self.JetUncertainty = None def correctAll(self,jets,rho,delta=0,metShift=[0,0]): """Applies 'correct' to all the jets, discard the ones that have bad corrections (corrected pt <= 0)""" badJets = [] @@ -44,16 +48,18 @@ def correct(self,jet,rho,delta=0,metShift=[0,0]): self.JetCorrector.setJetA(jet.jetArea()) self.JetCorrector.setRho(rho) corr = self.JetCorrector.getCorrection() - self.JetUncertainty.setJetEta(jet.eta()) - self.JetUncertainty.setJetPt(corr * jet.pt() * jet.rawFactor()) - try: - jet.jetEnergyCorrUncertainty = self.JetUncertainty.getUncertainty(True) - except RuntimeError, r: - print "Caught %s when getting uncertainty for jet of pt %.1f, eta %.2f\n" % (r,corr * jet.pt() * jet.rawFactor(),jet.eta()) - jet.jetEnergyCorrUncertainty = 0.5 - if jet.component(4).fraction() < 0.9 and jet.pt()*corr*jet.rawFactor() > 10: - metShift[0] -= jet.px()*(corr*jet.rawFactor() - 1)*(1-jet.component(3).fraction()) - metShift[1] -= jet.py()*(corr*jet.rawFactor() - 1)*(1-jet.component(3).fraction()) + if delta != 0: + if not self.JetUncertainty: raise RuntimeError, "Jet energy scale uncertainty shifts requested, but not available" + self.JetUncertainty.setJetEta(jet.eta()) + self.JetUncertainty.setJetPt(corr * jet.pt() * jet.rawFactor()) + try: + jet.jetEnergyCorrUncertainty = self.JetUncertainty.getUncertainty(True) + except RuntimeError, r: + print "Caught %s when getting uncertainty for jet of pt %.1f, eta %.2f\n" % (r,corr * jet.pt() * jet.rawFactor(),jet.eta()) + jet.jetEnergyCorrUncertainty = 0.5 + if jet.photonEnergyFraction() < 0.9 and jet.pt()*corr*jet.rawFactor() > 10: + metShift[0] -= jet.px()*(corr*jet.rawFactor() - 1)*(1-jet.muonEnergyFraction()) + metShift[1] -= jet.py()*(corr*jet.rawFactor() - 1)*(1-jet.muonEnergyFraction()) if delta != 0: #print " jet with corr pt %6.2f has an uncertainty %.2f " % (jet.pt()*jet.rawFactor()*corr, jet.jetEnergyCorrUncertainty) corr *= max(0, 1+delta*jet.jetEnergyCorrUncertainty) From f497042e77677c445d1dfc7632be5112378ca48a Mon Sep 17 00:00:00 2001 From: Giovanni Date: Thu, 18 Dec 2014 17:16:34 +0100 Subject: [PATCH 04/10] Heppy: Working points for the combinedMVABJetTags tagger (preliminary) --- PhysicsTools/Heppy/python/physicsobjects/classes/Jet.py | 3 +++ 1 file changed, 3 insertions(+) diff --git a/PhysicsTools/Heppy/python/physicsobjects/classes/Jet.py b/PhysicsTools/Heppy/python/physicsobjects/classes/Jet.py index 6aea5d8866883..22e6d89a74082 100644 --- a/PhysicsTools/Heppy/python/physicsobjects/classes/Jet.py +++ b/PhysicsTools/Heppy/python/physicsobjects/classes/Jet.py @@ -28,6 +28,9 @@ "CSVv2IVFL": ("combinedInclusiveSecondaryVertexV2BJetTags", 0.423), "CSVv2IVFM": ("combinedInclusiveSecondaryVertexV2BJetTags", 0.814), "CSVv2IVFT": ("combinedInclusiveSecondaryVertexV2BJetTags", 0.941), + "CMVAL": ("combinedMVABJetTags", 0.630), # for same b-jet efficiency of CSVv2IVFL on ttbar MC, jet pt > 30 + "CMVAM": ("combinedMVABJetTags", 0.732), # for same b-jet efficiency of CSVv2IVFM on ttbar MC, jet pt > 30 + "CMVAT": ("combinedMVABJetTags", 0.813), # for same b-jet efficiency of CSVv2IVFT on ttbar MC, jet pt > 30 } class Jet(PhysicsObject): From a9c42407470fef5935add13064113cb882734639 Mon Sep 17 00:00:00 2001 From: Giovanni Date: Thu, 18 Dec 2014 17:48:16 +0100 Subject: [PATCH 05/10] HeppyCore: allow splitting a file in multiple chunks This feature can be enabled by setting fineSplitFactor = N > 1 on a component. In that case, the scheduler will make N jobs for each file, and each job will process 1/N of the events in the file. To avoid code duplication, now heppy_batch takes the split function from heppy. Also, the print out of the time per event has been updated in order to work correctly with this feature, or in general when not processing events starting from the first one in the file. --- .../HeppyCore/python/framework/heppy.py | 11 +++++++++- .../HeppyCore/python/framework/looper.py | 17 +++++++++++--- PhysicsTools/HeppyCore/scripts/heppy_batch.py | 22 +------------------ 3 files changed, 25 insertions(+), 25 deletions(-) diff --git a/PhysicsTools/HeppyCore/python/framework/heppy.py b/PhysicsTools/HeppyCore/python/framework/heppy.py index d93c1e177edf5..c62cb44ef9659 100755 --- a/PhysicsTools/HeppyCore/python/framework/heppy.py +++ b/PhysicsTools/HeppyCore/python/framework/heppy.py @@ -95,7 +95,16 @@ def split(comps): # import pdb; pdb.set_trace() splitComps = [] for comp in comps: - if hasattr( comp, 'splitFactor') and comp.splitFactor>1: + if hasattr( comp, 'fineSplitFactor') and comp.fineSplitFactor>1: + subchunks = range(comp.fineSplitFactor) + for ichunk, chunk in enumerate([(f,i) for f in comp.files for i in subchunks]): + newComp = copy.deepcopy(comp) + newComp.files = [chunk[0]] + newComp.fineSplit = ( chunk[1], comp.fineSplitFactor ) + newComp.name = '{name}_Chunk{index}'.format(name=newComp.name, + index=ichunk) + splitComps.append( newComp ) + elif hasattr( comp, 'splitFactor') and comp.splitFactor>1: chunkSize = len(comp.files) / comp.splitFactor if len(comp.files) % comp.splitFactor: chunkSize += 1 diff --git a/PhysicsTools/HeppyCore/python/framework/looper.py b/PhysicsTools/HeppyCore/python/framework/looper.py index ac72609e1fb72..31802988e229d 100644 --- a/PhysicsTools/HeppyCore/python/framework/looper.py +++ b/PhysicsTools/HeppyCore/python/framework/looper.py @@ -8,6 +8,7 @@ import logging import pprint from platform import platform +from math import ceil from event import Event import time @@ -85,6 +86,15 @@ def __init__( self, name, errmsg = 'please provide at least an input file in the files attribute of this component\n' + str(self.cfg_comp) raise ValueError( errmsg ) self.events = config.events_class(self.cfg_comp.files, tree_name) + if hasattr(self.cfg_comp, 'fineSplit'): + fineSplitIndex, fineSplitFactor = self.cfg_comp.fineSplit + if fineSplitFactor > 1: + if len(self.cfg_comp.files) != 1: + raise RuntimeError, "Any component with fineSplit > 1 is supposed to have just a single file, while %s has %s" % (self.cfg_comp.name, self.cfg_comp.files) + totevents = min(len(self.events),int(nEvents)) if (nEvents and int(nEvents) not in [-1,0]) else len(self.events) + self.nEvents = int(ceil(totevents/float(fineSplitFactor))) + self.firstEvent = firstEvent + fineSplitIndex * self.nEvents + #print "For component %s will process %d events starting from the %d one" % (self.cfg_comp.name, self.nEvents, self.firstEvent) # self.event is set in self.process self.event = None services = dict() @@ -146,11 +156,12 @@ def loop(self): # break if iEv%100 ==0: # print 'event', iEv - if iEv == 100: + if not hasattr(self,'start_time'): print 'event', iEv self.start_time = time.time() - elif iEv > 100: - print 'event %d (%.1f ev/s)' % (iEv, (iEv-100)/float(time.time() - self.start_time)) + self.start_time_event = iEv + else: + print 'event %d (%.1f ev/s)' % (iEv, (iEv-self.start_time_event)/float(time.time() - self.start_time)) self.process( iEv ) if iEv1: - chunkSize = len(comp.files) / comp.splitFactor - if len(comp.files) % comp.splitFactor: - chunkSize += 1 - # print 'chunk size',chunkSize, len(comp.files), comp.splitFactor - for ichunk, chunk in enumerate( chunks( comp.files, chunkSize)): - newComp = copy.deepcopy(comp) - newComp.files = chunk - newComp.name = '{name}_Chunk{index}'.format(name=newComp.name, - index=ichunk) - splitComps.append( newComp ) - else: - splitComps.append( comp ) - return splitComps +from PhysicsTools.HeppyCore.framework.heppy import split def batchScriptPISA( index, remoteDir=''): From 065d46ebcec9b809777a9b642ceb178bd8f66b64 Mon Sep 17 00:00:00 2001 From: Andrea Date: Fri, 19 Dec 2014 11:35:01 +0100 Subject: [PATCH 06/10] remove fast objects --- .../Heppy/python/physicsobjects/Electron.py | 200 +++++++++++++++++- .../Heppy/python/physicsobjects/IsoTrack.py | 30 ++- .../Heppy/python/physicsobjects/Jet.py | 112 +++++++++- .../Heppy/python/physicsobjects/Muon.py | 95 ++++++++- .../Heppy/python/physicsobjects/Photon.py | 40 +++- .../Heppy/python/physicsobjects/Tau.py | 86 +++++++- .../python/physicsobjects/classes/Electron.py | 196 ----------------- .../python/physicsobjects/classes/IsoTrack.py | 29 --- .../python/physicsobjects/classes/Jet.py | 106 ---------- .../python/physicsobjects/classes/Muon.py | 89 -------- .../python/physicsobjects/classes/Photon.py | 39 ---- .../python/physicsobjects/classes/Tau.py | 85 -------- 12 files changed, 543 insertions(+), 564 deletions(-) delete mode 100644 PhysicsTools/Heppy/python/physicsobjects/classes/Electron.py delete mode 100644 PhysicsTools/Heppy/python/physicsobjects/classes/IsoTrack.py delete mode 100644 PhysicsTools/Heppy/python/physicsobjects/classes/Jet.py delete mode 100644 PhysicsTools/Heppy/python/physicsobjects/classes/Muon.py delete mode 100644 PhysicsTools/Heppy/python/physicsobjects/classes/Photon.py delete mode 100644 PhysicsTools/Heppy/python/physicsobjects/classes/Tau.py diff --git a/PhysicsTools/Heppy/python/physicsobjects/Electron.py b/PhysicsTools/Heppy/python/physicsobjects/Electron.py index e2594e3ad6785..0132385f1046d 100644 --- a/PhysicsTools/Heppy/python/physicsobjects/Electron.py +++ b/PhysicsTools/Heppy/python/physicsobjects/Electron.py @@ -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) + diff --git a/PhysicsTools/Heppy/python/physicsobjects/IsoTrack.py b/PhysicsTools/Heppy/python/physicsobjects/IsoTrack.py index 8b7d82b5122d6..e062f85349e73 100644 --- a/PhysicsTools/Heppy/python/physicsobjects/IsoTrack.py +++ b/PhysicsTools/Heppy/python/physicsobjects/IsoTrack.py @@ -1 +1,29 @@ -from classes.IsoTrack import IsoTrack +from PhysicsTools.Heppy.physicsobjects.PhysicsObject import * +import ROOT +import math + +class IsoTrack( PhysicsObject ): + + def __init__(self, isoTrack): + self.isoTrack = isoTrack + super(IsoTrack, self).__init__(isoTrack) + + def absIso(self, dummy1, dummy2): + '''Just making the tau behave as a lepton.''' + return -1 + + def relIso(self, dummy1, dummy2): + '''Just making the tau behave as a lepton.''' + return -1 + + def __str__(self): + lep = super(IsoTrack, self).__str__() + return lep + #spec = '\t\tTau: decay = {decMode:<15}, eOverP = {eOverP:4.2f}, isoMVALoose = {isoMVALoose}'.format( + # decMode = tauDecayModes.intToName( self.decayMode() ), + # eOverP = self.calcEOverP(), + # isoMVALoose = self.tauID('byLooseIsoMVA') + # ) + #return '\n'.join([lep, spec]) + + diff --git a/PhysicsTools/Heppy/python/physicsobjects/Jet.py b/PhysicsTools/Heppy/python/physicsobjects/Jet.py index e3366f65fa101..910b97e96d1ea 100644 --- a/PhysicsTools/Heppy/python/physicsobjects/Jet.py +++ b/PhysicsTools/Heppy/python/physicsobjects/Jet.py @@ -1,7 +1,107 @@ -from classes.Jet import Jet,GenJet +from PhysicsTools.Heppy.physicsobjects.PhysicsObject import * + +loose_WP = [ + (0, 2.5, -0.8), + (2.5, 2.75, -0.74), + (2.75, 3.0, -0.68), + (3.0, 5.0, -0.77), + ] + +# Working point 2 May 2013 (Phil via H2tau list) +loose_53X_WP = [ + (0, 2.5, -0.63), + (2.5, 2.75, -0.60), + (2.75, 3.0, -0.55), + (3.0, 5.2, -0.45), + ] + +_btagWPs = { + "TCHEL": ("trackCountingHighEffBJetTags", 1.7), + "TCHEM": ("trackCountingHighEffBJetTags", 3.3), + "TCHPT": ("trackCountingHighPurBJetTags", 3.41), + "JPL": ("jetProbabilityBJetTags", 0.275), + "JPM": ("jetProbabilityBJetTags", 0.545), + "JPT": ("jetProbabilityBJetTags", 0.790), + "CSVL": ("combinedSecondaryVertexBJetTags", 0.244), + "CSVM": ("combinedSecondaryVertexBJetTags", 0.679), + "CSVT": ("combinedSecondaryVertexBJetTags", 0.898), + "CSVv2IVFL": ("combinedInclusiveSecondaryVertexV2BJetTags", 0.423), + "CSVv2IVFM": ("combinedInclusiveSecondaryVertexV2BJetTags", 0.814), + "CSVv2IVFT": ("combinedInclusiveSecondaryVertexV2BJetTags", 0.941), + "CMVAL": ("combinedMVABJetTags", 0.630), # for same b-jet efficiency of CSVv2IVFL on ttbar MC, jet pt > 30 + "CMVAM": ("combinedMVABJetTags", 0.732), # for same b-jet efficiency of CSVv2IVFM on ttbar MC, jet pt > 30 + "CMVAT": ("combinedMVABJetTags", 0.813), # for same b-jet efficiency of CSVv2IVFT on ttbar MC, jet pt > 30 + +} + +class Jet(PhysicsObject): + def __init__(self, *args, **kwargs): + super(Jet, self).__init__(*args, **kwargs) + self._physObjInit() + + def _physObjInit(self): + self._rawFactorMultiplier = 1.0 + + def jetID(self,name=""): + if not self.isPFJet(): + raise RuntimeError, "jetID implemented only for PF Jets" + eta = abs(self.eta()); + chf = self.chargedHadronEnergyFraction(); + nhf = self.neutralHadronEnergyFraction(); + phf = self.neutralEmEnergyFraction(); + muf = self.muonEnergyFraction(); + elf = self.chargedEmEnergyFraction(); + chm = self.chargedHadronMultiplicity(); + npr = self.chargedMultiplicity() + self.neutralMultiplicity(); + #if npr != self.nConstituents(): + # import pdb; pdb.set_trace() + if name == "POG_PFID": + if self.jetID("POG_PFID_Tight") : return 3; + elif self.jetID("POG_PFID_Medium") : return 2; + elif self.jetID("POG_PFID_Loose") : return 1; + else : return 0; + + if name == "POG_PFID_Loose": return (npr>1 and phf<0.99 and nhf<0.99) and (eta>2.4 or (elf<0.99 and chf>0 and chm>0)); + if name == "POG_PFID_Medium": return (npr>1 and phf<0.95 and nhf<0.95) and (eta>2.4 or (elf<0.99 and chf>0 and chm>0)); + if name == "POG_PFID_Tight": return (npr>1 and phf<0.90 and nhf<0.90) and (eta>2.4 or (elf<0.99 and chf>0 and chm>0)); + if name == "VBFHBB_PFID_Loose": return (npr>1 and phf<0.99 and nhf<0.99); + if name == "VBFHBB_PFID_Medium": return (npr>1 and phf<0.99 and nhf<0.99) and ((eta<=2.4 and nhf<0.9 and phf<0.9 and elf<0.99 and muf<0.99 and chf>0 and chm>0) or eta>2.4); + if name == "VBFHBB_PFID_Tight": return (npr>1 and phf<0.99 and nhf<0.99) and ((eta<=2.4 and nhf<0.9 and phf<0.9 and elf<0.70 and muf<0.70 and chf>0 and chm>0) or eta>2.4); + raise RuntimeError, "jetID '%s' not supported" % name + + def looseJetId(self): + '''PF Jet ID (loose operation point) [method provided for convenience only]''' + return self.jetID("POG_PFID_Loose") + + def puMva(self, label="pileupJetId:fullDiscriminant"): + return self.userFloat(label) + + def puJetId(self, label="pileupJetId:fullDiscriminant"): + '''Full mva PU jet id''' + + puMva = self.puMva(label) + wp = loose_53X_WP + eta = abs(self.eta()) + + for etamin, etamax, cut in wp: + if not(eta>=etamin and etacut + + def rawFactor(self): + return self.jecFactor('Uncorrected') * self._rawFactorMultiplier + def setRawFactor(self, factor): + self._rawFactorMultiplier = factor/self.jecFactor('Uncorrected') + + def btag(self,name): + return self.bDiscriminator(name) + + def btagWP(self,name): + global _btagWPs + (disc,val) = _btagWPs[name] + return self.bDiscriminator(disc) > val + + +class GenJet( PhysicsObject): + pass -#comment for old -import ROOT -import FastObjects -FastObjects.decorate(ROOT.pat.Jet,Jet) -Jet=FastObjects.AddPhysObjAndCallInit diff --git a/PhysicsTools/Heppy/python/physicsobjects/Muon.py b/PhysicsTools/Heppy/python/physicsobjects/Muon.py index 9c4338d45f4d1..a4bed67b83e2d 100644 --- a/PhysicsTools/Heppy/python/physicsobjects/Muon.py +++ b/PhysicsTools/Heppy/python/physicsobjects/Muon.py @@ -1,6 +1,89 @@ -from classes.Muon import Muon -#comment next three lines for old objects -import ROOT -import FastObjects -FastObjects.decorate(ROOT.pat.Muon,Muon) -Muon=FastObjects.AddPhysObj +from PhysicsTools.Heppy.physicsobjects.Lepton import Lepton + +class Muon( Lepton ): + + def looseId( self ): + '''Loose ID as recommended by mu POG.''' + return self.physObj.isLooseMuon() + + def tightId( self ): + '''Tight ID as recommended by mu POG.''' + return self.muonID("POG_ID_Tight") + + def muonID(self, name, vertex=None): + if name == "" or name is None: + return True + if name.startswith("POG_"): + if name == "POG_ID_Loose": return self.physObj.isLooseMuon() + if vertex is None: + vertex = getattr(self, 'associatedVertex', None) + if name == "POG_ID_Tight": return self.physObj.isTightMuon(vertex) + if name == "POG_ID_HighPt": return self.physObj.isHighPtMuon(vertex) + if name == "POG_ID_Soft": return self.physObj.isSoftMuon(vertex) + if name == "POG_ID_TightNoVtx": return self.looseId() and \ + self.isGlobalMuon() and \ + self.globalTrack().normalizedChi2() < 10 and \ + self.globalTrack().hitPattern().numberOfValidMuonHits() > 0 and \ + self.numberOfMatchedStations()>1 and \ + self.innerTrack().hitPattern().numberOfValidPixelHits()>0 and \ + self.innerTrack().hitPattern().trackerLayersWithMeasurement() > 5 + return self.physObj.muonID(name) + + def mvaId(self): + '''For a transparent treatment of electrons and muons. Returns -99''' + return -99 + + + + def absEffAreaIso(self,rho,effectiveAreas): + return self.absIsoFromEA(rho,self.eta(),effectiveAreas.muon) + + + + def dxy(self, vertex=None): + '''either pass the vertex, or set associatedVertex before calling the function. + note: the function does not work with standalone muons as innerTrack + is not available. + ''' + if vertex is None: + vertex = self.associatedVertex + return self.innerTrack().dxy( vertex.position() ) + + + def dz(self, vertex=None): + '''either pass the vertex, or set associatedVertex before calling the function. + note: the function does not work with standalone muons as innerTrack + is not available. + ''' + if vertex is None: + vertex = self.associatedVertex + return self.innerTrack().dz( vertex.position() ) + + def chargedHadronIsoR(self,R=0.4): + if R == 0.3: return self.physObj.pfIsolationR03().sumChargedHadronPt + elif R == 0.4: return self.physObj.pfIsolationR04().sumChargedHadronPt + raise RuntimeError, "Muon chargedHadronIso missing for R=%s" % R + + def neutralHadronIsoR(self,R=0.4): + if R == 0.3: return self.physObj.pfIsolationR03().sumNeutralHadronEt + elif R == 0.4: return self.physObj.pfIsolationR04().sumNeutralHadronEt + raise RuntimeError, "Muon neutralHadronIso missing for R=%s" % R + + def photonIsoR(self,R=0.4): + if R == 0.3: return self.physObj.pfIsolationR03().sumPhotonEt + elif R == 0.4: return self.physObj.pfIsolationR04().sumPhotonEt + raise RuntimeError, "Muon photonIso missing for R=%s" % R + + def chargedAllIsoR(self,R=0.4): + if R == 0.3: return self.physObj.pfIsolationR03().sumChargedParticlePt + elif R == 0.4: return self.physObj.pfIsolationR04().sumChargedParticlePt + raise RuntimeError, "Muon chargedAllIso missing for R=%s" % R + + def puChargedHadronIsoR(self,R=0.4): + if R == 0.3: return self.physObj.pfIsolationR03().sumPUPt + elif R == 0.4: return self.physObj.pfIsolationR04().sumPUPt + raise RuntimeError, "Muon chargedHadronIso missing for R=%s" % R + + + + diff --git a/PhysicsTools/Heppy/python/physicsobjects/Photon.py b/PhysicsTools/Heppy/python/physicsobjects/Photon.py index 8013b9e0a7ee8..371b25a057375 100644 --- a/PhysicsTools/Heppy/python/physicsobjects/Photon.py +++ b/PhysicsTools/Heppy/python/physicsobjects/Photon.py @@ -1 +1,39 @@ -from classes.Photon import Photon +from PhysicsTools.Heppy.physicsobjects.PhysicsObject import * + +class Photon(PhysicsObject ): + + ''' return object from the photon + ''' + def hOVERe(self): +# return self.physObj.full5x5_hadTowOverEm() + hadTowDepth1O = self.physObj.hadTowDepth1OverEm() * (self.physObj.superCluster().energy()/self.physObj.full5x5_e5x5() if self.physObj.full5x5_e5x5() else 1) + hadTowDepth2O = self.physObj.hadTowDepth2OverEm() * (self.physObj.superCluster().energy()/self.physObj.full5x5_e5x5() if self.physObj.full5x5_e5x5() else 1) + return hadTowDepth1O + hadTowDepth2O + + def r9(self): + return self.physObj.full5x5_r9() + + def sigmaIetaIeta(self): + return self.physObj.full5x5_sigmaIetaIeta() + + def chargedHadronIso(self): + return self.physObj.chargedHadronIso() + + def neutralHadronIso(self): + return self.physObj.neutralHadronIso() + + def photonIso(self): + return self.physObj.photonIso() + + def photonIDCSA14(self, name): + keepThisPhoton = True + if name == "PhotonCutBasedIDLoose_CSA14": + if abs(self.physObj.eta())<1.479 : + if self.physObj.sigmaIetaIeta() > 0.012 : keepThisPhoton = False + if self.hOVERe() > 0.0559 : keepThisPhoton = False + else : + if self.physObj.sigmaIetaIeta() > 0.035 : keepThisPhoton = False + if self.hOVERe() > 0.049 : keepThisPhoton = False + return keepThisPhoton + + pass diff --git a/PhysicsTools/Heppy/python/physicsobjects/Tau.py b/PhysicsTools/Heppy/python/physicsobjects/Tau.py index 2c3578a653840..e20506148b06c 100644 --- a/PhysicsTools/Heppy/python/physicsobjects/Tau.py +++ b/PhysicsTools/Heppy/python/physicsobjects/Tau.py @@ -1 +1,85 @@ -from classes.Tau import Tau,isTau +from PhysicsTools.Heppy.physicsobjects.Lepton import Lepton +from PhysicsTools.Heppy.physicsutils.TauDecayModes import tauDecayModes +import math + +cutsElectronMVA3Medium = [0.933,0.921,0.944,0.945,0.918,0.941,0.981,0.943,0.956,0.947,0.951,0.95,0.897,0.958,0.955,0.942] + +class Tau( Lepton ): + + def __init__(self, tau): + self.tau = tau + super(Tau, self).__init__(tau) + self.eOverP = None + + # JAN FIXME - check what's available in miniAOD and if we need this at all + def calcEOverP(self): + if self.eOverP is not None: + return self.eOverP + self.leadChargedEnergy = self.tau.leadChargedHadrEcalEnergy() \ + + self.tau.leadChargedHadrHcalEnergy() + # self.leadChargedMomentum = self.tau.leadChargedHadrPt() / math.sin(self.tau.theta()) + self.leadChargedMomentum = self.tau.leadPFChargedHadrCand().energy() + self.eOverP = self.leadChargedEnergy / self.leadChargedMomentum + return self.eOverP + + def relIso(self, dummy1, dummy2): + '''Just making the tau behave as a lepton.''' + return -1 + + def mvaId(self): + '''For a transparent treatment of electrons, muons and taus. Returns -99''' + return -99 + + def dxy(self, vertex=None): + if vertex is None: + vertex = self.associatedVertex + vtx = self.vertex() # FIXME + p4 = self.p4() + return ( - (vtx.x()-vertex.position().x()) * p4.y() + + (vtx.y()-vertex.position().y()) * p4.x() ) / p4.pt() + + def dz(self, vertex=None): + if vertex is None: + vertex = self.associatedVertex + vtx = self.vertex() # FIXME + p4 = self.p4() + return (vtx.z()-vertex.position().z()) - ((vtx.x()-vertex.position().x())*p4.x()+(vtx.y()-vertex.position().y())*p4.y())/ p4.pt() * p4.z()/ p4.pt() + + def zImpact(self, vertex=None): + '''z impact at ECAL surface''' + if vertex is None: + vertex = self.associatedVertex + return vertex.z() + 130./math.tan(self.theta()) + + def __str__(self): + lep = super(Tau, self).__str__() + return lep + #spec = '\t\tTau: decay = {decMode:<15}, eOverP = {eOverP:4.2f}, isoMVALoose = {isoMVALoose}'.format( + # decMode = tauDecayModes.intToName( self.decayMode() ), + # eOverP = self.calcEOverP(), + # isoMVALoose = self.tauID('byLooseIsoMVA') + # ) + #return '\n'.join([lep, spec]) + + def electronMVA3Medium(self): + '''Custom electron MVA 3 medium WP used for H->tau tau''' + icat = int(round(self.tauID('againstElectronMVA3category'))) + if icat < 0: + return False + elif icat > 15: + return True + + rawMVA = self.tauID('againstElectronMVA3raw') + return rawMVA > cutsElectronMVA3Medium[icat] + + +def isTau(leg): + '''Duck-typing a tau''' + try: + # Taken from BaseTau to work for both PFTaus and PAT Taus + # (can maybe find a less expensive method) + leg.leadTrack() + except AttributeError: + return False + return True + diff --git a/PhysicsTools/Heppy/python/physicsobjects/classes/Electron.py b/PhysicsTools/Heppy/python/physicsobjects/classes/Electron.py deleted file mode 100644 index 0132385f1046d..0000000000000 --- a/PhysicsTools/Heppy/python/physicsobjects/classes/Electron.py +++ /dev/null @@ -1,196 +0,0 @@ -from PhysicsTools.Heppy.physicsobjects.Lepton import Lepton -from PhysicsTools.Heppy.physicsutils.ElectronMVAID import ElectronMVAID_Trig, ElectronMVAID_NonTrig, ElectronMVAID_TrigNoIP -import ROOT - -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) - diff --git a/PhysicsTools/Heppy/python/physicsobjects/classes/IsoTrack.py b/PhysicsTools/Heppy/python/physicsobjects/classes/IsoTrack.py deleted file mode 100644 index e062f85349e73..0000000000000 --- a/PhysicsTools/Heppy/python/physicsobjects/classes/IsoTrack.py +++ /dev/null @@ -1,29 +0,0 @@ -from PhysicsTools.Heppy.physicsobjects.PhysicsObject import * -import ROOT -import math - -class IsoTrack( PhysicsObject ): - - def __init__(self, isoTrack): - self.isoTrack = isoTrack - super(IsoTrack, self).__init__(isoTrack) - - def absIso(self, dummy1, dummy2): - '''Just making the tau behave as a lepton.''' - return -1 - - def relIso(self, dummy1, dummy2): - '''Just making the tau behave as a lepton.''' - return -1 - - def __str__(self): - lep = super(IsoTrack, self).__str__() - return lep - #spec = '\t\tTau: decay = {decMode:<15}, eOverP = {eOverP:4.2f}, isoMVALoose = {isoMVALoose}'.format( - # decMode = tauDecayModes.intToName( self.decayMode() ), - # eOverP = self.calcEOverP(), - # isoMVALoose = self.tauID('byLooseIsoMVA') - # ) - #return '\n'.join([lep, spec]) - - diff --git a/PhysicsTools/Heppy/python/physicsobjects/classes/Jet.py b/PhysicsTools/Heppy/python/physicsobjects/classes/Jet.py deleted file mode 100644 index 22e6d89a74082..0000000000000 --- a/PhysicsTools/Heppy/python/physicsobjects/classes/Jet.py +++ /dev/null @@ -1,106 +0,0 @@ -from PhysicsTools.Heppy.physicsobjects.PhysicsObject import * - -loose_WP = [ - (0, 2.5, -0.8), - (2.5, 2.75, -0.74), - (2.75, 3.0, -0.68), - (3.0, 5.0, -0.77), - ] - -# Working point 2 May 2013 (Phil via H2tau list) -loose_53X_WP = [ - (0, 2.5, -0.63), - (2.5, 2.75, -0.60), - (2.75, 3.0, -0.55), - (3.0, 5.2, -0.45), - ] - -_btagWPs = { - "TCHEL": ("trackCountingHighEffBJetTags", 1.7), - "TCHEM": ("trackCountingHighEffBJetTags", 3.3), - "TCHPT": ("trackCountingHighPurBJetTags", 3.41), - "JPL": ("jetProbabilityBJetTags", 0.275), - "JPM": ("jetProbabilityBJetTags", 0.545), - "JPT": ("jetProbabilityBJetTags", 0.790), - "CSVL": ("combinedSecondaryVertexBJetTags", 0.244), - "CSVM": ("combinedSecondaryVertexBJetTags", 0.679), - "CSVT": ("combinedSecondaryVertexBJetTags", 0.898), - "CSVv2IVFL": ("combinedInclusiveSecondaryVertexV2BJetTags", 0.423), - "CSVv2IVFM": ("combinedInclusiveSecondaryVertexV2BJetTags", 0.814), - "CSVv2IVFT": ("combinedInclusiveSecondaryVertexV2BJetTags", 0.941), - "CMVAL": ("combinedMVABJetTags", 0.630), # for same b-jet efficiency of CSVv2IVFL on ttbar MC, jet pt > 30 - "CMVAM": ("combinedMVABJetTags", 0.732), # for same b-jet efficiency of CSVv2IVFM on ttbar MC, jet pt > 30 - "CMVAT": ("combinedMVABJetTags", 0.813), # for same b-jet efficiency of CSVv2IVFT on ttbar MC, jet pt > 30 -} - -class Jet(PhysicsObject): - def __init__(self, *args, **kwargs): - super(Jet, self).__init__(*args, **kwargs) - self._physObjInit() - - def _physObjInit(self): - self._rawFactorMultiplier = 1.0 - - def jetID(self,name=""): - if not self.isPFJet(): - raise RuntimeError, "jetID implemented only for PF Jets" - eta = abs(self.eta()); - chf = self.chargedHadronEnergyFraction(); - nhf = self.neutralHadronEnergyFraction(); - phf = self.neutralEmEnergyFraction(); - muf = self.muonEnergyFraction(); - elf = self.chargedEmEnergyFraction(); - chm = self.chargedHadronMultiplicity(); - npr = self.chargedMultiplicity() + self.neutralMultiplicity(); - #if npr != self.nConstituents(): - # import pdb; pdb.set_trace() - if name == "POG_PFID": - if self.jetID("POG_PFID_Tight") : return 3; - elif self.jetID("POG_PFID_Medium") : return 2; - elif self.jetID("POG_PFID_Loose") : return 1; - else : return 0; - - if name == "POG_PFID_Loose": return (npr>1 and phf<0.99 and nhf<0.99) and (eta>2.4 or (elf<0.99 and chf>0 and chm>0)); - if name == "POG_PFID_Medium": return (npr>1 and phf<0.95 and nhf<0.95) and (eta>2.4 or (elf<0.99 and chf>0 and chm>0)); - if name == "POG_PFID_Tight": return (npr>1 and phf<0.90 and nhf<0.90) and (eta>2.4 or (elf<0.99 and chf>0 and chm>0)); - if name == "VBFHBB_PFID_Loose": return (npr>1 and phf<0.99 and nhf<0.99); - if name == "VBFHBB_PFID_Medium": return (npr>1 and phf<0.99 and nhf<0.99) and ((eta<=2.4 and nhf<0.9 and phf<0.9 and elf<0.99 and muf<0.99 and chf>0 and chm>0) or eta>2.4); - if name == "VBFHBB_PFID_Tight": return (npr>1 and phf<0.99 and nhf<0.99) and ((eta<=2.4 and nhf<0.9 and phf<0.9 and elf<0.70 and muf<0.70 and chf>0 and chm>0) or eta>2.4); - raise RuntimeError, "jetID '%s' not supported" % name - - def looseJetId(self): - '''PF Jet ID (loose operation point) [method provided for convenience only]''' - return self.jetID("POG_PFID_Loose") - - def puMva(self, label="pileupJetId:fullDiscriminant"): - return self.userFloat(label) - - def puJetId(self, label="pileupJetId:fullDiscriminant"): - '''Full mva PU jet id''' - - puMva = self.puMva(label) - wp = loose_53X_WP - eta = abs(self.eta()) - - for etamin, etamax, cut in wp: - if not(eta>=etamin and etacut - - def rawFactor(self): - return self.jecFactor('Uncorrected') * self._rawFactorMultiplier - def setRawFactor(self, factor): - self._rawFactorMultiplier = factor/self.jecFactor('Uncorrected') - - def btag(self,name): - return self.bDiscriminator(name) - - def btagWP(self,name): - global _btagWPs - (disc,val) = _btagWPs[name] - return self.bDiscriminator(disc) > val - - -class GenJet( PhysicsObject): - pass - diff --git a/PhysicsTools/Heppy/python/physicsobjects/classes/Muon.py b/PhysicsTools/Heppy/python/physicsobjects/classes/Muon.py deleted file mode 100644 index a4bed67b83e2d..0000000000000 --- a/PhysicsTools/Heppy/python/physicsobjects/classes/Muon.py +++ /dev/null @@ -1,89 +0,0 @@ -from PhysicsTools.Heppy.physicsobjects.Lepton import Lepton - -class Muon( Lepton ): - - def looseId( self ): - '''Loose ID as recommended by mu POG.''' - return self.physObj.isLooseMuon() - - def tightId( self ): - '''Tight ID as recommended by mu POG.''' - return self.muonID("POG_ID_Tight") - - def muonID(self, name, vertex=None): - if name == "" or name is None: - return True - if name.startswith("POG_"): - if name == "POG_ID_Loose": return self.physObj.isLooseMuon() - if vertex is None: - vertex = getattr(self, 'associatedVertex', None) - if name == "POG_ID_Tight": return self.physObj.isTightMuon(vertex) - if name == "POG_ID_HighPt": return self.physObj.isHighPtMuon(vertex) - if name == "POG_ID_Soft": return self.physObj.isSoftMuon(vertex) - if name == "POG_ID_TightNoVtx": return self.looseId() and \ - self.isGlobalMuon() and \ - self.globalTrack().normalizedChi2() < 10 and \ - self.globalTrack().hitPattern().numberOfValidMuonHits() > 0 and \ - self.numberOfMatchedStations()>1 and \ - self.innerTrack().hitPattern().numberOfValidPixelHits()>0 and \ - self.innerTrack().hitPattern().trackerLayersWithMeasurement() > 5 - return self.physObj.muonID(name) - - def mvaId(self): - '''For a transparent treatment of electrons and muons. Returns -99''' - return -99 - - - - def absEffAreaIso(self,rho,effectiveAreas): - return self.absIsoFromEA(rho,self.eta(),effectiveAreas.muon) - - - - def dxy(self, vertex=None): - '''either pass the vertex, or set associatedVertex before calling the function. - note: the function does not work with standalone muons as innerTrack - is not available. - ''' - if vertex is None: - vertex = self.associatedVertex - return self.innerTrack().dxy( vertex.position() ) - - - def dz(self, vertex=None): - '''either pass the vertex, or set associatedVertex before calling the function. - note: the function does not work with standalone muons as innerTrack - is not available. - ''' - if vertex is None: - vertex = self.associatedVertex - return self.innerTrack().dz( vertex.position() ) - - def chargedHadronIsoR(self,R=0.4): - if R == 0.3: return self.physObj.pfIsolationR03().sumChargedHadronPt - elif R == 0.4: return self.physObj.pfIsolationR04().sumChargedHadronPt - raise RuntimeError, "Muon chargedHadronIso missing for R=%s" % R - - def neutralHadronIsoR(self,R=0.4): - if R == 0.3: return self.physObj.pfIsolationR03().sumNeutralHadronEt - elif R == 0.4: return self.physObj.pfIsolationR04().sumNeutralHadronEt - raise RuntimeError, "Muon neutralHadronIso missing for R=%s" % R - - def photonIsoR(self,R=0.4): - if R == 0.3: return self.physObj.pfIsolationR03().sumPhotonEt - elif R == 0.4: return self.physObj.pfIsolationR04().sumPhotonEt - raise RuntimeError, "Muon photonIso missing for R=%s" % R - - def chargedAllIsoR(self,R=0.4): - if R == 0.3: return self.physObj.pfIsolationR03().sumChargedParticlePt - elif R == 0.4: return self.physObj.pfIsolationR04().sumChargedParticlePt - raise RuntimeError, "Muon chargedAllIso missing for R=%s" % R - - def puChargedHadronIsoR(self,R=0.4): - if R == 0.3: return self.physObj.pfIsolationR03().sumPUPt - elif R == 0.4: return self.physObj.pfIsolationR04().sumPUPt - raise RuntimeError, "Muon chargedHadronIso missing for R=%s" % R - - - - diff --git a/PhysicsTools/Heppy/python/physicsobjects/classes/Photon.py b/PhysicsTools/Heppy/python/physicsobjects/classes/Photon.py deleted file mode 100644 index 371b25a057375..0000000000000 --- a/PhysicsTools/Heppy/python/physicsobjects/classes/Photon.py +++ /dev/null @@ -1,39 +0,0 @@ -from PhysicsTools.Heppy.physicsobjects.PhysicsObject import * - -class Photon(PhysicsObject ): - - ''' return object from the photon - ''' - def hOVERe(self): -# return self.physObj.full5x5_hadTowOverEm() - hadTowDepth1O = self.physObj.hadTowDepth1OverEm() * (self.physObj.superCluster().energy()/self.physObj.full5x5_e5x5() if self.physObj.full5x5_e5x5() else 1) - hadTowDepth2O = self.physObj.hadTowDepth2OverEm() * (self.physObj.superCluster().energy()/self.physObj.full5x5_e5x5() if self.physObj.full5x5_e5x5() else 1) - return hadTowDepth1O + hadTowDepth2O - - def r9(self): - return self.physObj.full5x5_r9() - - def sigmaIetaIeta(self): - return self.physObj.full5x5_sigmaIetaIeta() - - def chargedHadronIso(self): - return self.physObj.chargedHadronIso() - - def neutralHadronIso(self): - return self.physObj.neutralHadronIso() - - def photonIso(self): - return self.physObj.photonIso() - - def photonIDCSA14(self, name): - keepThisPhoton = True - if name == "PhotonCutBasedIDLoose_CSA14": - if abs(self.physObj.eta())<1.479 : - if self.physObj.sigmaIetaIeta() > 0.012 : keepThisPhoton = False - if self.hOVERe() > 0.0559 : keepThisPhoton = False - else : - if self.physObj.sigmaIetaIeta() > 0.035 : keepThisPhoton = False - if self.hOVERe() > 0.049 : keepThisPhoton = False - return keepThisPhoton - - pass diff --git a/PhysicsTools/Heppy/python/physicsobjects/classes/Tau.py b/PhysicsTools/Heppy/python/physicsobjects/classes/Tau.py deleted file mode 100644 index e20506148b06c..0000000000000 --- a/PhysicsTools/Heppy/python/physicsobjects/classes/Tau.py +++ /dev/null @@ -1,85 +0,0 @@ -from PhysicsTools.Heppy.physicsobjects.Lepton import Lepton -from PhysicsTools.Heppy.physicsutils.TauDecayModes import tauDecayModes -import math - -cutsElectronMVA3Medium = [0.933,0.921,0.944,0.945,0.918,0.941,0.981,0.943,0.956,0.947,0.951,0.95,0.897,0.958,0.955,0.942] - -class Tau( Lepton ): - - def __init__(self, tau): - self.tau = tau - super(Tau, self).__init__(tau) - self.eOverP = None - - # JAN FIXME - check what's available in miniAOD and if we need this at all - def calcEOverP(self): - if self.eOverP is not None: - return self.eOverP - self.leadChargedEnergy = self.tau.leadChargedHadrEcalEnergy() \ - + self.tau.leadChargedHadrHcalEnergy() - # self.leadChargedMomentum = self.tau.leadChargedHadrPt() / math.sin(self.tau.theta()) - self.leadChargedMomentum = self.tau.leadPFChargedHadrCand().energy() - self.eOverP = self.leadChargedEnergy / self.leadChargedMomentum - return self.eOverP - - def relIso(self, dummy1, dummy2): - '''Just making the tau behave as a lepton.''' - return -1 - - def mvaId(self): - '''For a transparent treatment of electrons, muons and taus. Returns -99''' - return -99 - - def dxy(self, vertex=None): - if vertex is None: - vertex = self.associatedVertex - vtx = self.vertex() # FIXME - p4 = self.p4() - return ( - (vtx.x()-vertex.position().x()) * p4.y() - + (vtx.y()-vertex.position().y()) * p4.x() ) / p4.pt() - - def dz(self, vertex=None): - if vertex is None: - vertex = self.associatedVertex - vtx = self.vertex() # FIXME - p4 = self.p4() - return (vtx.z()-vertex.position().z()) - ((vtx.x()-vertex.position().x())*p4.x()+(vtx.y()-vertex.position().y())*p4.y())/ p4.pt() * p4.z()/ p4.pt() - - def zImpact(self, vertex=None): - '''z impact at ECAL surface''' - if vertex is None: - vertex = self.associatedVertex - return vertex.z() + 130./math.tan(self.theta()) - - def __str__(self): - lep = super(Tau, self).__str__() - return lep - #spec = '\t\tTau: decay = {decMode:<15}, eOverP = {eOverP:4.2f}, isoMVALoose = {isoMVALoose}'.format( - # decMode = tauDecayModes.intToName( self.decayMode() ), - # eOverP = self.calcEOverP(), - # isoMVALoose = self.tauID('byLooseIsoMVA') - # ) - #return '\n'.join([lep, spec]) - - def electronMVA3Medium(self): - '''Custom electron MVA 3 medium WP used for H->tau tau''' - icat = int(round(self.tauID('againstElectronMVA3category'))) - if icat < 0: - return False - elif icat > 15: - return True - - rawMVA = self.tauID('againstElectronMVA3raw') - return rawMVA > cutsElectronMVA3Medium[icat] - - -def isTau(leg): - '''Duck-typing a tau''' - try: - # Taken from BaseTau to work for both PFTaus and PAT Taus - # (can maybe find a less expensive method) - leg.leadTrack() - except AttributeError: - return False - return True - From 411274c2f5aea76ab85aaa3f9104d28676fa612a Mon Sep 17 00:00:00 2001 From: Andrea Date: Fri, 19 Dec 2014 11:35:29 +0100 Subject: [PATCH 07/10] add __str__ to reco::GenParticle and remove legacy GenParticle/GenLepton --- .../python/physicsobjects/GenParticle.py | 30 +++++-------------- .../python/physicsobjects/PhysicsObjects.py | 2 +- 2 files changed, 8 insertions(+), 24 deletions(-) diff --git a/PhysicsTools/Heppy/python/physicsobjects/GenParticle.py b/PhysicsTools/Heppy/python/physicsobjects/GenParticle.py index 54219fddd70d3..c5a26a85247a4 100644 --- a/PhysicsTools/Heppy/python/physicsobjects/GenParticle.py +++ b/PhysicsTools/Heppy/python/physicsobjects/GenParticle.py @@ -1,28 +1,12 @@ from PhysicsTools.Heppy.physicsobjects.PhysicsObject import * -class GenParticle( PhysicsObject ): - def __str__(self): - base = super(GenParticle, self).__str__() +#add __str__ to reco::GenParticle python wrapper +import ROOT +def printGenParticle(self): + base = basePrint(self) theStr = '{base}, status = {status:>2}'.format(base=base, status=self.status()) return theStr +setattr(ROOT.reco.GenParticle,"basePrint",Particle.__str__) +setattr(ROOT.reco.GenParticle,"__str__",printGenParticle) - -class GenLepton( GenParticle ): - def sip3D(self): - '''Just to make generic code work on GenParticles''' - return 0 - def relIso(self, dummy): - '''Just to make generic code work on GenParticles''' - return 0 - - def absIso(self, dummy): - '''Just to make generic code work on GenParticles''' - return 0 - - def absEffAreaIso(self,rho): - '''Just to make generic code work on GenParticles''' - return 0 - - def relEffAreaIso(self,rho): - '''Just to make generic code work on GenParticles''' - return 0 +from ROOT.reco import GenParticle diff --git a/PhysicsTools/Heppy/python/physicsobjects/PhysicsObjects.py b/PhysicsTools/Heppy/python/physicsobjects/PhysicsObjects.py index 99bd7024afe1c..68973a642f54a 100644 --- a/PhysicsTools/Heppy/python/physicsobjects/PhysicsObjects.py +++ b/PhysicsTools/Heppy/python/physicsobjects/PhysicsObjects.py @@ -17,5 +17,5 @@ def printOut(objects): # COLIN need to import MVA ID recipe. # from PhysicsTools.Heppy.physicsobjects.Electron import Electron from PhysicsTools.Heppy.physicsobjects.Tau import Tau, isTau -from PhysicsTools.Heppy.physicsobjects.GenParticle import GenParticle,GenLepton +from PhysicsTools.Heppy.physicsobjects.GenParticle import GenParticle From 191d2680444e4e78fd5dc5a0680a9057e4859eea Mon Sep 17 00:00:00 2001 From: Andrea Date: Fri, 19 Dec 2014 11:36:00 +0100 Subject: [PATCH 08/10] fix typo in autovars showing up when isMC = False --- PhysicsTools/Heppy/python/analyzers/core/autovars.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PhysicsTools/Heppy/python/analyzers/core/autovars.py b/PhysicsTools/Heppy/python/analyzers/core/autovars.py index 5126cd921c48e..5ba729e6aef2b 100644 --- a/PhysicsTools/Heppy/python/analyzers/core/autovars.py +++ b/PhysicsTools/Heppy/python/analyzers/core/autovars.py @@ -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 = {} From 5be700a78fe7b582b1b91f42bc16f7f6ba81aee5 Mon Sep 17 00:00:00 2001 From: Andrea Date: Fri, 19 Dec 2014 11:36:39 +0100 Subject: [PATCH 09/10] Add default configs and option to store individual bits for trigger studies --- .../analyzers/core/TriggerBitAnalyzer.py | 55 ++++++++++++++++++- 1 file changed, 54 insertions(+), 1 deletion(-) diff --git a/PhysicsTools/Heppy/python/analyzers/core/TriggerBitAnalyzer.py b/PhysicsTools/Heppy/python/analyzers/core/TriggerBitAnalyzer.py index 7cf022904f42a..5e9195fc5bdc3 100644 --- a/PhysicsTools/Heppy/python/analyzers/core/TriggerBitAnalyzer.py +++ b/PhysicsTools/Heppy/python/analyzers/core/TriggerBitAnalyzer.py @@ -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() @@ -17,15 +20,31 @@ 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 ) @@ -33,6 +52,40 @@ def process(self, event): 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 = { + # "" : [ 'HLT__v*', 'HLT__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" ], + } +) +) From 35228a2f2902fd19066f50e0f9144951b6654b09 Mon Sep 17 00:00:00 2001 From: Andrea Date: Fri, 19 Dec 2014 12:19:15 +0100 Subject: [PATCH 10/10] reusing Particle.__str__ seem not to work. better would be to have a free function baseP4Printing to be shared for different cases --- PhysicsTools/Heppy/python/physicsobjects/GenParticle.py | 9 +++++++-- 1 file changed, 7 insertions(+), 2 deletions(-) diff --git a/PhysicsTools/Heppy/python/physicsobjects/GenParticle.py b/PhysicsTools/Heppy/python/physicsobjects/GenParticle.py index c5a26a85247a4..984800b34a5bf 100644 --- a/PhysicsTools/Heppy/python/physicsobjects/GenParticle.py +++ b/PhysicsTools/Heppy/python/physicsobjects/GenParticle.py @@ -3,10 +3,15 @@ #add __str__ to reco::GenParticle python wrapper import ROOT def printGenParticle(self): - base = basePrint(self) + tmp = '{className} : {pdgId:>3}, pt = {pt:5.1f}, eta = {eta:5.2f}, phi = {phi:5.2f}, mass = {mass:5.2f}' + base= tmp.format( className = self.__class__.__name__, + pdgId = self.pdgId(), + pt = self.pt(), + eta = self.eta(), + phi = self.phi(), + mass = self.mass() ) theStr = '{base}, status = {status:>2}'.format(base=base, status=self.status()) return theStr -setattr(ROOT.reco.GenParticle,"basePrint",Particle.__str__) setattr(ROOT.reco.GenParticle,"__str__",printGenParticle) from ROOT.reco import GenParticle