Repository: sunpy/sunpy

Issue:
Add a helpful error message to `GenericMap.submap()` for helioprojective-related NaNs
It's a recurring source of confusion for users when they try to submap a helioprojective map from one observatory using 2D helioprojective coordinates from a different observatory (or the same observatory at a different time) and one of the corners is off-disk.  Even if the observer-location difference or the observation-time difference is small, formally there is no unambiguous conversion of an off-disk 2D coordinate to a different helioprojective frame because we don't know what assumption to make the coordinate 3D.  (For an on-disk 2D coordinate, we make the broadly safe assumption that the coordinate is on the solar surface in 3D space.)  The coordinate transformation of the corner point returns NaNs, and `submap()` errors ungracefully.

We should add an error message to `submap()` so that it checks if it is about to work with NaNs and tells the user that they likely want to use `Helioprojective.assume_spherical_screen()` so that the off-disk 2D coordinate can be converted to a 3D coordinate.

Discussion:
nabobalis: TODO:
- [ ] Rebase this PR when https://github.com/sunpy/sunpy/pull/7115 is merged in
ayshih: Here's an example of the common use case:
```python
>>> import astropy.units as u
>>> from astropy.coordinates import SkyCoord

>>> import sunpy.map
>>> from sunpy.data.sample import AIA_171_IMAGE

>>> aiamap = sunpy.map.Map(AIA_171_IMAGE)
```
The following submap that includes off-disk pixels works fine because it is in the native coordinate frame, so 2D coordinates do not need to be promoted:
```python
>>> coord_native = SkyCoord(0*u.arcsec, 0*u.arcsec, frame=aiamap.coordinate_frame)
>>> aiamap.submap(coord_native, width=1000*u.arcsec, height=1000*u.arcsec)
```
The following, very similar submap fails because the observer is slightly different, so 2D coordinates need to be transformed to the native coordinate frame, and thus the off-disk 2D corners get a NaN distance when internally promoted.
```python
>>> coord_other = SkyCoord(0*u.arcsec, 0*u.arcsec, frame='helioprojective', observer='earth', obstime=aiamap.date)
>>> aiamap.submap(coord_other, width=1000*u.arcsec, height=1000*u.arcsec)
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
...
ValueError: cannot convert float NaN to integer
```
The check needs to be on the output of `self._parse_submap_input()` right here:
https://github.com/sunpy/sunpy/blob/e7536b1ddd34eaa9a5f6990627a82f89ab64c6be/sunpy/map/mapbase.py#L1912-L1913

Patch:
diff --git a/changelog/7543.feature.rst b/changelog/7543.feature.rst
new file mode 100644
index 000000000..1d850c102
--- /dev/null
+++ b/changelog/7543.feature.rst
@@ -0,0 +1,1 @@
+Updated :meth:`sunpy.map.GenericMap.submap` to check if it is about to work on locations with NaNs now errors and informs the user that they likely want to use :meth:`~sunpy.coordinates.Helioprojective.assume_spherical_screen` so that the off-disk 2D coordinate can be converted to a 3D coordinate.
diff --git a/sunpy/map/mapbase.py b/sunpy/map/mapbase.py
index 1031872d2..fc602a703 100644
--- a/sunpy/map/mapbase.py
+++ b/sunpy/map/mapbase.py
@@ -1912,6 +1912,17 @@ class GenericMap(NDData):
         pixel_corners = u.Quantity(self._parse_submap_input(
             bottom_left, top_right, width, height)).T
 
+        msg = (
+            "The provided input coordinates to ``submap`` when transformed to the target "
+            "coordinate frame contain NaN values and cannot be used to crop the map. "
+            "The most common reason for NaN values is transforming off-disk 2D "
+            "coordinates without specifying an assumption (e.g., via the"
+            "`Helioprojective.assume_spherical_screen()` context manager) that allows "
+            "such coordinates to be interpreted as 3D coordinates."
+        )
+        if np.any(np.isnan(pixel_corners)):
+            raise ValueError(msg)
+
         # The pixel corners result is in Cartesian order, so the first index is
         # columns and the second is rows.
         bottom = np.min(pixel_corners[1]).to_value(u.pix)


Test patch:
diff --git a/sunpy/map/tests/test_mapbase.py b/sunpy/map/tests/test_mapbase.py
index db08e7ebe..67903dadc 100644
--- a/sunpy/map/tests/test_mapbase.py
+++ b/sunpy/map/tests/test_mapbase.py
@@ -1819,3 +1819,12 @@ def test_plot_deprecated_positional_args(aia171_test_map):
     with pytest.warns(SunpyDeprecationWarning, match=r"Pass annotate=interpolation, axes=True as keyword args."):
         with pytest.raises(TypeError, match="non-boolean value"):
             aia171_test_map.plot('interpolation', True)
+
+
+def test_submap_nan_error(aia171_test_map):
+    # See https://github.com/sunpy/sunpy/pull/7543#issuecomment-2167019208 for more context
+    coord_native = SkyCoord(0*u.arcsec, 0*u.arcsec, frame=aia171_test_map.coordinate_frame)
+    aia171_test_map.submap(coord_native, width=1000*u.arcsec, height=1000*u.arcsec)
+    coord_other = SkyCoord(0*u.arcsec, 0*u.arcsec, frame='helioprojective', observer='earth', obstime=aia171_test_map.date)
+    with pytest.raises(ValueError, match="The provided input coordinates to"):
+        aia171_test_map.submap(coord_other, width=1000*u.arcsec, height=1000*u.arcsec)