Getting Started | Documentation | Glish | Learn More | Programming | Contact Us |
Version 1.6 Build 363 |
|
M.H. Wieringa
This chapter describes how to reduce data from the ATNF's Compact Array (ATCA) in AIPS++. It should be read in conjunction with several other chapters in Getting Results, including the chapters which cover i) synthesis calibration, ii) synthesis imaging, including general imaging and mosaicing, iii) display of data, and iv) image analysis.
The data reduction categories covered in this chapter include those for continuum data, continuum polarimetry, spectral line data, and spectral line polarimetry. Individual steps in each data reduction sequence are linked to underlying Glish recipes, or illustrated using Glish script excerpts where appropriate. These data reduction steps can equivalently be executed from the Toolmanager GUI interface. The primary data reduction tools used in this chapter are Calibrater and Imager, which constitute the lowest level synthesis interface.
The reduction sequences described in this chapter require some familiarity with the tools listed here:
The User Reference Manual contains the primary documentation for all these tools and their functions. The aim of this document is to provide pointers to the relevant chapters and examples, and to describe their use in ATCA reduction.
Users are also referred to the Miriad-aips++ dictionary which is a separate chapter in Getting Results, which provides a mapping of some common Miriad tasks to their counterparts in AIPS++.
Visibility data in AIPS++ are stored in a particular type of AIPS++ Table known as a MeasurementSet (MS). An MS consists of a main table containing the visibility data, with associated sub-tables containing auxilliary or secondary information. The full MS data format is described fully in AIPS++ Note 229.
To load the basic tools likely to be used in ATCA reduction, when working from the command line, the following initialization is recommended at the start of each AIPS++ or Glish command-line session:
include 'synthesis.g' include 'flagger.g' include 'image.g'
This is not required when working from the Toolmanager GUI, in which case the tools are loaded automatically when they are selected.
ATCA data in on-line RPFITS format are read into an AIPS++ MS using the atcafiller tool.
The atcafiller tool can be used to read data from disk or CDROM. It allows selection on scan number, frequency band, frequency range and field name using the select function.
To run the atcafiller:
include 'atcafiller.g' # af := atcafiller('c999.ms', filenames="1999-01-01_0101.C999 1999-01-01_0202.C999", options=reweight) # Construct af.select(firstscan=3); # data selection af.fill() # fill the actual data af.done() # destroy filler tool # m := ms('c999.ms') # access the MS
In this example, we construct the filler tool from from two input RPFITS files, specifying that we wish to apply lag reweighting to suppress ringing and sources. Then the select function is used to select which data we wish to read from the RPFITS files. Then we call the fill function to actually read and convert the data to the MS.
Data which have been exported from AIPS or Miriad as UVFITS files can be read into an AIPS++ MS using the general UVFITS filler ms.fitstoms. For example:
include 'ms.g' m := fitstoms(msfile='atca.ms', fitsfile='/aips++/data/demo/ATCA.UV') m.done()
Once the data havd been filled to a MeasurementSet, they can be examined or displayed using various separate tools. The Catalog tool provides an overview of AIPS++ data files at the highest level, as individual AIPS++ tables of different types (e.g. images, MS, etc) in a given directory. The Catalog tool can be started from Glish as:
include 'catalog.g' dc.gui()
Here, dc is the default catalog tool. It can also be accessed via the Toolmanager interface from the Tools in Use selection. A mouse click on an individual MS listed in the catalog GUI followed by a click on the Summary button, will write a summary of the data in the MS (including frequency and polarization settings and fieldnames) to the logger. These will be handy to have ready during calibration.
When you click on the View button (to view the presently selected catalog item) you will launch the Tablebrowser, which provides facilities to browse and examine all data fields in the MS main and sub-tables. Note that indices in the MS, such as FIELD_ID, will appear as zero-relative when viewed at the low level provided by the Tablebrowser. All AIPS++ tools present selection in a one-relative basis.
The primary graphical visibility data visualization tool in AIPS++ at present is msplot. It can be started from Glish via:
include 'msplot.g' p := msplot('atca.ms', edit=T); # Allow editing as well as inspection
It can also be launched from the Toolmanager via the Synthesis package, and msplot module.
The msplot tool allows raster and line plots of the observed, corrected, and model data against a range of MS indices and parameters, as well as associated plots such as uv-coverage. Facilities are provided to support versatile data selection, range setting and plot iteration. Interactive editing is also supported, as described in the next section.
Command-based, interactive and automated editing capabilities are supported, and are considered separately below.
The Flagger tool offers several directed editing (flagging) operations. These include editing in specific time intervals, on specific MS indices, such as FIELD_ID, or for specific polarizations.
fg:= flagger('mydata.ms') fg.auto() # Flag all autocorrelations fg.setantennas(ants=5) # Select antenna 5 fg.setids(spectralwindowid=1,fieldid=2) # Select spectral id and field id fg.timerange(starttime="24-FEB-2000/10:00:05", # Flag selected data in specified time range endtime="24-FEB-2000/10:09:37");
A full range of command-based flagging operations are available, as described fully in the flagger documentation.
It can be extremely time consuming to make directed editing commands. Better is to use intelligent algorithms to automatically edit out bad data.
Automated editing methods are available in the Flagger and Autoflag tools. These include clipping, median flagging and related methods.
For example, to clip data outside of a specified range with the Flagger tool:
fg.setflagmode('flag'); fg.filter(column='DATA', operation='range', comparison='Amplitude', range=['0.00001Jy','1.34Jy']); fg.done();
The specified lower limit ensures that spurious data with zero amplitude are removed. To reverse all flags in the data:
fg.setflagmode('unflag'); fg.query(query='FIELD_ID >= 0');
The Autoflag tool has more sophisticated algorithms. Here are a few examples.
Interactive editing is possible on almost all graphical views of the data which are possible in msplot. To enable the editing option, it is necessary to set parameter edit=T when the msplot tool is created. The GUI panel controlling editin in msplot allows global flag masks to be set, such as that for all polarizations if any one is bad, as well as access to the list of accumulated edits, reversal of edits already applied, and other related utility operations. Raster displays of visibility data are possible also, and maybe useful for ATCA line data editing.
As noted above, AIPS++ visibility data are stored in a MeasurementSet (MS). The observed data are stored in the DATA column in the MAIN table of the MS; these are the raw data as loaded from the ATCA archive or from a UVFITS file. Associated with the DATA column are related columns to hold the most recent version of the calibrated data (CORRECTED_DATA), and the latest visibility plane representation of the source model (MODEL_DATA). The latter two columns are filled in by the calibrater and imager tools, respectively. The observed DATA do not change during reduction, but the related columns do, and should be viewed in this context. When plotting the data using a tool such as msplot, the data sub-types (observed, corrected and model) can be selected individually.
Calibration and imaging responsibilities are divided between the Calibrater and Imager tools. The calibrater tool handles visibility-plane calibration while the Imager tool deals with image formation, deconvolution and image-plane calibration. As such, they broadly deal with the left-hand and right-hand sides of the AIPS++ Measurement Equation (ME) respectively. This calibration formalism was defined by Hamaker, Bregman and Sault (1996a, 1996b) and is discussed more comprehensively in the Synthesis Calibration chapter of Getting Results.
When solving for a visibility-plane calibration component, the Calibrater tool requires that the calibrator model be set first. This is stored in the MODEL_DATA column, and used by the Calibrater tool to form when solving for calibration components. In general, this process requires transforming a source model from the image-plane to the visibility-plane, including any image-plane calibration effects. For this reason, the AIPS++ functions for converting source models to the visibility plane are attached to the Imager tool.
For the ATCA, the default source model is usually either 1934638 or point sources of unit flux density. If you have used 0823500 as your primary calibrator, or given 1934638 another name, you'll need to set the flux density for it by hand (see below).
To compute source models of this type:
imagr:= imager('atca.ms'); # Create an imager tool imagr.setjy(fieldid=1); # Compute source model for primary calibrator OR imagr.setjy(fieldid=1,spwid=1,fluxdensity=5.6); # Set flux density explicitly for source imagr.setjy(fieldid=1,spwid=2,fluxdensity=5.8); imagr.done()
In these examples, the source is specified through its fieldid argument. Running the Imager.setjy function without arguments will compute the default point source flux density for all recognized standard sources for all spectral windows in the data. The Imager.setjy function also allows selection by field id. and spectral window id. as required, as well as the specification of an arbitrary Stokes (I, Q, U, V) for individual sources where necessary. The standard sources recognized are [(3C286, 1328+307, 1331+305), (3C48, 0134+329, 0137+331), (3C147, 0538+498, 0542+498), (3C138, 0518+165, 0521+166), (1934-638), and (3C295, 1409+524, 1411+522)]. At the present time, the flux density is computed for the frequency mid-point of each spectral window using the Perley-Taylor (1995) spectral index coefficients by default. All unrecognized sources, will be set to unit Stokes I flux density. (Note that the MODEL_DATA for all sources is already set to unity when this column if first created, so you don't normally need to do this again here.)
Arbitrary source models can be set using the Imager.ft method, which will transform images and source model components, as required.
The individual calibration effects in the ME, which are expressed as Jones matrices, are described in the Synthesis Calibration chapter of Getting Results, and the reader is referred there for a more detailed discussion. Each Jones matrix has a specific type; the complex electronic gain is the G Jones component, the bandpass gain is the B Jones component and the polarization leakage is the D Jones component.
For ATCA data with 4 polarizations you will want to solve for G, D and B. If you have spectral-line data with only 1 or 2 polarizations, then you need solve for G and B only. You can solve for more than one spectral window at a time. Solving for D Jones with XX/YY only data will not do anything useful, but generally won't harm either.
First solve for G, D and B Jones using the primary calibrator:
cal := calibrater('atca.ms') # Create a calibrater tool cal.setdata(msselect='FIELD_ID==1') # Select the primary calibrator (in this case id=1) cal.setapply(type='P', t=10.0); # Correct for parallactic angle with 10 sec sampling. cal.setsolve(type='G', t=30.0, table='cal.G') # Solve for G Jones in 30 sec intervals cal.setsolve(type='B', t=1.e6, table='cal.B') # Solve for B Jones with an interval larger than your observation cal.setsolve(type='D', t=1.e6, table='cal.D') # Solve for D Jones with an interval larger than your observation cal.solve(); # Do the actual solution cal.plotcal(tablename='cal.G'); # Plot the calibration solution obtained cal.done() # Destroy the calibrater tool
Now solve for the G Jones term using the secondary calibrator(s), e.g.:
cal := calibrater('atca.ms') # Create a new calibrater tool cal.setdata(msselect='FIELD_ID==2') # Select the secondary calibrator (in this case id=2) cal.setapply(type='P', t=10.0) # Correct for parallactic angle with 10 sec sampling. cal.setapply(type='B', table='cal.B') # Apply B Jones solved for on primary cal.setapply(type='D', table='cal.D') # Apply D Jones solved for on primary cal.setsolve(type='G', t=60.0, table='cal.G',append=T) # Solve for G Jones in 60 sec intervals; append output cal.solve() # Do the actual solution cal.plotcal(tablename='cal.G',plottype='PHASE',fields=2); # Plot the calibration solution obtained
Generally you should use the above scheme, but if you want or need to (e.g. if you used 0823-500 as your primary calibrator), you can solve for D (and/or B) using the secondary instead of the primary, and just use the latter for the absolute fluxscale. If you want to solve for D using the secondary the process would look like this:
cal := calibrater('atca.ms') # Create a calibrater tool cal.setdata(msselect='FIELD_ID==1') # Select the primary calibrator (in this case id=1) cal.setapply(type='P', t=10.0); # Correct for parallactic angle with 10 sec sampling. cal.setsolve(type='G', t=30.0, table='cal.G') # Solve for G Jones in 30 sec intervals cal.setsolve(type='B', t=1.e6, table='cal.B') # Solve for B Jones with an interval larger than your observation cal.solve(); # Do the actual solution cal.done() # Destroy the calibrater tool
Now solve for the G Jones and D Jones terms using the secondary calibrator, e.g.:
cal := calibrater('atca.ms') # Create a new calibrater tool cal.setdata(msselect='FIELD_ID==2') # Select the secondary calibrator (in this case id=2) cal.setapply(type='P', t=10.0) # Correct for parallactic angle with 10 sec sampling. cal.setapply(type='B', table='cal.B') # Apply B Jones solved for on primary cal.setsolve(type='G', t=60.0, table='cal.G',append=T) # Solve for G Jones in 60 sec intervals; append output cal.setsolve(type='D', t=1.e6, table='cal.D') # Solve for D Jones using secondary cal.solve() # Do the actual solution cal.plotcal(tablename='cal.G',plottype='PHASE',fields=2); # Plot the calibration solution obtained
Note that the G Jones term, although designated as electronic gain, will include an atmospheric contribution, as solved for here. If separable, atmospheric corrections (e.g., WVR corrections) can be separately set or solved for as T Jones. The inclusion of a parallactic angle correction P Jones is required for all ATCA reduction, see the following section.
The gain solutions were generated with an unpolarized unit point source model for all calibrators unrecognized by Imager.setjy. In practice calibrators often have noticable linear polarization of up to 10% or so. This is not a problem for telescopes with circular feeds, but for telescopes with linear feeds, like the ATCA, this causes a small time variation in the gains proportional to the percentage polarization of the calibrator. The calibrator polarization can be determined from the variation of the X and Y gain solutions with changing feed position angle. After determining the polarization we can correct the gains for this effect. Note that this only works correctly if the calibrator has been observed over a sufficient range in position angle. Generally 5 or more scans with a span of at least 90 degrees in position angle is sufficient. The quality of the polarization solution can be checked by looking at the fit error, it should generally be smaller than 0.5% in Q and U. The correction is performed using the Calibrater.linpolcor function.
cal.linpolcor(tablein='cal.G', fields="0405-385 2254-367");
This corrects the calibration table entries for sources 0405-385 and 2254-367 and updates the table. If you add a tableout argument, it will write the updated table to a new file. (Note: using double quotes like in the example above produces a vector of values, split at each space, another way to specify the fields argument would be: fields=['0405-385','2254-367']).
For point sources of unknown flux density, but for which valid solutions have been derived against a unit point source model, the absolute flux density can be established by scaling the mean gain moduli against those of the known amplitude calibrators recognized by Imager.setjy. This scaling is achieved by adjusting the G Jones calibration table, after determining the amplitude calibration factors, as:
cal.fluxscale(tablein='cal.G', reference='1934-638', # scale the gains transfer="0405-385 2254-367"); cal.done();
In this case, the amplitude scale is set with 1934638 as the reference. The calibration table cal.G will be updated with the scaled G Jones calibration factors, and will be used going forward in the reduction.
Imaging is performed using the Imager tool, which supports a broad range of imaging modes and deconvolution algorithms, as described in the Synthesis Imaging chapter in Getting Results.
Before imaging, the data need to be corrected by applying the Jones matrices determined in the calibration:
cal:= calibrater('atca.ms') cal.setapply(type='P', t=10.0) # Apply parallactic angle cal.setapply(type='G', table='cal.G') # Apply the G Jones solutions cal.setapply(type='D', table='cal.D') # Apply the D Jones solutions cal.setapply(type='B', table='cal.B') # Apply the B Jones solutions cal.correct() # Correct the data cal.done()
This action fills the CORRECTED_DATA column, which can be used directly in further imaging.
The Imager tool, which can be used to make images for each field or source, requires initial selection of the visibility data and initialization of the imaging coordinates. Imaging, deconvolution and image restoration operations are then possible, as described in the Imaging chapter of Getting Results or in more detail in the Imager documentation. A single-field, continuum ATCA image can be made as follows:
imagr := imager('atca.ms') # Create an imager tool imagr.setdata(fieldid=3, mode='channel', nchan=25, # Select field id. 3, 25 channels start=5, step=1, spwid=2) # and spectral window id. 2 # imagr.setimage(nx=256, ny=256, cellx='0.4arcsec', # Set imaging parameters celly='0.4arcsec', stokes='IQUV', # use mfs mode to combine all selected channels mode='mfs', nchan=1) # imagr.clean(algorithm='clark', niter=500, gain=0.1, # Image, and deconvolve using the Clark CLEAN. threshold='0.0Jy', model='1331-g1.model', # Save model image='1331-g1.restored', # Save restored residual='1331-g1.residual') # Save residual imagr.done(); # im:= image('1331-g1.restored') # View restored image im.view() # # im2 := image('1331-g1.residual'); # View the residual image im2.view() # im.done() im2.done()
Note that the deconvolved image model is written to the image file 1331-g1.model. If non-zero on entry to the clean function, it will be subtracted before imaging. This is how you restart the deconvolution. In the above example we are not restarting so the model image is made for you.
The imaging parameters may be selected automatically using Imager.advise (not recommended for ATCA data at this time). For further information on any imaging parameters or options, please consult the full Imager documentation.
This chapter is based on the VLA chapter, with modifications where appropriate.