 News |
FAQ
|
| Search
|
Home
|
| Getting Started | Documentation | Glish | Learn More | Programming | Contact Us |
|
| Version 1.6 Build 363 |
|
A.J. Kemball & G.A. Moellenbrock
This chapter describes how to calibrate and image VLA data in AIPS++. For fairly standard VLA continuum and spectral line observations (including polarimetry), it can stand on its own, but for complicated or unconventional observing strategies, it will likely be necessary to consult several other chapters in Getting Results in AIPS++, including those which cover: i) synthesis calibration; ii) synthesis imaging, including general imaging, mosaicing and wide-field imaging; iii) use of the viewer, 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 illustrated using Glish script excerpts. These data reduction steps can equivalently be executed from the Toolmanager GUI interface, which is the recommended approach for new users.
The reduction sequences described in this chapter presume 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 the present document is to provide pointers to the relevant sections and examples, and to describe their use in VLA reduction.
Users are also referred to the AIPS-AIPS++ dictionary, which is a separate chapter in Getting Results in AIPS++, which provides a mapping of some common AIPS tasks and adverbs to their counterparts in AIPS++.
Since reduction of continuum, spectral line, and polarimetry data in AIPS++differ only the addition or omission of certain steps, and in the details of a few parameter settings, only one pass through the basic reduction path is described here, with labelled sequences of glish for both continuum polarimetry and spectral line reduction included in parallel. Some actions apply to only one type of reduction, but all methods are instructive for either type of observation.
The two sample datasets used in the examples are available from the AIPS++ data respository. The continuum polarimetry dataset is a full-polarization VLA CnB-configuration observation of ten VLA calibrators at 5 GHz from the VLA/VLBA polarization position angle calibration monitoring project. In this dataset, there is no target source per se, but the principles of standard calibration (described below) apply. The spectral line dataset is a short VLA D-configuration observation of HI (1413 MHz) in NGC5921. This dataset was observed in total intensity only (RR & LL), with 63 channels in a single spectral window. A single, near-by phase calibrater was observed, as well as 3C386 for flux density calibration.
The glish sequences listed below to describe calibration and imaging steps are available on the AIPS++ recipes page as full end-to-end scripts, easily adaptable to one's own needs. The glish sequences can also be used as a guide for operation using the tool manager GUI interface. All parameter names and settings correspond exactly to those listed in the GUI interface (only relevant parameters are listed in the glish sequences; others may be safely left set to default values). Consult the Getting Started in AIPS++ document for tips on use of the tool manager.
To load the basic tools likely to be used in VLA reduction, when working from the command line, the following initialization is necessary at the start of each AIPS++ or Glish command-line session:
- include 'vlafiller.g'; - include 'ms.g' - include 'flagger.g'; - include 'calibrater.g' - includ 'imager.g' - include 'image.g';
This is not required when working entirely from the GUI interface, in which case the tools are loaded automtically when they are selected.
Visibility data in AIPS++ are stored in AIPS++ Tables known as MeasurementSets (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 definition is described fully in AIPS++ Note 229.
AIPS++ data files, including MeasurementSets, are written into the current working directory by default, with each AIPS++ Table represented as a separate sub-directory. MS names therefore need only comply with UNIX file or directory naming conventions. The MS being filled can be created in a location other than the default directory by specifying an explicit path in the MS name.
There two likely paths for importing VLA data into AIPS++: directly from a VLA archive or distribution tape, or from a UVFITS disk file generated by the AIPS task FITTP. Both methods are described here.
- include 'vlafiller.g'
The vlafiller tool can be used to read data from disk or tape. Tape drives can be located on the local system or on a remote system accessible on the network. Use of a remote tape drive requires knowledge of the network host name (e.g. gibbon.cv.nrao.edu, or gibbon), and the tape device name (e.g. /dev/rmt/0ln). This can often be determined from the physical drive itself, or can otherwise be determined by logging on to the remote host and using a UNIX utility such as sysinfo or the UNIX command ls -lut /dev/r[ms]t* to list the device name.
Remote tapes need not be explictly mounted from AIPS++ as the vlafiller itself will be run on the remote system directly. If you encounter problems in using the filler on a remote host, make sure that the AIPS++ initialization command source /aips++/release/aipsinit.csh is included in your .cshrc file (or the sh equivalent).
To run the vlafiller using a local tape drive:
# Load project TESTT from tape to AIPS++ file cals00feb.ms:
- vlafillerfromtape(msname='cals00feb.ms', device='/dev/rmt/0ln',
project='TESTT');
To run the vlafiller using a remote tape drive:
# Load project TESTT from remote tape on gibbon to a
# local AIPS++ file cals00feb.ms:
- vlafillerfromtape(msname='/mydirectory/cals00feb.ms', device='/dev/rmt/0ln',
project='TESTT', host='gibbon');
Note that it is necessary to specify a full path name (recognizable by the remote host) for the output file when filling from a remote tape drive.
The vlafiller allows selection on frequency band, calibrator name, frequency range, project, time, subarray and source, and supports concatenation to an existing uv-data set.
It is usually desirable to obtain a summary of the dataset's properties after filling. This can be done by creating an ms tool, and running its summary function as follows:
# Create an ms tool:
- m := ms('cals00feb.ms')
# Obtain the summary:
- m.summary(verbose=T);
# Finish the ms tool:
- m.done();
Data which have been exported to a disk file from AIPS using the task FITTP can be read into AIPS++ using the general UVFITS filler ms, which itself creates an ms tool. This is achieved as follows:
# Create an AIPS++ uv-data file from a UVFITS file:
- m:=fitstoms(fitsfile='/aips++/data/demo/CALS00FEB.fits',
msfile='cals00feb.ms');
# Obtain a summary of the dataset:
- m.summary(verbose=T);
# Finish the MS tool:
- m.done();
# Create an AIPS++ uv-data file from a UVFITS file:
- m:=fitstoms(fitsfile='/aips++/data/demo/NGC5921.fits',
msfile='ngc5921.ms');
# Obtain a summary of the dataset:
- m.summary(verbose=T);
# Finish the MS tool:
- m.done();
Once the data are loaded as an AIPS++ MS, 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 your working directory. The catalog tool can be started and viewed from Glish as:
- dc.gui()
Here, dc is the default catalog tool. It can also be accessed via the toolmanager interface from the Tools in Use selection. Selecting an individual MS listed in the catalog GUI and pressing the View button 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 count from zero (not one) when viewed at the low level provided by the tablebrowser. However, all AIPS++ sythesis module tools hide this artifact of c++ and will present selection on indices counting from one. The ms.summary function used above is useful for determining the (one-based) indices used in the synthesis module tools.
The primary graphical uv-data visualization tool in AIPS++ is msplot. It can be started from Glish directly as:
# Initialize the msplot tool:
- include 'msplot.g'
# Create an msplot tool (editing capabilities turned on):
- mp:=msplot('cals00feb.ms', edit=T); # or 'ngc5921.ms'
It can also be launched from the Toolmanager via the general 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 primary tool for command-based editing is flagger, which allows several directed flagging operations. These include flagging in specific time intervals, on specific MS indices, such as FIELD_ID, or for specific polarizations. The flagger tool can be started from the Toolmanager via the synthesis package, and flagger module, or directly from Glish as:
- fg:= flagger('mydata.ms');
# Create a flagger tool for a given MS
Typical VLA flagging commands include:
# Flag all autocorrelations:
- fg.auto()
# Flag the first 10 sec of each scan (defined by gaps of 60s):
- fg.quack(scaninterval='60s', delat='10s');
# Flag all data involving antenna 5 in the specified time range:
- fg.setantennas(ants=5);
- fg.timerange(starttime="24-FEB-2000/10:00:05",
endtime="24-FEB-2000/10:09:37");
A full range of command-based flagging operations are available, as described fully in the flagger documentation.
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 flagging 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 flags, reversal of flags already applied, and other related utility operations. Raster displays of uv-data are possible also, and are recommended for VLA data editing.
Automated editing methods are available in the flagger and ms tools. These include clipping, median flagging and related methods. For example, to clip data outside of a specified range:
- 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');
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 (bounded at the Fourier transform operation) 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 in Getting Results. To orient the user, a brief summary is given here.
In the typical VLA observation, a number of calibrators of known structure (usually points, and at least one of known flux density) are observed in addition to the scientific target in order to obtain a measure of the systematic properties of the instrument (and perhaps the medium above it). The effects of these systematic properties are then removed from the data for the scientific target. This technique is often referred to as cross-calibration and it success depends implicitly on the relationship between three aspects of a radio astronomy observation: the true brightness distribution of the observed radio sources, the formal response of the instrument and its variety of imperfections of varying importance, and the data that is actually recorded. The Measurement Equation is a formally complete and explicit description of this relationship. For standard VLA continuum and spectral line observations, the unknown properties of the instrument can usually be handled entirely in the visibility plane (the ``left-hand'' side of the ME) in the calibrater tool. In simple terms, the Fourier transform of the source brightness distribution (as given by the formal theory of the instrument) is multiplied by a number of factorable, feed-based "gains" to yield the observed data. Knowledge of the observed data and brightness distribution models for the calibrators yield the gains, which may be applied to the data for the scientific target. The corrected data for the scientific target is then Fourier transformed and deconvoled to obtain its brightness distribution model. The method of self-calibration extends this process to iterative alternate improvement of the gains and brightness distribution model of a single source, with additional dependence on the fact that the true source brightness distribution is positive (in Stokes I) and of limited spatial extent.
As noted above, AIPS++ uv-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 VLA 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 or field model (MODEL_DATA). The latter two columns are filled in by the calibrater and imager tools. The actual calibration information itself is stored in separate calibration tables. The observed DATA column does not change during reduction, but the related columns do. When plotting the data using a tool such as msplot, the data sub-types (observed, corrected and model) can be selected individually.
In the following sections, the means of solving for and applying (in the visibility plane) the standard (instrumental/atmospheric) complex gain (G), the polarization leakage (D), and the relative bandpass (B) (all in the visibility plane) are described for typical VLA datasets. The individual calibration components in the ME (which are expressed formally as Jones matrices to fully describe their polarization properties) are described in detail 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 
and its gradients 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 available in the imager
tool.
For the VLA, the default source models are customarily point sources defined by the Baars or Perley-Taylor flux density scales, or point sources of unit flux density if the flux density is unknown. The first time that an imager tool is constructed for a given MS, the MODEL_DATA column is created and initialized with unit point source flux density visibilities for all sources. The imager.setjy function can then be used to set the proper flux density for flux density calibrators. For sources among recognized flux density calibrators [(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)], the Perley-Taylor (1995) flux densities will be calculated, transformed (trivially) and written to the MODEL_DATA column. For other sources, the user can specify an arbitrary Stokes (I, Q, U, V), or unit flux density (Stokes I) will be used.
Set 3C286 (field 2) model according to Perley-Taylor (1995):
# Make imager tool:
- include 'imager.g'
- imagr:=imager('cals00feb.ms');
# Initialize source model for 1331+305 (3C286):
- imagr.setjy(fieldid=2);
Set 3C286 (field 1) model according to user-specified value:
# Make imager tool:
- include 'imager.g'
- imagr:=imager('ngc5921.ms');
# Set 1331+305 (3C286) model:
- imagr.setjy(fieldid=1, fluxdensity=[14.8009,0,0,0]);
At this point, the source model and observed data exist for all calibrators at all uv-points, and it is possible to begin solving for the gains.
It is almost invariably the case that the instrumental and atmospheric amplitude and phase gains (hereafter ``complex gains'') are the dominant calibration error in VLA data. Therefore we solve for this factor first. In the ME formalism, in fact, two factors (of nearly identical functional form, but at different points in the ME) are available to describe the net complex gain; these are the electronic gain (G) and the atmospheric gain (T), the latter of which has no polarization dependence. The practical advantages this factorization enables are described in the Synthesis Calibration chapter of Getting Results in AIPS++. For the standard cross-calibration scheme described here, only one term (G) will be used to describe the net complex gain (and so it will include atmospheric effects). Note that one complex gain value will be determined for each spectral window, polarization, and antenna in the MS. For spectral line observations, the gains' function of frequency channel within a VLA IF will be calibrated later (B).
To solve for G:
# Start calibrater tool:
- cal:=calibrater('cals00feb.ms');
# Select calibrator sources (in this case, all sources)
- cal.setdata('FIELD_ID IN [1:10] && SPECTRAL_WINDOW_ID==1');
# Apply parallactic angle correction (for polarization only):
- cal.setapply(type='P', t=5.0);
# Arrange to solve for G on 60.0 second timescale
- cal.setsolve(type='G', t=60.0, refant=3, table='cals00feb.gcal1');
# Obtain the solution:
- cal.solve();
# Examine solutions:
- cal.plotcal(tablename='cals00feb.gcal1');
Note that all calibrater.setapply and calibrater.setsolve executions prior to the solve execution will be in effect when obtaining the solutions. Thus, the parallactic angle correction, P, is applied to the data before the G solutions are obtained. Note that the parallactic angle correction is only required for polarimetry reduction, and may be omitted otherwise (in which case the parallactic angle variation is irrelevant and subsumed by the G term). Note that data imported from AIPS may have already had the parallactic angle correction applied using the AIPS task CLCOR. If so, the calibrater.setapply step should be omitted. (This is rarely, if ever, necessary for VLA data, but it is the user's responsibility to keep track of this.)
# Start calibrater tool:
- include 'calibrater.g'
- cal:=calibrater('ngc5921.ms');
# Select data for calibrators (1 & 2)
- cal.setdata(msselect='FIELD_ID <= 2');
# Arrange to solve for G on 300-second timescale:
- cal.setsolve(type='G', t=300.0, preavg=300.0, refant=15,
table='ngc5921.gcal');
# Determine solutions:
- cal.solve();
# Examine solutions:
- cal.plotcal(tablename='ngc5921.gcal')
In this case, we have omitted application of 'P' since there is no polarization data.
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. This scaling is achieved using the calibrater.fluxscale function as follows:
# Transfer flux density scale from 1331+305 (3C286) to others:
- cal.fluxscale(tablein='cals00feb.gcal', tableout='cals00feb.fluxcal',
reference='1331+305',
transfer=['0927+390', '0713+438', '0854+201', '1310+323',
'1337-129', '1252-336', '1534-354', '1743-038',
'1751+096']);
# Examine solutions:
- cal.plotcal(tablename='ngc5921.fluxcal');
# Finish calibrater tool (for now):
- cal.done()
# Transfer flux density scale from 3c286 to others:
- cal.fluxscale(tablein='ngc5921.gcal',
tableout='ngc5921.fluxcal',
reference='1331+30500002',
transfer=['1445+09900002']);
# Examine solutions:
- cal.plotcal(tablename='ngc5921.fluxcal');
# Finish calibrater tool:
- cal.done();
In both cases, the amplitude scale is set with 1331+305 (3C286) as
the reference. The output calibration table contains the scaled G
calibration factors for all of the transfer sources, and will be used
going forward in the reduction. Proper scaling is indicated by
smoothly varying (ideally constant) amplitude gains for each antenna,
despite source changes. At high frequencies (
> 15 GHz), this
scaling step may be biased by variation of the atmospheric
contribution to the system noise temperature. Methods for treating
this case properly will be available in the near future.
In the case of spectral line data, it is necessary to remove gain fluctuations that are a function of frequency. This is achieved by using the calibrater tool to solve for B calibration, as follows:
# Start calibrater tool:
- cal:=calibrater('ngc5921.ms');
# Select data bandpass calibrator (1331+305):
- cal.setdata(msselect='FIELD_ID==1');
# Apply G solutions:
- cal.setapply(type='G', t=0.0, table='ngc5921.gcalout');
# Arrange to solve for bandpass
- cal.setsolve(type='B', t=86400.0, refant=15,
table='ngc5921.bcal');
# Determine solutions:
- cal.solve();
# Examine solutions:
- cal.plotcal(tablename='ngc5921.bcal');
# Finish calibrater tool:
- cal.done();
In this example, the G solutions determined above are applied, and relative B solutions are obtained on a timescale of one day, i.e., one solution per antenna per day.
Before imaging, the data need to be corrected by applying the calibrations determined above. In the context of applying calibration, the calibrater.setdata function selects which data are to be corrected, and the calibrater.setapply function selects which solutions are to be applied. The setapply function may be run several times (in any order) to specify simultaneous application of calibration tables of different Jones types. In the case, the different tables will be applied in the appropriate order according to the underlying ME formalism.
These examples fill the CORRECTED_DATA column, which will be used directly in imaging:
# Start calibrater tool:
- cal:=calibrater('cals00feb.ms');
# Select sources to which calibration is applied (all):
- cal.setdata(msselect='FIELD_ID IN [1:10] && SPECTRAL_WINDOW_ID==1');
# Arrange to apply parallactic angle correction:
- cal.setapply(type='P', t=5.0);
# Arrange to apply G solutions:
- cal.setapply(type='G', t=0.0, table='cals00feb.fluxcal');
# Apply solutions and write CORRECTED_DATA column:
- cal.correct();
# Finish calibrater tool:
- cal.done();
# Start calibrater tool:
- cal:=calibrater('ngc5921.ms');
# Select data to which calibration will be applied:
- cal.setdata(msselect='FIELD_ID IN [1:3]');
# Arrange to apply G solutions:
- cal.setapply(type='G', t=0.0, table='ngc5921.fluxcal',
select='FIELD_ID==2');
# Arrange to apply B solutions:
- cal.setapply(type='B', t=0.0, table='ngc5921.bcal');
# Apply solutions and write CORRECTED_DATA column:
- cal.correct();
# Finish calibrater tool:
- cal.done();
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 in AIPS++ and in the User Reference Manual. Here, we describe simple single-field imaging/deconvolution for continuum and spectral line datasets.
Forming an image requires selection of the data for the desired field and spectral window/channel(s) via the setdata function, and selection of the imaging coordinates via the setimage function. The selection of imaging coordinates includes not only the image scale (size and number of pixels), but also field, polarization, and spectral selections. The unique selection of field and spectral window/channel in both of these functions enables flexible control of the image process, e.g., to form an integrated continuum image from multi-channel data, or to form an image from data obtained for one field at the position of another (as in mosaicing). This dual selection mechanism reflects the distinction between the uv- and image- domains, as prescribed by the ME.
The image-plane scaling parameters may be selected automatically using advise. For more information on any imaging parameters or options, including alternative deconvolution methods, please consult the full imager documentation.
For the continuum dataset, we will first make a single-field, full-polarization (I, Q, U, & V), continuum VLA image of a calibrator (3C286), since this is useful for checking the initial calibration and for obtaining a full-Stokes model for polarization calibration (described below). This image is made with the Clark CLEAN as follows:
# Create an imager tool
- imagr:= imager('cals00feb.ms');
# Select data for field 2, continuum, and spectral window 1
- imagr.setdata(fieldid=2, mode='none', nchan=1, start=1, step=1,
spwid=1);
# Set the image coords for field 2, IQUV, spectral window 1
- imagr.setimage(nx=256, ny=256, cellx='0.4arcsec', celly='0.4arcsec',
stokes='IQUV', fieldid=2, spwid=1, start=1, step=1,
nchan=1, mode='none');
# Create an empty model image in these coordinates
- imagr.make('1331-g1.model');
# Image, and deconvolve using the Clark CLEAN.
- imagr.clean(algorithm='clark', niter=5000, gain=0.1,
threshold='0.0Jy', model='1331-g1.model',
image='1331-g1.restored',
residual='1331-g1.residual');
# Generate restored image over full image:
- imagr.restore()
# Finish the imager tool:
- imagr.done();
Note that the clean component model is written as an image file. If necessary or desired, cleaning may be restarted by specifying an existing (non-zero) model image, which will be subtracted before the deconvolution is resumed. This provides a simple means of monitoring progress of the deconvolution. (Note that the setdata and setimage functions need not be run again when restarting a clean.)
For the spectral line dataset, total intensity images must be generated for NGC5921 for each spectral channel to form a spectral line cube. This is achieved by specifying mode='channel' in setdata and setimage in the imager tool.
# Select data for field 3, all channels in spectral window 1
- imagr.setdata(fieldid=3, mode='channel', nchan=63, start=1, step=1,
spwid=1);
# Select image of field 3, I, all channels in spectral window 1
- imagr.setimage(nx=256, ny=256, cellx='10arcsec', celly='10arcsec',
stokes='I', fieldid=3, spwid=1, start=1, step=1, nchan=63,
mode='channel');
- # Select default imaging weighting (uniform)
- imagr.weight();
# Form empty model image:
- imagr.make('ngc5921.model');
# Image and deconvolve using the Clark CLEAN
- imagr.clean(algorithm='clark', niter=3000, gain=0.1, threshold='0.0Jy',
model='ngc5921.model', image='ngc5921.restored',
residual='ngc5921.residual');
# Generate restored image over full image plane
- imagr.restore(model='ngc5921.model',
image='ngc5921.restored2');
# Finish this imager tool:
- imagr.done();
In both cases the imager.restore has been used to transform the final residual data to the full image plane and add the restored model. The final ``restored'' image generated by the imager.clean is only restored over the region where clean components were searched for, in this case, the inner quarter.
The images generated in the last section can be view, analyzed, and rendered for publication using the AIPS++ viewer and image tools.
When launched, the viewer displays the sky plane of the image by default. Other axes may be selected and explored using the tapedeck buttons to the right of the display. To the left of the display, buttons to actuate mouse-driven zoom, colormap and brightness/constrast fiddling, and region and position setting are available. The region- and position-setting tools enable access to image tool functions that provide pixel value properities and statistics. Consult the documentation for these tools for more information.
# Create an image tool
- im:= image('1331-g1.restored');
# View the restored image
- im.view()
# Close this image, keeping image tool
- im.close()
# View the residual image
- im.open('1331-g1.residual');
# Close the image tool
- im.done()
Spectral line images can be viewed similarly.
Extensive image analysis tools are available in AIPS++ , as implemented in the image module. These include versatile image algebra, concatenation, convolution, moment calculation and some source component fitting capabilities. These are fully described in the image documentation.
As an example, the peak Stokes (I, Q, U, V) can be derived from an image in a specific region using the following script:
# Create an image tool
- im:= image('1331-g1.model');
# Array to hold min, max
- pmax:= array(0.0, 4);
# Loop over I,Q.U,V
- for (plane in 1:4) {
# Define image region
pregion:= drm.box([56,56,plane], [192,192,plane]);
# Retrieve statistics
im.statistics(statsout=pstats, region=pregion);
# Extract max(abs)
if (abs(pstats.min) > abs(pstats.max)) {
pmax[plane]:= pstats.min;
} else {
pmax[plane]:= pstats.max;
};
};
# Finish image tool
- im.done();
With the peak flux densities for 3C286 now available in the glish variable pmax, it is possible to proceed with instrumental polarization calibration (see below).
Another example of image analysis comes from the common need to subtract the continuum emission from a spectral line observation. In the image plane, this is done by first isolating those image channels which contain only continuum, and then subtracting their average from the full spectral line cube (note the image is 4D, with axes RA, Decl, Stokes, and channel):
# Start an image tool containing the restored image:
- restim:=image('ngc5921.restored2');
# Mask bad channels (first few, last few):
- restim.set(pixelmask=F,region=drm.box([1,1,1,1],[256,256,1,5]));
- restim.set(pixelmask=F,region=drm.box([1,1,1,60],[256,256,1,63]));
# Copy restored image to final image for in-place continuum subtraction:
- csubim:=imagefromimage(outfile='ngc5921.final',
infile='ngc5921.restored2');
# Define regions of image which is continuum (channels 6-9,55-58):
- cregion1:=drm.box([1,1,1,6],[256,256,1,9]);
- cont1:=restim.subimage(outfile='cont1.im',region=cregion1);
- cregion2:=drm.box([1,1,1,55],[256,256,1,58]);
- cont2:=restim.subimage(outfile='cont2.im',region=cregion2);
# Form a single image containing all continuum channels:
# (this is a virtual image referencing cont1.im & cont2.im)
- cont:=imageconcat(infiles="cont1.im cont2.im",relax=T);
# Form integrated continuum as moment=-1 of continuum channels
- cont.moments(outfile='ngc5921.cont',moments=-1,axis=4);
# Delete intermediate images/tools:
- cont.done();
- cont1.delete(); cont1.done();
- cont2.delete(); cont2.done();
# Subtract continuum from cube (in place):
- csubim.calc('$csubim - ngc5921.cont')
# View final image
- csubim.view();
This glish sequence has made use of image regions through the default regionmanager (drm), the image.moments function, and the image.calc function. Consult the image tool documentation for more information on the use of these and other image analysis functions.
With a model for the linear polarization of 3C286, and the G calibration obtained above, it is possible to proceed with instrumental polarization for the continuum dataset as follows:
# Use Stokes image peaks to set I, Q, & U models (assume V=0):
- imagr.setjy(fieldid=2, fluxdensity=[pmax[1], pmax[2], pmax[3], 0.0]);
# Start calibrater tool:
- cal:= calibrater('cals00feb.ms');
# Select data for polarization calibrater:
- cal.setdata(msselect='FIELD_ID==2 && SPECTRAL_WINDOW_ID==1');
# Apply parallactic angle:
- cal.setapply(type='P', t=5.0);
# Apply G solutions
- cal.setapply(type='G', t=0.0, table='cals00feb.gcalout1');
# Arrange to solve for D on day-long timescale, averaging data
# within the solution to no more than 600 sec per chunk:
- cal.setsolve(type='D', t=86400.0, preavg=600.0, table='cals00feb.dcal1');
# Obtain the D solutions
cal.solve();
Now, the data can be corrected for both G and D as follows:
# Select fields to which calibation is to be applied: - cal.setdata(msselect='FIELD_ID IN [1:10] && SPECTRAL_WINDOW_ID==1'); # Apply parallactic angle correction: - cal.setapply(type='P', t=5.0); # Apply G solutions: - cal.setapply(type='G', t=0.0, table='cals00feb.gcalout1'); # Apply D solutions: - cal.setapply(type='D', t=0.0, table='cals00feb.dcal1'); # Correct the data: - cal.correct(); # Finish calibrater tool: - cal.done();
All fields are now calibrated and can be imaged in full polarization as described above.
Hamaker, J.P., Bregman, J.D., and Sault, R.J., 1996a, A&AS, 117, 137.
Sault, R.J., Hamaker, J.P., and Bregman, J.D., 1996b, A&AS, 117, 149.