INVERT INVERT is the program used to make maps from calibrated visibility data by Fourier inversion (aperture synthesis). Both MERGE-formats 1 and 2 are accepted. Example $ INVERT INPUT "[TJP]PHASED.TST" MAPFILE "[SCRATCH]MAP.TST" BEAMFILE "[SCRATCH]BEAM.TST" PLOTFILE "INVERT.PLT/PR" LISTFILE "INVERT.LIS" IDIM = 128 MAXP = 2000 PRINTMAP PLOTBEAM = 6 PLOTMAP XYINT=0.35 MAPSIZE=100 FLUX=3.2 SHIFT 10,10 / PARAMETERS 1 FILES: 1. INPUT="file1","file2" [default "FROM"] 2. MAXP = n [default 5000] 3. MAPFILE="filename" [default "MAP"] 4. BEAMFILE="filename" [default "BEAM"] 5. PLOTFILE="device/type" [default "PLOTOUT"] 6. LISTFILE="filename" [default "invert.lis"] INPUT specifies the input file: this is a standard MERGE-format file (format 1 or 2), which must have valid phases. In normal hybrid mapping, it is the output file from the program AMPHI. For simulations, it is the output file from program FAKE. Up to two input files can be specified; the data from the two files are combined during gridding to make a single map. INVERT Page 2 MAXP specifies the maximum number of visibility data points (uv-points) that INVERT can handle. For large datasets, it will be necessary to increase MAXP from the default value. The maximum number of stations accepted is 40. If you do not want to use all the data to make the map, the unwanted data can be removed with VLBEDIT. The parameter UVCUT (see below) can also be used to reject long-baseline data. MAPFILE and BEAMFILE specify the (new) output files to contain the map and the beam respectively. Usually these files will be used as input to program CLEAN. PLOTFILE specifies the device and type of device to be used for plotting. Normally a hardcopy plot is generated, in which case PLOTFILE will be something like "INVERT.PLT/PRINTRONIX" or "INVERT.PLT/VPS". Any device supported by the PGPLOT graphics package can be used. LISTFILE specifies a file name for the listing file which records the progress of the program. The listing file can be printed or examined with an editor. 2 SWITCHES CONTROLLING PROGRAM ACTION: 1. NOMAP : if specified, do not calculate map 2. NOBEAM : if specified, do not calculate beam 3. FLUX = flux [default: do not include zero-baseline flux] The default action is to calculate both map and beam. In hybrid mapping, the beam usually remains the same from iteration to iteration, so time can be saved by not recalculating it. The parameter FLUX specifies the total flux-density of the source to be included at the origin of the u,v plane. If it is not specified, the zero-baseline flux is not included, i.e. it is given zero weight. When it is included, it is assigned a weight in the same way as each of the measured visibility points. 3 MAP SIZE AND SCALE: 1. IDIM = n [default = 256 ] 2. MAPSIZE = n [default = IDIM/2 ] 3. XYINT = xyint MILLIARCSEC | ARCSEC The size of the map and beam in gridpoints is fixed by the parameter IDIM: maps are IDIM x IDIM points square. IDIM must be a power of 2 INVERT Page 3 (16, 32, ..., 1024, 2048). This parameter should be as small as possible to minimize running time and memory requirements. Usually the dirty map made by INVERT will be cleaned, and CLEAN requires a beam array that is twice as big as the map array. For this reason, INVERT writes out only the central MAPSIZE gridpoints of the map, where MAPSIZE is usually one half of IDIM, although other values up to IDIM can be specified. The size of each cell of the map is specified with the XYINT parameter; the units of this parameter can be explicitly given as either MILLIARCSEC or ARCSEC, the default (if no units are specified) being MILLIARCSEC. The cellsize in the u,v-plane is determined by XYINT and by the arraysize (IDIM) as follows: UVINT = 1 / ( IDIM * XYINT ) where XYINT is in radians and UVINT in wavelengths. The array size for the u,v-plane is also IDIM, so the maximum baseline that can be coped with is UVMAX = UVINT*IDIM/2. The program inspects the data to determine the actual value of UVMAX, and evaluates the corresponding maximum allowed XYINT = 1/(2*UVMAX) radians; an error message is printed if the requested XYINT is larger than this. (Note: the actual maximum XYINT is slightly smaller than this formula implies, to allow for the effects of convolution.) The program also evaluates an "oversampling factor" = 1/(2*UVMAX*XYINT). If this factor is less than 1, spatial frequencies corresponding to UVMAX would be sampled at less than the Nyquist frequency on the map, and hence be attenuated; the program does not permit this. In practice the oversampling should be at least 2, and for CLEAN to run successfully it is sometimes advisable to have higher oversampling. Note that high oversampling (small XYINT) leads to large UVINT and possible undersampling of the u,v-plane. For this reason one should inspect the gray-scale plot of the u,v-plane to check that the elliptical tracks are reasonably well separated. For a given XYINT, UVINT can only be decreased by going to a larger array size (ie, by increasing IDIM). The default value of XYINT will give an oversampling factor of 2. 4 WEIGHTING AND TAPERING: 1. UVWEIGHT = "UNIFORM" | "NATURAL" | "RADIAL" [default = "UNIFORM"] 2. UVRAD = radius [default = 1.0] 3. WTPOWER = exponent for SNR weighting [default: 0.0] 4. UVTAPER = taper,baseline [default: no uv-taper] 5. UVCUT = cutoff-baseline [default: no cutoff] The weighting of data before convolution and gridding is controlled by these five parameters. With natural weighting (UVWEIGHT="N"), each measured visibility point is given equal weight in the transform. This generally gives rise to a very broad beamshape, because the measured INVERT Page 4 points are more densely crowded near the origin of the u,v-plane, and is not recommended. With uniform weighting (UVWEIGHT="U"), each point is weighted with the inverse of the number of points in a square of dimension UVRAD*UVINT. [This is the algorithm that was called "fast uniform weighting" in earlier versions of the program. The earlier scheme of calculating the distances between all pairs of points was very slow and has been abandoned.] UVRAD=2 is recommended with uniform weighting, although UVRAD=1 is the default. Note that UVRAD need not be an integer, but for reasons of memory space it may not be less than 1.0. Radial weighting (UVWEIGHT="R") is sometimes useful: each point is weighted as its radial distance from the origin of the u,v-plane. This compensates exactly for the density effect only if the u,v-tracks are circles around the u,v origin. Additional Gaussian tapering of the u,v-plane can be imposed with the UVTAPER parameter. The first value of UVTAPER is the fractional taper, and the second value is the baseline (in wavelengths) at which the gaussian attains this value; if only one value is given, it is the value of the taper function at the maximum baseline (UVMAX); eg if UVTAPER=0.3, the Gaussian taper function has a value 1 at the origin of the u,v-plane and a value 0.3 at the radius UVMAX. Tapering is sometimes used to get a cleaner beamshape - it reduces sidelobes at the expense of resolution, and it also increases the sensitivity to low surface brightness emission. If you want to make a map with a large beam, it is better to use a suitable UVTAPER than to simply specify a larger clean beam in CLEAN. The parameter WTPOWER allows the use of the error bars associated with the data to do weighting based on signal-to-noise ratio. After the weights for gridding are calculated as above (uniform, natural, or radial weighting), they are divided by the data error to the WTPOWER power. The default is 0 and leaves the weights unmodified. WTPOWER=2 causes the data to be weighted according to its SNR. For the best beam, uniform weighting and WTPOWER=0 should be used. For greatest sensitivity, natural weighting and WTPOWER=2 should be used. Intermediate values of WTPOWER might be useful in certain cases. The errors in typical VLB data have a tremendous range so, if full SNR weighting is used, the beam is likely to be very bad. It is not recommended that full SNR weighting be used until the calibration is very good - ie. the hybrid mapping iterations are nearly done. At that point, SNR weighting is likely to reduce the off-source noise significantly. Beware of short bits of very high SNR data (eg. Arecibo baselines) as they can badly distort the beam with full SNR weighting - it could be worth modifying the errors of such data to make them similar to that of other good data in the experiment. Note that natural weighting and WTPOWER=2 is equivalent to UVWTFN='NA' in AIPS. The parameter UVCUT is used to reject long-baseline data completely. The value of UVCUT is a baseline in wavelengths. Data from points in the u,v-plane outside a circle with this radius are ignored. INVERT Page 5 5 CONVOLUTION AND GRIDDING: 1. SIGMA = sigma [default 0.7] 2. PRINTCF : if specified, print the interpolation fn Gridding is the process by which the measured visibilities, at random points in the u,v-plane, are assigned to gridpoints in the u,v array. Essentially the measured points, regarded as a set of delta-functions, are convolved with an "interpolation function" to obtain a function which is continuous over the u,v-plane; this continuous function is then sampled at the gridpoints. The effect of the convolution is to multiply the map by the Fourier transform of the convolution function: INVERT removes this effect by multiplying by a "sensitivity correction function" that is simply the reciprocal of the Fourier transform of the convolution function. The effect of sampling is to introduce "aliasing": features which should appear outside the mapped area are translated into it. The interpolation function must be chosen to (a) minimize the numerical errors inherent in the sensitivity correction, and (b) reduce the response to aliases. These are conflicting requirements. The interpolation function used (at present) is a truncated Gaussian with FWHM = 2*SIGMA (in units of UVINT). The default SIGMA is chosen to give reasonable alias rejection without introducing large numerical errors. I do not recommend changing SIGMA, but if problems occur with aliasing, it may be necessary. An alternative solution to the alias problem is to increase IDIM. If the parameter PRINTCF is specified, the values of both the interpolation function and the sensitivity correction function are tabulated. 6 PLOTTING: 1. PRINTUV : if specified, plot the u,v-coverage (gray-scale) 2. PRINTMAP : if specified, plot the map (gray-scale) 3. PRINTBEAM : if specified, plot the beam (gray-scale) 4. PLOTMAP [=size] : if specified, plot the map (contours) 5. PLOTBEAM [=size] : if specified, plot the beam (contours) 6. CONTOURS=v1,v2,v3,v4,...,v20 The default is no plotting. Plotting, especially gray-scale plotting, uses a lot of cpu time; INVERT will run quite a bit faster if you do not request plots that you do not need. INVERT Page 6 The gray-scale and contour plots are sent to the device specified by parameter PLOTFILE (see above). Each gray-scale plot is scaled so that zero appears white (or the background color) and the maximum value appears black (or the default pen color); all negative values appear in the background color. PRINTUV produces a plot of the u,v-coverage after convolution and weighting: it represents the real part of the visibility function used for calculating the beam. Only one half of the u,v-plane is plotted. The imaginary part of the visibility function is zero for a beam. PRINTUV has no effect if NOBEAM is specified. If PLOTMAP (or PLOTBEAM) is specified with no value, the entire map (or beam) is plotted; note that the beam has twice the area of the map. If either parameter is given a value, a square area centered at the map center with side twice "size" is plotted. The units of "size" are the same as those of XYINT, ie milliarcsec by default. Up to 20 contour levels can be specified with CONTOURS. They are expressed as percentages of the maximum, and apply to both the map and the beam. The default contours are at -45,-35,...,-5,5,15,...,95 % of maximum. 7 SHIFTING AND ROTATING THE MAP: 1. SHIFT = dx,dy [default 0,0 ] 2. ROTATE = angle [default 0 ] SHIFT is used to change the phase-reference point (map center). A positive value for dx moves the map-center to the right (West) and a positive value for dy moves the map-center down (South); in other words, if both dx and dy are positive, the source will move up and to the left relative to the map frame. Shifting is done by multiplying each complex visibility value by a complex phase-shift factor. The units of dx and dy are the same as those of XYINT, ie milliarcsec by default. ROTATE rotates the map. A positive angle rotates the source anticlockwise relative to the map (or in other words, it rotates the map grid clockwise relative to the sky). Shifting is done before rotation, so the rotation is about the shifted map center. 8 SOURCE SUBTRACTION: 1. SUBSOURCE = flux,x,y A point source may be subtracted from the data before mapping. The values of the SUBSOURCE parameter are (1) flux density in Jy, (2) x-coordinate, (3) y-coordinate; the coordinates are in the same units as INVERT Page 7 XYINT (ie milliarcsec by default) and are referred to the phase center. x increases to the East and y increases to the North. Subtraction is done before shifting or rotating. History Version 8.0: 1983 Mar 15 - add PLOTFILE parameter, and change format of gray-scale plots (TJP). Version 9.0: 1983 Jul 28 - extend to read MERGE version 2 files (DLM). Version 9.1: 1983 Aug 4 - add W to MERGE version 2 (DLM). Version 9.2: 1983 Aug 15 - sort uvdata before weighting (TJP). Version 9.3: 1983 Aug 18 - add FAST weighting option; label uvplot in wavelengths (TJP). Version 9.4: 1984 Dec 10 - add noise calculations; make UVWEIGHT='U' equivalent to UVWEIGHT='F'. Version 9.5: 1985 Jun 26 - pass antennas and epoch through to CLEAN. Version 9.6: 1986 Apr 4 - pass map center through to CLEAN. Version 9.7: 1986 Aug 28 - increase to 20 stations. Version 9.8: 1988 Jan 11 - Add SNR weighting (WTPOWER) (RCW); minor changes for Convex (TJP). Version 9.9: 1988 Feb 25 - Add beam size calculation (TJP). Version 10.0: 1988 Apr 10 - Convex-compatible version (TJP). Version 10.1: 1988 Dec 19 - Minor changes to plotting; issue warning if oversampling is < 2.5 (TJP). Version 10.2: 1989 Jul 13 - Two input files (TJP). Version 10.3: 1991 Mar 21 - increase to 40 stations (TJP).