00001 import ROOT, array, os, re, random
00002 from math import *
00003 import time
00004 import pickle
00005
00006
00007 try:
00008 import json
00009 except ImportError:
00010 import simplejson as json
00011
00012
00013 from signConventions import *
00014
00015
00016 from mutypes import *
00017
00018 CPP_LOADED = False
00019
00020
00021 MAP_RESULTS_SAWTOOTH = {}
00022 MAP_RESULTS_FITSIN = {}
00023 MAP_RESULTS_BINS = {}
00024
00025
00026 TEST_RESULTS = {}
00027
00028
00029
00030
00031 def wheelm2only(dt, wheel, station, sector): return dt == "DT" and wheel == -2
00032 def wheelm1only(dt, wheel, station, sector): return dt == "DT" and wheel == -1
00033 def wheel0only(dt, wheel, station, sector): return dt == "DT" and wheel == 0
00034 def wheelp1only(dt, wheel, station, sector): return dt == "DT" and wheel == 1
00035 def wheelp2only(dt, wheel, station, sector): return dt == "DT" and wheel == 2
00036
00037 def wheelLetter(wheel):
00038 if wheel == -2: return "A"
00039 elif wheel == -1: return "B"
00040 elif wheel == 0: return "C"
00041 elif wheel == +1: return "D"
00042 elif wheel == +2: return "E"
00043 else: raise Exception
00044
00045 def wheelNumber(wheell):
00046 if wheell == "A": return -2
00047 elif wheell == "B": return -1
00048 elif wheell == "C": return 0
00049 elif wheell == "D": return 1
00050 elif wheell == "E": return 2
00051 else: raise Exception
00052
00053 def mean(xlist):
00054 s, n = 0., 0.
00055 for x in xlist:
00056 s += x
00057 n += 1.
00058 return s/n
00059
00060 def rms(xlist):
00061 s2, n = 0., 0.
00062 for x in xlist:
00063 s2 += x**2
00064 n += 1.
00065 return sqrt(s2/n)
00066
00067 def stdev(xlist):
00068 s, s2, n = 0., 0., 0.
00069 for x in xlist:
00070 s += x
00071 s2 += x**2
00072 n += 1.
00073 return sqrt(s2/n - (s/n)**2)
00074
00075 def wmean(xlist):
00076 s, w = 0., 0.
00077 for x, e in xlist:
00078 if e > 0.:
00079 wi = 1./e**2
00080 s += x*wi
00081 w += wi
00082 return s/w, sqrt(1./w)
00083
00084
00085
00086 tdrStyle = None
00087 def setTDRStyle():
00088 global tdrStyle
00089 tdrStyle = ROOT.TStyle("tdrStyle","Style for P-TDR")
00090
00091 tdrStyle.SetCanvasBorderMode(0)
00092 tdrStyle.SetCanvasColor(ROOT.kWhite)
00093 tdrStyle.SetCanvasDefH(600)
00094 tdrStyle.SetCanvasDefW(600)
00095 tdrStyle.SetCanvasDefX(0)
00096 tdrStyle.SetCanvasDefY(0)
00097
00098
00099 tdrStyle.SetPadBorderMode(0)
00100
00101 tdrStyle.SetPadColor(ROOT.kWhite)
00102 tdrStyle.SetPadGridX(False)
00103 tdrStyle.SetPadGridY(False)
00104 tdrStyle.SetGridColor(0)
00105 tdrStyle.SetGridStyle(3)
00106 tdrStyle.SetGridWidth(1)
00107
00108
00109 tdrStyle.SetFrameBorderMode(0)
00110 tdrStyle.SetFrameBorderSize(1)
00111 tdrStyle.SetFrameFillColor(0)
00112 tdrStyle.SetFrameFillStyle(0)
00113 tdrStyle.SetFrameLineColor(1)
00114 tdrStyle.SetFrameLineStyle(1)
00115 tdrStyle.SetFrameLineWidth(1)
00116
00117
00118
00119
00120 tdrStyle.SetHistLineColor(1)
00121 tdrStyle.SetHistLineStyle(0)
00122 tdrStyle.SetHistLineWidth(1)
00123
00124
00125
00126 tdrStyle.SetEndErrorSize(2)
00127
00128 tdrStyle.SetErrorX(0.)
00129
00130 tdrStyle.SetMarkerStyle(20)
00131
00132
00133 tdrStyle.SetOptFit(1)
00134 tdrStyle.SetFitFormat("5.4g")
00135 tdrStyle.SetFuncColor(2)
00136 tdrStyle.SetFuncStyle(1)
00137 tdrStyle.SetFuncWidth(1)
00138
00139
00140 tdrStyle.SetOptDate(0)
00141
00142
00143
00144
00145 tdrStyle.SetOptFile(0)
00146 tdrStyle.SetOptStat(0)
00147 tdrStyle.SetStatColor(ROOT.kWhite)
00148 tdrStyle.SetStatFont(42)
00149 tdrStyle.SetStatFontSize(0.025)
00150 tdrStyle.SetStatTextColor(1)
00151 tdrStyle.SetStatFormat("6.4g")
00152 tdrStyle.SetStatBorderSize(1)
00153 tdrStyle.SetStatH(0.1)
00154 tdrStyle.SetStatW(0.15)
00155
00156
00157
00158
00159
00160 tdrStyle.SetPadTopMargin(0.05)
00161 tdrStyle.SetPadBottomMargin(0.13)
00162 tdrStyle.SetPadLeftMargin(0.13)
00163 tdrStyle.SetPadRightMargin(0.05)
00164
00165
00166 tdrStyle.SetOptTitle(0)
00167 tdrStyle.SetTitleFont(42)
00168 tdrStyle.SetTitleColor(1)
00169 tdrStyle.SetTitleTextColor(1)
00170 tdrStyle.SetTitleFillColor(10)
00171 tdrStyle.SetTitleFontSize(0.05)
00172
00173
00174
00175
00176
00177
00178
00179
00180 tdrStyle.SetTitleColor(1, "XYZ")
00181 tdrStyle.SetTitleFont(42, "XYZ")
00182 tdrStyle.SetTitleSize(0.06, "XYZ")
00183
00184
00185 tdrStyle.SetTitleXOffset(0.9)
00186 tdrStyle.SetTitleYOffset(1.05)
00187
00188
00189
00190 tdrStyle.SetLabelColor(1, "XYZ")
00191 tdrStyle.SetLabelFont(42, "XYZ")
00192 tdrStyle.SetLabelOffset(0.007, "XYZ")
00193 tdrStyle.SetLabelSize(0.05, "XYZ")
00194
00195
00196 tdrStyle.SetAxisColor(1, "XYZ")
00197 tdrStyle.SetStripDecimals(True)
00198 tdrStyle.SetTickLength(0.03, "XYZ")
00199 tdrStyle.SetNdivisions(510, "XYZ")
00200 tdrStyle.SetPadTickX(1)
00201 tdrStyle.SetPadTickY(1)
00202
00203
00204 tdrStyle.SetOptLogx(0)
00205 tdrStyle.SetOptLogy(0)
00206 tdrStyle.SetOptLogz(0)
00207
00208
00209 tdrStyle.SetPaperSize(20.,20.)
00210
00211
00212
00213
00214
00215
00216
00217
00218
00219
00220
00221
00222 tdrStyle.cd()
00223
00224 setTDRStyle()
00225
00226 def set_palette(name=None, ncontours=999):
00227 """Set a color palette from a given RGB list
00228 stops, red, green and blue should all be lists of the same length
00229 see set_decent_colors for an example"""
00230
00231 if name == "halfgray":
00232 stops = [0.00, 0.34, 0.61, 0.84, 1.00]
00233 red = map(lambda x: 1. - (1.-x)/2., [1.00, 0.84, 0.61, 0.34, 0.00])
00234 green = map(lambda x: 1. - (1.-x)/2., [1.00, 0.84, 0.61, 0.34, 0.00])
00235 blue = map(lambda x: 1. - (1.-x)/2., [1.00, 0.84, 0.61, 0.34, 0.00])
00236 elif name == "gray":
00237 stops = [0.00, 0.34, 0.61, 0.84, 1.00]
00238 red = [1.00, 0.84, 0.61, 0.34, 0.00]
00239 green = [1.00, 0.84, 0.61, 0.34, 0.00]
00240 blue = [1.00, 0.84, 0.61, 0.34, 0.00]
00241 elif name == "blues":
00242 stops = [0.00, 0.34, 0.61, 0.84, 1.00]
00243 red = [1.00, 0.84, 0.61, 0.34, 0.00]
00244 green = [1.00, 0.84, 0.61, 0.34, 0.00]
00245 blue = [1.00, 1.00, 1.00, 1.00, 1.00]
00246 elif name == "reds":
00247 stops = [0.00, 0.34, 0.61, 0.84, 1.00]
00248 red = [1.00, 1.00, 1.00, 1.00, 1.00]
00249 green = [1.00, 0.84, 0.61, 0.34, 0.00]
00250 blue = [1.00, 0.84, 0.61, 0.34, 0.00]
00251 elif name == "antigray":
00252 stops = [0.00, 0.34, 0.61, 0.84, 1.00]
00253 red = [1.00, 0.84, 0.61, 0.34, 0.00]
00254 green = [1.00, 0.84, 0.61, 0.34, 0.00]
00255 blue = [1.00, 0.84, 0.61, 0.34, 0.00]
00256 red.reverse()
00257 green.reverse()
00258 blue.reverse()
00259 elif name == "fire":
00260 stops = [0.00, 0.20, 0.80, 1.00]
00261 red = [1.00, 1.00, 1.00, 0.50]
00262 green = [1.00, 1.00, 0.00, 0.00]
00263 blue = [0.20, 0.00, 0.00, 0.00]
00264 elif name == "antifire":
00265 stops = [0.00, 0.20, 0.80, 1.00]
00266 red = [0.50, 1.00, 1.00, 1.00]
00267 green = [0.00, 0.00, 1.00, 1.00]
00268 blue = [0.00, 0.00, 0.00, 0.20]
00269 else:
00270
00271 stops = [0.00, 0.34, 0.61, 0.84, 1.00]
00272 red = [0.00, 0.00, 0.87, 1.00, 0.51]
00273 green = [0.00, 0.81, 1.00, 0.20, 0.00]
00274 blue = [0.51, 1.00, 0.12, 0.00, 0.00]
00275
00276 s = array.array('d', stops)
00277 r = array.array('d', red)
00278 g = array.array('d', green)
00279 b = array.array('d', blue)
00280
00281 npoints = len(s)
00282 ROOT.TColor.CreateGradientColorTable(npoints, s, r, g, b, ncontours)
00283 ROOT.gStyle.SetNumberContours(ncontours)
00284
00285 set_palette()
00286
00287
00288
00289
00290
00291
00292
00293 phiedgesCSC36 = [pi/180.*(-5. + 10.*i) for i in range(36)]
00294 phiedgesCSC18 = [pi/180.*(-5. + 20.*i) for i in range(18)]
00295 phiedges = [
00296 phiedgesCSC36,
00297 phiedgesCSC36,
00298 phiedgesCSC36,
00299 phiedgesCSC36,
00300 phiedgesCSC18,
00301 phiedgesCSC36,
00302 phiedgesCSC18,
00303 phiedgesCSC36,
00304 phiedgesCSC18,
00305 phiedgesCSC36,
00306 [0.35228048120123945, 0.87587781482541827, 1.3994776462193192, 1.923076807996136, 2.4466741416203148, 2.970273973014216,
00307 -2.7893121723885534, -2.2657148387643748, -1.7421150073704739, -1.2185158455936571, -0.69491851196947851, -0.17131868057557731],
00308 [0.22000706229660855, 0.74360690430428489, 1.267204926935573, 1.7908033890915052, 2.3144032310991816, 2.8380012537304697,
00309 -2.9215855912931841, -2.3979857492855081, -1.8743877266542202, -1.3507892644982882, -0.82718942249061178, -0.30359139985932365],
00310 [0.29751957124275596, 0.82111826253905784, 1.3447162969496083, 1.8683158980376524, 2.3919145893339548, 2.915512623744505,
00311 -2.844073082347037, -2.3204743910507353, -1.7968763566401849, -1.2732767555521407, -0.74967806425583894, -0.22608002984528835],
00312 [3.0136655290752188, -2.7530905195097337, -2.2922883025568734, -1.9222915077192773, -1.5707963267948966, -1.2193011458705159,
00313 -0.84930435103291968, -0.38850213408005951, 0.127927124514574, 0.65152597487624719, 1.1322596819239259, 1.5707963267948966,
00314 2.0093329716658674, 2.4900666787135459]]
00315
00316 def phiedges2c():
00317 lines = []
00318 for ed in phiedges[:]:
00319 ed.sort()
00320
00321 ed.extend([999 for n in range(0,37-len(ed))])
00322 lines.append('{' + ', '.join(map(str, ed)) + '}')
00323
00324 res = ', '.join(lines)
00325 ff = open("phiedges_export.h",mode="w")
00326 print>>ff,'double phiedges[14][37] = {' + res + '};'
00327 ff.close()
00328
00329 class SawTeethFunction:
00330 def __init__(self, name):
00331 self.name = name
00332 self.edges = (phiedges[stationIndex(name)])[:]
00333 self.ed = sorted(self.edges)
00334
00335 self.ed.append(pi+1.)
00336 self.n = len(self.edges)
00337 def __call__(self, xx, par):
00338
00339 x = xx[0]
00340 if x < self.ed[0]: x += 2*pi
00341
00342 for i in range(0,self.n):
00343 if x <= self.ed[i]: continue
00344 if x > self.ed[i+1]: continue
00345 return par[i*2] + par[i*2+1]*(x - self.ed[i])
00346 return 0
00347 def pp(self):
00348 print self.name, self.n
00349 print self.edges
00350 print self.ed
00351
00352
00353 def stationIndex(name):
00354 if ("MB" in name or "ME" in name):
00355
00356 pa = idToPostalAddress(name)
00357 if pa is None: return None
00358 if pa[0]=="CSC":
00359 if pa[2]==1 and pa[3]==1: return 0
00360 if pa[2]==1 and pa[3]==2: return 1
00361 if pa[2]==1 and pa[3]==3: return 2
00362 if pa[2]==1 and pa[3]==4: return 3
00363 if pa[2]==2 and pa[3]==1: return 4
00364 if pa[2]==2 and pa[3]==2: return 5
00365 if pa[2]==3 and pa[3]==1: return 6
00366 if pa[2]==3 and pa[3]==2: return 7
00367 if pa[2]==4 and pa[3]==1: return 8
00368 if pa[2]==4 and pa[3]==2: return 9
00369 if pa[0]=="DT":
00370 if pa[2]==1: return 10
00371 if pa[2]==2: return 11
00372 if pa[2]==3: return 12
00373 if pa[2]==4: return 13
00374 else:
00375 if ("mem11" in name or "mep11" in name): return 0
00376 if ("mem12" in name or "mep12" in name): return 1
00377 if ("mem13" in name or "mep13" in name): return 2
00378 if ("mem14" in name or "mep14" in name): return 3
00379 if ("mem21" in name or "mep21" in name): return 4
00380 if ("mem22" in name or "mep22" in name): return 5
00381 if ("mem31" in name or "mep31" in name): return 6
00382 if ("mem32" in name or "mep32" in name): return 7
00383 if ("mem41" in name or "mep41" in name): return 8
00384 if ("mem42" in name or "mep42" in name): return 9
00385 if ("st1" in name): return 10
00386 if ("st2" in name): return 11
00387 if ("st3" in name): return 12
00388 if ("st4" in name): return 13
00389
00390
00391
00392 def philines(name, window, abscissa):
00393 global philine_tlines, philine_labels
00394 philine_tlines = []
00395 edges = phiedges[stationIndex(name)]
00396
00397 for phi in edges:
00398 if abscissa is None or abscissa[0] < phi < abscissa[1]:
00399 philine_tlines.append(ROOT.TLine(phi, -window, phi, window))
00400 philine_tlines[-1].SetLineStyle(2)
00401 philine_tlines[-1].Draw()
00402 if "st" in name:
00403 philine_labels = []
00404 edges = edges[:]
00405 edges.sort()
00406 if "st4" in name:
00407 labels = [" 7", " 8", " 9", "14", "10", "11", "12", " 1", " 2", " 3", "13", " 4", " 5", " 6"]
00408 else:
00409 labels = [" 8", " 9", "10", "11", "12", " 1", " 2", " 3", " 4", " 5", " 6"]
00410 edges = edges[1:]
00411 for phi, label in zip(edges, labels):
00412 littlebit = 0.
00413 if label in (" 7", " 9", "14", "10", "11"): littlebit = 0.05
00414 philine_labels.append(ROOT.TText(phi-0.35+littlebit, -0.9*window, label))
00415 philine_labels[-1].Draw()
00416 philine_labels.append(ROOT.TText(-2.9, -0.75*window, "Sector:"))
00417 philine_labels[-1].Draw()
00418 if "CSC" in name:
00419 philine_labels = []
00420 edges = edges[:]
00421 edges.sort()
00422 labels = [" 1", " 2", " 3", " 4", " 5", " 6", " 7", " 8", " 9", "10", "11", "12", "13", "14", "15", "16", "17", "18",
00423 "19", "20", "21", "22", "23", "24", "25", "26", "27", "28", "29", "30", "31", "32", "33", "34", "35", "36"]
00424
00425
00426
00427 for phi, label in zip(edges, labels):
00428 littlebit = 0.
00429
00430 philine_labels.append(ROOT.TText(phi+littlebit, -0.9*window, label))
00431 philine_labels[-1].SetTextFont(42)
00432 philine_labels[-1].SetTextSize(0.028)
00433 philine_labels[-1].Draw()
00434 philine_labels.append(ROOT.TText(0, -0.78*window, "Chamber:"))
00435 philine_labels[-1].SetTextSize(0.035)
00436 philine_labels[-1].Draw()
00437
00438 def zlines(window, abscissa):
00439 global zline_tlines
00440 zline_tlines = []
00441 for z in -401.625, -133.875, 133.875, 401.625:
00442 if abscissa is None or abscissa[0] < z < abscissa[1]:
00443 zline_tlines.append(ROOT.TLine(z, -window, z, window))
00444 zline_tlines[-1].SetLineStyle(2)
00445 zline_tlines[-1].Draw()
00446 zline_labels = []
00447 zline_labels.append(ROOT.TText(-550, -0.9*window, "-2"))
00448 zline_labels.append(ROOT.TText(-300, -0.9*window, "-1"))
00449 zline_labels.append(ROOT.TText(-10, -0.9*window, "0"))
00450 zline_labels.append(ROOT.TText(250, -0.9*window, "+1"))
00451 zline_labels.append(ROOT.TText(500, -0.9*window, "+2"))
00452 for z in zline_labels: z.Draw()
00453 zline_labels.append(ROOT.TText(-600, -0.75*window, "Wheel:")); zline_labels[-1].Draw()
00454
00455 def rlines(disk, window, abscissa):
00456 global rline_tlines
00457 rline_tlines = []
00458 if disk == 1: rl = [150., 270., 480.]
00459 else: rl = [350.]
00460 for r in rl:
00461 if abscissa is None or abscissa[0] < r < abscissa[1]:
00462 rline_tlines.append(ROOT.TLine(r, -window, r, window))
00463 rline_tlines[-1].SetLineStyle(2)
00464 rline_tlines[-1].Draw()
00465
00466
00467
00468 def getReportByPostalAddress(postal_address, report):
00469 for r in report:
00470 if postal_address == r.postal_address:
00471 return r
00472 return None
00473
00474
00475
00476
00477 def DBMC(database, reports, window=10., windows=None, selection=None, phi=False,
00478 color=ROOT.kBlue-8, style=1, bins=50, normalized=False, getvalues=False, name="", canvas=None, reportdiff=False, inlog=True):
00479 return DBdiff(database, None, reports, None, window, windows, selection, phi, color, style, bins, normalized, getvalues,
00480 name, canvas, reportdiff, inlog)
00481
00482
00483 def DBdiff(database1, database2, reports1, reports2,
00484 window=10., windows=None, selection=None, phi=False, color=ROOT.kBlue-8,
00485 style=1, bins=50, normalized=False, getvalues=False, name="tmp", canvas=None, reportdiff=False, inlog=False ):
00486
00487 tdrStyle.SetOptStat("emrou")
00488 tdrStyle.SetStatW(0.40)
00489
00490 wnd = [window]*6
00491 if windows is not None:
00492 i=0
00493 for w in windows:
00494 wnd[i] = windows[i]
00495 i+=1
00496
00497 global hx, hy, hz, hphix, hphiy, hphiz
00498
00499 if phi:
00500 hx = ROOT.TH1F("%s_phi" % name, "", bins, -wnd[0], wnd[0])
00501 else:
00502 hx = ROOT.TH1F("%s_x" % name, "", bins, -wnd[0], wnd[0])
00503 hy = ROOT.TH1F("%s_y" % name, "", bins, -wnd[1], wnd[1])
00504 hz = ROOT.TH1F("%s_z" % name, "", bins, -wnd[2], wnd[2])
00505 hphix = ROOT.TH1F("%s_phix" % name, "", bins, -wnd[3], wnd[3])
00506 hphiy = ROOT.TH1F("%s_phiy" % name, "", bins, -wnd[4], wnd[4])
00507 hphiz = ROOT.TH1F("%s_phiz" % name, "", bins, -wnd[5], wnd[5])
00508
00509 for r1 in reports1:
00510 if selection is None or (selection.func_code.co_argcount == len(r1.postal_address) and selection(*r1.postal_address)):
00511 if reports2 is None:
00512 r2 = Report(r1.chamberId, r1.postal_address, r1.name)
00513 r2.add_parameters(ValErr(0., 0., 0.), ValErr(0., 0., 0.), ValErr(0., 0., 0.),
00514 ValErr(0., 0., 0.), ValErr(0., 0., 0.), ValErr(0., 0., 0.), 0., 0., 0., 0.)
00515 else:
00516 r2 = getReportByPostalAddress(r1.postal_address, reports2)
00517 if r2 is None: continue
00518
00519 found = False
00520 if r1.postal_address[0] == "DT":
00521 if r1.postal_address[1:] in database1.dt:
00522 found = True
00523 db1 = database1.dt[r1.postal_address[1:]]
00524 if database2 is None:
00525 db2 = DTAlignable()
00526 db2.x = db2.y = db2.z = db2.phix = db2.phiy = db2.phiz = 0.
00527 db2.xx = db2.xy = db2.xz = db2.yx = db2.yy = db2.yz = db2.zx = db2.zy = db2.zz = 0.
00528 else:
00529 db2 = database2.dt[r1.postal_address[1:]]
00530
00531 else:
00532
00533 if r1.postal_address[2]==1 and r1.postal_address[3]==4: continue
00534 if r1.postal_address[1:] in database1.csc:
00535 found = True
00536 db1 = database1.csc[r1.postal_address[1:]]
00537 if database2 is None:
00538 db2 = CSCAlignable()
00539 db2.x = db2.y = db2.z = db2.phix = db2.phiy = db2.phiz = 0.
00540 db2.xx = db2.xy = db2.xz = db2.yx = db2.yy = db2.yz = db2.zx = db2.zy = db2.zz = 0.
00541 else:
00542 db2 = database2.csc[r1.postal_address[1:]]
00543
00544 if found and r1.status == "PASS" and r2.status == "PASS":
00545 if r1.deltax is not None and r2.deltax is not None and r1.deltax.error is not None and \
00546 r2.deltax.error is not None and (r1.deltax.error**2 + r2.deltax.error**2) > 0.:
00547 delta = db1.x - db2.x
00548 if reportdiff: delta -= r1.deltax.value
00549 if normalized:
00550 fill = delta/sqrt(r1.deltax.error**2 + r2.deltax.error**2) * signConventions[r1.postal_address][0]
00551 else:
00552 if phi:
00553 fill = delta/signConventions[r1.postal_address][3] * 1000. * signConventions[r1.postal_address][0]
00554 else:
00555 fill = delta * 10. * signConventions[r1.postal_address][0]
00556 hx.Fill(fill)
00557 if getvalues not in (False, None):
00558 getvalues["x"].append((fill, 10. * sqrt(r1.deltax.error**2 + r2.deltax.error**2)))
00559
00560 if r1.deltay is not None and r2.deltay is not None and r1.deltay.error is not None and \
00561 r2.deltay.error is not None and (r1.deltay.error**2 + r2.deltay.error**2) > 0.:
00562 delta = db1.y - db2.y
00563 if reportdiff: delta -= r1.deltay.value
00564 if normalized:
00565 fill = delta/sqrt(r1.deltay.error**2 + r2.deltay.error**2) * signConventions[r1.postal_address][1]
00566 else:
00567 fill = delta * 10. * signConventions[r1.postal_address][1]
00568 hy.Fill(fill)
00569 if getvalues not in (False, None):
00570 getvalues["y"].append((fill, 10. * sqrt(r1.deltay.error**2 + r2.deltay.error**2)))
00571
00572 if r1.deltaz is not None and r2.deltaz is not None and r1.deltaz.error is not None and \
00573 r2.deltaz.error is not None and (r1.deltaz.error**2 + r2.deltaz.error**2) > 0.:
00574 delta = db1.z - db2.z
00575 if reportdiff: delta -= r1.deltaz.value
00576 if normalized:
00577 fill = delta/sqrt(r1.deltaz.error**2 + r2.deltaz.error**2) * signConventions[r1.postal_address][2]
00578 else:
00579 fill = delta * 10. * signConventions[r1.postal_address][2]
00580 hz.Fill(fill)
00581 if getvalues not in (False, None):
00582 getvalues["z"].append((fill, 10. * sqrt(r1.deltaz.error**2 + r2.deltaz.error**2)))
00583
00584 if r1.deltaphix is not None and r2.deltaphix is not None and r1.deltaphix.error is not None and \
00585 r2.deltaphix.error is not None and (r1.deltaphix.error**2 + r2.deltaphix.error**2) > 0.:
00586 delta = db1.phix - db2.phix
00587 if reportdiff: delta -= r1.deltaphix.value
00588 if normalized:
00589 fill = delta/sqrt(r1.deltaphix.error**2 + r2.deltaphix.error**2)
00590 else:
00591 fill = delta * 1000.
00592 hphix.Fill(fill)
00593 if getvalues not in (False, None):
00594 getvalues["phix"].append((fill, 10. * sqrt(r1.deltaphix.error**2 + r2.deltaphix.error**2)))
00595
00596 if r1.deltaphiy is not None and r2.deltaphiy is not None and r1.deltaphiy.error is not None and \
00597 r2.deltaphiy.error is not None and (r1.deltaphiy.error**2 + r2.deltaphiy.error**2) > 0.:
00598 delta = db1.phiy - db2.phiy
00599 if reportdiff:
00600 delta -= r1.deltaphiy.value
00601 if abs(delta)>0.02/1000: print r1.postal_address, 1000*delta, "=", 1000*db1.phiy - 1000*db2.phiy, "-", 1000*r1.deltaphiy.value, "... ",1000*db1.phiy , 1000*db2.phiy
00602 if normalized:
00603 fill = delta/sqrt(r1.deltaphiy.error**2 + r2.deltaphiy.error**2)
00604 else:
00605 fill = delta * 1000.
00606 hphiy.Fill(fill)
00607 if getvalues not in (False, None):
00608 getvalues["phiy"].append((fill, 10. * sqrt(r1.deltaphiy.error**2 + r2.deltaphiy.error**2)))
00609
00610 if r1.deltaphiz is not None and r2.deltaphiz is not None and r1.deltaphiz.error is not None and \
00611 r2.deltaphiz.error is not None and (r1.deltaphiz.error**2 + r2.deltaphiz.error**2) > 0.:
00612 delta = db1.phiz - db2.phiz
00613 if reportdiff: delta -= r1.deltaphiz.value
00614 if normalized:
00615 fill = delta/sqrt(r1.deltaphiz.error**2 + r2.deltaphiz.error**2)
00616 else:
00617 fill = delta * 1000.
00618 hphiz.Fill(fill)
00619 if getvalues not in (False, None):
00620 getvalues["phiz"].append((fill, 10. * sqrt(r1.deltaphiz.error**2 + r2.deltaphiz.error**2)))
00621
00622 if not normalized:
00623 if phi:
00624 hx.SetXTitle("#delta_{#phi} position (mrad)")
00625 else:
00626 hx.SetXTitle("#delta_{x'} (mm)")
00627 hy.SetXTitle("#delta_{y'} (mm)")
00628 hz.SetXTitle("#delta_{z'} (mm)")
00629 hphix.SetXTitle("#delta_{#phi_{x}} (mrad)")
00630 hphiy.SetXTitle("#delta_{#phi_{y}} (mrad)")
00631 hphiz.SetXTitle("#delta_{#phi_{z}} (mrad)")
00632 if reportdiff:
00633 if phi:
00634 hx.SetXTitle("#delta_{#phi}(XML) - #delta_{#phi}(report) position (mrad)")
00635 else:
00636 hx.SetXTitle("#delta_{x'}(XML) - #delta_{x'}(report) (mm)")
00637 hy.SetXTitle("#delta_{y'}(XML) - #delta_{y'}(report) (mm)")
00638 hz.SetXTitle("#delta_{z'}(XML) - #delta_{z'}(report) (mm)")
00639 hphix.SetXTitle("#delta_{#phi_{x}}(XML) - #delta_{#phi_{x}}(report) (mrad)")
00640 hphiy.SetXTitle("#delta_{#phi_{y}}(XML) - #delta_{#phi_{y}}(report) (mrad)")
00641 hphiz.SetXTitle("#delta_{#phi_{z}}(XML) - #delta_{#phi_{z}}(report) (mrad)")
00642 else:
00643 if phi:
00644 hx.SetXTitle("#delta_{#phi}/#sigma_{#phi} position")
00645 else:
00646 hx.SetXTitle("#delta_{x'}/#sigma_{x'}")
00647 hy.SetXTitle("#delta_{y'}/#sigma_{y'}")
00648 hz.SetXTitle("#delta_{z'}/#sigma_{z'}")
00649 hphix.SetXTitle("#delta_{#phi_{x}}/#sigma_{#phi_{x}}")
00650 hphiy.SetXTitle("#delta_{#phi_{y}}/#sigma_{#phi_{y}}")
00651 hphiz.SetXTitle("#delta_{#phi_{z}}/#sigma_{#phi_{z}}")
00652
00653 for h in hx, hy, hz, hphix, hphiy, hphiz:
00654 h.GetXaxis().CenterTitle()
00655 h.GetYaxis().CenterTitle()
00656 h.SetFillColor(color)
00657 h.SetLineStyle(style)
00658
00659 if canvas is not None: c = canvas
00660 else: c = c1
00661
00662 if normalized:
00663 fx = ROOT.TF1("fx", "%g * exp(-x**2/2.)/sqrt(2.*3.1415926)" % (hx.GetEntries()*2.*window/bins), -window, window)
00664 fy = ROOT.TF1("fy", "%g * exp(-x**2/2.)/sqrt(2.*3.1415926)" % (hy.GetEntries()*2.*window/bins), -window, window)
00665 fz = ROOT.TF1("fz", "%g * exp(-x**2/2.)/sqrt(2.*3.1415926)" % (hz.GetEntries()*2.*window/bins), -window, window)
00666 fphix = ROOT.TF1("fphix", "%g * exp(-x**2/2.)/sqrt(2.*3.1415926)" % (hphix.GetEntries()*2.*window/bins), -window, window)
00667 fphiy = ROOT.TF1("fphiy", "%g * exp(-x**2/2.)/sqrt(2.*3.1415926)" % (hphiy.GetEntries()*2.*window/bins), -window, window)
00668 fphiz = ROOT.TF1("fphiz", "%g * exp(-x**2/2.)/sqrt(2.*3.1415926)" % (hphiz.GetEntries()*2.*window/bins), -window, window)
00669 for f in fx, fy, fz, fphix, fphiy, fphiz:
00670 f.SetLineWidth(2)
00671 f.SetLineColor(ROOT.kBlue)
00672 for h, f in (hx, fx), (hy, fy), (hz, fz), (hphix, fphix), (hphiy, fphiy), (hphiz, fphiz):
00673 h.SetAxisRange(0, 1.1*max(h.GetMaximum(), f.GetMaximum()), "Y")
00674
00675 c.Clear()
00676 c.Divide(3, 2)
00677 c.GetPad(1).cd(); hx.Draw(); fx.Draw("same")
00678 c.GetPad(2).cd(); hy.Draw(); fy.Draw("same")
00679 c.GetPad(3).cd(); hz.Draw(); fz.Draw("same")
00680 c.GetPad(4).cd(); hphix.Draw(); fphix.Draw("same")
00681 c.GetPad(5).cd(); hphiy.Draw(); fphiy.Draw("same")
00682 c.GetPad(6).cd(); hphiz.Draw(); fphiz.Draw("same")
00683 return hx, hy, hz, hphix, hphiy, hphiz, fx, fy, fz, fphix, fphiy, fphiz
00684 else:
00685 nvar = 6
00686
00687 c.Clear()
00688 if nvar == 4: c.Divide(2, 2)
00689 if nvar == 6: c.Divide(3, 2)
00690 c.GetPad(1).cd(); hx.Draw()
00691 c.GetPad(2).cd(); hy.Draw()
00692 if nvar == 4:
00693 c.GetPad(3).cd(); hphiy.Draw()
00694 c.GetPad(4).cd(); hphiz.Draw()
00695 if nvar == 6:
00696 c.GetPad(3).cd(); hz.Draw()
00697 c.GetPad(4).cd(); hphix.Draw()
00698 c.GetPad(5).cd(); hphiy.Draw()
00699 c.GetPad(6).cd(); hphiz.Draw()
00700
00701 if inlog:
00702 if hx.GetEntries()>0: c.GetPad(1).SetLogy(1)
00703 if hy.GetEntries()>0: c.GetPad(2).SetLogy(1)
00704 if nvar == 4:
00705 if hphiy.GetEntries()>0: c.GetPad(3).SetLogy(1)
00706 if hphiz.GetEntries()>0: c.GetPad(4).SetLogy(1)
00707 if nvar == 6:
00708 if hz.GetEntries()>0: c.GetPad(3).SetLogy(1)
00709 if hphix.GetEntries()>0: c.GetPad(4).SetLogy(1)
00710 if hphiy.GetEntries()>0: c.GetPad(5).SetLogy(1)
00711 if hphiz.GetEntries()>0: c.GetPad(6).SetLogy(1)
00712
00713 return hx, hy, hz, hphix, hphiy, hphiz
00714
00715
00716
00717 def DBMCVersus(quantity, versus, database, reports, window=10., selection=None, color=ROOT.kBlack):
00718 return DBdiffVersus(quantity, versus, database, None, reports, None, window, selection, color)
00719
00720 def DBdiffVersus(quantity, versus, database1, database2, reports1, reports2, windwselection=None, color=ROOT.kBlack):
00721 tdrStyle.SetOptStat("")
00722
00723 domain = []
00724 values = []
00725 errors = []
00726
00727 for r1 in reports1:
00728 if selection is None or (selection.func_code.co_argcount == len(r1.postal_address) and selection(*r1.postal_address)):
00729 if reports2 is None:
00730 r2 = Report(r1.chamberId, r1.postal_address, r1.name)
00731 r2.add_parameters(ValErr(0., 0., 0.), ValErr(0., 0., 0.), ValErr(0., 0., 0.),
00732 ValErr(0., 0., 0.), ValErr(0., 0., 0.), ValErr(0., 0., 0.), 0., 0., 0.)
00733 else:
00734 found = False
00735 for r2 in reports2:
00736 if r1.postal_address == r2.postal_address:
00737 found = True
00738 break
00739 if not found: continue
00740
00741 found = False
00742 if r1.postal_address[0] == "DT":
00743 if r1.postal_address[1:] in database1.dt:
00744 found = True
00745 db1 = database1.dt[r1.postal_address[1:]]
00746 if database2 is None:
00747 db2 = DTAlignable()
00748 db2.x = db2.y = db2.z = db2.phix = db2.phiy = db2.phiz = 0.
00749 db2.xx = db2.xy = db2.xz = db2.yx = db2.yy = db2.yz = db2.zx = db2.zy = db2.zz = 0.
00750 else:
00751 db2 = database2.dt[r1.postal_address[1:]]
00752 else:
00753 if r1.postal_address[1:] in database1.csc:
00754 found = True
00755 db1 = database1.csc[r1.postal_address[1:]]
00756 if database2 is None:
00757 db2 = CSCAlignable()
00758 db2.x = db2.y = db2.z = db2.phix = db2.phiy = db2.phiz = 0.
00759 db2.xx = db2.xy = db2.xz = db2.yx = db2.yy = db2.yz = db2.zx = db2.zy = db2.zz = 0.
00760 else:
00761 db2 = database2.csc[r1.postal_address[1:]]
00762
00763 if found and r1.status == "PASS" and r2.status == "PASS":
00764 okay = False
00765
00766 if quantity == "phi":
00767 if r1.deltax is not None and r2.deltax is not None and r1.deltax.error is not None and \
00768 r2.deltax.error is not None and (r1.deltax.error**2 + r2.deltax.error**2) > 0.:
00769 okay = True
00770 values.append((db1.x - db2.x)/
00771 signConventions[r1.postal_address][3] * 1000. * signConventions[r1.postal_address][0])
00772 errors.append((r1.deltax.error**2 + r2.deltax.error**2)/
00773 signConventions[r1.postal_address][3] * 1000. * signConventions[r1.postal_address][0])
00774
00775 elif quantity == "x":
00776 if r1.deltax is not None and r2.deltax is not None and r1.deltax.error is not None and \
00777 r2.deltax.error is not None and (r1.deltax.error**2 + r2.deltax.error**2) > 0.:
00778 okay = True
00779 values.append((db1.x - db2.x) * 10. * signConventions[r1.postal_address][0])
00780 errors.append((r1.deltax.error**2 + r2.deltax.error**2) * 10. * signConventions[r1.postal_address][0])
00781
00782 elif quantity == "y":
00783 if r1.deltay is not None and r2.deltay is not None and r1.deltay.error is not None and \
00784 r2.deltay.error is not None and (r1.deltay.error**2 + r2.deltay.error**2) > 0.:
00785 okay = True
00786 values.append((db1.y - db2.y) * 10. * signConventions[r1.postal_address][1])
00787 errors.append((r1.deltay.error**2 + r2.deltay.error**2) * 10. * signConventions[r1.postal_address][1])
00788
00789 elif quantity == "z":
00790 if r1.deltaz is not None and r2.deltaz is not None and r1.deltaz.error is not None and \
00791 r2.deltaz.error is not None and (r1.deltaz.error**2 + r2.deltaz.error**2) > 0.:
00792 okay = True
00793 values.append((db1.z - db2.z) * 10. * signConventions[r1.postal_address][2])
00794 errors.append((r1.deltaz.error**2 + r2.deltaz.error**2) * 10. * signConventions[r1.postal_address][2])
00795
00796 elif quantity == "phix":
00797 if r1.deltaphix is not None and r2.deltaphix is not None and r1.deltaphix.error is not None and \
00798 r2.deltaphix.error is not None and (r1.deltaphix.error**2 + r2.deltaphix.error**2) > 0.:
00799 okay = True
00800 values.append((db1.phix - db2.phix) * 1000.)
00801 errors.append((r1.deltaphix.error**2 + r2.deltaphix.error**2) * 1000.)
00802
00803 elif quantity == "phiy":
00804 if r1.deltaphiy is not None and r2.deltaphiy is not None and r1.deltaphiy.error is not None and \
00805 r2.deltaphiy.error is not None and (r1.deltaphiy.error**2 + r2.deltaphiy.error**2) > 0.:
00806 okay = True
00807 values.append((db1.phiy - db2.phiy) * 1000.)
00808 errors.append((r1.deltaphiy.error**2 + r2.deltaphiy.error**2) * 1000.)
00809
00810 elif quantity == "phiz":
00811 if r1.deltaphiz is not None and r2.deltaphiz is not None and r1.deltaphiz.error is not None and \
00812 r2.deltaphiz.error is not None and (r1.deltaphiz.error**2 + r2.deltaphiz.error**2) > 0.:
00813 okay = True
00814 values.append((db1.phiz - db2.phiz) * 1000.)
00815 errors.append((r1.deltaphiz.error**2 + r2.deltaphiz.error**2) * 1000.)
00816
00817 else: raise Exception
00818
00819 if okay:
00820 if versus == "r": domain.append(signConventions[r1.postal_address][3])
00821 elif versus == "phi": domain.append(signConventions[r1.postal_address][4])
00822 elif versus == "z": domain.append(signConventions[r1.postal_address][5])
00823 else: raise Exception
00824
00825 if versus == "r":
00826 bkgndhist = ROOT.TH1F("bkgndhist", "", 100, 0., 800.)
00827 bkgndhist.SetXTitle("R (cm)")
00828 elif versus == "phi":
00829 bkgndhist = ROOT.TH1F("bkgndhist", "", 100, -pi, pi)
00830 bkgndhist.SetXTitle("#phi (rad)")
00831 elif versus == "z":
00832 bkgndhist = ROOT.TH1F("bkgndhist", "", 100, -1100., 1100.)
00833 bkgndhist.SetXTitle("z (cm)")
00834 bkgndhist.GetXaxis().CenterTitle()
00835
00836 bkgndhist.SetAxisRange(-window, window, "Y")
00837 if quantity == "phi": bkgndhist.SetYTitle("#delta_{#phi} position (mrad)")
00838 elif quantity == "x": bkgndhist.SetYTitle("#delta_{x'} (mm)")
00839 elif quantity == "y": bkgndhist.SetYTitle("#delta_{y'} (mm)")
00840 elif quantity == "z": bkgndhist.SetYTitle("#delta_{z'} (mm)")
00841 elif quantity == "phix": bkgndhist.SetYTitle("#delta_{#phi_{x}} (mrad)")
00842 elif quantity == "phiy": bkgndhist.SetYTitle("#delta_{#phi_{y}} (mrad)")
00843 elif quantity == "phiz": bkgndhist.SetYTitle("#delta_{#phi_{z}} (mrad)")
00844 else: raise Exception
00845 bkgndhist.GetYaxis().CenterTitle()
00846
00847 if len(domain) == 0:
00848 tgraph = ROOT.TGraphErrors(0)
00849 else:
00850 tgraph = ROOT.TGraphErrors(len(domain), array.array("d", domain), array.array("d", values),
00851 array.array("d", [0.]*len(domain)), array.array("d", errors))
00852 tgraph.SetMarkerColor(color)
00853 tgraph.SetLineColor(color)
00854
00855 bkgndhist.Draw()
00856 if tgraph.GetN() > 0: tgraph.Draw("p")
00857 return bkgndhist, tgraph, domain, values, errors
00858
00859
00860
00861 def idToPostalAddress(id):
00862
00863 if len(id)!=9: return None
00864 if id[0:2]=="MB":
00865
00866 pa = ("DT", int(id[2:4]), int(id[5]), int(id[7:9]))
00867
00868 if pa[1]<-2 or pa[1]>2: return None
00869 if pa[2]>4: return None
00870 if pa[3]<1 or pa[3]>14 or (pa[3]==4 and pa[3]>12): return None
00871 return pa
00872 elif id[0:2]=="ME":
00873 if id[2]=="+": ec=1
00874 elif id[2]=="-": ec=2
00875 else: return None
00876 pa = ("CSC", ec, int(id[3]), int(id[5]), int(id[7:9]))
00877 if pa[2]<1 or pa[2]>4: return None
00878 if pa[3]<1 or pa[3]>4 or (pa[2]>1 and pa[3]>2): return None
00879 if pa[4]<1 or pa[4]>36 or (pa[2]>1 and pa[3]==1 and pa[4]>18): return None
00880 return pa
00881 else: return None
00882
00883
00884 def postalAddressToId(postal_address):
00885 if postal_address[0] == "DT":
00886 wheel, station, sector = postal_address[1:]
00887 w = "%+d"%wheel
00888 if w=="+0": w = "-0"
00889 return "MB%s/%d/%02d" % (w, station, sector)
00890 elif postal_address[0] == "CSC":
00891 endcap, station, ring, chamber = postal_address[1:]
00892 if endcap != 1: station = -1 * abs(station)
00893 return "ME%+d/%d/%02d" % (station, ring, chamber)
00894
00895
00896 def nameToId(name):
00897 if name[0:2] == "MB":
00898 wh = name[4]
00899 if wh == "A": w = "-2"
00900 elif wh == "B": w = "-1"
00901 elif wh == "C": w = "-0"
00902 elif wh == "D": w = "+1"
00903 elif wh == "E": w = "+2"
00904 else: return ""
00905 station = name[7]
00906 sector = name[11:13]
00907 return "MB%s/%s/%s" % (w, station, sector)
00908 elif name[0:2] == "ME":
00909 if name[2]=="p": endcap = "+"
00910 elif name[2]=="m": endcap = "-"
00911 else: return ""
00912 station = name[3]
00913 ring = name[4]
00914 chamber = name[6:8]
00915 return "ME%s%s/%s/%s" % (endcap, station, ring, chamber)
00916 return None
00917
00918
00919 def availableCellsDT(reports):
00920 dts = []
00921
00922 for iwheel in DT_TYPES:
00923 if iwheel[1]=="ALL": continue
00924 dts.append(iwheel[0])
00925
00926 for wheel in DT_TYPES:
00927 if wheel[1]=="ALL": continue
00928 for station in wheel[2]:
00929 dts.append(wheel[0]+'/'+station[1])
00930
00931 for wheel in DT_TYPES:
00932 if wheel[1]!="ALL": continue
00933 for station in wheel[2]:
00934 for sector in range(1,station[2]+1):
00935 ssector = "%02d" % sector
00936 dts.append(wheel[0]+'/'+station[1]+'/'+ssector)
00937
00938 for wheel in DT_TYPES:
00939 if wheel[1]!="ALL": continue
00940 for station in wheel[2]:
00941 dts.append(wheel[0]+'/'+station[1])
00942
00943 for wheel in DT_TYPES:
00944 if wheel[1]=="ALL": continue
00945 for station in wheel[2]:
00946 for sector in range(1,station[2]+1):
00947 ssector = "%02d" % sector
00948 label = "MBwh%sst%ssec%s" % (wheelLetter(int(wheel[1])),station[1],ssector)
00949 if len(reports)==0:
00950
00951
00952 continue
00953 found = False
00954 for r in reports:
00955 if r.name == label:
00956 found = True
00957 break
00958 if not found: continue
00959 if r.status == "TOOFEWHITS" and r.posNum+r.negNum==0: continue
00960 if r.status == "NOFIT": continue
00961 dts.append(wheel[0]+'/'+station[1]+'/'+ssector)
00962 return dts
00963
00964
00965 def availableCellsCSC(reports):
00966 cscs = []
00967
00968 for endcap in CSC_TYPES:
00969 for station in endcap[2]:
00970 cscs.append("%s%s" % (endcap[0], station[1]))
00971
00972 for endcap in CSC_TYPES:
00973 for station in endcap[2]:
00974 for ring in station[2]:
00975 if ring[1]=="ALL": continue
00976
00977 cscs.append("%s%s/%s" % (endcap[0], station[1],ring[1]))
00978
00979 for endcap in CSC_TYPES:
00980 for station in endcap[2]:
00981 for ring in station[2]:
00982 if ring[1]!="ALL": continue
00983 for chamber in range(1,ring[2]+1):
00984
00985 cscs.append("%s%s/ALL/%02d" % (endcap[0], station[1],chamber))
00986
00987 for endcap in CSC_TYPES:
00988 for station in endcap[2]:
00989 for ring in station[2]:
00990 if ring[1]!="ALL": continue
00991
00992 cscs.append("%s%s/ALL" % (endcap[0], station[1]))
00993
00994 for endcap in CSC_TYPES:
00995 for station in endcap[2]:
00996 for ring in station[2]:
00997 if ring[1]=="ALL": continue
00998 for chamber in range(1,ring[2]+1):
00999
01000 if station[1]=="4" and ring[1]=="2":
01001 if endcap[1]=="m": continue
01002 if chamber<9 or chamber>13: continue
01003 schamber = "%02d" % chamber
01004 label = "ME%s%s%s_%s" % (endcap[1], station[1], ring[1], schamber)
01005 if len(reports)==0:
01006
01007
01008 continue
01009 found = False
01010 for r in reports:
01011 if r.name == label:
01012 found = True
01013 break
01014 if not found: continue
01015 if r.status == "TOOFEWHITS" and r.posNum+r.negNum==0: continue
01016 if r.status == "NOFIT": continue
01017 cscs.append(endcap[0]+station[1]+'/'+ring[1]+'/'+schamber)
01018 return cscs
01019
01020
01021 DQM_SEVERITY = [
01022 {"idx":0, "name": "NONE", "color": "lightgreen", "hex":"#90EE90"},
01023 {"idx":1, "name": "LOWSTAT05", "color": "lightgreen", "hex":"#96D953"},
01024 {"idx":2, "name": "LOWSTAT075", "color": "lightgreen", "hex":"#94E26F"},
01025 {"idx":3, "name": "LOWSTAT1", "color": "yellowgreen", "hex":"#9ACD32"},
01026 {"idx":4, "name": "LOWSTAT", "color": "yellow", "hex":"#FFFF00"},
01027 {"idx":5, "name": "TOLERABLE", "color": "lightpink", "hex":"#FFB6C1"},
01028 {"idx":6, "name": "SEVERE", "color": "orange", "hex":"#FFA500"},
01029 {"idx":7, "name": "CRITICAL", "color": "red", "hex":"#FF0000"}];
01030
01031
01032 def addToTestResults(c,res):
01033 if len(res)>0:
01034 if TEST_RESULTS.has_key(c): TEST_RESULTS[c].extend(res)
01035 else: TEST_RESULTS[c] = res
01036
01037
01038 def testEntry(testID,scope,descr,severity):
01039 s = 0
01040 for sev in DQM_SEVERITY:
01041 if sev["name"]==severity: s = sev["idx"]
01042 return {"testID":testID,"scope":scope,"descr":descr,"severity":s}
01043
01044
01045 def testZeroWithin5Sigma(x):
01046 if abs(x[1])==0.: return 0.
01047 pull = abs(x[0])/abs(x[1])
01048 if pull <= 5: return 0.
01049 else: return pull
01050
01051
01052 def testDeltaWithin5Sigma(x,sx):
01053 n = len(x)
01054 res = []
01055 dr = []
01056
01057
01058 for i in range(1,n+1):
01059 x1 = x[i-1]
01060 sx1 = sx[i-1]
01061 x2 = x[0]
01062 sx2 = sx[0]
01063 if i<n:
01064 x2 = x[i]
01065 sx2 = sx[i]
01066 sig1 = sqrt( (sx1[0]-sx1[1])**2 + x1[1]**2 )
01067 sig2 = sqrt( (sx2[0]-sx2[1])**2 + x2[1]**2 )
01068 df = abs(x1[0]-x2[0]) - 3*( sig1 + sig2 )
01069
01070
01071 dr.append(df)
01072 if df > 0: res.append(i)
01073
01074
01075 return res
01076
01077
01078 def doTestsForReport(cells,reports):
01079 for c in cells:
01080
01081 postal_address = idToPostalAddress(c)
01082 if not postal_address: continue
01083
01084
01085 found = False
01086 for r in reports:
01087 if r.postal_address == postal_address:
01088 found = True
01089 break
01090 if not found: continue
01091
01092
01093 res = []
01094 scope = postal_address[0]
01095
01096
01097 if r.status == "FAIL" or r.status == "MINUITFAIL":
01098 res.append(testEntry("FAILURE",scope,r.status+" failure","CRITICAL"))
01099 addToTestResults(c,res)
01100 continue
01101
01102
01103 nseg = r.posNum + r.negNum
01104 if r.status == "TOOFEWHITS" and nseg>0:
01105 res.append(testEntry("LOW_STAT",scope,"low stat, #segments=%d"%nseg,"LOWSTAT"))
01106 addToTestResults(c,res)
01107 continue
01108
01109
01110 sdx = 10.*r.deltax.error
01111 if sdx>0.5:
01112 if sdx<0.75: res.append(testEntry("LOW_STAT_DDX05",scope,"low stat, delta(dx)=%f #segments=%d" % (sdx,nseg),"LOWSTAT05"))
01113 elif sdx<1.: res.append(testEntry("LOW_STAT_DDX075",scope,"low stat, delta(dx)=%f #segments=%d" % (sdx,nseg),"LOWSTAT075"))
01114 else: res.append(testEntry("LOW_STAT_DDX1",scope,"low stat, delta(dx)=%f #segments=%d" % (sdx,nseg),"LOWSTAT1"))
01115
01116
01117 if r.redchi2 > 20.:
01118 res.append(testEntry("BIG_CHI2",scope,"chi2=%f>20" % r.redchi2,"TOLERABLE"))
01119
01120
01121 medx, meddx = 10.*r.median_x, 1000.*r.median_dxdz
01122
01123 if medx>2: res.append(testEntry("BIG_MED_X",scope,"median dx=%f>2 mm"%medx,"SEVERE"))
01124
01125 if meddx>2: res.append(testEntry("BIG_MED_DXDZ",scope,"median d(dx/dz)=%f>2 mrad"%meddx,"SEVERE"))
01126
01127
01128
01129 isDTst4 = False
01130 if postal_address[0] == "DT" and postal_address[2]==4: isDTst4 = True
01131 dx, dy, dpy, dpz = 10.*r.deltax.value, 0., 1000.*r.deltaphiy.value, 1000.*r.deltaphiz.value
01132 if not isDTst4: dy = 10.*r.deltay.value
01133 if dx>0.2: res.append(testEntry("BIG_LAST_ITR_DX",scope,"dx=%f>0.2 mm"%dx,"CRITICAL"))
01134 if dy>0.2: res.append(testEntry("BIG_LAST_ITR_DY",scope,"dy=%f>0.2 mm"%dy,"CRITICAL"))
01135 if dpy>0.2: res.append(testEntry("BIG_LAST_ITR_DPHIY",scope,"dphiy=%f>0.2 mrad"%dpy,"CRITICAL"))
01136 if dpz>0.2: res.append(testEntry("BIG_LAST_ITR_DPHIZ",scope,"dphiz=%f>0.2 mrad"%dpz,"CRITICAL"))
01137
01138
01139 addToTestResults(c,res)
01140
01141
01142 def doTestsForMapPlots(cells):
01143 for c in cells:
01144 res = []
01145
01146 scope = "zzz"
01147 if c[0:2]=="MB": scope = "DT"
01148 if c[0:2]=="ME": scope = "CSC"
01149 if scope == "zzz":
01150 print "strange cell ID: ", c
01151 return None
01152
01153 if MAP_RESULTS_FITSIN.has_key(c):
01154 t = MAP_RESULTS_FITSIN[c]
01155 t_a = testZeroWithin5Sigma(t['a'])
01156 t_s = testZeroWithin5Sigma(t['sin'])
01157 t_c = testZeroWithin5Sigma(t['cos'])
01158 if t_a+t_s+t_c >0:
01159 descr = "map fitsin 5 sigma away from 0; pulls : a=%.2f sin=%.2f, cos=%.2f" % (t_a,t_s,t_c)
01160 res.append(testEntry("MAP_FITSIN",scope,descr,"SEVERE"))
01161
01162 if MAP_RESULTS_SAWTOOTH.has_key(c):
01163 t = MAP_RESULTS_SAWTOOTH[c]
01164
01165 t_a = testDeltaWithin5Sigma(t['a'],t['da'])
01166 if len(t_a)>0:
01167 descr = "map discontinuities: %s" % ",".join(map(str,t_a))
01168 res.append(testEntry("MAP_DISCONTIN",scope,descr,"SEVERE"))
01169
01170 t_b = map(testZeroWithin5Sigma, t['b'])
01171 t_bi = []
01172 for i in range(0,len(t_b)):
01173 if t_b[i]>0: t_bi.append(i+1)
01174 if len(t_bi)>0:
01175 descr = "map sawteeth: %s" % ",".join(map(str,t_bi))
01176 res.append(testEntry("MAP_SAWTEETH",scope,descr,"TOLERABLE"))
01177
01178 addToTestResults(c,res)
01179
01180
01181 def saveTestResultsMap(run_name):
01182 if len(MAP_RESULTS_SAWTOOTH)+len(MAP_RESULTS_FITSIN)==0: return None
01183 ff = open("tmp_test_results_map__%s.pkl" % run_name, "wb")
01184 pickle.dump(MAP_RESULTS_SAWTOOTH, ff)
01185 pickle.dump(MAP_RESULTS_FITSIN, ff)
01186 ff.close()
01187
01188
01189 def loadTestResultsMap(run_name):
01190 print "tmp_test_results_map__%s.pkl" % run_name, os.access("tmp_test_results_map__%s.pkl" % run_name,os.F_OK)
01191 if not os.access("tmp_test_results_map__%s.pkl" % run_name,os.F_OK): return None
01192 global MAP_RESULTS_FITSIN, MAP_RESULTS_SAWTOOTH
01193 ff = open("tmp_test_results_map__%s.pkl" % run_name, "rb")
01194 MAP_RESULTS_SAWTOOTH = pickle.load(ff)
01195 MAP_RESULTS_FITSIN = pickle.load(ff)
01196 ff.close()
01197
01198
01199 return True
01200
01201 def writeDQMReport(fname_dqm, run_name):
01202 tests = []
01203 for c in TEST_RESULTS:
01204 tests.append({"objID":c, "name":c, "list":TEST_RESULTS[c]})
01205 lt = time.localtime(time.time())
01206 lts = "%04d-%02d-%02d %02d:%02d:%02d %s" % (lt[0], lt[1], lt[2], lt[3], lt[4], lt[5], time.tzname[1])
01207 dqm_report = {"run":run_name, "genDate": lts, "report":tests}
01208 ff = open(fname_dqm,mode="w")
01209 print >>ff, "var DQM_REPORT = "
01210 json.dump(dqm_report,ff)
01211
01212 ff.close()
01213
01214
01215 def doTests(reports, pic_ids, fname_base, fname_dqm, run_name):
01216
01217 dts = []
01218 cscs = []
01219 if len(reports)>0:
01220 dts = availableCellsDT(reports)
01221 cscs = availableCellsCSC(reports)
01222 elif len(pic_ids)>0:
01223 dts = [id for id in pic_ids if 'MB' in id]
01224 cscs = [id for id in pic_ids if 'ME' in id]
01225 mulist = ['Run: '+run_name,['ALL',['MU']],['DT',dts],['CSC',cscs]]
01226 ff = open(fname_base,mode="w")
01227 print >>ff, "var MU_LIST = ["
01228 json.dump(mulist,ff)
01229 print >>ff, "];"
01230 ff.close()
01231
01232 doTestsForReport(dts,reports)
01233 doTestsForReport(cscs,reports)
01234
01235 loadTestResultsMap(run_name)
01236 doTestsForMapPlots(dts)
01237 doTestsForMapPlots(cscs)
01238
01239 writeDQMReport(fname_dqm, run_name)
01240
01241
01242
01243
01244 def plotmedians(reports1, reports2, selection=None, binsx=100, windowx=5., ceilingx=None, binsy=100, windowy=5.,
01245 ceilingy=None, binsdxdz=100, windowdxdz=5., ceilingdxdz=None, binsdydz=100, windowdydz=5., ceilingdydz=None,
01246 r1text=" before", r2text=" after", which="median"):
01247 tdrStyle.SetOptStat("emrou")
01248 tdrStyle.SetStatW(0.40)
01249 tdrStyle.SetStatFontSize(0.05)
01250
01251 global hmediandxdz_after, hmediandxdz_before, hmediandxdz_beforecopy, \
01252 hmediandydz_after, hmediandydz_before, hmediandydz_beforecopy, \
01253 hmedianx_after, hmedianx_before, hmedianx_beforecopy, \
01254 hmediany_after, hmediany_before, hmediany_beforecopy, tlegend
01255
01256 hmedianx_before = ROOT.TH1F("hmedianx_before", "", binsx, -windowx, windowx)
01257 hmediany_before = ROOT.TH1F("hmediany_before", "", binsy, -windowy, windowy)
01258 hmediandxdz_before = ROOT.TH1F("hmediandxdz_before", "", binsdxdz, -windowdxdz, windowdxdz)
01259 hmediandydz_before = ROOT.TH1F("hmediandydz_before", "", binsdydz, -windowdydz, windowdydz)
01260 hmedianx_after = ROOT.TH1F("hmedianx_after", "", binsx, -windowx, windowx)
01261 hmediany_after = ROOT.TH1F("hmediany_after", "", binsy, -windowy, windowy)
01262 hmediandxdz_after = ROOT.TH1F("hmediandxdz_after", "", binsdxdz, -windowdxdz, windowdxdz)
01263 hmediandydz_after = ROOT.TH1F("hmediandydz_after", "", binsdydz, -windowdydz, windowdydz)
01264
01265 if which == "median":
01266 whichx = whichy = whichdxdz = whichdydz = "median"
01267 elif which == "bigmean":
01268 whichx = "mean30"
01269 whichy = "mean30"
01270 whichdxdz = "mean20"
01271 whichdydz = "mean50"
01272 elif which == "mean":
01273 whichx = "mean15"
01274 whichy = "mean15"
01275 whichdxdz = "mean10"
01276 whichdydz = "mean25"
01277 elif which == "bigwmean":
01278 whichx = "wmean30"
01279 whichy = "wmean30"
01280 whichdxdz = "wmean20"
01281 whichdydz = "wmean50"
01282 elif which == "wmean":
01283 whichx = "wmean15"
01284 whichy = "wmean15"
01285 whichdxdz = "wmean10"
01286 whichdydz = "wmean25"
01287 elif which == "bigstdev":
01288 whichx = "stdev30"
01289 whichy = "stdev30"
01290 whichdxdz = "stdev20"
01291 whichdydz = "stdev50"
01292 elif which == "stdev":
01293 whichx = "stdev15"
01294 whichy = "stdev15"
01295 whichdxdz = "stdev10"
01296 whichdydz = "stdev25"
01297 else:
01298 raise Exception, which + " not recognized"
01299
01300 for r1 in reports1:
01301 if selection is None or (selection.func_code.co_argcount == len(r1.postal_address) and selection(*r1.postal_address)):
01302 found = False
01303 for r2 in reports2:
01304 if r1.postal_address == r2.postal_address:
01305 found = True
01306 break
01307 if not found: continue
01308
01309
01310 if r1.postal_address[0]=='CSC':
01311 if r1.postal_address[2]==1 and r1.postal_address[3]==4: continue
01312
01313 if r1.status == "PASS" and r2.status == "PASS":
01314 hmedianx_before.Fill(10.*eval("r1.%s_x" % whichx))
01315 hmediandxdz_before.Fill(1000.*eval("r1.%s_dxdz" % whichdxdz))
01316 hmedianx_after.Fill(10.*eval("r2.%s_x" % whichx))
01317 hmediandxdz_after.Fill(1000.*eval("r2.%s_dxdz" % whichdxdz))
01318
01319 if r1.median_y is not None:
01320 hmediany_before.Fill(10.*eval("r1.%s_y" % whichy))
01321 hmediandydz_before.Fill(1000.*eval("r1.%s_dydz" % whichdydz))
01322 hmediany_after.Fill(10.*eval("r2.%s_y" % whichy))
01323 hmediandydz_after.Fill(1000.*eval("r2.%s_dydz" % whichdydz))
01324
01325 hmedianx_beforecopy = hmedianx_before.Clone()
01326 hmediany_beforecopy = hmediany_before.Clone()
01327 hmediandxdz_beforecopy = hmediandxdz_before.Clone()
01328 hmediandydz_beforecopy = hmediandydz_before.Clone()
01329 hmedianx_beforecopy.SetLineStyle(2)
01330 hmediany_beforecopy.SetLineStyle(2)
01331 hmediandxdz_beforecopy.SetLineStyle(2)
01332 hmediandydz_beforecopy.SetLineStyle(2)
01333
01334 hmedianx_before.SetFillColor(ROOT.kMagenta+2)
01335 hmediany_before.SetFillColor(ROOT.kMagenta+2)
01336 hmediandxdz_before.SetFillColor(ROOT.kMagenta+2)
01337 hmediandydz_before.SetFillColor(ROOT.kMagenta+2)
01338 hmedianx_after.SetFillColor(ROOT.kYellow)
01339 hmediany_after.SetFillColor(ROOT.kYellow)
01340 hmediandxdz_after.SetFillColor(ROOT.kYellow)
01341 hmediandydz_after.SetFillColor(ROOT.kYellow)
01342
01343 hmedianx_after.SetXTitle("median(#Deltax) (mm)")
01344 hmediany_after.SetXTitle("median(#Deltay) (mm)")
01345 hmediandxdz_after.SetXTitle("median(#Deltadx/dz) (mrad)")
01346 hmediandydz_after.SetXTitle("median(#Deltadydz) (mrad)")
01347 hmedianx_after.GetXaxis().CenterTitle()
01348 hmediany_after.GetXaxis().CenterTitle()
01349 hmediandxdz_after.GetXaxis().CenterTitle()
01350 hmediandydz_after.GetXaxis().CenterTitle()
01351
01352 if ceilingx is not None: hmedianx_after.SetAxisRange(0., ceilingx, "Y")
01353 if ceilingy is not None: hmediany_after.SetAxisRange(0., ceilingy, "Y")
01354 if ceilingdxdz is not None: hmediandxdz_after.SetAxisRange(0., ceilingdxdz, "Y")
01355 if ceilingdydz is not None: hmediandydz_after.SetAxisRange(0., ceilingdydz, "Y")
01356
01357 c1.Clear()
01358 c1.Divide(2, 2)
01359
01360 c1.GetPad(1).cd()
01361 hmedianx_after.Draw()
01362 hmedianx_before.Draw("same")
01363 hmedianx_after.Draw("same")
01364 hmedianx_beforecopy.Draw("same")
01365 hmedianx_after.Draw("axissame")
01366
01367 tlegend = ROOT.TLegend(0.17, 0.75-0.05, 0.45+0.05, 0.9)
01368 tlegend.SetFillColor(ROOT.kWhite)
01369 tlegend.SetBorderSize(0)
01370 tlegend.AddEntry(hmedianx_after, r2text, "f")
01371 tlegend.AddEntry(hmedianx_before, r1text, "f")
01372 tlegend.Draw()
01373
01374 c1.GetPad(2).cd()
01375 hmediandxdz_after.Draw()
01376 hmediandxdz_before.Draw("same")
01377 hmediandxdz_after.Draw("same")
01378 hmediandxdz_beforecopy.Draw("same")
01379 hmediandxdz_after.Draw("axissame")
01380
01381 c1.GetPad(3).cd()
01382 hmediany_after.Draw()
01383 hmediany_before.Draw("same")
01384 hmediany_after.Draw("same")
01385 hmediany_beforecopy.Draw("same")
01386 hmediany_after.Draw("axissame")
01387
01388 c1.GetPad(4).cd()
01389 hmediandydz_after.Draw()
01390 hmediandydz_before.Draw("same")
01391 hmediandydz_after.Draw("same")
01392 hmediandydz_beforecopy.Draw("same")
01393 hmediandydz_after.Draw("axissame")
01394
01395 return hmediandxdz_after, hmediandxdz_before, hmediandxdz_beforecopy, \
01396 hmediandydz_after, hmediandydz_before, hmediandydz_beforecopy, \
01397 hmedianx_after, hmedianx_before, hmedianx_beforecopy, \
01398 hmediany_after, hmediany_before, hmediany_beforecopy, tlegend
01399
01400
01401
01402
01403 def createPeaksProfile(the2d, rebin=1):
01404 htmp = ROOT.gROOT.FindObject(the2d.GetName()+"_peaks")
01405 if htmp != None: htmp.Delete()
01406
01407 hpeaks = the2d.ProjectionX(the2d.GetName()+"_peaks")
01408 hpeaks.Reset()
01409 hpeaks.Rebin(rebin)
01410 bad_fit_bins = []
01411 for i in xrange(0, int(the2d.GetNbinsX()), rebin):
01412 tmp = the2d.ProjectionY("tmp", i+1, i + rebin)
01413 nn = tmp.GetEntries()
01414
01415 drange = tmp.GetRMS()
01416 drange = 2.*drange
01417 fgaus = ROOT.TF1("fgaus","gaus", tmp.GetMean() - drange, tmp.GetMean() + drange)
01418 fgaus.SetParameter(0,nn)
01419 fgaus.SetParameter(1,tmp.GetMean())
01420 fgaus.SetParameter(2,tmp.GetRMS())
01421
01422
01423 fitOk = False
01424 if nn > 10:
01425 fr = tmp.Fit("fgaus","RNSQ")
01426
01427 hpeaks.SetBinContent(i/rebin+1, fgaus.GetParameter(1))
01428 hpeaks.SetBinError(i/rebin+1, fgaus.GetParError(1))
01429 if fr.Status()==0 and fr.CovMatrixStatus()==3 : fitOk = True
01430 if not fitOk:
01431 bad_fit_bins.append(i/rebin+1)
01432 if nn > 1. and tmp.GetRMS() > 0:
01433 hpeaks.SetBinContent(i/rebin+1, tmp.GetMean())
01434 hpeaks.SetBinError(i/rebin+1, ROOT.TMath.StudentQuantile(0.841345,nn) * tmp.GetRMS() / sqrt(nn))
01435 else:
01436 hpeaks.SetBinContent(i/rebin+1, 0.)
01437 hpeaks.SetBinError(i/rebin+1, 0.)
01438 if len(bad_fit_bins): print "createPeaksProfile bad fit bins: ", bad_fit_bins
01439 return hpeaks
01440
01441
01442
01443
01444 def mapplot(tfiles, name, param, mode="from2d", window=10., abscissa=None, title="",
01445 widebins=False, fitsine=False, fitline=False, reset_palette=False, fitsawteeth=False, fitpeaks=False, peaksbins=1, fixfitpars={}, **args):
01446 tdrStyle.SetOptTitle(1)
01447 tdrStyle.SetTitleBorderSize(0)
01448 tdrStyle.SetOptStat(0)
01449
01450 tdrStyle.SetOptFit(0)
01451 tdrStyle.SetTitleFontSize(0.05)
01452 tdrStyle.SetPadRightMargin(0.1)
01453
01454 c1.Clear()
01455 c1.ResetAttPad()
01456
01457 if reset_palette: set_palette("blues")
01458 global hist, hist2d, hist2dweight, tline1, tline2, tline3
01459
01460 if fitsine or fitsawteeth:
01461 id = mapNameToId(name)
01462 if not id:
01463 print "bad id for ", name
01464 raise Exception
01465
01466 hdir = "AlignmentMonitorMuonSystemMap1D/iter1/"
01467 hpref= "%s_%s" % (name, param)
01468 hhh = hdir+hpref
01469
01470 combine_all = False
01471 if "ALL" in name and ("CSCvsr" in name or "DTvsz" in name): combine_all = True
01472
01473 add1d = ("vsphi" in name) and (param == "x")
01474
01475 if "h2d" in args:
01476 hist2d = args["h2d"].Clone(hpref+"_2d_")
01477 if "CSC" in name and add1d: hist1d = args["h1d"].Clone(hpref+"_1d_")
01478
01479 elif combine_all:
01480 nch = 12
01481 if "DT" in name and name[6:9]=='st4': nch = 14
01482 if "CSC" in name: nch = 36
01483 chambers = ["%02d" % ch for ch in range (2,nch+1)]
01484
01485 ch_hhh = hhh.replace('ALL','01')
01486 ch_hpref = hpref.replace('ALL','01')
01487 hist2d = tfiles[0].Get(ch_hhh+"_2d").Clone(ch_hpref+"_2d_")
01488 if "CSC" in name and add1d: hist1d = tfiles[0].Get(ch_hhh+"_1d").Clone(ch_hpref+"_1d_")
01489
01490 for ch in chambers:
01491 ch_hhh = hhh.replace('ALL',ch)
01492 ch_hpref = hpref.replace('ALL',ch)
01493 hist2d.Add(tfiles[0].Get(ch_hhh+"_2d"))
01494 if "CSC" in name and add1d: hist1d.Add(tfiles[0].Get(ch_hhh+"_1d"))
01495 for tfile in tfiles[1:]:
01496 hist2d.Add(tfile.Get(ch_hhh+"_2d"))
01497 if "CSC" in name and add1d: hist1d.Add(tfile.Get(ch_hhh+"_1d"))
01498
01499 else:
01500 hist2d = tfiles[0].Get(hhh+"_2d").Clone(hpref+"_2d_")
01501 if "CSC" in name and add1d: hist1d = tfiles[0].Get(hhh+"_1d").Clone(hpref+"_1d_")
01502 for tfile in tfiles[1:]:
01503 hist2d.Add(tfile.Get(hhh+"_2d"))
01504 if "CSC" in name and add1d: hist1d.Add(tfile.Get(hhh+"_1d"))
01505
01506
01507 if mode == "from2d":
01508 the2d = hist2d
01509
01510 hist = the2d.ProjectionX()
01511 hist.Reset()
01512
01513 skip = 1
01514 if widebins:
01515 hist.Rebin(10)
01516 skip = 10
01517
01518
01519 for i in xrange(0, int(the2d.GetNbinsX()), skip):
01520 tmp = the2d.ProjectionY("tmp", i+1, i + skip)
01521 if tmp.GetEntries() > 1:
01522
01523 hist.SetBinContent(i/skip+1, tmp.GetMean())
01524 hist.SetBinError(i/skip+1, ROOT.TMath.StudentQuantile(0.841345,tmp.GetEntries()) * tmp.GetRMS() / sqrt(tmp.GetEntries()))
01525
01526
01527 else:
01528
01529
01530 hist.SetBinContent(i/skip+1, 0.)
01531 hist.SetBinError(i/skip+1, 0.)
01532
01533 hpeaks = createPeaksProfile(the2d, peaksbins)
01534
01535 else:
01536 raise Exception
01537
01538 hist.SetAxisRange(-window, window, "Y")
01539 if abscissa is not None: hist.SetAxisRange(abscissa[0], abscissa[1], "X")
01540 hist.SetMarkerStyle(20)
01541 hist.SetMarkerSize(0.75)
01542 hist.GetXaxis().CenterTitle()
01543 hist.GetYaxis().CenterTitle()
01544 hist.GetYaxis().SetTitleOffset(0.75)
01545 hist.GetXaxis().SetTitleSize(0.05)
01546 hist.GetYaxis().SetTitleSize(0.05)
01547 hist.SetTitle(title)
01548 if "vsphi" in name: hist.SetXTitle("Global #phi position (rad)")
01549 elif "vsz" in name: hist.SetXTitle("Global z position (cm)")
01550 elif "vsr" in name: hist.SetXTitle("Global R position (cm)")
01551 if "DT" in name:
01552 if param == "x": hist.SetYTitle("x' residual (mm)")
01553 if param == "dxdz": hist.SetYTitle("dx'/dz residual (mrad)")
01554 if param == "y": hist.SetYTitle("y' residual (mm)")
01555 if param == "dydz": hist.SetYTitle("dy'/dz residual (mrad)")
01556 if "CSC" in name:
01557 if param == "x": hist.SetYTitle("r#phi residual (mm)")
01558 if param == "dxdz": hist.SetYTitle("d(r#phi)/dz residual (mrad)")
01559 hist.SetMarkerColor(ROOT.kBlack)
01560 hist.SetLineColor(ROOT.kBlack)
01561 hist.Draw()
01562 hist2d.Draw("colzsame")
01563 if widebins: hist.Draw("samee1")
01564 else: hist.Draw("same")
01565
01566 hpeaks.SetMarkerStyle(20)
01567 hpeaks.SetMarkerSize(0.9)
01568 hpeaks.SetMarkerColor(ROOT.kRed)
01569 hpeaks.SetLineColor(ROOT.kRed)
01570 hpeaks.SetLineWidth(2)
01571
01572 hpeaks.Draw("same")
01573
01574 if fitsine and "vsphi" in name:
01575 global fitsine_const, fitsine_sin, fitsine_cos, fitsine_chi2, fitsine_ndf
01576 if 'CSC' in name:
01577 f = ROOT.TF1("f", "[0] + [1]*sin(x) + [2]*cos(x)", -pi/180.*5., pi*(2.-5./180.))
01578 else:
01579 f = ROOT.TF1("f", "[0] + [1]*sin(x) + [2]*cos(x)", -pi, pi)
01580 f.SetLineColor(ROOT.kRed)
01581 f.SetLineWidth(2)
01582 if len(fixfitpars)>0:
01583 for fpar in fixfitpars.keys():
01584 f.FixParameter(fpar, fixfitpars[fpar])
01585
01586 if fitpeaks: hpeaks.Fit(f,"NQ")
01587 else: hist.Fit(f,"NEQ")
01588 if len(fixfitpars)>0:
01589 for fpar in fixfitpars.keys():
01590 f.ReleaseParameter(fpar)
01591 fitsine_const = f.GetParameter(0), f.GetParError(0)
01592 fitsine_sin = f.GetParameter(1), f.GetParError(1)
01593 fitsine_cos = f.GetParameter(2), f.GetParError(2)
01594 fitsine_chi2 = f.GetChisquare()
01595 fitsine_ndf = f.GetNDF()
01596 global MAP_RESULTS_FITSIN
01597
01598 MAP_RESULTS_FITSIN[id] = {'a':fitsine_const, 'phi':fitsine_const, 'sin': fitsine_sin, 'cos': fitsine_cos, 'chi2': fitsine_chi2, 'ndf': fitsine_ndf}
01599 f.Draw("same")
01600 global fitsine_ttext, fitsine_etext
01601 text_xposition = -1.
01602 if 'CSC' in name: text_xposition = 2.
01603 fitsine_ttext = ROOT.TLatex(text_xposition, 0.8*window,
01604 "%+.3f %+.3f sin#phi %+.3f cos#phi" % (fitsine_const[0], fitsine_sin[0], fitsine_cos[0]))
01605 fitsine_ttext.SetTextColor(ROOT.kRed)
01606 fitsine_ttext.SetTextSize(0.05)
01607 fitsine_ttext.Draw()
01608 fitsine_etext = ROOT.TLatex(text_xposition, 0.70*window,
01609 " #pm%.3f #pm%.3f #pm%.3f" % (fitsine_const[1], fitsine_sin[1], fitsine_cos[1]))
01610 fitsine_etext.SetTextColor(ROOT.kRed)
01611 fitsine_etext.SetTextSize(0.045)
01612 fitsine_etext.Draw()
01613
01614
01615 if 'CSC' in name and add1d:
01616
01617 f0 = ROOT.TF1("f0", "gaus", hist1d.GetBinLowEdge(1), -hist1d.GetBinLowEdge(1))
01618 fit = hist1d.Fit(f0,"NRQ")
01619 rangea, rangeb = hist1d.GetMean() - hist1d.GetRMS(), hist1d.GetMean() + hist1d.GetRMS()
01620 if fit==0: rangea, rangeb = f0.GetParameter(1) - f0.GetParameter(2), f0.GetParameter(1) + f0.GetParameter(2)
01621
01622
01623
01624 f1 = ROOT.TF1("f1", "gaus", rangea, rangeb)
01625 fit = hist1d.Fit(f1,"NRQ")
01626 nn = hist1d.GetEntries()
01627 dphiz, ephiz = 0, 0
01628 if nn>0: dphiz, ephiz = hist1d.GetMean(), ROOT.TMath.StudentQuantile(0.841345,nn) * hist1d.GetRMS() / sqrt(nn)
01629 if fit==0: dphiz, ephiz = f1.GetParameter(1), f1.GetParError(1)
01630
01631 MAP_RESULTS_FITSIN[id]['phi'] = (dphiz, ephiz)
01632
01633 global ttex_sine_, ttex_sine, ttex_1d_, ttex_1d
01634 postal_address = idToPostalAddress(id+'/01')
01635 ttex_sine_ = ROOT.TLatex(0, 0.8*window,"#Delta#phi_{z}^{sine} (mrad):")
01636 ttex_sine_.SetTextColor(ROOT.kGreen+2); ttex_sine_.SetTextSize(0.04); ttex_sine_.Draw()
01637 ttex_sine = ROOT.TLatex(0, 0.7*window," %+.3f#pm%.3f" %
01638 (-100*fitsine_const[0]/signConventions[postal_address][3],
01639 100*fitsine_const[1]/signConventions[postal_address][3]))
01640 ttex_sine.SetTextColor(ROOT.kGreen+2); ttex_sine.SetTextSize(0.04); ttex_sine.Draw()
01641 ttex_1d_ = ROOT.TLatex(0, 0.6*window,"#Delta#phi_{z}^{phi} (mrad):")
01642 ttex_1d_.SetTextColor(ROOT.kGreen+2); ttex_1d_.SetTextSize(0.04); ttex_1d_.Draw()
01643 ttex_1d = ROOT.TLatex(0, 0.5*window," %+.3f#pm%.3f" % (-dphiz, ephiz))
01644 ttex_1d.SetTextColor(ROOT.kGreen+2); ttex_1d.SetTextSize(0.04); ttex_1d.Draw()
01645 ROOT.gPad.Update()
01646
01647 if fitline:
01648 f = ROOT.TF1("f", "[0] + [1]*x", -1000., 1000.)
01649 hist2d.Fit(f, "q")
01650 hist2d.GetFunction("f").SetLineColor(ROOT.kRed)
01651 global fitline_const, fitline_linear, fitline_chi2, fitline_ndf
01652 fitline_const = hist2d.GetFunction("f").GetParameter(0), hist2d.GetFunction("f").GetParError(0)
01653 fitline_linear = hist2d.GetFunction("f").GetParameter(1), hist2d.GetFunction("f").GetParError(1)
01654 fitline_chi2 = hist2d.GetFunction("f").GetChisquare()
01655 fitline_ndf = hist2d.GetFunction("f").GetNDF()
01656 hist2d.GetFunction("f").Draw("same")
01657 global fitline_ttext
01658 if "vsz" in name: which = "Z"
01659 elif "vsr" in name: which = "R"
01660 fitline_ttext = ROOT.TText(hist.GetBinCenter(hist.GetNbinsX()/4),
01661 0.8*window, "%.3g %+.3g %s" % (fitline_const[0], fitline_linear[0], which))
01662 fitline_ttext.SetTextColor(ROOT.kRed)
01663 fitline_ttext.Draw()
01664
01665 ROOT.gPad.RedrawAxis()
01666
01667 if "vsphi" in name:
01668 if not widebins: philines(name, window, abscissa)
01669 if abscissa is None:
01670 if 'CSC' in name:
01671 tline1 = ROOT.TLine(-pi/180.*5., 0, pi*(2.-5./180.), 0); tline1.Draw()
01672 tline2 = ROOT.TLine(-pi/180.*5., -window, pi*(2.-5./180.), -window); tline2.SetLineWidth(2); tline2.Draw()
01673 tline3 = ROOT.TLine(-pi/180.*5., window, pi*(2.-5./180.), window); tline3.Draw()
01674 else:
01675 tline1 = ROOT.TLine(-pi, 0, pi, 0); tline1.Draw()
01676 tline2 = ROOT.TLine(-pi, -window, pi, -window); tline2.SetLineWidth(2); tline2.Draw()
01677 tline3 = ROOT.TLine(-pi, window, pi, window); tline3.Draw()
01678 else:
01679 tline1 = ROOT.TLine(abscissa[0], 0, abscissa[1], 0); tline1.Draw()
01680 tline2 = ROOT.TLine(abscissa[0], -window, abscissa[1], -window); tline2.SetLineWidth(2); tline2.Draw()
01681 tline3 = ROOT.TLine(abscissa[0], window, abscissa[1], window); tline3.Draw()
01682 elif "vsz" in name:
01683 if not widebins: zlines(window, abscissa)
01684 if abscissa is None:
01685 tline1 = ROOT.TLine(-660, 0, 660, 0); tline1.Draw()
01686 tline2 = ROOT.TLine(-660, -window, 660, -window); tline2.SetLineWidth(2); tline2.Draw()
01687 tline3 = ROOT.TLine(-660, window, 660, window); tline3.Draw()
01688 else:
01689 tline1 = ROOT.TLine(abscissa[0], 0, abscissa[1], 0); tline1.Draw()
01690 tline2 = ROOT.TLine(abscissa[0], -window, abscissa[1], -window); tline2.SetLineWidth(2); tline2.Draw()
01691 tline3 = ROOT.TLine(abscissa[0], window, abscissa[1], window); tline3.Draw()
01692 elif "vsr" in name:
01693 if "mem1" in name or "mep1" in name and not widebins: rlines(1, window, abscissa)
01694 if "mem2" in name or "mep2" in name and not widebins: rlines(2, window, abscissa)
01695 if "mem3" in name or "mep3" in name and not widebins: rlines(3, window, abscissa)
01696 if "mem4" in name or "mep4" in name and not widebins: rlines(4, window, abscissa)
01697 if abscissa is None:
01698 tline1 = ROOT.TLine(100, 0, 700, 0); tline1.Draw()
01699 tline2 = ROOT.TLine(100, -window, 700, -window); tline2.SetLineWidth(2); tline2.Draw()
01700 tline3 = ROOT.TLine(100, window, 700, window); tline3.Draw()
01701 else:
01702 tline1 = ROOT.TLine(abscissa[0], 0, abscissa[1], 0); tline1.Draw()
01703 tline2 = ROOT.TLine(abscissa[0], -window, abscissa[1], -window); tline2.SetLineWidth(2); tline2.Draw()
01704 tline3 = ROOT.TLine(abscissa[0], window, abscissa[1], window); tline3.Draw()
01705
01706 if "vsphi" in name and fitsawteeth:
01707 global CPP_LOADED
01708 if not CPP_LOADED:
01709 phiedges2c()
01710 ROOT.gROOT.ProcessLine(".L phiedges_fitfunctions.C++")
01711 CPP_LOADED = True
01712 fn={0: ROOT.fitf0,
01713 1: ROOT.fitf2,
01714 2: ROOT.fitf2,
01715 3: ROOT.fitf3,
01716 4: ROOT.fitf4,
01717 5: ROOT.fitf5,
01718 6: ROOT.fitf6,
01719 7: ROOT.fitf7,
01720 8: ROOT.fitf8,
01721 9: ROOT.fitf9,
01722 10: ROOT.fitf10,
01723 11: ROOT.fitf11,
01724 12: ROOT.fitf12,
01725 13: ROOT.fitf13
01726 } [stationIndex(name)]
01727 fn.SetNpx(5000)
01728 fn.SetLineColor(ROOT.kYellow)
01729 hist.Fit(fn,"N")
01730 fn.Draw("same")
01731
01732
01733 edges = (phiedges[stationIndex(name)])[:]
01734 ed = sorted(edges)
01735
01736 ed.append(pi+abs(ed[0]))
01737
01738 global sawtooth_a, sawtooth_b
01739 sawtooth_a = []
01740 sawtooth_da = []
01741 sawtooth_b = []
01742 for pr in range(0,fn.GetNpar(),2):
01743 sawtooth_a.append( (fn.GetParameter(pr), fn.GetParError(pr)) )
01744 sawtooth_b.append( (fn.GetParameter(pr+1), fn.GetParError(pr+1)) )
01745 sawtooth_da.append( (fn.Eval(ed[pr/2]+0.01), fn.Eval(ed[pr/2+1]-0.01)) )
01746 global MAP_RESULTS_SAWTOOTH
01747 MAP_RESULTS_SAWTOOTH[id] = {'a': sawtooth_a, 'da': sawtooth_da, 'b': sawtooth_b, 'chi2': fn.GetChisquare(), 'ndf': fn.GetNDF()}
01748
01749
01750
01751
01752
01753 ROOT.SetOwnership(hist2d,False)
01754 ROOT.SetOwnership(hist,False)
01755 ROOT.SetOwnership(tline1,False)
01756 ROOT.SetOwnership(tline2,False)
01757 ROOT.SetOwnership(tline3,False)
01758 return hist
01759
01760
01761 def mapNameToId(name):
01762 if "DT" in name:
01763 wh = "-ALL"
01764 if name.find('wh')>1: wh = name[name.find('wh')+2]
01765 if wh == "A": w = "-2"
01766 elif wh == "B": w = "-1"
01767 elif wh == "C": w = "-0"
01768 elif wh == "D": w = "+1"
01769 elif wh == "E": w = "+2"
01770 elif wh == "-ALL": w = "-ALL"
01771 else: return None
01772 station=''
01773 if wh == "-ALL":
01774 if name.find('sec')<0: return None
01775 station = name[name.find('sec')-1]
01776 sector = ''
01777 sector = name[name.find('sec')+3:name.find('sec')+5]
01778 return "MB%s/%s/%s" % (w, station, sector)
01779 if name.find('st')>1: station = name[name.find('st')+2]
01780 else: return None
01781 return "MB%s/%s" % (w, station)
01782 elif "CSC" in name:
01783 p = name.find('me')
01784 if p<0: return None
01785 if name[p+2]=="p": endcap = "+"
01786 elif name[p+2]=="m": endcap = "-"
01787 else: return None
01788 station = name[p+3]
01789 pch = name.find('ch')
01790 if pch<0:
01791 ring = name[p+4]
01792 return "ME%s%s/%s" % (endcap, station, ring)
01793 ring = 'ALL'
01794 chamber = name[pch+2:pch+4]
01795 return "ME%s%s/%s/%s" % (endcap, station, ring, chamber)
01796 return None
01797
01798
01799
01800
01801
01802
01803
01804 def curvatureplot(tfiles, name, param, mode="from2d", window=15., widebins=False, title="", fitgauss=False, fitconst=False, fitline=False, fitpeaks=True, reset_palette=False):
01805 tdrStyle.SetOptTitle(1)
01806 tdrStyle.SetTitleBorderSize(0)
01807 tdrStyle.SetOptStat(0)
01808 tdrStyle.SetOptFit(0)
01809 tdrStyle.SetTitleFontSize(0.05)
01810
01811 c1.Clear()
01812 if reset_palette: set_palette("blues")
01813 global hist, histCOPY, hist2d, tline1, tline2, tline3, tline4, tline5
01814
01815 hdir = "AlignmentMonitorMuonVsCurvature/iter1/"
01816
01817 if name not in ("all", "top", "bottom"):
01818 hsuffix = "_%s_%s" % (name, param)
01819 prof = tfiles[0].Get(hdir+"tprofile"+hsuffix).Clone("tprofile_"+hsuffix)
01820 hist2d = tfiles[0].Get(hdir+"th2f"+hsuffix).Clone("th2f_"+hsuffix)
01821 for tfile in tfiles[1:]:
01822 prof.Add(tfile.Get(hdir+"tprofile"+hsuffix))
01823 hist2d.Add(tfile.Get(hdir+"th2f"+hsuffix))
01824 else:
01825 prof = None
01826 hist2d = None
01827 for wheel in "m2", "m1", "z", "p1", "p2":
01828 if name == "all": sectors = "01", "02", "03", "04", "05", "06", "07", "08", "09", "10", "11", "12"
01829 elif name == "top": sectors = "01", "02", "03", "04", "05", "06"
01830 elif name == "bottom": sectors = "07", "08", "09", "10", "11", "12"
01831 else: raise Exception
01832
01833 for sector in sectors:
01834 hsuffix = "_%s_%s" % ("wheel%s_sector%s" % (wheel, sector), param)
01835 for tfile in tfiles:
01836 if prof is None:
01837 prof = tfiles[0].Get(hdir+"tprofile"+hsuffix).Clone("tprofile_"+hsuffix)
01838 hist2d = tfiles[0].Get(hdir+"th2f"+hsuffix).Clone("tprofile_"+hsuffix)
01839 else:
01840 prof.Add(tfile.Get(hdir+"tprofile"+hsuffix))
01841 hist2d.Add(tfile.Get(hdir+"th2f"+hsuffix))
01842
01843 hist = ROOT.TH1F("hist", "", prof.GetNbinsX(), prof.GetBinLowEdge(1), -prof.GetBinLowEdge(1))
01844 for i in xrange(1, prof.GetNbinsX()+1):
01845 hist.SetBinContent(i, prof.GetBinContent(i))
01846 hist.SetBinError(i, prof.GetBinError(i))
01847
01848 if mode == "plain":
01849 hist = prof
01850 elif mode == "from2d":
01851 skip = 1
01852 if widebins:
01853 hist.Rebin(5)
01854 skip = 5
01855 htmp = ROOT.gROOT.FindObject("tmp")
01856 if htmp != None: htmp.Delete()
01857
01858 for i in xrange(0, int(prof.GetNbinsX()), skip):
01859 tmp = hist2d.ProjectionY("tmp", i+1, i + skip)
01860 if tmp.GetEntries() > 1:
01861 hist.SetBinContent(i/skip+1, tmp.GetMean())
01862 hist.SetBinError(i/skip+1, ROOT.TMath.StudentQuantile(0.841345,tmp.GetEntries()) * tmp.GetRMS() / sqrt(tmp.GetEntries()))
01863
01864 else:
01865
01866
01867 hist.SetBinContent(i/skip+1, 0.)
01868 hist.SetBinError(i/skip+1, 0.)
01869
01870 hpeaks = createPeaksProfile(hist2d, skip)
01871
01872 else:
01873 raise Exception
01874
01875
01876 if fitgauss:
01877 f = ROOT.TF1("f", "[0] + [1]*exp(-x**2/2/0.01**2)", hist.GetBinLowEdge(1), -hist.GetBinLowEdge(1))
01878 f.SetParameters(0, 0., 0.01)
01879 if fitpeaks: hpeaks.Fit(f, "q")
01880 else: hist.Fit(f, "q")
01881 f.SetLineColor(ROOT.kRed)
01882 global fitgauss_diff, fitgauss_chi2, fitgauss_ndf
01883
01884
01885
01886 fitgauss_diff = f.GetParameter(1), f.GetParError(1)
01887 fitgauss_chi2 = f.GetChisquare()
01888 fitgauss_ndf = f.GetNDF()
01889
01890 global fitline_intercept, fitline_slope
01891 if fitconst:
01892 f = ROOT.TF1("f", "[0]", hist.GetBinLowEdge(1), -hist.GetBinLowEdge(1))
01893 if fitpeaks: hpeaks.Fit(f, "q")
01894 else: hist.Fit(f, "q")
01895 f.SetLineColor(ROOT.kRed)
01896 fitline_intercept = f.GetParameter(0), f.GetParError(0)
01897
01898 if fitline:
01899 f = ROOT.TF1("f", "[0] + [1]*x", hist.GetBinLowEdge(1), -hist.GetBinLowEdge(1))
01900 if fitpeaks: hpeaks.Fit(f, "qNE")
01901 else: hist.Fit(f, "qNE")
01902 f.SetLineColor(ROOT.kRed)
01903 global f2, f3
01904 f2 = ROOT.TF1("2", "[0] + [1]*x", hist.GetBinLowEdge(1), -hist.GetBinLowEdge(1))
01905 f3 = ROOT.TF1("2", "[0] + [1]*x", hist.GetBinLowEdge(1), -hist.GetBinLowEdge(1))
01906 f2.SetParameters(f.GetParameter(0), f.GetParameter(1) + f.GetParError(1))
01907 f3.SetParameters(f.GetParameter(0), f.GetParameter(1) - f.GetParError(1))
01908 f2.SetLineColor(ROOT.kRed)
01909 f3.SetLineColor(ROOT.kRed)
01910 f2.SetLineStyle(2)
01911 f3.SetLineStyle(2)
01912 fitline_intercept = f.GetParameter(0), f.GetParError(0)
01913 fitline_slope = f.GetParameter(1), f.GetParError(1)
01914
01915 hist2d.SetAxisRange(-window, window, "Y")
01916 hist2d.SetMarkerStyle(20)
01917 hist2d.SetMarkerSize(0.75)
01918 hist2d.GetXaxis().CenterTitle()
01919 hist2d.GetYaxis().CenterTitle()
01920 if param == "curverr":
01921 hist2d.GetYaxis().SetTitleOffset(1.35)
01922 else:
01923 hist2d.GetYaxis().SetTitleOffset(0.75)
01924 hist2d.GetXaxis().SetTitleOffset(1.2)
01925 hist2d.GetXaxis().SetTitleSize(0.05)
01926 hist2d.GetYaxis().SetTitleSize(0.05)
01927 hist2d.SetTitle(title)
01928 if param == "pterr": hist2d.SetXTitle("qp_{T} (GeV/c)")
01929 else: hist2d.SetXTitle("q/p_{T} (c/GeV)")
01930 if param == "deltax": hist2d.SetYTitle("#Deltax' (mm)")
01931 if param == "deltadxdz": hist2d.SetYTitle("#Deltadx'/dz (mrad)")
01932 if param == "pterr": hist2d.SetYTitle("#Deltap_{T}/p_{T} (%)")
01933 if param == "curverr": hist2d.SetYTitle("#Deltaq/p_{T} (c/GeV)")
01934 hist2d.Draw("colz")
01935 hist.SetMarkerColor(ROOT.kBlack)
01936 hist.SetLineColor(ROOT.kBlack)
01937 hist.Draw("same")
01938
01939
01940
01941
01942
01943
01944
01945
01946
01947
01948
01949 if fitline:
01950 f.Draw("same")
01951
01952
01953
01954 hpeaks.SetMarkerStyle(20)
01955 hpeaks.SetMarkerSize(0.9)
01956 hpeaks.SetMarkerColor(ROOT.kRed)
01957 hpeaks.SetLineColor(ROOT.kRed)
01958 hpeaks.SetLineWidth(2)
01959
01960 hpeaks.Draw("same")
01961
01962
01963
01964
01965
01966 tline5 = ROOT.TLine(-hist.GetBinLowEdge(1), 0., hist.GetBinLowEdge(1), 0.)
01967 tline5.Draw()
01968
01969
01970
01971 def curvatureDTsummary(tfiles, window=15., pdgSfactor=False):
01972 global h, gm2, gm1, gz, gp1, gp2, tlegend
01973
01974 set_palette("blues")
01975 phis = {-2: [], -1: [], 0: [], 1: [], 2: []}
01976 diffs = {-2: [], -1: [], 0: [], 1: [], 2: []}
01977 differrs = {-2: [], -1: [], 0: [], 1: [], 2: []}
01978 for wheelstr, wheel in ("m2", "-2"), ("m1", "-1"), ("z", "0"), ("p1", "+1"), ("p2", "+2"):
01979 for sector in "01", "02", "03", "04", "05", "06", "07", "08", "09", "10", "11", "12":
01980 curvatureplot(tfiles, "wheel%s_sector%s" % (wheelstr, sector), "deltax",
01981 title="Wheel %s, sector %s" % (wheel, sector), fitgauss=True, reset_palette=False)
01982 if fitgauss_diff[1] < window:
01983 uncertainty = fitgauss_diff[1]
01984 if pdgSfactor and (fitgauss_chi2/fitgauss_ndf) > 1.: uncertainty *= sqrt(fitgauss_chi2/fitgauss_ndf)
01985
01986 phis[int(wheel)].append(signConventions["DT", int(wheel), 1, int(sector)][4])
01987 diffs[int(wheel)].append(fitgauss_diff[0])
01988 differrs[int(wheel)].append(uncertainty)
01989
01990 h = ROOT.TH1F("h", "", 1, -pi, pi)
01991 h.SetAxisRange(-window, window, "Y")
01992 h.SetXTitle("#phi (rad)")
01993 h.SetYTitle("#Deltax(p_{T} #rightarrow #infty) - #Deltax(p_{T} #rightarrow 0) (mm)")
01994 h.GetXaxis().CenterTitle()
01995 h.GetYaxis().CenterTitle()
01996
01997 gm2 = ROOT.TGraphErrors(len(phis[-2]), array.array("d", phis[-2]), array.array("d", diffs[-2]),
01998 array.array("d", [0.]*len(phis[-2])), array.array("d", differrs[-2]))
01999 gm1 = ROOT.TGraphErrors(len(phis[-1]), array.array("d", phis[-1]), array.array("d", diffs[-1]),
02000 array.array("d", [0.]*len(phis[-1])), array.array("d", differrs[-1]))
02001 gz = ROOT.TGraphErrors(len(phis[0]), array.array("d", phis[0]), array.array("d", diffs[0]),
02002 array.array("d", [0.]*len(phis[0])), array.array("d", differrs[0]))
02003 gp1 = ROOT.TGraphErrors(len(phis[1]), array.array("d", phis[1]), array.array("d", diffs[1]),
02004 array.array("d", [0.]*len(phis[1])), array.array("d", differrs[1]))
02005 gp2 = ROOT.TGraphErrors(len(phis[2]), array.array("d", phis[2]), array.array("d", diffs[2]),
02006 array.array("d", [0.]*len(phis[2])), array.array("d", differrs[2]))
02007
02008 gm2.SetMarkerStyle(21); gm2.SetMarkerColor(ROOT.kRed); gm2.SetLineColor(ROOT.kRed)
02009 gm1.SetMarkerStyle(22); gm1.SetMarkerColor(ROOT.kBlue); gm1.SetLineColor(ROOT.kBlue)
02010 gz.SetMarkerStyle(3); gz.SetMarkerColor(ROOT.kBlack); gz.SetLineColor(ROOT.kBlack)
02011 gp1.SetMarkerStyle(26); gp1.SetMarkerColor(ROOT.kBlue); gp1.SetLineColor(ROOT.kBlue)
02012 gp2.SetMarkerStyle(25); gp2.SetMarkerColor(ROOT.kRed); gp2.SetLineColor(ROOT.kRed)
02013
02014 h.Draw()
02015 tlegend = ROOT.TLegend(0.25, 0.2, 0.85, 0.5)
02016 tlegend.SetFillColor(ROOT.kWhite)
02017 tlegend.SetBorderSize(0)
02018 tlegend.AddEntry(gm2, "Wheel -2", "p")
02019 tlegend.AddEntry(gm1, "Wheel -1", "p")
02020 tlegend.AddEntry(gz, "Wheel 0", "p")
02021 tlegend.AddEntry(gp1, "Wheel +1", "p")
02022 tlegend.AddEntry(gp2, "Wheel +2", "p")
02023 tlegend.Draw()
02024
02025 gm2.Draw("p")
02026 gm1.Draw("p")
02027 gz.Draw("p")
02028 gp1.Draw("p")
02029 gp2.Draw("p")
02030
02031
02032 def getname(r):
02033 if r.postal_address[0] == "DT":
02034 wheel, station, sector = r.postal_address[1:]
02035 return "DT wheel %d, station %d, sector %d" % (wheel, station, sector)
02036 elif r.postal_address[0] == "CSC":
02037 endcap, station, ring, chamber = r.postal_address[1:]
02038 if endcap != 1: station = -1 * abs(station)
02039 return "CSC ME%d/%d chamber %d" % (station, ring, chamber)
02040
02041 ddt=[0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.]
02042 def clearDDT():
02043 for i in range(0,15):
02044 ddt[i]=0.
02045
02046 def printDeltaTs():
02047 n = 0
02048 for t in ddt:
02049 if n==0 or n==7 or n==15: print "%d calls" % t
02050 else: print "%d : %0.3f ms" % (n,t*1000.0)
02051 n += 1
02052
02053 def bellcurves(tfile, reports, name, twobin=True, suppressblue=False):
02054 t1 = time.time()
02055 ddt[0] += 1
02056 tdrStyle.SetOptTitle(1)
02057 tdrStyle.SetTitleBorderSize(1)
02058 tdrStyle.SetTitleFontSize(0.1)
02059 tdrStyle.SetOptStat(0)
02060 tdrStyle.SetHistMinimumZero()
02061
02062 c1.Clear()
02063 c1.ResetAttPad()
02064
02065 found = False
02066 for r in reports:
02067 if r.name == name:
02068 found = True
02069 break
02070 if not found: raise Exception, "Not a valid name"
02071 if r.status == "FAIL":
02072
02073 print "Fit failed"
02074 c1.Clear()
02075 return
02076
02077 Pos = "Pos"; Neg = "Neg"
02078 if not twobin:
02079 Pos = ""; Neg = ""
02080
02081 pdirPos = "MuonAlignmentFromReference/%s%s" % (name, Pos)
02082 pdirNeg = "MuonAlignmentFromReference/%s%s" % (name, Neg)
02083
02084 t2 = time.time()
02085 ddt[1] = 1./ddt[0]*((ddt[0]-1)*ddt[1] + t2-t1)
02086
02087 chamber_x = tfile.Get(pdirPos+"_x")
02088 chamber_x_fit = tfile.Get(pdirPos+"_x_fit")
02089 chamber_y = tfile.Get(pdirPos+"_y")
02090 chamber_y_fit = tfile.Get(pdirPos+"_y_fit")
02091 chamber_dxdz = tfile.Get(pdirPos+"_dxdz")
02092 chamber_dxdz_fit = tfile.Get(pdirPos+"_dxdz_fit")
02093 chamber_dydz = tfile.Get(pdirPos+"_dydz")
02094 chamber_dydz_fit = tfile.Get(pdirPos+"_dydz_fit")
02095 chamber_alphax = tfile.Get(pdirPos+"_alphax")
02096 chamber_alphax_fit = tfile.Get(pdirPos+"_alphax_fit")
02097 chamber_alphay = tfile.Get(pdirPos+"_alphay")
02098 chamber_alphay_fit = tfile.Get(pdirPos+"_alphay_fit")
02099 if twobin:
02100 chamber_x_fit2 = tfile.Get(pdirNeg+"_x_fit")
02101 chamber_y_fit2 = tfile.Get(pdirNeg+"_y_fit")
02102 chamber_dxdz_fit2 = tfile.Get(pdirNeg+"_dxdz_fit")
02103 chamber_dydz_fit2 = tfile.Get(pdirNeg+"_dydz_fit")
02104 chamber_alphax_fit2 = tfile.Get(pdirNeg+"_alphax_fit")
02105 chamber_alphay_fit2 = tfile.Get(pdirNeg+"_alphay_fit")
02106
02107 if not chamber_x:
02108 chamber_x = tfile.Get(pdirPos+"_residual")
02109 chamber_x_fit = tfile.Get(pdirPos+"_residual_fit")
02110 chamber_dxdz = tfile.Get(pdirPos+"_resslope")
02111 chamber_dxdz_fit = tfile.Get(pdirPos+"_resslope_fit")
02112 chamber_alphax = tfile.Get(pdirPos+"_alpha")
02113 chamber_alphax_fit = tfile.Get(pdirPos+"_alpha_fit")
02114 if twobin:
02115 chamber_x_fit2 = tfile.Get(pdirNeg+"_residual_fit")
02116 chamber_dxdz_fit2 = tfile.Get(pdirNeg+"_resslope_fit")
02117 chamber_alphax_fit2 = tfile.Get(pdirNeg+"_alpha_fit")
02118
02119 if not chamber_x:
02120 print "Can't find neither "+pdirPos+"_x nor "+pdirPos+"_residual"
02121 return
02122
02123 t3 = time.time()
02124 ddt[2] = 1./ddt[0]*((ddt[0]-1)*ddt[2] + t3-t2)
02125
02126
02127 chamber_x.SetAxisRange(-50., 50., "X")
02128 if chamber_x.GetRMS()>15: chamber_x.SetAxisRange(-75., 75., "X")
02129 chamber_dxdz.SetAxisRange(-30., 30., "X")
02130 chamber_alphax.SetAxisRange(-50., 50., "X")
02131 if not not chamber_y:
02132 chamber_y.SetAxisRange(-75., 75., "X")
02133 chamber_dydz.SetAxisRange(-120., 120., "X")
02134 chamber_alphay.SetAxisRange(-120., 120., "X")
02135 chamber_alphay.SetAxisRange(-75., 75., "Y")
02136
02137
02138 chamber_x.SetXTitle("Local x residual (mm)")
02139 chamber_dxdz.SetXTitle("Local dx/dz residual (mrad)")
02140 chamber_alphax.SetXTitle("Local dx/dz residual (mrad)")
02141 chamber_alphax.SetYTitle("Local x residual (mm)")
02142 if not not chamber_y:
02143 chamber_y.SetXTitle("Local y residual (mm)")
02144 chamber_dydz.SetXTitle("Local dy/dz residual (mrad)")
02145 chamber_alphay.SetXTitle("Local dy/dz residual (mrad)")
02146 chamber_alphay.SetYTitle("Local y residual (mm)")
02147 if name[0:2] == "ME":
02148 chamber_x.SetXTitle("Local r#phi residual (mm)")
02149 chamber_dxdz.SetXTitle("Local d(r#phi)/dz residual (mrad)")
02150 chamber_alphax.SetXTitle("Local d(r#phi)/dz residual (mrad)")
02151 chamber_alphax.SetYTitle("Local r#phi residual (mm)")
02152
02153 t4 = time.time()
02154 ddt[3] = 1./ddt[0]*((ddt[0]-1)*ddt[3] + t4-t3)
02155
02156 for h in chamber_x, chamber_dxdz, chamber_alphax, chamber_alphax, \
02157 chamber_y, chamber_dydz, chamber_alphay, chamber_alphay:
02158 if not not h:
02159 h.GetXaxis().CenterTitle()
02160 h.GetYaxis().CenterTitle()
02161 h.GetXaxis().SetLabelSize(0.05)
02162 h.GetYaxis().SetLabelSize(0.05)
02163 h.GetXaxis().SetTitleSize(0.07)
02164 h.GetYaxis().SetTitleSize(0.07)
02165 h.GetXaxis().SetTitleOffset(0.9)
02166 h.GetYaxis().SetTitleOffset(0.9)
02167
02168 if twobin:
02169 for f in chamber_x_fit2, chamber_y_fit2, chamber_dxdz_fit2, chamber_dydz_fit2, \
02170 chamber_alphax_fit2, chamber_alphay_fit2:
02171 if not not f:
02172 f.SetLineColor(4)
02173 if not twobin:
02174 suppressblue = True
02175
02176 t5 = time.time()
02177 ddt[4] = 1./ddt[0]*((ddt[0]-1)*ddt[4] + t5-t4)
02178
02179 global l1, l2, l3, l4
02180 if not not chamber_y:
02181 c1.Clear()
02182 c1.Divide(3, 2)
02183 chamber_x.SetTitle(getname(r))
02184
02185 c1.GetPad(1).cd()
02186 chamber_x.Draw()
02187 if not suppressblue: chamber_x_fit2.Draw("same")
02188 chamber_x_fit.Draw("same")
02189 l1 = ROOT.TLatex(0.67,0.8,"#splitline{#mu: %0.2f#pm%0.2f}{#sigma: %0.1f#pm%0.1f}" % (
02190 chamber_x_fit.GetParameter(1), chamber_x_fit.GetParError(1),
02191 chamber_x_fit.GetParameter(2), chamber_x_fit.GetParError(2)))
02192 l1.Draw()
02193
02194 c1.GetPad(2).cd()
02195 chamber_dxdz.Draw()
02196 if not suppressblue: chamber_dxdz_fit2.Draw("same")
02197 chamber_dxdz_fit.Draw("same")
02198 l2 = ROOT.TLatex(0.67,0.8,"#splitline{#mu: %0.2f#pm%0.2f}{#sigma: %0.1f#pm%0.1f}" % (
02199 chamber_dxdz_fit.GetParameter(1), chamber_dxdz_fit.GetParError(1),
02200 chamber_dxdz_fit.GetParameter(2), chamber_dxdz_fit.GetParError(2)))
02201 l2.Draw()
02202
02203 c1.GetPad(3).cd()
02204 chamber_alphax.Draw("col")
02205 if not suppressblue: chamber_alphax_fit2.Draw("same")
02206 chamber_alphax_fit.Draw("same")
02207
02208 c1.GetPad(4).cd()
02209 chamber_y.Draw()
02210 if not suppressblue: chamber_y_fit2.Draw("same")
02211 chamber_y_fit.Draw("same")
02212 l3 = ROOT.TLatex(0.67,0.8,"#splitline{#mu: %0.2f#pm%0.2f}{#sigma: %0.1f#pm%0.1f}" % (
02213 chamber_y_fit.GetParameter(1), chamber_y_fit.GetParError(1),
02214 chamber_y_fit.GetParameter(2), chamber_y_fit.GetParError(2)))
02215 l3.Draw()
02216
02217 c1.GetPad(5).cd()
02218 chamber_dydz.Draw()
02219 if not suppressblue: chamber_dydz_fit2.Draw("same")
02220 chamber_dydz_fit.Draw("same")
02221 l4 = ROOT.TLatex(0.67,0.8,"#splitline{#mu: %0.2f#pm%0.2f}{#sigma: %0.1f#pm%0.1f}" % (
02222 chamber_dydz_fit.GetParameter(1), chamber_dydz_fit.GetParError(1),
02223 chamber_dydz_fit.GetParameter(2), chamber_dydz_fit.GetParError(2)))
02224 l4.Draw()
02225
02226 for lb in l1,l2,l3,l4:
02227 lb.SetNDC(1)
02228 lb.SetTextColor(ROOT.kRed)
02229
02230 c1.GetPad(6).cd()
02231 chamber_alphay.Draw("col")
02232 if not suppressblue: chamber_alphay_fit2.Draw("same")
02233 chamber_alphay_fit.Draw("same")
02234
02235 else:
02236 c1.Clear()
02237 c1.Divide(3, 1)
02238 chamber_x.SetTitle(getname(r))
02239
02240 c1.GetPad(1).cd()
02241 chamber_x.Draw()
02242 if not suppressblue: chamber_x_fit2.Draw("same")
02243 chamber_x_fit.Draw("same")
02244 l1 = ROOT.TLatex(0.67,0.8,"#splitline{#mu: %0.2f#pm%0.2f}{#sigma: %0.1f#pm%0.1f}" % (
02245 chamber_x_fit.GetParameter(1), chamber_x_fit.GetParError(1),
02246 chamber_x_fit.GetParameter(2), chamber_x_fit.GetParError(2)))
02247 l1.Draw()
02248
02249 c1.GetPad(2).cd()
02250 chamber_dxdz.Draw()
02251 if not suppressblue: chamber_dxdz_fit2.Draw("same")
02252 chamber_dxdz_fit.Draw("same")
02253 l2 = ROOT.TLatex(0.67,0.8,"#splitline{#mu: %0.2f#pm%0.2f}{#sigma: %0.1f#pm%0.1f}" % (
02254 chamber_dxdz_fit.GetParameter(1), chamber_dxdz_fit.GetParError(1),
02255 chamber_dxdz_fit.GetParameter(2), chamber_dxdz_fit.GetParError(2)))
02256 l2.Draw()
02257
02258 c1.GetPad(3).cd()
02259 chamber_alphax.Draw("col")
02260 if not suppressblue: chamber_alphax_fit2.Draw("same")
02261 chamber_alphax_fit.Draw("same")
02262
02263 for lb in l1,l2:
02264 lb.SetNDC(1)
02265 lb.SetTextColor(ROOT.kRed)
02266
02267 t6 = time.time()
02268 ddt[5] = 1./ddt[0]*((ddt[0]-1)*ddt[5] + t6-t5)
02269 ddt[6] = 1./ddt[0]*((ddt[0]-1)*ddt[6] + t6-t1)
02270
02271
02272 def polynomials(tfile, reports, name, twobin=True, suppressblue=False):
02273 t1 = time.time()
02274 ddt[7] += 1
02275 global label1, label2, label3, label4, label5, label6, label7, label8, label9
02276 plotDirectory = "MuonAlignmentFromReference"
02277 tdrStyle.SetOptTitle(1)
02278 tdrStyle.SetTitleBorderSize(1)
02279 tdrStyle.SetTitleFontSize(0.1)
02280 tdrStyle.SetOptStat(0)
02281
02282 c1.Clear()
02283 c1.ResetAttPad()
02284
02285 found = False
02286 for r in reports:
02287 if r.name == name:
02288 found = True
02289 break
02290 if not found: raise Exception, "Not a valid name"
02291
02292 if r.status == "FAIL":
02293
02294 print "Fit failed"
02295 c1.Clear()
02296 return
02297
02298 Pos = "Pos"; Neg = "Neg"
02299 if not twobin:
02300 Pos = ""; Neg = ""
02301
02302 pdirPos = "MuonAlignmentFromReference/%s%s" % (name, Pos)
02303 pdirNeg = "MuonAlignmentFromReference/%s%s" % (name, Neg)
02304
02305 global chamber_x_trackx, chamber_x_trackx_fit, chamber_y_trackx, chamber_y_trackx_fit, \
02306 chamber_dxdz_trackx, chamber_dxdz_trackx_fit, chamber_dydz_trackx, chamber_dydz_trackx_fit, \
02307 chamber_x_trackx_fit2, chamber_y_trackx_fit2, chamber_dxdz_trackx_fit2, chamber_dydz_trackx_fit2
02308 global chamber_x_tracky, chamber_x_tracky_fit, chamber_y_tracky, chamber_y_tracky_fit, \
02309 chamber_dxdz_tracky, chamber_dxdz_tracky_fit, chamber_dydz_tracky, chamber_dydz_tracky_fit, \
02310 chamber_x_tracky_fit2, chamber_y_tracky_fit2, chamber_dxdz_tracky_fit2, chamber_dydz_tracky_fit2
02311 global chamber_x_trackdxdz, chamber_x_trackdxdz_fit, chamber_y_trackdxdz, chamber_y_trackdxdz_fit, \
02312 chamber_dxdz_trackdxdz, chamber_dxdz_trackdxdz_fit, chamber_dydz_trackdxdz, chamber_dydz_trackdxdz_fit, \
02313 chamber_x_trackdxdz_fit2, chamber_y_trackdxdz_fit2, chamber_dxdz_trackdxdz_fit2, chamber_dydz_trackdxdz_fit2
02314 global chamber_x_trackdydz, chamber_x_trackdydz_fit, chamber_y_trackdydz, chamber_y_trackdydz_fit, \
02315 chamber_dxdz_trackdydz, chamber_dxdz_trackdydz_fit, chamber_dydz_trackdydz, chamber_dydz_trackdydz_fit, \
02316 chamber_x_trackdydz_fit2, chamber_y_trackdydz_fit2, chamber_dxdz_trackdydz_fit2, chamber_dydz_trackdydz_fit2
02317
02318 chamber_x_trackx = tfile.Get(pdirPos+"_x_trackx")
02319 chamber_x_trackx_fit = tfile.Get(pdirPos+"_x_trackx_fitline")
02320 chamber_y_trackx = tfile.Get(pdirPos+"_y_trackx")
02321 chamber_y_trackx_fit = tfile.Get(pdirPos+"_y_trackx_fitline")
02322 chamber_dxdz_trackx = tfile.Get(pdirPos+"_dxdz_trackx")
02323 chamber_dxdz_trackx_fit = tfile.Get(pdirPos+"_dxdz_trackx_fitline")
02324 chamber_dydz_trackx = tfile.Get(pdirPos+"_dydz_trackx")
02325 chamber_dydz_trackx_fit = tfile.Get(pdirPos+"_dydz_trackx_fitline")
02326 chamber_x_trackx_fit2 = tfile.Get(pdirNeg+"_x_trackx_fitline")
02327 chamber_y_trackx_fit2 = tfile.Get(pdirNeg+"_y_trackx_fitline")
02328 chamber_dxdz_trackx_fit2 = tfile.Get(pdirNeg+"_dxdz_trackx_fitline")
02329 chamber_dydz_trackx_fit2 = tfile.Get(pdirNeg+"_dydz_trackx_fitline")
02330
02331 chamber_x_tracky = tfile.Get(pdirPos+"_x_tracky")
02332 chamber_x_tracky_fit = tfile.Get(pdirPos+"_x_tracky_fitline")
02333 chamber_y_tracky = tfile.Get(pdirPos+"_y_tracky")
02334 chamber_y_tracky_fit = tfile.Get(pdirPos+"_y_tracky_fitline")
02335 chamber_dxdz_tracky = tfile.Get(pdirPos+"_dxdz_tracky")
02336 chamber_dxdz_tracky_fit = tfile.Get(pdirPos+"_dxdz_tracky_fitline")
02337 chamber_dydz_tracky = tfile.Get(pdirPos+"_dydz_tracky")
02338 chamber_dydz_tracky_fit = tfile.Get(pdirPos+"_dydz_tracky_fitline")
02339 chamber_x_tracky_fit2 = tfile.Get(pdirNeg+"_x_tracky_fitline")
02340 chamber_y_tracky_fit2 = tfile.Get(pdirNeg+"_y_tracky_fitline")
02341 chamber_dxdz_tracky_fit2 = tfile.Get(pdirNeg+"_dxdz_tracky_fitline")
02342 chamber_dydz_tracky_fit2 = tfile.Get(pdirNeg+"_dydz_tracky_fitline")
02343
02344 chamber_x_trackdxdz = tfile.Get(pdirPos+"_x_trackdxdz")
02345 chamber_x_trackdxdz_fit = tfile.Get(pdirPos+"_x_trackdxdz_fitline")
02346 chamber_y_trackdxdz = tfile.Get(pdirPos+"_y_trackdxdz")
02347 chamber_y_trackdxdz_fit = tfile.Get(pdirPos+"_y_trackdxdz_fitline")
02348 chamber_dxdz_trackdxdz = tfile.Get(pdirPos+"_dxdz_trackdxdz")
02349 chamber_dxdz_trackdxdz_fit = tfile.Get(pdirPos+"_dxdz_trackdxdz_fitline")
02350 chamber_dydz_trackdxdz = tfile.Get(pdirPos+"_dydz_trackdxdz")
02351 chamber_dydz_trackdxdz_fit = tfile.Get(pdirPos+"_dydz_trackdxdz_fitline")
02352 chamber_x_trackdxdz_fit2 = tfile.Get(pdirNeg+"_x_trackdxdz_fitline")
02353 chamber_y_trackdxdz_fit2 = tfile.Get(pdirNeg+"_y_trackdxdz_fitline")
02354 chamber_dxdz_trackdxdz_fit2 = tfile.Get(pdirNeg+"_dxdz_trackdxdz_fitline")
02355 chamber_dydz_trackdxdz_fit2 = tfile.Get(pdirNeg+"_dydz_trackdxdz_fitline")
02356
02357 chamber_x_trackdydz = tfile.Get(pdirPos+"_x_trackdydz")
02358 chamber_x_trackdydz_fit = tfile.Get(pdirPos+"_x_trackdydz_fitline")
02359 chamber_y_trackdydz = tfile.Get(pdirPos+"_y_trackdydz")
02360 chamber_y_trackdydz_fit = tfile.Get(pdirPos+"_y_trackdydz_fitline")
02361 chamber_dxdz_trackdydz = tfile.Get(pdirPos+"_dxdz_trackdydz")
02362 chamber_dxdz_trackdydz_fit = tfile.Get(pdirPos+"_dxdz_trackdydz_fitline")
02363 chamber_dydz_trackdydz = tfile.Get(pdirPos+"_dydz_trackdydz")
02364 chamber_dydz_trackdydz_fit = tfile.Get(pdirPos+"_dydz_trackdydz_fitline")
02365 chamber_x_trackdydz_fit2 = tfile.Get(pdirNeg+"_x_trackdydz_fitline")
02366 chamber_y_trackdydz_fit2 = tfile.Get(pdirNeg+"_y_trackdydz_fitline")
02367 chamber_dxdz_trackdydz_fit2 = tfile.Get(pdirNeg+"_dxdz_trackdydz_fitline")
02368 chamber_dydz_trackdydz_fit2 = tfile.Get(pdirNeg+"_dydz_trackdydz_fitline")
02369
02370 if not chamber_x_trackx:
02371 chamber_x_trackx = tfile.Get(pdirPos+"_residual_trackx")
02372 chamber_x_trackx_fit = tfile.Get(pdirPos+"_residual_trackx_fitline")
02373 chamber_dxdz_trackx = tfile.Get(pdirPos+"_resslope_trackx")
02374 chamber_dxdz_trackx_fit = tfile.Get(pdirPos+"_resslope_trackx_fitline")
02375 chamber_x_trackx_fit2 = tfile.Get(pdirNeg+"_residual_trackx_fitline")
02376 chamber_dxdz_trackx_fit2 = tfile.Get(pdirNeg+"_resslope_trackx_fitline")
02377
02378 chamber_x_tracky = tfile.Get(pdirPos+"_residual_tracky")
02379 chamber_x_tracky_fit = tfile.Get(pdirPos+"_residual_tracky_fitline")
02380 chamber_dxdz_tracky = tfile.Get(pdirPos+"_resslope_tracky")
02381 chamber_dxdz_tracky_fit = tfile.Get(pdirPos+"_resslope_tracky_fitline")
02382 chamber_x_tracky_fit2 = tfile.Get(pdirNeg+"_residual_tracky_fitline")
02383 chamber_dxdz_tracky_fit2 = tfile.Get(pdirNeg+"_resslope_tracky_fitline")
02384
02385 chamber_x_trackdxdz = tfile.Get(pdirPos+"_residual_trackdxdz")
02386 chamber_x_trackdxdz_fit = tfile.Get(pdirPos+"_residual_trackdxdz_fitline")
02387 chamber_dxdz_trackdxdz = tfile.Get(pdirPos+"_resslope_trackdxdz")
02388 chamber_dxdz_trackdxdz_fit = tfile.Get(pdirPos+"_resslope_trackdxdz_fitline")
02389 chamber_x_trackdxdz_fit2 = tfile.Get(pdirNeg+"_residual_trackdxdz_fitline")
02390 chamber_dxdz_trackdxdz_fit2 = tfile.Get(pdirNeg+"_resslope_trackdxdz_fitline")
02391
02392 chamber_x_trackdydz = tfile.Get(pdirPos+"_residual_trackdydz")
02393 chamber_x_trackdydz_fit = tfile.Get(pdirPos+"_residual_trackdydz_fitline")
02394 chamber_dxdz_trackdydz = tfile.Get(pdirPos+"_resslope_trackdydz")
02395 chamber_dxdz_trackdydz_fit = tfile.Get(pdirPos+"_resslope_trackdydz_fitline")
02396 chamber_x_trackdydz_fit2 = tfile.Get(pdirNeg+"_residual_trackdydz_fitline")
02397 chamber_dxdz_trackdydz_fit2 = tfile.Get(pdirNeg+"_resslope_trackdydz_fitline")
02398
02399 if not chamber_x_trackx:
02400 print "Can't find neither "+pdirPos+"_residual nor "+pdirPos+"_residual_trackx"
02401 return
02402
02403 chamber_x_trackx = chamber_x_trackx.Clone()
02404 chamber_dxdz_trackx = chamber_dxdz_trackx.Clone()
02405 chamber_x_tracky = chamber_x_tracky.Clone()
02406 chamber_dxdz_tracky = chamber_dxdz_tracky.Clone()
02407 chamber_x_trackdxdz = chamber_x_trackdxdz.Clone()
02408 chamber_dxdz_trackdxdz = chamber_dxdz_trackdxdz.Clone()
02409 chamber_x_trackdydz = chamber_x_trackdydz.Clone()
02410 chamber_dxdz_trackdydz = chamber_dxdz_trackdydz.Clone()
02411
02412 if not not chamber_y_trackx:
02413 chamber_y_trackx = chamber_y_trackx.Clone()
02414 chamber_dydz_trackx = chamber_dydz_trackx.Clone()
02415 chamber_y_tracky = chamber_y_tracky.Clone()
02416 chamber_dydz_tracky = chamber_dydz_tracky.Clone()
02417 chamber_y_trackdxdz = chamber_y_trackdxdz.Clone()
02418 chamber_dydz_trackdxdz = chamber_dydz_trackdxdz.Clone()
02419 chamber_y_trackdydz = chamber_y_trackdydz.Clone()
02420 chamber_dydz_trackdydz = chamber_dydz_trackdydz.Clone()
02421
02422 if not not chamber_y_trackx:
02423 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_x_trackx")); chamber_x_trackx.Merge(tlist)
02424 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_dxdz_trackx")); chamber_dxdz_trackx.Merge(tlist)
02425 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_x_tracky")); chamber_x_tracky.Merge(tlist)
02426 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_dxdz_tracky")); chamber_dxdz_tracky.Merge(tlist)
02427 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_x_trackdxdz")); chamber_x_trackdxdz.Merge(tlist)
02428 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_dxdz_trackdxdz")); chamber_dxdz_trackdxdz.Merge(tlist)
02429 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_x_trackdydz")); chamber_x_trackdydz.Merge(tlist)
02430 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_dxdz_trackdydz")); chamber_dxdz_trackdydz.Merge(tlist)
02431 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_y_trackx")); chamber_y_trackx.Merge(tlist)
02432 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_dydz_trackx")); chamber_dydz_trackx.Merge(tlist)
02433 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_y_tracky")); chamber_y_tracky.Merge(tlist)
02434 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_dydz_tracky")); chamber_dydz_tracky.Merge(tlist)
02435 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_y_trackdxdz")); chamber_y_trackdxdz.Merge(tlist)
02436 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_dydz_trackdxdz")); chamber_dydz_trackdxdz.Merge(tlist)
02437 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_y_trackdydz")); chamber_y_trackdydz.Merge(tlist)
02438 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_dydz_trackdydz")); chamber_dydz_trackdydz.Merge(tlist)
02439 else:
02440 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_residual_trackx")); chamber_x_trackx.Merge(tlist)
02441 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_resslope_trackx")); chamber_dxdz_trackx.Merge(tlist)
02442 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_residual_tracky")); chamber_x_tracky.Merge(tlist)
02443 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_resslope_tracky")); chamber_dxdz_tracky.Merge(tlist)
02444 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_residual_trackdxdz")); chamber_x_trackdxdz.Merge(tlist)
02445 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_resslope_trackdxdz")); chamber_dxdz_trackdxdz.Merge(tlist)
02446 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_residual_trackdydz")); chamber_x_trackdydz.Merge(tlist)
02447 tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_resslope_trackdydz")); chamber_dxdz_trackdydz.Merge(tlist)
02448
02449 rr1=10.
02450 rr2=10.
02451 chamber_x_trackx.SetAxisRange(-rr1, rr1, "Y")
02452 chamber_dxdz_trackx.SetAxisRange(-rr2, rr2, "Y")
02453 chamber_x_tracky.SetAxisRange(-rr1, rr1, "Y")
02454 chamber_dxdz_tracky.SetAxisRange(-rr2, rr2, "Y")
02455 chamber_x_trackdxdz.SetAxisRange(-rr1, rr1, "Y")
02456 chamber_dxdz_trackdxdz.SetAxisRange(-rr2, rr2, "Y")
02457 chamber_x_trackdydz.SetAxisRange(-rr1, rr1, "Y")
02458 chamber_dxdz_trackdydz.SetAxisRange(-rr2, rr2, "Y")
02459
02460 rr3=10.
02461 if not not chamber_y_trackx:
02462 chamber_y_trackx.SetAxisRange(-rr3, rr3, "Y")
02463 chamber_dydz_trackx.SetAxisRange(-rr3, rr3, "Y")
02464 chamber_y_tracky.SetAxisRange(-rr3, rr3, "Y")
02465 chamber_dydz_tracky.SetAxisRange(-rr3, rr3, "Y")
02466 chamber_y_trackdxdz.SetAxisRange(-rr3, rr3, "Y")
02467 chamber_dydz_trackdxdz.SetAxisRange(-rr3, rr3, "Y")
02468 chamber_y_trackdydz.SetAxisRange(-rr3, rr3, "Y")
02469 chamber_dydz_trackdydz.SetAxisRange(-rr3, rr3, "Y")
02470
02471 for h in chamber_x_trackx, chamber_y_trackx, chamber_dxdz_trackx, chamber_dydz_trackx, \
02472 chamber_x_tracky, chamber_y_tracky, chamber_dxdz_tracky, chamber_dydz_tracky, \
02473 chamber_x_trackdxdz, chamber_y_trackdxdz, chamber_dxdz_trackdxdz, chamber_dydz_trackdxdz, \
02474 chamber_x_trackdydz, chamber_y_trackdydz, chamber_dxdz_trackdydz, chamber_dydz_trackdydz:
02475 if not not h:
02476 h.SetMarkerStyle(20)
02477 h.SetMarkerSize(0.5)
02478 h.GetXaxis().SetLabelSize(0.12)
02479 h.GetYaxis().SetLabelSize(0.12)
02480 h.GetXaxis().SetNdivisions(505)
02481 h.GetYaxis().SetNdivisions(505)
02482 h.GetXaxis().SetLabelOffset(0.03)
02483 h.GetYaxis().SetLabelOffset(0.03)
02484
02485 trackdxdz_minimum, trackdxdz_maximum = None, None
02486 for h in chamber_x_trackdxdz, chamber_y_trackdxdz, chamber_dxdz_trackdxdz, chamber_dydz_trackdxdz:
02487 if not not h:
02488 for i in xrange(1, h.GetNbinsX()+1):
02489 if h.GetBinError(i) > 0.01 and h.GetBinContent(i) - h.GetBinError(i) < 10. and \
02490 h.GetBinContent(i) + h.GetBinError(i) > -10.:
02491 if not trackdxdz_minimum or trackdxdz_minimum > h.GetBinCenter(i):
02492 trackdxdz_minimum = h.GetBinCenter(i)
02493 if trackdxdz_maximum < h.GetBinCenter(i):
02494 trackdxdz_maximum = h.GetBinCenter(i)
02495 if not not trackdxdz_minimum and not not trackdxdz_maximum:
02496 for h in chamber_x_trackdxdz, chamber_y_trackdxdz, chamber_dxdz_trackdxdz, chamber_dydz_trackdxdz:
02497 if not not h:
02498 h.SetAxisRange(trackdxdz_minimum, trackdxdz_maximum, "X")
02499
02500 trackdydz_minimum, trackdydz_maximum = None, None
02501 for h in chamber_x_trackdydz, chamber_y_trackdydz, chamber_dxdz_trackdydz, chamber_dydz_trackdydz:
02502 if not not h:
02503 for i in xrange(1, h.GetNbinsX()+1):
02504 if h.GetBinError(i) > 0.01 and h.GetBinContent(i) - h.GetBinError(i) < 10. and \
02505 h.GetBinContent(i) + h.GetBinError(i) > -10.:
02506 if not trackdydz_minimum or trackdydz_minimum > h.GetBinCenter(i):
02507 trackdydz_minimum = h.GetBinCenter(i)
02508 if trackdydz_maximum < h.GetBinCenter(i):
02509 trackdydz_maximum = h.GetBinCenter(i)
02510 if not not trackdydz_minimum and not not trackdydz_maximum:
02511 for h in chamber_x_trackdydz, chamber_y_trackdydz, chamber_dxdz_trackdydz, chamber_dydz_trackdydz:
02512 if not not h:
02513 h.SetAxisRange(trackdydz_minimum, trackdydz_maximum, "X")
02514
02515 for f in chamber_x_trackx_fit2, chamber_y_trackx_fit2, chamber_dxdz_trackx_fit2, chamber_dydz_trackx_fit2, \
02516 chamber_x_tracky_fit2, chamber_y_tracky_fit2, chamber_dxdz_tracky_fit2, chamber_dydz_tracky_fit2, \
02517 chamber_x_trackdxdz_fit2, chamber_y_trackdxdz_fit2, chamber_dxdz_trackdxdz_fit2, chamber_dydz_trackdxdz_fit2, \
02518 chamber_x_trackdydz_fit2, chamber_y_trackdydz_fit2, chamber_dxdz_trackdydz_fit2, chamber_dydz_trackdydz_fit2:
02519 if not not f:
02520 f.SetLineColor(4)
02521
02522 if not not chamber_y_trackx:
02523 c1.Clear()
02524
02525 pads = [None]
02526 pads.append(ROOT.TPad("p1" ,"",0.00,0.78,0.07,1.00,0,0,0))
02527 pads.append(ROOT.TPad("p2" ,"",0.07,0.78,0.34,1.00,0,0,0))
02528 pads.append(ROOT.TPad("p3" ,"",0.34,0.78,0.56,1.00,0,0,0))
02529 pads.append(ROOT.TPad("p4" ,"",0.56,0.78,0.78,1.00,0,0,0))
02530 pads.append(ROOT.TPad("p5" ,"",0.78,0.78,1.00,1.00,0,0,0))
02531 pads.append(ROOT.TPad("p6" ,"",0.00,0.56,0.07,0.78,0,0,0))
02532 pads.append(ROOT.TPad("p7" ,"",0.07,0.56,0.34,0.78,0,0,0))
02533 pads.append(ROOT.TPad("p8" ,"",0.34,0.56,0.56,0.78,0,0,0))
02534 pads.append(ROOT.TPad("p9" ,"",0.56,0.56,0.78,0.78,0,0,0))
02535 pads.append(ROOT.TPad("p10","",0.78,0.56,1.00,0.78,0,0,0))
02536 pads.append(ROOT.TPad("p11","",0.00,0.34,0.07,0.56,0,0,0))
02537 pads.append(ROOT.TPad("p12","",0.07,0.34,0.34,0.56,0,0,0))
02538 pads.append(ROOT.TPad("p13","",0.34,0.34,0.56,0.56,0,0,0))
02539 pads.append(ROOT.TPad("p14","",0.56,0.34,0.78,0.56,0,0,0))
02540 pads.append(ROOT.TPad("p15","",0.78,0.34,1.00,0.56,0,0,0))
02541 pads.append(ROOT.TPad("p16","",0.00,0.07,0.07,0.34,0,0,0))
02542 pads.append(ROOT.TPad("p17","",0.07,0.07,0.34,0.34,0,0,0))
02543 pads.append(ROOT.TPad("p18","",0.34,0.07,0.56,0.34,0,0,0))
02544 pads.append(ROOT.TPad("p19","",0.56,0.07,0.78,0.34,0,0,0))
02545 pads.append(ROOT.TPad("p20","",0.78,0.07,1.00,0.34,0,0,0))
02546 pads.append(ROOT.TPad("p21","",0.00,0.00,0.07,0.07,0,0,0))
02547 pads.append(ROOT.TPad("p22","",0.07,0.00,0.34,0.07,0,0,0))
02548 pads.append(ROOT.TPad("p23","",0.34,0.00,0.56,0.07,0,0,0))
02549 pads.append(ROOT.TPad("p24","",0.56,0.00,0.78,0.07,0,0,0))
02550 pads.append(ROOT.TPad("p25","",0.78,0.00,1.00,0.07,0,0,0))
02551 for p in pads:
02552 if not not p:
02553 p.Draw()
02554 ROOT.SetOwnership(p,False)
02555
02556 label1 = ROOT.TPaveLabel(0, 0, 1, 1, "x residuals (mm)","")
02557 label2 = ROOT.TPaveLabel(0, 0, 1, 1, "y residuals (mm)","")
02558 label3 = ROOT.TPaveLabel(0, 0, 1, 1, "dx/dz residuals (mrad)","")
02559 label4 = ROOT.TPaveLabel(0, 0, 1, 1, "dy/dz residuals (mrad)","")
02560 label5 = ROOT.TPaveLabel(0, 0, 1, 1, "x position (cm)","")
02561 label6 = ROOT.TPaveLabel(0, 0, 1, 1, "y position (cm)","")
02562 label7 = ROOT.TPaveLabel(0, 0, 1, 1, "dx/dz angle (rad)","")
02563 label8 = ROOT.TPaveLabel(0, 0, 1, 1, "dy/dz angle (rad)","")
02564 label9 = ROOT.TPaveLabel(0, 0.85, 1, 1, getname(r),"NDC")
02565
02566 for l in label1, label2, label3, label4, label5, label6, label7, label8, label9:
02567 l.SetBorderSize(0)
02568 l.SetFillColor(ROOT.kWhite)
02569
02570 for l in label1, label2, label3, label4:
02571 l.SetTextAngle(90)
02572 l.SetTextSize(0.09)
02573
02574
02575 label9.SetTextSize(0.59)
02576
02577 pads[1].cd(); label1.Draw()
02578 pads[6].cd(); label2.Draw()
02579 pads[11].cd(); label3.Draw()
02580 pads[16].cd(); label4.Draw()
02581 pads[22].cd(); label5.Draw()
02582 pads[23].cd(); label6.Draw()
02583 pads[24].cd(); label7.Draw()
02584 pads[25].cd(); label8.Draw()
02585
02586 pads[2].SetRightMargin(1e-5)
02587 pads[2].SetBottomMargin(1e-5)
02588 pads[2].SetLeftMargin(0.17)
02589 pads[3].SetLeftMargin(1e-5)
02590 pads[3].SetRightMargin(1e-5)
02591 pads[3].SetBottomMargin(1e-5)
02592 pads[4].SetLeftMargin(1e-5)
02593 pads[4].SetRightMargin(1e-5)
02594 pads[4].SetBottomMargin(1e-5)
02595 pads[5].SetLeftMargin(1e-5)
02596 pads[5].SetBottomMargin(1e-5)
02597
02598 pads[7].SetRightMargin(1e-5)
02599 pads[7].SetBottomMargin(1e-5)
02600 pads[7].SetTopMargin(1e-5)
02601 pads[7].SetLeftMargin(0.17)
02602 pads[8].SetLeftMargin(1e-5)
02603 pads[8].SetRightMargin(1e-5)
02604 pads[8].SetBottomMargin(1e-5)
02605 pads[8].SetTopMargin(1e-5)
02606 pads[9].SetLeftMargin(1e-5)
02607 pads[9].SetRightMargin(1e-5)
02608 pads[9].SetBottomMargin(1e-5)
02609 pads[9].SetTopMargin(1e-5)
02610 pads[10].SetLeftMargin(1e-5)
02611 pads[10].SetBottomMargin(1e-5)
02612 pads[10].SetTopMargin(1e-5)
02613
02614 pads[12].SetRightMargin(1e-5)
02615 pads[12].SetBottomMargin(1e-5)
02616 pads[12].SetTopMargin(1e-5)
02617 pads[12].SetLeftMargin(0.17)
02618 pads[13].SetLeftMargin(1e-5)
02619 pads[13].SetRightMargin(1e-5)
02620 pads[13].SetBottomMargin(1e-5)
02621 pads[13].SetTopMargin(1e-5)
02622 pads[14].SetLeftMargin(1e-5)
02623 pads[14].SetRightMargin(1e-5)
02624 pads[14].SetBottomMargin(1e-5)
02625 pads[14].SetTopMargin(1e-5)
02626 pads[15].SetLeftMargin(1e-5)
02627 pads[15].SetBottomMargin(1e-5)
02628 pads[15].SetTopMargin(1e-5)
02629
02630 pads[17].SetRightMargin(1e-5)
02631 pads[17].SetTopMargin(1e-5)
02632 pads[17].SetLeftMargin(0.17)
02633 pads[18].SetLeftMargin(1e-5)
02634 pads[18].SetRightMargin(1e-5)
02635 pads[18].SetTopMargin(1e-5)
02636 pads[19].SetLeftMargin(1e-5)
02637 pads[19].SetRightMargin(1e-5)
02638 pads[19].SetTopMargin(1e-5)
02639 pads[20].SetLeftMargin(1e-5)
02640 pads[20].SetTopMargin(1e-5)
02641
02642 chamber_x_trackx.GetXaxis().SetLabelColor(ROOT.kWhite)
02643 chamber_x_tracky.GetXaxis().SetLabelColor(ROOT.kWhite)
02644 chamber_x_tracky.GetYaxis().SetLabelColor(ROOT.kWhite)
02645 chamber_x_trackdxdz.GetXaxis().SetLabelColor(ROOT.kWhite)
02646 chamber_x_trackdxdz.GetYaxis().SetLabelColor(ROOT.kWhite)
02647 chamber_x_trackdydz.GetXaxis().SetLabelColor(ROOT.kWhite)
02648 chamber_x_trackdydz.GetYaxis().SetLabelColor(ROOT.kWhite)
02649 chamber_y_trackx.GetXaxis().SetLabelColor(ROOT.kWhite)
02650 chamber_y_tracky.GetXaxis().SetLabelColor(ROOT.kWhite)
02651 chamber_y_tracky.GetYaxis().SetLabelColor(ROOT.kWhite)
02652 chamber_y_trackdxdz.GetXaxis().SetLabelColor(ROOT.kWhite)
02653 chamber_y_trackdxdz.GetYaxis().SetLabelColor(ROOT.kWhite)
02654 chamber_y_trackdydz.GetXaxis().SetLabelColor(ROOT.kWhite)
02655 chamber_y_trackdydz.GetYaxis().SetLabelColor(ROOT.kWhite)
02656 chamber_dxdz_trackx.GetXaxis().SetLabelColor(ROOT.kWhite)
02657 chamber_dxdz_tracky.GetXaxis().SetLabelColor(ROOT.kWhite)
02658 chamber_dxdz_tracky.GetYaxis().SetLabelColor(ROOT.kWhite)
02659 chamber_dxdz_trackdxdz.GetXaxis().SetLabelColor(ROOT.kWhite)
02660 chamber_dxdz_trackdxdz.GetYaxis().SetLabelColor(ROOT.kWhite)
02661 chamber_dxdz_trackdydz.GetXaxis().SetLabelColor(ROOT.kWhite)
02662 chamber_dxdz_trackdydz.GetYaxis().SetLabelColor(ROOT.kWhite)
02663
02664
02665 chamber_dydz_tracky.GetYaxis().SetLabelColor(ROOT.kWhite)
02666 chamber_dydz_trackdxdz.GetYaxis().SetLabelColor(ROOT.kWhite)
02667 chamber_dydz_trackdydz.GetYaxis().SetLabelColor(ROOT.kWhite)
02668
02669 pads[2].cd()
02670 chamber_x_trackx.Draw("e1")
02671 if not suppressblue: chamber_x_trackx_fit2.Draw("samel")
02672 chamber_x_trackx_fit.Draw("samel")
02673
02674 print getname(r)
02675
02676 label9.Draw()
02677
02678
02679 pads[3].cd()
02680 chamber_x_tracky.Draw("e1")
02681 if not suppressblue: chamber_x_tracky_fit2.Draw("samel")
02682 chamber_x_tracky_fit.Draw("samel")
02683
02684 pads[4].cd()
02685 chamber_x_trackdxdz.Draw("e1")
02686 if not suppressblue: chamber_x_trackdxdz_fit2.Draw("samel")
02687 chamber_x_trackdxdz_fit.Draw("samel")
02688
02689 pads[5].cd()
02690 chamber_x_trackdydz.Draw("e1")
02691 if not suppressblue: chamber_x_trackdydz_fit2.Draw("samel")
02692 chamber_x_trackdydz_fit.Draw("samel")
02693
02694 pads[7].cd()
02695 chamber_y_trackx.Draw("e1")
02696 if not suppressblue: chamber_y_trackx_fit2.Draw("samel")
02697 chamber_y_trackx_fit.Draw("samel")
02698
02699 pads[8].cd()
02700 chamber_y_tracky.Draw("e1")
02701 if not suppressblue: chamber_y_tracky_fit2.Draw("samel")
02702 chamber_y_tracky_fit.Draw("samel")
02703
02704 pads[9].cd()
02705 chamber_y_trackdxdz.Draw("e1")
02706 if not suppressblue: chamber_y_trackdxdz_fit2.Draw("samel")
02707 chamber_y_trackdxdz_fit.Draw("samel")
02708
02709 pads[10].cd()
02710 chamber_y_trackdydz.Draw("e1")
02711 if not suppressblue: chamber_y_trackdydz_fit2.Draw("samel")
02712 chamber_y_trackdydz_fit.Draw("samel")
02713
02714 pads[12].cd()
02715 chamber_dxdz_trackx.Draw("e1")
02716 if not suppressblue: chamber_dxdz_trackx_fit2.Draw("samel")
02717 chamber_dxdz_trackx_fit.Draw("samel")
02718
02719 pads[13].cd()
02720 chamber_dxdz_tracky.Draw("e1")
02721 if not suppressblue: chamber_dxdz_tracky_fit2.Draw("samel")
02722 chamber_dxdz_tracky_fit.Draw("samel")
02723
02724 pads[14].cd()
02725 chamber_dxdz_trackdxdz.Draw("e1")
02726 if not suppressblue: chamber_dxdz_trackdxdz_fit2.Draw("samel")
02727 chamber_dxdz_trackdxdz_fit.Draw("samel")
02728
02729 pads[15].cd()
02730 chamber_dxdz_trackdydz.Draw("e1")
02731 if not suppressblue: chamber_dxdz_trackdydz_fit2.Draw("samel")
02732 chamber_dxdz_trackdydz_fit.Draw("samel")
02733
02734 pads[17].cd()
02735 chamber_dydz_trackx.Draw("e1")
02736 if not suppressblue: chamber_dydz_trackx_fit2.Draw("samel")
02737 chamber_dydz_trackx_fit.Draw("samel")
02738
02739 pads[18].cd()
02740 chamber_dydz_tracky.Draw("e1")
02741 if not suppressblue: chamber_dydz_tracky_fit2.Draw("samel")
02742 chamber_dydz_tracky_fit.Draw("samel")
02743
02744 pads[19].cd()
02745 chamber_dydz_trackdxdz.Draw("e1")
02746 if not suppressblue: chamber_dydz_trackdxdz_fit2.Draw("samel")
02747 chamber_dydz_trackdxdz_fit.Draw("samel")
02748
02749 pads[20].cd()
02750 chamber_dydz_trackdydz.Draw("e1")
02751 if not suppressblue: chamber_dydz_trackdydz_fit2.Draw("samel")
02752 chamber_dydz_trackdydz_fit.Draw("samel")
02753
02754 else:
02755 c1.Clear()
02756
02757 pads = [None]
02758 pads.append(ROOT.TPad("p1" ,"",0.00,0.55,0.07,1.00,0,0,0))
02759 pads.append(ROOT.TPad("p2" ,"",0.07,0.55,0.34,1.00,0,0,0))
02760 pads.append(ROOT.TPad("p3" ,"",0.34,0.55,0.56,1.00,0,0,0))
02761 pads.append(ROOT.TPad("p4" ,"",0.56,0.55,0.78,1.00,0,0,0))
02762 pads.append(ROOT.TPad("p5" ,"",0.78,0.55,1.00,1.00,0,0,0))
02763 pads.append(ROOT.TPad("p6" ,"",0.00,0.1,0.07,0.55,0,0,0))
02764 pads.append(ROOT.TPad("p7" ,"",0.07,0.1,0.34,0.55,0,0,0))
02765 pads.append(ROOT.TPad("p8" ,"",0.34,0.1,0.56,0.55,0,0,0))
02766 pads.append(ROOT.TPad("p9" ,"",0.56,0.1,0.78,0.55,0,0,0))
02767 pads.append(ROOT.TPad("p10","",0.78,0.1,1.00,0.55,0,0,0))
02768 pads.append(ROOT.TPad("p11","",0.00,0.,0.07,0.1,0,0,0))
02769 pads.append(ROOT.TPad("p12","",0.07,0.,0.34,0.1,0,0,0))
02770 pads.append(ROOT.TPad("p13","",0.34,0.,0.56,0.1,0,0,0))
02771 pads.append(ROOT.TPad("p14","",0.56,0.,0.78,0.1,0,0,0))
02772 pads.append(ROOT.TPad("p15","",0.78,0.,1.00,0.1,0,0,0))
02773 for p in pads:
02774 if not not p:
02775 p.Draw()
02776 ROOT.SetOwnership(p,False)
02777
02778 label1 = ROOT.TPaveLabel(0, 0, 1, 1, "x residuals (mm)")
02779 label2 = ROOT.TPaveLabel(0, 0, 1, 1, "dx/dz residuals (mrad)")
02780 label3 = ROOT.TPaveLabel(0, 0.3, 1, 1, "x position (cm)")
02781 label4 = ROOT.TPaveLabel(0, 0.3, 1, 1, "y position (cm)")
02782 label5 = ROOT.TPaveLabel(0, 0.3, 1, 1, "dx/dz angle (rad)")
02783 label6 = ROOT.TPaveLabel(0, 0.3, 1, 1, "dy/dz angle (rad)")
02784 label9 = ROOT.TPaveLabel(0, 0.85, 1, 1, getname(r),"NDC")
02785
02786 if name[0:2] == "ME":
02787 label1 = ROOT.TPaveLabel(0, 0, 1, 1, "r#phi residuals (mm)")
02788 label2 = ROOT.TPaveLabel(0, 0, 1, 1, "d(r#phi)/dz residuals (mrad)")
02789
02790 for l in label1, label2, label3, label4, label5, label6, label9:
02791 l.SetBorderSize(0)
02792 l.SetFillColor(ROOT.kWhite)
02793
02794 for l in label1, label2:
02795 l.SetTextAngle(90)
02796 l.SetTextSize(0.09)
02797
02798
02799 label9.SetTextSize(0.29)
02800
02801 pads[1].cd(); label1.Draw()
02802 pads[6].cd(); label2.Draw()
02803 pads[12].cd(); label3.Draw()
02804 pads[13].cd(); label4.Draw()
02805 pads[14].cd(); label5.Draw()
02806 pads[15].cd(); label6.Draw()
02807
02808
02809 pads[2].SetRightMargin(1e-5)
02810 pads[2].SetBottomMargin(1e-5)
02811 pads[3].SetLeftMargin(1e-5)
02812 pads[3].SetRightMargin(1e-5)
02813 pads[3].SetBottomMargin(1e-5)
02814 pads[4].SetLeftMargin(1e-5)
02815 pads[4].SetRightMargin(1e-5)
02816 pads[4].SetBottomMargin(1e-5)
02817 pads[5].SetLeftMargin(1e-5)
02818 pads[5].SetBottomMargin(1e-5)
02819
02820 pads[7].SetRightMargin(1e-5)
02821 pads[7].SetTopMargin(1e-5)
02822 pads[8].SetLeftMargin(1e-5)
02823 pads[8].SetRightMargin(1e-5)
02824 pads[8].SetTopMargin(1e-5)
02825 pads[9].SetLeftMargin(1e-5)
02826 pads[9].SetRightMargin(1e-5)
02827 pads[9].SetTopMargin(1e-5)
02828 pads[10].SetLeftMargin(1e-5)
02829 pads[10].SetTopMargin(1e-5)
02830
02831 chamber_x_trackx.GetXaxis().SetLabelColor(ROOT.kWhite)
02832 chamber_x_tracky.GetXaxis().SetLabelColor(ROOT.kWhite)
02833 chamber_x_tracky.GetYaxis().SetLabelColor(ROOT.kWhite)
02834 chamber_x_trackdxdz.GetXaxis().SetLabelColor(ROOT.kWhite)
02835 chamber_x_trackdxdz.GetYaxis().SetLabelColor(ROOT.kWhite)
02836 chamber_x_trackdydz.GetXaxis().SetLabelColor(ROOT.kWhite)
02837 chamber_x_trackdydz.GetYaxis().SetLabelColor(ROOT.kWhite)
02838
02839 chamber_dxdz_tracky.GetYaxis().SetLabelColor(ROOT.kWhite)
02840 chamber_dxdz_trackdxdz.GetYaxis().SetLabelColor(ROOT.kWhite)
02841 chamber_dxdz_trackdydz.GetYaxis().SetLabelColor(ROOT.kWhite)
02842
02843 pads[2].cd()
02844 chamber_x_trackx.Draw("e1")
02845 if not suppressblue: chamber_x_trackx_fit2.Draw("samel")
02846 chamber_x_trackx_fit.Draw("samel")
02847 label9.Draw()
02848
02849 pads[3].cd()
02850 chamber_x_tracky.Draw("e1")
02851 if not suppressblue: chamber_x_tracky_fit2.Draw("samel")
02852 chamber_x_tracky_fit.Draw("samel")
02853
02854 pads[4].cd()
02855 chamber_x_trackdxdz.Draw("e1")
02856 if not suppressblue: chamber_x_trackdxdz_fit2.Draw("samel")
02857 chamber_x_trackdxdz_fit.Draw("samel")
02858
02859 pads[5].cd()
02860 chamber_x_trackdydz.Draw("e1")
02861 if not suppressblue: chamber_x_trackdydz_fit2.Draw("samel")
02862 chamber_x_trackdydz_fit.Draw("samel")
02863
02864 pads[7].cd()
02865 chamber_dxdz_trackx.Draw("e1")
02866 if not suppressblue: chamber_dxdz_trackx_fit2.Draw("samel")
02867 chamber_dxdz_trackx_fit.Draw("samel")
02868
02869 pads[8].cd()
02870 chamber_dxdz_tracky.Draw("e1")
02871 if not suppressblue: chamber_dxdz_tracky_fit2.Draw("samel")
02872 chamber_dxdz_tracky_fit.Draw("samel")
02873
02874 pads[9].cd()
02875 chamber_dxdz_trackdxdz.Draw("e1")
02876 if not suppressblue: chamber_dxdz_trackdxdz_fit2.Draw("samel")
02877 chamber_dxdz_trackdxdz_fit.Draw("samel")
02878
02879 pads[10].cd()
02880 chamber_dxdz_trackdydz.Draw("e1")
02881 if not suppressblue: chamber_dxdz_trackdydz_fit2.Draw("samel")
02882 chamber_dxdz_trackdydz_fit.Draw("samel")
02883
02884 tn = time.time()
02885 ddt[8] = 1./ddt[7]*((ddt[7]-1)*ddt[8] + tn-t1)
02886
02887
02888
02889 def segdiff(tfiles, component, pair, **args):
02890 tdrStyle.SetOptFit(1)
02891 tdrStyle.SetOptTitle(1)
02892 tdrStyle.SetTitleBorderSize(1)
02893 tdrStyle.SetTitleFontSize(0.05)
02894 tdrStyle.SetStatW(0.2)
02895 tdrStyle.SetStatY(0.9)
02896 tdrStyle.SetStatFontSize(0.06)
02897
02898 if component[0:2] == "dt":
02899 wheel = args["wheel"]
02900 wheelletter = wheelLetter(wheel)
02901 sector = args["sector"]
02902 profname = "%s_%s_%02d_%s" % (component, wheelletter, sector, str(pair))
02903 posname = "pos" + profname
02904 negname = "neg" + profname
02905
02906
02907 station1 = int(str(pair)[0])
02908 station2 = int(str(pair)[1])
02909 phi1 = signConventions["DT", wheel, station1, sector][4]
02910 phi2 = signConventions["DT", wheel, station2, sector][4]
02911 if abs(phi1 - phi2) > 1.:
02912 if phi1 > phi2: phi1 -= 2.*pi
02913 else: phi1 += 2.*pi
02914 phi = (phi1 + phi2) / 2.
02915 while (phi < -pi): phi += 2.*pi
02916 while (phi > pi): phi -= 2.*pi
02917
02918 elif component[0:3] == "csc":
02919 endcap = args["endcap"]
02920 if endcap=="m":
02921 endcapnum=2
02922 endcapsign="-"
02923 elif endcap=="p":
02924 endcapnum=1
02925 endcapsign="+"
02926 else: raise Exception
02927
02928 ring = args["ring"]
02929 if ring>2 or ring<1: raise Exception
02930 station1 = int(str(pair)[0])
02931 station2 = int(str(pair)[1])
02932 if ring==1: ringname="inner"
02933 elif ring==2: ringname="outer"
02934 else: raise Exception
02935
02936 chamber = args["chamber"]
02937 if (ring==1 and chamber>18) or (ring==2 and chamber>36): raise Exception
02938
02939 profname = "csc%s_%s_%s_%02d_%s" % (ringname,component[4:], endcap, chamber, str(pair))
02940 posname = "pos" + profname
02941 negname = "neg" + profname
02942
02943
02944 station1 = int(str(pair)[0])
02945 station2 = int(str(pair)[1])
02946 phi1 = signConventions["CSC", endcapnum, station1, ring, chamber][4]
02947 phi2 = signConventions["CSC", endcapnum, station2, ring, chamber][4]
02948 if abs(phi1 - phi2) > 1.:
02949 if phi1 > phi2: phi1 -= 2.*pi
02950 else: phi1 += 2.*pi
02951 phi = (phi1 + phi2) / 2.
02952 while (phi < -pi*5./180.): phi += 2.*pi
02953 while (phi > pi*(2.-5./180.)): phi -= 2.*pi
02954
02955 else: raise Exception
02956
02957 if "window" in args: window = args["window"]
02958 else: window = 5.
02959
02960 global tmpprof, tmppos, tmpneg
02961 pdir = "AlignmentMonitorSegmentDifferences/iter1/"
02962 tmpprof = tfiles[0].Get(pdir + profname).Clone()
02963 tmpprof.SetMarkerStyle(8)
02964 tmppos = tfiles[0].Get(pdir + posname).Clone()
02965 tmpneg = tfiles[0].Get(pdir + negname).Clone()
02966 for tfile in tfiles[1:]:
02967 tmpprof.Add(tfile.Get(pdir + profname))
02968 tmppos.Add(tfile.Get(pdir + posname))
02969 tmpneg.Add(tfile.Get(pdir + negname))
02970
02971 for i in xrange(1, tmpprof.GetNbinsX()+1):
02972 if tmpprof.GetBinError(i) < 1e-5:
02973 tmpprof.SetBinError(i, 100.)
02974 tmpprof.SetAxisRange(-window, window, "Y")
02975
02976 f = ROOT.TF1("p1", "[0] + [1]*x", tmpprof.GetBinLowEdge(1), -tmpprof.GetBinLowEdge(1))
02977 f.SetParameters((tmppos.GetMean() + tmpneg.GetMean())/2., 0.)
02978
02979 tmpprof.SetXTitle("q/p_{T} (c/GeV)")
02980 if component == "dt13_resid":
02981 tmpprof.SetYTitle("#Deltax^{local} (mm)")
02982 tmppos.SetXTitle("#Deltax^{local} (mm)")
02983 tmpneg.SetXTitle("#Deltax^{local} (mm)")
02984 f.SetParNames("#Deltax^{local}_{0}", "Slope")
02985 if component == "dt13_slope":
02986 tmpprof.SetYTitle("#Deltadx/dz^{local} (mrad)")
02987 tmppos.SetXTitle("#Deltadx/dz^{local} (mrad)")
02988 tmpneg.SetXTitle("#Deltadx/dz^{local} (mrad)")
02989 f.SetParNames("#Deltadx/dz^{local}_{0}", "Slope")
02990 if component == "dt2_resid":
02991 tmpprof.SetYTitle("#Deltay^{local} (mm)")
02992 tmppos.SetXTitle("#Deltay^{local} (mm)")
02993 tmpneg.SetXTitle("#Deltay^{local} (mm)")
02994 f.SetParNames("#Deltay^{local}_{0}", "Slope")
02995 if component == "dt2_slope":
02996 tmpprof.SetYTitle("#Deltady/dz^{local} (mrad)")
02997 tmppos.SetXTitle("#Deltady/dz^{local} (mrad)")
02998 tmpneg.SetXTitle("#Deltady/dz^{local} (mrad)")
02999 f.SetParNames("#Deltady/dz^{local}_{0}", "Slope")
03000 if component == "csc_resid":
03001 tmpprof.SetXTitle("q/p_{z} (c/GeV)")
03002 tmpprof.SetYTitle("#Delta(r#phi)^{local} (mm)")
03003 tmppos.SetXTitle("#Delta(r#phi)^{local} (mm)")
03004 tmpneg.SetXTitle("#Delta(r#phi)^{local} (mm)")
03005 f.SetParNames("#Delta(r#phi)^{local}_{0}", "Slope")
03006 if component == "csc_slope":
03007 tmpprof.SetXTitle("q/p_{z} (c/GeV)")
03008 tmpprof.SetYTitle("#Deltad(r#phi)/dz^{local} (mrad)")
03009 tmppos.SetXTitle("#Deltad(r#phi)/dz^{local} (mrad)")
03010 tmpneg.SetXTitle("#Deltad(r#phi)/dz^{local} (mrad)")
03011 f.SetParNames("#Deltad(r#phi)/dz^{local}_{0}", "Slope")
03012
03013 tmpprof.GetXaxis().CenterTitle()
03014 tmpprof.GetYaxis().CenterTitle()
03015 tmppos.GetXaxis().CenterTitle()
03016 tmpneg.GetXaxis().CenterTitle()
03017 if component[0:2] == "dt":
03018 tmpprof.SetTitle("MB%d - MB%d, wheel %d, sector %02d" % (station1, station2, int(wheel), int(sector)))
03019 elif component[0:3] == "csc":
03020 tmpprof.SetTitle("ME%d - ME%d, for ME%s%d/%d/%d" % (station1, station2, endcapsign, station2, ring, chamber))
03021 else: raise Exception
03022
03023 tmppos.SetTitle("Positive muons")
03024 tmpneg.SetTitle("Negative muons")
03025
03026 c1.Clear()
03027 c1.Divide(2, 1)
03028 c1.GetPad(1).cd()
03029 fit1 = tmpprof.Fit("p1", "qS")
03030 tmpprof.Draw("e1")
03031 c1.GetPad(2).cd()
03032 c1.GetPad(2).Divide(1, 2)
03033 c1.GetPad(2).GetPad(1).cd()
03034 tmppos.Draw()
03035 f = ROOT.TF1("gausR", "[0]*exp(-(x - [1])**2 / 2. / [2]**2) / sqrt(2.*3.1415926) / [2]",
03036 tmppos.GetMean() - tmppos.GetRMS(), tmppos.GetMean() + tmppos.GetRMS())
03037 f.SetParameters(tmppos.GetEntries() * ((10. - -10.)/100.), tmppos.GetMean(), tmppos.GetRMS())
03038 f.SetParNames("Constant", "Mean", "Sigma")
03039 fit2 = tmppos.Fit("gausR", "qRS")
03040 c1.GetPad(2).GetPad(2).cd()
03041 tmpneg.Draw()
03042 f = ROOT.TF1("gausR", "[0]*exp(-(x - [1])**2 / 2. / [2]**2) / sqrt(2.*3.1415926) / [2]",
03043 tmpneg.GetMean() - tmpneg.GetRMS(), tmpneg.GetMean() + tmpneg.GetRMS())
03044 f.SetParameters(tmpneg.GetEntries() * ((10. - -10.)/100.), tmpneg.GetMean(), tmpneg.GetRMS())
03045 f.SetParNames("Constant", "Mean", "Sigma")
03046 fit3 = tmpneg.Fit("gausR", "qRS")
03047
03048 fit1ok = fit1.Status()==0 and fit1.CovMatrixStatus()==3
03049 fit2ok = fit2.Status()==0 and fit2.CovMatrixStatus()==3
03050 fit3ok = fit3.Status()==0 and fit3.CovMatrixStatus()==3
03051
03052 fitresult1 = None, None
03053 if fit1ok:
03054 fitresult1 = tmpprof.GetFunction("p1").GetParameter(0), tmpprof.GetFunction("p1").GetParError(0)
03055 fitresult2 = None, None
03056 if fit2ok and fit3ok:
03057 fitresult2 = (tmppos.GetFunction("gausR").GetParameter(1) + tmpneg.GetFunction("gausR").GetParameter(1)) / 2., \
03058 sqrt(tmppos.GetFunction("gausR").GetParError(1)**2 + tmpneg.GetFunction("gausR").GetParError(1)**2) / 2.
03059 return phi, fitresult1[0], fitresult1[1], fitresult2[0], fitresult2[1], fit1ok, fit2ok, fit3ok
03060
03061
03062
03063
03064
03065 def segdiff_xalign(tfiles, component, **args):
03066 tdrStyle.SetOptFit(1)
03067 tdrStyle.SetOptTitle(1)
03068 tdrStyle.SetTitleBorderSize(1)
03069 tdrStyle.SetTitleFontSize(0.05)
03070 tdrStyle.SetStatW(0.2)
03071 tdrStyle.SetStatY(0.9)
03072 tdrStyle.SetStatFontSize(0.06)
03073
03074 if component[0:4] == "x_dt":
03075 wheel = int(args["wheel"])
03076 if int(wheel)<0:
03077 wheell = "m%d" % abs(wheel)
03078 endcapsign="-"
03079 else:
03080 wheell = "p%d" % abs(wheel)
03081 endcapsign="+"
03082 station_dt = component[4]
03083 station_csc_1 = args["cscstations"][0]
03084 if station_csc_1=='1': ring_1 = 3
03085 else: ring_1 = 2
03086 sector = args["sector"]
03087 profname = "%s%s_W%sS%02d" % (component, station_csc_1, wheell, sector)
03088 posname_1 = "pos_" + profname
03089 negname_1 = "neg_" + profname
03090 if len(args["cscstations"]) > 1:
03091 station_csc_2 = args["cscstations"][1]
03092 if station_csc_2=='1': ring_2 = 3
03093 else: ring_2 = 2
03094 profname = "%s%s_W%sS%02d" % (component, station_csc_2, wheell, sector)
03095 posname_2 = "pos_" + profname
03096 negname_2 = "neg_" + profname
03097
03098 phi = signConventions["DT", wheel, int(station_dt), sector][4]
03099 while (phi < -pi): phi += 2.*pi
03100 while (phi > pi): phi -= 2.*pi
03101
03102 else: raise Exception
03103
03104 if "window" in args: window = args["window"]
03105 else: window = 5.
03106
03107 global tmppos, tmpneg, tmppos_2, tmpneg_2
03108 pdir = "AlignmentMonitorSegmentDifferences/iter1/"
03109 tmppos = tfiles[0].Get(pdir + posname_1).Clone()
03110 tmpneg = tfiles[0].Get(pdir + negname_1).Clone()
03111 if len(args["cscstations"]) > 1:
03112 tmppos_2 = tfiles[0].Get(pdir + posname_2).Clone()
03113 tmpneg_2 = tfiles[0].Get(pdir + negname_2).Clone()
03114 tmpneg.Rebin(2); tmppos.Rebin(2)
03115 for tfile in tfiles[1:]:
03116 tmppos.Add(tfile.Get(pdir + posname_1))
03117 tmpneg.Add(tfile.Get(pdir + negname_1))
03118 if len(args["cscstations"]) > 1:
03119 tmppos_2.Add(tfile.Get(pdir + posname_2))
03120 tmpneg_2.Add(tfile.Get(pdir + negname_2))
03121 tmpneg_2.Rebin(2); tmppos_2.Rebin(2)
03122
03123 result = {}
03124 result["fit_ok"] = False
03125 result["phi"] = phi
03126 ntot = tmppos.GetEntries() + tmpneg.GetEntries()
03127 if ntot == 0:
03128 return result
03129
03130 tmppos.SetXTitle("#Deltax^{loc}_{MB} - r_{DT}/r_{CSC}#times#Deltax^{loc}_{ME} (mm)")
03131 tmpneg.SetXTitle("#Deltax^{loc}_{MB} - r_{DT}/r_{CSC}#times#Deltax^{loc}_{ME} (mm)")
03132 title1 = "MB(W%+d St%s Sec%d) - ME%s%s/%d" % (wheel, station_dt, int(sector), endcapsign, station_csc_1, ring_1)
03133 tmppos.SetTitle("Positive #mu: %s" % title1);
03134 tmpneg.SetTitle("Negative #mu: %s" % title1);
03135 tmppos.GetXaxis().CenterTitle()
03136 tmpneg.GetXaxis().CenterTitle()
03137 if len(args["cscstations"]) > 1:
03138 tmppos.SetXTitle("#Deltax^{loc}_{DT} - r_{DT}/r_{CSC}#times#Deltax^{loc}_{CSC} (mm)")
03139 tmpneg.SetXTitle("#Deltax^{loc}_{DT} - r_{DT}/r_{CSC}#times#Deltax^{loc}_{CSC} (mm)")
03140 title2 = "MB(W%+d St%s Sec%d) - ME%s%s/%d" % (wheel, station_dt, int(sector), endcapsign, station_csc_2, ring_2)
03141 tmppos_2.SetTitle("Positive #mu: %s" % title2);
03142 tmpneg_2.SetTitle("Negative #mu: %s" % title2);
03143 tmppos_2.GetXaxis().CenterTitle()
03144 tmpneg_2.GetXaxis().CenterTitle()
03145
03146 c1.Clear()
03147 c1.Divide(2, 2)
03148
03149 c1.GetPad(1).cd()
03150 tmppos.Draw()
03151 fpos = ROOT.TF1("gpos", "gaus", tmppos.GetMean() - tmppos.GetRMS(), tmppos.GetMean() + tmppos.GetRMS())
03152 fpos.SetParameters(tmppos.GetEntries() * 2.5, tmppos.GetMean(), tmppos.GetRMS())
03153 fit_pos = tmppos.Fit("gpos", "qRS")
03154
03155 c1.GetPad(3).cd()
03156 tmpneg.Draw()
03157 fneg = ROOT.TF1("gneg", "gaus", tmpneg.GetMean() - tmpneg.GetRMS(), tmpneg.GetMean() + tmpneg.GetRMS())
03158 fneg.SetParameters(tmpneg.GetEntries() * 2.5, tmpneg.GetMean(), tmpneg.GetRMS())
03159 fit_neg = tmpneg.Fit("gneg", "qRS")
03160
03161 result["fit_ok"] = (fit_pos.Status()==0 and fit_pos.CovMatrixStatus()==3 and fit_neg.Status()==0 and fit_neg.CovMatrixStatus()==3)
03162 result["fit_peak"] = (fpos.GetParameter(1)*tmppos.GetEntries() + fneg.GetParameter(1)*tmpneg.GetEntries()) / ntot
03163 result["fit_peak_error"] = sqrt( (fpos.GetParError(1)*tmppos.GetEntries())**2 + (fneg.GetParError(1)*tmpneg.GetEntries())**2) / ntot
03164
03165 if len(args["cscstations"]) > 1:
03166 c1.GetPad(2).cd()
03167 tmppos_2.Draw()
03168 fpos_2 = ROOT.TF1("gpos2", "gaus", tmppos_2.GetMean() - tmppos_2.GetRMS(), tmppos_2.GetMean() + tmppos_2.GetRMS())
03169 fpos_2.SetParameters(tmppos_2.GetEntries() * 2.5, tmppos_2.GetMean(), tmppos_2.GetRMS())
03170 fit_pos_2 = tmppos_2.Fit("gpos2", "qRS")
03171
03172 c1.GetPad(4).cd()
03173 tmpneg_2.Draw()
03174 fneg_2 = ROOT.TF1("gneg2", "gaus", tmpneg_2.GetMean() - tmpneg_2.GetRMS(), tmpneg_2.GetMean() + tmpneg_2.GetRMS())
03175 fneg_2.SetParameters(tmpneg_2.GetEntries() * 2.5, tmpneg_2.GetMean(), tmpneg_2.GetRMS())
03176 fit_neg_2 = tmpneg_2.Fit("gneg2", "qRS")
03177
03178 result["fit_ok_2"] = (fit_pos_2.Status()==0 and fit_pos_2.CovMatrixStatus()==3 and fit_neg_2.Status()==0 and fit_neg_2.CovMatrixStatus()==3)
03179 ntot = tmppos_2.GetEntries() + tmpneg_2.GetEntries()
03180 result["fit_peak_2"] = (fpos_2.GetParameter(1)*tmppos_2.GetEntries() + fneg_2.GetParameter(1)*tmpneg_2.GetEntries()) / ntot
03181 result["fit_peak_error_2"] = sqrt( (fpos_2.GetParError(1)*tmppos_2.GetEntries())**2 + (fneg_2.GetParError(1)*tmpneg_2.GetEntries())**2) / ntot
03182
03183 return result
03184
03185
03186
03187 def segdiffvsphi_xalign(tfiles, wheel, window=10.):
03188 tdrStyle.SetOptTitle(1)
03189 tdrStyle.SetTitleBorderSize(1)
03190 tdrStyle.SetTitleFontSize(0.05)
03191
03192 global htemp, gtemp_12, gtemp_21, gtemp_11, tlegend
03193 htemp = ROOT.TH1F("htemp", "", 1, -pi, pi)
03194 gtemp_11_phi, gtemp_11_val, gtemp_11_err = [], [], []
03195 gtemp_12_phi, gtemp_12_val, gtemp_12_err = [], [], []
03196 gtemp_21_phi, gtemp_21_val, gtemp_21_err = [], [], []
03197 for sector in xrange(1, 12+1):
03198
03199 r1 = segdiff_xalign(tfiles, "x_dt1_csc", wheel=wheel, sector=sector, cscstations = "12")
03200 r2 = segdiff_xalign(tfiles, "x_dt2_csc", wheel=wheel, sector=sector, cscstations = "1")
03201
03202 if r1["fit_ok"]:
03203 gtemp_11_phi.append(r1["phi"])
03204 gtemp_11_val.append(r1["fit_peak"])
03205 gtemp_11_err.append(r1["fit_peak_error"])
03206 if r1["fit_ok_2"]:
03207 gtemp_12_phi.append(r1["phi"])
03208 gtemp_12_val.append(r1["fit_peak_2"])
03209 gtemp_12_err.append(r1["fit_peak_error_2"])
03210 if r2["fit_ok"]:
03211 gtemp_21_phi.append(r2["phi"])
03212 gtemp_21_val.append(r2["fit_peak"])
03213 gtemp_21_err.append(r2["fit_peak_error"])
03214
03215
03216
03217
03218 if len(gtemp_11_phi) > 0:
03219 gtemp_11 = ROOT.TGraphErrors(len(gtemp_11_phi), array.array("d", gtemp_11_phi), array.array("d", gtemp_11_val),
03220 array.array("d", [0.] * len(gtemp_11_phi)), array.array("d", gtemp_11_err))
03221 if len(gtemp_12_phi) > 0:
03222 gtemp_12 = ROOT.TGraphErrors(len(gtemp_12_phi), array.array("d", gtemp_12_phi), array.array("d", gtemp_12_val),
03223 array.array("d", [0.] * len(gtemp_12_phi)), array.array("d", gtemp_12_err))
03224 if len(gtemp_11_phi) > 0:
03225 gtemp_21 = ROOT.TGraphErrors(len(gtemp_21_phi), array.array("d", gtemp_21_phi), array.array("d", gtemp_21_val),
03226 array.array("d", [0.] * len(gtemp_21_phi)), array.array("d", gtemp_21_err))
03227
03228 if len(gtemp_11_phi) > 0:
03229 gtemp_11.SetMarkerStyle(20); gtemp_11.SetMarkerSize(1.5);
03230 gtemp_11.SetMarkerColor(ROOT.kRed); gtemp_11.SetLineColor(ROOT.kRed)
03231 if len(gtemp_12_phi) > 0:
03232 gtemp_12.SetMarkerStyle(22); gtemp_12.SetMarkerSize(1.);
03233 gtemp_12.SetMarkerColor(ROOT.kGreen+2); gtemp_12.SetLineColor(ROOT.kGreen+2)
03234 if len(gtemp_21_phi) > 0:
03235 gtemp_21.SetMarkerStyle(21); gtemp_21.SetMarkerSize(1.5);
03236 gtemp_21.SetMarkerColor(ROOT.kBlue); gtemp_21.SetLineColor(ROOT.kBlue)
03237
03238 htemp.SetTitle("Wheel %+d" % wheel)
03239 htemp.SetAxisRange(-window, window, "Y")
03240 htemp.SetXTitle("#phi of MB")
03241 htemp.SetYTitle("#Deltax^{loc}_{DT} - r_{DT}/r_{CSC}#times#Deltax^{loc}_{CSC} (mm)")
03242 htemp.GetXaxis().CenterTitle()
03243 htemp.GetYaxis().CenterTitle()
03244 htemp.GetYaxis().SetTitleOffset(0.75)
03245
03246 c1.Clear()
03247 htemp.Draw()
03248 if len(gtemp_12_phi) > 0:
03249 gtemp_12.Draw("p")
03250 if len(gtemp_21_phi) > 0:
03251 gtemp_21.Draw("p")
03252 if len(gtemp_11_phi) > 0:
03253 gtemp_11.Draw("p")
03254
03255 tlegend = ROOT.TLegend(0.59, 0.75, 0.99, 0.92)
03256 tlegend.SetBorderSize(0)
03257 tlegend.SetFillColor(ROOT.kWhite)
03258 if len(gtemp_11_phi) > 0:
03259 tlegend.AddEntry(gtemp_11, "MB1 - ME1/3 (mean: %4.2f, RMS: %4.2f)" % (mean(gtemp_11_val), stdev(gtemp_11_val)), "pl")
03260 if len(gtemp_21_phi) > 0:
03261 tlegend.AddEntry(gtemp_21, "MB2 - ME1/3 (mean: %4.2f, RMS: %4.2f)" % (mean(gtemp_21_val), stdev(gtemp_21_val)), "pl")
03262 if len(gtemp_12_phi) > 0:
03263 tlegend.AddEntry(gtemp_12, "MB1 - ME2/2 (mean: %4.2f, RMS: %4.2f)" % (mean(gtemp_12_val), stdev(gtemp_12_val)), "pl")
03264
03265
03266
03267
03268 tlegend.Draw()
03269
03270 f_11 = ROOT.TF1("f11", "[0] + [1]*sin(x) + [2]*cos(x)", -pi, pi)
03271 f_11.SetLineColor(ROOT.kRed)
03272 f_11.SetLineWidth(2)
03273 f_21 = ROOT.TF1("f21", "[0] + [1]*sin(x) + [2]*cos(x)", -pi, pi)
03274 f_21.SetLineColor(ROOT.kBlue)
03275 f_21.SetLineWidth(2)
03276 if len(gtemp_11_phi) > 0:
03277 gtemp_11.Fit(f_11,"")
03278 if len(gtemp_21_phi) > 0:
03279 gtemp_21.Fit(f_21,"")
03280
03281 global f_txt,f_11_txt, f_21_txt
03282 f_txt = ROOT.TLatex(-2.9, -0.7*window, "ME1/3 ring corrections equivalents:")
03283 f_txt.SetTextSize(0.028)
03284 f_txt.Draw()
03285 if len(gtemp_11_phi) > 0:
03286 rdt = signConventions[("DT", 2, 1, 1)][3]*10
03287 f_11_txt = ROOT.TLatex(-2.9, -0.8*window, "#Deltax=%.2f#pm%.2f mm #Deltay=%.2f#pm%.2f mm #Delta#phi_{z}=%.2f#pm%.2f mrad" % (
03288 -f_11.GetParameter(1), f_11.GetParError(1), f_11.GetParameter(2), f_11.GetParError(2), -f_11.GetParameter(0)/rdt*1000, f_11.GetParError(0)/rdt*1000))
03289 f_11_txt.SetTextSize(0.028)
03290 f_11_txt.SetTextColor(ROOT.kRed)
03291 f_11_txt.Draw()
03292 if len(gtemp_11_phi) > 0:
03293 rdt = signConventions[("DT", 2, 2, 1)][3]*10
03294 f_21_txt = ROOT.TLatex(-2.9, -0.9*window, "#Deltax=%.2f#pm%.2f mm #Deltay=%.2f#pm%.2f mm #Delta#phi_{z}=%.2f#pm%.2f mrad" % (
03295 -f_21.GetParameter(1), f_21.GetParError(1), f_21.GetParameter(2), f_21.GetParError(2), -f_21.GetParameter(0)/rdt*1000, f_21.GetParError(0)/rdt*1000))
03296 f_21_txt.SetTextSize(0.028)
03297 f_21_txt.SetTextColor(ROOT.kBlue)
03298 f_21_txt.Draw()
03299
03300
03301
03302 def segdiffvsphi(tfiles, reports, component, wheel, window=5., excludesectors=()):
03303 tdrStyle.SetOptTitle(1)
03304 tdrStyle.SetTitleBorderSize(1)
03305 tdrStyle.SetTitleFontSize(0.05)
03306
03307 global htemp, gtemp_12, gtemp2_12, gtemp_23, gtemp2_23, gtemp_34, gtemp2_34, tlegend
03308 htemp = ROOT.TH1F("htemp", "", 1, -pi, pi)
03309 gtemp_12_phi, gtemp_12_val, gtemp_12_err, gtemp_12_val2, gtemp_12_err2 = [], [], [], [], []
03310 gtemp_23_phi, gtemp_23_val, gtemp_23_err, gtemp_23_val2, gtemp_23_err2 = [], [], [], [], []
03311 gtemp_34_phi, gtemp_34_val, gtemp_34_err, gtemp_34_val2, gtemp_34_err2 = [], [], [], [], []
03312 for sector in xrange(1, 12+1):
03313
03314 r1_found, r2_found, r3_found, r4_found = False, False, False, False
03315 for r1 in reports:
03316 if r1.postal_address == ("DT", wheel, 1, sector):
03317 r1_found = True
03318 break
03319 for r2 in reports:
03320 if r2.postal_address == ("DT", wheel, 2, sector):
03321 r2_found = True
03322 break
03323 for r3 in reports:
03324 if r3.postal_address == ("DT", wheel, 3, sector):
03325 r3_found = True
03326 break
03327 for r4 in reports:
03328 if r4.postal_address == ("DT", wheel, 4, sector):
03329 r4_found = True
03330 break
03331
03332
03333 if sector not in excludesectors:
03334 if r1_found and r2_found and r1.status == "PASS" and r2.status == "PASS":
03335 phi, val, err, val2, err2, fit1, fit2, fit3 = segdiff(tfiles, component, 12, wheel=wheel, sector=sector)
03336
03337 if fit1 and fit2 and fit3:
03338 gtemp_12_phi.append(phi)
03339 gtemp_12_val.append(val)
03340 gtemp_12_err.append(err)
03341 gtemp_12_val2.append(val2)
03342 gtemp_12_err2.append(err2)
03343 if r2_found and r3_found and r2.status == "PASS" and r3.status == "PASS":
03344 phi, val, err, val2, err2, fit1, fit2, fit3 = segdiff(tfiles, component, 23, wheel=wheel, sector=sector)
03345
03346 if fit1 and fit2 and fit3:
03347 gtemp_23_phi.append(phi)
03348 gtemp_23_val.append(val)
03349 gtemp_23_err.append(err)
03350 gtemp_23_val2.append(val2)
03351 gtemp_23_err2.append(err2)
03352 if component[:4] == "dt13":
03353 if r3_found and r4_found and r3.status == "PASS" and r4.status == "PASS":
03354 phi, val, err, val2, err2, fit1, fit2, fit3 = segdiff(tfiles, component, 34, wheel=wheel, sector=sector)
03355
03356 if fit1 and fit2 and fit3:
03357 gtemp_34_phi.append(phi)
03358 gtemp_34_val.append(val)
03359 gtemp_34_err.append(err)
03360 gtemp_34_val2.append(val2)
03361 gtemp_34_err2.append(err2)
03362
03363
03364
03365
03366 if len(gtemp_12_phi) > 0:
03367 gtemp_12 = ROOT.TGraphErrors(len(gtemp_12_phi), array.array("d", gtemp_12_phi), array.array("d", gtemp_12_val),
03368 array.array("d", [0.] * len(gtemp_12_phi)), array.array("d", gtemp_12_err))
03369 gtemp2_12 = ROOT.TGraphErrors(len(gtemp_12_phi), array.array("d", gtemp_12_phi), array.array("d", gtemp_12_val2),
03370 array.array("d", [0.] * len(gtemp_12_phi)), array.array("d", gtemp_12_err2))
03371 if len(gtemp_23_phi) > 0:
03372 gtemp_23 = ROOT.TGraphErrors(len(gtemp_23_phi), array.array("d", gtemp_23_phi), array.array("d", gtemp_23_val),
03373 array.array("d", [0.] * len(gtemp_23_phi)), array.array("d", gtemp_23_err))
03374 gtemp2_23 = ROOT.TGraphErrors(len(gtemp_23_phi), array.array("d", gtemp_23_phi), array.array("d", gtemp_23_val2),
03375 array.array("d", [0.] * len(gtemp_23_phi)), array.array("d", gtemp_23_err2))
03376 if len(gtemp_34_phi) > 0:
03377 gtemp_34 = ROOT.TGraphErrors(len(gtemp_34_phi), array.array("d", gtemp_34_phi), array.array("d", gtemp_34_val),
03378 array.array("d", [0.] * len(gtemp_34_phi)), array.array("d", gtemp_34_err))
03379 gtemp2_34 = ROOT.TGraphErrors(len(gtemp_34_phi), array.array("d", gtemp_34_phi), array.array("d", gtemp_34_val2),
03380 array.array("d", [0.] * len(gtemp_34_phi)), array.array("d", gtemp_34_err2))
03381
03382 if len(gtemp_12_phi) > 0:
03383 gtemp_12.SetMarkerStyle(20); gtemp_12.SetMarkerSize(1.);
03384 gtemp_12.SetMarkerColor(ROOT.kBlue); gtemp_12.SetLineColor(ROOT.kBlue)
03385 gtemp2_12.SetMarkerStyle(24); gtemp2_12.SetMarkerSize(1.);
03386 gtemp2_12.SetMarkerColor(ROOT.kBlue); gtemp2_12.SetLineColor(ROOT.kBlue)
03387 if len(gtemp_23_phi) > 0:
03388 gtemp_23.SetMarkerStyle(21); gtemp_23.SetMarkerSize(1.);
03389 gtemp_23.SetMarkerColor(ROOT.kRed); gtemp_23.SetLineColor(ROOT.kRed)
03390 gtemp2_23.SetMarkerStyle(25); gtemp2_23.SetMarkerSize(1.);
03391 gtemp2_23.SetMarkerColor(ROOT.kRed); gtemp2_23.SetLineColor(ROOT.kRed)
03392 if len(gtemp_34_phi) > 0 and component[:4] == "dt13":
03393 gtemp_34.SetMarkerStyle(22); gtemp_34.SetMarkerSize(1.25);
03394 gtemp_34.SetMarkerColor(ROOT.kGreen+2); gtemp_34.SetLineColor(ROOT.kGreen+2)
03395 gtemp2_34.SetMarkerStyle(26); gtemp2_34.SetMarkerSize(1.25);
03396 gtemp2_34.SetMarkerColor(ROOT.kGreen+2); gtemp2_34.SetLineColor(ROOT.kGreen+2)
03397
03398 if wheel == 0: htemp.SetTitle("Wheel %d" % wheel)
03399 else: htemp.SetTitle("Wheel %+d" % wheel)
03400 htemp.SetAxisRange(-window, window, "Y")
03401 htemp.SetXTitle("Average #phi of pair (rad)")
03402 if component == "dt13_resid": htemp.SetYTitle("#Deltax^{local} (mm)")
03403 if component == "dt13_slope": htemp.SetYTitle("#Deltadx/dz^{local} (mrad)")
03404 if component == "dt2_resid": htemp.SetYTitle("#Deltay^{local} (mm)")
03405 if component == "dt2_slope": htemp.SetYTitle("#Deltady/dz^{local} (mrad)")
03406 htemp.GetXaxis().CenterTitle()
03407 htemp.GetYaxis().CenterTitle()
03408 htemp.GetYaxis().SetTitleOffset(0.75)
03409
03410 c1.Clear()
03411 htemp.Draw()
03412 if len(gtemp_12_phi) > 0:
03413 gtemp_12.Draw("p")
03414 gtemp2_12.Draw("p")
03415 if len(gtemp_23_phi) > 0:
03416 gtemp_23.Draw("p")
03417 gtemp2_23.Draw("p")
03418 if len(gtemp_34_phi) > 0:
03419 gtemp_34.Draw("p")
03420 gtemp2_34.Draw("p")
03421
03422 tlegend = ROOT.TLegend(0.5, 0.72, 0.9, 0.92)
03423 tlegend.SetBorderSize(0)
03424 tlegend.SetFillColor(ROOT.kWhite)
03425 if len(gtemp_12_phi) > 0:
03426 tlegend.AddEntry(gtemp_12, "MB1 - MB2 (mean: %4.2f, RMS: %4.2f)" % (mean(gtemp_12_val), stdev(gtemp_12_val)), "pl")
03427 if len(gtemp_23_phi) > 0:
03428 tlegend.AddEntry(gtemp_23, "MB2 - MB3 (mean: %4.2f, RMS: %4.2f)" % (mean(gtemp_23_val), stdev(gtemp_23_val)), "pl")
03429 if len(gtemp_34_phi) > 0:
03430 tlegend.AddEntry(gtemp_34, "MB3 - MB4 (mean: %4.2f, RMS: %4.2f)" % (mean(gtemp_34_val), stdev(gtemp_34_val)), "pl")
03431 if len(gtemp_12_phi) > 0:
03432 tlegend.AddEntry(gtemp_12, "total mean: %4.2f, total RMS: %4.2f" % \
03433 (mean(gtemp_12_val + gtemp_23_val + gtemp_34_val),
03434 stdev(gtemp_12_val + gtemp_23_val + gtemp_34_val)), "")
03435 tlegend.Draw()
03436
03437
03438
03439
03440 def segdiffvsphicsc(tfiles, component, pair, window=5., **args):
03441 tdrStyle.SetOptTitle(1)
03442 tdrStyle.SetTitleBorderSize(1)
03443 tdrStyle.SetTitleFontSize(0.05)
03444
03445 if not component[0:3] == "csc": Exception
03446
03447 endcap = args["endcap"]
03448 if endcap=="m":
03449 endcapnum=2
03450 endcapsign="-"
03451 elif endcap=="p":
03452 endcapnum=1
03453 endcapsign="+"
03454 else: raise Exception
03455
03456 station1 = int(str(pair)[0])
03457 station2 = int(str(pair)[1])
03458 if not station2-station1==1: raise Exception
03459
03460 rings = [1,2]
03461 if station2==4: rings = [1]
03462
03463
03464 global htemp, gtemp_1, gtemp2_1, gtemp_2, gtemp2_2, tlegend
03465 htemp = ROOT.TH1F("htemp", "", 1, -pi*5./180., pi*(2.-5./180.))
03466 gtemp_1_phi, gtemp_1_val, gtemp_1_err, gtemp_1_val2, gtemp_1_err2 = [], [], [], [], []
03467 gtemp_2_phi, gtemp_2_val, gtemp_2_err, gtemp_2_val2, gtemp_2_err2 = [], [], [], [], []
03468
03469 for ring in rings:
03470 chambers = xrange(1,37)
03471 if ring == 1: chambers = xrange(1,19)
03472
03473 for chamber in chambers:
03474 phi, val, err, val2, err2, fit1, fit2, fit3 = segdiff(tfiles, component, pair, endcap=endcap, ring=ring, chamber=chamber)
03475 if fit1 and fit2 and fit3:
03476 if ring==1:
03477 gtemp_1_phi.append(phi)
03478 gtemp_1_val.append(val)
03479 gtemp_1_err.append(err)
03480 gtemp_1_val2.append(val2)
03481 gtemp_1_err2.append(err2)
03482 if ring==2:
03483 gtemp_2_phi.append(phi)
03484 gtemp_2_val.append(val)
03485 gtemp_2_err.append(err)
03486 gtemp_2_val2.append(val2)
03487 gtemp_2_err2.append(err2)
03488
03489
03490
03491
03492 if len(gtemp_1_phi) > 0:
03493 gtemp_1 = ROOT.TGraphErrors(len(gtemp_1_phi), array.array("d", gtemp_1_phi), array.array("d", gtemp_1_val),
03494 array.array("d", [0.] * len(gtemp_1_phi)), array.array("d", gtemp_1_err))
03495 gtemp2_1 = ROOT.TGraphErrors(len(gtemp_1_phi), array.array("d", gtemp_1_phi), array.array("d", gtemp_1_val2),
03496 array.array("d", [0.] * len(gtemp_1_phi)), array.array("d", gtemp_1_err2))
03497 if len(gtemp_2_phi) > 0:
03498 gtemp_2 = ROOT.TGraphErrors(len(gtemp_2_phi), array.array("d", gtemp_2_phi), array.array("d", gtemp_2_val),
03499 array.array("d", [0.] * len(gtemp_2_phi)), array.array("d", gtemp_2_err))
03500 gtemp2_2 = ROOT.TGraphErrors(len(gtemp_2_phi), array.array("d", gtemp_2_phi), array.array("d", gtemp_2_val2),
03501 array.array("d", [0.] * len(gtemp_2_phi)), array.array("d", gtemp_2_err2))
03502
03503 if len(gtemp_1_phi) > 0:
03504 gtemp_1.SetMarkerStyle(20); gtemp_1.SetMarkerSize(1.);
03505 gtemp_1.SetMarkerColor(ROOT.kBlue); gtemp_1.SetLineColor(ROOT.kBlue)
03506 gtemp2_1.SetMarkerStyle(24); gtemp2_1.SetMarkerSize(1.);
03507 gtemp2_1.SetMarkerColor(ROOT.kBlue); gtemp2_1.SetLineColor(ROOT.kBlue)
03508 if len(gtemp_2_phi) > 0:
03509 gtemp_2.SetMarkerStyle(21); gtemp_2.SetMarkerSize(1.);
03510 gtemp_2.SetMarkerColor(ROOT.kRed); gtemp_2.SetLineColor(ROOT.kRed)
03511 gtemp2_2.SetMarkerStyle(25); gtemp2_2.SetMarkerSize(1.);
03512 gtemp2_2.SetMarkerColor(ROOT.kRed); gtemp2_2.SetLineColor(ROOT.kRed)
03513
03514 htemp.SetTitle("ME%s%d - ME%s%d" % (endcapsign,station2,endcapsign,station1))
03515 htemp.SetAxisRange(-window, window, "Y")
03516 htemp.SetXTitle("Average #phi of pair (rad)")
03517 if component == "csc_resid": htemp.SetYTitle("#Delta(r#phi)^{local} (mm)")
03518 if component == "csc_slope": htemp.SetYTitle("#Deltad(r#phi)/dz^{local} (mrad)")
03519 htemp.GetXaxis().CenterTitle()
03520 htemp.GetYaxis().CenterTitle()
03521 htemp.GetYaxis().SetTitleOffset(0.75)
03522
03523 c1.Clear()
03524 htemp.Draw()
03525 if len(gtemp_1_phi) > 0:
03526 gtemp_1.Draw("p")
03527 gtemp2_1.Draw("p")
03528 if len(gtemp_2_phi) > 0:
03529 gtemp_2.Draw("p")
03530 gtemp2_2.Draw("p")
03531
03532 tlegend = ROOT.TLegend(0.5, 0.72, 0.9, 0.92)
03533 tlegend.SetBorderSize(0)
03534 tlegend.SetFillColor(ROOT.kWhite)
03535 if len(gtemp_1_phi) > 0:
03536 tlegend.AddEntry(gtemp_1, "ring 1 (mean: %4.2f, RMS: %4.2f)" % (mean(gtemp_1_val), stdev(gtemp_1_val)), "pl")
03537 if len(gtemp_2_phi) > 0:
03538 tlegend.AddEntry(gtemp_2, "ring 2 (mean: %4.2f, RMS: %4.2f)" % (mean(gtemp_2_val), stdev(gtemp_2_val)), "pl")
03539
03540
03541
03542
03543 tlegend.Draw()
03544
03545
03546
03547
03548
03549
03550
03551 def corrections2D(reportsX=None, reportsY=None, geometry0=None, geometryX=None, geometryY=None,
03552 window=25., selection=None, name="tmp", canvas=None, pre_title_x=None, pre_title_y=None,
03553 which="110011"):
03554
03555 tdrStyle.SetOptStat(0)
03556 tdrStyle.SetStatW(0.40)
03557
03558
03559 mode = None
03560 check_reports = False
03561 if reportsX is not None and reportsY is not None:
03562 mode = "reports"
03563 check_reports = True
03564 if geometry0 is not None and geometryX is not None and geometryY is not None:
03565 mode = "xmls"
03566 if mode is None:
03567 print "Either couple of reports or three geometries have to be given as input. Exiting..."
03568 return
03569
03570
03571 wnd = [window]*6
03572 wnd_adaptive = [.1]*6
03573
03574 global hx, hy, hz, hphix, hphiy, hphiz
03575 bins=2000
03576 hx = ROOT.TH2F("%s_x" % name, "", bins, -wnd[0], wnd[0], bins, -wnd[0], wnd[0])
03577 hy = ROOT.TH2F("%s_y" % name, "", bins, -wnd[1], wnd[1], bins, -wnd[1], wnd[1])
03578 hz = ROOT.TH2F("%s_z" % name, "", bins, -wnd[2], wnd[2], bins, -wnd[2], wnd[2])
03579 hphix = ROOT.TH2F("%s_phix" % name, "", bins, -wnd[3], wnd[3], bins, -wnd[3], wnd[3])
03580 hphiy = ROOT.TH2F("%s_phiy" % name, "", bins, -wnd[4], wnd[4], bins, -wnd[4], wnd[4])
03581 hphiz = ROOT.TH2F("%s_phiz" % name, "", bins, -wnd[5], wnd[5], bins, -wnd[5], wnd[5])
03582 hhh = [hx, hy, hz, hphix, hphiy, hphiz]
03583
03584
03585 global pca_x, pca_y, pca_z, pca_phix, pca_phiy, pca_phiz
03586 pca_x = ROOT.TPrincipal(2,"D")
03587 pca_y = ROOT.TPrincipal(2,"D")
03588 pca_z = ROOT.TPrincipal(2,"D")
03589 pca_phix = ROOT.TPrincipal(2,"D")
03590 pca_phiy = ROOT.TPrincipal(2,"D")
03591 pca_phiz = ROOT.TPrincipal(2,"D")
03592 pcas = [pca_x, pca_y, pca_z, pca_phix, pca_phiy, pca_phiz]
03593
03594
03595 ax=[]; ay=[]; az=[]; aphix=[]; aphiy=[]; aphiz=[]
03596 aaa = [ax, ay, az, aphix, aphiy, aphiz]
03597
03598
03599 postal_addresses = []
03600
03601
03602 if check_reports:
03603 for r1 in reportsX:
03604
03605 if r1.postal_address[0]=='CSC' and r1.postal_address[2]==1 and r1.postal_address[3]==4: continue
03606 if selection is None or (selection.func_code.co_argcount == len(r1.postal_address) and selection(*r1.postal_address)):
03607 r2 = getReportByPostalAddress(r1.postal_address, reportsY)
03608 if r2 is None:
03609 print "bad r2 in ",r1.postal_address
03610 continue
03611
03612 if r1.status != "PASS" or r2.status != "PASS":
03613 print "bad status", r1.postal_address, r1.status, r2.status
03614 continue
03615 postal_addresses.append(r1.postal_address)
03616
03617 else:
03618 for key in geometry0.dt.keys():
03619 if len(key)==3 and key in geometryX.dt and key in geometryY.dt:
03620 postal_addresses.append( tuple(['DT'] + list(key)) )
03621 for key in geometry0.csc.keys():
03622
03623 if key[2]==1 and key[3]==4: continue
03624 if len(key)==4 and key in geometryX.csc and key in geometryY.csc:
03625 postal_addresses.append( tuple(['CSC'] + list(key)) )
03626
03627
03628 for addr in postal_addresses:
03629
03630
03631 if not (selection is None or (selection.func_code.co_argcount == len(addr) and selection(*addr)) ): continue
03632
03633 factors = [10. * signConventions[addr][0], 10. * signConventions[addr][1], 10. * signConventions[addr][2],
03634 1000., 1000., 1000. ]
03635
03636 if check_reports:
03637 rX = getReportByPostalAddress(addr, reportsX)
03638 rY = getReportByPostalAddress(addr, reportsY)
03639 deltasX = [rX.deltax, rX.deltay, rX.deltaz, rX.deltaphix, rX.deltaphiy, rX.deltaphiz]
03640 deltasY = [rY.deltax, rY.deltay, rY.deltaz, rY.deltaphix, rY.deltaphiy, rY.deltaphiz]
03641
03642 if mode == "reports":
03643
03644 checks = map( lambda d1, d2: d1 is not None and d2 is not None and d1.error is not None \
03645 and d2.error is not None and (d1.error**2 + d2.error**2) > 0. , \
03646 deltasX, deltasY)
03647
03648 for i in range(len(checks)):
03649 if not checks[i]: continue
03650 fillX = deltasX[i].value * factors[i]
03651 fillY = deltasY[i].value * factors[i]
03652 aaa[i].append([fillX,fillY])
03653 pcas[i].AddRow(array.array('d',[fillX,fillY]))
03654 mx = max(abs(fillX), abs(fillY))
03655 if mx > wnd_adaptive[i]: wnd_adaptive[i] = mx
03656
03657 if mode == "xmls":
03658
03659 db0 = dbX = dbY = None
03660 if addr[0] == "DT":
03661 db0, dbX, dbY = geometry0.dt[addr[1:]], geometryX.dt[addr[1:]], geometryY.dt[addr[1:]]
03662 if addr[0] == 'CSC':
03663 db0, dbX, dbY = geometry0.csc[addr[1:]], geometryX.csc[addr[1:]], geometryY.csc[addr[1:]]
03664
03665 checks = [True]*6
03666 if check_reports:
03667 checks = map( lambda d1, d2: d1 is not None and d2 is not None , deltasX, deltasY)
03668
03669 gdeltas0 = [db0.x, db0.y, db0.z, db0.phix, db0.phiy, db0.phiz]
03670 gdeltasX = [dbX.x, dbX.y, dbX.z, dbX.phix, dbX.phiy, dbX.phiz]
03671 gdeltasY = [dbY.x, dbY.y, dbY.z, dbY.phix, dbY.phiy, dbY.phiz]
03672
03673 for i in range(len(checks)):
03674 if not checks[i]: continue
03675 fillX = (gdeltasX[i] - gdeltas0[i]) * factors[i]
03676 fillY = (gdeltasY[i] - gdeltas0[i]) * factors[i]
03677 aaa[i].append([fillX,fillY])
03678 pcas[i].AddRow(array.array('d',[fillX,fillY]))
03679 mx = max(abs(fillX), abs(fillY))
03680 if mx > wnd_adaptive[i]: wnd_adaptive[i] = mx
03681
03682
03683
03684
03685 if mode == "xmls":
03686 if pre_title_x is None: pre_title_x = "geometry 1 "
03687 if pre_title_y is None: pre_title_y = "geometry 2 "
03688 if mode == "reports":
03689 if pre_title_x is None: pre_title_x = "iteration's "
03690 if pre_title_y is None: pre_title_y = "other iteration's "
03691 tmptitles = ["#Deltax (mm)", "#Deltay (mm)", "#Deltaz (mm)",
03692 "#Delta#phi_{x} (mrad)", "#Delta#phi_{y} (mrad)", "#Delta#phi_{z} (mrad)"]
03693 htitles = []
03694 for t in tmptitles: htitles.append([pre_title_x + t, pre_title_y + t])
03695
03696 if canvas is not None: c = canvas
03697 else: c = c1
03698 c.Clear()
03699 ndraw = which.count('1')
03700 if ndraw > 4: c.Divide(3, 2)
03701 elif ndraw > 2: c.Divide(2, 2)
03702 elif ndraw > 1: c.Divide(2, 1)
03703
03704 global lines, graphs, texs
03705 lines = []; graphs = []; texs = []
03706
03707 ipad = 0
03708 for i in range(6):
03709
03710
03711 if ( int(which,2) & (1<<i) ) == 0: continue
03712
03713 ipad += 1
03714 c.GetPad(ipad).cd()
03715 c.GetPad(ipad).SetGridx(1)
03716 c.GetPad(ipad).SetGridy(1)
03717
03718 wn = 1.08 * wnd_adaptive[i]
03719 hhh[i].GetXaxis().SetRangeUser(-wn, wn)
03720 hhh[i].GetYaxis().SetRangeUser(-wn, wn)
03721 hhh[i].SetXTitle(htitles[i][0])
03722 hhh[i].SetYTitle(htitles[i][1])
03723 hhh[i].GetXaxis().CenterTitle()
03724 hhh[i].GetYaxis().CenterTitle()
03725 hhh[i].Draw()
03726
03727 if len(aaa[i]) == 0: continue
03728
03729 a1, a2 = map( lambda x: array.array('d',x), zip(*aaa[i]) )
03730 g = ROOT.TGraph(len(a1), a1, a2)
03731 g.SetMarkerStyle(5)
03732 g.SetMarkerSize(0.3)
03733 g.SetMarkerColor(ROOT.kBlue)
03734 graphs.append(g)
03735
03736 pcas[i].MakePrincipals()
03737
03738
03739 b = pcas[i].GetEigenVectors()(1,0) / pcas[i].GetEigenVectors()(0,0)
03740 a = pcas[i].GetMeanValues()[1] - b * pcas[i].GetMeanValues()[0]
03741
03742
03743 cov = pcas[i].GetCovarianceMatrix()
03744 r = cov(0,1)/sqrt(cov(1,1)*cov(0,0))
03745 print "r, RMSx, RMSy =", r, g.GetRMS(1), g.GetRMS(2)
03746 texrms = ROOT.TLatex(0.17,0.87, "RMS x,y = %.02g, %.02g" % (g.GetRMS(1),g.GetRMS(2)))
03747 texr = ROOT.TLatex(0.17,0.80, "r = %.02g" % r)
03748 for t in texr, texrms:
03749 t.SetNDC(1)
03750 t.SetTextColor(ROOT.kBlue)
03751 t.SetTextSize(0.053)
03752 t.Draw()
03753 texs.append(t)
03754
03755 g.Draw("p")
03756
03757 if not isnan(b):
03758 wn = wnd_adaptive[i]
03759 line = ROOT.TLine(-wn, a - b*wn, wn, a + b*wn)
03760 line.SetLineColor(ROOT.kRed)
03761 line.Draw()
03762 lines.append(line)
03763
03764
03765 return aaa
1.7.3