Back to home page

sPhenix code displayed by LXR
 
 

    

File indexing completed on 2025-08-05 08:11:19 

0001 #! /usr/bin/env python
0002 from optparse import OptionParser
0003 import sys
0004 import os
0005 import datetime
0006 from array import *
0007 from ROOT import *
0008 import numpy
0009 import math
0010 import glob
0011 from plotUtil import *
0012 
0013 gROOT.SetBatch(True)
0014 
0015 TickSize = 0.03
0016 AxisTitleSize = 0.05
0017 AxisLabelSize = 0.04
0018 LeftMargin = 0.15
0019 RightMargin = 0.08
0020 TopMargin = 0.08
0021 BottomMargin = 0.13
0022 
0023 # with ratio in the second pad        
0024 def Draw_1Dhist_datasimcomp(hdata, hsims, gpadmargin, norm, logy, ymaxscale, XaxisName, Ytitle_unit, prelim, simlegtex, evtseltexts, outname):
0025     hsimcolor = colorset(len(hsims))
0026 
0027     hdata.Sumw2()
0028     for hsim in hsims:
0029         hsim.Sumw2()
0030 
0031     binwidth = hdata.GetXaxis().GetBinWidth(1)
0032 
0033     if norm == 'unity':
0034         hdata.Scale(1. / hdata.Integral(-1, -1))
0035         for hsim in hsims:
0036             hsim.Scale(1. / hsim.Integral(-1, -1))
0037     elif norm == 'data':
0038         for hsim in hsims:
0039             hsim.Scale(hdata.Integral(-1, -1) / hsim.Integral(-1, -1))
0040     else:
0041         if norm != 'none':
0042             print('Invalid normalization option: {}'.format(norm))
0043             sys.exit(1)
0044     
0045     # Get the maximum bin content 
0046     maxbincontent = max(hdata.GetMaximum(), hsim.GetMaximum())
0047 
0048     c = TCanvas('c', 'c', 800, 700)
0049     pad1 = TPad( 'pad1', ' ', 0, 0.3, 1, 1)
0050     pad2 = TPad( 'pad2', ' ', 0, 0, 1, 0.3)
0051     pad1.SetRightMargin(gpadmargin[0])
0052     pad1.SetTopMargin(gpadmargin[1])
0053     pad1.SetLeftMargin(gpadmargin[2])
0054     pad1.SetBottomMargin(0.0)
0055     pad1.Draw()
0056     pad2.SetRightMargin(gpadmargin[0])
0057     pad2.SetTopMargin(0.0)
0058     pad2.SetLeftMargin(gpadmargin[2])
0059     pad2.SetBottomMargin(gpadmargin[3])
0060     pad2.Draw() # Draw the TPad on the TCanvas before plotting something on TPad
0061     # cd to the pad1
0062     pad1.cd()
0063     if logy:
0064         pad1.SetLogy()
0065     
0066     for i, hsim in enumerate(hsims):
0067         if i == 0:
0068             if norm == 'unity':
0069                 if Ytitle_unit == '':
0070                     hsim.GetYaxis().SetTitle(
0071                         'Normalized entries / ({:g})'.format(binwidth))
0072                 else:
0073                     hsim.GetYaxis().SetTitle(
0074                         'Normalized entries / ({:g} {unit})'.format(binwidth, unit=Ytitle_unit))
0075             else:
0076                 if Ytitle_unit == '':
0077                     hsim.GetYaxis().SetTitle('Entries / ({:g})'.format(binwidth))
0078                 else:
0079                     hsim.GetYaxis().SetTitle(
0080                         'Entries / ({:g} {unit})'.format(binwidth, unit=Ytitle_unit))
0081 
0082             if logy:
0083                 hsim.GetYaxis().SetRangeUser(hsim.GetMinimum(0)*0.5, (hsim.GetMaximum()) * ymaxscale)
0084             else:
0085                 hsim.GetYaxis().SetRangeUser(0., (hsim.GetMaximum()) * ymaxscale)
0086 
0087             hsim.GetXaxis().SetTitle(XaxisName)
0088             hsim.GetXaxis().SetTitleOffset(1.1)
0089             hsim.GetYaxis().SetTitleOffset(1.35)
0090             hsim.SetLineColor(TColor.GetColor(hsimcolor[i]))
0091             hsim.SetLineWidth(2)
0092             hsim.SetMarkerSize(0)
0093             hsim.Draw('histe')
0094         else:
0095             hsim.SetLineColor(TColor.GetColor(hsimcolor[i]))
0096             hsim.SetLineWidth(2)
0097             hsim.SetMarkerSize(0)
0098             hsim.Draw('histe same')
0099 
0100     hdata.SetMarkerStyle(20)
0101     hdata.SetMarkerSize(1)
0102     hdata.SetMarkerColor(1)
0103     hdata.SetLineColor(1)
0104     hdata.SetLineWidth(2)
0105     hdata.Draw('same PE1')
0106     shift = 0.45 if prelim else 0.75
0107     legylow = 0.2 + 0.04 * (len(hsims) - 1)
0108     leg = TLegend((1-RightMargin)-shift, (1-TopMargin)-legylow,
0109                   (1-RightMargin)-0.1, (1-TopMargin)-0.03)
0110     leg.SetTextSize(0.04)
0111     leg.SetFillStyle(0)
0112     # prelimtext = 'Preliminary' if prelim else 'Work-in-progress'
0113     prelimtext = 'Preliminary' if prelim else 'Internal'
0114     leg.AddEntry('', '#it{#bf{sPHENIX}} '+prelimtext, '')
0115     leg.AddEntry('', 'Au+Au #sqrt{s_{NN}}=200 GeV', '')
0116     leg.AddEntry(hdata, 'Data', "PE1")
0117     for i, lt in enumerate(simlegtex):
0118         leg.AddEntry(hsims[i], lt, "le")
0119     leg.Draw('same')
0120     
0121     # event selection text
0122     legylow_evtselshift = 0.04 * len(evtseltexts)
0123     leg2 = TLegend((1-RightMargin)-shift, (1-TopMargin)-legylow-legylow_evtselshift,
0124                     (1-RightMargin)-0.1, (1-TopMargin)-legylow)
0125     leg2.SetTextSize(0.035)
0126     leg2.SetFillStyle(0)
0127     for i, evtseltext in enumerate(evtseltexts):
0128         leg2.AddEntry('', evtseltext, '')
0129     leg2.Draw('same')
0130     c.RedrawAxis()
0131     
0132     c.Update()
0133     # cd to the pad2
0134     textscale = 2.7
0135     pad2.cd()
0136     pad2.SetLogy(0)
0137     pad2.SetGridy(0)
0138     # take the ratio = data/simulation
0139     l_hM_ratio = []
0140     for i, hsim in enumerate(hsims):
0141         hM_ratio = hdata.Clone('hM_ratio'+str(i))
0142         hM_ratio.Divide(hsim)
0143         hM_ratio.SetMarkerStyle(20)
0144         hM_ratio.SetMarkerSize(0.7)
0145         hM_ratio.SetMarkerColor(TColor.GetColor(hsimcolor[i]))
0146         hM_ratio.SetLineColor(TColor.GetColor(hsimcolor[i]))
0147         hM_ratio.SetLineWidth(2)
0148         if i == 0:
0149             hM_ratio.GetYaxis().SetRangeUser(0.001, 1.999)
0150             hM_ratio.GetYaxis().SetNdivisions(505)
0151             hM_ratio.GetYaxis().SetTitle('Data/Sim.')
0152             hM_ratio.GetYaxis().SetTitleOffset(0.5)
0153             hM_ratio.GetYaxis().SetTickSize(TickSize)
0154             hM_ratio.GetYaxis().SetTitleSize(AxisTitleSize*textscale)
0155             hM_ratio.GetYaxis().SetLabelSize(AxisLabelSize*textscale)
0156             hM_ratio.GetXaxis().SetTickSize(TickSize*textscale)
0157             hM_ratio.GetXaxis().SetTitleSize(AxisTitleSize*textscale)
0158             hM_ratio.GetXaxis().SetLabelSize(AxisLabelSize*textscale)
0159             hM_ratio.GetXaxis().SetTitleOffset(1.1)
0160             hM_ratio.GetXaxis().SetTitle(XaxisName)
0161             hM_ratio.Draw('PE1X0')
0162         else:
0163             hM_ratio.Draw('PE1X0 same')
0164             
0165         l_hM_ratio.append(hM_ratio)
0166     
0167     c.Draw()
0168     c.SaveAs(outname+'.pdf')
0169     c.SaveAs(outname+'.png')
0170     if(c):
0171         c.Close()
0172         gSystem.ProcessEvents()
0173         del c
0174         c = 0
0175 
0176 if __name__ == '__main__':
0177     plotpath = './CombinatorialBkg/'
0178     os.makedirs(plotpath, exist_ok=True)
0179 
0180     histbasepath = '/sphenix/user/hjheng/sPHENIXRepo/analysis/dNdEta_Run2023/analysis_INTT/plot/hists'
0181     datahistpath = 'Data_CombinedNtuple_Run20869_HotDead_BCO_ADC_Survey/dRcut0p5_NominalVtxZ_RandomClusSet0_clusAdcCutSet0'
0182     simhistspath_randclusset = ['Sim_Ntuple_HIJING_new_20240424/dRcut0p5_NominalVtxZ_RandomClusSet{}_clusAdcCutSet0'.format(i) for i in range(5)]
0183     
0184     if os.path.isfile("{}/{}/hists_merged.root".format(histbasepath, datahistpath)):
0185         os.system("rm {}/{}/hists_merged.root".format(histbasepath, datahistpath))
0186         os.system("hadd -f -j 20 {}/{}/hists_merged.root {}/{}/hists_*.root".format(histbasepath, datahistpath, histbasepath, datahistpath))
0187     else:
0188         os.system("hadd -f -j 20 {}/{}/hists_merged.root {}/{}/hists_*.root".format(histbasepath, datahistpath, histbasepath, datahistpath))
0189         
0190     for i, simhistspath in enumerate(simhistspath_randclusset):
0191         if os.path.isfile("{}/{}/hists_merged.root".format(histbasepath, simhistspath)):
0192             os.system("rm {}/{}/hists_merged.root".format(histbasepath, simhistspath))
0193             os.system("hadd -f -j 20 {}/{}/hists_merged.root {}/{}/hists_*.root".format(histbasepath, simhistspath, histbasepath, simhistspath))
0194         else:
0195             os.system("hadd -f -j 20 {}/{}/hists_merged.root {}/{}/hists_*.root".format(histbasepath, simhistspath, histbasepath, simhistspath))
0196         
0197     hM_dR_reco_Centrality0to70_MBDAsymLe0p75_VtxZm30tom10_data = GetHistogram(histbasepath+'/'+datahistpath+'/hists_merged.root', 'hM_dR_reco_Centrality0to70_MBDAsymLe0p75_VtxZm30tom10')
0198     l_hM_dR_reco_Centrality0to70_MBDAsymLe0p75_VtxZm30tom10_sim = [GetHistogram(histbasepath+'/'+simhistspath+'/hists_merged.root', 'hM_dR_reco_Centrality0to70_MBDAsymLe0p75_VtxZm30tom10') for simhistspath in simhistspath_randclusset]
0199     hM_dR_reco_altrange_Centrality0to70_MBDAsymLe0p75_VtxZm30tom10_data = GetHistogram(histbasepath+'/'+datahistpath+'/hists_merged.root', 'hM_dR_reco_altrange_Centrality0to70_MBDAsymLe0p75_VtxZm30tom10')
0200     l_hM_dR_reco_altrange_Centrality0to70_MBDAsymLe0p75_VtxZm30tom10_sim = [GetHistogram(histbasepath+'/'+simhistspath+'/hists_merged.root', 'hM_dR_reco_altrange_Centrality0to70_MBDAsymLe0p75_VtxZm30tom10') for simhistspath in simhistspath_randclusset]
0201 
0202     # Draw_1Dhist_datasimcomp(hdata, hsims, gpadmargin, norm, logy, ymaxscale, XaxisName, Ytitle_unit, prelim, simlegtex, evtseltexts, outname):
0203     Draw_1Dhist_datasimcomp(hM_dR_reco_Centrality0to70_MBDAsymLe0p75_VtxZm30tom10_data, l_hM_dR_reco_Centrality0to70_MBDAsymLe0p75_VtxZm30tom10_sim, [0.05,0.06,0.15,0.32], 'data', True, 100, '#DeltaR_{reco}', '', False, ['Simulation, baseline', 'Simulation, +1% random clusters', 'Simulation, +3% random clusters', 'Simulation, +5% random clusters', 'Simulation, +10% random clusters'], ['Centrality 0-70%, Cluster ADC > 35, |Asymm._{MBD}|#leq0.75, -30#leqZ_{vtx}#leq-10[cm]'], plotpath+'dR_reco_Centrality0to70_MBDAsymLe0p75_VtxZm30tom10_randcluscomp')
0204     Draw_1Dhist_datasimcomp(hM_dR_reco_altrange_Centrality0to70_MBDAsymLe0p75_VtxZm30tom10_data, l_hM_dR_reco_altrange_Centrality0to70_MBDAsymLe0p75_VtxZm30tom10_sim, [0.05,0.06,0.15,0.32], 'data', False, 2, '#DeltaR_{reco}', '', False, ['Simulation, baseline', 'Simulation, +1% random clusters', 'Simulation, +3% random clusters', 'Simulation, +5% random clusters', 'Simulation, +10% random clusters'], ['Centrality 0-70%, Cluster ADC > 35, |Asymm._{MBD}|#leq0.75, -30#leqZ_{vtx}#leq-10[cm]'], plotpath+'dR_reco_altrange_Centrality0to70_MBDAsymLe0p75_VtxZm30tom10_randcluscomp')
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!