Reference Manual

Go to the documentation of this file.
00001 import ROOT, array, os, re, random
00002 from math import *
00003 import time
00004 import pickle
00005 
00006 # python 2.6 has json modue; <2.6 could use simplejson
00007 try:
00008   import json
00009 except ImportError:
00010   import simplejson as json
00011 
00012 # my common muon types structures
00013 from mutypes import *
00014 
00015 CPP_LOADED = False
00016 
00017 # containers for test results for map plots
00018 MAP_RESULTS_SAWTOOTH = {}
00019 MAP_RESULTS_FITSIN = {}
00020 MAP_RESULTS_BINS = {}
00021 
00022 # general container for all test results
00023 TEST_RESULTS = {}
00024 
00025 #############################################################
00026 # Convenience functions
00027 
00028 def wheelm2only(dt, wheel, station, sector): return dt == "DT" and wheel == -2
00029 def wheelm1only(dt, wheel, station, sector): return dt == "DT" and wheel == -1
00030 def wheel0only(dt, wheel, station, sector): return dt == "DT" and wheel == 0
00031 def wheelp1only(dt, wheel, station, sector): return dt == "DT" and wheel == 1
00032 def wheelp2only(dt, wheel, station, sector): return dt == "DT" and wheel == 2
00033 
00034 def wheelLetter(wheel):
00035   if   wheel == -2: return "A"
00036   elif wheel == -1: return "B"
00037   elif wheel ==  0: return "C"
00038   elif wheel == +1: return "D"
00039   elif wheel == +2: return "E"
00040   else: raise Exception
00041 
00042 def wheelNumber(wheell):
00043   if   wheell == "A": return -2
00044   elif wheell == "B": return -1
00045   elif wheell == "C": return 0
00046   elif wheell == "D": return 1
00047   elif wheell == "E": return 2
00048   else: raise Exception
00049 
00050 def mean(xlist):
00051   s, n = 0., 0.
00052   for x in xlist:
00053     s += x
00054     n += 1.
00055   return s/n
00056 
00057 def rms(xlist):
00058   s2, n = 0., 0.
00059   for x in xlist:
00060     s2 += x**2
00061     n += 1.
00062   return sqrt(s2/n)
00063 
00064 def stdev(xlist):
00065   s, s2, n = 0., 0., 0.
00066   for x in xlist:
00067     s += x
00068     s2 += x**2
00069     n += 1.
00070   return sqrt(s2/n - (s/n)**2)
00071 
00072 def wmean(xlist):
00073   s, w = 0., 0.
00074   for x, e in xlist:
00075     if e > 0.:
00076       wi = 1./e**2
00077       s += x*wi
00078       w += wi
00079   return s/w, sqrt(1./w)
00080 
00081 #############################################################
00082 
00083 tdrStyle = None
00084 def setTDRStyle():
00085   global tdrStyle
00086   tdrStyle = ROOT.TStyle("tdrStyle","Style for P-TDR")
00087 # For the canvas:
00088   tdrStyle.SetCanvasBorderMode(0)
00089   tdrStyle.SetCanvasColor(ROOT.kWhite)
00090   tdrStyle.SetCanvasDefH(600) #Height of canvas
00091   tdrStyle.SetCanvasDefW(600) #Width of canvas
00092   tdrStyle.SetCanvasDefX(0)   #POsition on screen
00093   tdrStyle.SetCanvasDefY(0)
00094 
00095 # For the Pad:
00096   tdrStyle.SetPadBorderMode(0)
00097   # tdrStyle.SetPadBorderSize(Width_t size = 1)
00098   tdrStyle.SetPadColor(ROOT.kWhite)
00099   tdrStyle.SetPadGridX(False)
00100   tdrStyle.SetPadGridY(False)
00101   tdrStyle.SetGridColor(0)
00102   tdrStyle.SetGridStyle(3)
00103   tdrStyle.SetGridWidth(1)
00104 
00105 # For the frame:
00106   tdrStyle.SetFrameBorderMode(0)
00107   tdrStyle.SetFrameBorderSize(1)
00108   tdrStyle.SetFrameFillColor(0)
00109   tdrStyle.SetFrameFillStyle(0)
00110   tdrStyle.SetFrameLineColor(1)
00111   tdrStyle.SetFrameLineStyle(1)
00112   tdrStyle.SetFrameLineWidth(1)
00113 
00114 # For the histo:
00115   # tdrStyle.SetHistFillColor(1)
00116   # tdrStyle.SetHistFillStyle(0)
00117   tdrStyle.SetHistLineColor(1)
00118   tdrStyle.SetHistLineStyle(0)
00119   tdrStyle.SetHistLineWidth(1)
00120   # tdrStyle.SetLegoInnerR(Float_t rad = 0.5)
00121   # tdrStyle.SetNumberContours(Int_t number = 20)
00122 
00123   tdrStyle.SetEndErrorSize(2)
00124 #  tdrStyle.SetErrorMarker(20)
00125   tdrStyle.SetErrorX(0.)
00126 
00127   tdrStyle.SetMarkerStyle(20)
00128 
00129 #For the fit/function:
00130   tdrStyle.SetOptFit(1)
00131   tdrStyle.SetFitFormat("5.4g")
00132   tdrStyle.SetFuncColor(2)
00133   tdrStyle.SetFuncStyle(1)
00134   tdrStyle.SetFuncWidth(1)
00135 
00136 #For the date:
00137   tdrStyle.SetOptDate(0)
00138   # tdrStyle.SetDateX(Float_t x = 0.01)
00139   # tdrStyle.SetDateY(Float_t y = 0.01)
00140 
00141 # For the statistics box:
00142   tdrStyle.SetOptFile(0)
00143   tdrStyle.SetOptStat(0) # To display the mean and RMS:   SetOptStat("mr")
00144   tdrStyle.SetStatColor(ROOT.kWhite)
00145   tdrStyle.SetStatFont(42)
00146   tdrStyle.SetStatFontSize(0.025)
00147   tdrStyle.SetStatTextColor(1)
00148   tdrStyle.SetStatFormat("6.4g")
00149   tdrStyle.SetStatBorderSize(1)
00150   tdrStyle.SetStatH(0.1)
00151   tdrStyle.SetStatW(0.15)
00152   # tdrStyle.SetStatStyle(Style_t style = 1001)
00153   # tdrStyle.SetStatX(Float_t x = 0)
00154   # tdrStyle.SetStatY(Float_t y = 0)
00155 
00156 # Margins:
00157   tdrStyle.SetPadTopMargin(0.05)
00158   tdrStyle.SetPadBottomMargin(0.13)
00159   tdrStyle.SetPadLeftMargin(0.13)
00160   tdrStyle.SetPadRightMargin(0.05)
00161 
00162 # For the Global title:
00163   tdrStyle.SetOptTitle(0)
00164   tdrStyle.SetTitleFont(42)
00165   tdrStyle.SetTitleColor(1)
00166   tdrStyle.SetTitleTextColor(1)
00167   tdrStyle.SetTitleFillColor(10)
00168   tdrStyle.SetTitleFontSize(0.05)
00169   # tdrStyle.SetTitleH(0) # Set the height of the title box
00170   # tdrStyle.SetTitleW(0) # Set the width of the title box
00171   # tdrStyle.SetTitleX(0) # Set the position of the title box
00172   # tdrStyle.SetTitleY(0.985) # Set the position of the title box
00173   # tdrStyle.SetTitleStyle(Style_t style = 1001)
00174   # tdrStyle.SetTitleBorderSize(2)
00175 
00176 # For the axis titles:
00177   tdrStyle.SetTitleColor(1, "XYZ")
00178   tdrStyle.SetTitleFont(42, "XYZ")
00179   tdrStyle.SetTitleSize(0.06, "XYZ")
00180   # tdrStyle.SetTitleXSize(Float_t size = 0.02) # Another way to set the size?
00181   # tdrStyle.SetTitleYSize(Float_t size = 0.02)
00182   tdrStyle.SetTitleXOffset(0.9)
00183   tdrStyle.SetTitleYOffset(1.05)
00184   # tdrStyle.SetTitleOffset(1.1, "Y") # Another way to set the Offset
00185 
00186 # For the axis labels:
00187   tdrStyle.SetLabelColor(1, "XYZ")
00188   tdrStyle.SetLabelFont(42, "XYZ")
00189   tdrStyle.SetLabelOffset(0.007, "XYZ")
00190   tdrStyle.SetLabelSize(0.05, "XYZ")
00191 
00192 # For the axis:
00193   tdrStyle.SetAxisColor(1, "XYZ")
00194   tdrStyle.SetStripDecimals(True)
00195   tdrStyle.SetTickLength(0.03, "XYZ")
00196   tdrStyle.SetNdivisions(510, "XYZ")
00197   tdrStyle.SetPadTickX(1)  # To get tick marks on the opposite side of the frame
00198   tdrStyle.SetPadTickY(1)
00199 
00200 # Change for log plots:
00201   tdrStyle.SetOptLogx(0)
00202   tdrStyle.SetOptLogy(0)
00203   tdrStyle.SetOptLogz(0)
00204 
00205 # Postscript options:
00206   tdrStyle.SetPaperSize(20.,20.)
00207   # tdrStyle.SetLineScalePS(Float_t scale = 3)
00208   # tdrStyle.SetLineStyleString(Int_t i, const char* text)
00209   # tdrStyle.SetHeaderPS(const char* header)
00210   # tdrStyle.SetTitlePS(const char* pstitle)
00211 
00212   # tdrStyle.SetBarOffset(Float_t baroff = 0.5)
00213   # tdrStyle.SetBarWidth(Float_t barwidth = 0.5)
00214   # tdrStyle.SetPaintTextFormat(const char* format = "g")
00215   # tdrStyle.SetPalette(Int_t ncolors = 0, Int_t* colors = 0)
00216   # tdrStyle.SetTimeOffset(Double_t toffset)
00217   # tdrStyle.SetHistMinimumZero(True)
00218 
00219   tdrStyle.cd()
00220 
00221 setTDRStyle()
00222 
00223 def set_palette(name=None, ncontours=999):
00224     """Set a color palette from a given RGB list
00225     stops, red, green and blue should all be lists of the same length
00226     see set_decent_colors for an example"""
00227 
00228     if name == "halfgray":
00229         stops = [0.00, 0.34, 0.61, 0.84, 1.00]
00230         red   = map(lambda x: 1. - (1.-x)/2., [1.00, 0.84, 0.61, 0.34, 0.00])
00231         green = map(lambda x: 1. - (1.-x)/2., [1.00, 0.84, 0.61, 0.34, 0.00])
00232         blue  = map(lambda x: 1. - (1.-x)/2., [1.00, 0.84, 0.61, 0.34, 0.00])
00233     elif name == "gray":
00234         stops = [0.00, 0.34, 0.61, 0.84, 1.00]
00235         red   = [1.00, 0.84, 0.61, 0.34, 0.00]
00236         green = [1.00, 0.84, 0.61, 0.34, 0.00]
00237         blue  = [1.00, 0.84, 0.61, 0.34, 0.00]
00238     elif name == "blues":
00239         stops = [0.00, 0.34, 0.61, 0.84, 1.00]
00240         red   = [1.00, 0.84, 0.61, 0.34, 0.00]
00241         green = [1.00, 0.84, 0.61, 0.34, 0.00]
00242         blue  = [1.00, 1.00, 1.00, 1.00, 1.00]
00243     elif name == "reds":
00244         stops = [0.00, 0.34, 0.61, 0.84, 1.00]
00245         red   = [1.00, 1.00, 1.00, 1.00, 1.00]
00246         green = [1.00, 0.84, 0.61, 0.34, 0.00]
00247         blue  = [1.00, 0.84, 0.61, 0.34, 0.00]
00248     elif name == "antigray":
00249         stops = [0.00, 0.34, 0.61, 0.84, 1.00]
00250         red   = [1.00, 0.84, 0.61, 0.34, 0.00]
00251         green = [1.00, 0.84, 0.61, 0.34, 0.00]
00252         blue  = [1.00, 0.84, 0.61, 0.34, 0.00]
00253         red.reverse()
00254         green.reverse()
00255         blue.reverse()
00256     elif name == "fire":
00257         stops = [0.00, 0.20, 0.80, 1.00]
00258         red   = [1.00, 1.00, 1.00, 0.50]
00259         green = [1.00, 1.00, 0.00, 0.00]
00260         blue  = [0.20, 0.00, 0.00, 0.00]
00261     elif name == "antifire":
00262         stops = [0.00, 0.20, 0.80, 1.00]
00263         red   = [0.50, 1.00, 1.00, 1.00]
00264         green = [0.00, 0.00, 1.00, 1.00]
00265         blue  = [0.00, 0.00, 0.00, 0.20]
00266     else:
00267         # default palette, looks cool
00268         stops = [0.00, 0.34, 0.61, 0.84, 1.00]
00269         red   = [0.00, 0.00, 0.87, 1.00, 0.51]
00270         green = [0.00, 0.81, 1.00, 0.20, 0.00]
00271         blue  = [0.51, 1.00, 0.12, 0.00, 0.00]
00272 
00273     s = array.array('d', stops)
00274     r = array.array('d', red)
00275     g = array.array('d', green)
00276     b = array.array('d', blue)
00277 
00278     npoints = len(s)
00279     ROOT.TColor.CreateGradientColorTable(npoints, s, r, g, b, ncontours)
00280     ROOT.gStyle.SetNumberContours(ncontours)
00281 
00282 set_palette()
00283 
00284 ######################################################################################################
00285 ## sector phi edges in: me11 me12 me13 me14 me21 me22 me31 me32 me41 me42 mb1 mb2 mb3 mb4
00286 ## index:               0    1    2    3    4    5    6    7    8    9    10  11  12  13
00287 phiedges = [
00288    [0.087266462599716474, 0.26179938550504211, 0.43633230751381297, 0.61086524309298951, 0.78539818789089832, 0.95993106410343132,
00289     1.13446400890134, 1.3089969444805165, 1.4835298664892873, 1.6580627893946129, 1.8325957122999386, 2.0071286343087094,
00290     2.1816615698878858, 2.3561945146857948, 2.5307273908983277, 2.7052603356962366, 2.8797932712754131, 3.0543261932841839,
00291     -3.0543261909900767, -2.8797932680847511, -2.7052603460759803, -2.5307274104968038, -2.3561944656988949, -2.181661589486362,
00292     -2.0071286446884531, -1.8325957091092766, -1.6580627871005058, -1.4835298641951802, -1.3089969412898546, -1.1344640192810838,
00293     -0.95993108370190716, -0.78539813890399834, -0.61086526269146535, -0.43633231789355653, -0.26179938231437999, -0.087266460305609153],
00294    [0.087266462599716474, 0.26179938297741073, 0.43633231700542385, 0.61086526005981812, 0.78539815872971441, 0.95993109326461523,
00295     1.1344639919345114, 1.3089969349889057, 1.4835298690169187, 1.6580627893946129, 1.8325957097723073, 2.0071286438003204, 
00296     2.1816615868547147, 2.3561944855246111, 2.5307274200595118, 2.7052603187294082, 2.879793261783802, 3.0543261958118153, 
00297     -3.0543261909900767, -2.8797932706123825, -2.7052603365843693, -2.5307273935299754, -2.356194494860079, -2.1816615603251783,
00298     -2.0071286616552819, -1.8325957186008877, -1.6580627845728746, -1.4835298641951802, -1.308996943817486, -1.1344640097894729,
00299     -0.95993106673507855, -0.78539816806518226, -0.61086523353028144, -0.43633233486038514, -0.26179939180599088, -0.087266457777977771],
00300    [0.087266462599716474, 0.26179938235213535, 0.43633230952414037, 0.61086523916470359, 0.78539817763669606, 0.95993107435763347,
00301     1.1344640128296259, 1.3089969424701891, 1.4835298696421941, 1.6580627893946129, 1.832595709147032, 2.0071286363190368,
00302     2.1816615659596001, 2.3561945044315924, 2.53072740115253, 2.7052603396245227, 2.8797932692650856, 3.0543261964370907,
00303     -3.0543261909900767, -2.8797932712376579, -2.7052603440656529, -2.53072741442509, -2.3561944759530973, -2.1816615792321596,
00304     -2.0071286407601674, -1.8325957111196041, -1.6580627839475992, -1.4835298641951802, -1.3089969444427614, -1.1344640172707563,
00305     -0.95993108763019308, -0.7853981491582005, -0.61086525243726308, -0.43633231396527061, -0.2617993843247074, -0.087266457152702412],
00306    [0.087266462599716474, 0.26179938550504211, 0.43633230751381297, 0.61086524309298951, 0.78539818789089832, 0.95993106410343132,
00307     1.13446400890134, 1.3089969444805165, 1.4835298664892873, 1.6580627893946129, 1.8325957122999386, 2.0071286343087094,
00308     2.1816615698878858, 2.3561945146857948, 2.5307273908983277, 2.7052603356962366, 2.8797932712754131, 3.0543261932841839,
00309     -3.0543261909900767, -2.8797932680847511, -2.7052603460759803, -2.5307274104968038, -2.3561944656988949, -2.181661589486362,
00310     -2.0071286446884531, -1.8325957091092766, -1.6580627871005058, -1.4835298641951802, -1.3089969412898546, -1.1344640192810838,
00311     -0.95993108370190716, -0.78539813890399834, -0.61086526269146535, -0.43633231789355653, -0.26179938231437999, -0.087266460305609153],
00312    [0.26179938481428705, 0.6108652193791777, 0.95993108859688125, 1.3089969578145848, 1.6580627923794755, 2.0071286538798305,
00313     2.356194498693418, 2.7052603320901376, 3.0543261769037247, -2.8797932687755066, -2.5307274342106156, -2.1816615649929121,
00314     -1.8325956957752083, -1.4835298612103178, -1.1344639997099626, -0.78539815489637521, -0.43633232149965551, -0.087266476686068212],
00315    [0.087266462599716474, 0.26179938871066555, 0.43633231557670243, 0.61086524129631259, 0.785398172964478, 0.95993107902985153,
00316     1.1344640106980168, 1.308996936417627, 1.483529863283664, 1.6580627893946129, 1.8325957155055621, 2.0071286423715993,
00317     2.1816615680912093, 2.3561944997593747, 2.5307274058247482, 2.7052603374929136, 2.8797932632125236, 3.0543261900785605,
00318     -3.0543261909900767, -2.8797932648791278, -2.7052603380130908, -2.5307274122934809, -2.3561944806253154, -2.1816615745599419,
00319     -2.0071286428917765, -1.8325957171721663, -1.6580627903061294, -1.4835298641951802, -1.3089969380842312, -1.1344640112181943,
00320     -0.95993108549858397, -0.78539815383041856, -0.61086524776504503, -0.43633231609687961, -0.26179939037726946, -0.087266463511232586],
00321    [0.26179938498198485, 0.61086523665761272, 0.95993108859688125, 1.3089969405361499, 1.6580627922117777, 2.0071286313120122,
00322     2.3561944778405319, 2.7052603529430232, 3.0543261994715434, -2.8797932686078087, -2.530727416932181, -2.1816615649929121,
00323     -1.8325957130536434, -1.4835298613780155, -1.1344640222777811, -0.78539817574926085, -0.43633230064676976, -0.087266454118249653],
00324    [0.087266462599716474, 0.26179938871066555, 0.43633231557670243, 0.61086524129631259, 0.785398172964478, 0.95993107902985153,
00325     1.1344640106980168, 1.308996936417627, 1.483529863283664, 1.6580627893946129, 1.8325957155055621, 2.0071286423715993,
00326     2.1816615680912093, 2.3561944997593747, 2.5307274058247482, 2.7052603374929136, 2.8797932632125236, 3.0543261900785605,
00327     -3.0543261909900767, -2.8797932648791278, -2.7052603380130908, -2.5307274122934809, -2.3561944806253154, -2.1816615745599419,
00328     -2.0071286428917765, -1.8325957171721663, -1.6580627903061294, -1.4835298641951802, -1.3089969380842312, -1.1344640112181943,
00329     -0.95993108549858397, -0.78539815383041856, -0.61086524776504503, -0.43633231609687961, -0.26179939037726946, -0.087266463511232586],
00330    [0.26179938879942166, 0.61086525092924071, 0.95993108859688125, 1.3089969262645218, 1.6580627883943408, 2.0071286288299772, 
00331     2.3561945088997609, 2.7052603218837943, 3.0543262019535784, -2.8797932647903717, -2.5307274026605526, -2.1816615649929121, 
00332     -1.8325957273252713, -1.4835298651954525, -1.1344640247598159, -0.785398144690032, -0.43633233170599861, -0.087266451636214853],
00333    [0.087266462599716474, 0.26179938871066555, 0.43633231557670243, 0.61086524129631259, 0.785398172964478, 0.95993107902985153,
00334     1.1344640106980168, 1.308996936417627, 1.483529863283664, 1.6580627893946129, 1.8325957155055621, 2.0071286423715993,
00335     2.1816615680912093, 2.3561944997593747, 2.5307274058247482, 2.7052603374929136, 2.8797932632125236, 3.0543261900785605,
00336     -3.0543261909900767, -2.8797932648791278, -2.7052603380130908, -2.5307274122934809, -2.3561944806253154, -2.1816615745599419,
00337     -2.0071286428917765, -1.8325957171721663, -1.6580627903061294, -1.4835298641951802, -1.3089969380842312, -1.1344640112181943,
00338     -0.95993108549858397, -0.78539815383041856, -0.61086524776504503, -0.43633231609687961, -0.26179939037726946, -0.087266463511232586],
00339    [0.35228048120123945, 0.87587781482541827, 1.3994776462193192, 1.923076807996136, 2.4466741416203148, 2.970273973014216,
00340     -2.7893121723885534, -2.2657148387643748, -1.7421150073704739, -1.2185158455936571, -0.69491851196947851, -0.17131868057557731],
00341    [0.22000706229660855, 0.74360690430428489, 1.267204926935573, 1.7908033890915052, 2.3144032310991816, 2.8380012537304697,
00342     -2.9215855912931841, -2.3979857492855081, -1.8743877266542202, -1.3507892644982882, -0.82718942249061178, -0.30359139985932365],
00343    [0.29751957124275596, 0.82111826253905784, 1.3447162969496083, 1.8683158980376524, 2.3919145893339548, 2.915512623744505,
00344     -2.844073082347037, -2.3204743910507353, -1.7968763566401849, -1.2732767555521407, -0.74967806425583894, -0.22608002984528835],
00345    [3.0136655290752188, -2.7530905195097337, -2.2922883025568734, -1.9222915077192773, -1.5707963267948966, -1.2193011458705159,
00346     -0.84930435103291968, -0.38850213408005951, 0.127927124514574, 0.65152597487624719, 1.1322596819239259, 1.5707963267948966, 
00347     2.0093329716658674, 2.4900666787135459]]
00348 
00349 def phiedges2c():
00350   lines = []
00351   for ed in phiedges[:]:
00352     ed.sort()
00353     #print ed
00354     ed.extend([999 for n in range(0,37-len(ed))])
00355     lines.append('{' + ', '.join(map(str, ed)) + '}')
00356   #print lines
00357   res = ', '.join(lines)
00358   ff = open("phiedges_export.h",mode="w")
00359   print>>ff,'double phiedges[14][37] = {' + res + '};'
00360   ff.close()
00361 
00362 class SawTeethFunction:
00363   def __init__(self, name):
00364     self.name = name
00365     self.edges = (phiedges[stationIndex(name)])[:]
00366     self.ed = sorted(self.edges)
00367     # add some padding to the end
00368     self.ed.append(pi+1.)
00369     self.n = len(self.edges)
00370   def __call__(self, xx, par):
00371     # wrap x in the most negative phi sector into positive phi
00372     x = xx[0]
00373     if x < self.ed[0]: x += 2*pi
00374     # locate sector
00375     for i in range(0,self.n):
00376       if x <= self.ed[i]: continue
00377       if x > self.ed[i+1]: continue
00378       return par[i*2] + par[i*2+1]*(x - self.ed[i])
00379     return 0
00380   def pp(self):
00381     print self.name, self.n
00382     print self.edges
00383     print self.ed
00384 
00385 
00386 def stationIndex(name):
00387   if ("MB" in name or "ME" in name):
00388     # assume the name is ID
00389     pa = idToPostalAddress(name)
00390     if pa is None: return None
00391     if pa[0]=="CSC":
00392       if pa[2]==1 and pa[3]==1: return 0
00393       if pa[2]==1 and pa[3]==1: return 1
00394       if pa[2]==1 and pa[3]==1: return 2
00395       if pa[2]==1 and pa[3]==1: return 3
00396       if pa[2]==1 and pa[3]==1: return 4
00397       if pa[2]==1 and pa[3]==1: return 5
00398       if pa[2]==1 and pa[3]==1: return 6
00399       if pa[2]==1 and pa[3]==1: return 7
00400       if pa[2]==1 and pa[3]==1: return 8
00401       if pa[2]==1 and pa[3]==1: return 9
00402     if pa[0]=="DT":
00403       if pa[2]==1: return 10
00404       if pa[2]==2: return 11
00405       if pa[2]==3: return 12
00406       if pa[2]==2: return 13
00407   else:
00408     if ("mem11" in name or "mep11" in name): return 0
00409     if ("mem12" in name or "mep12" in name): return 1
00410     if ("mem13" in name or "mep13" in name): return 2
00411     if ("mem14" in name or "mep14" in name): return 3
00412     if ("mem21" in name or "mep21" in name): return 4
00413     if ("mem22" in name or "mep22" in name): return 5
00414     if ("mem31" in name or "mep31" in name): return 6
00415     if ("mem32" in name or "mep32" in name): return 7
00416     if ("mem41" in name or "mep41" in name): return 8
00417     if ("mem42" in name or "mep42" in name): return 9
00418     if ("st1" in name): return 10
00419     if ("st2" in name): return 11
00420     if ("st3" in name): return 12
00421     if ("st4" in name): return 13
00422 
00423 
00424 
00425 def philines(name, window, abscissa):
00426     global philine_tlines, philine_labels
00427     philine_tlines = []
00428     edges = phiedges[stationIndex(name)]
00429     #print name, len(edges)
00430     for phi in edges:
00431         if abscissa is None or abscissa[0] < phi < abscissa[1]:
00432             philine_tlines.append(ROOT.TLine(phi, -window, phi, window))
00433             philine_tlines[-1].SetLineStyle(2)
00434             philine_tlines[-1].Draw()
00435     if "st" in name: # DT labels
00436         philine_labels = []
00437         edges = edges[:]
00438         edges.sort()
00439         if "st4" in name:
00440             labels = [" 7", " 8", " 9", "14", "10", "11", "12", " 1", " 2", " 3", "13", " 4", " 5", " 6"]
00441         else: 
00442             labels = [" 8", " 9", "10", "11", "12", " 1", " 2", " 3", " 4", " 5", " 6"]
00443             edges = edges[1:]
00444         for phi, label in zip(edges, labels):
00445             littlebit = 0.
00446             if label in (" 7", " 9", "14", "10", "11"): littlebit = 0.05
00447             philine_labels.append(ROOT.TText(phi-0.35+littlebit, -0.9*window, label))
00448             philine_labels[-1].Draw()
00449         philine_labels.append(ROOT.TText(-2.9, -0.75*window, "Sector:"))
00450         philine_labels[-1].Draw()
00451 
00452 def zlines(window, abscissa):
00453     global zline_tlines
00454     zline_tlines = []
00455     for z in -401.625, -133.875, 133.875, 401.625:
00456         if abscissa is None or abscissa[0] < z < abscissa[1]:
00457             zline_tlines.append(ROOT.TLine(z, -window, z, window))
00458             zline_tlines[-1].SetLineStyle(2)
00459             zline_tlines[-1].Draw()
00460     zline_labels = []
00461     zline_labels.append(ROOT.TText(-550, -0.9*window, "-2"))
00462     zline_labels.append(ROOT.TText(-300, -0.9*window, "-1"))
00463     zline_labels.append(ROOT.TText(-10, -0.9*window, "0"))
00464     zline_labels.append(ROOT.TText(250, -0.9*window, "+1"))
00465     zline_labels.append(ROOT.TText(500, -0.9*window, "+2"))
00466     for z in zline_labels: z.Draw()
00467     zline_labels.append(ROOT.TText(-600, -0.75*window, "Wheel:")); zline_labels[-1].Draw()
00468 
00469 def rlines(disk, window, abscissa):
00470     global rline_tlines
00471     rline_tlines = []
00472     if disk == 1: rl = [150., 270., 480.]
00473     else: rl = [350.]
00474     for r in rl:
00475         if abscissa is None or abscissa[0] < r < abscissa[1]:
00476             rline_tlines.append(ROOT.TLine(r, -window, r, window))
00477             rline_tlines[-1].SetLineStyle(2)
00478             rline_tlines[-1].Draw()
00479 
00480 ######################################################################################################
00481 
00482 def DBMC(database, reports, window=10., selection=None, phi=False, 
00483          color=ROOT.kBlue-8, style=1, bins=50, normalized=False, getvalues=False, name=""):
00484     return DBdiff(database, None, reports, None, window, selection, phi, color, style, bins, normalized, getvalues, name)
00485 
00486 def DBdiff(database1, database2, reports1, reports2, window=10., selection=None, phi=False, color=ROOT.kBlue-8,
00487            style=1, bins=50, normalized=False, getvalues=False, name="tmp"):
00488     tdrStyle.SetOptStat("emrou")
00489     tdrStyle.SetStatW(0.40)
00490 
00491     if phi:
00492         hx = ROOT.TH1F("%s_phi" % name, "", bins, -window, window)
00493     else:
00494         hx = ROOT.TH1F("%s_x" % name, "", bins, -window, window)
00495     hy = ROOT.TH1F("%s_y" % name, "", bins, -window, window)
00496     hz = ROOT.TH1F("%s_z" % name, "", bins, -window, window)
00497     hphix = ROOT.TH1F("%s_phix" % name, "", bins, -window, window)
00498     hphiy = ROOT.TH1F("%s_phiy" % name, "", bins, -window, window)
00499     hphiz = ROOT.TH1F("%s_phiz" % name, "", bins, -window, window)
00500         
00501     for r1 in reports1:
00502         if selection is None or (selection.func_code.co_argcount == len(r1.postal_address) and selection(*r1.postal_address)):
00503             if reports2 is None:
00504                 r2 = Report(r1.chamberId, r1.postal_address, r1.name)
00505                 r2.add_parameters(ValErr(0., 0., 0.), ValErr(0., 0., 0.), ValErr(0., 0., 0.), 
00506                     ValErr(0., 0., 0.), ValErr(0., 0., 0.), ValErr(0., 0., 0.), 0., 0., 0., 0.)
00507             else:
00508                 found = False
00509                 for r2 in reports2:
00510                     if r1.postal_address == r2.postal_address:
00511                         found = True
00512                         break
00513                 if not found: continue
00514 
00515             found = False
00516             if r1.postal_address[0] == "DT":
00517                 if r1.postal_address[1:] in database1.dt:
00518                     found = True
00519                     db1 = database1.dt[r1.postal_address[1:]]
00520                     if database2 is None:
00521                         db2 = DTAlignable()
00522                         db2.x = db2.y = db2.z = db2.phix = db2.phiy = db2.phiz = 0.
00523                         db2.xx = db2.xy = db2.xz = db2.yx = db2.yy = db2.yz = db2.zx = db2.zy = db2.zz = 0.
00524                     else:
00525                         db2 = database2.dt[r1.postal_address[1:]]
00526             else:
00527                 if r1.postal_address[1:] in database1.csc:
00528                     found = True
00529                     db1 = database1.csc[r1.postal_address[1:]]
00530                     if database2 is None:
00531                         db2 = CSCAlignable()
00532                         db2.x = db2.y = db2.z = db2.phix = db2.phiy = db2.phiz = 0.
00533                         db2.xx = db2.xy = db2.xz = db2.yx = db2.yy = db2.yz = db2.zx = db2.zy = db2.zz = 0.
00534                     else:
00535                         db2 = database2.csc[r1.postal_address[1:]]
00536 
00537             if found and r1.status == "PASS" and r2.status == "PASS":
00538                 if r1.deltax is not None and r2.deltax is not None and r1.deltax.error is not None and \
00539                    r2.deltax.error is not None and (r1.deltax.error**2 + r2.deltax.error**2) > 0.:
00540                     if normalized:
00541                         fill = (db1.x - db2.x)/sqrt(r1.deltax.error**2 + r2.deltax.error**2) * signConventions[r1.postal_address][0]
00542                     else:
00543                         if phi:
00544                             fill = (db1.x - db2.x)/signConventions[r1.postal_address][3] * 1000. * signConventions[r1.postal_address][0]
00545                         else:
00546                             fill = (db1.x - db2.x) * 10. * signConventions[r1.postal_address][0]
00547                     hx.Fill(fill)
00548                     if getvalues not in (False, None):
00549                         getvalues["x"].append((fill, 10. * sqrt(r1.deltax.error**2 + r2.deltax.error**2)))
00550 
00551                 if r1.deltay is not None and r2.deltay is not None and r1.deltay.error is not None and \
00552                    r2.deltay.error is not None and (r1.deltay.error**2 + r2.deltay.error**2) > 0.:
00553                     if normalized:
00554                         fill = (db1.y - db2.y)/sqrt(r1.deltay.error**2 + r2.deltay.error**2) * signConventions[r1.postal_address][1]
00555                     else:
00556                         fill = (db1.y - db2.y) * 10. * signConventions[r1.postal_address][1]
00557                     hy.Fill(fill)
00558                     if getvalues not in (False, None):
00559                         getvalues["y"].append((fill, 10. * sqrt(r1.deltay.error**2 + r2.deltay.error**2)))
00560 
00561                 if r1.deltaz is not None and r2.deltaz is not None and r1.deltaz.error is not None and \
00562                    r2.deltaz.error is not None and (r1.deltaz.error**2 + r2.deltaz.error**2) > 0.:
00563                     if normalized:
00564                         fill = (db1.z - db2.z)/sqrt(r1.deltaz.error**2 + r2.deltaz.error**2) * signConventions[r1.postal_address][2]
00565                     else:
00566                         fill = (db1.z - db2.z) * 10. * signConventions[r1.postal_address][2]
00567                     hz.Fill(fill)
00568                     if getvalues not in (False, None):
00569                         getvalues["z"].append((fill, 10. * sqrt(r1.deltaz.error**2 + r2.deltaz.error**2)))
00570 
00571                 if r1.deltaphix is not None and r2.deltaphix is not None and r1.deltaphix.error is not None and \
00572                    r2.deltaphix.error is not None and (r1.deltaphix.error**2 + r2.deltaphix.error**2) > 0.:
00573                     if normalized:
00574                         fill = (db1.phix - db2.phix)/sqrt(r1.deltaphix.error**2 + r2.deltaphix.error**2)
00575                     else:
00576                         fill = (db1.phix - db2.phix) * 1000.
00577                     hphix.Fill(fill)
00578                     if getvalues not in (False, None):
00579                         getvalues["phix"].append((fill, 10. * sqrt(r1.deltaphix.error**2 + r2.deltaphix.error**2)))
00580 
00581                 if r1.deltaphiy is not None and r2.deltaphiy is not None and r1.deltaphiy.error is not None and \
00582                    r2.deltaphiy.error is not None and (r1.deltaphiy.error**2 + r2.deltaphiy.error**2) > 0.:
00583                     if normalized:
00584                         fill = (db1.phiy - db2.phiy)/sqrt(r1.deltaphiy.error**2 + r2.deltaphiy.error**2)
00585                     else:
00586                         fill = (db1.phiy - db2.phiy) * 1000.
00587                     hphiy.Fill(fill)
00588                     if getvalues not in (False, None):
00589                         getvalues["phiy"].append((fill, 10. * sqrt(r1.deltaphiy.error**2 + r2.deltaphiy.error**2)))
00590 
00591                 if r1.deltaphiz is not None and r2.deltaphiz is not None and r1.deltaphiz.error is not None and \
00592                    r2.deltaphiz.error is not None and (r1.deltaphiz.error**2 + r2.deltaphiz.error**2) > 0.:
00593                     if normalized:
00594                         fill = (db1.phiz - db2.phiz)/sqrt(r1.deltaphiz.error**2 + r2.deltaphiz.error**2)
00595                     else:
00596                         fill = (db1.phiz - db2.phiz) * 1000.
00597                     hphiz.Fill(fill)
00598                     if getvalues not in (False, None):
00599                         getvalues["phiz"].append((fill, 10. * sqrt(r1.deltaphiz.error**2 + r2.deltaphiz.error**2)))
00600 
00601     if not normalized:
00602         if phi:
00603             hx.SetXTitle("#delta_{#phi} position (mrad)")
00604         else:
00605             hx.SetXTitle("#delta_{x'} (mm)")
00606         hy.SetXTitle("#delta_{y'} (mm)")
00607         hz.SetXTitle("#delta_{z'} (mm)")
00608         hphix.SetXTitle("#delta_{#phi_{x}} (mrad)")
00609         hphiy.SetXTitle("#delta_{#phi_{y}} (mrad)")
00610         hphiz.SetXTitle("#delta_{#phi_{z}} (mrad)")
00611     else:
00612         if phi:
00613             hx.SetXTitle("#delta_{#phi}/#sigma_{#phi} position")
00614         else:
00615             hx.SetXTitle("#delta_{x'}/#sigma_{x'}")
00616         hy.SetXTitle("#delta_{y'}/#sigma_{y'}")
00617         hz.SetXTitle("#delta_{z'}/#sigma_{z'}")
00618         hphix.SetXTitle("#delta_{#phi_{x}}/#sigma_{#phi_{x}}")
00619         hphiy.SetXTitle("#delta_{#phi_{y}}/#sigma_{#phi_{y}}")
00620         hphiz.SetXTitle("#delta_{#phi_{z}}/#sigma_{#phi_{z}}")
00621 
00622     for h in hx, hy, hz, hphix, hphiy, hphiz:
00623         h.GetXaxis().CenterTitle()
00624         h.GetYaxis().CenterTitle()
00625         h.SetFillColor(color)
00626         h.SetLineStyle(style)
00627 
00628     if normalized:
00629         fx = ROOT.TF1("fx", "%g * exp(-x**2/2.)/sqrt(2.*3.1415926)" % (hx.GetEntries()*2.*window/bins), -window, window)
00630         fy = ROOT.TF1("fy", "%g * exp(-x**2/2.)/sqrt(2.*3.1415926)" % (hy.GetEntries()*2.*window/bins), -window, window)
00631         fz = ROOT.TF1("fz", "%g * exp(-x**2/2.)/sqrt(2.*3.1415926)" % (hz.GetEntries()*2.*window/bins), -window, window)
00632         fphix = ROOT.TF1("fphix", "%g * exp(-x**2/2.)/sqrt(2.*3.1415926)" % (hphix.GetEntries()*2.*window/bins), -window, window)
00633         fphiy = ROOT.TF1("fphiy", "%g * exp(-x**2/2.)/sqrt(2.*3.1415926)" % (hphiy.GetEntries()*2.*window/bins), -window, window)
00634         fphiz = ROOT.TF1("fphiz", "%g * exp(-x**2/2.)/sqrt(2.*3.1415926)" % (hphiz.GetEntries()*2.*window/bins), -window, window)
00635         for f in fx, fy, fz, fphix, fphiy, fphiz:
00636             f.SetLineWidth(2)
00637             f.SetLineColor(ROOT.kBlue)
00638         for h, f in (hx, fx), (hy, fy), (hz, fz), (hphix, fphix), (hphiy, fphiy), (hphiz, fphiz):
00639             h.SetAxisRange(0, 1.1*max(h.GetMaximum(), f.GetMaximum()), "Y")
00640         c1.Clear()
00641         c1.Divide(3, 2)
00642         c1.GetPad(1).cd(); hx.Draw(); fx.Draw("same")
00643         c1.GetPad(2).cd(); hy.Draw(); fy.Draw("same")
00644         c1.GetPad(3).cd(); hz.Draw(); fz.Draw("same")
00645         c1.GetPad(4).cd(); hphix.Draw(); fphix.Draw("same")
00646         c1.GetPad(5).cd(); hphiy.Draw(); fphiy.Draw("same")
00647         c1.GetPad(6).cd(); hphiz.Draw(); fphiz.Draw("same")
00648         return hx, hy, hz, hphix, hphiy, hphiz, fx, fy, fz, fphix, fphiy, fphiz
00649     else:
00650         c1.Clear()
00651         c1.Divide(3, 2)
00652         c1.GetPad(1).cd(); hx.Draw()
00653         c1.GetPad(2).cd(); hy.Draw()
00654         c1.GetPad(3).cd(); hz.Draw()
00655         c1.GetPad(4).cd(); hphix.Draw()
00656         c1.GetPad(5).cd(); hphiy.Draw()
00657         c1.GetPad(6).cd(); hphiz.Draw()
00658         return hx, hy, hz, hphix, hphiy, hphiz
00659 
00660 def DBMCVersus(quantity, versus, database, reports, window=10., selection=None, color=ROOT.kBlack):
00661     return DBdiffVersus(quantity, versus, database, None, reports, None, window, selection, color)
00662 
00663 def DBdiffVersus(quantity, versus, database1, database2, reports1, reports2, windwselection=None, color=ROOT.kBlack):
00664     tdrStyle.SetOptStat("")
00665 
00666     domain = []
00667     values = []
00668     errors = []
00669         
00670     for r1 in reports1:
00671         if selection is None or (selection.func_code.co_argcount == len(r1.postal_address) and selection(*r1.postal_address)):
00672             if reports2 is None:
00673                 r2 = Report(r1.chamberId, r1.postal_address, r1.name)
00674                 r2.add_parameters(ValErr(0., 0., 0.), ValErr(0., 0., 0.), ValErr(0., 0., 0.), 
00675                                   ValErr(0., 0., 0.), ValErr(0., 0., 0.), ValErr(0., 0., 0.), 0., 0., 0.)
00676             else:
00677                 found = False
00678                 for r2 in reports2:
00679                     if r1.postal_address == r2.postal_address:
00680                         found = True
00681                         break
00682                 if not found: continue
00683 
00684             found = False
00685             if r1.postal_address[0] == "DT":
00686                 if r1.postal_address[1:] in database1.dt:
00687                     found = True
00688                     db1 = database1.dt[r1.postal_address[1:]]
00689                     if database2 is None:
00690                         db2 = DTAlignable()
00691                         db2.x = db2.y = db2.z = db2.phix = db2.phiy = db2.phiz = 0.
00692                         db2.xx = db2.xy = db2.xz = db2.yx = db2.yy = db2.yz = db2.zx = db2.zy = db2.zz = 0.
00693                     else:
00694                         db2 = database2.dt[r1.postal_address[1:]]
00695             else:
00696                 if r1.postal_address[1:] in database1.csc:
00697                     found = True
00698                     db1 = database1.csc[r1.postal_address[1:]]
00699                     if database2 is None:
00700                         db2 = CSCAlignable()
00701                         db2.x = db2.y = db2.z = db2.phix = db2.phiy = db2.phiz = 0.
00702                         db2.xx = db2.xy = db2.xz = db2.yx = db2.yy = db2.yz = db2.zx = db2.zy = db2.zz = 0.
00703                     else:
00704                         db2 = database2.csc[r1.postal_address[1:]]
00705 
00706             if found and r1.status == "PASS" and r2.status == "PASS":
00707                 okay = False
00708 
00709                 if quantity == "phi":
00710                     if r1.deltax is not None and r2.deltax is not None and r1.deltax.error is not None and \
00711                        r2.deltax.error is not None and (r1.deltax.error**2 + r2.deltax.error**2) > 0.:
00712                         okay = True
00713                         values.append((db1.x - db2.x)/
00714                                       signConventions[r1.postal_address][3] * 1000. * signConventions[r1.postal_address][0])
00715                         errors.append((r1.deltax.error**2 + r2.deltax.error**2)/
00716                                       signConventions[r1.postal_address][3] * 1000. * signConventions[r1.postal_address][0])
00717 
00718                 elif quantity == "x":
00719                     if r1.deltax is not None and r2.deltax is not None and r1.deltax.error is not None and \
00720                        r2.deltax.error is not None and (r1.deltax.error**2 + r2.deltax.error**2) > 0.:
00721                         okay = True
00722                         values.append((db1.x - db2.x) * 10. * signConventions[r1.postal_address][0])
00723                         errors.append((r1.deltax.error**2 + r2.deltax.error**2) * 10. * signConventions[r1.postal_address][0])
00724 
00725                 elif quantity == "y":
00726                     if r1.deltay is not None and r2.deltay is not None and r1.deltay.error is not None and \
00727                        r2.deltay.error is not None and (r1.deltay.error**2 + r2.deltay.error**2) > 0.:
00728                         okay = True
00729                         values.append((db1.y - db2.y) * 10. * signConventions[r1.postal_address][1])
00730                         errors.append((r1.deltay.error**2 + r2.deltay.error**2) * 10. * signConventions[r1.postal_address][1])
00731 
00732                 elif quantity == "z":
00733                     if r1.deltaz is not None and r2.deltaz is not None and r1.deltaz.error is not None and \
00734                        r2.deltaz.error is not None and (r1.deltaz.error**2 + r2.deltaz.error**2) > 0.:
00735                         okay = True
00736                         values.append((db1.z - db2.z) * 10. * signConventions[r1.postal_address][2])
00737                         errors.append((r1.deltaz.error**2 + r2.deltaz.error**2) * 10. * signConventions[r1.postal_address][2])
00738 
00739                 elif quantity == "phix":
00740                     if r1.deltaphix is not None and r2.deltaphix is not None and r1.deltaphix.error is not None and \
00741                        r2.deltaphix.error is not None and (r1.deltaphix.error**2 + r2.deltaphix.error**2) > 0.:
00742                         okay = True
00743                         values.append((db1.phix - db2.phix) * 1000.)
00744                         errors.append((r1.deltaphix.error**2 + r2.deltaphix.error**2) * 1000.)
00745 
00746                 elif quantity == "phiy":
00747                     if r1.deltaphiy is not None and r2.deltaphiy is not None and r1.deltaphiy.error is not None and \
00748                        r2.deltaphiy.error is not None and (r1.deltaphiy.error**2 + r2.deltaphiy.error**2) > 0.:
00749                         okay = True
00750                         values.append((db1.phiy - db2.phiy) * 1000.)
00751                         errors.append((r1.deltaphiy.error**2 + r2.deltaphiy.error**2) * 1000.)
00752 
00753                 elif quantity == "phiz":
00754                     if r1.deltaphiz is not None and r2.deltaphiz is not None and r1.deltaphiz.error is not None and \
00755                        r2.deltaphiz.error is not None and (r1.deltaphiz.error**2 + r2.deltaphiz.error**2) > 0.:
00756                         okay = True
00757                         values.append((db1.phiz - db2.phiz) * 1000.)
00758                         errors.append((r1.deltaphiz.error**2 + r2.deltaphiz.error**2) * 1000.)
00759 
00760                 else: raise Exception
00761 
00762                 if okay:
00763                     if versus == "r": domain.append(signConventions[r1.postal_address][3])
00764                     elif versus == "phi": domain.append(signConventions[r1.postal_address][4])
00765                     elif versus == "z": domain.append(signConventions[r1.postal_address][5])
00766                     else: raise Exception
00767 
00768     if versus == "r":
00769         bkgndhist = ROOT.TH1F("bkgndhist", "", 100, 0., 800.)
00770         bkgndhist.SetXTitle("R (cm)")
00771     elif versus == "phi":
00772         bkgndhist = ROOT.TH1F("bkgndhist", "", 100, -pi, pi)
00773         bkgndhist.SetXTitle("#phi (rad)")
00774     elif versus == "z":
00775         bkgndhist = ROOT.TH1F("bkgndhist", "", 100, -1100., 1100.)
00776         bkgndhist.SetXTitle("z (cm)")
00777     bkgndhist.GetXaxis().CenterTitle()
00778 
00779     bkgndhist.SetAxisRange(-window, window, "Y")
00780     if quantity == "phi": bkgndhist.SetYTitle("#delta_{#phi} position (mrad)")
00781     elif quantity == "x": bkgndhist.SetYTitle("#delta_{x'} (mm)")
00782     elif quantity == "y": bkgndhist.SetYTitle("#delta_{y'} (mm)")
00783     elif quantity == "z": bkgndhist.SetYTitle("#delta_{z'} (mm)")
00784     elif quantity == "phix": bkgndhist.SetYTitle("#delta_{#phi_{x}} (mrad)")
00785     elif quantity == "phiy": bkgndhist.SetYTitle("#delta_{#phi_{y}} (mrad)")
00786     elif quantity == "phiz": bkgndhist.SetYTitle("#delta_{#phi_{z}} (mrad)")
00787     else: raise Exception
00788     bkgndhist.GetYaxis().CenterTitle()
00789 
00790     if len(domain) == 0:
00791         tgraph = ROOT.TGraphErrors(0)
00792     else:
00793         tgraph = ROOT.TGraphErrors(len(domain), array.array("d", domain), array.array("d", values), 
00794                                                 array.array("d", [0.]*len(domain)), array.array("d", errors))
00795     tgraph.SetMarkerColor(color)
00796     tgraph.SetLineColor(color)
00797 
00798     bkgndhist.Draw()
00799     if tgraph.GetN() > 0: tgraph.Draw("p")
00800     return bkgndhist, tgraph, domain, values, errors
00801 
00802 ######################################################################################################
00803 
00804 def idToPostalAddress(id):
00805   # only len==9 ids can correspond to valid postal address
00806   if len(id)!=9: return None
00807   if id[0:2]=="MB":
00808     #print id
00809     pa = ("DT", int(id[2:4]), int(id[5]), int(id[7:9]))
00810     #print pa
00811     if pa[1]<-2 or pa[1]>2: return None
00812     if pa[2]>4: return None
00813     if pa[3]<1 or pa[3]>14 or (pa[3]==4 and pa[3]>12): return None
00814     return pa
00815   elif id[0:2]=="ME":
00816     if id[2]=="+": ec=1
00817     elif id[2]=="-": ec=2
00818     else: return None
00819     pa = ("CSC", ec, int(id[3]), int(id[5]), int(id[7:9]))
00820     if pa[2]<1 or pa[2]>4: return None
00821     if pa[3]<1 or pa[3]>4 or (pa[2]>1 and pa[3]>2): return None
00822     if pa[4]<1 or pa[4]>36 or (pa[2]>1 and pa[3]==1 and pa[4]>18): return None
00823     return pa
00824   else: return None
00825 
00826 
00827 def postalAddressToId(postal_address):
00828   if postal_address[0] == "DT":
00829     wheel, station, sector = postal_address[1:]
00830     w = "%+d"%wheel
00831     if w=="+0": w = "-0"
00832     return "MB%s/%d/%02d" % (w, station, sector)
00833   elif postal_address[0] == "CSC":
00834     endcap, station, ring, chamber = postal_address[1:]
00835     if endcap != 1: station = -1 * abs(station)
00836     return "ME%+d/%d/%02d" % (station, ring, chamber)
00837 
00838 
00839 def nameToId(name):
00840   if name[0:2] == "MB":
00841     wh = name[4]
00842     if   wh == "A": w = "-2"
00843     elif wh == "B": w = "-1"
00844     elif wh == "C": w = "-0"
00845     elif wh == "D": w = "+1"
00846     elif wh == "E": w = "+2"
00847     else: return ""
00848     station = name[7]
00849     sector = name[11:13]
00850     return "MB%s/%s/%s" % (w, station, sector)
00851   elif name[0:2] == "ME":
00852     if name[2]=="p": endcap = "+"
00853     elif name[2]=="m": endcap = "-"
00854     else: return ""
00855     station = name[3]
00856     ring = name[4]
00857     chamber = name[6:8]
00858     return "ME%s%s/%s/%s" % (endcap, station, ring, chamber)
00859   return None
00860 
00861 
00862 def availableCellsDT(reports):
00863   dts = []
00864   # DT wheels
00865   for iwheel in DT_TYPES:
00866     if iwheel[1]=="ALL": continue
00867     dts.append(iwheel[0])
00868   # DT wheel & station
00869   for wheel in DT_TYPES:
00870     if wheel[1]=="ALL": continue
00871     for station in wheel[2]:
00872       dts.append(wheel[0]+'/'+station[1])
00873   # DT station & sector 
00874   for wheel in DT_TYPES:
00875     if wheel[1]!="ALL": continue
00876     for station in wheel[2]:
00877       for sector in range(1,station[2]+1):
00878         ssector = "%02d" % sector
00879         dts.append(wheel[0]+'/'+station[1]+'/'+ssector)
00880   # DT chambers
00881   for wheel in DT_TYPES:
00882     if wheel[1]=="ALL": continue
00883     for station in wheel[2]:
00884       for sector in range(1,station[2]+1):
00885         ssector = "%02d" % sector
00886         label = "MBwh%sst%ssec%s" % (wheelLetter(int(wheel[1])),station[1],ssector)
00887         if len(reports)==0:
00888           # no reports case: do not include chambers 
00889           #dts.append(wheel[0]+'/'+station[1]+'/'+ssector)
00890           continue
00891         found = False
00892         for r in reports:
00893           if r.name == label:
00894             found = True
00895             break
00896         if not found: continue
00897         if r.status == "TOOFEWHITS" and r.posNum+r.negNum==0: continue
00898         if r.status == "NOFIT": continue
00899         dts.append(wheel[0]+'/'+station[1]+'/'+ssector)
00900   return dts
00901 
00902 
00903 def availableCellsCSC(reports):
00904   cscs = []
00905   # CSC station & ring 
00906   for endcap in CSC_TYPES:
00907     for station in endcap[2]:
00908       for ring in station[2]:
00909         if ring[1]=="ALL": continue
00910         #label = "CSCvsphi_me%s%s%s" % (endcap[1], station[1], ring[1])
00911         cscs.append("%s%s/%s" % (endcap[0], station[1],ring[1]))
00912   # CSC station and chamber
00913   for endcap in CSC_TYPES:
00914     for station in endcap[2]:
00915       for ring in station[2]:
00916         if ring[1]!="ALL": continue
00917         for chamber in range(1,ring[2]+1):
00918           #label = "CSCvsr_me%s%sch%02d" % (endcap[1], station[1], chamber)
00919           cscs.append("%s%s/ALL/%02d" % (endcap[0], station[1],chamber))
00920   # CSC chambers
00921   for endcap in CSC_TYPES:
00922     for station in endcap[2]:
00923       for ring in station[2]:
00924         if ring[1]=="ALL": continue
00925         for chamber in range(1,ring[2]+1):
00926           # exclude non instrumented ME4/2 
00927           if station[1]=="4" and ring[1]=="2":
00928             if endcap[1]=="m": continue
00929             if chamber<9 or chamber>13: continue
00930           schamber = "%02d" % chamber
00931           label = "ME%s%s%s_%s" % (endcap[1], station[1], ring[1], schamber)
00932           if len(reports)==0:
00933             # no reports case: do not include chambers 
00934             #cscs.append(endcap[0]+station[1]+'/'+ring[1]+'/'+schamber)
00935             continue
00936           found = False
00937           for r in reports:
00938             if r.name == label:
00939               found = True
00940               break
00941           if not found: continue
00942           if r.status == "TOOFEWHITS" and r.posNum+r.negNum==0: continue
00943           if r.status == "NOFIT": continue
00944           cscs.append(endcap[0]+station[1]+'/'+ring[1]+'/'+schamber)
00945   return cscs
00946 
00947 
00948 DQM_SEVERITY = [
00949   {"idx":0, "name": "NONE", "color": "lightgreen", "hex":"#90EE90"},
00950   {"idx":1, "name": "LOWSTAT05", "color": "lightgreen", "hex":"#96D953"},
00951   {"idx":2, "name": "LOWSTAT075", "color": "lightgreen", "hex":"#94E26F"},
00952   {"idx":3, "name": "LOWSTAT1", "color": "yellowgreen", "hex":"#9ACD32"},
00953   {"idx":4, "name": "LOWSTAT", "color": "yellow", "hex":"#FFFF00"},
00954   {"idx":5, "name": "TOLERABLE", "color": "lightpink", "hex":"#FFB6C1"},
00955   {"idx":6, "name": "SEVERE", "color": "orange", "hex":"#FFA500"},
00956   {"idx":7, "name": "CRITICAL", "color": "red", "hex":"#FF0000"}];
00957 
00958 
00959 def addToTestResults(c,res):
00960   if len(res)>0:
00961     if TEST_RESULTS.has_key(c): TEST_RESULTS[c].extend(res)
00962     else: TEST_RESULTS[c] = res
00963 
00964 
00965 def testEntry(testID,scope,descr,severity):
00966   s = 0
00967   for sev in DQM_SEVERITY:
00968     if sev["name"]==severity: s = sev["idx"]
00969   return {"testID":testID,"scope":scope,"descr":descr,"severity":s}
00970 
00971 
00972 def testZeroWithin5Sigma(x):
00973   if abs(x[1])==0.: return 0.
00974   pull = abs(x[0])/abs(x[1])
00975   if pull <= 5: return 0.
00976   else: return pull
00977 
00978 
00979 def testDeltaWithin5Sigma(x,sx):
00980   n = len(x)
00981   res = []
00982   dr = []
00983   #print x
00984   #print sx
00985   for i in range(1,n+1):
00986     x1 = x[i-1]
00987     sx1 = sx[i-1]
00988     x2 = x[0]
00989     sx2 = sx[0]
00990     if i<n: 
00991       x2 = x[i]
00992       sx2 = sx[i]
00993     sig1 = sqrt( (sx1[0]-sx1[1])**2 + x1[1]**2 )
00994     sig2 = sqrt( (sx2[0]-sx2[1])**2 + x2[1]**2 )
00995     df = abs(x1[0]-x2[0]) - 3*( sig1 + sig2 )
00996     #df = abs(sx1[1]-sx2[0]) - 5*(abs(x1[1]) + abs(x2[1]))
00997     #print i, df, '= abs(',sx1[1],'-',sx2[0],')-5*(abs(',x1[1],')+abs(',x2[1],'))'
00998     dr.append(df)
00999     if df > 0: res.append(i)
01000   #print dr
01001   #print res
01002   return res
01003 
01004 
01005 def doTestsForReport(cells,reports):
01006   for c in cells:
01007     # can a cell be converted to a chamber postal address?
01008     postal_address = idToPostalAddress(c)
01009     if not postal_address: continue
01010     
01011     # is this chamber in _report?
01012     found = False
01013     for r in reports:
01014       if r.postal_address == postal_address:
01015         found = True
01016         break
01017     if not found: continue
01018 
01019     # chamber's tests result
01020     res = []
01021     scope = postal_address[0]
01022     
01023     # noting could be done if fitting fails
01024     if r.status == "FAIL" or r.status == "MINUITFAIL":
01025       res.append(testEntry("FAILURE",scope,r.status+" failure","CRITICAL"))
01026       addToTestResults(c,res)
01027       continue
01028 
01029     # noting could be done if TOOFEWHITS
01030     nseg = r.posNum + r.negNum
01031     if r.status == "TOOFEWHITS" and nseg>0:
01032       res.append(testEntry("LOW_STAT",scope,"low stat, #segments=%d"%nseg,"LOWSTAT"))
01033       addToTestResults(c,res)
01034       continue
01035 
01036     # set shades of light green according to sidma(dx)
01037     sdx = 10.*r.deltax.error
01038     if sdx>0.5:
01039       if sdx<0.75: res.append(testEntry("LOW_STAT_DDX05",scope,"low stat, delta(dx)=%f #segments=%d" % (sdx,nseg),"LOWSTAT05"))
01040       elif sdx<1.: res.append(testEntry("LOW_STAT_DDX075",scope,"low stat, delta(dx)=%f #segments=%d" % (sdx,nseg),"LOWSTAT075"))
01041       else: res.append(testEntry("LOW_STAT_DDX1",scope,"low stat, delta(dx)=%f #segments=%d" % (sdx,nseg),"LOWSTAT1"))
01042 
01043     # check chi2
01044     if r.redchi2 > 2.5:
01045       res.append(testEntry("BIG_CHI2",scope,"chi2=%f>2.5" % r.redchi2,"TOLERABLE"))
01046 
01047     # check medians
01048     medx, meddx = 10.*r.median_x, 1000.*r.median_dxdz
01049     #medy, meddy = 10.*r.median_y, 1000.*r.median_dydz
01050     if medx>2:  res.append(testEntry("BIG_MED_X",scope,"median dx=%f>2 mm"%medx,"SEVERE"))
01051     #if medy>3: res.append(testEntry("BIG_MED_Y",scope,"median dy=%f>3 mm"%medy,"SEVERE"))
01052     if meddx>2: res.append(testEntry("BIG_MED_DXDZ",scope,"median d(dx/dz)=%f>2 mrad"%meddx,"SEVERE"))
01053     #if meddy>3: res.append(testEntry("BIG_MED_DYDZ",scope,"median d(dy/dz)=%f>3 mrad"%meddy,"SEVERE"))
01054 
01055     # check residuals far from zero
01056     isDTst4 = False
01057     if postal_address[0] == "DT" and postal_address[2]==4: isDTst4 = True
01058     dx, dy, dpy, dpz = 10.*r.deltax.value, 0., 1000.*r.deltaphiy.value, 1000.*r.deltaphiz.value
01059     if not isDTst4: dy = 10.*r.deltay.value
01060     if dx>0.2:   res.append(testEntry("BIG_LAST_ITR_DX",scope,"dx=%f>0.2 mm"%dx,"CRITICAL"))
01061     if dy>0.2:   res.append(testEntry("BIG_LAST_ITR_DY",scope,"dy=%f>0.2 mm"%dy,"CRITICAL"))
01062     if dpy>0.2:   res.append(testEntry("BIG_LAST_ITR_DPHIY",scope,"dphiy=%f>0.2 mrad"%dpy,"CRITICAL"))
01063     if dpz>0.2:   res.append(testEntry("BIG_LAST_ITR_DPHIZ",scope,"dphiz=%f>0.2 mrad"%dpz,"CRITICAL"))
01064     #if ddx>0.03: res.append(testEntry("BIG_DX",scope,"dphix=%f>0.03 mrad"%ddx,"CRITICAL"))
01065 
01066     addToTestResults(c,res)
01067 
01068 
01069 def doTestsForMapPlots(cells):
01070   for c in cells:
01071     res = []
01072     
01073     scope = "zzz"
01074     if c[0:2]=="MB": scope = "DT"
01075     if c[0:2]=="ME": scope = "CSC"
01076     if scope == "zzz":
01077       print "strange cell ID: ", c
01078       return None
01079     
01080     if MAP_RESULTS_FITSIN.has_key(c):
01081       t = MAP_RESULTS_FITSIN[c]
01082       t_a = testZeroWithin5Sigma(t['a'])
01083       t_s = testZeroWithin5Sigma(t['sin'])
01084       t_c = testZeroWithin5Sigma(t['cos'])
01085       if t_a+t_s+t_c >0:
01086         descr = "map fitsin 5 sigma away from 0; pulls : a=%.2f sin=%.2f, cos=%.2f" % (t_a,t_s,t_c)
01087         res.append(testEntry("MAP_FITSIN",scope,descr,"SEVERE"))
01088     
01089     if MAP_RESULTS_SAWTOOTH.has_key(c):
01090       t = MAP_RESULTS_SAWTOOTH[c]
01091       
01092       t_a = testDeltaWithin5Sigma(t['a'],t['da'])
01093       if len(t_a)>0:
01094         descr = "map discontinuities: %s" % ",".join(map(str,t_a))
01095         res.append(testEntry("MAP_DISCONTIN",scope,descr,"SEVERE"))
01096 
01097       t_b = map(testZeroWithin5Sigma, t['b'])
01098       t_bi = []
01099       for i in range(0,len(t_b)):
01100         if t_b[i]>0: t_bi.append(i+1)
01101       if len(t_bi)>0:
01102         descr = "map sawteeth: %s" % ",".join(map(str,t_bi))
01103         res.append(testEntry("MAP_SAWTEETH",scope,descr,"TOLERABLE"))
01104 
01105     addToTestResults(c,res)
01106 
01107 
01108 def saveTestResultsMap(run_name):
01109   if len(MAP_RESULTS_SAWTOOTH)+len(MAP_RESULTS_FITSIN)==0: return None
01110   ff = open("tmp_test_results_map__%s.pkl" % run_name, "wb")
01111   pickle.dump(MAP_RESULTS_SAWTOOTH, ff)
01112   pickle.dump(MAP_RESULTS_FITSIN, ff)
01113   ff.close()
01114 
01115 
01116 def loadTestResultsMap(run_name):
01117   print "tmp_test_results_map__%s.pkl" % run_name, os.access("tmp_test_results_map__%s.pkl" % run_name,os.F_OK)
01118   if not os.access("tmp_test_results_map__%s.pkl" % run_name,os.F_OK): return None
01119   global MAP_RESULTS_FITSIN, MAP_RESULTS_SAWTOOTH
01120   ff = open("tmp_test_results_map__%s.pkl" % run_name, "rb")
01121   MAP_RESULTS_SAWTOOTH = pickle.load(ff)
01122   MAP_RESULTS_FITSIN = pickle.load(ff)
01123   ff.close()
01124   #execfile("tmp_test_results_map__%s.py" % run_name)
01125   #print 'asasas', MAP_RESULTS_FITSIN
01126   return True
01127 
01128 def writeDQMReport(fname_dqm, run_name):
01129   tests = []
01130   for c in TEST_RESULTS:
01131     tests.append({"objID":c, "name":c, "list":TEST_RESULTS[c]})
01132   lt = time.localtime(time.time())
01133   lts = "%04d-%02d-%02d %02d:%02d:%02d %s" % (lt[0], lt[1], lt[2], lt[3], lt[4], lt[5], time.tzname[1])
01134   dqm_report = {"run":run_name, "genDate": lts, "report":tests}
01135   ff = open(fname_dqm,mode="w")
01136   print >>ff, "var DQM_REPORT = "
01137   json.dump(dqm_report,ff)
01138   #print >>ff, "];"
01139   ff.close()
01140 
01141 
01142 def doTests(reports,fname_base, fname_dqm, run_name):
01143   # find available baseline
01144   dts = availableCellsDT(reports)
01145   cscs = availableCellsCSC(reports)
01146   mulist = ['Run: '+run_name,['ALL',['MU']],['DT',dts],['CSC',cscs]]
01147   ff = open(fname_base,mode="w")
01148   print >>ff, "var MU_LIST = ["
01149   json.dump(mulist,ff)
01150   print >>ff, "];"
01151   ff.close()
01152   
01153   doTestsForReport(dts,reports)
01154   doTestsForReport(cscs,reports)
01155   
01156   loadTestResultsMap(run_name)
01157   doTestsForMapPlots(dts)
01158   doTestsForMapPlots(cscs)
01159   
01160   writeDQMReport(fname_dqm, run_name)
01161 
01162 
01163 ######################################################################################################
01164 
01165 def plotmedians(reports1, reports2, selection=None, binsx=50, windowx=3., ceilingx=None, binsy=50, windowy=3., 
01166                 ceilingy=None, binsdxdz=50, windowdxdz=3., ceilingdxdz=None, binsdydz=50, windowdydz=3., ceilingdydz=None, 
01167                 r1text=" before", r2text=" after", which="median"):
01168     tdrStyle.SetOptStat("emrou")
01169     tdrStyle.SetStatW(0.40)
01170     tdrStyle.SetStatFontSize(0.05)
01171 
01172     global hmediandxdz_after, hmediandxdz_aftercopy, hmediandxdz_before, hmediandxdz_beforecopy, \
01173            hmediandydz_after, hmediandydz_aftercopy, hmediandydz_before, hmediandydz_beforecopy, \
01174            hmedianx_after, hmedianx_aftercopy, hmedianx_before, hmedianx_beforecopy, \
01175            hmediany_after, hmediany_aftercopy, hmediany_before, hmediany_beforecopy, tlegend 
01176 
01177     hmedianx_before = ROOT.TH1F("hmedianx_before", "", binsx, -windowx, windowx)
01178     hmediany_before = ROOT.TH1F("hmediany_before", "", binsy, -windowy, windowy)
01179     hmediandxdz_before = ROOT.TH1F("hmediandxdz_before", "", binsdxdz, -windowdxdz, windowdxdz)
01180     hmediandydz_before = ROOT.TH1F("hmediandydz_before", "", binsdydz, -windowdydz, windowdydz)
01181     hmedianx_after = ROOT.TH1F("hmedianx_after", "", binsx, -windowx, windowx)
01182     hmediany_after = ROOT.TH1F("hmediany_after", "", binsy, -windowy, windowy)
01183     hmediandxdz_after = ROOT.TH1F("hmediandxdz_after", "", binsdxdz, -windowdxdz, windowdxdz)
01184     hmediandydz_after = ROOT.TH1F("hmediandydz_after", "", binsdydz, -windowdydz, windowdydz)
01185 
01186     if which == "median":
01187         whichx = whichy = whichdxdz = whichdydz = "median"
01188     elif which == "bigmean":
01189         whichx = "mean30"
01190         whichy = "mean30"
01191         whichdxdz = "mean20"
01192         whichdydz = "mean50"
01193     elif which == "mean":
01194         whichx = "mean15"
01195         whichy = "mean15"
01196         whichdxdz = "mean10"
01197         whichdydz = "mean25"
01198     elif which == "bigwmean":
01199         whichx = "wmean30"
01200         whichy = "wmean30"
01201         whichdxdz = "wmean20"
01202         whichdydz = "wmean50"
01203     elif which == "wmean":
01204         whichx = "wmean15"
01205         whichy = "wmean15"
01206         whichdxdz = "wmean10"
01207         whichdydz = "wmean25"
01208     elif which == "bigstdev":
01209         whichx = "stdev30"
01210         whichy = "stdev30"
01211         whichdxdz = "stdev20"
01212         whichdydz = "stdev50"
01213     elif which == "stdev":
01214         whichx = "stdev15"
01215         whichy = "stdev15"
01216         whichdxdz = "stdev10"
01217         whichdydz = "stdev25"
01218     else:
01219         raise Exception, which + " not recognized"
01220 
01221     for r1 in reports1:
01222         if selection is None or (selection.func_code.co_argcount == len(r1.postal_address) and selection(*r1.postal_address)):
01223             found = False
01224             for r2 in reports2:
01225                 if r1.postal_address == r2.postal_address:
01226                     found = True
01227                     break
01228             if not found: continue
01229 
01230             if r1.status == "PASS" and r2.status == "PASS":
01231                 hmedianx_before.Fill(10.*eval("r1.%s_x" % whichx))
01232                 hmediandxdz_before.Fill(1000.*eval("r1.%s_dxdz" % whichdxdz))
01233                 hmedianx_after.Fill(10.*eval("r2.%s_x" % whichx))
01234                 hmediandxdz_after.Fill(1000.*eval("r2.%s_dxdz" % whichdxdz))
01235 
01236                 if r1.median_y is not None:
01237                     hmediany_before.Fill(10.*eval("r1.%s_y" % whichy))
01238                     hmediandydz_before.Fill(1000.*eval("r1.%s_dydz" % whichdydz))
01239                     hmediany_after.Fill(10.*eval("r2.%s_y" % whichy))
01240                     hmediandydz_after.Fill(1000.*eval("r2.%s_dydz" % whichdydz))
01241 
01242     hmedianx_beforecopy = hmedianx_before.Clone()
01243     hmediany_beforecopy = hmediany_before.Clone()
01244     hmediandxdz_beforecopy = hmediandxdz_before.Clone()
01245     hmediandydz_beforecopy = hmediandydz_before.Clone()
01246     hmedianx_beforecopy.SetLineStyle(2)
01247     hmediany_beforecopy.SetLineStyle(2)
01248     hmediandxdz_beforecopy.SetLineStyle(2)
01249     hmediandydz_beforecopy.SetLineStyle(2)
01250 
01251     hmedianx_before.SetFillColor(ROOT.kMagenta+2)
01252     hmediany_before.SetFillColor(ROOT.kMagenta+2)
01253     hmediandxdz_before.SetFillColor(ROOT.kMagenta+2)
01254     hmediandydz_before.SetFillColor(ROOT.kMagenta+2)
01255     hmedianx_after.SetFillColor(ROOT.kYellow)
01256     hmediany_after.SetFillColor(ROOT.kYellow)
01257     hmediandxdz_after.SetFillColor(ROOT.kYellow)
01258     hmediandydz_after.SetFillColor(ROOT.kYellow)
01259 
01260     hmedianx_aftercopy = hmedianx_after.Clone()
01261     hmediany_aftercopy = hmediany_after.Clone()
01262     hmediandxdz_aftercopy = hmediandxdz_after.Clone()
01263     hmediandydz_aftercopy = hmediandydz_after.Clone()
01264     hmedianx_aftercopy.GetXaxis().SetLabelColor(ROOT.kWhite)
01265     hmediany_aftercopy.GetXaxis().SetLabelColor(ROOT.kWhite)
01266     hmediandxdz_aftercopy.GetXaxis().SetLabelColor(ROOT.kWhite)
01267     hmediandydz_aftercopy.GetXaxis().SetLabelColor(ROOT.kWhite)
01268     hmedianx_aftercopy.GetYaxis().SetLabelColor(ROOT.kWhite)
01269     hmediany_aftercopy.GetYaxis().SetLabelColor(ROOT.kWhite)
01270     hmediandxdz_aftercopy.GetYaxis().SetLabelColor(ROOT.kWhite)
01271     hmediandydz_aftercopy.GetYaxis().SetLabelColor(ROOT.kWhite)
01272 
01273     hmedianx_after.SetXTitle("median(#Deltax) (mm)")
01274     hmediany_after.SetXTitle("median(#Deltay) (mm)")
01275     hmediandxdz_after.SetXTitle("median(#Deltadx/dz) (mrad)")
01276     hmediandydz_after.SetXTitle("median(#Deltadydz) (mrad)")
01277     hmedianx_after.GetXaxis().CenterTitle()
01278     hmediany_after.GetXaxis().CenterTitle()
01279     hmediandxdz_after.GetXaxis().CenterTitle()
01280     hmediandydz_after.GetXaxis().CenterTitle()
01281 
01282     if ceilingx is not None: hmedianx_aftercopy.SetAxisRange(0., ceilingx, "Y")
01283     if ceilingy is not None: hmediany_aftercopy.SetAxisRange(0., ceilingy, "Y")
01284     if ceilingdxdz is not None: hmediandxdz_aftercopy.SetAxisRange(0., ceilingdxdz, "Y")
01285     if ceilingdydz is not None: hmediandydz_aftercopy.SetAxisRange(0., ceilingdydz, "Y")
01286 
01287     c1.Clear()
01288     c1.Divide(2, 2)
01289 
01290     c1.GetPad(1).cd()
01291     hmedianx_aftercopy.Draw()
01292     hmedianx_before.Draw("same")
01293     hmedianx_after.Draw("same")
01294     hmedianx_beforecopy.Draw("same")
01295     hmedianx_after.Draw("axissame")
01296 
01297     tlegend = ROOT.TLegend(0.17, 0.75-0.05, 0.45+0.05, 0.9)
01298     tlegend.SetFillColor(ROOT.kWhite)
01299     tlegend.SetBorderSize(0)
01300     tlegend.AddEntry(hmedianx_after, r2text, "f")
01301     tlegend.AddEntry(hmedianx_before, r1text, "f")
01302     tlegend.Draw()
01303 
01304     c1.GetPad(2).cd()
01305     hmediandxdz_aftercopy.Draw()
01306     hmediandxdz_before.Draw("same")
01307     hmediandxdz_after.Draw("same")
01308     hmediandxdz_beforecopy.Draw("same")
01309     hmediandxdz_after.Draw("axissame")
01310 
01311     c1.GetPad(3).cd()
01312     hmediany_aftercopy.Draw()
01313     hmediany_before.Draw("same")
01314     hmediany_after.Draw("same")
01315     hmediany_beforecopy.Draw("same")
01316     hmediany_after.Draw("axissame")
01317 
01318     c1.GetPad(4).cd()
01319     hmediandydz_aftercopy.Draw()
01320     hmediandydz_before.Draw("same")
01321     hmediandydz_after.Draw("same")
01322     hmediandydz_beforecopy.Draw("same")
01323     hmediandydz_after.Draw("axissame")
01324 
01325     return hmediandxdz_after, hmediandxdz_aftercopy, hmediandxdz_before, hmediandxdz_beforecopy, \
01326            hmediandydz_after, hmediandydz_aftercopy, hmediandydz_before, hmediandydz_beforecopy, \
01327            hmedianx_after, hmedianx_aftercopy, hmedianx_before, hmedianx_beforecopy, \
01328            hmediany_after, hmediany_aftercopy, hmediany_before, hmediany_beforecopy, tlegend 
01329 
01330 ######################################################################################################
01331 
01332 def mapplot(tfiles, name, param, mode="from2d", window=40., abscissa=None, title="", 
01333             widebins=False, fitsine=False, fitline=False, reset_palette=False, fitsawteeth=False):
01334     tdrStyle.SetOptTitle(1)
01335     tdrStyle.SetTitleBorderSize(0)
01336     tdrStyle.SetOptStat(0)
01337     #tdrStyle.SetOptStat("emrou")
01338     tdrStyle.SetOptFit(0)
01339     tdrStyle.SetTitleFontSize(0.05)
01340     tdrStyle.SetPadRightMargin(0.1) # to see the pallete labels on the left
01341     
01342     c1.Clear()
01343     c1.ResetAttPad()
01344     
01345     if reset_palette: set_palette("blues")
01346     global hist, hist2d, hist2dweight, tline1, tline2, tline3
01347     
01348     if fitsine or fitsawteeth:
01349         id = mapNameToId(name)
01350         if not id:
01351             print "bad id for ", name
01352             raise Exception
01353     
01354     hdir = "AlignmentMonitorMuonSystemMap1D/iter1/"
01355     hpref= "%s_%s" % (name, param)
01356     hhh  = hdir+hpref
01357 
01358     prof = tfiles[0].Get(hhh+"_prof").Clone(hpref+"_prof_")
01359     profPos = tfiles[0].Get(hhh+"_profPos").Clone(hpref+"_profPos_")
01360     profNeg = tfiles[0].Get(hhh+"_profNeg").Clone(hpref+"_profNeg_")
01361     weights = tfiles[0].Get(hhh+"_weights").Clone(hpref+"_weights_")
01362     valweights = tfiles[0].Get(hhh+"_valweights").Clone(hpref+"_valweights_")
01363     hist2d = tfiles[0].Get(hhh+"_2d").Clone(hpref+"_2d_")
01364     hist2dweight = tfiles[0].Get(hhh+"_2dweight").Clone(hpref+"_2dweight_")
01365     
01366     for tfile in tfiles[1:]:
01367         prof.Add(tfile.Get(hhh+"_prof"))
01368         profPos.Add(tfile.Get(hhh+"_profPos"))
01369         profNeg.Add(tfile.Get(hhh+"_profNeg"))
01370         weights.Add(tfile.Get(hhh+"_weights"))
01371         valweights.Add(tfile.Get(hhh+"_valweights"))
01372         hist2d.Add(tfile.Get(hhh+"_2d"))
01373         hist2dweight.Add(tfile.Get(hhh+"_2dweight"))
01374 
01375     if mode == "plain":
01376         hist = prof
01377 
01378     elif mode in ("from2d", "from2dweight"):
01379         if mode == "from2d": the2d = hist2d
01380         else: the2d = hist2dweight
01381 
01382         hist = weights.Clone()
01383         hist.Reset()
01384         skip = 1
01385         if widebins:
01386             hist.Rebin(10)
01387             skip = 10
01388 
01389         #f = ROOT.TF1("g", "gaus", -40., 40)
01390         for i in xrange(0, int(weights.GetNbinsX()), skip):
01391             tmp = the2d.ProjectionY("tmp", i+1, i + skip)
01392             if tmp.GetEntries() > 1:
01393                 #tmp.Fit("g","LNq")
01394                 hist.SetBinContent(i/skip+1, tmp.GetMean())
01395                 hist.SetBinError(i/skip+1, ROOT.TMath.StudentQuantile(0.841345,tmp.GetEntries()) * tmp.GetRMS() / sqrt(tmp.GetEntries()))
01396                 #hist.SetBinError(i/skip+1, tmp.GetRMS() / sqrt(tmp.GetEntries()))
01397                 #hist.SetBinError(i/skip+1, f.GetParameter(2))
01398             else:
01399                 #hist.SetBinContent(i/skip+1, 2000.)
01400                 #hist.SetBinError(i/skip+1, 1000.)
01401                 hist.SetBinContent(i/skip+1, 0.)
01402                 hist.SetBinError(i/skip+1, 0.)
01403         
01404     elif mode == "weighted":
01405         if weights.GetEntries() == 0:
01406             averageweight = 0.
01407         else:
01408             sumofweights = 0.
01409             for i in xrange(0, int(weights.GetNbinsX())+2):
01410                 sumofweights += weights.GetBinContent(i)
01411             averageweight = sumofweights / weights.GetEntries()
01412         hist = weights.Clone()
01413         for i in xrange(1, int(weights.GetNbinsX())+1):
01414             if weights.GetBinContent(i) > 0:
01415                 thisweight = weights.GetBinContent(i) / averageweight
01416                 hist.SetBinContent(i, valweights.GetBinContent(i) / thisweight)
01417                 hist.SetBinError(i, sqrt(1. / thisweight))
01418             else:
01419                 hist.SetBinContent(i, 2000.)
01420                 hist.SetBinError(i, 1000.)
01421 
01422     else:
01423         raise Exception
01424 
01425     hist.SetAxisRange(-window, window, "Y")
01426     if abscissa is not None: hist.SetAxisRange(abscissa[0], abscissa[1], "X")
01427     hist.SetMarkerStyle(20)
01428     hist.SetMarkerSize(0.75)
01429     hist.GetXaxis().CenterTitle()
01430     hist.GetYaxis().CenterTitle()
01431     hist.GetYaxis().SetTitleOffset(0.75)
01432     hist.GetXaxis().SetTitleSize(0.05)
01433     hist.GetYaxis().SetTitleSize(0.05)
01434     hist.SetTitle(title)
01435     if "vsphi" in name: hist.SetXTitle("Global #phi position (rad)")
01436     elif "vsz" in name: hist.SetXTitle("Global z position (cm)")
01437     elif "vsr" in name: hist.SetXTitle("Global R position (cm)")
01438     if "DT" in name:
01439         if param == "x": hist.SetYTitle("x' residual (mm)")
01440         if param == "dxdz": hist.SetYTitle("dx'/dz residual (mrad)")
01441         if param == "y": hist.SetYTitle("y' residual (mm)")
01442         if param == "dydz": hist.SetYTitle("dy'/dz residual (mm)")
01443     if "CSC" in name:
01444         if param == "x": hist.SetYTitle("r#phi residual (mm)")
01445         if param == "dxdz": hist.SetYTitle("d(r#phi)/dz residual (mrad)")
01446     hist.SetMarkerColor(ROOT.kBlack)
01447     hist.SetLineColor(ROOT.kBlack)
01448     hist.Draw()
01449     hist2d.Draw("colzsame")
01450     if widebins: hist.Draw("samee1")
01451     else: hist.Draw("same")
01452 
01453     if fitsine and "vsphi" in name:
01454         global fitsine_const, fitsine_sin, fitsine_cos, fitsine_chi2, fitsine_ndf
01455         f = ROOT.TF1("f", "[0] + [1]*sin(x) + [2]*cos(x)", -pi, pi)
01456         #hist.Fit(f,"N")
01457         hist.Fit(f,"NE")
01458         f.SetLineColor(ROOT.kRed)
01459         fitsine_const = f.GetParameter(0), f.GetParError(0)
01460         fitsine_sin = f.GetParameter(1), f.GetParError(1)
01461         fitsine_cos = f.GetParameter(2), f.GetParError(2)
01462         fitsine_chi2 = f.GetChisquare()
01463         fitsine_ndf = f.GetNDF()
01464         global MAP_RESULTS_FITSIN
01465         MAP_RESULTS_FITSIN[id] = {'a':fitsine_const, 'sin': fitsine_sin, 'cos': fitsine_cos, 'chi2': fitsine_chi2, 'ndf': fitsine_ndf}
01466         f.Draw("same")
01467         global fitsine_ttext
01468         fitsine_ttext = ROOT.TLatex(-1., 0.8*window, 
01469                 "%.3g %+.3g sin(#phi) %+.3g cos(#phi)" % (fitsine_const[0], fitsine_sin[0], fitsine_cos[0]))
01470         fitsine_ttext.SetTextColor(ROOT.kRed)
01471         fitsine_ttext.Draw()
01472 
01473     if fitline:
01474         f = ROOT.TF1("f", "[0] + [1]*x", -1000., 1000.)
01475         hist2d.Fit(f, "q")
01476         hist2d.GetFunction("f").SetLineColor(ROOT.kRed)
01477         global fitline_const, fitline_linear, fitline_chi2, fitline_ndf
01478         fitline_const = hist2d.GetFunction("f").GetParameter(0), hist2d.GetFunction("f").GetParError(0)
01479         fitline_linear = hist2d.GetFunction("f").GetParameter(1), hist2d.GetFunction("f").GetParError(1)
01480         fitline_chi2 = hist2d.GetFunction("f").GetChisquare()
01481         fitline_ndf = hist2d.GetFunction("f").GetNDF()
01482         hist2d.GetFunction("f").Draw("same")
01483         global fitline_ttext
01484         if "vsz" in name: which = "Z"
01485         elif "vsr" in name: which = "R"
01486         fitline_ttext = ROOT.TText(hist.GetBinCenter(hist.GetNbinsX()/4), 
01487                 0.8*window, "%.3g %+.3g %s" % (fitline_const[0], fitline_linear[0], which))
01488         fitline_ttext.SetTextColor(ROOT.kRed)
01489         fitline_ttext.Draw()
01490 
01491     ROOT.gPad.RedrawAxis()
01492 
01493     if "vsphi" in name: 
01494         if not widebins: philines(name, window, abscissa)
01495         if abscissa is None:
01496             tline1 = ROOT.TLine(-pi, 0, pi, 0); tline1.Draw()
01497             tline2 = ROOT.TLine(-pi, -window, pi, -window); tline2.SetLineWidth(2); tline2.Draw()
01498             tline3 = ROOT.TLine(-pi, window, pi, window); tline3.Draw()
01499         else:
01500             tline1 = ROOT.TLine(abscissa[0], 0, abscissa[1], 0); tline1.Draw()
01501             tline2 = ROOT.TLine(abscissa[0], -window, abscissa[1], -window); tline2.SetLineWidth(2); tline2.Draw()
01502             tline3 = ROOT.TLine(abscissa[0], window, abscissa[1], window); tline3.Draw()
01503     elif "vsz" in name:
01504         if not widebins: zlines(window, abscissa)
01505         if abscissa is None:
01506             tline1 = ROOT.TLine(-660, 0, 660, 0); tline1.Draw()
01507             tline2 = ROOT.TLine(-660, -window, 660, -window); tline2.SetLineWidth(2); tline2.Draw()
01508             tline3 = ROOT.TLine(-660, window, 660, window); tline3.Draw()
01509         else:
01510             tline1 = ROOT.TLine(abscissa[0], 0, abscissa[1], 0); tline1.Draw()
01511             tline2 = ROOT.TLine(abscissa[0], -window, abscissa[1], -window); tline2.SetLineWidth(2); tline2.Draw()
01512             tline3 = ROOT.TLine(abscissa[0], window, abscissa[1], window); tline3.Draw()
01513     elif "vsr" in name:
01514         if "mem1" in name or "mep1" in name and not widebins: rlines(1, window, abscissa)
01515         if "mem2" in name or "mep2" in name and not widebins: rlines(2, window, abscissa)
01516         if "mem3" in name or "mep3" in name and not widebins: rlines(3, window, abscissa)
01517         if "mem4" in name or "mep4" in name and not widebins: rlines(4, window, abscissa)
01518         if abscissa is None:
01519             tline1 = ROOT.TLine(100, 0, 700, 0); tline1.Draw()
01520             tline2 = ROOT.TLine(100, -window, 700, -window); tline2.SetLineWidth(2); tline2.Draw()
01521             tline3 = ROOT.TLine(100, window, 700, window); tline3.Draw()
01522         else:
01523             tline1 = ROOT.TLine(abscissa[0], 0, abscissa[1], 0); tline1.Draw()
01524             tline2 = ROOT.TLine(abscissa[0], -window, abscissa[1], -window); tline2.SetLineWidth(2); tline2.Draw()
01525             tline3 = ROOT.TLine(abscissa[0], window, abscissa[1], window); tline3.Draw()
01526 
01527     if "vsphi" in name and fitsawteeth:
01528         global CPP_LOADED
01529         if not CPP_LOADED:
01530             phiedges2c()
01531             ROOT.gROOT.ProcessLine(".L phiedges_fitfunctions.C++")
01532             CPP_LOADED = True
01533         fn={0: ROOT.fitf0,
01534             1: ROOT.fitf2,
01535             2: ROOT.fitf2,
01536             3: ROOT.fitf3,
01537             4: ROOT.fitf4,
01538             5: ROOT.fitf5,
01539             6: ROOT.fitf6,
01540             7: ROOT.fitf7,
01541             8: ROOT.fitf8,
01542             9: ROOT.fitf9,
01543             10: ROOT.fitf10,
01544             11: ROOT.fitf11,
01545             12: ROOT.fitf12,
01546             13: ROOT.fitf13
01547         } [stationIndex(name)]
01548         fn.SetNpx(5000)
01549         fn.SetLineColor(ROOT.kYellow)
01550         hist.Fit(fn,"N")
01551         fn.Draw("same")
01552 
01553         # get properly arranged phi edges
01554         edges = (phiedges[stationIndex(name)])[:]
01555         ed = sorted(edges)
01556         # add some padding to the end
01557         ed.append(pi+abs(ed[0]))
01558 
01559         global sawtooth_a, sawtooth_b
01560         sawtooth_a = []
01561         sawtooth_da = []
01562         sawtooth_b = []
01563         for pr in range(0,fn.GetNpar(),2):
01564             sawtooth_a.append( (fn.GetParameter(pr), fn.GetParError(pr)) )
01565             sawtooth_b.append( (fn.GetParameter(pr+1), fn.GetParError(pr+1)) )
01566             sawtooth_da.append( (fn.Eval(ed[pr/2]+0.01), fn.Eval(ed[pr/2+1]-0.01)) )
01567         global MAP_RESULTS_SAWTOOTH
01568         MAP_RESULTS_SAWTOOTH[id] = {'a': sawtooth_a, 'da': sawtooth_da, 'b': sawtooth_b, 'chi2': fn.GetChisquare(), 'ndf': fn.GetNDF()}
01569 
01570     # fill number of contributiong bins
01571     
01572     
01573     #ROOT.SetOwnership(hist,False)
01574     ROOT.SetOwnership(hist2d,False)
01575     ROOT.SetOwnership(hist2dweight,False)
01576     ROOT.SetOwnership(tline1,False)
01577     ROOT.SetOwnership(tline2,False)
01578     ROOT.SetOwnership(tline3,False)
01579     return hist
01580 
01581 
01582 def mapNameToId(name):
01583   if "DT" in name:
01584     wh = "-ALL"
01585     if name.find('wh')>1: wh = name[name.find('wh')+2]
01586     if   wh == "A": w = "-2"
01587     elif wh == "B": w = "-1"
01588     elif wh == "C": w = "-0"
01589     elif wh == "D": w = "+1"
01590     elif wh == "E": w = "+2"
01591     elif wh == "-ALL": w = "-ALL"
01592     else: return None
01593     station=''
01594     if wh == "-ALL": 
01595         if name.find('sec')<0: return None
01596         station = name[name.find('sec')-1]
01597         sector = ''
01598         sector = name[name.find('sec')+3:name.find('sec')+5]
01599         return "MB%s/%s/%s" % (w, station, sector)
01600     if name.find('st')>1: station = name[name.find('st')+2]
01601     else: return None
01602     return "MB%s/%s" % (w, station)
01603   elif "CSC" in name:
01604     p = name.find('me')
01605     if p<0: return None
01606     if name[p+2]=="p": endcap = "+"
01607     elif name[p+2]=="m": endcap = "-"
01608     else: return None
01609     station = name[p+3]
01610     pch = name.find('ch')
01611     if pch<0:
01612         ring = name[p+4]
01613         return "ME%s%s/%s" % (endcap, station, ring)
01614     ring = 'ALL'
01615     chamber = name[pch+2:pch+4]
01616     return "ME%s%s/%s/%s" % (endcap, station, ring, chamber)
01617   return None
01618 
01619 
01620 def curvatureplot(tfiles, name, param, mode="from2d", window=15., widebins=False, title="", fitgauss=False, fitconst=False, fitline=False, reset_palette=False):
01621 # "param" may be one of "deltax" (Delta x position residuals), "deltadxdz" (Delta (dx/dz) angular residuals), "curverr" (Delta x * d(Delta q/pT)/d(Delta x) = Delta q/pT in the absence of misalignment)
01622     tdrStyle.SetOptTitle(1)
01623     tdrStyle.SetTitleBorderSize(0)
01624     tdrStyle.SetOptStat(0)
01625     tdrStyle.SetOptFit(0)
01626     tdrStyle.SetTitleFontSize(0.05)
01627 
01628     c1.Clear()
01629     if reset_palette: set_palette("blues")
01630     global hist, histCOPY, hist2d, tline1, tline2, tline3, tline4, tline5
01631 
01632     hdir = "AlignmentMonitorMuonVsCurvature/iter1/"
01633 
01634     if name not in ("all", "top", "bottom"):
01635         hsuffix = "_%s_%s" % (name, param)
01636         prof = tfiles[0].Get(hdir+"tprofile"+hsuffix).Clone("tprofile_"+hsuffix)
01637         hist2d = tfiles[0].Get(hdir+"th2f"+hsuffix).Clone("th2f_"+hsuffix)
01638         for tfile in tfiles[1:]:
01639             prof.Add(tfile.Get(hdir+"tprofile"+hsuffix))
01640             hist2d.Add(tfile.Get(hdir+"th2f"+hsuffix))
01641     else:
01642         prof = None
01643         hist2d = None
01644         for wheel in "m2", "m1", "z", "p1", "p2":
01645             if name == "all": sectors = "01", "02", "03", "04", "05", "06", "07", "08", "09", "10", "11", "12"
01646             elif name == "top": sectors = "01", "02", "03", "04", "05", "06"
01647             elif name == "bottom": sectors = "07", "08", "09", "10", "11", "12"
01648             else: raise Exception
01649 
01650             for sector in sectors:
01651                 hsuffix = "_%s_%s" % ("wheel%s_sector%s" % (wheel, sector), param)
01652                 for tfile in tfiles:
01653                     if prof is None:
01654                         prof = tfiles[0].Get(hdir+"tprofile"+hsuffix).Clone("tprofile_"+hsuffix)
01655                         hist2d = tfiles[0].Get(hdir+"th2f"+hsuffix).Clone("tprofile_"+hsuffix)
01656                     else:
01657                         prof.Add(tfile.Get(hdir+"tprofile"+hsuffix))
01658                         hist2d.Add(tfile.Get(hdir+"th2f"+hsuffix))
01659 
01660     hist = ROOT.TH1F("hist", "", prof.GetNbinsX(), prof.GetBinLowEdge(1), -prof.GetBinLowEdge(1))
01661     for i in xrange(1, prof.GetNbinsX()+1):
01662         hist.SetBinContent(i, prof.GetBinContent(i))
01663         hist.SetBinError(i, prof.GetBinError(i))
01664 
01665     if mode == "plain":
01666         hist = prof
01667     elif mode == "from2d":
01668         skip = 1
01669         if widebins:
01670             hist.Rebin(5)
01671             skip = 5
01672         for i in xrange(0, int(prof.GetNbinsX()), skip):
01673             tmp = hist2d.ProjectionY("tmp", i+1, i + skip)
01674             if tmp.GetEntries() > 1:
01675                 hist.SetBinContent(i/skip+1, tmp.GetMean())
01676                 hist.SetBinError(i/skip+1, ROOT.TMath.StudentQuantile(0.841345,tmp.GetEntries()) * tmp.GetRMS() / sqrt(tmp.GetEntries()))
01677                 #hist.SetBinError(i/skip+1, tmp.GetRMS() / sqrt(tmp.GetEntries()))
01678             else:
01679                 #hist.SetBinContent(i/skip+1, 2000.)
01680                 #hist.SetBinError(i/skip+1, 1000.)
01681                 hist.SetBinContent(i/skip+1, 0.)
01682                 hist.SetBinError(i/skip+1, 0.)
01683     else:
01684         raise Exception
01685 
01686 
01687     if fitgauss:
01688         f = ROOT.TF1("f", "[0] + [1]*exp(-x**2/2/0.01**2)", hist.GetBinLowEdge(1), -hist.GetBinLowEdge(1))
01689         f.SetParameters(0, 0., 0.01)
01690         hist.Fit(f, "q")
01691         f.SetLineColor(ROOT.kRed)
01692         global fitgauss_diff, fitgauss_chi2, fitgauss_ndf
01693 #         fitter = ROOT.TVirtualFitter.GetFitter()
01694 #         fitgauss_diff = f.GetParameter(0) - f.GetParameter(1), \
01695 #                         sqrt(f.GetParError(0)**2 + f.GetParError(1)**2 + 2.*fitter.GetCovarianceMatrixElement(0, 1))
01696         fitgauss_diff = f.GetParameter(1), f.GetParError(1)
01697         fitgauss_chi2 = f.GetChisquare()
01698         fitgauss_ndf = f.GetNDF()
01699 
01700     global fitline_intercept, fitline_slope
01701     if fitconst:
01702         f = ROOT.TF1("f", "[0]", hist.GetBinLowEdge(1), -hist.GetBinLowEdge(1))
01703         hist.Fit(f, "q")
01704         f.SetLineColor(ROOT.kRed)
01705         fitline_intercept = f.GetParameter(0), f.GetParError(0)
01706 
01707     if fitline:
01708         f = ROOT.TF1("f", "[0] + [1]*x", hist.GetBinLowEdge(1), -hist.GetBinLowEdge(1))
01709         hist.Fit(f, "qNE")
01710         f.SetLineColor(ROOT.kRed)
01711         global f2, f3
01712         f2 = ROOT.TF1("2", "[0] + [1]*x", hist.GetBinLowEdge(1), -hist.GetBinLowEdge(1))
01713         f3 = ROOT.TF1("2", "[0] + [1]*x", hist.GetBinLowEdge(1), -hist.GetBinLowEdge(1))
01714         f2.SetParameters(f.GetParameter(0), f.GetParameter(1) + f.GetParError(1))
01715         f3.SetParameters(f.GetParameter(0), f.GetParameter(1) - f.GetParError(1))
01716         f2.SetLineColor(ROOT.kRed)
01717         f3.SetLineColor(ROOT.kRed)
01718         f2.SetLineStyle(2)
01719         f3.SetLineStyle(2)
01720         fitline_intercept = f.GetParameter(0), f.GetParError(0)
01721         fitline_slope = f.GetParameter(1), f.GetParError(1)
01722 
01723     hist.SetAxisRange(-window, window, "Y")
01724     hist.SetMarkerStyle(20)
01725     hist.SetMarkerSize(0.75)
01726     hist.GetXaxis().CenterTitle()
01727     hist.GetYaxis().CenterTitle()
01728     if param == "curverr":
01729         hist.GetYaxis().SetTitleOffset(1.35)
01730     else:
01731         hist.GetYaxis().SetTitleOffset(0.75)
01732     hist.GetXaxis().SetTitleOffset(1.2)
01733     hist.GetXaxis().SetTitleSize(0.05)
01734     hist.GetYaxis().SetTitleSize(0.05)
01735     hist.SetTitle(title)
01736     if param == "pterr": hist.SetXTitle("qp_{T} (GeV/c)")
01737     else: hist.SetXTitle("q/p_{T} (c/GeV)")
01738     if param == "deltax": hist.SetYTitle("#Deltax' (mm)")
01739     if param == "deltadxdz": hist.SetYTitle("#Deltadx'/dz (mrad)")
01740     if param == "pterr": hist.SetYTitle("#Deltap_{T}/p_{T} (%)")
01741     if param == "curverr": hist.SetYTitle("#Deltaq/p_{T} (c/GeV)")
01742     hist.SetMarkerColor(ROOT.kBlack)
01743     hist.SetLineColor(ROOT.kBlack)
01744     hist.Draw()
01745     hist2d.Draw("colzsame")
01746     histCOPY = hist.Clone()
01747     histCOPY.SetXTitle("")
01748     histCOPY.SetYTitle("")
01749 
01750     if widebins:
01751         histCOPY.Draw("samee1")
01752         histCOPY.Draw("sameaxis")
01753     else:
01754         histCOPY.Draw("same")
01755         histCOPY.Draw("sameaxis")
01756 
01757     if fitline:
01758         f.Draw("same")
01759         #f2.Draw("same")
01760         #f3.Draw("same")
01761 
01762     tline1 = ROOT.TLine(hist.GetBinLowEdge(1), -window, hist.GetBinLowEdge(1), window)
01763     tline2 = ROOT.TLine(hist.GetBinLowEdge(1), window, -hist.GetBinLowEdge(1), window)
01764     tline3 = ROOT.TLine(-hist.GetBinLowEdge(1), window, -hist.GetBinLowEdge(1), -window)
01765     tline4 = ROOT.TLine(-hist.GetBinLowEdge(1), -window, hist.GetBinLowEdge(1), -window)
01766     tline5 = ROOT.TLine(-hist.GetBinLowEdge(1), 0., hist.GetBinLowEdge(1), 0.)
01767     for t in tline1, tline2, tline3, tline4, tline5: t.Draw()
01768 
01769 
01770 def curvatureDTsummary(tfiles, window=15., pdgSfactor=False):
01771     global h, gm2, gm1, gz, gp1, gp2, tlegend
01772 
01773     set_palette("blues")
01774     phis = {-2: [], -1: [], 0: [], 1: [], 2: []}
01775     diffs = {-2: [], -1: [], 0: [], 1: [], 2: []}
01776     differrs = {-2: [], -1: [], 0: [], 1: [], 2: []}
01777     for wheelstr, wheel in ("m2", "-2"), ("m1", "-1"), ("z", "0"), ("p1", "+1"), ("p2", "+2"):
01778         for sector in "01", "02", "03", "04", "05", "06", "07", "08", "09", "10", "11", "12":
01779             curvatureplot(tfiles, "wheel%s_sector%s" % (wheelstr, sector), "deltax", 
01780                     title="Wheel %s, sector %s" % (wheel, sector), fitgauss=True, reset_palette=False)
01781             if fitgauss_diff[1] < window:
01782                 uncertainty = fitgauss_diff[1]
01783                 if pdgSfactor and (fitgauss_chi2/fitgauss_ndf) > 1.: uncertainty *= sqrt(fitgauss_chi2/fitgauss_ndf)
01784 
01785                 phis[int(wheel)].append(signConventions["DT", int(wheel), 1, int(sector)][4])
01786                 diffs[int(wheel)].append(fitgauss_diff[0])
01787                 differrs[int(wheel)].append(uncertainty)
01788 
01789     h = ROOT.TH1F("h", "", 1, -pi, pi)
01790     h.SetAxisRange(-window, window, "Y")
01791     h.SetXTitle("#phi (rad)")
01792     h.SetYTitle("#Deltax(p_{T} #rightarrow #infty) - #Deltax(p_{T} #rightarrow 0) (mm)")
01793     h.GetXaxis().CenterTitle()
01794     h.GetYaxis().CenterTitle()
01795 
01796     gm2 = ROOT.TGraphErrors(len(phis[-2]), array.array("d", phis[-2]), array.array("d", diffs[-2]), 
01797             array.array("d", [0.]*len(phis[-2])), array.array("d", differrs[-2]))
01798     gm1 = ROOT.TGraphErrors(len(phis[-1]), array.array("d", phis[-1]), array.array("d", diffs[-1]), 
01799             array.array("d", [0.]*len(phis[-1])), array.array("d", differrs[-1]))
01800     gz = ROOT.TGraphErrors(len(phis[0]), array.array("d", phis[0]), array.array("d", diffs[0]), 
01801             array.array("d", [0.]*len(phis[0])), array.array("d", differrs[0]))
01802     gp1 = ROOT.TGraphErrors(len(phis[1]), array.array("d", phis[1]), array.array("d", diffs[1]), 
01803             array.array("d", [0.]*len(phis[1])), array.array("d", differrs[1]))
01804     gp2 = ROOT.TGraphErrors(len(phis[2]), array.array("d", phis[2]), array.array("d", diffs[2]), 
01805             array.array("d", [0.]*len(phis[2])), array.array("d", differrs[2]))
01806 
01807     gm2.SetMarkerStyle(21); gm2.SetMarkerColor(ROOT.kRed); gm2.SetLineColor(ROOT.kRed)
01808     gm1.SetMarkerStyle(22); gm1.SetMarkerColor(ROOT.kBlue); gm1.SetLineColor(ROOT.kBlue)
01809     gz.SetMarkerStyle(3); gz.SetMarkerColor(ROOT.kBlack); gz.SetLineColor(ROOT.kBlack)
01810     gp1.SetMarkerStyle(26); gp1.SetMarkerColor(ROOT.kBlue); gp1.SetLineColor(ROOT.kBlue)
01811     gp2.SetMarkerStyle(25); gp2.SetMarkerColor(ROOT.kRed); gp2.SetLineColor(ROOT.kRed)
01812 
01813     h.Draw()
01814     tlegend = ROOT.TLegend(0.25, 0.2, 0.85, 0.5)
01815     tlegend.SetFillColor(ROOT.kWhite)
01816     tlegend.SetBorderSize(0)
01817     tlegend.AddEntry(gm2, "Wheel -2", "p")
01818     tlegend.AddEntry(gm1, "Wheel -1", "p")
01819     tlegend.AddEntry(gz, "Wheel 0", "p")
01820     tlegend.AddEntry(gp1, "Wheel +1", "p")
01821     tlegend.AddEntry(gp2, "Wheel +2", "p")
01822     tlegend.Draw()
01823 
01824     gm2.Draw("p")
01825     gm1.Draw("p")
01826     gz.Draw("p")
01827     gp1.Draw("p")
01828     gp2.Draw("p")
01829 
01830 
01831 def getname(r):
01832     if r.postal_address[0] == "DT":
01833         wheel, station, sector = r.postal_address[1:]
01834         return "DT wheel %d, station %d, sector %d" % (wheel, station, sector)
01835     elif r.postal_address[0] == "CSC":
01836         endcap, station, ring, chamber = r.postal_address[1:]
01837         if endcap != 1: station = -1 * abs(station)
01838         return "CSC ME%d/%d chamber %d" % (station, ring, chamber)
01839 
01840 ddt=[0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.]
01841 def clearDDT():
01842     for i in range(0,15):
01843         ddt[i]=0.
01844 
01845 def printDeltaTs():
01846     n = 0
01847     for t in ddt:
01848         if n==0 or n==7 or n==15: print "%d calls" % t
01849         else: print "%d : %0.3f ms" % (n,t*1000.0)
01850         n += 1
01851 
01852 def bellcurves(tfile, reports, name, twobin=True, suppressblue=False):
01853     t1 = time.time()
01854     ddt[0] += 1
01855     tdrStyle.SetOptTitle(1)
01856     tdrStyle.SetTitleBorderSize(1)
01857     tdrStyle.SetTitleFontSize(0.1)
01858     tdrStyle.SetOptStat(0)
01859     tdrStyle.SetHistMinimumZero()
01860 
01861     c1.Clear()
01862     c1.ResetAttPad()
01863 
01864     found = False
01865     for r in reports:
01866         if r.name == name:
01867             found = True
01868             break
01869     if not found: raise Exception, "Not a valid name"
01870     if r.status == "FAIL":
01871         #raise Exception, "Fit failed"
01872         print "Fit failed"
01873         c1.Clear()
01874         return
01875     
01876     Pos = "Pos"; Neg = "Neg"
01877     if not twobin:
01878         Pos = ""; Neg = ""
01879 
01880     pdirPos = "MuonAlignmentFromReference/%s%s" % (name, Pos)
01881     pdirNeg = "MuonAlignmentFromReference/%s%s" % (name, Neg)
01882 
01883     t2 = time.time()
01884     ddt[1] = 1./ddt[0]*((ddt[0]-1)*ddt[1] + t2-t1)
01885 
01886     chamber_x = tfile.Get(pdirPos+"_x")
01887     chamber_x_fit = tfile.Get(pdirPos+"_x_fit")
01888     chamber_y = tfile.Get(pdirPos+"_y")
01889     chamber_y_fit = tfile.Get(pdirPos+"_y_fit")
01890     chamber_dxdz = tfile.Get(pdirPos+"_dxdz")
01891     chamber_dxdz_fit = tfile.Get(pdirPos+"_dxdz_fit")
01892     chamber_dydz = tfile.Get(pdirPos+"_dydz")
01893     chamber_dydz_fit = tfile.Get(pdirPos+"_dydz_fit")
01894     chamber_alphax = tfile.Get(pdirPos+"_alphax")
01895     chamber_alphax_fit = tfile.Get(pdirPos+"_alphax_fit")
01896     chamber_alphay = tfile.Get(pdirPos+"_alphay")
01897     chamber_alphay_fit = tfile.Get(pdirPos+"_alphay_fit")
01898     chamber_x_fit2 = tfile.Get(pdirNeg+"_x_fit")
01899     chamber_y_fit2 = tfile.Get(pdirNeg+"_y_fit")
01900     chamber_dxdz_fit2 = tfile.Get(pdirNeg+"_dxdz_fit")
01901     chamber_dydz_fit2 = tfile.Get(pdirNeg+"_dydz_fit")
01902     chamber_alphax_fit2 = tfile.Get(pdirNeg+"_alphax_fit")
01903     chamber_alphay_fit2 = tfile.Get(pdirNeg+"_alphay_fit")
01904 
01905     if not chamber_x:
01906         chamber_x = tfile.Get(pdirPos+"_residual")
01907         chamber_x_fit = tfile.Get(pdirPos+"_residual_fit")
01908         chamber_dxdz = tfile.Get(pdirPos+"_resslope")
01909         chamber_dxdz_fit = tfile.Get(pdirPos+"_resslope_fit")
01910         chamber_alphax = tfile.Get(pdirPos+"_alpha")
01911         chamber_alphax_fit = tfile.Get(pdirPos+"_alpha_fit")
01912         chamber_x_fit2 = tfile.Get(pdirNeg+"_residual_fit")
01913         chamber_dxdz_fit2 = tfile.Get(pdirNeg+"_resslope_fit")
01914         chamber_alphax_fit2 = tfile.Get(pdirNeg+"_alpha_fit")
01915 
01916     if not chamber_x:
01917         print "Can't find neither "+pdirPos+"_x  nor "+pdirPos+"_residual"
01918         return
01919 
01920     t3 = time.time()
01921     ddt[2] = 1./ddt[0]*((ddt[0]-1)*ddt[2] + t3-t2)
01922 
01923     ####
01924     chamber_x.SetAxisRange(-30., 30., "X")
01925     chamber_dxdz.SetAxisRange(-50., 50., "X")
01926     if not not chamber_y:
01927         chamber_y.SetAxisRange(-150., 150., "X")
01928         chamber_dydz.SetAxisRange(-200., 200., "X")
01929     ####
01930 
01931     chamber_x.SetXTitle("Local x residual (mm)")
01932     chamber_dxdz.SetXTitle("Local dx/dz residual (mrad)")
01933     chamber_alphax.SetXTitle("Local dx/dz residual (mrad)")
01934     chamber_alphax.SetYTitle("Local x residual (mm)")
01935     if not not chamber_y:
01936         chamber_y.SetXTitle("Local y residual (mm)")
01937         chamber_dydz.SetXTitle("Local dy/dz residual (mrad)")
01938         chamber_alphay.SetXTitle("Local dy/dz residual (mrad)")
01939         chamber_alphay.SetYTitle("Local y residual (mm)")
01940     if name[0:2] == "ME":
01941         chamber_x.SetXTitle("Local r#phi residual (mm)")
01942         chamber_dxdz.SetXTitle("Local d(r#phi)/dz residual (mrad)")
01943         chamber_alphax.SetXTitle("Local d(r#phi)/dz residual (mrad)")
01944         chamber_alphax.SetYTitle("Local r#phi residual (mm)")
01945 
01946     t4 = time.time()
01947     ddt[3] = 1./ddt[0]*((ddt[0]-1)*ddt[3] + t4-t3)
01948 
01949     for h in chamber_x, chamber_dxdz, chamber_alphax, chamber_alphax, \
01950              chamber_y, chamber_dydz, chamber_alphay, chamber_alphay:
01951         if not not h:
01952             h.GetXaxis().CenterTitle()
01953             h.GetYaxis().CenterTitle()
01954             h.GetXaxis().SetLabelSize(0.05)
01955             h.GetYaxis().SetLabelSize(0.05)
01956             h.GetXaxis().SetTitleSize(0.07)
01957             h.GetYaxis().SetTitleSize(0.07)
01958             h.GetXaxis().SetTitleOffset(0.9)
01959             h.GetYaxis().SetTitleOffset(0.9)
01960 
01961     for f in chamber_x_fit2, chamber_y_fit2, chamber_dxdz_fit2, chamber_dydz_fit2, \
01962              chamber_alphax_fit2, chamber_alphay_fit2:
01963         if not not f:
01964             f.SetLineColor(4)
01965 
01966     t5 = time.time()
01967     ddt[4] = 1./ddt[0]*((ddt[0]-1)*ddt[4] + t5-t4)
01968 
01969     if not not chamber_y:
01970         c1.Clear()
01971         c1.Divide(3, 2)
01972         chamber_x.SetTitle(getname(r))
01973 
01974         c1.GetPad(1).cd()
01975         chamber_x.Draw()
01976         if not suppressblue: chamber_x_fit2.Draw("same")
01977         chamber_x_fit.Draw("same")
01978         
01979         c1.GetPad(2).cd()
01980         chamber_dxdz.Draw()
01981         if not suppressblue: chamber_dxdz_fit2.Draw("same")
01982         chamber_dxdz_fit.Draw("same")
01983         
01984         c1.GetPad(3).cd()
01985         chamber_alphax.Draw()
01986         if not suppressblue: chamber_alphax_fit2.Draw("same")
01987         chamber_alphax_fit.Draw("same")
01988         
01989         c1.GetPad(4).cd()
01990         chamber_y.Draw()
01991         if not suppressblue: chamber_y_fit2.Draw("same")
01992         chamber_y_fit.Draw("same")
01993         
01994         c1.GetPad(5).cd()
01995         chamber_dydz.Draw()
01996         if not suppressblue: chamber_dydz_fit2.Draw("same")
01997         chamber_dydz_fit.Draw("same")
01998         
01999         c1.GetPad(6).cd()
02000         chamber_alphay.Draw()
02001         if not suppressblue: chamber_alphay_fit2.Draw("same")
02002         chamber_alphay_fit.Draw("same")
02003 
02004     else:
02005         c1.Clear()
02006         c1.Divide(3, 1)
02007         chamber_x.SetTitle(getname(r))
02008 
02009         c1.GetPad(1).cd()
02010         chamber_x.Draw()
02011         if not suppressblue: chamber_x_fit2.Draw("same")
02012         chamber_x_fit.Draw("same")
02013         
02014         c1.GetPad(2).cd()
02015         chamber_dxdz.Draw()
02016         if not suppressblue: chamber_dxdz_fit2.Draw("same")
02017         chamber_dxdz_fit.Draw("same")
02018         
02019         c1.GetPad(3).cd()
02020         chamber_alphax.Draw()
02021         if not suppressblue: chamber_alphax_fit2.Draw("same")
02022         chamber_alphax_fit.Draw("same")
02023 
02024     t6 = time.time()
02025     ddt[5] = 1./ddt[0]*((ddt[0]-1)*ddt[5] + t6-t5)
02026     ddt[6] = 1./ddt[0]*((ddt[0]-1)*ddt[6] + t6-t1)
02027 
02028 
02029 def polynomials(tfile, reports, name, twobin=True, suppressblue=False):
02030     t1 = time.time()
02031     ddt[7] += 1
02032     global label1, label2, label3, label4, label5, label6, label7, label8, label9
02033     plotDirectory = "MuonAlignmentFromReference"
02034     tdrStyle.SetOptTitle(1)
02035     tdrStyle.SetTitleBorderSize(1)
02036     tdrStyle.SetTitleFontSize(0.1)
02037     tdrStyle.SetOptStat(0)
02038 
02039     c1.Clear()
02040     c1.ResetAttPad()
02041 
02042     found = False
02043     for r in reports:
02044         if r.name == name:
02045             found = True
02046             break
02047     if not found: raise Exception, "Not a valid name"
02048 
02049     if r.status == "FAIL":
02050         #raise Exception, "Fit failed"
02051         print "Fit failed"
02052         c1.Clear()
02053         return
02054 
02055     Pos = "Pos"; Neg = "Neg"
02056     if not twobin:
02057         Pos = ""; Neg = ""
02058 
02059     pdirPos = "MuonAlignmentFromReference/%s%s" % (name, Pos)
02060     pdirNeg = "MuonAlignmentFromReference/%s%s" % (name, Neg)
02061 
02062     global chamber_x_trackx, chamber_x_trackx_fit, chamber_y_trackx, chamber_y_trackx_fit, \
02063         chamber_dxdz_trackx, chamber_dxdz_trackx_fit, chamber_dydz_trackx, chamber_dydz_trackx_fit, \
02064         chamber_x_trackx_fit2, chamber_y_trackx_fit2, chamber_dxdz_trackx_fit2, chamber_dydz_trackx_fit2
02065     global chamber_x_tracky, chamber_x_tracky_fit, chamber_y_tracky, chamber_y_tracky_fit, \
02066         chamber_dxdz_tracky, chamber_dxdz_tracky_fit, chamber_dydz_tracky, chamber_dydz_tracky_fit, \
02067         chamber_x_tracky_fit2, chamber_y_tracky_fit2, chamber_dxdz_tracky_fit2, chamber_dydz_tracky_fit2
02068     global chamber_x_trackdxdz, chamber_x_trackdxdz_fit, chamber_y_trackdxdz, chamber_y_trackdxdz_fit, \
02069         chamber_dxdz_trackdxdz, chamber_dxdz_trackdxdz_fit, chamber_dydz_trackdxdz, chamber_dydz_trackdxdz_fit, \
02070         chamber_x_trackdxdz_fit2, chamber_y_trackdxdz_fit2, chamber_dxdz_trackdxdz_fit2, chamber_dydz_trackdxdz_fit2
02071     global chamber_x_trackdydz, chamber_x_trackdydz_fit, chamber_y_trackdydz, chamber_y_trackdydz_fit, \
02072         chamber_dxdz_trackdydz, chamber_dxdz_trackdydz_fit, chamber_dydz_trackdydz, chamber_dydz_trackdydz_fit, \
02073         chamber_x_trackdydz_fit2, chamber_y_trackdydz_fit2, chamber_dxdz_trackdydz_fit2, chamber_dydz_trackdydz_fit2
02074 
02075     chamber_x_trackx = tfile.Get(pdirPos+"_x_trackx")
02076     chamber_x_trackx_fit = tfile.Get(pdirPos+"_x_trackx_fitline")
02077     chamber_y_trackx = tfile.Get(pdirPos+"_y_trackx")
02078     chamber_y_trackx_fit = tfile.Get(pdirPos+"_y_trackx_fitline")
02079     chamber_dxdz_trackx = tfile.Get(pdirPos+"_dxdz_trackx")
02080     chamber_dxdz_trackx_fit = tfile.Get(pdirPos+"_dxdz_trackx_fitline")
02081     chamber_dydz_trackx = tfile.Get(pdirPos+"_dydz_trackx")
02082     chamber_dydz_trackx_fit = tfile.Get(pdirPos+"_dydz_trackx_fitline")
02083     chamber_x_trackx_fit2 = tfile.Get(pdirNeg+"_x_trackx_fitline")
02084     chamber_y_trackx_fit2 = tfile.Get(pdirNeg+"_y_trackx_fitline")
02085     chamber_dxdz_trackx_fit2 = tfile.Get(pdirNeg+"_dxdz_trackx_fitline")
02086     chamber_dydz_trackx_fit2 = tfile.Get(pdirNeg+"_dydz_trackx_fitline")
02087 
02088     chamber_x_tracky = tfile.Get(pdirPos+"_x_tracky")
02089     chamber_x_tracky_fit = tfile.Get(pdirPos+"_x_tracky_fitline")
02090     chamber_y_tracky = tfile.Get(pdirPos+"_y_tracky")
02091     chamber_y_tracky_fit = tfile.Get(pdirPos+"_y_tracky_fitline")
02092     chamber_dxdz_tracky = tfile.Get(pdirPos+"_dxdz_tracky")
02093     chamber_dxdz_tracky_fit = tfile.Get(pdirPos+"_dxdz_tracky_fitline")
02094     chamber_dydz_tracky = tfile.Get(pdirPos+"_dydz_tracky")
02095     chamber_dydz_tracky_fit = tfile.Get(pdirPos+"_dydz_tracky_fitline")
02096     chamber_x_tracky_fit2 = tfile.Get(pdirNeg+"_x_tracky_fitline")
02097     chamber_y_tracky_fit2 = tfile.Get(pdirNeg+"_y_tracky_fitline")
02098     chamber_dxdz_tracky_fit2 = tfile.Get(pdirNeg+"_dxdz_tracky_fitline")
02099     chamber_dydz_tracky_fit2 = tfile.Get(pdirNeg+"_dydz_tracky_fitline")
02100 
02101     chamber_x_trackdxdz = tfile.Get(pdirPos+"_x_trackdxdz")
02102     chamber_x_trackdxdz_fit = tfile.Get(pdirPos+"_x_trackdxdz_fitline")
02103     chamber_y_trackdxdz = tfile.Get(pdirPos+"_y_trackdxdz")
02104     chamber_y_trackdxdz_fit = tfile.Get(pdirPos+"_y_trackdxdz_fitline")
02105     chamber_dxdz_trackdxdz = tfile.Get(pdirPos+"_dxdz_trackdxdz")
02106     chamber_dxdz_trackdxdz_fit = tfile.Get(pdirPos+"_dxdz_trackdxdz_fitline")
02107     chamber_dydz_trackdxdz = tfile.Get(pdirPos+"_dydz_trackdxdz")
02108     chamber_dydz_trackdxdz_fit = tfile.Get(pdirPos+"_dydz_trackdxdz_fitline")
02109     chamber_x_trackdxdz_fit2 = tfile.Get(pdirNeg+"_x_trackdxdz_fitline")
02110     chamber_y_trackdxdz_fit2 = tfile.Get(pdirNeg+"_y_trackdxdz_fitline")
02111     chamber_dxdz_trackdxdz_fit2 = tfile.Get(pdirNeg+"_dxdz_trackdxdz_fitline")
02112     chamber_dydz_trackdxdz_fit2 = tfile.Get(pdirNeg+"_dydz_trackdxdz_fitline")
02113 
02114     chamber_x_trackdydz = tfile.Get(pdirPos+"_x_trackdydz")
02115     chamber_x_trackdydz_fit = tfile.Get(pdirPos+"_x_trackdydz_fitline")
02116     chamber_y_trackdydz = tfile.Get(pdirPos+"_y_trackdydz")
02117     chamber_y_trackdydz_fit = tfile.Get(pdirPos+"_y_trackdydz_fitline")
02118     chamber_dxdz_trackdydz = tfile.Get(pdirPos+"_dxdz_trackdydz")
02119     chamber_dxdz_trackdydz_fit = tfile.Get(pdirPos+"_dxdz_trackdydz_fitline")
02120     chamber_dydz_trackdydz = tfile.Get(pdirPos+"_dydz_trackdydz")
02121     chamber_dydz_trackdydz_fit = tfile.Get(pdirPos+"_dydz_trackdydz_fitline")
02122     chamber_x_trackdydz_fit2 = tfile.Get(pdirNeg+"_x_trackdydz_fitline")
02123     chamber_y_trackdydz_fit2 = tfile.Get(pdirNeg+"_y_trackdydz_fitline")
02124     chamber_dxdz_trackdydz_fit2 = tfile.Get(pdirNeg+"_dxdz_trackdydz_fitline")
02125     chamber_dydz_trackdydz_fit2 = tfile.Get(pdirNeg+"_dydz_trackdydz_fitline")
02126 
02127     if not chamber_x_trackx:
02128         chamber_x_trackx = tfile.Get(pdirPos+"_residual_trackx")
02129         chamber_x_trackx_fit = tfile.Get(pdirPos+"_residual_trackx_fitline")
02130         chamber_dxdz_trackx = tfile.Get(pdirPos+"_resslope_trackx")
02131         chamber_dxdz_trackx_fit = tfile.Get(pdirPos+"_resslope_trackx_fitline")
02132         chamber_x_trackx_fit2 = tfile.Get(pdirNeg+"_residual_trackx_fitline")
02133         chamber_dxdz_trackx_fit2 = tfile.Get(pdirNeg+"_resslope_trackx_fitline")
02134 
02135         chamber_x_tracky = tfile.Get(pdirPos+"_residual_tracky")
02136         chamber_x_tracky_fit = tfile.Get(pdirPos+"_residual_tracky_fitline")
02137         chamber_dxdz_tracky = tfile.Get(pdirPos+"_resslope_tracky")
02138         chamber_dxdz_tracky_fit = tfile.Get(pdirPos+"_resslope_tracky_fitline")
02139         chamber_x_tracky_fit2 = tfile.Get(pdirNeg+"_residual_tracky_fitline")
02140         chamber_dxdz_tracky_fit2 = tfile.Get(pdirNeg+"_resslope_tracky_fitline")
02141 
02142         chamber_x_trackdxdz = tfile.Get(pdirPos+"_residual_trackdxdz")
02143         chamber_x_trackdxdz_fit = tfile.Get(pdirPos+"_residual_trackdxdz_fitline")
02144         chamber_dxdz_trackdxdz = tfile.Get(pdirPos+"_resslope_trackdxdz")
02145         chamber_dxdz_trackdxdz_fit = tfile.Get(pdirPos+"_resslope_trackdxdz_fitline")
02146         chamber_x_trackdxdz_fit2 = tfile.Get(pdirNeg+"_residual_trackdxdz_fitline")
02147         chamber_dxdz_trackdxdz_fit2 = tfile.Get(pdirNeg+"_resslope_trackdxdz_fitline")
02148 
02149         chamber_x_trackdydz = tfile.Get(pdirPos+"_residual_trackdydz")
02150         chamber_x_trackdydz_fit = tfile.Get(pdirPos+"_residual_trackdydz_fitline")
02151         chamber_dxdz_trackdydz = tfile.Get(pdirPos+"_resslope_trackdydz")
02152         chamber_dxdz_trackdydz_fit = tfile.Get(pdirPos+"_resslope_trackdydz_fitline")
02153         chamber_x_trackdydz_fit2 = tfile.Get(pdirNeg+"_residual_trackdydz_fitline")
02154         chamber_dxdz_trackdydz_fit2 = tfile.Get(pdirNeg+"_resslope_trackdydz_fitline")
02155 
02156     if not chamber_x_trackx:
02157         print "Can't find neither "+pdirPos+"_residual  nor "+pdirPos+"_residual_trackx"
02158         return
02159 
02160     chamber_x_trackx = chamber_x_trackx.Clone()
02161     chamber_dxdz_trackx = chamber_dxdz_trackx.Clone()
02162     chamber_x_tracky = chamber_x_tracky.Clone()
02163     chamber_dxdz_tracky = chamber_dxdz_tracky.Clone()
02164     chamber_x_trackdxdz = chamber_x_trackdxdz.Clone()
02165     chamber_dxdz_trackdxdz = chamber_dxdz_trackdxdz.Clone()
02166     chamber_x_trackdydz = chamber_x_trackdydz.Clone()
02167     chamber_dxdz_trackdydz = chamber_dxdz_trackdydz.Clone()
02168 
02169     if not not chamber_y_trackx:
02170         chamber_y_trackx = chamber_y_trackx.Clone()
02171         chamber_dydz_trackx = chamber_dydz_trackx.Clone()
02172         chamber_y_tracky = chamber_y_tracky.Clone()
02173         chamber_dydz_tracky = chamber_dydz_tracky.Clone()
02174         chamber_y_trackdxdz = chamber_y_trackdxdz.Clone()
02175         chamber_dydz_trackdxdz = chamber_dydz_trackdxdz.Clone()
02176         chamber_y_trackdydz = chamber_y_trackdydz.Clone()
02177         chamber_dydz_trackdydz = chamber_dydz_trackdydz.Clone()
02178 
02179     if not not chamber_y_trackx:
02180         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_x_trackx")); chamber_x_trackx.Merge(tlist)
02181         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_dxdz_trackx")); chamber_dxdz_trackx.Merge(tlist)
02182         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_x_tracky")); chamber_x_tracky.Merge(tlist)
02183         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_dxdz_tracky")); chamber_dxdz_tracky.Merge(tlist)
02184         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_x_trackdxdz")); chamber_x_trackdxdz.Merge(tlist)
02185         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_dxdz_trackdxdz")); chamber_dxdz_trackdxdz.Merge(tlist)
02186         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_x_trackdydz")); chamber_x_trackdydz.Merge(tlist)
02187         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_dxdz_trackdydz")); chamber_dxdz_trackdydz.Merge(tlist)
02188         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_y_trackx")); chamber_y_trackx.Merge(tlist)
02189         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_dydz_trackx")); chamber_dydz_trackx.Merge(tlist)
02190         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_y_tracky")); chamber_y_tracky.Merge(tlist)
02191         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_dydz_tracky")); chamber_dydz_tracky.Merge(tlist)
02192         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_y_trackdxdz")); chamber_y_trackdxdz.Merge(tlist)
02193         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_dydz_trackdxdz")); chamber_dydz_trackdxdz.Merge(tlist)
02194         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_y_trackdydz")); chamber_y_trackdydz.Merge(tlist)
02195         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_dydz_trackdydz")); chamber_dydz_trackdydz.Merge(tlist)
02196     else:
02197         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_residual_trackx")); chamber_x_trackx.Merge(tlist)
02198         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_resslope_trackx")); chamber_dxdz_trackx.Merge(tlist)
02199         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_residual_tracky")); chamber_x_tracky.Merge(tlist)
02200         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_resslope_tracky")); chamber_dxdz_tracky.Merge(tlist)
02201         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_residual_trackdxdz")); chamber_x_trackdxdz.Merge(tlist)
02202         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_resslope_trackdxdz")); chamber_dxdz_trackdxdz.Merge(tlist)
02203         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_residual_trackdydz")); chamber_x_trackdydz.Merge(tlist)
02204         tlist = ROOT.TList(); tlist.Add(tfile.Get(pdirNeg+"_resslope_trackdydz")); chamber_dxdz_trackdydz.Merge(tlist)
02205 
02206     rr1=19.99
02207     rr2=9.99
02208     chamber_x_trackx.SetAxisRange(-rr1, rr1, "Y")
02209     chamber_dxdz_trackx.SetAxisRange(-rr2, rr2, "Y")
02210     chamber_x_tracky.SetAxisRange(-rr1, rr1, "Y")
02211     chamber_dxdz_tracky.SetAxisRange(-rr2, rr2, "Y")
02212     chamber_x_trackdxdz.SetAxisRange(-rr1, rr1, "Y")
02213     chamber_dxdz_trackdxdz.SetAxisRange(-rr2, rr2, "Y")
02214     chamber_x_trackdydz.SetAxisRange(-rr1, rr1, "Y")
02215     chamber_dxdz_trackdydz.SetAxisRange(-rr2, rr2, "Y")
02216 
02217     rr3=49.99
02218     if not not chamber_y_trackx:
02219         chamber_y_trackx.SetAxisRange(-rr3, rr3, "Y")
02220         chamber_dydz_trackx.SetAxisRange(-rr3, rr3, "Y")
02221         chamber_y_tracky.SetAxisRange(-rr3, rr3, "Y")
02222         chamber_dydz_tracky.SetAxisRange(-rr3, rr3, "Y")
02223         chamber_y_trackdxdz.SetAxisRange(-rr3, rr3, "Y")
02224         chamber_dydz_trackdxdz.SetAxisRange(-rr3, rr3, "Y")
02225         chamber_y_trackdydz.SetAxisRange(-rr3, rr3, "Y")
02226         chamber_dydz_trackdydz.SetAxisRange(-rr3, rr3, "Y")
02227 
02228     for h in chamber_x_trackx, chamber_y_trackx, chamber_dxdz_trackx, chamber_dydz_trackx, \
02229              chamber_x_tracky, chamber_y_tracky, chamber_dxdz_tracky, chamber_dydz_tracky, \
02230              chamber_x_trackdxdz, chamber_y_trackdxdz, chamber_dxdz_trackdxdz, chamber_dydz_trackdxdz, \
02231              chamber_x_trackdydz, chamber_y_trackdydz, chamber_dxdz_trackdydz, chamber_dydz_trackdydz:
02232         if not not h:
02233             h.SetMarkerStyle(20)
02234             h.SetMarkerSize(0.5)
02235             h.GetXaxis().SetLabelSize(0.12)
02236             h.GetYaxis().SetLabelSize(0.12)
02237             h.GetXaxis().SetNdivisions(505)
02238             h.GetYaxis().SetNdivisions(505)
02239             h.GetXaxis().SetLabelOffset(0.03)
02240             h.GetYaxis().SetLabelOffset(0.03)
02241 
02242     trackdxdz_minimum, trackdxdz_maximum = None, None
02243     for h in chamber_x_trackdxdz, chamber_y_trackdxdz, chamber_dxdz_trackdxdz, chamber_dydz_trackdxdz:
02244         if not not h:
02245             for i in xrange(1, h.GetNbinsX()+1):
02246                 if h.GetBinError(i) > 0.01 and h.GetBinContent(i) - h.GetBinError(i) < 10. and \
02247                    h.GetBinContent(i) + h.GetBinError(i) > -10.:
02248                     if not trackdxdz_minimum or trackdxdz_minimum > h.GetBinCenter(i): 
02249                         trackdxdz_minimum = h.GetBinCenter(i)
02250                     if trackdxdz_maximum < h.GetBinCenter(i): 
02251                         trackdxdz_maximum = h.GetBinCenter(i)
02252     if not not trackdxdz_minimum and not not trackdxdz_maximum:
02253         for h in chamber_x_trackdxdz, chamber_y_trackdxdz, chamber_dxdz_trackdxdz, chamber_dydz_trackdxdz:
02254             if not not h:
02255                 h.SetAxisRange(trackdxdz_minimum, trackdxdz_maximum, "X")
02256 
02257     trackdydz_minimum, trackdydz_maximum = None, None
02258     for h in chamber_x_trackdydz, chamber_y_trackdydz, chamber_dxdz_trackdydz, chamber_dydz_trackdydz:
02259         if not not h:
02260             for i in xrange(1, h.GetNbinsX()+1):
02261                 if h.GetBinError(i) > 0.01 and h.GetBinContent(i) - h.GetBinError(i) < 10. and \
02262                    h.GetBinContent(i) + h.GetBinError(i) > -10.:
02263                     if not trackdydz_minimum or trackdydz_minimum > h.GetBinCenter(i): 
02264                         trackdydz_minimum = h.GetBinCenter(i)
02265                     if trackdydz_maximum < h.GetBinCenter(i): 
02266                         trackdydz_maximum = h.GetBinCenter(i)
02267     if not not trackdydz_minimum and not not trackdydz_maximum:
02268         for h in chamber_x_trackdydz, chamber_y_trackdydz, chamber_dxdz_trackdydz, chamber_dydz_trackdydz:
02269             if not not h:
02270                 h.SetAxisRange(trackdydz_minimum, trackdydz_maximum, "X")
02271 
02272     for f in chamber_x_trackx_fit2, chamber_y_trackx_fit2, chamber_dxdz_trackx_fit2, chamber_dydz_trackx_fit2, \
02273              chamber_x_tracky_fit2, chamber_y_tracky_fit2, chamber_dxdz_tracky_fit2, chamber_dydz_tracky_fit2, \
02274              chamber_x_trackdxdz_fit2, chamber_y_trackdxdz_fit2, chamber_dxdz_trackdxdz_fit2, chamber_dydz_trackdxdz_fit2, \
02275              chamber_x_trackdydz_fit2, chamber_y_trackdydz_fit2, chamber_dxdz_trackdydz_fit2, chamber_dydz_trackdydz_fit2:
02276         if not not f:
02277             f.SetLineColor(4)
02278 
02279     if not not chamber_y_trackx:
02280         c1.Clear()
02281         #c1.Divide(5, 5, 1e-5, 1e-5)
02282         pads = [None]
02283         pads.append(ROOT.TPad("p1" ,"",0.00,0.78,0.07,1.00,0,0,0))
02284         pads.append(ROOT.TPad("p2" ,"",0.07,0.78,0.34,1.00,0,0,0))
02285         pads.append(ROOT.TPad("p3" ,"",0.34,0.78,0.56,1.00,0,0,0))
02286         pads.append(ROOT.TPad("p4" ,"",0.56,0.78,0.78,1.00,0,0,0))
02287         pads.append(ROOT.TPad("p5" ,"",0.78,0.78,1.00,1.00,0,0,0))
02288         pads.append(ROOT.TPad("p6" ,"",0.00,0.56,0.07,0.78,0,0,0))
02289         pads.append(ROOT.TPad("p7" ,"",0.07,0.56,0.34,0.78,0,0,0))
02290         pads.append(ROOT.TPad("p8" ,"",0.34,0.56,0.56,0.78,0,0,0))
02291         pads.append(ROOT.TPad("p9" ,"",0.56,0.56,0.78,0.78,0,0,0))
02292         pads.append(ROOT.TPad("p10","",0.78,0.56,1.00,0.78,0,0,0))
02293         pads.append(ROOT.TPad("p11","",0.00,0.34,0.07,0.56,0,0,0))
02294         pads.append(ROOT.TPad("p12","",0.07,0.34,0.34,0.56,0,0,0))
02295         pads.append(ROOT.TPad("p13","",0.34,0.34,0.56,0.56,0,0,0))
02296         pads.append(ROOT.TPad("p14","",0.56,0.34,0.78,0.56,0,0,0))
02297         pads.append(ROOT.TPad("p15","",0.78,0.34,1.00,0.56,0,0,0))
02298         pads.append(ROOT.TPad("p16","",0.00,0.07,0.07,0.34,0,0,0))
02299         pads.append(ROOT.TPad("p17","",0.07,0.07,0.34,0.34,0,0,0))
02300         pads.append(ROOT.TPad("p18","",0.34,0.07,0.56,0.34,0,0,0))
02301         pads.append(ROOT.TPad("p19","",0.56,0.07,0.78,0.34,0,0,0))
02302         pads.append(ROOT.TPad("p20","",0.78,0.07,1.00,0.34,0,0,0))
02303         pads.append(ROOT.TPad("p21","",0.00,0.00,0.07,0.07,0,0,0))
02304         pads.append(ROOT.TPad("p22","",0.07,0.00,0.34,0.07,0,0,0))
02305         pads.append(ROOT.TPad("p23","",0.34,0.00,0.56,0.07,0,0,0))
02306         pads.append(ROOT.TPad("p24","",0.56,0.00,0.78,0.07,0,0,0))
02307         pads.append(ROOT.TPad("p25","",0.78,0.00,1.00,0.07,0,0,0))
02308         for p in pads:
02309           if not not p:
02310             p.Draw()
02311             ROOT.SetOwnership(p,False)
02312 
02313         label1 = ROOT.TPaveLabel(0, 0, 1, 1, "x residuals (mm)","")
02314         label2 = ROOT.TPaveLabel(0, 0, 1, 1, "y residuals (mm)","")
02315         label3 = ROOT.TPaveLabel(0, 0, 1, 1, "dx/dz residuals (mrad)","")
02316         label4 = ROOT.TPaveLabel(0, 0, 1, 1, "dy/dz residuals (mrad)","")
02317         label5 = ROOT.TPaveLabel(0, 0, 1, 1, "x position (cm)","")
02318         label6 = ROOT.TPaveLabel(0, 0, 1, 1, "y position (cm)","")
02319         label7 = ROOT.TPaveLabel(0, 0, 1, 1, "dx/dz angle (rad)","")
02320         label8 = ROOT.TPaveLabel(0, 0, 1, 1, "dy/dz angle (rad)","")
02321         label9 = ROOT.TPaveLabel(0, 0.85, 1, 1, getname(r),"NDC")
02322 
02323         for l in label1, label2, label3, label4, label5, label6, label7, label8, label9:
02324             l.SetBorderSize(0)
02325             l.SetFillColor(ROOT.kWhite)
02326 
02327         for l in label1, label2, label3, label4:
02328             l.SetTextAngle(90)
02329             l.SetTextSize(0.09)
02330         
02331         #label9.SetTextAngle(30)
02332         label9.SetTextSize(0.59)
02333 
02334         pads[1].cd(); label1.Draw()
02335         pads[6].cd(); label2.Draw()
02336         pads[11].cd(); label3.Draw()
02337         pads[16].cd(); label4.Draw()
02338         pads[22].cd(); label5.Draw()
02339         pads[23].cd(); label6.Draw()
02340         pads[24].cd(); label7.Draw()
02341         pads[25].cd(); label8.Draw()
02342 
02343         pads[2].SetRightMargin(1e-5)
02344         pads[2].SetBottomMargin(1e-5)
02345         pads[2].SetLeftMargin(0.17)
02346         pads[3].SetLeftMargin(1e-5)
02347         pads[3].SetRightMargin(1e-5)
02348         pads[3].SetBottomMargin(1e-5)
02349         pads[4].SetLeftMargin(1e-5)
02350         pads[4].SetRightMargin(1e-5)
02351         pads[4].SetBottomMargin(1e-5)
02352         pads[5].SetLeftMargin(1e-5)
02353         pads[5].SetBottomMargin(1e-5)
02354 
02355         pads[7].SetRightMargin(1e-5)
02356         pads[7].SetBottomMargin(1e-5)
02357         pads[7].SetTopMargin(1e-5)
02358         pads[7].SetLeftMargin(0.17)
02359         pads[8].SetLeftMargin(1e-5)
02360         pads[8].SetRightMargin(1e-5)
02361         pads[8].SetBottomMargin(1e-5)
02362         pads[8].SetTopMargin(1e-5)
02363         pads[9].SetLeftMargin(1e-5)
02364         pads[9].SetRightMargin(1e-5)
02365         pads[9].SetBottomMargin(1e-5)
02366         pads[9].SetTopMargin(1e-5)
02367         pads[10].SetLeftMargin(1e-5)
02368         pads[10].SetBottomMargin(1e-5)
02369         pads[10].SetTopMargin(1e-5)
02370 
02371         pads[12].SetRightMargin(1e-5)
02372         pads[12].SetBottomMargin(1e-5)
02373         pads[12].SetTopMargin(1e-5)
02374         pads[12].SetLeftMargin(0.17)
02375         pads[13].SetLeftMargin(1e-5)
02376         pads[13].SetRightMargin(1e-5)
02377         pads[13].SetBottomMargin(1e-5)
02378         pads[13].SetTopMargin(1e-5)
02379         pads[14].SetLeftMargin(1e-5)
02380         pads[14].SetRightMargin(1e-5)
02381         pads[14].SetBottomMargin(1e-5)
02382         pads[14].SetTopMargin(1e-5)
02383         pads[15].SetLeftMargin(1e-5)
02384         pads[15].SetBottomMargin(1e-5)
02385         pads[15].SetTopMargin(1e-5)
02386 
02387         pads[17].SetRightMargin(1e-5)
02388         pads[17].SetTopMargin(1e-5)
02389         pads[17].SetLeftMargin(0.17)
02390         pads[18].SetLeftMargin(1e-5)
02391         pads[18].SetRightMargin(1e-5)
02392         pads[18].SetTopMargin(1e-5)
02393         pads[19].SetLeftMargin(1e-5)
02394         pads[19].SetRightMargin(1e-5)
02395         pads[19].SetTopMargin(1e-5)
02396         pads[20].SetLeftMargin(1e-5)
02397         pads[20].SetTopMargin(1e-5)
02398         
02399         chamber_x_trackx.GetXaxis().SetLabelColor(ROOT.kWhite)
02400         chamber_x_tracky.GetXaxis().SetLabelColor(ROOT.kWhite)
02401         chamber_x_tracky.GetYaxis().SetLabelColor(ROOT.kWhite)
02402         chamber_x_trackdxdz.GetXaxis().SetLabelColor(ROOT.kWhite)
02403         chamber_x_trackdxdz.GetYaxis().SetLabelColor(ROOT.kWhite)
02404         chamber_x_trackdydz.GetXaxis().SetLabelColor(ROOT.kWhite)
02405         chamber_x_trackdydz.GetYaxis().SetLabelColor(ROOT.kWhite)
02406         chamber_y_trackx.GetXaxis().SetLabelColor(ROOT.kWhite)
02407         chamber_y_tracky.GetXaxis().SetLabelColor(ROOT.kWhite)
02408         chamber_y_tracky.GetYaxis().SetLabelColor(ROOT.kWhite)
02409         chamber_y_trackdxdz.GetXaxis().SetLabelColor(ROOT.kWhite)
02410         chamber_y_trackdxdz.GetYaxis().SetLabelColor(ROOT.kWhite)
02411         chamber_y_trackdydz.GetXaxis().SetLabelColor(ROOT.kWhite)
02412         chamber_y_trackdydz.GetYaxis().SetLabelColor(ROOT.kWhite)
02413         chamber_dxdz_trackx.GetXaxis().SetLabelColor(ROOT.kWhite)
02414         chamber_dxdz_tracky.GetXaxis().SetLabelColor(ROOT.kWhite)
02415         chamber_dxdz_tracky.GetYaxis().SetLabelColor(ROOT.kWhite)
02416         chamber_dxdz_trackdxdz.GetXaxis().SetLabelColor(ROOT.kWhite)
02417         chamber_dxdz_trackdxdz.GetYaxis().SetLabelColor(ROOT.kWhite)
02418         chamber_dxdz_trackdydz.GetXaxis().SetLabelColor(ROOT.kWhite)
02419         chamber_dxdz_trackdydz.GetYaxis().SetLabelColor(ROOT.kWhite)
02420         
02421         # chamber_dydz_trackx
02422         chamber_dydz_tracky.GetYaxis().SetLabelColor(ROOT.kWhite)
02423         chamber_dydz_trackdxdz.GetYaxis().SetLabelColor(ROOT.kWhite)
02424         chamber_dydz_trackdydz.GetYaxis().SetLabelColor(ROOT.kWhite)
02425 
02426         pads[2].cd()
02427         chamber_x_trackx.Draw("e1")
02428         if not suppressblue: chamber_x_trackx_fit2.Draw("samel")
02429         chamber_x_trackx_fit.Draw("samel")
02430         #label99 = ROOT.TPaveLabel(0, 0.8, 1, 1, getname(r),"NDC")
02431         print getname(r)
02432         #label99 = ROOT.TPaveLabel(0, 0.8, 1, 1, "aaa","NDC")
02433         label9.Draw()
02434         #pads[2].Modified()
02435         
02436         pads[3].cd()
02437         chamber_x_tracky.Draw("e1")
02438         if not suppressblue: chamber_x_tracky_fit2.Draw("samel")
02439         chamber_x_tracky_fit.Draw("samel")
02440         
02441         pads[4].cd()
02442         chamber_x_trackdxdz.Draw("e1")
02443         if not suppressblue: chamber_x_trackdxdz_fit2.Draw("samel")
02444         chamber_x_trackdxdz_fit.Draw("samel")
02445         
02446         pads[5].cd()
02447         chamber_x_trackdydz.Draw("e1")
02448         if not suppressblue: chamber_x_trackdydz_fit2.Draw("samel")
02449         chamber_x_trackdydz_fit.Draw("samel")
02450         
02451         pads[7].cd()
02452         chamber_y_trackx.Draw("e1")
02453         if not suppressblue: chamber_y_trackx_fit2.Draw("samel")
02454         chamber_y_trackx_fit.Draw("samel")
02455         
02456         pads[8].cd()
02457         chamber_y_tracky.Draw("e1")
02458         if not suppressblue: chamber_y_tracky_fit2.Draw("samel")
02459         chamber_y_tracky_fit.Draw("samel")
02460         
02461         pads[9].cd()
02462         chamber_y_trackdxdz.Draw("e1")
02463         if not suppressblue: chamber_y_trackdxdz_fit2.Draw("samel")
02464         chamber_y_trackdxdz_fit.Draw("samel")
02465         
02466         pads[10].cd()
02467         chamber_y_trackdydz.Draw("e1")
02468         if not suppressblue: chamber_y_trackdydz_fit2.Draw("samel")
02469         chamber_y_trackdydz_fit.Draw("samel")
02470         
02471         pads[12].cd()
02472         chamber_dxdz_trackx.Draw("e1")
02473         if not suppressblue: chamber_dxdz_trackx_fit2.Draw("samel")
02474         chamber_dxdz_trackx_fit.Draw("samel")
02475         
02476         pads[13].cd()
02477         chamber_dxdz_tracky.Draw("e1")
02478         if not suppressblue: chamber_dxdz_tracky_fit2.Draw("samel")
02479         chamber_dxdz_tracky_fit.Draw("samel")
02480         
02481         pads[14].cd()
02482         chamber_dxdz_trackdxdz.Draw("e1")
02483         if not suppressblue: chamber_dxdz_trackdxdz_fit2.Draw("samel")
02484         chamber_dxdz_trackdxdz_fit.Draw("samel")
02485         
02486         pads[15].cd()
02487         chamber_dxdz_trackdydz.Draw("e1")
02488         if not suppressblue: chamber_dxdz_trackdydz_fit2.Draw("samel")
02489         chamber_dxdz_trackdydz_fit.Draw("samel")
02490         
02491         pads[17].cd()
02492         chamber_dydz_trackx.Draw("e1")
02493         if not suppressblue: chamber_dydz_trackx_fit2.Draw("samel")
02494         chamber_dydz_trackx_fit.Draw("samel")
02495         
02496         pads[18].cd()
02497         chamber_dydz_tracky.Draw("e1")
02498         if not suppressblue: chamber_dydz_tracky_fit2.Draw("samel")
02499         chamber_dydz_tracky_fit.Draw("samel")
02500         
02501         pads[19].cd()
02502         chamber_dydz_trackdxdz.Draw("e1")
02503         if not suppressblue: chamber_dydz_trackdxdz_fit2.Draw("samel")
02504         chamber_dydz_trackdxdz_fit.Draw("samel")
02505         
02506         pads[20].cd()
02507         chamber_dydz_trackdydz.Draw("e1")
02508         if not suppressblue: chamber_dydz_trackdydz_fit2.Draw("samel")
02509         chamber_dydz_trackdydz_fit.Draw("samel")
02510 
02511     else:
02512         c1.Clear()
02513         #c1.Divide(5, 3, 1e-5, 1e-5)
02514         pads = [None]
02515         pads.append(ROOT.TPad("p1" ,"",0.00,0.55,0.07,1.00,0,0,0))
02516         pads.append(ROOT.TPad("p2" ,"",0.07,0.55,0.34,1.00,0,0,0))
02517         pads.append(ROOT.TPad("p3" ,"",0.34,0.55,0.56,1.00,0,0,0))
02518         pads.append(ROOT.TPad("p4" ,"",0.56,0.55,0.78,1.00,0,0,0))
02519         pads.append(ROOT.TPad("p5" ,"",0.78,0.55,1.00,1.00,0,0,0))
02520         pads.append(ROOT.TPad("p6" ,"",0.00,0.1,0.07,0.55,0,0,0))
02521         pads.append(ROOT.TPad("p7" ,"",0.07,0.1,0.34,0.55,0,0,0))
02522         pads.append(ROOT.TPad("p8" ,"",0.34,0.1,0.56,0.55,0,0,0))
02523         pads.append(ROOT.TPad("p9" ,"",0.56,0.1,0.78,0.55,0,0,0))
02524         pads.append(ROOT.TPad("p10","",0.78,0.1,1.00,0.55,0,0,0))
02525         pads.append(ROOT.TPad("p11","",0.00,0.,0.07,0.1,0,0,0))
02526         pads.append(ROOT.TPad("p12","",0.07,0.,0.34,0.1,0,0,0))
02527         pads.append(ROOT.TPad("p13","",0.34,0.,0.56,0.1,0,0,0))
02528         pads.append(ROOT.TPad("p14","",0.56,0.,0.78,0.1,0,0,0))
02529         pads.append(ROOT.TPad("p15","",0.78,0.,1.00,0.1,0,0,0))
02530         for p in pads:
02531           if not not p:
02532             p.Draw()
02533             ROOT.SetOwnership(p,False)
02534 
02535         label1 = ROOT.TPaveLabel(0, 0, 1, 1, "x residuals (mm)")
02536         label2 = ROOT.TPaveLabel(0, 0, 1, 1, "dx/dz residuals (mrad)")
02537         label3 = ROOT.TPaveLabel(0, 0.3, 1, 1, "x position (cm)")
02538         label4 = ROOT.TPaveLabel(0, 0.3, 1, 1, "y position (cm)")
02539         label5 = ROOT.TPaveLabel(0, 0.3, 1, 1, "dx/dz angle (rad)")
02540         label6 = ROOT.TPaveLabel(0, 0.3, 1, 1, "dy/dz angle (rad)")
02541         label9 = ROOT.TPaveLabel(0, 0.85, 1, 1, getname(r),"NDC")
02542 
02543         if name[0:2] == "ME":
02544             label1 = ROOT.TPaveLabel(0, 0, 1, 1, "r#phi residuals (mm)")
02545             label2 = ROOT.TPaveLabel(0, 0, 1, 1, "d(r#phi)/dz residuals (mrad)")
02546 
02547         for l in label1, label2, label3, label4, label5, label6, label9:
02548             l.SetBorderSize(0)
02549             l.SetFillColor(ROOT.kWhite)
02550 
02551         for l in label1, label2:
02552             l.SetTextAngle(90)
02553             l.SetTextSize(0.09)
02554 
02555         #label9.SetTextAngle(30)
02556         label9.SetTextSize(0.29)
02557 
02558         pads[1].cd(); label1.Draw()
02559         pads[6].cd(); label2.Draw()
02560         pads[12].cd(); label3.Draw()
02561         pads[13].cd(); label4.Draw()
02562         pads[14].cd(); label5.Draw()
02563         pads[15].cd(); label6.Draw()
02564         #pads[11].cd(); label9.Draw()
02565 
02566         pads[2].SetRightMargin(1e-5)
02567         pads[2].SetBottomMargin(1e-5)
02568         pads[3].SetLeftMargin(1e-5)
02569         pads[3].SetRightMargin(1e-5)
02570         pads[3].SetBottomMargin(1e-5)
02571         pads[4].SetLeftMargin(1e-5)
02572         pads[4].SetRightMargin(1e-5)
02573         pads[4].SetBottomMargin(1e-5)
02574         pads[5].SetLeftMargin(1e-5)
02575         pads[5].SetBottomMargin(1e-5)
02576 
02577         pads[7].SetRightMargin(1e-5)
02578         pads[7].SetTopMargin(1e-5)
02579         pads[8].SetLeftMargin(1e-5)
02580         pads[8].SetRightMargin(1e-5)
02581         pads[8].SetTopMargin(1e-5)
02582         pads[9].SetLeftMargin(1e-5)
02583         pads[9].SetRightMargin(1e-5)
02584         pads[9].SetTopMargin(1e-5)
02585         pads[10].SetLeftMargin(1e-5)
02586         pads[10].SetTopMargin(1e-5)
02587 
02588         chamber_x_trackx.GetXaxis().SetLabelColor(ROOT.kWhite)
02589         chamber_x_tracky.GetXaxis().SetLabelColor(ROOT.kWhite)
02590         chamber_x_tracky.GetYaxis().SetLabelColor(ROOT.kWhite)
02591         chamber_x_trackdxdz.GetXaxis().SetLabelColor(ROOT.kWhite)
02592         chamber_x_trackdxdz.GetYaxis().SetLabelColor(ROOT.kWhite)
02593         chamber_x_trackdydz.GetXaxis().SetLabelColor(ROOT.kWhite)
02594         chamber_x_trackdydz.GetYaxis().SetLabelColor(ROOT.kWhite)
02595         # chamber_dxdz_trackx
02596         chamber_dxdz_tracky.GetYaxis().SetLabelColor(ROOT.kWhite)
02597         chamber_dxdz_trackdxdz.GetYaxis().SetLabelColor(ROOT.kWhite)
02598         chamber_dxdz_trackdydz.GetYaxis().SetLabelColor(ROOT.kWhite)
02599 
02600         pads[2].cd()
02601         chamber_x_trackx.Draw("e1")
02602         if not suppressblue: chamber_x_trackx_fit2.Draw("samel")
02603         chamber_x_trackx_fit.Draw("samel")
02604         label9.Draw()
02605         
02606         pads[3].cd()
02607         chamber_x_tracky.Draw("e1")
02608         if not suppressblue: chamber_x_tracky_fit2.Draw("samel")
02609         chamber_x_tracky_fit.Draw("samel")
02610         
02611         pads[4].cd()
02612         chamber_x_trackdxdz.Draw("e1")
02613         if not suppressblue: chamber_x_trackdxdz_fit2.Draw("samel")
02614         chamber_x_trackdxdz_fit.Draw("samel")
02615         
02616         pads[5].cd()
02617         chamber_x_trackdydz.Draw("e1")
02618         if not suppressblue: chamber_x_trackdydz_fit2.Draw("samel")
02619         chamber_x_trackdydz_fit.Draw("samel")
02620         
02621         pads[7].cd()
02622         chamber_dxdz_trackx.Draw("e1")
02623         if not suppressblue: chamber_dxdz_trackx_fit2.Draw("samel")
02624         chamber_dxdz_trackx_fit.Draw("samel")
02625         
02626         pads[8].cd()
02627         chamber_dxdz_tracky.Draw("e1")
02628         if not suppressblue: chamber_dxdz_tracky_fit2.Draw("samel")
02629         chamber_dxdz_tracky_fit.Draw("samel")
02630         
02631         pads[9].cd()
02632         chamber_dxdz_trackdxdz.Draw("e1")
02633         if not suppressblue: chamber_dxdz_trackdxdz_fit2.Draw("samel")
02634         chamber_dxdz_trackdxdz_fit.Draw("samel")
02635         
02636         pads[10].cd()
02637         chamber_dxdz_trackdydz.Draw("e1")
02638         if not suppressblue: chamber_dxdz_trackdydz_fit2.Draw("samel")
02639         chamber_dxdz_trackdydz_fit.Draw("samel")
02640 
02641     tn = time.time()
02642     ddt[8] = 1./ddt[7]*((ddt[7]-1)*ddt[8] + tn-t1)
02643 
02644 ##################################################################################
02645 
02646 def segdiff(tfiles, component, pair, **args):
02647     tdrStyle.SetOptFit(1)
02648     tdrStyle.SetOptTitle(1)
02649     tdrStyle.SetTitleBorderSize(1)
02650     tdrStyle.SetTitleFontSize(0.05)
02651     tdrStyle.SetStatW(0.2)
02652     tdrStyle.SetStatY(0.9)
02653     tdrStyle.SetStatFontSize(0.06)
02654 
02655     if component[0:2] == "dt":
02656         wheel = args["wheel"]
02657         wheelletter = wheelLetter(wheel)
02658         sector = args["sector"]
02659         profname = "%s_%s_%02d_%s" % (component, wheelletter, sector, str(pair))
02660         posname = "pos" + profname
02661         negname = "neg" + profname
02662         #print profname
02663 
02664         station1 = int(str(pair)[0])
02665         station2 = int(str(pair)[1])
02666         phi1 = signConventions["DT", wheel, station1, sector][4]
02667         phi2 = signConventions["DT", wheel, station2, sector][4]
02668         if abs(phi1 - phi2) > 1.:
02669             if phi1 > phi2: phi1 -= 2.*pi
02670             else: phi1 += 2.*pi
02671         phi = (phi1 + phi2) / 2.
02672         while (phi < -pi): phi += 2.*pi
02673         while (phi > pi): phi -= 2.*pi
02674 
02675     elif component[0:3] == "csc":
02676         endcap = args["endcap"]
02677         if endcap=="m":
02678             endcapnum=2
02679             endcapsign="-"
02680         elif endcap=="p":
02681             endcapnum=1
02682             endcapsign="+"
02683         else: raise Exception
02684         
02685         ring = args["ring"]
02686         if ring>2 or ring<1: raise Exception
02687         station1 = int(str(pair)[0])
02688         station2 = int(str(pair)[1])
02689         if   ring==1: ringname="inner"
02690         elif ring==2: ringname="outer"
02691         else: raise Exception
02692         
02693         chamber = args["chamber"]
02694         if (ring==1 and chamber>18) or (ring==2 and chamber>36): raise Exception
02695         
02696         profname = "csc%s_%s_%s_%02d_%s" % (ringname,component[4:], endcap, chamber, str(pair))
02697         posname = "pos" + profname
02698         negname = "neg" + profname
02699         #print profname
02700 
02701         station1 = int(str(pair)[0])
02702         station2 = int(str(pair)[1])
02703         phi1 = signConventions["CSC", endcapnum, station1, ring, chamber][4]
02704         phi2 = signConventions["CSC", endcapnum, station1, ring, chamber][4]
02705         if abs(phi1 - phi2) > 1.:
02706             if phi1 > phi2: phi1 -= 2.*pi
02707             else: phi1 += 2.*pi
02708         phi = (phi1 + phi2) / 2.
02709         while (phi < -pi): phi += 2.*pi
02710         while (phi > pi): phi -= 2.*pi
02711 
02712     else: raise Exception
02713 
02714     if "window" in args: window = args["window"]
02715     else: window = 5.
02716 
02717     global tmpprof, tmppos, tmpneg
02718     pdir = "AlignmentMonitorSegmentDifferences/iter1/"
02719     tmpprof = tfiles[0].Get(pdir + profname).Clone()
02720     tmpprof.SetMarkerStyle(8)
02721     tmppos = tfiles[0].Get(pdir + posname).Clone()
02722     tmpneg = tfiles[0].Get(pdir + negname).Clone()
02723     for tfile in tfiles[1:]:
02724         tmpprof.Add(tfile.Get(pdir + profname))
02725         tmppos.Add(tfile.Get(pdir + posname))
02726         tmpneg.Add(tfile.Get(pdir + negname))
02727 
02728     for i in xrange(1, tmpprof.GetNbinsX()+1):
02729         if tmpprof.GetBinError(i) < 1e-5:
02730             tmpprof.SetBinError(i, 100.)
02731     tmpprof.SetAxisRange(-window, window, "Y")
02732 
02733     f = ROOT.TF1("p1", "[0] + [1]*x", tmpprof.GetBinLowEdge(1), -tmpprof.GetBinLowEdge(1))
02734     f.SetParameters((tmppos.GetMean() + tmpneg.GetMean())/2., 0.)
02735 
02736     tmpprof.SetXTitle("q/p_{T} (c/GeV)")
02737     if component == "dt13_resid":
02738         tmpprof.SetYTitle("#Deltax^{local} (mm)")
02739         tmppos.SetXTitle("#Deltax^{local} (mm)")
02740         tmpneg.SetXTitle("#Deltax^{local} (mm)")
02741         f.SetParNames("#Deltax^{local}_{0}", "Slope")
02742     if component == "dt13_slope":
02743         tmpprof.SetYTitle("#Deltadx/dz^{local} (mrad)")
02744         tmppos.SetXTitle("#Deltadx/dz^{local} (mrad)")
02745         tmpneg.SetXTitle("#Deltadx/dz^{local} (mrad)")
02746         f.SetParNames("#Deltadx/dz^{local}_{0}", "Slope")
02747     if component == "dt2_resid":
02748         tmpprof.SetYTitle("#Deltay^{local} (mm)")
02749         tmppos.SetXTitle("#Deltay^{local} (mm)")
02750         tmpneg.SetXTitle("#Deltay^{local} (mm)")
02751         f.SetParNames("#Deltay^{local}_{0}", "Slope")
02752     if component == "dt2_slope":
02753         tmpprof.SetYTitle("#Deltady/dz^{local} (mrad)")
02754         tmppos.SetXTitle("#Deltady/dz^{local} (mrad)")
02755         tmpneg.SetXTitle("#Deltady/dz^{local} (mrad)")
02756         f.SetParNames("#Deltady/dz^{local}_{0}", "Slope")
02757     if component == "csc_resid":
02758         tmpprof.SetYTitle("#Delta(r#phi)^{local} (mm)")
02759         tmppos.SetXTitle("#Delta(r#phi)^{local} (mm)")
02760         tmpneg.SetXTitle("#Delta(r#phi)^{local} (mm)")
02761         f.SetParNames("#Delta(r#phi)^{local}_{0}", "Slope")
02762     if component == "csc_slope":
02763         tmpprof.SetYTitle("#Deltad(r#phi)/dz^{local} (mrad)")
02764         tmppos.SetXTitle("#Deltad(r#phi)/dz^{local} (mrad)")
02765         tmpneg.SetXTitle("#Deltad(r#phi)/dz^{local} (mrad)")
02766         f.SetParNames("#Deltad(r#phi)/dz^{local}_{0}", "Slope")
02767     
02768     tmpprof.GetXaxis().CenterTitle()
02769     tmpprof.GetYaxis().CenterTitle()
02770     tmppos.GetXaxis().CenterTitle()
02771     tmpneg.GetXaxis().CenterTitle()
02772     if component[0:2] == "dt":
02773         tmpprof.SetTitle("MB%d - MB%d, wheel %d, sector %02d" % (station1, station2, int(wheel), int(sector)))
02774     elif component[0:3] == "csc":
02775         tmpprof.SetTitle("ME%d - ME%d, for ME%s%d/%d/%d" % (station1, station2, endcapsign, station2, ring, chamber))
02776     else: raise Exception
02777 
02778     tmppos.SetTitle("Positive muons")
02779     tmpneg.SetTitle("Negative muons")
02780 
02781     c1.Clear()
02782     c1.Divide(2, 1)
02783     c1.GetPad(1).cd()
02784     fit1 = tmpprof.Fit("p1", "q")
02785     tmpprof.Draw("e1")
02786     c1.GetPad(2).cd()
02787     c1.GetPad(2).Divide(1, 2)
02788     c1.GetPad(2).GetPad(1).cd()
02789     tmppos.Draw()
02790     f = ROOT.TF1("gausR", "[0]*exp(-(x - [1])**2 / 2. / [2]**2) / sqrt(2.*3.1415926) / [2]", 
02791                  tmppos.GetMean() - tmppos.GetRMS(), tmppos.GetMean() + tmppos.GetRMS())
02792     f.SetParameters(tmppos.GetEntries() * ((10. - -10.)/100.), tmppos.GetMean(), tmppos.GetRMS())
02793     f.SetParNames("Constant", "Mean", "Sigma")
02794     fit2 = tmppos.Fit("gausR", "qR")
02795     c1.GetPad(2).GetPad(2).cd()
02796     tmpneg.Draw()
02797     f = ROOT.TF1("gausR", "[0]*exp(-(x - [1])**2 / 2. / [2]**2) / sqrt(2.*3.1415926) / [2]", 
02798                  tmpneg.GetMean() - tmpneg.GetRMS(), tmpneg.GetMean() + tmpneg.GetRMS())
02799     f.SetParameters(tmpneg.GetEntries() * ((10. - -10.)/100.), tmpneg.GetMean(), tmpneg.GetRMS())
02800     f.SetParNames("Constant", "Mean", "Sigma")
02801     fit3 = tmpneg.Fit("gausR", "qR")
02802 
02803     fitresult1 = None, None
02804     if fit1 == 0:
02805         fitresult1 = tmpprof.GetFunction("p1").GetParameter(0), tmpprof.GetFunction("p1").GetParError(0)
02806     fitresult2 = None, None
02807     if fit2 == 0 and fit3 == 0:
02808         fitresult2 = (tmppos.GetFunction("gausR").GetParameter(1) + tmpneg.GetFunction("gausR").GetParameter(1)) / 2., \
02809                      sqrt(tmppos.GetFunction("gausR").GetParError(1)**2 + tmpneg.GetFunction("gausR").GetParError(1)**2) / 2.
02810     return phi, fitresult1[0], fitresult1[1], fitresult2[0], fitresult2[1], fit1, fit2, fit3
02811 
02812 def segdiffvsphi(tfiles, reports, component, wheel, window=5., excludesectors=()):
02813     tdrStyle.SetOptTitle(1)
02814     tdrStyle.SetTitleBorderSize(1)
02815     tdrStyle.SetTitleFontSize(0.05)
02816 
02817     global htemp, gtemp_12, gtemp2_12, gtemp_23, gtemp2_23, gtemp_34, gtemp2_34, tlegend
02818     htemp = ROOT.TH1F("htemp", "", 1, -pi, pi)
02819     gtemp_12_phi, gtemp_12_val, gtemp_12_err, gtemp_12_val2, gtemp_12_err2 = [], [], [], [], []
02820     gtemp_23_phi, gtemp_23_val, gtemp_23_err, gtemp_23_val2, gtemp_23_err2 = [], [], [], [], []
02821     gtemp_34_phi, gtemp_34_val, gtemp_34_err, gtemp_34_val2, gtemp_34_err2 = [], [], [], [], []
02822     for sector in xrange(1, 12+1):
02823         #print sector
02824         r1_found, r2_found, r3_found, r4_found = False, False, False, False
02825         for r1 in reports:
02826             if r1.postal_address == ("DT", wheel, 1, sector):
02827                 r1_found = True
02828                 break
02829         for r2 in reports:
02830             if r2.postal_address == ("DT", wheel, 2, sector):
02831                 r2_found = True
02832                 break
02833         for r3 in reports:
02834             if r3.postal_address == ("DT", wheel, 3, sector):
02835                 r3_found = True
02836                 break
02837         for r4 in reports:
02838             if r4.postal_address == ("DT", wheel, 4, sector):
02839                 r4_found = True
02840                 break
02841         #print "rfounds: ", r1_found, r2_found, r3_found, r4_found
02842         
02843         if sector not in excludesectors:
02844             if r1_found and r2_found and r1.status == "PASS" and r2.status == "PASS":
02845                 phi, val, err, val2, err2, fit1, fit2, fit3 = segdiff(tfiles, component, 12, wheel=wheel, sector=sector)
02846                 #print "segdif 12", phi, val, err, val2, err2, fit1, fit2, fit3
02847                 if fit1 == 0 and fit2 == 0 and fit3 == 0:
02848                 #if fit1 == 0 and fit2 == 0:
02849                     gtemp_12_phi.append(phi)
02850                     gtemp_12_val.append(val)
02851                     gtemp_12_err.append(err)
02852                     gtemp_12_val2.append(val2)
02853                     gtemp_12_err2.append(err2)
02854             if r2_found and r3_found and r2.status == "PASS" and r3.status == "PASS":
02855                 phi, val, err, val2, err2, fit1, fit2, fit3 = segdiff(tfiles, component, 23, wheel=wheel, sector=sector)
02856                 #print "segdif 23", phi, val, err, val2, err2, fit1, fit2, fit3
02857                 if fit1 == 0 and fit2 == 0 and fit3 == 0:
02858                 #if fit1 == 0 and fit2 == 0:
02859                     gtemp_23_phi.append(phi)
02860                     gtemp_23_val.append(val)
02861                     gtemp_23_err.append(err)
02862                     gtemp_23_val2.append(val2)
02863                     gtemp_23_err2.append(err2)
02864             if component[:4] == "dt13":
02865                 if r3_found and r4_found and r3.status == "PASS" and r4.status == "PASS":
02866                     phi, val, err, val2, err2, fit1, fit2, fit3 = segdiff(tfiles, component, 34, wheel=wheel, sector=sector)
02867                     #print "segdif 34", phi, val, err, val2, err2, fit1, fit2, fit3
02868                     if fit1 == 0 and fit2 == 0 and fit3 == 0:
02869                     #if fit1 == 0 and fit2 == 0:
02870                         gtemp_34_phi.append(phi)
02871                         gtemp_34_val.append(val)
02872                         gtemp_34_err.append(err)
02873                         gtemp_34_val2.append(val2)
02874                         gtemp_34_err2.append(err2)
02875 
02876     #print "len(gtemp_12_phi) ", len(gtemp_12_phi)
02877     #print "len(gtemp_23_phi) ",len(gtemp_23_phi)
02878     #print "len(gtemp_34_phi) ",len(gtemp_34_phi)
02879     if len(gtemp_12_phi) > 0:
02880         gtemp_12 = ROOT.TGraphErrors(len(gtemp_12_phi), array.array("d", gtemp_12_phi), array.array("d", gtemp_12_val), 
02881                                      array.array("d", [0.] * len(gtemp_12_phi)), array.array("d", gtemp_12_err))
02882         gtemp2_12 = ROOT.TGraphErrors(len(gtemp_12_phi), array.array("d", gtemp_12_phi), array.array("d", gtemp_12_val2), 
02883                                       array.array("d", [0.] * len(gtemp_12_phi)), array.array("d", gtemp_12_err2))
02884     if len(gtemp_23_phi) > 0:
02885         gtemp_23 = ROOT.TGraphErrors(len(gtemp_23_phi), array.array("d", gtemp_23_phi), array.array("d", gtemp_23_val), 
02886                                      array.array("d", [0.] * len(gtemp_23_phi)), array.array("d", gtemp_23_err))
02887         gtemp2_23 = ROOT.TGraphErrors(len(gtemp_23_phi), array.array("d", gtemp_23_phi), array.array("d", gtemp_23_val2), 
02888                                       array.array("d", [0.] * len(gtemp_23_phi)), array.array("d", gtemp_23_err2))
02889     if len(gtemp_34_phi) > 0:
02890         gtemp_34 = ROOT.TGraphErrors(len(gtemp_34_phi), array.array("d", gtemp_34_phi), array.array("d", gtemp_34_val), 
02891                                      array.array("d", [0.] * len(gtemp_34_phi)), array.array("d", gtemp_34_err))
02892         gtemp2_34 = ROOT.TGraphErrors(len(gtemp_34_phi), array.array("d", gtemp_34_phi), array.array("d", gtemp_34_val2), 
02893                                       array.array("d", [0.] * len(gtemp_34_phi)), array.array("d", gtemp_34_err2))
02894 
02895     if len(gtemp_12_phi) > 0:
02896         gtemp_12.SetMarkerStyle(20);  gtemp_12.SetMarkerSize(1.);  
02897         gtemp_12.SetMarkerColor(ROOT.kBlue);  gtemp_12.SetLineColor(ROOT.kBlue)
02898         gtemp2_12.SetMarkerStyle(24); gtemp2_12.SetMarkerSize(1.); 
02899         gtemp2_12.SetMarkerColor(ROOT.kBlue); gtemp2_12.SetLineColor(ROOT.kBlue)
02900     if len(gtemp_23_phi) > 0:
02901         gtemp_23.SetMarkerStyle(21);  gtemp_23.SetMarkerSize(1.);  
02902         gtemp_23.SetMarkerColor(ROOT.kRed);   gtemp_23.SetLineColor(ROOT.kRed)
02903         gtemp2_23.SetMarkerStyle(25); gtemp2_23.SetMarkerSize(1.); 
02904         gtemp2_23.SetMarkerColor(ROOT.kRed);  gtemp2_23.SetLineColor(ROOT.kRed)
02905     if len(gtemp_34_phi) > 0 and component[:4] == "dt13":
02906         gtemp_34.SetMarkerStyle(22);  gtemp_34.SetMarkerSize(1.25);  
02907         gtemp_34.SetMarkerColor(ROOT.kGreen+2);  gtemp_34.SetLineColor(ROOT.kGreen+2)
02908         gtemp2_34.SetMarkerStyle(26); gtemp2_34.SetMarkerSize(1.25); 
02909         gtemp2_34.SetMarkerColor(ROOT.kGreen+2); gtemp2_34.SetLineColor(ROOT.kGreen+2)
02910 
02911     if wheel == 0: htemp.SetTitle("Wheel %d" % wheel)
02912     else: htemp.SetTitle("Wheel %+d" % wheel)
02913     htemp.SetAxisRange(-window, window, "Y")
02914     htemp.SetXTitle("Average #phi of pair (rad)")
02915     if component == "dt13_resid": htemp.SetYTitle("#Deltax^{local} (mm)")
02916     if component == "dt13_slope": htemp.SetYTitle("#Deltadx/dz^{local} (mrad)")
02917     if component == "dt2_resid": htemp.SetYTitle("#Deltay^{local} (mm)")
02918     if component == "dt2_slope": htemp.SetYTitle("#Deltady/dz^{local} (mrad)")
02919     htemp.GetXaxis().CenterTitle()
02920     htemp.GetYaxis().CenterTitle()
02921     htemp.GetYaxis().SetTitleOffset(0.75)
02922 
02923     c1.Clear()
02924     htemp.Draw()
02925     if len(gtemp_12_phi) > 0:
02926         gtemp_12.Draw("p")
02927         gtemp2_12.Draw("p")
02928     if len(gtemp_23_phi) > 0:
02929         gtemp_23.Draw("p")
02930         gtemp2_23.Draw("p")
02931     if len(gtemp_34_phi) > 0:
02932         gtemp_34.Draw("p")
02933         gtemp2_34.Draw("p")
02934 
02935     tlegend = ROOT.TLegend(0.5, 0.72, 0.9, 0.92)
02936     tlegend.SetBorderSize(0)
02937     tlegend.SetFillColor(ROOT.kWhite)
02938     if len(gtemp_12_phi) > 0:
02939         tlegend.AddEntry(gtemp_12, "MB1 - MB2 (mean: %4.2f, RMS: %4.2f)" % (mean(gtemp_12_val), stdev(gtemp_12_val)), "pl")
02940     if len(gtemp_23_phi) > 0:
02941         tlegend.AddEntry(gtemp_23, "MB2 - MB3 (mean: %4.2f, RMS: %4.2f)" % (mean(gtemp_23_val), stdev(gtemp_23_val)), "pl")
02942     if len(gtemp_34_phi) > 0:
02943         tlegend.AddEntry(gtemp_34, "MB3 - MB4 (mean: %4.2f, RMS: %4.2f)" % (mean(gtemp_34_val), stdev(gtemp_34_val)), "pl")
02944     if len(gtemp_12_phi) > 0:
02945         tlegend.AddEntry(gtemp_12, "total mean: %4.2f, total RMS: %4.2f" % \
02946                                    (mean(gtemp_12_val + gtemp_23_val + gtemp_34_val), 
02947                                    stdev(gtemp_12_val + gtemp_23_val + gtemp_34_val)), "")
02948     tlegend.Draw()
02949 
02950
/data/refman/pasoursint/CMSSW_4_1_8_patch13/src/Alignment/MuonAlignmentAlgorithms/scripts/plotscripts.py
