Reference Manual

Go to the documentation of this file.
00001 #!/usr/bin/env python
00002 #____________________________________________________________
00003 #
00004 #
00005 # A very simple way to make plots with ROOT via an XML file
00006 #
00007 # Francisco Yumiceva
00008 # yumiceva@fnal.gov
00009 #
00010 # Fermilab, 2010
00011 #
00012 #____________________________________________________________
00013 
00014 """
00015    plotBeamSpotDB
00016 
00017    A very simple script to plot the beam spot data stored in condDB
00018 
00019    usage: %prog -t <tag name>
00020    -a, --auth    = AUTH: DB authorization path. online(/nfshome0/popcondev/conddb).
00021    -b, --batch : Run ROOT in batch mode.
00022    -c, --create  = CREATE: name for beam spot data file.
00023    -d, --data    = DATA: input beam spot data file.
00024    -D, --destDB  = DESTDB: destination DB string. online(oracle://cms_orcon_prod/CMS_COND_31X_BEAMSPOT).
00025    -i, --initial = INITIAL: First IOV. Options: run number, or run:lumi, eg. \"133200:21\"
00026    -f, --final   = FINAL: Last IOV. Options: run number, or run:lumi
00027    -g, --graph : create a TGraphError instead of a TH1 object
00028    -n, --noplot : Only extract beam spot data, plots are not created..
00029    -o, --output  = OUTPUT: filename of ROOT file with plots.
00030    -p, --payload = PAYLOAD: filename of output text file. Combine and splits lumi IOVs.
00031    -P, --Print : create PNG plots from canvas.
00032    -s, --suffix = SUFFIX: Suffix will be added to plots filename.
00033    -t, --tag     = TAG: Database tag name.
00034    -T, --Time : create plots with time axis.
00035    -I, --IOVbase = IOVBASE: options: runbase(default), lumibase, timebase
00036    -w, --wait : Pause script after plotting a new histograms.
00037    -W, --weighted : Create a weighted result for a range of lumi IOVs, skip lumi IOV combination and splitting.
00038    -x, --xcrossing = XCROSSING : Bunch crossing number.
00039    
00040    Francisco Yumiceva (yumiceva@fnal.gov)
00041    Fermilab 2010
00042    
00043 """
00044 
00045 
00046 import os, string, re, sys, math
00047 import commands, time
00048 from BeamSpotObj import BeamSpot
00049 from IOVObj import IOV
00050 from CommonMethods import *
00051 
00052 try:
00053     import ROOT
00054 except:
00055     print "\nCannot load PYROOT, make sure you have setup ROOT in the path"
00056     print "and pyroot library is also defined in the variable PYTHONPATH, try:\n"
00057     if (os.getenv("PYTHONPATH")):
00058         print " setenv PYTHONPATH ${PYTHONPATH}:$ROOTSYS/lib\n"
00059     else:
00060         print " setenv PYTHONPATH $ROOTSYS/lib\n"
00061         sys.exit()
00062 
00063 from ROOT import TFile, TGraphErrors, TGaxis, TDatime
00064 from ROOT import TCanvas, TH1F
00065 
00066 # ROOT STYLE
00067 #############################
00068 def SetStyle():
00069 
00070     # canvas
00071     ROOT.gStyle.SetCanvasBorderMode(0)
00072     ROOT.gStyle.SetCanvasColor(0)
00073     ROOT.gStyle.SetCanvasDefH(600)
00074     ROOT.gStyle.SetCanvasDefW(600)
00075     ROOT.gStyle.SetCanvasDefX(0)
00076     ROOT.gStyle.SetCanvasDefY(0)
00077     # pad
00078     ROOT.gStyle.SetPadBorderMode(0)
00079     ROOT.gStyle.SetPadColor(0)
00080     ROOT.gStyle.SetPadGridX(False)
00081     ROOT.gStyle.SetPadGridY(False)
00082     ROOT.gStyle.SetGridColor(0)
00083     ROOT.gStyle.SetGridStyle(3)
00084     ROOT.gStyle.SetGridWidth(1)
00085                   
00086     ROOT.gStyle.SetFrameBorderMode(0)
00087     ROOT.gStyle.SetFrameFillColor(0)
00088     ROOT.gStyle.SetTitleColor(1)
00089     ROOT.gStyle.SetStatColor(0)
00090 
00091     # set the paper & margin sizes
00092     ROOT.gStyle.SetPaperSize(20,26)
00093     ROOT.gStyle.SetPadTopMargin(0.04)
00094     ROOT.gStyle.SetPadRightMargin(0.04)
00095     ROOT.gStyle.SetPadBottomMargin(0.14)
00096     ROOT.gStyle.SetPadLeftMargin(0.11)
00097     ROOT.gStyle.SetPadTickX(1)
00098     ROOT.gStyle.SetPadTickY(1)
00099 
00100     ROOT.gStyle.SetTextFont(42) #132
00101     ROOT.gStyle.SetTextSize(0.09)
00102     ROOT.gStyle.SetLabelFont(42,"xyz")
00103     ROOT.gStyle.SetTitleFont(42,"xyz")
00104     ROOT.gStyle.SetLabelSize(0.035,"xyz")
00105     ROOT.gStyle.SetTitleSize(0.045,"xyz")
00106     ROOT.gStyle.SetTitleOffset(1.1,"y")
00107 
00108     # use bold lines and markers
00109     ROOT.gStyle.SetMarkerStyle(8)
00110     ROOT.gStyle.SetHistLineWidth(2)
00111     ROOT.gStyle.SetLineWidth(1)
00112     #ROOT.gStyle.SetLineStyleString(2,"[12 12]") // postscript dashes
00113 
00114     ROOT.gStyle.SetMarkerSize(0.6)
00115     
00116     # do not display any of the standard histogram decorations
00117     ROOT.gStyle.SetOptTitle(0)
00118     ROOT.gStyle.SetOptStat(0) #("m")
00119     ROOT.gStyle.SetOptFit(0)
00120     
00121     #ROOT.gStyle.SetPalette(1,0)
00122     ROOT.gStyle.cd()
00123     ROOT.gROOT.ForceStyle()
00124 #########################################
00125 
00126 
00127 if __name__ == '__main__':
00128 
00129     # style
00130     SetStyle()
00131     
00132     printCanvas = False
00133     printFormat = "png"
00134     printBanner = False
00135     Banner = "CMS Preliminary"
00136 
00137     # COMMAND LINE OPTIONS
00138     #################################
00139     option,args = parse(__doc__)
00140     if not args and not option: exit()
00141     
00142     tag = ''
00143     if not option.tag and not option.data: 
00144         print " need to provide DB tag name or beam spot data file"
00145         exit()
00146     else:
00147         tag = option.tag
00148 
00149     if option.batch:
00150         ROOT.gROOT.SetBatch()
00151         
00152     datafilename = "tmp_beamspot.dat"
00153     if option.create:
00154         datafilename = option.create
00155     
00156     getDBdata = True
00157     if option.data:
00158         getDBdata = False
00159    
00160     IOVbase = 'runbase'
00161     if option.IOVbase:
00162         if option.IOVbase != "runbase" and option.IOVbase != "lumibase" and option.IOVbase != "timebase":
00163             print "\n\n unknown iov base option: "+ option.IOVbase +" \n\n\n"
00164             exit()
00165         IOVbase = option.IOVbase
00166     
00167     firstRun = "0"
00168     lastRun  = "4999999999"
00169     if IOVbase == "lumibase":
00170         firstRun = "0:0"
00171         lastRun = "4999999999:4999999999"
00172     
00173     if option.initial:
00174         firstRun = option.initial
00175     if option.final:
00176         lastRun = option.final
00177         
00178     # GET IOVs
00179     ################################
00180 
00181     if getDBdata:
00182 
00183         print " read DB to get list of IOVs for the given tag"
00184         mydestdb = 'frontier://PromptProd/CMS_COND_31X_BEAMSPOT'
00185         if option.destDB:
00186                 mydestdb = option.destDB
00187         acommand = 'cmscond_list_iov -c '+mydestdb+' -P /afs/cern.ch/cms/DB/conddb -t '+ tag
00188         tmpstatus = commands.getstatusoutput( acommand )
00189         tmplistiov = tmpstatus[1].split('\n')
00190         #print tmplistiov
00191 
00192         iovlist = []
00193         passline = False
00194         iline = jline = 0
00195         totlines = len(tmplistiov)
00196         for line in tmplistiov:
00197         
00198             if line.find('since') != -1:
00199                 passline = True
00200                 jline = iline
00201             if passline and iline > jline and iline < totlines-1:
00202                 linedata = line.split()
00203                 #print linedata
00204                 aIOV = IOV()
00205                 aIOV.since = int(linedata[0])
00206                 aIOV.till = int(linedata[1])
00207                 iovlist.append( aIOV )
00208             iline += 1
00209     
00210         print " total number of IOVs = " + str(len(iovlist))
00211 
00212 
00213         #  GET DATA
00214         ################################
00215 
00216         otherArgs = ''
00217         if option.destDB:
00218             otherArgs = " -d " + option.destDB
00219             if option.auth:
00220                 otherArgs = otherArgs + " -a "+ option.auth
00221         
00222         print " get beam spot data from DB for IOVs. This can take a few minutes ..."
00223 
00224         tmpfile = open(datafilename,'w')
00225 
00226         for iIOV in iovlist:
00227             passiov = False
00228             tmprunfirst = firstRun
00229             tmprunlast = lastRun
00230             tmplumifirst = 1
00231             tmplumilast = 9999999
00232             if IOVbase=="lumibase":
00233                 #tmprunfirst = int(firstRun.split(":")[0])
00234                 #tmprunlast  = int(lastRun.split(":")[0])
00235                 #tmplumifirst = int(firstRun.split(":")[1])
00236                 #tmplumilast  = int(lastRun.split(":")[1])
00237                 tmprunfirst = pack( int(firstRun.split(":")[0]) , int(firstRun.split(":")[1]) )
00238                 tmprunlast  = pack( int(lastRun.split(":")[0]) , int(lastRun.split(":")[1]) )
00239             #print "since = " + str(iIOV.since) + " till = "+ str(iIOV.till)
00240             if iIOV.since >= int(tmprunfirst) and int(tmprunlast) < 0 and iIOV.since <= int(tmprunfirst):
00241                 print " IOV: " + str(iIOV.since)
00242                 passiov = True
00243             if iIOV.since >= int(tmprunfirst) and int(tmprunlast) > 0 and iIOV.till <= int(tmprunlast):
00244                 print " a IOV: " + str(iIOV.since) + " to " + str(iIOV.till)
00245                 passiov = True
00246             #if iIOV.since >= int(tmprunlast) and iIOV.till >= 4294967295:
00247             #    print " b IOV: " + str(iIOV.since) + " to " + str(iIOV.till)
00248             #    passiov = True                
00249             if passiov:
00250                 acommand = 'getBeamSpotDB.py -t '+ tag + " -r " + str(iIOV.since) +otherArgs
00251                 if IOVbase=="lumibase":
00252                     tmprun = unpack(iIOV.since)[0]
00253                     tmplumi = unpack(iIOV.since)[1]
00254                     acommand = 'getBeamSpotDB.py -t '+ tag + " -r " + str(tmprun) +" -l "+str(tmplumi) +otherArgs
00255                     print acommand
00256                 status = commands.getstatusoutput( acommand )
00257                 tmpfile.write(status[1])
00258     
00259         print " beam spot data collected and stored in file " + datafilename
00260     
00261         tmpfile.close()
00262 
00263     
00264     # PROCESS DATA
00265     ###################################
00266 
00267     # check if input data exists if given
00268     if option.data:
00269         if os.path.isdir(option.data):
00270             tmp = commands.getstatusoutput("ls "+option.data)
00271             files = tmp[1].split()
00272             datafilename = "combined_all.txt"
00273             output = open(datafilename,"w")
00274             
00275             for f in files:
00276                 input = open(option.data +"/"+f)
00277                 output.writelines(input.readlines())
00278             output.close()
00279             print " data files have been collected in "+datafilename
00280             
00281         elif os.path.exists(option.data):
00282             datafilename = option.data
00283         else:
00284             print " input beam spot data DOES NOT exist, file " + option.data
00285             exit()
00286 
00287     listbeam = []
00288 
00289     if option.xcrossing:
00290         listmap = readBeamSpotFile(datafilename,listbeam,IOVbase,firstRun,lastRun)
00291         # bx
00292         print "List of bunch crossings in the file:"
00293         print listmap.keys()
00294         listbeam = listmap[option.Xrossing]
00295     else:
00296         readBeamSpotFile(datafilename,listbeam,IOVbase,firstRun,lastRun)
00297     
00298     sortAndCleanBeamList(listbeam,IOVbase)
00299             
00300     if IOVbase == "lumibase" and option.payload:
00301         weighted = True;
00302         if not option.weighted:
00303             weighted = False
00304         createWeightedPayloads(option.payload,listbeam,weighted)
00305     if option.noplot:
00306         print " no plots requested, exit now."
00307         sys.exit()
00308     # MAKE PLOTS
00309     ###################################    
00310     TGaxis.SetMaxDigits(8)
00311 
00312     graphlist = []
00313     graphnamelist = ['X','Y','Z','SigmaZ','dxdz','dydz','beamWidthX', 'beamWidthY']
00314     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']
00315     graphXaxis = 'Run number'
00316     if IOVbase == 'runbase':
00317         graphXaxis = "Run number"
00318     if IOVbase == 'lumibase':
00319         graphXaxis = 'Lumi section'
00320     if IOVbase == 'timebase' or option.Time:
00321         graphXaxis = "Time"
00322         #dh = ROOT.TDatime(2010,06,01,00,00,00)
00323         ROOT.gStyle.SetTimeOffset(0) #dh.Convert())
00324                     
00325     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]']
00326 
00327     cvlist = []
00328 
00329     for ig in range(0,8):
00330         cvlist.append( TCanvas(graphnamelist[ig],graphtitlelist[ig], 1200, 600) )
00331         if option.graph:
00332             graphlist.append( TGraphErrors( len(listbeam) ) )
00333         else:
00334             graphlist.append( TH1F("name","title",len(listbeam),0,len(listbeam)) )
00335         
00336         graphlist[ig].SetName(graphnamelist[ig])
00337         graphlist[ig].SetTitle(graphtitlelist[ig])
00338         ipoint = 0
00339         for ii in range(0,len(listbeam)):
00340             
00341             ibeam = listbeam[ii]
00342             datax = dataxerr = 0.
00343             datay = datayerr = 0.
00344             if graphnamelist[ig] == 'X':
00345                 datay = ibeam.X
00346                 datayerr = ibeam.Xerr
00347             if graphnamelist[ig] == 'Y':
00348                 datay = ibeam.Y
00349                 datayerr = ibeam.Yerr
00350             if graphnamelist[ig] == 'Z':
00351                 datay = ibeam.Z
00352                 datayerr = ibeam.Zerr
00353             if graphnamelist[ig] == 'SigmaZ':
00354                 datay = ibeam.sigmaZ
00355                 datayerr = ibeam.sigmaZerr
00356             if graphnamelist[ig] == 'dxdz':
00357                 datay = ibeam.dxdz
00358                 datayerr = ibeam.dxdzerr
00359             if graphnamelist[ig] == 'dydz':
00360                 datay = ibeam.dydz
00361                 datayerr = ibeam.dydzerr
00362             if graphnamelist[ig] == 'beamWidthX':
00363                 datay = ibeam.beamWidthX
00364                 datayerr = ibeam.beamWidthXerr
00365             if graphnamelist[ig] == 'beamWidthY':
00366                 datay = ibeam.beamWidthY
00367                 datayerr = ibeam.beamWidthYerr
00368 
00369             datax = ibeam.IOVfirst
00370             if IOVbase=="lumibase":
00371                 datax = str(ibeam.Run) + ":" + str(ibeam.IOVfirst)
00372                 if ibeam.IOVfirst != ibeam.IOVlast:
00373                     datax = str(ibeam.Run) + ":" + str(ibeam.IOVfirst)+"-"+str(ibeam.IOVlast)
00374             #print datax
00375             if option.graph:
00376                 if IOVbase=="lumibase":
00377                     #first = int( pack( int(ibeam.Run) , int(ibeam.IOVfirst) ) )
00378                     #last = int( pack( int(ibeam.Run) , int(ibeam.IOVlast) ) )
00379                     first = ibeam.IOVfirst
00380                     last = ibeam.IOVlast
00381                     if option.Time:
00382                         atime = ibeam.IOVBeginTime
00383                         first = time.mktime( time.strptime(atime.split()[0] +  " " + atime.split()[1] + " " + atime.split()[2],"%Y.%m.%d %H:%M:%S %Z") )
00384                         atime = ibeam.IOVEndTime
00385                         last = time.mktime( time.strptime(atime.split()[0] +  " " + atime.split()[1] + " " + atime.split()[2],"%Y.%m.%d %H:%M:%S %Z") )
00386                         da_first = TDatime(time.strftime('%Y-%m-%d %H:%M:%S',time.localtime(first - time.timezone)))
00387                         da_last = TDatime(time.strftime('%Y-%m-%d %H:%M:%S',time.localtime(last - time.timezone)))
00388                         if ipoint == 0:
00389                             ## print local time
00390                             da_first.Print()
00391                             ## print gmt time
00392                             print "GMT = " + str(time.strftime('%Y-%m-%d %H:%M:%S',time.gmtime(first - time.timezone)))
00393                             reftime = first
00394                             ptm = time.localtime(reftime)
00395                             da = TDatime(time.strftime('%Y-%m-%d %H:%M:%S',ptm))
00396                             if time.daylight and ptm.tm_isdst:
00397                                 offset_daylight = time.timezone - time.altzone
00398                             ROOT.gStyle.SetTimeOffset(da.Convert(1) - 3600)
00399 
00400                     datax = (float(last) - float(first))/2 + float(first) - da.Convert() + 3600
00401                     ## Comment out this block if running on Mac ##
00402                     if time.daylight and ptm.tm_isdst:
00403                         datax += offset_daylight
00404                     ##################################
00405                     
00406                     dataxerr =  (float(last) - float(first))/2
00407                 graphlist[ig].SetPoint(ipoint, float(datax), float(datay) )
00408                 graphlist[ig].SetPointError(ipoint, float(dataxerr), float(datayerr) )
00409             else:
00410                 graphlist[ig].GetXaxis().SetBinLabel(ipoint +1 , str(datax) )
00411                 graphlist[ig].SetBinContent(ipoint +1, float(datay) )
00412                 graphlist[ig].SetBinError(ipoint +1, float(datayerr) )
00413 
00414             ipoint += 1
00415         if option.Time:
00416             ## print local time
00417             da_last.Print()
00418             print "GMT = " + str(time.strftime('%Y-%m-%d %H:%M:%S',time.gmtime(last - time.timezone)))
00419             graphlist[ig].GetXaxis().SetTimeDisplay(1);
00420             graphlist[ig].GetXaxis().SetTimeFormat("#splitline{%Y/%m/%d}{%H:%M}")           
00421         if option.graph:
00422             graphlist[ig].Draw('AP')
00423         else:
00424             graphlist[ig].Draw('P E1 X0')
00425             graphlist[ig].GetXaxis().SetTitle(graphXaxis)
00426             graphlist[ig].GetYaxis().SetTitle(graphYaxis[ig])
00427             #graphlist[ig].Fit('pol1')
00428         cvlist[ig].Update()
00429         cvlist[ig].SetGrid()
00430         if option.Print:
00431             suffix = ''
00432             if option.suffix:
00433                 suffix = option.suffix
00434             cvlist[ig].Print(graphnamelist[ig]+"_"+suffix+".png")
00435         if option.wait:
00436             raw_input( 'Press ENTER to continue\n ' )
00437         #graphlist[0].Print('all')
00438 
00439     if option.output:
00440         outroot = TFile(option.output,"RECREATE")
00441         for ig in graphlist:
00442             ig.Write()
00443 
00444         outroot.Close()
00445         print " plots have been written to "+option.output
00446         
00447     
00448 
00449     # CLEAN temporal files
00450     ###################################
00451     #os.system('rm tmp_beamspotdata.log')
/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_5_3_13_patch3/src/RecoVertex/BeamSpotProducer/scripts/plotBeamSpotDB.py
