CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
plotResidualsPerLayer.py
Go to the documentation of this file.
1 import ROOT
2 from fitResidual import fitResidual
3 from drawHistoAllChambers import drawHisto
4 
5 layerCorrectionFactors = {'SL1':(1.17,1.16,1.15,1.14),
6  'SL2':(1.83,1.20,1.20,1.83),
7  'SL3':(1.14,1.15,1.16,1.17)}
8 
9 def plotResLayer(fileName,sl,layer,
10  dir='DQMData/Run 1/DT/Run summary/DTCalibValidation',
11  option="HISTOPE1",draw=True):
12 
13  mean_ymin = -0.02
14  mean_ymax = 0.02
15  sig_ymin = 0.
16  sig_ymax = 0.1
17 
18  slType = sl
19  slStr = "SL%d" % slType
20  layerType = layer
21  layerStr = "Layer%d" % layerType
22  verbose = False
23  nSigmas = 2
24 
25  ROOT.TH1.AddDirectory(False)
26 
27  file = ROOT.TFile(fileName,'read')
28 
29  wheels = (-2,-1,0,1,2)
30  stations = (1,2,3,4)
31 
32  # (Wh-2 MB1 Sec1 ... Wh-2 MB1 Sec12 ... Wh-1 MB1 Sec1 ... Wh-1 MB1 Sec12 ...)
33  # (Wh-2 MB2 Sec1 ... Wh-2 MB2 Sec12 ... Wh-1 MB2 Sec1 ... Wh-1 MB1 Sec12 ...) ...
34  nBins = 250
35  if slType == 2: nBins = 180
36  histoMean = ROOT.TH1F("h_ResMeanAll_%s_%s" % (slStr,layerStr),"Mean of residuals",nBins,0,nBins)
37  histoSigma = ROOT.TH1F("h_ResSigmaAll_%s_%s" % (slStr,layerStr),"Sigma of residuals",nBins,0,nBins)
38  for st in stations:
39  nSectors = 12
40  if st == 4: nSectors = 14
41  if st == 4 and slType == 2: continue
42  if verbose: print "Station",st
43  for wh in wheels:
44  if verbose: print "Wheel",wh
45  for sec in range(1,nSectors+1):
46  if verbose: print "Sector",sec
47  # Get histogram
48  histoName = "%s/Wheel%d/Station%d/Sector%d/%s/hResDist_STEP3_W%d_St%d_Sec%d_%s_%s" % (dir,wh,st,sec,slStr,wh,st,sec,slStr,layerStr)
49  print "Accessing",histoName
50  histo = file.Get(histoName)
51  (histo,fitFunc) = fitResidual(histo,nSigmas,verbose)
52  fitMean = fitFunc.GetParameter(1)
53  fitMeanErr = fitFunc.GetParError(1)
54  fitSigma = fitFunc.GetParameter(2)
55  fitSigmaErr = fitFunc.GetParError(2)
56 
57  layerIdx = (layer - 1)
58  corrFactor = layerCorrectionFactors[slStr][layerIdx]
59 
60  binHistoNew = (st - 1)*60 + (wh + 2)*nSectors + sec
61  if verbose: print "Bin in summary histo",binHistoNew
62  histoMean.SetBinContent(binHistoNew,fitMean)
63  histoMean.SetBinError(binHistoNew,fitMeanErr)
64  histoSigma.SetBinContent(binHistoNew,fitSigma*corrFactor)
65  histoSigma.SetBinError(binHistoNew,fitSigmaErr*corrFactor)
66 
67  if sec == 1:
68  label = "Wheel %d" % wh
69  if wh == -2: label += " MB%d" % st
70  histoMean.GetXaxis().SetBinLabel(binHistoNew,label)
71  histoSigma.GetXaxis().SetBinLabel(binHistoNew,label)
72 
73  objectsMean = drawHisto(histoMean,title="Mean of residuals (cm)",
74  ymin=mean_ymin,ymax=mean_ymax,option=option,draw=draw)
75  objectsSigma = drawHisto(histoSigma,title="Resolution (cm)",
76  ymin=sig_ymin,ymax=sig_ymax,option=option,draw=draw)
77 
78  return (objectsMean,objectsSigma)
79 
80 def plot(fileName,sl,
81  dir='DQMData/Run 1/DT/Run summary/DTCalibValidation',type='mean',option='HISTOPE1'):
82  colors = (2,4,12,44,55,38,27,46)
83  markers = (24,25,26,27,28,30,32,5)
84  labels=['Layer 1','Layer 2','Layer 3','Layer 4']
85  idx_type = None
86  if type == 'mean': idx_type = 0
87  elif type == 'sigma': idx_type = 1
88  else: raise RuntimeError("Wrong option: %s" % type)
89 
90  idx = 0
91  canvas = None
92  objects = None
93  histos = []
94  for layer in range(1,5):
95  draw = False
96  if not idx: draw = True
97 
98  objs = plotResLayer(fileName,sl,layer,dir,option,draw)
99  histos.append(objs[idx_type][1])
100  histos[-1].SetName( "%s_%d" % (histos[-1].GetName(),idx) )
101  if not idx:
102  canvas = objs[idx_type][0]
103  objects = objs[idx_type][2]
104 
105  canvas.cd()
106  if idx:
107  histos[-1].SetLineColor(colors[ (idx - 1) % len(colors) ])
108  histos[-1].SetMarkerColor(colors[ (idx - 1) % len(colors) ])
109  histos[-1].SetMarkerStyle(markers[ (idx - 1) % len(markers) ])
110 
111  histos[-1].Draw(option + "SAME")
112 
113  idx += 1
114 
115  legend = ROOT.TLegend(0.4,0.7,0.95,0.8)
116  for idx in range( len(histos) ):
117  histo = histos[idx]
118  label = histo.GetName()
119  if len(labels): label = labels[idx]
120  legend.AddEntry(histo,label,"LP")
121 
122  idx += 1
123 
124  canvas.cd()
125  legend.SetFillColor( canvas.GetFillColor() )
126  legend.Draw("SAME")
127  objects.append(legend)
128 
129  # Compute averages
130  # (Wh-2 MB1 Sec1 ... Wh-2 MB1 Sec12 ... Wh-1 MB1 Sec1 ... Wh-1 MB1 Sec12 ...)
131  # (Wh-2 MB2 Sec1 ... Wh-2 MB2 Sec12 ... Wh-1 MB2 Sec1 ... Wh-1 MB1 Sec12 ...) ...
132  import math
133  wheels = (-2,-1,0,1,2)
134  stations = (1,2,3,4)
135  slType = sl
136  slStr = "SL%d" % slType
137 
138  nBinsAve = len(stations)*len(wheels)
139  histoAverage = ROOT.TH1F("h_AverageAll_" + slStr,"",nBinsAve,0,nBinsAve)
140  averages = {}
141  averagesErr = {}
142  averagesSumw = {}
143  print "Averages:"
144  for st in stations:
145  nSectors = 12
146  if st == 4: nSectors = 14
147  if st == 4 and slType == 2: continue
148  for wh in wheels:
149  binHistoAve = (st - 1)*5 + (wh + 2) + 1
150  label = "Wheel %d" % wh
151  if wh == -2: label += " MB%d" % st
152  histoAverage.GetXaxis().SetBinLabel(binHistoAve,label)
153 
154  averages[(st,wh)] = 0.
155  averagesSumw[(st,wh)] = 0.
156  for sec in range(1,nSectors+1):
157  binHisto = (st - 1)*60 + (wh + 2)*nSectors + sec
158  for idx in range( len(histos) ):
159  histo = histos[idx]
160  value = histo.GetBinContent( binHisto )
161  error = histo.GetBinError( binHisto )
162  averages[(st,wh)] += value/( error*error )
163  averagesSumw[(st,wh)] += 1./( error*error )
164  # Average per (st,wh)
165  averages[(st,wh)] = averages[(st,wh)]/averagesSumw[(st,wh)]
166  averagesErr[(st,wh)] = math.sqrt( 1./averagesSumw[(st,wh)] )
167  histoAverage.SetBinContent(binHistoAve,averages[(st,wh)])
168  histoAverage.SetBinError(binHistoAve,averagesErr[(st,wh)])
169  print "Station %d, Wheel %d: %.4f +/- %.6f" % (st,wh,averages[(st,wh)],averagesErr[(st,wh)])
170 
171  canvasAverage = ROOT.TCanvas("c_" + histoAverage.GetName())
172  canvasAverage.SetGridx()
173  canvasAverage.SetGridy()
174  canvasAverage.SetFillColor( 0 )
175  canvasAverage.cd()
176  mean_ymin = -0.02
177  mean_ymax = 0.02
178  sig_ymin = 0.
179  sig_ymax = 0.1
180  if type == 'mean':
181  histoAverage.GetYaxis().SetTitle("Mean of residuals (cm)")
182  histoAverage.GetYaxis().SetRangeUser(mean_ymin,mean_ymax)
183  elif type == 'sigma':
184  histoAverage.GetYaxis().SetTitle("Resolution (cm)")
185  histoAverage.GetYaxis().SetRangeUser(sig_ymin,sig_ymax)
186 
187  histoAverage.SetStats(0)
188  histoAverage.SetLineWidth(2)
189  histoAverage.SetMarkerStyle( 27 )
190  histoAverage.SetMarkerSize( 1.5 )
191  histoAverage.LabelsOption("d","X")
192  histoAverage.Draw("E2")
193 
194  return ( (canvas,canvasAverage),(histos,histoAverage),objects )
195 
196 def plotMean(fileName,sl,dir='DQMData/Run 1/DT/Run summary/DTCalibValidation',option='HISTOPE1'):
197  type = 'mean'
198  objs = plot(fileName,sl,dir,type,option)
199  return objs
200 
201 def plotSigma(fileName,sl,dir='DQMData/Run 1/DT/Run summary/DTCalibValidation',option='HISTOPE1'):
202  type = 'sigma'
203  objs = plot(fileName,sl,dir,type,option)
204  return objs
205 
206 def plotSigmaAll(fileName,dir='DQMData/Run 1/DT/Run summary/DTCalibValidation',option='HISTOPE1',outputFileName=''):
207  colors = (2,4,12,44,55,38,27,46)
208  markers = (24,25,26,27,28,30,32,5)
209 
210  slList = (1,2,3)
211  labels = ('R-#phi SL1','R-z SL2','R-#phi SL3')
212  canvas = None
213  objects = None
214  histos = []
215  idx = 0
216  for sl in slList:
217  draw = False
218  if not idx: draw = True
219 
220  objs = plotSigma(fileName,sl,dir,option)
221  histos.append(objs[1][1])
222  histos[-1].SetName( "%s_%d" % (histos[-1].GetName(),idx) )
223  if not idx:
224  canvas = objs[0][1]
225  #objects = objs[2][1]
226 
227  canvas.cd()
228  if idx:
229  histos[-1].SetLineColor(colors[ (idx - 1) % len(colors) ])
230  histos[-1].SetMarkerColor(colors[ (idx - 1) % len(colors) ])
231  histos[-1].SetMarkerStyle(markers[ (idx - 1) % len(markers) ])
232 
233  histos[-1].Draw(option + "SAME")
234 
235  idx += 1
236 
237  legend = ROOT.TLegend(0.4,0.7,0.95,0.8)
238  for idx in range( len(histos) ):
239  histo = histos[idx]
240  label = histo.GetName()
241  if len(labels): label = labels[idx]
242  legend.AddEntry(histo,label,"LP")
243 
244  idx += 1
245 
246  canvas.cd()
247  legend.SetFillColor( canvas.GetFillColor() )
248  legend.Draw("SAME")
249  if not objects: objects = [legend]
250  else: objects.append(legend)
251 
252  if outputFileName:
253  outputFile = ROOT.TFile(outputFileName,'recreate')
254  outputFile.cd()
255  for histo in histos: histo.Write()
256  outputFile.Close()
257  return 0
258  else:
259  return (canvas,histos,objects)
260 
261 #def plotDataVsMCFromFile(fileNameData,fileNameMC,labels=[]):
262 def plotFromFile(fileNames,labels=[]):
263 
264  AddDirectoryStatus_ = ROOT.TH1.AddDirectoryStatus()
265  ROOT.TH1.AddDirectory(False)
266 
267  #fileData = ROOT.TFile(fileNameData,'read')
268  #fileMC = ROOT.TFile(fileNameMC,'read')
269  rootFiles = []
270  for file in fileNames: rootFiles.append( ROOT.TFile(file,'read') )
271 
272  variables = ['h_AverageAll_SL1_0',
273  'h_AverageAll_SL2_1',
274  'h_AverageAll_SL3_2']
275 
276  colors = (1,2,4,12,44,55,38,27,46)
277  markers = (20,24,25,26,27,28,30,32,5)
278  objects = None
279  canvases = []
280  legends = []
281  histos = []
282  idx_var = 0
283  for var in variables:
284  print "Accessing",var
285  #histoData = fileData.Get(var)
286  #histoData.SetName(histoData.GetName() + "_data")
287  #histoMC = fileMC.Get(var)
288  #histoMC.SetName(histoMC.GetName() + "_mc")
289  #histoData.SetLineColor(1)
290  #histoData.SetMarkerStyle(20)
291  #histoData.SetMarkerSize(1.4)
292  #histoData.SetMarkerColor(1)
293 
294  #histoMC.SetLineColor(2)
295  #histoMC.SetMarkerStyle(24)
296  #histoMC.SetMarkerSize(1.4)
297  #histoMC.SetMarkerColor(2)
298 
299  histos_tmp = []
300  idx = 0
301  for file in rootFiles:
302  histos_tmp.append( file.Get(var) )
303  histos_tmp[-1].SetName( "%s_%d" % (histos_tmp[-1].GetName(),idx) )
304  print "Created",histos_tmp[-1].GetName()
305  histos_tmp[-1].SetLineColor(colors[ idx % len(colors) ])
306  histos_tmp[-1].SetMarkerColor(colors[ idx % len(colors) ])
307  histos_tmp[-1].SetMarkerStyle(markers[ idx % len(markers) ])
308  histos_tmp[-1].SetMarkerSize(1.4)
309  idx += 1
310  histos.append( histos_tmp )
311 
312  canvases.append( ROOT.TCanvas("c_" + var,var) )
313  canvases[-1].SetGridx()
314  canvases[-1].SetGridy()
315  canvases[-1].SetFillColor(0)
316  canvases[-1].cd()
317  #histoData.Draw()
318  #histoMC.Draw("SAME")
319  #histos.append( (histoData,histoMC) )
320  histos[-1][0].Draw()
321  for histo in histos[-1][1:]: histo.Draw("SAME")
322 
323  if len(labels):
324  #labelData = labels[0]
325  #labelMC = labels[1]
326  legends.append( ROOT.TLegend(0.4,0.7,0.95,0.8) )
327  idx = 0
328  for histo in histos[-1]:
329  legends[-1].AddEntry(histo,labels[idx],"LP")
330  idx += 1
331 
332  legends[-1].SetFillColor( canvases[-1].GetFillColor() )
333  legends[-1].Draw("SAME")
334 
335  idx_var += 1
336 
337  if not objects: objects = [legends]
338  else: objects.append(legends)
339 
340  ROOT.TH1.AddDirectory(AddDirectoryStatus_)
341 
342  return (canvases,histos,objects)
Generated for CMSSW Reference Manual by doxygen 1.8.5