# Path: skyfield/constants.py
# AU_M = 149597870700             # per IAU 2012 Resolution B2
# AU_KM = 149597870.700
# ASEC360 = 1296000.0
# DAY_S = 86400.0
# ASEC2RAD = 4.848136811095359935899141e-6
# DEG2RAD = 0.017453292519943296
# RAD2DEG = 57.295779513082321
# C = 299792458.0                            # m/s
# ANGVEL = 7.2921150e-5           # radians/s
# ERAD = 6378136.6                # meters
# IERS_2010_INVERSE_EARTH_FLATTENING = 298.25642
# GS = 1.32712440017987e+20
# T0 = 2451545.0
# B1950 = 2433282.4235
# C_AUDAY = C * DAY_S / AU_M
# 
# Path: skyfield/functions.py
# def dots(v, u):
#     """Given one or more vectors in `v` and `u`, return their dot products.
# 
#     This works whether `v` and `u` each have the shape ``(3,)``, or
#     whether they are each whole arrays of corresponding x, y, and z
#     coordinates and have shape ``(3, N)``.
# 
#     """
#     return (v * u).sum(axis=0)


from numpy import (abs, arcsin, arccos, arctan2, array, clip, cos,
                   minimum, pi, sin, sqrt, tan, where, zeros_like)
from .constants import (AU_M, ANGVEL, DAY_S, DEG2RAD, ERAD,
                        IERS_2010_INVERSE_EARTH_FLATTENING, RAD2DEG, T0, tau)
from .functions import dots

    a = ERAD / AU_M
    f = 1.0 / IERS_2010_INVERSE_EARTH_FLATTENING
    e2 = 2.0*f - f*f
    i = 0
    C = 1.0
    while i < iterations:
        i += 1
        C = 1.0 / sqrt(1.0 - e2 * (sin(lat) ** 2.0))
        lat = arctan2(z + a * C * e2 * sin(lat), R)
    elevation_m = ((R / cos(lat)) - a * C) * AU_M
    return lat, lon, elevation_m

def compute_limb_angle(position_au, observer_au):
    """Determine the angle of an object above or below the Earth's limb.

    Given an object's GCRS `position_au` |xyz| vector and the position
    of an `observer_au` as a vector in the same coordinate system,
    return a tuple that provides `(limb_ang, nadir_ang)`:

    limb_angle
        Angle of observed object above (+) or below (-) limb in degrees.
    nadir_angle
        Nadir angle of observed object as a fraction of apparent radius
        of limb: <1.0 means below the limb, =1.0 means on the limb, and
        >1.0 means above the limb.

    """
    # Compute the distance to the object and the distance to the observer.


# Path: skyfield/constants.py
# AU_M = 149597870700             # per IAU 2012 Resolution B2
# AU_KM = 149597870.700
# ASEC360 = 1296000.0
# DAY_S = 86400.0
# ASEC2RAD = 4.848136811095359935899141e-6
# DEG2RAD = 0.017453292519943296
# RAD2DEG = 57.295779513082321
# C = 299792458.0                            # m/s
# ANGVEL = 7.2921150e-5           # radians/s
# ERAD = 6378136.6                # meters
# IERS_2010_INVERSE_EARTH_FLATTENING = 298.25642
# GS = 1.32712440017987e+20
# T0 = 2451545.0
# B1950 = 2433282.4235
# C_AUDAY = C * DAY_S / AU_M
# 
# Path: skyfield/functions.py
# def dots(v, u):
#     """Given one or more vectors in `v` and `u`, return their dot products.
# 
#     This works whether `v` and `u` each have the shape ``(3,)``, or
#     whether they are each whole arrays of corresponding x, y, and z
#     coordinates and have shape ``(3, N)``.
# 
#     """
#     return (v * u).sum(axis=0)


# Path: skyfield/earthlib.py
from numpy import (abs, arcsin, arccos, arctan2, array, clip, cos,
                   minimum, pi, sin, sqrt, tan, where, zeros_like)
from .constants import (AU_M, ANGVEL, DAY_S, DEG2RAD, ERAD,
                        IERS_2010_INVERSE_EARTH_FLATTENING, RAD2DEG, T0, tau)
from .functions import dots


    a = ERAD / AU_M
    f = 1.0 / IERS_2010_INVERSE_EARTH_FLATTENING
    e2 = 2.0*f - f*f
    i = 0
    C = 1.0
    while i < iterations:
        i += 1
        C = 1.0 / sqrt(1.0 - e2 * (sin(lat) ** 2.0))
        lat = arctan2(z + a * C * e2 * sin(lat), R)
    elevation_m = ((R / cos(lat)) - a * C) * AU_M
    return lat, lon, elevation_m

def compute_limb_angle(position_au, observer_au):
    """Determine the angle of an object above or below the Earth's limb.

    Given an object's GCRS `position_au` |xyz| vector and the position
    of an `observer_au` as a vector in the same coordinate system,
    return a tuple that provides `(limb_ang, nadir_ang)`:

    limb_angle
        Angle of observed object above (+) or below (-) limb in degrees.
    nadir_angle
        Nadir angle of observed object as a fraction of apparent radius
        of limb: <1.0 means below the limb, =1.0 means on the limb, and
        >1.0 means above the limb.

    """
    # Compute the distance to the object and the distance to the observer.


    disobj = sqrt(dots(position_au, position_au))