fallback for non-astropy and python 2.7

Author: scivision

README.rst

@@ -19,13 +19,19 @@ Python 3-D coordinate conversions
 ==================================
 
 Python coordinate conversions, following convention of several popular Matlab routines.
+Tested from Python 2.7 through Python 3.7+
+
+For those not having AstroPy, lower accuracy fallback functions are included for some functions.
 
 .. contents::
 
+
+:prereqs: AstroPy, Numpy, Python 3 or Python 2.7
+
 Install
 =======
 Development::
-
+    
     python -m pip install -e .
 
 
@@ -36,12 +42,16 @@ simple::
 
 Usage
 =====
-a few quick examples::
+a few quick examples
+
+.. code:: python
 
-   from pymap3d import *
+   import pymap3d as pm
 
-   lat,lon,alt = eci2geodetic(eci,t)
-   az,el,range = geodetic2aer(lat,lon,alt,42,-82,0)
+   lat,lon,alt = pm.eci2geodetic(eci, t)
+   
+   az,el,range = pm.geodetic2aer(lat, lon, alt, 42, -82, 0)
+   
 
 Functions
 ==========
@@ -60,5 +70,6 @@ Popular mapping toolbox functions ported to Python include::
 
 Caveats
 =======
-Atmospheric effects neglected in all functions not invoking AstroPy.
-Planetary perturbations and nutation etc. not fully considered.
+
+* Atmospheric effects neglected in all functions not invoking AstroPy. Need to update code to add these input parameters (just start a GitHub Issue to request).
+* Planetary perturbations and nutation etc. not fully considered.

pymap3d/__init__.py

@@ -13,10 +13,8 @@
 
 see test.py for example uses.
 """
-
 from __future__ import division
 from six import string_types,PY2
-from dateutil.parser import parse
 from datetime import datetime
 import numpy as np
 from numpy import sin, cos, tan, sqrt, radians, arctan2, hypot, degrees
@@ -26,6 +24,10 @@
     from astropy.coordinates import Angle,SkyCoord, EarthLocation, AltAz, ICRS
 except ImportError:
     Time = None
+#
+from .vallado import vazel2radec, vradec2azel
+from .timeconv import str2dt
+
 
 class EarthEllipsoid:
 
@@ -418,20 +420,6 @@ def get_radius_normal(lat_radians, ell):
     return a**2 / sqrt(
         a**2 * (cos(lat_radians))**2 + b**2 *
         (sin(lat_radians))**2)
-
-
-def str2dt(t):
-    """
-    output: datetime
-    """
-
-    t = np.atleast_1d(t)
-    if isinstance(t[0], string_types):
-        t = [parse(T) for T in t]
-
-    assert isinstance(t[0], datetime), 'did not convert {} to datetime'.format(type(t[0]))
-
-    return t
 #%% internal use
 def _depack(x0):
     m, n = x0.shape
@@ -490,12 +478,9 @@ def _uvw2enu(u, v, w, lat0, lon0, deg):
 
 #%% azel radec
 def azel2radec(az_deg, el_deg, lat_deg, lon_deg, t):
-    if PY2:
-        raise RuntimeError('Python 2 is not supported for this AstroPy operation. It is time to upgrade to Python 3 (released 2008)')
-
-    if Time is None:
-        raise ImportError('You need to install AstroPy')
 
+    if PY2 or Time is None: # non-AstroPy method, less accurate
+        return vazel2radec(az_deg, el_deg, lat_deg, lon_deg, t)
 
     t = str2dt(t)
 
@@ -511,7 +496,7 @@ def azel2radec(az_deg, el_deg, lat_deg, lon_deg, t):
 
 def radec2azel(ra_deg, dec_deg, lat_deg, lon_deg, t):
     if Time is None:
-        raise ImportError('You need to install AstroPy')
+        return vradec2azel(ra_deg, dec_deg, lat_deg, lon_deg, t)
 #%% input trapping
     t = str2dt(t)
     lat_deg = np.atleast_1d(lat_deg)

pymap3d/datetime2hourangle.py

@@ -1,40 +1,54 @@
 #!/usr/bin/env python
+from __future__ import division
 from numpy import atleast_1d, empty_like, pi, nan
 from datetime import datetime
-from astropy.time import Time
-import astropy.units as u
-from astropy.coordinates import Longitude
+try:
+    from astropy.time import Time
+    import astropy.units as u
+    from astropy.coordinates import Longitude
+except ImportError:
+    Time = None
+
+#
+from .timeconv import str2dt
 """
-input longitude and output are in RADIANS!
+The "usevallado" datetime to julian runs 4 times faster than astropy.
+However, AstroPy is more accurate.
 
-Note: The "usevallado" datetime to julian runs 4 times faster than astropy.
-AstroPy is more accurate, but version 1.0.0 doesn't auto-download new IERS data
-for current/future times, giving an
-IndexError: (some) times are outside of range covered by IERS table.
 """
 
 
 def datetime2sidereal(t, lon_radians, usevallado=True):
     """
+    algorithm by D. Vallado Fundamentals of Astrodynamics and Applications
      lon: longitude in RADIANS
     """
     if usevallado:
         jd = datetime2julian(t)
-        gst = julian2sidereal(jd)  # Greenwich Sidereal time RADIANS
-        return gst + lon_radians  # radians # Algorithm 15 p. 188 rotate to LOCAL SIDEREAL TIME
+# %% Greenwich Sidereal time RADIANS
+        gst = julian2sidereal(jd)
+# %% Algorithm 15 p. 188 rotate GST to LOCAL SIDEREAL TIME
+        tsr = gst + lon_radians  # radians
     else:  # astropy
-        return Time(t).sidereal_time(kind='apparent',
-                                     longitude=Longitude(lon_radians, unit=u.radian)).radian
+        if Time is not None:
+            tsr = Time(t).sidereal_time(kind='apparent',
+                                    longitude=Longitude(lon_radians, unit=u.radian)).radian
+        else:
+            raise ImportError('AstroPy required, or use "usevallado=True"')
+
+    return tsr
 
 
 def datetime2julian(t):
     """
     from D.Vallado Fundamentals of Astrodynamics and Applications p.187
      and J. Meeus Astronomical Algorithms 1991 Eqn. 7.1 pg. 61
     """
-    t= atleast_1d(t)
 
-    assert isinstance(t[0],datetime)
+    t = str2dt(t)
+    t = atleast_1d(t)
+
+    assert isinstance(t[0], datetime)
 
     jDate = empty_like(t, dtype=float)  # yes we need the dtype!
 
@@ -53,21 +67,25 @@ def datetime2julian(t):
         A = int(year / 100.0)
         B = 2 - A + int(A / 4.)
         C = ((d.second / 60. + d.minute) / 60. + d.hour) / 24.
-        jDate[i] = (int(365.25 * (year + 4716)) + 
+        jDate[i] = (int(365.25 * (year + 4716)) +
                     int(30.6001 * (month + 1)) + d.day + B - 1524.5 + C)
 
     return jDate
 
 
 def julian2sidereal(juliandate):
-    # D. Vallado Ed. 4
-    # TODO needs unit testing
-    # Julian centuries from J2000.0
-    
-     # Vallado Eq. 3-42 p. 184, Seidelmann 3.311-1
-    tUT1 = (juliandate - 2451545.0) / 36525. 
-
-    gmst_sec = (67310.54841 + (876600 * 3600 + 8640184.812866) * 
+    """
+    D. Vallado Ed. 4
+    TODO needs unit testing
+
+    input:
+    juliandate: Julian centuries from J2000.0
+    """
+
+# %% Vallado Eq. 3-42 p. 184, Seidelmann 3.311-1
+    tUT1 = (juliandate - 2451545.0) / 36525.
+
+    gmst_sec = (67310.54841 + (876600 * 3600 + 8640184.812866) *
         tUT1 + 0.093104 * tUT1**2 - 6.2e-6 * tUT1**3) # Eqn. 3-47 p. 188
 
     # 1/86400 and %(2*pi) implied by units of radians

pymap3d/vallado.py

@@ -1,8 +1,8 @@
 #!/usr/bin/env python
 """
-converts right ascension, declination to azimuth, elevation and vice versa
-
-
+converts right ascension, declination to azimuth, elevation and vice versa.
+Normally do this via AstroPy.
+These functions are fallbacks for those who don't wish to use AstroPy (perhaps Python 2.7 users).
 
 inputs:
 ra_deg: numpy ndarray of right ascension values [degrees]
@@ -17,13 +17,13 @@
 
 Michael Hirsch implementation of algorithms from D. Vallado
 """
+from __future__ import division
 from numpy import sin, cos, degrees, radians, arcsin, arctan2, atleast_1d
-
 #
 from .datetime2hourangle import datetime2sidereal
 
 
-def azel2radec(az_deg, el_deg, lat_deg, lon_deg, t):
+def vazel2radec(az_deg, el_deg, lat_deg, lon_deg, t):
     """
     from D.Vallado Fundamentals of Astrodynamics and Applications
     p.258-259
@@ -40,7 +40,7 @@ def azel2radec(az_deg, el_deg, lat_deg, lon_deg, t):
     el = radians(el_deg)
     lat = radians(lat_deg)
     lon = radians(lon_deg)
-    # Vallado "algorithm 28" p 268
+# %% Vallado "algorithm 28" p 268
     dec = arcsin(sin(el) * sin(lat) + cos(el) * cos(lat) * cos(az))
 
     lha = arctan2(-(sin(az) * cos(el)) / cos(dec),
@@ -52,7 +52,7 @@ def azel2radec(az_deg, el_deg, lat_deg, lon_deg, t):
     return degrees(lst - lha) % 360, degrees(dec)
 
 
-def radec2azel(dtime,ra_deg,dec_deg,lat_deg,lon_deg):
+def vradec2azel(ra_deg,dec_deg,lat_deg,lon_deg,t):
     """
     from D. Vallado "Fundamentals of Astrodynamics and Applications "
        4th Edition Ch. 4.4 pg. 266-268
@@ -62,11 +62,12 @@ def radec2azel(dtime,ra_deg,dec_deg,lat_deg,lon_deg):
     lat = radians(lat_deg)
     lon = radians(lon_deg)
 
-    lst = datetime2sidereal(dtime,lon) #RADIANS
-    lha = lst - ra #Eq. 4-11 p. 267 LOCAL HOUR ANGLE
-
-    el = arcsin(sin(lat) * sin(dec) + cos(lat) * cos(dec) * cos(lha) ) #Eq. 4-12 p. 267
-    #combine Eq. 4-13 and 4-14 p. 268
+    lst = datetime2sidereal(t, lon) #RADIANS
+# %% Eq. 4-11 p. 267 LOCAL HOUR ANGLE
+    lha = lst - ra
+# %% #Eq. 4-12 p. 267
+    el = arcsin(sin(lat) * sin(dec) + cos(lat) * cos(dec) * cos(lha) )
+# %% combine Eq. 4-13 and 4-14 p. 268
     az = arctan2( -sin(lha) * cos(dec) / cos(el),
                    ( sin(dec) - sin(el) * sin(lat) ) / (cos(el)*cos(lat))  )
 

setup.py

@@ -1,12 +1,12 @@
 #!/usr/bin/env python
-install_requires = ['python-dateutil','pytz','numpy','astropy']
+install_requires = ['six','python-dateutil','pytz','numpy','astropy']
 tests_require = ['nose','coveralls']
 # %%
 from setuptools import setup,find_packages
 
 setup(name='pymap3d',
       packages=find_packages(),
-      version = '1.2.6',
+      version = '1.3.0',
       description='Python coordinate conversions, following convention of several popular Matlab routines.',
       long_description=open('README.rst').read(),
       author = 'Michael Hirsch, Ph.D.',

tests/test.py

@@ -3,43 +3,64 @@
 runs tests
 """
 from datetime import datetime
-import six
-import logging
 from numpy import asarray,radians
 from numpy.testing import assert_allclose, assert_almost_equal,run_module_suite
 from pytz import UTC
-from dateutil.parser import parse
 #
 import pymap3d as pm
 from pymap3d.haversine import angledist,angledist_meeus
 from pymap3d.vincenty import vreckon,vdist
 from pymap3d.datetime2hourangle import datetime2sidereal
+from pymap3d.vallado import vazel2radec, vradec2azel
+from pymap3d.timeconv import str2dt
+
+def test_str2dt():
+
+    assert str2dt('2014-04-06T08:00:00Z') == datetime(2014,4,6,8, tzinfo=UTC)
+    ti = [str2dt('2014-04-06T08:00:00Z'), str2dt('2014-04-06T08:01:02Z')]
+    to = [datetime(2014,4,6,8, tzinfo=UTC),  datetime(2014,4,6,8,1,2, tzinfo=UTC)]
+    assert  ti == to   # even though ti is numpy array of datetime and to is list of datetime
+
+# %%
+t = '2014-04-06T08:00:00Z'
+rdlat, rdlon = (65, -148)
+ra, dec = (166.5032081149338,55.000011165405752)
+ha = 45.482789587392013
+azi, eli = (180.1,80)
 
 def test_datetime2sidereal():
-    sdrtest = datetime2sidereal(datetime(2014,4,6,8), radians(-148), False)
-    assert_allclose(sdrtest,2.9065780550600806,rtol=1e-5)
+    sdrapy = datetime2sidereal(t, radians(rdlon), False)
+    assert_allclose(sdrapy, 2.9065780550600806, rtol=1e-5)
+
+    sdrvallado = datetime2sidereal(t, radians(rdlon), True)
+    assert_allclose(sdrvallado, 2.9065780550600806, rtol=1e-5)
+
 
 def test_azel2radec():
-    if six.PY2:
-        logging.error('AstroPy does not support Python 2 for this function. Python 3 was released in 2008.')
-        return
+    R,D = pm.azel2radec(azi, eli, rdlat, rdlon, t)
+    assert_allclose(R, ra, rtol=1e-2)
+    assert_allclose(D, dec, rtol=1e-2)
+
+    Rv, Dv = vazel2radec(azi, eli, rdlat, rdlon, t)
+    assert_allclose(Rv, ra)
+    assert_allclose(Dv, dec)
 
-    ra,dec = pm.azel2radec(180.1, 80,
-                           65, -148,
-                           '2014-04-06T08:00:00Z')
-    assert_allclose(ra,166.5032081149338,rtol=1e-2)
-    assert_allclose(dec,55.000011165405752,rtol=1e-2)
 
 def test_radec2azel():
-    aztest, eltest = pm.radec2azel(166.5032081149338,55.000011165405752, 65, -148,
-                            parse('2014-04-06T08:00:00'))
-    assert_allclose(aztest, 179.40287898,rtol=1e-2) #180.1 from vallado, but his is less precise
-    assert_allclose(eltest, 79.92186504,rtol=1e-2)    #80.0 from vallado, but  his is less precise
+    azapy, elapy = pm.radec2azel(ra,dec,  rdlat, rdlon, t)
+    assert_allclose(azapy, azi, rtol=1e-2)
+    assert_allclose(elapy, eli, rtol=1e-2)
+
+    azvallado, elvallado = vradec2azel(ra, dec, rdlat, rdlon, t)
+    assert_allclose(azvallado, azi, rtol=1e-2)
+    assert_allclose(elvallado, eli, rtol=1e-2)
+
 
 def test_haversine():
-    assert_almost_equal(angledist(35,23,84,20),45.482789587392013)
+    assert_almost_equal(angledist(35,23, 84,20), ha)
     #%% compare with astropy
-    assert_almost_equal(45.482789587392013, angledist_meeus(35,23, 84,20))
+    assert_almost_equal(ha, angledist_meeus(35,23, 84,20))
+
 
 def test_vreckon():
     lat2,lon2,a21 = vreckon(10,20,3000,38)