Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-11-27 03:18:06

0001 #!/usr/bin/env python3
0002 #____________________________________________________________
0003 #
0004 #
0005 # A very simple way to make plots with ROOT via an XML file
0006 #
0007 # Francisco Yumiceva
0008 # yumiceva@fnal.gov
0009 #
0010 # Fermilab, 2010
0011 #
0012 #____________________________________________________________
0013 
0014 """
0015    plotBeamSpotDB
0016 
0017    A very simple script to plot the beam spot data stored in condDB
0018 
0019    usage: %prog -t <tag name>
0020    -a, --auth    = AUTH: DB authorization path. online(/nfshome0/popcondev/conddb).
0021    -b, --batch : Run ROOT in batch mode.
0022    -c, --create  = CREATE: name for beam spot data file.
0023    -d, --data    = DATA: input beam spot data file.
0024    -D, --destDB  = DESTDB: destination DB string. online(oracle://cms_orcon_prod/CMS_COND_31X_BEAMSPOT).
0025    -i, --initial = INITIAL: First IOV. Options: run number, or run:lumi, eg. \"133200:21\"
0026    -f, --final   = FINAL: Last IOV. Options: run number, or run:lumi
0027    -g, --graph : create a TGraphError instead of a TH1 object
0028    -n, --noplot : Only extract beam spot data, plots are not created..
0029    -o, --output  = OUTPUT: filename of ROOT file with plots.
0030    -p, --payload = PAYLOAD: filename of output text file. Combine and splits lumi IOVs.
0031    -P, --Print : create PNG plots from canvas.
0032    -s, --suffix = SUFFIX: Suffix will be added to plots filename.
0033    -t, --tag     = TAG: Database tag name.
0034    -T, --Time : create plots with time axis.
0035    -I, --IOVbase = IOVBASE: options: runbase(default), lumibase, timebase
0036    -w, --wait : Pause script after plotting a new histograms.
0037    -W, --weighted : Create a weighted result for a range of lumi IOVs, skip lumi IOV combination and splitting.
0038    -x, --xcrossing = XCROSSING : Bunch crossing number.
0039    
0040    Francisco Yumiceva (yumiceva@fnal.gov)
0041    Fermilab 2010
0042    
0043 """
0044 
0045 
0046 from builtins import range
0047 import os, string, re, sys, math
0048 import subprocess, time
0049 from BeamSpotObj import BeamSpot
0050 from IOVObj import IOV
0051 from CommonMethods import *
0052 
0053 try:
0054     import ROOT
0055 except:
0056     print("\nCannot load PYROOT, make sure you have setup ROOT in the path")
0057     print("and pyroot library is also defined in the variable PYTHONPATH, try:\n")
0058     if (os.getenv("PYTHONPATH")):
0059         print(" setenv PYTHONPATH ${PYTHONPATH}:$ROOTSYS/lib\n")
0060     else:
0061         print(" setenv PYTHONPATH $ROOTSYS/lib\n")
0062         sys.exit()
0063 
0064 from ROOT import TFile, TGraphErrors, TGaxis, TDatime
0065 from ROOT import TCanvas, TH1F
0066 
0067 # ROOT STYLE
0068 #############################
0069 def SetStyle():
0070 
0071     # canvas
0072     ROOT.gStyle.SetCanvasBorderMode(0)
0073     ROOT.gStyle.SetCanvasColor(0)
0074     ROOT.gStyle.SetCanvasDefH(600)
0075     ROOT.gStyle.SetCanvasDefW(600)
0076     ROOT.gStyle.SetCanvasDefX(0)
0077     ROOT.gStyle.SetCanvasDefY(0)
0078     # pad
0079     ROOT.gStyle.SetPadBorderMode(0)
0080     ROOT.gStyle.SetPadColor(0)
0081     ROOT.gStyle.SetPadGridX(False)
0082     ROOT.gStyle.SetPadGridY(False)
0083     ROOT.gStyle.SetGridColor(0)
0084     ROOT.gStyle.SetGridStyle(3)
0085     ROOT.gStyle.SetGridWidth(1)
0086 
0087     ROOT.gStyle.SetFrameBorderMode(0)
0088     ROOT.gStyle.SetFrameFillColor(0)
0089     ROOT.gStyle.SetTitleColor(1)
0090     ROOT.gStyle.SetStatColor(0)
0091 
0092     # set the paper & margin sizes
0093     ROOT.gStyle.SetPaperSize(20,26)
0094     ROOT.gStyle.SetPadTopMargin(0.04)
0095     ROOT.gStyle.SetPadRightMargin(0.04)
0096     ROOT.gStyle.SetPadBottomMargin(0.14)
0097     ROOT.gStyle.SetPadLeftMargin(0.11)
0098     ROOT.gStyle.SetPadTickX(1)
0099     ROOT.gStyle.SetPadTickY(1)
0100 
0101     ROOT.gStyle.SetTextFont(42) #132
0102     ROOT.gStyle.SetTextSize(0.09)
0103     ROOT.gStyle.SetLabelFont(42,"xyz")
0104     ROOT.gStyle.SetTitleFont(42,"xyz")
0105     ROOT.gStyle.SetLabelSize(0.035,"xyz")
0106     ROOT.gStyle.SetTitleSize(0.045,"xyz")
0107     ROOT.gStyle.SetTitleOffset(1.1,"y")
0108 
0109     # use bold lines and markers
0110     ROOT.gStyle.SetMarkerStyle(8)
0111     ROOT.gStyle.SetHistLineWidth(2)
0112     ROOT.gStyle.SetLineWidth(1)
0113     #ROOT.gStyle.SetLineStyleString(2,"[12 12]") // postscript dashes
0114 
0115     ROOT.gStyle.SetMarkerSize(0.6)
0116 
0117     # do not display any of the standard histogram decorations
0118     ROOT.gStyle.SetOptTitle(0)
0119     ROOT.gStyle.SetOptStat(0) #("m")
0120     ROOT.gStyle.SetOptFit(0)
0121 
0122     #ROOT.gStyle.SetPalette(1,0)
0123     ROOT.gStyle.cd()
0124     ROOT.gROOT.ForceStyle()
0125 #########################################
0126 
0127 
0128 if __name__ == '__main__':
0129 
0130     # style
0131     SetStyle()
0132 
0133     printCanvas = False
0134     printFormat = "png"
0135     printBanner = False
0136     Banner = "CMS Preliminary"
0137 
0138     # COMMAND LINE OPTIONS
0139     #################################
0140     option,args = parse(__doc__)
0141     if not args and not option: exit()
0142 
0143     tag = ''
0144     if not option.tag and not option.data: 
0145         print(" need to provide DB tag name or beam spot data file")
0146         exit()
0147     else:
0148         tag = option.tag
0149 
0150     if option.batch:
0151         ROOT.gROOT.SetBatch()
0152 
0153     datafilename = "tmp_beamspot.dat"
0154     if option.create:
0155         datafilename = option.create
0156 
0157     getDBdata = True
0158     if option.data:
0159         getDBdata = False
0160 
0161     IOVbase = 'runbase'
0162     if option.IOVbase:
0163         if option.IOVbase != "runbase" and option.IOVbase != "lumibase" and option.IOVbase != "timebase":
0164             print("\n\n unknown iov base option: "+ option.IOVbase +" \n\n\n")
0165             exit()
0166         IOVbase = option.IOVbase
0167 
0168     firstRun = "0"
0169     lastRun  = "4999999999"
0170     if IOVbase == "lumibase":
0171         firstRun = "0:0"
0172         lastRun = "4999999999:4999999999"
0173 
0174     if option.initial:
0175         firstRun = option.initial
0176     if option.final:
0177         lastRun = option.final
0178 
0179     # GET IOVs
0180     ################################
0181 
0182     if getDBdata:
0183 
0184         print(" read DB to get list of IOVs for the given tag")
0185         mydestdb = 'frontier://PromptProd/CMS_COND_31X_BEAMSPOT'
0186         if option.destDB:
0187             mydestdb = option.destDB
0188         acommand = 'cmscond_list_iov -c '+mydestdb+' -P /afs/cern.ch/cms/DB/conddb -t '+ tag
0189         tmpstatus = subprocess.getstatusoutput( acommand )
0190         tmplistiov = tmpstatus[1].split('\n')
0191         #print tmplistiov
0192 
0193         iovlist = []
0194         passline = False
0195         iline = jline = 0
0196         totlines = len(tmplistiov)
0197         for line in tmplistiov:
0198 
0199             if line.find('since') != -1:
0200                 passline = True
0201                 jline = iline
0202             if passline and iline > jline and iline < totlines-1:
0203                 linedata = line.split()
0204                 #print linedata
0205                 aIOV = IOV()
0206                 aIOV.since = int(linedata[0])
0207                 aIOV.till = int(linedata[1])
0208                 iovlist.append( aIOV )
0209             iline += 1
0210 
0211         print(" total number of IOVs = " + str(len(iovlist)))
0212 
0213 
0214         #  GET DATA
0215         ################################
0216 
0217         otherArgs = ''
0218         if option.destDB:
0219             otherArgs = " -d " + option.destDB
0220             if option.auth:
0221                 otherArgs = otherArgs + " -a "+ option.auth
0222 
0223         print(" get beam spot data from DB for IOVs. This can take a few minutes ...")
0224 
0225         tmpfile = open(datafilename,'w')
0226 
0227         for iIOV in iovlist:
0228             passiov = False
0229             tmprunfirst = firstRun
0230             tmprunlast = lastRun
0231             tmplumifirst = 1
0232             tmplumilast = 9999999
0233             if IOVbase=="lumibase":
0234                 #tmprunfirst = int(firstRun.split(":")[0])
0235                 #tmprunlast  = int(lastRun.split(":")[0])
0236                 #tmplumifirst = int(firstRun.split(":")[1])
0237                 #tmplumilast  = int(lastRun.split(":")[1])
0238                 tmprunfirst = pack( int(firstRun.split(":")[0]) , int(firstRun.split(":")[1]) )
0239                 tmprunlast  = pack( int(lastRun.split(":")[0]) , int(lastRun.split(":")[1]) )
0240             #print "since = " + str(iIOV.since) + " till = "+ str(iIOV.till)
0241             if iIOV.since >= int(tmprunfirst) and int(tmprunlast) < 0 and iIOV.since <= int(tmprunfirst):
0242                 print(" IOV: " + str(iIOV.since))
0243                 passiov = True
0244             if iIOV.since >= int(tmprunfirst) and int(tmprunlast) > 0 and iIOV.till <= int(tmprunlast):
0245                 print(" a IOV: " + str(iIOV.since) + " to " + str(iIOV.till))
0246                 passiov = True
0247             #if iIOV.since >= int(tmprunlast) and iIOV.till >= 4294967295:
0248             #    print " b IOV: " + str(iIOV.since) + " to " + str(iIOV.till)
0249             #    passiov = True                
0250             if passiov:
0251                 acommand = 'getBeamSpotDB.py -t '+ tag + " -r " + str(iIOV.since) +otherArgs
0252                 if IOVbase=="lumibase":
0253                     tmprun = unpack(iIOV.since)[0]
0254                     tmplumi = unpack(iIOV.since)[1]
0255                     acommand = 'getBeamSpotDB.py -t '+ tag + " -r " + str(tmprun) +" -l "+str(tmplumi) +otherArgs
0256                     print(acommand)
0257                 status = subprocess.getstatusoutput( acommand )
0258                 tmpfile.write(status[1])
0259 
0260         print(" beam spot data collected and stored in file " + datafilename)
0261 
0262         tmpfile.close()
0263 
0264 
0265     # PROCESS DATA
0266     ###################################
0267 
0268     # check if input data exists if given
0269     if option.data:
0270         if os.path.isdir(option.data):
0271             tmp = subprocess.getstatusoutput("ls "+option.data)
0272             files = tmp[1].split()
0273             datafilename = "combined_all.txt"
0274             output = open(datafilename,"w")
0275 
0276             for f in files:
0277                 input = open(option.data +"/"+f)
0278                 output.writelines(input.readlines())
0279             output.close()
0280             print(" data files have been collected in "+datafilename)
0281 
0282         elif os.path.exists(option.data):
0283             datafilename = option.data
0284         else:
0285             print(" input beam spot data DOES NOT exist, file " + option.data)
0286             exit()
0287 
0288     listbeam = []
0289 
0290     if option.xcrossing:
0291         listmap = readBeamSpotFile(datafilename,listbeam,IOVbase,firstRun,lastRun)
0292         # bx
0293         print("List of bunch crossings in the file:")
0294         print(listmap.keys())
0295         listbeam = listmap[option.Xrossing]
0296     else:
0297         readBeamSpotFile(datafilename,listbeam,IOVbase,firstRun,lastRun)
0298 
0299     sortAndCleanBeamList(listbeam,IOVbase)
0300 
0301     if IOVbase == "lumibase" and option.payload:
0302         weighted = True;
0303         if not option.weighted:
0304             weighted = False
0305         createWeightedPayloads(option.payload,listbeam,weighted)
0306     if option.noplot:
0307         print(" no plots requested, exit now.")
0308         sys.exit()
0309     # MAKE PLOTS
0310     ###################################    
0311     TGaxis.SetMaxDigits(8)
0312 
0313     graphlist = []
0314     graphnamelist = ['X','Y','Z','SigmaZ','dxdz','dydz','beamWidthX', 'beamWidthY']
0315     graphtitlelist = ['beam spot X','beam spot Y','beam spot Z','beam spot #sigma_Z','beam spot dX/dZ','beam spot dY/dZ','beam width X','beam width Y']
0316     graphXaxis = 'Run number'
0317     if IOVbase == 'runbase':
0318         graphXaxis = "Run number"
0319     if IOVbase == 'lumibase':
0320         graphXaxis = 'Lumi section'
0321     if IOVbase == 'timebase' or option.Time:
0322         graphXaxis = "Time"
0323         #dh = ROOT.TDatime(2010,06,01,00,00,00)
0324         ROOT.gStyle.SetTimeOffset(0) #dh.Convert())
0325 
0326     graphYaxis = ['beam spot X [cm]','beam spot Y [cm]','beam spot Z [cm]', 'beam spot #sigma_{Z} [cm]', 'beam spot dX/dZ', 'beam spot dY/dZ','beam width X [cm]', 'beam width Y [cm]']
0327 
0328     cvlist = []
0329 
0330     for ig in range(0,8):
0331         cvlist.append( TCanvas(graphnamelist[ig],graphtitlelist[ig], 1200, 600) )
0332         if option.graph:
0333             graphlist.append( TGraphErrors( len(listbeam) ) )
0334         else:
0335             graphlist.append( TH1F("name","title",len(listbeam),0,len(listbeam)) )
0336 
0337         graphlist[ig].SetName(graphnamelist[ig])
0338         graphlist[ig].SetTitle(graphtitlelist[ig])
0339         ipoint = 0
0340         for ii in range(0,len(listbeam)):
0341 
0342             ibeam = listbeam[ii]
0343             datax = dataxerr = 0.
0344             datay = datayerr = 0.
0345             if graphnamelist[ig] == 'X':
0346                 datay = ibeam.X
0347                 datayerr = ibeam.Xerr
0348             if graphnamelist[ig] == 'Y':
0349                 datay = ibeam.Y
0350                 datayerr = ibeam.Yerr
0351             if graphnamelist[ig] == 'Z':
0352                 datay = ibeam.Z
0353                 datayerr = ibeam.Zerr
0354             if graphnamelist[ig] == 'SigmaZ':
0355                 datay = ibeam.sigmaZ
0356                 datayerr = ibeam.sigmaZerr
0357             if graphnamelist[ig] == 'dxdz':
0358                 datay = ibeam.dxdz
0359                 datayerr = ibeam.dxdzerr
0360             if graphnamelist[ig] == 'dydz':
0361                 datay = ibeam.dydz
0362                 datayerr = ibeam.dydzerr
0363             if graphnamelist[ig] == 'beamWidthX':
0364                 datay = ibeam.beamWidthX
0365                 datayerr = ibeam.beamWidthXerr
0366             if graphnamelist[ig] == 'beamWidthY':
0367                 datay = ibeam.beamWidthY
0368                 datayerr = ibeam.beamWidthYerr
0369 
0370             datax = ibeam.IOVfirst
0371             if IOVbase=="lumibase":
0372                 datax = str(ibeam.Run) + ":" + str(ibeam.IOVfirst)
0373                 if ibeam.IOVfirst != ibeam.IOVlast:
0374                     datax = str(ibeam.Run) + ":" + str(ibeam.IOVfirst)+"-"+str(ibeam.IOVlast)
0375             #print datax
0376             if option.graph:
0377                 if IOVbase=="lumibase":
0378                     #first = int( pack( int(ibeam.Run) , int(ibeam.IOVfirst) ) )
0379                     #last = int( pack( int(ibeam.Run) , int(ibeam.IOVlast) ) )
0380                     first = ibeam.IOVfirst
0381                     last = ibeam.IOVlast
0382                     if option.Time:
0383                         atime = ibeam.IOVBeginTime
0384                         first = time.mktime( time.strptime(atime.split()[0] +  " " + atime.split()[1] + " " + atime.split()[2],"%Y.%m.%d %H:%M:%S %Z") )
0385                         atime = ibeam.IOVEndTime
0386                         last = time.mktime( time.strptime(atime.split()[0] +  " " + atime.split()[1] + " " + atime.split()[2],"%Y.%m.%d %H:%M:%S %Z") )
0387                         da_first = TDatime(time.strftime('%Y-%m-%d %H:%M:%S',time.localtime(first - time.timezone)))
0388                         da_last = TDatime(time.strftime('%Y-%m-%d %H:%M:%S',time.localtime(last - time.timezone)))
0389                         if ipoint == 0:
0390                             ## print local time
0391                             da_first.Print()
0392                             ## print gmt time
0393                             print("GMT = " + str(time.strftime('%Y-%m-%d %H:%M:%S',time.gmtime(first - time.timezone))))
0394                             reftime = first
0395                             ptm = time.localtime(reftime)
0396                             da = TDatime(time.strftime('%Y-%m-%d %H:%M:%S',ptm))
0397                             if time.daylight and ptm.tm_isdst:
0398                                 offset_daylight = time.timezone - time.altzone
0399                             ROOT.gStyle.SetTimeOffset(da.Convert(1) - 3600)
0400 
0401                     datax = (float(last) - float(first))/2 + float(first) - da.Convert() + 3600
0402                     ## Comment out this block if running on Mac ##
0403                     if time.daylight and ptm.tm_isdst:
0404                         datax += offset_daylight
0405                     ##################################
0406 
0407                     dataxerr =  (float(last) - float(first))/2
0408                 graphlist[ig].SetPoint(ipoint, float(datax), float(datay) )
0409                 graphlist[ig].SetPointError(ipoint, float(dataxerr), float(datayerr) )
0410             else:
0411                 graphlist[ig].GetXaxis().SetBinLabel(ipoint +1 , str(datax) )
0412                 graphlist[ig].SetBinContent(ipoint +1, float(datay) )
0413                 graphlist[ig].SetBinError(ipoint +1, float(datayerr) )
0414 
0415             ipoint += 1
0416         if option.Time:
0417             ## print local time
0418             da_last.Print()
0419             print("GMT = " + str(time.strftime('%Y-%m-%d %H:%M:%S',time.gmtime(last - time.timezone))))
0420             graphlist[ig].GetXaxis().SetTimeDisplay(1);
0421             graphlist[ig].GetXaxis().SetTimeFormat("#splitline{%Y/%m/%d}{%H:%M}")           
0422         if option.graph:
0423             graphlist[ig].Draw('AP')
0424         else:
0425             graphlist[ig].Draw('P E1 X0')
0426             graphlist[ig].GetXaxis().SetTitle(graphXaxis)
0427             graphlist[ig].GetYaxis().SetTitle(graphYaxis[ig])
0428             #graphlist[ig].Fit('pol1')
0429         cvlist[ig].Update()
0430         cvlist[ig].SetGrid()
0431         if option.Print:
0432             suffix = ''
0433             if option.suffix:
0434                 suffix = option.suffix
0435             cvlist[ig].Print(graphnamelist[ig]+"_"+suffix+".png")
0436         if option.wait:
0437             raw_input( 'Press ENTER to continue\n ' )
0438         #graphlist[0].Print('all')
0439 
0440     if option.output:
0441         outroot = TFile(option.output,"RECREATE")
0442         for ig in graphlist:
0443             ig.Write()
0444 
0445         outroot.Close()
0446         print(" plots have been written to "+option.output)
0447 
0448 
0449 
0450     # CLEAN temporal files
0451     ###################################
0452     #os.system('rm tmp_beamspotdata.log')