Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-11-26 02:34:10

0001 #!/usr/bin/env python3
0002 
0003 import ROOT
0004 from ROOT import TBufferFile, TH1F, TProfile, TProfile2D, TH2F, TFile, TH1D, TH2D
0005 import re
0006 import os
0007 
0008 
0009 def draw_line(lineList,x1,x2,y1,y2,width=1,style=1,color=1):
0010     from ROOT import TLine
0011     l=TLine(x1,y1,x2,y2)
0012     l.SetBit(ROOT.kCanDelete)
0013     l.SetLineWidth(width)
0014     l.SetLineStyle(style)
0015     l.SetLineColor(color)
0016     l.Draw()
0017     lineList.append(l)
0018 
0019 def draw_box(boxList,xl,xr,yl,yr,opacity=1,color=1,style=1001,lstyle=1,lw=3):
0020     from ROOT import TBox
0021     b=TBox(xl,yl,xr,yr)
0022     b.SetBit(ROOT.kCanDelete)
0023     b.SetFillStyle(style)
0024     b.SetFillColorAlpha(color, opacity)
0025     b.SetLineColor(color)
0026     b.SetLineWidth(lw)
0027     b.SetLineStyle(lstyle)
0028     b.Draw()
0029     boxList.append(b)
0030 
0031 def renderPluginBPIX(lineList,layer) :
0032     from ROOT import TCanvas,TLine
0033     nlad=[6,14,22,32]
0034     coordSign=[(-1,-1),(-1,1),(1,-1),(1,1)]
0035     for xsign,ysign in coordSign:
0036         xlow = xsign*0.5
0037         xhigh= xsign*(0.5+4)
0038         ylow = ysign*0.5
0039         yhigh= ysign*(0.5 + nlad[layer-1])
0040         # Outside Box
0041         draw_line(lineList,xlow,  xhigh,  ylow,  ylow) # bottom
0042         draw_line(lineList,xlow,  xhigh, yhigh, yhigh) # top
0043         draw_line(lineList,xlow,   xlow,  ylow, yhigh) # left
0044         draw_line(lineList,xhigh, xhigh,  ylow, yhigh) # right
0045         # Inner Horizontal lines
0046         for lad in range(nlad[layer-1]):
0047             lad+=1
0048             if lad != nlad[layer-1]:
0049                 y = ysign * (lad+0.5)
0050                 draw_line(lineList,xlow, xhigh,  y,  y)
0051             y = ysign * (lad);
0052             draw_line(lineList,xlow, xhigh,  y,  y, 1, 3);
0053         # Inner Vertical lines
0054         for mod in range(3) : 
0055             mod+=1
0056             x = xsign * (mod + 0.5);
0057             draw_line(lineList,x, x,  ylow,  yhigh);
0058 
0059         # Draw ROC0
0060         for mod in range(4):
0061             mod+=1
0062             for lad in range(nlad[layer-1]):
0063                 lad+=1
0064                 if ysign==1:
0065                     flipped = not(lad%2==0)
0066                 else :
0067                     flipped = not(lad%2==1)
0068                 if flipped : roc0_orientation = -1
0069                 else : roc0_orientation = 1
0070                 if xsign==-1 : roc0_orientation *= -1
0071                 if ysign==-1 : roc0_orientation *= -1
0072                 x1 = xsign * (mod+0.5)
0073                 x2 = xsign * (mod+0.5 - 1./8);
0074                 y1 = ysign * (lad)
0075                 y2 = ysign * (lad + roc0_orientation*1./2)
0076                 if layer == 1 and xsign == -1 :
0077                     x1 = xsign * (mod-0.5)
0078                     x2 = xsign * (mod-0.5 + 1./8)
0079                     y1 = ysign * (lad)
0080                     y2 = ysign * (lad - roc0_orientation*1./2)
0081 
0082                     draw_line(lineList,x1, x2, y1, y1, 1)
0083                     draw_line(lineList,x2, x2, y1, y2, 1)
0084                   
0085                 else:
0086                     draw_line(lineList,x1, x2, y1, y1, 1)
0087                     draw_line(lineList,x2, x2, y1, y2, 1)
0088 
0089 def maskBPixROC(boxList,xsign,ysign,layer,lad,mod,roc):
0090     if roc<8 : 
0091         rocShiftX=roc*1./8
0092         rocShiftY=0
0093     else : 
0094         rocShiftX=(15-roc)*1./8
0095         rocShiftY=1./2
0096     if ysign==1:
0097         flipped = not(lad%2==0)
0098     else :
0099         flipped = not(lad%2==1)
0100     if flipped : roc0_orientation = -1
0101     else : roc0_orientation = 1
0102     if xsign==-1 : roc0_orientation *= -1
0103     if ysign==-1 : roc0_orientation *= -1
0104     x1 = xsign * (mod+0.5-rocShiftX)
0105     x2 = xsign * (mod+0.5 - 1./8-rocShiftX);
0106     y1 = ysign * (lad-roc0_orientation*rocShiftY)
0107     y2 = ysign * (lad + roc0_orientation*1./2-roc0_orientation*rocShiftY)
0108     if layer == 1 and xsign == -1 :
0109         x1 = xsign * (mod-0.5)-rocShiftX
0110         x2 = xsign * (mod-0.5 + 1./8)-rocShiftX
0111         y1 = ysign * (lad +rocShiftY)
0112         y2 = ysign * (lad - roc0_orientation*1./2+rocShiftY)
0113     draw_box(boxList,min(x1,x2),max(x1,x2),min(y1, y2),max(y1,y2),0.75)
0114 
0115 
0116                   
0117 def renderPluginFPIX(lineList,ring) :
0118     from ROOT import TCanvas,TLine
0119     coordSign=[(-1,-1),(-1,1),(1,-1),(1,1)]
0120     for dsk in range(3) :
0121         dsk+=1
0122         for xsign,ysign in coordSign:            
0123             for bld in range(5+ring*6):
0124                 bld+=1
0125                 # Panel 2 has dashed mid-plane
0126                 x1      = xsign * (0.5 + dsk - 1)
0127                 x2      = xsign * (0.5 + dsk)
0128                 sign = ysign
0129                 y1      = ysign * (bld + sign*0.5)
0130                 y2      = ysign * (bld)
0131                 yp2_mid = ysign * (bld - sign*0.25)
0132                 y3      = ysign * (bld - sign*0.5)
0133                 draw_line(lineList,x1, x2, y1, y1)
0134                 draw_line(lineList,x1, x2, y2, y2)
0135                 draw_line(lineList,x1, x2, yp2_mid, yp2_mid,1,2)
0136                 draw_line(lineList,x1, x2, y3, y3)
0137                 # Vertical lines
0138                 x = xsign * (0.5 + dsk - 1)
0139                 draw_line(lineList,x,  x,  y1,  y2)
0140                 draw_line(lineList,x,  x,  y2,  y3)
0141                 if ring==2 :
0142                     x = xsign * (0.5 + dsk)
0143                     draw_line(lineList,x,  x,  y1,  y2)
0144                     draw_line(lineList,x,  x,  y2,  y3)
0145                 #Make a BOX around ROC 0
0146                 x1 = xsign * (0.5 + dsk - 1/8.)
0147                 x2 = xsign * (0.5 + dsk)
0148                 y1_p1 = ysign * (bld + sign*0.25)
0149                 y2_p1 = ysign * (bld + sign*0.25 + xsign*ysign*0.25)
0150                 draw_line(lineList,x1, x2, y1_p1, y1_p1, 1)
0151                 draw_line(lineList,x1, x1, y1_p1, y2_p1, 1)
0152                 y1_p2 = ysign * (bld - sign*0.25)
0153                 y2_p2 = ysign * (bld - sign*0.25 - xsign*ysign*0.25)
0154                 draw_line(lineList,x1, x2, y1_p2, y1_p2)
0155                 draw_line(lineList,x1, x1, y1_p2, y2_p2)
0156 
0157 def maskFPixROC(boxList,xsign,ysign,dsk,bld,pnl,roc) :
0158     from ROOT import TCanvas,TLine       
0159     if roc<8 : 
0160         rocShiftX=roc*1./8
0161         rocShiftY=0
0162     else : 
0163         rocShiftX=(15-roc)*1./8
0164         rocShiftY=1./4
0165     sign=ysign
0166     x1 = xsign * (0.5 + dsk - 1/8.-rocShiftX)
0167     x2 = xsign * (0.5 + dsk-rocShiftX)
0168     if pnl==1:
0169         y1 = ysign * (bld + sign*0.25)-xsign*rocShiftY
0170         y2 = ysign * (bld + sign*0.25 + xsign*ysign*0.25)-xsign*rocShiftY
0171     else:
0172         y1 = ysign * (bld - sign*0.25)+xsign*rocShiftY
0173         y2 = ysign * (bld - sign*0.25 - xsign*ysign*0.25)+xsign*rocShiftY
0174     draw_box(boxList,min(x1,x2),max(x1,x2),min(y1,y2),max(y1,y2),0.75)
0175 
0176 
0177 def dqm_get_dataset(server, match, run, type="offline_data"):
0178     datareq = urllib2.Request(('%s/data/json/samples?match=%s') % (server, match))
0179     datareq.add_header('User-agent', ident)
0180     # Get data                                                                                                                              
0181     data = eval(re.sub(r"\bnan\b", "0", urllib2.build_opener(X509CertOpen()).open(datareq).read()),
0182                        { "__builtins__": None }, {})
0183     ret = ""
0184     for l in data['samples']:
0185         if l['type'] == type:
0186             for x in l['items']:
0187                 if int(x['run']) == int(run):
0188                     ret=x['dataset']
0189                     break
0190     print(ret)
0191     return ret
0192 
0193 
0194 
0195 
0196 def main():
0197     import sys
0198     import os
0199     import ROOT
0200 
0201     if len(sys.argv) != 2:
0202         print("input file needed!")
0203         return
0204     else:
0205         filename=sys.argv[1]
0206     print(filename+" -- "+filename[19:25])
0207     runNum=filename[19:25]
0208 
0209     dir="DQMData/Run " + runNum + "/PixelPhase1/Run summary/Pahse1_MechanicalView/"
0210     dirFED="DQMData/Run " + runNum + "/PixelPhase1/Run summary/FED/"
0211 
0212 
0213     dirBPix=dir + "PXBarrel/"
0214     dirFPix=dir + "PXForward/"
0215 
0216     hoccB="digi_occupancy_per_SignedModuleCoord_per_SignedLadderCoord_PXLayer_"
0217     hoccF="digi_occupancy_per_SignedDiskCoord_per_SignedBladePanelCoord_PXRing_"
0218     hdeadB="Dead Channels per ROC_per_SignedModuleCoord_per_SignedLadderCoord_PXLayer_"
0219     hdeadF="Dead Channels per ROC_per_SignedDiskCoord_per_SignedBladePanelCoord_PXRing_"
0220         
0221     ROOT.gROOT.SetBatch(1)    
0222     ROOT.gStyle.SetOptStat(0)
0223     ROOT.gStyle.SetPalette(1) #104 kTemperatureMap // 55 kRainBow
0224     ROOT.gStyle.SetNumberContours(128)
0225     rootf=ROOT.TFile(filename)
0226     
0227 
0228 
0229     c=ROOT.TCanvas("c","c",1250,1000)
0230     #BPIX
0231     print("----> Build maps for BPix")
0232     histOccList=[]
0233     histDeadList=[]
0234     for lyr in range(1,5):
0235         histOccList.append(rootf.FindObjectAny(hoccB+str(lyr)))
0236         histDeadList.append(rootf.FindObjectAny(hdeadB+str(lyr)))
0237     for hist1, hist2 in zip(histOccList, histDeadList):
0238         if hist1 != None or hist2 !=None:
0239             hist1.Draw("colz")
0240             match=re.search('(?<=PXLayer_)[0-9]',hist1.GetName())
0241             if match != None and "per_SignedModuleCoord_per_SignedLadderCoord" in hist1.GetName():
0242                 lyr=int(match.group(0))
0243                 hist1.SetTitle("Digi Occupancy Layer {0}".format(lyr))
0244                 boxList=[]
0245                 lineList=[]
0246                 renderPluginBPIX(lineList,lyr)
0247                 lineWd=3
0248                 if lyr==4 :
0249                     lineWd=2
0250                 if lyr==1: 
0251                     tbmRoc=4
0252                 else:
0253                     tbmRoc=8
0254                 binTBM=[]
0255                 for biny in range(1,hist2.GetNbinsY()+1):
0256                     if len(binTBM)!=0:
0257                         x1=hist2.GetXaxis().GetBinLowEdge(binTBM[0])
0258                         x2=hist2.GetXaxis().GetBinUpEdge(binTBM[len(binTBM)-1])
0259                         y1=hist2.GetYaxis().GetBinLowEdge(biny-1)
0260                         y2=hist2.GetYaxis().GetBinUpEdge(biny-1)
0261                         if hist2.GetBinContent(binTBM[0],biny-1) >= 0.97*hist2.GetMaximum():
0262                             draw_box(boxList,x1,x2,y1,y2,0.2,1,0,1,lineWd)
0263                         else:
0264                             draw_box(boxList,x1,x2,y1,y2,0.2,633,0,1,lineWd)
0265                     binTBM=[]
0266                     for binx in range(1,hist2.GetNbinsX()+1):
0267                         if hist2.GetBinContent(binx,biny)!=0:
0268                             if len(binTBM)==0:
0269                                 binTBM.append(binx)
0270                             else:
0271                                 if len(binTBM)==tbmRoc: 
0272                                     x1=hist2.GetXaxis().GetBinLowEdge(binTBM[0])
0273                                     x2=hist2.GetXaxis().GetBinUpEdge(binTBM[len(binTBM)-1])
0274                                     y1=hist2.GetYaxis().GetBinLowEdge(biny)
0275                                     y2=hist2.GetYaxis().GetBinUpEdge(biny)
0276                                     if hist2.GetBinContent(binTBM[0],biny) >= 0.97*hist2.GetMaximum():
0277                                         draw_box(boxList,x1,x2,y1,y2,0.2,1,0,1,lineWd)
0278                                     else:
0279                                         draw_box(boxList,x1,x2,y1,y2,0.2,633,0,1,lineWd)
0280                                     binTBM=[]
0281                                     binTBM.append(binx)
0282                                 else:
0283                                     binTBM.append(binx)
0284                         else:
0285                             if len(binTBM)!=0:
0286                                 x1=hist2.GetXaxis().GetBinLowEdge(binTBM[0])
0287                                 x2=hist2.GetXaxis().GetBinUpEdge(binTBM[len(binTBM)-1])
0288                                 y1=hist2.GetYaxis().GetBinLowEdge(biny)
0289                                 y2=hist2.GetYaxis().GetBinUpEdge(biny)
0290                                 if hist2.GetBinContent(binTBM[0],biny) >= 0.97*hist2.GetMaximum():
0291                                     draw_box(boxList,x1,x2,y1,y2,0.2,1,0,1,lineWd)
0292                                 else:
0293                                     draw_box(boxList,x1,x2,y1,y2,0.2,633,0,1,lineWd)
0294                 c.SaveAs('MergedOccupancyDeadROC_BPix_Layer{0}_TBM.pdf'.format(lyr))
0295                 os.system('gs -dBATCH -dNOPAUSE -sDEVICE=png16m -dUseCropBox -sOutputFile=MergedOccupancyDeadROC_BPix_Layer{0}_TBM.png -r144 -q MergedOccupancyDeadROC_BPix_Layer{0}_TBM.pdf'.format(lyr))
0296                 os.system('rm -f MergedOccupancyDeadROC_BPix_Layer{0}_TBM.pdf'.format(lyr))
0297         else :
0298             print("Some Error in get the histograms for FPIX")
0299     #FPIX
0300     print("----> Build maps for FPix")
0301     for rng in range(1,3):
0302         histOccList.append(rootf.FindObjectAny(hoccF+str(rng)))
0303         histDeadList.append(rootf.FindObjectAny(hdeadF+str(rng)))
0304     for hist1, hist2 in zip(histOccList, histDeadList):
0305         if hist1 != None or hist2 !=None:
0306             hist1.Draw("colz")
0307             match=re.search('(?<=PXRing_)[0-9]',hist1.GetName())
0308             if match != None and "per_SignedDiskCoord_per_SignedBladePanelCoord" in hist1.GetName():
0309                 ring=int(match.group(0))
0310                 hist1.SetTitle("Digi Occupancy Ring {0}".format(ring))
0311                 boxList=[]
0312                 lineList=[]
0313                 renderPluginFPIX(lineList,ring)
0314                 lineWd=3
0315                 if ring==2 :
0316                     lineWd=2
0317                 tbmRoc=8
0318                 binTBM=[]
0319                 for biny in range(1,hist2.GetNbinsY()+1):
0320                     if len(binTBM)!=0:
0321                         x1=hist2.GetXaxis().GetBinLowEdge(binTBM[0])
0322                         x2=hist2.GetXaxis().GetBinUpEdge(binTBM[len(binTBM)-1])
0323                         y1=hist2.GetYaxis().GetBinLowEdge(biny-1)
0324                         y2=hist2.GetYaxis().GetBinUpEdge(biny-1)
0325                         if hist2.GetBinContent(binTBM[0],biny-1) >= 0.97*hist2.GetMaximum():
0326                             draw_box(boxList,x1,x2,y1,y2,0.2,1,0,1,lineWd)
0327                         else:
0328                             draw_box(boxList,x1,x2,y1,y2,0.2,633,0,1,lineWd)
0329                     binTBM=[]
0330                     for binx in range(1,hist2.GetNbinsX()+1):
0331                         if hist2.GetBinContent(binx,biny)!=0:
0332                             if len(binTBM)==0:
0333                                 binTBM.append(binx)
0334                             else:
0335                                 if len(binTBM)==tbmRoc: 
0336                                     x1=hist2.GetXaxis().GetBinLowEdge(binTBM[0])
0337                                     x2=hist2.GetXaxis().GetBinUpEdge(binTBM[len(binTBM)-1])
0338                                     y1=hist2.GetYaxis().GetBinLowEdge(biny)
0339                                     y2=hist2.GetYaxis().GetBinUpEdge(biny)
0340                                     if hist2.GetBinContent(binTBM[0],biny) >= 0.97*hist2.GetMaximum():
0341                                         draw_box(boxList,x1,x2,y1,y2,0.2,1,0,1,lineWd)
0342                                     else:
0343                                         draw_box(boxList,x1,x2,y1,y2,0.2,633,0,1,lineWd)
0344                                     binTBM=[]
0345                                     binTBM.append(binx)
0346                                 else:
0347                                     binTBM.append(binx)
0348                         else:
0349                             if len(binTBM)!=0:
0350                                 x1=hist2.GetXaxis().GetBinLowEdge(binTBM[0])
0351                                 x2=hist2.GetXaxis().GetBinUpEdge(binTBM[len(binTBM)-1])
0352                                 y1=hist2.GetYaxis().GetBinLowEdge(biny)
0353                                 y2=hist2.GetYaxis().GetBinUpEdge(biny)
0354                                 if hist2.GetBinContent(binTBM[0],biny) >= 0.97*hist2.GetMaximum():
0355                                     draw_box(boxList,x1,x2,y1,y2,0.2,1,0,1,lineWd)
0356                                 else:
0357                                     draw_box(boxList,x1,x2,y1,y2,0.2,633,0,1,lineWd)
0358                 c.SaveAs('MergedOccupancyDeadROC_FPix_Ring{0}_TBM.pdf'.format(ring))
0359                 os.system('gs -dBATCH -dNOPAUSE -sDEVICE=png16m -dUseCropBox -sOutputFile=MergedOccupancyDeadROC_FPix_Ring{0}_TBM.png -r144 -q MergedOccupancyDeadROC_FPix_Ring{0}_TBM.pdf'.format(ring))
0360                 os.system('rm -f MergedOccupancyDeadROC_FPix_Ring{0}_TBM.pdf'.format(ring))
0361 
0362         else :
0363             print("Some Error in get the histograms for FPIX")
0364 
0365 
0366 
0367 
0368 if __name__ == '__main__':
0369     main()