Repository: astropy/reproject

Issue:
Failure w/astropy-dev from spectral-cube: APE14 issue?
I haven't tracked this all the way down yet, but this looks to me like a problem b/w reproject and astropy-dev?

```
   def test_mosaic_cubes(use_memmap, data_adv, use_dask, spectral_block_size):

        pytest.importorskip('reproject')

        # Read in data to use
        cube, data = cube_and_raw(data_adv, use_dask=use_dask)

        # cube is doppler-optical by default, which uses the rest wavelength,
        # which isn't auto-computed, resulting in nan pixels in the WCS transform
        cube._wcs.wcs.restwav = constants.c.to(u.m/u.s).value / cube.wcs.wcs.restfrq

        expected_wcs = WCS(combine_headers(cube.header, cube.header)).celestial

        # Make two overlapping cubes of the data
        part1 = cube[:, :round(cube.shape[1]*2./3.), :]
        part2 = cube[:, round(cube.shape[1]/3.):, :]

        assert part1.wcs.wcs.restwav != 0
        assert part2.wcs.wcs.restwav != 0

>       result = mosaic_cubes([part1, part2], order='nearest-neighbor',
                              roundtrip_coords=False,
                              spectral_block_size=spectral_block_size)

spectral_cube/tests/test_regrid.py:623:
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _
spectral_cube/cube_utils.py:832: in mosaic_cubes
    cube_repr = cube.reproject(header,
spectral_cube/utils.py:49: in wrapper
    return function(self, *args, **kwargs)
spectral_cube/spectral_cube.py:2697: in reproject
    newcube, newcube_valid = reproject_interp((data,
../astropy/astropy/utils/decorators.py:604: in wrapper
    return function(*args, **kwargs)
../python-reprojection/reproject/interpolation/high_level.py:121: in reproject_interp
    return _reproject_blocked(
../python-reprojection/reproject/utils.py:323: in _reproject_blocked
    da.store(
../../mambaforge/envs/py39forge/lib/python3.9/site-packages/dask/array/core.py:1236: in store
    compute_as_if_collection(Array, store_dsk, map_keys, **kwargs)
../../mambaforge/envs/py39forge/lib/python3.9/site-packages/dask/base.py:342: in compute_as_if_collection
    return schedule(dsk2, keys, **kwargs)
../../mambaforge/envs/py39forge/lib/python3.9/site-packages/dask/local.py:557: in get_sync
    return get_async(
../../mambaforge/envs/py39forge/lib/python3.9/site-packages/dask/local.py:500: in get_async
    for key, res_info, failed in queue_get(queue).result():
../../mambaforge/envs/py39forge/lib/python3.9/concurrent/futures/_base.py:439: in result
    return self.__get_result()
../../mambaforge/envs/py39forge/lib/python3.9/concurrent/futures/_base.py:391: in __get_result
    raise self._exception
../../mambaforge/envs/py39forge/lib/python3.9/site-packages/dask/local.py:542: in submit
    fut.set_result(fn(*args, **kwargs))
../../mambaforge/envs/py39forge/lib/python3.9/site-packages/dask/local.py:238: in batch_execute_tasks
    return [execute_task(*a) for a in it]
../../mambaforge/envs/py39forge/lib/python3.9/site-packages/dask/local.py:238: in <listcomp>
    return [execute_task(*a) for a in it]
../../mambaforge/envs/py39forge/lib/python3.9/site-packages/dask/local.py:229: in execute_task
    result = pack_exception(e, dumps)
../../mambaforge/envs/py39forge/lib/python3.9/site-packages/dask/local.py:224: in execute_task
    result = _execute_task(task, data)
../../mambaforge/envs/py39forge/lib/python3.9/site-packages/dask/core.py:119: in _execute_task
    return func(*(_execute_task(a, cache) for a in args))
../../mambaforge/envs/py39forge/lib/python3.9/site-packages/dask/optimization.py:990: in __call__
    return core.get(self.dsk, self.outkey, dict(zip(self.inkeys, args)))
../../mambaforge/envs/py39forge/lib/python3.9/site-packages/dask/core.py:149: in get
    result = _execute_task(task, cache)
../../mambaforge/envs/py39forge/lib/python3.9/site-packages/dask/core.py:119: in _execute_task
    return func(*(_execute_task(a, cache) for a in args))
../../mambaforge/envs/py39forge/lib/python3.9/site-packages/dask/array/core.py:523: in _pass_extra_kwargs
    return func(*args[len(keys) :], **kwargs)
../python-reprojection/reproject/utils.py:285: in reproject_single_block
    array, footprint = reproject_func(
../python-reprojection/reproject/interpolation/core.py:98: in _reproject_full
    _validate_wcs(wcs_in, wcs_out, array.shape, shape_out)
../python-reprojection/reproject/interpolation/core.py:28: in _validate_wcs
    if has_celestial(wcs_in) and not has_celestial(wcs_out):
../python-reprojection/reproject/wcs_utils.py:22: in has_celestial
    for world_axis_class in wcs.low_level_wcs.world_axis_object_classes.values():
../astropy/astropy/wcs/wcsapi/wrappers/sliced_wcs.py:295: in world_axis_object_classes
    keys_keep = [item[0] for item in self.world_axis_object_components]
../astropy/astropy/wcs/wcsapi/wrappers/sliced_wcs.py:291: in world_axis_object_components
    return [self._wcs.world_axis_object_components[idx] for idx in self._world_keep]
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _

.0 = <iterator object at 0x17957bf40>

>   return [self._wcs.world_axis_object_components[idx] for idx in self._world_keep]
E   AttributeError: 'HighLevelWCSWrapper' object has no attribute 'world_axis_object_components'

../astropy/astropy/wcs/wcsapi/wrappers/sliced_wcs.py:291: AttributeError
```

Discussion:
codecov[bot]: ## [Codecov](https://codecov.io/gh/astropy/reproject/pull/349?src=pr&el=h1&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=astropy) Report
> Merging [#349](https://codecov.io/gh/astropy/reproject/pull/349?src=pr&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=astropy) (b7edfa6) into [main](https://codecov.io/gh/astropy/reproject/commit/214544874e838d03ce9b44dacc13621a375df15f?el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=astropy) (2145448) will **decrease** coverage by `1.79%`.
> The diff coverage is `80.00%`.

```diff
@@            Coverage Diff             @@
##             main     #349      +/-   ##
==========================================
- Coverage   92.44%   90.66%   -1.79%     
==========================================
  Files          24       24              
  Lines         821      825       +4     
==========================================
- Hits          759      748      -11     
- Misses         62       77      +15     
```


| [Impacted Files](https://codecov.io/gh/astropy/reproject/pull/349?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=astropy) | Coverage Δ | |
|---|---|---|
| [reproject/interpolation/core.py](https://codecov.io/gh/astropy/reproject/pull/349?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=astropy#diff-cmVwcm9qZWN0L2ludGVycG9sYXRpb24vY29yZS5weQ==) | `82.08% <ø> (ø)` | |
| [reproject/utils.py](https://codecov.io/gh/astropy/reproject/pull/349?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=astropy#diff-cmVwcm9qZWN0L3V0aWxzLnB5) | `84.27% <75.00%> (-0.34%)` | :arrow_down: |
| [reproject/wcs\_utils.py](https://codecov.io/gh/astropy/reproject/pull/349?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=astropy#diff-cmVwcm9qZWN0L3djc191dGlscy5weQ==) | `96.00% <100.00%> (+0.16%)` | :arrow_up: |

... and [2 files with indirect coverage changes](https://codecov.io/gh/astropy/reproject/pull/349/indirect-changes?src=pr&el=tree-more&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=astropy)

:mega: We’re building smart automated test selection to slash your CI/CD build times. [Learn more](https://about.codecov.io/iterative-testing/?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=astropy)

astrofrog: I think the best fix might be to actually make it so that ``has_celestial`` can take a non-FITS WCS high level WCS object - basically we should check in the else clause that if it is a high level WCS object we should extract ``.low_level_wcs`` before proceeding with the check.
keflavich: I'm happy with that approach, I just wasn't sure it was intended.  I'd already prepped a version that worked like that; I'll push it shortly

keflavich: Ah, the weird thing is, what `has_celestial` is getting passed right now is something like `HighLevelWCS(SlicedWCS(HighLevelWCS)))`, which isn't really right
keflavich: My spammed changes will need to be squashed of course
keflavich: @astrofrog I'm going to revert the change you suggested; it is invalid.  If we use `wcs_out = <something>`, it becomes a local variable in the scope of the function `reproject_single_block`, which isn't what we want.

Patch:
diff --git a/reproject/utils.py b/reproject/utils.py
index d79b13b1..a098d095 100644
--- a/reproject/utils.py
+++ b/reproject/utils.py
@@ -277,7 +277,12 @@ def _reproject_blocked(
         if a.ndim == 0 or block_info is None or block_info == []:
             return np.array([a, a])
         slices = [slice(*x) for x in block_info[None]["array-location"][-wcs_out.pixel_n_dim :]]
-        wcs_out_sub = HighLevelWCSWrapper(SlicedLowLevelWCS(wcs_out, slices=slices))
+
+        if isinstance(wcs_out, BaseHighLevelWCS):
+            low_level_wcs = SlicedLowLevelWCS(wcs_out.low_level_wcs, slices=slices)
+        else:
+            low_level_wcs = SlicedLowLevelWCS(wcs_out, slices=slices)
+        wcs_out_sub = HighLevelWCSWrapper(low_level_wcs)
         if isinstance(array_in_or_path, str):
             array_in = np.memmap(array_in_or_path, dtype=float, shape=shape_in)
         else:
diff --git a/reproject/wcs_utils.py b/reproject/wcs_utils.py
index 1a6baf76..9fe71862 100644
--- a/reproject/wcs_utils.py
+++ b/reproject/wcs_utils.py
@@ -8,6 +8,7 @@ import numpy as np
 from astropy.coordinates import SkyCoord
 from astropy.wcs import WCS
 from astropy.wcs.utils import pixel_to_pixel
+from astropy.wcs.wcsapi.high_level_api import BaseHighLevelWCS
 
 __all__ = ["has_celestial", "pixel_to_pixel_with_roundtrip"]
 


Test patch:
diff --git a/reproject/interpolation/tests/test_core.py b/reproject/interpolation/tests/test_core.py
index 2dac25c5..1941eefb 100644
--- a/reproject/interpolation/tests/test_core.py
+++ b/reproject/interpolation/tests/test_core.py
@@ -672,8 +672,9 @@ def test_broadcast_reprojection(input_extra_dims, output_shape, input_as_wcs, ou
 @pytest.mark.parametrize("input_extra_dims", (1, 2))
 @pytest.mark.parametrize("output_shape", (None, "single", "full"))
 @pytest.mark.parametrize("parallel", [True, False])
+@pytest.mark.parametrize("header_or_wcs", (lambda x: x, WCS))
 @pytest.mark.filterwarnings("ignore::astropy.wcs.wcs.FITSFixedWarning")
-def test_blocked_broadcast_reprojection(input_extra_dims, output_shape, parallel):
+def test_blocked_broadcast_reprojection(input_extra_dims, output_shape, parallel, header_or_wcs):
     image_stack, array_ref, footprint_ref, header_in, header_out = _setup_for_broadcast_test()
     # Test both single and multiple dimensions being broadcast
     if input_extra_dims == 2:
@@ -689,6 +690,9 @@ def test_blocked_broadcast_reprojection(input_extra_dims, output_shape, parallel
         # Provide the broadcast dimensions as part of the output shape
         output_shape = image_stack.shape
 
+    # test different behavior when the output projection is a WCS
+    header_out = header_or_wcs(header_out)
+
     array_broadcast, footprint_broadcast = reproject_interp(
         (image_stack, header_in), header_out, output_shape, parallel=parallel, block_size=[5, 5]
     )
@@ -701,9 +705,12 @@ def test_blocked_broadcast_reprojection(input_extra_dims, output_shape, parallel
 @pytest.mark.parametrize("block_size", [[500, 500], [500, 100], None])
 @pytest.mark.parametrize("return_footprint", [False, True])
 @pytest.mark.parametrize("existing_outputs", [False, True])
+@pytest.mark.parametrize("header_or_wcs", (lambda x: x, WCS))
 @pytest.mark.remote_data
 @pytest.mark.filterwarnings("ignore::astropy.wcs.wcs.FITSFixedWarning")
-def test_blocked_against_single(parallel, block_size, return_footprint, existing_outputs):
+def test_blocked_against_single(
+    parallel, block_size, return_footprint, existing_outputs, header_or_wcs
+):
     # Ensure when we break a reprojection down into multiple discrete blocks
     # it has the same result as if all pixels where reprejcted at once
 
@@ -727,7 +734,7 @@ def test_blocked_against_single(parallel, block_size, return_footprint, existing
 
     result_test = reproject_interp(
         hdu2,
-        hdu1.header,
+        header_or_wcs(hdu1.header),
         parallel=parallel,
         block_size=block_size,
         return_footprint=return_footprint,
@@ -737,7 +744,7 @@ def test_blocked_against_single(parallel, block_size, return_footprint, existing
 
     result_reference = reproject_interp(
         hdu2,
-        hdu1.header,
+        header_or_wcs(hdu1.header),
         parallel=False,
         block_size=None,
         return_footprint=return_footprint,
diff --git a/reproject/tests/test_utils.py b/reproject/tests/test_utils.py
index 88c2510d..d842cb94 100644
--- a/reproject/tests/test_utils.py
+++ b/reproject/tests/test_utils.py
@@ -7,6 +7,7 @@ from astropy.wcs import WCS
 
 from reproject.tests.helpers import assert_wcs_allclose
 from reproject.utils import parse_input_data, parse_input_shape, parse_output_projection
+from reproject.wcs_utils import has_celestial
 
 
 @pytest.mark.filterwarnings("ignore:unclosed file:ResourceWarning")
@@ -89,3 +90,22 @@ def test_parse_output_projection_invalid_header(simple_celestial_fits_wcs):
 def test_parse_output_projection_invalid_wcs(simple_celestial_fits_wcs):
     with pytest.raises(ValueError, match="Need to specify shape"):
         parse_output_projection(simple_celestial_fits_wcs)
+
+
+@pytest.mark.filterwarnings("ignore::astropy.utils.exceptions.AstropyUserWarning")
+@pytest.mark.filterwarnings("ignore::astropy.wcs.wcs.FITSFixedWarning")
+def test_has_celestial():
+    from .test_high_level import INPUT_HDR
+
+    hdr = fits.Header.fromstring(INPUT_HDR)
+    ww = WCS(hdr)
+    assert ww.has_celestial
+    assert has_celestial(ww)
+
+    from astropy.wcs.wcsapi import HighLevelWCSWrapper, SlicedLowLevelWCS
+
+    wwh = HighLevelWCSWrapper(SlicedLowLevelWCS(ww, Ellipsis))
+    assert has_celestial(wwh)
+
+    wwh2 = HighLevelWCSWrapper(SlicedLowLevelWCS(ww, [slice(0, 1), slice(0, 1)]))
+    assert has_celestial(wwh2)