casa6

###############################################

## To plot stuff in weather tables, saved to MSname+.plotWX.png

## and estimate zenith opacity per spw, returned as a list named myTau

##

##

## J. Marvil 2.6.12

## revised 4.27.12 to add support for missing/empty weather table

## revised 11.05.12 to address CASA 4.0 changes

###############################################

​

​

import casac

from taskinit import *

import pylab as pl

from math import pi,floor

import os.path as osp

​

​

###############

## hides the extreme Y-axis ticks, helps stack plots close together without labels overlaping

def jm_clip_Yticks():

    xa=pl.gca()

    nlabels=0

    for label in xa.yaxis.get_ticklabels():

        nlabels+=1

    thislabel=0

    if nlabels>3:

        for label in xa.yaxis.get_ticklabels():

            if thislabel==0: label.set_alpha(0)

            if thislabel==nlabels-1: label.set_alpha(0)

            thislabel+=1

​

##############

##  sets the position of the y-axis label to the right side of the plot, can also move up/down

def jm_set_Ylabel_pos(pos=(0.5,0.5)):

    ax=pl.gca();

    ax.yaxis.set_label_position('right')

    ax.yaxis.label.set_rotation(270)

    ax.yaxis.label.set_position(pos)

​

​

###############

## fixed y-ticks, from myMin to myMax

def jm_set_Ylim_ticks(myMin=-1,myMax=1):

    myYlocs=pl.linspace(round(myMin,1),round(myMax,1),5)

    myLocator = pl.FixedLocator(myYlocs)

    ax=pl.gca()

    ax.yaxis.set_major_locator( myLocator )

    pl.ylim(myMin,myMax)

    jm_clip_Yticks()

​

###############

## variable y-ticks, but not more than 1+ this argument

def jm_set_Yvar_ticks(myScale=4):

    xa=pl.gca()

    xa.yaxis.set_major_locator(pl.MaxNLocator(myScale))

    jm_clip_Yticks()

​

###############

## calculates K-band zenith opacity from temperature and dewpoint

def Tau_K_Calc(D,T,day, weights=(.5,.5)):

    P = pl.exp(1.81+(17.27*D)/(D+237.3)) # water vapor partial pressure

    h = 324.7*P/(T+273.15) # PWV in mm

    tau_w = 3.8 + 0.23*h + 0.065*h**2 # tau from weather, in %, at 22GHz

    if day > 199: day = day - 365.

    m = day + 165. # modified day of the year

    tau_d = 22.1 - 0.178*m + 0.00044*m**2 # tau from seaonal model, in %

    tau_k = weights[0]*tau_w + weights[1]*tau_d # the average, with equal weights (as in the AIPS default)

    return tau_k, h

​

################

## calculates elevation of the sun

def jm_sunEL(mytime):

    me.doframe(me.observatory('VLA'))

    me.doframe(me.epoch('utc',mytime))

    mysun=me.measure(me.direction('SUN'),'AZELGEO')

    return mysun['m1']['value']

​

################

## gets and plots data from the weather table of the given MS

def plotweather(vis='', seasonal_weight=0.5, doPlot=True, plotName = ''):

    myMS=vis

    if plotName == '':

        if myMS.endswith("/"):

            plotName = myMS + myMS.rstrip("/") + '.plotweather.png'

        else:

            plotName = myMS + '.plotweather.png'

​

    # check for weather table

​

    if osp.isdir(myMS+'/WEATHER'):