python/mpsvalidate/timeStructure.py

0001 ##########################################################################
0002 # Creates a histogram where the the names of the structures are present
0003 # as humanreadable text. Multiple MillePedeUser TTrees are used to
0004 # get a time dependent plot.
0005 ##
0006
0007 from builtins import range
0008 import logging
0009
0010 import ROOT
0011 ROOT.PyConfig.IgnoreCommandLineOptions = True
0012 ROOT.gROOT.SetBatch()
0013
0014 import Alignment.MillePedeAlignmentAlgorithm.mpsvalidate.style as mpsv_style
0015 import Alignment.MillePedeAlignmentAlgorithm.mpsvalidate.classes as mpsv_classes
0016
0017
0018 def plot(treeFile, alignables, config):
0019     logger = logging.getLogger("mpsvalidate")
0020
0021     for mode in ["xyz", "rot"]:
0022
0023         time = mpsv_classes.PlotData(mode)
0024
0025         # list of all avaible TTrees
0026         listMillePedeUser = []
0027         MillePedeUser = []
0028         for i in range(config.firsttree, 101):
0029             if (treeFile.GetListOfKeys().Contains("MillePedeUser_{0}".format(i))):
0030                 listMillePedeUser.append(i)
0031
0032         # load MillePedeUser_X TTrees
0033         for i in listMillePedeUser:
0034             MillePedeUser.append(treeFile.Get("MillePedeUser_{0}".format(i)))
0035
0036         ######################################################################
0037         # remove TTrees without results
0038         #
0039
0040         # check if there is a TTree without any results
0041         # therefor search for the first alignable
0042         first = 0
0043         newlistMillePedeUser = []
0044         # find first alignable
0045         for line in MillePedeUser[0]:
0046             if (line.ObjId != 1 and any(abs(line.Par[time.data[i]]) != 999999 for i in [0, 1, 2])):
0047                 first = line.Id
0048                 newlistMillePedeUser.append(config.firsttree)
0049                 break
0050
0051         # check the following TTrees
0052         for ttreeNumber, ttree in enumerate(MillePedeUser[1:]):
0053             for line in ttree:
0054                 if (line.Id == first):
0055                     if (any(abs(line.Par[time.data[i]]) != 999999 for i in [0, 1, 2])):
0056                         # note that the first tree was checked
0057                         newlistMillePedeUser.append(
0058                             ttreeNumber + config.firsttree + 1)
0059                     break
0060
0061         listMillePedeUser = newlistMillePedeUser
0062
0063         # reload MillePedeUser_X TTrees
0064         MillePedeUser = []
0065         for i in listMillePedeUser:
0066             MillePedeUser.append(treeFile.Get("MillePedeUser_{0}".format(i)))
0067
0068         if not listMillePedeUser:
0069             logger.error("Timeplots: no TTrees found")
0070             return
0071
0072         if not MillePedeUser:
0073             logger.error("Timeplots: no TTree could be opened")
0074             return
0075
0076         ######################################################################
0077         # initialize data hierarchy
0078         #
0079
0080         plots = []
0081         # objids which were found in the TTree
0082         objids = []
0083         obj_names = []
0084
0085         # loop over first tree to initialize
0086         for line in MillePedeUser[0]:
0087             if (line.ObjId != 1 and any(abs(line.Par[time.data[i]]) != 999999 for i in [0, 1, 2])):
0088                 plots.append(mpsv_classes.PlotData(mode))
0089
0090                 # new objid?
0091                 if (line.ObjId not in objids):
0092                     objids.append(line.ObjId)
0093                     obj_names.append(str(line.Name))
0094
0095                 # initialize histograms
0096                 for i in range(3):
0097                     plots[-1].histo.append(ROOT.TH1F("Time Structure {0} {1} {2} {3}".format(mode, str(line.Name),
0098                         len(plots), i), "", len(listMillePedeUser), 0, len(listMillePedeUser)))
0099                     plots[-1].label = line.Id
0100                     plots[-1].objid = line.ObjId
0101
0102                     if (time.unit!=""):
0103                         plots[-1].histo[i].SetYTitle("#Delta"+time.xyz[i]+" ["+time.unit+"]")
0104                     else:
0105                         plots[-1].histo[i].SetYTitle("#Delta"+time.xyz[i])
0106                     plots[-1].histo[i].SetXTitle("IOV")
0107                     plots[-1].histo[i].SetStats(0)
0108                     plots[-1].histo[i].SetMarkerStyle(21)
0109                     # bigger labels for the text
0110                     plots[-1].histo[i].GetXaxis().SetLabelSize(0.08)
0111                     plots[-1].histo[i].GetYaxis().SetTitleOffset(1.6)
0112
0113         ######################################################################
0114         # fill histogram
0115         #
0116
0117         # loop over TTrees
0118         for treeNumber, tree in enumerate(MillePedeUser):
0119             for line in tree:
0120                 if (line.ObjId != 1 and any(abs(line.Par[time.data[i]]) != 999999 for i in [0, 1, 2])):
0121                     # find the right plot
0122                     for plot in plots:
0123                         if (plot.label == line.Id):
0124                             for i in range(3):
0125                                 # note that the first bin is referenced by 1
0126                                 plot.histo[i].GetXaxis().SetBinLabel(
0127                                     treeNumber + 1, str(listMillePedeUser[treeNumber]))
0128                                 # transform xyz data from cm to #mu m
0129                                 if (mode == "xyz"):
0130                                     plot.histo[i].SetBinContent(
0131                                         treeNumber + 1, 10000 * line.Par[plot.data[i]])
0132                                 else:
0133                                     plot.histo[i].SetBinContent(
0134                                         treeNumber + 1, line.Par[plot.data[i]])
0135
0136         ######################################################################
0137         # find maximum/minimum
0138         #
0139
0140         maximum = [[0, 0, 0] for x in range(len(objids))]
0141         minimum = [[0, 0, 0] for x in range(len(objids))]
0142
0143         for index, objid in enumerate(objids):
0144             for plot in plots:
0145                 if (plot.objid == objid):
0146                     for i in range(3):
0147                         # maximum
0148                         if (plot.histo[i].GetMaximum() > maximum[index][i]):
0149                             maximum[index][i] = plot.histo[i].GetMaximum()
0150                         # minimum
0151                         if (plot.histo[i].GetMinimum() < minimum[index][i]):
0152                             minimum[index][i] = plot.histo[i].GetMinimum()
0153
0154         ######################################################################
0155         # make the plots
0156         #
0157
0158         # loop over all objids
0159         for index, objid in enumerate(objids):
0160
0161             canvas = ROOT.TCanvas("canvasTimeBigStrucutres_{0}_{1}".format(
0162                 mode, obj_names[index]), "Parameter", 300, 0, 800, 600)
0163             canvas.Divide(2, 2)
0164
0165             # add text
0166             title = ROOT.TPaveLabel(0.1, 0.8, 0.9, 0.9, "{0} over time {1}".format(
0167                 obj_names[index], mode))
0168
0169             legend = ROOT.TLegend(0.05, 0.1, 0.95, 0.75)
0170
0171             # draw on canvas
0172             canvas.cd(1)
0173             title.Draw()
0174
0175             # draw identification
0176             ident = mpsv_style.identification(config)
0177             ident.Draw()
0178
0179             # TGraph copies to hide outlier
0180             copy = []
0181
0182             # reset y range of first plot
0183             # two types of ranges
0184             for i in range(3):
0185                 for plot in plots:
0186                     if (plot.objid == objid):
0187                         # 1. show all
0188                         if config.rangemodeHL == "all":
0189                             plot.usedRange[i] = max(
0190                                 abs(maximum[index][i]), abs(minimum[index][i]))
0191
0192                         # 2. use given values
0193                         if (config.rangemodeHL == "given"):
0194                             # loop over coordinates
0195                             if mode == "xyz":
0196                                 valuelist = config.rangexyzHL
0197                             if mode == "rot":
0198                                 valuelist = config.rangerotHL
0199                             # loop over given values
0200                             # without last value
0201                             for value in valuelist:
0202                                 # maximum smaller than given value
0203                                 if max(abs(maximum[index][i]), abs(minimum[index][i])) < value:
0204                                     plot.usedRange[i] = value
0205                                     break
0206                             # if not possible, force highest
0207                             if (max(abs(maximum[index][i]), abs(minimum[index][i])) > valuelist[-1]):
0208                                 plot.usedRange[i] = valuelist[-1]
0209
0210             # draw plots on canvas
0211             for i in range(3):
0212                 canvas.cd(2 + i)
0213
0214                 number = 1
0215
0216                 for plot in plots:
0217                     if (plot.objid == objid):
0218                         # all the same range
0219                         if (config.samerangeHL == 1):
0220                             plot.usedRange[i] = max(map(abs, plot.usedRange))
0221
0222                         # set new range
0223                         plot.histo[i].GetYaxis(
0224                         ).SetRangeUser(-1.2 * abs(plot.usedRange[i]), 1.2 * abs(plot.usedRange[i]))
0225
0226                         plot.histo[i].SetLineColorAlpha(number + 2, 0.5)
0227                         plot.histo[i].SetMarkerColorAlpha(number + 2, 1)
0228
0229                         # option "AXIS" to only draw the axis
0230                         plot.histo[i].SetLineColor(0)
0231                         plot.histo[i].Draw("PSAME")
0232
0233                         # TGraph object to hide outlier
0234                         copy.append(ROOT.TGraph(plot.histo[i]))
0235                         # set the new range
0236                         copy[-1].SetMaximum(1.2 * abs(plot.usedRange[i]))
0237                         copy[-1].SetMinimum(-1.2 * abs(plot.usedRange[i]))
0238                         # draw the data
0239                         copy[-1].SetLineColorAlpha(number + 2, 0.5)
0240                         copy[-1].Draw("LPSAME")
0241
0242                         if (i == 0):
0243                             legend.AddEntry(
0244                                 plot.histo[i], "{0}".format(number))
0245                         number += 1
0246
0247             canvas.cd(1)
0248             legend.Draw()
0249
0250             canvas.Update()
0251
0252             # save as pdf
0253             canvas.Print("{0}/plots/pdf/timeStructures_{1}_{2}.pdf".format(
0254                 config.outputPath, mode, obj_names[index]))
0255
0256             # export as png
0257             image = ROOT.TImage.Create()
0258             image.FromPad(canvas)
0259             image.WriteImage("{0}/plots/png/timeStructures_{1}_{2}.png".format(
0260                 config.outputPath, mode, obj_names[index]))
0261
0262             # add to output list
0263             output = mpsv_classes.OutputData(plottype="time", name=obj_names[index],
0264                                              parameter=mode, filename="timeStructures_{0}_{1}".format(mode, obj_names[index]))
0265             config.outputList.append(output)