Repository: astropy/astropy

Issue:
add_index on multi-dimensional masked column failing with index error
<!-- This comments are hidden when you submit the issue,
so you do not need to remove them! -->

<!-- Please be sure to check out our contributing guidelines,
https://github.com/astropy/astropy/blob/main/CONTRIBUTING.md .
Please be sure to check out our code of conduct,
https://github.com/astropy/astropy/blob/main/CODE_OF_CONDUCT.md . -->

<!-- Please have a search on our GitHub repository to see if a similar
issue has already been posted.
If a similar issue is closed, have a quick look to see if you are satisfied
by the resolution.
If not please go ahead and open an issue! -->

<!-- Please check that the development version still produces the same bug.
You can install development version with
pip install git+https://github.com/astropy/astropy
command. -->

### Description
<!-- Provide a general description of the bug. -->
Masked columns with multi-dimensional arrays do not work correctly:

```python
In [127]: t=Table(); t.add_column(np.ma.masked_inside(np.arange(0,6),2,4).reshape((3,-1) ), name='arr')

In [128]: t
Out[128]: 
<Table length=3>
  arr   
int64[2]
--------
  0 .. 1
-- .. --
 -- .. 5

In [129]: t.add_index('arr')
/Users/jdsmith/code/python/astropy/astropy/table/index.py:92: MaskedArrayFutureWarning: In the future the default for argsort will be axis=-1, not the current None, to match its documentation and np.argsort. Explicitly pass -1 or None to silence this warning.
  row_index = Column(col.argsort())
Traceback (most recent call last):
  Input In [129] in <cell line: 1>
    t.add_index('arr')
  File ~/code/python/astropy/astropy/table/table.py:1016 in add_index
    index = Index(columns, engine=engine, unique=unique)
  File ~/code/python/astropy/astropy/table/index.py:93 in __init__
    data = Table([col[row_index]])
  File astropy/table/_column_mixins.pyx:87 in astropy.table._column_mixins._MaskedColumnGetitemShim.__getitem__
  File astropy/table/_column_mixins.pyx:57 in astropy.table._column_mixins.base_getitem
  File astropy/table/_column_mixins.pyx:81 in astropy.table._column_mixins.masked_column_getitem
  File /usr/local/lib/python3.10/site-packages/numpy/ma/core.py:3222 in __getitem__
    dout = self.data[indx]
  File /usr/local/lib/python3.10/site-packages/numpy/ma/core.py:3222 in __getitem__
    dout = self.data[indx]
IndexError: index 5 is out of bounds for axis 0 with size 3
```




### System Details
<!-- Even if you do not think this is necessary, it is useful information for the maintainers.
Please run the following snippet and paste the output below:
import platform; print(platform.platform())
import sys; print("Python", sys.version)
import numpy; print("Numpy", numpy.__version__)
import erfa; print("pyerfa", erfa.__version__)
import astropy; print("astropy", astropy.__version__)
import scipy; print("Scipy", scipy.__version__)
import matplotlib; print("Matplotlib", matplotlib.__version__)
-->
```
macOS-12.3.1-x86_64-i386-64bit
Python 3.10.4 (main, Apr 26 2022, 19:42:59) [Clang 13.1.6 (clang-1316.0.21.2)]
Numpy 1.22.3
pyerfa 2.0.0.1
astropy 5.2.dev92+gf0e2129aa
Scipy 1.7.3
Matplotlib 3.5.2
```

Discussion:
github-actions[bot]: 👋 Thank you for your draft pull request! Do you know that you can use `[ci skip]` or `[skip ci]` in your commit messages to skip running continuous integration tests until you are ready?
taldcroft: This stops the crash but I'm not sure it is creating a functional index. This snippet sets up the data like in the test (but as np.ma.array not MaskedColumn) and then runs some of the index code near the patch. The outputs don't entirely make sense to me.
```
In [20]: data = np.arange(0, 6)
    ...: data_ma = np.ma.masked_inside(data, 2, 4)
    ...: data_ma3 = data_ma.reshape(3, -1)

# I'm not quite sure what this is telling us, but if I understand this should
# be a 1-d array of indices from 0..2.
In [21]: data_ma3.argsort(kind="stable", axis=-1)
Out[21]: 
array([[0, 1],
       [0, 1],
       [1, 0]])

In [22]: row_index = data_ma3.argsort(kind="stable", axis=-1)
In [24]: print(data_ma3[row_index])  # ??
[[[0 1]
  [-- --]]

 [[0 1]
  [-- --]]

 [[-- --]
  [0 1]]]
```

taldcroft: Even for the unmasked case I don't think a multi-dim column can be used for the index. When I tried the index gets created but it seems nonsensical and I could not use it via `.loc`. I think the fix here may be raising a useful exception for the user in that case.
neutrinoceros: Indeed, my patch makes the behavior consistent between masked and unmasked multi-dim indexes, so if the result is nonsensical in one case, it must be the case for the other one too ! I'll admit I didn't pay attention to anything other than consistency since the report was specifically about the masked case and my brain quickly overheats from multi-dim thinking 😅 .

> I think the fix here may be raising a useful exception for the user in that case.

Agreed, let's go with that !


neutrinoceros: @taldcroft do you like it now ? 😸 
neutrinoceros: Needed to rebase to refresh CI. @taldcroft is this good to go ?

Patch:
diff --git a/.github/workflows/ci_workflows.yml b/.github/workflows/ci_workflows.yml
index 5d0287c488..cd57ab5627 100644
--- a/.github/workflows/ci_workflows.yml
+++ b/.github/workflows/ci_workflows.yml
@@ -62,10 +62,8 @@ jobs:
       envs: |
         # NOTE: this coverage test is needed for tests and code that
         #       run only with minimal dependencies.
-        # 2024-05-30: added predeps to ensure new PRs do not break numpy 2.0
-        # docs.  Fine to remove once numpy 2.0 is released (if problematic).
         - name: Python 3.12 with minimal dependencies and full coverage
-          linux: py312-test-cov-predeps
+          linux: py312-test-cov
           coverage: codecov
 
         - name: Python 3.11 in Parallel with all optional dependencies
diff --git a/astropy/io/fits/hdu/hdulist.py b/astropy/io/fits/hdu/hdulist.py
index 5b1f1eec2f..92f53915ed 100644
--- a/astropy/io/fits/hdu/hdulist.py
+++ b/astropy/io/fits/hdu/hdulist.py
@@ -24,6 +24,7 @@ from astropy.io.fits.util import (
 )
 from astropy.io.fits.verify import VerifyError, VerifyWarning, _ErrList, _Verify
 from astropy.utils import indent
+from astropy.utils.compat.numpycompat import NUMPY_LT_2_0
 
 # NOTE: Python can be built without bz2.
 from astropy.utils.compat.optional_deps import HAS_BZ2
@@ -1506,7 +1507,7 @@ class HDUList(list, _Verify):
                     del hdu.data
                 hdu._file = ffo
 
-            if sys.platform.startswith("win"):
+            if sys.platform.startswith("win") and NUMPY_LT_2_0:
                 # On Windows, all the original data mmaps were closed above.
                 # However, it's possible that the user still has references to
                 # the old data which would no longer work (possibly even cause
@@ -1517,6 +1518,8 @@ class HDUList(list, _Verify):
                 # lead to odd behavior in practice.  Better to just not keep
                 # references to data from files that had to be resized upon
                 # flushing (on Windows--again, this is no problem on Linux).
+                # Note that this hack is only possible on numpy 1.x:
+                # in 2.x, we cannot write directly to the data attribute
                 for idx, mmap, arr in mmaps:
                     if mmap is not None:
                         # https://github.com/numpy/numpy/issues/8628
diff --git a/astropy/table/table.py b/astropy/table/table.py
index 02c1bf3d9f..c91f692de9 100644
--- a/astropy/table/table.py
+++ b/astropy/table/table.py
@@ -1086,6 +1086,11 @@ class Table:
                     f'Cannot create an index on column "{col.info.name}", '
                     f'of type "{type(col)}"'
                 )
+            if col.ndim > 1:
+                raise ValueError(
+                    f"Multi-dimensional column {col.info.name!r} "
+                    "cannot be used as an index."
+                )
 
         is_primary = not self.indices
         index = Index(columns, engine=engine, unique=unique)
diff --git a/astropy/time/core.py b/astropy/time/core.py
index c10d110c6e..4d827e19f6 100644
--- a/astropy/time/core.py
+++ b/astropy/time/core.py
@@ -1659,9 +1659,9 @@ class TimeBase(ShapedLikeNDArray):
     def ptp(self, axis=None, out=None, keepdims=False):
         """Peak to peak (maximum - minimum) along a given axis.
 
-        This is similar to :meth:`~numpy.ndarray.ptp`, but adapted to ensure
-        that the full precision given by the two doubles ``jd1`` and ``jd2``
-        is used.
+        This method is similar to the :func:`numpy.ptp` function, but
+        adapted to ensure that the full precision given by the two doubles
+        ``jd1`` and ``jd2`` is used.
 
         Note that the ``out`` argument is present only for compatibility with
         `~numpy.ptp`; since `Time` instances are immutable, it is not possible
diff --git a/astropy/units/format/generic.py b/astropy/units/format/generic.py
index ad735be30d..c7175787aa 100644
--- a/astropy/units/format/generic.py
+++ b/astropy/units/format/generic.py
@@ -19,6 +19,7 @@ from __future__ import annotations
 import re
 import unicodedata
 import warnings
+from contextlib import suppress
 from copy import copy
 from fractions import Fraction
 from typing import TYPE_CHECKING
@@ -28,7 +29,6 @@ from astropy.utils.misc import did_you_mean
 
 from . import core
 from .base import Base
-from .utils import did_you_mean_units, unit_deprecation_warning
 
 if TYPE_CHECKING:
     from astropy.units import NamedUnit, UnitBase
@@ -604,18 +604,55 @@ class Generic(Base):
             if detailed_exception:
                 raise ValueError(
                     f"Unit '{unit}' not supported by the {cls.__name__} standard. "
-                    + did_you_mean_units(
-                        unit,
-                        cls._units,
-                        cls._deprecated_units,
-                        cls._to_decomposed_alternative,
-                    ),
+                    + cls._did_you_mean_units(unit)
                 )
             raise ValueError()
         if unit in cls._deprecated_units:
-            unit_deprecation_warning(
-                unit, cls._units[unit], cls.__name__, cls._to_decomposed_alternative
+            from astropy.units.core import UnitsWarning
+
+            message = (
+                f"The unit '{unit}' has been deprecated in the {cls.__name__} standard."
             )
+            if (decomposed := cls._try_decomposed(cls._units[unit])) is not None:
+                message += f" Suggested: {decomposed}."
+            warnings.warn(message, UnitsWarning)
+
+    @classmethod
+    def _did_you_mean_units(cls, unit: str) -> str:
+        """
+        A wrapper around `astropy.utils.misc.did_you_mean` that deals with
+        the display of deprecated units.
+
+        Parameters
+        ----------
+        unit : str
+            The invalid unit string
+
+        Returns
+        -------
+        msg : str
+            A message with alternatives, or the empty string.
+        """
+
+        def fix_deprecated(x: str) -> list[str]:
+            if x not in cls._deprecated_units:
+                return [x]
+            results = [x + " (deprecated)"]
+            if (decomposed := cls._try_decomposed(cls._units[x])) is not None:
+                results.append(decomposed)
+            return results
+
+        return did_you_mean(unit, cls._units, fix=fix_deprecated)
+
+    @classmethod
+    def _try_decomposed(cls, unit: UnitBase) -> str | None:
+        if (represents := getattr(unit, "_represents", None)) is not None:
+            with suppress(ValueError):
+                return cls._to_decomposed_alternative(represents)
+        if (decomposed := unit.decompose()) is not unit:
+            with suppress(ValueError):
+                return cls._to_decomposed_alternative(decomposed)
+        return None
 
     @classmethod
     def _to_decomposed_alternative(cls, unit: UnitBase) -> str:
diff --git a/astropy/units/format/utils.py b/astropy/units/format/utils.py
index a8384b3d75..e65f61d839 100644
--- a/astropy/units/format/utils.py
+++ b/astropy/units/format/utils.py
@@ -6,22 +6,13 @@ Utilities shared by the different formats.
 
 from __future__ import annotations
 
-import warnings
 from keyword import iskeyword
 from typing import TYPE_CHECKING
 
 from astropy.units.utils import maybe_simple_fraction
-from astropy.utils.misc import did_you_mean
 
 if TYPE_CHECKING:
-    from collections.abc import (
-        Callable,
-        Container,
-        Generator,
-        Iterable,
-        Mapping,
-        Sequence,
-    )
+    from collections.abc import Callable, Generator, Iterable, Sequence
     from numbers import Real
     from typing import TypeVar
 
@@ -153,108 +144,6 @@ def format_power(power: Real) -> str:
     return str(power)
 
 
-def _try_decomposed(
-    unit: UnitBase, format_decomposed: Callable[[UnitBase], str | None]
-) -> str | None:
-    represents = getattr(unit, "_represents", None)
-    if represents is not None:
-        try:
-            represents_string = format_decomposed(represents)
-        except ValueError:
-            pass
-        else:
-            return represents_string
-
-    decomposed = unit.decompose()
-    if decomposed is not unit:
-        try:
-            decompose_string = format_decomposed(decomposed)
-        except ValueError:
-            pass
-        else:
-            return decompose_string
-
-    return None
-
-
-def did_you_mean_units(
-    s: str,
-    all_units: Mapping[str, UnitBase],
-    deprecated_units: Container[str],
-    format_decomposed: Callable[[UnitBase], str | None],
-) -> str:
-    """
-    A wrapper around `astropy.utils.misc.did_you_mean` that deals with
-    the display of deprecated units.
-
-    Parameters
-    ----------
-    s : str
-        The invalid unit string
-
-    all_units : dict
-        A mapping from valid unit names to unit objects.
-
-    deprecated_units : sequence
-        The deprecated unit names
-
-    format_decomposed : callable
-        A function to turn a decomposed version of the unit into a
-        string.  Should return `None` if not possible
-
-    Returns
-    -------
-    msg : str
-        A string message with a list of alternatives, or the empty
-        string.
-    """
-
-    def fix_deprecated(x: str) -> list[str] | tuple[str]:
-        if x in deprecated_units:
-            results = [x + " (deprecated)"]
-            decomposed = _try_decomposed(all_units[x], format_decomposed)
-            if decomposed is not None:
-                results.append(decomposed)
-            return results
-        return (x,)
-
-    return did_you_mean(s, all_units, fix=fix_deprecated)
-
-
-def unit_deprecation_warning(
-    s: str,
-    unit: UnitBase,
-    standard_name: str,
-    format_decomposed: Callable[[UnitBase], str | None],
-) -> None:
-    """
-    Raises a UnitsWarning about a deprecated unit in a given format.
-    Suggests a decomposed alternative if one is available.
-
-    Parameters
-    ----------
-    s : str
-        The deprecated unit name.
-
-    unit : astropy.units.core.UnitBase
-        The unit object.
-
-    standard_name : str
-        The name of the format for which the unit is deprecated.
-
-    format_decomposed : callable
-        A function to turn a decomposed version of the unit into a
-        string.  Should return `None` if not possible
-    """
-    from astropy.units.core import UnitsWarning
-
-    message = f"The unit '{s}' has been deprecated in the {standard_name} standard."
-    decomposed = _try_decomposed(unit, format_decomposed)
-    if decomposed is not None:
-        message += f" Suggested: {decomposed}."
-    warnings.warn(message, UnitsWarning)
-
-
 def get_non_keyword_units(
     bases: Iterable[str], prefixes: Sequence[str]
 ) -> Generator[tuple[str, str], None, None]:
diff --git a/astropy/units/format/vounit.py b/astropy/units/format/vounit.py
index ff0fd63bfb..e806e49de2 100644
--- a/astropy/units/format/vounit.py
+++ b/astropy/units/format/vounit.py
@@ -92,12 +92,9 @@ class VOUnit(generic.Generic):
 
             if cls._custom_unit_regex.match(t.value):
                 warnings.warn(
-                    f"Unit {t.value!r} not supported by the VOUnit standard. "
-                    + utils.did_you_mean_units(
-                        t.value,
-                        cls._units,
-                        cls._deprecated_units,
-                        cls._to_decomposed_alternative,
+                    (
+                        f"Unit {t.value!r} not supported by the VOUnit standard. "
+                        + cls._did_you_mean_units(t.value)
                     ),
                     core.UnitsWarning,
                 )
@@ -108,14 +105,7 @@ class VOUnit(generic.Generic):
 
     @classmethod
     def _parse_unit(cls, unit, detailed_exception=True):
-        if unit not in cls._units:
-            raise ValueError()
-
-        if unit in cls._deprecated_units:
-            utils.unit_deprecation_warning(
-                unit, cls._units[unit], "VOUnit", cls._to_decomposed_alternative
-            )
-
+        super()._validate_unit(unit, detailed_exception=False)
         return cls._units[unit]
 
     @classmethod
diff --git a/docs/changes/io.fits/16581.bugfix.rst b/docs/changes/io.fits/16581.bugfix.rst
new file mode 100644
index 0000000000..41edd54473
--- /dev/null
+++ b/docs/changes/io.fits/16581.bugfix.rst
@@ -0,0 +1,5 @@
+Fixed a crash that occurred for files opened via
+``fits.open(..., mode='update')``, on Windows, and with numpy 2.0 installed.
+It is possible, though unlikely, that this patch reveals memory issues
+for downstream code run under these exact conditions. Please report any
+regression found.
diff --git a/docs/changes/table/16360.api.rst b/docs/changes/table/16360.api.rst
new file mode 100644
index 0000000000..047dba4756
--- /dev/null
+++ b/docs/changes/table/16360.api.rst
@@ -0,0 +1,2 @@
+An exception is now raised when trying to add a multi-dimensional column as an
+index via ``Table.add_index``.
diff --git a/pyproject.toml b/pyproject.toml
index d34063a256..6ab7739f54 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -84,7 +84,8 @@ all = [
     "jplephem",
     "mpmath",
     "asdf-astropy>=0.3",
-    "bottleneck",
+    # https://github.com/astropy/astropy/issues/16574
+    "bottleneck>=1.4.0rc5",
     "fsspec[http]>=2023.4.0",
     "s3fs>=2023.4.0",
 ]
@@ -156,7 +157,7 @@ wcslint = "astropy.wcs.wcslint:main"
 requires = ["setuptools",
             "setuptools_scm>=6.2",
             "cython>=3.0.0,<3.1.0",
-            "numpy>=2.0.0rc1", # see https://github.com/astropy/astropy/issues/16257
+            "numpy>=2.0.0",
             "extension-helpers==1.*"]
 build-backend = "setuptools.build_meta"
 
diff --git a/tox.ini b/tox.ini
index 74677936a4..a4e941ddd7 100644
--- a/tox.ini
+++ b/tox.ini
@@ -61,6 +61,8 @@ deps =
     numpy125: numpy==1.25.*
 
     mpl360: matplotlib==3.6.0
+    # mpl 3.6.0 isn't compatible with numpy 2.0 but didn't include a pin
+    mpl360: numpy<2.0
 
     image: latex
     image: scipy


Test patch:
diff --git a/astropy/coordinates/tests/test_name_resolve.py b/astropy/coordinates/tests/test_name_resolve.py
index 033c5ca26e..a4b4e72e56 100644
--- a/astropy/coordinates/tests/test_name_resolve.py
+++ b/astropy/coordinates/tests/test_name_resolve.py
@@ -4,7 +4,6 @@
 This module contains tests for the name resolve convenience module.
 """
 
-import re
 import time
 import urllib.request
 
@@ -114,10 +113,8 @@ def test_names():
     # "all" choice should ask for SIMBAD, then NED, then Vizier: "SNV" in the url
     with pytest.raises(
         NameResolveError,
-        match=re.escape(
-            "Unable to find coordinates for name 'm87h34hhh' "
-            "using https://cds.unistra.fr/cgi-bin/nph-sesame/SNV?m87h34hhh"
-        ),
+        # avoid hard-coding an exact url as it might depend on external state
+        match=r"Unable to find coordinates for name 'm87h34hhh' using http.*SNV.*",
     ):
         get_icrs_coordinates("m87h34hhh")
 
diff --git a/astropy/table/tests/test_index.py b/astropy/table/tests/test_index.py
index 9153bbff72..ab884fc50e 100644
--- a/astropy/table/tests/test_index.py
+++ b/astropy/table/tests/test_index.py
@@ -604,3 +604,21 @@ def test_hstack_qtable_table():
 def test_index_slice_exception():
     with pytest.raises(TypeError, match="index_slice must be tuple or slice"):
         SlicedIndex(None, None)
+
+
+@pytest.mark.parametrize(
+    "masked",
+    [pytest.param(False, id="raw-array"), pytest.param(True, id="masked array")],
+)
+def test_nd_columun_as_index(masked):
+    # see https://github.com/astropy/astropy/issues/13292
+    # and https://github.com/astropy/astropy/pull/16360
+    t = Table()
+    data = np.arange(0, 6)
+    if masked:
+        data = np.ma.masked_inside(data, 2, 4)
+    t.add_column(data.reshape(3, -1), name="arr")
+    with pytest.raises(
+        ValueError, match="Multi-dimensional column 'arr' cannot be used as an index."
+    ):
+        t.add_index("arr")