Initial commit

15 files changed, 4171 insertions, 0 deletions

diff --git a/noao/astutil/doc/airmass.hlp b/noao/astutil/doc/airmass.hlp
new file mode 100644
index 00000000..b35c461e
--- /dev/null
+++ b/noao/astutil/doc/airmass.hlp
@@ -0,0 +1,31 @@
+.help airmass Mar84 noao.astutil
+.ih
+NAME
+airmass -- compute the airmass at a given elevation above horizon
+.ih
+USAGE
+airmass elevation
+.ih
+PARAMETERS
+.ls elevation
+Elevation above horizon in either degrees or radians.
+.le
+.ls scale = 750.0
+Scale factor of the Earth's atmosphere.
+.le
+.ls radians = no
+Input elevation in radians instead of degrees.
+.le
+.ls airmass
+On output, contains the computed airmass.
+.le
+.ih
+EXAMPLE
+Compute the airmass at an elevation of 30 degrees above the horizon
+
+.nf
+	cl> airmass 30
+	airmass 1.996 at an elevation of 30. degrees (0.5236 radians)
+	above horizon
+.fi
+.endhelp
diff --git a/noao/astutil/doc/astcalc.hlp b/noao/astutil/doc/astcalc.hlp
new file mode 100644
index 00000000..f138e8ea
--- /dev/null
+++ b/noao/astutil/doc/astcalc.hlp
@@ -0,0 +1,654 @@
+.help astcalc Jan96 astutil
+.ih
+NAME
+astcalc -- astronomical calculator
+.ih
+USAGE
+astcalc
+.ih
+PARAMETERS
+.ls commands
+A file of commands using the simple syntax given in the DESCRIPTION.  If no
+file name is given then the commands are read interactively from the
+standard input with a prompt given by the \fIprompt\fR parameter.  The
+command input ends with either EOF or "quit".  If a list of images and/or a
+table is specified the commands are repeated for each image or until the
+end of the table is reached.  Comments beginning with '#', blank lines, and
+escaped newlines are allowed.
+.le
+.ls images = ""
+Optional list of images.  The command input is repeated for each image.
+Image header keyword values may be read and used as variables and
+keywords may be created, modified, or deleted provided the image has
+read-write permissions.
+.le
+.ls table = ""
+Optional text file containing columns of values.  The table consists of
+one or more lines of whitespace separated columns of values.  Note that a
+string with whitespace needs to be quoted.  The lines may be scanned and
+the values assigned to variables for use in expressions.  If the command
+input includes reading from the table then the commands will be repeated
+until all the lines in the table have been read.
+.le
+.ls prompt = "astcalc> "
+When no command file is specified the input commands are read from the
+standard input (the terminal) and the value of the \fIprompt\fR string is
+printed as a prompt.  Note that if the input command file is specified as
+"STDIN" there will be no prompt even though commands will also be read from
+the standard input.
+.le
+.ls verbose = no
+Print each variable assignment?  This is useful for debugging command
+files.
+.le
+.ih
+DESCRIPTION
+\fBAstcalc\fR evaluates statements using variables, constants, and
+functions.  Of special interest are many astronomical functions and the
+ability to read and write image header keywords (\fIimages\fR), read
+from a text file (\fItable\fR), and read and write CL parameters.
+
+This task may be used interactively or with input from a command file
+(\fIcommands\fR).  If no command file is specified a prompt (\fIprompt\fR)
+is printed and commands are entered interactively.  The input is terminated
+with either the end-of-file character (EOF) or the command "quit".  Input
+command files simply contain the same input in a file and end with the end
+of the file or "quit".  The input commands, either those entered
+interactively or from a file, are repeated for each image in the image list
+and until the end of the input text table is reached, whichever comes
+first.  The image list and the table are optional and if neither is
+specified the commands are executed just once.
+
+The command input consists of statements with each statement on a
+line by itself.  However long statements may be broken up with
+escaped newlines using the back-slash as the escape character;
+i.e. \<newline>.  Comments beginning with '#', blank lines,
+and whitespace are ignored.
+
+There are three types of statements: assignment, expressions, and
+conditional.  Each statement is on a line by itself though long statements
+may be broken up with escaped newlines (\<newline>).  Assignment statements
+have a variable name, an equal sign, and an expression.  Expression
+statements consist of only the expression with the value of the expression
+being ignored.  Expression statements are generally used with certain
+functions.  Conditional statements are blocks of if-endif and if-else-endif
+with assignment and expression statements between the if-else-endif
+statements.  These may not be nested.
+
+A variable is an arbitrary identifier which must begin with an alphabetic
+character or '$' followed by an alphabetic character and may use alphabetic
+characters, digits, or the characters '_', '$', or '.'.  Other special
+characters may be used but they must be set and referenced with the
+special '@' operator described below.  Lower and upper
+case characters may be used and are considered different characters; i.e.
+identifiers are case sensitive (as are function names).
+
+There are a few special predefined variables: "$D" contains the current
+local date (in new FITS YYYY-MM-DD), "$T" contains the current local
+time, "$GMD" contains the current Greenwich meridian date (in FITS
+YYYY-MM-DD format), "$GMT" contains the current Greenwich meridian time,
+and "$GMDT" contains the current date and time in FITS YYYY-MM-DDTHH:MM:SS
+format.
+
+The expression syntax is described below.  Expressions may use previously
+define variable names, constants (both quoted strings and numeric values),
+and functions.  The functions are given below.  Input from image headers,
+and text files, and CL parameters, and output to image headers is performed
+by I/O functions.
+
+In \fBastcalc\fR variables are maintained internally and input and output
+are performed explicitly by functions.  A related task is \fBasthedit\fR.
+In that task variables are image header keywords and references to keywords
+(assignments, use in expressions, and by themselves with no expression)
+read and write to the image headers.  Updating of the image headers,
+however, can be suppressed.  Also a line of a text table is read
+automatically at the beginning of the command input so that column values
+can be referenced directly.
+
+STATEMENTS
+
+The following gives a more formal description of the statement syntax
+and the special words "if", "else", "endif", and "quit".
+
+.nf
+        <variable> = <expression>
+        <expression>
+        if (<expression>)
+            <statements>
+        endif
+        if (<expression>)
+            <statements>
+        else
+            <statements>
+        endif
+        quit
+.fi
+
+The result of the expression in the "if" statement is normally a logical
+value.  However, a numeric value of 0 is false while any other value is
+true and any string beginning with either "y" or "Y" is true with
+any other value being false; i.e. string values of yes and no may be used.
+
+VARIABLES
+
+Variables may formally be defined as:
+
+.nf
+        [$]{a-zA-Z}[{a-zA-Z0-9._$}]*
+.fi
+
+where [] indicate optional, {} indicates a class, - indicates an
+ASCII range of characters, and * indicates zero or more occurrences.
+Stated in words, a variable must begin with an alphabetic character (ignoring
+an option leading $) and may be followed by any combinations of
+alphabetic, digit, or '.', '_', and '$' characters.
+
+There are a few predefined variables which may be referenced in
+expressions.
+
+.nf
+        $I      The name of the current image (if used)
+        $D      The current date in the YYYY-MM-DD format
+        $T      The current (local) time as a sexagesimal string
+.fi
+
+The date and time are set once at the beginning of execution.
+
+Though not recommended it is possible to use any set of characters
+for a variable provided the variable is referenced as @"<name>".
+For example one could use @"date-obs" to include the character '-'.
+
+EXPRESSIONS
+
+Expressions consist of operands and operators.  The operands may be any
+PREVIOUSLY DEFINED variables, quoted string constants, numeric constants,
+and functions.  Values given as sexagesimal strings are automatically
+converted to decimal numbers.  The operators are arithmetic, logical, and
+string.  The expression syntax is equivalent to that used in the CL and SPP
+languages.
+
+Additional information may be found in the help for \fBhedit\fR except that
+all unquoted nonnumeric strings are considered to be variables and so the
+'(', ')' operators are not used.  The "field" references are not needed so
+the references "." and  "$" are not used and are not legal variable
+names in this task.
+
+operators:
+
+The following operators are recognized in expressions.  With the exception
+of the operators "?" and "?=", the operator set is equivalent to that
+available in the CL and SPP languages.
+
+
+.nf
+        +  -  *  /              arithmetic operators
+        **                      exponentiation
+        //                      string concatenation
+        !  -                    boolean not, unary negation
+        <  <= >  >=             order comparison (works for strings)
+        == != && ||             equals, not equals, and, or
+        ?=                      string equals pattern
+        ? :                     conditional expression
+	@			reference a variable
+.fi
+
+
+The operators "==", "&&", and "||" may be abbreviated as "=", "&", and "|"
+if desired.  The ?= operator performs pattern matching upon strings.
+
+A point to be aware of is that in the ?: conditional expression both
+possible result values are evaluated though the result of the expression
+is only one of them.  This means that one should not use this to
+call I/O functions that one wants to be executed only if a certain
+condition holds.
+
+intrinsic functions:
+
+A number of standard intrinsic functions are recognized within expressions.
+The set of functions currently supported is shown below.
+
+
+.nf
+	abs     atan2   deg     log     min     real    sqrt
+	acos    bool    double  log10   mod     short   str
+	asin    cos     exp     long    nint    sin     tan
+	atan    cosh    int     max     rad     sinh    tanh
+.fi
+
+
+The trigonometric functions operate in units of radians.
+The \fImin\fR and \fImax\fR functions may have any number of arguments up
+to a maximum of sixteen or so (configurable).  The arguments need not all
+be of the same datatype.
+
+A function call may take either of the following forms:
+
+.nf
+        <identifier> '(' arglist ')'
+or
+        <string_expr> '(' arglist ')'
+.fi
+
+The first form is the conventional form found in all programming languages.
+The second permits the generation of function names by string valued
+expressions and might be useful on rare occasions.
+
+astronomical functions:
+
+In addition to the above intrinsic functions there are a number of
+astronomical functions.  More will be added in time.  These are:
+
+.nf
+     sexstr - convert a number to a sexagesimal string (xx:mm:ss.ss)
+      epoch - compute an epoch given a date and time
+     julday - compute a Julian day given a date and time
+        mst - compute a mean sidereal time w/ date, time, and longitude
+ ra_precess - precess ra from one epoch to another
+dec_precess - precess dec from one epoch to another
+    airmass - compute airmass w/ ra, dec, sidereal time, and latitude
+   eairmass - compute effective airmass given
+                ra, dec, sidereal time, exposure time, and latitude
+      obsdb - get parameters from the observatory database
+.fi
+
+.ls sexstr (number), sexstr (number, digits)
+Convert a number to a sexagesimal string in the format X:MM:SS.SS.  There
+is an optional second argument (the default is 0) which is the number of
+decimal digits in the seconds field.
+.le
+.ls epoch (date[, ut])
+Compute an epoch given a date and time.  The date is a string in the
+format DD/MM/YY, YYYY-MM-DD, or YYYY-MM-DDTHH:MM:SS.
+Typically this argument will be the standard FITS
+keyword DATE-OBS.  Because of possible confusion of the hyphen with
+subtraction this keyword would be specified as @"date-obs".  The time
+argument is optional.  If it is not given the time from the date
+string will be used and if absent a time of 0h is used.
+.le
+.ls julday (date[, ut])
+Compute a Julian day given a date and time.  The date and time are
+specified as described previously.
+.le
+.ls mst (date[, ut], longitude)
+Compute a mean sidereal time given a date, time, and longitude in degrees.  The
+date and (optional) time are specified as described previously.  The longitude
+may be given as a constant or using the observatory database function
+as shown in the examples.  The returned value is a sexagesimal
+string with two decimals in the seconds.
+.le
+.ls precess (ra, dec, epoch1, epoch2)
+Precess coordinates from one epoch to another.  The ra is the
+right ascension in hours, the dec in the declination in degrees,
+and the epochs are in years.  This function returns a formatted string with
+the precessed right ascension, declination, and epoch.  Numerical
+values for the right ascension and declination are obtained with the
+functions ra_precess and dec_precess.
+.le
+.ls ra_precess (ra, dec, epoch1, epoch2)
+Precess a right ascension from one epoch to another.  The ra is the
+input right ascension in hours, the dec is the declination in degrees,
+and the epochs are in years.  Because a function can return only one
+value there is a second function to return the precessed declination.
+The returned value is a sexagesimal string with two decimals in the seconds.
+.le
+.ls dec_precess (ra1, dec1, epoch1, epoch2)
+Precess a declination from one epoch to another.  The ra is the
+input right ascension in hours, the dec is the declination in degrees,
+and the epochs are in years.  Because a function can return only one
+value there is a second function to return the precessed right ascension.
+The returned value is a sexagesimal string with two decimals in the seconds.
+.le
+.ls arcsep (ra1, dec1, ra2, dec2)
+Compute the separation between two spherical coordinates.  The parameters
+ra1 and ra2 are coordinates in hours (right ascension, longitude, etc.)
+and the dec1 and dec2 parameters are coordinates in degrees (declination,
+latitude, etc.).  The computed value is returned in seconds of arc.
+.le
+.ls airmass (ra, dec, st, latitude)
+Compute an airmass given right ascension in hours, declination in
+degrees, sidereal time in hours, and latitude in degrees.  The latitude
+is often specified using the observatory database function as shown
+in the examples.
+.le
+.ls eairmass (ra, dec, st, exptime, latitude)
+Compute an "effective" airmass given right ascension in hours, declination
+in degrees, beginning sidereal time in hours, exposure time in seconds, and
+latitude in degrees.  The The latitude is often specified using the
+observatory database function as shown in the examples.  The effective
+airmass is based on a Simpson's rule weighting of the beginning, middle,
+and ending airmass (with no provision for paused exposure).  The weights
+are:
+
+.nf
+    effective = beginning + 4 * middle + ending
+.fi
+.le
+.ls obsdb (observatory, parameter)
+Return a value from the observatory database.  The observatory parameter is
+a observatory identification string as defined in the database.  Another
+special value is "observatory" which then follows a name resolution
+scheme.  The observatory database mechanism is described by the help topic
+\fBobservatory\fR.  The parameter is a string given the quantity desired.
+Typically this would be "longitude" or "latitude" but there are other
+possible parameters.
+.le
+
+input/output functions:
+
+There are special functions for formatting, printing, error aborts,
+reading, writing, and deleting image header keywords, reading a text file,
+and reading and writing CL parameters.
+
+.nf
+     print  - print a set of arguments with default format
+     printf - print a set arguments with specified format
+     format - format a string
+     error  - print an error message and abort
+     clget  - get a value from a CL parameter
+     clput  - put a value to a CL parameter
+      scan  - scan a string and parse into variables
+     fscan  - scan a line of a text file
+     imget  - get the value of an image header keyword
+     imput  - put (add or modify) the value of an image header keyword
+     imdel  - delete an image header keyword
+.fi
+
+.ls print ([argument, ...])
+Print the arguments with default formats based on the type of value ending
+with a newline.  There may be zero or more arguments.  With zero arguments
+only a newline will be printed.
+.le
+.ls printf (fmt [, argument, ...])
+Print a list of arguments using the formatting syntax described later.
+Parameters to be formatted are given by the % fields and the values are
+passed as further arguments in the order in which they are referenced.
+There is no automatic newline so the format must include "\n" to
+produce newlines.
+.le
+.ls error (message)
+Print the "message", which can be any string variable such as might
+be produced by "format", and abort the task.  This is useful in
+conjunction with the conditional operator to abort if a variable
+takes an inappropriate value.
+.le
+.ls clget (parameter)
+Get the value of a CL parameter.  The argument must be a string.  The
+function value is the value of the parameter.
+.le
+.ls clput (parameter, value)
+Put a value into a CL parameter.  The parameter argument must be a
+string and the value can be anything.  The function returns a string
+of the form "clput: parameter = value" where parameter and value are
+the actual values.
+.le
+.ls scan (string, var, ...)
+Parse a string of whitespace separated words into a list of
+variables.  The number of variables assigned is
+the returned value of the function.
+.le
+.ls fscan (var, ...)
+Scan a line of a text file into a list of variables.  The arguments
+are zero or more variable names to which to assign the values of
+the whitespace separated fields.  The number of variables assigned
+is the returned value of the function.
+.le
+.ls imget (parameter)
+Get the value of an image header keyword from the current image.  The
+argument must be a string.  The function value is the value of the keyword.
+.le
+.ls imput (parameter, value)
+Put a value into an image header keyword for the current image.  The
+parameter argument must be a string and the value can be anything.  If the
+keyword exists it will be modified and if it does not exist it will be
+added.  The function returns a string of the form "imput: parameter =
+value" for new keywords or "imput: parameter = old_value -> value" for
+modified keywords where parameter and value are the actual values.
+.le
+.ls imdel (parameter)
+Delete an image header keyword.  The parameter argument must be a string.
+The returned values are the strings "imdel: parameter not found"
+or "imdel: parameter = value (DELETED)" where parameter is the parameter
+name and value is the old value.
+.le
+
+.ih
+FORMATS
+A  format  specification has the form "%w.dCn", where w is the field
+width, d is the number of decimal places or the number of digits  of
+precision,  C  is  the  format  code,  and  n is radix character for
+format code "r" only.  The w and d fields are optional.  The  format
+codes C are as follows:
+    
+.nf
+b       boolean (YES or NO)
+c       single character (c or '\c' or '\0nnn')
+d       decimal integer
+e       exponential format (D specifies the precision)
+f       fixed format (D specifies the number of decimal places)
+g       general format (D specifies the precision)
+h       hms format (hh:mm:ss.ss, D = no. decimal places)
+m       minutes, seconds (or hours, minutes) (mm:ss.ss)
+o       octal integer
+rN      convert integer in any radix N
+s       string (D field specifies max chars to print)
+t       advance To column given as field W
+u       unsigned decimal integer 
+w       output the number of spaces given by field W
+x       hexadecimal integer
+z       complex format (r,r) (D = precision)
+    
+    
+Conventions for w (field width) specification:
+    
+    W =  n      right justify in field of N characters, blank fill
+        -n      left justify in field of N characters, blank fill
+        0n      zero fill at left (only if right justified)
+absent, 0       use as much space as needed (D field sets precision)
+    
+    
+Escape sequences (e.g. "\n" for newline):
+    
+\b      backspace   (not implemented)
+\f      formfeed
+\n      newline (crlf)
+\r      carriage return
+\t      tab
+\"      string delimiter character
+\'      character constant delimiter character
+\\      backslash character
+\nnn    octal value of character
+    
+Examples
+    
+%s          format a string using as much space as required
+%-10s       left justify a string in a field of 10 characters
+%-10.10s    left justify and truncate a string in a field of 10 characters
+%10s        right justify a string in a field of 10 characters
+%10.10s     right justify and truncate a string in a field of 10 characters
+    
+%7.3f       print a real number right justified in floating point format
+%-7.3f      same as above but left justified
+%15.7e      print a real number right justified in exponential format
+%-15.7e     same as above but left justified
+%12.5g      print a real number right justified in general format
+%-12.5g     same as above but left justified
+
+%h          format as nn:nn:nn.n
+%15h        right justify nn:nn:nn.n in field of 15 characters
+%-15h       left justify nn:nn:nn.n in a field of 15 characters
+%12.2h      right justify nn:nn:nn.nn
+%-12.2h     left justify nn:nn:nn.nn
+    
+%H          / by 15 and format as nn:nn:nn.n
+%15H        / by 15 and right justify nn:nn:nn.n in field of 15 characters
+%-15H       / by 15 and left justify nn:nn:nn.n in field of 15 characters
+%12.2H      / by 15 and right justify nn:nn:nn.nn
+%-12.2H     / by 15 and left justify nn:nn:nn.nn
+
+\n          insert a newline
+.fi
+
+.ih
+EXAMPLES
+1.  This example shows interactive use.
+
+.nf
+    cl> astcalc
+    astcalc> print ((1 + 2 + 3) / 2 - 2 * 2)
+    -1
+    astcalc> observatory = "kpno"
+    astcalc> date = "05/04/87"
+    astcalc> ut = 9:27:27
+    astcalc> ra = 13:29:24
+    astcalc> dec = 47:15:34
+    astcalc> epoch = epoch (date, ut)
+    astcalc> mst = mst (date, ut, obsdb (observatory, "longitude"))
+    astcalc> print (epoch)
+    1987.257752395672
+    astcalc> print (mst)
+    14:53:39.81
+    astcalc> print (julday (date, ut))
+    2446890.894062519
+    astcalc> print (ra_precess (ra, dec, epoch, 1950))
+    13:27:49.84
+    astcalc> print (dec_precess (ra, dec, epoch, 1950))
+    47:27:05.72
+    astcalc> print (airmass (ra, dec, mst, obsdb (observatory, "latitude")))
+    1.07968417231416
+    astcalc> printf ("Hello World: %s\n", precess (ra, dec, epoch, 1950))
+    Hello World: 13:27:49.84  47:27:05.7   1950.
+    astcalc> quit
+.fi
+
+2.  This example shows the same commands as in the previous example
+read from a file.
+
+.nf
+    cl> type example2.dat
+    # Define variables.
+    observatory = "kpno"
+    date = "05/04/87"
+    ut = 9:27:27
+    ra = 13:29:24
+    dec = 47:15:34
+    epoch = epoch (date, ut)
+    mst = mst (date, ut, obsdb (observatory, "longitude"))
+
+    # Print results of some expressions.
+    print ((1 + 2 + 3) / 2 - 2 * 2)       # Calculation with constants
+    print (epoch)                         # Print variable
+    print (mst)                           # Print variable
+    print (julday (date, ut))             # Print result of function
+    print (ra_precess (ra, dec, epoch, 1950))
+    print (dec_precess (ra, dec, epoch, 1950))
+    print (airmass (ra, dec, mst, obsdb (observatory, "latitude")))
+
+    # Formatted print with arguments.  Note newline.
+    printf ("Hello World: %s\n", precess (ra, dec, epoch, 1950))
+    cl> astcalc commands=example2.dat
+    -1
+    1987.257752395672
+    14:53:39.81
+    2446890.894062519
+    13:27:49.84
+    47:27:05.72
+    1.07968417231416
+    Hello World: 13:27:49.84  47:27:05.7   1950.
+.fi
+
+3.  This example precesses coordinates given in a text file.
+
+.nf
+    cl> type example3.dat,table.dat
+    ===> example3.dat <===
+    # Read table of RA, DEC, and optional EPOCH and precess to 2000.
+
+    epoch = 1900            # Default input epoch
+    epoch1 = 2000           # Precession epoch
+
+    # Scan table and precess coordinates.
+    if (fscan ("ra", "dec", "epoch") >= 2)
+	ra1 = ra_precess (ra, dec, epoch, epoch1)
+	dec1 = dec_precess (ra, dec, epoch, epoch1)
+	printf ("%h %h %d -> %h %h %d\n", ra, dec, epoch, ra1, dec1, epoch1)
+    else
+	printf ("Missing coordinates\n")
+    endif
+
+    ===> table.dat <===
+    12:22:31        31:10:15        1950
+    13:52:44        10:21:32        1996.1
+    14:52:44        11:21:32
+    10:20:30
+
+    cl> astcalc commands=example3.dat table=table.dat
+    12:22:31.0 31:10:15.0 1950 -> 12:25:00.56 30:53:38.13 2000
+    13:52:44.0 10:21:32.0 1996 -> 13:52:55.54 10:20:23.11 2000
+    14:52:44.0 11:21:32.0 1900 -> 14:57:33.16 10:57:24.74 2000
+    Missing coordinates
+.fi
+
+4.  This complex example illustrates reading from CL parameters and
+image header keywords.  It precesses coordinates to a standard epoch
+and computes the arc separation between the coordinates and a center
+coordinate.  If the separation is less than a specified amount it
+prints the image name and additional information.  This is the
+data file for the \fBastradius\fR script task.
+
+.nf
+    cl> type astutil$astradius.dat
+    # Print images which are within a given radius in the sky.
+
+    # Get parameters.
+    racenter = clget ("astradius.racenter")
+    deccenter = clget ("astradius.deccenter")
+    epcenter = clget ("astradius.epcenter")
+    radius = clget ("astradius.radius")
+    ra = imget(clget("keywpars.ra"))
+    dec = imget(clget("keywpars.dec"))
+
+    epoch = imget(clget("keywpars.epoch"))
+    if (str(epoch) == "" || real(epoch) == 0.)
+	date = imget(clget("keywpars.date_obs"))
+	ut = imget(clget("keywpars.ut"))
+	epoch = epoch (date, ut)
+    endif
+
+    # Precess image coordinates to center epoch and compute separation.
+    radec = precess (ra, dec, epoch, epcenter)
+    ra1 = ra_precess (ra, dec, epoch, epcenter)
+    dec1 = dec_precess (ra, dec, epoch, epcenter)
+    sep = arcsep (racenter, deccenter, ra1, dec1)
+
+    # Print result if within radius.
+    if (sep < real (radius))
+	printf ("%-15s %s %4d %s\n", $I, radec, sep, imget ("title"))
+    endif
+    cl> astcalc commands=astutil$astradius.dat images=dev$pix
+    RA center (hours) (13:31): 
+    DEC center (degrees) (47:00): 
+    Epoch of center (2000.): 
+    Radius in arc seconds (3600.): 
+    dev$pix         13:29:56.16  47:11:37.9   2000.  955 m51  B  600s
+.fi
+
+.ih
+REVISIONS
+.ls ASTCALC V2.15
+The $D variable was changed from the old MM/DD/YY format to the post-Y2K
+YYYY-MM-DD format.
+.le
+.ls ASTCALC V2.11.2
+Y2K update:  The epoch, julday, and mst functions now take either the old
+or new FITS style date strings.  The time argument is optional and if
+it is not specified the time from the date string is used and if neither
+time is present a value of 0h is used.  New internal variables $GMD,
+$GMT, and $GMDT for the current time Greenwich time are defined.
+.le
+.ls ASTCALC V2.11
+This task is new in this release.
+.le
+.ih
+SEE ALSO
+astradius, asthedit, setairmass, setjd, asttimes, precess, observatory, hedit
+.endhelp
diff --git a/noao/astutil/doc/asthedit.hlp b/noao/astutil/doc/asthedit.hlp
new file mode 100644
index 00000000..2806cf58
--- /dev/null
+++ b/noao/astutil/doc/asthedit.hlp
@@ -0,0 +1,789 @@
+.help asthedit Jan96 astutil
+.ih
+NAME
+asthedit -- astronomical header editor
+.ih
+USAGE
+asthedit images commands
+.ih
+PARAMETERS
+.ls images
+List of images to be used.  The image header keywords are used in this task
+as variables which are read, modified, created, or deleted.  If the images
+do not have write permission or the \fIupdate\fR parameter is "no" then the
+image headers will not be modified.  If no images are specified then this
+task can be used as a calculator (though see \fBastcalc\fR).
+.le
+.ls commands
+A file of commands using the simple syntax given in the DESCRIPTION.  If no
+file name is given then the commands are read interactively from the
+standard input with a prompt given by the \fIprompt\fR parameter.  The
+command input ends with either EOF or "quit".  If a list of images and/or a
+table is specified the commands are repeated for each image or until the
+end of the table is reached.  Comments beginning with '#', blank lines, and
+escaped newlines are allowed.
+.le
+.ls table = ""
+Optional text file containing columns of values.  The table consists of
+one or more lines of whitespace separated columns of values.  Note that a
+string with whitespace needs to be quoted.  One line of the table is
+scanned for each image.  There must be at lest as many fields as are
+defined by the column names.
+.le
+.ls colnames = ""
+List of whitespace separated column names.  These are the names referenced
+in the command file by $<name>.  The leading '$' is not included in the
+column name specification.  There may be fewer columns than the number of
+columns in the table.  Dummy names must be used if some columns occur
+before a column to be referenced.
+.le
+.ls prompt = "asthedit> "
+When no command file is specified the input commands are read from the
+standard input (the terminal) and the value of the \fIprompt\fR string is
+printed as a prompt.  Note that if the input command file is specified as
+"STDIN" there will be no prompt even though commands will also be read from
+the standard input.
+.le
+.ls update = yes
+Update the image headers?  If no then any new, modified, or deleted
+keywords will not be recorded in the image headers.  This allows using the
+task only for computing and printing quantities.  Also this allows
+accessing read-only images.
+.le
+.ls verbose = no
+Print each keyword added or modified?
+.le
+.ls oldstyle = no
+Use the old style syntax of this task from versions prior to V2.11.  This
+parameter allows backward compatibility for command files previously
+developed.  Some aspects of the new syntax are still available.
+.le
+.ih
+DESCRIPTION
+\fBAsthedit\fR evaluates expressions using image header keywords, column
+names from a text table, CL parameters, internal variables, constants, and
+functions to create or modify image header keywords.  This task is
+particularly useful for adding keywords from a table and deriving keywords
+used by IRAF tasks which are not present in the images.  It differs from
+\fBhedit\fR in that it includes astronomical functions, operates from a
+command file which may perform many edits, and references columns from a
+text table.  The command file may be omitted in which case commands may be
+entered interactively for the first image and then the same commands will
+be repeated for any subsequent images.
+
+This task may be used interactively or with input from a command file
+(\fIcommands\fR).  If no command file is specified a prompt (\fIprompt\fR)
+is printed and commands are entered interactively.  The input is terminated
+with either the end-of-file character (EOF) or the command "quit".  Input
+command files simply contain the same input in a file and end with the end
+of the file or "quit".  The input commands, either those entered
+interactively or from a file, are repeated for each image in the image list
+and until the end of the input text table is reached, whichever comes
+first.  Generally this task is used on one or more images but if no
+image is specified the commands are executed just once and task behaves
+like an calculator.
+
+The command input consists of statements with each statement on a
+line by itself.  However long statements may be broken up with
+escaped newlines using the back-slash as the escape character;
+i.e. \<newline>.  Comments beginning with '#', blank lines,
+and whitespace are ignored.
+
+There are three types of statements: assignment, expressions, and
+conditional.  Each statement is on a line by itself though long statements
+may be broken up with escaped newlines (\<newline>).  Assignment statements
+have an image header keyword name (or variable name beginning with $), an
+equal sign (but see the \fIoldstyle\fR parameter), and an expression.
+Expression statements consist of only the expression with the value of the
+expression being ignored.  Expression statements are generally used with
+certain functions.  Conditional statements are blocks of if-endif and
+if-else-endif with assignment and expression statements between the
+if-else-endif statements.  These may not be nested.
+
+In earlier versions of this task there were only assignment statements
+and these did not use an equal sign; i.e. all statements consisted
+of an image header keyword and an expression separated by whitespace
+except that a keyword name by itself indicates deletion of a keyword.
+In order to interpret old command files the \fIoldstyle\fR parameter
+may be set to yes.  This will insert an equal sign internally.  It
+also only allows a subset of statements to not begin with a keyword
+or variable.  These are if, else, endif, print, printf, and quit.
+Note that with the old style syntax one may still include an equal
+sign.  It is recommended that the old style syntax not be used because
+of the greater flexibility in the new syntax.
+
+An image header keyword name is an arbitrary identifier which must begin
+with an alphabetic character or '$' followed by an alphabetic character and
+may use alphabetic characters, digits, or the characters '_', '$', or '.'.
+Keyword names are case insensitive.  Because some additional characters are
+allowed in the FITS definition of keyword names, such names may be
+referenced with the special '@' operator described below.
+
+One may also use internal variables which have the same identifier rules
+but begin with '$'.  Note that these variables are case sensitive (as are
+function names).  There are a few special predefined variables: "$I"
+contains the current image name, "$D" contains the current local date (in
+old FITS DD/MM/YY format), "$T" contains the current local time, "$GMD"
+contains the current Greenwich meridian date (in FITS YYYY-MM-DD format),
+"$GMT" contains the current Greenwich meridian time, and "$GMDT" contains
+the current date and time in FITS YYYY-MM-DDTHH:MM:SS format.
+
+Before the commands are interpreted for each image a line of a text
+file may be read.  This occurs when a file is specified by the
+\fItable\fR parameter.  The line is scanned and the values of each
+column are stored in the variable names specified by the \fIcolnames\fR
+parameter.  The values may be referenced in expressions by the
+specified column name preceded with '$'.  Note that additional lines
+may be scanned with the "fscan" function.  The user is then responsible
+for the table containing the correct sequence of lines when there
+are multiple images.
+
+In \fBasthedit\fR identifiers are image header keywords and lines
+for the table file are read automatically.  A related task is \fBastcalc\fR.
+In this task all variables are maintained internally and input and output
+are performed explicitly by functions.  There are functions to read,
+write, and delete image header keywords from a list of images.
+
+STATEMENTS
+
+The following gives a more formal description of the statement syntax
+and the special words "if", "else", "endif", and "quit".
+
+.nf
+	<keyword>
+        <keyword> = <expression>
+	$<variable> = <expression> 
+        <expression>
+        if (<expression>)
+            <statements>
+        endif
+        if (<expression>)
+            <statements>
+        else
+            <statements>
+        endif
+        quit
+.fi
+
+The result of the expression in the "if" statement is normally a logical
+value.  However, a numeric value of 0 is false while any other value is
+true and any string beginning with either "y" or "Y" is true with
+any other value being false; i.e. string values of yes and no may be used.
+
+The old style syntax allows the following statements.
+
+.nf
+	<keyword>
+        <keyword>     <expression>
+	$<variable>   <expression> 
+        <keyword> = <expression>
+	$<variable> = <expression> 
+        print (...)
+        printf (...)
+        if (<expression>)
+            <statements>
+        endif
+        if (<expression>)
+            <statements>
+        else
+            <statements>
+        endif
+        quit
+.fi
+
+Old style command files would only use the first two statements.
+
+KEYWORD NAMES AND VARIABLES
+
+Keyword names and variables may formally be defined as:
+
+.nf
+        [$]{a-zA-Z}[{a-zA-Z0-9._$}]*
+.fi
+
+where [] indicate optional, {} indicates a class, - indicates an ASCII
+range of characters, and * indicates zero or more occurrences.  In words, a
+keyword must begin with an alphabetic character, a variable or text file
+column name begins with '$' and an alphabetic character, and both may be
+followed by any combinations of alphabetic, digit, or '.', '_', and '$'
+characters.
+
+There are a few predefined variables which may be referenced in
+expressions.
+
+.nf
+        $I      The name of the current image (if used)
+        $D      The current date in the DD/MM/YY format
+        $T      The current (local) time as a sexagesimal string
+.fi
+
+The date and time are set once at the beginning of execution.
+
+Though not recommended it is possible to use any set of characters
+for a variable provided the variable is referenced as @"<name>".
+For example one could use @"date-obs" to include the character '-'.
+This option is primarily used for FITS keywords that use '-' as
+a hyphen character and must be escaped from interpretation as the
+an arithmetic subtraction operator.
+
+EXPRESSIONS
+
+Expressions consist of operands and operators.  The operands may be any
+image header keyword, previously defined variable, column name, quoted
+string constants, numeric constants, and functions.  Values given as
+sexagesimal strings are automatically converted to decimal numbers.  The
+operators are arithmetic, logical, and string.  The expression syntax is
+equivalent to that used in the CL and SPP languages.
+
+Additional information may be found in the help for \fBhedit\fR except that
+all unquoted nonnumeric strings are considered to be keywords or variables
+ and so the '(', ')' operators are not used.  The "field" references are
+not needed so the references "." and  "$" are not used and are not legal
+variable names in this task.
+
+operators:
+
+The following operators are recognized in expressions.  With the exception
+of the operators "?", "?=", and "@", the operator set is equivalent to that
+available in the CL and SPP languages.
+
+
+.nf
+        +  -  *  /              arithmetic operators
+        **                      exponentiation
+        //                      string concatenation
+        !  -                    boolean not, unary negation
+        <  <= >  >=             order comparison (works for strings)
+        == != && ||             equals, not equals, and, or
+        ?=                      string equals pattern
+        ? :                     conditional expression
+        @                       reference a variable
+.fi
+
+
+The operators "==", "&&", and "||" may be abbreviated as "=", "&", and "|"
+if desired.  The ?= operator performs pattern matching upon strings.
+The @ operator is required to reference keywords with
+one of the operator characters.  This is most like to be used as:
+
+        @"date-obs"
+
+A point to be aware of is that in the ?: conditional expression both
+possible result values are evaluated though the result of the expression
+is only one of them.  This means that one should not use this to
+call I/O functions that one wants to be executed only if a certain
+condition holds.
+
+intrinsic functions:
+
+A number of standard intrinsic functions are recognized within expressions.
+The set of functions currently supported is shown below.
+
+
+.nf
+	abs     atan2   deg     log     min     real    sqrt
+	acos    bool    double  log10   mod     short   str
+	asin    cos     exp     long    nint    sin     tan
+	atan    cosh    int     max     rad     sinh    tanh
+.fi
+
+
+The trigonometric functions operate in units of radians.
+The \fImin\fR and \fImax\fR functions may have any number of arguments up
+to a maximum of sixteen or so (configurable).  The arguments need not all
+be of the same datatype.
+
+A function call may take either of the following forms:
+
+.nf
+        <identifier> '(' arglist ')'
+or
+        <string_expr> '(' arglist ')'
+.fi
+
+The first form is the conventional form found in all programming languages.
+The second permits the generation of function names by string valued
+expressions and might be useful on rare occasions.
+
+special functions:
+
+In addition to the above intrinsic functions there are a number of
+astronomical functions. More will be added in time.  These are:
+
+.nf
+     sexstr - convert a number to a sexagesimal string (xx:mm:ss.ss)
+      epoch - compute an epoch given a date and time
+     julday - compute a Julian day given a date and time
+        mst - compute a mean sidereal time given a date, time, and longitude
+ ra_precess - precess ra from one epoch to another
+dec_precess - precess dec from one epoch to another
+    airmass - compute airmass given ra, dec, sidereal time, and latitude
+   eairmass - compute effective airmass given
+		ra, dec, sidereal time, exposure time, and latitude
+      obsdb - get parameters from the observatory database
+.fi
+
+.ls sexstr (number), sexstr (number, digits)
+Convert a number to a sexagesimal string in the format X:MM:SS.SS.  There
+is an optional second argument (the default is 0) which is the number of
+decimal digits in the seconds field.
+.le
+.ls epoch (date[, ut])
+Compute an epoch given a date and time.  The date is a string in the
+format DD/MM/YY, YYYY-MM-DD, or YYYY-MM-DDTHH:MM:SS.
+Typically this argument will be the standard FITS
+keyword DATE-OBS.  Because of possible confusion of the hyphen with
+subtraction this keyword would be specified as @"date-obs".  The time
+argument is optional.  If it is not given the time from the date
+string will be used and if absent a time of 0h is used.
+.le
+.ls julday (date[, ut])
+Compute a Julian day given a date and time.  The date and time are
+specified as described previously.
+.le
+.ls mst (date[, ut], longitude)
+Compute a mean sidereal time given a date, time, and longitude in degrees.  The
+date and (optional) time are specified as described previously.  The longitude
+may be given as a constant or using the observatory database function
+as shown in the examples.  The returned value is a sexagesimal
+string with two decimals in the seconds.
+.le
+.ls precess (ra, dec, epoch1, epoch2)
+Precess coordinates from one epoch to another.  The ra is the
+right ascension in hours, the dec in the declination in degrees,
+and the epochs are in years.  This function returns a formatted string with
+the precessed right ascension, declination, and epoch.  Numerical
+values for the right ascension and declination are obtained with the
+functions ra_precess and dec_precess.
+.le
+.ls ra_precess (ra, dec, epoch1, epoch2)
+Precess a right ascension from one epoch to another.  The ra is the
+input right ascension in hours, the dec is the declination in degrees,
+and the epochs are in years.  Because a function can return only one
+value there is a second function to return the precessed declination.
+The returned value is a sexagesimal string with two decimals in the seconds.
+.le
+.ls dec_precess (ra1, dec1, epoch1, epoch2)
+Precess a declination from one epoch to another.  The ra is the
+input right ascension in hours, the dec is the declination in degrees,
+and the epochs are in years.  Because a function can return only one
+value there is a second function to return the precessed right ascension.
+The returned value is a sexagesimal string with two decimals in the seconds.
+.le
+.ls arcsep (ra1, dec1, ra2, dec2)
+Compute the separation between two spherical coordinates.  The parameters
+ra1 and ra2 are coordinates in hours (right ascension, longitude, etc.)
+and the dec1 and dec2 parameters are coordinates in degrees (declination,
+latitude, etc.).  The computed value is returned in seconds of arc.
+.le
+.ls airmass (ra, dec, st, latitude)
+Compute an airmass given right ascension in hours, declination in
+degrees, sidereal time in hours, and latitude in degrees.  The latitude
+is often specified using the observatory database function as shown
+in the examples.
+.le
+.ls eairmass (ra, dec, st, exptime, latitude)
+Compute an "effective" airmass given right ascension in hours, declination
+in degrees, beginning sidereal time in hours, exposure time in seconds, and
+latitude in degrees.  The The latitude is often specified using the
+observatory database function as shown in the examples.  The effective
+airmass is based on a Simpson's rule weighting of the beginning, middle,
+and ending airmass (with no provision for paused exposure).  The weights
+are:
+
+.nf
+    effective = beginning + 4 * middle + ending
+.fi
+.le
+.ls obsdb (observatory, parameter)
+Return a value from the observatory database.  The observatory parameter is
+a observatory identification string as defined in the database.  Often this
+is the value stored in the OBSERVAT keyword.  Another special value is
+"observatory" which then follows a name resolution scheme.  The observatory
+database mechanism is described by the help topic \fBobservatory\fR.  The
+parameter is a string given the quantity desired.  Typically this would be
+"longitude" or "latitude" but there are other possible parameters.
+.le
+
+input/output functions:
+
+There are special functions for formatting, printing, error aborts,
+reading, writing, and deleting image header keywords, reading a text file,
+and reading and writing CL parameters.  Note that in \fBasthedit\fR
+one would not normally use the image input/output functions or
+the text file scanning function since any keyword reference reads or
+writes to the image header and one line of the text file is scanned
+automatically for each image.
+
+.nf
+     print  - print a set of arguments with default format
+     printf - print a set arguments with specified format
+     format - format a string
+     error  - print an error message and abort
+     clget  - get a value from a CL parameter
+     clput  - put a value to a CL parameter
+     scan   - scan a string and parse into keywords or variables
+     fscan  - scan a line of a text file
+     imget  - get the value of an image header keyword
+     imput  - put (add or modify) the value of an image header keyword
+     imdel  - delete an image header keyword
+.fi
+
+.ls print ([argument, ...])
+Print the arguments with default formats based on the type of value ending
+with a newline.  There may be zero or more arguments.  With zero arguments
+only a newline will be printed.
+.le
+.ls printf (fmt [, argument, ...])
+Print a list of arguments using the formatting syntax described later.
+Parameters to be formatted are given by the % fields and the values are
+passed as further arguments in the order in which they are referenced.
+There is no automatic newline so the format must include "\n" to
+produce newlines.
+.le
+.ls error (message)
+Print the "message", which can be any string variable such as might
+be produced by "format", and abort the task.  This is useful in
+conjunction with the conditional operator to abort if a variable
+takes an inappropriate value.
+.le
+.ls clget (parameter)
+Get the value of a CL parameter.  The argument must be a string.  The
+function value is the value of the parameter.
+.le
+.ls clput (parameter, value)
+Put a value into a CL parameter.  The parameter argument must be a
+string and the value can be anything.  The function returns a string
+of the form "clput: parameter = value" where parameter and value are
+the actual values.
+.le
+.ls scan (string, var, ...)
+Parse a string of whitespace separated words into a list of
+keywords or variables.  The number of variables assigned is
+the returned value of the function.
+.le
+.ls fscan (var, ...)
+Scan a line of a text file into a list of keywords or variables.  The arguments
+are zero or more variable names to which to assign the values of
+the whitespace separated fields.  The number of variables assigned
+is the returned value of the function.
+.le
+.ls imget (parameter)
+Get the value of an image header keyword from the current image.  The
+argument must be a string.  The function value is the value of the keyword.
+.le
+.ls imput (parameter, value)
+Put a value into an image header keyword for the current image.  The
+parameter argument must be a string and the value can be anything.  If the
+keyword exists it will be modified and if it does not exist it will be
+added.  The function returns a string of the form "imput: parameter =
+value" for new keywords or "imput: parameter = old_value -> value" for
+modified keywords where parameter and value are the actual values.
+.le
+.ls imdel (parameter)
+Delete an image header keyword.  The parameter argument must be a string.
+The returned values are the strings "imdel: parameter not found"
+or "imdel: parameter = value (DELETED)" where parameter is the parameter
+name and value is the old value.
+.le
+
+.ih
+FORMATS
+A  format  specification has the form "%w.dCn", where w is the field
+width, d is the number of decimal places or the number of digits  of
+precision,  C  is  the  format  code,  and  n is radix character for
+format code "r" only.  The w and d fields are optional.  The  format
+codes C are as follows:
+    
+.nf
+b       boolean (YES or NO)
+c       single character (c or '\c' or '\0nnn')
+d       decimal integer
+e       exponential format (D specifies the precision)
+f       fixed format (D specifies the number of decimal places)
+g       general format (D specifies the precision)
+h       hms format (hh:mm:ss.ss, D = no. decimal places)
+m       minutes, seconds (or hours, minutes) (mm:ss.ss)
+o       octal integer
+rN      convert integer in any radix N
+s       string (D field specifies max chars to print)
+t       advance To column given as field W
+u       unsigned decimal integer 
+w       output the number of spaces given by field W
+x       hexadecimal integer
+z       complex format (r,r) (D = precision)
+    
+    
+Conventions for w (field width) specification:
+    
+    W =  n      right justify in field of N characters, blank fill
+        -n      left justify in field of N characters, blank fill
+        0n      zero fill at left (only if right justified)
+absent, 0       use as much space as needed (D field sets precision)
+    
+    
+Escape sequences (e.g. "\n" for newline):
+    
+\b      backspace   (not implemented)
+\f      formfeed
+\n      newline (crlf)
+\r      carriage return
+\t      tab
+\"      string delimiter character
+\'      character constant delimiter character
+\\      backslash character
+\nnn    octal value of character
+    
+Examples
+    
+%s          format a string using as much space as required
+%-10s       left justify a string in a field of 10 characters
+%-10.10s    left justify and truncate a string in a field of 10 characters
+%10s        right justify a string in a field of 10 characters
+%10.10s     right justify and truncate a string in a field of 10 characters
+    
+%7.3f       print a real number right justified in floating point format
+%-7.3f      same as above but left justified
+%15.7e      print a real number right justified in exponential format
+%-15.7e     same as above but left justified
+%12.5g      print a real number right justified in general format
+%-12.5g     same as above but left justified
+
+%h	    format as nn:nn:nn.n
+%15h	    right justify nn:nn:nn.n in field of 15 characters
+%-15h	    left justify nn:nn:nn.n in a field of 15 characters
+%12.2h	    right justify nn:nn:nn.nn
+%-12.2h	    left justify nn:nn:nn.nn
+    
+%H	    / by 15 and format as nn:nn:nn.n
+%15H	    / by 15 and right justify nn:nn:nn.n in field of 15 characters
+%-15H	    / by 15 and left justify nn:nn:nn.n in field of 15 characters
+%12.2H	    / by 15 and right justify nn:nn:nn.nn
+%-12.2H	    / by 15 and left justify nn:nn:nn.nn
+
+\n          insert a newline
+.fi
+
+.ih
+EXAMPLES
+1.  The following command file exercises the astronomical functions:
+
+.nf
+    cl> type cmds
+    observat = "kpno"
+    time = sexstr (1.2345)
+    epoch = epoch (@'date-obs', ut)
+    jd = julday (@'date-obs', ut)
+    mst = mst (@'date-obs', ut, obsdb (observat, "longitude"))
+    rap = ra_precess (ra, dec, epoch, 1950)
+    dap = dec_precess (ra, dec, epoch, 1950)
+    airmass =  airmass (ra, dec, mst, obsdb (observat, "latitude"))
+    airmass
+    airmass = " "
+    airmass = eairmass (ra, dec, mst, itime, obsdb (observat, "latitude"))
+    cl> imhead obj001 l+
+        ...
+	DATE-OBS= '05/04/87'            /  DATE DD/MM/YY
+	RA      = '13:29:24.00'         /  RIGHT ASCENSION
+	DEC     = '47:15:34.00'         /  DECLINATION
+	UT      = ' 9:27:27.00'         /  UNIVERSAL TIME
+	ITIME   =                  600  /  REQUESTED INTEGRATION TIME (SECS)
+        ...
+    cl> asthedit obj001 cmds table="" verbose+
+    obj001:
+      $I = pix
+      $D = 22/01/96
+      $T = 19:14:41
+      observat = kpno
+      time = 1:14:04
+      epoch = 1987.257752395672
+      jd = 2446890.894062519
+      mst = 14:53:39.81
+      rap = 13:27:49.84
+      dap = 47:27:05.72
+      airmass = 1.079684154511483
+      airmass = 1.07968415451148 -> DELETED
+      airmass =  
+      airmass =  -> 1.08519059292424
+.fi
+
+Note the use of the keyword deletion and syntax for adding an empty
+value.
+
+2.  The following command file shows computing a mid-ut and using a table
+of values.
+
+.nf
+    cl> type cmds
+    midut = sexstr ($ut + $itime/3600./2.)
+    imagetyp = $imagetyp
+    cl> type table
+    object	9:27:27		600
+    comp	9:48:00		10
+    object	9:49:00		600
+    flat	12:00:00	2
+    cl> asthedit obj* cmds table=table colnames="imagetyp ut itime" verbose+
+    obj001.imh:
+      $I = obj001.imh
+      $D = 22/01/96
+      $T = 20:38:39
+      midut = 9:32:27
+      imagetyp = object
+    obj002.imh:
+      $I = obj002.imh
+      midut = 9:48:05
+      imagetyp = comp
+    ...
+.fi
+
+3.  The following example computes quantities used by some NOAO tasks from
+a minimal ESO/IHAP header.
+
+.nf
+    cl> type eso.dat
+    observat = "eso"
+    ut = sexstr ((@'tm-start'+0.1) / 3600.)
+    utend = sexstr ((@'tm-end'+0.1) / 3600.)
+    epoch = epoch (@'date-obs', ut)
+    st = mst (@'date-obs', ut, obsdb (observat, "longitude"))
+    exptime = (utend>ut)?(utend-ut)*3600.:(utend+24-ut)*3600.
+    ra = sexstr (@'postn-ra' / 15)
+    dec = sexstr (@'postn-dec')
+    airmass = airmass (ra, dec, st, obsdb (observat, "latitude"))
+    imagetyp =  $imagetyp
+    filter = $filter
+    cl> type table.dat
+    object	V
+    as> imhead eso
+        ....
+	DATE-OBS= '12/12/92'            / Date this data created dd/mm/yy
+	TM-START=             84854.    / '23:34:14' measurement start time
+	TM-END  =             84974.    / '23:36:14' measurement end time (U
+	TIME-SID=                 1.    / '00:00:01' sidereal start time
+	POSTN-RA=           354.0709    / '23:36:17' tel. position right-asc
+	POSTN-DE=           6.556945    /'+06:33:25' tel. position declinati
+        ....
+    as> asthedit eso eso.dat table=table.dat col="imagetyp filter" verbose+
+    eso:
+      $I = eso
+      $D = 23/01/96
+      $T = 09:02:55
+      observat = eso
+      ut = 23:34:14
+      utend = 23:36:14
+      epoch = 1992.948616307863
+      st = 0:18:56.76
+      exptime = 120.000000000006
+      ra = 23:36:17
+      dec = 6:33:25
+      airmass = 1.255875187126549
+      imagetyp = object
+      filter = V
+    as> imhead eso
+	...
+	DATE-OBS= '12/12/92'            / Date this data created dd/mm/yy
+	TM-START=             84854.    / '23:34:14' measurement start time
+	TM-END  =             84974.    / '23:36:14' measurement end time (U
+	TIME-SID=                 1.    / '00:00:01' sidereal start time
+	POSTN-RA=           354.0709    / '23:36:17' tel. position right-asc
+	POSTN-DE=           6.556945    /'+06:33:25' tel. position declinati
+	OBSERVAT= 'eso     '
+	UT      = '23:34:14'
+	UTEND   = '23:36:14'
+	EPOCH   =     1992.94861630786
+	ST      = '0:18:56.76'
+	EXPTIME =     120.000000000006
+	RA      = '23:36:17'
+	DEC     = '6:33:25 '
+	AIRMASS =     1.25587518712655
+	IMAGETYP= 'object  '
+	FILTER  = 'V       '
+	...
+.fi
+
+The 0.1 in the UT calculation are to account for round-off.
+Note the use of the conditional expression for the exposure time.
+
+4.  The following example is for a case where there was no telescope
+information but there is date and time information.  This example is
+relevant to data from the Kitt Peak Schmidt telescope circa 1993.
+A table is prepared with the RA, Dec, and Epoch of each observation
+and all other information is derived from the date, ut, and observatory
+database. 
+
+.nf
+    cl> type table.dat
+    12:45:32  +49:34:12   1950
+    13:12:02  -01:12:05   1950
+    cl> type cmds.hast
+    epoch   = epoch (@'date-obs', ut)
+    ra      = ra_precess ($ra, $dec, $epoch, epoch)
+    dec     = dec_precess ($ra, $dec, $epoch, epoch)
+    st      = mst (@'date-obs', ut, obsdb (observat, "longitude"))
+    airmass = eairmass (ra, dec, st, exptime, obsdb (observat, "latitude"))
+    midut   = sexstr (ut + exptime/3600./2.)
+    cl> asthedit *.imh cmds.hast table=table.dat colnames="ra dec epoch" ver+
+    sbs0119.imh:
+      $I = sbs0119.imh
+      $D = 23/01/96
+      $T = 10:38:32
+      epoch = 1987.257752395672
+      ra = 12:47:14.84
+      dec = 49:22:00.39
+      st = 14:53:39.81
+      airmass = 1.154765212092646
+      midut = 9:32:27
+    sbs0120.imh:
+      $I = sbs0120.imh
+      epoch = 1987.257752395672
+      ra = 13:13:56.90
+      dec = -1:23:54.30
+      st = 14:53:39.81
+      airmass = 1.336016291162518
+      midut = 9:32:27
+.fi
+
+Note the use of the table and image header epochs in the precession.
+
+5.  The following example shows the use of the printf function,
+and a null image name, and interactive command input.
+
+.nf
+    cl> asthedit "" ""
+    astcalc> ra = 12:20:30
+    astcalc> dec = 45:00:10
+    astcalc> ep1 = 1950
+    astcalc> ep2 = 2000
+    astcalc> ra1 = ra_precess (ra, dec, ep1, ep2)
+    astcalc> printf ("ra=%h dec=%h\n", ra1, dec_precess (ra, dec, ep1, ep2))
+    ra=12:22:57.4 dec=44:43:32.25
+.fi
+
+.ih
+REVISIONS
+.ls ASTHEDIT V2.11.2
+Y2K update:  The epoch, julday, and mst functions now take either the old
+or new FITS style date strings.  The time argument is optional and if
+it is not specified the time from the date string is used and if neither
+time is present a value of 0h is used.  New internal variables $GMD,
+$GMT, and $GMDT for the current time Greenwich time are defined.
+.le
+.ls ASTHEDIT V2.11
+There are new astronomical functions and input/output functions.
+
+The command syntax may now use "=" as a delimiter as well as the whitespace.
+
+A new parameter "update" allows protecting images and accessing read-only
+images for the purpose of calculating and printing quantities.
+
+The special variable name "$I" has the value of the image name, $D
+the current date, and $T the current time.
+
+The case of no image name creates and deletes a temporary image so the
+task can be used purely as a calculator (but see \fBastcalc\fR).
+.le
+.ih
+SEE ALSO
+astcalc, hedit, hfix, mkheader, setairmass, setjd, asttimes, precess,
+observatory
+.endhelp
diff --git a/noao/astutil/doc/astradius.hlp b/noao/astutil/doc/astradius.hlp
new file mode 100644
index 00000000..f19a47ec
--- /dev/null
+++ b/noao/astutil/doc/astradius.hlp
@@ -0,0 +1,138 @@
+.help astradius Jan96 astutil
+.ih
+NAME
+astradius -- find images within a circle on the sky
+.ih
+USAGE
+astradius images racenter deccenter epcenter radius
+.ih
+PARAMETERS
+.ls images
+List of images for which the radius to a point on the sky is to be
+determined.
+.le
+.ls racenter, deccenter, epcenter
+Right ascension in hours, declination in degrees, and epoch of a position
+on the sky to use as the center of a circle.
+.le
+.ls radius
+Search radius in arc seconds about the center position.
+.le
+.ls keywpars = "" (pset)
+Parameter set defining the image header keywords.  This task requires
+keywords for the right ascension, declination, and epoch.  If
+there is no epoch in the image header keywords for the date of observation
+and the universal time are used for the epoch.  The default parameter
+set (specified by the empty string) is \fBkeywpars\fR.
+.le
+.ls commands = "astutil$astradius.dat"
+Command file used to compute the distance from the coordinate center
+and print a result if the distance is less than the specified radius.
+The command file uses the syntax described for \fBastcalc\fR.
+Users may copy and modify this file if desired.
+.le
+.ih
+DESCRIPTION
+\fBAstradius\fR computes the spherical distance from a specified point on
+the sky for each image in a list of images (\fIimages\fR).  The point on
+the sky is specified by the parameters \fIracenter\fR, \fIdeccenter\fR, and
+\fIepcenter\fR which give a right ascension in hours, a declination in
+degrees, and an epoch.  Each image is required to have keywords for the
+right ascension (hours), declination (degrees), and epoch.  However, if no
+epoch is defined in the image header then an epoch is computed from the
+observation date and universal time.  The spherical distance is compared to
+a specified radius (\fIradius\fR) in arc seconds.  If the distance is less
+than the radius the image name and title are printed.
+
+The image header keywords giving the observation coordinates are defined
+by the parameter set selected with the \fIkeywpars\fR parameter.
+If no value is given then the parameters from the \fBkeywpars\fR
+parameter set task are used.  The keywords required are those
+select by the \fIkeywpars.ra\fR, \fIkeywpars.dec\fR, and
+\fIkeywpars.epoch\fR.  If the epoch is absent or zero then the
+keywords selected by \fIkeywpars.date_obs\fR and \fIkeywpars.ut\fR
+are used to compute an epoch.
+
+\fBAstradius\fR is a simple script which calls \fBastcalc\fR.  The
+command file is specified by the parameter \fIcommands\fR.  The
+default file precesses the observation coordinates to the epoch
+of the search center coordinates and then computes the spherical
+distance between the search center and the observation.  Finally
+it tests the distance against the specified radius and prints
+the image name and title if the observation is within the radius.
+Users may copy the default command file and modify it.  The
+command syntax is described in the help for \fBastcalc\fR.
+.ih
+EXAMPLES
+1.  Page the script task and the command file.
+
+.nf
+    cl> page astutil$astradius.cl,astutil$astradius.dat
+    # ASTRADIUS -- Find images within a radius.
+
+    procedure astradius (images, racenter, deccenter, epcenter, radius)
+
+    string  images = ""             {prompt="List of images"}
+    string  racenter = ""           {prompt="RA center (hours)"}
+    string  deccenter = ""          {prompt="DEC center (degrees)"}
+    real    epcenter = 2000.        {prompt="Epoch of center"}
+    real    radius = 60.            {prompt="Radius in arc seconds"}
+    pset    keywpars = ""           {prompt="Keywords for RA, DEC, EPOCH\n"}
+
+    file    commands = "astutil$astradius.dat"      {prompt="ASTCALC file"}
+
+    begin
+	    astcalc (commands=commands, images=images, table="", verbose=no)
+    end
+
+     Print images which are within a given radius in the sky.
+
+    # Get parameters.
+    racenter = clget ("astradius.racenter")
+    deccenter = clget ("astradius.deccenter")
+    epcenter = clget ("astradius.epcenter")
+    radius = clget ("astradius.radius")
+    ra = imget(clget("keywpars.ra"))
+    dec = imget(clget("keywpars.dec"))
+
+    epoch = imget(clget("keywpars.epoch"))
+    if (str(epoch) == "" || real(epoch) == 0.)
+	date = imget(clget("keywpars.date_obs"))
+	ut = imget(clget("keywpars.ut"))
+	epoch = epoch (date, ut)
+    endif
+
+    # Precess image coordinates to center epoch and compute separation.
+    radec = precess (ra, dec, epoch, epcenter)
+    ra1 = ra_precess (ra, dec, epoch, epcenter)
+    dec1 = dec_precess (ra, dec, epoch, epcenter)
+    sep = arcsep (racenter, deccenter, ra1, dec1)
+
+    # Print result if within radius.
+    if (sep < real (radius))
+	printf ("%-15s %s\n", $I, imget ("title"))
+    endif
+.fi
+
+2. Find images within an arc minute of a particular position.
+
+.nf
+cl> astradius
+List of images: *.imh
+RA center (hours): 13:31
+DEC center (degrees): 47:00
+Epoch of center (2000.):
+Radius in arc seconds (60.):
+obj0020.imh         m51 B 600s
+obj0021.imh         m51 V 600s
+obj0022.imh         m51 R 600s
+.fi
+.ih
+REVISIONS
+.ls ASTRADIUS V2.11
+This task is new in this release.
+.le
+.ih
+SEE ALSO
+astcalc, hselect
+.endhelp
diff --git a/noao/astutil/doc/asttimes.hlp b/noao/astutil/doc/asttimes.hlp
new file mode 100644
index 00000000..eefc17df
--- /dev/null
+++ b/noao/astutil/doc/asttimes.hlp
@@ -0,0 +1,128 @@
+.help asttimes May93 astutil
+.ih
+NAME
+asttimes -- Compute UT, Julian day, epoch, and sidereal time
+.ih
+USAGE
+asttimes
+.ih
+PARAMETERS
+.ls files = ""
+List of files containing local dates and times for which the astronomical
+dates and times are desired.  If no input files are specified then task
+parameters are used.
+.le
+.ls header = yes
+Print header and observatory information to output?
+.le
+.ls observatory = ")_.observatory"
+Observatory for  which times are to be computed.  The default is a
+redirection to look in the parameters for the parent package for a value.
+The final value of this parameter may be one of the
+observatories in the observatory database, "observatory" to select the
+observatory defined by the environment variable "observatory" or the
+parameter \fBobservatory.observatory\fR, or "obspars" to select the
+current parameters set in the \fBobservatory\fR task.  See help for
+\fBobservatory\fR for additional information.
+.le
+.ls year, month, day, time
+If no input files are specified then the date and time for which the
+astronomical date and time is computed are given by these parameters.
+If the year is less than 100 then the century is assumed to be 1900.
+The month is specified as an integer between 1 and 12, and the local
+time for the specified time zone is in hours (sexagesimal format is
+acceptable).
+.le
+.ls ut, epoch, jd, lmst
+If no input files are specified then the universal time, J2000 Julian epoch,
+Julian day, and local mean sidereal time (at the specified longitude)
+are recorded in these parameters for possible reference as CL
+variables.  This is in addition to the usual printed output.
+.le
+.ih
+DESCRIPTION
+The astronomical quantities of universal time, J2000 Julian epoch, Julian day,
+and local mean sidereal time at the specified observatory are computed and
+printed for the given dates and times.  To compute parameters for a
+location not specified in the observatory database use the observatory name
+"obspars" which will use the values defined by the parameters
+\fIobservatory.longitude\fR and \fIobservatory.timezone\fR.  The input
+dates and times may be taken from files containing the year, month (as an
+integer between 1 and 12), day, and local time (sexagesimal notation is
+acceptable) in the specified time zone.  If no files are specified then task
+parameters are used.  The output consists of a printed table with optional
+header and the input data and derived astronomical data.  In addition, if
+the input date and time is from the task parameters then the astronomical
+times are recorded in the user's parameter file (provided the task is not
+run as a background job).  These parameters may then be used as CL
+parameters.
+.ih
+EXAMPLES
+1. For use directly without data files set the date and time using
+the parameter editor, with explicit assignments, or on the command line:
+
+.nf
+    cl> asttimes year=1987 month=10 day=28 time=15:30 obs=kpno
+    # ASTTIMES: Observatory parameters for Kitt Peak National Observatory
+    #       timezone = 7
+    #       longitude = 111:36.0
+    ##YR MON   DAY          ZT         UT      EPOCH           JD       LMST
+    1987  10 28 WED 15:30:00.0 22:30:00.0 1987.82324 2447097.4375 17:30:31.8
+    cl> =asttimes.lmst
+    17.508823973881
+.fi
+
+2. To make a table using a CL loop:
+
+.nf
+    cl> asttimes.observatory="kpno"
+    cl> asttimes.year=1987
+    cl> asttimes.month=10
+    cl> asttimes.time=0
+    cl> for (i=10; i<16; i+=1) {
+    >>> asttimes (day=i, header=no)
+    >>> }
+    1987  10 10 SAT  0:00:00.0  7:00:00.0 1987.77219 2447078.7917  0:47:01.0
+    1987  10 11 SUN  0:00:00.0  7:00:00.0 1987.77493 2447079.7917  0:50:57.5
+    1987  10 12 MON  0:00:00.0  7:00:00.0 1987.77766 2447080.7917  0:54:54.1
+    1987  10 13 TUE  0:00:00.0  7:00:00.0 1987.78040 2447081.7917  0:58:50.7
+    1987  10 14 WED  0:00:00.0  7:00:00.0 1987.78314 2447082.7917  1:02:47.2
+    1987  10 15 THU  0:00:00.0  7:00:00.0 1987.78588 2447083.7917  1:06:43.8
+.fi
+
+In practice the output would be directed to a file:
+
+    >>> asttimes (day=i, header=no, >>"table")
+
+3. To use an input file:
+
+.nf
+    cl> asttimes f=dates > table
+    cl> type table
+    # ASTTIMES: Observatory parameters for Kitt Peak National Observatory
+    #       timezone = 7
+    #       longitude = 111:36.0
+    ##YR MON   DAY          ZT         UT      EPOCH           JD       LMST
+    1987  10 28 WED 22:00:00.0  5:00:00.0 1987.82398 2447097.7083  0:01:35.8
+    1987  10 28 WED 23:00:00.0  6:00:00.0 1987.82409 2447097.7500  1:01:45.7
+    1987  10 29 THU  0:00:00.0  7:00:00.0 1987.82421 2447097.7917  2:01:55.5
+    1987  10 29 THU  1:00:00.0  8:00:00.0 1987.82432 2447097.8333  3:02:05.4
+.fi
+.ih
+REVISIONS
+.ls ASTTIMES V2.10.3
+The epoch was changed from day of the year divided by 365.25 to the
+precise J2000 Julian epoch definition.  In addition to changing
+the output value this fixes incorrect values JD and LMST around the
+new year.
+
+The times are now always printed in the proper 24 hour interval instead
+of using negative or values greater than 24 to indicate the day difference
+with Greenwich.
+
+The header parameter now suppress printing the observatory information.
+.le
+.ih
+SEE ALSO
+observatory
+.endhelp
diff --git a/noao/astutil/doc/ccdtime.hlp b/noao/astutil/doc/ccdtime.hlp
new file mode 100644
index 00000000..318c89b6
--- /dev/null
+++ b/noao/astutil/doc/ccdtime.hlp
@@ -0,0 +1,354 @@
+.help ccdtime Aug98 noao.astutil
+.ih
+NAME
+ccdtime -- compute time, magnitude, and signal-to-noise for CCDs
+.ih
+USAGE
+ccdtime
+.ih
+PARAMETERS
+.ls time = INDEF
+Time in seconds for output of magnitude at the specified signal-to-noise and
+signal-to-noise at the specified magnitude.  This time applies to all
+filters.  If specified as INDEF then no output at fixed exposure time will
+be produced.  If the value is not greater than zero or less than 100000
+an error is reported.
+.le
+.ls magnitude = 20.
+Magnitude for output of time at the specified signal-to-noise and
+signal-to-noise at the specified time.  This magnitude applies to all
+filters. If specified as INDEF then no output at fixed magnitude will
+be produced.  If the absolute value of the magnitude is greater than 40
+an error will be reported.
+.le
+.ls snr = 20.
+Signal-to-noise ratio for output of time at the specified magnitude and
+magnitude at the specified time.  This signal-to-noise ratio applies to all
+filters. If specified as INDEF then no output at fixed signal-to-noise
+ratio will be produced.  If the value is not greater than zero or less than
+100000 an error is reported.
+.le
+
+.ls database = "ccdtime$kpno.dat"
+Database file for telescope, filter, and detector information.  The format
+of this file is described elsewhere.  This file is typically a standard
+file from the logical directory "ccdtime$" or a personal copy in a
+user's directory.
+.le
+.ls telescope = "?"
+Telescope entry from the database.  If "?" a list of telescopes in the
+database is produced.  The name must match the entry name in the database
+but ignoring case.  If the same telescope has multiple focal ratios then
+there must be multiple entries in the database.
+.le
+.ls detector = ""
+Detector entry from the database.  If "?" a list of detectors in the
+database is produced.  The name must match the entry name in the database
+but ignoring case.
+.le
+.ls sum = 1
+CCD on-chip summing or binning factor.
+.le
+.ls seeing = 1.5
+Expected seeing (FWHM) in arc seconds.  The number of pixels used for computing
+the total star counts and the signal-to-noise is given by 1.4 times the square
+of the seeing converted to pixels and rounded up.
+.le
+.ls airmass = 1.2
+Airmass for observation.
+.le
+.ls phase = 0.
+Moon phase in days (0-28) for the estimation of sky brightness.  A
+phase of zero is new moon or dark sky conditions and a phase of 14
+is full moon.
+.le
+
+.ls f1 = "U", f2 = "B", f3 = "V", f4 = "R", f5 = "I"
+Filters for which to compute the CCD information.  If given as "?"
+a list of filters in the database is produced.  If the name (ignoring
+case) is not found then it is ignored.  A null name, that is "",
+is used to eliminate listing of a filter.  There may be many filters
+in the database but the task is currently limited to displaying no
+more than five.
+.le
+.ih
+DESCRIPTION
+A telescope, CCD detector, and list of filters is selected from a database
+to define the expected photon/electron count rates.  These rates along with
+a specified seeing and airmass are used to estimate the signal-to-noise
+ratio (SNR) for a stellar observation in each filter.  The output provides
+three results per filter; the exposure time to achieve a desired SNR for a
+given magnitude, the magnitude to achieve a desired SNR in a given time, and
+the SNR at a specified magnitude and exposure time.  With each of these,
+the number of star photons (or CCD electrons) in an area 1.4 times the
+square of the seeing, the number of sky photons per pixel, and the RMS noise
+contributions from photon noise in the star, the sky, and the detector
+noise from dark current and read out noise are given.  Note that least two
+of the time, magnitude, and signal-to-noise ratio must be specified but if
+one is INDEF then output with that quantity fixed will be skipped or, in
+other words, only the output where the quantity is computed is produced.
+
+The calibration information needed to define the count rates are
+taken from a database file.  This file may be standard ones given in
+the logical directory "ccdtime$" or the user may create their own.
+The database contains entries organized by telescope name (which may
+include a focal ratio if there are multiple ones), detector name,
+and filter name.  One of the standard files may be used as a template.
+
+The file is actually in free format with whitespace and comments ignored.
+However, following the template formatting makes it easy to see the logical
+structure.  All lines, except the "end" line which separates the different
+categories of entries, consist of a keyword an equal sign, and a value
+separated by whitespace.  An entry begins with one of the keywords
+"telescope", "detector", or "filter" and ends with the beginning of
+a new entry or the "end" separator.
+
+A keyword is one of the words shown in the example below.  These keywords
+can also be indexed by the name of a telescope, filter, and/or detector
+entry.  This allows having different transmissions in different filters
+due to correctors, different scales for different detectors which may
+have fore-optics, etc.
+
+Specifically a keyword in the telescope section may have arguments
+from the filter or detector entries, a keyword in the filter section may
+have arguments from the telescope and detector entries, and a keyword
+in the detector section may have arguments from the telescope and filter
+entries.  The formats are keyword, keyword(arg), and keyword(arg,arg).
+The arg fields must match an entry name exactly (without the quotes)
+and there can be no whitespace between the keyword and (, between (
+and the argument, between the arguments and the comma, and between the
+last argument and the closing ).  The software will first look for
+keywords with both arguments in either order, then for keywords with
+one argument, and then for keywords with no arguments.
+
+Below is an example of each type of entry:
+
+.nf
+    telescope = "0.9m"
+	    aperture = 0.91
+	    scale = 30.2
+	    transmission = 1.0
+	    transmission(U) = 0.8
+	    transmission(U,T1KA) = 0.7
+
+    filter = "U"
+	    mag = 20
+	    star = 18.0
+	    extinction = 0.2
+	    sky0 = 22.0
+	    sky1 = -0.2666
+	    sky2 = -.00760
+
+    detector = "T1KA"
+	    rdnoise = 3.5
+	    dark = 0.001
+	    pixsize = 24
+	    U = 0.36
+	    B = 0.61
+	    V = 0.71
+	    R = 0.78
+	    I = 0.60
+.fi
+
+In the example, a transmission of 0.7 will be used if the filter is U
+and the detector is T1KA, a value of 0.8 if the filter is U and the
+detector is not T1KA, and a value of 1 for all other cases.
+
+The telescope entry contains the aperture diameter in meters, the
+scale in arcsec/mm, and a transmission factor.  The transmission factor is
+mostly a fudge factor but may be useful if a telescope has various
+configurations with additional mirrors and optics.
+
+The filter entry contains a fiducial magnitude and the total photon count
+rate for a star of that magnitude.  The units are photons per second
+per square meter of aperture.  An effective extinction in magnitudes/airmass is
+given here.  The sky is defined by a quadratic
+function of lunar phase in days:
+
+.nf
+	if (phase < 14)
+	    sky = sky0 + sky1 * phase + sky2 * phase**2
+	else
+	    sky = sky0 + sky1 * (14 - phase) + sky2 * (14 - phase)**2
+.fi
+
+One may set the higher order terms to zero if the moon contribution
+is to be ignored.  The units are magnitudes per square arc second.
+
+The detector entry contains the read out noise in electrons, the
+dark current rate in electrons per second, the pixel size in
+microns, and the detective quantum efficiency (DQE); the fraction of
+detected photons converted to electrons.  Note that the actual
+values used are the DQE times the rates given by the filter entries.
+Thus, one may set the DQE values to 1 and adjust the filter values
+or set the star count rates to 1 in the filter and set the actual
+count rates in the DQE values.
+
+The computed quantities are formally given as follows.  The
+star count rates for the specified telescope/detector/filter are:
+
+.nf
+	r(star) = star * aperture**2 * transmission *
+	    10**(0.4*(1-airmass)*extinction) * dqe
+.fi
+
+where the "star", "aperture", "transmission", "extinction", are those
+in the database and the "dqe" is the appropriate filter value.  The sky
+rate per pixel is:
+
+.nf
+	r(sky) = r(star) * 10 ** (0.4 * (mag - sky)) * pixel**2
+	pixel = pixsize * scale * sum
+.fi
+
+where mag is the fiducial magnitude, sky is the value computed using
+the quadratic formula for the specified moon phase and the database
+coefficients, the "pixel" size is computed using the CCD pixel size and
+the telescope scale from the database, and sum is the
+specified CCD binning factor.
+
+The number of pixels per star is computed from the seeing as:
+
+.nf
+	npix = 1.4 * (seeing / pixel) ** 2
+.fi
+
+where the number is rounded up to the next integer and a minimum of 9
+pixels is enforced.  This number is a compromise between a large aperture
+for high SNR stars and a smaller aperture for fainter stars.
+
+The number of star photons/electrons per star of magnitude m,
+the number of sky photons per pixel, and the number of dark current
+electrons, all in exposure time t, are given by:
+
+.nf
+	nstar = r(star) * 10 ** (0.4 * (mag - m)) * t
+	nsky = r(sky) * t
+	ndark = dark * t
+.fi
+
+where dark is taken from the detector database entry.
+
+Finally the noise contributions, total noise, and signal-to-noise are
+given by:
+
+.nf
+	noise star = nstar ** 1/2
+	noise sky = (npix * nsky) ** 1/2
+	noise ccd = (npix * (ndark + rdnoise**2)) ** 1/2
+	noise total = (nstar+npix*(nsky+ndark+rdnoise**2)) ** 1/2
+	SNR = nstar / noise total
+.fi
+.ih
+EXAMPLES
+1.  To get a list of the telescopes, filters, and detectors in a database:
+
+.nf
+    cl> ccdtime telescope=? detector=? f1=?
+    Entries for telescope in database ccdtime$kpno.dat:
+	    0.9m
+	    ...
+	    4m
+    Entries for detector in database ccdtime$kpno.dat:
+	    T1KA
+	    T2KA
+	    T2KB
+	    TI2
+	    TI3
+	    T5HA
+	    S2KA
+    Entries for filter in database ccdtime$kpno.dat:
+	    U
+	    B
+	    V
+	    R
+	    I
+.fi
+
+2.  The following is for the default magnitude and SNR and with
+a 1 second exposure time specified.  The output has some
+whitespace removed to fit on this page.
+
+.nf
+    cl> ccdtime time=1
+    Telescope: 0.9m
+    Detector: t1ka
+    Database: ccdtime$kpno.dat Telescope: 0.9m    Detector: t1ka
+      Sum: 1 Arcsec/pixel: 0.72  Pixels/star: 6.0
+      Seeing: 1.5  Airmass: 1.20  Phase: 0.0
+
+
+     Filter  Time   Mag   SNR   Star Sky/pix Noise contributions
+					      Star    Sky    CCD
+
+	  U  70.2  20.0  10.0  196.6    8.8  14.02   8.90  10.53
+	  B  13.0  20.0  10.0  208.8   13.0  14.45  10.82  10.51
+	  V  13.2  20.0  10.0  250.7   29.8  15.83  16.37  10.51
+	  R  17.3  20.0  10.0  365.8   95.9  19.13  29.38  10.51
+	  I 126.4  20.0  10.0 1259.2 1609.8  35.49 120.37  10.55
+
+	  U   1.0  15.6  10.0  166.6    0.1  12.91   1.06  10.50
+	  B   1.0  17.4  10.0  170.0    1.0  13.04   3.00  10.50
+	  V   1.0  17.6  10.0  174.6    2.3  13.21   4.50  10.50
+	  R   1.0  17.6  10.0  186.0    5.5  13.64   7.06  10.50
+	  I   1.0  16.7  10.0  207.9   12.7  14.42  10.71  10.50
+
+	  U   1.0  20.0   0.3    2.8    0.1   1.67   1.06  10.50
+	  B   1.0  20.0   1.4   16.0    1.0   4.00   3.00  10.50
+	  V   1.0  20.0   1.6   19.0    2.3   4.36   4.50  10.50
+	  R   1.0  20.0   1.6   21.1    5.5   4.59   7.06  10.50
+	  I   1.0  20.0   0.7   10.0   12.7   3.16  10.71  10.50
+
+.fi
+
+Note that the default of 1 second in the last section
+gives the count rates per second for star and sky.
+
+3.  Sometimes one may want to vary one parameter easily on the command
+line or query.  This can be done by changing the parameter to query
+mode.  In the following example we want to change the magnitude.
+
+.nf
+    cl> ccdtime.magnitude.p_mode=query
+    cl> ccdtime.telescope="0.9m"
+    cl> ccdtime.detector="t1ka"
+    cl> ccdtime.f1=""; ccdtime.f5=""
+    cl> ccdtime
+    Magnitude (20.):
+    Database: ccdtime$kpno.dat   Telescope: 0.9m     Detector: t1ka
+      Sum: 1 Arcsec/pixel: 0.72  Pixels/star: 6.0
+      Seeing: 1.5  Airmass: 1.20  Phase: 0.0
+
+     Filter  Time   Mag   SNR  Star Sky/pix  Noise contributions
+					       Star   Sky    CCD
+
+	  B  13.0  20.0  10.0 208.8    13.0  14.45  10.82  10.51
+	  V  13.2  20.0  10.0 250.7    29.8  15.83  16.37  10.51
+	  R  17.3  20.0  10.0 365.8    95.9  19.13  29.38  10.51
+
+    cl> ccdtime 21
+    ...
+    cl> ccdtime 22
+    ...
+.fi
+.ih
+REVISIONS
+.ls CCDTIME V2.11.4
+A error will be reported if the requested time or SNR is not greater
+than zero and less than 100000., or if the absolute value
+of the magnitude is greater than 40.
+.le
+.ls CCDTIME V2.11.2
+The incorrect usage of a 1 mag/airmass extinction was fixed by adding an
+expected "extinction" entry in the filter entries.  Note that old files
+will still give the same result by using an extinction of 1 if the keyword
+is not found.
+
+The database keywords can not be indexed by telescope, filter, and/or
+detector.
+
+The number of pixels per aperture now has a minimum of 9 pixels.
+.le
+.ih
+SEE ALSO
+.endhelp
diff --git a/noao/astutil/doc/galactic.hlp b/noao/astutil/doc/galactic.hlp
new file mode 100644
index 00000000..f6beedc9
--- /dev/null
+++ b/noao/astutil/doc/galactic.hlp
@@ -0,0 +1,68 @@
+.help galactic Oct87 noao.astutil
+.ih
+NAME
+galactic -- convert between equatorial and galactic coordinates
+.ih
+USAGE
+galactic files
+.ih
+PARAMETERS
+.ls files
+The name of a file (or a file list or template) containing the coordinates
+to be converted.
+.le
+.ls in_coords = "equatorial"
+Type of input coordinates.  May be either "equatorial" (RA and DEC) or
+"galactic" (l and b).
+.le
+.ls print_coords = yes
+If \fBprint_coords\fR = yes, the RA, DEC and epoch (as well as lII and bII) 
+will be listed on the output file.
+.le
+.ih
+DESCRIPTION
+Program \fBgalactic\fR is used to convert between equatorial and
+galactic coordinates.  It converts in either direction based on the
+specified input coordinates.  Coordinates are read from the input file
+as RA and DEC or galactic longitude and latitude pairs, one pair per
+input line.  Each coordinate pair may optionally be followed by the
+epoch of the equatorial coordinates, in which case the coordinates are
+precessed to 1950.0 (the epoch of definition for the galactic center)
+before conversion for equatorial to galactic or to the specified epoch
+for galactic to equatorial.  Coordinates may be entered in either
+decimal or sexagesimal notation.
+.ih
+EXAMPLES
+1. Convert the given RA and DEC coordinates to galactic coordinates.  When
+the epoch is specified as other than 1950.0, precess before converting.
+The lines input by the user are marked:
+
+.nf
+	cl> galactic STDIN 	         		[input]
+	12:30:10.12 10:18:27.5 1930.			[input]
+  	12:30:10.12   10:18:27.5  1930.00     288.4695   72.2884
+	12:30 10:18					[input]
+  	12:30:00.00   10:18:00.0  1950.00     287.4598   72.3202
+	12.5  10:18                                     [input]
+  	12:30:00.00   10:18:00.0  1950.00     287.4598   72.3202
+	(eof=<ctrl/z>)					[input]
+.fi
+
+2. The following is equivalent, except that coordinate input is taken from
+the file "coords", rather than from the terminal:
+
+.nf
+	cl> galactic coords 				[input]
+  	12:30:10.12   10:18:27.5  1930.00     288.4695   72.2884
+  	12:30:00.00   10:18:00.0  1950.00     287.4598   72.3202
+  	12:30:00.00   10:18:00.0  1950.00     287.4598   72.3202
+.fi
+
+3. If image headers contain the coordinates, in this case RA, DEC, and EPOCH,
+then one can get the galactic coordinates for the image by:
+
+	cl> hselect *.imh ra,dec,epoch yes | galactic STDIN
+
+(Consult the help for the task \fBhselect\fR for information about selecting
+fields from image headers.)
+.endhelp
diff --git a/noao/astutil/doc/gratings.hlp b/noao/astutil/doc/gratings.hlp
new file mode 100644
index 00000000..593ddd62
--- /dev/null
+++ b/noao/astutil/doc/gratings.hlp
@@ -0,0 +1,252 @@
+.help gratings Mar91 noao.astutil
+.ih
+NAME
+gratings -- Compute and print grating parameters
+.ih
+USAGE
+gratings
+.ih
+PARAMETERS
+.ls echelle = no
+Is the grating an echelle grating?  This selects whether the angle of
+incidence is greater or less than blaze angle when the angle of incidence
+or blaze angle are not specified.  For an echelle the angle of incidence
+is generally greater than the blaze angle.
+.le
+.ls f = 590.
+Focal length in millimeters.  Technically it is defined by the equation x =
+f * tan (theta) where x is distance from the optical axis on the detector
+and theta is the diffraction angle; i.e. it converts angular measures to
+millimeters on the detector.  If the focal length is specified as INDEF it
+is computed from the dispersion, which is required in this case, and the
+other parameters.
+.le
+.ls gmm = 226.
+Grating grooves per millimeter.  If specified as INDEF it is computed
+from the order, which is required in this case, and the other parameters.
+.le
+.ls blaze = 4.5
+Blaze angle in degrees.  It is always specified or printed as a positive
+angle relative to the grating normal.  If specified as INDEF it is
+computed from the other parameters.
+.le
+.ls theta = -10.5
+Angle of incidence in degrees.  The angle of incidence must be in the plane
+perpendicular to face of the grating.  The angle of incidence may be
+specified relative to the grating normal or the blaze angle though it is
+always printed relative to the grating normal.  To specify it relative to
+the blaze angle add 360 degrees; for example to have an angle of 15 degrees
+less than the blaze angle specify 360 - 15 = 345.  If the angle of
+incidence is specified as INDEF it is computed from the other parameters.
+.le
+.ls order = 1
+Order for which the wavelength and dispersion are specified.  If specified
+as INDEF it will be computed from the grooves per mm, which is required in
+this case, and the other parameters.
+.le
+.ls wavelength = INDEF
+Blaze wavelength in Angstroms.  If specified as INDEF it will be computed
+from the other parameters.
+.le
+.ls dispersion = INDEF
+Blaze dispersion in Angstroms per millimeter.  If specified as INDEF it
+will be computed from the focal length, which is required in this case,
+and the other parameters.
+.le
+.ih
+DESCRIPTION
+This task computes the grating parameters specified as INDEF from the other
+grating parameters and prints the final set of self-consistent parameters.
+The parameters are the focal length to the detector, the grooves per
+millimeter of the grating, the blaze angle of the grating, the angle of
+incidence of the incoming light to the grating (which is required to be in
+the plane perpendicular to the face of the grating), the diffraction order,
+and the blaze wavelength and dispersion at the blaze wavelength on the
+detector for that order.  There must be five of these parameters specified
+to compute the remaining two with the exceptions that the combinations
+of the grooves per millimeter and the order or the focal length and
+dispersion must not be simultaneously unspecified.  There are two cases in
+which the computation will not succeed, if not enough parameters are
+specified or when the combination of parameters is not possible.  In these
+cases a warning is printed and the input parameters, including INDEF
+values, are printed.
+
+If more than the minimum number of parameters are specified then some of
+the specified parameters will be adjusted to give a self-consistent set.
+In particular, if all parameters are specified the input wavelength and
+dispersion are ignored and new values are calculated.  If only one
+parameter is not specified then the dispersion is adjusted if it is not the
+dispersion the wavelength is adjusted if it is the dispersion.
+
+When the order is not specified, the nearest integer order is computed from
+the other non-integer parameters and then the wavelength and dispersion are
+recomputed based on the integer order.
+
+The basic grating equation used is
+
+.nf
+(1)	m * lambda = (sin(theta) + sin(beta)) / g
+.fi
+
+where m is the order, lambda the wavelength, g the grooves per wavelength unit,
+theta the angle of incidence to the grating normal, and beta the angle of
+diffraction to the normal.  The diffraction angle relative to that
+of the blaze maximum, psi, is given by
+
+.nf
+(2)	beta = psi + 2 * blaze - theta
+.fi
+
+where blaze is the blaze angle.  The diffraction angle psi is related to
+position on the detector, again measured from the blaze peak, by
+
+.nf
+(3)	x = f * tan(psi)
+.fi
+
+where f is the effective focal length (as defined by this equation).
+At the blaze maximum psi = x = 0 and the wavelength and dispersion
+per millimeter on the detector are given by (1) and the derivative of (1)
+with respect to x:
+
+.nf
+(4)	wavelength = 1E7*(sin(theta)+sin(2*blaze-theta))/(gmm*order)
+(5)	dispersion = 1E7*cos(2*blaze-theta)/(gmm*order*f)
+.fi
+
+where the variable names are the same as the program parameters and
+the factor of 1E7 is the conversion between millimeters and Angstroms.
+
+Equations (4) and (5) are the ones solved by this task.  There are a some
+interesting points to note about the angle of incidence.  There are two
+solutions of these equations one with the angle of incidence less than the
+blaze angle and one greater than the blaze angle.  For an echelle the angle
+of incidence is generally set greater than the blaze angle to avoid light
+lost by reflections back along the angle of incidence.  The \fIechelle\fR
+parameter is used to determine which side of the blaze angle the angle of
+incidence will be computed in the cases in which it is not specified;
+greater than the blaze angle when yes and less than the blaze angle when
+no.
+
+In spectrographs it is often the case that the angle between the
+incoming beam and center of the diffracted beam, delta, is fixed where
+
+.nf
+(6)	delta = 2 * |theta - blaze|
+.fi
+
+This fixes the angle between the blaze angle and the angle of incidence
+needed to center the blaze function on the detector.  If one wants to
+solve (4) and (5) for the blaze angle with this difference fixed the
+angle of incidence may be specified relative to the blaze angle by
+adding 360 degrees to the difference.  An example best describes this.
+The Kitt Peak 4m Echelle Spectrograph has a 12 degree angle
+between the incoming beam to the echelle grating and the beam to the
+crossdisperser.  Then |theta - blaze| = 6 degrees.  For an echelle the
+angle of incidence is greater than the blaze angle (relative to the
+grating normal) so if we set the angle of incidence to 6 + 360
+and the blaze angle to INDEF the resulting computation will
+determine blaze and theta with a fixed 6 degree angle.
+.ih
+EXAMPLES
+1.  The default values are for a grating of 226 grooves per millimeter
+in a 590 mm focal length camera.  For a blaze angle of 4.5 degrees
+and an angle of incidence of -10.5 degrees (the angle is on the
+other side of the grating normal relative to the blaze angle) the
+first order wavelength and dispersion at the blaze peak is:
+
+.nf
+	cl> gratings
+	Grating parameters:
+	  Focal length = 590. mm
+	  Grating = 226. grooves/mm
+	  Blaze angle = 4.5 degrees
+	  Incidence angle = -10.5 degrees
+	  Order = 1
+	  Blaze wavelength = 6706.696 Angstroms
+	  Blaze dispersion = 70.69458 Angstroms/mm
+.fi
+
+2.  To find nearest order and the dispersion for a wavelength of 3400
+Angstroms:
+
+.nf
+	cl> gratings order=INDEF wave=3400
+	Grating parameters:
+	  Focal length = 590. mm
+	  Grating = 226. grooves/mm
+	  Blaze angle = 4.5 degrees
+	  Incidence angle = -10.5 degrees
+	  Order = 2
+	  Blaze wavelength = 3353.348 Angstroms
+	  Blaze dispersion = 35.34729 Angstroms/mm
+.fi
+
+3.  To find the grating parameters need to center 8000 Angstroms with
+a dispersion of 90 Angstroms per millimeter:
+
+.nf
+	cl> gratings gmm=INDEF blaze=INDEF theta=345 wave=8000 disp=90
+	Grating parameters:
+	  Focal length = 590. mm
+	  Grating = 177.8237 grooves/mm
+	  Blaze angle = 4.223008 degrees
+	  Incidence angle = -10.77702 degrees
+	  Order = 1
+	  Blaze wavelength = 8000. Angstroms
+	  Blaze dispersion = 90. Angstroms/mm
+.fi
+
+4.  What focal length should be used to get a dispersion of 20 Angstroms/mm
+at 6700 Angstroms:
+
+.nf
+	cl> gratings f=INDEF wave=6700 disp=20
+	Grating parameters:
+	  Focal length = 2085.49 mm
+	  Grating = 226. grooves/mm
+	  Blaze angle = 4.5 degrees
+	  Incidence angle = -10.5 degrees
+	  Order = 1
+	  Blaze wavelength = 6706.696 Angstroms
+	  Blaze dispersion = 20. Angstroms/mm
+.fi
+
+5.  What are the first order wavelength parameters for an echelle of
+31.6 grooves per millimeter with a 63 degree blaze, and a 6 degree
+angle of incidence relative to the blaze angle.  Then what are
+the wavelength parameters in 80th order and what order is 6563 in.
+
+.nf
+	cl> gratings gmm=31.6 blaze=63 theta=69
+	Grating parameters:
+	  Focal length = 590. mm
+	  Grating = 31.6 grooves/mm
+	  Blaze angle = 63. degrees
+	  Incidence angle = 69. degrees
+	  Order = 1
+	  Blaze wavelength = 560838.9 Angstroms
+	  Blaze dispersion = 292.1256 Angstroms/mm
+	cl> gratings gmm=31.6 blaze=63 theta=69 order=80
+	Grating parameters:
+	  Focal length = 590. mm
+	  Grating = 31.6 grooves/mm
+	  Blaze angle = 63. degrees
+	  Incidence angle = 69. degrees
+	  Order = 80
+	  Blaze wavelength = 7010.487 Angstroms
+	  Blaze dispersion = 3.651571 Angstroms/mm
+	cl> gratings gmm=31.6 blaze=63 theta=69 order=INDEF wave=6563
+	Grating parameters:
+	  Focal length = 590. mm
+	  Grating = 31.6 grooves/mm
+	  Blaze angle = 63. degrees
+	  Incidence angle = 69. degrees
+	  Order = 85
+	  Blaze wavelength = 6598.105 Angstroms
+	  Blaze dispersion = 3.436772 Angstroms/mm
+.fi
+.ih
+SEE ALSO
+artdata.mkechelle
+.endhelp
diff --git a/noao/astutil/doc/keywpars.hlp b/noao/astutil/doc/keywpars.hlp
new file mode 100644
index 00000000..e5b50d4a
--- /dev/null
+++ b/noao/astutil/doc/keywpars.hlp
@@ -0,0 +1,94 @@
+.help keywpars May93 noao.astutil
+.ih
+NAME
+keywpars -- edit the image header keywords used by the package
+.ih
+USAGE
+keywpars
+.ih
+PARAMETERS
+.ls ra = "RA"
+Right Ascension keyword. (Value in HMS format).
+.le
+.ls dec = "DEC"
+Declination keyword. (Value in HMS format).
+.le
+.ls ut = "UT"
+UT of observation keyword.  This field is the UT start of the observation.
+(Value in HMS Format).
+.le
+.ls utmiddle = "UTMIDDLE"
+UT mid-point of observation keyword.  This field is the UT mid-point of 
+the observation.  (Value in HMS Format).
+.le
+.ls exptime = "EXPTIME"
+Exposure time keyword. (Value in Seconds).
+.le
+.ls epoch = "EPOCH"
+Epoch of coordinates keyword. (Value in Years).
+.le
+.ls date_obs = "DATE-OBS"
+Date of observation keyword.  Format for this field should be
+dd/mm/yy (old FITS format), yyyy-mm-dd (new FITS format), or
+yyyy-mm-ddThh:mm:ss.sss (new FITS format with time).
+.le
+
+.ce 
+OUTPUT KEYWORDS
+.ls hjd = "HJD"
+Heliocentric Julian date keyword. (Value in Days).
+.le
+.ls mjd_obs = "MJD-OBS"
+Modified Julian Data keyword.  The MJD is defined as the Julian date of
+the mid-point of the observation - 2440000.5.  (Value in Days).
+.le
+.ls vobs = "VOBS"
+Observed radial velocity keyword.  (Value in Km/sec).
+.le
+.ls vrel = "VREL"
+Observed radial velocity keyword. (Value in Km/sec).
+.le
+.ls vhelio = "VHELIO"
+Corrected heliocentric radial velocity keyword.  (Value in Km/sec).
+.le
+.ls vlsr = "VLSR"
+Local Standard of Rest velocity keyword.  (Value in Km/sec).
+.le
+.ls vsun = "VSUN"
+Epoch of solar motion.  (Character string with four real valued fields 
+describing the solar velocity (km/sec), the RA of the solar velocity (hours),
+the declination of the solar velocity (degrees), and the epoch of solar
+coordinates (years)).
+.le
+.ih
+DESCRIPTION
+The image header keywords used by the \fIfxcor\fR task can be 
+edited if they differ
+from the NOAO standard keywords.  For example, if the image header keyword
+giving the exposure time for the image is written out as "EXP-TIME" instead
+of the standard "OTIME" at a given site, the keyword accessed for 
+that information
+may be changed based on the value of the \fIexptime\fR parameter.
+
+.ih
+EXAMPLES
+1. List the image header keywords.
+
+.nf
+	as> lpar keywpars
+.fi
+
+2. Edit the image header keywords
+
+.nf
+	as> keywpars
+.fi
+.ih
+REVISIONS
+.ls KEYPARS V2.10.3
+First version.  Currently only used by the \fIRVCORRECT\fR task.
+.le
+.ih
+SEE ALSO
+fxcor, rvcorrect
+.endhelp
diff --git a/noao/astutil/doc/obs.hlp b/noao/astutil/doc/obs.hlp
new file mode 100644
index 00000000..4b5a8a88
--- /dev/null
+++ b/noao/astutil/doc/obs.hlp
@@ -0,0 +1,390 @@
+.help observatory Jan92 noao
+.ih
+NAME
+observatory -- Examine and set observatory parameters
+.ih
+USAGE
+observatory command obsid [images]
+.ih
+PARAMETERS
+.ls command = "list" (set|list|images)
+Command option which is one of "set", "list", or "images".  The set command
+sets the default observatory task parameters for the specified
+observatory.  The list command lists the observatory parameters for the
+specified observatory but does not modify the task parameters.  The images
+command lists the observatory parameters for a list of images.  The list
+and images commands examine and verify the observatory parameters applied
+by other tasks using the observatory database facility.
+.le
+.ls obsid = "?"
+Observatory identification to be set, listed, or used as the default for
+images without the OBSERVAT keyword.  The observatory ID is one of those in
+the database (case ignored), the special string "observatory" to default to
+the environment variable "observatory" or the \fIobservatory.observatory\fR
+parameter, "obspars" to select the parameters in the \fBobservatory\fR
+task, or "?" to list the observatories defined in the database.
+.le
+.ls images
+List of images to be examined with the "images" command.  The images are
+checked for the OBSERVAT keyword to determine the observatory parameters
+to be listed, otherwise the observatory given by \fIobsid\fR is used.
+.le
+.ls verbose = no
+Verbose output?  Because there are a number of different ways in which
+observatory information is determine this option prints detailed
+information on how the observatory database and parameters are
+ultimately selected.
+.le
+
+.ls observatory
+The default observatory used by tasks which use the special
+observatory identification "observatory".  The value is one of the
+observatory names in the observatory database (case ignored)
+or the special value "obspars" to select the parameters defined in this
+task.  There is no default to force users to set it at least once.
+.le
+.ls name
+Observatory name.
+.le
+.ls longitude
+Observatory longitude given in degrees west.
+.le
+.ls latitude
+Observatory latitude in degrees.  Positive latitudes are north and negative
+latitudes are south.
+.le
+.ls altitude
+Observatory altitude in meters above sea level.
+.le
+.ls timezone
+Observatory time zone.  The time zone is the number of hours west of
+Greenwich or the number of hours to be added to local time to obtain
+Greenwich time.
+.le
+.ih
+ENVIRONMENT VARIABLES
+.ls obsdb
+This variable selects the observatory database.  If not defined it defaults
+to noao$lib/obsdb.dat.
+.le
+.ls observatory
+This variable selects the observatory entry whenever a task uses the
+observatory name "observatory".  If not defined the value of the task
+parameter \fIobservatory.observatory\fR is used.
+.le
+.ih
+IMAGE HEADER KEYWORDS
+The observatory identification for images is first sought under the
+image header keyword OBSERVAT.  This always takes precedence over any
+other means of defining the observatory.
+.ih
+DESCRIPTION
+
+OBSERVATORY PARAMETERS IN THE NOAO PACKAGE
+
+Some astronomical data reduction and analysis tasks perform
+computations requiring information about where the data was observed.
+For example a number of \fBnoao\fR tasks make corrections for the
+airmass.  Generally they look for an airmass in the image header and
+if it is not present they attempt to compute it from other image header
+parameters.  The information about time and telescope coordinates
+of the observation are often in the image header but the observatory
+latitude is not.  The task must get this information somehow.
+
+Prior to IRAF V2.10 tasks generally had explicit parameters, such as
+latitude, with default values pointing (using parameter redirection) to
+the parameter of the same name in the \fBobservatory\fR task.  The
+user was required to know the values of the observatory parameters and
+manually change them for data from different observatories.  In V2.10
+an observatory database has been implemented.  Observatory parameters
+are stored in a simple text file and tasks obtain observatory related
+parameters by specifying an observatory identification.
+
+In general the information about the observatory should be directly
+associated with the image data.  Unless stated otherwise in the
+description of a task,  tasks which require observatory information
+will first look for the image header keyword OBSERVAT.  The value of
+this keyword is the observatory identification used to index the
+observatory database.  The task will then look up any observatory
+parameters it needs in the observatory database.  Data from
+observatories that support this keyword will, therefore, always use the
+correct observatory parameters without user intervention.  All
+observatories which export FITS image data are urged to adopt the
+OBSERVAT keyword (a keyword recommended by the FITS standard).
+
+For image data which do not identify the observatory in this way
+and in tasks which do not operate on images (such as astronomical
+calculator tools), the observatory must be specified by the user.
+Most tasks provide an "observatory" parameter which either directly
+selects the observatory or use special values for defining the
+observatory with an environment variable or the parameters
+from the \fBobservatory\fR task.
+
+An observatory is specified by the identification name used in the
+observatory database.  The names in the database may be listed using
+the \fBobservatory\fR task as described below.  If the desired observatory
+is not in the database a user may copy/create their own database and
+select it with the environment variable "obsdb", modify the standard
+database if allowed (any changes to the distributed version should
+be forwarded to iraf$noao.edu), or use the special observatory name
+"obspars".  The last option directly uses the parameters in the
+\fBobservatory\fR task which can be set to any values using the normal
+parameter editing mechanism.
+
+The default value for the observatory parameter in a task is generally
+"observatory".  This special name directs the task to look first
+for the environment variable of the same name and then at the
+\fIobservatory\fR parameter of the \fBobservatory\fR task.  The environment
+variable allows users or sites to set the default observatory in their
+login files and site defaults.  Also it is simple to change the
+default observatory either with a \fBreset\fR command or the
+\fBobservatory\fR command.
+
+The observatory database is selected by the environment variable
+"obsdb".  The default when the variable is not defined is the
+\fBnoao\fR package library database file "noao$lib/obsdb.dat".  The use
+of an environment variable allows users to permanently change the
+default database in the OS environment (when IRAF has access to it such
+as in UNIX systems) or in the startup IRAF environment as set in the
+"login.cl" or "loginuser.cl" files.  One can, of course, change it
+during a session with the set or reset commands.  For sites which want
+to customize the observatory mechanism the environment variables can
+also be set and changed in the files "hlib$zzsetenv.def",
+"noao$lib/zzsetenv.def", and the template login file "hlib$login.cl".
+
+An observatory database file consist of a simple list of keyword=value
+pairs with arbitrary whitespace allowed.  An observatory entry begins
+with the observatory keyword and extends to the next observatory
+keyword or the end of the file.  The observatory identification should
+be the same as the string used in the OBSERVAT image header parameter
+for data from that observatory.  The default file noao$lib/obsdb.dat
+begins as follows:
+
+.nf
+# Observatory Parameters.  Taken from the Almanac.
+#
+# Observatories wishing to be added or make changes in the default
+# distributed database should send information to iraf@noao.edu.
+
+observatory = "kpno"
+	name = "Kitt Peak National Observatory"
+	longitude = 111:36.0
+	latitude = 31:58.8
+	altitude = 2120.
+	timezone = 7
+
+observatory = "ctio"
+	<etc>
+.fi
+
+In summary, access to observatory parameters is now done by referencing
+the image header keyword OBSERVAT and, if not defined, determine the
+observatory name from a task parameter.  The environment variables
+"observatory" and "obsdb" can be set by the user to select alternate
+observatories and observatory database files.  For data without an
+observatory entry the observatory can be set to "obspars" or the user
+may make their own observatory database.
+
+THE OBSERVATORY TASK
+
+The \fBobservatory\fR task serves a number of functions.  It may be used to
+examine the observatory database, verify the observatory parameters which
+will be used by other tasks, particularly those operating on images, set
+the default observatory if not defined by other means, set observatory
+parameters explicitly, especially when there is no observatory database
+entry, and as a parameter set for tasks which explicitly reference
+observatory parameters.  The \fBverbose\fR parameter also provides a
+detailed check of the steps used to determine the observatory database,
+observatory identification, and observatory parameters.
+
+The \fIcommand\fR parameter takes the values "set", "list", or "images".
+The \fIobsid\fR parameter supplies the observatory identification and the
+\fIimages\fR parameter is used to specify a list of images for the "images"
+command.  The parameters are query parameters and so may be either queried
+or simply typed on the command line.
+
+The "set" command prints the observatory parameters for the specified
+observatory and sets many of these in the \fBobservatory\fR task
+parameters.  This command is used to set the default observatory parameters
+for tasks where images are not used, the images do not contain the
+observatory id, or direct references to specific parameters with parameter
+redirection (for example ")observatory.latitude") are used.
+
+The "list" command is similar to the "set" command except the task parameters
+are not modified.  It is used to list observatory parameters.  It is also
+use with the special observatory identifications to list the entries in
+an observatory database and verify the observatory to be used by
+tasks which do not operate on images.  The special value "?" lists
+the entries in the database.  The special value "observatory" lists
+the observatory defined by the "observatory" environment variable or
+that given by the \fIobservatory.observatory\fR parameter.  The special
+value "obspars" simply lists the observatory task parameters.
+
+The "images" command lists the observatory information applicable to
+one or more images.  In particular, the observatory identification is
+first sought in OBSERVAT image header keyword and, if not found, the
+\fIobsid\fR parameter is used.  Often the default observatory is
+"observatory" to follow the same search path used by other tasks.
+
+The \fIverbose\fR parameter prints additional detailed information.  It
+prints the database used and whether it is selected by default
+(noao$lib/obsdb.dat) or by the "obsdb" environment variable.  When the
+observatory is defined as "observatory" it indicates whether the
+observatory is defined by the environment variable "observatory" or by the
+observatory task.  When listing images it prints the OBSERVAT keyword or
+the default observatory assigned.
+
+For observatories not in a database the name, latitude, longitude,
+altitude, and time zone parameters may be set using \fBeparam\fR.
+The observatory id must be set to "obspars" in this case.
+These parameters will then be referenced by other tasks in which
+the observatory is specified as "obspars".  This allows arbitrary
+observatory parameters to be set without creating or modifying
+an observatory database.  However, it is advisable to create a
+local database and also send the observatory information to the
+IRAF group at NOAO for inclusion in the default database.
+.ih
+EXAMPLES
+1.  List the observatory entries in the database:
+
+.nf
+	cl> observatory list ? v+
+	Using default observatory database: noao$lib/obsdb.dat
+
+	default: Kitt Peak National Observatory
+	kpno: Kitt Peak National Observatory
+	ctio: Cerro Tololo Interamerican Observatory
+	eso: European Southern Observatory
+	lick: Lick Observatory
+	mmt: Whipple Observatory
+	cfht: Canada-France-Hawaii Telescope
+	lapalma: Roque de los Mucachos, La Palma
+.fi
+
+2.  Set the observatory parameters for Cerro Tololo:
+
+.nf
+	cl> observatory set ctio
+	Observatory parameters for Cerro Tololo...
+		observatory = ctio
+		timezone = 5
+		altitude = 2215.
+		latitude = -30:09.9
+		longitude = 70:48.9
+	         name = 'Cerro Tololo Interamerican Observatory'
+	cl> lpar observatory
+	      command = "set"		Command (set|list|images)
+	     argument = ctio		Observatory or images
+	 (observatory = "ctio")         Observatory identification
+	        (name = "Cerro Tololo...") Observatory name
+	   (longitude = 70.815)         Observatory longitude (degrees)
+	    (latitude = -30.165)        Observatory latitude (degrees)
+	    (altitude = 2215.)          Observatory altitude (meters)
+	    (timezone = 4)              Observatory time zone
+	     (verbose = no)             Verbose output?
+	        (mode = "q")            
+.fi
+
+3.  Set the observatory parameters to use the environment variable
+"observatory" and verify it.
+
+.nf
+	cl> set observatory=cfht
+	cl> observatory list observatory
+	Observatory parameters for Canada-France-Hawaii Telescope
+		observatory = cfht
+		timezone = 10
+		altitude = 4215
+		latitude = 19:49.6
+		longitude = 155:28.3
+		name = 'Canada-France-Hawaii Telescope'
+.fi
+
+4.  Change the default observatory database and verify verbosely:
+
+.nf
+	cl> set observatory="sco"
+	cl> set obsdb="/local/iraf/obsdb.dat"
+	cl> type obsdb$
+	# Local Observatory Parameters.
+
+	observatory = "sco"
+		name = "Small College Observatory"
+		longitude = 100:20.0
+		latitude = 35:58.8
+		altitude = 212.
+		timezone = 6
+	cl> observ set observatory v+
+	Using database defined by 'obsdb' environment variable:
+		/tmp/test/obsdb.dat
+	Using obs... defined by 'obs...' environment variable: sco
+	Using observatory parameters for database entry: sco
+	Observatory parameters for Small College Observatory
+		observatory = sco
+		timezone = 6
+		altitude = 212.
+		latitude = 35:58.8
+		longitude = 100:20.0
+		name = 'Small College Observatory'
+.fi
+
+5.  List the observatory assigned to some images with a default observatory
+determined either by the "observatory" environment variable or that set
+in the observatory task.
+
+.nf
+	cl> observ images observatory dev$pix,demoobj1
+	Observatory parameters for Small College Observatory
+		observatory = sco
+		timezone = 6
+		altitude = 212.
+		latitude = 35:58.8
+		longitude = 100:20.0
+		name = 'Small College Observatory'
+		Images: dev$pix (default observatory)
+	Observatory parameters for Kitt Peak National Observatory
+		observatory = kpno
+		timezone = 7
+		altitude = 2120.
+		latitude = 31:58.8
+		longitude = 111:36.0
+		name = 'Kitt Peak National Observatory'
+		Images: demoobj1 (OBSERVAT keyword)
+
+.fi
+
+6.  Set explicit observatory parameters:
+
+.nf
+	cl> epar observatory
+	<set observatory parameters>
+	cl> observ list obspars
+	Observatory parameters for North Pole
+		observatory = obspars
+		timezone = 0
+		altitude = 0.
+		latitude = 90.
+		longitude = 0.
+		name = 'North Pole'
+.fi
+
+7.  Use observatory parameters in expressions:
+
+.nf
+	cl> observ set kpno
+	Observatory parameters for Kitt Peak National Observatory
+		observatory = kpno
+		timezone = 7
+		altitude = 2120.
+		latitude = 31:58.8
+		longitude = 111:36.0
+		name = 'Kitt Peak National Observatory'
+	cl> = observ.lat
+	31.98
+	cl> = sin (3.14159/180 * observ.lat)
+	0.52962280742153
+.fi
+.ih
+SEE ALSO
+Tasks in astutil, imred, onedspec, and twodspec.
+.endhelp
diff --git a/noao/astutil/doc/pdm.hlp b/noao/astutil/doc/pdm.hlp
new file mode 100644
index 00000000..51b24d5f
--- /dev/null
+++ b/noao/astutil/doc/pdm.hlp
@@ -0,0 +1,372 @@
+.help pdm May87 noao.astutil
+.ih
+NAME
+pdm -- Find periods in lightcurve data.
+.ih
+USAGE
+pdm infile
+.ih
+PARAMETERS
+
+.ls infiles
+Input file template.  If more than one file matches the template, data
+from all the files will be concatenated to produce the working data set.
+.le
+.ls metafile = "pdmmeta"
+File in which to store metacode when running in batch mode.  All of the
+plots saved will be put here with formfeeds between them.
+.le
+.ls batchfile = "pdmbatch"
+File in which to store information about the period, amplitude, epoch
+and fit function when running in batch mode.
+.le
+.ls device = "stdgraph"
+The output device for interactive graphics.
+.le
+.ls interactive = yes
+Interactive flag.  If set to no, the analysis is done in batch mode with output
+to a file and no interactive graphics.  Metacode will be saved for the data
+plot, the theta plot, and the phase plot.  If set to yes, various types of
+plots can be made on the user's terminal and cursor commands are available.
+.le
+.ls flip = no
+Flag to tell the program to flip the y-axis.  This is useful when the y-scale
+is in magnitudes (decreasing numbers mean increasing brightness).
+.le
+.ls minp = 0
+Minimum period to be searched.  This parameter, if set, tells the program
+the bottom end of the period window to be searched.   If not set, the
+program uses as a value the smallest chronological distance between
+any two adjacent data points.  When the program is run, it writes a value
+into this parameter as stored in the uparm directory.  This means the
+next time the program is run, the default minp will be the value assigned
+or calculated the last time the program was run by this user.  We say the
+program 'remembers' the last value used.
+.le
+.ls maxp = 0
+Maximum period to be searched.  This parameter, if set, tells the program
+the top end of the period window to be searched.  If not set, the program
+uses as a value 4 times the distance between the first and last data
+point.  This parameter is remembered as minp is.
+.le
+.ls ntheta = 200
+Resolution of the theta plot.  This parameter tells how many points in
+the period window should have their theta statistic calculated.  The points
+are spaced equidistant from one another in frequency space.
+.le
+.ls pluspoint = 50
+Maximum number of data points for which to use plus symbols.  If there
+are more data points then points are plotted.
+.le
+.ls autoranges = no
+This flag, when set, instructs the program to look for gaps in
+the data and, if large gaps are found, divide the data up into ranges
+discarding the gaps and doing the analysis only on the ranges.  This
+helps remove side lobes from the spectra.
+.le
+.ls nsigma = 3
+Number of standard deviations for autorange break.  If ranges are to 
+be automatically calculated, this parameter tells how large a gap in
+the data should constitute a division between ranges.  The mean
+and standard deviation of the distribution of chronological spacing
+of input points are calculated.  Then the points are scanned in
+increasing order and if an inter-data gap bigger than nsigma
+standard deviations is found, a new range is started.
+.le
+.ls cursor = "stdgcur"
+The source of graphics cursor input.
+.le
+.ih
+DESCRIPTION
+Pdm applies a phase dispersion minimization algorithm (R. F. Stellingwerf,
+"Period Determination by Phase Dispersion Minimization", ApJ 224, 1978,
+953) to lightcurve data to determine periodicities in the data.  It also
+calculates amplitude and epoch information.
+
+Pdm can be used in batch or interactive mode.  In batch
+mode the
+output is period, amplitude, and epoch for the minimum found within
+the period window.  Metacode will be produced for the data plot,
+the theta statistic plot, and the phasecurve plot.
+The metacode will be saved in the metafile.  In interactive mode the user
+can plot the data at different stages in the analysis, fit and remove
+curves from the data, reject points, set data ranges, plot and fit
+phasecurves, etc.
+
+Pdm guesses at the period/frequency window to be searched unless
+the minimum
+and maximum period for the window are specified using minp and maxp.  The
+minimum period is taken as twice the chronological distance between the closest
+two points in the data.  The maximum period is taken as 4 times the distance
+between the first and last data points.
+
+Pdm will work on one object at a time and the input data may
+be contained in multiple input files if desired.  The program will
+concatenate data in all the files which match the input template.
+The input files are text files containing one (x,y) pair per line or
+just a (y) value per line.  If only one value per line is found the
+program will number x sequentially (1,2,3,4,...).  If a third value
+is included on each line it will be read as the error in that
+measurement.   (The 'e' key is used to toggle error bars on the phase
+plot.)
+
+At startup, if the interactive flag is set, the user will be presented
+with a plot of the data and the cursor will be turned on.
+
+When the user plots a phasecurve, points that are deleted or undeleted from
+the phasecurve plot will be deleted or undeleted from the working data set.
+
+The ICFIT keystrokes are described elsewhere. (see help for icfit)
+
+
+Phase Dispersion Minimization User Interface (keystrokes)
+
+When the program starts up it reads the data file(s) and displays
+the data on the screen as a standard mark plot.  The user is
+then placed in a graphics cursor loop with the following options
+available in addition to the standard graphics commands:
+
+Note:
+The remembered period is for the last minimum found.  This
+minimum calculation is done whenever a new theta plot is graphed
+and whenever the "m" key is used.
+
+.ls ? -- list options
+
+Print out the menu.
+.le
+.ls h -- graph data
+
+Make a plot on the screen, using marks, of observation time vs observed
+value. If there are more than 50 points, use dots, else use pluses.  If
+points have been deleted, draw an x through them on the plot.  If ranges
+are in effect, draw range bars along the abscissa of the plot marking
+the ranges.
+.le
+.ls e -- toggle error bars on or off
+
+When the phase plot is on the screen and error data has been supplied,
+this key will toggle the drawing of error bars on the phase plot so that
+the user can determine how well the period found works with the data
+including this error information.
+.le
+.ls i,k -- graph frequency or period respectively versus theta
+
+Calculate the theta statistic in the period/frequency range specified
+previously.  If no period/frequency range has been specified,
+pdm guesses one.  The minimum period is taken as twice the chronological
+distance between the closest two points in the data.  The maximum
+period is taken as 4 times the distance between the first and last
+data points.  The number of theta points in this range is also a
+parameter which can be specified.
+
+Next, plot theta on the screen using line drawing mode.  Plot
+either period vs theta or frequency vs theta.  Calculate the minimum
+value of theta displayed, turn the cursor back on (clgcur) and put
+the cursor x position at that minimum.
+.le
+.ls p -- graph phase curve for period/frequency at cursor position
+
+Calculate the phase curve for the period/frequency under the
+cursor.  This assumes the user has a theta plot on the screen and
+an error message will be given otherwise.
+
+The phase curve will be plotted in mark mode with two copies displayed
+and placed end to end to clarify the plot by providing continuity at
+all phases.  The amplitude and epoch values for this period are calculated
+and the phases are plotted relative to this epoch.
+.le
+.ls d,u -- delete/undelete respectively point nearest the cursor
+
+Points deleted will have an x drawn through them.  The x will be
+erased when the point is undeleted.
+.le
+.ls f -- call ICFIT on displayed data
+
+ICFIT is used for interactive curve fitting.
+It is called with either the data values or the phase values,
+depending on which type of plot is on the screen at the time.
+Any point deleted in ICFIT will be removed from consideration in
+all subsequent calculations until restored.
+
+The fit curve is retained by PDM after the return from ICFIT and
+may be subsequently subtracted from the data using the j command.
+
+Note: The user must exit ICFIT using the q key before he is placed
+back into PDM.
+.le
+.ls j -- subtract fit from data, use residuals
+
+Just as it says. The original data is retained and can be reinstated
+with the :origdata command.  This command only applies to the data.
+The user cannot subtract a fit from the phase plot.
+.le
+.ls s -- set sample range for calculations
+
+This command is used to set ranges of data to be used.  The cursor is
+first positioned to the beginning of the range of interest, an s is
+struck, the program prints the prompt again:, the cursor is
+repositioned to the end of the range and a second s is struck
+completing the command.  Multiple ranges may be set and all the data
+inside the union of the ranges will be used.  Data points outside the
+ranges will be ignored until the data is reset with an :alldata
+or an :origdata command.
+
+This also forces the boolean flag segments to be set true.
+.le
+.ls ,, -- Set minp or minf to cursor x position
+
+When the theta plot is on the screen, this keystroke can be used
+to set the minimum period (frequency) to the current cursor position.
+.le
+.ls . -- Set maxp or maxf to cursor x position
+
+When the theta plot is on the screen, this keystroke can be used
+to set the maximum period (frequency) to the current cursor position.
+.le
+.ls g -- significance of theta at cursor x position
+
+The statistical significance of the period/frequency under the
+cursor is calculated by Fisher's method of randomization.
+This value is printed at the bottom of the screen.
+
+This assumes that a theta plot is on the screen.
+.le
+.ls a -- amplitude and epoch at cursor x position
+
+For the period/frequency under the cursor or of the plot, the amplitude
+and epoch are calculated and returned to the user.
+
+This assumes that a theta plot is on the screen.
+.le
+.ls m -- mark range and find minimum in this range
+
+This command is used exactly like the s command but has a different
+effect.  After the user has positioned the cursor and struck the m
+key twice, defining the range, the minimum value of theta is found
+in this range and its associated period/frequency is returned.
+.le
+.ls r -- replot
+
+Redraw the plot on the screen.
+.le
+.ls x -- remove a trend from the data by removing a bestfit line
+
+This command calls the curfit package to fit a straight line to the
+data and then subtracts it point by point from the data.
+.le
+.ls z -- flip the y-axis scale
+
+This command toggles a y-axis flip for the plots.  This is useful when
+the user is plotting magnitudes where the smaller the ordinate value the
+larger the intensity.
+.le
+.ls q -- quit
+
+Exit PDM.
+
+.le
+The following commands may be abbreviated.  If entered without an
+argument; :minp, :maxp, :minf, :maxf, and :ntheta will display the named
+parameter; :show, :vshow will print to STDOUT; :signif, :ampep, and :phase,
+will do the calculation at the remembered period.
+
+.ls :show [file]		show parameter settings
+
+Print on the screen the min/max period, the remembered minimum,
+the range if it is in effect, the start and end of the ranges
+if they are defined, the mean and variance of the data in each
+range. If file is specified, the output will go there.
+.le
+.ls :vshow [file]		show verbose information
+
+This command will display all the information displayed by the :show
+command plus curfit information if the any curves have been fit.  Also,
+the residual data will be shown if residuals have been calculated. If
+file is specified, the output will go there.
+.le
+.nf
+
+:minp :maxp [period]		set min/max search period
+:minf :maxf [frequency]		set min/max search frequency
+.fi
+.ls
+These commands are self explanatory.  Whichever value is set,
+period or frequency, the corresponding frequency or period is
+automatically calculated and remembered.
+.le
+.ls :ntheta [num]		set number of points for theta
+
+Set the number of equally spaced points in the period window for
+which theta should be calculated.  This is really a setting of
+the resolution of the theta plot and defaults to 200 since
+the calculation time for 200 points is only a few seconds.  Very
+large numbers entered here will cause the program to warn the user
+that the theta calculation may take some time.
+.le
+.ls :sample [value]		set/show the sample ranges
+
+The start and end values for the ranges will be printed on the screen.
+If value is present, it has the form begin:end where begin
+and end are real numbers specifying a new range.
+.le
+.ls :signif [period]		find theta significance
+
+Same as the g key.  The colon command allows the user to 
+set the period exactly, instead of using the cursor.  If no period
+is entered, the calculation will be done using the remembered period.
+.le
+.ls :ampep [period]		amplitude and epoch
+
+Same as the e key.  Without an argument, use remembered minima.
+.le
+.ls :phase [period]		graph phase curve
+
+Same as the h key.  Without an argument, use remembered minima.
+.le
+.ls :unreject			unreject all points
+
+This tells the program to use all of the data points. If a fit
+has been subtracted from a subset of the data points then this command
+causes the original data set to be restored since, otherwise, we would
+restore a mixture of data and residuals.
+.le
+.ls :alldata			reset range to entire dataset
+
+The effect of this command is to turn off the range settings.  All
+of the data will be used if the ranges settings are off.  Rejected
+points remain rejected though.  Again, if these data are residuals,
+the original data are restored.
+.le
+.ls :origdata			reset data to original dataset
+
+Copy the original data vector into the working data vector.
+.le
+.ih
+EXAMPLES
+1. To find the main period in the data contained in the file 'vstar645',
+whose period is within the bounds (200., 800.) interactively
+the command might be:
+
+	cl> pdm vstar645 minp=200. maxp=800.
+
+2. To do the same thing in batch mode, allowing the program to guess the 
+period window, with no lightcurve analysis, and saving the metacode
+in vstar645.m, the command might be:
+
+	cl> pdm vstar645 inter=no meta="vstar645.m"
+
+.ih
+BUGS
+Pdm has some problems with data sets containing a small number (<20)
+points.  Generally, it will do fairly well but the theta curve may look
+strange.
+
+The amplitude and epoch calculation might be improved by fitting a parabola
+to the phase curve near the minimum and near the maximum and using points
+on these parabolas for the min and max points instead of actual data points.
+
+.ih
+SEE ALSO
+icfit
+.endhelp
diff --git a/noao/astutil/doc/precess.hlp b/noao/astutil/doc/precess.hlp
new file mode 100644
index 00000000..64475c56
--- /dev/null
+++ b/noao/astutil/doc/precess.hlp
@@ -0,0 +1,63 @@
+.help precess Oct87 noao.astutil
+.ih
+NAME
+precess -- general astronomical coordinate precession
+.ih
+USAGE
+precess files startyear endyear
+.ih
+PARAMETERS
+.ls files
+The name of a file (or a file list or template) containing the coordinates
+to be precessed.
+.le
+.ls startyear
+The default equinox of the input coordinates.
+.le
+.ls endyear
+The default target year to which the coordinates will be precessed.
+.le
+.ls stdepoch = 0
+If nonzero, coordinates will be output precessed to both \fBendyear\fR
+and the specified standard epoch.
+.le
+.ih
+DESCRIPTION
+Coordinates are read from the input file as RA and DEC pairs,
+one pair per input line.  Each coordinate pair may optionally be followed
+by the equinox of the input coordinates (if different from the default)
+and the epoch of the output coordinates.
+Coordinates may be entered in either decimal or sexagesimal notation.
+The given coordinates are rotated according to the
+precession rates to the requested year and printed on the standard output.
+Basic data is taken from the Explanation to the American Ephemeris.
+.ih
+EXAMPLES
+Precess coordinate entered interactively from 1950 to 1990, except where
+the dates are specified otherwise on the command line (lines input by the
+user are marked:
+
+.nf
+	cl> precess STDIN 1950 1990			[input]
+	12:30:10.12 10:18:27.5				[input]
+	  12:32:11.79   10:05:13.09  1990.0
+	12:30 10:18					[input]
+	  12:32:01.68   10:04:45.51  1990.0
+	12:30 -20 1900					[input]
+	  12:34:42.89  -20:29:46.29  1990.0
+	12:30 -20 1900 2000				[input]
+	  12:35:14.40  -20:33:04.40  2000.0
+	(eof=<ctrl/z>)					[input]
+.fi
+
+The following is equivalent, except that coordinate input is taken from
+the file "coords", rather than from the terminal:
+
+.nf
+	cl> precess coords 1950 1990			[input]
+	  12:32:11.79   10:05:13.09  1990.0
+	  12:32:01.68   10:04:45.51  1990.0
+	  12:34:42.89  -20:29:46.29  1990.0
+	  12:35:14.40  -20:33:04.40  2000.0
+.fi
+.endhelp
diff --git a/noao/astutil/doc/rvcorrect.hlp b/noao/astutil/doc/rvcorrect.hlp
new file mode 100644
index 00000000..0609e725
--- /dev/null
+++ b/noao/astutil/doc/rvcorrect.hlp
@@ -0,0 +1,373 @@
+.help rvcorrect Nov90 astutil
+.ih
+NAME
+rvcorrect -- Compute radial velocity corrections
+.ih
+USAGE
+rvcorrect
+.ih
+PARAMETERS
+.ls files = ""
+List of files containing date, time, coordinates of observation, and possibly
+an observed radial velocity.
+.le
+.ls images = ""
+List of images containing date, time, coordinates of observation, and possibly
+an observed radial velocity.
+.le
+.ls header = yes
+Print header for output?
+.le
+.ls input = no
+Print input data in output?
+.le
+.ls imupdate = no
+Update the image header with the computed values of heliocentric correction
+(in the \fIVHELIO\fR keyword), Heliocentric Julian Date (in the \fIHJD\fR
+keyword), Local Standard of Rest velocity (in the \fIVLSR\fR keyword), and
+information describing the solar motion with respect to the desired standard
+of rest (in the \fIVSUN\fR keyword).
+.le
+
+.ls epoch = INDEF
+Epoch of observation coordinates in Julian years. If zero or INDEF then the
+epoch is assumed to be the same as the date of observation.
+.le
+.ls observatory = ")_.observatory"
+Observatory for  which corrections are to be computed.  The default is a
+redirection to look in the parameters for the parent package for a value.
+This may be one of the observatories in the observatory database,
+"observatory" to select the observatory defined by the environment variable
+"observatory" or the parameter \fBobservatory.observatory\fR, or "obspars"
+to select the current parameters set in the \fBobservatory\fR task.  See
+help for \fBobservatory\fR for additional information.  If the input
+consists of images then the observatory is defined by the OBSERVAT keyword
+if present.
+.le
+.ls vsun = 20.
+Velocity in km/s of the sun relative to the desired standard of rest.  The
+default is for the Local Standard of Rest (LSR).
+.le
+.ls ra_vsun = 18:00:00
+Right ascension in hours of the solar motion relative to the desired standard
+of rest.  The default is for the Local Standard of Rest (LSR).
+.le
+.ls dec_vsun = 30:00:00
+Declination in degrees of the solar motion relative to the desired standard
+of rest.  The default is for the Local Standard of Rest (LSR).
+.le
+.ls epoch_vsun = 1900.
+Epoch in years for the solar motion components.
+.le
+
+If no input files or images are specified then the following parameters
+are used for input.
+.ls year, month, day, ut
+Date and time of observation.  If the year is less than 100 then the century is
+assumed to be 1900.  The month is specified as an integer between 1 and 12.
+The date of observation is the Greenwich date; i.e. the new day begins at
+0 hours universal time.  Universal time of observation in hours.
+.le
+.ls ra , dec 
+Right ascension (hours) and declination (degrees) of observation.
+.le
+.ls vobs = 0.
+Observed velocity (km/s) to be corrected.
+.le
+.ls keywpars = ""
+The image header keyword translation table as described in
+the \fIkeywpars\fR named pset.
+.le
+
+If no input files or images are specified the following parameters are
+set by the task.
+.ls hjd
+Heliocentric Julian date.  The date of observation is corrected for
+light travel difference to the sun.
+.le
+.ls vhelio
+Heliocentric radial velocity in km/s.  The observed velocity is corrected
+for the rotation of the Earth, the motion of the Earth about the Earth-Moon
+barycenter, and the orbit of the barycenter about the Sun.
+.le
+.ls vlsr
+Local standard of rest radial velocity in km/s.
+The heliocentric radial velocity is corrected for the motion of the Sun
+relative to the specified standard of rest.
+.le
+.ih
+DESCRIPTION
+The observed radial velocity is corrected for the motion of the
+observer in the direction of the observation.  The components of the
+observer's motion corrected are those due to the Earth's rotation
+(diurnal velocity), the motion of the Earth's center about the
+Earth-Moon barycenter (lunar velocity), the motion of the Earth-Moon
+barycenter about the center of the Sun (annual velocity), and the
+motion of the Sun (solar velocity) relative to some specified standard
+of rest.
+
+The input parameters consist of the date and time of the observation, the
+direction of observation, the location of the observation, the direction
+and magnitude of the solar motion relative to some standard of rest, and
+the observed radial velocity.  In all cases years between 0 and 99 are
+treated as 20th century years.  The observatory for the observations
+defaults to that specified by the environment variable "observatory" if
+defined or that set for the task \fBobservatory\fR.  If the input consists
+of images the observatory is defined by the OBSERVAT image header parameter
+if present.  See \fBovservatory\fR for additional information.  The solar
+motion defaults to that relative to the galactic local standard of rest
+(LSR).  Note that one can make the local standard of rest velocity be
+equivalent to the heliocentric velocity by setting the velocity of the
+solar motion to zero.
+
+The observed velocity, date, time, and direction of observation may be
+specified in three ways; from files, images, or the task parameters.  If a
+list of files is given then the files are read for the observation
+parameters.  The format of the files is lines containing the year, month
+(as an integer), day, universal time, right ascension, declination,
+(optional) coordinate epoch, and (optional) observed radial velocity.  If
+no file list is specified but a list of images is given then the
+observation parameters are determined from the image header parameters
+specified through the keywpars parameters.  If the observation date
+includes the time then it is used in preference to universal time keyword.
+Finally, if no list of files or images is given then the task parameters
+are used.  If no observed radial velocity is given in the file list or
+found in the image header then a value of zero is assumed.  In this case
+the corrected velocities are interpreted as the corrections to be added to
+a measured velocity to correct to the desired standard of rest.
+
+The results of the radial velocity calculations are output in three
+ways.  The velocities are always printed on the standard output with an
+optional header.  If the observation parameters are set with the task
+parameters (no file or image list) then the results are also stored in
+the parameter file.  This mechanism allows the task to be used easily
+in a script and to obtain greater precision.  If the observation
+parameters are taken from the image headers and the \fIimupdate\fR
+parameter is set, then the heliocentric
+Julian day is recorded as HJD, the heliocentric velocity as VHELIO,
+the LSR velocity as VLSR, and the velocity, ra and dec, and epoch
+of the solar motion used in VLSR is recorded as VSUN.
+
+The printed output may include the input data if desired.  This produces two
+lines per observation, one for the input data and one for the output
+velocities.  The calculated data consists of the heliocentric Julian
+date, the observed velocity, the observed heliocentric velocity, and
+the observed local standard of rest velocity.  Following this are
+component corrections for the diurnal, lunar, annual, and solar
+velocities.
+.ih
+DIURNAL VELOCITY
+The geodetic latitude to geocentric latitude correction is given by
+
+.nf
+	dlat = -(11. * 60. + 32.743000) * sin (2*lat) +
+		1.163300 * sin (4*lat) - 0.002600 * sin (6*lat)
+.fi
+
+where lat is the geodetic latitude and dlat is the additive correction.
+The distance, r, of the observer from the Earth's center in meters is given by
+
+.nf
+	r = 6378160.0 * (0.998327073 + 0.00167643800 * cos(2*lat) -
+	    0.00000351 * cos(4*lat) + 0.000000008 * cos(6*lat)) +
+	    altitude
+.fi
+
+where lat is the corrected latitude and altitude is the altitude above
+sea level.  The rotational velocity (perpendicular to the radius vector)
+in km/s is given by
+
+	v = TWOPI * (r / 1000.)  / (23.934469591229 * 3600.)
+
+where 23.934469591229 is the sidereal day in hours for 1986 and TWOPI is the
+ratio of the circumference to the radius of a circle.  The projection of
+this velocity along the line of sight is
+
+	vdiurnal = v * cos (lat) * cos (dec) * sin (ra-lmst)
+
+where lmst is the local mean sidereal time.
+.ih
+BARYCENTRIC VELOCITY
+The orbital elements of the lunar orbit are computed from the following
+interpolation formulas
+
+.nf
+	t = (JD - 2415020) / 36525.
+
+	oblq = 23.452294-t*(0.0130125+t*(0.00000164-t*0.000000503))
+	omega = 259.183275-t*(1934.142008+t*(0.002078+t*0.000002))
+	mlong = 270.434164+t*(481267.88315+t*(-0.001133+t*0.0000019))-
+	    omega
+	lperi = 334.329556+t*(4069.034029-t*(0.010325+t*0.000012))-
+	    omega
+	em = 0.054900489
+	inclin = 5.1453964
+.fi
+
+where t is the time from the Julian day 2415020 (~J1900) in Julian centuries,
+oblq is the mean obliquity of the ecliptic, omega is the longitude of the mean
+ascending node, mlong is the mean lunar longitude, lperi is the mean lunar
+longitude of perigee, em is the eccentricity of the lunar orbit, and inclin
+is the inclination of the orbit to the ecliptic.  The true lunar longitude,
+tlong, is computed from the mean longitude using the correction for the mean
+anomaly to the true anomaly (radians)
+
+.nf
+	manom = mlong - lperi
+	tanom = manom + (2 * em - 0.25 * em**3) * sin (manom) +
+	    1.25 * em**2 * sin (2 * manom) + 13/12 * em**3 *
+	    sin (3 * manom)
+	tlong = tanom + lperi
+.fi
+
+The velocity of the Moon around the Earth's center in the plane of the orbit
+in km/s is
+
+.nf
+	vmoon = (TWOPI * 384403.12040) / (27.321661 * 86400) /
+	    sqrt (1. - em**2)
+.fi
+
+where 384403.12040 is the mean lunar distance (km) and 27.321661 is the mean
+lunar month (days).  The component along the line to the observation is
+
+	v = vmoon * cos (bm) * (sin (tlong-lm) - em*sin (lperi-lm))
+
+where lm and bm are the longitude and latitude of the observation
+along the lunar orbital plane relative to the ascending node using a standard
+coordinate transformation.  The barycentric velocity is that reduced by
+the ratio of the Earth's mass to the Moon's mass.
+
+	vlunar = v / 81.53
+.ih
+ANNUAL VELOCITY
+The orbital elements of the Earth's orbit are computed from the following
+interpolation formulas
+
+.nf
+	t = (ast_julday (epoch) - 2415020) / 36525.
+
+	manom = 358.47583+t*(35999.04975-t*(0.000150+t*0.000003))
+	oblq = 23.452294-t*(0.0130125+t*(0.00000164-t*0.000000503))
+	lperi = 101.22083+t*(1.7191733+t*(0.000453+t*0.000003))
+	eccen = 0.01675104-t*(0.00004180+t*0.000000126)
+.fi
+
+where t is the time from the Julian day 2415020 (~J1900) in Julian centuries,
+manom is the mean anomaly (degrees), oblq is the mean obliquity of the ecliptic
+(degrees), lperi is the mean longitude of perihelion (degrees), and
+eccen is the eccentricity of the orbit.  The true anomaly (radians) is 
+obtained from the mean anomaly (radians) by
+
+.nf
+	tanom = manom + (2 * eccen - 0.25 * eccen**3) * sin (manom) +
+	    1.25 * eccen**2 * sin (2 * manom) +
+	    13./12. * eccen**3 * sin (3 * manom)
+.fi
+
+The orbital velocity of the Earth-Moon barycenter perpendicular to
+the radius vector is given by
+
+.nf
+	v = ((TWOPI * 149598500.) / (365.2564 * 86400.)) /
+	    sqrt (1. - eccen**2)
+.fi
+
+where the semi-major axis is 149598500 km and the year is 365.2564 days.
+To compute the projection of this velocity along the line of observation
+the direction of observation (precessed to the epoch of observation)
+is converted into ecliptic latitude and
+longitude, l and b, measured from the point of the ascending node using
+a standard spherical coordinate transformation.  The component is then
+
+	vannual = v * cos(b) * (sin(slong-l) - eccen*sin(lperi-l))
+
+where the longitude of the Sun as seen from the Earth, slong, is given by
+
+	slong = lperi + tanom + 180
+.ih
+SOLAR MOTION
+The solar motion is computed by precessing the coordinates of the solar
+motion to the observation epoch and taking the appropriate component
+along the line of sight.
+.ih
+ACCURACY
+The calculations are done using IRAF double precision.
+No correction is made for the perturbation of the other planets.  The
+precession does not include nutation.  The interpolation formulas are
+only approximations.  The accuracy of the heliocentric
+velocity are better than a 0.005 of a kilometer per second.
+Relative velocities over short intervals are even better.
+.ih
+EXAMPLES
+1. For use directly without data files or images there are two common modes.
+Because of the large number of parameters the parameter values are often
+set using the task \fBeparam\fR.  Then simply execute the command
+
+	cl> rvcorrect
+
+2. To set some of the parameters on the command line
+
+	cl> rvcorrect ra=12:22:1.116 dec=15:55:16.244 ut=5:30
+
+3. To use a text file generate a file containing the year, month, day, ut,
+ra, and dec with one observation per line.
+
+.nf
+cl> type rv.obs
+1987 10 21 11:00:24  3:36:15   0:22:04
+1987 10 21 11:08:00  8:19:35  -0:51:35
+1987 10 21 11:15:47  8:35:12   6:40:29
+1987 10 21 12:12:10  9:13:20  61:28:49
+1987 10 21 12:16:03  9:27:48   9:07:08
+1987 10 21 12:20:43  9:50:45  -6:06:58
+1979  3 25 11:22:59 16:07:28 -23:37:49 0 -67.5
+cl> rvcorrect f=rv.obs > rv.dat
+cl> type rv.dat
+##   HJD          VOBS   VHELIO     VLSR   VDIURNAL   VLUNAR  VANNUAL   VSOLAR
+2447089.96358     0.00    11.07    -2.74     -0.189    0.008   11.246  -13.808
+2447089.96296     0.00    28.05    13.56      0.253    0.010   27.790  -14.498
+2447089.96813     0.00    29.04    16.64      0.262    0.011   28.770  -12.401
+2447090.00834     0.00    22.06    25.26      0.114    0.010   21.940    3.200
+2447090.00884     0.00    27.70    18.55      0.250    0.009   27.438   -9.152
+2447090.01129     0.00    23.99    13.50      0.275    0.007   23.704  -10.484
+2443957.97716   -67.50   -41.37   -31.48      0.002    0.012   26.117    9.884
+.fi
+
+4. To use observation parameters from a set of images the command is
+
+	cl> rvcorrect images=hz44.001,aboo.001 > rv.dat
+
+5. A CL loop can be used to compute a table in which one parameter varies.
+
+.nf
+	cl> for (x=0.; x<=12.; x=x+1)
+	>>> rvcorrect (ut=x, header=no)
+.fi
+
+6. To get the total velocity correction in a script the following may be done.
+
+.nf
+	rvcorrect (vobs=12.3, ra=12:33, dec=30:22, ut=5:30, > "dev$null")
+	vlsr = rvcorrect.vlsr
+.fi
+
+Note that this does not work when the task is run as a background job!
+.ih
+REVISIONS
+.ls RVCORRECT V2.11.4
+The ut keyword can be in either date plus time or hours.
+.le
+.ls RVCORRECT V2.11
+Y2K update: The date keyword can be in the full format with full
+year and time.  The time takes precedence over a time keyword.
+.le
+.ih
+ACKNOWLEDGMENTS
+Some of the formulas used were obtained by inspection of the code
+for the subroutine DOP in the program DOPSET written by R. N. Manchester
+and M. A. Gordon of NRAO dated January 1970.
+.ih
+SEE ALSO
+observatory, asttimes
+.endhelp
diff --git a/noao/astutil/doc/setairmass.hlp b/noao/astutil/doc/setairmass.hlp
new file mode 100644
index 00000000..bbb1cbb8
--- /dev/null
+++ b/noao/astutil/doc/setairmass.hlp
@@ -0,0 +1,243 @@
+.help setairmass Nov00 astutil
+.ih
+NAME
+setairmass -- update image headers with the effective airmass 
+.ih
+USAGE
+setairmass images
+.ih
+PARAMETERS
+.ls images
+The list of images for which to calculate the airmass.  The image headers may
+optionally be updated with the airmass and the mid-UT of the exposure.
+.le
+.ls observatory = ")_.observatory"
+Observatory for which airmasses are to be computed if the observatory is not
+specified in the image header by the keyword OBSERVAT. The default is a
+redirection to look in the parameters for the parent package for a value. The
+observatory may be one of the observatories in the observatory database,
+"observatory" to select the observatory defined by the environment variable
+"observatory" or the parameter \fBobservatory.observatory\fR, or "obspars" to
+select the current parameters set in the \fBobservatory\fR task. See help for
+\fBobservatory\fR for additional information. If the input consists of images
+then the observatory is defined by the OBSERVAT keyword if present.
+.le
+.ls intype = "beginning"
+The time stamp of the observation as recorded at the telescope for the time
+dependent header keywords.  The choices are the "beginning", "middle" or "end"
+of the observation.
+.le
+.ls outtype = "effective"
+The output time stamp desired for the airmass. The choices are the "effective"
+airmass, or the airmass at the "beginning", "middle" or "end" of the
+observation.
+.le
+.ls ra = "ra"
+The name of the keyword that contains the right ascension. The right ascension
+is assumed to be in hours unless ra is one of the standard CRVALn keywords in
+which case it is assumed to be in degrees.
+.le
+.ls dec = "dec"
+The name of the keyword that contains the declination in degrees.
+.le
+.ls equinox = "epoch"
+The name of the keyword that contains the equinox of the right ascension and
+declination coordinates in years.
+.le
+.ls st = "st"
+The name of the keyword containing the sidereal time in hours. 
+.le
+.ls ut = "ut"
+The name of the keyword containing the ut time.  This keyword can either
+be in date plus time format or in hours.  Note that this allows setting
+both the "date-obs" and "ut".  If no time is found then
+a time of 0hrs is used.
+.le
+.ls date = "date-obs"
+The name of the keyword that contains the UT date of the observation. The
+format should be `DD/MM/YY' (old FITS format), YYYY-MM-DD (new FITS format),
+or YYYY-MM-DDTHH:MM:SS (new FITS format with time).  If there is a time
+and no time is found in the ut keyword then it is used for the ut.
+.le
+.ls exposure = "exptime"
+The name of the keyword that contains the exposure time (in seconds) of the
+image.
+.le
+.ls airmass = "airmass"
+The name of the output keyword that will receive the computed airmass.
+.le
+.ls utmiddle = "utmiddle"
+The name of the output keyword that will receive the universal time for
+the middle of the observation.  The format of the keyword will be the same
+as that specifying the universal time.
+.le
+.ls scale = 750.0
+The atmospheric scale height.
+.le
+.ls show = yes
+Print the airmasses and mid-UT's for each image?
+.le
+.ls update = yes
+Update the image headers with the airmasses and the mid-UT's?
+.le
+.ls override = yes
+If updating the image headers, override values that were previously recorded ?
+.le
+
+.ih
+DESCRIPTION
+
+SETAIRMASS will calculate the effective airmass of an astronomical image, as
+described below under "ALGORITHMS".  The task requires the instantaneous
+zenith distance at the beginning, middle and end of the exposure. These are
+calculated using the right ascension, declination, and equinox as well as the
+sidereal time, exposure time, UT date, and observatory from the image header.
+If the observatory is not available in the image header under the keyword
+OBSERVAT, the observatory is defined by the \fIobservatory\fR parameter. See
+help for \fIobservatory\fR for further information.
+
+The right ascension and declination will be precessed from the given equinox to
+the date of observation. The name of the right ascension, declination, equinox,
+sidereal time, ut time, exposure time, and date keywords can be specified as
+parameters. These keywords should express the right ascension in hours,
+the declination in degrees, the equinox in years, the sidereal time in hours,
+the universal time in hours, the exposure time in seconds, and the date in
+FITS format. If any of the required keywords are missing from the image
+headers, they can be added using the hedit or asthedit tasks.  Note that
+the universal time keyword may be in either a date plus time format or
+in hours and any output middle universal time will be in the same format.
+
+Before using this task, you will need to know the "time stamp" of the time
+varying header quantities (e.g. sidereal time).  Do the recorded values
+represent the beginning, the middle or the end of the exposure ? This should
+be set in the \fBintype\fR parameter.
+
+If for some reason the effective airmass is not desired, the value of the
+airmass at the beginning, middle or end of the exposure can be recorded in the
+header keyword specified by the \fIairmass\fR parameter. The \fBshow\fR
+parameter can be used to control the output to the terminal. The \fBupdate\fR
+and \fBoverride\fR parameters control the header keyword output.
+
+SETAIRMASS will also calculate the universal time of the middle of the exposure
+and place the value in the header keyword specified by the \fIutmiddle\fR
+parameter.  This assumes that the value for the UT is in the date keyword
+or ut keyword, with the same time stamp as the sidereal time. The
+mid-observation UT is useful for interpolating calibration arc dispersion
+solutions using REFSPECTRA, especially when the exposure time is
+long.
+
+.ih
+ALGORITHMS
+The mean airmass is calculated uses the formula described in "Some
+Factors Affecting the Accuracy of Stellar Photometry with CCDs" by P.
+Stetson, DAO preprint, September 1988.  This simple formula is:
+
+.nf
+	    AM (eff) = [AM (beginning) + 4*AM (middle) + AM (end)] / 6
+.fi
+
+and is derived by using Simpson's 1/3 rule to approximate the integral
+that represents the mean airmass.
+
+The beginning, middle and end airmasses are calculated using the
+relation between airmass and elevation (or zenith distance) in John
+Ball's book on Algorithms for the HP-45:
+
+.nf
+	    AM = sqrt (x**2 + 2*scale + 1) - x, where
+
+	     x = scale * sin(elevation) = scale * cos(ZD)
+.fi
+
+The atmospheric scaling parameter is \fIscale\fR (see "Astrophysical
+Quantities" by Allen, 1973 p.125,133).
+
+.ih
+KEYWORDS
+The input keywords are:
+.ls OBSERVAT
+Observatory at which the data was taken.  If absent the observatory is
+determined using the \fIobservatory\fR parameter.
+.le
+.ls \fIra\fR
+Right ascension in hours at the beginning, middle, or end of the observation.
+If ra is one of the CRVALn keywords it is assumed to be in degrees.
+.le
+.ls \fIdec\fR
+Declination in degrees at the beginning, middle, or end of the observation.
+.le
+.ls \fIequinox\fR
+The equinox of the coordinates.  The right ascension and declination will
+be precessed from this epoch to the date of the observation before being
+used.
+.le
+.ls \fIst\fR
+Sidereal time in hours at the beginning, middle, or end of the observation.
+.le
+.ls \fIut\fR
+Universal time in hours at the beginning, middle, or end of the observation.
+This may be in either date plus time format or just in hours.  
+.le
+.ls \fIdate\fR
+The value of the date parameter is the keyword name to be used for the date of
+the observation.  The date must be in either the old or new FITS format.
+.le
+.ls \fIexposure\fR
+The value of the exposure parameter is the keyword name to be used for the
+exposure time in seconds.
+.le
+
+The output keywords are:
+.ls \fIairmass\fR
+The value of the airmass parameter is the keyword name to be used for
+the computed airmass at either the beginning, middle, or end of the
+exposure, or for the weighted effective value over the exposure.
+.le
+.ls \fIutmiddle\fR
+The value of the utmiddle parameter is the keyword name to be used for
+the universal time at the middle of the exposure.
+.le
+
+.ih
+EXAMPLES
+
+1. Calculate the effective airmass of the IRAF test picture, dev$pix.
+
+.nf
+    cl> setairmass dev$pix exposure=itime update-
+.fi
+
+Note that the test picture does not have the correct coordinate epoch
+listed in its header, so no precession will be performed. 
+
+2. Calculate the effective airmass of the IRAF test picture dev$ypix in two
+ways.
+
+.nf
+    cl> setairmass dev$ypix exposure=itime update-
+
+    cl> setairmass dev$ypix ra=crval1 dec=crval2 equinox=equinox \
+        exposure=itime update-
+.fi
+
+Note the first way gives the same results as example 1. The second way
+uses the J2000 equatorial system rather then the ra and dec at the time
+of observation.
+
+.ih
+REVISIONS
+.ls SETAIRMASS V2.11.4
+The ut keyword now has precedence over any time in the date keyword
+and it can be either date plus time or hours.
+.le
+.ls SETAIRMASS V2.11.3
+The right ascension, declination, equinox, st, and ut keywords were made 
+parameters rather than being hard wired.
+.le
+.ls SETAIRMASS V2.11.2
+Y2K update: This task was updated to use the new FITS date format.
+.le
+.ih
+SEE ALSO
+airmass, hedit, refspectra, observatory
+.endhelp
diff --git a/noao/astutil/doc/setjd.hlp b/noao/astutil/doc/setjd.hlp
new file mode 100644
index 00000000..c010d533
--- /dev/null
+++ b/noao/astutil/doc/setjd.hlp
@@ -0,0 +1,222 @@
+.help setjd Jan92 astutil
+.ih
+NAME
+setjd -- set various Julian dates in image headers
+.ih
+USAGE
+setjd images
+.ih
+PARAMETERS
+.ls images
+The list of images for which to compute Julian dates.  If the \fIlistonly\fR
+parameter is not set the image headers will be modified to add or update
+the calculated Julian date values.
+.le
+.ls observatory = ")_.observatory"
+Observatory of observation, used to define the time zone relative to
+Greenwich, if not specified in the image header by the keyword OBSERVAT.
+The default is a redirection to look in the parameters for the parent
+package for a value.  The observatory may be one of the observatories in
+the observatory database, "observatory" to select the observatory defined
+by the environment variable "observatory" or the parameter
+\fBobservatory.observatory\fR, or "obspars" to select the current
+parameters set in the \fBobservatory\fR task.  See \fBobservatory\fR for
+additional information.
+.le
+.ls date = "date-obs"
+Date of observation keyword.  The value must be in FITS format.
+This is one of DD/MM/YY (old FITS format), YYYY-MM-DD (new FITS format),
+or YYYY-MM-DDTHH:MM:SS (new FITS format with time).  The date should be
+in universal time though the \fIutdate\fR parameter can be used if
+this is not the case.  If a time is included it is used in preference
+to the \fItime\fR value.
+.le
+.ls time = "ut"
+Time of observation keyword with value given in decimal hours or HH:MM:SS.S
+format.  The date may be a local time or universal time as selected by the
+\fIuttime\fR parameter.  The time is used as is
+unless an exposure time keyword is specified, in which case
+the time will be corrected to the midpoint of the exposure from the
+beginning or end of the exposure.  This time is not used if a time
+is given in the date keyword.
+.le
+.ls exposure = "exptime"
+Exposure time keyword with value in seconds.  If specified the time
+is corrected to the midpoint of the exposure.  The time is assumed
+to be the beginning of the exposure unless
+the exposure time keyword name begins with a minus sign, for example
+"-exptime", in which case the time is used as the end of the exposure.
+.le
+.ls ra = "ra", dec = "dec", epoch = "epoch"
+If the heliocentric Julian date is requested the right ascension (in hours)
+and declination (in degrees) of the observation is determined from these
+keywords.  The values may be in either decimal or sexagesimal notation.
+An epoch keyword is optional and if given is used to precess
+the coordinates from the specified epoch to the observation epoch.
+If an epoch keyword is given but is not found in the header or can't
+be interpreted then it is an error.  The epoch keyword value may begin
+with 'B' or 'J'.  If the value is before 1800 or after 2100 a warning
+will be printed though the task will still compute the values.
+.le
+
+.ls jd = "jd"
+If specified compute the geocentric Julian day (Greenwich) at the
+midpoint of the exposure and record the value in the specified
+header keyword.
+.le
+.ls hjd = "hjd"
+If specified compute the heliocentric Julian day (Greenwich) at the
+midpoint of the exposure and record the value in the specified
+header keyword.
+.le
+.ls ljd = "ljd"
+If specified compute the local Julian day number.  This is an integer
+number which is constant for all observations made during the same night.
+It may be used to group observations by night in such tasks as
+\fBrefspectra\fR.
+.le
+
+.ls utdate = yes, uttime = yes
+Define whether the date and time of observation are in local standard
+time or in universal time.
+.le
+.ls listonly = no
+List the computed values only and do not modify the image headers.
+When simply listing the images need not have write permission.
+.le
+.ih
+DESCRIPTION
+\fBSetjd\fR computes the geocentric, heliocentric, and integer
+local Julian dates from information given in the headers of
+the input list of images.  This information may simply be listed or
+the values may be added or modified in the image headers.  Only
+those values which have a keyword specified are computed, printed,
+and entered in the images.  Thus, one need not compute all values
+and the dependent image header parameters required for computing them
+need not be present.  For example, if the coordinates of the
+observation are not available one should set the \fIhjd\fR parameter
+to an empty string.
+
+Often the date and time of observation are recorded either at the
+beginning or the end of an exposure.  To compute the Julian dates
+at the midpoint of the exposure the exposure keyword is specified.
+A negative sign preceding the keyword name defines correcting from
+the end of the exposure otherwise the correction is from the
+beginning of the exposure.  The exposure time must be in seconds and
+there is no allowance made for exposures which are interrupted.
+See also the task \fBsetairmass\fR which may be used to compute a
+universal time midexposure value.
+
+The date and time of observations should be given either in universal
+time.  However, if they are given in local standard time (there is no
+provisions for daylight savings times) the \fIutdate\fR and \fIuttime\fR
+parameters may be used.  Conversion between local and universal times, as
+well as the computation of the local integer date, requires the time zone
+in (positive) hours behind Greenwich or (negative) hours ahead of
+Greenwich.  This information is determined from the observatory at which
+the observations were made.  If the observatory is specified in the image
+header under the keyword OBSERVAT with a value which has an entry in the
+NOAO, local, or user observatory database then the value from the database
+is used.  This is the safest way since the observatory is tied to the
+actual image.  Otherwise, the \fIobservatory\fR parameter defines the
+observatory.  The special value "observatory" allows defining a default
+observatory with an environment variable or the \fBobservatory\fR task.
+Explicitly use the parameter \fIobservatory.timezone\fR use the value
+"obspars".  For more information see help under \fBobservatory\fR.
+
+The heliocentric Julian date is computed by defining a keyword for
+this value and also defining the keywords for the right ascension (in hours)
+and declination (in degrees).  An optional epoch keyword may be
+used if the RA and DEC are not for the observation epoch.
+
+The local integer Julian day number is the Julian date which begins at
+local noon.  Thus, all observations made during a night will have the
+same day number.  This day number may be useful in grouping
+observations by nights.  Note that in some time zones the UT
+date of observation may also be constant over a night.
+
+Among the uses for this task is to define keywords to be used by the task
+\fBrefspectra\fR.  In particular, the exposure midpoint geocentric Julian
+date makes a good sort parameter and the local Julian day number makes a
+good group parameter.
+.ih
+EXAMPLES
+1.  Compute all the Julian date quantities for 4 arc exposures with
+header parameters given below.
+
+.nf
+    demoarc1:
+	OBSERVAT= 'KPNO              '  /  observatory
+	EXPTIME =                  60.  /  actual integration time
+	DATE-OBS= '26/11/91          '  /  date (dd/mm/yy) of obs.
+	UT      = '12:11:30.00       '  /  universal time
+	RA      = '06:37:02.00       '  /  right ascension
+	DEC     = '06:09:03.00       '  /  declination
+	EPOCH   =               1991.9  /  epoch of ra and dec
+
+    demoarc2:
+	OBSERVAT= 'KPNO              '  /  observatory
+	EXPTIME =                  60.  /  actual integration time
+	DATE-OBS= '26/11/91          '  /  date (dd/mm/yy) of obs.
+	UT      = '12:41:30.00       '  /  universal time
+	RA      = '06:37:02.00       '  /  right ascension
+	DEC     = '06:09:03.00       '  /  declination
+	EPOCH   =               1991.9  /  epoch of ra and dec
+
+    demoarc3:
+	OBSERVAT= 'CTIO              '  /  observatory
+	EXPTIME =                  60.  /  actual integration time
+	DATE-OBS= '27/11/91          '  /  date (dd/mm/yy) of obs.
+	UT      = '11:11:30.00       '  /  universal time
+	RA      = '06:37:02.00       '  /  right ascension
+	DEC     = '06:09:03.00       '  /  declination
+	EPOCH   =               1991.9  /  epoch of ra and dec
+
+    demoarc4:
+	OBSERVAT= 'CTIO              '  /  observatory
+	EXPTIME =                  60.  /  actual integration time
+	DATE-OBS= '27/11/91          '  /  date (dd/mm/yy) of obs.
+	UT      = '12:21:30.00       '  /  universal time
+	RA      = '06:37:02.00       '  /  right ascension
+	DEC     = '06:09:03.00       '  /  declination
+	EPOCH   =               1991.9  /  epoch of ra and dec
+
+    cl> setjd demoarc?.imh
+    # SETJD: Observatory parameters for Kitt Peak ...
+    #              Image            JD           HJD   LOCALJD
+	    demoarc1.imh  2448587.0083  2448587.0127   2448586
+	    demoarc2.imh  2448587.0292  2448587.0336   2448586
+    # SETJD: Observatory parameters for Cerro Tololo ...
+	    demoarc3.imh  2448587.9667  2448587.9711   2448587
+	    demoarc4.imh  2448588.0153  2448588.0197   2448587
+.fi
+
+Note the use of the observatory parameter to switch observatories and
+the local Julian day number which is constant over a night even though
+the Julian date may change during the observations.
+
+2.  To compute only the geocentric Julian date from the "DATE" and
+"TIME" keywords in an image,
+
+.nf
+    cl> setjd obs1 date=date time=time exp="" hjd="" ljd=""
+.fi
+.ih
+REVISIONS
+.ls SETJD V2.11.2
+Y2K update:  Updated to use the new FITS format for the date.  If the
+time is given in the date keyword it is used in preference to the
+time keyword.
+.le
+.ls SETJD V2.11
+The checking of the epoch keyword value was improved.  Previously if
+there was a problem with the keyword value (missing or malformed) the
+task would use the epoch of the observation.  Now it is an error
+if an epoch keyword is specified but the epoch value can't be determined.
+Also a leading 'B' or 'J' is allowed and a warning will be given if
+the epoch value is unlikely.
+.le
+.ih
+SEE ALSO
+setairmass, hedit, refspectra, observatory
+.endhelp