Original: Juan Uson, 8/99
Updated: J. Hibbard, 7/05
VLA HI AIPS Data Reduction Guide
This is a step-to-step instruction how to calibrate, map, and reduce a
typical observation of a small galaxy in the HI line using the VLA in
the D-configuration. Some familiarity with AIPS fundamentals is
assumed. Instructions are not necessarily complete; it is strongly
advised to carefully read the HELP file for most of the tasks
mentioned here. See also the
AIPS COOKBOOK, especially
Chapter 4: Calibrating Interferometer Data (ps.gz) and
Chapter 8: Spectral-Line Software (ps.gz). For a more streamlined guide, see
Appendix B: A Ste-p-by-step guide to Spectral Line Data Analysis in AIPS (ps.gz).
The two files needed are called
SUM_SPEC_0.gz and
SUM_SPEC_LINE.gz.
Using ftp, transfer these files to your own workstation. Size of the
files is 7.4 and 43 MB, respectively, so you may want to do this
during a time of the day with low net traffic. Once on your own disk,
these files have to be 'unzipped':
%gunzip *.gz
After this, the sizes are 15 MB and 70 MB.
Set up a variable pointing to the directory where you
put the relevant files to be read into aips.
%setenv MYFITS '/users/myname/fitsdirectory/'
%aips
pick a printer, then enter your aips number
Next, reset aips parameters and load in the procedures you'll
need for the default calibration:
>RESTORE 0
>RUN VLAPROCS
Make it so that the command line comes back after starting
each aips task.
>DOWAIT=0
see how much free disk space there is:
>FREE
Chose the disk with the largest free space which isn't marked
as RESERVED. If this is disk 'n' then designate it by:
>INDISK n
>OUTDIS n
Now read your data into aips:
FITLD or FILLM
FITLD is used to read in the UV data from the FITS files
SUM_SPEC_0 and SUM_SPEC_LINE.
FILLM is used for archival data or native
For the FITLD case:
>TASK 'FITLD'
>UVCOMP=0
>CLRONA
This makes OUTNAME, OUTCLASS etc. blank
> GO
This runs the task with the given inputs. In this case,
The datasets will both have NAME='13/04/95';
one has CLASS 'CH0' and one class 'LINE'. In subsequent steps,
the flux and phase calibration will be done on the CH0 data
and then transferred to the LINE data.
For the FILLM case:
>TASK 'FILLM'
>CLRONA
. >UVCOMP=0
>DOUVCOM=-1
>DOWEIGHT 2
>CPARM 0;CPARM(7)=150
>GO
>UC
Gives you a listing of the UV data in your user area.
For VLA Archival data:
Example:
want to read in 2 files VLAOBS:
AV220_1
AV220_2
>INFILE 'T:AV220_'
(where T is you export defined directory that contains both files)
for the first file, AV220_1
>NFILES 0
>NCOUNT 1
(this starts at 0 and increments to file _1 so to read in AV220_1)
for the second file:
>NFILES 1
>NCOUNT 1
(to read in AV220_2 ..and so on )
>GETN <ch0>;zap
We don't want to use the on-line "channel 0". We will make
a new version below using AVSPC.
>GETN <line UV data>
>IMH
Look carefully at the header. Check the number of
frequency channels (127). RA and Dec are 0 since there are
several sources in the dataset (multi source file). A list of
the sources is in the SU (source) table. You can look at any
table using PRTAB. Note the number of visibilities in the CH0
data and the LINE data: they should be identical. In CH0, each
visibility contains just one spectral line channel; in LINE,
there are 127 channels.
Also note the listing of tables attached to the dataset. See
CH1 of the AIPS cookbook for a description of these tables.
In the following calibration steps, AIPS simply writes values
into tables. It does not modify the raw UV data at all until
the SPLIT step. To undo calibration or flagging, simply delete
the relevant calibration (SN or CL;2 tables) or flagging (FG)
tables as follows:
>INEXT 'CL'
>INVER 2
>EXTD
BE VERY CAREFUL TO NOT EVER DELETE CL TABLE VERSION 1! IF BY ACCIDENT YOU DO, YOU WILL HAVE TO READ IN YOUR DATA FROM SCRATCH
The calibration procedure is:
- SETJY and GETJY write flux values into the SU table
- VLACALIB writes calibration coefficients into the SN table for times when the array was observing calibration sources
- VLACLCAL interpolates calibration coefficients from an existing SN table at regularly spaced time intervals, writing them into CL table version 2 (CL table version 1 contains 1's for all antenna gains and 0's for phases)
- SPLIT applies the CL table version 2 to the UV dataset
PRTAN
Get a printout of antennae positions, for use in
flagging and chosing a reference antenna:
>TASK 'PRTAN';DOCRT 0;GO
LISTR
Get a printout of the observing sequence, frequencies, etc.
>TASK 'LISTR';
>OPTYPE='SCAN'; DOCRT=-1; OUTFIL=''; GO
Pick up the output and look over the observing sequence.
Identify the flux calibrator, phase calibrator, source, and
observing sequence. Also the total time on each source.
If the LISTR output shows that you have more than one
FREQID in the observation, copy each out to its own file. If
only one FREQID is present (as it is in the example dataset)
you can skip the next step.
In the practice dataset, the flux calibrator is 1331+30500002,
the phase calibrator is 1445+09900002, and the source is
N5921. There is only 1 freqid (corresponding too...)
The flux calibrator was observed... for ...
**
and the phase calibrator was observed... for ...
**
The source was observed...
**
UVCOP [may not be necessary]
Make a dataset containing only a single frequency
>TASK 'UVCOP';
>FREQID=n;OUTNA=<give name>;OUTCL 'UVLINE';GO
>GETN <uvcop output>;GO INDXR
>GO 'LISTR'
This copies the data out for FREQID=n, re-indexes the observe
sequence to include only sources observed with that FREQID,
and then re-runs listr to get a new scan list output. Use
that from now on. Repeat for each FREQID.
SETJY
The flux of the flux calibrator is well tabulated and
its frequency dependence is known. SETJY will calculate this
flux based on the frequency it finds in the header, and puts
the result in the SU table.
>TASK 'SETJY'
>SOURCE '<flux calibrator>' ' ';OPTYPE 'CALC';GO
You can check that it was put into the SU table via:
>INEXT 'SU';DOCRT 132;GO PRTAB
AVSPC
Create a data file that averages data over all channels
except those at the band edges. For calibrators, this file
will have a higher signal to noise than the channel data.
It will be used to derive calibration and flag table entries,
which will then be copied to the line data.
The range of channels should be as many as you expect to use
in your final datacube. They are specified by ICHANSEL, where
the syntax is ICHANSEL [start ch] [stop ch] [increment] 1.
In practice, these should be as follows:
- 127 channel dataset: ICHANSEL 5 124 1 1
- 63 channel dataset: ICHANSEL 4 59 1 1
- 31 channel dataset: ICHANSEL 3 29 1 1
>TASK 'AVSPC'
>ICHANSEL 5 124 1 1
>OUTNA INNA; OUTCL 'UVCH0';GO
TVFLG [calibrators, no calibration]
Run TVFLG on the CH0 data, only looking at the calibrator
scans (SOUR='' and CALCO '*'). In this case, the data are
of sufficient quality that flagging of bad data is not
necessary, but TVFLG is an excellent way to get acquainted
with the data.
>TASK 'TVFLG'
>GETN <ch0>
>DOCAL 0; DOCAT 0;SOUR ' ';
>CALCO '*'
>DPARM 0;DPARM(6) <integration time>;GO
You should know the intgegration time from the observe file.
If it is archival data, you may have to figure it out via
trial and error. If the integration time is set too short,
blank lines will appear in between each line of real data.
Change DPARM(6) until the data are displayed without blanks,
except for the times between calibrator scans.
Most useful ways to dislay data here: delta amplitude;
delta vector amplitude. Only flag grossly incorrect things:
bad first scans; obviously bad times (horizontal features);
obviously bad baselines (vertical features).
When flagging, make sure that the "all channels" flag is set
(see writing at bottom of screen). This way, when you copy
the FG table from your CH0 data to the LINE data, the flags
will apply to your channels.
You also probably want to make sure the "all source" flag is
set, so that when you flag an antenna or baseline it will apply
also to your source data. But then you must be careful if you
select the "flag timerange" option to not flag in the area
between calibrator scans (unless you really mean to).
You will only want to flag the stokes data that you are looking
at. To flag RR only, set stokes flag to 1000. To flag LL, set
stokes flag to 0100. To flag both RR and LL, set stokes flag
to 1100.
VLACALIB
Begin of the actual calibration. First, need to know if there
are any UVLIMITS for either of your calibration sources.
These are the range of spatial frequencies over which the
calibrators look like unresolved ("point") sources.
You get this information from the VLA calibrator manual, which
can be found online via the following navigation:
www.nrao.edu->VLA->Astronomers->Under "VLA Calibration": VLA
Calibration Manual->Item 4.2 List of Calibrators
(URL in 2005 is
http://www.aoc.nrao.edu/~gtaylor/csource.html).
Each source has two names, one based on B1950 coordinates, the
other based on J2000.
If the UVLIMITS do not overlap, you can run VLACALIB once with
one set of limits and calibrate both sources. If they differ,
you will need to run VLACALIB separately for each source.
In the practice dataset, the flux calibrator 1331+305 has the
following entries in the calibration for L-band:
BAND A B C D FLUX(Jy) UVMIN(kL) UVMAX(kL)
=====================================================
20cm L S S P P 15.00 visplot
So it is a Prefered calibrator when the VLA is in D-array,
with a flux at L-band of ~15 Jy, and no UV limits. The visplot
link will plot its visibilities - print that out for your notes.
Similarly, the phase calibrator 1445+099 is a Prefered
with no uvlimits and a flux of ~2.6 Jy.
In this case, there will be no UVLIMITs, and a single pass of
VLACALIB with CALSOU ' ' and CALCO='*' will work.
RUN VLAPROCS
>TASK "VLACALIB'
>CALSOU, UVRA, CALCO=<values based on calibration manual>
>MINAMP 10; MINPHS 10; (only print out closure errors that
deviate by more than 10% or 10deg)
>REFANT n (chose an antenna near the center of the array, as
determined from the PRTAN output)
>DOCAL 0;FLAGVER 1;DOPRINT 1
>GO VLACALIB
Three things will be printed out: closure errors, table of
antenna gains for times observing each calibrator, and table
of antenna phases (relative to the REFANT).
There will be separate solutions for each polarization and
each IF.
The closure errors should not be too far from MINAMP, MINPHS.
From the tables, check that the amplitudes are reasonably
constant and that there are no large phase jumps. If they
do, identify and flag the bad data using a combination of
TVFLG, UVPLT, UVFND, UVFLG. It is often easier to find bad
data after calibration has been applied, so at this step I
usually only look for obviously bad baselines, times or
antennas.
GETJY
Based on the flux of the flux calibrator and on the
values in the two SN tables, GETJY will calculate the flux of
the phase calibrator and put it in the SU table. This process
is necessary since the fluxes of phase calibrators are not as
constant and well monitored as are the few flux calibrators in
the sky. Again, check with PRTAB on the SU table that the
value found has been put in the table.
>TASK 'GETJY'
>SOUR '<phase cal>'
>CALSO '<flux cal>'
>GO
Keep note of the reported value.
VLACLCAL
Interpolates SN table(s); puts results in CL table version 2.
A 'skeleton' CL table was already present in version 1.
How you interpolate the values in the SN table depends on
the observing sequence. If each time on the source is preceded
and followed by a scan on the phase calibrator, then use
SIMP interpolation. Otherwise, use a BOXcar.
First, calibrate your source by interpolating the SN
solutions from the phase calibrator:
>TASK 'VLACLCAL'
>REFANT n
>SOUR '<source>'
>CALSO '<phase cal>'
>INTERPOL 'SIMP' (or INTERPOL 'box')
>BPARM 0 (or BPARM 1,1,0 for INTERPOL 'box' and a 1hr boxcar)
>GAINVER 0;GAINUSE 2;GO VLACLCAL
Next, calibrate your calibrators by using themselves:
>SOUR ' ';CALSO ' ';CALCO '*';GAINUSE 0
>GO
UVPLT
Check calibration on phase and flux calibrator
>TASK 'UVPLT'
>UVRA 0; source <phasecal>
>XINC 5;
>DOCAL 2;GAINUSE 2;
>APARM 0; DOTV 1;GO
If this takes too long, then decrease XINC. If the plot looks too sparse, increase it. Make a hardcopy of the results for your records:
>DOTV 0;GO;WAIT
>DPARM 0;GO LWPLA
All visibilities should have uniform value of level specified
by SETJY or GETJY, at least over applicable uvlimits. If not,
something is seriously wrong. You will need to track down
the bad data using tasks like TVFLG, UVPLT, UVFND, UVFLG,
then delete the current calibration (SN table and CL table
version 2: INEXT 'SN';INVER 0;EXTD;INEXT 'CL';INVER 2;EXTD)
reset the flux values in the SU table (SOUR ' ';OPCO 'REJY';
GO SETJY), and redo all steps from SETJY to UVPLT.
TVFLG [calibrators, with calibration]
Rerun TVFLG on calibrators (CALCO '*'), this time
with DOCAL=2. All baselines should have about the same
amplitude now. Besides the displays you used last time, its
also good to look at amplitude and delta phase.
IMAGR [phase cal ch0]
If UVPLT looks good, make sure there are no major imaging
errors (see chapter on Error Recognition in Synthesis Imaging
book) by mapping your phase and/or bandpass calibrator.
IMAGR is an all-in-one task which applies the complex
gains in the CL and BP tables (i.e. performs the actual
calibration), does the Fourier transform to make maps in the
sky plane, and optionally "cleans" the data as well.
Map CH0 of the phase calibrator first; the result should be
a point source in the field center, just like the beam.
For cellsize, you want ~3 pixels across a beam width.
The beam width is printed to the message
window when you run IMAGR. Rough cellsizes: 12", 4",
1", 0.25" for data taken with the VLA in the D-, C-,
B-, or A-arrays at L-band. Here, the image size doesn't
need to be large, and should be a power of two.
Try IMSIZ 64 64. Use ROBUST=-1 for maximum resolution.
>TASK 'IMAGR'
>OUTNA 'PHCAL.DR-1'
>DOCALIB=2; FLAGV=1
>ROBUST=-1
>CELLSI 12 12
>IMSIZE 64 64
>NITER 500
>DOTV 0
>GO
>MC
MC gives a listing of maps in your catalog. You should see
the map you made (extension IIM), and the dirty beam
(extension IBM). Get the map and display it to the tv:
>GETN <phasecal map>
>TVALL
FFT (includes COMB & UVNOU)
If the above has artefacts, you have to track them down. A very
good diagnostic is looking at the baseline amplitudes and phases
in the UV plane. You do this by doing a Fast Fourier Transform (FFT)
of your map.
Running the task FFT on an image produces a Real & Imaginary FFT, which is
basically all the antenna baselines as seen from the source. Short
baselines are in the middle (note the hole in the center, since there
are no baselines with zero spacing), and long baselines are in the
outer regions.
To make an amplitude and phase FFT:
>TASK 'FFT'
>OPCO ' '
>OUTNA INNA;GO
>TASK 'COMB'
>GETN <real part of FFT>
>GET2N <imaginary part of FFT>
>OPCOD 'POLI';GO
>OPCOD 'POLA';GO
>GETN <poli output>
>TVALL
Bad baselines should show up very clearly. You need to find
them and flag them. Use either TVFLG, UVPLT or UVFND. If you flag
significant data, it may be necessary to redo your calibration.
(See steps outlined above under UVPLT notes).
Look at R. Bracewell's "The Fourier Transform and its applications" to
help you understand what you might expect in the Fourier plan for various
types of sources (gaussians, point sources, etc). Or there is a great applet
at
http://www.narrabri.atnf.csiro.au/astronomy/vri.html. So, e.g. a
point source at the phase center should the same amplitude on all
baselines, and no phase structure (random phases). A point source
offset from the phase center should have a uniform amplitude and a
sinusoidal phase pattern. Try some at the above applet.
This is a really powerful diagnostic. The reason is that real
sources, or weak sources, NEVER look like the array or are confined to
a discrete number of points in the "UV" plane (which is basically what
the FFT is showing you). So any high points that are confined to a few
number of points is showing you problems with a correlator (if only 2
points - the UV point and its complex conjugate), or a baseline (if
an arc), or a set of points (e.g. cross talk, which affects adjacent
antennas). It also shows you when you have problems with short
baselines or long baselines.
Sometimes its good to FFT the continuum, since the S/N is highest
there, but often its best to FFT a continuum subtracted channel with
little or no signal, since that should have NO structure in either
amplitude or phase. But you don't have spectral line data, so I just
throw this out there FYI. If you expect the source to be unpolarized,
then you can use the stoke V, U maps as good diagnostics.
The trouble is tracking down the offending data in the UV visibilities
based on the FFT map - there is no easy way to go between them in
AIPS. Use TVFLG, UVFND, UVPLT. For example, if you see e.g. a
strong SW-NE stripe in your map, this suggests a
problem with a NW-SE baseline. After an FFT, you might see
are a couple of whopping points, which should be one baseline at one
time with a high amplitude. This should be fairly easy to find using
TVFLG (see below), or even by doing a UVPLT using the defaults (plots
amplitude vs. baseline length; set XINC equal to something like 2-5 and
DOTV=1; if not getting enough points then decrease XINC) - the vast
majority of points will cluster at some reasonable flux level, and
there will be a few at high amplitude. These need to be clipped
(either CLIP or TVFLG), and then re-image, re-FFT and see if there are
still problems at this position angle. Another example is if you see an
arc of high amplitudes in the FFT - this means a single baseline is bad
for a certain period of time. If you see a set of arcs, it means an
antenna is bad for a certain period of time.
A special case is when you find the high points all lying along the x-axis
in the FFT. This is a sign of RFI (see the EXPLAIN for the task UVNOU).
To get rid of it, use the task UVNOU. The format for UVRANGE is different
for this task than most others; instead of it being an annulus in
UVdistance, its the distance away from the U and/or the V planes. You
want to use UVRANGE= [a small number],0. I determine the small number
by trial and error - start with something like 100 and increase or
decrease based on what you see in the FFT. You want to flag just the
points very close to the axis. NOTE: unlike some other flagging tasks,
UVNOU flags the actually UVdata instead of producing a flag table.
So if you run UVNOU on your ch0 data, you need to run it with the same
parameters on your line data.
When satisfied:
TVFLG [source, with calibration]
Now run TVFLG on your calibrated source data
>TGET TVFLG
>CALCO ' ';SOURCE 'N5921';DOCAL=2;GAINUS 2;FLAGVER=1;GO
In general, this will not look like your calibrator
data (similar amplitudes on all baselines). Set
"Sort by length" - usually, your shortest baselines
(now displayed to the left side of the tv) will have
more flux than your longer baselines. Flag stuff that
sticks out in time or baseline. After that, "Enter
Scan Time" of about 10 integration times (performs
a rolling boxcar average), and display amplitude then
phase differences, and clip large deviations ("clip
interactively).
Note that no matter how much you flag here (and it
shouldnt be a lot), it does not affect your calibration,
so you will not need to delete tables and redo anything.
IMAGR [source ch0]
Now map the continuum sources in the field of your source.
Map the entire primary beam (30' FWHM at L-band).
Set DOTV 1 to run TVFLG in interactive mode. Set boxes
tightly around any sources, update them as cleaning progresses.
>TASK 'IMAGR'
>OUTNA 'N5921D.CH0-1'
>DOCALIB=2;FLAGV=1;ROBUST=-1
>ROBUST=-1
>CELLSI 12 12
>IMSIZE 512 512
>NITER 500
>DOTV 1
>GO
If problems, try FFT/COMB, UVPLT, etc.
If ok, then spatial calibration is done; copy CL;2 and FG
tables to the line data, and proceed with Bandpass
calibration.
TACOP
Copy table CL version 2 and FG table from CH0 to LINE data.
>TASK 'TACOP'
>GETN <line uvdata>
>outna inna;outcl incl;
>GETN <ch0 uvdata>
>INEXT 'fg';inver 1;outver 1;go
>INEXT 'CL';inver 2;outver 2;go
This assumes that the complex gains can be split in a time
dependent part and a frequency dependent part. CALIB and CLCAL
take care of the time dependent - but frequency independent -
part; in the next step we will determine the frequency dependent
complex gains.
BPASS
Creates BP table; gains across band. Use 1331+30500002
as bandpass calibrator (it is the only calibrator which is
bright enough that it has sufficient S/N
in each channel to
do a proper bandpass calibration).
There are many ways to do this. The following calculates
the bandpass response by dividing each visibility by its
average over the many channels (as specified by ICHANSEL).
>TASK 'BPASS'
>GETN <line uvdata>
>CALSO '<bpass calibrator>'
>UVRA 0; DOCAL 0;FLAGVER 1;SOLINT 0;DOBAN 0;
>REFANT n (same as vlacalib)
>BPASSPRM 0,1,0,0,-2,0
>ICHANSEL <same as AVSPC>
>GO
POSSM
Plot the bandpass response for each antenna, and make sure no
sudden jumps, which would indicate RFI.
>TASK 'POSSM'
>APARM 0,1,0.9,1.1,-10,10,0,2,0
>NPLOTS 9
>DOTV 1
Inspect the bp table for each antenna, showing 9
plots per page. If ok, make hardcopies for posterity.
>DOTV -1;go
Now check how flat your bandpass solution makes your phase
calibrator calculating the global spectrum of your phase
calibrator with both bandpass and calibration applied
>APARM 1,0
>NPLOTS 0
>SOUR '<phase cal>'
>DOCAL 2; DOBAND 1
>DOTV 1
Should be flat, and at the flux level given by GETJY.
If ok, make hardcopies for posterity.
>DOTV -1;go
Print out copies of the above plot files. Do an >IMH on your
uvdata to find the number of plot files (call this n), then:
>task 'lwpla';dparm 0
>for i=1:n;plver i;wait;go;end
If you only observed your BP calibrator once, you are done.
If you observed it several times, you have another decision
to make - how to apply the bandpass. Type ">help doband" for
the options. To help decide, repeat the last POSSM option
with DOBAND=2, and again with DOBAND=3. Print out each.
Pick the option which makes the phasecal spectrum look flattest.
Quantify the flatness of the bandpass by looking at the
maximum deviation across the band. Call the devaition delta
over N channels, then the bandpass is flat to Delta/Flux/Nchan.
This is a good quantity to know.
ALTDEF
For older data, it is necessary to enter the velocity which
corresponds to the observed sky frequency relationship into
the header. First do an IMH on the UVData. If it does not
have lines that look like the following:
AIPS 1: Rest freq 1420.406 Vel type: OPTICAL wrt SUN
AIPS 1: Alt ref. value 1.4800E+06 wrt pixel 32.00
then you will need to do this. You need to know the observed
velocity either from the original LISTR 'scan' output or the
observe file. For the practice data:
>AXTYPE 'OPTHEL'
>AXVAL=<the radial velocity of source in m/s>
>AXREF=<central channel>
>RESTFREQ=1420E9,40575
>ALTDEF
IMAGR [phasecal channels]
It is good practice to make a channel cube of the phase cal
first; the result should be a point source in the field
center, just like the beam. I always use BCHAN=1,
so I can easily go from a plan in the image cube to the
correspinding channel in the dataset. ECHAN=125 to avoid
the band edges. Clean the image lightly, i.e. NITER=200.
>TASK 'IMAGR'
>OUTNA 'phcal<array>-1'
>DOCALIB=2; GAINUSE 2; FLAGV 1;
>DOBAND=n (from POSSM)
>BCHAN 1; ECHAN 125
>CELLSI 12 12
>IMSIZE 64 64;
>UVWT=''; ROBUST=-1;
>NITER 200
>DOTV 0
View using TVMOVIE, LTYPE=6 so that the
channel number is displayed in each frame. You might have
to change PIXR in order to see the signal, e.g.
>GETN <phasecal cube>
>PIXR -0.01,0.03
>TVMOV
ISPEC
Check that spectra are flat using ISPEC:
>TASK 'ISPEC'
>GETN <phasecal data cube>
>OPTY 'FLUX';LTYPE 6
>BLC 0 0 4;TRC 0 0 125;DOTV 1;GO
This better be flat with the flux reported by GETJY
Make a hardcopy (dotv -1;go;plver 0;go lwpla).
IMAGR [source line+continuum data]
If that looks good, re-run IMAGR on your source with the
continuum still in it. We will not keep this; its to get an
idea of which channels contain line emission.
Note that IMAGR has adverbs BCHAN and ECHAN, so it needs
to be run only once to do all channels.
You may want to first run IMAGR on a single test channel
(e.g BCHAN=64; ECHAN=64); this also allows you to test
your mapping and cleaning parameters.
Use a IMSIZE large enough to map the primary beam (30').
CELLSIZE as before. Shallow clean to remove large sidelobes.
>TGET IMAGR
>SOURCE='N5921'
>OUTNA 'N5921D.r=1c'
>ROBUST=1
>IMSIZ 512 512
>NITER 500;DOTV 0;
>GO
The result will be a 'cube' of maps and a 'cube'of beams,
with frequency along the third axis. Check this with
IMHEAD. Make note of the "Conv size" & "Position angle"
from the header - this is the size (FWHM) of the Gaussian
used to replace dirty point sources in the original map. It
is therefore the resolution of your map. You will need it
below. This map contains signal from both continuum and
spectra line sources. If it looks good, then its finally time
to split the source data into its own dataset, with all
calibration and flags from the tables applied to the data.
TVMOV with LTYPE=6 a few times and try to identify the
first and last channels where line emission appear.
SPLIT
Use SPLIT to convert your multisource file into datafiles
containing just a single source. SPLIT applies the complex
gains in the CL and BP tables (i.e. it performs the actual
calibration).
>TASK 'SPLIT'
>SOUR 'N5921'
>CALCO ' ';DOCAL 2;GAINUS 2;FLAGVE 1
>APARM 0;NCHAV 1;ICHANSEl 0
>DOBAND=<from POSSM>;
>GO
I'll generally rerun the last IMAGR again, using the newly
split dataset (but now setting DOCAL=-2; DOBAND=-1), to make
sure everything got applied correctly.
Give the SPLIT file a descriptive name:
>GETN <split dataset>
>OUTNA 'N5291D';OUTCL 'UVLINE';RENA
OPTIONAL: Self Calibration
- AVSPC line free channels from split data
- IMAGR Robust-1, UVRANG, DOTV 1, use tight boxes
- CCMRG on resulting image
- PRTCC to find number of clean components to use
- CALIB getna [AVSPC data], get2n [AVSPC map], ncomp=from prtcc, uvra as in IMAGR, calso [source], solmo 'p' for phase only selfcal, refant [as in original vlacalib], solint 10, aparm 4,0,0,0,2,0, cparm 0, go
- IMAGR new uvdataset. If better, EXTD 'sn' file from avspc data, ccmrg,prtcc,calib, etc.
- TACOP SN file from AVSPC data to SPLIT UVLINE data
- SPLIT UVLINE data with SN table applied. OUTNA .SC
TASAV
Make a dummy file which contains no data, but all the
calibration tables. You will keep the split data and
resulting products, but there is no need to keep the
raw data, since you can always get this back from the
archive. But if you keep a file with your calibration
tables attached, then if you do need to go back to the
raw data, you can easily copy the calibration and flags
from this file.
> Task 'TASAV'
> GETN <uvdata before split>
> OUTNA 'N5921D.MS';OUTCL 'TASAV';GO
ZAP
This is a good time to get rid of files you no longer need.
Basically, at this stage you only really need the SPLIT
and TASAV uvdata, and the best cleaned ch0 map.
We will make a new ch0 data later. We don't need any phase
calibrator maps, nor source maps with continuum in them.
So for each file you don't want:
>GETN <i>
>ZAP
To do a whole sequence of files, set up a loop:
>FOR i=4:18;getn i;zap;end
When this is cleaned up, get remaining data into adjacent
slots on disk:
>recat
UVLIN
Now its time to remove the continuum from the line data.
Subtracting a continuum from UV data can be done using the
task UVLIN. This is an iterative procedure: UVLIN it, map it,
review it; if negatives appear or ISPEC doesn't look right,
try a new set of channels and repeat.
From the shallowly cleaned source datacube, you identified
line-free channels. Lets say that the line emission appears
from channel A to B. We'll guess that the line free channels
are b=A-2, c=B+2. For the test data, these are ????
We also need the first and last channels used in our original
AVSPC (used to make ch0). Call these a and d (e.g., for the
test dataset these were 5, 124. So the line-free channels
that we want to fit are a->b and c->d. These are specified
in UVLIN using ICHANSEL=a b 1 0 c d 1 0
But before running, we note that UVLIN also flags data that
deviate in some channel by more than FLUX from the mean over
the included channels. See EXPLAIN UVLIN for more on this.
In practice, I determin FLUX interactively - I start with
FLUX~3, and increase or decrease it based on how many
correlators it flags - on otherwise well behaved data,
I want to flag ~<1% of the correlators.
>TASK 'UVLIN'
>GETN <UVLINE DATA>
>OUTNA INNA
>OUTCL 'UVLIN'
>DOCONT 0;ORDER 1;PRTLEV 1
>ICHANSEL=5 ?? 1 0 ?? 125 1 0
>FLUX 3
>GO
IMAGR,ISPEC,UVLIN, repeat
Make shallowly cleaned datacube, inspect it, use ISPEC to
take the spectra through the location of continuum sources
(get the pixel location of continuum sources from your ch0
map made with the same CELLSI, IMSIZE as your channel map
using the verbs CURV or IMPOS, and use these for BLC,TRC
in ISPEC) and your line source. When you are satisfied
that UVLIN only included line-free data, make a deeper
cleaned datacube.
IMAGR [UVLIN'd source, R=+1]
Use the same inputs as before, but a larger NITER
Also set BMAJ,BMIN,BPA so that each channel has the
exact same resolution. Use Robust=1 to improve sensitivity
to extended resolution.
??(NBOX? CLBOX?)
>TGET IMAGR
>SOURCE='N5921'
>OUTNA 'N5921D.r=1c'
>ROBUST=1.
>BMAJ <major axis FWHM from previous IMAGR run, rounded off>
>BMIN <minor axis FWHM from previous IMAGR run, rounded off>
>BPA <position angle from previous IMAGR run>
>IMSIZ 512 512
>NITER 1000;DOTV 0;
>GO
Get the cube and change the 3rd axis from frequency to
velocity, using ALTSW
Also, I always rename cleaned cubes with the OUTCL 'LMV',
which designates the cube order as being RA, DEC, VEL.
>GETN <r+1 cube>
>ALTSW
>IMH (notice units are now m/s on 3rd axis, and increment is
your channel spacing in m/s)
>OUTNA INNA,OUTCL 'LMV';RENAM
IMAGR [UVLIN'd source, R=-1]
I always also make a R=-1 version of the data. You'll
have to run it on the central channel first to get
BMAJ,BMIN,BPA
>TGET IMAGR
>ROBUST=-1
>BCHAN 64
>ECHAN 64
>OUTNA 'ch64.r-1'
>GO
>GETN <ch64.r-1>
>IMH (make note of beam parameters)
>TGET IMAGR
>BCHAN 1
>ECHAN 125
>OUTNA 'N5921D.r=1c'
>GO
ALTSW and rename as above
>GETN <r+1 cube>
>ALTSW
>OUTNA INNA,OUTCL 'LMV';RENAM
When satisfied, we want to make two new ch0-like datasets from
the split (i.e. not UVLIN'd) datacube.
_1 A new ch0, but this one should be nicer b/c the bandpass
shape has been removed and we'll use UVLIN to flag outliers
_2 A line-free continuum map, made from line free channels.
UVLIN (flagging only)
First run UVLINE only to flag outliers; retain the continuum
>TGET UVLIN
>DOCONT 1
>OUTCLA 'LINE+C';GO
AVSPC
now make the new ch0 and line-free continuum
>TASK 'AVSPC'
>GETN <uvlin docont=1 output>
>ICHANSEL 5 125
>OUTNA INNA
>OUTCLA 'UV CH0'
>GO
>INCHANSEL 5 *A* *B* 125
>OUTCL 'UVCONT'
>GO
IMAGR [new Ch0, line-free continuum]
Use IMAGR to map both of those datasets.
Use same CELLSI, IMSIZE so that you can correlate
continuum sources with line data.
>TGET IMAGR
>ROBUST -1
>GETN <ch0>
>OUTNA 'N5921Dr-1.ch0'
>NITER <large number - depends on number of sources. ~5000>
>DOTV 1
>GO
>GETN <uvcont>
>OUTNA 'N5921Dr-1.nol'
>GO
IMEAN
Measure noise in each map using TVWIN, GO IMEAN
Make note of this value. Also note number of pixels
per beam. These are needed later.
TRANS
Many routines that work on cubes and velocity profiles
require the data with frequency - or velocity - on the first
axis. This is done using the task TRANS with TRANSCOD='312' on
the cleaned cube. Use IMHEAD on the result to check whether
TRANS did what you wanted.
>TASK 'TRANS'
>GETN <r+1 datacube>
>BLC 0; TRC 0
>TRANSCO '312'
>OUTNA INNA
>OUTCL 'VLM'
>GO
>GETN <r-1 datacube>
>OUTNA INNA
>GO
JEH...I'M HERE...
MOMNT
MOMNT creates maps containing the various moments of
the velocity profiles. It works on the cube you just created
with TRANS. Read the HELP file for MOMNT carefully, and play
around with the values of CELLSIZE, ICUT, and FLUX.
>TASK 'MOMNT'
>CELLSI 3 5
>BLC a 0 0
>TRC b 0 0
>OUTCL '0'
>FLUX -1
>ICUT (1sigma)
>GO
# ICUT~ 1-1.5* σ
#repeat for ZMOM0:
>FLUX -1
>ICUT~ 1-1.5* σ (Jy/b)
>OUTN INNA
>OUTCL '0'
>GO
MOM0,1,2
>FLUX 0
>OUTN INNA
>OUTCL '012'
>GO
--
AmandaH - 25 Jul 2005
NOTES:
Naming convention:
I find it very useful to give each
file a descriptive name so that you can tell what it is
without having to do an IMH. My naming convention is:
UVDATA:
OUTNA='<source><array><MS if dataset contains multisource data>'
OUTCL=one of the following:
'UVLINE': uv line data.
'UVCH0': average over many channels.
'UVCONT': average over line-free channels.
'UVLIN': uv line data after UVLIN.
'TASAV': uv dummy file with all calibration tables attached.
MAPS:
OUTNA='<source><array>.<Robust value><ch0 or cont or nothing>'
OUTCL=one of the following:
'ICLN': cleaned single channel map.
'IMAP': uncleaned map or cube.
'LMV': cleaned cube in [x,y,z]=RA,DEC,Vel order
'VLM': cleaned cube in [x,y,z]=Vel,RA,DEC order
'ZMOM0': mom0 made with flux=-1
'MOM0,MOM1,MOM2': moment maps made with flux=0
To change the name of a file:
>getn <file no.>
>outna 'n5921d.ms'
>outcl 'uvline'
>renam
Final Products:
Print Outs:
- Antenna positions from PRTAN
- Original LISTR OPTY='SCAN'
- Final LISTR OPTY='SCAN' (after UVCOP)
- printout from original VLACALIB
- printout from final VLACALIB
Plots:
- UVPLT for each of your calibrators and source
- BP table for each antenna (POSSM, APARM(8)=2)
- Spectrum of Phase Cal with BP applied (POSSM, APARM(1)=1)
- ISPEC of source
- ISPEC through cube at location of brightest continuum source
- KNTR of line-free continuum greys with Ch0 contours
- KNTR of channel maps (contours on optical greys if possible)
- KNTR of mom1 contours on greys of mom0
- KNTR of mom1 contours on greys of mom2
UV data:
- [name][array].MS.TASAV (dummy uvfile with BP,CL,SN,FG tables)
- [name][array].UVLINE (uv line data)
- [name][array].UVCH0 (average over many channels)
- [name][array].UVCONT (average over line-free channels)
- [name][array].UVLIN (uv line data after UVLIN)
MAPS:
- [name][array]-1.CH0.ICLN (final ch0)
- [name][array]-1.-L.ICLN (line-free continuum)
- [name][array].R-1.LMV
- [name][array].R-1.ZMOM0
- [name][array].R-1.MOM0
- [name][array].R-1.MOM1
- [name][array].R-1.MOM2
- [name][array].R+1.LMV
- [name][array].R+1.ZMOM0
- [name][array].R+1.MOM0
- [name][array].R+1.MOM1
- [name][array].R+1.MOM2
Things to document for each map:
- theoretical and measured noise (mJy/b)
- clean beam size FWHMxFWHM at P.A.
- For continuum maps: channels which went into making them
- For moment maps: icut, flux, cellsi,
- for moment 0 maps: total flux (Jy km/s)
Calculating Theoretical Noise:
S_rms = k*fw/sqrt(N_ANT*(N_ANT-1)*N_IFs*Δt_hr*Δν_MHz) [mJy/beam]
k~8.0 @ L-band (from VLA "observational status summary")
fw ~ 1.03 (Robust = +1)
~ 1.2-1.4 (Robust =-1 (see task 'imagr' output)
N_ANT ~27
N_IFs = # of polirization pairs (LL & RR)
Δt_hr = time on source from observing log in hours
Δν_MHz = bandwidth (from imh of data cube)
--
AmandaHeiderman - 25 Jul 2005
Measuring Noise:
>GETN <final image cube R=+1>
>TVWIN (select area near source)
>GO IMEAN
write down rms, # pixels/beam, and # pixels in chosen area
repeat for R=-1 cube
--
AmandaHeiderman - 25 Jul 2005
Measuring Flux:
Flagging based on UVSUB residuals
--
JohnHibbard - 07 Sep 2005
--
AmandaHeiderman - 08 Sep 2005 (updated)