Index: /trunk/doc/pslib/psLibSDRS.tex
===================================================================
--- /trunk/doc/pslib/psLibSDRS.tex	(revision 626)
+++ /trunk/doc/pslib/psLibSDRS.tex	(revision 627)
@@ -1,3 +1,3 @@
-%%% $Id: psLibSDRS.tex,v 1.40 2004-05-04 04:13:29 price Exp $
+%%% $Id: psLibSDRS.tex,v 1.41 2004-05-08 03:15:51 price Exp $
 \documentclass[panstarrs]{panstarrs}
 
@@ -82,8 +82,9 @@
 The installed base of code for PSLib consists of header files, the
 binary library code, \code{libpslib.a} and the shared-library
-equivalent, \code{libpslib.so}.  Assuming these components have been
-installed into the library and search path, PSLib may be used within a
-program by including the line \code{#include <pslib.h>} into the C
-code and linking with \code{-lpslib}.
+equivalent, \code{libpslib.so} (or \code{libpslib.dylib} in the case
+of OS/X).  Assuming these components have been installed into the
+library and search path, PSLib may be used within a program by
+including the line \code{#include <pslib.h>} into the C code and
+linking with \code{-lpslib}.
 
 This document describes the data structures and details the functions
@@ -113,8 +114,9 @@
 (\href{heasarc.gsfc.nasa.gov/docs/software/fitsio/}{\tt
 heasarc.gsfc.nasa.gov/docs/software/fitsio/}); and
-\item Many of the astronomy routines will wrap the StarLink Astronomy
-libraries (SLALib;
+\item Many of the astronomy routines will wrap the StarLink Positional
+Astronomy libraries (SLALib;
 \href{star-www.rl.ac.uk/star/docs/sun67.htx/sun67.html}{\tt
 star-www.rl.ac.uk/star/docs/sun67.htx/sun67.html}.
+\item \tbd{Some graphics library, possibly the SM library.}
 \end{itemize}
 
@@ -486,6 +488,6 @@
 management routines respect the use of the \code{refCounter} field:
 \code{psFree} will not free a block for which \code{refCounter != 1},
-and \code{psAlloc} will initialize the field to 1.  \tbd{\code{psFree}
-should generate an error if \code{refCounter != 1}.}  However, they do
+and \code{psAlloc} will initialize the field to 1.  \code{psFree}
+must generate an error if \code{refCounter != 1}.  However, they do
 not (and cannot practically) enforce the use of the counters; this is
 a requirement of external APIs which intend to use the feature.
@@ -512,4 +514,6 @@
 (\S\ref{sec:psDlist}).
 
+\paragraph{Destructors}
+
 \subsubsection{Relation of Memory Management to Structures}
 \label{sec:free}
@@ -517,14 +521,37 @@
 In this document, we specify several C \code{struct}s.  It is expected
 that instances of, for example, \code{struct psMyType} will be
-constructed using \code{psMyTypeAlloc()} calls, and destructed using
+constructed using \code{psMyTypeAlloc()} calls, and destroyed using
 \code{psMyTypeFree()} calls.  The allocator will allocate the required
 memory with \code{psAlloc} and increment the appropriate
-\code{refCounter}.  The destructor will decrement the appropriate
-\code{refCounter} and free the memory with \code{psFree} \emph{if and
-only if} the new value of the \code{refCounter} is zero.
-
-This allows the \code{psMyType} to be imported into the metadata
-(\S\ref{sec:metadata}) without the user worrying about the details of
-the memory allocation/deallocation.  For example:
+\code{refCounter}.
+
+The existence of complicated structures which include pointers to
+other structures require that we lay out a rule regarding destructors
+(i.e., \code{psMyTypeFree}) and reference counters.  Simply put,
+\textit{the destructor for every structure should only call
+\code{psFree} if \code{refCounter == 1}; otherwise, it decrements the
+reference counter and returns.}  An example destructor is shown below:
+
+\begin{verbatim}
+void psMyTypeFree(psMyType *myType ///< Object to destroy
+    )
+{
+    /* No operation if object is NULL */
+    if (myType == NULL) {
+	return;
+    }
+    /* Only call psFree if reference counter is 1 */
+    if (psMemGetRefCounter(myType) == 1) {
+	psFree(myType);
+	return;
+    }
+    /* Otherwise, decrement the reference counter only */
+    psMemDecrRefCounter(myType);
+}
+\end{verbatim}
+
+This allows, for example, the \code{psMyType} to be imported into the
+metadata (\S\ref{sec:metadata}) without the user worrying about the
+details of the memory allocation/deallocation:
 
 \begin{verbatim}
@@ -596,10 +623,11 @@
 relationship is seen in two ways.  First, a facility inherits the
 trace level of its parent unless explicitly specified.  Second, the
-hierarchy is used to format the listing of the trace facilities.  The
-first of these rules provides a mechanism to define the default trace
-levels for any facility even if it has not been registered explicitly
-since all named facilities are implicitly children of the top level
-facility (\code{.}).  The second rule is simply an organizational
-technique to make the listing of facility information clear.
+hierarchy is used to format the listing of the trace facilities so
+that they are easy to read.  The first of these rules provides a
+mechanism to define the default trace levels for any facility even if
+it has not been registered explicitly since all named facilities are
+implicitly children of the top level facility (\code{.}).  The second
+rule is simply an organizational technique to make the listing of
+facility information clear.
 
 The API to place a trace message in the code, and simultaneously set
@@ -621,5 +649,5 @@
 %
 \begin{verbatim}
-int psSetTraceLevel(char *facil, int level);
+int psTraceSetLevel(char *facil, int level);
 \end{verbatim}
 % 
@@ -629,5 +657,5 @@
 %
 \begin{verbatim}
-int psGetTraceLevel(char *facil);
+int psTraceGetLevel(char *facil);
 \end{verbatim}
 % 
@@ -635,5 +663,5 @@
 rules specified above.  A specified trace message (identified by
 \code{psTrace}) shall be printed if and only if
-\code{psGetTraceLevel(facil)} returns a value greater than or equal to
+\code{psTraceGetLevel(facil)} returns a value greater than or equal to
 the value of \code{myLevel} for that message.  That is, a larger
 number for the trace level corresponds to lower-level statements, and
@@ -641,17 +669,18 @@
 
 PSLib will include a utility function for examining the current
-tracing levels of all facilities: \code{void psPrintTraceLevels(void);}.
+tracing levels of all facilities: \code{void psTracePrintLevels(void);}.
 This function will print the hierarchy of trace facilities along with
 the current trace level for each facility.  For example, a particular
 program may have a few facilities defined, along with their trace
-levels.  A call to \code{psPrintTraceLevels} may produce a listing
+levels.  A call to \code{psTracePrintLevels} may produce a listing
 which appears as:
 \begin{verbatim}
-.                    0
-.debias              1
-.flatten             1
-.flatten.divide      2
-.object.findpeak     1
-.object.getsky       1
+.                        0
+ .IPP                    0
+ .IPP.debias             1
+ .IPP.flatten            1
+  .IPP.flatten.divide    2
+  .IPP.object.findpeak   1
+  .IPP.object.getsky     1
 \end{verbatim}
 
@@ -661,10 +690,10 @@
 %
 \begin{verbatim}
-psTrace("flatten", 2, "starting flatten function\n");
-psTrace("flatten", 0, "ERROR: flat-field image \%s is invalid\n", filename);
-psTrace("flatten.divide", 2, "doing the divide\n");
-psTrace("flatten.divide", 3, "trying the loop, i = \%d\n", i);
-psTrace("flatten.divide", 1, "got an invalid pixel value (NaN) at \%d,\%d\n");
-psTrace("flatten.divide", 2, "divide is done\n");
+psTrace("IPP.flatten", 2, "starting flatten function\n");
+psTrace("IPP.flatten", 0, "flat-field image is \%s\n", filename);
+psTrace("IPP.flatten.divide", 2, "doing the divide\n");
+psTrace("IPP.flatten.divide", 3, "trying the loop, i = \%d\n", i);
+psTrace("IPP.flatten.divide", 1, "got an invalid pixel value (NaN) at \%d,\%d\n");
+psTrace("IPP.flatten.divide", 2, "divide is done\n");
 \end{verbatim}
 %
@@ -673,8 +702,8 @@
 %
 \begin{verbatim}
-ERROR: flat-field image foo.fits is invalid
-doing the divide
-got an invalid pixel value (NaN) at 500,20
-divide is done
+flat-field image is foo.fits
+  doing the divide
+ got an invalid pixel value (NaN) at 500,20
+  divide is done
 \end{verbatim}
 %
@@ -683,8 +712,8 @@
 %
 \begin{verbatim}
-starting flatten function
-trying the loop, i = 0   
-trying the loop, i = 1   
-trying the loop, i = 2   
+  starting flatten function
+   trying the loop, i = 0   
+   trying the loop, i = 1   
+   trying the loop, i = 2   
 ...
 \end{verbatim}
@@ -698,5 +727,13 @@
 \code{void psTraceReset(void)} will set all trace levels to 0.
 
-\tbd{ability to send trace messages to output locations other than stdout?}
+
+The trace may optionally be written to a file on specifying a filename
+with \code{psTraceSetDestination}:
+\begin{verbatim}
+void psTraceSetDestination(const char *filename // File to write to
+    );
+\end{verbatim}
+If the \code{filename} is \code{NULL}, then the trace is sent to
+standard output.
 
 \subsubsection{Message Logging}
@@ -722,5 +759,5 @@
 \begin{verbatim}
 void psLogMsg(char *name, int myLevel, char *fmt, ...); 
-void psVLogMsg(char *name, int myLevel, char *fmt, va_list ap); 
+void psLogVMsg(char *name, int myLevel, char *fmt, va_list ap); 
 \end{verbatim}
 where \code{name} is a word to describe the source of the message,
@@ -728,5 +765,5 @@
 is a printf-style formatting statement defining the actual message,
 and is followed by the arguments to the formatting code.  The second
-form, \code{psVLogMsg} is an equivalent command which takes a
+form, \code{psLogVMsg} is an equivalent command which takes a
 \code{va_list} argument.
 
@@ -744,5 +781,5 @@
 %
 \begin{verbatim}
-int psSetLogLevel(int level);           
+int psLogSetLevel(int level);           
 \end{verbatim}
 %
@@ -750,21 +787,23 @@
 invoked with \code{psLogMsg} is only printed if its value of
 \code{myLevel} is less than the current value set by
-\code{psSetLogLevel}.
+\code{psLogSetLevel}.
 
 Log messages are sent to the destination most recently set using:
 %
 \begin{verbatim}
-int psSetLogDestination(int dest);      
-\end{verbatim}
-%
-The only two values that are initially defined are
-\code{PS_LOG_TO_STDERR} and \code{PS_LOG_TO_STDOUT} to write to
-\code{stderr} and \code{stdout} respectively.  \tbd{define other
-destinations? why not use file descriptors?, etc?}.
+int psLogSetDestination(int dest);      
+\end{verbatim}
+%
+The only values that are initially defined are \code{PS_LOG_TO_STDERR}
+and \code{PS_LOG_TO_STDOUT} to write to \code{stderr} and
+\code{stdout} respectively, and \code{PS_LOG_TO_NONE} to turn off
+logging.  \tbd{Log to a server via TCP/IP; will possibly set these
+values to negative integers and allow psLogSetDestination to take
+a positive-integer file descriptor.}.
 
 The output format is controlled with the function:
 %
 \begin{verbatim}
-void psSetLogFormat(const char *fmt);   
+void psLogSetFormat(const char *fmt);   
 \end{verbatim}
 %
@@ -795,7 +834,11 @@
 %
 The possible order of the format entries is fixed and not determined
-by the order of the letters used in \code{psSetLogFormat}.  Selecting
+by the order of the letters used in \code{psLogSetFormat}.  Selecting
 an output format with fewer than the complete set of 5 entries simply
 removes those entries from the output messages.
+
+Specifying a \code{fmt} of \code{NULL} turns off logging (equivalent
+to calling \code{psLogSetDestination(PS_LOG_TO_NONE)}, whereas if the
+\code{fmt} is \code{""}, then the format reverts to the default.
 
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@@ -833,5 +876,5 @@
 error that has already resulted in a call to \code{psError}.  The
 final required argument, \code{fmt}, is a \code{printf}-style format
-that is passed to \code{psLogMsg} with code \code{PS_LOG_ERROR}.
+that is passed to \code{psLogVMsg} with code \code{PS_LOG_ERROR}.
 
 The result of a call to \code{psError} shall be to push an error onto
@@ -855,6 +898,6 @@
 
 \begin{verbatim}
-const psErr *psGetError(int which);     ///< return specified error (or an "error" with code PS_ERR_NONE)
-const psErr *psLastError(void);         ///< return last error (or an "error" with code PS_ERR_NONE)
+const psErr *psErrorGet(int which);     ///< return specified error (or an "error" with code PS_ERR_NONE)
+const psErr *psErrorLast(void);         ///< return last error (or an "error" with code PS_ERR_NONE)
 \end{verbatim}
 
@@ -866,5 +909,5 @@
 
 The entire error stack may be printed to an open file descriptor by
-calling \code{psErrorStackPrint} (or \code{psVErrorStackPrint}); if
+calling \code{psErrorStackPrint} (or \code{psErrorVStackPrint}); if
 and only if there are current errors, the printf-style string
 \code{fmt} is first printed to the file descriptor \code{fd}. In this
@@ -872,10 +915,10 @@
 defined in the next section.  Note that these are also available in
 the \code{psErr} structure. The successive lines of the traceback
-should be indented by an additional space (see
-example). \code{psVErrorStackPrint} shall not invoke \code{va_end}.
+should be indented by an additional space (see example).
+\code{psErrorVStackPrint} shall not invoke \code{va_end}.
 
 \begin{verbatim}
 void psErrorStackPrint(FILE *fd, const char *fmt, ...); ///< print the errorstack to this file descriptor
-void psVErrorStackPrint(FILE *fd, const char *fmt, va_list va); ///< print the errorstack to this file
+void psErrorVStackPrint(FILE *fd, const char *fmt, va_list va); ///< print the errorstack to this file
                                                                 ///< descriptor
 \end{verbatim}
@@ -916,6 +959,6 @@
 \begin{verbatim}
 typedef struct {
-    psErrorCode code;                  // An error code
-    char *descrip;                     // the associated description
+    psErrorCode code;                  ///< An error code
+    const char *descrip;               ///< the associated description
 } psErrorDescription;
 \end{verbatim}
@@ -998,6 +1041,4 @@
 to be valid values of \code{errno}.
 
-The implementation may add extra fields (e.g., \code{PS_ERR_N_ERR_CLASSES}).
-
 \file{psErrorCodes.c} shall define the necessary function to register
 the error codes.
@@ -1042,8 +1083,8 @@
         psErrorStackPrint(stdout, "Traceback:\n");
 
-        if (psLastError()->code == PS_ERR_UNKNOWN) {
+        if (psErrorLast()->code == PS_ERR_UNKNOWN) {
             fprintf(stderr, "Last error is of unknown type\n");
         }
-        if (psGetError(2)->code == PS_ERR_IO) {
+        if (psErrorGet(2)->code == PS_ERR_IO) {
             fprintf(stderr, "Third oldest error is of type IO\n");
         }
@@ -1053,5 +1094,5 @@
     psErrorStackPrint(stdout, "Traceback:\n");
 
-    if (psLastError()->code == PS_ERR_NONE) {
+    if (psErrorLast()->code == PS_ERR_NONE) {
         fprintf(stderr, "No errors. Hurrah\n");
     }
@@ -1178,5 +1219,5 @@
 \end{verbatim}
 We discuss the application of \code{psType} in more detail in
-section~\ref{arithmetic}.  
+section~\ref{sec:arithmetic}.  
 
 \subsection{Arrays of Simple Types}
@@ -1269,6 +1310,8 @@
 array.\footnote{\eg{} \code{va->arr[i] = psMemIncrRefCounter(ptr);}}
 If \code{psVoidPtrArrayFree}'s argument \code{elemFree} is NULL, the
-list should be deleted, but not the elements on it (although their
-\code{refcounter}'s should be decremented).
+list should be deleted, but not the elements on it (an error should be
+raised if the \code{refCounter} is 1; otherwise their
+\code{refCounter}'s should be decremented) --- this is to account for
+an array of heterogeneous types.
 
 \subsection{Doubly-linked lists}
@@ -1379,16 +1422,24 @@
 \code{iter} is provided by the functions:
 \begin{verbatim}
-void psDlistSetIterator(psDlist *list, int where);
-void *psDlistGetNext(psDlist *list);
-void *psDlistGetPrev(psDlist *list);
+void psDlistSetIterator(psDlist *list, int where, int which);
+void *psDlistGetNext(psDlist *list, int which);
+void *psDlistGetPrev(psDlist *list, int which);
 \end{verbatim}
 The first of these functions uses the value of \code{where} to set the
 iteration cursor for the given list to the beginning
-\code{PS_DLIST_HEAD} or the end \code{PS_DLIST_TAIL}.  The next two
-functions move the iteration cursor forward or backwards, returning
-the data item from the resulting list entry, or returning \code{NULL}
-at the end of the list.  Explicit traversal of the list using the
-\code{psDlistElem}s \code{prev} and \code{next} pointers is also
-supported.
+\code{PS_DLIST_HEAD} or the end \code{PS_DLIST_TAIL} for the iterator
+specified by \code{which}.  The next two functions move the iteration
+cursor forward or backwards, returning the data item from the
+resulting list entry, or returning \code{NULL} at the end of the list.
+Explicit traversal of the list using the \code{psDlistElem}s
+\code{prev} and \code{next} pointers is also supported.
+
+A list may be sorted using \code{psDlistSort}, which requires the
+specification of a comparison function to specify how the objects on
+the list should be sorted.  The motivation is primarily to be able to
+iterate on a sorted list of keys from a hash.
+\begin{verbatim}
+psDlist *psDlistSort(psDlist *list, int (*compare)(const void *a, const void *b) );
+\end{verbatim}
 
 \subsection{Hash Tables}
@@ -1411,5 +1462,5 @@
 
 For PSLib, we define a hash table and hash buckets as follows:
-\footnote{ We choose not to use the posix function \code{hcreate},
+\footnote{ We choose not to use the POSIX function \code{hcreate},
 \code{hdestroy}, and \code{hsearch} as they only support a single hash
 table at any one time.}
@@ -1550,5 +1601,5 @@
 psBitset *
 psBitsetSet(psBitset *restrict myBits, ///< Input bitset
-             int bit                     ///< Bit to set
+            int bit                    ///< Bit to set
     );
 \end{verbatim}
@@ -1582,4 +1633,14 @@
 performing the specified \code{operator} (one of \code{"AND"},
 \code{"OR"}, or \code{"XOR"}) on \code{inBits1} and \code{inBits2}.
+
+Finally, \code{psBitsetNot} applies a unary \code{NOT} to a bitset:
+
+\begin{verbatim}
+/** Apply unary NOT to a bitset */
+psBitset *
+psBitsetNot(psBitset *out,		///< Output bitset, or NULL
+	    psBitset *in		///< Input bitset to be NOT-ed
+    );
+\end{verbatim}
 
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@@ -1749,4 +1810,5 @@
 
 \subsubsection{Histograms}
+\label{sec:histograms}
 
 We also require to be able to generate histograms, given a list of
@@ -1754,11 +1816,11 @@
 
 \begin{verbatim}
-/** Histograms */
+/** Histograms  */
 typedef struct {
-    const psFloatArray *restrict lower; ///< Lower bounds for the bins
-    const psFloatArray *restrict upper; ///< Upper bounds for the bins
-    psIntArray *nums;                   ///< Number in each of the bins
-    const float minVal, maxVal;         ///< Minimum and maximum values
-    int minNum, maxNum;                 ///< Number below the minimum and above the maximum
+    const psFloatArray *restrict bounds; ///< Bounds for the bins
+    psIntArray *nums;			///< Number in each of the bins
+    const float minVal, maxVal;		///< Minimum and maximum values
+    int minNum, maxNum;			///< Number below the minimum and above the maximum
+    int uniform;			///< Is it a uniform distribution?
 } psHistogram;
 \end{verbatim}
@@ -1772,7 +1834,7 @@
 /** Constructor */
 psHistogram *
-psHistogramAlloc(float lower,           ///< Lower limit for the bins
-                 float upper,           ///< Upper limit for the bins
-                 float size             ///< Size of the bins
+psHistogramAlloc(float lower,		///< Lower limit for the bins
+		 float upper,		///< Upper limit for the bins
+		 int n			///< Number of the bins
     );
 \end{verbatim}
@@ -1781,6 +1843,5 @@
 /** Generic constructor */
 psHistogram *
-psHistogramAllocGeneric(const psFloatArray *restrict lower, ///< Lower bounds for the bins
-                        const psFloatArray *restrict upper, ///< Upper bounds for the bins
+psHistogramAllocGeneric(const psFloatArray *restrict bounds, ///< Bounds for the bins
                         float minVal,   ///< Minimum value
                         float maxVal    ///< Maximum value
@@ -1825,5 +1886,6 @@
 \item Calculate a matrix determinant;
 \item Perform matrix addition, subtraction and multiplication;
-\item Transpose a matrix; and
+\item Transpose a matrix;
+\item Get eigenvectors for a matrix; and
 \item Convert a matrix to a vector.
 \end{itemize}
@@ -1905,4 +1967,13 @@
 \end{verbatim}
 
+Eigenvectors of a matrix are calculated by \code{psMatrixEigenvectors}:
+
+\begin{verbatim}
+/** Eigenvectors of a matrix */
+psVector *
+psMatrixEigenvectors(psImage *myMatrix	///< Matrix to get eigenvectors for
+    );
+\end{verbatim}
+
 Finally, we specify two functions to convert between matrices and
 vectors.  This allows the use of vectors in matrix arithmetic, since
@@ -1930,148 +2001,87 @@
 
 We require the ability to calculate the (fast) fourier transforms of
-floating-point two-dimensional arrays (including one-dimensional
-arrays as a limiting case), and use these to perform filtering, and
-calculate the power spectrum and cross-correlation.  We expect that
-these will be implemented through wrapping an external library.  We
-define a FFT type:
-
-\begin{verbatim}
-/** Fast Fourier Transform */
-typedef struct {
-    p_psFFTDetails *details;            ///< Details on FFT implementation (private)
-    int nx, ny;                         ///< Size in x and y
-    float *real;                        ///< Data in real space: a 2D array using the [nx*y + x] stuff
-    void *fourier;                      ///< Data in fourier space; implementation dependent
-} psFFT;
-\end{verbatim}
-
-The \code{details} entry allows us to wrap the external library in the
-event that it carries information around from one transform to the
-next (as is the case for FFTW), and hence it is
-implementation-dependent.  The size of the data in real space are
-stored in \code{nx} and \code{ny}, and the data is stored as a
-Fortran-style array, \code{real}.  Since different libraries store the
-data in the Fourier space differently, \code{fourier}, we have
-specified it as of type \code{void*}.
-
-The intent is that every function that employs a FFT will act on a
-\code{psFFT}.  A user will create a \code{psFFT} by stuffing it with
-the real data, perform operations in Fourier space using the
-\code{psFFT}, and then extract the output as an image when done.
-
-We specify two constructors --- one for use with images, and the other
-for use with a one-dimensional array.  The destructor returns an image
-that includes the real space data.
-
-\begin{verbatim}
-/** Constructor */
-psFFT *
-psFFTAlloc(psImage *image               //!< Image to transform
-           );
-
-/** Constructor for 1D case */
-psFFT *
-psFFTAlloc1D(const psFloatArray *arr    //!< Array to transform
-             );
-
-/** Destructor. Returns the data in the real space as an image. */
+floating-point two-dimensional images and one-dimensional arrays.  We
+expect that these will be implemented through wrapping an external
+library.
+
+The forward and reverse Fourier transforms may be calculated using the
+following APIs:
+
+\begin{verbatim}
+/** Forward FFT an image.  Returns a complex float image. */
 psImage *
-psFFTFree(psImage *out,                 //!< Image to write the data to, or NULL
-          psFFT *restrict fft           //!< FFT to destroy
-          );
-\end{verbatim}
-
-The forward and reverse Fourier transforms may be calculated,
-returning the altered \code{psFFT}:
-
-\begin{verbatim}
-/** Forward FFT: from real to fourier space */
-psFFT *
-psFFTForward(psFFT *fft                 ///< FFT to apply (input and output)
-             );
-
-/** Reverse FFT: from fourier to real space */
-psFFT *
-psFFTReverse(psFFT *fft                 ///< FFT to apply (input and output)
-             );
-\end{verbatim}
-
-The data in Fourier space may be filtered using functions that return
-a multiplicative factor for a given position in Fourier space.  These
-function return the altered \code{psFFT}.  If the filter function
-specified for \code{psFFTFilter()} returns a real value, $r$, then the
-corresponding value in the Fourier plane should be multiplied by $r$.
-If the real and imaginary filter functions specified for
-\code{psFFTFilterComplex()} return the values $r$ and $s$,
-respectively, then the corresponding value in the Fourier plane should
-be multiplied by the complex number $r + si$.
-
-\begin{verbatim}
-/** Apply filter function in fourier space */
-psFFT *
-psFFTFilter(psFFT *fft,                 //!< FFT to use (input and output)
-            float (*filterFunc)(int kx, int ky) //!< External filter function
-            );
-
-/** Apply complex filter function */
-psFFT *
-psFFTFilterComplex(psFFT *fft,          //!< FFT to use (input and output)
-                   float (*realFilterFunc)(int kx, int ky), //!< External filter function, real part
-                   float (*imagFilterFunc)(int kx, int ky) //!< External filter function, imaginary part
-                   );
-\end{verbatim}
-
-The cross-correlation function may be calculated from two Fourier transforms,
-and returns a new \code{psFFT}:
-
-\begin{verbatim}
-/** Calculate FFT of the cross-correlation.  Multiplication of the first FFT with the complex conjugate of
- *  the second.
- */
-psFFT *
-psFFTCrossCorrelate(psFFT *out          //!< Output FFT (or NULL)
-                    psFFT *fft1, psFFT *fft2 //!< FFTs to use in cross-correlation
-                    );
-\end{verbatim}
-
-In convolution, two Fourier transforms are multiplied, and the resultant
-\code{psFFT} is returned.
-
-\begin{verbatim}
-/** Calculate FFT of the convolution.  Straight multiplication of the FFTs */
-psFFT *
-psFFTConvolve(psFFT *out,               ///< Output FFT (or NULL)
-              const psFFT *fft1, const psFFT *fft2 ///< FFTs to multiply
-              );
-\end{verbatim}
-
-The power spectrum is calculated from the Fourier transform, returning
-an array of floating-point values containing the power spectrum.
-
-\begin{verbatim}
-/** Calculate power spectrum */
+psImageFFT(const psImage *image		///< Image to transform
+    );
+
+/** Forward FFT an array of floating-point numbers. */
 psFloatArray *
-psFFTPowerSpec(psFFT *fft               ///< FFT to use (input and output)
-               );
-\end{verbatim}
-
-And finally, the user may extract both the real data and the complex
-Fourier space data as \code{psImage}s (without needing to destroy the
-\code{psFFT}).  For \code{psFFTGetFT}, note that the output image will be of
-type \code{complex float}, since it extracts the data in the Fourier plane,
-which is complex.
-
-\begin{verbatim}
-/* Convert the real data in the FFT struct to an image again */
+psFloatArrayFFT(const psFloatArray *arr	///< Array to transform
+    );
+
+/** Inverse FFT an image. */
 psImage *
-psFFTGetImage(psImage *out,             //!< Image to write to (or NULL)
-              const psFFT *fft          //!< FFT to get image from
-              );
-
-/** Convert the Fourier transform data in the FFT struct to an image of complex numbers */
+psImageInverseFFT(const psImage *fftImage ///< FFT-ed image to inverse-transform
+    );
+
+/** Inverse FFT an array of complex floating-point numbers */
+psComplexArray *
+psComplexArrayInverseFFT(const psComplexArray *fftArray ///< FFT-ed array to inverse-transform
+    );
+\end{verbatim}
+
+Neither the forward or inverse transforms shall multiply by $1/N$ (or
+$1/N^{1/2}$), and so it falls to the responsibility of the user to
+multiply an image that has been forward- and reverse-transformed by
+$1/N$.
+
+The \code{psImage}s output by \code{psImageFFT} and
+\code{psImageInverseFFT} should be of complex floating-point type.  Following
+an inverse transform, therefore, the user will likely want to extract the
+real parts:
+
+\begin{verbatim}
+/** Get the real part of an image */
 psImage *
-psFFTGetFT(psImage *out,                //!< Image to write to (or NULL)
-           const psFFT *fft             //!< FFT to get Fourier transform from
-           );
+psImageReal(psImage *out,		///< Image for output (or NULL)
+	    const psImage *in		///< Image to get the real part of
+	    );
+
+/** Get the real part of an array of complex floating-point numbers */
+psFloatArray *
+psComplexArrayReal(psFloatArray *out,	///< Array for output (or NULL)
+		   const psComplexArray *in ///< Array to get the real part of
+    );
+\end{verbatim}
+
+To aid operations on FFT-ed images or arrays, we also specify
+functions to return the complex conjugate of an image or array:
+
+\begin{verbatim}
+/** Get the complex conjugate of an image */
+psImage *
+psImageConjugate(psImage *out,		///< Image for output (or NULL)
+		 const psImage *in	///< Image to get the complex conjugate of
+    );
+
+/** Get the complex conjugate of an array of complex floating-point numbers */
+psComplexArray *
+psComplexArrayConjugate(psComplexArray *out, ///< Array for output (or NULL)
+			const psComplexArray *in ///< Array to get the complex conjugate of
+    );
+\end{verbatim}
+
+Finally, the power spectrum may be calculated directly from the
+original image or array:
+
+\begin{verbatim}
+/** Calculate power spectrum of an image */
+psImage *
+psImagePowerSpectrum(const psImage *image ///< Image to obtain power spectrum of
+    );
+
+/** Calculate power spectrum of an array of floating-point numbers */
+psFloatArray *
+psFloatArrayPowerSpectrum(const psFloatArray *arr ///< Array to obtain power spectrum of
+    );
 \end{verbatim}
 
@@ -2085,5 +2095,7 @@
 \begin{verbatim}
 /** Evaluate a non-normalized Gaussian with the given mean and sigma at the given coordianate.  Note that this
-    is not a Gaussian deviate.  The evaluated Gaussian is: \f[ exp(-\frac{(x-mean)^2}{2\sigma^2}) \f] */
+ *  is not a Gaussian deviate.  The evaluated Gaussian is:
+ *  \f[ 1/(\sqrt(2\pi)\sigma) exp(-\frac{(x-mean)^2}{2\sigma^2}) \f]
+ */
 float
 psGaussian(float x,                     ///< Value at which to evaluate
@@ -2180,10 +2192,11 @@
 
 \begin{verbatim}
-/** Minimize a particular non-linear function */
+/** Find the minimum of a particular non-linear function */
 psFloatArray *
-psMinimize(float (*myFunction)(const psFloatArray *restrict), ///< Function to minimize
-           psFloatArray *restrict initialGuess, ///< Initial guess
-           psIntArray *restrict paramMask //!< 1 = fit for parameter, 0 = hold parameter constant
-           );
+psMinimize(psFloatArray *restrict initialGuess,	///< Initial guess and answer
+	   float (*myFunction)(const psFloatArray *restrict), ///< Function to minimize
+	   float (*myFuncDeriv)(const psFloatArray *restrict), ///< Derivatives of function, or NULL
+	   const psIntArray *restrict paramMask //!< 1 = fit for parameter, 0 = hold parameter constant
+	   );
 \end{verbatim}
 
@@ -2191,17 +2204,18 @@
 the specified function.  It takes as input a function,
 \code{myFunction}, an initial guess for the parameters that minimize
-the function, \code{initialGuess}, and an optional mask for the
-parameters to minimize, \code{paramMask}.
+the function, \code{initialGuess} which is returned with the best-fit
+parameters, and an optional mask for the parameters to minimize,
+\code{paramMask} (all parameters are fit if \code{NULL}).
 
 \begin{verbatim}
 /** Minimize chi^2 for input data */
 psFloatArray *
-psMinimizeChi2(float (*evalModel)(const psFloatArray *restrict,
-                                  const psFloatArray *restrict), ///< Model to fit; (domain and params)
-               const psFloatArray *restrict domain, ///< The domain values for the corresponding measurements
-               const psFloatArray *restrict data, ///< Data to fit
-               const psFloatArray *restrict errors, ///< Errors in the data
-               psFloatArray *restrict initialGuess, ///< Initial guess
-               const psIntArray *restrict paramMask ///< 1 = fit for parameter, 0 = hold parameter constant
+psMinimizeChi2(psFloatArray *restrict initialGuess, ///< Initial guess and answer
+	       float (*evalModel)(const psFloatArray *restrict,
+				  const psFloatArray *restrict), ///< Model to fit; (domain and params)
+	       const psFloatArray *restrict domain, ///< The domain values for the corresponding measurements
+	       const psFloatArray *restrict data, ///< Data to fit
+	       const psFloatArray *restrict errors, ///< Errors in the data
+	       const psIntArray *restrict paramMask ///< 1 = fit for parameter, 0 = hold parameter constant
     );
 \end{verbatim}
@@ -2212,6 +2226,7 @@
 parameters, \code{evalModel}; a list of observations, (\code{domain},
 \code{data}, and \code{errors}); an initial guess at the best-fit
-parameters, \code{initialGuess}, and an optional mask specifying which
-parameters are to be fit, \code{paramMask}.
+parameters, \code{initialGuess} which is returned with the best-fit
+parameters, and an optional mask specifying which parameters are to be
+fit, \code{paramMask} (all parameters are fit if \code{NULL}).
 
 \begin{verbatim}
@@ -2243,22 +2258,22 @@
 \begin{verbatim}
 typedef struct psImage {
-    psType type;                        ///< image data type and dimension
-    int nx, ny;                         ///< size of image 
-    int x0, y0;                         ///< data region relative to parent 
+    psType type; 			///< image data type and dimension
+    const int nx, ny;			///< size of image 
+    const int x0, y0;			///< data region relative to parent 
     union {
-        float **rows;                   ///< Pointers to floating-point data (default)
-        short int **rows_si;            ///< Pointers to short-integer data
-        int **rows_i;                   ///< Pointers to integer data
-        long int **rows_li;             ///< Pointers to long-integer data
-        unsigned short int **rows_usi;  ///< Pointers to unsigned-short-integer data
-        unsigned int **rows_ui;         ///< Pointers to unsigned-integer data
-        unsigned long int **rows_uli;   ///< Pointers to unsigned-long-integer data
-        float **rows_f;                 ///< Pointers to floating-point data
-        double **rows_d;                ///< Pointers to double-precision data
-        complex float **rows_complex;   ///< Pointers to complex floating-point data
+	float **rows;			///< Pointers to floating-point data (default)
+	short int **rows_si;		///< Pointers to short-integer data
+	int **rows_i;			///< Pointers to integer data
+	long int **rows_li;		///< Pointers to long-integer data
+	unsigned short int **rows_usi;	///< Pointers to unsigned-short-integer data
+	unsigned int **rows_ui;		///< Pointers to unsigned-integer data
+	unsigned long int **rows_uli;	///< Pointers to unsigned-long-integer data
+	float **rows_f;			///< Pointers to floating-point data
+	double **rows_d;		///< Pointers to double-precision data
+	complex float **rows_complex;	///< Pointers to complex floating-point data
     } rows;
-    struct psImage *parent;             ///< parent, if a subimage 
-    int Nchildren;                      ///< number of subimages 
-    struct psImage *children;           ///< children of this region; array of Nchildren pointers
+    const struct psImage *parent;	///< parent, if a subimage 
+    int Nchildren;			///< number of subimages 
+    struct psImage *children;		///< children of this region; array of Nchildren pointers
 } psImage;
 \end{verbatim}
@@ -2266,17 +2281,17 @@
 This structure represents an image consisting of a 2-D array of
 pixels.  The size of this array is given by the elements \code{(nx,
-ny)}.  The data type of the pixel is defined by the \code{psType type}
-entry (see \ref{arithmetic}).  (n.b. that for FITS images, these
-values are restricted to the datatypes equivalent to the valid BITPIX
-values 8, 16, 32, -32, -64).  The image represented in the data
-structure may represent a subset of the pixels in a complete array, in
-which case the image is considered to be the child of that parent
-array.  The offset of the \code{(0,0)} pixel in this array relative to
-the parent array is given by the elements \code{(x0,y0)}.  The
-structure may include references to subrasters (\code{children,
-Nchildren}) and/or to a containing array (\code{parent}).  Unless this
-is image is a child of another image (represents a subset of the
-pixels of another image), the image data is allocated in a contiguous
-block.
+ny)}.  The data type of the pixel is defined in the \code{psType type}
+entry (specifically, the \code{psElemType} member, \code{type}; see
+\ref{sec:arithmetic}).  (n.b. that for FITS images, these values are
+restricted to the datatypes equivalent to the valid BITPIX values 8,
+16, 32, -32, -64).  The image represented in the data structure may
+represent a subset of the pixels in a complete array, in which case
+the image is considered to be the child of that parent array.  The
+offset of the \code{(0,0)} pixel in this array relative to the parent
+array is given by the elements \code{(x0,y0)}.  The structure may
+include references to subrasters (\code{children, Nchildren}) and/or
+to a containing array (\code{parent}).  Unless this is image is a
+child of another image (represents a subset of the pixels of another
+image), the image data is allocated in a contiguous block.
 
 We require a variety of functions to manipulate these image
@@ -2291,5 +2306,5 @@
 to the valid FITS BITPIX types.
 \begin{verbatim}
-psImage *psImageAlloc (int nx, int ny, psType type);
+psImage *psImageAlloc (int nx, int ny, psElemType type);
 \end{verbatim}
 where \code{nx} and \code{ny} specify the size of the image and
@@ -2301,12 +2316,11 @@
 outside of the parent image.
 \begin{verbatim}
-psImage *psImageSubset(psImage *out, psImage *image, int nx, int ny, int x0, int y0);
+psImage *psImageSubset(const psImage *image, int nx, int ny, int x0, int y0);
 \end{verbatim}
 where \code{image} is the parent image, \code{nx,ny} specify the
 dimensions of the desired subraster, and \code{x0, y0} specify the
 starting pixel of the subraster.  The entire subraster must be
-contained within the raster of the parent image.  An image structure
-to populate may be specified by \code{out}, or if this is \code{NULL},
-a new structure is allocated.  
+contained within the raster of the parent image.  Note that the
+\code{refCounter} for the parent should be incremented.
 
 Free the memory associated with a specific image, including the pixel
@@ -2321,16 +2335,11 @@
 \end{verbatim}
 
-Free the image structure memory associated with the children of a
-specific image, returning the number of children freed.  Do not free
-the pixel data:
-\begin{verbatim}
-int psImageFreeChildren(psImage *image);
-\end{verbatim}
-
-Create a copy of the specified image.  If the output target pointer is
-not NULL, place the result in the specified structure.  The output
-image data must be allocated as a single, contiguous block of memory.
-\begin{verbatim}
-psImage *psImageCopy(psImage *output, psImage *input);
+Create a copy of the specified image, converting the type in the
+process.  If the output target pointer is not NULL, place the result
+in the specified structure.  The output image data must be allocated
+as a single, contiguous block of memory.  The output image may not be
+the input image.
+\begin{verbatim}
+psImage *psImageCopy(psImage *output, const psImage *input, psElemType type);
 \end{verbatim}
 
@@ -2346,10 +2355,10 @@
 specified by \code{stats}.
 \begin{verbatim}
-psFloatArray *psImageSlice(psFloatArray *out, psImage *input, int x, int y, int nx, int ny, 
+psFloatArray *psImageSlice(psFloatArray *out, const psImage *input, int x, int y, int nx, int ny, 
                            int direction, const psStats *stats);
 \end{verbatim}
 
 Extract pixels from an image along a line to a vector (array of
-floats).  The vector \code{xs,ys} - \code{xe,ye} forms the basis of
+floats).  The vector \code{(xs,ys)} - \code{(xe,ye)} forms the basis of
 the output vector.  Pixels are considered in a rectangular region of
 width \code{dw} about this vector.  The input region is collapsed in
@@ -2360,6 +2369,6 @@
 specified by \code{stats}.
 \begin{verbatim}
-psFloatArray *psImageCut(psFloatArray *out, psImage *input, float xs, float ys, float xe, float ye, 
-                         float dw, psStats *stats);
+psFloatArray *psImageCut(psFloatArray *out, const psImage *input, float xs, float ys, float xe, float ye, 
+                         float dw, const psStats *stats);
 \end{verbatim}
 
@@ -2371,6 +2380,6 @@
 output vector value is specified by \code{stats}
 \begin{verbatim}
-psFloatArray *psImageRadialCut(psFloatArray *out, psImage *input, float x, float y, float radius, float dr, 
-                               psStats *stats);
+psFloatArray *psImageRadialCut(psFloatArray *out, const psImage *input, float x, float y,
+                               const psFloatArray *radii, const psStats *stats);
 \end{verbatim}
 
@@ -2379,10 +2388,11 @@
 Rebin image to new scale.  A new image is constructed in which the
 dimensions are reduced by a factor of \code{scale} $\le 1$ (it is an
-error for \code{scale} $> 1$).  The output image is generated from all
-input image pixels.  Each pixel in the output image is derived from
-the statistics of the corresponding set of input image pixels based on
-the statistics specified by \code{stats}.
-\begin{verbatim}
-psImage *psImageRebin(psImage *out, psImage *input, float scale, psStats *stats);
+error for \code{scale} $> 1$).  The \code{scale} is equal in each
+dimension.  The output image is generated from all input image pixels.
+Each pixel in the output image is derived from the statistics of the
+corresponding set of input image pixels based on the statistics
+specified by \code{stats}.
+\begin{verbatim}
+psImage *psImageRebin(psImage *out, const psImage *input, float scale, const psStats *stats);
 \end{verbatim}
 
@@ -2394,5 +2404,5 @@
 that a positive angle is an anti-clockwise rotation.  
 \begin{verbatim}
-psImage *psImageRotate(psImage *out, psImage *input, float angle);
+psImage *psImageRotate(psImage *out, const psImage *input, float angle);
 \end{verbatim}
 
@@ -2404,5 +2414,5 @@
 set to the value given by \code{exposed}.  
 \begin{verbatim}
-psImage *psImageShift(psImage *out, psImage *input, float dx, float dy, float exposed);
+psImage *psImageShift(psImage *out, const psImage *input, float dx, float dy, float exposed);
 \end{verbatim}
 
@@ -2411,5 +2421,5 @@
 image.  Edge pixels wrap to the other side (no values are lost).
 \begin{verbatim}
-psImage *psImageRoll(psImage *out, psImage *input, int dx, int dy);
+psImage *psImageRoll(psImage *out, const psImage *input, int dx, int dy);
 \end{verbatim}
 
@@ -2419,12 +2429,12 @@
 determined are specified by \code{stats}.
 \begin{verbatim}
-psStats *psImageGetStats(psImage *input, psStats *stats);
+psStats *psImageGetStats(psStats *stats, const psImage *input);
 \end{verbatim}
 
 Construct a histogram from an image (or subimage).  The histogram to
 generate is specified by \code{psHistogram hist} (see
-section~\ref{histogram}).
-\begin{verbatim}
-psHistogram *psImageHistogram(psHistogram *hist, psImage *input);
+section~\ref{sec:histograms}).
+\begin{verbatim}
+psHistogram *psImageHistogram(psHistogram *hist, const psImage *input);
 \end{verbatim}
 
@@ -2433,12 +2443,12 @@
 interest.
 \begin{verbatim}
-psPolynomial2D *psImageFitPolynomial(psImage *input, psPolynomial2D *coeffs);
+psPolynomial2D *psImageFitPolynomial(psPolynomial2D *coeffs, const psImage *input);
 \end{verbatim}
 
 Evaluate a 2-D polynomial surface for the image pixels.  Given the
 input polynomial coefficients, set the image pixel values on the basis
-of the polynomial function.  \tbd{Return value is not yet specified.}
-\begin{verbatim}
-int psImageEvalPolynomial(psImage *input, psPolynomial2D *coeffs);
+of the polynomial function.
+\begin{verbatim}
+psImage *psImageEvalPolynomial(psImage *input, const psPolynomial2D *coeffs);
 \end{verbatim}
 
@@ -2449,17 +2459,17 @@
 full image or a subimage to be read.  The starting pixel of the region
 is specified by \code{x,y}, while the dimensions of the requested
-region are specified by \code{nx,ny}.  A value of -1 for these two
-parameters specifies the full array of the requested image.  If the
-native image is a cube, the value of z specifies the requested slice
-of the image.  The data is read from the extension specified by
-extname (matching the EXTNAME keyword) or by the extnum value (with 0
-representing the PHU, 1 the first extension, etc).  This function must
-generate an error and return \code{NULL} if any of the specified
-parameters are out of range for the data in the image file, if the
-specified image file does not exist, or the image on disk is zero- or
-one-dimensional.
-\begin{verbatim}
-psImage *psImageReadSection(psImage *output, int x, int y, int nx, int ny, int z, 
-                            char *extname, int extnum, char *filename);
+region are specified by \code{nx,ny}.  A negative value of -1 for
+these two parameters specifies the full array of the requested image,
+less absolute value of the variable.  If the native image is a cube,
+the value of z specifies the requested slice of the image.  The data
+is read from the extension specified by extname (matching the EXTNAME
+keyword) or by the extnum value (with 0 representing the PHU, 1 the
+first extension, etc).  This function must call \code{psError} and
+return \code{NULL} if any of the specified parameters are out of range
+for the data in the image file, if the specified image file does not
+exist, or the image on disk is zero- or one-dimensional.
+\begin{verbatim}
+psImage *psImageReadSection(psImage *output, int x0, int y0, int nx, int ny, int z, 
+                            const char *extname, int extnum, const char *filename);
 \end{verbatim}
  
@@ -2469,5 +2479,5 @@
 \begin{verbatim}
 psImage *psImageFReadSection(psImage *output, int x, int y, int nx, int ny, int z, 
-                             char *extname, int extnum, FILE *f);
+                             const char *extname, int extnum, FILE *f);
 \end{verbatim}
 \tbd{The use of \code{FILE*} to carry around the file descriptor is to be reviewed.}
@@ -2480,22 +2490,22 @@
 created.
 \begin{verbatim}
-psImage *psImageWriteSection(psImage *input, int x, int y, int z, 
-                             char *extname, int extnum, char *filename);
+psImage *psImageWriteSection(const psImage *input, int x, int y, int z, 
+                             const char *extname, int extnum, const char *filename);
 \end{verbatim}
 
 Write an image section to file descriptor as above:
 \begin{verbatim}
-psImage *psImageFWriteSection(psImage *input, int x, int y, int z, 
-                              char *extname, int extnum, FILE *f);
+psImage *psImageFWriteSection(const psImage *input, int x, int y, int z, 
+                              const char *extname, int extnum, FILE *f);
 \end{verbatim}
 \tbd{The use of \code{FILE*} to carry around the file descriptor is to be reviewed.}
 
 Read header data from a FITS image file into a \code{psMetaData}
-structure (see section~\ref{metadata}.  The \code{extname} and
+structure (see section~\ref{sec:metadata}.  The \code{extname} and
 \code{extnum} parameters specify the extension of interest as above.
 If the named extension does not exist, the function should return an
 error.
 \begin{verbatim}
-psMetadata *psImageReadHeader(psMetadata *output, char *extname, int extnum, char *filename);
+psMetadata *psImageReadHeader(psMetadata *output, const char *extname, int extnum, const char *filename);
 \end{verbatim}
 
@@ -2503,5 +2513,5 @@
 structure.  
 \begin{verbatim}
-psMetadata *psImageFReadHeader(psMetadata *output, char *extname, int extnum, FILE *f);
+psMetadata *psImageFReadHeader(psMetadata *output, const char *extname, int extnum, FILE *f);
 \end{verbatim}
 
@@ -2531,5 +2541,5 @@
 (multiply overlay times image), \code{/} (divide image by overlay).
 \begin{verbatim}
-int psImageOverlaySection(psImage *image, psImage *overlay, int x0, int y0, char *op);
+int psImageOverlaySection(psImage *image, const psImage *overlay, int x0, int y0, const char *op);
 \end{verbatim}
 
@@ -2928,6 +2938,6 @@
 appropriate destructor for the \code{val} (recall that it is the duty
 of the \code{psMyTypeFree}s to decrement the \code{refCounter} and
-free the memory if and only if the new value of the \code{refCounter}
-is zero --- see \S\ref{sec:free}).
+free the memory if and only if the \code{refCounter == 1} --- see
+\S\ref{sec:free}).
 
 \begin{verbatim}
@@ -3010,5 +3020,6 @@
 key, e.g., if the user only wants to iterate through
 \code{IPP.machines.sky} and doesn't want to bother with
-\code{IPP.machines.detector}.
+\code{IPP.machines.detector}.  The iterator should iterate over every
+item of metadata --- even those that are non-unique.
 
 \begin{verbatim}
@@ -3017,7 +3028,8 @@
     );
 
-/// get the next entry in the sequence
+/// get the next item in the sequence
 psMetadataItem *psMetadataGetNext(psMetadata *restrict md, ///< metadata to get from
-                                  const char *restrict match ///< Match this
+				  const char *restrict match, ///< Match this
+				  int which ///< Which iterator to use
     );
 \end{verbatim}
@@ -3046,34 +3058,27 @@
 \subsection{Detector and sky positions}
 
-\tbd{the definition of psCoord must be review: split into psSphere and psPlane?}
-
 Both detector and sky positions will be used extensively in the IPP.
 The first are linear coordinates which conform to Euclidean geometry
 while the second are angular coordinates for which additional care
-must often be taken.  Since these both consist of two coordinates with
-their associated errors, we bundle these into a single generic
-structure, \code{psCoord}, containing \code{union}s to handle the
-semantic differences.  A variety of functions should operate
-equivalently on both types of structures (\tbd{list}), while other
-functions must behave differently when the operand is an angular vs a
-linear coordinate system.
-
-\begin{verbatim}
-/** A point in 2-D space, with errors.
- */
-typedef union {
-    struct {
-        double x;                       ///< x position
-        double y;                       ///< y position
-        double xErr;                    ///< Error in x position
-        double yErr;                    ///< Error in y position
-    } xy;
-    struct {
-        double r;                       ///< RA
-        double d;                       ///< Dec
-        double rErr;                    ///< Error in RA
-        double dErr;                    ///< Error in Dec
-    } rd;
-} psCoord;
+must often be taken.  We put these into two structures,
+\code{psPlaneCoord} and \code{psSphereCoord}, respectively.
+Partitioning these two will enable error-checking.
+
+\begin{verbatim}
+/** A point in 2-D space, with errors. */
+typedef struct {
+    double x;				///< x position
+    double y;				///< y position
+    double xErr;			///< Error in x position
+    double yErr;			///< Error in y position
+} psPlaneCoord;
+
+/** A point on the surface of a sphere, with errors */
+typedef struct {
+    double r;				///< RA
+    double d;				///< Dec
+    double rErr;			///< Error in RA
+    double dErr;			///< Error in Dec
+} psSphereCoord;
 \end{verbatim}
 
@@ -3120,5 +3125,5 @@
     psDPolynomial2D *x;
     psDPolynomial2D *y;
-} psCoordXform;
+} psPlaneCoordXform;
 \end{verbatim}
 
@@ -3129,5 +3134,5 @@
 frame \code{x,y} (say a single CCD) to a second frame \code{p,q} (say,
 a second CCD image). If we have only first order terms in the
-transformation \code{psCoordXform T}, the new coordinates would be
+transformation \code{psPlaneCoordXform T}, the new coordinates would be
 represented by the terms:
 %
@@ -3147,8 +3152,8 @@
     psDPolynomial4D *x;
     psDPolynomial4D *y;
-} psDistortion;
-\end{verbatim}
-
-Like \code{psCoordXform}, \code{psDistortion} contains two
+} psPlaneDistortion;
+\end{verbatim}
+
+Like \code{psPlaneCoordXform}, \code{psPlaneDistortion} contains two
 \code{psDPolynomial4D} structures representing polynomials of
 arbitrary order as a function of four, rather than two dimensions.
@@ -3159,5 +3164,5 @@
 single CCD) of an object with magnitude and color \code{m,c} to a
 second frame \code{p,q} (say, a second CCD image). If we have only
-first order terms in the transformation \code{psCoordXform T}, the new
+first order terms in the transformation \code{psPlaneCoordXform T}, the new
 coordinates would be represented by the terms:
 %
@@ -3176,17 +3181,15 @@
 \begin{verbatim}
 /** apply the coordinate transformation to the given coordinate */
-psCoord *psCoordXformApply (psCoord *out, ///< Output coordinates, or NULL
-                            const psCoordXform *frame, ///< coordinate transformation
-                            const psCoord *coords ///< input coordiate
-    );
-\end{verbatim}
-%
-\begin{verbatim}
+psPlaneCoord *psPlaneCoordXformApply (psPlaneCoord *out, ///< Output coordinates, or NULL
+				      const psPlaneCoordXform *frame, ///< coordinate transformation
+				      const psPlaneCoord *coords ///< input coordiate
+    );
+
 /** apply the optical distortion to the given coordinate, magnitude, color */
-psCoord *psDistortionApply (psCoord *out, ///< Output coordinates, or NULL
-                            const psdistortion *pattern, ///< optical distortion pattern
-                            const psCoord *coords, ///< input coordinate
-                            float mag,  ///< magnitude of object
-                            float color ///< color of object
+psPlaneCoord *psPlaneDistortionApply (psPlaneCoord *out, ///< Output coordinates, or NULL
+				      const psPlaneDistortion *pattern, ///< optical distortion pattern
+				      const psPlaneCoord *coords, ///< input coordinate
+				      float mag, ///< magnitude of object
+				      float color ///< color of object
     );
 \end{verbatim}
@@ -3199,53 +3202,64 @@
 systems.  All of these basic spherical transformations represent
 rotations of the spherical coordinate reference.  We specify a general
-transformation function which takes a structure \code{psCoordSphere}
-defining the relationship between two spherical coordinate systems.
-We also define a function to generate \code{psCoordSphere}, which
-contains the appropriate sines and cosines of the angles, based on the
-three angles describing the location of the pole and the relative
-equatorial rotations of the two systems.  We also specify special
-coordinate transformations between standard coordinate systems which
-define the value of psCoordSphere internally and automatically.
-
-\begin{verbatim}
-/** Convert to new spherical coordinate system */
-psCoord *
-psCoordTransform (psCoord *coord, 
-                  psCoordSphere *sphere) ///< ICRS coordinates to convert
-;
-\end{verbatim}
-
-\begin{verbatim}
-/** Construct  to new spherical coordinate system */
-psCoordSphere *
-psCoordSphereDefine (double pole1, double pole2, double rotation);
-\end{verbatim}
-
-\begin{verbatim}
-/** Convert ICRS to Ecliptic */
-psCoord *
-psCoordinatesItoE(const psCoord *restrict coordinates ///< ICRS coordinates to convert
-    );
-\end{verbatim}
-
-\begin{verbatim}
-/** Convert Ecliptic to ICRS */
-psCoord *
-psCoordinatesEtoI(const psCoord *restrict coordinates ///< Ecliptic coordinates to convert
-    );
-\end{verbatim}
-
-\begin{verbatim}
-/** Convert ICRS to Galactic */
-psCoord *
-psCoordinatesItoG(const psCoord *restrict coordinates ///< ICRS coordinates to convert
-    );
-\end{verbatim}
-
-\begin{verbatim}
-/** Convert Galactic to ICRS */
-psCoord *
-psCoordinatesGtoI(const psCoord *restrict coordinates ///< Galactic coordinates to convert
-    );
+transformation function which takes a structure,
+\code{psSphereCoordTransformation}, defining the transformation
+between two spherical coordinate systems (the structure contains the
+sines and cosines of the angles involved so as to minimize computation
+time for repeated transformations).  We also define a function to
+generate \code{psSphereCoordTransformation}, based on the three angles
+describing the location of the pole and the relative equatorial
+rotations of the two systems.  We also specify special functions to
+return the \code{psSphereCoordTransformation} for transformations
+between standard coordinate systems.
+
+\begin{verbatim}
+/** General spherical transformation */
+typedef struct {
+    double sin1, sin2, sin3, cos1, cos2, cos3; ///< Sines and cosines for transformation
+} psSphereCoordTransformation;
+\end{verbatim}
+
+The constructor and destructor are defined as follows:
+
+\begin{verbatim}
+/** Constructor */
+psSphereCoordTransformation *
+psSphereCoordTransformationAlloc(double pole1, ///< First location of pole
+				 double pole2, ///< Second location of pole
+				 double rotation ///< Rotation between systems
+    );
+
+/** Destructor */
+void psSphereCoordTranformationFree(psSphereCoordTransformation *trans ///< Transformation to destroy
+    );
+\end{verbatim}
+
+Spherical coordinates may be transformed by providing the transformation to
+\code{psSphereCoordTransform}:
+
+\begin{verbatim}
+/** Apply general spherical transformation */
+psSphereCoord *
+psSphereCoordTransform(const psSphereCoord *coord, ///< Coordinates to convert
+		       psSphereCoordSystem *sys	///< System to use to convert
+    );
+\end{verbatim}
+
+The following functions simply return the appropriate
+\code{psSphereCoordTransformation} to convert between predefined
+spherical coordinate systems (i.e., ICRS, Ecliptic and Galactic).
+
+\begin{verbatim}
+/** Return transformation structure to convert ICRS to Ecliptic */
+psSphereCoordTransformation *psSphereCoordTransformationItoE(void);
+
+/** Return transformation structure to convert Ecliptic to ICRS */
+psSphereCoordTransformation *psSphereCoordTransformationEtoI(void);
+
+/** Return transformation structure to convert ICRS to Galactic */
+psSphereCoordTransformation *psSphereCoordTransformationItoG(void);
+
+/** Return transformation structure to convert Galactic to ICRS */
+psSphereCoordTransformation *psSphereCoordTransformationGtoI(void);
 \end{verbatim}
 
@@ -3264,6 +3278,19 @@
 \end{itemize}
 
-\begin{verbatim}
-psCoord *psCoordProject (psCoord *coord, char *projection);
+The following functions will project and deproject (respectively)
+spherical coordinates:
+
+\begin{verbatim}
+/** Project spherical system onto a plane */
+psPlaneCoord *
+psCoordProject(const psSphereCoord *coord, ///< Spherical coordinates to project
+	       const char *projection	///< Projection to use
+    );
+
+/** Deproject plane onto spherical system */
+psSphereCoord *
+psCoordDeproject(const psPlaneCoord *coord, ///< Plane coordinates to deproject
+		 const char *projection	///< Projection to use
+    );
 \end{verbatim}
 
@@ -3275,19 +3302,21 @@
 \begin{verbatim}
 /** Get offset (RA,Dec) on the sky between two positions position1 and position2 may not be identical */
-psCoord *
-psGetOffset(const psCoord *restrict position1, ///< Position 1
-            const psCoord *restrict position2, ///< Position 2
-            const char *type            ///< Type of offset: Linear, Spherical/Arcsec, Spherical/Degreees etc
-    );
-\end{verbatim}
-
-\begin{verbatim}
+psSphereCoord *
+psSphereCoordGetOffset(const psSphereCoord *restrict position1, ///< Position 1
+		       const psSphereCoord *restrict position2, ///< Position 2
+		       const char *type		///< Type of offset: Linear, Spherical/Arcsec,
+						///< Spherical/Degreees etc
+    );
+
 /** Apply an offset to a position */
-psCoord *
-psApplyOffset(const psCoord *restrict position, ///< Position
-              const psCoord *restrict offset, ///< Offset
-              const char *type          ///< Type of offset: Linear, Spherical/Arcsec, Spherical/Degreees etc
-    );
-\end{verbatim}
+psSphereCoord *
+psSphereCoordApplyOffset(const psSphereCoord *restrict position, ///< Position
+			 const psSphereCoord *restrict offset, ///< Offset
+			 const char *type		///< Type of offset: Linear, Spherical/Arcsec,
+							///< Spherical/Degreees etc
+    );
+\end{verbatim}
+
+Note that these should propagate the errors appropriately.
 
 
@@ -3318,7 +3347,6 @@
 single readout of the detector, along with metadata needed to define
 that readout: the origin and binning of the image relative to the
-original detector pixels explicitly in the structure, a pointer to the
-appropriate overscan region, and pointers to the general metadata and
-derived objects, if any.
+original detector pixels explicitly in the structure, and pointers to
+the general metadata and derived objects, if any.
 
 A single detector may produce more than one read which is associated.
@@ -3382,18 +3410,16 @@
 
 A readout is the result of a single read of a cell (or a portion
-thereof).  It contains a pointer to the pixel data, a separate pointer
-to the overscan pixels, and additional pointers to the objects found
-in the readout, and the readout metadata.  It also contains the offset
-from the lower-left corner of the chip, in the case that the CCD was
-windowed.
+thereof).  It contains a pointer to the pixel data, and additional
+pointers to the objects found in the readout, and the readout
+metadata.  It also contains the offset from the lower-left corner of
+the chip, in the case that the CCD was windowed.
 
 \begin{verbatim}
 /** a Readout: a collection of pixels */
 typedef struct {
-    int x0, y0;                         ///< Offset from the lower-left corner
-    int nx, ny;                         ///< Image binning
-    psImage *image;                     ///< imaging area of cell 
-    psDlist *objects;                   ///< objects derived from cell
-    psImage *overscan;                  ///< bias region (subimage) of cell
+    const int x0, y0;                ///< Offset from the lower-left corner
+    const int nx, ny;                ///< Image binning
+    psImage *image;                  ///< imaging area of cell 
+    psDlist *objects;                ///< objects derived from cell
     psMetadata *md;                  ///< Readout-level metadata
 } psReadout;
@@ -3421,14 +3447,14 @@
  */
 typedef struct {
-    int nReadouts;                      ///< number of readouts in this cell realization; each may have its own
-                                        ///< image, objects and overscan.
-    struct psReadout *readouts;         ///< Readouts from the cell
-    psMetadata *md;                  ///< Cell-level metadata
-
-    psCoordXform *cellToChip;           ///< Transformations from cell coordinates to chip coordinates
-    psCoordXform *cellToFPA;            ///< Transformations from cell coordinates to FPA coordinates
-    psCoordXform *cellToSky;            ///< Quick and Dirty transformations from cell coordinates to sky
-
-    struct psChip  *parentChip;         ///< chip which contains this cell
+    int nReadouts;			///< number of readouts in this cell realization; each may have its
+					///< own image, objects and overscan.
+    struct psReadout *readouts;		///< Readouts from the cell
+    psMetadata *md;			///< Cell-level metadata
+
+    psPlaneCoordXform *cellToChip;	///< Transformations from cell coordinates to chip coordinates
+    psPlaneCoordXform *cellToFPA;	///< Transformations from cell coordinates to FPA coordinates
+    psPlaneCoordXform *cellToSky;	///< Quick and Dirty transformations from cell coordinates to sky
+
+    struct psChip  *parentChip;		///< chip which contains this cell
 } psCell;
 \end{verbatim}
@@ -3442,5 +3468,5 @@
 a coordinate transform from the chip to the focal plane.  It is
 expected that this transforms will consist of two second-order 2D
-polynomials; hence we expect that it is prudent to include a reverse
+polynomials; hence we think that it is prudent to include a reverse
 transformation which will be derived from numerically inverting the
 forward transformation.
@@ -3452,11 +3478,11 @@
 typedef struct {
     int nCells;                         ///< Number of Cells assigned
-    psCell *cells;                      ///< Cells in the Chip
-
-    psMetadata *md;                  ///< Chip-level metadata
-    psCoordXform *chipToFPA;            ///< Transformations from chip coordinates to FPA coordinates
-    psCoordXform *FPAtoChip;            ///< Transformations from FPA coordinates to chip
-
-    struct psFPA *parentFPA;            ///< FPA which contains this chip
+    struct psCell *cells;		///< Cells in the Chip
+
+    psMetadata *md;			///< Chip-level metadata
+    psPlaneCoordXform *chipToFPA;	///< Transformations from chip coordinates to FPA coordinates
+    psPlaneCoordXform *FPAtoChip;	///< Transformations from FPA coordinates to chip
+
+    struct psFPA *parentFPA;		///< FPA which contains this chip
 } psChip;
 \end{verbatim}
@@ -3472,5 +3498,5 @@
 two coordinates in position, the magnitude of the object, and the
 color of the object) in order to correct for optical distortions and
-the effects of the atmosphere; hence we expect that it is prudent to
+the effects of the atmosphere; hence we think that it is prudent to
 include a reverse transformation which will be derived from
 numerically inverting the forward transformation.  Since colors are
@@ -3488,11 +3514,11 @@
     int nChips;                         ///< Number of Cells assigned
     int nAlloc;                         ///< Number of Cells available
-    psChip *chips;                      ///< Chips in the Focal Plane Array
-
-    psMetadata *md;                     ///< FPA-level metadata 
-    psDistortion *TPtoFP;               ///< Transformation term from 
-    psDistortion *FPtoTP;               ///< Transformation term from 
-    psFixedPattern *pattern;            ///< Fixed pattern residual offsets
-    psExposure *exp;                    ///< information about this exposure
+    struct psChip *chips;		///< Chips in the Focal Plane Array
+
+    psMetadata *md;			///< FPA-level metadata 
+    psPlaneDistortion *TPtoFP;		///< Transformation term from 
+    psPlaneDistortion *FPtoTP;		///< Transformation term from 
+    psFixedPattern *pattern;		///< Fixed pattern residual offsets
+    const psExposure *exp;		///< information about this exposure
     psPhotSystem colorPlus, colorMinus; ///< Colour reference
     float rmsX, rmsY;                   ///< Dispersion in astrometric solution
@@ -3511,23 +3537,22 @@
 typedef struct {
     // Telescope longitude, latitude and height are stored separately, since they don't change with pointing
-    double ra, dec;                     ///< Telescope boresight
-    double ha;                          ///< Hour angle
-    double zd;                          ///< Zenith distance
-    double az;                          ///< Azimuth
-    double lst;                         ///< Local Sidereal Time
-    float mjd;                          ///< MJD of observation
-    float rotAngle;                     ///< Rotator position angle
-    float temp;                         ///< Air temperature, for estimating refraction
-    float pressure;                     ///< Air pressure, for calculating refraction
-    float humidity;                     ///< Relative humidity, for calculating refraction
-    float exptime;                      ///< Exposure time
+    const double ra, dec;		///< Telescope boresight
+    const double ha;			///< Hour angle
+    const double zd;			///< Zenith distance
+    const double az;			///< Azimuth
+    const double lst;			///< Local Sidereal Time
+    const float mjd;			///< MJD of observation
+    const float rotAngle;		///< Rotator position angle
+    const float temp;			///< Air temperature, for estimating refraction
+    const float pressure;		///< Air pressure, for calculating refraction
+    const float humidity;		///< Relative humidity, for calculating refraction
+    const float exptime;		///< Exposure time
     /* Derived quantities */
-    float posAngle;                     ///< Position angle
-    float parallactic;                  ///< Parallactic angle
-    float airmass;                      ///< Airmass, calculated from zenith distance
-    float pf;                           ///< Parallactic factor
-    char *cameraName;                   ///< name of camera which provided exposure
-    char *telescopeName;                ///< name of telescope which provided exposure
-    psGrommit *grommit;                 ///< Data needed to convert from the sky to the tangent plane
+    const float posAngle;		///< Position angle
+    const float parallactic;		///< Parallactic angle
+    const float airmass;		///< Airmass, calculated from zenith distance
+    const float pf;			///< Parallactic factor
+    const char *cameraName;		///< name of camera which provided exposure
+    const char *telescopeName;		///< name of telescope which provided exposure
 } psExposure;
 \end{verbatim}
@@ -3583,5 +3608,5 @@
 
 The following steps are required to convert from the cell coordinates to
-the sky \tbd{convert Cell->Chip and Chip->FP to psDistortion}:
+the sky \tbd{convert Cell-$>$Chip and Chip-$>$FP to psDistortion}:
 \begin{itemize}
 \item Cell $\longleftrightarrow$ Chip: two 2D polynomials, $(X,Y) = f(x,y)$;
@@ -3613,19 +3638,34 @@
 
 \begin{verbatim}
-/** Information needed (by SLALib) to convert Apparent to Observed Position */
+/** Information needed (by SLALIB) to convert Apparent to Observed Position */
 typedef struct {
-    double latitude;                    ///< geodetic latitude (radians)
-    double sinLat, cosLat;              ///< sine and cosine of geodetic latitude
-    double abberationMag;               ///< magnitude of diurnal aberration vector
-    double height;                      ///< height (HM)
-    double temperature;                 ///< ambient temperature (TDK)
-    double pressure;                    ///< pressure (PMB)
-    double humidity;                    ///< relative humidity (RH)
-    double wavelength;                  ///< wavelength (WL)
-    double lapseRate;                   ///< lapse rate (TLR)
-    double refractA, refractB;          ///< refraction constants A and B (radians)
-    double longitudeOffset;             ///< longitude + eqn of equinoxes + ``sidereal UT'' (radians)
-    double siderealTime;                ///< local apparent sidereal time (radians)
+    const double latitude;		///< geodetic latitude (radians)
+    const double sinLat, cosLat;	///< sine and cosine of geodetic latitude
+    const double abberationMag;		///< magnitude of diurnal aberration vector
+    const double height;		///< height (HM)
+    const double temperature;		///< ambient temperature (TDK)
+    const double pressure;		///< pressure (PMB)
+    const double humidity;		///< relative humidity (RH)
+    const double wavelength;		///< wavelength (WL)
+    const double lapseRate;		///< lapse rate (TLR)
+    const double refractA, refractB;	///< refraction constants A and B (radians)
+    const double longitudeOffset;	///< longitude + eqn of equinoxes + ``sidereal UT'' (radians)
+    const double siderealTime;		///< local apparent sidereal time (radians)
 } psGrommit;
+\end{verbatim}
+
+The \code{psGrommit} is calculated from telescope information for the
+particular exposure:
+
+\begin{verbatim}
+/** Constructor */
+psGrommit *
+psGrommitAlloc(const psExposure *exp	///< Relevant exposure
+    );
+
+/** Destructor */
+void
+psGrommitFree(psGrommit *grommit	///< Grommit to destroy
+    );
 \end{verbatim}
 
@@ -3664,30 +3704,22 @@
 /** returns Chip in FPA which contains the given FPA coordinate */
 psChip *
-psChipInFPA (psChip *out,               ///< Chip to return, or NULL
-             const psFPA *fpa,          ///< FPA description
-             const psCoord *coord       ///< coordinate in FPA
-             );
-\end{verbatim}
-
-\begin{verbatim}
+psChipInFPA (psChip *out,		///< Chip to return, or NULL
+	     const psPlaneCoord *coord	///< coordinate in FPA
+	     const psFPA *fpa, 		///< FPA description
+	     );
+
 /** returns Cell in Chip which contains the given chip coordinate */
 psCell *
-psCellInChip(psCell *out,               ///< Cell to return, or NULL
-             const psChip *chip,        ///< chip description
-             const psCoord *coord       ///< coordinate in chip
-             );
-\end{verbatim}
-
-Usually we will want to go directly from the FPA to the Cell, so we
-also specify the following function, which performs the above two
-functions in order.
-
-\begin{verbatim}
+psCellInChip(psCell *out,		///< Cell to return, or NULL
+	     const psPlaneCoord *coord	///< coordinate in chip
+	     const psChip *chip,	///< chip description
+	     );
+
 /** Return the cell in FPA which contains the given FPA coordinates */
 psCell *
-psCellInFPA(psCell *out,                ///< Cell to return, or NULL
-            const psFPA *fpa,           ///< FPA description
-            const psCoord *coord        ///< Coordinate in FPA
-            );
+psCellInFPA(psCell *out,		///< Cell to return, or NULL
+	    const psPlaneCoord *coord	///< Coordinate in FPA
+	    const psFPA *fpa,		///< FPA description
+	    );
 \end{verbatim}
 
@@ -3724,108 +3756,87 @@
 \begin{verbatim}
 /** Convert (RA,Dec) to cell and cell coordinates */
-psCoord *
-psCoordSkyToCell(psCoord *out,          ///< Coordinates to return, or NULL
-                 psCell *cell,          ///< Cell to return
-                 const psFPA *fpa       ///< FPA description
-                 );
-\end{verbatim}
-
-\begin{verbatim}
+psPlaneCoord *
+psCoordSkyToCell(psPlaneCoord *out,	///< Coordinates to return, or NULL
+		 psCell *cell,		///< Cell to return
+		 const psSphereCoord *in, ///< Input coordinates
+		 const psFPA *fpa	///< FPA description
+		 );
+
 /** Convert cell and cell coordinate to (RA,Dec) */
-psCoord *
-psCoordCellToSky(psCoord *out,          ///< Coordinates to return, or NULL
-                 const psCell *cell,    ///< Cell to get coordinates for
-                 const psCoord *coord   ///< cell coordinates to transform
-                 );
-\end{verbatim}
-
-\begin{verbatim}
+psSphereCoord *
+psCoordCellToSky(psSphereCoord *out,	///< Coordinates to return, or NULL
+		 const psPlaneCoord *coord ///< cell coordinates to transform
+		 const psCell *cell,	///< Cell to get coordinates for
+		 );
+
 /** Quick and dirty cell to (RA,Dec) --- employs cellToSky transformation */
-psCoord *
-psCoordCellToSkyQuick(psCoord *out,     ///< Coordinates to return, or NULL
-                      const psCell *cell, ///< Cell description
-                      const psCoord *coord ///< cell coordinates to transform
-                      );
-\end{verbatim}
-
-\begin{verbatim}
+psSphereCoord *
+psCoordCellToSkyQuick(psSphereCoord *out, ///< Coordinates to return, or NULL
+		      const psPlaneCoord *coord	///< cell coordinates to transform
+		      const psCell *cell, ///< Cell description
+		      );
+
 /** Convert (RA,Dec) to tangent plane coords */
-psCoord *
-psCoordSkyToTP(psCoord *out,            ///< Coordinates to return, or NULL
-               const psExposure *exp,   ///< Exposure description
-               const psCoord *coord     ///< input Sky coordinate
-               );
-\end{verbatim}
-
-\begin{verbatim}
+psPlaneCoord *
+psCoordSkyToTP(psPlaneCoord *out,	///< Coordinates to return, or NULL
+	       const psSphereCoord *coord ///< input Sky coordinate
+	       const psGrommit *grommit, ///< Grommit for fast conversion
+	       );
+
 /** Convert tangent plane coords to focal plane coordinates */
-psCoord *
-psCoordTPtoFPA(psCoord *out,            ///< Coordinates to return, or NULL
-               const psFPA *fpa,        ///< FPA description
-               const psCoord *coord     ///< input TP coordinate
-               );
-\end{verbatim}
-
-\begin{verbatim}
+psPlaneCoord *
+psCoordTPtoFPA(psPlaneCoord *out,	///< Coordinates to return, or NULL
+	       const psPlaneCoord *coord ///< input TP coordinate
+	       const psFPA *fpa,	///< FPA description
+	       );
+
 /** converts the specified FPA coord to the coord on the given Chip */
-psCoord *
-psCoordFPAtoChip (psCoord *out,         ///< Coordinates to return, or NULL
-                  const psChip *chip,   ///< Chip of interest
-                  const psCoord *coord  ///< input FPA coordinate
-                  ); 
-\end{verbatim}
-
-\begin{verbatim}
+psPlaneCoord *
+psCoordFPAtoChip (psPlaneCoord *out,	///< Coordinates to return, or NULL
+		  const psPlaneCoord *coord ///< input FPA coordinate
+		  const psChip *chip,	///< Chip of interest
+		  ); 
+
 /** converts the specified Chip coord to the coord on the given Cell */
-psCoord *
-psCoordChiptoCell (psCoord *out,        ///< Coordinates to return, or NULL
-                   const psCell *cell,  ///< Cell of interest
-                   const psCoord *coord ///< input Chip coordinate
-                   );
-\end{verbatim}
-
-\begin{verbatim}
+psPlaneCoord *
+psCoordChiptoCell (psPlaneCoord *out,	///< Coordinates to return, or NULL
+		   const psPlaneCoord *coord ///< input Chip coordinate
+		   const psCell *cell, 	///< Cell of interest
+		   );
+
 /** converts the specified Cell coord to the coord on the parent Chip */
-psCoord *
-psCoordCelltoChip (psCoord *out,        ///< Coordinates to return, or NULL
-                   const psCell *cell,  ///< Cell description
-                   const psCoord *coord ///< input Cell coordinate
-                   );
-\end{verbatim}
-
-\begin{verbatim}
+psPlaneCoord *
+psCoordCelltoChip (psPlaneCoord *out,	///< Coordinates to return, or NULL
+		   const psPlaneCoord *coord ///< input Cell coordinate
+		   const psCell *cell, 	///< Cell description
+		   );
+
 /** converts the specified Chip coord to the coord on the parent FPA */
-psCoord *
-psCoordChiptoFPA (psCoord *out,         ///< Coordinates to return, or NULL
-                  const psChip *chip,   ///< Chip description
-                  const psCoord *coord  ///< input Chip coordinate
-                  );
-\end{verbatim}
-
-\begin{verbatim}
+psPlaneCoord *
+psCoordChiptoFPA (psPlaneCoord *out,		///< Coordinates to return, or NULL
+		  const psPlaneCoord *coord	///< input Chip coordinate
+		  const psChip *chip, 	///< Chip description
+		  );
+
 /** Convert focal plane coords to tangent plane coordinates */
-psCoord *
-psCoordFPAToTP(psCoord *out,            ///< Coordinates to return, or NULL
-               const psFPA *fpa,        ///< FPA description
-               const psCoord *coord     ///< input FPA coordinate
-               );
-\end{verbatim}
-
-\begin{verbatim}
+psPlaneCoord *
+psCoordFPAToTP(psPlaneCoord *out,		///< Coordinates to return, or NULL
+	       const psPlaneCoord *coord ///< input FPA coordinate
+	       const psFPA *fpa,	///< FPA description
+	       );
+
 /** Convert tangent plane coords to (RA,Dec) */
-psCoord *
-psCoordTPtoSky(psCoord *out,            ///< Coordinates to return, or NULL
-               const psExposure *exp,   ///< Exposure description
-               const psCoord *coord     ///< input TP coordinate
-               );
-\end{verbatim}
-
-\begin{verbatim}
+psSphereCoord *
+psCoordTPtoSky(psSphereCoord *out,	///< Coordinates to return, or NULL
+	       const psPlaneCoord *coord ///< input TP coordinate
+	       const psGrommit *grommit, ///< Grommit for fast conversion
+	       );
+
 /** Convert Cell coords to FPA coordinates */
-psCoord *
-psCoordCellToFPA(psCoord *out,          ///< Coordinates to return, or NULL
-                 const psCell *cell,    ///< Cell description
-                 const psCoord *coord   ///< Input cell coordinates
-                 );
+psPlaneCoord *
+psCoordCellToFPA(psPlaneCoord *out,	///< Coordinates to return, or NULL
+		 const psPlaneCoord *coord ///< Input cell coordinates
+		 const psCell *cell,	///< Cell description
+		 );
 \end{verbatim}
 
@@ -3840,7 +3851,7 @@
 /** Get the airmass for a given position and sidereal time */
 float
-psGetAirmass(const psCoord *coord,      ///< Position on the sky
-             double siderealTime,       ///< Sidereal time
-             float height               ///< Height above sea level
+psGetAirmass(const psSphereCoord *coord, ///< Position on the sky
+             double siderealTime,	///< Sidereal time
+	     float height		///< Height above sea level
              );
 \end{verbatim}
@@ -3849,5 +3860,5 @@
 /** Get the parallactic angle for a given position and sidereal time */
 float
-psGetParallactic(const psCoord *coord,  ///< Position on the sky
+psGetParallactic(const psSphereCoord *coord, ///< Position on the sky
                  double siderealTime    ///< Sidereal time
                  );
@@ -3893,9 +3904,9 @@
 \begin{verbatim}
 typedef struct {
-    int ID;                             ///< ID number for this photometric system
-    char *name;                         ///< Name of photometric system
-    char *camera;                       ///< Camera for photometric system
-    char *filter;                       ///< Filter used for photometric system
-    char *detector;                     ///< Detector used for photometric system
+    const int ID;                       ///< ID number for this photometric system
+    const char *name;                   ///< Name of photometric system
+    const char *camera;                 ///< Camera for photometric system
+    const char *filter;                 ///< Filter used for photometric system
+    const char *detector;               ///< Detector used for photometric system
 } psPhotSystem;
 \end{verbatim}
@@ -3939,4 +3950,7 @@
 M_{\rm pM} - pA, M_{\rm sP} - M_{\rm sM} - sA)$.
 
+\TBD{Really want a set of polynomials defined for specific colour
+ranges.}
+
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@@ -3952,28 +3966,22 @@
 \begin{verbatim}
 /** Get Sun Position */
-psCoord *
-psGetSunPos(float mjd                   ///< MJD to get position for
-    );
-\end{verbatim}
-
-\begin{verbatim}
+psSphereCoord *
+psGetSunPos(float mjd			///< MJD to get position for
+    );
+
 /** Get Moon position */
-psCoord *
-psGetMoonPos(float mjd,                 ///< MJD to get position for
-             double latitude,           ///< Latitude for apparent position
-             double longitude           ///< Longitude for apparent position
-    );
-\end{verbatim}
-
-\begin{verbatim}
+psSphereCoord *
+psGetMoonPos(float mjd,			///< MJD to get position for
+	     double latitude,		///< Latitude for apparent position
+	     double longitude		///< Longitude for apparent position
+    );
+
 /** Get Moon phase */
 float
 psGetMoonPhase(float mjd                ///< MJD to get phase for
     );
-\end{verbatim}
-
-\begin{verbatim}
+
 /** Get Planet positions */
-psCoord *
+psSphereCoord *
 psGetSolarSystemPos(const char *solarSystemObject, ///< Named S.S. object
                     float mjd           ///< MJD to get position for
@@ -3994,24 +4002,24 @@
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
-\appendix
-
-\pagebreak
-\section{API Summary: all functions}
-
-\subsection{System Utilities}
-\input{psSystemGroup.tex}
-
-\subsection{Data Containers}
-\input{psDataGroup.tex}
-
-\subsection{Math Utilities}
-\input{psMathGroup.tex}
-
-\subsection{Astronomy Functions}
-\input{psAstroGroup.tex}
-
-\pagebreak
-\section{API Summary: all structures}
-\input{psStructures.tex}
+%\appendix
+%
+%\pagebreak
+%\section{API Summary: all functions}
+%
+%\subsection{System Utilities}
+%\input{psSystemGroup.tex}
+%
+%\subsection{Data Containers}
+%\input{psDataGroup.tex}
+%
+%\subsection{Math Utilities}
+%\input{psMathGroup.tex}
+%
+%\subsection{Astronomy Functions}
+%\input{psAstroGroup.tex}
+%
+%\pagebreak
+%\section{API Summary: all structures}
+%\input{psStructures.tex}
 
 \bibliographystyle{plain} \bibliography{panstarrs}
