CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
METAnalyzer.py
Go to the documentation of this file.
1 import random
2 from PhysicsTools.Heppy.analyzers.core.Analyzer import Analyzer
3 from PhysicsTools.Heppy.analyzers.core.AutoHandle import AutoHandle
4 from PhysicsTools.Heppy.physicsobjects.PhysicsObjects import Jet
6 from PhysicsTools.HeppyCore.statistics.counter import Counter, Counters
7 from PhysicsTools.Heppy.physicsutils.JetReCalibrator import JetReCalibrator
9 
10 import operator
11 import itertools
12 import copy
13 from ROOT import TLorentzVector, TVectorD
14 import ROOT
15 import math
16 
17 from copy import deepcopy
18 
19 class METAnalyzer( Analyzer ):
20  def __init__(self, cfg_ana, cfg_comp, looperName ):
21  super(METAnalyzer,self).__init__(cfg_ana,cfg_comp,looperName)
22 
23  def declareHandles(self):
24  super(METAnalyzer, self).declareHandles()
25  self.handles['met'] = AutoHandle( self.cfg_ana.metCollection, 'std::vector<pat::MET>' )
26  self.handles['nopumet'] = AutoHandle( self.cfg_ana.noPUMetCollection, 'std::vector<pat::MET>' )
27  self.handles['cmgCand'] = AutoHandle( self.cfg_ana.candidates, self.cfg_ana.candidatesTypes )
28  self.handles['vertices'] = AutoHandle( "offlineSlimmedPrimaryVertices", 'std::vector<reco::Vertex>', fallbackLabel="offlinePrimaryVertices" )
29  self.mchandles['packedGen'] = AutoHandle( 'packedGenParticles', 'std::vector<pat::PackedGenParticle>' )
30 
31  def beginLoop(self, setup):
32  super(METAnalyzer,self).beginLoop(setup)
33  self.counters.addCounter('events')
34  count = self.counters.counter('events')
35  count.register('all events')
36 
37  def applyDeltaMet(self, met, deltaMet):
38  px,py = self.met.px()+deltaMet[0], self.met.py()+deltaMet[1]
39  met.setP4(ROOT.reco.Particle.LorentzVector(px,py, 0, math.hypot(px,py)))
40 
41  def makeTkMETs(self, event):
42 
43  charged = []
44  chargedchs = []
45  chargedPVLoose = []
46  chargedPVTight = []
47 
48  pfcands = self.handles['cmgCand'].product()
49 
50  for i in xrange(pfcands.size()):
51 
52 ## ===> require the Track Candidate charge and with a minimum dz
53 
54  if (pfcands.at(i).charge()!=0):
55 
56  if abs(pfcands.at(i).dz())<=self.cfg_ana.dzMax:
57  charged.append(pfcands.at(i))
58 
59  if pfcands.at(i).fromPV()>0:
60  chargedchs.append(pfcands.at(i))
61 
62  if pfcands.at(i).fromPV()>1:
63  chargedPVLoose.append(pfcands.at(i))
64 
65  if pfcands.at(i).fromPV()>2:
66  chargedPVTight.append(pfcands.at(i))
67 
68  import ROOT
69  setattr(event, "tkMet"+self.cfg_ana.collectionPostFix, \
70  ROOT.reco.Particle.LorentzVector(-1.*(sum([x.px() for x in charged])) , -1.*(sum([x.py() for x in charged])), 0, math.hypot((sum([x.px() for x in charged])),(sum([x.py() for x in charged]))) ))
71  setattr(event, "tkMetPVchs"+self.cfg_ana.collectionPostFix, \
72  ROOT.reco.Particle.LorentzVector(-1.*(sum([x.px() for x in chargedchs])) , -1.*(sum([x.py() for x in chargedchs])), 0, math.hypot((sum([x.px() for x in chargedchs])),(sum([x.py() for x in chargedchs]))) ))
73  setattr(event, "tkMetPVLoose"+self.cfg_ana.collectionPostFix, \
74  ROOT.reco.Particle.LorentzVector(-1.*(sum([x.px() for x in chargedPVLoose])) , -1.*(sum([x.py() for x in chargedPVLoose])), 0, math.hypot((sum([x.px() for x in chargedPVLoose])),(sum([x.py() for x in chargedPVLoose]))) ))
75  setattr(event, "tkMetPVTight"+self.cfg_ana.collectionPostFix, \
76  ROOT.reco.Particle.LorentzVector(-1.*(sum([x.px() for x in chargedPVTight])) , -1.*(sum([x.py() for x in chargedPVTight])), 0, math.hypot((sum([x.px() for x in chargedPVTight])),(sum([x.py() for x in chargedPVTight]))) ))
77 ## print 'tkmet',self.tkMet.pt(),'tkmetphi',self.tkMet.phi()
78 
79 
80  event.tkSumEt = sum([x.pt() for x in charged])
81  event.tkPVchsSumEt = sum([x.pt() for x in chargedchs])
82  event.tkPVLooseSumEt = sum([x.pt() for x in chargedPVLoose])
83  event.tkPVTightSumEt = sum([x.pt() for x in chargedPVTight])
84 
85 
86  def makeGenTkMet(self, event):
87  genCharged = [ x for x in self.mchandles['packedGen'].product() if x.charge() != 0 and abs(x.eta()) < 2.4 ]
88  event.tkGenMet = ROOT.reco.Particle.LorentzVector(-1.*(sum([x.px() for x in genCharged])) , -1.*(sum([x.py() for x in genCharged])), 0, math.hypot((sum([x.px() for x in genCharged])),(sum([x.py() for x in genCharged]))) )
89 
90  def makeMETNoMu(self, event):
91  self.metNoMu = copy.deepcopy(self.met)
92  if self.cfg_ana.doMetNoPU: self.metNoMuNoPU = copy.deepcopy(self.metNoPU)
93 
94  mupx = 0
95  mupy = 0
96  #sum muon momentum
97  for mu in event.selectedMuons:
98  mupx += mu.px()
99  mupy += mu.py()
100 
101  #subtract muon momentum and construct met
102  px,py = self.metNoMu.px()+mupx, self.metNoMu.py()+mupy
103  self.metNoMu.setP4(ROOT.reco.Particle.LorentzVector(px,py, 0, math.hypot(px,py)))
104  px,py = self.metNoMuNoPU.px()+mupx, self.metNoMuNoPU.py()+mupy
105  self.metNoMuNoPU.setP4(ROOT.reco.Particle.LorentzVector(px,py, 0, math.hypot(px,py)))
106  setattr(event, "metNoMu"+self.cfg_ana.collectionPostFix, self.metNoMu)
107  if self.cfg_ana.doMetNoPU: setattr(event, "metNoMuNoPU"+self.cfg_ana.collectionPostFix, self.metNoMuNoPU)
108 
109 
110  def makeMETNoEle(self, event):
111  self.metNoEle = copy.deepcopy(self.met)
112  if self.cfg_ana.doMetNoPU: self.metNoEleNoPU = copy.deepcopy(self.metNoPU)
113 
114  elepx = 0
115  elepy = 0
116  #sum electron momentum
117  for ele in event.selectedElectrons:
118  elepx += ele.px()
119  elepy += ele.py()
120 
121  #subtract electron momentum and construct met
122  px,py = self.metNoEle.px()+elepx, self.metNoEle.py()+elepy
123  self.metNoEle.setP4(ROOT.reco.Particle.LorentzVector(px,py, 0, math.hypot(px,py)))
124 
125  px,py = self.metNoEleNoPU.px()+elepx, self.metNoEleNoPU.py()+elepy
126  self.metNoEleNoPU.setP4(ROOT.reco.Particle.LorentzVector(px,py, 0, math.hypot(px,py)))
127  setattr(event, "metNoEle"+self.cfg_ana.collectionPostFix, self.metNoEle)
128  if self.cfg_ana.doMetNoPU: setattr(event, "metNoEleNoPU"+self.cfg_ana.collectionPostFix, self.metNoEleNoPU)
129 
130  def makeMETNoPhoton(self, event):
131  self.metNoPhoton = copy.deepcopy(self.met)
132 
133  phopx = 0
134  phopy = 0
135  #sum photon momentum
136  for pho in event.selectedPhotons:
137  phopx += pho.px()
138  phopy += pho.py()
139 
140  #subtract photon momentum and construct met
141  px,py = self.metNoPhoton.px()+phopx, self.metNoPhoton.py()+phopy
142  self.metNoPhoton.setP4(ROOT.reco.Particle.LorentzVector(px,py, 0, math.hypot(px,py)))
143  setattr(event, "metNoPhoton"+self.cfg_ana.collectionPostFix, self.metNoPhoton)
144  if self.cfg_ana.doMetNoPU:
145  self.metNoPhotonNoPU = copy.deepcopy(self.metNoPU)
146  px,py = self.metNoPhotonNoPU.px()+phopx, self.metNoPhotonNoPU.py()+phopy
147  self.metNoPhotonNoPU.setP4(ROOT.reco.Particle.LorentzVector(px,py, 0, math.hypot(px,py)))
148  setattr(event, "metNoPhotonNoPU"+self.cfg_ana.collectionPostFix, self.metNoPhotonNoPU)
149 
150 
151  def makeMETs(self, event):
152  import ROOT
153  if self.cfg_ana.copyMETsByValue:
154  self.met = ROOT.pat.MET(self.handles['met'].product()[0])
155  if self.cfg_ana.doMetNoPU: self.metNoPU = ROOT.pat.MET(self.handles['nopumet'].product()[0])
156  else:
157  self.met = self.handles['met'].product()[0]
158  if self.cfg_ana.doMetNoPU: self.metNoPU = self.handles['nopumet'].product()[0]
159 
160  #Shifted METs
161  #Uncertainties defined in https://github.com/cms-sw/cmssw/blob/CMSSW_7_2_X/DataFormats/PatCandidates/interface/MET.h#L168
162  #event.met_shifted = []
163  if not self.cfg_ana.copyMETsByValue:
164  for i in range(self.met.METUncertaintySize):
165  m = ROOT.pat.MET(self.met)
166  px = m.shiftedPx(i);
167  py = m.shiftedPy(i);
168  m.setP4(ROOT.reco.Particle.LorentzVector(px,py, 0, math.hypot(px,py)))
169  #event.met_shifted += [m]
170  setattr(event, "met{0}_shifted_{1}".format(self.cfg_ana.collectionPostFix, i), m)
171 
172  self.met_sig = self.met.significance()
173  self.met_sumet = self.met.sumEt()
174 
175 
176  ###https://github.com/cms-sw/cmssw/blob/CMSSW_7_2_X/DataFormats/PatCandidates/interface/MET.h
177  if not self.cfg_ana.copyMETsByValue:
178  self.metraw = self.met.shiftedPt(12, 0)
179  self.metType1chs = self.met.shiftedPt(12, 1)
180  setattr(event, "metraw"+self.cfg_ana.collectionPostFix, self.metraw)
181  setattr(event, "metType1chs"+self.cfg_ana.collectionPostFix, self.metType1chs)
182 
183  if self.cfg_ana.recalibrate and hasattr(event, 'deltaMetFromJetSmearing'+self.cfg_ana.jetAnalyzerCalibrationPostFix):
184  deltaMetSmear = getattr(event, 'deltaMetFromJetSmearing'+self.cfg_ana.jetAnalyzerCalibrationPostFix)
185  self.applyDeltaMet(self.met, deltaMetSmear)
186  if self.cfg_ana.doMetNoPU:
187  self.applyDeltaMet(self.metNoPU, deltaMetSmear)
188  if self.cfg_ana.recalibrate and hasattr(event, 'deltaMetFromJEC'+self.cfg_ana.jetAnalyzerCalibrationPostFix):
189  deltaMetJEC = getattr(event, 'deltaMetFromJEC'+self.cfg_ana.jetAnalyzerCalibrationPostFix)
190 # print 'before JEC', self.cfg_ana.collectionPostFix, self.met.px(),self.met.py(), 'deltaMetFromJEC'+self.cfg_ana.jetAnalyzerCalibrationPostFix, deltaMetJEC
191  self.applyDeltaMet(self.met, deltaMetJEC)
192  if self.cfg_ana.doMetNoPU:
193  self.applyDeltaMet(self.metNoPU, deltaMetJEC)
194 # print 'after JEC', self.cfg_ana.collectionPostFix, self.met.px(),self.met.py(), 'deltaMetFromJEC'+self.cfg_ana.jetAnalyzerCalibrationPostFix, deltaMetJEC
195 
196  setattr(event, "met"+self.cfg_ana.collectionPostFix, self.met)
197  if self.cfg_ana.doMetNoPU: setattr(event, "metNoPU"+self.cfg_ana.collectionPostFix, self.metNoPU)
198  setattr(event, "met_sig"+self.cfg_ana.collectionPostFix, self.met_sig)
199  setattr(event, "met_sumet"+self.cfg_ana.collectionPostFix, self.met_sumet)
200 
201  if self.cfg_ana.doMetNoMu and hasattr(event, 'selectedMuons'):
202  self.makeMETNoMu(event)
203 
204  if self.cfg_ana.doMetNoEle and hasattr(event, 'selectedElectrons'):
205  self.makeMETNoEle(event)
206 
207  if self.cfg_ana.doMetNoPhoton and hasattr(event, 'selectedPhotons'):
208  self.makeMETNoPhoton(event)
209 
210  def process(self, event):
211  self.readCollections( event.input)
212  self.counters.counter('events').inc('all events')
213 
214  self.makeMETs(event)
215 
216  if self.cfg_ana.doTkMet:
217  self.makeTkMETs(event);
218 
219  if self.cfg_comp.isMC and hasattr(event, 'genParticles'):
220  self.makeGenTkMet(event)
221 
222  return True
223 
224 
225 setattr(METAnalyzer,"defaultConfig", cfg.Analyzer(
226  class_object = METAnalyzer,
227  metCollection = "slimmedMETs",
228  noPUMetCollection = "slimmedMETs",
229  copyMETsByValue = False,
230  recalibrate = True,
231  jetAnalyzerCalibrationPostFix = "",
232  doTkMet = False,
233  doMetNoPU = True,
234  doMetNoMu = False,
235  doMetNoEle = False,
236  doMetNoPhoton = False,
237  candidates='packedPFCandidates',
238  candidatesTypes='std::vector<pat::PackedCandidate>',
239  dzMax = 0.1,
240  collectionPostFix = "",
241  )
242 )
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
def makeGenTkMet
print &#39;tkmet&#39;,self.tkMet.pt(),&#39;tkmetphi&#39;,self.tkMet.phi()
Definition: METAnalyzer.py:86