Merge pull request #13 from bianchini/Lorenzo_fixCrab

Lorenzo fix crab

crab/crab_cfg.py

@@ -12,14 +12,14 @@
 config.JobType.psetName = 'PSet.py'
 config.JobType.scriptExe = 'crab_script.sh'
 config.JobType.inputFiles = ['crab_script.py','../scripts/haddnano.py','../python/postprocessing/wmass/keep_and_drop_MC.txt', '../python/postprocessing/wmass/keep_and_drop_Data.txt']
-config.JobType.scriptArgs = ['isMC=1','passall=0']
+config.JobType.scriptArgs = ['isMC=1','passall=0', 'dataYear=2016']
 config.JobType.sendPythonFolder	 = True
 config.section_("Data")
 config.Data.inputDataset = '/WJetsToLNu_TuneCUETP8M1_13TeV-amcatnloFXFX-pythia8/RunIISummer16NanoAODv3-PUMoriond17_94X_mcRun2_asymptotic_v3-v1/NANOAODSIM'
 config.Data.inputDBS = 'global'
 #config.Data.splitting = 'Automatic'
 config.Data.splitting = 'FileBased'
-config.Data.unitsPerJob = 50
+config.Data.unitsPerJob = 10
 #config.Data.totalUnits = 10
 config.Data.outLFNDirBase = '/store/user/%s/NanoAODv3-TEST' % (getUsernameFromSiteDB())
 config.Data.publication = False
@@ -42,5 +42,5 @@
         config.General.requestName = dataset.split('/')[1]+'_sub'+str(n)+'_v1'
         config.Data.outputDatasetTag = dataset.split('/')[2]
         print config.General.requestName, '==>', config.Data.outLFNDirBase+'/'+dataset.split('/')[1]+'/'+config.Data.outputDatasetTag
-        crabCommand('submit', '--dryrun', config = config)
+        crabCommand('submit', config=config)
         n += 1

crab/crab_script.py

@@ -18,18 +18,20 @@
 parser.add_argument('-jobNum', '--jobNum',   type=int, default=1,      help="")
 parser.add_argument('-passall', '--passall', type=int, default=0,      help="")
 parser.add_argument('-isMC', '--isMC',       type=int, default=1,      help="")
+parser.add_argument('-dataYear', '--dataYear',type=int, default=2016, help="")
 args = parser.parse_args()
 isMC    = args.isMC
 passall = args.passall
-print "isMC =", isMC, ", passall=", passall
+dataYear = args.dataYear
+print "isMC =", isMC, ", passall =", passall, ", dataYear =", dataYear
 
 from PhysicsTools.NanoAODTools.postprocessing.framework.crabhelper import inputFiles,runsAndLumis
 
 modules = []
 if isMC:
-    modules = [puAutoWeight(), preselection(isMC=isMC, passall=passall), additionalVariables(), leptonSelectModule(), CSAngleModule(), WproducerModule()]
+    modules = [puAutoWeight(), preselection(isMC=isMC, passall=passall, dataYear=dataYear), additionalVariables(isMC=isMC), leptonSelectModule(), CSAngleModule(), WproducerModule()]
 else:
-    modules = [preselection(isMC=isMC, passall=passall), additionalVariables()]
+    modules = [preselection(isMC=isMC, passall=passall, dataYear=dataYear), additionalVariables(isMC=isMC)]
 
 p = PostProcessor(outputDir=".",        
                   inputFiles=inputFiles(),

crab/crab_script.sh

@@ -24,5 +24,5 @@ mv module $CMSSW_BASE/module
 mv python $CMSSW_BASE/python
 
 echo Found Proxy in: $X509_USER_PROXY
-python crab_script.py --jobNum=$1 --$2 --$3
+python crab_script.py --jobNum=$1 --$2 --$3 --$4
 fi

python/postprocessing/wmass/CSVariables.py

@@ -15,7 +15,7 @@ def getCSangles(muon, neutrino):
 
     w = m + n
 
-    if(w.Z()==0):
+    if(w.Z()==0.):
 	sign=1
     else :
 	sign  = abs(w.Z())/w.Z()
@@ -79,15 +79,15 @@ def analyze(self, event):
             preFSRMuonIdx = event.GenPart_preFSRMuonIdx
             dressMuonIdx = event.GenDressedLepton_dressMuonIdx
             if(bareMuonIdx>=0) :
-                CStheta_bare, CSphi_bare     = getCSangles(genParticles[bareMuonIdx], genParticles[NeutrinoIdx])
+                CStheta_bare, CSphi_bare  = getCSangles(genParticles[bareMuonIdx], genParticles[NeutrinoIdx])
             else :
                 CStheta_bare, CSphi_bare = 0.0, 0.0
             if(preFSRMuonIdx>=0) :
                 CStheta_preFSR, CSphi_preFSR = getCSangles(genParticles[preFSRMuonIdx], genParticles[NeutrinoIdx])
             else :
                 CStheta_preFSR, CSphi_preFSR = 0.0, 0.0
             if(dressMuonIdx>=0) :
-                CStheta_dress, CSphi_dress   = getCSangles(genDressedLeptons[dressMuonIdx], genParticles[NeutrinoIdx])
+                CStheta_dress, CSphi_dress = getCSangles(genDressedLeptons[dressMuonIdx], genParticles[NeutrinoIdx])
             else :
                 CStheta_dress, CSphi_dress = 0.0, 0.0
 

python/postprocessing/wmass/additionalVariables.py

@@ -43,13 +43,15 @@ def evaluateMtRecoBased(mu_pt, mu_phi, met_pt, met_phi) :
 
 
 class additionalVariables(Module):
-    def __init__(self):
+    def __init__(self, isMC=True):
         self.metdict = {
-        "pf" : "",
-        "gen" : "Gen",
-        "tk" : "Tk",
-        "puppi" : "Puppi"
-        }
+            "pf" : "",
+            #"gen" : "Gen",
+            "tk" : "Tk",
+            "puppi" : "Puppi"
+            }
+        if isMC:
+            self.metdict["gen"] = "Gen"
         pass
 
     def beginJob(self):

python/postprocessing/wmass/keep_and_drop_MC.txt

@@ -17,3 +17,5 @@ keep PV*
 keep Muon*
 keep MET*
 keep TkMET*
+keep GenMET*
+keep PuppiMET*

python/postprocessing/wmass/lepSelection.py

@@ -39,6 +39,8 @@ def analyze(self, event):
 
         #look at the flavour of the highest pt neutrino to decide if accept or not the event
         neutrini.sort(key = lambda x: x[1].pt, reverse=True )
+
+        # return if there are no neutrini or the highest-pt nu is not of type mu
         if len(neutrini)==0 or abs(neutrini[0][1].pdgId) != 14:
             self.out.fillBranch("genVtype", -1)
             self.out.fillBranch("GenPart_bareMuonIdx", -1)
@@ -53,7 +55,6 @@ def analyze(self, event):
         myIdx = baremuons[0][0]
         myNuIdx = neutrini[0][0]
 
-
         self.out.fillBranch("GenPart_bareMuonIdx",myIdx)
         self.out.fillBranch("GenPart_NeutrinoIdx", myNuIdx)
 

python/postprocessing/wmass/postproc.py

@@ -15,13 +15,15 @@
 from Wproducer import *
 
 parser = argparse.ArgumentParser("")
-parser.add_argument('-jobNum', '--jobNum',   type=int, default=1,      help="")
-parser.add_argument('-passall', '--passall', type=int, default=0,      help="")
-parser.add_argument('-isMC', '--isMC',       type=int, default=1,      help="")
+parser.add_argument('-jobNum', '--jobNum',    type=int, default=1,      help="")
+parser.add_argument('-passall', '--passall',  type=int, default=0,      help="")
+parser.add_argument('-isMC', '--isMC',        type=int, default=1,      help="")
+parser.add_argument('-dataYear', '--dataYear',type=int, default=2016, help="")
 args = parser.parse_args()
 isMC    = args.isMC
 passall = args.passall
-print "isMC =", isMC, ", passall=", passall
+dataYear = args.dataYear
+print "isMC =", isMC, ", passall =", passall, ", dataYear =", dataYear
 
 input_dir = "/gpfs/ddn/srm/cms/store/"
 
@@ -33,12 +35,12 @@
         #input_dir+"mc/RunIISummer16NanoAODv3/DYJetsToLL_M-50_TuneCUETP8M1_13TeV-madgraphMLM-pythia8/NANOAODSIM/PUMoriond17_94X_mcRun2_asymptotic_v3_ext1-v2/120000/06544C90-EFDF-E811-80E6-842B2B6F5D5C.root",
         input_dir+"mc/RunIISummer16NanoAODv3/DYJetsToLL_Pt-50To100_TuneCUETP8M1_13TeV-amcatnloFXFX-pythia8/NANOAODSIM/PUMoriond17_94X_mcRun2_asymptotic_v3-v2/280000/26DE6A2F-9329-E911-8766-002590DE6E8A.root",
         )
-    modules = [puAutoWeight(), preselection(isMC=isMC, passall=passall), additionalVariables(), leptonSelectModule(), CSAngleModule(), WproducerModule()]
+    modules = [puAutoWeight(), preselection(isMC=isMC, passall=passall, dataYear=dataYear), additionalVariables(isMC=isMC), leptonSelectModule(), CSAngleModule(), WproducerModule()]
 else:
     input_files.append(
         input_dir+"data/Run2016D/DoubleEG/NANOAOD/Nano14Dec2018-v1/280000/481DA5C0-DF96-5640-B5D1-208F52CAC829.root"
         )
-    modules = [preselection(isMC=isMC, passall=passall), additionalVariables()]
+    modules = [preselection(isMC=isMC, passall=passall, dataYear=dataYear), additionalVariables(isMC=isMC)]
 
 p = PostProcessor(outputDir=".",        
                   inputFiles=input_files,

python/postprocessing/wmass/preselection.py

@@ -11,11 +11,49 @@ def medium_muon_id(mu):
     return (abs(mu.eta)<2.4 and mu.pt>20 and abs(mu.dxy)<0.05 and abs(mu.dz)<0.2 and mu.mediumId and mu.pfRelIso04_all<=0.10)
 def medium_aiso_muon_id(mu):
     return (abs(mu.eta)<2.4 and mu.pt>20 and abs(mu.dxy)<0.05 and abs(mu.dz)<0.2 and mu.mediumId and mu.pfRelIso04_all >0.10 and mu.pfRelIso04_all<0.30)
+def veto_electron_id(ele):
+    etaSC = ele.eta + ele.deltaEtaSC
+    if (abs(etaSC) <= 1.479):
+        return (ele.pt>10 and abs(ele.dxy)<0.05 and abs(ele.dz)<0.1 and ele.cutBased>=1 and ele.pfRelIso03_all< 0.30)
+    else:
+        return (ele.pt>10 and abs(ele.dxy)<0.10 and abs(ele.dz)<0.2 and ele.cutBased>=1 and ele.pfRelIso03_all< 0.30)
 
 class preselection(Module):
-    def __init__(self, isMC=True, passall=False):
+    def __init__(self, isMC=True, passall=False, dataYear=2016):
         self.isMC = isMC
         self.passall = passall
+        self.dataYear = str(dataYear)
+
+        # baded on https://twiki.cern.ch/twiki/bin/view/CMS/MissingETOptionalFiltersRun2
+        # val = 2*(IS FOR DATA) + 1*(IS FOR MC)
+        self.met_filters = {
+            "2016": { "goodVertices" : 3,
+                      "globalTightHalo2016Filter": 3,
+                      "HBHENoiseFilter": 3,
+                      "HBHENoiseIsoFilter" : 3,
+                      "EcalDeadCellTriggerPrimitiveFilter": 3
+                      },
+            "2017": { "goodVertices": 3,
+                      "globalTightHalo2016Filter": 3,
+                      "HBHENoiseFilter": 3,
+                      "HBHENoiseIsoFilter": 3,
+                      "EcalDeadCellTriggerPrimitiveFilter": 3,
+                      "BadPFMuonFilter": 3,
+                      "BadChargedCandidateFilter": 3,
+                      "eeBadScFilter" : 2,
+                      #"ecalBadCalibFilter": 3, # needs to be rerun
+                      },
+            "2018": { "goodVertices": 3,
+                      "globalTightHalo2016Filter": 3,
+                      "HBHENoiseFilter": 3,
+                      "HBHENoiseIsoFilter": 3,
+                      "EcalDeadCellTriggerPrimitiveFilter": 3,
+                      "BadPFMuonFilter" : 3,
+                      "BadChargedCandidateFilter" : 3,
+                      "eeBadScFilter" : 2,
+                      #"ecalBadCalibFilter": 3,  # needs to be rerun 
+                      },
+            }
         pass
     def beginJob(self):
         pass
@@ -27,6 +65,8 @@ def beginFile(self, inputFile, outputFile, inputTree, wrappedOutputTree):
         self.out.branch("Muon_idx1", "I", title="index of W-like muon / index of Z-like 1st muon")
         self.out.branch("Muon_idx2", "I", title="index of Z-like 2nd muon")
         self.out.branch("Vtype", "I", title="0:W-like; 1:Fake-like; 2:Z-like; 3:SS-dimuon")
+        self.out.branch("MET_filters", "I", title="AND of all MET filters")
+        self.out.branch("nVetoElectrons", "I", title="Number of veto electrons")
 
     def endFile(self, inputFile, outputFile, inputTree, wrappedOutputTree):
         pass
@@ -40,14 +80,17 @@ def analyze(self, event):
             for hlt in triggers_OR:
                 HLT_pass |= getattr(event, "HLT_"+hlt, False)
             self.out.fillBranch("Muon_HLT", int(HLT_pass))
-            if not HLT_pass and not self.passall:
-                return False
         else:
             self.out.fillBranch("Muon_HLT", 1)
 
+        # MET filters
+        met_filters_AND = True
+        for key,val in self.met_filters[self.dataYear].items():
+            met_filters_AND &=  (not (val & (1 << (1-int(self.isMC)) )) or getattr(event, "Flag_"+key))
+        self.out.fillBranch("MET_filters", int(met_filters_AND))
+
         # Muon selection
         all_muons = Collection(event, "Muon")
-
         loose_muons       = [ [mu,imu] for imu,mu in enumerate(all_muons) if loose_muon_id(mu) ]
         medium_muons      = [ [mu,imu] for imu,mu in enumerate(all_muons) if medium_muon_id(mu)]
         medium_aiso_muons = [ [mu,imu] for imu,mu in enumerate(all_muons) if medium_aiso_muon_id(mu)]
@@ -56,6 +99,11 @@ def analyze(self, event):
         medium_muons.sort(key = lambda x: x[0].pt, reverse=True )
         medium_aiso_muons.sort(key = lambda x: x[0].pt, reverse=True )
 
+        # Electron selection
+        all_electrons  = Collection(event, "Electron")    
+        loose_electrons = [ [ele,iele] for iele,ele in enumerate(all_electrons) if veto_electron_id(ele) ]
+        self.out.fillBranch("nVetoElectrons", len(loose_electrons))
+
         event_flag = -1        
         (idx1, idx2) = (-1, -1)
         # W-like event: 1 loose, 1 medium
@@ -70,7 +118,7 @@ def analyze(self, event):
         elif len(loose_muons)==2:
             (idx1, idx2) = (loose_muons[0][1], loose_muons[1][1])
             event_flag = 2 if (loose_muons[0][0].charge+loose_muons[1][0].charge)==0 else 3
-        # anything else
+        # anything else        
         else:   
             event_flag = -1