Getting Started | Documentation | Glish | Learn More | Programming | Contact Us |
Version 1.8 Build 235 |
|
Kumar Golap and Tim Cornwell
The relationship between the visibility domain to the image plane is given by the Van Cittert-Zernicke theorem. This can be reduced to a Fourier relationship for planar arrays and also for small field of views for non-planar arrays (Thompson, Moran and Swenson). When this happens the Fourier inversion of observed visibilities gives the true image convolved with the point spread function (PSF). Both the true image and the PSF are functions of 2 variables only (l,m). Deconvolving is conceptually simple.
When visibilities are sampled by a non-planar array (u, v and w components), the relationship between them and the image plane can no longer be represented by a 2-D Fourier transform. The relationship is now a 3-D Fourier transform. Of course our image of the sky is truly 2-D and when deconvolving it, one should ideally use a different PSF for each point in the image. This can be very difficult to achieve if the image has a few million pixels.
It becomes worse when one has to consider that the field of view for deconvolution is much more than the field of view of interest (at most the primary beam of the antennas). Deconvolution may need to consider the sidelobes and at low frequencies the whole visible sky may have to be considered, especially if strong sources like Cygnus-A or Cas-A are above the horizon. The PSF in the sidelobes suffers both distortion from the 'w' term and the amplitude and gain difference difference of the sidelobes of each antenna. Wide-field imaging is different from mosaicing as the whole observation can be made with a single pointing. Wide-field imaging is needed when the term w varies significantly over the field of view. Mosaicing is the technique used to combine visibilities (or images) obtained from different pointings; in effect it increases the field of view beyond the primary beam pattern.
An important question is how do you know that you need wide-field imaging techniques? If your array is planar (e.g. an E-W array) you do not need them. Otherwise, one simple way is to make an image via basic 2-D Fourier techniqes. If point-source responses on the edge of the image have the same shape as those in the center then it is not needed.
A simple criterion (due to Barry Clark) is that:
(Field of view in radians) * (Field of view in beams) > 1
There are in fact a number of different possible approaches to solving the wide-field imaging problem. These are reviewed by Cornwell and Perley. Here we give a brief account:
Calculation of the facet planes can be handled in two ways:
The basic tools for calibration and imaging are Calibrater and Imager. Combining these is a tool called Dragon. This tool allows the user to loop through many cycles of wide-field imaging and self-calibration successively with little or no inputs from the user in the process.
Imager handles wide-field imaging via the faceting described above. For various reasons, one may like to change the number of facets being used when processing an image. For example, one may wish to have few facets when deconvolving bright sources and have many more facets to have deep CLEANing.
During the self-calibration part of the imaging process (also sometimes called difference mapping), the flux levels where phase self-cal and amplitude self are needed can be specified.
The imaging process can be helped by specifying the bright sources via a Componentlist tool.
If there are bright sources outside the field of interest, their sidelobes can be deconvolved from the image of interest by specifying an outlier region around that source. Inherently there are no limits on the number of facets or number of outlier fields to be specified.
Imager and Dragon have a function which advises on the optimal number of facets, tangent point of image, to be used for given desired signal-to-noise ratio and size of image to be made.
The number of facets needed to image a given field of view depends on several factors. Among them are the uv distance, the amount of amplitude loss allowable, and the positional error allowable. The number of facets used ultimately depends on what the user wants. If flux information only is needed the number of facets can be much less than if, say, very accurate positions are needed.
The advise function gives both the number of facets needed to get the accurate positions or just accurate flux measurements.
For accurate flux, a gradient in the w term phase is not important. The advise function's least squares fit to uvw gives the best fitting plane, dispersion in w' gives the beam size in w. Max amplitude loss plus simple geometry gives the facet size. But for accurate positioning the facet size should be such that even a gradient in the w term is not important for a given amplitude loss.
The advice for accurate flux may sometimes be too small to the extent of causing artifacts. The positional error might be large enough that sources close to the edge of the facets may end up appearing on 2 facets. In general if one requires accurate flux values, then it would be advisable to use a few more facets than the minimum suggested by advise. Make sure that there is not any source splitting over the image.
If very accurate positions are needed then either one has to image with the large number of facets suggested or try an alternative solution. That alternative solution is to image with a number of facets which is practical then get the best fit the position of each of the component from the UV-data. The fitting of components in the UV data is yet to be implemented in AIPS++.
Deconvolution in the wide-field context is presently available only via CLEAN and implemented only for the Hogbom and Clark CLEAN variants (see Imager.clean and Dragon.image).
If no self-calibration is needed, either Imager can be used directly or else Dragon can be used but with the arguments levels and amplitudelevel in Dragon.image left at '0Jy'.
If you do need self-calibration, then it is set to start at an amplitude level where it is expected that the phase or amplitude errors start to become significant. Usually the user will do phase self-calibration only initially, followed later by phase and amplitude self-calibration at a much lower level of residuals.
drgn := dragon('data.ms') drgn.setimage(...) # Setup image parameters # drgn.image(levels='1Jy 0.3Jy 0.1Jy',amplitudelevel='0.2Jy', timescales='200s 200s 100s', niter=50000, gain=0.05, threshold='0.05Jy', plot=F, display=F, algorithm='wfclark')
In this example, CLEANing will be done and self-calibration will be executed at residual levels of 1Jy, 0.3Jy and 0.1Jy. As amplitudelevel is set to 0.2 Jy, the calibration steps at 1Jy and 0.3Jy will be phase calibration and correction only. However, the step at 0.1Jy will be a phase and amplitude phase calibration and correction. The timescales argument defines the averaging time period for the calibration solutions.
Deconvolution converges faster when the algorithm is told where to expect flux. This is especially the case for CLEAN in the presence of bright sources. A description of how to handle this within Imager can be found in the basic Imaging chapter of Getting Results.
In Dragon, the argument maskmodification in Dragon.image can be used to define masks at each self-calibration step. The argument can take the following: 'none' for no masking, 'auto' for automatic masking and 'interactive' for user defined ones.
When maskmodification='auto', the restored image at the level CLEAN reached is used to make a threshold mask. All regions which are above the flux limit defined in levels is set to be used for searching for CLEAN components.
When maskmodification='interactive', a Viewer display panel is launched when CLEAN reaches the flux level defined and the user can define CLEAN regions interactively.
Outlier fields are regions outside the main image where there are sources that whose sidelobes fall into the field of view. Such regions would be around bright radio sources such as Cyg-A, Vir-A, the Sun etc. Such sources have to be included in the CLEAN and self-calibration loops so as to achieve good dynamic range. A special outlier field is the pole regions where interference tend to pile up in rotation synthesis arrays. This can be treated as an outlier field and the effect of such interference can thus be somewhat reduced from field of view.
The position of the center of the outlier fields can be given manually or can be looked up from a catalog by using the Measures module.
First let us see the steps that one would undertake to handle outlier fields ddirectly in Imager:
dir1:=dm.direction('J2000', '19h00m00', '50d00m00') # Outlier 1 imgr.setimage(nx=200,ny=200,cellx='10arcsec', celly='10arcsec', doshift=T, phasecenter=dir1, mode='mfs', spwid=[1,2], facets=1) imgr.make('outlier1'); # dir2:=dm.direction('J2000', '19h05m00', '55d00m00') # Outlier 2 imgr.setimage(nx=200,ny=200,cellx='10arcsec', celly='10arcsec', doshift=T, phasecenter=dir2, mode='mfs', spwid=[1,2], facets=1) imgr.make('outlier2'); # # Main field imgr.setimage(nx=2000,ny=2000,cellx='10arcsec', celly='10arcsec', doshift=F, mode='mfs', spwid=[1,2] facets=5) imgr.make('mainfield') # Deconvolve all together imgr.clean(algorithm='wfclark', niter=1000, gain=0.1, model="mainfield outlier1 outlier2", image="mainfield.restored outlier1.restored outlier2.restored", residual="mainfield.residual outlier1.residual outlier2.residual")
In Dragon to image with outlier fields one would proceed as follows:
# # Main image drgn.setimage(name='mainfield' nx=2000,ny=2000,cellx='10arcsec', celly='10arcsec', doshift=F, mode='mfs', spwid=[1,2], facets=5) # # Outlier 1 dir1:=dm.direction('J2000', '19h00m00', '50d00m00') drgn.setoutlier(name='outlier1',nx=200,ny=200,cellx='10arcsec', celly='10arcsec', doshift=T, phasecenter=dir1, mode='mfs', spwid=[1,2], facets=1) # # Outlier 2 dir2:=dm.direction('J2000', '19h05m00', '55d00m00') drgn.setoutlier(name='outlier2',nx=200,ny=200,cellx='10arcsec', celly='10arcsec', doshift=T, phasecenter=dir2, mode='mfs', spwid=[1,2], facets=1) # # Image and deconvolve drgn.image(levels='1Jy 0.3Jy 0.1Jy',amplitudelevel='0.2Jy', timescales='200s 200s 100s', niter=50000, gain=0.05, threshold='0.05Jy')
The process of imaging can be helped and made to converge much faster if a model is passed to Dragon/Imager to start with. The models can either be position, flux and shape of sources (via a Componentlist tool) or a model image (or both).
Model images are usually images made of the region but with poorer sensitivity and dynamic range. The model(s) can be put in the model data column of the MeasurementSet by using the Imager.ft function. Otherwise these can be passed as the model argument for the Dragon.image or Imager.clean function.
A Componentlist tool, as its name suggests, allows one to make a list of source models. Componentlists can be constructed from from standard catalogues (e.g the NVSS, VLA FIRST etc.) via the asciitocomponentlist constructor. Using the catalogues can be very helpful to Dragon or Imager (helps convergence especially in complex regions).
A Componentlist can also be created from an image using the interactive Imagefitter tool. This tool allows the user to interactively select regions where there are components and estimate its parameters.
Here is an example of a Glish script which was used to image the Coma cluster at 74MHz with data from the VLA in the B and C configurations.
include 'dragon.g' myms := fitstoms('coma.ms','COMA-4CB-CUT.FITS',T,F) # Convert from FITS to MS myms.done() # drgn := dragon('coma.ms') # Create Dragon tool drgn.setimage(name='coma',nx=1800, ny=1800, # Set imaging parameters cellx='30arcsec', celly='30arcsec', doshift=F, phasecenter=dm.direction('J2000','0deg','0deg'), mode='mfs', facets=25) # drgn.setoutlier(name='CenA' , nx=400, ny=400, # CenA outlier field cellx='37arcsec', celly='37arcsec' , doshift=T, mode='mfs', phasecenter=dm.direction('J2000','201.3deg','-42.6deg'), nchan=1, start=1, step=1, spwid=1, fieldid=1) # #setting an outlier field towards Virgo-A; note for direction we have # made use of the inbuilt radio catalog drgn.setoutlier(name='VirA' , nx=400, ny=400, # VirgoA outlier from builtin catalog cellx='37arcsec', celly='37arcsec' , doshift=T, mode='mfs', phasecenter=dm.source('VIRGOA'), nchan=1, start=1, step=1, spwid=1, fieldid=1) # drgn.uvrange(uvmin=0, uvmax=10000000) # uv range drgn.weight(type="uniform") # uniform weighting drgn.setoptions(padding=1.5, cache=0) # padding for the Fourier Transform # drgn.image(levels='1Jy 0.3Jy 0.1Jy', # image and selfcalibrate amplitudelevel='0.2Jy', timescales='200s 200s 100s', niter=50000, gain=0.05, threshold='0.05Jy', plot=F, display=F) # drgn.done()
This problem is detected when all the sources away from the center of the image are broadened. This should not be confused with the effect of bandwidth smearing (decorrelation). Another effect of too few facets is that sources near the edge of facets appear as double. This is beacause the positional error is large enough to cause sources near the edge of the facets appear on 2 facets.
The effect of non-isoplanatism is that the phase/gain of the longer baselines gets distorted while the shorter ones are not affected. In the image one would notice broadening of all known point sources and they will have similar structures. This will affect all the sources in the field.
The solution to this is to image with self-calibration corrections starting at much higher flux levels. One should be aware of the time scale of the non-isoplanar changes. If one passes a too low time scale to self-cal, the signal-to-noise ratio of the self-calibration estimates will be worse. If it is much larger than necessary then it's equivalent to no self-calibration.
The capability to correct for directional dependent part of phase and amplitude errors is going to be implemented soon. This should be able to take care of assymetric primary beam and effect of the atmosphere/ionosphere which are direction dependent.
The problem is most probably the effect of ringing on the edge of the facets which can be reduced by padding. So in the Dragon.setoptions or Imager.setoptions function the padding factor can be given a value of more than 1 (usually 1.2 to 1.5) depending on how severe the edge get picked up by clean. The more the padding the longer will FFT take.