VHbb/python/plot_systematics.py

#!/usr/bin/env python
import ROOT 
ROOT.gROOT.SetBatch(True)
from ROOT import TFile
from optparse import OptionParser
import sys
from myutils import BetterConfigParser, TdrStyles, getRatio


argv = sys.argv
parser = OptionParser()
parser.add_option("-C", "--config", dest="config", default=[], action="append",
                      help="configuration file")
(opts, args) = parser.parse_args(argv)
config = BetterConfigParser()
config.read(opts.config)


#---------- yes, this is not in the config yet---------
mode = 'BDT'
xMin=-1
xMax=1
masses = ['125']
Abins = ['HighPt','LowPt']
channels= ['Zee','Zmm']
#------------------------------------------------------
#---------- Mjj ---------------------------------------
#mode = 'Mjj'
#xMin=0
#xMax=255
#masses = ['125']
#Abins = ['highPt','lowPt','medPt']
#channels= ['Zee','Zmm']
#------------------------------------------------------

path = config.get('Directories','limits')
outpath = config.get('Directories','plotpath')

setup = eval(config.get('LimitGeneral','setup'))
Dict = eval(config.get('LimitGeneral','Dict'))
MCs = [Dict[s] for s in setup]

sys_BDT= eval(config.get('LimitGeneral','sys_BDT'))
systematicsnaming = eval(config.get('LimitGeneral','systematicsnaming'))
systs=[systematicsnaming[s] for s in sys_BDT]

if eval(config.get('LimitGeneral','weightF_sys')): systs.append('UEPS')

def myText(txt="CMS Preliminary",ndcX=0,ndcY=0,size=0.8):
    ROOT.gPad.Update()
    text = ROOT.TLatex()
    text.SetNDC()
    text.SetTextColor(ROOT.kBlack)
    text.SetTextSize(text.GetTextSize()*size)
    text.DrawLatex(ndcX,ndcY,txt)
    return text


#for mass in ['110','115','120','125','130','135']:
for mass in masses:
    for Abin in Abins:
        for channel in channels:

            if mode == 'BDT':
                input = TFile.Open(path+'/vhbb_TH_BDT_M'+mass+'_'+channel+Abin+'_8TeV.root','read')
            if mode == 'Mjj':
                input = TFile.Open(path+'/vhbb_TH_Mjj_'+Abin+'_M'+mass+'_'+channel+'.root','read')

            for MC in MCs:
                print MC
                for syst in systs:
                    print syst
                #['CMS_res_j','CMS_scale_j','CMS_eff_b','CMS_fake_b_8TeV','UEPS']:
                #for syst in ['CMS_vhbb_stats_']:


                    TdrStyles.tdrStyle()

                    c = ROOT.TCanvas('','', 600, 600)
                    c.SetFillStyle(4000)
                    c.SetFrameFillStyle(1000)
                    c.SetFrameFillColor(0)
                    oben = ROOT.TPad('oben','oben',0,0.3 ,1.0,1.0)
                    oben.SetBottomMargin(0)
                    oben.SetFillStyle(4000)
                    oben.SetFrameFillStyle(1000)
                    oben.SetFrameFillColor(0)
                    unten = ROOT.TPad('unten','unten',0,0.0,1.0,0.3)
                    unten.SetTopMargin(0.)
                    unten.SetBottomMargin(0.35)
                    unten.SetFillStyle(4000)
                    unten.SetFrameFillStyle(1000)
                    unten.SetFrameFillColor(0)
                    oben.Draw()
                    unten.Draw()
                    oben.cd()

                    ROOT.gPad.SetTicks(1,1)


                    input.cd(channel+Abin+'_8TeV')
                    Ntotal=ROOT.gDirectory.Get(MC)
                    Utotal=ROOT.gDirectory.Get(MC+syst+'Up')
                    #Utotal=input.Get(MC+syst+MC+'_'+channel+'Up')
                    Dtotal=ROOT.gDirectory.Get(MC+syst+'Down')
                    #Dtotal=input.Get(MC+syst+MC+'_'+channel+'Down')
                    l = ROOT.TLegend(0.17, 0.8, 0.37, 0.65)
                    
                    l.SetLineWidth(2)
                    l.SetBorderSize(0)
                    l.SetFillColor(0)
                    l.SetFillStyle(4000)
                    l.SetTextFont(62)
                    l.SetTextSize(0.035)

                    
                    l.AddEntry(Ntotal,'nominal','PL')
                    l.AddEntry(Utotal,'up','PL')
                    l.AddEntry(Dtotal,'down','PL')
                    Ntotal.GetYaxis().SetRangeUser(0,1.5*Ntotal.GetBinContent(Ntotal.GetMaximumBin()))
                    Ntotal.SetMarkerStyle(8)
                    Ntotal.SetLineColor(1)
                    Ntotal.SetStats(0)
                    Ntotal.SetTitle(MC +' '+syst)
                    Ntotal.Draw("P0")
                    Ntotal.Draw("same")
                    Utotal.SetLineColor(4)    
                    Utotal.SetLineStyle(4)
                    Utotal.SetLineWidth(2)        
                    Utotal.Draw("same hist")
                    Dtotal.SetLineColor(2)
                    Dtotal.SetLineStyle(3)
                    Dtotal.SetLineWidth(2)  
                    Dtotal.Draw("same hist")
                    l.SetFillColor(0)
                    l.SetBorderSize(0)
                    l.Draw()
                    
                    title=myText('Shape Systematic %s in %s'%(syst,MC),0.17,0.85)
                    
                    print 'Shape Systematic %s in %s'%(syst,MC)
                    print 'Up:     \t%s'%Utotal.Integral()
                    print 'Nominal:\t%s'%Ntotal.Integral()
                    print 'Down:   \t%s'%Dtotal.Integral()
                    
                    unten.cd()
                    ROOT.gPad.SetTicks(1,1)

                    ratioU, errorU  = getRatio(Utotal,Ntotal,xMin,xMax)
                    ratioD, errorD  = getRatio(Dtotal,Ntotal,xMin,xMax)

                    ksScoreU = Ntotal.KolmogorovTest( Utotal )
                    chiScoreU = Ntotal.Chi2Test( Utotal , "WWCHI2/NDF")
                    ksScoreD = Ntotal.KolmogorovTest( Dtotal )
                    chiScoreD = Ntotal.Chi2Test( Dtotal , "WWCHI2/NDF")


                    ratioU.SetStats(0)
                    ratioU.SetMinimum(0.01)
                    ratioU.SetMaximum(2.49)
                    ratioU.GetYaxis().SetNdivisions(505)
                    #ratioU.GetYaxis().SetLabelSize(0.2)
                    #ratioU.GetYaxis().SetTitleSize(0.2)
                    #ratioU.GetYaxis().SetTitleOffset(0.2)
                    #ratioU.GetXaxis().SetLabelColor(10)
                    ratioU.SetLineColor(4)    
                    ratioU.SetLineStyle(4)
                    ratioU.SetLineWidth(2)
                    ratioU.GetXaxis().SetTitle('BDT output')
                    ratioU.GetYaxis().SetTitle('Ratio') 
                    ratioU.Draw("hist")
                    fitRatioU = ratioU.Fit("pol2","S")
                    fitRatioU.Draw("SAME")
                    ratioU.SetTitle("")
                    ratioD.SetStats(0)
                    ratioD.GetYaxis().SetRangeUser(0.5,1.5)
                    ratioD.GetYaxis().SetNdivisions(502,0)
                    #ratioD.GetYaxis().SetLabelSize(0.2)
                    #ratioD.GetYaxis().SetTitleSize(0.2)
                    #ratioD.GetYaxis().SetTitleOffset(0.2)
                    #ratioD.GetXaxis().SetLabelColor(10)
                    ratioD.SetLineColor(2)
                    ratioD.SetLineStyle(3)
                    ratioD.SetLineWidth(2)  
                    ratioD.Draw("hist same")
                    fitRatioD = ratioD.Fit("pol2","S")
                    fitRatioD.Draw("SAME")
                    ratioD.SetTitle("")
                    m_one_line = ROOT.TLine(xMin,1,xMax,1)
                    m_one_line.SetLineStyle(7)
                    m_one_line.SetLineColor(4)
                    m_one_line.Draw("Same")


                    #name = outpath+Abin+'_M'+mass+'_'+channel+'_'+MC+syst+'.png'
                    #c.Print(name)
                    name = outpath+'systPlot_'+Abin+'_M'+mass+'_'+channel+'_'+MC+syst+'.pdf'
                    c.Print(name)


            input.Close()
Revision:	1.9
Committed:	Fri Jul 19 10:13:50 2013 UTC (11 years, 9 months ago) by nmohr
Content type:	text/x-python
Branch:	MAIN
CVS Tags:	HEAD
Changes since 1.8:	+6 -3 lines
Log Message:	For git migration
#	User	Rev	Content
1	peller	1.1	#!/usr/bin/env python
2			import ROOT
3	nmohr	1.9	ROOT.gROOT.SetBatch(True)
4	peller	1.1	from ROOT import TFile
5	peller	1.3	from optparse import OptionParser
6			import sys
7	peller	1.8	from myutils import BetterConfigParser, TdrStyles, getRatio
8	peller	1.3
9	nmohr	1.6
10	peller	1.3	argv = sys.argv
11			parser = OptionParser()
12			parser.add_option("-C", "--config", dest="config", default=[], action="append",
13			help="configuration file")
14			(opts, args) = parser.parse_args(argv)
15			config = BetterConfigParser()
16			config.read(opts.config)
17	peller	1.1
18
19	peller	1.3	#---------- yes, this is not in the config yet---------
20	peller	1.8	mode = 'BDT'
21			xMin=-1
22			xMax=1
23			masses = ['125']
24			Abins = ['HighPt','LowPt']
25			channels= ['Zee','Zmm']
26			#------------------------------------------------------
27			#---------- Mjj ---------------------------------------
28			#mode = 'Mjj'
29			#xMin=0
30			#xMax=255
31	peller	1.5	#masses = ['125']
32	peller	1.8	#Abins = ['highPt','lowPt','medPt']
33	peller	1.5	#channels= ['Zee','Zmm']
34			#------------------------------------------------------
35	peller	1.3
36	peller	1.8	path = config.get('Directories','limits')
37			outpath = config.get('Directories','plotpath')
38	peller	1.3
39			setup = eval(config.get('LimitGeneral','setup'))
40			Dict = eval(config.get('LimitGeneral','Dict'))
41			MCs = [Dict[s] for s in setup]
42
43			sys_BDT= eval(config.get('LimitGeneral','sys_BDT'))
44	peller	1.8	systematicsnaming = eval(config.get('LimitGeneral','systematicsnaming'))
45			systs=[systematicsnaming[s] for s in sys_BDT]
46	peller	1.1
47	peller	1.3	if eval(config.get('LimitGeneral','weightF_sys')): systs.append('UEPS')
48	peller	1.2
49	peller	1.3	def myText(txt="CMS Preliminary",ndcX=0,ndcY=0,size=0.8):
50			ROOT.gPad.Update()
51			text = ROOT.TLatex()
52			text.SetNDC()
53			text.SetTextColor(ROOT.kBlack)
54			text.SetTextSize(text.GetTextSize()*size)
55			text.DrawLatex(ndcX,ndcY,txt)
56			return text
57	peller	1.2
58
59	peller	1.3	#for mass in ['110','115','120','125','130','135']:
60			for mass in masses:
61			for Abin in Abins:
62			for channel in channels:
63
64	peller	1.5	if mode == 'BDT':
65			input = TFile.Open(path+'/vhbb_TH_BDT_M'+mass+'_'+channel+Abin+'_8TeV.root','read')
66			if mode == 'Mjj':
67			input = TFile.Open(path+'/vhbb_TH_Mjj_'+Abin+'_M'+mass+'_'+channel+'.root','read')
68	peller	1.3
69			for MC in MCs:
70	peller	1.8	print MC
71	peller	1.3	for syst in systs:
72	peller	1.8	print syst
73	peller	1.3	#['CMS_res_j','CMS_scale_j','CMS_eff_b','CMS_fake_b_8TeV','UEPS']:
74	peller	1.2	#for syst in ['CMS_vhbb_stats_']:
75
76	peller	1.3
77	peller	1.8	TdrStyles.tdrStyle()
78	peller	1.3
79			c = ROOT.TCanvas('','', 600, 600)
80			c.SetFillStyle(4000)
81			c.SetFrameFillStyle(1000)
82			c.SetFrameFillColor(0)
83			oben = ROOT.TPad('oben','oben',0,0.3 ,1.0,1.0)
84			oben.SetBottomMargin(0)
85			oben.SetFillStyle(4000)
86			oben.SetFrameFillStyle(1000)
87			oben.SetFrameFillColor(0)
88			unten = ROOT.TPad('unten','unten',0,0.0,1.0,0.3)
89			unten.SetTopMargin(0.)
90			unten.SetBottomMargin(0.35)
91			unten.SetFillStyle(4000)
92			unten.SetFrameFillStyle(1000)
93			unten.SetFrameFillColor(0)
94	peller	1.2	oben.Draw()
95			unten.Draw()
96			oben.cd()
97	peller	1.3
98	peller	1.2	ROOT.gPad.SetTicks(1,1)
99
100
101	peller	1.8	input.cd(channel+Abin+'_8TeV')
102			Ntotal=ROOT.gDirectory.Get(MC)
103			Utotal=ROOT.gDirectory.Get(MC+syst+'Up')
104	peller	1.2	#Utotal=input.Get(MC+syst+MC+'_'+channel+'Up')
105	peller	1.8	Dtotal=ROOT.gDirectory.Get(MC+syst+'Down')
106	peller	1.2	#Dtotal=input.Get(MC+syst+MC+'_'+channel+'Down')
107	peller	1.3	l = ROOT.TLegend(0.17, 0.8, 0.37, 0.65)
108
109			l.SetLineWidth(2)
110			l.SetBorderSize(0)
111			l.SetFillColor(0)
112			l.SetFillStyle(4000)
113			l.SetTextFont(62)
114			l.SetTextSize(0.035)
115
116
117	peller	1.2	l.AddEntry(Ntotal,'nominal','PL')
118			l.AddEntry(Utotal,'up','PL')
119			l.AddEntry(Dtotal,'down','PL')
120	peller	1.3	Ntotal.GetYaxis().SetRangeUser(0,1.5*Ntotal.GetBinContent(Ntotal.GetMaximumBin()))
121	peller	1.2	Ntotal.SetMarkerStyle(8)
122			Ntotal.SetLineColor(1)
123			Ntotal.SetStats(0)
124			Ntotal.SetTitle(MC +' '+syst)
125			Ntotal.Draw("P0")
126			Ntotal.Draw("same")
127			Utotal.SetLineColor(4)
128			Utotal.SetLineStyle(4)
129			Utotal.SetLineWidth(2)
130			Utotal.Draw("same hist")
131			Dtotal.SetLineColor(2)
132			Dtotal.SetLineStyle(3)
133			Dtotal.SetLineWidth(2)
134			Dtotal.Draw("same hist")
135			l.SetFillColor(0)
136			l.SetBorderSize(0)
137			l.Draw()
138
139	peller	1.3	title=myText('Shape Systematic %s in %s'%(syst,MC),0.17,0.85)
140
141			print 'Shape Systematic %s in %s'%(syst,MC)
142			print 'Up: \t%s'%Utotal.Integral()
143			print 'Nominal:\t%s'%Ntotal.Integral()
144			print 'Down: \t%s'%Dtotal.Integral()
145	peller	1.2
146			unten.cd()
147			ROOT.gPad.SetTicks(1,1)
148
149	peller	1.3	ratioU, errorU = getRatio(Utotal,Ntotal,xMin,xMax)
150			ratioD, errorD = getRatio(Dtotal,Ntotal,xMin,xMax)
151
152			ksScoreU = Ntotal.KolmogorovTest( Utotal )
153			chiScoreU = Ntotal.Chi2Test( Utotal , "WWCHI2/NDF")
154			ksScoreD = Ntotal.KolmogorovTest( Dtotal )
155			chiScoreD = Ntotal.Chi2Test( Dtotal , "WWCHI2/NDF")
156
157
158
159	peller	1.2	ratioU.SetStats(0)
160	nmohr	1.4	ratioU.SetMinimum(0.01)
161			ratioU.SetMaximum(2.49)
162			ratioU.GetYaxis().SetNdivisions(505)
163	peller	1.3	#ratioU.GetYaxis().SetLabelSize(0.2)
164			#ratioU.GetYaxis().SetTitleSize(0.2)
165			#ratioU.GetYaxis().SetTitleOffset(0.2)
166			#ratioU.GetXaxis().SetLabelColor(10)
167	peller	1.2	ratioU.SetLineColor(4)
168			ratioU.SetLineStyle(4)
169	peller	1.3	ratioU.SetLineWidth(2)
170			ratioU.GetXaxis().SetTitle('BDT output')
171			ratioU.GetYaxis().SetTitle('Ratio')
172	peller	1.2	ratioU.Draw("hist")
173	nmohr	1.9	fitRatioU = ratioU.Fit("pol2","S")
174			fitRatioU.Draw("SAME")
175	peller	1.2	ratioU.SetTitle("")
176			ratioD.SetStats(0)
177			ratioD.GetYaxis().SetRangeUser(0.5,1.5)
178			ratioD.GetYaxis().SetNdivisions(502,0)
179	peller	1.3	#ratioD.GetYaxis().SetLabelSize(0.2)
180			#ratioD.GetYaxis().SetTitleSize(0.2)
181			#ratioD.GetYaxis().SetTitleOffset(0.2)
182			#ratioD.GetXaxis().SetLabelColor(10)
183	peller	1.2	ratioD.SetLineColor(2)
184			ratioD.SetLineStyle(3)
185			ratioD.SetLineWidth(2)
186			ratioD.Draw("hist same")
187	nmohr	1.9	fitRatioD = ratioD.Fit("pol2","S")
188			fitRatioD.Draw("SAME")
189	peller	1.2	ratioD.SetTitle("")
190			m_one_line = ROOT.TLine(xMin,1,xMax,1)
191			m_one_line.SetLineStyle(7)
192			m_one_line.SetLineColor(4)
193			m_one_line.Draw("Same")
194
195
196			#name = outpath+Abin+'_M'+mass+'_'+channel+'_'+MC+syst+'.png'
197			#c.Print(name)
198	nmohr	1.9	name = outpath+'systPlot_'+Abin+'_M'+mass+'_'+channel+'_'+MC+syst+'.pdf'
199	peller	1.2	c.Print(name)
200	peller	1.1
201
202	peller	1.3	input.Close()