def compareSubTables(input,reference,order=None,excluded_cols=[]):
    
    tbinput = tbtool()
    tbinput.open(input)
    if order is not None:
        tbinput_sorted = tbinput.taql("SELECT * from " + input + " order by " + order)
    else:
        tbinput_sorted = tbinput
    
    tbreference = tbtool()
    tbreference.open(reference)
    if order is not None:
        tbreference_sorted = tbreference.taql("SELECT * from " + reference + " order by " + order)
    else:
        tbreference_sorted = tbreference
    
    columns = tbinput.colnames()
    for col in columns:
        if not col in excluded_cols:
            col_input = tbinput_sorted.getcol(col)
            col_reference = tbreference_sorted.getcol(col)
            if not (col_input == col_reference).all():
                tbinput.close()
                tbreference.close()
                del tbinput
                del tbreference
                return (False,col)
    
    tbinput.close()
    tbreference.close()
    del tbinput
    del tbreference
    
    return (True,"OK")

    def test_vla_mixed_polarizations_mms2(self):
        
        self.outputms = 'test_vla_mixed_polarizations_2.mms'
        
        mstransform(vis=self.vis,outputvis=self.outputms,scan='16',datacolumn='DATA', createmms=True,
                    separationaxis='spw',spw='16~18',correlation='XX')
        
        # Check that DDI sub-table is consistent with POLARIZATION sub-table
        mytb = tbtool()
        mytb.open(self.outputms + '/POLARIZATION')
        npols = mytb.nrows()
        mytb.close()
        
        mytb = tbtool()
        mytb.open(self.outputms + '/DATA_DESCRIPTION')
        polIds = mytb.getcol('POLARIZATION_ID')
        mytb.close()    
        for id in polIds:
            self.assertTrue(id in range(npols),'PolarizationId in DATA_DESCRIPTION not consistent with POLARIZATION table')
        
#        self.assertTrue(all(polIds < npols), 'PolarizationId in DATA_DESCRIPTION not consistent with POLARIZATION table') 
        
        # Check that flagdata can run properly with output MS
        summary = flagdata(vis=self.outputms,mode='summary')
        self.assertTrue(summary.has_key('correlation'), 'Flagdata failure due to missformated MS') 

def fixfeedpa(vis,defband='',forceband=''):

    mytb=taskinit.tbtool()

    mytb.open(vis+'/SPECTRAL_WINDOW')
    spwnames=mytb.getcol('NAME')
    mytb.close()
    if len(forceband)>0:
        print 'Forcing band = ',forceband
        spwnames[:]=forceband
        defband=forceband
    mytb.open(vis+'/FEED',nomodify=False)
    spwids=mytb.getcol('SPECTRAL_WINDOW_ID')
    ra=mytb.getcol('RECEPTOR_ANGLE')
    ra[:,:]=0.0
    spwmask=(spwids>-1)
    ra[0,spwmask]=[bandpa(spwnames[ispw]) for ispw in spwids[spwmask]]
    spwmask=pl.logical_not(spwmask)
    if (sum(spwmask)>0):
        if (len(defband)>0):
            print 'NB: Setting spwid=-1 rows in FEED table to RECEPTOR_ANGLE for band='+str(defband)
            ra[0,spwmask]=bandpa(defband)
        else:
            print 'NB: Setting spwid=-1 rows in FEED table to RECEPTOR_ANGLE=(0,pi/2)'
    ra[1,:]=ra[0,:]+(pi/2.)
    mytb.putcol('RECEPTOR_ANGLE',ra)
    mytb.close()

def abschanwidth(vis="", spw=""):

    """
    Usage: abschanwidth(vis, spw)
           Get the absolute channel width for the given spw.
           Returns 0 upon error.
    """

    if vis == "" or spw == "":
        print "Usage: abschanwidth(vis, spw)"
        return 0

    myvis = vis
    myspw = spw
    mytb = taskinit.tbtool()

    mytb.open(myvis + "/SPECTRAL_WINDOW")
    if spw >= mytb.nrows() or spw < 0:
        print "Error: spw out of range. Min is 0. Max is ", mytb.nrows() - 1
        return 0

    mychw = mytb.getcell("CHAN_WIDTH", spw)[0]
    mytb.close()

    return numpy.fabs(mychw)

def scanbystate(vis,undo=False):

    mytb=taskinit.tbtool()

    mytb.open(vis,nomodify=False)
    scans=mytb.getcol('SCAN_NUMBER')
    states=mytb.getcol('STATE_ID')
    print 'Unique STATE_IDs = ',str(pl.unique(states))
    maxstate=states.max()

    if undo:
        d=10**int(floor(log10(scans.min())))
        if d<10:
            mytb.close()
            raise Exception, 'Apparently, nothing to undo'
        scans-=states
        scans/=d
        print 'New SCAN_NUMBER = (SCAN_NUMBER - STATE_ID) / '+str(d)
    else:
        m=10**int(floor(log10(states.max())+1.0))
        scans*=m
        scans+=states
        print 'New SCAN_NUMBER = SCAN_NUMBER * '+str(m)+' + STATE_ID'

    mytb.putcol('SCAN_NUMBER',scans)
    mytb.close()

def getnch(vis="", spw=""):

    """
    Usage: getnch(vis, spw)
           Get the nchan for the given spw.
           Returns 0 upon error.
    """

    if vis == "" or spw == "":
        print "Usage: abschanwidth(vis, spw)"
        return 0

    myvis = vis
    myspw = spw
    mytb = taskinit.tbtool()

    mytb.open(myvis + "/SPECTRAL_WINDOW")
    if spw >= mytb.nrows() or spw < 0:
        print "Error: spw out of range. Min is 0. Max is ", mytb.nrows() - 1
        return 0

    mynch = mytb.getcell("NUM_CHAN", spw)
    mytb.close()

    return numpy.abs(mynch)

def dxy(dtab,xytab,dout):

    mytb=taskinit.tbtool()

    os.system('cp -r '+dtab+' '+dout)

    # How many spws
    mytb.open(dtab+'/SPECTRAL_WINDOW')
    nspw=mytb.nrows()
    mytb.close()


    for ispw in range(nspw):
        mytb.open(xytab)
        st=mytb.query('SPECTRAL_WINDOW_ID=='+str(ispw))
        x=st.getcol('CPARAM')
        st.close()
        mytb.close()

        mytb.open(dout,nomodify=False)
        st=mytb.query('SPECTRAL_WINDOW_ID=='+str(ispw))
        d=st.getcol('CPARAM')

        # the following assumes all antennas and chans same in both tables.

        # Xinv.D.X:
        d[0,:,:]*=pl.conj(x[0,:,:])
        d[1,:,:]*=x[0,:,:]

        st.putcol('CPARAM',d)
        st.close()
        mytb.close()

def getColShape(table,col,start_row=0,nrow=1,row_inc=1):
    """ Get the shape of the given column.
    Keyword arguments:
        table      --    input table or MS
        col        --    column to get the shape
        start_row  --    start row (default 0)
        nrow       --    number of rows to read (default 1)
        row_inc    --    increment of rows to read (default 1)
        
        Return a list of strings with the shape of each row in the column.
    
    """

    col_shape = []
    try:
        try:
            tblocal = tbtool()
            tblocal.open(table)
            col_shape = tblocal.getcolshapestring(col,start_row,nrow,row_inc)
        except:
            print 'Cannot get shape of col %s from table %s '%(col,table)

    finally:
        tblocal.close()
            
    return col_shape

    def test_unapply_clip_and_unset_flagrow(self):
        '''flagcmd: Check that FLAG_ROW is unset after un-applying a clip agent'''
        # Remove any cmd from table
        flagcmd(vis=self.vis, action='clear', clearall=True)

        # Flag using manual agent
        myinput = "scan='4'"
        filename = create_input(myinput)
        flagcmd(vis=self.vis, inpmode='list', inpfile=filename, action='apply', savepars=False)
        
        # Check FLAG_ROW is all set to true
        mytb = tbtool()
        mytb.open(self.vis)
        selectedtb = mytb.query('SCAN_NUMBER in [4]')
        FLAG_ROW = selectedtb.getcol('FLAG_ROW')
        mytb.close()        
        selectedtb.close()
        self.assertEqual(FLAG_ROW.sum(), FLAG_ROW.size)
        
        # Flag using tfcrop agent from file
        myinput = "scan='4' mode=clip "
        filename = create_input(myinput)
        flagcmd(vis=self.vis, inpmode='list', inpfile=filename, action='apply', savepars=True,
                flagbackup=False)
        
        # Check FLAG_ROW is all set to true
        mytb = tbtool()
        mytb.open(self.vis)
        selectedtb = mytb.query('SCAN_NUMBER in [4]')
        FLAG_ROW = selectedtb.getcol('FLAG_ROW')
        mytb.close()           
        selectedtb.close()
        self.assertEqual(FLAG_ROW.sum(), FLAG_ROW.size)
        
        # Unapply only the tfcrop line
        flagcmd(vis=self.vis, action='unapply', useapplied=True, tablerows=0, savepars=False)
       
        # Check FLAG_ROW is now all set to false
        mytb = tbtool()
        mytb.open(self.vis)
        selectedtb = mytb.query('SCAN_NUMBER in [4]')
        FLAG_ROW = selectedtb.getcol('FLAG_ROW')
        mytb.close()        
        selectedtb.close()
        self.assertEqual(FLAG_ROW.sum(), 0)

def xyamb(xytab,qu,xyout=''):

    mytb=taskinit.tbtool()

    if not isinstance(qu,tuple):
        raise Exception,'qu must be a tuple: (Q,U)'

    if xyout=='':
        xyout=xytab
    if xyout!=xytab:
        os.system('cp -r '+xytab+' '+xyout)

    QUexp=complex(qu[0],qu[1])
    print 'Expected QU = ',qu   # , '  (',pl.angle(QUexp)*180/pi,')'

    mytb.open(xyout,nomodify=False)

    QU=mytb.getkeyword('QU')['QU']
    P=pl.sqrt(QU[0,:]**2+QU[1,:]**2)

    nspw=P.shape[0]
    for ispw in range(nspw):
        st=mytb.query('SPECTRAL_WINDOW_ID=='+str(ispw))
        if (st.nrows()>0):
            q=QU[0,ispw]
            u=QU[1,ispw]
            qufound=complex(q,u)
            c=st.getcol('CPARAM')
            fl=st.getcol('FLAG')
            xyph0=pl.angle(pl.mean(c[0,:,:][pl.logical_not(fl[0,:,:])]),True)
            print 'Spw = '+str(ispw)+': Found QU = '+str(QU[:,ispw])  # +'   ('+str(pl.angle(qufound)*180/pi)+')'
            #if ( (abs(q)>0.0 and abs(qu[0])>0.0 and (q/qu[0])<0.0) or
            #     (abs(u)>0.0 and abs(qu[1])>0.0 and (u/qu[1])<0.0) ):
            if ( pl.absolute(pl.angle(qufound/QUexp)*180/pi)>90.0 ):
                c[0,:,:]*=-1.0
                xyph1=pl.angle(pl.mean(c[0,:,:][pl.logical_not(fl[0,:,:])]),True)
                st.putcol('CPARAM',c)
                QU[:,ispw]*=-1
                print '   ...CONVERTING X-Y phase from '+str(xyph0)+' to '+str(xyph1)+' deg'
            else:
                print '      ...KEEPING X-Y phase '+str(xyph0)+' deg'
            st.close()
    QUr={}
    QUr['QU']=QU
    mytb.putkeyword('QU',QUr)
    mytb.close()
    QUm=pl.mean(QU[:,P>0],1)
    QUe=pl.std(QU[:,P>0],1)
    Pm=pl.sqrt(QUm[0]**2+QUm[1]**2)
    Xm=0.5*atan2(QUm[1],QUm[0])*180/pi

    print 'Ambiguity resolved (spw mean): Q=',QUm[0],'U=',QUm[1],'(rms=',QUe[0],QUe[1],')','P=',Pm,'X=',Xm

    stokes=[1.0,QUm[0],QUm[1],0.0]
    print 'Returning the following Stokes vector: '+str(stokes)
    
    return stokes

def clearcal(
    vis=None,
    field=None,
    spw=None,
    intent=None,
    addmodel=None,
    ):

    casalog.origin('clearcal')

    # Do the trivial parallelization
    if ParallelTaskHelper.isParallelMS(vis):
        helper = ParallelTaskHelper('clearcal', locals())
        helper.go()
        return

    # Local versions of the tools
    tblocal = tbtool()
    cblocal = cbtool()
    mslocal = mstool()

    try:

        # we will initialize scr cols only if we don't create them
        doinit = False

        if (type(vis) == str) & os.path.exists(vis):
            tblocal.open(vis)
            doinit = tblocal.colnames().count('CORRECTED_DATA') > 0
            tblocal.close()

            # We ignore selection if creating the scratch columns
            if not doinit:
                casalog.post('Need to create scratch columns; ignoring selection.'
                             )

            cblocal.open(vis, addmodel=addmodel)
        else:
            raise Exception, \
                'Visibility data set not found - please verify the name'

        # If necessary (scr col not just created), initialize scr cols
        if doinit:
            cblocal.selectvis(field=field, spw=spw, intent=intent)
            cblocal.initcalset(1)
        cblocal.close()

        # Write history to the MS
        param_names = clearcal.func_code.co_varnames[:clearcal.func_code.co_argcount]
        param_vals = [eval(p) for p in param_names]
        casalog.post('Updating the history in the output', 'DEBUG1')
        write_history(mslocal, vis, 'clearcal', param_names,
                      param_vals, casalog)
        
    except Exception, instance:

        print '*** Error ***', instance

def copy_model_RRtoLL(vis):
    # copy RR column of model_data column to LL
    tb = tbtool()
    tb.open(vis,nomodify=False)
    model_vis = tb.getcol('MODEL_DATA')
    model_vis[3,:,:] = model_vis[0,:,:] # copy RR model column to LL model column
    tb.putcol('MODEL_DATA',model_vis)
    tb.unlock()
    tb.close()

def getavweight(vis="", field=[], spw=""):

    """
    Usage: getavweight(vis, field, spw)
           Get the average weight for the given field and spw.
           The field parameter takes a list of fields.
    """

    if vis == "" or spw == "" or field == [] or not type(field) == list:
        print "Usage: getavweight(vis, field, spw)"
        print "       The field parameter takes a list of fields."
        return False

    myvis = vis
    myspw = spw
    myfields = field
    mytb = taskinit.tbtool()

    mytb.open(myvis)
    w = mytb.getcol("WEIGHT")
    dd = mytb.getcol("DATA_DESC_ID")
    ff = mytb.getcol("FIELD_ID")
    mytb.close()

    mytb.open(myvis + "/DATA_DESCRIPTION")
    mydds = []

    for i in range(0, mytb.nrows()):
        if mytb.getcell("SPECTRAL_WINDOW_ID", i) != myspw:
            continue
        else:
            mydds.append(i)

    mytb.close()

    mynrows = 0
    mysumw = 0

    npol = len(w)

    if len(mydds) > 0:
        for row in range(0, len(dd)):
            if (dd[row] in mydds) and (ff[row] in myfields):
                mynrows += 1
                for i in range(0, npol):
                    mysumw += w[i][row]

    rval = 0.0

    if mynrows > 0:
        rval = mysumw / float(npol) / float(mynrows)
        print "Average weight is ", rval
    else:
        print "No rows selected."

    return rval

def copy_data_RRtoLL(vis):
    # copy RR column of data column to LL
    # doesn't account for different flags on RR and LL
    tb = tbtool()
    tb.open(vis,nomodify=False)
    data_vis = tb.getcol('DATA')
    data_vis[3,:,:] = data_vis[0,:,:] # copy RR model column to LL model column
    tb.putcol('DATA',data_vis)
    tb.unlock()
    tb.close()

 def test_model_keys(self):
     '''partition: CAS-4398, handle the MODEL keywords correctly'''
     
     print '*** Check that MODEL_DATA is not present in MS first'
     mytb = tbtool()
     try:
         mytb.open(self.msfile+'/MODEL_DATA')
     except Exception, instance:
         print '*** Expected exception. \"%s\"'%instance
         mytb.close()

def adjustweights2(vis="", field="", spws=[]):

    """
       Usage: adjustweights2(vis, field, spws)
              Scale the weights in specified spws by a factor
              2*df*dt/nchan, where df is the channel bandwidth,
              dt is the integration time, and nchan is the number
              of channels in the spw.  This enables imaging of
              mixed mode spws in CASA v4.2 and earlier.  (Note
              that there will be no net effect if the spws
              share the same df, dt, and nchan.

              spws are of type list,
              field should be given as field id.
    """

    myvis = vis
    myfield = int(field)
    mytb = taskinit.tbtool()

    if vis == "" or myfield == "" or spws == [] or not type(spws) == list:
        print "Usage: adjustweights2(vis, field, spws)"
        print "       spws are of type list,"
        print "       field should be given as field id"
        return False

    # get avweight and chanwidth from spws
    chw = []
    nch = []
    for spw in spws:
        cw = abschanwidth(myvis, spw)
        if cw == 0:
            print "Error reading channel width of spw ", spw
            return False
        chw.append(cw)

        nc = getnch(myvis, spw)
        if nc == 0:
            print "Error: nch of spw ", spw, " is zero (could also mean no data)."
            return False
        nch.append(nc)

        print "Spw ", spw, ", channelwidth", cw, ", nchan ", nc

    # calculate scale factor and apply scaling to the spws

    for i in range(0, len(spws)):
        myscale = chw[i] / nch[i]
        print "Scale factor for weights in spw ", spws[i], " is ", myscale
        scaleweights(myvis, [myfield], spws[i], myscale, True)  # include integ time

    print "Done."

    return True

def scan_reindex(vis,gap=50.):
   # replace scan numbers to count up from 1, changing every time there is a gap in time
   #   greater than gap (in seconds) between entries
   # assumes ms is already chronological; does not account for multiple fields
   # (but should be fine w/ multiple fields if gap < slew time)
   t = tbtool()
   t.open(vis, nomodify=False)
   scanlist0 = t.getcol('SCAN_NUMBER')
   times = t.getcol('TIME')
   scanlist = numpy.ones(scanlist0.shape)
   dt = times[1:]-times[:-1]   # should be in seconds
   scan_inc = numpy.cumsum(dt>gap)
   scanlist[1:] += scan_inc
   scanlist = scanlist.astype(int)
   t.putcol('SCAN_NUMBER',scanlist)
   t.unlock()
   t.close()

 def test_alma_wvr_correlation_products_mms1(self):
     
     self.outputms = 'test_alma_wvr_correlation_products_1.mms'
     # Only spw=2 exist in MS
     mstransform(vis=self.vis,outputvis=self.outputms,spw='0,1,2',datacolumn='DATA',createmms=True)
     
     # Check that POLARIZATION sub-table is properly sorted
     mytb = tbtool()
     mytb.open(self.outputms + '/POLARIZATION')
     numCorr = mytb.getcol('NUM_CORR')
     mytb.close()    
     
     self.assertEqual(numCorr[0],2,'POLARIZATION table miss-sorted')         
     self.assertEqual(numCorr[1],1, 'POLARIZATION table miss-sorted')         
     
     # Check that flagdata can run properly with output MS
     summary = flagdata(vis=self.outputms,mode='summary')
     self.assertTrue(summary.has_key('correlation'), 'Flagdata failure due to missformated MS')   

    def test_shape3(self):
        '''mstransform: DATA and FLAG tileshapes should be the same'''
        self.outputms = "shape3.ms"
        inptsh = [4,10,1024]
        mstransform(vis=self.vis, outputvis=self.outputms, createmms=True, tileshape=inptsh)

        self.assertTrue(os.path.exists(self.outputms))

        # Get the tile shape for the DATA output
        tblocal = tbtool()
        tblocal.open(self.outputms)
        outdm = tblocal.getdminfo()
        tblocal.close()
        outtsh = th.getTileShape(outdm)
        # And for the FLAG column
        flagtsh = th.getTileShape(outdm, 'FLAG')

        self.assertTrue((outtsh==flagtsh).all(), 'Tile shapes are different')

 def test_CAS6206(self):
     '''mstransform: verify that all columns are re-indexed in SPW sub-table'''
     self.outputmms='test.mms'
     self.outputms='assoc.ms'
     self.setUp_CAS_5013()
     mstransform(vis=self.vis, outputvis=self.outputmms,createmms=True, datacolumn='corrected')
     
     # Check that optional ASSOC_SPW_ID is the same in input and output
     tblocal = tbtool()
     tblocal.open(self.vis+'/SPECTRAL_WINDOW',nomodify=True)
     in_assoc = tblocal.iscelldefined('ASSOC_SPW_ID',0)
     tblocal.close()
     tblocal.open(self.outputmms+'/SPECTRAL_WINDOW',nomodify=True)
     out_assoc = tblocal.iscelldefined('ASSOC_SPW_ID',0)
     tblocal.close()
     self.assertEqual(in_assoc, out_assoc, 'Error in SPW sub-table creation; ASSOC_SPW_ID is different')
     
     # if SPW sub-table is not correct, the next step might fail
     self.assertTrue(mstransform(vis=self.outputmms, outputvis=self.outputms, hanning=True, datacolumn='data'))