Index: trunk/doc/release.2015/ps1.analysis/analysis.tex
===================================================================
--- trunk/doc/release.2015/ps1.analysis/analysis.tex	(revision 39947)
+++ trunk/doc/release.2015/ps1.analysis/analysis.tex	(revision 39948)
@@ -237,5 +237,5 @@
 \code{psLib}.  Components of the photometry code were integrated into
 the IPP's mid-level astronomy data analysis toolkit called
-\code{psModules}.  The resulting software, `PSPhot', can be used either
+\code{psModules}.  The resulting software, `\code{psphot}', can be used either
 as a stand-alone C program, or as a set of library functions which may
 be integrated into other programs
@@ -243,5 +243,5 @@
 \note{add refs to the psLib and psModules ADDs}
 
-The main version of PSPhot is a stand-alone program which is run on a
+The main version of \code{psphot} is a stand-alone program which is run on a
 single image, or a group of related images representing the data read
 from a camera in a single exposure.  The images are expected to have
@@ -253,10 +253,10 @@
   integrated library call}
 
-The version called PSPhotStack accepts a set of images, each
+The version called \code{psphotStack} accepts a set of images, each
 representing the same patch of sky in a different filter, nominally
 the full $grizy$ filter set for the analysis of the PS1 PV3 stack
 images, though where insufficient data were available in a given
 filter, a subset of these filters was processed as a group.  As
-discussed in detail below, the PSPhotStack analysis includes the
+discussed in detail below, the \code{psphotStack} analysis includes the
 capability of measuring forced PSF photometry in some filter images
 based on the position of sources detected in the other filters.  It
@@ -265,6 +265,6 @@
 photometry.
 
-Another version of PSPhot used in the PV3 analysis is called
-PSPhotFullForce.  In this version, a set of image all representing the
+Another version of \code{psphot} used in the PV3 analysis is called
+\code{psphotFullForce}.  In this version, a set of image all representing the
 same pixels are processed together, with the positions of sources to
 be analysed loaded from a supplied file.  In this version of the
@@ -276,18 +276,18 @@
 supplied guess model.  
 
-\section{PSPhot Design Goals}
-
-PSPhot has a number of important requirements that it must meet, and a
+\section{\code{psphot} Design Goals}
+
+\code{psphot} has a number of important requirements that it must meet, and a
 number of design goals which we believe will help to make usable in a
 wide range of circumstances.  The critical requirements of the
-Pan-STARRS IPP which drive the requirements for PSPhot:
+Pan-STARRS IPP which drive the requirements for \code{psphot}:
 
 \begin{itemize}
-\item {\bf 10 millimagnitude photometric accuracy}.  For PSPhot, this
+\item {\bf 10 millimagnitude photometric accuracy}.  For \code{psphot}, this
   implies that the measured photometry of stellar sources must be
   substantially better than this 10 mmag since the photometry error
   per image is combined with an error in the flat-field calibration
   and an error in measuring the atmospheric effects.  We have set a
-  goal for PSPhot of 3mmag photometric consistency for bright stars
+  goal for \code{psphot} of 3mmag photometric consistency for bright stars
   between pairs of images obtained in photometric conditions at the
   same pointing, ie to remove sensitivity to flat-field errors.  This
@@ -298,14 +298,14 @@
 \item {\bf 10 milliarcsecond astrometric accuracy}. Relative
   astrometric calibration depends on the consistency of the individual
-  measurements.  The measurements from PSPhot must be sufficiently
+  measurements.  The measurements from \code{psphot} must be sufficiently
   representative of the true source position to enable astrometric
   calibration at the 10mas level.  The error in the individual
   measurements will be folded together with the errors introduced by
   the optical system, the effects of seeing, and by the available
-  reference catalogs.  We have set a goal for PSPhot of 5mas
+  reference catalogs.  We have set a goal for \code{psphot} of 5mas
   consistency between the true source postion and the measured
   position given reasonable PSF variations under simulations.  This
   level must be reached for images with 250 mas pixels, implying
-  PSPhot must introduce measurement errors less than 1/50th of a
+  \code{psphot} must introduce measurement errors less than 1/50th of a
   pixel. The choice of 32 bit floating point data values for the
   source centroids places a numerical limit of 1e-7 on the accuracy of
@@ -315,5 +315,5 @@
 \end{itemize}
 
-The design goals for PSPhot are chosen to make the program flexible,
+The design goals for \code{psphot} are chosen to make the program flexible,
 general, and able to meet the unknown usages cases future projects may
 require:
@@ -328,31 +328,31 @@
   naturally incorporate 2-D variations.
 
-\item {\bf Flexible non-PSF models} PSPhot must be able to represent
+\item {\bf Flexible non-PSF models} \code{psphot} must be able to represent
   PSF-like sources as well as non-PSF sources (e.g., galaxies).  It
   must be easy to add new source models as interesting representations
   of sources are invented.
 
-\item {\bf Clean code base} PSPhot should incorporate a high-degree of
+\item {\bf Clean code base} \code{psphot} should incorporate a high-degree of
   abstraction and encapsulation so that changes to the code structure
   can be performed without pulling the code apart and starting from scratch.
 
-\item {\bf PSF validity tests} PSPhot should include the ability to
+\item {\bf PSF validity tests} \code{psphot} should include the ability to
   choose different types of PSF models for diffent situations, or to
   provide the user with methods for assessing the different PSF models.
 
-\item {\bf Careful systematic corrections} PSPhot must carefully
+\item {\bf Careful systematic corrections} \code{psphot} must carefully
   measure and correct for the photometric and astrometric trends
   introduced by using analytical PSF models.
 
-\item {\bf User Configurable} PSPhot should allow users to change the
+\item {\bf User Configurable} \code{psphot} should allow users to change the
   options easily and to allow different approaches to the analysis.
 
 \end{itemize}
 
-\section{PSPhot Analysis Process}
+\section{\code{psphot} Analysis Process}
 
 \subsection{Overview}
 
-The PSPhot analysis is divided into several major stages:
+The \code{psphot} analysis is divided into several major stages:
 
 \begin{enumerate}
@@ -383,5 +383,5 @@
 \end{enumerate}
 
-PSPhot is highly configurable.  Users may choose via the configuration
+\code{psphot} is highly configurable.  Users may choose via the configuration
 system which of the above analyses are performed.  This is useful for
 testing, but also allows for specialized use cases.  For example, the
@@ -405,5 +405,5 @@
 references to the mask and variance are provided in the configuration
 information.  As in the stand-alone C-program, the variance and mask may
-be constructed automatically by PSPhot.
+be constructed automatically by \code{psphot}.
 
 The mask is represented as a 16-bit integer image in which a value of
@@ -445,5 +445,5 @@
 
 \begin{table*}
-\caption{\label{tab:mask_values} PSPhot / GPC1 Mask Image Pixel Values}\vspace{-0.5cm}
+\caption{\label{tab:mask_values} \code{psphot} / GPC1 Mask Image Pixel Values}\vspace{-0.5cm}
 \begin{center}
 \begin{tabular}{lcl}
@@ -650,5 +650,5 @@
 The peaks detected in the image may correspond to real sources, but
 they may also correspond to noise fluctuations, especially in the
-wings of bright stars.  PSPhot attempts to identify peaks which may be
+wings of bright stars.  \code{psphot} attempts to identify peaks which may be
 formally significant, but are not locally significant.  It first
 generates a set of ``footprints'', contiguous collections of pixels in
@@ -830,5 +830,11 @@
 \subsubsection{PSF Model vs Source Model}
 
-The PSF model used by \code{psphot} consists of an analytical function
+The point-spread-function (PSF) of an image describes the shape of all
+unresolved sources in the image.  In a typical wide-field image, the
+shape of unresolved sources varies as a function of position in the
+image.  The full PSF thus needs to include a model with parameters
+which vary across the image.
+
+The PSF used by \code{psphot} consists of an analytical function
 combined with a pixelized representation of the residual differences
 between the analytical model and the true PSF.  Both the shape
@@ -837,8 +843,7 @@
 
 Within \code{psphot}, several analytical models may be used to
-describe the PSF, but all share a few common characteristics.
-
-Any source in an image may be represented by some analytical model,
-for example, a 2-D elliptical Gaussian:
+describe the smooth portion of the PSF, but all share a few common
+characteristics.  As an example, a simple model consists of a 2-D
+elliptical Gaussian:
 \begin{eqnarray}
 f(x,y) & = & I_o e^{-z} + S  \\
@@ -847,31 +852,24 @@
     y  & = & y_{\rm ccd} - y_o 
 \end{eqnarray}
-The source model will have a variety of model parameters, in this case
-the centroid coordinates ($x_o, y_o$), the elliptical shape parameters
-($\sigma_x, \sigma_y, \sigma_{\rm xy}$), the model normalization
-($I_o$) and the local value of the background ($S$).  A specific
-source will have a particular set of values for these different
-parameters.
-
-The point-spread-function (PSF) of an image describes the shape of all
-unresolved sources in the image.  In a typical image, the shape of
-point sources is not well described by a single function.  Instead,
-the shape will vary as a function of position in the image.  The PSF
-model therefore must describe the parameter variation as a function of
-the position of the source on the image.  Note that the source model
-consists of a certain number of parameters which are defined by the
-PSF model, and another set of parameters which are independent from
-source to source.  For the case of the elliptical Gaussian model, the
-PSF parameters would be the shape terms ($\sigma_x, \sigma_y,
-\sigma_{\rm xy}$) while the independent parameters would be the
-centroid, normalization and local sky values ($x_o, y_o, I_o, S$).
-Thus these parameters are each a function of the source centroid
+Here the model parameters consist of the centroid coordinates ($x_o,
+y_o$), the elliptical shape parameters ($\sigma_x, \sigma_y,
+\sigma_{\rm xy}$), the model normalization ($I_o$) and the local value
+of the background ($S$).  
+
+A specific source will have a particular set of values for the model
+parameters, some of which depend on the PSF model and the position of
+the source in the image, while the rest are unique to the individual
+source.  For the case of the elliptical Gaussian model, the PSF
+parameters would be the shape terms ($\sigma_x, \sigma_y, \sigma_{\rm
+  xy}$) while the independent parameters would be the centroid,
+normalization and local sky values ($x_o, y_o, I_o, S$).  Thus the
+shape parameters are each a function of the source centroid
 coordinates:
 \begin{eqnarray}
-\sigma_x    & = & f_1(x,y) \\
-\sigma_y    & = & f_2(x,y) \\
-\sigma_{xy} & = & f_3(x,y) \\
+\sigma_x    & = & f_1(x_{\rm ccd},y_{\rm ccd}) \\
+\sigma_y    & = & f_2(x_{\rm ccd},y_{\rm ccd}) \\
+\sigma_{xy} & = & f_3(x_{\rm ccd},y_{\rm ccd}) \\
 \end{eqnarray}
-PSPhot represents the variation in the PSF parameters as a function of
+\code{psphot} represents the variation in the PSF parameters as a function of
 position in the image in two possible ways, specified by the
 configuration.  The first option is to use a 2-D polynomial which is
@@ -890,63 +888,56 @@
 some of the observed PSF variations in the images
 
-% XXX specify the rule for the polynomial order and grid scale
-% XXX discuss the improvements in the astrometric modeling PV1 - PV3
-
-PSPhot uses a single structure to represent the source model and
-another structure to represent the PSF model.  The source model
-structure consists of the collection of measured source model
-parameters, carried as a \code{psLib} vector (\code{psVector}) along
-with an equal-length vector with the parameter errors.  The structure
-also includes an integer giving the identifier of the model used in
-the particular case, as well as model fit statistics such as the
-Chi-Square of the fit and the magnitude representation of the ratio
-between the model flux and an aperture flux (see below for more
-details on this value).
-
-The PSPhot representation of the PSF consists of an array of
-polynomials, each representing the variation in the source model PSF
-parameters (\code{psArray} of \code{psPolynomial2D}).  The PSF model
-structure also includes the same integer used to identify which model
-corresponds to particular instance of the PSF.  At the moment, the
-number of PSF parameters is a fixed number (4) fewer than the number
-of parameters of the corresponding source model.  For example, the
-elliptical Gaussian model uses 7 parameters to represent the source and
-3 for the PSF model.  
-
-PSPhot is written so that the source detection, measurement, and
-classification code does not depend on the specific form of the
-available source model functions.  Access to the characteristics of
-the models is provided through a simple function abstraction method.
-Throughout PSPhot, there are many places where it is necessary for the
-code to refer to an aspect of the source or PSF model.  Often, these
-quantities are needed deep within other parts of the code.  For
-example, when attempting to fit the pixel flux values for a source,
-it is necessary to generate a guess for the model parameters.  Or, in
-order to limit the domain of the fit, it is necessary to determine an
-isophotal radius for a model.  
-
-In order to avoid having the code depend on the specific form of a
-model, the function calls needed in these types of circumstances are
-abstracted, and a method is provided to return the necessary function
-to the higher-level software.  For example, each model type has its
-own function to define an initial guess for the model, or a function
-to determine the radius for a given flux level.  These are then
-registered as part of the model function code.  Another function is
-then used to return the appropriate function for a specific model
-type.  For example, the \code{psModelLookup_GetFunction} will return
-the \code{psModelLookup} function for a given model type.  This
-mechanism makes it very easy to add new model functions into the
-PSPhot code base.  To add a new model function, the programmer simply
-defines a new model name (a string), the set of all necessary model
-lookup functions, and places the reference to the model code at the
-appropriate location in the psModelInit.c routine.
-
-When a new model is provided to PSPhot, it is not necessary to specify
-the intended use of the source model function (ie, PSF-like source,
-galaxy, comet, etc).  Any model can be used for the PSF model, or to
-describe the flux distributions of the non-PSF sources.  The code
-currently uses a fixed translation between the source model parameters
-and the PSF model parameters.  It also defines a specific order for
-the 4 independent parameters.  
+\note{need to describe fitting the pixel residual image}
+
+\note{write up the fitting process to define the grid?}
+
+\notespecify the rule for the polynomial order and grid scale}
+
+\note{discuss the improvements in the astrometric modeling PV1 - PV3}
+
+Several analytical functions which are likely candidates to describe
+the smooth portion of the PSF are available in \code{psphot}:
+\begin{itemize}
+\item Gaussian : $f = I_0 e^{-z}$
+\item Pseudo-Gaussian : $f = I_0 (1 + z + \frac{1}{2} z^2 + \frac{1}{6} z^3)^{-1}$ \code{[PGAUSS]}
+\item Variable Power-Law : $f = I_0 (1 + z + z^{\alpha})^{-1}$ \code{[RGAUSS]}
+\item Steep Power-Law : $f = I_0 (1 + \kappa z + z^{2.25})^{-1}$ \code{[QGAUSS]}
+\item PS1 Power-Law : $f = I_0 (1 + \kappa z + z^{1.67})^{-1}$ \code{[PS1_V1]}
+\end{itemize}
+The Pseudo-Gaussian is a Taylor expansion of the Gaussian and is used
+by Dophot \citep{1993PASP..105.1342S}.  The latter profiles are
+similar to the Moffat profile form
+\citep{1969AA.....3..455M,1983AA...126..278B}, with small differences.
+A user may choose to try more than one analytical function for a given
+image.  As discussed below (Section~\ref{sec:psf.model.choice}),
+\code{psphot} can automatically choose the best model based on the
+quality of the PSF fits.
+
+For the PS1 GPC1 analysis, we used the \code{PS1_V1} model, which we
+found by experimentation to match well to the observed profiles
+generated by PS1.  Figure~\ref{fig:radial.profiles} shows example
+radial profiles for moderately bright stars in fairly good (0.9
+arcsec) and poor (2.2 arcsec) seeing.  Using a fixed power-law
+exponent results in somewhat faster profile fitting compared to the
+variable power-law exponent model.
+
+The analytical models in \code{psphot} are written with a high degree
+of code abstraction making it relatively easy to add different
+analytical models to the software.  The same portion of code used to
+describe the analytical portion of the PSF sources is also used to for
+galaxy models. 
+
+% moffat : 1969A&A.....3..455M
+% buonanno : 1983A&AS...51...83B
+
+\begin{figure}[htbp]
+  \begin{center}
+  \includegraphics[width=\hsize]{{pics/radial.profiles}.\plotext}
+  \caption{\label{fig:radial.profiles} Radial profiles of stellar images from PS1.  These two
+    profiles illustrate the radial trend of the PS1 PSFs for a star
+    with FWHM 0.9 arcsec (red) and 2.2 arcsec (blue).  The black line
+    shows the PSF model with radial trend of the form $(1 + \kappa r^2 + r^{3.33})^{-1}$.}
+  \end{center}
+\end{figure}
 
 \subsubsection{Candidate PSF Source Selection}
@@ -955,29 +946,26 @@
 The first stage of determining the PSF model for an image is to
 identify a collection of sources in the image which are {\em likely}
-to be PSF-like.  PSPhot uses the source moments to make the initial
-guess at a collection of PSF-like sources.  At this point, the program
-has measured the second order moments for all sources identified by
-their peaks, as well as an approximate signal-to-noise ratio.  All
-sources with a S/N ratio greater than a user-defined parameter
-(\code{PSF_SHAPE_NSIGMA} = 20.0) are selected by PSPhot, though
-sources which have more than a certain number of saturated pixels are
-excluded at this stage.  PSPhot then examines the 2-D plane of
-$\sigma_x, \sigma_y$ in search of a concentrated clump of sources (see
-Figure~\ref{fig:moment.class}).  To
-do this, it constructs an artificial image with pixels representing
-the value of $\sigma_x, \sigma_y$, using a user-defined scale for the
-size of a pixel in this artificial image (note that the units of the
-$\sigma_x, \sigma_y$ plane are the size of the second-moment in pixels
-in the original image).  A typical value for the bin size is
-approximately 0.1 image pixels.  The binned $\sigma_x, \sigma_y$ plane
-is then examined to find a peak which has a significance greater than
-XXX.  Unless the image is extremely sparse, such a peak will be
-well-defined and should represent the sources which are all very
-similar in shape.  Other sources in the image will tend to land in
-very different locations, failing to produce a single peak.  To avoid
-detecting a peak from the unresolved cosmic rays, sources which have
-second-moments very close to 0 are ignored.  The only danger is if the
-PSF is very small and too many of these sources are rejected as cosmic
-rays.
+to be unresolved (i.e., stars).  \code{psphot} uses the source sizes as
+estimated from the second moments to make the initial guess at a
+collection of unresolved sources.  At this point, the program has
+measured the second order moments for all sources identified by their
+peaks, as well as an approximate signal-to-noise ratio, above the
+bright threshold.  All sources with a S/N ratio greater than a
+user-defined parameter (\code{PSF_SN_LIM} = 20.0 for PS1 PV3) are
+selected by \code{psphot}, though sources which have more than a
+certain number of saturated pixels are excluded at this stage.  The
+program then examines the 2-D plane of $\sigma_x, \sigma_y$ in search
+of a concentrated clump of sources (see
+Figure~\ref{fig:moment.class}).  To do this, it constructs an
+artificial image with pixels representing the value of $\sigma_x,
+\sigma_y$, using $0.1 \sigma_w$ as the size of a pixel in this
+artificial image.  The binned $\sigma_x, \sigma_y$ plane is then
+examined to find a significant peak.  Unless the image is extremely
+sparse, such a peak will be well-defined and should represent the
+sources which are all very similar in shape.  Other sources in the
+image will tend to land in very different locations, failing to
+produce a single peak.  To avoid detecting a peak from the unresolved
+cosmic rays, sources which have second-moments very close to 0 are
+ignored.
 
 Once a peak has been detected in this plane, the centroid and second
@@ -985,4 +973,6 @@
 $\sigma$ of this centroid are selected as likely PSF-like sources in
 the image.  
+
+\note{work out the logic for selecting the PSF stars}
 
 \begin{figure}[htbp]
@@ -1000,14 +990,18 @@
 
 \subsubsection{Candidate PSF Source Model Fits}
-
+\label{sec:psf.model.choice}
 % \note{link to psLibADD}
+
+% Madsen:
+%% http://www2.imm.dtu.dk/pubdb/views/edoc_download.php/3215/pdf/imm3215.pdf
+% Press
 
 All candidate PSF sources are then fitted with the selected source
 model, allowing all of the parameters (PSF and independent) to vary in
-the fit.  PSPhot uses the Levenberg-Marquardt minimization technique
-for the non-linear fitting.  Non-linear
+the fit.  The software uses the Levenberg-Marquardt minimization
+technique \citep{Press,Madsen} for the non-linear fitting.  Non-linear
 fitting can be very computationally intensive, particularly for if the
-starting parameters are far from the minimization values.  PSPhot uses
-the first and second moments to make a good guess for the centroid and
+starting parameters are far from the minimization values.  The first
+and second moments are used to make a good guess for the centroid and
 shape parameters for the PSF models.  Any sources which fail to
 converge in the fit are flagged as invalid.
@@ -1015,32 +1009,43 @@
 For the resulting collection of source model parameters, the
 PSF-dependent parameters of the models are all fitted as a function of
-position to a 2-D polynomial.  The order of this polynomial is a
-user-defined parameter.  The fitting process for these polynomials is
-iterative, and rejects the $3-\sigma$ outliers in each of three
-passes.  This fitting technique results in a robust measurement of the
-variation of the PSF model parameters as a function of position
-without being excessively biased by individual sources which fail
-drastically.  Sources whose model parameters are rejected by this
-iterative fitting technique are also marked as invalid and ignored in
-the later PSF model fitting stages.
+position using either the 2-D polynomial or the gridded superpixel
+representation.  The maximum order of these fits depends on the number
+of PSF sources (see Table~\ref{tab:order}).  The fitting process for
+these polynomials is iterative, and rejects the $3\sigma$ outliers in
+each of three passes.  This fitting technique results in a robust
+measurement of the variation of the PSF model parameters as a function
+of position without being excessively biased by individual sources
+which are not well described by the PSF model (e.g., galaxies which
+snuck into the sample).  Sources whose model parameters are rejected
+by this iterative fitting technique are also marked as invalid PSF
+sources and ignored in the later PSF model fitting stages.
+
+%% table of orders:
+%% N stars | max order | max Ncells
+%%  16   |  1; //  4 cells, 4 per cell
+%%  54   |  2; //  9 cells, 6 per cell
+%% 128   |  3; // 16 cells, 8 per cell
+%% 300   |  4; // 25 cells, 12 per cell
+%% 576   |  5; // 36 cells, 16 per cell
 
 All of the PSF-candidate sources are then re-fitted using the PSF
-model to specify the dependent model parameter values for each source.
-For example, in the case of the elliptical Gaussian model, the shape
-parameters ($\sigma_x, \sigma_y, \sigma_{xy}$) for each source are
-set by the coordinates of the source centroid and fixed (not allowed
-to vary) in the fitting procedure.  The resulting fitted models are
-then used to determine a metric which tests the quality of the PSF
-model for this particular image.  
-
-The metric used by PSPhot to assess the PSF model is the scatter in
+model to specify the PSF-dependent model parameter values for each
+source.  For example, in the case of the elliptical Gaussian model,
+the shape parameters ($\sigma_x, \sigma_y, \sigma_{xy}$) for each
+source are set by the coordinates of the source centroid and fixed
+(not allowed to vary) in the fitting procedure.  The resulting fitted
+models are then used to determine a metric which tests the quality of
+the PSF model for this particular image.
+
+The metric used by \code{psphot} to assess the PSF model is the scatter in
 the differences between the aperture and fit magnitudes for the PSF
-sources.  The difference between the aperture and fit magnitudes ({\em
-ApResid}) is a critical parameter for any PSF modeling software which
-uses an analytical model to represent the flux distribution of the
-sources in an image.  An approximate correction is measured here, with
-a more detailed correction measured after all source analysis is
-performed.  The PSF model with the best consistency of the aperture
-correction is judged to be the best model.
+sources.  This difference is a critical parameter for any PSF modeling
+software as it is a measurement of how well the PSF model captures the
+flux of the star.  An approximate correction is measured here, with a
+more detailed correction measured after all source analysis is
+performed (see Section~\ref{sec:aperture.correction}).  The PSF model
+with the best consistency of the aperture correction is judged to be
+the best model.  \note{are we making a decision on the order or
+  anything based on apresid?}
 
 \subsection{Bright Source Analysis}
@@ -1076,10 +1081,15 @@
 \subsubsection{Very Bright Stars}
 
-The PSF modeling code fails to fit the wings of highly saturated stars
-if the core of the star is too contaminated by saturated pixels. For
-stars with estimated instrumental magnitudes brighter than XXX, we fit
-and subtract a radial profile modeled with a spline (?).
-
-\note{more here}
+The standard \code{psphot} PSF modeling code fails to fit the wings of
+highly saturated stars, especially if the core of the star is too
+contaminated by saturated pixels.  For stars with more than a single
+saturated pixel, we model the radial profile of the logarithmic
+instrumental flux in logarithmically spaced radial bins.  For each
+radial bin, we determine the median of the log-flux.  This median
+profile is then interpolated to generate the full radial flux
+distribution.
+
+% logRdel = 0.1
+% logRmax = log(320)
 
 \subsubsection{Fast Ensemble PSF Fitting}
@@ -1113,5 +1123,5 @@
 diagonal square matrix.  The dimension is the number of sources,
 likely to be 1000s or 10,000s.  Direct inversion of the matrix would
-be computationally very slow.  However, an interative solution quickly
+be computationally very slow.  However, an iterative solution quickly
 yields a result with sufficient accuracy.  In the iterative solution,
 a guess at the solution $\bar{A}$ is made assuming $M_{i,j}$ is purely
@@ -1119,8 +1129,12 @@
 compared with the observed vector $\bar{F_j}$.  The difference is used
 to modify the initial guess.  This proces is repeated several times to
-achieve a good convergence.
+achieve a good convergence.  Convergence is quick (a few iterations)
+because of the highly diagonal matrix with small off-diagonal terms:
+the dot product of source $i$ and source $j$ is 1 where $i = j$ and
+much less than 1 where $i \noteq j$.
 
 Once a solution set for $A_i$ is found, all of the sources are
-subtracted from the by applying these values to the unit-flux PSF.
+subtracted from the image by applying these values to the unit-flux
+PSF.
 
 \subsubsection{Full PSF Model Fitting}
@@ -1128,14 +1142,17 @@
 % \note{review the discussion below}
 
-Once a PSF model has been selected for an image, PSPhot attempts to
+Once a PSF model has been selected for an image, \code{psphot} attempts to
 fit all of the detected sources, above a user-defined signal-to-noise
 ratio with the PSF model.  For these fits, the dependent parameters
 are fixed by the PSF model and only the 4 independent source model
-parameters are allowed to vary in the fit.  PSPhot again uses
+parameters are allowed to vary in the fit.  \code{psphot} again uses
 Levenberg-Marquardt minimization for the non-linear fitting.  The
 sources are fitted in their S/N order, starting with the brightest and
-working down to the user-specified limit.
-
-Once a solution has been achieved for a source, PSPhot attempts to
+working down to the user-specified limit, with the other sources
+subtracted as discussed above.
+
+\node{code review for the next bit}
+
+Once a solution has been achieved for a source, \code{psphot} attempts to
 judge the quality of the PSF model as a representation of the source
 shape.  To do this, it calculates the next step of the minimization
@@ -1149,5 +1166,5 @@
 $\sigma_y$.  For a generic model, the shape parameters may be defined
 differently, but there should always be two parameters which scale the
-source size in two dimensions.  Currently, PSPhot requires the two
+source size in two dimensions.  Currently, \code{psphot} requires the two
 relevant shape parameters to be the first two dependent parameters in
 the list of model parameters (ie, parameters 4 \& 5).
@@ -1172,5 +1189,5 @@
 as a likely defect.  
 
-At this stage of the analysis, PSPhot uses two additional indicators
+At this stage of the analysis, \code{psphot} uses two additional indicators
 to identify good and poor PSF fits.  The first of these is the
 signal-to-noise ratio.  It is possible for the peak finding algorithm
@@ -1183,5 +1200,5 @@
 smoothed image).  The fit can either fail to converge or it can
 converge on a fit with very low or negative peak flux / flux
-normalization.  PSPhot will flag any non-convergent PSF fit and any
+normalization.  \code{psphot} will flag any non-convergent PSF fit and any
 source with PSF S/N ratio lower than a user-defined cutoff.  It is
 also useful to identify very poor fits by setting a maximum Chi-Square
@@ -1224,6 +1241,8 @@
 \label{sec:source.size}
 
+\note{is this in the right place?}
+
 After the PSF model has been fitted to all sources, and the Kron flux
-has been measured for all sources, PSPhot uses these two measurements,
+has been measured for all sources, \code{psphot} uses these two measurements,
 along with some additional pixel-level analysis, to determine the size class
 of the source.  If the source is large compared to a PSF, it is
@@ -1233,13 +1252,14 @@
 PSF, it is considered to be a {\em cosmic ray} and masked. 
 
-Extended sources are identified as those for which the Kron magnitude is
-significantly brighter than the PSF magnitude when compared to a PSF
-star.  The value $dMagKP = m_{\rm Kron} - m_{\rm PSF}$, the difference between the PSF 
-and Kron magnitudes, is calculated for each source.  The median of
-$dMagKP$ is calculated for the PSF stars.  This median is subtracted
-from $dMagKP$ for each star.  The result is divided by the quadrature
-error of the PSF and Kron magnitudes and called \code{extNsigma}.  If
-\code{extNsigma} is larger than \code{PSPHOT.EXT.NSIGMA.LIMIT} (3.0),
-the source is considered to be extended.
+Extended sources are identified as those for which the Kron magnitude
+is significantly brighter than the PSF magnitude when compared to a
+PSF star.  The value $\delta M_{rm KP} = m_{\rm Kron} - m_{\rm PSF}$,
+the difference between the PSF and Kron magnitudes, is calculated for
+each source.  The median of $\delta M_{rm KP}$ is calculated for the
+PSF stars.  This median is subtracted from $\delta M_{rm KP}$ for each
+star.  The result is divided by the quadrature error of the PSF and
+Kron magnitudes and called \code{extNsigma}.  If \code{extNsigma} is
+larger than \code{PSPHOT.EXT.NSIGMA.LIMIT} (3.0), the source is
+considered to be extended.
 
 Cosmic Rays are identified by a combination of the Kron magnitude and
@@ -1254,4 +1274,6 @@
 and the associated pixels are masked.
 
+\note{how are / were these parameters set?}
+
 \subsubsection{Non-PSF Sources}
 
@@ -1263,5 +1285,5 @@
 aperture) and working to a user defined S/N limit.  
 
-PSPhot will use the user-selected galaxy model to attempt the galaxy
+\code{psphot} will use the user-selected galaxy model to attempt the galaxy
 model fits.  In the configuration system, the keyword \code{GAL_MODEL}
 is set to the model of interest.  All suspected extended sources are
@@ -1294,5 +1316,5 @@
 After a first pass through the image, in which the brighter sources
 above a high threshold level have been detected, measured, and
-subtracted, PSPhot optionally begins a second pass at the image.  In
+subtracted, \code{psphot} optionally begins a second pass at the image.  In
 this stage, the new peaks are detected on the image with the bright
 sources subtracted.  In this pass, the peak detection process uses the
@@ -1310,9 +1332,8 @@
 \subsection{Extended Source Analysis}
 
-\note{intro paragraph: After the initial, fast analysis of the image
-  relying primarily on the PSF model, a complete analysis of the
-  extended source properties may be performed.  For PS1 processing,
-  this step is the nightly (PV0) analysis of individual exposures and
-  only performed for the stacks. }
+After the initial, fast analysis of the image relying primarily on the
+PSF model, a complete analysis of the extended source properties may
+be performed.  For PS1 processing, this step is the nightly (PV0)
+analysis of individual exposures and only performed for the stacks.
 
 \subsubsection{Radial Profiles}
@@ -1326,5 +1347,5 @@
 contour and radial profile.  
 
-In order to facilitate the Petrosian photometry analysis below, PSPhot
+In order to facilitate the Petrosian photometry analysis below, \code{psphot}
 generates a radial profile for each suspected galaxy.  This analysis
 starts by generating a radial profile in 24 azimuthal segments.  Near
@@ -1341,10 +1362,10 @@
 profiles, pairs of radial profiles from opposite sides of the source
 are compared.  Any masked values are replaced by the corresponding
-value in the other profile.  The minimum of both profiles is the kept
+value in the other profile.  The minimum of both profiles is then kept
 for both profiles.  The result of this analysis is a set of profiles
 of the form $f_i(r_i)$.  In this case, $f_i$ is effectively the
 surface brightness for each radius in instrumental counts per pixel.
 
-The surface brightness profiles are then used to define the radial
+The surface brightness profiles are then used to define the azimuthal
 contour at a specific isophotal level.  This contour will be used to
 rescale the radial profiles into a single set of profiles normalized
@@ -1366,11 +1387,10 @@
 (\code{RADIAL.ANNULAR.BINS.LOWER} \&
 \code{RADIAL.ANNULAR.BINS.UPPER}).  For each source, the resulting
-surface brightness profile is saved in the output cmf-file as an
-N-element value in the FITS table (\code{PROF_SB}).  The flux at each
-radial position and the fill-factor (fraction of pixels used to the
-total possible) as also saved as equal-length vectors in the FITS
-table (\code{PROF_FLUX} and \code{PROF_FILL}).  The values of the
-radial bins are saved in the cmf header (\code{RMIN_NN},
-\code{RMAX_NN}).
+surface brightness profile is saved in the output FITS table as a
+vector (\code{PROF_SB}).  The flux at each radial position and the
+fill-factor (fraction of pixels used to the total possible) are also
+saved as equal-length vectors in the FITS table (\code{PROF_FLUX} and
+\code{PROF_FILL}).  The values of the radial bins are saved in the
+output file FITS header (\code{RMIN_NN}, \code{RMAX_NN}).
 
 % \note{these profiles are not saved in PSPS}
@@ -1382,7 +1402,7 @@
 aperture which can be determined for galaxies without significant
 biases as a function of distance.  Since surface brightness in a
-resolved source is conserved, using a ratio of surface brightness to
-define a spatial scale results in a spatial scale which is constant
-regardless of galaxy distance.
+resolved source is conserved as a function of distance, using a ratio
+of surface brightness to define a spatial scale results in a spatial
+scale which is constant regardless of galaxy distance.
 
 To measure the Petrosian radius and flux, we start by defining a
@@ -1440,21 +1460,22 @@
 
 Preliminary Kron radius and flux values \citep{1980ApJS...43..305K}
-are calculated soon after sources are detected (Section~\ref{sec:moments}).
-However, these preliminary values are not accurate due to the
-window-functions applied.  After sources have been characterized and
-the PSF model is well-determined, the Kron parameters are
-re-calculated more carefully.  In this version of the calculation, the
-image is first smoothed by Gaussian kernel with $\sigma = 1.7$ pixels,
-corresponding to a FWHM of 1.0\arcsec\ in the PS1 stack images.  Next,
-the Kron radius is determined in an iterative process: the first
-radial moment is measured using the pixels in an aperture 6$\times$
-the first radial moment from the previous iteration.  On the first
-iteration, the sky radius is used in place of the first radial moment.
-By default, 2 iterations are performed.  The Kron radius is defined
-the be 2.5$\times$ the first radial moment.  The Kron flux is the sum
-of pixel fluxes within the Kron radius.  We also calculate the flux in
-two related annular apertures: the Kron inner flux is the sum of pixel
-values for the annulus $R_1 < r < 2.5 R_1$, while the Kron outer flux
-is the sum of pixel values for $2.5 R_1 < r < 4 R_1$.
+are calculated soon after sources are detected
+(Section~\ref{sec:moments}).  However, these preliminary values are
+not accurate due to the window-functions applied.  After sources have
+been characterized and the PSF model is well-determined, the Kron
+parameters are re-calculated more carefully.  In this version of the
+calculation, following the algorithm described by \cite{sextractor},
+the image is first smoothed by Gaussian kernel with $\sigma = 1.7$
+pixels, corresponding to a FWHM of 1.0\arcsec\ in the PS1 stack
+images.  Next, the Kron radius is determined in an iterative process:
+the first radial moment is measured using the pixels in an aperture
+6$\times$ the first radial moment from the previous iteration.  On the
+first iteration, the sky radius is used in place of the first radial
+moment.  By default, 2 iterations are performed.  The Kron radius is
+defined the be 2.5$\times$ the first radial moment.  The Kron flux is
+the sum of pixel fluxes within the Kron radius.  We also calculate the
+flux in two related annular apertures: the Kron inner flux is the sum
+of pixel values for the annulus $R_1 < r < 2.5 R_1$, while the Kron
+outer flux is the sum of pixel values for $2.5 R_1 < r < 4 R_1$.
 
 Two details in the calculation above should be noted.  First, for
@@ -1479,4 +1500,6 @@
 effect of reducing the impact of pixels which include flux from near
 neighbors.
+
+\note{give a test example}
 
 \subsubsection{Convolved Galaxy Model Fits}
@@ -1546,5 +1569,7 @@
 a function of the Sersic index.
 
-The PSF-convolved galaxy model fitting analysys uses the
+\note{special handling for central pixel}
+
+The PSF-convolved galaxy model fitting analysis uses the
 Levenberg-Marquardt minimization method to determine the best fit.  In this
 process, the $\chi^2$ value to be minimized is:
@@ -1587,5 +1612,5 @@
 The gradient vector and Hessian matrix are used in the
 Levenberg-Marquardt minimization analysis using the standard
-techinique of determining a step from the current set of model
+technique of determining a step from the current set of model
 parameters to a new set by solving the matrix equation:
 \[
@@ -1611,10 +1636,10 @@
 (galaxy model) image to be convolved before multiplying by the PSF
 model profile at that radial coordinate.  This approximation reduces
-the number of multiplications by a factor of near 8 for larger radii.
+the number of multiplications by a factor of \approx 8 for larger radii.
 For the small size of the PSF model used to convolve the galaxy model
 images, it was found that this direct convolution was faster than
 using an FFT-based convolution.
 
-% \note{(examples?)}
+\note{examples?  show timing comparisons?}
 
 For the Exponential and DeVaucouleur fits, all parameters are fitted
@@ -1639,6 +1664,6 @@
 previous, values in the grid above.  E.g., if the minimum fitted index
 value is 3.0, then the LMM fits are performed using $n$ = 2.5, 3.0, 3.5.
-The resulting $\chi^2$ values are then used to perform quadratid
-interpolation to find the index $n$ which produces the locally minium
+The resulting $\chi^2$ values are then used to perform quadratic
+interpolation to find the index $n$ which produces the locally minimum
 $\chi^2$ value.  Finally, this best-fit index value is held constant
 while Levenberg-Marquardt minimization is used to find the best fit
@@ -1653,9 +1678,9 @@
 important than an accurate total magnitude.  In the case of PS1, the image
 quality variations for stacks of different filters presents a serious
-challenge for the determination of precise colors.  PSPhot determines
+challenge for the determination of precise colors.  \code{psphot} determines
 a set of PSF-matched radial aperture flux measurements in order to
 minimize the impact of the stack image quality variations.
 
-In PSPhotStack, the stack analysis version of PSPhot, the 5 filter
+In \code{psphotStack}, the stack analysis version of \code{psphot}, the 5 filter
 images are processed together.  After the PSF models have been fitted
 and a best set of galaxy models have been determined, three sets of
@@ -1683,51 +1708,51 @@
 
 \subsection{Aperture Correction}
-
-The important concept here is that an analytical model will always
-fail to describe the flux of the sources at some level.  In the end,
-all astronomical photometry is in some sense a relative measurement
-between two images.  Whether the goal is calibration of a science
-image taken at one location to a standard star image at another
-location, or the goal is simply the repetitive photometry of the same
-star at the same location in the image, it is always necessary to
-compare the photometry between two images.  If this measurement is to
-be consistent, then the measurement must represent the flux of the
-stars in the same way regardless of the conditions under which the
-images were taken, at least within some range of normal image
-conditions.  So, for example, two images with different image quality,
-or with different tracking and focus errors, will have different PSF
-models.  Since an analytical model will always fail to represent the
-flux of the star at some level, the measured flux of the same source
-in the two images will be different (even assuming all other
-atmospheric and instrumental effects have been corrected).  The
-amplitude of the error will by determined by how inconsistently the
-models represent the actual source flux.  For example, if the first
-image PSF model flux is consistently 10\% too low and the second is 5\%
-too high, then the comparison between the two images will be in error
-by 15\%.  
-
-Aperture photometry avoids these problems, by trading for other
-difficulties.  In aperture photometry, if a large enough aperture is
-chosen, the amount of flux which is lost will be a small fraction of
-the total source flux.  Even more importantly, as the image conditions
-change, the amount lost will change by an even smaller fraction, at
-least for a large aperture.  This can be seen by the fact that the
-dominant variations in the image quality are in the focus, tracking
-and seeing.  All of these errors initially affect the cores of the
-stellar images, rather than the wide wings.  The wide wings are
-largely dominated by scattering in the optics and scattering in the
-atmosphere.  The amplitude and distribution of these two scattering
-functions do not change significantly or quickly for a single
-telescope and site.  
+\label{sec:aperture.correction}
+
+A PSF model will always fail to describe the flux of the stellar
+sources at some level.  For high-precision photometry, we need to be
+able to correct for the difference between the PSF model fluxes and
+the total flux of the sources.  In the end, all astronomical
+photometry is in some sense a relative measurement between two images.
+Whether the goal is calibration of a science image taken at one
+location to a standard star image at another location, or the goal is
+simply the repetitive photometry of the same star at the same location
+in the image, it is always necessary to compare the photometry between
+two images.  If this measurement is to be consistent, then the
+measurement must represent the flux of the stars in the same way
+regardless of the conditions under which the images were taken, at
+least within some range of normal image conditions.  So, for example,
+two images with different image quality, or with different tracking
+and focus errors, will have different PSF models.  Since an analytical
+model will always fail to represent the flux of the star at some
+level, the measured flux of the same source in the two images will be
+different (even assuming all other atmospheric and instrumental
+effects have been corrected).  The amplitude of the error will by
+determined by how inconsistently the models represent the actual
+source flux.  
+
+Aperture photometry attempts to avoid these problems, but introduces
+other difficulties.  In aperture photometry, if a large enough
+aperture is chosen, the amount of flux which is lost will be a small
+fraction of the total source flux.  Even more importantly, as the
+image conditions change, the amount lost will change by an even
+smaller fraction, at least for a large aperture.  This can be seen by
+the fact that the dominant variations in the image quality are in the
+focus, tracking and seeing.  All of these errors initially affect the
+cores of the stellar images, rather than the wide wings.  The wide
+wings are largely dominated by scattering in the optics and scattering
+in the atmosphere.  The amplitude and distribution of these two
+scattering functions do not change significantly or quickly for a
+single telescope and site.
 
 The difficulty for aperture photometry is the need to make an accurate
 measurement of the local background for each source.  As the aperture
 grows, errors in the measurement of the sky flux start to become
-dominant.  If the aperture is too small, then variation in the image
+dominant.  If the aperture is too small, then variations in the image
 quality are dominant.  The brighter is the source, the smaller is the
 error introduced by the large size of the aperture.  However, the
 number of very bright stars is limited in any image, and of course the
 brighter stars are more likely to suffer from non-linearity or
-saturation.  PSPhot measures the aperture correction ({\em ApResid})
+saturation.  \code{psphot} measures the aperture correction ({\em ApResid})
 for every PSF candidate source and applies this correction to the PSF
 model photometry.
@@ -1747,5 +1772,5 @@
 % magnitude}.
 
-%%% PSPhot measures the aperture correction ({\em ApResid}) for every PSF
+%%% \code{psphot} measures the aperture correction ({\em ApResid}) for every PSF
 %%% candidate source, then calculates the trend of this correction as a
 %%% function of the magnitude.  This trend is fitted with a line.  The
@@ -1762,8 +1787,8 @@
 %%% term.
 
-PSPhot allows a collection of PSF model functions to be tried on all
+\code{psphot} allows a collection of PSF model functions to be tried on all
 PSF candidate sources.  For each model test, the above corrected
 ApResid scatter is measured.  The PSF model function with the smallest
-value for the ApResid scatter is then used by PSPhot as the best PSF
+value for the ApResid scatter is then used by \code{psphot} as the best PSF
 model for this image.  The number of models to be tested is specified
 by the configuration keyword \code{PSF_MODEL_N}.  The configuration
@@ -1772,39 +1797,4 @@
 tested.
 
-Several likely PSF model classes are available within \code{psphot}:
-\begin{itemize}
-\item Gaussian : $f = I_0 e^{-z}$
-\item Pseudo-Gaussian : $f = I_0 (1 + z + \frac{1}{2} z^2 + \frac{1}{6} z^3)^{-1}$ \code{[PGAUSS]}
-\item Variable Power-Law : $f = I_0 (1 + z + z^{\alpha})^{-1}$ \code{[RGAUSS]}
-\item Steep Power-Law : $f = I_0 (1 + \kappa z + z^{2.25})^{-1}$ \code{[QGAUSS]}
-\item PS1 Power-Law : $f = I_0 (1 + \kappa z + z^{1.67})^{-1}$ \code{[PS1_V1]}
-\end{itemize}
-where $z \propto r^2$ ($z = \frac{x^2}{2\sigma_x^2} +
-\frac{y^2}{2\sigma_y^2} + \sigma_{\rm xy} x y $).  The Pseudo-Gaussian
-is a Taylor expansion of the Gaussian and is used by Dophot
-\citep{1993PASP..105.1342S}.  The latter profiles are similar to the
-Moffat profile form \citep{1969AA.....3..455M,1983AA...126..278B},
-with small differences.  For the PS1 GPC1 analysis, we used the
-\code{PS1_V1} model, which we found by experimentation to match well
-to the observed profiles generated by PS1.
-Figure~\ref{fig:radial.profiles} shows example radial profiles for
-moderately bright stars in fairly good (0.9 arcsec) and poor (2.2
-arcsec) seeing.  Using a fixed power-law exponent results in somewhat
-faster profile fitting compared to the variable power-law exponent
-model.
-
-% moffat : 1969A&A.....3..455M
-% buonanno : 1983A&AS...51...83B
-
-\begin{figure}[htbp]
-  \begin{center}
-  \includegraphics[width=\hsize]{{pics/radial.profiles}.\plotext}
-  \caption{\label{fig:radial.profiles} Radial profiles of stellar images from PS1.  These two
-    profiles illustrate the radial trend of the PS1 PSFs for a star
-    with FWHM 0.9 arcsec (red) and 2.2 arcsec (blue).  The black line
-    shows the PSF model with radial trend of the form $(1 + \kappa r^2 + r^{3.33})^{-1}$.}
-  \end{center}
-\end{figure}
-
 \subsection{Output Formats}
 
@@ -1823,5 +1813,5 @@
 For a difference image, both positive and negative sources will be
 present.  The basic peak detection algorithm will only trigger for the
-positive sources.  One solution is to simply apply PSPhot to both the
+positive sources.  One solution is to simply apply \code{psphot} to both the
 difference image and its negative value.  \note{do we want to code in
 an automatic switch to get both positive and negative excursions in
