lim added files to sphere branch

git-svn-id: http://svn.stsci.edu/svn/ssb/stsci_python/stsci_python/branches/sphere@15878 fe389314-cf27-0410-b35b-8c050e845b92

Former-commit-id: b51412b5b545598ae101152eeaed470b08616176

1 files changed, 213 insertions, 0 deletions

diff --git a/lib/vector.py b/lib/vector.py
new file mode 100644
index 0000000..7323c92
--- /dev/null
+++ b/lib/vector.py
@@ -0,0 +1,213 @@
+# -*- coding: utf-8 -*-
+
+# Copyright (C) 2011 Association of Universities for Research in
+# Astronomy (AURA)
+#
+# Redistribution and use in source and binary forms, with or without
+# modification, are permitted provided that the following conditions
+# are met:
+#
+#     1. Redistributions of source code must retain the above
+#       copyright notice, this list of conditions and the following
+#       disclaimer.
+#
+#     2. Redistributions in binary form must reproduce the above
+#       copyright notice, this list of conditions and the following
+#       disclaimer in the documentation and/or other materials
+#       provided with the distribution.
+#
+#     3. The name of AURA and its representatives may not be used to
+#       endorse or promote products derived from this software without
+#       specific prior written permission.
+#
+# THIS SOFTWARE IS PROVIDED BY AURA ``AS IS'' AND ANY EXPRESS OR
+# IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
+# WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
+# ARE DISCLAIMED. IN NO EVENT SHALL AURA BE LIABLE FOR ANY DIRECT,
+# INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
+# (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
+# SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
+# HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT,
+# STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
+# ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED
+# OF THE POSSIBILITY OF SUCH DAMAGE.
+
+"""
+The `sphere.vector` module contains the basic operations for handling
+vectors and converting them to and from other representations.
+"""
+
+# THIRD-PARTY
+import numpy as np
+
+
+__all__ = ['radec_to_vector', 'vector_to_radec', 'normalize_vector',
+           'rotate_around']
+
+
+def radec_to_vector(ra, dec, degrees=True):
+    ur"""
+    Converts a location on the unit sphere from right-ascension and
+    declination to an *x*, *y*, *z* vector.
+
+    Parameters
+    ----------
+    ra, dec : scalars or 1-D arrays
+
+    degrees : bool, optional
+
+       If `True`, (default) *ra* and *dec* are in decimal degrees,
+       otherwise in radians.
+
+    Returns
+    -------
+    x, y, z : tuple of scalars or 1-D arrays of the same length
+
+    Notes
+    -----
+    Where right-ascension is *α* and declination is *δ*:
+
+    .. math::
+        x = \cos\alpha \cos\delta
+
+        y = \sin\alpha \cos\delta
+
+        z = \sin\delta
+    """
+    ra = np.asanyarray(ra)
+    dec = np.asanyarray(dec)
+
+    if degrees:
+        ra_rad = np.deg2rad(ra)
+        dec_rad = np.deg2rad(dec)
+    else:
+        ra_rad = ra
+        dec_rad = dec
+
+    cos_dec = np.cos(dec_rad)
+
+    return (
+        np.cos(ra_rad) * cos_dec,
+        np.sin(ra_rad) * cos_dec,
+        np.sin(dec_rad))
+
+
+def vector_to_radec(x, y, z, degrees=True):
+    ur"""
+    Converts a vector to right-ascension and declination.
+
+    Parameters
+    ----------
+    x, y, z : scalars or 1-D arrays
+        The input vectors
+
+    degrees : bool, optional
+        If `True` (default) the result is returned in decimal degrees,
+        otherwise radians.
+
+    Returns
+    -------
+    ra, dec : tuple of scalars or arrays of the same length
+
+    Notes
+    -----
+    Where right-ascension is *α* and declination is
+    *δ*:
+
+    .. math::
+        \alpha = \arctan2(y, x)
+
+        \delta = \arctan2(z, \sqrt{x^2 + y^2})
+    """
+    x = np.asanyarray(x)
+    y = np.asanyarray(y)
+    z = np.asanyarray(z)
+
+    result = (
+        np.arctan2(y, x),
+        np.arctan2(z, np.sqrt(x ** 2 + y ** 2)))
+
+    if degrees:
+        return np.rad2deg(result[0]), np.rad2deg(result[1])
+    else:
+        return result
+
+
+def normalize_vector(x, y, z, inplace=False):
+    ur"""
+    Normalizes a vector so it falls on the unit sphere.
+
+    *x*, *y*, *z* may be scalars or 1-D arrays
+
+    Parameters
+    ----------
+    x, y, z : scalars or 1-D arrays of the same length
+        The input vectors
+
+    inplace : bool, optional
+        When `True`, the original arrays will be normalized in place,
+        otherwise a normalized copy is returned.
+
+    Returns
+    -------
+    X, Y, Z : scalars of 1-D arrays of the same length
+        The normalized output vectors
+    """
+    x = np.asanyarray(x)
+    y = np.asanyarray(y)
+    z = np.asanyarray(z)
+
+    l = np.sqrt(x ** 2 + y ** 2 + z ** 2)
+
+    if inplace:
+        x /= l
+        y /= l
+        z /= l
+        return (x, y, z)
+    else:
+        return (x / l, y / l, z / l)
+
+
+def rotate_around(x, y, z, u, v, w, theta, degrees=True):
+    r"""
+    Rotates the vector (*x*, *y*, *z*) around the arbitrary axis defined by
+    vector (*u*, *v*, *w*) by *theta*.
+
+    It is assumed that both (*x*, *y*, *z*) and (*u*, *v*, *w*) are
+    already normalized.
+
+    Parameters
+    ----------
+    x, y, z : doubles
+        The normalized vector to rotate
+
+    u, v, w : doubles
+        The normalized vector to rotate around
+
+    theta : double, or array of doubles
+        The amount to rotate
+
+    degrees : bool, optional
+        When `True`, *theta* is given in degrees, otherwise radians.
+
+    Returns
+    -------
+    X, Y, Z : doubles
+        The rotated vector
+    """
+    if degrees:
+        theta = np.deg2rad(theta)
+
+    u2 = u**2
+    v2 = v**2
+    w2 = w**2
+    costheta = np.cos(theta)
+    sintheta = np.sin(theta)
+    icostheta = 1.0 - costheta
+
+    det = (-u*x - v*y - w*z)
+    X = (-u*det)*icostheta + x*costheta + (-w*y + v*z)*sintheta
+    Y = (-v*det)*icostheta + y*costheta + ( w*x - u*z)*sintheta
+    Z = (-w*det)*icostheta + z*costheta + (-v*x + u*y)*sintheta
+
+    return X, Y, Z