File indexing completed on 2024-04-06 12:23:28
0001 import math
0002 from ROOT import TLorentzVector
0003 from PhysicsTools.HeppyCore.utils.deltar import deltaPhi, deltaR2
0004
0005 class VBF( object ):
0006 '''Computes and holds VBF quantities'''
0007 def __init__(self, jets, diLepton, vbfMvaCalc, cjvPtCut):
0008 '''jets: jets cleaned from the diLepton legs.
0009 diLepton: the di-tau, for example. Necessary to compute input variables for MVA selection
0010 '''
0011 self.cjvPtCut = cjvPtCut
0012 self.vbfMvaCalc = vbfMvaCalc
0013 self.jets = jets
0014
0015
0016 self.met = diLepton.met()
0017 self.leadJets = jets[:2]
0018 self.otherJets = jets[2:]
0019 self.centralJets = self.findCentralJets( self.leadJets, self.otherJets )
0020
0021
0022 self.deta = self.leadJets[0].eta() - self.leadJets[1].eta()
0023
0024
0025
0026 self.dphi = deltaPhi(self.leadJets[0].phi(), self.leadJets[1].phi())
0027 dijetp4 = self.leadJets[0].p4() + self.leadJets[1].p4()
0028
0029 self.mjj = dijetp4.M()
0030
0031 self.dijetpt = dijetp4.pt()
0032
0033 self.dijetphi = dijetp4.phi()
0034
0035
0036 self.higgsp4 = diLepton.p4() + self.met.p4()
0037
0038 self.dphidijethiggs = deltaPhi( self.dijetphi, self.higgsp4.phi() )
0039
0040 visDiLepton = diLepton.leg1 ().p4 () + diLepton.leg2 ().p4 ()
0041 self.visjeteta = min (
0042 abs (self.leadJets[0].eta () - visDiLepton.eta ()),
0043 abs (self.leadJets[1].eta () - visDiLepton.eta ()))
0044
0045 self.ptvis = visDiLepton.pt()
0046
0047
0048 if self.vbfMvaCalc is not None:
0049 self.mva = self.vbfMvaCalc.val( self.mjj,
0050 abs(self.deta),
0051 self.visjeteta,
0052 self.ptvis )
0053 else:
0054 self.mva = -99.
0055
0056
0057
0058
0059
0060
0061
0062
0063
0064
0065
0066 def findCentralJets( self, leadJets, otherJets ):
0067 '''Finds all jets between the 2 leading jets, for central jet veto.'''
0068 if not len(otherJets):
0069 return []
0070 etamin = leadJets[0].eta()
0071 etamax = leadJets[1].eta()
0072 if etamin > etamax:
0073 etamin, etamax = etamax, etamin
0074 def isCentral( jet ):
0075 if jet.pt()<self.cjvPtCut:
0076 return False
0077 eta = jet.eta()
0078 if etamin < eta and eta < etamax:
0079 return True
0080 else:
0081 return False
0082 centralJets = list(filter( isCentral, otherJets ))
0083 return centralJets
0084
0085 def calcP4(self, jets):
0086 '''returns the sum p4 of a collection of objects.
0087 FIXME: remove this function, which is a bit stupid
0088 '''
0089 p4 = TLorentzVector()
0090 for jet in jets:
0091 p4 += TLorentzVector(jet.px(), jet.py(), jet.pz(), jet.energy())
0092 return p4
0093
0094 def __str__(self):
0095 header = 'VBF : deta={deta:4.2f}, Mjj={mjj:4.2f}, #centjets={ncentjets}'
0096 header = header.format( deta=self.deta, mjj=self.mjj, ncentjets=len(self.centralJets))
0097 leadJets = map( str, self.leadJets )
0098 centralJets = map( str, self.centralJets)
0099 tmp = [header]
0100 tmp.append('MVA input variables: dphi={dphi:4.2f}, dijetpt={dijetpt:4.2f}, dijetphi={dijetphi:4.2f}, dphidijethiggs={dphidijethiggs:4.2f}, visjeteta={visjeteta:4.2f}, ptvis={ptvis:4.2f}'.format(
0101 dphi = self.dphi,
0102 dijetpt = self.dijetpt,
0103 dijetphi = self.dijetphi,
0104 dphidijethiggs = self.dphidijethiggs,
0105 visjeteta = self.visjeteta,
0106 ptvis = self.ptvis
0107 ))
0108 tmp.append('Leading Jets:')
0109 tmp.extend( leadJets )
0110 tmp.append('Central Jets:')
0111 tmp.extend( centralJets )
0112 return '\n'.join( tmp )