Merge pull request matplotlib#1624 from lpsinger/mollweide_inverse_transform

implemented inverse transform for Mollweide axes

Author: dmcdougall

lib/matplotlib/projections/geo.py

@@ -476,8 +476,16 @@ def __init__(self, resolution):
             self._resolution = resolution
 
         def transform_non_affine(self, xy):
-            # MGDTODO: Math is hard ;(
-            return xy
+            x = xy[:, 0:1]
+            y = xy[:, 1:2]
+
+            # from Equations (7, 8) of
+            # http://mathworld.wolfram.com/MollweideProjection.html
+            theta = np.arcsin(y / np.sqrt(2))
+            lon = (np.pi / (2 * np.sqrt(2))) * x / np.cos(theta)
+            lat = np.arcsin((2 * theta + np.sin(2 * theta)) / np.pi)
+
+            return np.concatenate((lon, lat), 1)
         transform_non_affine.__doc__ = Transform.transform_non_affine.__doc__
 
         def inverted(self):

lib/matplotlib/tests/test_axes.py

@@ -947,6 +947,49 @@ def test_transparent_markers():
     ax = fig.add_subplot(111)
     ax.plot(data, 'D', mfc='none', markersize=100)
 
+@cleanup
+def test_mollweide_forward_inverse_closure():
+    # test that the round-trip Mollweide forward->inverse transformation is an
+    # approximate identity
+    fig = plt.figure()
+    ax = fig.add_subplot(111, projection='mollweide')
+
+    # set up 1-degree grid in longitude, latitude
+    lon = np.linspace(-np.pi, np.pi, 360)
+    lat = np.linspace(-np.pi / 2.0, np.pi / 2.0, 180)
+    lon, lat = np.meshgrid(lon, lat)
+    ll = np.vstack((lon.flatten(), lat.flatten())).T
+
+    # perform forward transform
+    xy = ax.transProjection.transform(ll)
+
+    # perform inverse transform
+    ll2 = ax.transProjection.inverted().transform(xy)
+
+    # compare
+    np.testing.assert_array_almost_equal(ll, ll2, 3)
+
+@cleanup
+def test_mollweide_inverse_forward_closure():
+    # test that the round-trip Mollweide inverse->forward transformation is an
+    # approximate identity
+    fig = plt.figure()
+    ax = fig.add_subplot(111, projection='mollweide')
+
+    # set up grid in x, y
+    x = np.linspace(0, 1, 500)
+    x, y = np.meshgrid(x, x)
+    xy = np.vstack((x.flatten(), y.flatten())).T
+
+    # perform inverse transform
+    ll = ax.transProjection.inverted().transform(xy)
+
+    # perform forward transform
+    xy2 = ax.transProjection.transform(ll)
+
+    # compare
+    np.testing.assert_array_almost_equal(xy, xy2, 3)
+
 if __name__=='__main__':
     import nose
     nose.runmodule(argv=['-s','--with-doctest'], exit=False)