IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 9785


Ignore:
Timestamp:
Oct 30, 2006, 11:45:46 AM (20 years ago)
Author:
eugene
Message:

added full derivation of ellipse rotations

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/doc/pslib/psLibADD.tex

    r6009 r9785  
    1 %%% $Id: psLibADD.tex,v 1.93 2006-01-16 01:11:40 eugene Exp $
     1%%% $Id: psLibADD.tex,v 1.94 2006-10-30 21:45:46 eugene Exp $
    22\documentclass[panstarrs]{panstarrs}
    33
     
    366366size is choosen to keep about 25 points per bin.  Thus, the bin size
    367367($d\sigma$) is set to about:
    368 \[
     368\begin{equation}
    369369d\sigma = 2.35 \sigma (N_{\rm bin}/N_{50}) = (25 \times 2.35) (\sigma/N_{50})
    370 \]
     370\end{equation}
    371371With the limitation that $\sigma/d\sigma$ should be limited on one end
    372372to the value 1, and the other to the value 5.  The easiest way to set
    373373this limit is to define dN to be:
    374 \[
     374\begin{equation}
    375375dN = (\sigma / d\sigma) = (N_{50} / N_{\rm bin}) / 2.35 = 0.017 * N_{50}
    376 \]
     376\end{equation}
    377377\subsubsection{Histograms}
    378378
     
    426426T_2(x) & = & 2x^2 - 1 \\
    427427T_3(x) & = & 4x^3 - 3x \\
    428 T_4(x) & = & 8x^4 - 8x^2 + 1 \\
     428T_4(x) & = & 8x^4 - 8x^2 + 1
    429429\end{eqnarray}
    430430Chebyshev polynomials follow the recurrence relation:
     
    437437\begin{eqnarray}
    438438d_j  & = & 2xd_{j+1} - d_{j+2} + c_j \\
    439 f(x) & = & x*d_1 - d_2 + 1/2 c_0 \\
     439f(x) & = & x*d_1 - d_2 + 1/2 c_0
    440440\end{eqnarray}
    441441
     
    721721$i+1,j$, $i,j+1$, and $i+1,j+1$ with pixel values $V_{0,0}$,
    722722$V_{1,0}$, $V_{0,1}$, $V_{1,1}$.  The value at $x,y$ is given by:
    723 \[ V = (V_{0,0}(1 - f_x) + V_{1,0}f_x)(1 - f_y) + (V_{0,1}(1-f_x) + V_{1,1}f_x)f_y \]
     723\begin{equation} V = (V_{0,0}(1 - f_x) + V_{1,0}f_x)(1 - f_y) + (V_{0,1}(1-f_x) + V_{1,1}f_x)f_y \end{equation}
    724724This expression is more efficiently evaluated by factoring and
    725725calculating the expresion as:
    726 \[ r_x = V_{0,0} + (V_{1,0} - V_{0,0})f_x \]
    727 \[ V = r_x + (V_{0,1} + (V_{1,1} - V_{0,1})f_x - r_x)f_y \]
     726\begin{equation} r_x = V_{0,0} + (V_{1,0} - V_{0,0})f_x \end{equation}
     727\begin{equation} V = r_x + (V_{0,1} + (V_{1,1} - V_{0,1})f_x - r_x)f_y \end{equation}
    728728
    729729Note that the values of $f_x$ and $f_y$ require some care.  Given a
     
    876876the output image ($i,j$) corresponds to a fractional pixel coordinate
    877877($x,y$) in the input image according to:
    878 \[ x = (i - i_o)*\cos\theta + (j - j_o)*\sin\theta \]
    879 \[ y = (i_o - i)*\sin\theta + (j - j_o)*\cos\theta \]
     878\begin{equation} x = (i - i_o)*\cos\theta + (j - j_o)*\sin\theta \end{equation}
     879\begin{equation} y = (i_o - i)*\sin\theta + (j - j_o)*\cos\theta \end{equation}
    880880where the offset coordinate ($i_o,j_o$) depends on the sign of the
    881881sine of the angle $\theta$.  If the sign of that sine is positive, the
     
    10921092P_0     & = & \left| C_0 \right|^2 / N^2 \\
    10931093P_j     & = & \left( \left| C_j \right|^2 + \left| C_{N-j} \right|^2 \right)/ N^2 \\
    1094 P_{N/2} & = & \left| C_{N/2} \right|^2 / N^2 \\
     1094P_{N/2} & = & \left| C_{N/2} \right|^2 / N^2
    10951095\end{eqnarray}
    10961096where $j = 1, 2, \ldots, (N/2 - 1)$.
     
    10981098Note that we leave the issue of ``windowing'' the data up to the
    10991099caller, and choose to normalise by $1/N^2$.
     1100
     1101\subsection{Ellipse Representations}
     1102
     1103Images of astronomical objects may often be represented using a model
     1104consisting of a radial profile combined with an elliptical contour.
     1105Two common ways to measure such a shape are to fit a model to the
     1106light distribution or to measure second-order moments, perhaps with
     1107some weighting profile.  In the special case of a 2D Gaussian with an
     1108elliptical contour, these representations are equivalent.  The
     1109following discussion shows how to relate the fitted parameters and
     1110second-order moments of a elliptical Gaussian of an arbitrary
     1111orientation with the parameters of an unrotate elliptical Gaussian.
     1112
     1113Consider a 2D Gaussian with an elliptical contour.  If the ellipse is
     1114oriented with the major axis along the x-axis, then the formula for
     1115such a Gaussian may be written $f = exp(-z)$ where
     1116\begin{equation}
     1117\label{aligned-ellipse}
     1118z = \frac{x^2}{2\sigma_a^2} + \frac{y^2}{2\sigma_b^2}
     1119\end{equation}
     1120with $\sigma_a$ the semi-major axis and $\sigma_b$ the semi-minor axis
     1121of the 1$\sigma$ contour.  Given such a Gaussian, we may measure its
     1122second moments by integration, and find that the second moment tensor
     1123is
     1124\begin{equation} \left| \begin{array}{cc}
     1125\sigma_a^2 & 0 \\
     11260 & \sigma_b^2 \\
     1127\end{array} \right| \end{equation}
     1128
     1129Now consider the same ellipse rotated to an arbitrary angle $\theta$.
     1130The formula for such a Gaussian may be written $f = exp(-z)$ where
     1131\begin{equation}
     1132\label{rotated-ellipse}
     1133z = \frac{x^2}{2\sigma_x^2} + \frac{y^2}{2\sigma_y^2} + \sigma_{xy}xy
     1134\end{equation}
     1135Note that, in the above form of the equation, $\sigma_{xy}$ goes to 0
     1136as the ellipse is rotated to be aligned with the x (or y) axis.  Thus,
     1137in this representation, $\sigma_{xy}$ is well behaved, but does not
     1138have the same units of length that $\sigma_x$ or $\sigma_y$ have.  Our
     1139goal is to determine the relationships between the rotated and
     1140unrotated components of the Gaussian formula as well as the second
     1141moments.
     1142
     1143To determine the behavior of $\sigma_x$, etc, under rotation, we start
     1144with the aligned ellipse (\ref{aligned-ellipse}) and rotate the
     1145coordinate frame by an angle $-\theta$:
     1146\begin{equation}
     1147\left( \begin{array}{c} x^\prime \\ y^\prime \end{array} \right) =
     1148\left| \begin{array}{cc} \cos \theta & \sin \theta \\
     1149                        -\sin \theta & \cos \theta
     1150\end{array} \right|
     1151\left( \begin{array}{c} x \\ y \end{array} \right)
     1152\end{equation}
     1153where $x^\prime$ and $y^\prime$ are the coordinates for the unrotated
     1154(aligned) ellipse.  Applying this rotation to (\ref{aligned-ellipse}) yields:
     1155\begin{equation}
     1156z = \frac{x^2 \cos^2 \theta + y^2 \sin^2 \theta + 2 x y \sin \theta \cos \theta}{2\sigma_a^2} +
     1157    \frac{x^2 \sin^2 \theta + y^2 \cos^2 \theta - 2 x y \sin \theta \cos \theta}{2\sigma_b^2}
     1158\end{equation}
     1159Grouping these terms together, we find:
     1160\begin{equation}
     1161z = \frac{x^2}{2}(\sigma_a^{-2} \cos^2 \theta + \sigma_b^{-2}\sin^2 \theta) +
     1162    \frac{y^2}{2}(\sigma_b^{-2} \cos^2 \theta + \sigma_a^{-2}\sin^2 \theta) -
     1163    \frac{xy}{2} \sin (2 \theta) (\sigma_b^{-2} - \sigma_a^{-2})
     1164\end{equation}
     1165We then associate the components of this equation with those of (\ref{rotated-ellipse}) and find:
     1166\begin{eqnarray}
     1167\sigma_x^{-2} & = & \sigma_a^{-2} \cos^2 \theta + \sigma_b^{-2}\sin^2 \theta \\
     1168\sigma_y^{-2} & = & \sigma_b^{-2} \cos^2 \theta + \sigma_a^{-2}\sin^2 \theta \\
     1169\sigma_{xy}   & = & \frac{1}{2} \sin (2 \theta) (\sigma_b^{-2} - \sigma_a^{-2})
     1170\end{eqnarray}
     1171Replacing $\cos^2$ and $\sin^2$ with the double-angle relationships, we find:
     1172\begin{eqnarray}
     1173\sigma_x^{-2} & = & \frac{1}{2}(\sigma_a^{-2} + \sigma_b^{-2}) - \frac{1}{2}(\sigma_b^{-2} - \sigma_a^{-2}) \cos (2 \theta) \\
     1174\sigma_y^{-2} & = & \frac{1}{2}(\sigma_a^{-2} + \sigma_b^{-2}) + \frac{1}{2}(\sigma_b^{-2} - \sigma_a^{-2}) \cos (2 \theta) \\
     1175\sigma_{xy}   & = & \frac{1}{2} \sin (2 \theta) (\sigma_b^{-2} - \sigma_a^{-2})
     1176\end{eqnarray}
     1177These formulae thus define the values of $\sigma_x$, $\sigma_y$, and
     1178$\sigma_{xy}$ given $\sigma_a$, $\sigma_b$, and $\theta$.  Note that
     1179in this equation and the one above, we represent the quantities in
     1180terms of $\sigma_b^{-2} - \sigma_a^{-2}$ which is always greater than
     11810, thus attributing the sign of the equation to the $\sin$ or $\cos$
     1182term.  This is necessary to determine the angle in the proper quadrant
     1183using the arctangent below.
     1184
     1185With the above relationships, we may now form combinations that help
     1186us to solve for $\sigma_a$, $\sigma_b$, and $\theta$:
     1187\begin{eqnarray}
     1188f_1 = \sigma_y^{-2} + \sigma_x^{-2} & = & \sigma_b^{-2} + \sigma_a^{-2} \\
     1189f_2 = \sigma_y^{-2} - \sigma_x^{-2} & = & (\sigma_b^{-2} - \sigma_a^{-2}) \cos (2 \theta) \\
     1190f_3 = \sqrt{f_2^2 + 4\sigma_{xy}^2} & = & \sigma_b^{-2} - \sigma_a^{-2}
     1191\end{eqnarray}
     1192From these, we may derive the equations for $\sigma_a$, $\sigma_b$, and $\theta$:
     1193\begin{eqnarray}
     1194\theta & = & \frac{1}{2} \arg (-2 \sigma_{xy}, f_2) \\
     1195\sigma_a & = & \sqrt{\frac{2}{f_1 - f_3}} \\
     1196\sigma_b & = & \sqrt{\frac{2}{f_1 + f_3}}
     1197\end{eqnarray}
     1198
     1199The relationship between the rotated ($m_{i,j}$) and unrotated
     1200($M_{i,j}$) second moments, the latter being equal to $\sigma_a^2$ and
     1201$\sigma_b^2$, is derived in a similar fashion.  We start with the
     1202point that the second moment is rotated as a tensor:
     1203\begin{equation}
     1204\left|
     1205\begin{array}{cc}
     1206m_{x,x} & m_{x,y} \\
     1207m_{y,x} & m_{y,y} \\
     1208\end{array} \right|
     1209=
     1210\left|
     1211\begin{array}{cc}
     1212+\cos \theta & +\sin \theta \\
     1213-\sin \theta & +\cos \theta \\
     1214\end{array} \right|
     1215\left|
     1216\begin{array}{cc}
     1217M_{x,x} & 0 \\
     12180       & M_{y,y} \\
     1219\end{array} \right|
     1220\left|
     1221\begin{array}{cc}
     1222+\cos \theta & -\sin \theta \\
     1223+\sin \theta & +\cos \theta \\
     1224\end{array} \right|
     1225\end{equation}
     1226Multiplying this out and substituting $\sigma_a^2$, $\sigma_b^2$ for $M_{x,x}$, $M_{y,y}$, we find:
     1227\begin{eqnarray}
     1228m_{x,x} & = & \sigma_a^{2} \cos^2 \theta + \sigma_b^{2}\sin^2 \theta \\
     1229m_{y,y} & = & \sigma_b^{2} \cos^2 \theta + \sigma_a^{2}\sin^2 \theta \\
     1230m_{x,y} & = & -\frac{1}{2} \sin (2 \theta) (\sigma_a^2 - \sigma_b^2)
     1231\end{eqnarray}
     1232Using the double-angle relationships, these become:
     1233\begin{eqnarray}
     1234m_{x,x} & = & \frac{1}{2}(\sigma_a^{2} + \sigma_b^{2}) + \frac{1}{2}(\sigma_a^{2} - \sigma_b^{2}) \cos (2 \theta) \\
     1235m_{y,y} & = & \frac{1}{2}(\sigma_a^{2} + \sigma_b^{2}) - \frac{1}{2}(\sigma_a^{2} - \sigma_b^{2}) \cos (2 \theta) \\
     1236m_{x,y} & = & -\frac{1}{2} \sin (2 \theta) (\sigma_a^{-2} - \sigma_b^{-2})
     1237\end{eqnarray}
     1238These three formulae define the second moments in terms of $\sigma_a$, $\sigma_b$, and $\theta$.
     1239
     1240We define equivalent intermediate products to the above:
     1241\begin{eqnarray}
     1242g_1 = m_{x,x} + m_{y,y}          & = & \sigma_a^{2} + \sigma_b^{2} \\
     1243g_2 = m_{x,x} - m_{y,y}          & = & (\sigma_a^{2} - \sigma_b^{2}) \cos (2 \theta) \\
     1244g_3 = \sqrt{f_2^2 + 4 m_{x,y}^2} & = & \sigma_a^{2} - \sigma_b^{-2}
     1245\end{eqnarray}
     1246From these, we may derive the equations for $\sigma_a$, $\sigma_b$, and $\theta$:
     1247\begin{eqnarray}
     1248\theta   & = & \frac{1}{2} \arg (-2 m_{x,y}, g_2) \\
     1249\sigma_a & = & \sqrt{\frac{g_1 - g_3}{2}} \\
     1250\sigma_b & = & \sqrt{\frac{g_1 + g_3}{2}}
     1251\end{eqnarray}
    11001252
    11011253\section{PSLib Astronomy Utilities}
     
    17931945\alpha_p & = & 180^\circ - 192.85948^\circ \\
    17941946\delta_p & = & 90^\circ - 62.87175^\circ \\
    1795 \phi_p & = & 90^\circ + 32.93192^\circ \\
     1947\phi_p & = & 90^\circ + 32.93192^\circ
    17961948\end{eqnarray}
    17971949
     
    27682920\begin{eqnarray}
    27692921x_p & = & \rho_x x \\
    2770 y_p & = & \rho_y y \\
     2922y_p & = & \rho_y y
    27712923\end{eqnarray}
    27722924%
     
    28302982\begin{eqnarray}
    28312983\alpha - \alpha_p & = & \arctan (\sin \alpha, \cos \alpha) \\
    2832 \delta            & = & \arcsin (\sin \delta) \\
     2984\delta            & = & \arcsin (\sin \delta)
    28332985\end{eqnarray}
    28342986%
     
    28452997y           & = & \frac{-\cos \theta \cos \phi}{\sin \theta} \\
    28462998\sin \theta & = & \zeta / \sqrt{1 + \zeta^2} \\
    2847 \cos \theta & = & 1 / \sqrt{1 + \zeta^2} \\
     2999\cos \theta & = & 1 / \sqrt{1 + \zeta^2}
    28483000\end{eqnarray}
    28493001
     
    28603012y           & = & -\cos \theta \cos \phi \\
    28613013\sin \theta & = & \sqrt{1 - R_\theta^2} \\
    2862 \cos \theta & = & R_\theta \\
     3014\cos \theta & = & R_\theta
    28633015\end{eqnarray}
    28643016
     
    29223074\begin{eqnarray}
    29233075x & = & \phi \left( 2 \cos \frac{2 \theta}{3} - 1 \right) \\
    2924 y & = & \pi \sin \frac{\theta}{3} \\
     3076y & = & \pi \sin \frac{\theta}{3}
    29253077\end{eqnarray}
    29263078
     
    29303082\theta & = & 3 \sin^{-1} \rho \\
    29313083\phi   & = & \frac{x}{1 - 4\rho^2} \\
    2932 {\rm where}\hspace{1cm} \rho & \equiv & y/\pi \\
     3084{\rm where}\hspace{1cm} \rho & \equiv & y/\pi
    29333085\end{eqnarray}
    29343086
     
    30313183This function is a two-dimensional Gaussian with an elliptical
    30323184cross-section and a constant local background:
    3033 \[
     3185\begin{equation}
    30343186f(x,y) = Z_o e^{-z} + S_o
    3035 \]
     3187\end{equation}
    30363188where
    3037 \[
     3189\begin{equation}
    30383190z = \frac{(x - x_o)^2}{2\sigma_x^2} + \frac{(y-y_o)^2}{2\sigma_y^2} + (x-x_o) (y - y_o) \sigma_{xy}
    3039 \]
     3191\end{equation}
    30403192
    30413193Below is the relationship between the \code{psModel} parameters and
     
    30803232This function is a polynomial approximation of a 2D Gaussian.  The
    30813233function is very similar to the real Gaussian:
    3082 \[
     3234\begin{equation}
    30833235f(x,y) = Z_o (1 + z + z^2/2 + z^3/6)^{-1} + S_o
    3084 \]
     3236\end{equation}
    30853237where
    3086 \[
     3238\begin{equation}
    30873239z = \frac{(x - x_o)^2}{2\sigma_x^2} + \frac{(y-y_o)^2}{2\sigma_y^2} + (x-x_o) (y - y_o) \sigma_{xy}
    3088 \]
     3240\end{equation}
    30893241
    30903242Below is the relationship between the \code{psModel} parameters and
     
    31333285than the Taylor series values of 1/2 and 1/6.  The
    31343286function is very similar to the pseudo-Gaussian:
    3135 \[
     3287\begin{equation}
    31363288f(x,y) = Z_o (1 + z + B_2 (z^2/2 + B_3 z^3/6))^{-1} + S_o
    3137 \]
     3289\end{equation}
    31383290where
    3139 \[
     3291\begin{equation}
    31403292z = \frac{(x - x_o)^2}{2\sigma_x^2} + \frac{(y-y_o)^2}{2\sigma_y^2} + (x-x_o) (y - y_o) \sigma_{xy}
    3141 \]
     3293\end{equation}
    31423294
    31433295Below is the relationship between the \code{psModel} parameters and
     
    31893341core, where the core has a different contour from the wings. 
    31903342
    3191 \[
     3343\begin{equation}
    31923344f(x,y) = Z_{\rm pk} (1 + z_1 + z_2^M)^{-1} + Sky
    3193 \]
     3345\end{equation}
    31943346where
    3195 \[
     3347\begin{equation}
    31963348z_1 = \frac{x^2}{2\sigma_{x,in}^2} + \frac{y^2}{2\sigma_{y,in}^2} + x y \sigma_{xy,in}
    31973349z_2 = \frac{x^2}{2\sigma_{x,out}^2} + \frac{y^2}{2\sigma_{y,out}^2} + x y \sigma_{xy,out}
    3198 \]
     3350\end{equation}
    31993351
    32003352\begin{verbatim}
Note: See TracChangeset for help on using the changeset viewer.