Reference Manual

Go to the documentation of this file.
00001 #!/usr/bin/env python
00002 #
00003 # dump all fills into files.
00004 # allfills.txt all the existing fills.
00005 # fill_num.txt all the runs in the fill
00006 # dumpFill -o outputdir
00007 # dumpFill -f fillnum generate runlist for the given fill
00008 #
00009 import os,os.path,sys,math,array,datetime,time,re
00010 import coral
00011 
00012 from RecoLuminosity.LumiDB import argparse,lumiTime,CommonUtil,lumiQueryAPI
00013 
00014 allfillname='allfills.txt'
00015         
00016 def calculateSpecificLumi(lumi,lumierr,beam1intensity,beam1intensityerr,beam2intensity,beam2intensityerr):
00017     '''
00018     '''
00019     specificlumi=0.0
00020     specificlumierr=0.0
00021     if lumi!=0.0 and beam1intensity!=0.0 and  beam2intensity!=0.0:
00022         specificlumi=float(lumi)/(float(beam1intensity)*float(beam2intensity))
00023         specificlumierr=specificlumi*math.sqrt(lumierr**2/lumi**2+beam1intensityerr**2/beam1intensity**2+beam2intensityerr**2/beam2intensity**2)
00024     return (specificlumi,specificlumierr)
00025 
00026 def getFillFromDB(dbsession,parameters,fillnum):
00027     runtimesInFill={}
00028     q=dbsession.nominalSchema().newQuery()
00029     fillrundict=lumiQueryAPI.runsByfillrange(q,fillnum,fillnum)
00030     del q
00031     if len(fillrundict)>0:
00032         for fill,runs in  fillrundict.items():
00033             for run in runs:
00034                 q=dbsession.nominalSchema().newQuery()
00035                 rresult=lumiQueryAPI.runsummaryByrun(q,run)
00036                 del q
00037                 if len(rresult)==0: continue
00038                 runtimesInFill[run]=rresult[3]
00039     return runtimesInFill
00040 
00041 def getFillFromFile(fillnum,inputdir):
00042     runtimesInFill={}
00043     #look for files 'fill_num.txt' in inputdir
00044     for filename in os.listdir(inputdir):
00045         filename=filename.strip()
00046         if filename.find('.')==-1: continue            
00047         basename,extension=filename.split('.')        
00048         if not extension or extension!='txt':
00049             continue
00050         if basename.find('_')==-1: continue
00051         prefix,number=basename.split('_')
00052         if not number : continue
00053         if fillnum!=int(number):continue
00054         f=open(os.path.join(inputdir,'fill_'+number+'.txt'),'r')
00055         for line in f:
00056             l=line.strip()
00057             fields=l.split(',')
00058             if len(fields)<2 : continue
00059             runtimesInFill[int(fields[0])]=fields[1]
00060         f.close()
00061     return runtimesInFill
00062 
00063 def getSpecificLumi(dbsession,parameters,fillnum,inputdir):
00064     '''
00065     specific lumi in 1e-30 (ub-1s-1) unit
00066     lumidetail occlumi in 1e-27
00067     1309_lumi_401_CMS.txt
00068     time(in seconds since January 1,2011,00:00:00 UTC) stab(fraction of time spent in stable beams for this time bin) l(lumi in Hz/ub) dl(point-to-point error on lumi in Hz/ub) sl(specific lumi in Hz/ub) dsl(error on specific lumi)
00069     20800119.0 1 -0.889948 0.00475996848729 0.249009 0.005583287562 -0.68359 6.24140208607 0.0 0.0 0.0 0.0 0.0 0.0 0.0383576 0.00430892097862 0.0479095 0.00430892097862 66.6447 4.41269758764 0.0 0.0 0.0
00070     result [(time,beamstatusfrac,lumi,lumierror,speclumi,speclumierror)]
00071     '''
00072     #result=[]
00073     runtimesInFill=getFillFromFile(fillnum,inputdir)#{runnum:starttimestr}
00074     beamstatusDict={}#{runnum:{(startorbit,cmslsnum):beamstatus}}
00075     t=lumiTime.lumiTime()
00076     fillbypos={}#{bxidx:(lstime,beamstatusfrac,lumi,lumierror,specificlumi,specificlumierror)}
00077     #referencetime=time.mktime(datetime.datetime(2010,1,1,0,0,0).timetuple())
00078     referencetime=0
00079     if fillnum and len(runtimesInFill)==0:
00080         runtimesInFill=getFillFromDB(dbsession,parameters,fillnum)#{runnum:starttimestr}
00081     #precheck
00082     totalstablebeamLS=0
00083     for runnum in runtimesInFill.keys():
00084         q=dbsession.nominalSchema().newQuery()
00085         runinfo=lumiQueryAPI.lumisummaryByrun(q,runnum,'0001',beamstatus=None)
00086         del q
00087         lsbeamstatusdict={}
00088         for perlsdata in runinfo:
00089             cmslsnum=perlsdata[0]
00090             startorbit=perlsdata[3]
00091             beamstatus=perlsdata[4]
00092             lsbeamstatusdict[(startorbit,cmslsnum)]=beamstatus            
00093             #print (startorbit,cmslsnum),beamstatus
00094             if beamstatus=='STABLE BEAMS':
00095                 totalstablebeamLS+=1
00096         beamstatusDict[runnum]=lsbeamstatusdict
00097     if totalstablebeamLS<10:#less than 10 LS in a fill has 'stable beam', it's no a good fill
00098         print 'fill ',fillnum,' , having less than 10 stable beam lS, is not good, skip'
00099         return fillbypos
00100     for runnum,starttime in runtimesInFill.items():
00101         if not runtimesInFill.has_key(runnum):
00102             print 'run '+str(runnum)+' does not exist'
00103             continue
00104         q=dbsession.nominalSchema().newQuery()
00105         occlumidata=lumiQueryAPI.calibratedDetailForRunLimitresult(q,parameters,runnum)#{(startorbit,cmslsnum):[(bxidx,lumivalue,lumierr)]} #values after cut
00106         del q
00107         #print occlumidata
00108         q=dbsession.nominalSchema().newQuery()
00109         beamintensitydata=lumiQueryAPI.beamIntensityForRun(q,parameters,runnum)#{startorbit:[(bxidx,beam1intensity,beam2intensity)]}
00110         #print 'beamintensity for run ',runnum,beamintensitydata
00111         del q
00112         for (startorbit,cmslsnum),lumilist in occlumidata.items():
00113             if len(lumilist)==0: continue
00114             beamstatusflag=beamstatusDict[runnum][(startorbit,cmslsnum)]
00115             beamstatusfrac=0.0
00116             if beamstatusflag=='STABLE BEAMS':
00117                 beamstatusfrac=1.0
00118             lstimestamp=t.OrbitToUTCTimestamp(starttime,startorbit)
00119             if beamintensitydata.has_key(startorbit) and len(beamintensitydata[startorbit])>0:
00120                 for lumidata in lumilist:
00121                     bxidx=lumidata[0]
00122                     lumi=lumidata[1]
00123                     lumierror=lumidata[2]                    
00124                     for beamintensitybx in beamintensitydata[startorbit]:
00125                         if beamintensitybx[0]==bxidx:
00126                             if not fillbypos.has_key(bxidx):
00127                                 fillbypos[bxidx]=[]
00128                             beam1intensity=beamintensitybx[1]
00129                             beam2intensity=beamintensitybx[2]
00130                             speclumi=calculateSpecificLumi(lumi,lumierror,beam1intensity,0.0,beam2intensity,0.0)
00131                             fillbypos[bxidx].append([lstimestamp-referencetime,beamstatusfrac,lumi,lumierror,beam1intensity,beam2intensity,speclumi[0],speclumi[1]])
00132     #print fillbypos
00133     return fillbypos
00134 
00135 #####output methods####
00136 def filltofiles(allfills,runsperfill,runtimes,dirname):
00137     f=open(os.path.join(dirname,allfillname),'w')
00138     for fill in allfills:
00139         print >>f,'%d'%(fill)
00140     f.close()
00141     for fill,runs in runsperfill.items():
00142         filename='fill_'+str(fill)+'.txt'
00143         if len(runs)!=0:
00144             f=open(os.path.join(dirname,filename),'w')
00145             for run in runs:
00146                 print >>f,'%d,%s'%(run,runtimes[run])
00147             f.close()
00148             
00149 def specificlumiTofile(fillnum,filldata,outdir):
00150     timedict={}#{lstime:[[stablebeamfrac,lumi,lumierr,speclumi,speclumierr]]}
00151     #
00152     #check outdir/fillnum subdir exists; if not, create it; else outdir=outdir/fillnum
00153     #
00154     filloutdir=os.path.join(outdir,str(fillnum))
00155     print 'filloutdir ',filloutdir
00156     if not os.path.exists(filloutdir):
00157         os.mkdir(filloutdir)
00158     for cmsbxidx,perbxdata in filldata.items():
00159         lhcbucket=0
00160         if cmsbxidx!=0:
00161             lhcbucket=(cmsbxidx-1)*10+1
00162         a=sorted(perbxdata,key=lambda x:x[0])
00163         filename=str(fillnum)+'_lumi_'+str(lhcbucket)+'_CMS.txt'
00164         linedata=[]
00165         for perlsdata in a:
00166             ts=int(perlsdata[0])
00167             beamstatusfrac=perlsdata[1]
00168             lumi=perlsdata[2]
00169             lumierror=perlsdata[3]
00170             beam1intensity=perlsdata[4]
00171             beam2intensity=perlsdata[5]
00172             speclumi=perlsdata[6]
00173             speclumierror= perlsdata[7]
00174             if lumi>0 and lumierror>0 and speclumi>0:
00175                 linedata.append([ts,beamstatusfrac,lumi,lumierror,speclumi,speclumierror])
00176             if not timedict.has_key(ts):
00177                   timedict[ts]=[]
00178             timedict[ts].append([beamstatusfrac,lumi,lumierror,speclumi,speclumierror])
00179         if len(linedata)>10:#at least 10 good ls
00180             f=open(os.path.join(filloutdir,filename),'w')
00181             for line in linedata:
00182                 print >>f, '%d\t%e\t%e\t%e\t%e\t%e'%(line[0],line[1],line[2],line[3],line[4],line[5])
00183             f.close()
00184     #print 'writing avg file'
00185     summaryfilename=str(fillnum)+'_lumi_CMS.txt'
00186     f=None
00187     lstimes=timedict.keys()
00188     lstimes.sort()
00189     fillseg=[]
00190     lscounter=0
00191     for lstime in lstimes:
00192         allvalues=timedict[lstime]
00193         transposedvalues=CommonUtil.transposed(allvalues,0.0)
00194         bstatfrac=transposedvalues[0][0]#beamstatus does not change with bx position
00195         lumivals=transposedvalues[1]
00196         lumitot=sum(lumivals)
00197         if bstatfrac==1.0 :
00198             fillseg.append([lstime,lumitot])
00199         lumierrs=transposedvalues[2]
00200         lumierrortot=math.sqrt(sum(map(lambda x:x**2,lumierrs)))
00201         specificvals=transposedvalues[3]
00202         specificavg=sum(specificvals)/float(len(specificvals))#avg spec lumi
00203         specificerrs=transposedvalues[4]
00204         specifictoterr=math.sqrt(sum(map(lambda x:x**2,specificerrs)))
00205         specificerravg=specifictoterr/float(len(specificvals))
00206         if lscounter==0:
00207             f=open(os.path.join(filloutdir,summaryfilename),'w')
00208         lscounter+=1
00209         print >>f,'%d\t%e\t%e\t%e\t%e\t%e'%(lstime,bstatfrac,lumitot,lumierrortot,specificavg,specificerravg)
00210     if f is not None:
00211         f.close()
00212     #print 'writing summary file'
00213     fillsummaryfilename=str(fillnum)+'_summary_CMS.txt'
00214     f=open(os.path.join(filloutdir,fillsummaryfilename),'w')    
00215     if len(fillseg)==0:
00216         print >>f,'%s'%('#no stable beams')
00217         f.close()
00218         return
00219     previoustime=fillseg[0][0]
00220     boundarytime=fillseg[0][0]
00221     #print 'boundary time ',boundarytime
00222     summaryls={}
00223     summaryls[boundarytime]=[]
00224     for [lstime,lumitot] in fillseg:#fillseg is everything with stable beam flag
00225         if lstime-previoustime>50.0:
00226             boundarytime=lstime
00227             #print 'found new boundary ',boundarytime
00228             summaryls[boundarytime]=[]
00229      #   print 'appending ',boundarytime,lstime,lumitot
00230         summaryls[boundarytime].append([lstime,lumitot])
00231         previoustime=lstime
00232     #print summaryls
00233    
00234     summarylstimes=summaryls.keys()
00235     summarylstimes.sort()
00236     for bts in summarylstimes:
00237         startts=bts
00238         tsdatainseg=summaryls[bts]
00239         #print 'tsdatainseg ',tsdatainseg
00240         stopts=tsdatainseg[-1][0]
00241         plu=max(CommonUtil.transposed(tsdatainseg,0.0)[1])
00242         lui=sum(CommonUtil.transposed(tsdatainseg,0.0)[1])*23.357
00243         print >>f,'%d\t%d\t%e\t%e'%(startts,stopts,plu,lui)
00244     f.close()
00245         
00246 if __name__ == '__main__':
00247     parser = argparse.ArgumentParser(prog=os.path.basename(sys.argv[0]),description = "Dump Fill",formatter_class=argparse.ArgumentDefaultsHelpFormatter)
00248     # parse arguments
00249     parser.add_argument('-c',dest='connect',action='store',required=False,help='connect string to lumiDB,optional',default='frontier://LumiCalc/CMS_LUMI_PROD')
00250     parser.add_argument('-P',dest='authpath',action='store',help='path to authentication file,optional')
00251     parser.add_argument('-i',dest='inputdir',action='store',required=False,help='output dir',default='.')
00252     parser.add_argument('-o',dest='outputdir',action='store',required=False,help='output dir',default='.')
00253     parser.add_argument('-f',dest='fillnum',action='store',required=False,help='specific fill',default=None)
00254     parser.add_argument('-norm',dest='norm',action='store',required=False,help='norm',default=None)
00255     parser.add_argument('-siteconfpath',dest='siteconfpath',action='store',help='specific path to site-local-config.xml file, optional. If path undefined, fallback to cern proxy&server')
00256     parser.add_argument('--debug',dest='debug',action='store_true',help='debug')
00257     parser.add_argument('--toscreen',dest='toscreen',action='store_true',help='dump to screen')
00258     options=parser.parse_args()
00259     if options.authpath:
00260         os.environ['CORAL_AUTH_PATH'] = options.authpath
00261     parameters = lumiQueryAPI.ParametersObject()
00262     if options.norm!=None:
00263         parameters.normFactor=float(options.norm)
00264     session,svc =  lumiQueryAPI.setupSession (options.connect or \
00265                                               'frontier://LumiCalc/CMS_LUMI_PROD',
00266                                                options.siteconfpath,parameters,options.debug)
00267 
00268     ##
00269     #query DB for all fills and compare with allfills.txt
00270     #if found newer fills, store  in mem fill number
00271     #reprocess anyway the last 1 fill in the dir
00272     #redo specific lumi for all marked fills
00273     ##
00274  
00275     allfillsFromFile=[]
00276     fillstoprocess=[]
00277     session.transaction().start(True)
00278     if options.fillnum: #if process a specific single fill
00279         fillstoprocess.append(int(options.fillnum))
00280     else: #if process fills automatically
00281         q=session.nominalSchema().newQuery()    
00282         allfillsFromDB=lumiQueryAPI.allfills(q)
00283         del q
00284         if os.path.exists(os.path.join(options.inputdir,allfillname)):
00285             allfillF=open(os.path.join(options.inputdir,allfillname),'r')
00286             for line in allfillF:
00287                 l=line.strip()
00288                 allfillsFromFile.append(int(l))
00289             allfillF.close()
00290             if len(allfillsFromDB)==0:
00291                 print 'no fill found in DB, exit'
00292                 sys.exit(-1)
00293             if len(allfillsFromDB)!=0:
00294                 allfillsFromDB.sort()
00295             if len(allfillsFromFile) != 0:
00296                 allfillsFromFile.sort()
00297             #print 'allfillsFromFile ',allfillsFromFile
00298             if max(allfillsFromDB)>max(allfillsFromFile) : #need not to be one to one match because data can be deleted in DB
00299                 print 'found new fill '
00300                 for fill in allfillsFromDB:
00301                     if fill>max(allfillsFromFile):
00302                         fillstoprocess.append(fill)
00303             else:
00304                 print 'no new fill '
00305                 fillstoprocess+=allfillsFromFile[-1:]
00306             #if len(allfillsFromFile)>5: #reprocess anyway old fills
00307             #    fillstoprocess+=allfillsFromFile[-5:]
00308         else:
00309             fillstoprocess=allfillsFromDB #process everything from scratch
00310     #print 'fillstoprocess ',fillstoprocess
00311     runsperfillFromDB={}
00312     q=session.nominalSchema().newQuery()
00313     runsperfillFromDB=lumiQueryAPI.runsByfillrange(q,int(min(fillstoprocess)),int(max(fillstoprocess)))
00314     del q
00315     #print 'runsperfillFromDB ',runsperfillFromDB
00316     runtimes={}
00317     runs=runsperfillFromDB.values()#list of lists
00318     allruns=[item for sublist in runs for item in sublist]
00319     allruns.sort()
00320     #print 'allruns ',allruns
00321     for run in allruns:
00322         q=session.nominalSchema().newQuery()
00323         runtimes[run]=lumiQueryAPI.runsummaryByrun(q,run)[3]
00324         del q
00325     #write specificlumi to outputdir
00326     #update inputdir
00327     if len(fillstoprocess)!=0 and options.fillnum is None:
00328         filltofiles(allfillsFromDB,runsperfillFromDB,runtimes,options.inputdir)
00329     print '===== Start Processing Fills',fillstoprocess
00330     print '====='
00331     for fillnum in fillstoprocess:
00332         filldata=getSpecificLumi(session,parameters,fillnum,options.inputdir)
00333         specificlumiTofile(fillnum,filldata,options.outputdir)
00334     session.transaction().commit()
/data/refman/pasoursint/CMSSW_4_1_8_patch9/src/RecoLuminosity/LumiDB/scripts/specificLumi-2011.py
