CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | Public Attributes
JetReCalibrator.JetReCalibrator Class Reference

Public Member Functions

def __init__
 
def correct
 
def correctAll
 

Public Attributes

 JetCorrector
 
 JetUncertainty
 
 L1JetPar
 
 L2JetPar
 
 L3JetPar
 
 ResJetPar
 
 vPar
 

Detailed Description

Definition at line 6 of file JetReCalibrator.py.

Constructor & Destructor Documentation

def JetReCalibrator.JetReCalibrator.__init__ (   self,
  globalTag,
  jetFlavour,
  doResidualJECs,
  jecPath,
  upToLevel = 3 
)
Create a corrector object that reads the payloads from the text dumps of a global tag under
    CMGTools/RootTools/data/jec  (see the getJec.py there to make the dumps).
   It will apply the L1,L2,L3 and possibly the residual corrections to the jets.

Definition at line 7 of file JetReCalibrator.py.

7 
8  def __init__(self,globalTag,jetFlavour,doResidualJECs,jecPath,upToLevel=3):
9  """Create a corrector object that reads the payloads from the text dumps of a global tag under
10  CMGTools/RootTools/data/jec (see the getJec.py there to make the dumps).
11  It will apply the L1,L2,L3 and possibly the residual corrections to the jets."""
12  # Make base corrections
13  path = jecPath #"%s/src/CMGTools/RootTools/data/jec" % os.environ['CMSSW_BASE'];
14  self.L1JetPar = ROOT.JetCorrectorParameters("%s/%s_L1FastJet_%s.txt" % (path,globalTag,jetFlavour),"");
15  self.L2JetPar = ROOT.JetCorrectorParameters("%s/%s_L2Relative_%s.txt" % (path,globalTag,jetFlavour),"");
16  self.L3JetPar = ROOT.JetCorrectorParameters("%s/%s_L3Absolute_%s.txt" % (path,globalTag,jetFlavour),"");
17  self.vPar = ROOT.vector(ROOT.JetCorrectorParameters)()
18  self.vPar.push_back(self.L1JetPar);
19  if upToLevel >= 2: self.vPar.push_back(self.L2JetPar);
20  if upToLevel >= 3: self.vPar.push_back(self.L3JetPar);
21  # Add residuals if needed
22  if doResidualJECs :
23  self.ResJetPar = ROOT.JetCorrectorParameters("%s/%s_L2L3Residual_%s.txt" % (path,globalTag,jetFlavour))
24  self.vPar.push_back(self.ResJetPar);
25  #Step3 (Construct a FactorizedJetCorrector object)
26  self.JetCorrector = ROOT.FactorizedJetCorrector(self.vPar)
27  if os.path.exists("%s/%s_Uncertainty_%s.txt" % (path,globalTag,jetFlavour)):
28  self.JetUncertainty = ROOT.JetCorrectionUncertainty("%s/%s_Uncertainty_%s.txt" % (path,globalTag,jetFlavour));
29  else:
30  print 'Missing JEC uncertainty file "%s/%s_Uncertainty_%s.txt", so jet energy uncertainties will not be available' % (path,globalTag,jetFlavour)
self.JetUncertainty = None

Member Function Documentation

def JetReCalibrator.JetReCalibrator.correct (   self,
  jet,
  rho,
  delta = 0,
  metShift = [0 
)
Corrects a jet energy (optionally shifting it also by delta times the JEC uncertainty)
   If a two-component list is passes as 'metShift', it will be modified adding to the first and second
   component the change to the MET along x and y due to the JEC, defined as the negative difference 
   between the new and old jet 4-vectors, for jets with corrected pt > 10.

Definition at line 41 of file JetReCalibrator.py.

References JetReCalibrator.JetReCalibrator.JetUncertainty, and bookConverter.max.

Referenced by JetReCalibrator.JetReCalibrator.correctAll().

41 
42  def correct(self,jet,rho,delta=0,metShift=[0,0]):
43  """Corrects a jet energy (optionally shifting it also by delta times the JEC uncertainty)
44  If a two-component list is passes as 'metShift', it will be modified adding to the first and second
45  component the change to the MET along x and y due to the JEC, defined as the negative difference
46  between the new and old jet 4-vectors, for jets with corrected pt > 10."""
47  self.JetCorrector.setJetEta(jet.eta())
48  self.JetCorrector.setJetPt(jet.pt() * jet.rawFactor())
49  self.JetCorrector.setJetA(jet.jetArea())
50  self.JetCorrector.setRho(rho)
51  corr = self.JetCorrector.getCorrection()
52  if delta != 0:
53  if not self.JetUncertainty: raise RuntimeError, "Jet energy scale uncertainty shifts requested, but not available"
54  self.JetUncertainty.setJetEta(jet.eta())
55  self.JetUncertainty.setJetPt(corr * jet.pt() * jet.rawFactor())
56  try:
57  jet.jetEnergyCorrUncertainty = self.JetUncertainty.getUncertainty(True)
58  except RuntimeError, r:
59  print "Caught %s when getting uncertainty for jet of pt %.1f, eta %.2f\n" % (r,corr * jet.pt() * jet.rawFactor(),jet.eta())
60  jet.jetEnergyCorrUncertainty = 0.5
61  if jet.photonEnergyFraction() < 0.9 and jet.pt()*corr*jet.rawFactor() > 10:
62  metShift[0] -= jet.px()*(corr*jet.rawFactor() - 1)*(1-jet.muonEnergyFraction())
63  metShift[1] -= jet.py()*(corr*jet.rawFactor() - 1)*(1-jet.muonEnergyFraction())
64  if delta != 0:
65  #print " jet with corr pt %6.2f has an uncertainty %.2f " % (jet.pt()*jet.rawFactor()*corr, jet.jetEnergyCorrUncertainty)
66  corr *= max(0, 1+delta*jet.jetEnergyCorrUncertainty)
67  if jet.pt()*jet.rawFactor()*corr > 10:
68  metShift[0] -= jet.px()*jet.rawFactor()*corr*delta*jet.jetEnergyCorrUncertainty
69  metShift[1] -= jet.py()*jet.rawFactor()*corr*delta*jet.jetEnergyCorrUncertainty
70  #print " jet with raw pt %6.2f eta %+5.3f phi %+5.3f: previous corr %.4f, my corr %.4f " % (jet.pt()*jet.rawFactor(), jet.eta(), jet.phi(), 1./jet.rawFactor(), corr)
71  if corr <= 0:
72  return False
73  jet.setP4(jet.p4() * (corr * jet.rawFactor()))
74  jet.setRawFactor(1.0/corr)
75  return True
def JetReCalibrator.JetReCalibrator.correctAll (   self,
  jets,
  rho,
  delta = 0,
  metShift = [0 
)
Applies 'correct' to all the jets, discard the ones that have bad corrections (corrected pt <= 0)

Definition at line 31 of file JetReCalibrator.py.

References ElectronCalibrator.EmbeddedElectronCalibrator.correct(), FFTJetCorrector< Jet, Adjustable >.correct(), JetReCalibrator.JetReCalibrator.correct(), FFTJetCorrectorSequence< Jet, InitialConverter, FinalConverter >.correct(), and MomentumScaleCorrector.correct().

31 
32  def correctAll(self,jets,rho,delta=0,metShift=[0,0]):
33  """Applies 'correct' to all the jets, discard the ones that have bad corrections (corrected pt <= 0)"""
34  badJets = []
35  for j in jets:
36  ok = self.correct(j,rho,delta,metShift)
37  if not ok: badJets.append(j)
38  if len(badJets) > 0:
39  print "Warning: %d out of %d jets flagged bad by JEC." % (len(badJets), len(jets))
40  for bj in badJets:
jets.remove(bj)

Member Data Documentation

JetReCalibrator.JetReCalibrator.JetCorrector

Definition at line 25 of file JetReCalibrator.py.

JetReCalibrator.JetReCalibrator.JetUncertainty

Definition at line 27 of file JetReCalibrator.py.

Referenced by JetReCalibrator.JetReCalibrator.correct().

JetReCalibrator.JetReCalibrator.L1JetPar

Definition at line 13 of file JetReCalibrator.py.

JetReCalibrator.JetReCalibrator.L2JetPar

Definition at line 14 of file JetReCalibrator.py.

JetReCalibrator.JetReCalibrator.L3JetPar

Definition at line 15 of file JetReCalibrator.py.

JetReCalibrator.JetReCalibrator.ResJetPar

Definition at line 22 of file JetReCalibrator.py.

JetReCalibrator.JetReCalibrator.vPar

Definition at line 16 of file JetReCalibrator.py.