Getting ATCA data into CASA

This recipe describes a method to convert ATCA spectral line data to a MS so that the frequencies are scaled correctly to velocity. The data are taken at fixed frequency but the centre of the band may be at a frequency shifted from the line rest frequency to allow for the source velocity.

The example used for illustration has:

Pre-processing in MIRIAD

Raw data in RPFITS format needs to be passed through MIRIAD: 

> atlod in="Data/*" restfequency=6.6685142 options=birdie,noif out=atca.uv ifsel=1
The resulting MIRIAD .uv file has the following header: 
Prthd: version 03-Jan-05
****************************************************************
Filename: ../testingcasa2/atca.uv
Instrument: ATCA            Telescope: ATCA
Object: 299p782             Observer: mmb_group
First time: 08JAN22:12:29:50.0
Number of antennae: 6
Polarisations Present: XX,YY
Type of correlations present: crosscorrelation
----------------------------------------------------------------
Spectral Correlations:
  Spectrum  Channels  Freq(chan=1)  Increment  Restfreq
      1       1025       6.66600     0.000004   6.66851 GHz
  Total number of correlations: 66331850
  Correlations are stored in 32-bit form
----------------------------------------------------------------
This is then converted into UVFITS using the MIRIAD fits command 
 >fits in=atca.uv op=uvout out=atca.fits line=channel,1024
This .fits file has the following frequency information in its header (you can simply use unix more to look at this, it helps to set the xterm width so that the line breaks come in the right place). The lines below have many other lines in between. 
NAXIS4  =                 1024  /   # Number of channels
CRVAL4  =    6.66599981147E+09  /   # Frequency in the reference channel
CRPIX4  =    1.00000000000E+00  /   # Number of reference channel (here, 1-indexed)
CDELT4  =    3.90624988952E+03  /   # Channel width
If you find it easier you can calculate the channel width from the bandwidth and the number of channels or v.v.

Importing and inspecting data into CASA

Use importuvfits in CASA to import the data by conversion to a Measurement Set and inspect it using listobs. The MS is called atca.ms. 

listobs    Spectral Windows:  (1 unique spectral windows and 1 unique polarization setups)
listobs	    SpwID  #Chans Frame Ch1(MHz)    ChanWid(kHz) TotBW(kHz)  Ref(MHz)    Corrs   
listobs      0       1024 LSRK  6665.77661  3.90624989    3999.99989 6665.77661  XX  YY
with a source rest frequency of: 
listobs      20   13p667       0     6668.51902252  0
One or more of these may be incorrect by comparison with the input FITS header values and known rest frequency. This is related to a problem identified by Crystal Brogan in CAS-1992, which is that the velocity reference frame should be TOPO but the FITS code output by MIRIAD is interpreted by CASA as LSRK. Check all the values; the recipes below explain how to fix all the parameters we have found to be wrong but you may not need to change all of them.

Correcting the frequencies

Set up variables for this data set. Note that all frequencies are in Hz. Note that the input format determines the type, e.g. 3 is an integer, 3. is a float.
rest_freq is the line rest frequency.
ref_freq is CRVAL4 from the FITS header, the expected frequency in the reference channel at the source velocity used.
chan_width is CDELT4, the channel width.
num_chan is the number of channels, from listobs
bw is the total bandwidth 

myvis = "atca.ms"
rest_freq = 6668514200.0
ref_freq   = 6.66599981e+09
chan_width = 3.90624988952E+03
num_source = 3
num_chan = 1024
bw = 4.0E+06
Set up an array of rest frequencies, one for each source (assuming that they are all the same). 
freq_arr=[]
for x in range(num_source):             
    freq_arr.append(rest_freq)
Define the frequency in the first channel, in this case it is also the reference channel: 
chan1_freq = ref_freq
Create arrays of channel frequencies and channel widths: 
chan_freq_arr = [chan1_freq]
chan_width_arr = [chan_width]
chan_freq = chan1_freq
for y in range(num_chan-1):
    chan_freq = chan_freq+chan_width
    chan_freq_arr.append(chan_freq)
    chan_width_arr.append(chan_width)

Fix the rest frequency in the SOURCE table: 

tb.open(myvis+"/SOURCE",nomodify=F)               
tb.putcol("REST_FREQUENCY",[freq_arr])                       
tb.close()

Fix the SPECTRAL_WINDOW table. Comment out anything which does not need correction.
Setting MEAS_FREQ_REF to '5' sets the reference frame to TOPO. 

tb.open("atca.ms/SPECTRAL_WINDOW",nomodify=F)
tb.putcell('REF_FREQUENCY',0,ref_freq)
tb.putcell("CHAN_FREQ",0,chan_freq_arr)
tb.putcell("CHAN_WIDTH",0,chan_width_arr)
tb.putcell("TOTAL_BANDWIDTH",0,bw)
tb.putcell("MEAS_FREQ_REF",0,5)
tb.close()

You can plot the data using plotms (use the Trans tab to set the rest frequency and the velocity convention and reference frame you wish to see, this does not have to be the one used in observing).

<> Inspect the CASA tables using 
browsetable('atca.ms')
Use the table keywords tab on the left to list the tables and click on the one you want to view.
NB CASA uses python, which is zero-indexed, so what in FITS was channel 1 is now channel 0 etc.
Close the MS in one tool before operating on it in a different way; if you get messages about locks (or if something is taking very long) check that you remembered to close the table you were editing and any visualisation and type 
clearstat()

A Avison, A M S Richards, Jan 2011