PerformanceDB/test/plotPerformanceDB.py

0001 #! /usr/bin/env python3
0002 #
0003 #
0004 # Francisco Yumiceva
0005 # yumiceva@fnal.gov
0006 #
0007 # Fermilab, 2010
0008 #
0009 #____________________________________________________________
0010
0011 from __future__ import print_function
0012 from builtins import range
0013 import sys
0014 import math
0015 import subprocess
0016 from array import array
0017 from ROOT import *
0018
0019 def main():
0020
0021     executable = "getbtagPerformance"
0022     rootFile = "rootFile=btag_performance_CSVM.root"
0023     payload  = "payload=BTAGCSVM" #"payload=MISTAGCSVM"
0024     OP = "CSVM"
0025     flavor   = "flavor=b" #l
0026     typeSF   = "type=SF"
0027     typeEff  = "type=eff"
0028     jetpt    = "pt="
0029     jeteta   = "eta="
0030
0031     sp = " "
0032
0033     ptmin = 25.
0034     ptmax = 240. #240.
0035     etamin = 0.
0036     etamax = 2.4
0037     ptNbins = 30
0038     etaNbins = 3
0039
0040     ptbinwidth = (ptmax - ptmin)/ptNbins
0041     etabinwidth = (etamax - etamin)/etaNbins
0042     ptarray = array('d')
0043     pterrarray = array('d')
0044     etaarray = array('d')
0045     etaerrarray = array('d')
0046     SFarray = array('d')
0047     SFerrarray = array('d')
0048     effarray = array('d')
0049     efferrarray = array('d')
0050     MCeffarray = array('d')
0051     MCefferrarray = array('d')
0052     SFarray_alleta = array('d')
0053     SFerrarray_alleta = array('d')
0054     effarray_alleta = array('d')
0055     efferrarray_alleta = array('d')
0056     MCeffarray_alleta = array('d')
0057     MCefferrarray_alleta = array('d')
0058
0059
0060     for ipt in range(0,ptNbins):
0061         SF_alleta = 0.
0062         SFerr_alleta = 0.
0063         Eff_alleta = 0.
0064         Efferr_alleta = 0.
0065         MCEff_alleta = 0.
0066         MCEfferr_alleta = 0.
0067
0068         apt = ptmin + ipt*ptbinwidth
0069         ptarray.append( apt )
0070         pterrarray.append( ptbinwidth*0.5 )
0071
0072         for ieta in range(0,etaNbins):
0073
0074             aeta = etamin = ieta*etabinwidth
0075
0076             print("pt = "+str(apt) + " eta = "+str(aeta))
0077             tmpjetpt = jetpt+str(apt)
0078             tmpjeteta = jeteta+str(aeta)
0079
0080             allcmd = executable+ sp +rootFile+ sp +payload+ sp +flavor+ sp +typeSF+ sp +tmpjetpt+ sp +tmpjeteta
0081             print(allcmd)
0082             output = subprocess.getstatusoutput(allcmd)
0083             print(output[1])
0084
0085             if output[0]!=0:
0086                 print(" Error retrieving data from file DB.")
0087
0088             aSF = float( output[1].split()[2] )
0089             aSFerr = float( output[1].split()[4] )
0090
0091             SFarray.append( aSF )
0092             SFerrarray.append( aSFerr )
0093
0094             allcmd = executable+ sp +rootFile+ sp +payload+ sp +flavor+ sp +typeEff+ sp +tmpjetpt+ sp +tmpjeteta
0095             #print allcmd
0096             output = subprocess.getstatusoutput(allcmd)
0097             #print output[1]
0098
0099             if output[0]!=0:
0100                 print(" Error retrieving data from file DB.")
0101
0102             aEff = float( output[1].split()[2] )
0103             aEfferr = float( output[1].split()[4] )
0104
0105             effarray.append( aEff )
0106             efferrarray.append( aEfferr )
0107
0108             MCeffarray.append( aEff/aSF )
0109             aMCEfferr = math.sqrt( (aSF*aSF*aEfferr*aEfferr + aEff*aEff*aSFerr*aSFerr)/(aSF*aSF*aSF*aSF) )
0110             MCefferrarray.append( aMCEfferr )
0111
0112             SF_alleta += aSF
0113             SFerr_alleta += aSF*aSF
0114             Eff_alleta += aEff
0115             Efferr_alleta += aEfferr*aEfferr
0116             MCEff_alleta += aEff/aSF
0117             MCEfferr_alleta += aMCEfferr*aMCEfferr
0118
0119         SF_alleta = SF_alleta/etaNbins
0120         SFerr_alleta = math.sqrt(SFerr_alleta)/etaNbins
0121         Eff_alleta = Eff_alleta/etaNbins
0122         Efferr_alleta = math.sqrt(Efferr_alleta)/etaNbins
0123         MCEff_alleta = MCEff_alleta/etaNbins
0124         MCEfferr_alleta = math.sqrt(MCEfferr_alleta)/etaNbins
0125
0126         SFarray_alleta.append( SF_alleta )
0127         SFerrarray_alleta.append( SFerr_alleta )
0128         effarray_alleta.append( Eff_alleta )
0129         efferrarray_alleta.append( Efferr_alleta )
0130         MCeffarray_alleta.append( MCEff_alleta )
0131         MCefferrarray_alleta.append( MCEfferr_alleta )
0132
0133     histogram = {}
0134     histogram["SF_jet_pt"] = TGraphErrors(len(SFarray_alleta), ptarray, SFarray_alleta, pterrarray, SFerrarray_alleta)
0135     histogram["Eff_jet_pt"] = TGraphErrors(len(effarray_alleta), ptarray, effarray_alleta, pterrarray, efferrarray_alleta)
0136     histogram["MC_Eff_jet_pt"] = TGraphErrors(len(MCeffarray_alleta), ptarray, MCeffarray_alleta, pterrarray, MCefferrarray_alleta)
0137
0138     cv1 = TCanvas("SF_jet_pt","SF_jet_pt",700,700)
0139     histogram["SF_jet_pt"].Draw("ap")
0140     histogram["SF_jet_pt"].SetName("SF_jet_pt")
0141     histogram["SF_jet_pt"].GetHistogram().SetXTitle("Jet p_{T} [GeV/c]")
0142     histogram["SF_jet_pt"].GetHistogram().SetYTitle("SF_{b}=#epsilon^{"+OP+"}_{data}/#epsilon^{"+OP+"}_{MC}")
0143
0144     cv2 = TCanvas("Eff_jet_pt","Eff_jet_pt",700,700)
0145     histogram["Eff_jet_pt"].Draw("ap")
0146     histogram["Eff_jet_pt"].SetName("Eff_jet_pt")
0147     histogram["Eff_jet_pt"].GetHistogram().SetXTitle("Jet p_{T} [GeV/c]")
0148     histogram["Eff_jet_pt"].GetHistogram().SetYTitle("b-tag efficiency #epsilon^{"+OP+"}_{data}")
0149
0150     cv3 = TCanvas("MC_Eff_jet_pt","MC_Eff_jet_pt",700,700)
0151     histogram["MC_Eff_jet_pt"].Draw("ap")
0152     histogram["MC_Eff_jet_pt"].SetName("MC_Eff_jet_pt")
0153     histogram["MC_Eff_jet_pt"].GetHistogram().SetXTitle("Jet p_{T} [GeV/c]")
0154     histogram["MC_Eff_jet_pt"].GetHistogram().SetYTitle("b-tag efficiency #epsilon^{"+OP+"}_{MC}")
0155
0156
0157     outfilename = "plotPerformanceDB.root"
0158     outputroot = TFile( outfilename, "RECREATE")
0159
0160     for key in histogram.keys():
0161         histogram[key].Write()
0162
0163     outputroot.Close()
0164
0165
0166     raw_input ("Enter to quit:")
0167
0168
0169 if __name__ == '__main__':
0170     main()
0171