Index: /trunk/doc/release.2015/ps1.analysis/analysis.tex
===================================================================
--- /trunk/doc/release.2015/ps1.analysis/analysis.tex	(revision 40589)
+++ /trunk/doc/release.2015/ps1.analysis/analysis.tex	(revision 40590)
@@ -121,5 +121,5 @@
 \note{the beginning section needs to be updated to mention the DR1 and
   DR2 releases, the PV0-PV3 analysis versions, and the basic idea of
-  the IPP stages).
+  the IPP stages}.
 
 This is the fourth in a series of seven papers describing the
@@ -1734,5 +1734,5 @@
 \item Exponential profile : $f = I_0 e^{-\rho}$
 \item DeVaucouleur profile \citep{1948AnAp...11..247D}: $f = I_0 e^{-\rho^{1/4}}$
-\item Sersic \citep{1963BAAA....6...41S} : $f = I_0 e^{-\rho^{1/n}}$
+\item S\'ersic \citep{1963BAAA....6...41S} : $f = I_0 e^{-\rho^{1/n}}$
 \end{itemize}
 where $\rho$ is a normalized radial term: $\rho =
@@ -1743,6 +1743,6 @@
 x_0, Y_{\rm chip} - y_0$).  Including the normalization ($I_0$) and a
 local sky value, the Exponential and DeVaucouleur profiles have 7 free
-parameters and the Sersic profile has the additional free parameter of
-the Sersic index $n$.  In this stage, the galaxy model is convolved
+parameters and the S\'ersic profile has the additional free parameter of
+the S\'ersic index $n$.  In this stage, the galaxy model is convolved
 with an approximation to our best guess for the PSF model at the
 location of the galaxy.
@@ -1772,6 +1772,6 @@
 quantive the relationships between the first radial moment used to
 calculated a Kron Magnitude and the effective radius for different
-Sersic index values, $n$.  Since the Exponential and DeVaucouleur
-models are equivalent to Sersic models with $n$ = 1 and 4,
+S\'ersic index values, $n$.  Since the Exponential and DeVaucouleur
+models are equivalent to S\'ersic models with $n$ = 1 and 4,
 respectively, this work can be used to generate the initial effective
 radius values for all 3 model types.  Once the effective radius is
@@ -1780,9 +1780,7 @@
 generate a guess for the normalization, applying an appropriate scale
 factor based on the ($R_{xx}$, $R_{yy}$ , $R_{xy}$) values, generated
-by integrating normalized Sersic models and determining the
+by integrating normalized S\'ersic models and determining the
 relationship between the central intensity and the integrated flux as
-a function of the Sersic index.
-
-\note{special handling for central pixel}
+a function of the S\'ersic index.
 
 The PSF-convolved galaxy model fitting analysis uses the
@@ -1859,10 +1857,10 @@
 
 For the Exponential and DeVaucouleur fits, all parameters are fitted
-in the non-linear minimization stage.  For the Sersic model, we do not
+in the non-linear minimization stage.  For the S\'ersic model, we do not
 fit the index within the Levenberg-Marquardt analysis.  Instead, we
 start with a coarse grid search over a range of possible index values,
 ($n = 0.5, 1.0, 1.5, 2.0, 3.0, 4.0, 5.0, 6.0$) and a range of possible
 values for $R_{\rm eff}$ based on the value of $R_1$, the first radial
-moment.  For a given value of the Sersic index, the $R_{\rm eff}$ is
+moment.  For a given value of the S\'ersic index, the $R_{\rm eff}$ is
 related to the 1st radial moment by the scale factor specificy by
 Graham \& Driver.  We use the observed value of the 1st radial moment
@@ -1874,5 +1872,5 @@
 We next perform 3 Levenberg-Marquardt minimization fits allowing the
 shape parameters ($R_{xx}$, $R_{yy}$ , $R_{xy}$) and the normalization
-to be fitted, holding the centroid ($x_0, y_0$), Sersic index $n$, and
+to be fitted, holding the centroid ($x_0, y_0$), S\'ersic index $n$, and
 sky constant.  In these fits, the index $n$ is set to the minimum
 value previously calculated as well as values halfway to the next, and
@@ -1887,4 +1885,42 @@
 % Graham & Driver : Graham A. W., Driver S. P.  2005, PASA 22, 118
 % DOI: https://doi.org/10.1071/AS05001
+
+The central pixel of the S\'ersic, DeVaucouleur, and Exponential
+models require special handling.  When comparing an analytical model
+to the pixelized image recorded by a CCD, one normally treats the
+value in a pixel as equivalent to the value of the model at the center
+of the pixel.  However, in reality, the number of counts observed in a
+pixel represents the integral of the surface brightness across the
+area of the pixel.  This average will be equal to the central surface
+brightness times the area of a pixel as long as the second and higher
+derivatives of the analytical model are zero.  However, if the first
+and second derivatives are non-zero, the curvature of the function
+within the pixel will make the integral differ from the central
+surface brightness times a fixed pixel area.  If the curvature of the
+model function is sufficiently large, this difference will have a
+significant impact on the evaluation of the model.   This situation is
+particularly true for the central portion of the S\'ersic-like model
+functions. 
+
+%% this can be seen by writing the taylor expansion of the function
+%% about the center of the pixel.  do this?
+
+In order to accurately compare the observed galaxy flux distribution
+to a model, it is necessary to integrate the pixel flux for a given
+set of model parameter values.  This could be done in a numerical
+fashion, but in practice brute-force evaluation of the numerical
+integral is computationally very expensive, requiring many evaluations
+of the model function.  Within \ippprog{psphot}, we bypass this
+problem by defining a set of pre-calculated images for the central 9
+pixels (the $3 \times 3$ grid of pixels centered on the peak).  These
+pixel images are defined at higher resolution, with 11 subpixels per
+real CCD pixel.  The pre-calculated images are generated for a series
+of values for the following parameters: S\'ersic index, effective
+radius, axial ratio.  We then select the closest image to our specific
+case, and integrate over the true sub-pixels relevant for our position
+and model.  We have thus turned the problem from thousands of
+evaluations of the full galaxy model to \approx 100 straight
+additions, or up to $6 \times$ that number if we interpolate between
+any of the parameters.
 
 \subsubsection{Convolved Radial Aperture Photometry}
@@ -2109,11 +2145,9 @@
 \section{Forced Photometry Modes}
 
-\note{edit this section to remove references to the IPP stages; just refer to the psphot concepts}
-
 Traditionally, projects which use multiple exposures to increase the
 depth and sensitivity of the observations have generated something
-equivalent to the \ippstage{stack} images produced by the IPP analysis
+equivalent to the stack images produced by the IPP analysis
 (c.f, CFHT Legacy survey, COSMOS, etc).  In theory, the photometry of
-the \ippstage{stack} images produces the ``best'' photometry catalog,
+the stack images produces the ``best'' photometry catalog,
 with best sensitivity and the best data quality at all magnitudes.  In
 practice, these images have some significant limitations due to the
@@ -2130,25 +2164,25 @@
 that point.  Because of the high mask fraction, the exposures which
 contributed to pixels at one location may be somewhat different just a
-few tens of pixels away.  In the end, the \ippstage{stack} images have
+few tens of pixels away.  In the end, the stack images have
 a effective point spread function which is not just variable, but
 changing significantly on small scales in a highly textured fashion.
 
 Any measurement which relies on a good knowledge of the PSF at the
-location of an object either needs to determine the PSF variations
-present in the \ippstage{stack} image or the measurement will be
-somewhat degraded.  The highly textured PSF variations make this a
-very challenging problem: not only would such a PSF model require an
-unusually fine-grained PSF model, there would likely not be enough PSF
-stars in a given \ippstage{stack} image to determine the model at the
-resolution required.  The IPP photometry analysis code uses a PSF
-model with 2D variations using a grid of at most $6\times 6$ samples
-per skycell, a number reasonably well-matched to the density of stars
-at most moderate Galactic latitudes.  This scale is far too large to
-track the fine-grained changes apparent in the stack images.
-
-Thus PSF photometry as well as convolved galaxy models in the stack
-are degraded by the PSF variations.  Aperture-like measurements are in
-general not as affected by the PSF variations, as long as the aperture
-in question is large compared to the FWHM of the PSF.
+location of an object needs to determine the PSF variations present in
+the stack image or the measurement will be somewhat degraded.  The
+highly textured PSF variations make this a very challenging problem:
+not only would such a PSF model require an unusually fine-grained PSF
+model, there would likely not be enough PSF stars in a given stack
+image to determine the model at the resolution required.  The IPP
+photometry analysis code uses a PSF model with 2D variations using a
+grid of at most $6\times 6$ samples per skycell, a number reasonably
+well-matched to the density of stars at most moderate Galactic
+latitudes.  This scale is far too large to track the fine-grained
+changes apparent in the stack images.
+
+As a result, PSF photometry as well as convolved galaxy models in the
+stack are degraded by the PSF variations.  Aperture-like measurements
+are in general not as affected by the PSF variations, as long as the
+aperture in question is large compared to the FWHM of the PSF.
 
 %% The IPP team initially explored the option of convolving each input
@@ -2158,6 +2192,7 @@
 The IPP analysis solves this problem by starting with the sources
 detected in the stack images and performing forced photometry on the
-individual warp images used to generate the stack.  This
-forced-photometry analysis is performed using the
+individual warp images used to generate the stack, and then combining
+the resulting measurements to determine a high-quality average value.
+This forced-photometry analysis is performed using the
 \ippprog{psphotFullForce} variant of \ippprog{psphot}.
 
@@ -2176,6 +2211,6 @@
 image; the measured flux may even be negative due to statistical
 fluctuations.  When combined together, these low-significance
-measurements will result in a signficant measurement as the
-signal-to-noise increases with the combination of more data.
+measurements result in a signficant measurement as the signal-to-noise
+increases with the combination of more data.
 
 Individual warp images are processed independently with separate
@@ -2192,44 +2227,43 @@
 \label{sec:galaxy.forced.fit}
 
-The convolved galaxy models are also re-measured on the
-\ippstage{warp} images by the \ippstage{fullforce} stage analysis.  In
-this analysis, the galaxy models determined by the
-\ippstage{staticsky} photometry analysis are used to seed the analysis
-in the individual \ippstage{warp} images.  The motivation of this
-analysis is the same as the \ippstage{fullforce} PSF photometry: the
-PSF of the \ippstage{stack} image is poorly determined due to the
-masking and PSF variations in the inputs.  Without a good PSF model,
-the PSF-convolved galaxy models are of limited accuracy.
-
-In the \ippstage{fullforce} galaxy model analysis, we assume that the
-galaxy position and position angle, along with the Sersic index if
-appropriate, have been sufficiently well determined in the
-\ippstage{staticsky} analysis.  In this case, the goal is to determine
-the best values for the major and minor axis of the elliptical contour
-and at the same time the best normalization corresponding to the best
+The convolved galaxy models are also re-measured on the warp images by
+the \ippprog{psphotFullForce} analysis.  In this analysis, the galaxy
+models determined from the stack image analysis are used to seed the
+analysis in the individual warp images.  The motivation of this
+analysis is the same as the forced PSF photometry: the PSF of the
+stack image is poorly determined due to the masking and PSF variations
+in the inputs.  Without a good PSF model, the PSF-convolved galaxy
+models are of limited accuracy.
+
+In the forced galaxy model analysis, we assume that the galaxy
+position and position angle, along with the S\'ersic index if
+appropriate, have been sufficiently well determined in the analysis of
+the stack image.  In this case, the goal is to determine the best
+values for the major and minor axis of the elliptical contour and at
+the same time the best normalization corresponding to the best
 elliptical shape, and thus the best galaxy magnitude value.
 
-For each \ippstage{warp} image, the \ippstage{staticsky} values for
-the major and minor axis are used as the center of a $5 \times 5$ grid
-search of the major and minor axis parameter values.  The grid spacing
-is defined as a function of the signal-to-noise of the galaxy in the
-stack image so that bright galaxies are measured with a much finer
-grid spacing than faint galaxies.  For both the major and minor axis
-directions, values of ($1 - \frac{3.0}{S/N}$, $1 - \frac{1.5}{S/N}$,
-1.0, $1 + \frac{1.5}{S/N}$, $1 + \frac{3.0}{S/N}$) times the dimension
-are tested.  For each grid point, the major and minor axis values at
-that point are used to generate the model.  The model is then
-convolved with the PSF model for the \ippstage{warp} image at that
-point.  The resulting convolved model is then compared to the
-\ippstage{warp} pixel data values and the best fit normalization value
-is determined.  The integrated flux, flux error, and the $\chi^2$
-value for each grid point are recorded.
+For each warp image, the stack values for the major and minor axis are
+used as the center of a grid search of the major and minor axis
+parameter values.  The grid spacing is defined as a function of the
+signal-to-noise of the galaxy in the stack image so that bright
+galaxies are measured with a much finer grid spacing than faint
+galaxies.  For the PV3 $3\pi$ analysis, a $5 \times 5$ grid was used;
+values in both the major and minor axis directions of ($1 -
+\frac{3.0}{S/N}$, $1 - \frac{1.5}{S/N}$, 1.0, $1 + \frac{1.5}{S/N}$,
+$1 + \frac{3.0}{S/N}$) times the dimension are tested.  For each grid
+point, the major and minor axis values at that point are used to
+generate the model.  The model is then convolved with the PSF model
+for the warp image at that point.  The resulting convolved model is
+then compared to the warp pixel data values and the best fit
+normalization value is determined.  The integrated flux, flux error,
+and the $\chi^2$ value for each grid point are recorded.
 
 For a given galaxy, the result is a collection of $\chi^2$ values,
 fluxes, and flux errors for each of the grid points spanning all
-\ippstage{warp} images.  A single $\chi^2$ grid can then be made by
+warp images.  A single $\chi^2$ grid can then be made by
 combining each grid point across the inputs.  The combined $\chi^2$
 for a single grid point is simply the sum of all $\chi^2$ values at
-that point.  If, for a single \ippstage{warp} image, the galaxy model
+that point.  If, for a single warp image, the galaxy model
 is excessively masked, then that image will be dropped for all grid
 points for that galaxy.  The reduced $\chi^2$ values can be determined
@@ -2240,66 +2274,96 @@
 axis values for the interpolated minimum $\chi^2$ value.  The errors
 on these two parameters is then found by determining the contour at
-which the \note{reduced?} $\chi^2$ increases by 1.
-
-In this way, the \ippstage{fullforce} galaxy analysis uses the PSF
-information from each \ippstage{warp} to determine a best set of
-convolved galaxy models for each object in the \ippstage{skycal}
-catalog.
+which the $\chi^2$ increases by 1.
+
+In this way, the forced galaxy model analysis uses the PSF information
+from each warp image to determine a best set of convolved galaxy
+models for each galaxy model measured for the stack image.
 
 \section{Difference Image Photometry}
 
-\note{need an intro paragraph or so}
-
-The variance map for a difference image must be generated from the two
-images used to construct the difference.  Otherwise, the low sky level
-will automatically result in inconsistent interpretation of the variance.
+Among the primary science drivers for Pan-STARRS are the detection of
+moving objects (e.g., asteroids) and explosive transient sources
+(e.g., supernovae).  For both of these situations, difference images
+are commonly used to remove the clutter of the static stars and
+galaxies.  In the Pan-STARRS system, difference images are generated
+using the PSF-matching technique described by
+\citep[e.g.,][]{1998ApJ...503..325A}.  The description of the
+Pan-STARRS implementation is given by \cite{price2017}.  The analysis
+of the sources detected in these difference images uses a portion of
+the \ippprog{psphot} code embedded in the program, \ippprog{ppSub},
+which generates those image.  
+
+The analysis of the difference image follows the same basic steps as
+other \ippprog{psphot} versions with some minor modifictions (see
+Table~\ref{tab:measurements}), as follows.  The background subtraction
+is performed before the PSF matching and image subtraction is
+performed.  The PSF model construction stage is not possible in the
+difference image due to the lack of valid sources.  Instead, the PSF
+model from is generated from the positive image, after PSF-matching
+but before the subtraction is performed.  Because we do not expect to
+have a large number of sources, only a single source detection pass is
+performed, and at the lowest signal-to-noise threshold.  Only linear
+PSF model fitting is performed using the centroid determined from the
+analysis of the source moments.  
+
+For the difference images, the galaxy model analysis is not relevant.
+In a properly-constructed difference image, galaxies are unlikely to
+remain behind as significant sources.  Most real sources in the
+difference image will be PSF-like and will consist of photometrically
+variable sources (flare stars, supernovae, etc) or astrometrically
+variable sources (high-proper motion stars or solar-system bodies).
+There are three likely classes of sources which will not be well
+represented by the PSF model, as discussed below.
+
+Fast-moving solar-system objects will appear as short streaks.  For
+example, a fast solar system object may have an apparent rate of 0.5
+degrees per hour, translating to 15 arcseconds in a 30 second
+exposure.  Even a main belt asteroid at roughly 1 AU has reflex motion
+of approximately 1 degree per day, equivalent to 1.25 arcsec in a 30
+second exposure, and may be noticeably smeared and non-PSF-like.  In
+\ippprog{psphot}, we use a trailed-star model to characterize these
+types of sources.  This model is fitted in the same portion of the
+code which performs the unconvolved galaxy model analysis.
+\note{describe the trailed analytical model}.
+
+In some cases, the stars in the two images may be somewhat offset.
+For specific stars, this offset may be due to differential chromatic
+aberration from the atmosphere or the optics, or from modest proper
+motion.  If the astrometric solution for one of the two images is
+insufficiently accurate, all stars in large portions of the images may
+be noticably displaced.  In both of these situations, the stars will
+appear as PSF dipoles in the difference images.  The positive and the
+negative images will have stellar profiles, but they will be offset
+and will not subtract well.  The two components may not have the same
+amplitude.  In theory, a PSF-dipole model could be used to fit these types of
+sources, with free parameters of the two centroids and the two
+fluxes.  In practice in \ippprog{psphot}, we use a number of non-parametric
+pixel-level statistics in an attempt to detect these cases.  
+
+\note{list the parameters}
+
+Comets appear in the difference images as a non-PSF sources.  Their
+2-D structure includes both the flux from the coma (with a typical
+power-law profile) and flux from the tail (with a more complex flux
+distribution).  We use the Kron magnitudes to identify possibly
+extended objects which may be cometary in nature.  \note{need some
+  info from MOPS folks on what is used}
 
 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 \code{psphot} to
-both the difference image and its negative value.
-
-In the case of a difference image, the PSF model construction stage
-will probably fail for lack of valid sources.  It is better in these
-cases to provide PSF model from some other source.  For example, the
-two images which are combined to generate the difference image
-represent the PSF.  Presumably, one or both have been convolved with a
-PSF-matching kernel.  The images which result from the convolution
-should be used to measure the PSF model.  \note{this is what we
-  actually do, so remove hypothetical wording.}
-
-The source classification scheme defaults to the galaxy models for
-sources which are not well represented by the PSF model.  In a
-properly-constructed difference image, galaxies are unlikely to remain
-behind as significant sources.  Most real sources in the difference
-image will be PSF-like and will consist of photometrically variable
-sources (flare stars, supernovae, etc) or astrometrically variable
-sources (high-proper motion stars or solar-system bodies).  There are
-three likely classes of sources which will not be well represented by
-the PSF model.  1) Fast-moving solar-system objects will appear as
-short streaks.  For example, a fast solar system object would have an
-apparent rate of 0.5 degrees per hour, translating to 15 arcseconds in
-a 30 second exposure.  Even a main belt asteroid at roughly 1 AU would
-have reflect motion of approximately 1 degree per day, equivalent to
-1.25 arcsec in a 30 second exposure, and could be noticeably smeared
-and non-PSF-like.  A trailed-star model can be used to characterize
-these types of sourcess.  2) Small offset stars, either due to
-atmospheric / color effects or modest proper motion will appear as PSF
-dipoles in the difference images.  The positive and the negative
-images will have stellar profiles, but they will be significantly
-offset and will not subtract well.  The two components may not have
-the same amplitude.  A PSF-dipole model can be used to fit these types
-of sources, with free parameters of the two centroids and the two
-fluxes.  3) Comets will appear in the difference images as a non-PSF
-sources.  Their 2-D structure includes both the flux from the coma
-(with a typical power-law profile) and flux from the tail (with a more
-complex flux distribution).  A comet flux model can be used to
-characterize these sources in difference images.  A major difficulty
-in applying these three types of models is in making a robust test of
-which model should be used.  This problem is akin to the issue of
-selecting and distinguishing between multiple galaxy models, as
-discussed in the section on Galaxy models.
+positive sources.  In the \ippprog{ppSub} program, both the $A - B$
+and the $B - A$ images are sent to the \ippprog{psphot} routine for
+source detection and characterization.
+
+Note that the variance image for a difference image must be generated
+from the two positive images used to construct the difference.  It is
+possible to run \ippprog{psphot} as an external program on a
+difference image generated previously.  In this case, the variance
+image and the PSF model must be supplied as well as the difference
+image.
 
 \section{Examples and Tests}
+
+\note{to be added}
 
 \acknowledgments