ATLAS Offline Software
beamspotCutsPlots.py
Go to the documentation of this file.
1 #!/usr/bin/env python
2 
3 # Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
4 
5 from optparse import OptionParser
6 import os, glob
7 import string
8 
9 __author__ = 'James Walder'
10 __version__ = '$Id$'
11 __usage__ = '%prog [options] '
12 
13 def makeListFromString(s,token=','):
14  names = s.split(token)
15  return names
16 
17 parser = OptionParser(usage=__usage__, version=__version__)
18 parser.add_option('-f','--files',dest='files', default='bs.root', help='Root files for data plots')
19 parser.add_option('-d','--dir',dest='dir', default='Beamspot', help='Root directory for beamspots')
20 parser.add_option('-b','--batch', dest='batch', action='store_true', default=False, help='Run in batch mode')
21 
22 (options,args) = parser.parse_args()
23 
24 if options.batch:
25  os.unsetenv('DISPLAY')
26 
27 
28 from ROOT import TH1D, TCanvas
29 from ROOT import gSystem, gROOT, TFile, TH1D, TH2D, TCanvas, TTree, TChain
30 from ROOT import TGraphErrors
31 from ROOT import TMultiGraph
32 from ROOT import EColor
33 from ROOT import TLegend,TLine
34 from ROOT import TLatex,TGaxis
35 
36 legends = [] # holder
37 pdglines= [] # holder
38 labels = []
39 ghists = [] # additional hists
40 failed = [] # list of failed plots
41 graphs = []
42 canvases = []
43 
44 PointSize = 1.6
45 MaxHeightFactor = 1.5
46 
47 plotTextSize = 0.04
48 gROOT.ProcessLine(".L AtlasStyle.C")
49 gROOT.ProcessLine(".L AtlasUtils.C")
50 gROOT.ProcessLine("SetAtlasStyle();")
51 
52 def findBeamspots(file,dir, pattern):
53  import re
54  from ROOT import gDirectory
55  file.cd(dir)
56  p = re.compile(pattern)
57  print p.pattern
58  n = gDirectory.GetListOfKeys().GetSize()
59  print "N: ", n
60  matched = []
61  for i in range(n):
62  name = gDirectory.GetListOfKeys().At(i).GetName()
63  if not p.match(name):
64  continue
65  print name
66  #print p.match(name)
67  h = gDirectory.Get(gDirectory.GetListOfKeys().At(i).GetName())
68  #print h, h.GetName()
69  matched.append(h)
70 
71  file.cd()
72  return matched
73 
74 
75 def makeFloatFromString(s):
76  subs = { 'p':'.', 'm':'-' }
77  n = ""
78  for i in s:
79  c = i
80  if i in subs.keys():
81  c = subs[i]
82  n = n+c
83  return float(n)
84 
85 
86 def splitBeamspotTreeTitle(title, cutToken="CUT"):
87  t = title.split(cutToken)
88  namestring = t[0].strip('_')
89  cutstring = t[-1].strip('_')
90  cuts = cutstring.split('_')
91  output = {}
92  for i in range(0,len(cuts),2):
93  output[cuts[i]] = makeFloatFromString(cuts[i+1])
94 
95  return output,namestring
96 
97 def getTree(name, file, dir):
98  tr = file.Get( dir +"/" + name)
99  print type(tr)
100  return tr
101 
102 # beamspot tree-type stuff
103 def makePull( tree, xval='x0',xcov='x0x0', xTrue=0., isCovErr=True, xlow=-8, xhigh=8, xbins=100, doGausFit=False, cuts=""):
104  h = TH1D("hPull_"+tree.GetName()+"_"+xval, "hPull_"+tree.GetName()+"_"+xval, xbins, xlow,xhigh)
105  if isCovErr:
106  p = '(%(x)s - %(t)s)/sqrt(%(e)s) >> %(h)s' % { 'x':xval, 't':xTrue, 'e':xcov, 'h':h.GetName()}
107  else:
108  p = '(%(x)s - %(t)s)/(%(e)s) >> %(h)s' % { 'x':xval, 't':xTrue, 'e':xcov, 'h':h.GetName()}
109  print p
110  tree.Draw( p , cuts);
111  h.Print()
112 
113  mean = h.GetMean()
114  meanErr = h.GetMeanError()
115  rms = h.GetRMS()
116  rmsErr = h.GetRMSError()
117  if doGausFit:
118  h.Fit("gaus","LL")
119  f = h.GetFunction("gaus")
120  mean = f.GetParameter(1)
121  meanErr = f.GetParError(1)
122  rms = f.GetParameter(2)
123  rmsErr = f.GetParError(2)
124 
125  return (mean, meanErr, rms, rmsErr)
126 
127 
128 # Start of main programs
129 
130 #Open the files
131 
132 files = makeListFromString(options.files)
133 if len(files) == 1:
134  fii = TFile(files[0])
135 else:
136  print "Multiple files not yet supported"
137  sys.exit()
138 
139 
140 # matched = findBeamspots(fii, options.dir, "Beamspots_.*MinVt.*")
141 
142 def makeBeamspotCutPlots(name, points):
143  nPoints = len(points)
144 
145  x = TGraphErrors(nPoints)
146  x.SetName(name)
147  x.SetTitle(name + " x")
148  y = TGraphErrors(nPoints)
149  y.SetName(name)
150  y.SetTitle(name + " y")
151  z = TGraphErrors(nPoints)
152  z.SetName(name)
153  z.SetTitle(name + " z")
154  sx = TGraphErrors(nPoints)
155  sx.SetName(name)
156  sx.SetTitle(name + " sx")
157  sy = TGraphErrors(nPoints)
158  sy.SetName(name)
159  sy.SetTitle(name + " sy")
160  sz = TGraphErrors(nPoints)
161  sz.SetName(name)
162  sz.SetTitle(name + " sz")
163 
164  xrms = TGraphErrors(nPoints)
165  xrms.SetName(name+"_rms")
166  xrms.SetTitle(name+"_rms" + " x")
167  yrms = TGraphErrors(nPoints)
168  yrms.SetName(name+"_rms")
169  yrms.SetTitle(name+"_rms" + " y")
170  zrms = TGraphErrors(nPoints)
171  zrms.SetName(name+"_rms")
172  zrms.SetTitle(name+"_rms" + " z")
173  sxrms = TGraphErrors(nPoints)
174  sxrms.SetName(name+"_rms")
175  sxrms.SetTitle(name+"_rms" + " sx")
176  syrms = TGraphErrors(nPoints)
177  syrms.SetName(name+"_rms")
178  syrms.SetTitle(name+"_rms" + " sy")
179  szrms = TGraphErrors(nPoints)
180  szrms.SetName(name+"_rms")
181  szrms.SetTitle(name+"_rms" + " sz")
182 
183  graphs.extend([x,y,z,sx,sy,sz,xrms,yrms,zrms,sxrms,syrms,szrms])
184  print len(points)
185  for i in range(len(points)):
186  xval = points[i][0]
187  yvals = points[i][1]
188 
189  xerr =0;
190 
191  x.SetPoint(i, xval, yvals[0][0])
192  y.SetPoint(i, xval, yvals[1][0])
193  z.SetPoint(i, xval, yvals[2][0])
194  sx.SetPoint(i, xval, yvals[3][0])
195  sy.SetPoint(i, xval, yvals[4][0])
196  sz.SetPoint(i, xval, yvals[5][0])
197 
198  x.SetPointError(i, xerr, yvals[0][1])
199  y.SetPointError(i, xerr, yvals[1][1])
200  z.SetPointError(i, xerr, yvals[2][1])
201  sx.SetPointError(i, xerr, yvals[3][1])
202  sy.SetPointError(i, xerr, yvals[4][1])
203  sz.SetPointError(i, xerr, yvals[5][1])
204 
205 
206  xrms.SetPoint(i, xval, yvals[0][2])
207  yrms.SetPoint(i, xval, yvals[1][2])
208  zrms.SetPoint(i, xval, yvals[2][2])
209  sxrms.SetPoint(i, xval, yvals[3][2])
210  syrms.SetPoint(i, xval, yvals[4][2])
211  szrms.SetPoint(i, xval, yvals[5][2])
212 
213  xrms.SetPointError(i, xerr, yvals[0][3])
214  yrms.SetPointError(i, xerr, yvals[1][3])
215  zrms.SetPointError(i, xerr, yvals[2][3])
216  sxrms.SetPointError(i, xerr, yvals[3][3])
217  syrms.SetPointError(i, xerr, yvals[4][3])
218  szrms.SetPointError(i, xerr, yvals[5][3])
219 
220  c = TCanvas("cx")
221  c.Divide(3,1)
222  c.cd(1)
223  x.Draw("ap")
224  c.cd(2)
225  y.Draw("ap+Y+")
226  c.cd(3)
227  yaxis_xrms = TGaxis( -1, 0.2, 1 ,0.2, -1, 2,510,"+R")
228  yaxis_xrms.ImportAxisAttributes( y.GetHistogram().GetYaxis() )
229  #g = TMultiGraph()
230  #g.Add(x,"aXp")
231  #g.Add(y,"aY*")
232  #g.Draw("a")
233  # graphs.append(g)
234  # aa = y.GetHistogram()
235 
236  #x.Draw("ap")
237  yaxis_xrms.Draw()
238  #axis.PaintAxis(0,0.5,0.1,0.6,0.4,1.4,510,"+R")
239  #y.Draw("pY+sames")
240  #b.Draw("sames")
241  c.Modified()
242  c.Update()
243  #c.Print()
244  canvases.append(c)
245 
246 def makeCutPlotPulls( file, dir, pattern, name):
247  trueVals = { 'x0':-0.15, 'y0':1., 'z':-9., 'sx':'0.72', 'sy':'0.42', 'sz':44, 'ax':'0.', 'ay':0., 'k':1, 'rhoxy':0.0}
248 
249  # list of trees
250  matched = findBeamspots(file, dir, pattern)
251  points= []
252 
253  for i in matched:
254  # loop over trees and get cut values
255  o, n = splitBeamspotTreeTitle( i.GetName(), "CUT" )
256  # make assumption of only one cut value
257  xval = o[o.keys()[0]]
258  print "JW: ", xval
259  xpos = makePull( getTree( i.GetName(), file, dir), xval='x0', xcov='x0x0', xTrue=trueVals['x0'])
260  ypos = makePull( getTree( i.GetName(), file, dir), xval='y0', xcov='y0y0', xTrue=trueVals['y0'])
261  zpos = makePull( getTree( i.GetName(), file, dir), xval='z', xcov='zz', xTrue=trueVals['z'])
262 
263  sx = makePull( getTree( i.GetName(), file, dir), xval='sx', xcov='sxsx', xTrue=trueVals['sx'])
264  sy = makePull( getTree( i.GetName(), file, dir), xval='sy', xcov='sysy', xTrue=trueVals['sy'])
265  sz = makePull( getTree( i.GetName(), file, dir), xval='sz', xcov='szsz', xTrue=trueVals['sz'])
266  ax = makePull( getTree( i.GetName(), file, dir), xval='ax', xcov='axax', xTrue=trueVals['ax'])
267  ay = makePull( getTree( i.GetName(), file, dir), xval='ay', xcov='ayay', xTrue=trueVals['ay'])
268  k = makePull( getTree( i.GetName(), file, dir), xval='k', xcov='kk', xTrue=trueVals['k'])
269  rhoxy = makePull( getTree( i.GetName(), file, dir), xval='rhoxy', xcov='rhoxyrhoxy', xTrue=trueVals['rhoxy'])
270  points.append( (xval, [xpos, ypos, zpos, sx, sy, sz, ax, ay, k, rhoxy ]) )
271 
272 
273  plots = makeBeamspotCutPlots(name, points)
274 
275 
276 makeCutPlotPulls(fii, options.dir, "Beamspots_.*Min.*","MinProb")
277 
278 
279 #vals = []
280 #for i in matched:
281 # vals.append( [i.GetName(), splitBeamspotTreeTitle(i.GetName())] )##
282 #
283 #for i in vals:
284 # i.append( makePull(getTree( i[0] ,fii,options.dir),doGausFit=True) )#
285 
289 
290 
291 s = raw_input('--> ')
292 exit()
293 
beamspotCutsPlots.makeCutPlotPulls
def makeCutPlotPulls(file, dir, pattern, name)
Definition: beamspotCutsPlots.py:246
beamspotCutsPlots.makeBeamspotCutPlots
def makeBeamspotCutPlots(name, points)
Definition: beamspotCutsPlots.py:142
beamspotCutsPlots.makeFloatFromString
def makeFloatFromString(s)
Definition: beamspotCutsPlots.py:75
beamspotCutsPlots.findBeamspots
def findBeamspots(file, dir, pattern)
Definition: beamspotCutsPlots.py:52
beamspotCutsPlots.getTree
def getTree(name, file, dir)
Definition: beamspotCutsPlots.py:97
beamspotCutsPlots.splitBeamspotTreeTitle
def splitBeamspotTreeTitle(title, cutToken="CUT")
Definition: beamspotCutsPlots.py:86
beamspotCutsPlots.makePull
def makePull(tree, xval='x0', xcov='x0x0', xTrue=0., isCovErr=True, xlow=-8, xhigh=8, xbins=100, doGausFit=False, cuts="")
Definition: beamspotCutsPlots.py:103
beamspotCutsPlots.makeListFromString
def makeListFromString(s, token=',')
Definition: beamspotCutsPlots.py:13
plotBeamSpotVxVal.range
range
Definition: plotBeamSpotVxVal.py:195
calibdata.exit
exit
Definition: calibdata.py:236
python.CaloScaleNoiseConfig.type
type
Definition: CaloScaleNoiseConfig.py:78
readCCLHist.float
float
Definition: readCCLHist.py:83
Generated on Fri May 17 2024 21:10:54 for ATLAS Offline Software by   doxygen 1.8.18