File indexing completed on 2023-03-17 10:56:42
0001
0002
0003 from __future__ import print_function
0004 import sys
0005 import os
0006 from ROOT import *
0007 from copy import deepcopy
0008 from array import array
0009
0010 gROOT.SetBatch()
0011
0012
0013
0014 def getFileInPath(rfile):
0015 import os
0016 for dir in os.environ['CMSSW_SEARCH_PATH'].split(":"):
0017 if os.path.exists(os.path.join(dir,rfile)): return os.path.join(dir,rfile)
0018 return None
0019
0020
0021
0022 inputFileName = "DQM_V0013_R000292154__StreamExpressCosmics__Commissioning2017-Express-v1__DQMIO.root"
0023 limitsFileName = "limits.dat"
0024 outputDirectoryName = "OUT/"
0025 minMaxFileName = "minmax.out"
0026
0027 detIDsFileName = getFileInPath('DQM/SiStripMonitorClient/data/detids.dat')
0028
0029 baseRootDirs = ["DQMData/Run 292154/PixelPhase1/Run summary/Phase1_MechanicalView"
0030 ,"DQMData/Run 292154/PixelPhase1/Run summary/Tracks"
0031 ]
0032
0033
0034 maxPxBarrel = 4
0035 maxPxForward = 3
0036 barrelLadderShift = [0, 14, 44, 90]
0037
0038 forwardDiskXShift = [25, 75, 125]
0039 forwardDiskYShift = 45;
0040
0041 plotWidth, plotHeight = 3000, 2000
0042 extremeBinsNum = 20
0043
0044 limits = ["num_digis 0.01 90 1 0",
0045 "num_clusters 0.01 25 1 0",
0046 "Trechitsize_y 0.01 10 0 0",
0047 "Trechitsize_x 0.01 10 0 0",
0048 "Tresidual_y 0.0000001 0.004 0 1",
0049 "Tresidual_x 0.0000001 0.004 0 1",
0050 "Tcharge 2000 80000 0 0",
0051 "Thitefficiency 0.95 1 0 0",
0052
0053 "Tnum_clusters_ontrack 0.01 15 1 0",
0054 "Tsize 0.01 15 0 0",
0055
0056 "adc 0.01 256 0 0",
0057 "charge 2000 80000 0 0",
0058 "size 0.01 15 0 0",]
0059
0060 class TH2PolyOfflineMaps:
0061
0062
0063
0064
0065
0066
0067
0068 def __TraverseDirTree(self, dir):
0069
0070 try:
0071 currPath = (dir.GetPath().split(":/"))[1]
0072 except:
0073 print("Exception raised: Path not found in the input file")
0074 return
0075
0076 for obj in dir.GetListOfKeys():
0077 if not obj.IsFolder():
0078 if obj.ReadObjectAny(TClass.GetClass("TH2")):
0079 th2 = deepcopy(obj.ReadObj())
0080 name = th2.GetName()
0081 if 6 < th2.GetNbinsX() < 10 and name.find("per") != -1 and name.find("Lumisection") == -1:
0082 print(''.join([dir.GetPath(), '/', name]))
0083
0084
0085 prefix = ""
0086 for i in self.dirs:
0087 if currPath.startswith(i):
0088 prefix = self.dirsAliases[i]
0089 break
0090
0091 th2.SetName(prefix + th2.GetName())
0092 self.listOfNumHistograms.append(th2)
0093 else:
0094 self.__TraverseDirTree(obj.ReadObj())
0095
0096 def __GetPartStr(self, isXlowerThanZero, isYlowerThanZero):
0097 if isXlowerThanZero and isYlowerThanZero:
0098 return "mO"
0099 if isXlowerThanZero and isYlowerThanZero == False:
0100 return "mI"
0101 if isXlowerThanZero == False and isYlowerThanZero:
0102 return "pO"
0103 if isXlowerThanZero == False and isYlowerThanZero == False:
0104 return "pI"
0105
0106 def __GetBarrelSector(self, layer, signedLadder, signedModule):
0107 theLadder = abs(signedLadder)
0108 theModule = abs(signedModule)
0109
0110 sector = 0
0111
0112 if layer == 1:
0113
0114 if theLadder == 1:
0115 if theModule >= 2:
0116 return 1
0117 else:
0118 return 2
0119 if theLadder == 2:
0120 if theModule >= 3:
0121 return 2
0122 else:
0123 return 3
0124 if theLadder == 3:
0125 if theModule >= 4:
0126 return 3
0127 else:
0128 return 4
0129 if theLadder == 4:
0130 if theModule >= 2:
0131 return 5
0132 else:
0133 return 6
0134 if theLadder == 5:
0135 if theModule >= 3:
0136 return 6
0137 else:
0138 return 7
0139 if theLadder == 6:
0140 if theModule >= 4:
0141 return 7
0142 else:
0143 return 8
0144
0145 elif layer == 2:
0146 i = theLadder // 5
0147 sector = i * 3
0148 shortLadder = theLadder - 5 * i
0149 for i in range(0, shortLadder, 2):
0150 sector = sector + 1
0151 return sector
0152 elif layer == 3:
0153 sector = 1
0154 for i in range(2, theLadder, 3):
0155 if (i + 1) % 3 == 0:
0156 sector = sector + 1
0157 return sector
0158 elif layer == 4:
0159 sector = (theLadder + 3) // 4
0160 return sector
0161
0162 def __BuildOnlineBarrelName(self, signedModule, signedLadder, layer):
0163 thePart = self.__GetPartStr(signedModule < 0, signedLadder < 0)
0164 theSector = str(self.__GetBarrelSector(layer, signedLadder, signedModule))
0165 return "BPix_B" + thePart + "_SEC" + theSector + "_LYR" + str(layer) + "_LDR" + str(abs(signedLadder)) + "F_MOD" + str(abs(signedModule))
0166
0167 def __BuildOnlineDiskName(self, signedDisk, signedBlade, panel, ring):
0168 thePart = self.__GetPartStr(signedDisk < 0, signedBlade < 0)
0169 return "FPix_B" + thePart + "_D" + str(abs(signedDisk)) + "_BLD" + str(abs(signedBlade)) + "_PNL" + str(panel) + "_RNG" + str(ring)
0170
0171 def __GroupHistograms(self):
0172 currentGroupName = ""
0173 groupOfHists = []
0174 self.groupedHistograms = []
0175
0176
0177 for obj in self.listOfNumHistograms:
0178 objName = obj.GetName()
0179 objNameSplit = objName.split("_")
0180 objNameCollected = ''.join(objNameSplit[0:-1])
0181 if objNameCollected != currentGroupName:
0182 if len(groupOfHists):
0183 self.groupedHistograms.append(groupOfHists)
0184 groupOfHists = []
0185
0186 currentGroupName = objNameCollected
0187 groupOfHists.append(obj)
0188 self.groupedHistograms.append(groupOfHists)
0189
0190 def __AddNamedBins(self, geoFile, tX, tY, sX, sY, applyModuleRotation = False):
0191
0192 for line in geoFile:
0193 lineSpl = line.strip().split("\"")
0194
0195 detId = lineSpl[0].split(" ")[0]
0196
0197 vertices = lineSpl[1]
0198 xy = vertices.split(" ")
0199 x, y = array('d'), array('d')
0200 verNum = 1
0201 for coord in xy:
0202 coordSpl = coord.split(",")
0203 if applyModuleRotation:
0204 x.append(-(float(coordSpl[0]) * sX + tX))
0205 y.append((float(coordSpl[1]) * sY + tY))
0206 else:
0207 x.append(float(coordSpl[0]) * sX + tX)
0208 y.append(float(coordSpl[1]) * sY + tY)
0209 verNum = verNum + 1
0210
0211 x.append(x[0])
0212 y.append(y[0])
0213
0214
0215
0216
0217 if applyModuleRotation:
0218 bin = TGraph(verNum, y, x)
0219 else:
0220 bin = TGraph(verNum, x, y)
0221
0222 bin.SetName(detId)
0223
0224 self.__BaseTrackerMap.AddBin(bin)
0225
0226 def __CreateTrackerBaseMap(self):
0227
0228 self.__BaseTrackerMap = TH2Poly("Summary", "", -10, 160, -70, 70)
0229
0230 self.__BaseTrackerMap.SetFloat(1)
0231 self.__BaseTrackerMap.GetXaxis().SetTitle("")
0232 self.__BaseTrackerMap.GetYaxis().SetTitle("")
0233 self.__BaseTrackerMap.SetOption("COLZ L")
0234 self.__BaseTrackerMap.SetStats(0)
0235
0236
0237 for i in range(maxPxBarrel):
0238 with open(self.geometryFilenames[i], "r") as geoFile:
0239 currBarrelTranslateX = 0
0240 currBarrelTranslateY = barrelLadderShift[i]
0241
0242 self.__AddNamedBins(geoFile, currBarrelTranslateX, currBarrelTranslateY, 1, 1, True)
0243
0244
0245
0246
0247 for i in range(-maxPxForward, 0):
0248 with open(self.geometryFilenames[maxPxBarrel + maxPxForward + i], "r") as geoFile:
0249 currForwardTranslateX = forwardDiskXShift[-i - 1]
0250 currForwardTranslateY = -forwardDiskYShift
0251
0252 self.__AddNamedBins(geoFile, currForwardTranslateX, currForwardTranslateY, 1, 1)
0253
0254
0255 for i in range(maxPxForward):
0256 with open(self.geometryFilenames[maxPxBarrel + maxPxForward + i], "r") as geoFile:
0257 currForwardTranslateX = forwardDiskXShift[i]
0258 currForwardTranslateY = forwardDiskYShift
0259
0260 self.__AddNamedBins(geoFile, currForwardTranslateX, currForwardTranslateY, 1, 1)
0261
0262
0263
0264 print("Base Tracker Map: constructed")
0265
0266
0267 def __init__(self, inputDQMName, outputDirName, minMaxFileName, limits, modDicName, runNumber, dirs, dirsAliases):
0268
0269 self.inputFileName = inputDQMName
0270 self.outputDirName = outputDirName
0271 self.minMaxFileName = minMaxFileName
0272
0273 self.detIDsFileName = modDicName
0274 self.limits = limits
0275
0276 self.runNumber = runNumber
0277 self.dirs = dirs
0278 self.dirsAliases = dirsAliases
0279
0280 self.inputFile = TFile(self.inputFileName)
0281 self.listOfNumHistograms = []
0282 self.availableNames = []
0283
0284 self.maxLadderToLayer = {6:1, 14:2, 22:3, 32:4}
0285 self.maxBladeToRing = {11:1, 17:2}
0286
0287 self.geometryFilenames = []
0288 for i in range(maxPxBarrel):
0289 self.geometryFilenames.append(getFileInPath("DQM/SiStripMonitorClient/data/Geometry/vertices_barrel_" + str(i + 1)))
0290
0291 for i in range(-maxPxForward, maxPxForward + 1):
0292 if i == 0:
0293 continue
0294 self.geometryFilenames.append(getFileInPath("DQM/SiStripMonitorClient/data/Geometry/vertices_forward_" + str(i)))
0295
0296
0297 self.internalData = {}
0298
0299 if self.inputFile.IsOpen():
0300 print("%s opened successfully!" % (self.inputFileName))
0301
0302 for dir in self.dirs:
0303 self.__TraverseDirTree(self.inputFile.Get(dir))
0304
0305
0306 self.detDict = {}
0307
0308 with open(self.detIDsFileName, "r") as detIDs:
0309 for entry in detIDs:
0310 items = entry.replace("\n", " ").split(" ")
0311 self.detDict.update({items[1] : int(items[0])})
0312
0313 self.internalData.update({int(items[0]) : {}})
0314
0315 self.rawToOnlineDict = dict((v,k) for k,v in self.detDict.items())
0316
0317 self.__GroupHistograms()
0318
0319 self.__CreateTrackerBaseMap()
0320
0321 else:
0322 print("Unable to open file %s" % (self.inputFileName))
0323
0324
0325
0326 self.limitsDic = {}
0327 for y in limits:
0328
0329 lineSpl = y.strip().split(" ")
0330
0331 if len(lineSpl) < 5:
0332 continue
0333
0334 currName = lineSpl[0]
0335 zMin = float(lineSpl[1])
0336 zMax = float(lineSpl[2])
0337 isLog = False if lineSpl[3] == "0" else True
0338 isAbs = False if lineSpl[4] == "0" else True
0339
0340 self.limitsDic.update({currName : {"zMin" : zMin, "zMax" : zMax, "isLog" : isLog, "isAbs" : isAbs}})
0341
0342
0343 def ReadHistograms(self):
0344 if self.inputFile.IsOpen():
0345 for group in self.groupedHistograms:
0346
0347 if len(group) == 0:
0348 return
0349 print(group[0].GetName())
0350 name = ''.join(group[0].GetName().split("_per_")[0])
0351 self.availableNames.append(name)
0352
0353 for obj in group:
0354 nbinsX = obj.GetNbinsX()
0355 nbinsY = obj.GetNbinsY()
0356
0357 if nbinsX == 9:
0358 maxX = nbinsX // 2
0359 maxY = nbinsY // 2
0360
0361 for x in range(-maxX, maxX + 1):
0362 if x == 0:
0363 continue
0364 for y in range(-maxY, maxY + 1, 1):
0365 if y == 0:
0366 continue
0367 onlineName = self.__BuildOnlineBarrelName(x, y, self.maxLadderToLayer[maxY])
0368 self.internalData[self.detDict[onlineName]].update({name : obj.GetBinContent(x + maxX + 1, y + maxY + 1)})
0369
0370 elif nbinsX == 7:
0371 maxX = nbinsX // 2
0372 maxY = nbinsY // 4
0373
0374 for x in range(-maxX, maxX + 1):
0375 if x == 0:
0376 continue
0377 for y in range(-maxY, maxY + 1):
0378 if int(y) == 0:
0379 continue
0380 for panel in range(1, 3):
0381 onlineName = self.__BuildOnlineDiskName(x, y, panel, self.maxBladeToRing[maxY])
0382 self.internalData[self.detDict[onlineName]].update({name : obj.GetBinContent(x + maxX + 1, (y + maxY) * 2 + (3-panel))})
0383 else:
0384 print("Unrecognized plot")
0385 else:
0386 print("Histograms saved to internal data structure")
0387
0388 def DumpData(self):
0389 for key in self.internalData:
0390 print("#"*20)
0391 print(key)
0392 module = self.internalData[key]
0393 for d in module:
0394 print((d, module[d]))
0395
0396 print(len(self.internalData))
0397
0398 for i in self.availableNames:
0399 print(i)
0400 print(len(self.availableNames))
0401
0402 def PrintTrackerMaps(self):
0403 monitoredValues = []
0404 gStyle.SetPalette(1)
0405 for key in self.internalData:
0406 monitoredValues = self.internalData[key].keys()
0407
0408 break
0409
0410 if os.path.exists(self.outputDirName) == False:
0411 os.system("mkdir " + self.outputDirName)
0412
0413 with open(self.outputDirName + self.minMaxFileName, "w") as minMaxFile:
0414
0415 for mv in monitoredValues:
0416 currentHist = deepcopy(self.__BaseTrackerMap)
0417
0418 histoTitle = "Run " + self.runNumber + ": Tracker Map for " + mv
0419
0420 applyLogScale = False
0421 applyAbsValue = False
0422 if mv in self.limitsDic:
0423 limitsElem = self.limitsDic[mv]
0424
0425 print(mv + " found in limits dictionary - applying custom limits...")
0426
0427 currentHist.SetMinimum(limitsElem["zMin"])
0428 currentHist.SetMaximum(limitsElem["zMax"])
0429 applyLogScale = limitsElem["isLog"]
0430 applyAbsValue = limitsElem["isAbs"]
0431
0432 listOfVals = []
0433 onlineName = ""
0434 for detId in self.internalData:
0435 val = (self.internalData[detId])[mv]
0436 onlineName = self.rawToOnlineDict[detId]
0437 listOfVals.append([val, detId, onlineName])
0438 if applyAbsValue:
0439 currentHist.Fill(str(detId), abs(val))
0440 else:
0441 currentHist.Fill(str(detId), val)
0442
0443 listOfVals = sorted(listOfVals, key = lambda item: item[0])
0444
0445 minMaxFile.write("\n" + mv + "\n\n")
0446
0447 minMaxFile.write("MIN:\n")
0448 for i in range(extremeBinsNum):
0449 minMaxFile.write("\t" + str(listOfVals[i][1]) + " " + str(listOfVals[i][2]) + " " + str(listOfVals[i][0]) + "\n")
0450
0451 minMaxFile.write("MAX:\n")
0452 for i in range(extremeBinsNum):
0453 minMaxFile.write("\t" + str(listOfVals[-i - 1][1]) + " " + str(listOfVals[-i - 1][2]) + " " + str(listOfVals[-i - 1][0]) + "\n")
0454
0455 c1 = TCanvas(mv, mv, plotWidth , plotHeight)
0456
0457 if applyLogScale:
0458 c1.SetLogz()
0459
0460
0461
0462 currentHist.Draw("AC COLZ L")
0463
0464 gPad.Update()
0465 palette = currentHist.FindObject("palette");
0466 palette.SetX1NDC(0.89);
0467 palette.SetX2NDC(0.91);
0468 palette.SetLabelSize(0.05);
0469 gPad.Update()
0470
0471
0472
0473
0474
0475
0476
0477
0478
0479
0480
0481 arrow = TArrow(0.05, 27.0, 0.05, -30.0, 0.02, "|>")
0482 arrow.SetLineWidth(4)
0483 arrow.Draw()
0484
0485 phiArrow = TArrow(0.0, 27.0, 30.0, 27.0, 0.02, "|>")
0486 phiArrow.SetLineWidth(4)
0487 phiArrow.Draw()
0488
0489 xArrow = TArrow(25.0, 44.5, 50.0, 44.5, 0.02, "|>")
0490 xArrow.SetLineWidth(4)
0491 xArrow.Draw()
0492
0493 yArrow = TArrow(25.0, 44.5, 25.0, 69.5, 0.02, "|>")
0494 yArrow.SetLineWidth(4)
0495 yArrow.Draw()
0496
0497
0498
0499
0500 txt = TLatex()
0501 txt.SetNDC()
0502 txt.SetTextFont(1)
0503 txt.SetTextColor(1)
0504 txt.SetTextAlign(22)
0505 txt.SetTextAngle(0)
0506
0507
0508 txt.SetTextSize(0.05)
0509 txt.DrawLatex(0.5, 0.95, histoTitle)
0510
0511 txt.SetTextSize(0.03)
0512
0513 txt.DrawLatex(0.5, 0.125, "-DISK")
0514 txt.DrawLatex(0.5, 0.075, "NUMBER ->")
0515 txt.DrawLatex(0.5, 0.875, "+DISK")
0516
0517 txt.DrawLatex(0.17, 0.35, "+z")
0518 txt.DrawLatexNDC(0.36, 0.685, "+phi")
0519 txt.DrawLatex(0.38, 0.73, "+x")
0520 txt.DrawLatex(0.26, 0.875, "+y")
0521
0522 txt.SetTextAngle(90)
0523 txt.DrawLatex(0.17, 0.5, "BARREL")
0524
0525
0526 c1.Print(self.outputDirName + mv + ".png")
0527
0528 def __del__(self):
0529 if self.inputFile :
0530 if self.inputFile.IsOpen():
0531 self.inputFile.Close()
0532
0533
0534 for i in range(1, len(sys.argv), 1):
0535 if i == 1:
0536 inputFileName = sys.argv[i]
0537 elif i == 2:
0538 plotWidth = int(sys.argv[i])
0539 elif i == 3:
0540 plotHeight = int(sys.argv[i])
0541
0542
0543
0544 elif i == 4:
0545 detIDsFileName = sys.argv[i]
0546
0547 deductedRunNumber = inputFileName.split("_R000")[1][0:6]
0548 print(deductedRunNumber)
0549
0550 baseRootDirs = ["DQMData/Run " + deductedRunNumber + "/PixelPhase1/Run summary/Phase1_MechanicalView"
0551 ,"DQMData/Run " + deductedRunNumber + "/PixelPhase1/Run summary/Tracks"
0552 ]
0553
0554 baseRootDirsAliases = {baseRootDirs[0]:""
0555 , baseRootDirs[1]:"T"
0556 }
0557
0558 readerObj = TH2PolyOfflineMaps(inputFileName, outputDirectoryName, minMaxFileName, limits, detIDsFileName, deductedRunNumber, baseRootDirs, baseRootDirsAliases)
0559
0560 readerObj.ReadHistograms()
0561
0562 readerObj.PrintTrackerMaps()