VHbb/python/gethistofromtree.py

from samplesclass import sample
from printcolor import printc
import pickle
import ROOT 
from ROOT import TFile, TTree
import ROOT
from array import array
from BetterConfigParser import BetterConfigParser
import sys


def getScale(job,path,config,rescale,subsample=-1):
    anaTag=config.get('Analysis','tag')
    input = TFile.Open(path+'/'+job.getpath())
    CountWithPU = input.Get("CountWithPU")
    CountWithPU2011B = input.Get("CountWithPU2011B")
    #print lumi*xsecs[i]/hist.GetBinContent(1)
    
    if subsample>-1:
        if type(job.xsec[subsample]) == str: xsec=float(eval(job.xsec[subsample]))
        else: xsec=float(job.xsec[subsample])
        sf=float(job.sf[subsample])
    else:
        if type(job.xsec) == str: xsec=float(eval(job.xsec))
        else: xsec=float(job.xsec)
        sf=float(job.sf)
    
    
    theScale = 1.
    if anaTag == '7TeV':
        theScale = float(job.lumi)*xsec*sf/(0.46502*CountWithPU.GetBinContent(1)+0.53498*CountWithPU2011B.GetBinContent(1))*rescale/float(job.split)
    elif anaTag == '8TeV':
        theScale = float(job.lumi)*xsec*sf/(CountWithPU.GetBinContent(1))*rescale/float(job.split)
    return theScale 

def getHistoFromTree(job,path,config,options,rescale=1,subsample=-1,which_weightF='weightF'):

    #print job.getpath()
    #print options
    treeVar=options[0]
    if subsample>-1:
        name=job.subnames[subsample]
        group=job.group[subsample]
    else:
        name=job.name
        group=job.group

    #title=job.plotname()
    nBins=int(options[3])
    xMin=float(options[4])
    xMax=float(options[5])
    #addOverFlow=eval(config.get('Plot_general','addOverFlow'))
    addOverFlow = False

    TrainFlag = eval(config.get('Analysis','TrainFlag'))
    if TrainFlag: traincut = " & EventForTraining == 0"
    if not TrainFlag: traincut=""

    if job.type != 'DATA':
    
        if type(options[7])==str:
            cutcut=config.get('Cuts',options[7])
        elif type(options[7])==list:
            cutcut=config.get('Cuts',options[7][0])
            cutcut=cutcut.replace(options[7][1],options[7][2])
            print cutcut
        if subsample>-1:
            treeCut='%s & %s%s'%(cutcut,job.subcuts[subsample],traincut)        
        else:
            treeCut='%s%s'%(cutcut,traincut)

    elif job.type == 'DATA':
        cutcut=config.get('Cuts',options[8])
        treeCut='%s'%(cutcut)


    input = TFile.Open(path+'/'+job.getpath(),'read')

    Tree = input.Get(job.tree)
    #Tree=tmpTree.CloneTree()
    #Tree.SetDirectory(0)
    
    #Tree=tmpTree.Clone()
    weightF=config.get('Weights',which_weightF)
    #hTree = ROOT.TH1F('%s'%name,'%s'%title,nBins,xMin,xMax)
    #hTree.SetDirectory(0)
    #hTree.Sumw2()
    #print 'drawing...'
    if job.type != 'DATA':
        #print treeCut
        #print job.name
        if Tree.GetEntries():
            Tree.Draw('%s>>%s(%s,%s,%s)' %(treeVar,name,nBins,xMin,xMax),'(%s)*(%s)' %(treeCut,weightF), "goff,e")
            full=True
        else:
            full=False
    elif job.type == 'DATA':
    
        if len(options)>10:
            if options[11] == 'blind':
                treeCut = treeCut + '&'+treeVar+'<0'
    
    
        Tree.Draw('%s>>%s(%s,%s,%s)' %(treeVar,name,nBins,xMin,xMax),treeCut, "goff,e")
        full = True
    if full:
        hTree = ROOT.gDirectory.Get(name)
    else:
        hTree = ROOT.TH1F('%s'%name,'%s'%name,nBins,xMin,xMax)
        hTree.Sumw2()
    #print job.name + ' Sumw2', hTree.GetEntries()

    if job.type != 'DATA':
        ScaleFactor = getScale(job,path,config,rescale,subsample)
        if ScaleFactor != 0:
            hTree.Scale(ScaleFactor)
    
    if addOverFlow:
            print 'Adding overflow'
            uFlow = hTree.GetBinContent(0)+hTree.GetBinContent(1)
            oFlow = hTree.GetBinContent(hTree.GetNbinsX()+1)+hTree.GetBinContent(hTree.GetNbinsX())
            uFlowErr = ROOT.TMath.Sqrt(ROOT.TMath.Power(hTree.GetBinError(0),2)+ROOT.TMath.Power(hTree.GetBinError(1),2))
            oFlowErr = ROOT.TMath.Sqrt(ROOT.TMath.Power(hTree.GetBinError(hTree.GetNbinsX()),2)+ROOT.TMath.Power(hTree.GetBinError(hTree.GetNbinsX()+1),2))
            hTree.SetBinContent(1,uFlow)
            hTree.SetBinContent(hTree.GetNbinsX(),oFlow)
            hTree.SetBinError(1,uFlowErr)
            hTree.SetBinError(hTree.GetNbinsX(),oFlowErr)
              
            
    print '\t-->import %s\t Integral: %s'%(job.name,hTree.Integral())
            
    hTree.SetDirectory(0)
    input.Close()  
    
    return hTree, group
    

######################


def orderandadd(histos,typs,setup):
#ORDER AND ADD TOGETHER

    ordnung=[]
    ordnungtyp=[]
    num=[0]*len(setup)
    for i in range(0,len(setup)):
        for j in range(0,len(histos)):
            if typs[j] in setup[i]:
                num[i]+=1
                ordnung.append(histos[j])
                ordnungtyp.append(typs[j])

    del histos
    del typs

    histos=ordnung
    typs=ordnungtyp

    print typs

    for k in range(0,len(num)):
        for m in range(0,num[k]):
            if m > 0:
            
                #add
                histos[k].Add(histos[k+1],1)
                printc('magenta','','\t--> added %s to %s'%(typs[k],typs[k+1]))
                del histos[k+1]
                del typs[k+1]

    del histos[len(setup):]
    del typs[len(setup):]

    return histos, typs


Revision:	1.18
Committed:	Fri Oct 19 15:44:07 2012 UTC (12 years, 6 months ago) by peller
Content type:	text/x-python
Branch:	MAIN
CVS Tags:	hcpApproval, HCP_unblinding
Changes since 1.17:	+5 -3 lines
Log Message:	BDT plots from TH files
#	User	Rev	Content
1	peller	1.1	from samplesclass import sample
2			from printcolor import printc
3			import pickle
4			import ROOT
5			from ROOT import TFile, TTree
6			import ROOT
7			from array import array
8	nmohr	1.3	from BetterConfigParser import BetterConfigParser
9	peller	1.1	import sys
10
11
12	nmohr	1.10	def getScale(job,path,config,rescale,subsample=-1):
13			anaTag=config.get('Analysis','tag')
14	nmohr	1.9	input = TFile.Open(path+'/'+job.getpath())
15	peller	1.1	CountWithPU = input.Get("CountWithPU")
16			CountWithPU2011B = input.Get("CountWithPU2011B")
17			#print lumi*xsecs[i]/hist.GetBinContent(1)
18	peller	1.7
19			if subsample>-1:
20	peller	1.18	if type(job.xsec[subsample]) == str: xsec=float(eval(job.xsec[subsample]))
21			else: xsec=float(job.xsec[subsample])
22	peller	1.7	sf=float(job.sf[subsample])
23			else:
24	peller	1.18	if type(job.xsec) == str: xsec=float(eval(job.xsec))
25			else: xsec=float(job.xsec)
26	peller	1.7	sf=float(job.sf)
27
28
29	nmohr	1.6	theScale = 1.
30			if anaTag == '7TeV':
31	peller	1.7	theScale = float(job.lumi)xsecsf/(0.46502CountWithPU.GetBinContent(1)+0.53498CountWithPU2011B.GetBinContent(1))*rescale/float(job.split)
32	nmohr	1.6	elif anaTag == '8TeV':
33	peller	1.7	theScale = float(job.lumi)xsecsf/(CountWithPU.GetBinContent(1))*rescale/float(job.split)
34	nmohr	1.6	return theScale
35	peller	1.1
36	peller	1.11	def getHistoFromTree(job,path,config,options,rescale=1,subsample=-1,which_weightF='weightF'):
37	peller	1.8
38			#print job.getpath()
39			#print options
40	peller	1.1	treeVar=options[0]
41	peller	1.7	if subsample>-1:
42			name=job.subnames[subsample]
43			group=job.group[subsample]
44			else:
45			name=job.name
46			group=job.group
47
48	peller	1.1	#title=job.plotname()
49			nBins=int(options[3])
50			xMin=float(options[4])
51			xMax=float(options[5])
52	peller	1.16	#addOverFlow=eval(config.get('Plot_general','addOverFlow'))
53			addOverFlow = False
54	peller	1.1
55	peller	1.17	TrainFlag = eval(config.get('Analysis','TrainFlag'))
56			if TrainFlag: traincut = " & EventForTraining == 0"
57			if not TrainFlag: traincut=""
58
59	peller	1.1	if job.type != 'DATA':
60	peller	1.12
61			if type(options[7])==str:
62			cutcut=config.get('Cuts',options[7])
63			elif type(options[7])==list:
64			cutcut=config.get('Cuts',options[7][0])
65			cutcut=cutcut.replace(options[7][1],options[7][2])
66	peller	1.18	print cutcut
67	peller	1.8	if subsample>-1:
68	peller	1.17	treeCut='%s & %s%s'%(cutcut,job.subcuts[subsample],traincut)
69	peller	1.7	else:
70	peller	1.17	treeCut='%s%s'%(cutcut,traincut)
71	peller	1.1
72			elif job.type == 'DATA':
73	peller	1.5	cutcut=config.get('Cuts',options[8])
74	peller	1.8	treeCut='%s'%(cutcut)
75	peller	1.5
76	peller	1.1
77	nmohr	1.9	input = TFile.Open(path+'/'+job.getpath(),'read')
78	peller	1.1
79			Tree = input.Get(job.tree)
80			#Tree=tmpTree.CloneTree()
81			#Tree.SetDirectory(0)
82
83			#Tree=tmpTree.Clone()
84	peller	1.11	weightF=config.get('Weights',which_weightF)
85	peller	1.1	#hTree = ROOT.TH1F('%s'%name,'%s'%title,nBins,xMin,xMax)
86			#hTree.SetDirectory(0)
87			#hTree.Sumw2()
88			#print 'drawing...'
89			if job.type != 'DATA':
90			#print treeCut
91			#print job.name
92			if Tree.GetEntries():
93			Tree.Draw('%s>>%s(%s,%s,%s)' %(treeVar,name,nBins,xMin,xMax),'(%s)*(%s)' %(treeCut,weightF), "goff,e")
94			full=True
95			else:
96			full=False
97			elif job.type == 'DATA':
98
99			if len(options)>10:
100			if options[11] == 'blind':
101			treeCut = treeCut + '&'+treeVar+'<0'
102
103
104			Tree.Draw('%s>>%s(%s,%s,%s)' %(treeVar,name,nBins,xMin,xMax),treeCut, "goff,e")
105			full = True
106			if full:
107			hTree = ROOT.gDirectory.Get(name)
108			else:
109			hTree = ROOT.TH1F('%s'%name,'%s'%name,nBins,xMin,xMax)
110			hTree.Sumw2()
111			#print job.name + ' Sumw2', hTree.GetEntries()
112
113			if job.type != 'DATA':
114	nmohr	1.10	ScaleFactor = getScale(job,path,config,rescale,subsample)
115	peller	1.1	if ScaleFactor != 0:
116			hTree.Scale(ScaleFactor)
117	nmohr	1.15
118			if addOverFlow:
119			print 'Adding overflow'
120			uFlow = hTree.GetBinContent(0)+hTree.GetBinContent(1)
121			oFlow = hTree.GetBinContent(hTree.GetNbinsX()+1)+hTree.GetBinContent(hTree.GetNbinsX())
122			uFlowErr = ROOT.TMath.Sqrt(ROOT.TMath.Power(hTree.GetBinError(0),2)+ROOT.TMath.Power(hTree.GetBinError(1),2))
123			oFlowErr = ROOT.TMath.Sqrt(ROOT.TMath.Power(hTree.GetBinError(hTree.GetNbinsX()),2)+ROOT.TMath.Power(hTree.GetBinError(hTree.GetNbinsX()+1),2))
124			hTree.SetBinContent(1,uFlow)
125			hTree.SetBinContent(hTree.GetNbinsX(),oFlow)
126			hTree.SetBinError(1,uFlowErr)
127			hTree.SetBinError(hTree.GetNbinsX(),oFlowErr)
128
129	peller	1.1
130			print '\t-->import %s\t Integral: %s'%(job.name,hTree.Integral())
131
132			hTree.SetDirectory(0)
133	peller	1.7	input.Close()
134
135			return hTree, group
136	peller	1.1
137
138			######################
139
140
141
142			def orderandadd(histos,typs,setup):
143			#ORDER AND ADD TOGETHER
144
145			ordnung=[]
146			ordnungtyp=[]
147			num=[0]*len(setup)
148			for i in range(0,len(setup)):
149			for j in range(0,len(histos)):
150	peller	1.2	if typs[j] in setup[i]:
151	peller	1.1	num[i]+=1
152			ordnung.append(histos[j])
153			ordnungtyp.append(typs[j])
154
155			del histos
156			del typs
157
158			histos=ordnung
159			typs=ordnungtyp
160
161	peller	1.14	print typs
162
163	peller	1.1	for k in range(0,len(num)):
164			for m in range(0,num[k]):
165			if m > 0:
166
167			#add
168			histos[k].Add(histos[k+1],1)
169	peller	1.2	printc('magenta','','\t--> added %s to %s'%(typs[k],typs[k+1]))
170	peller	1.1	del histos[k+1]
171			del typs[k+1]
172
173			del histos[len(setup):]
174			del typs[len(setup):]
175
176			return histos, typs
177
178