7 from ROOT
import TLorentzVector, TVectorD
9 from PhysicsTools.Heppy.analyzers.core.Analyzer
import Analyzer
10 from PhysicsTools.HeppyCore.framework.event
import Event
11 from PhysicsTools.HeppyCore.statistics.counter
import Counter, Counters
12 from PhysicsTools.Heppy.analyzers.core.AutoHandle
import AutoHandle
14 import PhysicsTools.HeppyCore.framework.config
as cfg
18 from ROOT.heppy
import Hemisphere
19 from ROOT.heppy
import ReclusterJets
21 from ROOT.heppy
import Davismt2
24 from ROOT.heppy
import mt2w_bisect
25 mt2wSNT = mt2w_bisect.mt2w()
33 def __init__(self, cfg_ana, cfg_comp, looperName ):
34 super(MT2Analyzer,self).
__init__(cfg_ana,cfg_comp,looperName)
40 self.handles[
'genJets'] = AutoHandle(
'slimmedGenJets',
'std::vector<reco::GenJet>')
41 self.handles[
'met'] = AutoHandle( self.cfg_ana.metCollection,
'std::vector<pat::MET>' )
45 self.counters.addCounter(
'pairs')
46 count = self.counters.counter(
'pairs')
47 count.register(
'all events')
54 metVector = array.array(
'd',[0.,metVec.px(), metVec.py()])
55 visaVector = array.array(
'd',[0.,visaVec.px(), visaVec.py()])
56 visbVector = array.array(
'd',[0.,visbVec.px(), visbVec.py()])
58 davismt2.set_momenta(visaVector,visbVector,metVector);
61 return davismt2.get_mt2()
63 def getMT2AKT(self, event, TMPobjects40jc, met , collectionPostFix, postFix):
66 if len(TMPobjects40jc)>=2:
68 objects = ROOT.std.vector(ROOT.reco.Particle.LorentzVector)()
69 for jet
in TMPobjects40jc:
70 objects.push_back(jet.p4())
72 hemisphereViaKt = ReclusterJets(objects, 1.,50.0)
73 groupingViaKt=hemisphereViaKt.getGroupingExclusive(2)
75 if len(groupingViaKt)>=2:
76 setattr(event,
"pseudoViaKtJet1"+collectionPostFix+postFix, ROOT.reco.Particle.LorentzVector(groupingViaKt[0]) )
77 setattr(event,
"pseudoViaKtJet2"+collectionPostFix+postFix, ROOT.reco.Particle.LorentzVector(groupingViaKt[1]) )
78 setattr(event,
"mt2ViaAKt"+collectionPostFix+postFix, self.
computeMT2(getattr(event,
'pseudoViaKtJet1'+collectionPostFix+postFix), getattr(event,
'pseudoViaKtJet2'+collectionPostFix+postFix), met) )
79 return self.
computeMT2(getattr(event,
'pseudoViaKtJet1'+collectionPostFix+postFix), getattr(event,
'pseudoViaKtJet2'+collectionPostFix+postFix), met)
81 if not self.cfg_ana.doOnlyDefault:
82 hemisphereViaAKt = ReclusterJets(objects, -1.,50.0)
83 groupingViaAKt=hemisphereViaAKt.getGroupingExclusive(2)
85 if len(groupingViaAKt)>=2:
86 setattr(event,
"pseudoViaAKtJet1"+collectionPostFix+postFix, ROOT.reco.Particle.LorentzVector(groupingViaAKt[0]) )
87 setattr(event,
"pseudoViaAKtJet2"+collectionPostFix+postFix, ROOT.reco.Particle.LorentzVector(groupingViaAKt[1]) )
88 setattr(event,
"mt2ViaAKt"+collectionPostFix+postFix, self.
computeMT2(getattr(event,
'pseudoViaAKtJet1'+collectionPostFix+postFix), getattr(event,
'pseudoViaAKtJet2'+collectionPostFix+postFix), met) )
89 return self.
computeMT2(getattr(event,
'pseudoViaAKtJet1'+collectionPostFix+postFix), getattr(event,
'pseudoViaAKtJet2'+collectionPostFix+postFix), met)
91 def getMT2Hemi(self, event, TMPobjects40jc, met, collectionPostFix, postFix):
93 if len(TMPobjects40jc)>=2:
95 pxvec = ROOT.std.vector(float)()
96 pyvec = ROOT.std.vector(float)()
97 pzvec = ROOT.std.vector(float)()
98 Evec = ROOT.std.vector(float)()
99 grouping = ROOT.std.vector(int)()
101 for jet
in TMPobjects40jc:
102 pxvec.push_back(jet.px())
103 pyvec.push_back(jet.py())
104 pzvec.push_back(jet.pz())
105 Evec.push_back(jet.energy())
107 hemisphere = Hemisphere(pxvec, pyvec, pzvec, Evec, 2, 3)
108 grouping=hemisphere.getGrouping()
122 for index
in range(0, len(pxvec)):
123 if(grouping[index]==1):
124 pseudoJet1px += pxvec[index]
125 pseudoJet1py += pyvec[index]
126 pseudoJet1pz += pzvec[index]
127 pseudoJet1energy += Evec[index]
129 if(grouping[index]==2):
130 pseudoJet2px += pxvec[index]
131 pseudoJet2py += pyvec[index]
132 pseudoJet2pz += pzvec[index]
133 pseudoJet2energy += Evec[index]
136 pseudoJet1pt2 = pseudoJet1px*pseudoJet1px + pseudoJet1py*pseudoJet1py
137 pseudoJet2pt2 = pseudoJet2px*pseudoJet2px + pseudoJet2py*pseudoJet2py
139 if pseudoJet1pt2 >= pseudoJet2pt2:
140 setattr(event,
"pseudoJet1"+collectionPostFix+postFix, ROOT.reco.Particle.LorentzVector( pseudoJet1px, pseudoJet1py, pseudoJet1pz, pseudoJet1energy ))
141 setattr(event,
"pseudoJet2"+collectionPostFix+postFix, ROOT.reco.Particle.LorentzVector( pseudoJet2px, pseudoJet2py, pseudoJet2pz, pseudoJet2energy ))
142 setattr(event,
"multPseudoJet1"+collectionPostFix+postFix, multPSJ1 )
143 setattr(event,
"multPseudoJet2"+collectionPostFix+postFix, multPSJ2 )
145 setattr(event,
"pseudoJet2"+collectionPostFix+postFix, ROOT.reco.Particle.LorentzVector( pseudoJet1px, pseudoJet1py, pseudoJet1pz, pseudoJet1energy ))
146 setattr(event,
"pseudoJet1"+collectionPostFix+postFix, ROOT.reco.Particle.LorentzVector( pseudoJet2px, pseudoJet2py, pseudoJet2pz, pseudoJet2energy ))
147 setattr(event,
"multPseudoJet1"+collectionPostFix+postFix, multPSJ2 )
148 setattr(event,
"multPseudoJet2"+collectionPostFix+postFix, multPSJ1 )
150 setattr(event,
"mt2"+collectionPostFix+postFix, self.
computeMT2(getattr(event,
'pseudoJet1'+collectionPostFix+postFix), getattr(event,
'pseudoJet2'+collectionPostFix+postFix), met) )
151 return self.
computeMT2(getattr(event,
'pseudoJet1'+collectionPostFix+postFix), getattr(event,
'pseudoJet2'+collectionPostFix+postFix), met)
158 self.
met = ROOT.pat.MET(self.handles[
'met'].product()[0])
163 objects40jc = [ j
for j
in event.cleanJets
if j.pt() > 40
and abs(j.eta())<2.5 ]
164 objectsXjc = [ j
for j
in event.cleanJets
if j.pt() > self.
jetPt and abs(j.eta())<2.5 ]
166 setattr(event,
"mt2ViaKt"+self.cfg_ana.collectionPostFix+
"had", -999)
167 setattr(event,
"mt2ViaKt"+self.cfg_ana.collectionPostFix+
"_Xj_had", -999)
168 setattr(event,
"pseudoViaKtJet1"+self.cfg_ana.collectionPostFix+
"_had", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
169 setattr(event,
"pseudoViaKtJet2"+self.cfg_ana.collectionPostFix+
"_had", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
170 setattr(event,
"pseudoViaKtJet1"+self.cfg_ana.collectionPostFix+
"_Xj_had", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
171 setattr(event,
"pseudoViaKtJet2"+self.cfg_ana.collectionPostFix+
"_Xj_had", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
173 setattr(event,
"mt2ViaAKt"+self.cfg_ana.collectionPostFix+
"had", -999)
174 setattr(event,
"mt2ViaAKt"+self.cfg_ana.collectionPostFix+
"_Xj_had", -999)
175 setattr(event,
"pseudoViaAKtJet1"+self.cfg_ana.collectionPostFix+
"_had", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
176 setattr(event,
"pseudoViaAKtJet2"+self.cfg_ana.collectionPostFix+
"_had", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
177 setattr(event,
"pseudoViaAKtJet1"+self.cfg_ana.collectionPostFix+
"_Xj_had", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
178 setattr(event,
"pseudoViaAKtJet2"+self.cfg_ana.collectionPostFix+
"_Xj_had", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
181 if len(objects40jc)>=2:
183 self.mt2ViaKt_had=self.
getMT2AKT(event, objects40jc, self.
met, self.cfg_ana.collectionPostFix,
"_had")
185 if len(objectsXjc)>=2:
187 self.mt2ViaKt_Xj_had=self.
getMT2AKT(event, objectsXjc, self.
met, self.cfg_ana.collectionPostFix,
"_Xj_had")
192 setattr(event,
"mt2"+self.cfg_ana.collectionPostFix+
"_had", -999)
193 setattr(event,
"mt2"+self.cfg_ana.collectionPostFix+
"_Xj_had", -999)
195 setattr(event,
"pseudoJet1"+self.cfg_ana.collectionPostFix+
"_had", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
196 setattr(event,
"pseudoJet2"+self.cfg_ana.collectionPostFix+
"_had", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
197 setattr(event,
"pseudoJet1"+self.cfg_ana.collectionPostFix+
"_Xj_had", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
198 setattr(event,
"pseudoJet2"+self.cfg_ana.collectionPostFix+
"_Xj_had", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
200 if len(objects40jc)>=2:
202 self.mt2_had = self.
getMT2Hemi(event,objects40jc, self.
met, self.cfg_ana.collectionPostFix,
"_had")
204 if len(objectsXjc)>=2:
206 self.mt2_Xj_had = self.
getMT2Hemi(event,objectsXjc, self.
met, self.cfg_ana.collectionPostFix,
"_Xj_had")
210 setattr(event,
"mt2"+self.cfg_ana.collectionPostFix+
"_gen", -999)
211 setattr(event,
"mt2"+self.cfg_ana.collectionPostFix+
"_Xj_gen", -999)
213 if self.cfg_comp.isMC
and self.met.genMET():
214 allGenJets = [ x
for x
in self.handles[
'genJets'].product() ]
215 objects40jc_Gen = [ j
for j
in allGenJets
if j.pt() > 40
and abs(j.eta())<2.5 ]
216 objectsXjc_Gen = [ j
for j
in allGenJets
if j.pt() > self.jetPt
and abs(j.eta())<2.5 ]
218 if len(objects40jc_Gen)>=2:
219 self.mt2_gen = self.getMT2Hemi(event,objects40jc_Gen, self.met.genMET(), self.cfg_ana.collectionPostFix,
"_gen")
221 if len(objectsXjc_Gen)>=2:
222 self.mt2_Xj_gen = self.getMT2Hemi(event,objectsXjc_Gen, self.met.genMET(), self.cfg_ana.collectionPostFix,
"_Xj_gen")
230 objects10lc = [ l
for l
in event.selectedLeptons
if l.pt() > 10
and abs(l.eta())<2.5 ]
231 if hasattr(event,
'selectedIsoCleanTrack'):
232 objects10lc = [ l
for l
in event.selectedLeptons
if l.pt() > 10
and abs(l.eta())<2.5 ] + [ t
for t
in event.selectedIsoCleanTrack ]
234 objects40j10lc = objects40jc + objects10lc
235 objects40j10lc.sort(key =
lambda obj : obj.pt(), reverse =
True)
237 objectsXj10lc = objectsXjc + objects10lc
238 objectsXj10lc.sort(key =
lambda obj : obj.pt(), reverse =
True)
240 setattr(event,
"mt2"+self.cfg_ana.collectionPostFix+
"", -999)
241 setattr(event,
"mt2"+self.cfg_ana.collectionPostFix+
"_Xj", -999)
243 setattr(event,
"pseudoJet1"+self.cfg_ana.collectionPostFix+
"", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
244 setattr(event,
"pseudoJet2"+self.cfg_ana.collectionPostFix+
"", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
245 setattr(event,
"pseudoJet1"+self.cfg_ana.collectionPostFix+
"_Xj", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
246 setattr(event,
"pseudoJet2"+self.cfg_ana.collectionPostFix+
"_Xj", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
248 if len(objects40j10lc)>=2:
250 self.mt2 = self.getMT2Hemi(event,objects40j10lc,self.met,self.cfg_ana.collectionPostFix,
"")
252 if len(objectsXj10lc)>=2:
254 self.mt2_Xj = self.getMT2Hemi(event,objectsXj10lc,self.met,self.cfg_ana.collectionPostFix,
"_Xj")
258 setattr(event,
"mt2"+self.cfg_ana.collectionPostFix+
"_gamma", -999)
260 setattr(event,
"pseudoJet1"+self.cfg_ana.collectionPostFix+
"_gamma", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
261 setattr(event,
"pseudoJet2"+self.cfg_ana.collectionPostFix+
"_gamma", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
263 if hasattr(event,
'gamma_met'):
265 gamma_objects40jc = [ j
for j
in event.gamma_cleanJets
if j.pt() > 40
and abs(j.eta())<2.5 ]
267 gamma_objects40j10lc = gamma_objects40jc + objects10lc
269 gamma_objects40j10lc.sort(key =
lambda obj : obj.pt(), reverse =
True)
272 if len(gamma_objects40jc)>=2:
274 self.gamma_mt2 = self.getMT2Hemi(event,gamma_objects40jc,event.gamma_met,self.cfg_ana.collectionPostFix,
"_gamma")
276 setattr(event,
"mt2"+self.cfg_ana.collectionPostFix+
"_Xj_gamma", -999)
277 setattr(event,
"pseudoJet1"+self.cfg_ana.collectionPostFix+
"_Xj_gamma", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
278 setattr(event,
"pseudoJet2"+self.cfg_ana.collectionPostFix+
"_Xj_gamma", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
280 if hasattr(event,
'gamma_met'):
282 gamma_objectsXjc = [ j
for j
in event.gamma_cleanJets
if j.pt() > self.jetPt
and abs(j.eta())<2.5 ]
284 gamma_objectsXj10lc = gamma_objectsXjc + objects10lc
286 gamma_objectsXj10lc.sort(key =
lambda obj : obj.pt(), reverse =
True)
288 if len(gamma_objectsXjc)>=2:
290 self.gamma_mt2_Xj = self.getMT2Hemi(event,gamma_objectsXjc,event.gamma_met,self.cfg_ana.collectionPostFix,
"_Xj_gamma")
296 setattr(event,
"mt2"+self.cfg_ana.collectionPostFix+
"_zll", -999)
297 setattr(event,
"pseudoJet1"+self.cfg_ana.collectionPostFix+
"_zll", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
298 setattr(event,
"pseudoJet2"+self.cfg_ana.collectionPostFix+
"_zll", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
300 if hasattr(event,
'zll_met'):
302 csLeptons = [ l
for l
in event.selectedLeptons
if l.pt() > 10
and abs(l.eta()) < 2.5 ]
304 if len(csLeptons)==2
and len(objects40jc)>=2:
306 self.zll_mt2 = self.getMT2Hemi(event,objects40jc,event.zll_met,self.cfg_ana.collectionPostFix,
"_zll")
308 setattr(event,
"mt2"+self.cfg_ana.collectionPostFix+
"_Xj_zll", -999)
309 setattr(event,
"pseudoJet1"+self.cfg_ana.collectionPostFix+
"_Xj_zll", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
310 setattr(event,
"pseudoJet2"+self.cfg_ana.collectionPostFix+
"_Xj_zll", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
312 if hasattr(event,
'zll_met'):
314 csLeptons = [ l
for l
in event.selectedLeptons
if l.pt() > 10
and abs(l.eta()) < 2.5 ]
316 if len(csLeptons)==2
and len(objectsXjc)>=2:
318 self.zll_mt2_Xj = self.getMT2Hemi(event,objectsXjc,event.zll_met,self.cfg_ana.collectionPostFix,
"_Xj_zll")
349 setattr(event,
"mt2"+self.cfg_ana.collectionPostFix+
"_rl", -999)
350 setattr(event,
"pseudoJet1"+self.cfg_ana.collectionPostFix+
"_rl", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
351 setattr(event,
"pseudoJet2"+self.cfg_ana.collectionPostFix+
"_rl", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
353 if hasattr(event,
'rl_met'):
355 csLeptons = [ l
for l
in event.selectedLeptons
if l.pt() > 10
and abs(l.eta()) < 2.5 ]
357 if len(csLeptons)==1
and len(objects40jc)>=2:
359 self.rl_mt2 = self.getMT2Hemi(event,objects40jc,event.rl_met,self.cfg_ana.collectionPostFix,
"_rl")
361 setattr(event,
"mt2"+self.cfg_ana.collectionPostFix+
"_Xj_rl", -999)
362 setattr(event,
"pseudoJet1"+self.cfg_ana.collectionPostFix+
"_Xj_rl", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
363 setattr(event,
"pseudoJet2"+self.cfg_ana.collectionPostFix+
"_Xj_rl", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
365 if hasattr(event,
'rl_met'):
367 csLeptons = [ l
for l
in event.selectedLeptons
if l.pt() > 10
and abs(l.eta()) < 2.5 ]
369 if len(csLeptons)==1
and len(objectsXjc)>=2:
371 self.rl_mt2_Xj = self.getMT2Hemi(event,objectsXjc,event.rl_met,self.cfg_ana.collectionPostFix,
"_Xj_rl")
376 setattr(event,
"mt2"+self.cfg_ana.collectionPostFix+
"_zllmt", -999)
377 setattr(event,
"pseudoJet1"+self.cfg_ana.collectionPostFix+
"_zllmt", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
378 setattr(event,
"pseudoJet2"+self.cfg_ana.collectionPostFix+
"_zllmt", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
380 if hasattr(event,
'zllmt_met'):
382 csLeptons = [ l
for l
in event.selectedLeptons
if l.pt() > 10
and abs(l.eta()) < 2.5 ]
384 if len(csLeptons)==2
and len(objects40jc)>=2:
386 self.zllmt_mt2 = self.getMT2Hemi(event,objects40jc,event.zllmt_met,self.cfg_ana.collectionPostFix,
"_zllmt")
388 setattr(event,
"mt2"+self.cfg_ana.collectionPostFix+
"_Xj_zllmt", -999)
389 setattr(event,
"pseudoJet1"+self.cfg_ana.collectionPostFix+
"_Xj_zllmt", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
390 setattr(event,
"pseudoJet2"+self.cfg_ana.collectionPostFix+
"_Xj_zllmt", ROOT.reco.Particle.LorentzVector( 0, 0, 0, 0 ))
392 if hasattr(event,
'zllmt_met'):
394 csLeptons = [ l
for l
in event.selectedLeptons
if l.pt() > 10
and abs(l.eta()) < 2.5 ]
396 if len(csLeptons)==2
and len(objectsXjc)>=1:
399 if (event.eventId%2):
400 csLeptons_mt.append(csLeptons[1])
402 csLeptons_mt.append(csLeptons[0])
404 self.zllmt_mt2_Xj = self.getMT2Hemi(event,objectsXjc+[l
for l
in csLeptons
if l
not in csLeptons_mt],event.zllmt_met,self.cfg_ana.collectionPostFix,
"_Xj_zllmt")
408 if len(event.bjetsMedium)>=2:
410 event.mt2bb = self.computeMT2(event.bjetsMedium[0], event.bjetsMedium[1], self.met)
411 event.mt2bb_Xj = self.computeMT2(event.bjetsMedium[0], event.bjetsMedium[1], self.met)
413 if len(event.bjetsMedium)==1:
415 objects40jcCSV = [ j
for j
in event.cleanJets
if j.pt() > 40
and abs(j.eta())<2.5
and j.p4()!=event.bjetsMedium[0].
p4() ]
416 objects40jcCSV.sort(key =
lambda l : l.btag(
'pfCombinedInclusiveSecondaryVertexV2BJetTags'), reverse =
True)
418 objectsXjcCSV = [ j
for j
in event.cleanJets
if j.pt() > self.jetPt
and abs(j.eta())<2.5
and j.p4()!=event.bjetsMedium[0].
p4() ]
419 objectsXjcCSV.sort(key =
lambda l : l.btag(
'pfCombinedInclusiveSecondaryVertexV2BJetTags'), reverse =
True)
421 if len(objects40jcCSV)>0:
422 self.mt2bb = self.computeMT2(event.bjetsMedium[0], objects40jcCSV[0], self.met)
423 setattr(event,
"mt2bb"+self.cfg_ana.collectionPostFix, self.mt2bb)
425 if len(objectsXjcCSV)>0:
426 self.mt2bb_Xj = self.computeMT2(event.bjetsMedium[0], objectsXjcCSV[0], self.met)
427 setattr(event,
"mt2bb_Xj"+self.cfg_ana.collectionPostFix, self.mt2bb_Xj)
432 if not self.cfg_ana.doOnlyDefault:
433 if len(event.selectedLeptons)>=2:
434 self.mt2lep = self.computeMT2(event.selectedLeptons[0], event.selectedLeptons[1], self.met)
435 setattr(event,
"mt2lep"+self.cfg_ana.collectionPostFix, self.mt2lep)
440 self.readCollections( event.input )
446 event.multPseudoJet1_had=0
447 event.multPseudoJet2_had=0
449 event.multPseudoJet1_Xj_had=0
450 event.multPseudoJet2_Xj_had=0
464 setattr(MT2Analyzer,
"defaultConfig", cfg.Analyzer(
465 class_object = MT2Analyzer,
466 metCollection =
"slimmedMETs",
467 collectionPostFix =
"",
468 doOnlyDefault =
True,
Abs< T >::type abs(const T &t)