diff --git a/plotting/HLTNominalvsHLTBest.py b/plotting/HLTNominalvsHLTBest.py new file mode 100644 index 0000000..7953144 --- /dev/null +++ b/plotting/HLTNominalvsHLTBest.py @@ -0,0 +1,550 @@ +# Compare Inputs between HLT Nominal and HLT Best +# --input1 ROOTFileWithHists_Nominal +# --input2 ROOTFileWithHists_Best +# --output outDirectory +import ROOT + + +ROOT.gROOT.SetBatch(True) + +import ROOTHelp.FancyROOTStyle + +from optparse import OptionParser +p = OptionParser() +p.add_option('--input1', type = 'string', default = "outBTag.FTKBtagging.ttbar.mwt2.All.root", dest = 'inFile1', help = 'intput File1' ) +p.add_option('--input2', type = 'string', default = "outBTag.FTKBtagging.ttbar.mwt2.All.root", dest = 'inFile2', help = 'intput File2' ) #another input forHLT2 +p.add_option('--output', type = 'string', default = "HLTNominalvsBest", dest = 'outDir', help = 'output dir' ) +p.add_option('--doCaloJets', action="store_true",help = '' ) +p.add_option('--cmsText', type = 'string', default = "Work in Progress", help = '' ) +#p.add_option('--lumiText', default = "", help = '' ) +(o,a) = p.parse_args() + + + +from ROOTHelp.Utils import do_variable_rebinning, makeCanvas +from ROOTHelp.Plotting import makeRatio +#from rocCurveUtils import drawWaterMarks +#import rebinning +from Rebinning import rebinningDB + +from JetLevelPlotUtils import getCMSText + +inFile1 = ROOT.TFile(o.inFile1, "READ") +inFile2 = ROOT.TFile(o.inFile2, "READ") + +import os +if not os.path.exists(o.outDir): + os.mkdir(o.outDir) + + +maxDict = {"jetNSelectedTracks":20, + "jetNTracks":30, + } + + +def getHist(inFile,dir,var,binning,color): + hist = inFile.Get(dir+"/"+var) + print dir+"/"+var + if type(binning ) == type(list()): + hist = do_variable_rebinning(hist,binning) + else: + hist.Rebin(binning) + + hist.SetLineColor(color) + hist.SetMarkerColor(color) + hist.Sumw2() + if hist.Integral(): + hist.Scale(1./hist.Integral()) + return hist + +#def doVar(var,binning,xTitle,setLogy=1): +# #NomiLF = getHist(inFile,"caloJets_matchedL", var,binning,ROOT.kBlack) +# NomiLF = getHist(inFile,"offJets_matched_L", var,binning,ROOT.kBlack) +# BestLF = getHist(inFile,"offJets_matchedJet_L",var,binning,ROOT.kBlack) +# #BestLF = getHist(inFile,"caloJets_matchedL",var,binning,ROOT.kBlack) +# #NomiBQ = getHist(inFile,"caloJets_matchedB", var,binning,ROOT.kRed) +# NomiBQ = getHist(inFile,"offJets_matched_B", var,binning,ROOT.kRed) +# BestBQ = getHist(inFile,"offJets_matchedJet_B",var,binning,ROOT.kRed) +# #BestBQ = getHist(inFile,"caloJets_matchedB",var,binning,ROOT.kRed) +# +# xpos = 0.18 +# ypos = 0.6 +# xwidth = 0.3 +# ywidth = 0.3 +# +# leg = ROOT.TLegend(xpos, ypos, xpos+xwidth, ypos+ywidth) +# leg.AddEntry(NomiBQ,"Offline tracks b-quark jets","L") +# leg.AddEntry(BestBQ,"HLT PF b-quark jets" ,"PEL") +# leg.AddEntry(NomiLF,"Offline tracks light-flavor jets","L") +# leg.AddEntry(BestLF,"HLT PF light-flavor jets" ,"PEL") +# +# +# can = ROOT.TCanvas(var,var) +# can.cd().SetLogy(setLogy) +# maxY = max(NomiLF.GetMaximum(),NomiBQ.GetMaximum(), +# BestLF.GetMaximum(),BestBQ.GetMaximum()) +# if setLogy: +# NomiLF.SetMaximum(4e0*maxY) +# NomiLF.SetMinimum(1e-6) +# else: +# NomiLF.SetMaximum(1.2*maxY) +# +# NomiLF.GetYaxis().SetTitle("Simulated Tracks") +# NomiLF.GetXaxis().SetTitle(xTitle ) +# NomiLF.Draw("hist") +# #BestLF.SetMarkerSize(0.75) +# #BestLF.SetMarkerStyle(21) +# BestLF.Draw("same pe") +# NomiBQ.Draw("hist same") +# #BestBQ.SetMarkerSize(0.75) +# #BestBQ.SetMarkerStyle(21) +# BestBQ.Draw("same pe") +# #leg.Draw("same") +# +# cmsLine1, cmsLine2 = getCMSText(xStart=0.2,yStart=0.87) +# cmsLine1.Draw("same") +# cmsLine2.Draw("same") +# +# +# #xatlas, yatlas = 0.18, 0.88 +# #atlas = ROOT.TLatex(xatlas+0.01, yatlas, "ATLAS") +# ##simulation = ROOT.TLatex(xatlas+0.11, yatlas, "Simulation Internal") +# #simulation = ROOT.TLatex(xatlas+0.11, yatlas, "Simulation Preliminary") +# #lumi = ROOT.TLatex(xatlas+0.01, yatlas-0.05, "#sqrt{s}=13 TeV, t#bar{t}") +# #jetText = ROOT.TLatex(xatlas+0.02, yatlas-0.1, "p_{T}^{jet} > 40 GeV, |#eta^{jet}| < 2.5" ) +# #wm = [atlas, simulation, lumi, jetText] +# +# #watermarks = drawWaterMarks(wm) +# +# varName = var.replace("tracks/","track_").replace("btags/","btag_") +# can.SaveAs(o.outDir+"/"+varName+".pdf") +# can.SaveAs(o.outDir+"/"+varName+".eps") +# can.SaveAs(o.outDir+"/"+varName+".png") + + +def doVarRatio(var,binning,xTitle,setLogy=1,minX=None,maxX=None,minY=None): + if o.doCaloJets: + NomiLF = getHist(inFile1,"offJets_matchedCaloJet_L", var,binning,ROOT.kBlack) + BestLF = getHist(inFile2,"offJets_matchedCaloJet_L",var,binning,ROOT.kBlack) + NomiBQ = getHist(inFile1,"offJets_matchedCaloJet_B", var,binning,ROOT.kRed) + BestBQ = getHist(inFile2,"offJets_matchedCaloJet_B",var,binning,ROOT.kRed) + # **matchedJet_L & B + # from Input1 and Input2 + else: + NomiLF = getHist(inFile1,"offJets_matchedJet_L", var,binning,ROOT.kBlack) + BestLF = getHist(inFile2,"offJets_matchedJet_L",var,binning,ROOT.kBlack) + NomiBQ = getHist(inFile1,"offJets_matchedJet_B", var,binning,ROOT.kRed) + BestBQ = getHist(inFile2,"offJets_matchedJet_B",var,binning,ROOT.kRed) + + + maxY = max(NomiLF.GetMaximum(),NomiBQ.GetMaximum(), + BestLF.GetMaximum(),BestBQ.GetMaximum()) + if setLogy: + NomiLF.SetMaximum(4e0*maxY) + NomiLF.SetMinimum(1.01e-5) + else: + NomiLF.SetMaximum(1.2*maxY) + + if not minY == None : + NomiLF.SetMinimum(minY) + + + NomiLF.GetYaxis().SetTitle("Simulated Tracks") + NomiLF.GetXaxis().SetTitle(xTitle ) + + if maxX: + NomiLF.GetXaxis().SetRangeUser(minX,maxX) + NomiBQ.GetXaxis().SetRangeUser(minX,maxX) + BestLF.GetXaxis().SetRangeUser(minX,maxX) + BestBQ.GetXaxis().SetRangeUser(minX,maxX) + + + LFRatio = makeRatio(num = BestLF.Clone(), den = NomiLF.Clone()) + BQRatio = makeRatio(num = BestBQ.Clone(), den = NomiBQ.Clone()) + + xpos = 0.2 + ypos = 0.07 + xwidth = 0.7 + ywidth = 0.1 + + leg = ROOT.TLegend(xpos, ypos, xpos+xwidth, ypos+ywidth) + leg.SetNColumns(2) + leg.AddEntry(NomiBQ,"Nominal tracks b-quark jets","L") + leg.AddEntry(NomiLF,"Nominal tracks light-flavor jets","L") + leg.AddEntry(BestBQ,"Best cut's tracks b-quark jets" ,"PEL") + leg.AddEntry(BestLF,"Best cut's tracks light-flavor jets" ,"PEL") + + canvas = makeCanvas(var, var, width=600, height=600) + split=0.3 + top_pad = ROOT.TPad("pad1", "The pad 80% of the height",0,split,1,1,0) + bottom_pad = ROOT.TPad("pad2", "The pad 20% of the height",0,0,1,split,0) + top_pad.Draw() + bottom_pad.Draw() + + axissep = 0.02 + top_pad.cd() + top_pad.SetLogy(setLogy) + top_pad.SetTopMargin(canvas.GetTopMargin()*1.0/(1.0-split)) + top_pad.SetBottomMargin(0.5*axissep) + top_pad.SetRightMargin(canvas.GetRightMargin()) + top_pad.SetLeftMargin(canvas.GetLeftMargin()); + top_pad.SetFillStyle(0) # transparent + top_pad.SetBorderSize(0) + + if var in maxDict.keys(): + NomiLF.GetXaxis().SetRangeUser(NomiLF.GetXaxis().GetXmin(),maxDict[var] ) + + NomiLF.Draw("hist") + #BestLF.SetMarkerSize(0.75) + #BestLF.SetMarkerStyle(21) + BestLF.Draw("same pe") + NomiBQ.Draw("hist same") + #BestBQ.SetMarkerSize(0.75) + #BestBQ.SetMarkerStyle(21) + BestBQ.Draw("same pe") + leg.Draw("same") + + + + + + + cmsLines = getCMSText(xStart=0.2,yStart=0.85,subtext=o.cmsText) + for cmsl in cmsLines: + cmsl.Draw("same") + + + + bottom_pad.cd() + bottom_pad.SetTopMargin(2*axissep) + bottom_pad.SetBottomMargin(canvas.GetBottomMargin()*1.0/split) + bottom_pad.SetRightMargin(canvas.GetRightMargin()) + bottom_pad.SetLeftMargin(canvas.GetLeftMargin()); + bottom_pad.SetFillStyle(0) # transparent + bottom_pad.SetBorderSize(0) + ratio_axis = NomiLF.Clone() + #ratio_axis.GetYaxis().SetTitle("PF to Calo") + ratio_axis.GetYaxis().SetTitle("Best to Nominal") + ratio_axis.GetXaxis().SetTitle(NomiLF.GetXaxis().GetTitle()) + ratio_axis.GetYaxis().SetNdivisions(507) + rMin = 0 + rMax = 2 + + + if var in maxDict.keys(): + LFRatio.GetXaxis().SetRangeUser(LFRatio.GetXaxis().GetXmin(),maxDict[var] ) + + ratio_axis.GetYaxis().SetRangeUser(rMin, rMax) + LFRatio.GetYaxis().SetRangeUser(rMin, rMax) + #LFRatio.GetYaxis().SetTitle("Calo to Offline") + #LFRatio.GetYaxis().SetTitle("PF to Calo") + LFRatio.GetYaxis().SetTitle("Best to Nominal") + + + LFRatio.Draw("PE") + LFRatio.Draw("PE same") + oldSize = LFRatio.GetMarkerSize() + LFRatio.SetMarkerSize(0) + LFRatio.DrawCopy("same e0") + LFRatio.SetMarkerSize(oldSize) + LFRatio.Draw("PE same") + + + BQRatio.Draw("PE same") + BQRatio.Draw("PE same") + oldSize = BQRatio.GetMarkerSize() + BQRatio.SetMarkerSize(0) + BQRatio.DrawCopy("same e0") + BQRatio.SetMarkerSize(oldSize) + BQRatio.Draw("PE same") + + + + line = ROOT.TLine() + line.DrawLine(NomiLF.GetXaxis().GetXmin(), 1.0, NomiLF.GetXaxis().GetXmax(), 1.0) + + ndivs=[505,503] + + pads = [top_pad, bottom_pad] + factors = [0.8/(1.0-split), 0.7/split] + for i_pad, pad in enumerate(pads): + + factor = factors[i_pad] + ndiv = ndivs[i_pad] + + prims = [ p.GetName() for p in pad.GetListOfPrimitives() ] + + # + # Protection for scaling hists multiple times + # + procedHist = [] + + for name in prims: + + if name in procedHist: continue + procedHist.append(name) + + h = pad.GetPrimitive(name) + if isinstance(h, ROOT.TH1) or isinstance(h, ROOT.THStack) or isinstance(h, ROOT.TGraph) or isinstance(h, ROOT.TGraphErrors) or isinstance(h, ROOT.TGraphAsymmErrors): + if isinstance(h, ROOT.TGraph) or isinstance(h, ROOT.THStack) or isinstance(h, ROOT.TGraphErrors) or isinstance(h, ROOT.TGraphAsymmErrors): + h = h.GetHistogram() + #print "factor is",factor,h.GetName(),split + + if i_pad == 1: + h.SetLabelSize(h.GetLabelSize('Y')*factor, 'Y') + h.SetTitleSize(h.GetTitleSize('X')*factor, 'X') + h.SetTitleSize(h.GetTitleSize('Y')*factor, 'Y') + h.SetTitleOffset(h.GetTitleOffset('Y')/factor, 'Y') + + if i_pad == 1: + h.GetYaxis().SetNdivisions(ndiv) + h.GetXaxis().SetNdivisions() + if i_pad == 0: + h.SetLabelSize(0.0, 'X') + h.GetXaxis().SetTitle("") + else: + h.SetLabelSize(h.GetLabelSize('X')*factor, 'X') + ## Trying to remove overlapping y-axis labels. Doesn't work. + # h.GetYaxis().Set(4, h.GetYaxis().GetXmin(), h.GetYaxis().GetXmax()) + # h.GetYaxis().SetBinLabel( h.GetYaxis().GetLast(), '') + + + + + + ##xatlas, yatlas = 0.18, 0.88 + ##atlas = ROOT.TLatex(xatlas+0.01, yatlas, "ATLAS") + ###simulation = ROOT.TLatex(xatlas+0.11, yatlas, "Simulation Internal") + ##simulation = ROOT.TLatex(xatlas+0.11, yatlas, "Simulation Preliminary") + ##lumi = ROOT.TLatex(xatlas+0.01, yatlas-0.05, "#sqrt{s}=13 TeV, t#bar{t}") + ##jetText = ROOT.TLatex(xatlas+0.02, yatlas-0.1, "p_{T}^{jet} > 40 GeV, |#eta^{jet}| < 2.5" ) + ##wm = [atlas, simulation, lumi, jetText] + + #watermarks = drawWaterMarks(wm) + + #varName = var.replace("tracks/","").replace("btags/","btags_") + varName = var.replace("tracks/","track_").replace("btags/","btag_").replace("btags_noV0/","btag_noV0_") + canvas.SaveAs(o.outDir+"/"+varName+".pdf") + #canvas.SaveAs(o.outDir+"/"+var+".eps") + #canvas.SaveAs(o.outDir+"/"+var+".png") + + + +#doVar("ip3d_sig_l", +# xTitle = "ip3d significance", +# #binning = [-20,-18,-16,-14,-12,-11,-10,-8,-6,-5,-4,-3,-2.5,-2,-1.5,-1,-0.5,0,0.5,1,1.5,2,2.5,3,4,5,6,8,10,11,12,14,16,18,20,22,24,28,32,36,40] +# binning = [-100 , -90,-80 , -70 , -60 , -50 , -40 , -34 , -32 , -30 , -28 , -26 , -24 , -22 , -20 , -18 , -16 , -14 , -12 , -10 , -8 , -6 , -4 , -2 , 0 , 2 , 4 , 6 , 8 , 10 , 12 , 14 , 16 , 18 , 20 , 22 , 24 , 26 , 28 , 30 , 32 , 34 , 40 , 50 , 60 , 70 , 80, 90 , 100] +# #binning = [-100,-90,-80,-70,-60,-50,-40,-30,-20,-10,-5,0,5,10,15,20,30,40,50,60,70,80,90,100] +# ) + +for v in ["tracks/ip3d_sig", + "tracks/ip2d_sig", + "CSVv2_l", + "DeepCSV_l", + #"deepcsv_bb", + "btags/sv_Flight2D", + "btags/sv_FlightSig2D", + "btags/sv_FlightSig", + "btags/sv_Flight", + "tracks/ip2d", + "tracks/ip2d_l", + "tracks/ip2d_sig", + "tracks/ip2d_sig_l", + "tracks/ip3d", + "tracks/ip3d_l", + "tracks/ip3d_sig", + "tracks/ip3d_sig_l", + "tracks/pt_s", + + "btags/ip2d", + "btags/ip2d_l", + "btags/ip2d_sig", + "btags/ip2d_sig_l", + "btags/ip3d", + "btags/ip3d_l", + "btags/ip3d_sig", + "btags/ip3d_sig_l", + + "pt_s", + "pt_m", + #"trackJetPt", + #"trackSip2dSigAboveCharm", + #"trackSip2dValAboveCharm", + #"trackSip3dSigAboveCharm", + #"trackSip3dValAboveCharm", + #"trackSumJetDeltaR", + #"vertexFitProb", + + "tracks/PtRel" , + "tracks/PtRatio" , + "tracks/PPar" , + "tracks/PParRatio" , + "tracks/Momentum" , + "tracks/DecayLenVal_l" , + "tracks/DecayLenVal" , + "tracks/algo", + "tracks/origAlgo", + + "btags/sv_Pt", + ]: + + vName = v.split("/")[-1] + if vName in rebinningDB: + binning = rebinningDB[vName] + else: + binning = 2 + + doVarRatio(v, + xTitle = v, + #binning = [-20,-18,-16,-14,-12,-11,-10,-8,-6,-5,-4,-3,-2.5,-2,-1.5,-1,-0.5,0,0.5,1,1.5,2,2.5,3,4,5,6,8,10,11,12,14,16,18,20,22,24,28,32,36,40] + binning = binning, + #binning = [-100 , -90,-80 , -70 , -60 , -50 , -40 , -34 , -32 , -30 , -28 , -26 , -24 , -22 , -20 , -18 , -16 , -14 , -12 , -10 , -8 , -6 , -4 , -2 , 0 , 2 , 4 , 6 , 8 , 10 , 12 , 14 , 16 , 18 , 20 , 22 , 24 , 26 , 28 , 30 , 32 , 34 , 40 , 50 , 60 , 70 , 80, 90 , 100] + #binning = [-100,-90,-80,-70,-60,-50,-40,-30,-20,-10,-5,0,5,10,15,20,30,40,50,60,70,80,90,100] + ) + + + + +for v in [ "tracks/eta", + "tracks/ip2d_err", + "tracks/ip2d_err_l", + "tracks/ip3d_err", + "tracks/ip3d_err_l", +# "neMult", + "phi", + "btags/sv_BoostOverSqrtJetPt", + "btags/sv_EnergyRatio", + "btags/sv_Eta", + "btags/sv_NDF", + "btags/sv_Phi", + "btags/sv_R", + "btags/sv_Z", + "btags/sv_massVertexEnergyFraction", + "btags/sv_Chi2", + "btags/sv_JetDeltaR", + "btags/sv_DistJetAxis", + "tracks/JetDistVal" , + "tracks/eta" , + "tracks/phi" , + "tracks/DeltaR" , + #"trackEtaRel" , + #"jetNSelectedTracks", + #"mult", + #"nTrk", + #"jetNTracksEtaRel", + + "btags/chargedHadronMultiplicity", + "btags/chargedMultiplicity", + "btags/elecMultiplicity", + "btags/muonMultiplicity", + "btags/neutralHadronMultiplicity", + "btags/neutralMultiplicity", + "btags/photonMultiplicity", + "btags/totalMultiplicity", + + #"vertexCategory", + #"ip2d", + "tracks/Chi2", + ]: + + if o.doCaloJets: + if v in [ "btags/chargedHadronMultiplicity", + "btags/chargedMultiplicity", + "btags/elecMultiplicity", + "btags/muonMultiplicity", + "btags/neutralHadronMultiplicity", + "btags/neutralMultiplicity", + "btags/photonMultiplicity", + "btags/totalMultiplicity", + ]: + continue + + vName = v.split("/")[-1] + if vName in rebinningDB: + binning = rebinningDB[vName] + else: + binning = 2 + + + doVarRatio(v, + xTitle = v, + binning = binning, + setLogy = 0, + ) + + + + +for v in [ + "tracks/HasInnerPixHit", + "tracks/NPixelHits", + "tracks/NTotalHits", + "tracks/NStripHits", + "btags/sv_nSVs", + "tracks/nTracks", + "btags_noV0/nTracks", + "btags/sv_NTracks", + ]: + + doVarRatio(v, + xTitle = v, + binning = 1, + setLogy = 0, + ) + + +for v in [ +"btags/chargedEmEnergyFraction", +"btags/chargedHadronEnergyFraction", +"btags/elecEnergyFraction", +"btags/muonEnergyFraction", +"btags/neutralEmEnergyFraction", +"btags/neutralHadronEnergyFraction", +"btags/photonEnergyFraction", + + + "tracks/IsFromV0", + "tracks/IsFromSV", +# "neutralHadronEnergyFraction", +# "trackSumJetEtRatio", + ]: + doVarRatio(v, + xTitle = v, + binning = 1, + setLogy = 0, + minY = 0, + ) + +# +# +#for v in [ +# "muonEnergyFraction", +# "muEF", +# +# ]: +# doVarRatio(v, +# xTitle = v, +# binning = 1, +# setLogy = 1, +# ) +# +for v in [ + "btags/sv_Mass", + "m", + ]: + doVarRatio(v, + xTitle = v, + binning = 1, + setLogy = 0, + minX = 0, + maxX = 15 + ) + + + +#doVar("trk_z0Sig_signed", +# xTitle = "z_{0} Significance", +# #binning = [-25,-20,-18,-14,-10,-8,-6,-5,-4,-3,-2.5,-2,-1.5,-1,-0.5,0,0.5,1,1.5,2,2.5,3,4,5,6,8,10,12,14,16,18,20,24,25] +# binning = [-12,-10,-8,-6,-5,-4,-3,-2.5,-2,-1.5,-1,-0.5,0,0.5,1,1.5,2,2.5,3,4,5,6,8,10,12,15,20] +# ) + diff --git a/plotting/joinROCs.py b/plotting/joinROCs.py new file mode 100644 index 0000000..d0ff5e2 --- /dev/null +++ b/plotting/joinROCs.py @@ -0,0 +1,52 @@ +### inFile = the output of saveROCs.py where ROC curve froms different cuts are saved in one ROOT File +### join the ROC curves together in one graph in the output ROOT file + +import ROOT + +inFile = ROOT.TFile("PFDeepCSV_ROCs_low.root","READ") + +inFile.ls() + +ptCuts = [ + "0p4", + "0p9", + "1p5", + "2p0", + "2p5", + "5p0", + "7p5", + "10p0", + "20p0", + ] + +graphs = [] +for cut in ptCuts: + #### need to edit depending on the names of inFiles + graphs.append(inFile.Get("PtScan_low_"+cut+"_0720;2")) + #### need to edit depending on the names of inFiles + +can = ROOT.TCanvas() +#can.cd().SetLogy(1) +#can.cd().SetLogx(1) +# draw axis +hist = ROOT.TH1F("Axis","PFDeepCSV_ROCs_LP;B-Jet Efficiency; Light Flavor Efficiency",2,0.6,1) +hist.GetYaxis().SetRangeUser(1e-3,0.6) +hist.SetStats(0) +hist.Draw() +# draw legend +leg = ROOT.TLegend(0.5, 0.6, 0.65, 0.89) +leg.SetHeader("Track Pt Cut at: ") + +for iG, graph in enumerate(graphs): + graph.SetLineColor(iG+1) + graph.SetMarkerColor(iG+1) + graph.Draw("PL") + if ptCuts[iG] == "7p5" or ptCuts[iG] == "10p0" or ptCuts[iG] == "20p0": + leg.AddEntry(graph, "("+ptCuts[iG]+"GeV )", "l") # add parentheses + else: + leg.AddEntry(graph, ptCuts[iG]+"GeV", "l") + + +leg.Draw() +can.SaveAs("0810low.root") + diff --git a/plotting/saveROCs.py b/plotting/saveROCs.py new file mode 100644 index 0000000..2861ff3 --- /dev/null +++ b/plotting/saveROCs.py @@ -0,0 +1,181 @@ +### extract multiple ROC Plots of the same variable from different root files in a folder and store them in one output root file +# --inFolder aFolderOfROOTFiles(HistsOutputFromRunningTriggerStudies) +# --outroot NameOfOutputRootFileWithIndividualGraphs +# --var name of the variable about wich the ROC plots are to be extracted, e.g., "PF_deepcsv" +import ROOT +import os +ROOT.gROOT.SetBatch(True) +import ROOTHelp.FancyROOTStyle + +# parser +from optparse import OptionParser +p = OptionParser() +p.add_option('--input', type = 'string', default = "outBTag.FTKBtagging.ttbar.mwt2.All.root", dest = 'inFile', help = 'intput file' ) +p.add_option('--inFolder', type = 'string', default = "", dest = 'inFolder', help = 'intput folder' ) +p.add_option('--output', type = 'string', default = "makeRocCurves", dest = 'outDir', help = 'output dir' ) +p.add_option('--outroot', type = 'string', default ='outputroot', dest = 'outrootName', help = 'output root file name') +p.add_option('--var', type = 'string', default ='PF_deepcsv', dest = 'varName', help = 'hist to be extracted') +p.add_option('--cmsText', type = 'string', default = "Work in Progress", help = '' ) +p.add_option('--doCaloJets', action="store_true", help = '' ) +(o,a) = p.parse_args() + +# store all files from folder in []files +path = o.inFolder +files = [] +for r, d, f in os.walk(path): # r=root, d=directories, f = files + for file in f: + files.append(os.path.join(r, file)) +print(files) + +if not os.path.exists(o.outDir): + os.mkdir(o.outDir) + +from rocCurveUtils import makeRoc +from JetLevelPlotUtils import getCMSText, getText + +outrootName = o.outrootName +varName = o.varName +hfile = ROOT.TFile(outrootName, "UPDATE") + +def makeRocPlot(inFile, name, var, bkg, sig, dir, varNorm=None,debug=False): + sigHist = inFile.Get(dir+"_"+sig+"/"+var) + bkgHist = inFile.Get(dir+"_"+bkg+"/"+var) + + if varNorm: + sigNormHist = inFile.Get(dir+"_"+sig+"/"+varNorm) + bkgNormHist = inFile.Get(dir+"_"+bkg+"/"+varNorm) + else : + sigNormHist = sigHist + bkgNormHist = bkgHist + + rocPlots = [] + for config in [("Rej",1,5e4),("Eff",5e-4,1)]: + temproc = makeRoc(sigHist, sigNormHist, bkgHist, bkgNormHist,doErr=False,bkgMode=config[0],cleanNoCut=True,debug=debug) + rocPlots.append(temproc) + can = ROOT.TCanvas(name+"_"+config[0], name+"_"+config[0]) + can.cd().SetLogy(1) + rocPlots[-1].SetLineWidth(5) + rocPlots[-1].GetXaxis().SetTitle("B-Jet Efficiency") + rocPlots[-1].GetXaxis().SetRangeUser(0.4,1) + if config[0] == "Rej": yTitle ="Light Flavor Rejection" + elif config[0] == "Eff": yTitle ="Light Flavor Efficiency" + rocPlots[-1].GetYaxis().SetTitle(yTitle) + rocPlots[-1].GetYaxis().SetRangeUser(config[1],config[2]) + rocPlots[-1].Draw("AL") + can.SaveAs(o.outDir+"/roc_"+name+"_"+config[0]+".pdf") + if name == varName: + hfile.cd() + ### the next line might need to change depending on the name of inFile + rocPlots[-1].Write(inFile.GetName()[:inFile.GetName().find(".")].replace("/","_")) # write the current rocPlot into + return rocPlots + +def plotSame(inFile, name,graphs,colors,styles, plotCaloJet=False, plotPFJet=False, plotOffJet=False,plotCSV=False,plotDeepCSV=False,workingPts= None,rocType=None): + + can = ROOT.TCanvas(name,name) + can.cd().SetLogy(1) + for gItr, g in enumerate(graphs): + g.SetLineColor(colors[gItr]) + g.SetLineStyle(styles[gItr]) + if not gItr: + g.Draw("AL") + else: + g.Draw("L") + + if not workingPts == None: + g_wrkPts = ROOT.TGraph(len(workingPts)) + g_wrkPts.SetMarkerSize(2) + g_wrkPts.SetMarkerColor(colors[1]) + g_wrkPts.SetMarkerStyle(34) + for wpItr, wp in enumerate(workingPts): + print wpItr,wp + + g_wrkPts.SetPoint(wpItr, wp[0],wp[1]) + + g_wrkPts.Draw("P") + + cmsLine1, cmsLine2 = getCMSText(xStart=0.2,yStart=0.875,subtext=o.cmsText) + cmsLine1.Draw("same") + cmsLine2.Draw("same") + + yStart = 0.75 + xStart = 0.2 + if rocType == "Rej": + xStart = 0.5 + yStart = 0.875 + + if plotOffJet: + offJetText = getText("Offline Jets ",xStart=xStart,yStart=yStart,size=0.04,color=ROOT.kBlack) + yStart = yStart - 0.05 + offJetText.Draw("same") + + if plotCaloJet: + caloJetText = getText("HLT Calo Jets",xStart=xStart,yStart=yStart,size=0.04,color=ROOT.kRed) + yStart = yStart - 0.05 + caloJetText.Draw("same") + if plotPFJet: + pfJetText = getText("HLT PF Jets ",xStart=xStart,yStart=yStart,size=0.04,color=ROOT.kBlue) + pfJetText.Draw("same") + + #offJetTextDeep = getText("Offline DeepCSV (Dashed) ",xStart=0.6,yStart=0.36,size=0.03,color=ROOT.kBlack) + + #offJetText = getText("Offline Jet ",xStart=0.6,yStart=0.4,size=0.03,color=ROOT.kBlack) + + yStart = 0.3 + xStart = 0.6 + if rocType == "Rej": + xStart = 0.2 + + + if plotDeepCSV: + if plotCSV: + deepCSVText = getText("DeepCSV (solid) ",xStart=xStart,yStart=yStart,size=0.04,color=ROOT.kBlack) + else: + deepCSVText = getText("DeepCSV",xStart=xStart,yStart=yStart,size=0.04,color=ROOT.kBlack) + deepCSVText.Draw("same") + yStart = yStart - 0.05 + + if plotCSV: + if plotDeepCSV: + CSVText = getText("CSV (dashed) ",xStart=xStart,yStart=yStart,size=0.04,color=ROOT.kBlack) + else: + CSVText = getText("CSV",xStart=xStart,yStart=yStart,size=0.04,color=ROOT.kBlack) + CSVText.Draw("same") + + + + + #offJetTextDeep.Draw("same") + + can.SaveAs(o.outDir+"/roc_"+name+".pdf") + + +# +# +# +def main(): + for f in files: + inFile = ROOT.TFile(f, "READ") + print(inFile) + print(inFile.GetName()) + off_deepcsv_roc = makeRocPlot(inFile, "Offline_deepcsv", "DeepCSV_l", bkg="matched_L",sig="matched_B",dir="offJets") + off_csv_roc = makeRocPlot(inFile, "Offline_csv", "CSVv2_l", bkg="matched_L",sig="matched_B",dir="offJets") + + pf_csv_roc = makeRocPlot(inFile, "PF_csv", "CSVv2_l", bkg="matchedJet_L",sig="matchedJet_B",dir="offJets") + pf_deepcsv_roc = makeRocPlot(inFile, "PF_deepcsv", "DeepCSV_l", bkg="matchedJet_L",sig="matchedJet_B",dir="offJets") + + + for i, rocType in enumerate(["Rej","Eff"]): + plotSame(inFile, "Off_vs_HLT_All_"+rocType, + [off_deepcsv_roc[i], off_csv_roc[i], pf_deepcsv_roc[i], pf_csv_roc[i]], + [ROOT.kBlack, ROOT.kBlack, ROOT.kBlue, ROOT.kBlue], + [ROOT.kSolid, ROOT.kDashed, ROOT.kSolid, ROOT.kDashed], + plotCaloJet = False, + plotPFJet = True, + plotOffJet = True, + plotCSV = True, + plotDeepCSV = True, + rocType = rocType + ) + +if __name__ == "__main__": + main() diff --git a/plotting/workingpoints_FixB.py b/plotting/workingpoints_FixB.py new file mode 100644 index 0000000..eaa62cd --- /dev/null +++ b/plotting/workingpoints_FixB.py @@ -0,0 +1,110 @@ +### take in a folder of root files with histograms, output pdfs with graph showing +### light flavor efficiency vs cut given certain B-jet efficiency +# --inFolder aFolderWithHistsOutputFromRunningTriggerStudies +# --output outPutDitractory +import ROOT +import os +ROOT.gROOT.SetBatch(True) +import ROOTHelp.FancyROOTStyle +from JetLevelPlotUtils import getCMSText, getText + +# parser +from optparse import OptionParser +p = OptionParser() +p.add_option('--inFolder', type = 'string', default = "", dest = 'inFolder', help = 'intput folder' ) +p.add_option('--output', type = 'string', default = "WorkingPoints", dest = 'outDir', help = 'output dir' ) +(o,a) = p.parse_args() + +# store all files from folder in []files +path = o.inFolder +files = [] +for r, d, f in os.walk(path): # r=root, d=directories, f = files + for file in f: + files.append(os.path.join(r, file)) +print "Read in files: ", files + +if not os.path.exists(o.outDir): + os.mkdir(o.outDir) + +def computebkgEff(inFile, var, bkg, sig, dir,workingpoint = 0.7, debug=False): + # input signal and background + sigNum = inFile.Get(dir+"_"+sig+"/"+var) + bkgNum = inFile.Get(dir+"_"+bkg+"/"+var) + sigDen = sigNum + bkgDen = bkgNum + nbins = sigNum.GetNbinsX() # number of bins + + sigTot = sigDen.Integral(0,nbins+1) # total signal + #print "total sig",sigTot + bkgTot = bkgDen.Integral(0,nbins+1) # total background + #print "total bkg",bkgTot + + for bin in range(nbins, 0, -1): + thisSig = abs(sigNum.Integral(bin,nbins+1)) + thisBkg = abs(bkgNum.Integral(bin,nbins+1)) + sigEff = thisSig/sigTot + bkgEff = thisBkg/bkgTot + if bkgEff: + bkgRej = 1.0/bkgEff + if sigEff > workingpoint: + cut = inFile.GetName()[:inFile.GetName().find(".")].replace("/","_") + #print "sigeff: ",sigEff, "bkgeff: ",bkgEff + return sigEff, bkgEff + +def makeWorkingPointsHists(workingpoint, cuts, bkgEffs): + c1 = ROOT.TCanvas("c1", "Light Flavor Efficiency as B-Jet Efficiency fixed to "+str(workingpoint)) + c1.SetGrid() + gr = ROOT.TGraph(len(cuts)) + #print("in make working points hists:", cuts, bkgEffs) + for i in range(len(cuts)): + gr.SetPoint(i,cuts[i],bkgEffs[i]) + gr.GetXaxis().SetTitle("TrackPtCut/GeV") + gr.GetYaxis().SetTitle("Light-Flavor Efficiency") + gr.SetLineColor( 2 ) + gr.SetLineWidth( 4 ) + gr.SetMarkerColor( 4 ) + gr.SetMarkerStyle( 21 ) + gr.SetTitle("Light Flavor Efficiency as B-Jet Efficiency fixed to "+str(workingpoint)) + gr.Draw("ALP") + cmsLine1, cmsLine2 = getCMSText(xStart=0.2,yStart=0.875,subtext="B-Jet Efficiency fixed to "+str(workingpoint)) + cmsLine1.Draw("same") + cmsLine2.Draw("same") + c1.Update() + c1.GetFrame().SetFillColor( 21 ) + c1.GetFrame().SetBorderSize( 12 ) + c1.Modified() + c1.Update() + c1.SaveAs(o.outDir+"/"+str(workingpoint*100)+"_high.pdf") + +def nametoCut(name): + if "0p4" in name: + return 0.4 + if "0p9" in name: + return 0.9 + if "1p5" in name: + return 1.5 + if "2p0" in name: + return 2.0 + if "2p5" in name: + return 2.5 + +def main(): + workingpoints = [0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95] + for workingpoint in workingpoints: + cuts = [] + sigEffs = [] + bkgEffs = [] + for f in files: + inFile = ROOT.TFile(f, "READ") + #print(inFile) + #print(inFile.GetName()) + sigEff, bkgEff = computebkgEff(inFile, "DeepCSV_l", bkg="matchedJet_L", sig="matchedJet_B", dir="offJets", workingpoint = workingpoint) + cut = nametoCut(f) + cuts.append(cut) + sigEffs.append(sigEff) + bkgEffs.append(bkgEff) + #print("at"+str(workingpoint)+" bkgEffs[]:", bkgEffs) + makeWorkingPointsHists(workingpoint, cuts, bkgEffs) + +if __name__ == "__main__": + main() diff --git a/plotting/workingpoints_FixLF.py b/plotting/workingpoints_FixLF.py new file mode 100644 index 0000000..d03fb65 --- /dev/null +++ b/plotting/workingpoints_FixLF.py @@ -0,0 +1,113 @@ +### take in a folder of root files with histograms, output pdfs with graph showing +### B-Jet efficiency vs cut given certain LightFlavor efficiency +# --inFolder aFolderWithHistsOutputFromRunningTriggerStudies +# --output outPutDitractory +import ROOT +import os +ROOT.gROOT.SetBatch(True) +import ROOTHelp.FancyROOTStyle +from JetLevelPlotUtils import getCMSText, getText + +# parser +from optparse import OptionParser +p = OptionParser() +p.add_option('--inFolder', type = 'string', default = "", dest = 'inFolder', help = 'intput folder' ) +p.add_option('--output', type = 'string', default = "makeRocCurves", dest = 'outDir', help = 'output dir' ) +(o,a) = p.parse_args() + +# store all files from folder in []files +path = o.inFolder +files = [] +for r, d, f in os.walk(path): # r=root, d=directories, f = files + for file in f: + files.append(os.path.join(r, file)) +print "Read in files: ", files + +if not os.path.exists(o.outDir): + os.mkdir(o.outDir) + +def computesigEff(inFile, var, bkg, sig, dir,workingpoint = 0.7, debug=False): + # input signal and background + sigNum = inFile.Get(dir+"_"+sig+"/"+var) + bkgNum = inFile.Get(dir+"_"+bkg+"/"+var) + sigDen = sigNum + bkgDen = bkgNum + nbins = sigNum.GetNbinsX() # number of bins + + sigTot = sigDen.Integral(0,nbins+1) # total signal + #print "total sig",sigTot + bkgTot = bkgDen.Integral(0,nbins+1) # total background + #print "total bkg",bkgTot + + for bin in range(nbins, 0, -1): + thisSig = abs(sigNum.Integral(bin,nbins+1)) + thisBkg = abs(bkgNum.Integral(bin,nbins+1)) + sigEff = thisSig/sigTot + bkgEff = thisBkg/bkgTot + if bkgEff: + bkgRej = 1.0/bkgEff + if bkgEff > workingpoint: + cut = inFile.GetName()[:inFile.GetName().find(".")].replace("/","_") + #print "sigeff: ",sigEff, "bkgeff: ",bkgEff + return sigEff, bkgEff + +def makeWorkingPointsHists(workingpoint, cuts, sigEffs): + c1 = ROOT.TCanvas("c1", "B-Jet Efficiency as Light-Flavor Efficiency fixed to "+str(workingpoint)) + c1.SetGrid() + gr = ROOT.TGraph(len(cuts)) + #print("in make working points hists:", cuts, bkgEffs) + for i in range(len(cuts)): + gr.SetPoint(i,cuts[i],sigEffs[i]) + gr.GetXaxis().SetTitle("TrackPtCut/GeV") + gr.GetYaxis().SetTitle("B-Jet Efficiency") + gr.SetLineColor( 2 ) + gr.SetLineWidth( 4 ) + gr.SetMarkerColor( 4 ) + gr.SetMarkerStyle( 21 ) + gr.SetTitle("B-Jet Efficiency as Light Flavor Efficiency fixed to "+str(workingpoint)) + gr.Draw("ALP") + cmsLine1, cmsLine2 = getCMSText(xStart=0.2,yStart=0.875,subtext="Light Flavor Efficiency fixed to "+str(workingpoint)) + cmsLine1.Draw("same") + cmsLine2.Draw("same") + c1.Update() + c1.GetFrame().SetFillColor( 21 ) + c1.GetFrame().SetBorderSize( 12 ) + c1.Modified() + c1.Update() + c1.SaveAs(o.outDir+"/"+str(workingpoint*100)+".pdf") + +def nametoCut(name): + if "0p4" in name: + return 0.4 + if "0p9" in name: + return 0.9 + if "1p5" in name: + return 1.5 + if "2p0" in name: + return 2.0 + if "2p5" in name: + return 2.5 + if "5p0" in name: + return 5.0 + +def main(): + workingpoints = [0.2, 0.1, 0.05, 0.01] + for workingpoint in workingpoints: + cuts = [] + sigEffs = [] + bkgEffs = [] + for f in files: + inFile = ROOT.TFile(f, "READ") + #print(inFile) + #print(inFile.GetName()) + sigEff, bkgEff = computesigEff(inFile, "DeepCSV_l", bkg="matchedJet_L", sig="matchedJet_B", dir="offJets", workingpoint = workingpoint) + cut = nametoCut(f) + cuts.append(cut) + sigEffs.append(sigEff) + bkgEffs.append(bkgEff) + #print("at"+str(workingpoint)+" bkgEffs[]:", bkgEffs) + makeWorkingPointsHists(workingpoint, cuts, sigEffs) + +if __name__ == "__main__": + main() +