convert limits on spectrum y-unit change (spacetelescope#3335)

* remove call to reset limits during unit conversion
* custom logic for mapping y-limits on unit change
* remove delay_callback

---------

Co-authored-by: gibsongreen <gibsongreen99@outlook.com>

CHANGES.rst

@@ -10,6 +10,8 @@ New Features
 
 * API method to toggle API hints. [#3336]
 
+* Changing flux/SB display units no longer resets viewer zoom levels. [#3335]
+
 Cubeviz
 ^^^^^^^
 

jdaviz/configs/cubeviz/plugins/tests/test_cubeviz_unit_conversion.py

@@ -356,3 +356,30 @@ def test_cubeviz_flux_sb_translation_counts(cubeviz_helper, angle_unit):
             f"Pixel x=00010.0 y=00008.0 Value +1.00000e+00 ct / {angle_str}",
             "World 13h39m59.7037s +27d00m03.2400s (ICRS)",
             "204.9987654313 27.0008999946 (deg)")
+
+
+@pytest.mark.parametrize("start_unit, end_unit, end_spectral_y_type, expected_limits",
+                         [('Jy', 'MJy', 'Flux', (5e-07, 6e-07, 1e-4, 1.05e-4)),
+                          ('MJy', 'ph / (Angstrom s cm2)', 'Surface Brightness', (5e-07, 6e-07, 25153169.66070254, 31692993.772485193))  # noqa
+                          ])
+def test_limits_on_unit_change(cubeviz_helper, start_unit, end_unit,
+                               end_spectral_y_type, expected_limits):
+    """
+    Test that the limits are reset when changing units
+    """
+
+    w, wcs_dict = cubeviz_wcs_dict()
+    flux = np.zeros((30, 20, 3001), dtype=np.float32)
+    flux[5:15, 1:11, :] = 1
+    cube = Spectrum1D(flux=flux * u.Unit(start_unit), wcs=w, meta=wcs_dict)
+    cubeviz_helper.load_data(cube, data_label="test")
+
+    uc_plg = cubeviz_helper.plugins['Unit Conversion']
+    sv = cubeviz_helper.viewers['spectrum-viewer']
+    sv.set_limits(x_min=5e-7, x_max=6e-7, y_min=100, y_max=105)
+
+    uc_plg.flux_unit = end_unit
+    uc_plg.spectral_y_type = end_spectral_y_type
+
+    new_limits = sv._obj.get_limits()
+    assert np.allclose(new_limits, expected_limits)

jdaviz/configs/specviz/plugins/unit_conversion/unit_conversion.py

@@ -362,9 +362,6 @@ def _handle_spectral_y_unit(self, *args):
             pass
         else:
             self.spectrum_viewer.set_plot_axes()
-            # until we can have upstream automatic limit updating on change
-            # in display units with equivalencies, we'll reset the limits
-            self.spectrum_viewer.reset_limits()
 
             # broadcast that there has been a change in the spectrum viewer y axis,
             self.hub.broadcast(GlobalDisplayUnitChanged('spectral_y',

jdaviz/core/freezable_state.py

@@ -2,11 +2,12 @@
 from echo import delay_callback, CallbackProperty
 import numpy as np
 
+from astropy import units as u
 from glue.viewers.profile.state import ProfileViewerState
 from glue_jupyter.bqplot.image.state import BqplotImageViewerState
 from glue.viewers.matplotlib.state import DeferredDrawCallbackProperty as DDCProperty
 
-from jdaviz.utils import get_reference_image_data
+from jdaviz.utils import get_reference_image_data, flux_conversion
 
 __all__ = ['FreezableState', 'FreezableProfileViewerState', 'FreezableBqplotImageViewerState']
 
@@ -53,6 +54,37 @@ def _reset_x_limits(self, *event):
             self.x_min = x_min
             self.x_max = x_max
 
+    def _convert_units_y_limits(self, old_unit, new_unit):
+        # override glue's _convert_units_y_limits to account
+        # for equivalencies.  This converts all four corners
+        # of the limits to set new limits that contain those
+        # same corners
+
+        if old_unit != new_unit:
+            if old_unit is None or new_unit is None:
+                self._reset_y_limits()
+                return
+            x_corners = np.array([self.x_min, self.x_min, self.x_max, self.x_max])
+            y_corners = np.array([self.y_min, self.y_max, self.y_min, self.y_max])
+            spectral_axis = x_corners * u.Unit(self.x_display_unit)
+
+            # NOTE: this uses the scale factor from the first found layer.  We may want to
+            # generalize this to iterate over all scale factors if/when we support multiple
+            # flux cubes (with potentially different pixel scale factors).
+            for layer in self.layers:
+                if psc := getattr(layer.layer, 'meta', {}).get('_pixel_scale_factor', None):  # noqa
+                    spectral_axis.info.meta = {'_pixel_scale_factor',
+                                               psc}
+                    break
+            else:
+                spectral_axis.info.meta = {}
+
+            y_corners_new = flux_conversion(y_corners, old_unit, new_unit, spectral_axis=spectral_axis)  # noqa
+
+            with delay_callback(self, 'y_min', 'y_max'):
+                self.y_min = np.nanmin(y_corners_new)
+                self.y_max = np.nanmax(y_corners_new)
+
 
 class FreezableBqplotImageViewerState(BqplotImageViewerState, FreezableState):
     linked_by_wcs = False

jdaviz/tests/test_utils.py

@@ -58,7 +58,7 @@ def test_spec_sb_flux_conversion():
                                                         targ_units=(u.ph / (u.s * u.cm**2 * u.Hz * solid_angle)),  # noqa
                                                         eqv=eqv,
                                                         spec_unit=spec_unit,
-                                                        image_data=None
+                                                        indirect_needs_spec_axis=None
                                                     )
         assert_allclose(returned_values, target_values)
         assert (return_units == expected_units)
@@ -71,7 +71,7 @@ def test_spec_sb_flux_conversion():
                                                     values=values,
                                                     orig_units=(u.ph / (u.Angstrom * u.s * u.cm**2 * u.sr)),  # noqa
                                                     targ_units=(u.MJy), eqv=eqv,
-                                                    spec_unit=spec_unit, image_data=None
+                                                    spec_unit=spec_unit, indirect_needs_spec_axis=None  # noqa
                                                 )
     assert_allclose(returned_values, target_values)
     assert (return_units == expected_units)
@@ -83,7 +83,7 @@ def test_spec_sb_flux_conversion():
     returned_values, return_units, unit_flag = _indirect_conversion(
                                                     values=values, orig_units=(u.Jy/u.sr),
                                                     targ_units=(u.MJy), eqv=eqv,
-                                                    spec_unit=spec_unit, image_data=None
+                                                    spec_unit=spec_unit, indirect_needs_spec_axis=None  # noqa
                                                 )
     assert_allclose(returned_values, target_values)
     assert (return_units == expected_units)
@@ -95,7 +95,7 @@ def test_spec_sb_flux_conversion():
     returned_values, return_units = _indirect_conversion(
                                         values=1, orig_units=(u.MJy/u.sr),
                                         targ_units=(u.erg / (u.s * u.cm**2 * u.Hz * u.sr)),
-                                        eqv=eqv, spec_unit=None, image_data=True
+                                        eqv=eqv, spec_unit=None, indirect_needs_spec_axis=True
                                     )
     assert_allclose(returned_values, target_value)
     assert return_units == expected_units
@@ -105,7 +105,7 @@ def test_spec_sb_flux_conversion():
     expected_units = (u.MJy / u.sr)
     returned_values, return_units = _indirect_conversion(
                                         values=10, orig_units=(u.MJy/u.sr), targ_units=(u.Jy/u.sr),
-                                        eqv=eqv, spec_unit=None, image_data=True
+                                        eqv=eqv, spec_unit=None, indirect_needs_spec_axis=True
                                     )
     assert_allclose(returned_values, target_value)
     assert return_units == expected_units

jdaviz/utils.py

@@ -342,7 +342,7 @@ def indirect_units():
     return units
 
 
-def flux_conversion(values, original_units, target_units, spec=None, eqv=None, slice=None):
+def flux_conversion(values, original_units, target_units, spec=None, eqv=None, spectral_axis=None):
     """
     Convert flux or surface brightness values from original units to target units.
 
@@ -367,8 +367,8 @@ def flux_conversion(values, original_units, target_units, spec=None, eqv=None, s
     eqv : list of :ref:`astropy:unit_equivalencies`, optional
         A list of Astropy equivalencies necessary for complex unit conversions/translations.
 
-    slice : `astropy.units.Quantity`, optional
-        The current slice of a data cube, with units. Necessary for complex unit
+    spectral_axis : `astropy.units.Quantity`, optional
+        Provide spectral axis associated with values array. Necessary for simple and complex unit
         conversions/translations that require spectral density equivalencies.
 
     Returns
@@ -383,7 +383,7 @@ def flux_conversion(values, original_units, target_units, spec=None, eqv=None, s
     # If there are only two values, this is likely the limits being converted, so then
     # in case we need to use the spectral density equivalency, we need to provide only
     # to spectral axis values. If there is only one value
-    image_data = False
+    indirect_needs_spec_axis = False
     if spec:
         if not np.isscalar(values) and len(values) == 2:
             spectral_values = spec.spectral_axis[0]
@@ -398,10 +398,18 @@ def flux_conversion(values, original_units, target_units, spec=None, eqv=None, s
             eqv = u.spectral_density(spectral_values)
         else:
             eqv = u.spectral_density(spec.spectral_axis[0])
-    elif slice is not None and eqv:
-        image_data = True
-        # Need this to convert Flux to Flux for complex conversions/translations of cube image data
-        eqv += u.spectral_density(slice)
+    elif spectral_axis is not None:
+        if not isinstance(spectral_axis, u.Quantity):
+            raise TypeError("spectral_axis must be an instance of astropy.units.Quantity.")
+
+        if eqv:
+            indirect_needs_spec_axis = True
+            # Needed to convert Flux to Flux for complex conversion/translation of non-spectrum data
+            eqv += u.spectral_density(spectral_axis)
+        elif not np.isscalar(values) and len(values) != spectral_axis.size:
+            raise ValueError(f"Size mismatch: values size : ({values.size}) must match spectral_axis size : ({spectral_axis.size}).")  # noqa
+        else:
+            eqv = u.spectral_density(spectral_axis)
 
     orig_units = u.Unit(original_units)
     targ_units = u.Unit(target_units)
@@ -410,7 +418,9 @@ def flux_conversion(values, original_units, target_units, spec=None, eqv=None, s
     solid_angle_in_targ = check_if_unit_is_per_solid_angle(targ_units, return_unit=True)
 
     # Ensure a spectrum passed through Spectral Extraction plugin
-    if (((spec and ('_pixel_scale_factor' in spec.meta))) and
+    if ((((spec and ('_pixel_scale_factor' in spec.meta))) or
+        (spectral_axis is not None and ('_pixel_scale_factor' in spectral_axis.info.meta)))
+        and
             (((solid_angle_in_orig) and (not solid_angle_in_targ)) or
              ((not solid_angle_in_orig) and (solid_angle_in_targ)))):
         # Data item in data collection does not update from conversion/translation.
@@ -419,7 +429,13 @@ def flux_conversion(values, original_units, target_units, spec=None, eqv=None, s
         n_values = len(values)
 
         # Make sure they are float (can be Quantity).
-        fac = spec.meta['_pixel_scale_factor']
+        if spec:
+            fac = spec.meta['_pixel_scale_factor']
+        elif spectral_axis is not None:
+            # Quantity.info.meta gets casted to set
+            fac = next((item for item in spectral_axis.info.meta if isinstance(item, float)), None)
+            spec_unit = orig_units
+
         if isinstance(fac, Quantity):
             fac = fac.value
 
@@ -454,10 +470,10 @@ def flux_conversion(values, original_units, target_units, spec=None, eqv=None, s
             elif selected_unit_updated == 'orig':
                 orig_units = updated_units
 
-    elif image_data:
+    elif indirect_needs_spec_axis:
         values, orig_units = _indirect_conversion(
             values=values, orig_units=orig_units, targ_units=targ_units,
-            eqv=eqv, image_data=image_data
+            eqv=eqv, indirect_needs_spec_axis=indirect_needs_spec_axis
             )
     elif solid_angle_in_orig == solid_angle_in_targ == PIX2:
         # in the case where we have 2 SBs per solid pixel that need
@@ -473,7 +489,7 @@ def flux_conversion(values, original_units, target_units, spec=None, eqv=None, s
 
 
 def _indirect_conversion(values, orig_units, targ_units, eqv,
-                         spec_unit=None, image_data=None):
+                         spec_unit=None, indirect_needs_spec_axis=None):
 
     # Note: is there a way we could write this to not require 'spec_unit'? It
     # seems like it falls back on this to get a solid angle unit, but can we
@@ -502,7 +518,7 @@ def _indirect_conversion(values, orig_units, targ_units, eqv,
 
         return values, targ_units, 'targ'
 
-    elif image_data or (spec_unit and solid_angle_in_spec):
+    elif indirect_needs_spec_axis or (spec_unit and solid_angle_in_spec):
         if not solid_angle_in_targ:
             targ_units /= solid_angle_in_spec
             solid_angle_in_targ = solid_angle_in_spec