Open in Colab: https://colab.research.google.com/github/casangi/examples/blob/master/community/ALMAephemimagingVenusMfsMosaicTrackfield.ipynb

---

ALMA ephemeris object imaging: continuum (mosaic gridder - trackfield)

Use ‘TRACKFIELD’ in phasecenter parameter

Original Author(s): ttsutsum@nrao.edu

Description

The following demo describes continuum imaging of ALMA observation of Venus. For imaigng of the sources that move significantly during the observations, proper ephemeris data need to be used. For the current ALMA data, the epehemeris tables used to track the source during the observation generally included in the ASDM. The task, importasdm attaches the ephemeris table(s) in FIELD sub-table in the generated MS. Here we starts from such MS data and also uses a short cut option called ‘TRACKFIELD’ in tclean to use the attached ephemeris table(s).

Data

The data, venus_ephem_test.ms (34MB) is available in the casatestdata repository, under https://open-bitbucket.nrao.edu/projects/CASA/repos/casatestdata/browse/measurementset/alma/.

Or at NRAO, casatestdata is accessible at /home/casa/data/casatestdata

Installation

Install CASA

Skip this step if CASA is already installed

import os
print("Installing CASA ...")
os.system('pip install casaconfig')
os.system('pip install casatools==6.7.0.31')
os.system('pip install casatasks==6.7.0.31')

# for Colab
mydatapath = '/content/.casa/data'
#
#mydatapath='data'

import pathlib
from casaconfig import config
if not pathlib.Path(mydatapath).exists():
  pathlib.Path(mydatapath).mkdir(parents=True)
config.measurespath=mydatapath

Install additional modules used in this notebook

import os
os.system("pip install aplpy")
os.system("pip install ipympl")
print("complete")

This notebook uses functions from EphemerisObjectImagingDemoFunctions.py so that file must be present in the directory running this notebook and need to import them.

# if you haven't downloaded EphemerisObjectImagingDemoFuntions.py
!wget https://raw.githubusercontent.com/casangi/examples/refs/heads/CAS-13662/community/EphemerisObjectImagingDemoFunctions.py

from EphemerisObjectImagingDemoFunctions import *

Import required tools and tasks

from casatasks import getephemtable, tclean, exportfits, listobs, imstat, imhead
from casatools import measures, quanta, table, image
import os
from astropy.io import fits
from astropy.wcs import WCS
from matplotlib import pylab as pl
import numpy as np
import aplpy
me = measures()
qa = quanta()
tb = table()
ia = image()

Imaging with epehemeris data attached to the MS

In order to track the ephemeris source during the observations, the proper ephemeris data is used. Ususally such ephemeris data table(s) are attached to the Measurement Set automatically during importadm when they are available. There may be multiple tables per MS but only single table is linked each field. The tables are used to calculate phases for gridding visibility. The ephemeis table is also used to stop the moving source. Setting phasecenter=’TRACKFIELD’ allows to use the same table for phase calculation and stopping of the moving source.

Obtain data

msfile = 'venus_ephem_test.ms'

# Use the ALMA Venus data in casatestdata
# At NRAO, point to the casatestdata directory
#rootdatapath = '/home/casa/data/casatestdata'
# or set a proper path to the local copy of the casatestdata repository
# e.g.
#rootdatapath = '/Volumes/ssd1/casatestdata'
#inputdatapath=f'{rootdatapath}/measurementset/alma/'

For Colab or for downloading the input data locally execute next cell, otherwise skip it

# If running on Colab

# sparse-checkout of the data repository

# install git lfs
os.system('pip install git-lfs -q')
os.system('git-lfs install')

os.system('git clone --filter=blob:none --sparse https://open-bitbucket.nrao.edu/scm/casa/casatestdata.git mytestdata')
%cd mytestdata
os.system('git sparse-checkout set measurementset/alma/venus_ephem_test.ms')
%cd ..

# If running on Colab
inputdatapath = 'mytestdata/measurementset/alma/'

# Check the data exists
from casatasks import listobs
listobs(inputdatapath+msfile)

Imaging

# set msfile to include the full path
msfile = inputdatapath+msfile

imgname='Venus-mosaic-mfs-trackfield'

delete_tcleanimages(imgname)

Run mosaic with usepointing = True to correctly shift the beam

ret=tclean(vis=msfile, field='Venus', imagename=imgname, imsize=[480, 432], cell=['0.14arcsec'],
           phasecenter='TRACKFIELD', specmode='mfs', gridder='mosaic', usepointing=True, niter=100)
ret['iterdone']

100

Examination of the resultant image

Convert a CASA image to FITS image for displaying the image via aplpy

# Define functions to use in displaying the image
def displayImage(imgname, markers={}, title=''):
    ''' convert casa image to fits and display the fits image'''
    # export to a FITS image (overwrite if the output exists)
    from casatools import quanta
    _qa = quanta()
    fitsimagename = imgname+'fits'
    exportfits(imagename=imgname, fitsimage=fitsimagename, overwrite=True)
    fits.setval(fitsimagename, 'TIMESYS', value='utc')
    fig = pl.figure(figsize=(15,15))
    img = aplpy.FITSFigure(fitsimagename, subplot=[0.1,0.1,0.5,0.5])
    img.show_colorscale()
    if markers!={}:
       for ky, pos in markers.items():
          if 'ra' in pos and 'dec' in pos:
             print(f'ky={ky}')
             if 'ext' in ky:
                 color = 'red'
             else:
                 color = 'yellow'
             img.show_markers(pos['ra'], pos['dec'], edgecolor=color, marker='o',s=10, alpha=0.5)
             marker_ra = _qa.time(_qa.quantity(pos['ra'],'deg'),prec=9)[0]
             marker_dec = _qa.angle(_qa.quantity(pos['dec'],'deg'),prec=9)[0]
             print(f'{ky} at: {marker_ra}, {marker_dec}')
    pl.title(title)
    fig.canvas.draw()
    print('image center:',printImageCenter(imgname))


def printImageCenter(imgname):
    from casatools import image, quanta
    _ia = image()
    _qa = quanta()
    _ia.open(imgname)
    shape = _ia.shape()
    csys = _ia.coordsys()
    center = csys.toworld([shape[0]/2.,shape[1]/2.])
    _ia.done()
    return _qa.time(_qa.quantity(center['numeric'][0], 'rad'),prec=9)[0], _qa.angle(_qa.quantity(center['numeric'][1],'rad'),prec=9)[0]

# Returns direction of the ephemeris object at a given time (parallax collected)
def ephem_dir(ephemtab, refep, observatory):
    from casatools import measures
    _me = measures()
    _me.framecomet(ephemtab)
    _me.doframe(_me.observatory(observatory))
    _me.doframe(_me.epoch('utc',refep))
    return _me.measure(_me.measure(_me.direction('COMET'), 'AZELGEO'),'ICRS')

Get the direction of Venus at the first timestamp of the select data from the ephemeris table

# Find the path of the attached ephem. table for fieldid=0
ephemtab = get_attachedEphemtablepath(msfile,'Venus')
print(ephemtab)

mytestdata/measurementset/alma/venus_ephem_test.ms/FIELD/EPHEM0_Venus_58491.4.tab/

from casatasks import imhead
iminfo=imhead(imgname+'.image', mode='list')
iminfo['date-obs']

'2019/01/08/11:21:50.208000'

ineph_dir = ephem_dir(ephemtab, iminfo['date-obs'],'ALMA')
print(ineph_dir)
ineph_rad = qa.convert(ineph_dir['m0'],'deg')
ineph_decd = qa.convert(ineph_dir['m1'],'deg')
inephmarker = {'ineph':{'ra':ineph_rad['value'], 'dec':ineph_decd['value']}}

{'m0': {'unit': 'rad', 'value': -2.1053946095161686}, 'm1': {'unit': 'rad', 'value': -0.29607600119188454}, 'refer': 'ICRS', 'type': 'direction'}

displayImage(imgname+'.image', markers=inephmarker, title='Venus mosaic image with the internal ephem table (TRACKFIELD)')

INFO: Setting slices=[0, 0] [aplpy.core]

WARNING: FITSFixedWarning: 'datfix' made the change 'Set MJD-OBS to 58491.473498 from DATE-OBS'. [astropy.wcs.wcs]
WARNING: FITSFixedWarning: 'obsfix' made the change 'Set OBSGEO-L to   -67.754929 from OBSGEO-[XYZ].
Set OBSGEO-B to   -23.022886 from OBSGEO-[XYZ].
Set OBSGEO-H to     5053.796 from OBSGEO-[XYZ]'. [astropy.wcs.wcs]

INFO: Auto-setting vmin to -3.801e+00 [aplpy.core]
INFO: Auto-setting vmax to  3.664e+00 [aplpy.core]
ky=ineph
ineph at: 15:57:28.746, -016.57.50.059
image center: ('15:57:28.746', '-016.57.50.059')

<Figure size 1500x1500 with 0 Axes>

Note: the marker in yellow is set to the direction based on the internal ephemeris table at the first time of the selected data and its coordintes agrees with the coordinates of the center of the image.

print("internal ephem dir = ",qa.time(ineph_dir['m0'],prec=9),qa.angle(ineph_dir['m1'],prec=9))
printImageCenter(imgname+'.image')

internal ephem dir =  ['15:57:28.746'] ['-016.57.50.059']

('15:57:28.746', '-016.57.50.059')

stats = imstat(imgname+'.image')
print(f"{stats['max']}@{stats['maxpos']}")

[3.33058977]@[177 172   0   0]