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 TCanvas, TFile, TGraphErrors
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