Readflags.py

text/x-python Readflags.py — 10.4 KB

File contents

######################################################################
#
# readflags.py
#
#    readflags: reads Antenna.xml and Flag.xml SDM tables and parses
#               flag commands into returned dictionary
#
#    listflags: prints flags from dictionary returned by readflags
#
#    useflags:  applies flags to ms
#
#    plotflags: plots flags from dictionary returned by readflags
#
# Created S.T. Myers 2010-05-04  v1.0 from flagant.py
# Updated S.T. Myers 2010-05-24  v2.0 include reasons
# Updated S.T. Myers 2010-08-19  v3.0 restructure (read/list/useflags)
# Updated J. Marvil  2010-08-23  v4.0 plotflags
#           
# Usage:
#         from readflags import *
#         sdmfile = '62_1_sb166824_1_000.55257.37544340278'
#         myflags = readflags(sdmfile)
#    or
#         myflags = readflags(sdmfile,tbuff=2.0)
#
#    then
#         listflags(myflags)
#         listflags(myflags,antenna='ea01')
#         listflags(myflags,reason='FOCUS_ERROR')
#
#    and/or
#         msfile = '62_1_sb166824_1_000.55257.37544340278.ms'
#         useflags(msfile,myflags)
#         useflags(msfile,myflags,antenna='ea01')
#         useflags(msfile,myflags,reason='FOCUS_ERROR')
#
#    then
#         plotflags(myflags)
#
# Notes: readflags
#    Returns a dictionary with the 'antenna' and 'timerange' parameters
#    used from the SDM (or MS) Antenna.xml and Flag.xml tables.  Also
#    individual flags in 'flaglist' structure.
#
#    The tbuff parameter pads the timeranges (default by 1sec)
#
# (I)Python Notes: (courtesy A.Deller)
#
#    Nominally, the flagant.py needs to be in the same directory you
#    are running casapy from.  If you want to put the source in a
#    different directory, then:
#
#    1) Set up a place where you will put your modules and point a 
#       variable at it, e.g.:
#         export CASAPYUTILS="~/src/python_mods/"
#    2) Then make sure ipython knows about it by putting:
#         ip.ex("sys.path.append(os.environ['CASAPYUTILS'])")
#       in your ~/.ipython/ipy_user_conf.py (under the main: block).
#
######################################################################

def readflags(sdmfile, tbuff=1.0):
    try:
        import casac
    except ImportError, e:
        print "failed to load casa:\n", e
        exit(1)
    qatool = casac.homefinder.find_home_by_name('quantaHome')
    qa = casac.qa = qatool.create()
    try:
        from xml.dom import minidom
    except ImportError, e:
        print "failed to load xml.dom.minidom:\n", e
        exit(1)
    try:
        from flagdata_cli import  flagdata_cli as flagdata
    except ImportError, e:
        print "failed to load flagdata:\n", e
        exit(1)
    # construct look-up dictionary of name vs. id from Antenna.xml
    xmlants = minidom.parse(sdmfile+'/Antenna.xml')
    antdict = {}
    rowlist = xmlants.getElementsByTagName("row")
    for rownode in rowlist:
        rowname = rownode.getElementsByTagName("name")
        ant = str(rowname[0].childNodes[0].nodeValue)
        rowid = rownode.getElementsByTagName("antennaId")
        id = str(rowid[0].childNodes[0].nodeValue)
        antdict[id] = ant
    print '  Found ',rowlist.length,' antennas in Antenna.xml'

    # now read Flag.xml into dictionary also and make a list
    xmlflags = minidom.parse(sdmfile+'/Flag.xml')
    flagdict = {}
    flaglist = []
    flagants = {}
    flagreas = {}
    rowlist = xmlflags.getElementsByTagName("row")
    for rownode in rowlist:
        rowfid = rownode.getElementsByTagName("flagId")
        fid = str(rowfid[0].childNodes[0].nodeValue)
        flagdict[fid] = {}
        flaglist.append(fid)
        rowid = rownode.getElementsByTagName("antennaId")
        id = str(rowid[0].childNodes[0].nodeValue)
        # start and end times in mjd ns
        rowstart = rownode.getElementsByTagName("startTime")
        start = int(rowstart[0].childNodes[0].nodeValue)
        startmjd = (float(start)*1.0E-9 - tbuff)/86400.0
        t = qa.quantity(startmjd,'d')
        starttime = qa.time(t,form="ymd")
        rowend = rownode.getElementsByTagName("endTime")
        end = int(rowend[0].childNodes[0].nodeValue)
        endmjd = (float(end)*1.0E-9 + tbuff)/86400.0
        t = qa.quantity(endmjd,'d')
        endtime = qa.time(t,form="ymd")
        # reasons
        rowreason = rownode.getElementsByTagName("reason")
        reas = str(rowreason[0].childNodes[0].nodeValue)
        # Construct antenna name and timerange and reason strings
        antname = antdict[id]
        flagdict[fid]['antenna'] = antname
        cmd = starttime+'~'+endtime
        flagdict[fid]['timerange'] = cmd
        flagdict[fid]['reason'] = reas
        # Also build per antenna strings
        if flagants.has_key(antname):
            flagants[antname] += ','+cmd
            flagreas[antname] += ','+reas
        else:
            flagants[antname] = cmd
            flagreas[antname] = reas
    print '  Found ',rowlist.length,' flags in Flag.xml'

    flags = {}
    flags['antenna'] = []
    flags['timerange'] = []
    flags['reason'] = []
    flags['flaglist'] = flagdict
    if rowlist.length > 0:
        #
        keylist = flagants.keys()
        keylist.sort()
        antenna = []
        timerange = []
        reason = []
        for ant in keylist:
            antenna.append(ant)
            timerange.append(flagants[ant])
            reason.append(flagreas[ant])
        # save the params in a dictionary
        flags['antenna'] = antenna
        flags['timerange'] = timerange
        flags['reason'] = reason
    else:
        print 'No flagging found in Flag.xml'

    # return the dictionary for later use
    return flags
# Done

def listflags(myflags,antenna='',reason=''):
    # Loop over flags
    for key in myflags['flaglist'].keys():
        myf = myflags['flaglist'][key]
        ant = myf['antenna']
        tim = myf['timerange']
        reas = myf['reason']
        if antenna=='' or ant==antenna:
            if reason=='' or reas==reason:
                print '%16s %8s %32s %48s' % (key, ant, reas, tim)
    return
# Done

def useflags(msfile,myflags,antenna='',reason=''):
    try:
        import casac
    except ImportError, e:
        print "failed to load casa:\n", e
        exit(1)
    try:
        from flagdata_cli import  flagdata_cli as flagdata
    except ImportError, e:
        print "failed to load flagdata:\n", e
        exit(1)

    keylist = myflags['flaglist'].keys()
    if keylist.__len__() > 0:
        #
        # Flag the data
        #
        keylist.sort()
        flagants = {}
        #
        # Construct flag list for selected ant,reason
        #
        for key in keylist:
            ant = myflags['flaglist'][key]['antenna']
            if antenna=='' or ant==antenna:
                reas = myflags['flaglist'][key]['reason']
                if reason=='' or reas==reason:
                    # print '%16s %8s %32s %48s' % (key, ant, reas, tim)
                    cmd = myflags['flaglist'][key]['timerange']
                    if flagants.has_key(ant):
                        flagants[ant] += ','+cmd
                    else:
                        flagants[ant] = cmd
        #
        # Form and apply flags if any
        flagkeys = flagants.keys()
        if flagkeys.__len__() > 0:
            antlist = []
            timelist = []
            flagkeys.sort()
            for key in flagkeys:
                cmd = flagants[key]
                antlist.append(key)
                timelist.append(cmd)
            print 'Flagging data with '+str(antlist.__len__())+' entries'
            flagdata(vis=msfile,mode='manualflag',antenna=antlist,timerange=timelist)
        else:
            print 'No flags found meeting criteria'
    else:
       print 'No flags found in dictionary'

    return
# Done

def plotflags(myflags):
    try: 
        import casac 
    except ImportError, e: 
        print "failed to load casa:\n", e 
        exit(1) 
    qatool = casac.homefinder.find_home_by_name('quantaHome') 
    qa = casac.qa = qatool.create()
    
    try:
        import pylab as pl
    except ImportError, e:
        print "failed to load pylab:\n", e
        exit(1)

    myants=myflags['antenna']

    antind=0
    pl.ioff()
    f1=pl.figure()
    ax1=f1.add_axes([.15,.1,.75,.85])
#    ax1.set_ylabel('antenna')
#    ax1.set_xlabel('time')
    badflags=[]
    for thisant in myflags['antenna']:
        antind+=1
	for thisflag in myflags['flaglist']:
	    if myflags['flaglist'][thisflag]['antenna']==thisant:
#	        print thisant, myflags['flaglist'][thisflag]['reason'], myflags['flaglist'][thisflag]['timerange']
                thisReason=myflags['flaglist'][thisflag]['reason']
		if thisReason=='FOCUS_ERROR': thisColor='red'
		elif thisReason=='SUBREFLECTOR_ERROR': thisColor='blue'
		elif thisReason=='ANTENNA_NOT_ON_SOURCE': thisColor='green'
		elif thisReason=='ANTENNA_NOT_IN_SUBARRAY': thisColor='black'
		mytimerange=myflags['flaglist'][thisflag]['timerange']
		t1=mytimerange[:mytimerange.find('~')]
		t2=mytimerange[mytimerange.find('~')+1:]
		t1s, t2s=qa.convert(t1,'s')['value'], qa.convert(t2,'s')['value']
                myTimeSpan=t2s-t1s
		if myTimeSpan < 10000: ax1.plot([t1s,t2s],[antind,antind], color=thisColor, lw=2, alpha=.4)
		else: badflags.append((thisant, myTimeSpan, thisReason))

    ##badflags are ones which span a time longer than that used above
    ##they can be so long that including them compresses the time axis so that none of the other flags are visible    
#    print 'badflags', badflags
    myXlim=ax1.get_xlim()
    myXrange=myXlim[1]-myXlim[0]
    ax1.text(myXlim[0]+.05*myXrange, 29, 'FOCUS_ERROR', color='red')
    ax1.text(myXlim[0]+.25*myXrange, 29, 'SUBREFLECTOR', color='blue')
    ax1.text(myXlim[0]+.45*myXrange, 29, 'OFF_SOURCE', color='green')
    ax1.text(myXlim[0]+.65*myXrange, 29, 'NOT_IN_SUBARRAY', color='black')
    ax1.set_ylim([0,30])   

    ax1.set_yticks(range(1,len(myants)+1))
    ax1.set_yticklabels(myants)
    ax1.set_xticks(pl.linspace(myXlim[0],myXlim[1],3))

    mytime=[myXlim[0],(myXlim[1]-myXlim[0])/2.0,myXlim[1]]
    myTimestr = []
    for time in mytime:
        q1=qa.quantity(time,'s')
        time1=qa.time(q1,form='ymd')
        myTimestr.append(time1)

    ax1.set_xticklabels([myTimestr[0],myTimestr[1][11:], myTimestr[2][11:]])
    f1.show()    

    return
# Done

Connect with NRAO

The NSF National Radio Astronomy Observatory and NSF Green Bank Observatory are facilities of the U.S. National Science Foundation operated under cooperative agreement by Associated Universities, Inc.

Founded in 1956, the NSF NRAO provides state-of-the-art radio telescope facilities for use by the international scientific community. NSF NRAO telescopes are open to all astronomers regardless of institutional or national affiliation. Observing time on NSF NRAO telescopes is available on a competitive basis to qualified scientists after evaluation of research proposals on the basis of scientific merit, the capability of the instruments to do the work, and the availability of the telescope during the requested time. NSF NRAO also provides both formal and informal programs in education and public outreach for teachers, students, the general public, and the media.

The U.S. National Science Foundation (U.S. NSF) is an independent federal agency created by Congress in 1950 "to promote the progress of science; to advance the national health, prosperity, and welfare; to secure the national defense..."

AUI collaborates with the scientific community and research sponsors to plan, build, and operate cutting-edge facilities. We cultivate excellence, deliver value, enhance education, and engage the public.