BUG: fix wrong masked median for some special cases

the masked nans which are equivalent to valid infs must be replaced
with infs earlier otherwise the inf is lost in the masked sum of the low
and high part.
Closes gh-8487

3 files changed, 47 insertions, 6 deletions

diff --git a/numpy/lib/tests/test_nanfunctions.py b/numpy/lib/tests/test_nanfunctions.py
index e1942338b..2b310457b 100644
--- a/numpy/lib/tests/test_nanfunctions.py
+++ b/numpy/lib/tests/test_nanfunctions.py
@@ -710,6 +710,18 @@ class TestNanFunctions_Median(TestCase):
                 a = np.array([[inf, inf], [inf, inf]])
                 assert_equal(np.nanmedian(a, axis=1), inf)
 
+                a = np.array([[inf, 7, -inf, -9],
+                              [-10, np.nan, np.nan, 5],
+                              [4, np.nan, np.nan, inf]],
+                              dtype=np.float32)
+                if inf > 0:
+                    assert_equal(np.nanmedian(a, axis=0), [4., 7., -inf, 5.])
+                    assert_equal(np.nanmedian(a), 4.5)
+                else:
+                    assert_equal(np.nanmedian(a, axis=0), [-10., 7., -inf, -9.])
+                    assert_equal(np.nanmedian(a), -2.5)
+                assert_equal(np.nanmedian(a, axis=-1), [-1., -2.5, inf])
+
                 for i in range(0, 10):
                     for j in range(1, 10):
                         a = np.array([([np.nan] * i) + ([inf] * j)] * 2)
diff --git a/numpy/ma/extras.py b/numpy/ma/extras.py
index 1774ece30..0a60ea331 100644
--- a/numpy/ma/extras.py
+++ b/numpy/ma/extras.py
@@ -758,6 +758,19 @@ def _median(a, axis=None, out=None, overwrite_input=False):
     ind[axis] = np.minimum(h, asorted.shape[axis] - 1)
     high = asorted[tuple(ind)]
 
+    def replace_masked(s):
+        # Replace masked entries with minimum_full_value unless it all values
+        # are masked. This is required as the sort order of values equal or
+        # larger than the fill value is undefined and a valid value placed
+        # elsewhere, e.g. [4, --, inf].
+        if np.ma.is_masked(s):
+            rep = (~np.all(asorted.mask, axis=axis)) & s.mask
+            s.data[rep] = np.ma.minimum_fill_value(asorted)
+            s.mask[rep] = False
+
+    replace_masked(low)
+    replace_masked(high)
+
     # duplicate high if odd number of elements so mean does nothing
     odd = counts % 2 == 1
     np.copyto(low, high, where=odd)
@@ -776,12 +789,6 @@ def _median(a, axis=None, out=None, overwrite_input=False):
     else:
         s = np.ma.mean([low, high], axis=0, out=out)
 
-    # if result is masked either the input contained enough minimum_fill_value
-    # so that it would be the median or all values masked
-    if np.ma.is_masked(s):
-        rep = (~np.all(asorted.mask, axis=axis)) & s.mask
-        s.data[rep] = np.ma.minimum_fill_value(asorted)
-        s.mask[rep] = False
     return s
 
 
diff --git a/numpy/ma/tests/test_extras.py b/numpy/ma/tests/test_extras.py
index fb16d92ce..58ac46f53 100644
--- a/numpy/ma/tests/test_extras.py
+++ b/numpy/ma/tests/test_extras.py
@@ -667,6 +667,15 @@ class TestMedian(TestCase):
         r = np.ma.median(np.ma.masked_array([[np.inf, np.inf],
                                              [np.inf, np.inf]]), axis=None)
         assert_equal(r, np.inf)
+        # all masked
+        r = np.ma.median(np.ma.masked_array([[np.inf, np.inf],
+                                             [np.inf, np.inf]], mask=True),
+                         axis=-1)
+        assert_equal(r.mask, True)
+        r = np.ma.median(np.ma.masked_array([[np.inf, np.inf],
+                                             [np.inf, np.inf]], mask=True),
+                         axis=None)
+        assert_equal(r.mask, True)
 
     def test_non_masked(self):
         x = np.arange(9)
@@ -992,6 +1001,19 @@ class TestMedian(TestCase):
             assert_equal(np.ma.median(a, axis=0), inf)
             assert_equal(np.ma.median(a, axis=1), inf)
 
+            a = np.array([[inf, 7, -inf, -9],
+                          [-10, np.nan, np.nan, 5],
+                          [4, np.nan, np.nan, inf]],
+                          dtype=np.float32)
+            a = np.ma.masked_array(a, mask=np.isnan(a))
+            if inf > 0:
+                assert_equal(np.ma.median(a, axis=0), [4., 7., -inf, 5.])
+                assert_equal(np.ma.median(a), 4.5)
+            else:
+                assert_equal(np.ma.median(a, axis=0), [-10., 7., -inf, -9.])
+                assert_equal(np.ma.median(a), -2.5)
+            assert_equal(np.ma.median(a, axis=1), [-1., -2.5, inf])
+
             for i in range(0, 10):
                 for j in range(1, 10):
                     a = np.array([([np.nan] * i) + ([inf] * j)] * 2)