add features and datafiles to run 5TeV

python/postprocessing/framework/eventloop.py

@@ -47,7 +47,7 @@ def eventLoop(modules, inputFile, outputFile, inputTree, wrappedOutputTree, maxE
     for m in modules: 
         m.beginFile(inputFile, outputFile, inputTree, wrappedOutputTree)
 
-    t0 = time.clock(); tlast = t0; doneEvents = 0; acceptedEvents = 0
+    t0 = time.time(); tlast = t0; doneEvents = 0; acceptedEvents = 0
     entries = inputTree.entries
     if eventRange: entries = len(eventRange)
     if maxEvents > 0: entries = min(entries, maxEvents)
@@ -67,11 +67,11 @@ def eventLoop(modules, inputFile, outputFile, inputTree, wrappedOutputTree, maxE
             wrappedOutputTree.fill()
         if progress:
             if ie > 0 and ie % progress[0] == 0:
-                t1 = time.clock()
+                t1 = time.time()
                 progress[1].write("Processed %8d/%8d entries, %5.2f%% (elapsed time %7.1fs, curr speed %8.3f kHz, avg speed %8.3f kHz), accepted %8d/%8d events (%5.2f%%)\n" % (
                         ie,entries, ie/float(0.01*entries), t1-t0, (progress[0]/1000.)/(max(t1-tlast,1e-9)), ie/1000./(max(t1-t0,1e-9)), acceptedEvents, doneEvents, acceptedEvents/(0.01*doneEvents) ))
                 tlast = t1
     for m in modules: 
         m.endFile(inputFile, outputFile, inputTree, wrappedOutputTree)
 
-    return (doneEvents, acceptedEvents, time.clock() - t0)
+    return (doneEvents, acceptedEvents, time.time() - t0)

python/postprocessing/modules/btv/btagSFProducer.py

@@ -1,45 +1,45 @@
 import ROOT
 import os
+from itertools import chain
+
 ROOT.PyConfig.IgnoreCommandLineOptions = True
 
 from PhysicsTools.NanoAODTools.postprocessing.framework.datamodel import Collection
 from PhysicsTools.NanoAODTools.postprocessing.framework.eventloop import Module
 
-def is_relevant_syst_for_shape_corr(flavor_btv, syst):
+def is_relevant_syst_for_shape_corr(flavor_btv, syst, jesSystsForShape=["jes"]):
     """Returns true if a flavor/syst combination is relevant"""
+    jesSysts = list(chain(*[ ("up_" + j, "down_" + j) for j in jesSystsForShape ]))
+
     if flavor_btv == 0:
         return syst in [ "central",
-                         "up_jes", "down_jes",
                          "up_lf", "down_lf",
                          "up_hfstats1", "down_hfstats1",
-                         "up_hfstats2", "down_hfstats2" ]
+                         "up_hfstats2", "down_hfstats2" ] + jesSysts
     elif flavor_btv == 1:
         return syst in [ "central",
                          "up_cferr1", "down_cferr1",
                          "up_cferr2", "down_cferr2" ]
     elif flavor_btv == 2:
         return syst in [ "central",
-                         "up_jes", "down_jes",
                          "up_hf", "down_hf",
                          "up_lfstats1", "down_lfstats1",
-                         "up_lfstats2", "down_lfstats2" ]
+                         "up_lfstats2", "down_lfstats2" ] + jesSysts
     else:
         raise ValueError("ERROR: Undefined flavor = %i!!" % flavor_btv)
     return True
 
 class btagSFProducer(Module):
     """Calculate btagging scale factors
-        algo has to be either 'csvv2' or 'cmva'
     """
-    def __init__(self, era, algo = 'csvv2', sfFileName = None, verbose = 0, collName='Jet',storeOutput=True):
+    def __init__(self, era, algo='csvv2', selectedWPs=['M', 'shape_corr'], sfFileName=None, verbose=0, jesSystsForShape=["jes"]):
 
         self.era = era
-
         self.algo = algo.lower()
-
+        self.selectedWPs = selectedWPs
         self.verbose = verbose
-        self.collName = collName
-        self.storeOutput = storeOutput
+        self.jesSystsForShape = jesSystsForShape
+
         # CV: Return value of BTagCalibrationReader::eval_auto_bounds() is zero
         #     in case jet abs(eta) > 2.4 !!
         self.max_abs_eta = 2.4
@@ -96,8 +96,8 @@ def __init__(self, era, algo = 'csvv2', sfFileName = None, verbose = 0, collName
                         1 : "comb",  # c
                         2 : "incl"   # light
                     },
-                    'supported_wp' : [ "L", "M", "T", "shape_corr"] 
-                }    
+                    'supported_wp' : [ "L", "M", "T", "shape_corr"]
+                }
             },
             'deepjet' : {
                 'Legacy2016' : {
@@ -126,7 +126,7 @@ def __init__(self, era, algo = 'csvv2', sfFileName = None, verbose = 0, collName
                         2 : "incl"   # light
                     },
                     'supported_wp' : [ "L", "M", "T", "shape_corr"]
-                }    
+                }
             },
             'cmva' : {
                 '2016' : {
@@ -143,7 +143,7 @@ def __init__(self, era, algo = 'csvv2', sfFileName = None, verbose = 0, collName
 
         supported_algos = []
         for algo in supported_btagSF.keys():
-            if era in supported_btagSF[algo].keys():
+            if self.era in supported_btagSF[algo].keys():
                 supported_algos.append(algo)
         if self.algo in supported_btagSF.keys():
             if self.era in supported_btagSF[self.algo].keys():
@@ -155,19 +155,9 @@ def __init__(self, era, algo = 'csvv2', sfFileName = None, verbose = 0, collName
                 raise ValueError("ERROR: Algorithm '%s' not supported for era = '%s'! Please choose among { %s }." % (self.algo, self.era, supported_algos))
         else:
             raise ValueError("ERROR: Algorithm '%s' not supported for era = '%s'! Please choose among { %s }." % (self.algo, self.era, supported_algos))
-
-        algoLabel = None
-        if self.algo == "csvv2":
-            algoLabel = "CSV (v2)"
-        elif self.algo == "deepcsv":
-            algoLabel = "deep-CSV (b)"
-        elif self.algo == "cmva":
-            algoLabel = "cMVA"
-        elif self.algo == "deepjet":
-            algoLabel = "DeepJet"
-        else:
-            raise ValueError("ERROR: Algorithm '%s' not supported for era = '%s'! Please choose among { %s }." % (self.algo, self.era, supported_algos))
-        print("Loading btagSF weights for %s algorithm from file '%s'" % (algoLabel, os.path.join(self.inputFilePath, self.inputFileName)))
+        for wp in self.selectedWPs:
+            if wp not in self.supported_wp:
+                raise ValueError("ERROR: Working point '%s' not supported for algo = '%s' and era = '%s'! Please choose among { %s }." % (wp, self.algo, self.era, self.supported_wp))
 
         # load libraries for accessing b-tag scale factors (SFs) from conditions database
         for library in [ "libCondFormatsBTauObjects", "libCondToolsBTau" ]:
@@ -183,36 +173,38 @@ def __init__(self, era, algo = 'csvv2', sfFileName = None, verbose = 0, collName
         self.central_and_systs.extend(self.systs)
 
         self.systs_shape_corr = []
-        for syst in [ 'jes',
-                      'lf', 'hf',
+        for syst in [ 'lf', 'hf',
                       'hfstats1', 'hfstats2',
                       'lfstats1', 'lfstats2',
-                      'cferr1', 'cferr2' ]:
+                      'cferr1', 'cferr2' ] + self.jesSystsForShape:
             self.systs_shape_corr.append("up_%s" % syst)
             self.systs_shape_corr.append("down_%s" % syst)
         self.central_and_systs_shape_corr = [ "central" ]
         self.central_and_systs_shape_corr.extend(self.systs_shape_corr)
 
         self.branchNames_central_and_systs = {}
-        for central_or_syst in self.central_and_systs:
-            if central_or_syst == "central":
-                self.branchNames_central_and_systs[central_or_syst] = "%s_btagSF"%self.collName
+        for wp in self.selectedWPs:
+            branchNames = {}
+            if wp == 'shape_corr':
+                central_and_systs = self.central_and_systs_shape_corr
+                baseBranchName = 'Jet_btagSF_{}_shape'.format(self.algo)
             else:
-                self.branchNames_central_and_systs[central_or_syst] = "%s_btagSF_%s" %(self.collName, central_or_syst)
+                central_and_systs = self.central_and_systs
+                baseBranchName = 'Jet_btagSF_{}_{}'.format(self.algo, wp)
+            for central_or_syst in central_and_systs:
+                if central_or_syst == "central":
+                    branchNames[central_or_syst] = baseBranchName
+                else:
+                    branchNames[central_or_syst] = baseBranchName + '_' + central_or_syst
+            self.branchNames_central_and_systs[wp] = branchNames
 
-        self.branchNames_central_and_systs_shape_corr = {}
-        for central_or_syst in self.central_and_systs_shape_corr:
-            if central_or_syst == "central":
-                self.branchNames_central_and_systs_shape_corr[central_or_syst] = "%s_btagSF_shape"%self.collName
-            else:
-                self.branchNames_central_and_systs_shape_corr[central_or_syst] = "%s_btagSF_shape_%s" % (self.collName,central_or_syst)
 
     def beginJob(self):
         # initialize BTagCalibrationReader
         # (cf. https://twiki.cern.ch/twiki/bin/viewauth/CMS/BTagCalibration )
         self.calibration = ROOT.BTagCalibration(self.algo, os.path.join(self.inputFilePath, self.inputFileName))
         self.readers = {}
-        for wp in self.supported_wp:
+        for wp in self.selectedWPs:
             wp_btv = { "l" : 0, "m" : 1, "t" : 2, "shape_corr" : 3 }.get(wp.lower(), None)
             syts = None
             if wp_btv in [ 0, 1, 2 ]:
@@ -235,22 +227,18 @@ def endJob(self):
 
     def beginFile(self, inputFile, outputFile, inputTree, wrappedOutputTree):
         self.out = wrappedOutputTree
-        if self.storeOutput:
-            for central_or_syst in self.central_and_systs:
-                self.out.branch(self.branchNames_central_and_systs[central_or_syst], "F", lenVar="n%s"%self.collName)
-            for central_or_syst in self.central_and_systs_shape_corr:
-                self.out.branch(self.branchNames_central_and_systs_shape_corr[central_or_syst], "F", lenVar="n%s"%self.collName)
+        for central_or_syst in self.branchNames_central_and_systs.values():
+            for branch in central_or_syst.values():
+                self.out.branch(branch, "F", lenVar="nJet")
 
     def endFile(self, inputFile, outputFile, inputTree, wrappedOutputTree):
         pass
 
-    def getReader(self, wp, shape_corr = False):
+    def getReader(self, wp):
         """
             Get btag scale factor reader.
             Convert working points: input is 'L', 'M', 'T', 'shape_corr' to 0, 1, 2, 3
         """
-        if shape_corr:
-            wp = "shape_corr"
         wp_btv = { "l" : 0, "m" : 1, "t" : 2, "shape_corr" : 3 }.get(wp.lower(), None)
         if wp_btv == None or not wp_btv in self.readers.keys():
             if self.verbose > 0:
@@ -294,7 +282,7 @@ def getSFs(self, jet_data, syst, reader,  measurement_type = 'auto', shape_corr
             # evaluate SF
             sf = None
             if shape_corr:
-                if is_relevant_syst_for_shape_corr(flavor_btv, syst):
+                if is_relevant_syst_for_shape_corr(flavor_btv, syst, self.jesSystsForShape):
                     sf = reader.eval_auto_bounds(syst, flavor_btv, eta, pt, discr)
                 else:
                     sf = reader.eval_auto_bounds('central', flavor_btv, eta, pt, discr)
@@ -322,26 +310,16 @@ def analyze(self, event):
         elif self.algo == "deepjet":
             discr = "btagDeepFlavB"
         else:
-            raise ValueError("ERROR: Invalid algorithm '%s'! Please choose either 'csvv2' or 'cmva'." % self.algo)
+            raise ValueError("ERROR: Invalid algorithm '%s'!" % self.algo)
 
         preloaded_jets = [(jet.pt, jet.eta, self.getFlavorBTV(jet.hadronFlavour), getattr(jet, discr)) for jet in jets]
-        reader = self.getReader('M', False)
-        for central_or_syst in self.central_and_systs:
-            central_or_syst = central_or_syst.lower()
-            scale_factors = list(self.getSFs(preloaded_jets, central_or_syst, reader, 'auto', False))
-            if self.storeOutput: 
-                self.out.fillBranch(self.branchNames_central_and_systs[central_or_syst], scale_factors)
-            else: 
-                setattr(event, self.branchNames_central_and_systs[central_or_syst], scale_factors)
-        # shape corrections
-        reader = self.getReader('shape_corr', True)
-        for central_or_syst in self.central_and_systs_shape_corr:
-            central_or_syst = central_or_syst.lower()
-            scale_factors = list(self.getSFs(preloaded_jets, central_or_syst, reader, 'auto', True))
-            if self.storeOutput:
-                self.out.fillBranch(self.branchNames_central_and_systs_shape_corr[central_or_syst], scale_factors)
-            else: 
-                setattr(event, self.branchNames_central_and_systs_shape_corr[central_or_syst], scale_factors)
+        for wp in self.selectedWPs:
+            reader = self.getReader(wp)
+            isShape = ( wp == 'shape_corr' )
+            central_and_systs = ( self.central_and_systs_shape_corr if isShape else self.central_and_systs )
+            for central_or_syst in central_and_systs:
+                scale_factors = list(self.getSFs(preloaded_jets, central_or_syst, reader, 'auto', isShape))
+                self.out.fillBranch(self.branchNames_central_and_systs[wp][central_or_syst], scale_factors)
 
         return True
 

python/postprocessing/modules/common/PrefireCorr.py

@@ -4,12 +4,29 @@
 import re
 import os
 
-from PhysicsTools.NanoAODTools.postprocessing.framework.datamodel import Collection 
+from PhysicsTools.NanoAODTools.postprocessing.framework.datamodel import Collection
 from PhysicsTools.NanoAODTools.postprocessing.framework.eventloop import Module
 
 class PrefCorr(Module):
-    def __init__(self, jetroot="L1prefiring_jetpt_2017BtoF.root", jetmapname="L1prefiring_jetpt_2017BtoF",
-                 photonroot="L1prefiring_photonpt_2017BtoF.root", photonmapname="L1prefiring_photonpt_2017BtoF"):
+    def __init__(self,
+                 jetroot="L1prefiring_jetpt_2017BtoF.root",
+                 jetmapname="L1prefiring_jetpt_2017BtoF",
+                 photonroot="L1prefiring_photonpt_2017BtoF.root",
+                 photonmapname="L1prefiring_photonpt_2017BtoF",
+                 branchnames=["PrefireWeight","PrefireWeight_Up", "PrefireWeight_Down"]):
+        """Module to compute prefiring weights
+
+        :param jetroot: Root file containing prefiring map for jets, defaults to "L1prefiring_jetpt_2017BtoF.root"
+        :type jetroot: str, optional
+        :param jetmapname: Name of jet prefiring map in ROOT file, defaults to "L1prefiring_jetpt_2017BtoF"
+        :type jetmapname: str, optional
+        :param photonroot: ROOT file containing prefiring map for photons, defaults to "L1prefiring_photonpt_2017BtoF.root"
+        :type photonroot: str, optional
+        :param photonmapname: Name of photon prefiring map in ROOT file, defaults to "L1prefiring_photonpt_2017BtoF"
+        :type photonmapname: str, optional
+        :param branchnames: Output branch names for nominal, up, down variations, defaults to ["PrefireWeight","PrefireWeight_Up", "PrefireWeight_Down"]
+        :type branchnames: list, optional
+        """
 
         cmssw_base = os.getenv('CMSSW_BASE')
 
@@ -20,7 +37,7 @@ def __init__(self, jetroot="L1prefiring_jetpt_2017BtoF.root", jetmapname="L1pref
         self.jet_map = self.get_root_obj(self.jet_file, jetmapname)
 
         self.UseEMpT = ("jetempt" in jetroot)
-
+        self.branchnames = branchnames
     def open_root(self, path):
         r_file = ROOT.TFile.Open(path)
         if not r_file.__nonzero__() or not r_file.IsOpen(): raise NameError('File ' + path + ' not open')
@@ -39,7 +56,6 @@ def endJob(self):
 
     def beginFile(self, inputFile, outputFile, inputTree, wrappedOutputTree):
         self.out = wrappedOutputTree
-        self.branchnames = ["PrefireWeight", "PrefireWeight_Up", "PrefireWeight_Down"]
         for bname in self.branchnames:
           self.out.branch(bname, "F")
 
@@ -115,9 +131,11 @@ def GetPrefireProbability(self, Map, eta, pt, maxpt):
       stat = Map.GetBinError(bin) # bin statistical uncertainty
       syst = 0.2*pref_prob # 20% of prefire rate
 
-      if self.variation == 1: 
+      if self.variation == 1:
         pref_prob = min(pref_prob + math.sqrt(stat*stat + syst*syst), 1.0)
       if self.variation == -1:
         pref_prob = max(pref_prob - math.sqrt(stat*stat + syst*syst), 0.0)
       return pref_prob
 
+prefCorr_2017G = lambda : PrefCorr(jetroot="L1prefiring_jetpt_2017G.root", jetmapname="L1prefiring_jetpt_2017G",
+                                   photonroot="L1prefiring_photonpt_2017G.root", photonmapname="L1prefiring_photonpt_2017G")

python/postprocessing/modules/common/muonScaleResProducer.py

@@ -83,6 +83,7 @@ def analyze(self, event):
         return True
 
 
-muonScaleRes2016 = lambda : muonScaleResProducer('roccor.Run2.v3', 'RoccoR2016.txt', 2016)
-muonScaleRes2017 = lambda : muonScaleResProducer('roccor.Run2.v3', 'RoccoR2017.txt', 2017)
-muonScaleRes2018 = lambda : muonScaleResProducer('roccor.Run2.v3', 'RoccoR2018.txt', 2018)
+muonScaleRes2016  = lambda : muonScaleResProducer('roccor.Run2.v3', 'RoccoR2016.txt', 2016)
+muonScaleRes2017  = lambda : muonScaleResProducer('roccor.Run2.v3', 'RoccoR2017.txt', 2017)
+muonScaleRes2017G = lambda : muonScaleResProducer('roccor.Run2.v3', 'RoccoR2017.txt', 2017)
+muonScaleRes2018  = lambda : muonScaleResProducer('roccor.Run2.v3', 'RoccoR2018.txt', 2018)

python/postprocessing/modules/common/puWeightProducer.py

@@ -100,3 +100,8 @@ def analyze(self, event):
 pufile_mc2018="%s/src/PhysicsTools/NanoAODTools/python/postprocessing/data/pileup/mcPileup2018.root" % os.environ['CMSSW_BASE']
 puWeight_2018 = lambda : puWeightProducer(pufile_mc2018,pufile_data2018,"pu_mc","pileup",verbose=False, doSysVar=True)
 puAutoWeight_2018 = lambda : puWeightProducer("auto",pufile_data2018,"pu_mc","pileup",verbose=False)
+
+pufile_data2017G="%s/src/PhysicsTools/NanoAODTools/python/postprocessing/data/pileup/pileup_5TeV.root" % os.environ['CMSSW_BASE']
+pufile_mc2017G="%s/src/PhysicsTools/NanoAODTools/python/postprocessing/data/pileup/PileupMC_5TeV.root" % os.environ['CMSSW_BASE']
+puWeight_2017G = lambda : puWeightProducer(pufile_mc2017G,pufile_data2017G,"pu_mc","pileup",verbose=False, doSysVar=True, nvtx_var='Pileup_nPU')
+puAutoWeight_2017G = lambda : puWeightProducer("auto",pufile_data2017G,"pu_mc","pileup",verbose=False, nvtx_var='Pileup_nPU')

python/postprocessing/modules/jme/jetmetHelperRun2.py

@@ -7,7 +7,8 @@
 
 # JEC dict
 jecTagsMC = {'2016' : 'Summer16_07Aug2017_V11_MC', 
-             '2017' : 'Fall17_17Nov2017_V32_MC', 
+#             '2017' : 'Fall17_17Nov2017_V32_MC', 
+             '2017' : 'Spring18_ppRef5TeV_V4_MC', 
              '2018' : 'Autumn18_V19_MC'}
 
 jecTagsFastSim = {'2016' : 'Summer16_FastSimV1_MC',
@@ -31,14 +32,16 @@
                 '2017D' : 'Fall17_17Nov2017DE_V32_DATA', 
                 '2017E' : 'Fall17_17Nov2017DE_V32_DATA', 
                 '2017F' : 'Fall17_17Nov2017F_V32_DATA', 
+                '2017G' : 'Spring18_ppRef5TeV_V4_DATA', 
                 '2018A' : 'Autumn18_RunA_V19_DATA',
                 '2018B' : 'Autumn18_RunB_V19_DATA',
                 '2018C' : 'Autumn18_RunC_V19_DATA',
                 '2018D' : 'Autumn18_RunD_V19_DATA',
                 } 
 
 jerTagsMC = {'2016' : 'Summer16_25nsV1_MC',
-             '2017' : 'Fall17_V3_MC',
+#             '2017' : 'Fall17_V3_MC',
+             '2017' : 'Spring18_ppRef5TeV_V4_MC',
              '2018' : 'Autumn18_V7_MC'
             }
 

python/postprocessing/modules/jme/jetmetUncertainties.py

@@ -522,6 +522,8 @@ def resolution_matching(jet, genjet):
 jetmetUncertainties2017 = lambda : jetmetUncertaintiesProducer("2017", "Fall17_17Nov2017_V32_MC", [ "Total" ])
 jetmetUncertainties2017METv2 = lambda : jetmetUncertaintiesProducer("2017", "Fall17_17Nov2017_V32_MC", metBranchName='METFixEE2017')
 jetmetUncertainties2017All = lambda : jetmetUncertaintiesProducer("2017", "Fall17_17Nov2017_V32_MC", [ "All" ])
+jetmetUncertainties20175TeV = lambda : jetmetUncertaintiesProducer("2017", "Spring18_ppRef5TeV_V4_MC", [ "Total" ])
+jetmetUncertainties20175TeVData = lambda : jetmetUncertaintiesProducer("2017", "Spring18_ppRef5TeV_V4_DATA", [ "Total" ], isData=True)
 
 jetmetUncertainties2018 = lambda : jetmetUncertaintiesProducer("2018", "Autumn18_V8_MC", [ "Total" ])
 jetmetUncertainties2018Data = lambda : jetmetUncertaintiesProducer("2018", "Autumn18_RunB_V8_DATA", archive="Autumn18_V8_DATA", isData=True)