Repository: sunpy/ndcube

Issue:
Detect Invariant Axes as part of NDCube.reproject_to
<!--
We know asking good questions takes effort, and we appreciate your time.
Thank you.

Please be aware that everyone has to follow our code of conduct:
https://sunpy.org/coc

These comments are hidden when you submit this github issue.

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!
-->

### Description

<!--
Provide a general description of the feature you would like.
If you prefer, you can also suggest a draft design or API.
-->
After an initial minimal implementation of `NDCube.reproject_to` has been made to work, the method needs to be able to detect invariant axes in the WCSs, i.e. axes over which no reprojection is to be performed.  This will open up a range of possibilities, such as more efficient reproject, supporting more reprojection algorithms, etc.

Discussion:

Patch:
diff --git a/changelog/459.trivial.rst b/changelog/459.trivial.rst
new file mode 100644
index 0000000..5c767af
--- /dev/null
+++ b/changelog/459.trivial.rst
@@ -0,0 +1,1 @@
+Adds a function to identify invariant axes between two WCS objects.
diff --git a/ndcube/utils/wcs.py b/ndcube/utils/wcs.py
index 228fc24..9154cdd 100644
--- a/ndcube/utils/wcs.py
+++ b/ndcube/utils/wcs.py
@@ -8,6 +8,7 @@ import numbers
 from collections import UserDict
 
 import numpy as np
+from astropy.wcs.utils import pixel_to_pixel
 from astropy.wcs.wcsapi import BaseHighLevelWCS, BaseLowLevelWCS, low_level_api
 
 __all__ = ['array_indices_for_world_objects', 'convert_between_array_and_pixel_axes',
@@ -477,3 +478,38 @@ def compare_wcs_physical_types(source_wcs, target_wcs):
     target_wcs = get_low_level_wcs(target_wcs, 'target_wcs')
 
     return source_wcs.world_axis_physical_types == target_wcs.world_axis_physical_types
+
+
+def identify_invariant_axes(source_wcs, target_wcs, input_shape, atol=1e-6, rtol=1e-6):
+    """
+    Performs a pixel to pixel transformation to identify if there are any invariant axes
+    between the given source and target WCS objects.
+
+    Parameters
+    ----------
+    source_wcs: `astropy.wcs.wcsapi.BaseHighLevelWCS` or `astropy.wcs.wcsapi.BaseLowLevelWCS`
+
+    target_wcs: `astropy.wcs.wcsapi.BaseHighLevelWCS` or `astropy.wcs.wcsapi.BaseLowLevelWCS`
+
+    input_shape: `tuple`
+        The array shape of the data.
+
+    atol: `float`
+        The absolute tolerance parameter for comparison.
+
+    rtol: `float`
+        The relative tolerance parameter for comparison.
+
+    Returns
+    -------
+    result: `list`
+        A list of booleans denoting whether the axis is invariant or not.
+        Follows the WCS ordering.
+    """
+
+    input_pixel_coords = np.meshgrid(*[np.arange(n) for n in input_shape])
+
+    output_pixel_coords = pixel_to_pixel(source_wcs, target_wcs, *input_pixel_coords)
+
+    return [np.allclose(input_coord, output_coord, atol=atol, rtol=rtol)
+            for input_coord, output_coord in zip(input_pixel_coords, output_pixel_coords)]


Test patch:
diff --git a/ndcube/utils/tests/test_utils_wcs.py b/ndcube/utils/tests/test_utils_wcs.py
index 6d811a8..d9730b9 100644
--- a/ndcube/utils/tests/test_utils_wcs.py
+++ b/ndcube/utils/tests/test_utils_wcs.py
@@ -162,3 +162,19 @@ def test_array_indices_for_world_objects_2(wcs_4d_lt_t_l_ln):
 def test_compare_wcs_physical_types(wcs_4d_t_l_lt_ln, wcs_3d_l_lt_ln):
     assert utils.wcs.compare_wcs_physical_types(wcs_4d_t_l_lt_ln, wcs_4d_t_l_lt_ln) is True
     assert utils.wcs.compare_wcs_physical_types(wcs_4d_t_l_lt_ln, wcs_3d_l_lt_ln) is False
+
+
+def test_identify_invariant_axes(wcs_3d_l_lt_ln):
+    source_wcs = wcs_3d_l_lt_ln
+
+    target_wcs_header = wcs_3d_l_lt_ln.low_level_wcs.to_header().copy()
+    target_wcs_header['CDELT2'] = 10
+    target_wcs_header['CDELT3'] = 20
+    target_wcs = WCS(header=target_wcs_header)
+
+    invariant_axes = utils.wcs.identify_invariant_axes(source_wcs, target_wcs, (4, 4, 4))
+    assert invariant_axes == [True, False, False]
+
+    invariant_axes = utils.wcs.identify_invariant_axes(source_wcs, target_wcs, (4, 4, 4),
+                                                       atol=1e-20, rtol=1e-20)
+    assert invariant_axes == [False, False, False]