Development of the EVN pipeline in a Jupyter notebook environment, based on the CASA data processing.
You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
 

1045 lines
32 KiB

{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Generic EVN spectral line pipeline"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Transforming the continuum pipeline into a notebook that can handle spectral line observations. At first we start with MASER line experiments. The continuum calibration is based on the tutorial developed by Minnie Mao, Des Small, Jack Radcliffe and many others. The spectral line calibraion is based on the tutorial by Ross Burns."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Essential settings and derived parameters"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Most of the calibration parameters will be taken from the listobs file. This requires rerunning of the first part of the notebook, which will give some warnings and errors, but is totally benign. If there is already information available from previous calibration runs, fill in the details as much as possible."
]
},
{
"cell_type": "markdown",
"metadata": {
"tags": []
},
"source": [
"### Settings for the observation\n",
"Prior to running the notebook, make sure to set the parameters below and run these cells so CASA knows what to do. To run the pre-calibration steps, only the 'Experiment information' is required. That gives you the listobs file you need to set the remaining parameters. For the imaging step you need to derive the image size and cell size from the interferometer settings. \n",
"\n",
"Make sure the archival EVN data are stored with the same filenames as they come out of the archive: \n",
"* <obsid\\>.uvflg\n",
"* <obsid\\>.antab\n",
"* all the original FITS-IDI files, e.g. n14c3_1_1.IDI1\n",
"\n",
"Since this notebook handles only continuum data, only the FITS-IDI files with \\_1_1 are required. Any other file will be ignored. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Experiment information\n",
"Note that the path setting is the path *inside* the virtual container. The default path is `/home/jupyter/work/` for the Jupyterhub environment. When launching the Jupyter-CASA environment locally from Docker you can set this with the -v option. If you choose a different path, make sure the adjust the entire `mypath` variable below.\n",
"\n",
"The values listed in the cells below serve as an example only, and are based on experiment N14C3. Make sure you include the proper values for your own experiment here, using the same formatting."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# The experiment is in upper-case font for the archive directory and download task\n",
"expname = \"EB063C\"\n",
"# the observation-ID is used in lower-case font in the remainder of the notebook\n",
"obsid = expname.lower()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Point the mypath variable to where your data live\n",
"mypath='/home/jupyter/work/Spectral_line/EVN_CASA_pipeline/'\n",
"sys.path.append(mypath)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Reference antenna"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#This can be a list, e.g. 'EF,WB,O8,JB'\n",
"#The fringe fitting step will take the list in order from left to right for scans \n",
"# where the reference antenna is missing.\n",
"refant = 'EF'"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Calibration sources and scans"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Scan on a bright fringe finder\n",
"sbdscan = '37'\n",
"\n",
"# Timerange in the scan above where all stations have good data (1-2 minutes)\n",
"sbdtimer = '13:18:00~13:20:00'\n",
"\n",
"# Phase calibrator field or source name\n",
"mbdfield = '1848+283,J1640+3946'\n",
"\n",
"# Bandpass calibrator field or source name\n",
"bpassfield = '1848+283'"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Flag settings"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#Flag file with manual flags derived from data inspection\n",
"myflags=None\n",
"\n",
"#Fraction of edge channels to flag on either side: 0=no flagging, 0.15 = 15% flagging \n",
"#For most continuum experiments 10% on either side is fine\n",
"#For line experiments make sure to identify where the line is before flagging\n",
"edgefraction = 0.0"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Imaging settings"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#Pair of target and associated phase calibrator source\n",
"target1 = 'J2202+4216'\n",
"\n",
"# Values below not adjusted for EB063C\n",
"phasecal1 = '1848+283'\n",
"\n",
"#Calibrated visibility data for the pair above\n",
"tp1vis = 'J1849.ms'\n",
"\n",
"#Imaging cell size (note the formatting!)\n",
"cell=['0.3mas']\n",
"\n",
"#Image size in pixels (x,y)\n",
"imsize=[500,500]\n",
"\n",
"#Imaging number of iterations to clean\n",
"# 500 iterations is OK for most well calibrated EVN data\n",
"# set to 0 for a dirty image\n",
"niter=500 "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Date-time stamp\n",
"To ensure that each run of the notebook results in a new set of calibration tables, comment the second line and uncomment the first line below to add the date and time to the calibration tables. By default this is off, to ensure the directory doesn't overflow."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#stamp = datetime.datetime.utcnow().strftime('%y%m%d_%H%M')\n",
"stamp='debug'"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Derived parameters & download\n",
"No need to edit anything below, but make sure to also always run this cell."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Global imports\n",
"import glob\n",
"import shutil\n",
"import os\n",
"from evn_tools.plotcalng import plotcalng\n",
"from evn_tools.evn_import_exp import evn_import_exp\n",
"\n",
"#Visibility data (needed for pre-calibration steps)\n",
"# Note that the glob sorting is alphabetical and breaks for more than 10 IDI files, it puts file \n",
"# 10 between 1 and 2. This requires natural sorting. That should be implemented.\n",
"idifiles = sorted(glob.glob(mypath+'*_1_1.IDI*'))\n",
"vis = mypath+obsid+'_cont.ms'\n",
"vis_line = mypath+obsid+'_line.ms'\n",
"\n",
"# For spectral line data we need to verify if there are 1 or 2 correlator passes in the data.\n",
"# For multiple passes the _N_1 is the line pass and 1_1 the continuum pass, for single pass _1_1 is the line pass.\n",
"# This can only handle the ANTAB file for the continuum pass in case of multiple correlator passes. This \n",
"# may be sufficient, the ANTAB file for the line has essentially the same data but with different frequency\n",
"# setup for the SPW and channels. This can be adjusted when applying the calibration table.\n",
"if os.path.isfile(obsid+'_2_1.IDI1')==True:\n",
" dualpass = True\n",
" antabfile = mypath+obsid+'_1.antab'\n",
" print(\"Found two correlator passes\")\n",
"else:\n",
" dualpass = False\n",
" antabfile = mypath+obsid+'.antab'\n",
" print(\"Found one correlator pass\")\n",
"\n",
"#EVN archival calibration tables (needed for pre-calibration steps)\n",
"AIPSflag = mypath+obsid+'.uvflg'\n",
"gcfile = mypath+'EVN.gc'\n",
"\n",
"#CASA calibration tables\n",
"gcaltab = mypath+obsid+'_'+stamp+'.gcal'\n",
"tsystab = mypath+obsid+'_'+stamp+'.tsys'\n",
"sbdtab = mypath+obsid+'_'+stamp+'.sbd'\n",
"mbdtab = mypath+obsid+'_'+stamp+'.mbd'\n",
"bpasstab = mypath+obsid+'_'+stamp+'.bpass'\n",
"\n",
"flagfile = mypath+obsid+'.flag'\n",
"listobsfile = mypath+obsid+'.listobs'"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Download the data from the EVN archive\n",
"For use in the Jupyterhub environment, this cell can be run to automatically download the above listed (meta-)data for the selected experiment."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Check if the IDI files are present, should probably check for all IDI files AND meta-data\n",
"if os.path.isfile(idifiles[0])==False:\n",
" evn_import_exp(expname)\n",
"else:\n",
" print(\"Skipping import step\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Precalibration steps\n",
"These steps are done to ensure the meta-data associated with the observation is properly handled by CASA. The pre-calibration needs to be run only once. The cells check if the step has been performed. In Step 4 any existing measurement set is cleaned up, and the flags are restored to the flags from the archive. This can be useful when applycal has flagged too much data. \n",
"\n",
"### Step 1. Append the system temperature to the FITS-IDI files."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Append the system temperature to the FITS-IDI files\n",
"# Catch error when IDI files cannot be found because the path has not been set in first cell\n",
"\n",
"# Local imports \n",
"import astropy.io.fits as pyfits\n",
"from casavlbitools import fitsidi\n",
"\n",
"# Check if the IDI files have a TSYS column, if not, add one\n",
"try:\n",
" hdulist = pyfits.open(idifiles[0])\n",
" hdu = hdulist['SYSTEM_TEMPERATURE']\n",
" print('TSYS table already present, skipping the append step')\n",
"except KeyError:\n",
" print('Appending TSYS, this takes some time, go for a walk')\n",
" fitsidi.append_tsys(antabfile, idifiles)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Step 2. Convert UVFLG file to CASA readable flagfile"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"from casavlbitools import fitsidi\n",
"fitsidi.convert_flags(AIPSflag,idifiles,outfile=flagfile)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Step 3. Convert the gain curve to a CASA readable gain curve \n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"from casavlbitools import casa as mdconvert\n",
"if os.path.isdir(gcfile)==False:\n",
" mdconvert.convert_gaincurve(antabfile,gcfile,min_elevation=0.0,max_elevation=90.0)\n",
"else:\n",
" print(\"GC conversion already done.\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Step 4. Import the FITS-IDI files and generate the CASA Measurement Set\n",
"If the MS file exists, this step cleans the previous calibration and flag settings. It identifies the number of correlator passes, and for spectral line data it creates a channel-averaged continuum MS and a line MS at full spectral resolution."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Handle single pass individually:\n",
"if dualpass==False:\n",
" if os.path.isdir(vis_line)==False:\n",
" print(\"Importing data\")\n",
" idifiles = sorted(glob.glob(mypath+'*_1_1.IDI*'))\n",
" importfitsidi(fitsidifile= idifiles, vis = vis_line, constobsid=True, scanreindexgap_s=15.0, specframe='GEO')\n",
" # If the MS already exists, clean up the old calibration and restore the flags\n",
" else:\n",
" print(\"Import already done, clearing calibration and flags\")\n",
" clearcal(vis_line)\n",
" flagmanager(vis=vis_line,mode='restore',versionname='precal_flags')\n",
" # Generate a continuum pass in the next cell\n",
"\n",
"# Handle dual pass with hard-coded numbers, but such that this can be modified to have N>2\n",
"if dualpass==True:\n",
" #Verify if the continuum MS exist, otherwise do importstep\n",
" if os.path.isdir(vis)==False:\n",
" print(\"Importing continuum data\")\n",
" idifiles_cont = sorted(glob.glob(mypath+'*_1_1.IDI*'))\n",
" importfitsidi(fitsidifile= idifiles_cont, vis = vis, constobsid=True, scanreindexgap_s=15.0, specframe='GEO')\n",
" else:\n",
" print(\"Import for continuum data already done, clearing calibration and flags\")\n",
" clearcal(vis)\n",
" flagmanager(vis=vis,mode='restore',versionname='precal_flags')\n",
" #Verify if the line MS exist, otherwise do importstep\n",
" if os.path.isdir(vis_line)==False:\n",
" print(\"Importing line data\")\n",
" idifiles_line = sorted(glob.glob(mypath+'*_2_1.IDI*'))\n",
" importfitsidi(fitsidifile= idifiles_line, vis = vis_line, constobsid=True, scanreindexgap_s=15.0, specframe='GEO')\n",
" # If the MS already exists, clean up the old calibration and restore the flags\n",
" else:\n",
" print(\"Import for line data already done, clearing calibration and flags\")\n",
" clearcal(vis_line)\n",
" flagmanager(vis=vis_line,mode='restore',versionname='precal_flags')"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# For single pass line data, generate a continuum dataset by averaging down the full resolution line data\n",
"# Implement for user to set averaging, for now average down to 256 (32 channels and 8 spw)\n",
"# For faster processing this can be averaged further\n",
"\n",
"# Determine the number of spectral windows and channels in the MS\n",
"# (also needed for edge-channel flagging and for spectral window mapping)\n",
"msmd.open(vis_line)\n",
"nspw=msmd.nspw()\n",
"nchan=msmd.nchan(0)\n",
"msmd.done()\n",
"print(f'Number of channels per spectral window: {nchan}')\n",
"print(f'Number of spectral windows: {nspw}')\n",
"\n",
"\n",
"# Calculate the bin size for channel averaging down to 256 channels:\n",
"# (assume that continuum data typically has 8 spectral windows of 32 channels each)\n",
"nbins= int(nchan/(256))\n",
"print(f'Averaging bin size {nbins}')\n",
"\n",
"if dualpass==False:\n",
" mstransform(vis=vis_line,outputvis=vis,datacolumn='data',chanaverage=True,chanbin=nbins)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Step 5. Apply the flags from the archival UVFLG data"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"flagdata(vis=vis, mode='list', inpfile=flagfile,\n",
" reason='any', action='apply', flagbackup=False, savepars=False)\n",
"flagmanager(vis=vis,mode='save',versionname='archival_flags',\n",
" comment='Archival flags prior to any calibration')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Step 6. Generate a list of scans"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"if os.path.isfile(listobsfile)==False:\n",
" listobs(vis, listfile=listobsfile)\n",
"else:\n",
" print(\"A file with listobs information is already present\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Quality control of raw auto-correlations"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Make plots of the uncalibrated auto-correlations, phase and amplitude versus time and frequency, for the target source. Use the line pass for this. Plot for all stations."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Need to have bigger points in the plots, for some reason this makes the plotting even slower...\n",
"\n",
"#Amp vs frequency\n",
"plotms(vis=vis_line,xaxis='frequency',yaxis='amp',ydatacolumn='data',field=target1,\n",
" correlation='LL,RR',coloraxis='antenna1',antenna='*&&&')"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Amp vs time\n",
"plotms(vis=vis_line,xaxis='time',yaxis='amp',ydatacolumn='data',field=target1,\n",
" correlation='LL,RR',coloraxis='antenna1',antenna='*&&&')"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Phase vs frequency\n",
"plotms(vis=vis_line,xaxis='frequency',yaxis='phase',ydatacolumn='data',field=target1,\n",
" correlation='LL,RR',coloraxis='antenna1',antenna='*&&&')"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Phase vs time\n",
"plotms(vis=vis_line,xaxis='time',yaxis='phase',ydatacolumn='data',field=target1,\n",
" correlation='LL,RR',coloraxis='antenna1',antenna='*&&&')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# A priori calibration"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Generate the calibration tables for gain and system temperature. Older versions of the calibration tables will be removed first."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# gencal generates a cumulative table if the table exists, remove old tables first\n",
"if os.path.isdir(tsystab)==True:\n",
" shutil.rmtree(tsystab)\n",
"\n",
"gencal(vis,caltable=tsystab,caltype='tsys',uniform=False)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"if os.path.isdir(gcaltab)==True:\n",
" shutil.rmtree(gcaltab)\n",
"\n",
"gencal(vis,caltable=gcaltab,caltype='gc',infile=gcfile)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"After these steps, it can be useful to do a first inspection of the data, and write up the most obvious bad data in a flag command file. E.g.: \n",
"`mode='manual' antenna='SH&*' spw='0,1' timerange='13:15:15.0~13:23:46.0'` \n",
"This file can be added to the notebook in the `myflags` variable in the Flag settings at the start of the notebook. The flags in this file will be applied in the Flagging section below.\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Adjust plots below to include LL and RR (as above)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"scrolled": true
},
"outputs": [],
"source": [
"applycal(vis=vis,gaintable=[tsystab,gcaltab], flagbackup=False, parang=True)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"plotms(vis=vis, xaxis='frequency', yaxis='amp',ydatacolumn='data', \n",
" antenna=refant, correlation='ll',coloraxis='baseline', \n",
" scan=sbdscan,averagedata=True, avgtime='600')"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"plotms(vis=vis, xaxis='frequency', yaxis='amp',ydatacolumn='corrected', \n",
" antenna=refant, correlation='ll',coloraxis='baseline', \n",
" scan=sbdscan,averagedata=True, avgtime='600')"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"plotms(vis=vis, xaxis='frequency', yaxis='phase',ydatacolumn='data', \n",
" antenna=refant, correlation='ll',coloraxis='baseline', \n",
" scan=sbdscan,averagedata=True, avgtime='600')"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"plotms(vis=vis, xaxis='frequency', yaxis='phase',ydatacolumn='corrected', \n",
" antenna=refant, correlation='ll',coloraxis='baseline', \n",
" scan=sbdscan,averagedata=True, avgtime='600')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Flagging"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Before flagging: identify the frequency range of the signal and other lines of interest. How is this done?"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Flag autocorrelations: this isusually done post-correlation, but not always. When using ACCOR, do not flag the autocorrelations."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"flagdata(vis,mode='manual',autocorr=True,flagbackup=False)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Flag bad antenna and poor scans. The user needs to add the commands here or make a flag file. "
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"flagdata(vis,mode='list',inpfile=myflags,flagbackup=False)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Flag the edge channels: typically 10% on either side. The fraction is set in the parameter box above."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Flag the edge channels\n",
"flagfraction = int(nchan/(100*edgefraction)) \n",
"start = str(flagfraction-1)\n",
"end = str(nchan-flagfraction)\n",
"spwflag = '*:0~'+start+';'+end+'~'+str(nchan-1)\n",
"\n",
"flagdata(vis,mode='manual',spw=spwflag,flagbackup=False)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Flag the first 5 seconds of every scan"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"flagdata(vis,mode='quack',quackinterval=5,flagbackup=False)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Now that all the flagging is done, save the flags to restore them prior to reruns of the calibration steps below."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"flagmanager(vis,mode='save',versionname='precal_flags',\n",
" comment='Flags from Tsys, gaincal, bad data and edge channels')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"This is where the data can be taken through AOFlagger or additional manual flagging is done"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Fringe fitting"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Instrumental delay \n",
"The initial step of fringe fitting includes only the instrumental delay. This requires a bright target source as input, and only a single scan. Within the scan a solution is determined for each spectral window. Rates cannot be determined due to the lack of a longer time baseline in a single scan. In this notebook the central few minutes of the scan on the fringe finder source are used. The full scan can also be used, but can result in somewhat poorer solutions due to missing information at the edges of the scan. Typically 1-2 minutes of good data is better than a full scan."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"fringefit(vis, caltable=sbdtab, timerange=sbdtimer, \n",
" solint='inf', zerorates=True, refant=refant,minsnr=50,\n",
" gaintable=[gcaltab,tsystab],\n",
" interp=['nearest','nearest,nearest'],\n",
" parang=True)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Check the output in the log to see how many solutions are found by fringefit. This is a good initial indication if this worked well or not. You should have a solution per spectral window."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"applycal(vis=vis,gaintable=[tsystab,gcaltab,sbdtab], flagbackup=False, parang=True)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"plotms(vis=vis, xaxis='frequency',yaxis='phase', \n",
" ydatacolumn='corrected', antenna=refant, correlation='ll',\n",
" coloraxis='baseline', scan=sbdscan, timerange=sbdtimer,\n",
" averagedata=True, avgtime='600')"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"flagmanager(vis=vis,mode='restore', versionname='precal_flags')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Multi-band delay\n",
"The next step is the multi-band delay, which includes phases as a function of time and frequency for all the phase-reference scans (delay and rate is determined). Solutions are determined per scan by setting `solint = 'inf'`. A lower value for the SNR is set as well. "
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"fringefit(vis, caltable=mbdtab, field=mbdfield,\n",
" solint='inf', zerorates=False, refant=refant,\n",
" combine='spw', minsnr=5,\n",
" gaintable=[gcaltab, tsystab, sbdtab],\n",
" parang=True)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Check the output in the log to see how many solutions are found by fringefit. This is a good initial indication if this worked well or not. Note that the second fringe fit step determines a single solution for all spectral-windows."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Some inspection of the solutions and corrected data after the fringe fitting step"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"plotcalng(caltable=mbdtab,xaxis='time',yaxis='delay')"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"tags": []
},
"outputs": [],
"source": [
"plotcalng(caltable=mbdtab,xaxis='time',yaxis='rate')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Bandpass"
]
},
{
"cell_type": "markdown",
"metadata": {
"tags": []
},
"source": [
"The final step in the calibration is the bandpass. CASA performs a complex bandpass, which means it fits both amplitude and phase. As a consequence the bandpass correction can in principle change the phase solutions which would then require a second round of fringe fitting to be corrected. It is important to verify the solutions after the bandpass step.\n",
"\n",
"The bandpass can be obtained from the combined data of the phase calibrator or on the specific bandpass calibrator scan. Note that the multi-band delay corrections do not apply for the direction of the bandpass calibrator, and its phases are therefore not calibrated. This should not affect the amplitude bandpass. As the phases have been corrected for instrumental delay, the phase bandpass should also be representative (VERIFY)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"bandpass(vis, caltable=bpasstab, field=bpassfield,\n",
" gaintable=[gcaltab, tsystab, sbdtab, mbdtab],\n",
" interp=['nearest','nearest,nearest','nearest','linear'],\n",
" spwmap=[[],[],[], int(nspw)*[0]],\n",
" solnorm=True, solint='inf', refant=refant, \n",
" bandtype='B', parang=True)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"plotcalng(bpasstab,iteration='antenna')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Apply calibration and split sources"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The last step is to apply all calibration tables to the data. Note the specific settings of the swpmap parameter. This step can take a while to run. The splitting is done per pair of target - phase reference source, which is a matter of personal taste."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"applycal(vis,\n",
" gaintable=[gcaltab, tsystab, sbdtab, mbdtab, bpasstab],\n",
" interp=['nearest','nearest,nearest','nearest','linear','linear,linear'],\n",
" spwmap=[[], [], [], int(nspw)*[0],[]],\n",
" parang=True)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"As a final calibration quality check, plot the amplitudes for the timerange in which the single-band delay (instrumental delay) was fitted. Then also plot the phases for the phase calibrator, averaged over the entire scan."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"plotms(vis=vis, xaxis='frequency', yaxis='amplitude',\n",
" ydatacolumn='corrected',antenna=refant, correlation='ll,rr',\n",
" coloraxis='baseline', timerange=sbdtimer)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"plotms(vis=vis, xaxis='frequency', yaxis='phase',\n",
" ydatacolumn='corrected', antenna=refant, correlation='ll',\n",
" coloraxis='baseline', field=phasecal1,\n",
" averagedata=True, avgscan=True)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"To do the splitting properly, pairs of phase-cal and target need to be specified by the user in the parameter window at the start. This is not yet implemented."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"field1=''.join((target1,',',phasecal1))\n",
"split(vis, outputvis=tp1vis, field=field1, datacolumn='corrected')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Visualization\n",
"Make an initial cleaned image of phase calibrator and target. Inspect the imaging settings above. For additional pairs of sources, copy and paste the four cells below and add the information to the paramter settings above."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"tclean(vis=tp1vis, field=phasecal1, imagename='clean_p1',\n",
" imsize=imsize, cell=cell,niter=niter)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"imview('clean_p1.image')"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"tclean(vis=tp1vis, field=target1, imagename='clean_t1',\n",
" imsize=imsize, cell=cell,niter=niter)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"imview('clean_t1.image')"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Casa",
"language": "casa",
"name": "casapy"
},
"language_info": {
"mimetype": "text/plain",
"name": "Casa"
}
},
"nbformat": 4,
"nbformat_minor": 4
}