mangadap.util.geometry module¶
Provides a set of utility functions dealing with computational geometry.
Revision history¶
2015: Original implementation by K. Westfall (KBW)20 May 2015: (KBW) Documentation and Sphinx tests28 Mar 2016: (KBW) Copied old projected_disk_plane toSemiMajorAxisCoo
. Coordinate projection functions (e.g.,SemiMajorAxisCoo.polar()
) can now take array-like objects as arguments.08 Sep 2016: (KBW) Allowpoint_inside_polygon()
to accept multiple coordinates.
Copyright © 2019, SDSS-IV/MaNGA Pipeline Group
-
class
mangadap.util.geometry.
SemiMajorAxisCoo
(xc=None, yc=None, rot=None, pa=None, ell=None)[source]¶ Bases:
object
Calculate the semi-major axis coordinates given a set of input parameters following \({\mathbf x} = {\mathbf A}^{-1}\ {\mathbf b}\), where
\[ \begin{align}\begin{aligned}\begin{split}{\mathbf A} = \left[ \begin{array}{rrrrrr} 1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 & 0 & 0 \\ \cos\psi & \sin\psi & -1 & 0 & 0 & 0 \\ -\sin\psi & \cos\psi & 0 & -1 & 0 & 0 \\ 0 & 0 & \sin\phi_0 & \cos\phi_0 & -1 & 0 \\ 0 & 0 & -\cos\phi_0 & \sin\phi_0 & 0 & \varepsilon-1 \end{array} \right]\end{split}\\\begin{split}{\mathbf b} = \left[ \begin{array}{r} x_f \\ y_f \\ -x_0 \\ -y_0 \\ 0 \\ 0 \end{array} \right]\end{split}\end{aligned}\end{align} \]such that
\[\begin{split}{\mathbf x} = \left[ \begin{array}{r} x_f \\ y_f \\ x_s \\ y_s \\ x_a \\ y_a \end{array} \right]\end{split}\]- and:
- \(\psi\) is the Cartesian rotation of the focal-plane relative to the sky-plane (+x toward East; +y toward North),
- \(\phi_0\) is the on-sky position angle of the major axis of the ellipse, defined as the angle from North through East
- \(\varepsilon=1-b/a\) is the ellipticity based on the the semi-minor to semi-major axis ratio (\(b/a\)).
- \((x_f,y_f)\) is the sky-right, focal-plane position relative to a reference on-sky position \((x_0,y_0)\) relative to the center of the ellipse (galaxy center),
- \((x_s,y_s)\) is the on-sky position of \((x_f,y_f)\) relative to the center of the ellipse, and
- \((x_a,y_a)\) is the Cartesian position of \((x_f,y_f)\) in units of the semi-major axis.
This form is used such that \({\mathbf A}\) need only be defined once per class instance.
The class also allows for inverse calculations, i.e., calculating the focal-plane positions provide the semi-major axis coordinates. In this case,
\[ \begin{align}\begin{aligned}\begin{split}{\mathbf C} = \left[ \begin{array}{rrrr} \cos\psi & \sin\psi & -1 & 0 \\ -\sin\psi & \cos\psi & 0 & -1 \\ 0 & 0 & \sin\phi_0 & \cos\phi_0 \\ 0 & 0 & -\cos\phi_0 & \sin\phi_0 \end{array} \right]\end{split}\\\begin{split}{\mathbf d} = \left[ \begin{array}{r} -x_0 \\ -y_0 \\ x_a \\ y_a (1-\varepsilon) \end{array} \right]\end{split}\end{aligned}\end{align} \]such that
\[\begin{split}{\mathbf f} = \left[ \begin{array}{r} x_f \\ y_f \\ x_s \\ y_s \end{array} \right]\end{split}\]and \({\mathbf f} = {\mathbf C}^{-1}\ {\mathbf d}\).
Parameters: - xc (float) – Same as \(x_0\), defined above
- yc (float) – Same as \(y_0\), defined above
- rot (float) – Same as \(\psi\), defined above
- pa (float) – Same as \(\phi_0\), defined above
- ell (float) – Same as \(\varepsilon\), defined above
-
xc,yc
a reference on-sky position relative to the center of the ellipse (galaxy center); same as \((x_0,y_0)\) defined above
Type: float,float
-
rot
¶ Cartesian rotation of the focal-plane relative to the sky-plane (+x toward East; +y toward North); same as \(\psi\) defined above
Type: float
-
pa
¶ On-sky position angle of the major axis of the ellipse, defined as the angle from North through East and is the same as \(\phi_0\) defined above
Type: float
-
ell
¶ Ellipticity define as \(\varepsilon=1-b/a\), based on the semi-minor to semi-major axis ratio (\(b/a\)) of the ellipse.
Type: float
-
A
¶ The coordinate transformation matrix
Type: numpy.ndarray
-
Alu
¶ The lu array returned by scipy.linalg.lu_factor, which is used to calculate the LU decomposition of \({\mathbf A}\)
Type: numpy.ndarray
-
Apiv
¶ The piv array returned by scipy.linalg.lu_factor, which is used to calculate the LU decomposition of \({\mathbf A}\)
Type: numpy.ndarray
-
B
¶ The vector \({\mathbf b}\), as defined above, used to calculate \({\mathbf x} = {\mathbf A}^{-1}\ {\mathbf b}\)
Type: numpy.ndarray
-
C
¶ The coordinate transformation matrix use for the inverse operations
Type: numpy.ndarray
-
Clu
¶ The lu array returned by scipy.linalg.lu_factor, which is used to calculate the LU decomposition of \({\mathbf C}\)
Type: numpy.ndarray
-
Cpiv
¶ The piv array returned by scipy.linalg.lu_factor, which is used to calculate the LU decomposition of \({\mathbf C}\)
Type: numpy.ndarray
-
D
¶ The vector \({\mathbf d}\), as defined above, used to calculate \({\mathbf f} = {\mathbf C}^{-1}\ {\mathbf d}\)
Type: numpy.ndarray
-
_calculate_cartesian
(r, theta)[source]¶ Invert the calculation of the semi-major-axis polar coordinates to calculate the semi-major-axis Cartesian coordinates \((x_a,y_a)\) using
\[\begin{split}x_a &= \pm R / \sqrt{1 + \tan^2\theta}\\ y_a &= -x_a\ \tan\theta\end{split}\]where \(x_a\) is negative when \(\pi/2 \leq \theta < 3\pi/2\).
Parameters: r,theta (array-like) – The semi-major-axis polar coordinates \((R,\theta)\). Returns: The semi-major-axis Cartesian coordinates: \(x_a, y_a\). Return type: numpy.ndarray
-
_calculate_polar
(x, y)[source]¶ Calculate the polar coordinates (radius and azimuth) provided the Cartesian semi-major-axis coordinates \((x_a,y_a)\) using
\[\begin{split}R &= \sqrt{x_a^2 + y_a^2} \\ \theta &= \tan^{-1}\left(\frac{-y_a}{x_a}\right)\end{split}\]Parameters: x,y (array-like) – The semi-major-axis Cartesian coordinates \((x_a,y_a)\). Returns: The semi-major-axis polar coordinates: \(R, \theta\). Return type: numpy.ndarray
-
_defined
()[source]¶ Determine if the object is defined such that its methods can be used to convert between coordinate systems.
-
_setA
()[source]¶ Set the transformation matrix and calculate its LU decomposition for forward operations.
-
_setB
(x, y)[source]¶ Set the on-sky coordinate vector for forward operations.
Parameters: y (x,) – Single values for use in calculating the semi-major-axis coordinates.
-
_setC
()[source]¶ Set the transformation matrix and calculate its LU decomposition for inverse operations.
-
_setD
(x, y)[source]¶ Set the semi-major-axis coordinate vector for inverse operations.
Parameters: y (x,) – Single values for use in calculating the on-sky focal plane coordinates.
-
cartesian
(x, y)[source]¶ Calculate \({\mathbf x}\) using
solve()
for the provided \((x_f,y_f)\) and return the semi-major-axis Cartesian and coordinates, \((x_a,y_a)\).Parameters: x,y (array-like) – The coordinate \((x_f,y_f)\), which is the sky-right, focal-plane position relative to a reference on-sky position \((x_0,y_0)\) relative to the center of the ellipse (galaxy center), Returns: Two arrays with the semi-major-axis Cartesian coordinates, \(x_a, y_a\). Return type: numpy.ndarray
-
cartesian_invert
(x, y)[source]¶ Calculate \({\mathbf f}\) using
solve()
for the provided \((x_a,y_a)\) and return focal-plane cartesian coordinates \((x_f,y_f)\).Parameters: x,y (array-like) – The semi-major-axis Cartesian coordinates \((x_a,y_a)\). Returns: The focal-plane Cartesian coordinates \((x_f,y_f)\). Return type: numpy.ndarray
-
coo
(x, y)[source]¶ Calculate \({\mathbf x}\) using
solve()
for the provided \((x_f,y_f)\) and return the semi-major-axis Cartesian and polar coordinates, \((x_a,y_a)\) and \((R,\theta)\). This combines the functionality ofcartesian()
andpolar()
, and so is more efficient than using these both separately.Parameters: x,y (array-like) – The coordinates \((x_f,y_f)\), which are the sky-right, focal-plane position relative to a reference on-sky position \((x_0,y_0)\) relative to the center of the ellipse (galaxy center), Returns: Four arrays with the semi-major-axis Cartesian and polar coordinates: \(x_a, y_a, R, \theta\). Return type: numpy.ndarray
-
polar
(x, y)[source]¶ Calculate \({\mathbf x}\) using
solve()
for the provided \((x_f,y_f)\) and return the semi-major-axis polar coordinates, \((R,\theta)\), where\[\begin{split}R &= \sqrt{x_a^2 + y_a^2} \\ \theta &= \tan^{-1}\left(\frac{-y_a}{x_a}\right)\end{split}\]Parameters: x,y (array-like) – The coordinate \((x_f,y_f)\), which is the sky-right, focal-plane position relative to a reference on-sky position \((x_0,y_0)\) relative to the center of the ellipse (galaxy center), Returns: Two arrays with the semi-major-axis polar coordinates: \(R, \theta\). Return type: numpy.ndarray
-
polar_invert
(r, theta)[source]¶ Calculate \({\mathbf f}\) using
solve()
for the provided \((R,\theta)\) and return focal-plane cartesian coordinates \((x_f,y_f)\).Parameters: r,theta (array-like) – The semi-major-axis polar coordinates \((R,\theta)\). Returns: Two arrays with the focal-plane Cartesian coordinates \((x_f,y_f)\). Return type: numpy.ndarray
-
solve
(x, y)[source]¶ Use scipy.linalg.lu_solve to solve \({\mathbf x} = {\mathbf A}^{-1}\ {\mathbf b}\).
Parameters: x,y (array-like) – The coordinates \((x_f,y_f)\), which are the sky-right, focal-plane Cartesian coordinates relative to a reference on-sky position \((x_0,y_0)\), which is relative to the center of the ellipse (galaxy center). Returns: The \({\mathbf x}\) vectors (separated by rows) as defined by the solution to \({\mathbf A}^{-1}\ {\mathbf b}\) Return type: numpy.ndarray Raises: ValueError
– Raised if object was not properly defined or if the X and Y arrays do not have the same size.
-
solve_inverse
(x, y)[source]¶ Use scipy.linalg.lu_solve to solve \({\mathbf f} = {\mathbf C}^{-1}\ {\mathbf d}\).
Parameters: x,y (array-like) – The semi-major-axis Cartesian coordinates \((x_a,y_a)\). Returns: The \({\mathbf f}\) vector as defined by the solution to \({\mathbf C}^{-1}\ {\mathbf d}\) Return type: numpy.ndarray Raises: ValueError
– Raised if object was not properly defined or if the X and Y arrays do not have the same size.
-
mangadap.util.geometry.
point_inside_polygon
(polygon, point)[source]¶ Determine if a point is inside a polygon using the winding number.
Parameters: - polygon (numpy.ndarray) – An Nx2 array containing the x,y coordinates of a polygon. The points should be ordered either counter-clockwise or clockwise.
- point (numpy.ndarray) – A 2-element array defining the x,y position of the point to use as a reference for the winding number.
Returns: True if the point is inside the polygon.
Return type: bool
Warning
If the point is on the polygon (or very close to it w.r.t. the machine precision), the returned value is False.
-
mangadap.util.geometry.
polygon_winding_number
(polygon, point)[source]¶ Determine the winding number of a 2D polygon about a point. The code does not check if the polygon is simple (no interesecting line segments). Algorithm taken from Numerical Recipies Section 21.4.
Parameters: - polygon (numpy.ndarray) – An Nx2 array containing the x,y coordinates of a polygon. The points should be ordered either counter-clockwise or clockwise.
- point (numpy.ndarray) – A 2-element array defining the x,y position of the point to use as a reference for the winding number.
Returns: Winding number of polygon w.r.t. point
Return type: int
Raises: ValueError
– Raised if polygon is not 2D, if polygon does not have two columns, or if point is not a 2-element array.