Index: /trunk/doc/pslib/ChangeLogADD.tex =================================================================== --- /trunk/doc/pslib/ChangeLogADD.tex (revision 3771) +++ /trunk/doc/pslib/ChangeLogADD.tex (revision 3772) @@ -52,2 +52,10 @@ \item moved some sections to reflect order in SDRS (matrix, fftw) \end{itemize} + +\subsection{Changes from version 10 (19 April 2005) to version 11 (27 April 2005)} + +\begin{itemize} +\item fixed some typos in the definition of the rotation from CEO to GCRS (Eqn~\ref{CEOtoGCRS}). +\item added references to the SDRS APIs for the Earth Orientation section + +\end{itemize} Index: /trunk/doc/pslib/ChangeLogSDRS.tex =================================================================== --- /trunk/doc/pslib/ChangeLogSDRS.tex (revision 3771) +++ /trunk/doc/pslib/ChangeLogSDRS.tex (revision 3772) @@ -1,3 +1,3 @@ -%%% $Id: ChangeLogSDRS.tex,v 1.90 2005-04-25 21:20:41 price Exp $ +%%% $Id: ChangeLogSDRS.tex,v 1.91 2005-04-27 19:59:03 eugene Exp $ \subsection{Changes from version 00 to version 01} @@ -520,5 +520,9 @@ \item Moved Fixed Pattern out of Astronomical Images \end{itemize} - + \end{itemize} + +\subsection{Changes from Revision 13 (30 March 2005) to Revision 14 (27 April 2005)} + +\begin{itemize} \item Restrictions on the use of \code{malloc}, \code{calloc}, \code{realloc}, and \code{free} should not be unintentionaly imposed on 3rd party code. \item Add database support for ``auto-incrementing'' @@ -547,3 +551,12 @@ \item Removed pre-defined LM minimization functions; these will be defined in the Modules SDRS. \end{itemize} -\end{itemize} + +\item defined \code{psEarthPole}, re-cast Earth Orientation + Calculations to use it for inputs and outputs. +\item minor name changes in Earth Orientation to match ADD changes. +\item dropped TBD for \code{psFitsUpdateImage} +\item \code{psAberration} return value defined. +\item added \code{psArrayRemove} function (already exists in psArray.c) +\item added \code{psVectorExtend} function +\item changed inputs to \code{psImageSlice} to use \code{psRegion} +\end{itemize} Index: /trunk/doc/pslib/psLibADD.tex =================================================================== --- /trunk/doc/pslib/psLibADD.tex (revision 3771) +++ /trunk/doc/pslib/psLibADD.tex (revision 3772) @@ -1,3 +1,3 @@ -%%% $Id: psLibADD.tex,v 1.72 2005-04-19 23:44:43 eugene Exp $ +%%% $Id: psLibADD.tex,v 1.73 2005-04-27 19:59:04 eugene Exp $ \documentclass[panstarrs]{panstarrs} @@ -14,5 +14,5 @@ \project{Pan-STARRS Image Processing Pipeline} \organization{Institute for Astronomy} -\version{10} +\version{11} \docnumber{PSDC-430-006} @@ -41,4 +41,5 @@ 09 & 2005 Feb 14 & Frozen for Cycle 5 \\ \hline 10 & 2005 Apr 19 & Frozen for Cycle 6 \\ \hline +11 & 2005 Apr 27 & Update for Cycle 6 \\ \hline \RevisionsEnd @@ -1553,6 +1554,4 @@ the quarternion for this transformation. -\tbd{can we drop this, since we do this with the quaternion?} - The relevant trigonometric relationships are: % @@ -1611,7 +1610,7 @@ \phi_p & = & 90^\circ + 0^\circ.6406161\, T + 0^\circ.0003041\, T^2 + 0^\circ.0000051\, T^3 \end{eqnarray} -where $T$ is $($MJD$_{\rm out} -$ MJD$_{\rm in})/36525$ is the difference -between the two epochs, in Julian centuries. - +where $T$ is $($MJD$_{\rm out} -$ MJD$_{\rm in})/36525$ is the +difference between the two epochs, in Julian centuries. This +precession form shall be used to implement \code{PS_PRECESS_ROUGH}. \subsubsection{Suggested test cases} @@ -1648,11 +1647,10 @@ There are two reference implementatins for the code to account for the motion of the Earth in space. The first are the sample routines -provided by the IERS to accompany chaper 5 of IERS Bulletin 32. This -document and the code can be downloaded from -http://maia.usno.navy.mil/conv2003.html . The second reference -implementation is the SOFA software package managed by the IAU and -available at http://www.iau-sofa.rl.ac.uk Only the 2003-04-29 version -of SOFA should be used. The IERS code requires a few of the rotation -matrix utility routines from SOFA. +provided by the IERS to accompany chaper 5 of IERS Bulletin +32.\footnote{http://maia.usno.navy.mil/conv2003.html} The second +reference implementation is the SOFA software package managed by the +IAU.\footnote{http://www.iau-sofa.rl.ac.uk} Only the 2003-04-29 +version of SOFA should be considered. The IERS code requires a few of +the rotation matrix utility routines from SOFA. Both implementations are in FORTRAN 77. The SOFA code has a more @@ -1663,5 +1661,6 @@ reference for psLib should be the IERS code. Note that the IERS code calculates the transform from terrestrial to celestial coordinates, -while the SOFA code calculates its inverse. +while the SOFA code calculates its inverse. This code may be using as +a comparison for testing purposes. \subsubsection{Coordinate Systems} @@ -1711,13 +1710,13 @@ The X axes of the intermediate coordinate systems are known as the -Celestial and Terrestrial Ephemeris Origins. (CEO and TEO). Both are defined -to be non-rotating origins. A non-rotating origin is a point on the equator -whose instantaneous motion is always orthogonal to the equator -(Kaplan 2003 IAU XXV Joint Discussion 16 -\footnote{http://aa.usno.navy.mil/kaplan/NROs\%5BJD16proc\%5D.pdf}). -Thus the CEO is defined by its position in the GCRS at some epoch and by the -motion of the CIP in the GCRS since that date. Similarly the TEO is -defined by its position in the ITRS at some epoch and the motion of the -CIP in the ITRS since that date. +Celestial and Terrestrial Ephemeris Origins. (CEO and TEO). Both are +defined to be non-rotating origins. A non-rotating origin is a point +on the equator whose instantaneous motion is always orthogonal to the +equator (Kaplan 2003 IAU XXV Joint Discussion +16\footnote{http://aa.usno.navy.mil/kaplan/NROs\%5BJD16proc\%5D.pdf}). +Thus the CEO is defined by its position in the GCRS at some epoch and +by the motion of the CIP in the GCRS since that date. Similarly the +TEO is defined by its position in the ITRS at some epoch and the +motion of the CIP in the ITRS since that date. \subsubsection{ICRS - GCRS} @@ -1793,5 +1792,5 @@ This section is largely a summary of Chapter 5 of IERS Technical Note -32 \footnote{http://maia.usno.navy.mil/conv2003.html} (hereafter +32\footnote{http://maia.usno.navy.mil/conv2003.html} (hereafter IERS32), which is a description of the implementation of the Resoltions of the XXIVth General Assembly of the IAU, available from @@ -1807,4 +1806,6 @@ accurate to the 0.2 mas level. For higher accuracy the user must apply corrections to the model, which are tabulated by the IERS. + +\subparagraph{IAU 200A Precession/Nutation Model : {\tt psEOC\_PrecessionModel}} The IAU 2000A precession-nutation model may be calculated in the @@ -1855,7 +1856,7 @@ The constants $p_j$, $w_{i,j,k}$, $(a_{{\rm s},j})_i$, and $(a_{{\rm c},j})_i$ are given in the ASCII files: -tab5.2a.txt \footnote{http://maia.usno.navy.mil/conv2000/chapter5/tab5.2a.txt} (for $X$), -tab5.2b.txt \footnote{http://maia.usno.navy.mil/conv2000/chapter5/tab5.2b.txt} (for $Y$), and -tab5.2c.txt \footnote{http://maia.usno.navy.mil/conv2000/chapter5/tab5.2c.txt} (for $s+XY/2$). +tab5.2a.txt\footnote{http://maia.usno.navy.mil/conv2000/chapter5/tab5.2a.txt} (for $X$), +tab5.2b.txt\footnote{http://maia.usno.navy.mil/conv2000/chapter5/tab5.2b.txt} (for $Y$), and +tab5.2c.txt\footnote{http://maia.usno.navy.mil/conv2000/chapter5/tab5.2c.txt} (for $s+XY/2$). Note that the expansion is given for $s+XY/2$, since this series converges more rapidly than the one for $s$ alone. @@ -1876,10 +1877,12 @@ A FORTRAN reference implementation for the precession/nutation model -is available from the IERS -\footnote{http://maia.usno.navy.mil/conv2000/chapter5/XYS2000A.f}. -The psLib results should agree with the reference implementation to within -the limits of numerical precision. - -Next, corrections to $X$, and $Y$ may be obtained from the IERS as +is available from the +IERS.\footnote{http://maia.usno.navy.mil/conv2000/chapter5/XYS2000A.f} +The psLib results should agree with the reference implementation to +within the limits of numerical precision. + +\subparagraph{Corrections to the Model : {\tt psEOC\_PrecessionCorr}} + +Corrections to $X$, and $Y$ may be obtained from the IERS as part of Bulletin A, or B. It is recommended to use the values published daily by USNO in the table @@ -1895,24 +1898,28 @@ the result as instantaneous values. -The final step is to use $X$, $Y$, and $s$ to calculate the rotation -matrix from the CIP/CEO system to the GCRS using IERS32 equation (10), -reproduced below: - -\begin{equation} +\subparagraph{Spherical Rotation from Polar Coordinates : {\tt psSphereRot\_CEOtoGCRS}} + +In order to relate the values $X$, $Y$, and $s$ to the rotation +components, the rotation matrix below must be used. The definitions +of $X$, $Y$, and $s$ transform from the CIP/CEO system to the GCRS +using IERS32 equation (10), reproduced below: + +\begin{equation} +\label{CEOtoGCRS} \begin{pmatrix}1-aX^2& -aXY& X\cr -aXY& 1-aY^2& Y\cr -X& -Y& 1-a(X^2+Y^2)\cr \end{pmatrix} \cdot R_3(s), \end{equation} -where $R_3$ denotes a rotation about the Z axis, -$a = 1/(1+\sqrt{1 - X^2 + Y^2})$, -and $X$ and $Y$ are expressed in radians. -A FORTRAN reference implementation for this calculation is given -by the IERS \footnote{http://maia.usno.navy.mil/conv2000/chapter5/BPN2000.f}. - -Note that above we gave the expression for the transform toward celestial -coordinates (upward in figure X), in order to match the IERS reference code. -The inverse transform may be found by inverting the resulting rotation. - -\paragraph{Rotation of the Earth} +where $R_3$ denotes a rotation about the Z axis, $a = 1/(1+\sqrt{1 - +(X^2 + Y^2})$, and $X$ and $Y$ are expressed in radians. A FORTRAN +reference implementation for this calculation is given by the +IERS.\footnote{http://maia.usno.navy.mil/conv2000/chapter5/BPN2000.f} + +Note that above we gave the expression for the transform toward +celestial coordinates (upward in Figure~\ref{earthrot}), in order to +match the IERS reference code. The inverse transform may be found by +inverting the resulting rotation. + +\paragraph{Earth Rotation} The transform from the CIP/CEO to CIP/TEO coordinate systems is a @@ -1931,11 +1938,14 @@ motion''. Similarly to precession/nutation, the instantaneous position of the CIP in the ITRS is specified by the quantites $x_p$, and $y_p$, -and a third quantity, $s'$, gives the position of the TEO with respect -to the ITRS. The values of $x_p$ and $y_p$ are published daily by the -IERS\footnote{http://maia.usno.navy.mil/ser7/finals2000A.daily}, with +and a third quantity, $s'$, which give the position of the TEO with +respect to the ITRS. The values of $x_p$ and $y_p$ are published +daily by the +IERS,\footnote{http://maia.usno.navy.mil/ser7/finals2000A.daily} with a format described by their \code{readme.finals2000A}\footnote{http://maia.usno.navy.mil/ser7/readme.finals2000A}. The UT1$-$UTC, and the precession/nutation corrections (discussed elsewhere in this document) come from this same source. + +\subparagraph{Polar Motion from Bulletin : {\tt psEOC\_GetPolarMotion}} The polar motion coordinates should be interpolated using a third @@ -1953,5 +1963,7 @@ The tidal effects should be included by using the Ray tidal model given in IERS Gazette \#13. The definition of this correction is -provided below. +provided below (Section~\ref{Raymodel}). + +\subparagraph{Polar Motion Nutation Correction : {\tt psEOC\_NutationCorr}} By definition of the CIP, nutation terms with periods less than 2 days @@ -1966,4 +1978,6 @@ over this century by $s' = -4.7 \times 10^{-5} t$ in arcseconds. There is no need to apply short timescale corrections to $s'$. + +\subparagraph{Spherical Rotation from Polar Motion : {\tt psSphereRot\_ITRStoTEO}} The transform from the ITRS to the CIP/TEO frame can be constructed by @@ -2009,5 +2023,5 @@ correction from the Ray Tidal Model applied. -\subsubsection{Ray Tidal Model} +\subsubsection{Ray Tidal Model : {\tt psEOC\_PolarTideCorr}} The Ray Model tidal corrections to X, Y, and dT are given by the the Index: /trunk/doc/pslib/psLibSDRS.tex =================================================================== --- /trunk/doc/pslib/psLibSDRS.tex (revision 3771) +++ /trunk/doc/pslib/psLibSDRS.tex (revision 3772) @@ -1,3 +1,3 @@ -%%% $Id: psLibSDRS.tex,v 1.207 2005-04-25 21:20:41 price Exp $ +%%% $Id: psLibSDRS.tex,v 1.208 2005-04-27 19:59:04 eugene Exp $ \documentclass[panstarrs,spec]{panstarrs} @@ -11,5 +11,5 @@ \project{Pan-STARRS Image Processing Pipeline} \organization{Institute for Astronomy} -\version{13} +\version{14} \docnumber{PSDC-430-007} @@ -44,4 +44,6 @@ 11 & 2005 Jan 21 & draft for cycle 5 \\ \hline 12 & 2005 Feb 09 & final for cycle 5 \\ +13 & 2005 Mar 30 & draft for cycle 6 \\ +14 & 2005 Apr 27 & final for cycle 6 \\ \RevisionsEnd @@ -1374,4 +1376,33 @@ If the value of \code{vector} is \code{NULL}, then \code{psVectorRealloc} must return an error. + +\begin{verbatim} +psVector *psVectorExtend(psVector *vector, int delta, int nExtend); +\end{verbatim} + +This function increments \code{psVector.n}, the number of elements in +the vector by \code{nExtend}. If the current length of the vector +plus {\em twice} the number of new elements is greater than the +allocated space, an additional \code{delta} elements are allocated. +If the value of \code{delta} is less than 1, 10 shall be used. + +Here is an example of how \code{psVectorExtend} is used to +automatically increment the vector length. +\begin{verbatim} + // create data vector + psVector *y = psVectorAlloc (100); + y->n = 0; + for (int i = 0; i < 1000; i++) { + y->data.F32[y->n + 0] = 2*i; + y->data.F32[y->n + 1] = 2*i; + y->data.F32[y->n + 2] = 2*i; + psVectorExtend (y, 100, 3); + // increments n by 1, extends length if needed by 100 + } +\end{verbatim} +Note that the specification that the allocation always be greater than +the number of elements by twice the number of new elements implies +that there will be room on the next loop for \code{nExtend} new +elements, as in this example. \subsection{Simple Images} @@ -1492,4 +1523,13 @@ \code{delta} defines how many elements to add on each pass (if this value is less than 1, 10 shall be used). + +\begin{verbatim} +psBool psArrayRemove(psArray *array, psPtr value); +\end{verbatim} + +This function removes all entries of \code{value} in the \code{array}, +reducing the total number of elements of \code{array} as needed. +Returns \code{TRUE} if any elements were removed, otherwise +\code{FALSE}. \begin{verbatim} @@ -2679,7 +2719,12 @@ } psImageCutDirection; -psVector *psImageSlice(psVector *out, psVector *coords, const psImage *input, - const psImage *mask, unsigned int maskVal, int x0, int y0, - int x1, int y1, psImageCutDirection direction, const psStats *stats); +psVector *psImageSlice(psVector *out, + psVector *coords, + const psImage *input, + const psImage *mask, + unsigned int maskVal, + int x0, int y0, int x1, int y1, + psImageCutDirection direction, + const psStats *stats); \end{verbatim} Extract pixels from rectlinear region to a vector (array of floats). @@ -3916,8 +3961,8 @@ the conventions of the \code{psList} iterators. \begin{verbatim} -psListIterator *psMetadataIteratorAlloc(psMetadata *md, int location, bool mutable); +psListIterator *psMetadataIteratorAlloc(psMetadata *md, int location, const char *regex); bool psMetadataIteratorSet(psListIterator *iterator, int location); -psMetadataItem *psMetadataGetAndIncrement(psListIterator *iterator, const char *regex); -psMetadataItem *psMetadataGetAndDecrement(psListIterator *iterator, const char *regex); +psMetadataItem *psMetadataGetAndIncrement(psListIterator *iterator); +psMetadataItem *psMetadataGetAndDecrement(psListIterator *iterator); \end{verbatim} @@ -4509,5 +4554,4 @@ bool psFitsUpdateImage(psFits *fits, const psImage *input, psRegion region, int z); \end{verbatim} -\tbd{we have discussed this as the alternate name} Write an image section to the open \code{psFits} file pointer. This operation may write a portion of an image over the existing bytes of @@ -4609,5 +4653,5 @@ \begin{verbatim} -bool psFitsUpdateTable(psFits* fits, psMetadata *header, psMetadata* data, int row); +bool psFitsUpdateTable(psFits* fits, psMetadata* data, int row); \end{verbatim} Writes the \code{psMetadata} data to a FITS table at the specified row @@ -5364,4 +5408,9 @@ \tbd{supply the velocity as an un-normalized 3 vector?} +\tbd{MHPCC: please code this section as currently specified. We will + define a function, and algorithm, to return the current velocity + vector given a time and position, which can be fed to this + function}. + \paragraph{Aberration} The following function calculates the \code{apparent} position of a @@ -5369,12 +5418,25 @@ observer, represented as a speed and a direction: \begin{verbatim} -psAberration(psSphere *apparent, psSphere *actual, psSphere direction, double speed); +psSphere *psAberration(psSphere *apparent, psSphere *actual, psSphere direction, double speed); \end{verbatim} The \code{actual} and \code{apparent} positions are represented as \code{psSphere} entries, as is the \code{direction} of motion. The -speed in that direction is given in units of the speed of light. +speed in that direction is given in units of the speed of light. If +the value of \code{apparent} is NULL, a new \code{psSphere} is +allocated, otherwise the point to \code{apparent} is used for the +result. \paragraph{Gravitational Deflection} +The following function calculates the \code{apparent} position of a +star, given its \code{actual} position and the position of the sun: +\begin{verbatim} +psSphere *psGravityDeflection(psSphere *apparent, psSphere *actual, psSphere *sun); +\end{verbatim} +The \code{actual} and \code{apparent} positions are represented as +\code{psSphere} entries, as is position of the sun. If the value of +\code{apparent} is NULL, a new \code{psSphere} is allocated, otherwise +the point to \code{apparent} is used for the result. + \paragraph{Parallax} @@ -5385,4 +5447,18 @@ \subsubsection{Transformation from GCRS to ITRS} + +The following functions calculate the components, $X$, $Y$, and $s$, +representing the location of the earth's pole at any moment, or they +determine the velocity of the pole $X'$, $Y'$, $s'$. We use the +following structure to carry the polar coordinate information. This +representation may be converted to a rotation between the frames. + +\begin{verbatim} +typedef struct { + double x; + double y; + double s; +} psEarthPole; +\end{verbatim} \paragraph{Precession/Nutation} @@ -5393,5 +5469,5 @@ % \begin{verbatim} -psSphere *psEOC_PrecessionModel(double *s, const psTime *time) +psEarthPole *psEOC_PrecessionModel(const psTime *time) \end{verbatim} % @@ -5401,20 +5477,20 @@ machine accuracy. -The following function provides interpolated corrections to $X$ and -$Y$ from the tables provided by the IERS, just as it does for UT1 and -polar motion. - -\begin{verbatim} -psSphere *psEOC_GetPolarCorr(const psTime *time, psTimeBulletin bulletin); +The following function provides interpolated corrections to the $X$ +and $Y$ components of the polar coordinates from the tables provided +by the IERS, just as it does for UT1 and polar motion. + +\begin{verbatim} +psEarthPole *psEOC_PrecessionCorr(const psTime *time, psTimeBulletin bulletin); \end{verbatim} The polar correction is applied to the $X$ and $Y$ elements of the rotation to provide higher accuracy. The spherical rotation term is -generated by providing the three elements of the rotation to the -following function: -\begin{verbatim} -psSphereRot *psSphereRot_CEOtoGCRS(double s, const psSphere *pole) -\end{verbatim} -The retulting \code{psSphereRot} may be used to determine the rotation +generated by providing the polar coordinate to the following function: +\begin{verbatim} +psSphereRot *psSphereRot_CEOtoGCRS(const psEarthPole *pole) +\end{verbatim} +This function constructs the rotation element as described in the ADD ( +The resulting \code{psSphereRot} may be used to determine the rotation from CIP/CEO to GCRS. This function must give results identical to the IERS BPN2000, within the limits of machine accuracy. @@ -5434,5 +5510,5 @@ motion components, $x_p$ and $y_p$, extracted from the IERS tables. \begin{verbatim} -psSphere *psEOC_GetPoleCoords(const psTime *time, psTimeBulletin bulletin); +psEarthPole *psEOC_GetPolarMotion(const psTime *time, psTimeBulletin bulletin); \end{verbatim} @@ -5441,20 +5517,21 @@ ADD). \begin{verbatim} -psSphere *psEOC_TidePolarCorr(const psTime *time); +psEarthPole *psEOC_PolarTideCorr(const psTime *time); \end{verbatim} The following function provides the additional corrections due to nutation -terms with periods less than or equal to two days: -\begin{verbatim} -psSphere *psEOC_NutationCorr(psTime *time); -\end{verbatim} - -The following function should generate the \code{psSphereRot} transform from -ITRS to CIP/TEO: -\begin{verbatim} -psSphereRot *psSphereRot_ITRStoTEO(psSphere pole, psTime *time); -\end{verbatim} -The time argument should be used to internally calculate $s'$. -This function should give identical results to the IERS POM2000 subroutine. +terms with periods less than or equal to two days, as well as the +correction to the $s'$ component of the polar motion: +\begin{verbatim} +psEarthPole *psEOC_NutationCorr(psTime *time); +\end{verbatim} + +The following function converts the polar motion corrections to a +spherical rotation using the prescription in the ADD: +\begin{verbatim} +psSphereRot *psSphereRot_ITRStoTEO(const psEarthPole *motion); +\end{verbatim} +This function should give identical results to the IERS POM2000 +subroutine. \subsubsection{Earth Orientation Wrappers}