eteq · Jun 15, 2016 · Jun 13, 2016 · Jun 13, 2016 · Jun 13, 2016 · Jun 13, 2016
diff --git a/astropy/coordinates/builtin_frames/ecliptic_transforms.py b/astropy/coordinates/builtin_frames/ecliptic_transforms.py
@@ -9,6 +9,7 @@
 import numpy as np
 
 from ... import units as u
+from ...utils.compat import NUMPY_LT_1_10
 from ..baseframe import frame_transform_graph
 from ..transformations import FunctionTransform, DynamicMatrixTransform
 from ..angles import rotation_matrix
@@ -26,7 +27,22 @@ def _ecliptic_rotation_matrix(equinox):
     jd1, jd2 = get_jd12(equinox, 'tt')
     rnpb = erfa.pnm06a(jd1, jd2)
     obl = erfa.obl06(jd1, jd2)*u.radian
-    return np.asarray(np.dot(rotation_matrix(obl, 'x'), rnpb))
+    """
+    The following code is the equivalent of looping over obl and rnpb,
+    creating a rotation matrix from obl and then taking
+    the dot product of the resulting matrices, finally combining
+    into a new array.
+    """
+    try:
+        rmat = np.array([rotation_matrix(this_obl, 'x') for this_obl in obl])
+        if NUMPY_LT_1_10:
+            result = np.einsum('...ij,...jk->...ik', rmat, rnpb)
+        else:
+            result = np.matmul(rmat, rnpb)
+    except:
+        # must be a scalar obliquity
+        result = np.asarray(np.dot(rotation_matrix(obl, 'x'), rnpb))
+    return result
 
 
 @frame_transform_graph.transform(FunctionTransform, GCRS, GeocentricTrueEcliptic)
@@ -45,7 +61,7 @@ def geoecliptic_to_gcrs(from_coo, gcrs_frame):
     newrepr = cartrepr_from_matmul(rmat, from_coo, transpose=True)
     gcrs = GCRS(newrepr, obstime=from_coo.equinox)
 
-    #now do any needed offsets (no-op if same obstime and 0 pos/vel)
+    # now do any needed offsets (no-op if same obstime and 0 pos/vel)
     return gcrs.transform_to(gcrs_frame)
 
 
@@ -64,6 +80,7 @@ def baryecliptic_to_icrs(from_coo, to_frame):
                      "no distance.  Heliocentric<->ICRS transforms cannot "
                      "function in this case because there is an origin shift.")
 
+
 @frame_transform_graph.transform(FunctionTransform, ICRS, HeliocentricTrueEcliptic)
 def icrs_to_helioecliptic(from_coo, to_frame):
     if not u.m.is_equivalent(from_coo.cartesian.x.unit):
@@ -72,7 +89,7 @@ def icrs_to_helioecliptic(from_coo, to_frame):
     pvh, pvb = erfa.epv00(*get_jd12(to_frame.equinox, 'tdb'))
     delta_bary_to_helio = pvh[..., 0, :] - pvb[..., 0, :]
 
-    #first offset to heliocentric
+    # first offset to heliocentric
     heliocart = (from_coo.cartesian.xyz).T + delta_bary_to_helio * u.au
 
     # now compute the matrix to precess to the right orientation

diff --git a/astropy/coordinates/solar_system.py b/astropy/coordinates/solar_system.py
@@ -214,6 +214,8 @@ def get_body_barycentric(body, time, ephemeris=None):
             cartesian_position_body = earth_pv_bary[..., 0, :]
         elif body == 'moon':
             cartesian_position_body = calc_moon(time).cartesian.xyz.to(u.au).value
+            cartesian_position_body = np.rollaxis(cartesian_position_body, 0,
+                                                  cartesian_position_body.ndim)
         else:
             sun_bary = earth_pv_bary[..., 0, :] - earth_pv_helio[..., 0, :]
             if body == 'sun':

diff --git a/astropy/coordinates/tests/test_regression.py b/astropy/coordinates/tests/test_regression.py
@@ -12,14 +12,40 @@
 import numpy as np
 
 from ... import units as u
-from .. import AltAz, EarthLocation, SkyCoord, get_sun, ICRS, CIRS, ITRS
+from .. import (AltAz, EarthLocation, SkyCoord, get_sun, ICRS, CIRS, ITRS,
+                GeocentricTrueEcliptic, Longitude, Latitude, GCRS)
 from ...time import Time
 from ...utils import iers
 
 from ...tests.helper import pytest, assert_quantity_allclose
 from .test_matching import HAS_SCIPY, OLDER_SCIPY
 
 
+def test_regression_5085():
+    """
+    PR #5085 was put in place to fix the following issue.
+
+    Issue: https://github.com/astropy/astropy/issues/5069
+    At root was the transformation of Ecliptic coordinates with
+    non-scalar times.
+    """
+    times = Time(["2015-08-28 03:30", "2015-09-05 10:30", "2015-09-15 18:35"])
+    latitudes = Latitude([3.9807075, -5.00733806, 1.69539491]*u.deg)
+    longitudes = Longitude([311.79678613,  72.86626741, 199.58698226]*u.deg)
+    distances = u.Quantity([0.00243266, 0.0025424, 0.00271296]*u.au)
+    coo = GeocentricTrueEcliptic(lat=latitudes,
+                                 lon=longitudes,
+                                 distance=distances, equinox=times)
+    # expected result
+    ras = Longitude([310.50095387, 314.67109863, 319.56507471]*u.deg)
+    decs = Latitude([-18.25190707, -17.1556641, -15.71616651]*u.deg)
+    distances = u.Quantity([1.78309902, 1.710874, 1.61326648]*u.au)
+    expected_result = GCRS(ra=ras, dec=decs,
+                           distance=distances, obstime="J2000").cartesian.xyz
+    actual_result = coo.transform_to(GCRS(obstime="J2000")).cartesian.xyz
+    assert_quantity_allclose(expected_result, actual_result)
+
+
 def test_regression_3920():
     """
     Issue: https://github.com/astropy/astropy/issues/3920
@@ -121,12 +147,12 @@ def test_regression_4082():
     Issue: https://github.com/astropy/astropy/issues/4082
     """
     from .. import search_around_sky, search_around_3d
-    cat = SkyCoord([10.076,10.00455], [18.54746, 18.54896], unit='deg')
+    cat = SkyCoord([10.076, 10.00455], [18.54746, 18.54896], unit='deg')
     search_around_sky(cat[0:1], cat, seplimit=u.arcsec * 60, storekdtree=False)
     # in the issue, this raises a TypeError
 
-    #also check 3d for good measure, although it's not really affected by this bug directly
-    cat3d = SkyCoord([10.076,10.00455]*u.deg, [18.54746, 18.54896]*u.deg, distance=[0.1,1.5]*u.kpc)
+    # also check 3d for good measure, although it's not really affected by this bug directly
+    cat3d = SkyCoord([10.076, 10.00455]*u.deg, [18.54746, 18.54896]*u.deg, distance=[0.1, 1.5]*u.kpc)
     search_around_3d(cat3d[0:1], cat3d, 1*u.kpc, storekdtree=False)
 
 

diff --git a/astropy/coordinates/tests/test_solar_system.py b/astropy/coordinates/tests/test_solar_system.py
@@ -291,3 +291,13 @@ def test_get_sun_consistency(time):
     builtin_get_sun = get_sun(time)
     sep = builtin_get_sun.separation(sun_jpl_gcrs)
     assert sep < 0.1*u.arcsec
+
+def test_get_moon_nonscalar_regression():
+    """
+    Test that the builtin ephemeris works with non-scalar times.
+
+    See Issue #5069.
+    """
+    times = Time(["2015-08-28 03:30", "2015-09-05 10:30"])
+    # the following line will raise an Exception if the bug recurs.
+    get_moon(times, ephemeris='builtin')