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    ntuplemaker
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    -o, --output  = OUTPUT: filename of ROOT file with plots.
0028    -x, --xcrossing = XCROSSING : Bunch crossing number.
0029 
0030    Francisco Yumiceva (yumiceva@fnal.gov)
0031    Fermilab 2010
0032 
0033 """
0034 
0035 
0036 from builtins import range
0037 import os, string, re, sys, math
0038 import subprocess, time
0039 from BeamSpotObj import BeamSpot
0040 from IOVObj import IOV
0041 from CommonMethods import *
0042 
0043 try:
0044     import ROOT
0045 except:
0046     print("\nCannot load PYROOT, make sure you have setup ROOT in the path")
0047     print("and pyroot library is also defined in the variable PYTHONPATH, try:\n")
0048     if (os.getenv("PYTHONPATH")):
0049         print(" setenv PYTHONPATH ${PYTHONPATH}:$ROOTSYS/lib\n")
0050     else:
0051         print(" setenv PYTHONPATH $ROOTSYS/lib\n")
0052         sys.exit()
0053 
0054 from ROOT import gROOT, TFile, TTree
0055 from array import array
0056 
0057 def getFill( json, run ):
0058 
0059     thefill = 0
0060     run = int(run)
0061     keys = json.keys()
0062 
0063     for i in keys:
0064 
0065         run0 = int(json[i][0])
0066         run1 = int(json[i][1])
0067         if run>= run0 and run<=run1:
0068             thefill = i
0069 
0070     return int(thefill)
0071 
0072 if __name__ == '__main__':
0073 
0074 
0075 
0076     # fill and runs
0077     FillList = {}
0078     runsfile = open("FillandRuns.txt")
0079     for line in runsfile:
0080         if line.find('fill:') != -1:
0081             aline = line.split()
0082             afill = aline[1]
0083             run0 = aline[3]
0084             run1 = aline[5]
0085             FillList[int(afill)] = [int(run0),int(run1)]
0086 
0087     #print FillList
0088 
0089     # create ntuple
0090     gROOT.ProcessLine(
0091         "struct spot {\
0092         Float_t   position[3];\
0093         Float_t   posError[3];\
0094         Float_t   width[3];\
0095         Float_t   widthError[3];\
0096         Float_t   slope[2];\
0097         Float_t   slopeError[2];\
0098         Float_t   time[2];\
0099         Int_t     run;\
0100         Int_t     lumi[2];\
0101         Int_t     fill;\
0102         };" );
0103 
0104     bntuple = spot()
0105     fntuple = TFile( 'bntuple.root', 'RECREATE' )
0106     tbylumi = TTree( 'bylumi', 'beam spot data lumi by lumi' )
0107     tbylumi.Branch('fill', AddressOf( bntuple, 'fill'), 'fill/I' )
0108     tbylumi.Branch('run', AddressOf( bntuple, 'run'), 'run/I' )
0109     tbylumi.Branch('lumi', AddressOf( bntuple, 'lumi'), 'lumi[2]/I' )
0110     tbylumi.Branch('position', AddressOf( bntuple, 'position'),'position[3]/F')
0111     tbylumi.Branch('posErr', AddressOf( bntuple, 'posError'),'posError[3]/F')
0112     tbylumi.Branch('width', AddressOf( bntuple, 'width'),'width[3]/F')
0113     tbylumi.Branch('widthErr', AddressOf( bntuple, 'widthError'),'widthError[3]/F')
0114     tbylumi.Branch('slope', AddressOf( bntuple, 'slope'),'slope[2]/F')
0115     tbylumi.Branch('slopeErr', AddressOf( bntuple, 'slopeError'),'slopeError[2]/F')
0116     tbylumi.Branch('time', AddressOf( bntuple, 'time'),'time[2]/F')
0117 
0118     tbyIOV = TTree( 'byIOV', 'beam spot data by IOV' )
0119     tbyIOV.Branch('fill', AddressOf( bntuple, 'fill'), 'fill/I' )
0120     tbyIOV.Branch('run', AddressOf( bntuple, 'run'), 'run/I' )
0121     tbyIOV.Branch('lumi', AddressOf( bntuple, 'lumi'), 'lumi[2]/I' )
0122     tbyIOV.Branch('position', AddressOf( bntuple, 'position'),'position[3]/F')
0123     tbyIOV.Branch('posErr', AddressOf( bntuple, 'posError'),'posError[3]/F')
0124     tbyIOV.Branch('width', AddressOf( bntuple, 'width'),'width[3]/F')
0125     tbyIOV.Branch('widthErr', AddressOf( bntuple, 'widthError'),'widthError[3]/F')
0126     tbyIOV.Branch('slope', AddressOf( bntuple, 'slope'),'slope[2]/F')
0127     tbyIOV.Branch('slopeErr', AddressOf( bntuple, 'slopeError'),'slopeError[2]/F')
0128     tbyIOV.Branch('time', AddressOf( bntuple, 'time'),'time[2]/F')
0129 
0130     tbyrun = TTree( 'byrun', 'beam spot data by run' )
0131     tbyrun.Branch('fill', AddressOf( bntuple, 'fill'), 'fill/I' )
0132     tbyrun.Branch('run', AddressOf( bntuple, 'run'), 'run/I' )
0133     tbyrun.Branch('lumi', AddressOf( bntuple, 'lumi'), 'lumi[2]/I' )
0134     tbyrun.Branch('position', AddressOf( bntuple, 'position'),'position[3]/F')
0135     tbyrun.Branch('posErr', AddressOf( bntuple, 'posError'),'posError[3]/F')
0136     tbyrun.Branch('width', AddressOf( bntuple, 'width'),'width[3]/F')
0137     tbyrun.Branch('widthErr', AddressOf( bntuple, 'widthError'),'widthError[3]/F')
0138     tbyrun.Branch('slope', AddressOf( bntuple, 'slope'),'slope[2]/F')
0139     tbyrun.Branch('slopeErr', AddressOf( bntuple, 'slopeError'),'slopeError[2]/F')
0140     tbyrun.Branch('time', AddressOf( bntuple, 'time'),'time[2]/F')
0141 
0142 
0143     # COMMAND LINE OPTIONS
0144     #################################
0145     option,args = parse(__doc__)
0146     if not args and not option: exit()
0147 
0148     if not option.data:
0149         print(" need to provide beam spot data file")
0150         exit()
0151 
0152     if option.batch:
0153         ROOT.gROOT.SetBatch()
0154 
0155     datafilename = "tmp_beamspot.dat"
0156     if option.create:
0157         datafilename = option.create
0158 
0159     getDBdata = True
0160     if option.data:
0161         getDBdata = False
0162 
0163     IOVbase = 'lumibase'
0164     firstRun = "0:0"
0165     lastRun = "4999999999:4999999999"
0166 
0167     if option.initial:
0168         firstRun = option.initial
0169     if option.final:
0170         lastRun = option.final
0171 
0172     # GET IOVs
0173     ################################
0174 
0175     if getDBdata:
0176 
0177         print(" read DB to get list of IOVs for the given tag")
0178         acommand = 'cmscond_list_iov -c frontier://PromptProd/CMS_COND_31X_BEAMSPOT -P /afs/cern.ch/cms/DB/conddb -t '+ tag
0179         tmpstatus = subprocess.getstatusoutput( acommand )
0180         tmplistiov = tmpstatus[1].split('\n')
0181         #print tmplistiov
0182 
0183         iovlist = []
0184         passline = False
0185         iline = jline = 0
0186         totlines = len(tmplistiov)
0187         for line in tmplistiov:
0188 
0189             if line.find('since') != -1:
0190                 passline = True
0191                 jline = iline
0192             if passline and iline > jline and iline < totlines-1:
0193                 linedata = line.split()
0194                 #print linedata
0195                 aIOV = IOV()
0196                 aIOV.since = int(linedata[0])
0197                 aIOV.till = int(linedata[1])
0198                 iovlist.append( aIOV )
0199             iline += 1
0200 
0201         print(" total number of IOVs = " + str(len(iovlist)))
0202 
0203 
0204         #  GET DATA
0205         ################################
0206 
0207         otherArgs = ''
0208         if option.destDB:
0209             otherArgs = " -d " + option.destDB
0210             if option.auth:
0211                 otherArgs = otherArgs + " -a "+ option.auth
0212 
0213         print(" get beam spot data from DB for IOVs. This can take a few minutes ...")
0214 
0215         tmpfile = open(datafilename,'w')
0216 
0217         for iIOV in iovlist:
0218             passiov = False
0219             tmprunfirst = firstRun
0220             tmprunlast = lastRun
0221             tmplumifirst = 1
0222             tmplumilast = 9999999
0223             if IOVbase=="lumibase":
0224                 #tmprunfirst = int(firstRun.split(":")[0])
0225                 #tmprunlast  = int(lastRun.split(":")[0])
0226                 #tmplumifirst = int(firstRun.split(":")[1])
0227                 #tmplumilast  = int(lastRun.split(":")[1])
0228                 tmprunfirst = pack( int(firstRun.split(":")[0]) , int(firstRun.split(":")[1]) )
0229                 tmprunlast  = pack( int(lastRun.split(":")[0]) , int(lasstRun.split(":")[1]) )
0230             #print "since = " + str(iIOV.since) + " till = "+ str(iIOV.till)
0231             if iIOV.since >= int(tmprunfirst) and int(tmprunlast) < 0 and iIOV.since <= int(tmprunfirst):
0232                 print(" IOV: " + str(iIOV.since))
0233                 passiov = True
0234             if iIOV.since >= int(tmprunfirst) and int(tmprunlast) > 0 and iIOV.till <= int(tmprunlast):
0235                 print(" IOV: " + str(iIOV.since) + " to " + str(iIOV.till))
0236                 passiov = True
0237             if iIOV.since >= int(tmprunlast) and iIOV.till >= 4294967295:
0238                 print(" IOV: " + str(iIOV.since) + " to " + str(iIOV.till))
0239                 passiov = True
0240             if passiov:
0241                 acommand = 'getBeamSpotDB.py -t '+ tag + " -r " + str(iIOV.since) +otherArgs
0242                 if IOVbase=="lumibase":
0243                     tmprun = unpack(iIOV.since)[0]
0244                     tmplumi = unpack(iIOV.since)[1]
0245                     acommand = 'getBeamSpotDB.py -t '+ tag + " -r " + str(tmprun) +" -l "+tmplumi +otherArgs
0246                 status = subprocess.getstatusoutput( acommand )
0247                 tmpfile.write(status[1])
0248 
0249         print(" beam spot data collected and stored in file " + datafilename)
0250 
0251         tmpfile.close()
0252 
0253 
0254     # PROCESS DATA
0255     ###################################
0256 
0257     # check if input data exists if given
0258     if option.data:
0259         if os.path.isdir(option.data):
0260             tmp = subprocess.getstatusoutput("ls "+option.data)
0261             files = tmp[1].split()
0262             datafilename = "combined_all.txt"
0263             output = open(datafilename,"w")
0264 
0265             for f in files:
0266                 if os.path.isdir(option.data+"/"+f) is False:
0267                     input = open(option.data +"/"+f)
0268                     output.writelines(input.readlines())
0269             output.close()
0270             print(" data files have been collected in "+datafilename)
0271 
0272         elif os.path.exists(option.data):
0273             datafilename = option.data
0274         else:
0275             print(" input beam spot data DOES NOT exist, file " + option.data)
0276             exit()
0277 
0278     listbeam = []
0279 
0280     if option.xcrossing:
0281         listmap = readBeamSpotFile(datafilename,listbeam,IOVbase,firstRun,lastRun)
0282         # bx
0283         print("List of bunch crossings in the file:")
0284         print(listmap.keys())
0285         listbeam = listmap[option.Xrossing]
0286     else:
0287         readBeamSpotFile(datafilename,listbeam,IOVbase,firstRun,lastRun)
0288 
0289     sortAndCleanBeamList(listbeam,IOVbase)
0290 
0291 
0292     ###################################
0293 
0294     for ii in range(0,len(listbeam)):
0295 
0296         ibeam = listbeam[ii]
0297 
0298         bntuple.position = array('f', [float(ibeam.X), float(ibeam.Y), float(ibeam.Z)])
0299         bntuple.posError = array('f', [float(ibeam.Xerr),float(ibeam.Yerr),float(ibeam.Zerr)])
0300         bntuple.width = array('f', [float(ibeam.beamWidthX), float(ibeam.beamWidthY), float(ibeam.sigmaZ)])
0301         bntuple.widthError = array('f',[float(ibeam.beamWidthXerr),float(ibeam.beamWidthYerr),float(ibeam.sigmaZerr)])
0302         bntuple.run = int(ibeam.Run)
0303         bntuple.fill = int( getFill( FillList, int(ibeam.Run) ) )
0304         bntuple.lumi = array('i', [int(ibeam.IOVfirst),int(ibeam.IOVlast)])
0305         line = ibeam.IOVBeginTime
0306         begintime = time.mktime( time.strptime(line.split()[0] +  " " + line.split()[1] + " " + line.split()[2],"%Y.%m.%d %H:%M:%S %Z") )
0307         line = ibeam.IOVEndTime
0308         endtime = time.mktime( time.strptime(line.split()[0] +  " " + line.split()[1] + " " + line.split()[2],"%Y.%m.%d %H:%M:%S %Z") )
0309         bntuple.time = array('f', [begintime, endtime])
0310         tbylumi.Fill()
0311 
0312 
0313     iovlist = listbeam
0314     iovlist = createWeightedPayloads("tmp.txt",iovlist,False)
0315 
0316     for ii in range(0,len(iovlist)):
0317 
0318         ibeam = iovlist[ii]
0319 
0320         bntuple.position = array('f', [float(ibeam.X), float(ibeam.Y), float(ibeam.Z)])
0321         bntuple.posError = array('f', [float(ibeam.Xerr),float(ibeam.Yerr),float(ibeam.Zerr)])
0322         bntuple.width = array('f', [float(ibeam.beamWidthX), float(ibeam.beamWidthY), float(ibeam.sigmaZ)])
0323         bntuple.widthError = array('f',[float(ibeam.beamWidthXerr),float(ibeam.beamWidthYerr),float(ibeam.sigmaZerr)])
0324         bntuple.run = int(ibeam.Run)
0325         bntuple.fill = int( getFill( FillList, int(ibeam.Run) ) )
0326         bntuple.lumi = array('i', [int(ibeam.IOVfirst),int(ibeam.IOVlast)])
0327         line = ibeam.IOVBeginTime
0328         begintime = time.mktime( time.strptime(line.split()[0] +  " " + line.split()[1] + " " + line.split()[2],"%Y.%m.%d %H:%M:%S %Z") )
0329         line = ibeam.IOVEndTime
0330         endtime = time.mktime( time.strptime(line.split()[0] +  " " + line.split()[1] + " " + line.split()[2],"%Y.%m.%d %H:%M:%S %Z") )
0331         bntuple.time = array('f', [begintime, endtime])
0332 
0333         tbyIOV.Fill()
0334 
0335     weightedlist = listbeam
0336     weightedlist = createWeightedPayloads("tmp.txt",weightedlist,True)
0337 
0338     for ii in range(0,len(weightedlist)):
0339 
0340         ibeam = weightedlist[ii]
0341 
0342         bntuple.position = array('f', [float(ibeam.X), float(ibeam.Y), float(ibeam.Z)])
0343         bntuple.posError = array('f', [float(ibeam.Xerr),float(ibeam.Yerr),float(ibeam.Zerr)])
0344         bntuple.width = array('f', [float(ibeam.beamWidthX), float(ibeam.beamWidthY), float(ibeam.sigmaZ)])
0345         bntuple.widthError = array('f',[float(ibeam.beamWidthXerr),float(ibeam.beamWidthYerr),float(ibeam.sigmaZerr)])
0346         bntuple.run = int(ibeam.Run)
0347         bntuple.fill = int( getFill( FillList, int(ibeam.Run) ) )
0348         bntuple.lumi = array('i', [int(ibeam.IOVfirst),int(ibeam.IOVlast)])
0349         line = ibeam.IOVBeginTime
0350         begintime = time.mktime( time.strptime(line.split()[0] +  " " + line.split()[1] + " " + line.split()[2],"%Y.%m.%d %H:%M:%S %Z") )
0351         line = ibeam.IOVEndTime
0352         endtime = time.mktime( time.strptime(line.split()[0] +  " " + line.split()[1] + " " + line.split()[2],"%Y.%m.%d %H:%M:%S %Z") )
0353         bntuple.time = array('f', [begintime, endtime])
0354 
0355         tbyrun.Fill()
0356 
0357 
0358     os.system('rm tmp.txt')
0359     fntuple.cd()
0360     tbylumi.Write()
0361     tbyIOV.Write()
0362     tbyrun.Write()
0363     fntuple.Close()
0364 
0365     # CLEAN temporal files
0366     ###################################
0367     #os.system('rm tmp_beamspotdata.log')