%-----------------------------------------------------------------------
%;  Copyright (C) 2004, 2009
%;  Associated Universities, Inc. Washington DC, USA.
%;
%;  This program is free software; you can redistribute it and/or
%;  modify it under the terms of the GNU General Public License as
%;  published by the Free Software Foundation; either version 2 of
%;  the License, or (at your option) any later version.
%;
%;  This program is distributed in the hope that it will be useful,
%;  but WITHOUT ANY WARRANTY; without even the implied warranty of
%;  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
%;  GNU General Public License for more details.
%;
%;  You should have received a copy of the GNU General Public
%;  License along with this program; if not, write to the Free
%;  Software Foundation, Inc., 675 Massachusetts Ave, Cambridge,
%;  MA 02139, USA.
%;
%;  Correspondence concerning AIPS should be addressed as follows:
%;         Internet email: aipsmail@nrao.edu.
%;         Postal address: AIPS Project Office
%;                         National Radio Astronomy Observatory
%;                         520 Edgemont Road
%;                         Charlottesville, VA 22903-2475 USA
%-----------------------------------------------------------------------
\documentclass[preprint]{aastex}
%
% probably don't need all these...
%
\newcommand{\AIPS}{{$\cal AIPS\/$}}
\newcommand{\POPS}{{$\cal POPS\/$}}
\newcommand{\ttaips}{{\tt AIPS}}
\newcommand{\DELZN}{{\tt DELZN}}
\newcommand{\eg}{{\it e.g.},}
\newcommand{\ie}{{\it i.e.},}
\newcommand{\daemon}{d\ae mon}

\newcommand{\boxit}[3]{\vbox{\hrule height#1\hbox{\vrule width#1\kern#2%
\vbox{\kern#2{#3}\kern#2}\kern#2\vrule width#1}\hrule height#1}}
%
\newcommand{\memnum}{110}
\newcommand{\whatmem}{\AIPS\ Memo \memnum}
\newcommand{\memtit}{Strategy for Removing Tropospheric and Clock Errors
using \DELZN \\ Version 2.0}
\title{
%   \hphantom{Hello World} \\
   \vskip -35pt
%   \fbox{AIPS Memo \memnum} \\
   \fbox{{\large\whatmem}} \\
   \vskip 28pt
   \memtit \\}
\author{ Amy J. Mioduszewski \& Leonia Kogan\\
National Radio Astronomy Observatory \\
\vskip 14pt
Revised: October 21, 2009}

%
\parskip 4mm
\linewidth 6.5in                     % was 6.5
\textwidth 6.5in                     % text width excluding margin 6.5
\textheight 8.91 in                  % was 8.81
\marginparsep 0in
\oddsidemargin .25in                 % EWG from -.25
\evensidemargin -.25in
\topmargin -.5in
\headsep 0.25in
\headheight 0.25in
\parindent 0in
\newcommand{\normalstyle}{\baselineskip 4mm \parskip 2mm \normalsize}
\newcommand{\tablestyle}{\baselineskip 2mm \parskip 1mm \small }
%
%
\begin{document}

\pagestyle{myheadings}
%\thispagestyle{empty}

\newcommand{\Rheading}{\whatmem \hfill \memtit \hfill Page~~}
\newcommand{\Lheading}{~~Page \hfill \memtit \hfill \whatmem}
\markboth{\Lheading}{\Rheading}
%
%

\vskip -.5cm
\pretolerance 10000
\listparindent 0cm
\labelsep 0cm
%
%

%\vskip -30pt
%\maketitle
%\vskip -30pt


% \normalstyle
%
\vspace{5mm}

\begin{abstract}
This memo provides a guide to scheduling and reducing phase
reference experiments in order to improve the astrometric accuracy and
image quality of a target source.  The recommended procedure includes
occasional short periods of observations of calibrator sources around
the sky, interspersed with the desired phase reference observations,
from which the troposphere modeling errors can be determined using a
new AIPS task DELZN.  The recommended observing procedure, data
reduction, running of DELZN and application of corrections to the
phase reference data are covered.
\end{abstract}

\section{Basic Strategy}
In order to get accurate relative positions and
good image quality with VLBI, two things are needed; a close
calibrator to the target source, and a good correlator model.
The VLBA correlator uses a model for the troposphere
computed from seasonal models as well as the best estimates
of the clocks, positions of each antenna and earth orientation parameters.
Since none of these estimates
are perfect, errors are introduced which effect the positional accuracy
and the success of phase transfer for phase referencing.
At high frequencies ($\gtrsim 5$\,GHz) the largest errors are from the
clocks and the troposphere.
At lower frequencies the major
source of error is the un-modeled ionosphere, see VLBA Scientific Memo
23 for a discussion of this.
This major source of error produces a systematic phase
difference between the calibrator and target which limits the target
positional accuracy and image quality.  The troposphere model error,
however, can be estimated by observing about ten calibrators over the
sky over a 30-minute period, and then correcting the model error in
the calibrator and target observations.  However, the model errors are
generally sufficiently large so that the measured phase of the all-sky
calibrator observations spans many turns and cannot be used to
determine the troposphere.

    The phase ambiguity problem is solved by the measurement of a
related term, the multi-band (or group) delay, which is the rate of
phase change with frequency.  For the troposphere (or any other
non-dispersive delay), the phase equals the multi-band delay (MBD)
times the frequency.  In order to accurately measure the MBD,
observations over a wide spanned-bandwidth must be made.  An example
of the measurement of the MBD is is given in Figure\,\ref{predelzn}.
Figure\,\ref{predelzn} shows the MBD as a phase slope over
a set of widely spaced IFs.
The MBD is a function of a variable clock delay, and a troposphere
term which has a well-defined function of elevation (called a mapping
function), times an unknown constant.  Thus, from observations of many
calibrators over a wide range of elevations, a smoothly changing clock
term and the troposphere offset (zenith-path delay error) can be
determined for each telescope.  The recent \AIPS\ task, written by L. Kogan,
\DELZN\ does such fitting.  Here we
describe the observation and data reduction method needed to use \DELZN\
to remove the tropospheric and clock delays.  Ideally, one would use a
combination of \DELZN\ and a method using
multiple calibrators near the target to take out local atmospheric gradients
between the calibrator and the target.  \AIPS\,Memo\,111 will describe this
method.

There are two basic issues to take into account while scheduling an experiment
in which these errors will be solved for: 1) calibrators selection and;
2) observing frequency setup.  The next two sections will discuss these issues
while \S 4 will discuss the data reduction method.

\begin{figure}[t!]
\epsscale{0.9}
\plotone{predelzn.ps}
\caption{{\tt POSSM} plot of the cross-correlation spectrum of a set of
widely spaced IFs in the 4cm band. Multi-band delays are evident in the
phase slopes across all the IFs.}
\label{predelzn}
\end{figure}

\section{Observing Calibrators}

This section is not intended to describe the choosing and scheduling of
calibrators for phase referencing.  For more information on this see VLBA
Scientific Memo 24.

First, include the time needed to observe
calibrators in the amount of time proposed for; approximately an extra
hour is needed for every 3-4 hours of normal experiment.
Since the tropospheric delay depends on the mapping function and
the clock delay does not, one wants to measure the MBD at
a wide range of elevations.  The sources should be strong and, to avoid delays
caused by errors in
the source position, use sources with most accurate possible positions.
For this reason,  International Celestial Reference Frame (ICRF) sources should be
preferentially chosen, as they will be both strong and have good positions, a
list of these can be found at the U.S. Navel
Observatory website: http://rorf.usno.navy.mil/ICRF/.  The VLBA calibrator
list also contains ICRF sources, however some of the objects labeled as
ICRF are {\it candidate} ICRF sources and therefore their positions
are not as accurate as true ICRF sources.

The basic idea is to observe 8-12 sources over a
large range of elevations every 4 hours.  These sources should be observed
in blocks, either 8-12 sources in a half hour block or 4-6 sources in a 15
minute block twice as often (\ie\ every 1-2 hours).  Note that since the
frequency setup will change (see next section) between these blocks and the
normal observing, the schedule should not switch to these kinds of blocks more
frequently than hourly because of limitations of the VLBA correlator.
You should have at least two
all-sky calibrator block, before and after the regular observing.  If the
experiment is longer than 4 hours then add at least one more block in the
middle of the run.
Here are the basic steps
to schedule the observation (here we assume 8-12 sources every 3-4 hours):

\begin{enumerate}
\item Choose $\sim 20$ strong calibrators, preferably ICRF sources, at various
declinations that are up during your
observation.  Be sure to choose some that will be at very low elevations
($< 20^\circ $).  You
are choosing more sources than you will observe because some will give you better
elevation coverage and shorter slewing times.
\item Observe $\sim 10$ sources each for about a minute, scheduling
them in such an order to minimize slew time and maximize elevation coverage
at all stations.  Unfortunately these are opposing goals, the most
important thing is to maximize the elevation coverage, but slew times $> 3$\,minutes
should probably be avoided.  Remember the goal is to fit about 10 sources in about
30 minutes, with one minute on source.  Also, it will probably be impossible
for all stations to observe all sources, particularly since it is important
to get low elevation observations at all stations.
\begin{itemize}
\item {\bf Check elevation coverage with {\tt SCHED} for each antenna}.
You can do this one of two ways,
the *.sum listing that {\tt SCHED} produces for your schedule has a listing of
the elevation in degrees at each station for each scan.  The second option
is to run {\tt SCHED} in it plot mode, and plot the source elevation verses time.
\item {\bf Check slew times with {\tt SCHED}}.  {\tt SCHED} will estimate slew times.
If you use {\tt DWELL} to schedule the observing scans, then the gaps between the
scans in the *.sum file will reflect the time it will take to get on source (or
at least {\tt SCHED}'s estimate).
\end{itemize}
\end{enumerate}

This scheduling is the most time consuming part of the entire process.

\begin{figure}[t!]
\epsscale{0.6}
\plotone{dlyfn.ps}
\caption{Delay function for various IF spacings.  Solid line: {\it Golomb ruler}
spacing of 8 IFs across 500\,MHz; dashed line: contiguous 64\,MHz; and dotted
line: 8 IFs evenly spaced across 500\,MHz.
}
\label{dlyfn}
\end{figure}

\section{Frequency Setup}

In order to estimate the MBD most accurately, it is important
to spread the observing IFs across the frequency band in the same band
as the regular observing blocks.  However this must
be done in such a way to minimize the ambiguity in the delay solution.  In
other words, there should be both narrow and wide spacings between the
IFs to rule out very large and very small delays.
For a more thorough discussion of this see the Bandwidth Synthesis section
of Thompson, Moran and Swenson (1994).  Figure\,\ref{dlyfn} shows the
delay function for 3 possible frequency setups.  The dashed line is for a
standard setup, like those which come with the VLBI scheduling software
{\tt SCHED}, where the 8 IFs are joined in one continuous 64\,MHz wide band.
It is obvious that a MBD fit to this data would likely
have an error of one nanosecond or more, which is unacceptable since the real
MBD should be at most one nanosecond.  The dotted line is the delay function
for eight 8\,MHz IFs spread evenly over a 500\,MHz band.  Even though the
peak is narrower there are fairly high secondary peaks at $\sim 4$\,nsec.  These
peaks are caused by the fact that a very high delay can be easily fit to
evenly spaced IFs.  The solid line in is the delay function
for IFs at  42.900, 42.928, 43.068, 43.166, 43.250, 43.320,
43.362 and 43.376\,GHz. Figure\,\ref{dlyfn} shows that this frequency spacing minimizes
secondary peaks in the delay function and provides a narrow central maximum, which
is exactly what is needed.
This last/best option is computed by taking the maximum frequency (43.376\,GHz) and
subtracting $14\,{\rm MHz} \times(0, 1, 4, 9, 15, 22, 32, 34)$.
The $0, 1, 4 ... 34$ is the unique spacing on a {\it Golomb ruler} with
eight marks.  This set of numbers has no pair of marks the same distance apart, so
this measures a unique set of gaps between 14 and 476\,MHz.  Although this
is not necessarily the perfect frequency setup, figure\,\ref{dlyfn} shows that
it is sufficient and it is a nice rule of thumb.

However, there is an important issue that
has been ignored, radio frequency interference (RFI).  It is all well and good
to find the perfect frequency setup, but the IFs must also avoid any known RFI.
Plots of the VLBA RFI survey can be found at
{\tt http://www.aoc.nrao.edu/vlba/html/rfi.shtml\#RFISurvey}.  The frequencies of
the IFs should be checked against these plots to make sure there is no significant
RFI on or near them.

\section{Reducing the Data}
We have assumed a familiarity with \AIPS\ and VLBI data reduction techniques.  For
more information on these topics see the \AIPS\ Cookbook, particularly
Appendix C (A Step--by--Step Recipe for VLBA Data Calibration in \AIPS).

If two different frequency setups were used for observing the calibrator blocks
and the standard phase referencing part of your experiment then {\tt FITLD} will
automatically separate these into two uv-data files.  Calibrate both data sets
using the standard method up to, but not including, global fringe fitting.

\subsection{Obtaining Multi-band Delays}

The next step is to find the MBD for the sources in the calibrator blocks.
There are two ways to get MBD.  One is to use {\tt FRING} with
{\tt APARM(5)=2}
the other is to use {\tt FRING} with {\tt APARM(5)=0} and the separate \AIPS\ task
{\tt MBDLY}\@.
The rule of thumb is to use {\tt FRING} on weak calibrators and {\tt MBDLY} on strong
calibrators.  Professional geodisists etc. use {\tt MBDLY}\@.  For most cases both
methods should produce similar results.  Although, running {\tt FRING} with
{\tt APARM(5)=2} has a significantly longer running time than the combination
of {\tt FRING (APARM(5)=0)} and {\tt MBDLY}.

For the {\tt FRING} only method, it is as simple as running {\tt FRING} with
{\tt APARM(5)=2}\,
applying the previous calibration.  To use {\tt MBDLY}\, run {\tt FRING} with {\tt APARM(5)=0} to get
the single band delays (\ie\ as if you were doing the global fringe fit).  Then run {\tt MBDLY}\,
with {\tt INVERS} set to the {\tt SN} table from {\tt FRING}\@.
{\tt MBDLY} fits a linear phase versus to the single band phases from {\tt FRING}\@.
For most cases the rest of
the inputs for {\tt MBDLY} can be left as the default.  Although you may want to set {\tt APARM(4)}
to a number 1 or 2 less than the total number of {\tt IF}s in the data, otherwise {\tt MBDLY}
will not fit a solution if all the {\tt IF}s are not present.

Use {\tt SNPLT (OPTYPE 'MDEL')} on the output {\tt SN} table to check the
solutions.  Look for outliers, the MBD should less than than a nanosecond.
Figure~\ref{mbdlyvstim} shows a {\tt SNPLT} of the MBD found
by {\tt MBDLY} versus time.  Most of the delays are significantly less than 1 nanosecond, however there are
some bigger delays for {\tt HN} and {\tt SC}.  An examination of the MBD verses
elevations ({\tt SNPLT; XASIS 2}), figure\,\ref{mbdlyvselv},
show that the {\tt HN} delays are probably fine, since they seem consistent with elevation, while
the $\sim -2$\,nsec delay at {\tt SC} is
probably wrong and should be deleted. Figure~\ref{mbdlyvselv} also shows that there may be
a problem with
the solutions for {\tt BR} as well since they do not seem to smoothly vary from low to high
elevation.  At this time the only way to flag bad MBDs is with
tasks {\tt TABED} and {\tt TAFLG}.  If you just have outlier MDBs,
which is the most common problem, you can use {\tt TABED} to clip them.  Set
{\tt OPTYPE = 'CLIP'} and set {\tt APARM(1)} to the column number for the MBD (usually
column 12) and {\tt KEYVAL} to the range of good MBD (\eg\ {\tt KEYVAL=-0.5E-09,
0.5E-09}).
For {\tt TAFLG}, first find which rows numbers to flag by running {\tt PRTAB}
this will print the {\tt SN} table.  To flag specific rows in the SN table
simply set: {\tt INEXT='SN'}; {\tt INVERS=}{\it input {\tt SN} table};
{\tt BCOUNT=}{\it number of the row in the {\tt SN} table to be
flagged}; {\tt ECOUNT=BCOUNT}; {\tt OPTYPE='<>'}; and {\tt APARM(1)=col\#}.
If there are significant problems try the other
method to obtain the MBD and compare the results, with modern computers
neither take much time.

\begin{figure}[t!]
\epsscale{0.75}
\plotone{mbdlyvstim.ps}
\caption{{\tt SNPLT} of MBD verses time.  The MBD were found
by the \AIPS\ task {\tt MBDLY}\@. Note that most of the delays are significantly less than
1 nanosecond.  {\tt FD} has zero delay because it is the reference antenna.}
\label{mbdlyvstim}
\end{figure}

\begin{figure}[t!]
\epsscale{0.75}
\plotone{mbdlyvselv.ps}
\caption{{\tt SNPLT} of MBD verses source elevation.  The MBDs were found
by the \AIPS\ task {\tt MBDLY}\@. {\tt FD} has zero delay because it is the reference antenna.}
\label{mbdlyvselv}
\end{figure}

\subsection{Running \DELZN }

Once an {\tt SN} table with reasonable MBDs has been obtained \DELZN\
can be run.
As mentioned in the first section \DELZN\ takes the MBDs and estimates
the zenith atmospheric delay and, if requested, the clock error.  It does this by fitting a
polynomial function.  If two different frequency setups were not used then
\DELZN\ can be used to correct a {\tt CL} table.  However, for the observing method
we describe here the {\tt CL} table that needs to be corrected is attached to
a separate uv-dataset.  For this case, \DELZN\ will produce a file on disk which
contains the
atmospheric and clock corrections,
which is in the right format to be read in by {\tt CLCOR} and applied to a
{\tt CL} table attached to the normal phase referencing dataset.

However, the best way to proceed is to run \DELZN\
{\it correcting the calibrators first.}
This allows you to examine the correction and do additional editing or
adjust the inputs to \DELZN\, if necessary.

\vfil
\eject

\begin{figure}[t!]
\epsscale{0.9}
\plotone{postdelzn.ps}
\caption{{\tt POSSM} plot of the same scan seen in figure\,\ref{predelzn}
after the tropospheric and clock delays were removed by \DELZN.}
\label{postdelzn}
\end{figure}

Inputs for \DELZN\ should be:
\begin{quote}
{\tt SNVER {\it snin\/} } {$\longrightarrow$ {\tt SN} table containing MBDs.}

{\tt GAINVER {\it clin\/} } {$\longrightarrow$ {\tt CL} table to correct.}

{\tt APARM(2) 3 } {$\longrightarrow$ use 3 term polynomial to fit zenith delay}

{\tt APARM(3) 2 }{$\longrightarrow$ use 2 term polynomial to fit clock}

{\tt APARM(4) 1 }{$\longrightarrow$ create a new CL table}

{\tt APARM(5) 1 }{$\longrightarrow$ to solve for atmosphere and clocks}

{\tt OPTYPE 'MDEL' }{$\longrightarrow$ to use MBD}

{\tt DOTV -1 }{$\longrightarrow$ to make {\tt PL} files}
\end{quote}
The rest of the inputs can be left as the defaults, {\tt CALSOUR} and {\tt
SOURCES} can be left blank because all the sources in this dataset are
calibrators.  \DELZN\ will produce a set of plots that show the fits it
obtained.  You should examine these.  Although, the best way to determine
the goodness of the fit is by applying the resulting {\tt CL} table to
the calibrator data in {\tt POSSM} and seeing if the the corrections
have removed the MBD, \ie\ flattened
the phases.  Figure\,\ref{postdelzn} shows the same scan as
figure\,\ref{predelzn}, but with troposphere and clock delays fitted and
removed.  The phases are much flatter, although in some cases not
completely flat.  If you are not satisfied
with the results you may want to play with {\tt APARM(6-7)} which does some
automatic editing of the {\tt SN} table or do more editing with {\tt TAFLG}.

To correct the target dataset simply set {\tt APARM(4)=0}, {\tt APARM(8)=0},
and {\tt OUTFILE={\it filename}}.  \DELZN\ will then produce an output
file with the zenith atmosphere and clock delays as well as their
derivatives.  This file is designed to be used by {\tt CLCOR
(OPCODE='ATMO')}. Run {\tt CLCOR} on the normal phase referencing data
set, setting {\tt OPCODE='ATMO', INFILE {\it filename}} and {\tt
GAINVER} to the {\tt CL} table with all the calibration except for
global fringe fitting.

After {\tt CLCOR} is run, finish processing as usual, \ie\
run {\tt FRING} ({\it applying the corrected {\tt CL} table}), {\tt CLCAL},
{\tt SPLIT}, and {\tt IMAGR}.  Figure\,\ref{images} shows
images of a weak source ($< 2$\,mJy) in the same data set which was
observed as described in \S 2 and \S 3.  The image on the left was reduced
in the standard manner (\ie\ no \DELZN ), while the image on the right
was reduced as described in this section.  It is obvious that using
\DELZN\ makes a significant improvement, the source on the right
is more compact and has a higher peak than the one on the left.  Also, the
peak is shifted by $\sim 0.4$\,milliarcsec, a significant shift if one wants
to do careful astrometry with errors $< 0.2$\,milliarcsec.  This
dataset was well behaved to begin with, the source was generally at high
elevation ($> 40^\circ$) and the observing frequency was 8.5\,GHz, so most
of the corrections where to the easternmost and westernmost antennas (St. Croix
and Mauna Kea) where the source was observed at much lower elevations. It
is clear that for a dataset at a higher frequency and/or with sources at lower
elevations a weak target source could well be undetectable without correcting
for the atmosphere.

\begin{figure}[t!]
%\epsscale{0.85}
\plottwo{nodelzn.ps}{delzn.ps}
\caption{Images of same source from the same observation at 8.5\,GHz.
The image of the left
was reduced in the standard phase referencing manner (\ie\ with no
atmospheric/clock correction), the image on the right was reduced as
described in this memo.}
\label{images}
\end{figure}

\bigskip

The authors acknowledge Mark Reid for initiating the creation of \DELZN\
and for useful comments and discussion on both the \AIPS\ task and this memo.
We would also like to thank Ed Fomalont, Greg Taylor and Lorant Sjouwerman
for additional comments and discussion of this memo.

\section{References and Further Reading}

\begin{description}
\item[For discussion of multi-band (group) delay:]
\item[]Thompson, A.R., Moran, J.M., \& Swenson, G.W., 1994, {\it Interferometry and
Synthesis in Radio Astronomy}, Krieger Publishing.
\end{description}
\begin{description}
\item[For general explanation of VLBI data reduction:]
\item[]The \AIPS\ Cookbook, Chapter 9 and/or Appendix C.
\end{description}
\begin{description}
\item[For phase referencing strategies including low frequency:]
\item[]Chatterjee, S., 1999, {\it How accurate is phase referencing at
L-Band, An assessment}, VLBA Scientific Memo 18
\item[]Wrobel, J.M., Walker, R.C. \& Benson, J.M., 2000 {\it Strategies for
Phase Referencing with the VLBA}, VLBA Scientific Memo 20
\item[]Ulvestad, J., 1999, {\it Phase-Referencing Cycle Times},
VLBA Scientific Memo 20
\end{description}
\begin{description}
\item[For atmospheric corrections including the ionosphere:]
\item[]Chatterjee, S., 1999, {\it Recipes for low frequency VLBI
Phase-referencing and GPS Ionospheric Correction}, VLBA Scientific Memo 22
\item[]Sovers, O.J., Fanselow, J.L. \& Jacobs, C.S., 1998, {\it Astrometry
and geodesy with radio interferometry: experiments, models, results},
Reviews of Modern Physics, 70, 1393
\item[]Walker, C. \& Chatterjee, S., 1999, {\it Ionospheric Corrections Using
GPS Based Models}, VLBA Scientific Memo 23
\end{description}

\end{document}
