VHbb/python/gethistofromtree.py

import pickle
import ROOT 
from ROOT import TFile, TTree
import ROOT
from array import array
import sys
from myutils import sample, printc

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