Fix "Unit conversion when equivalency needed breaks contours in cubev…

…iz" (spacetelescope#3207)

* TST: Unskip test_basic_unit_conversions
to make sure it still fails before I fix it.

* Fix bug

Data might not be selected it seems

Handle multiselect case

* BUG: Fix spectral extraction wrong for PIX2.

TST: Validate actual values in Clare tests.

jdaviz/configs/cubeviz/plugins/spectral_extraction/spectral_extraction.py

@@ -1,12 +1,11 @@
 import os
+from functools import cached_property
 from pathlib import Path
 
 import numpy as np
 import astropy
-from astropy.nddata import (
-    NDDataArray, StdDevUncertainty
-)
-from functools import cached_property
+from astropy import units as u
+from astropy.nddata import NDDataArray, StdDevUncertainty
 from traitlets import Any, Bool, Dict, Float, List, Unicode, observe
 
 from jdaviz.core.custom_traitlets import FloatHandleEmpty
@@ -522,7 +521,10 @@ def _extract_from_aperture(self, cube, uncert_cube, aperture,
                 # NOTE: just forcing these units for now!! this is in steradians and
                 # needs to be converted to the selected square angle unit but for now just
                 # force to correct units
-                aperture_area = self.cube.meta.get('PIXAR_SR', 1.0) * sq_angle_unit
+                if sq_angle_unit == u.sr:
+                    aperture_area = self.cube.meta.get('PIXAR_SR', 1.0) * sq_angle_unit
+                else:
+                    aperture_area = 1 * sq_angle_unit
                 collapsed_nddata = collapsed_nddata.multiply(aperture_area,
                                                              propagate_uncertainties=True)
         else:

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

@@ -2,7 +2,8 @@
 import pytest
 from astropy import units as u
 from astropy.wcs import WCS
-from regions import PixCoord, CirclePixelRegion
+from numpy.testing import assert_allclose
+from regions import PixCoord, CirclePixelRegion, RectanglePixelRegion
 from specutils import Spectrum1D
 
 from jdaviz.core.custom_units import PIX2, SPEC_PHOTON_FLUX_DENSITY_UNITS
@@ -18,7 +19,6 @@ def cubeviz_wcs_dict():
     return w, wcs_dict
 
 
-@pytest.mark.skip(reason="Unskip after JDAT 4785 resolved.")
 @pytest.mark.parametrize("angle_unit", [u.sr, PIX2])
 def test_basic_unit_conversions(cubeviz_helper, angle_unit):
     """
@@ -35,19 +35,92 @@ def test_basic_unit_conversions(cubeviz_helper, angle_unit):
 
     # load cube with flux units of MJy
     w, wcs_dict = cubeviz_wcs_dict()
-    flux = np.zeros((3, 4, 5), dtype=np.float32)
-    cube = Spectrum1D(flux=flux * u.MJy / angle_unit, wcs=w, meta=wcs_dict)
+    flux = np.ones((3, 4, 5), dtype=np.float32)
+    cube = Spectrum1D(flux=flux * (u.MJy / angle_unit), wcs=w, meta=wcs_dict)
     cubeviz_helper.load_data(cube, data_label="test")
+    viewer = cubeviz_helper.app.get_viewer("spectrum-viewer")
 
     # get all available flux units for translation. Since cube is loaded
     # in Jy, this will be all items in 'spectral_and_photon_flux_density_units'
 
     uc_plg = cubeviz_helper.plugins['Unit Conversion']
+    ap_plg = cubeviz_helper.plugins["Aperture Photometry"]._obj
+    label_mouseover = cubeviz_helper.app.session.application._tools['g-coords-info']
+
+    cubeviz_helper.load_regions(RectanglePixelRegion(PixCoord(1, 1), 1, 1))
+    ap_plg.background_selected = "Subset 1"
 
     for flux_unit in SPEC_PHOTON_FLUX_DENSITY_UNITS:
+        if flux_unit == 'MJy':
+            if angle_unit == u.sr:
+                ans = 9.6e-10  # 12 * 8e-11
+            else:
+                ans = 12  # 12 * 1
+            bg_ans = 1  # MJy / angle_unit
+        elif flux_unit == 'Jy':
+            if angle_unit == u.sr:
+                ans = 9.6e-04
+            else:
+                ans = 1.2e+07
+            bg_ans = 1e6
+        elif flux_unit == 'mJy':
+            if angle_unit == u.sr:
+                ans = 0.96
+            else:
+                ans = 1.2e+10
+            bg_ans = 1e9
+        elif flux_unit == 'uJy':
+            if angle_unit == u.sr:
+                ans = 960
+            else:
+                ans = 1.2e+13
+            bg_ans = 1e12
+        elif flux_unit == 'W / (Hz m2)':
+            if angle_unit == u.sr:
+                ans = 9.6e-30
+            else:
+                ans = 1.2e-19
+            bg_ans = 1e-20
+        elif flux_unit == 'eV / (Hz s m2)':
+            if angle_unit == u.sr:
+                ans = 5.99185e-11
+            else:
+                ans = 7.48981e-01
+            bg_ans = 0.06241509074460764
+        elif flux_unit == 'erg / (Angstrom s cm2)':
+            if angle_unit == u.sr:
+                ans = 1.34580e-15
+            else:
+                ans = 1.68225e-05
+            bg_ans = 1.4018163962192981e-06
+        elif flux_unit == 'erg / (Hz s cm2)':
+            if angle_unit == u.sr:
+                ans = 9.6e-27
+            else:
+                ans = 1.2e-16
+            bg_ans = 1e-17
+        elif flux_unit == 'ph / (Angstrom s cm2)':
+            if angle_unit == u.sr:
+                ans = 3.133e-04
+            else:
+                ans = 3.91624e+06
+            bg_ans = 326346.6709140777
+        elif flux_unit == 'ph / (Hz s cm2)':
+            if angle_unit == u.sr:
+                ans = 2.23486e-15
+            else:
+                ans = 2.79357e-05
+            bg_ans = 2.328027206660126e-06
+
         uc_plg.flux_unit = flux_unit
         assert cubeviz_helper.app._get_display_unit('spectral_y') == flux_unit
 
+        label_mouseover._viewer_mouse_event(viewer, {
+            'event': 'mousemove', 'domain': {'x': 4.6245e-7, 'y': ans}})
+        assert_allclose(label_mouseover._dict['value'], ans, rtol=1e-4)
+
+        assert_allclose(ap_plg.background_value, bg_ans, rtol=1e-4)
+
 
 @pytest.mark.parametrize("flux_unit, expected_choices", [(u.count, ['ct']),
                                                          (u.Jy, SPEC_PHOTON_FLUX_DENSITY_UNITS),

jdaviz/configs/imviz/plugins/aper_phot_simple/aper_phot_simple.py

@@ -24,7 +24,7 @@
                                         SubsetSelect, ApertureSubsetSelectMixin,
                                         TableMixin, PlotMixin, MultiselectMixin, with_spinner)
 from jdaviz.core.validunits import check_if_unit_is_per_solid_angle
-from jdaviz.utils import PRIHDR_KEY
+from jdaviz.utils import PRIHDR_KEY, _eqv_flux_to_sb_pixel, _eqv_pixar_sr
 
 __all__ = ['SimpleAperturePhotometry']
 
@@ -205,15 +205,32 @@ def _on_display_units_changed(self, event={}):
                 # convert background to new unit
                 if self.background_value is not None:
 
-                    prev_unit = u.Unit(prev_display_unit)
-                    new_unit = u.Unit(self.display_unit)
+                    # FIXME: Kinda hacky, better solution welcomed.
+                    # Basically we want to forget about PIX2 and do normal conversion.
+                    # Since traitlet does not keep unit, we do not need to reapply PIX2.
+                    if "pix2" in prev_display_unit and "pix2" in self.display_unit:
+                        prev_unit = u.Unit(prev_display_unit) * PIX2
+                        new_unit = u.Unit(self.display_unit) * PIX2
+                    else:
+                        prev_unit = u.Unit(prev_display_unit)
+                        new_unit = u.Unit(self.display_unit)
 
                     bg = self.background_value * prev_unit
 
-                    # will need to add additonal custom equiv here once
-                    # pix2<>angle is enabled
-                    self.background_value = bg.to_value(
-                        new_unit, u.spectral_density(self._cube_wave))
+                    if self.multiselect:
+                        if len(self.dataset.selected) == 1:
+                            data = self.dataset.selected_dc_item[0]
+                        else:
+                            raise ValueError("cannot calculate background median in multiselect")
+                    else:
+                        data = self.dataset.selected_dc_item
+
+                    with u.set_enabled_equivalencies(
+                            u.spectral() + u.spectral_density(self._cube_wave) +
+                            _eqv_flux_to_sb_pixel() +
+                            _eqv_pixar_sr(data.meta.get('_pixel_scale_factor', 1)
+                                          if data else 1)):
+                        self.background_value = bg.to_value(new_unit)
 
                 # convert flux scaling to new unit
                 if self.flux_scaling is not None:

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

@@ -1,12 +1,14 @@
-from astropy import units as u
+from contextlib import nullcontext
 from functools import cached_property
+
+from astropy import units as u
 from glue.core.subset_group import GroupedSubset
 from glue_jupyter.bqplot.image import BqplotImageView
 from specutils import Spectrum1D
 from traitlets import List, Unicode, observe, Bool
 
 from jdaviz.configs.default.plugins.viewers import JdavizProfileView
-from jdaviz.core.events import GlobalDisplayUnitChanged, AddDataMessage
+from jdaviz.core.events import GlobalDisplayUnitChanged, AddDataMessage, SliceValueUpdatedMessage
 from jdaviz.core.registries import tray_registry
 from jdaviz.core.template_mixin import (PluginTemplateMixin, UnitSelectPluginComponent,
                                         SelectPluginComponent, PluginUserApi)
@@ -15,6 +17,7 @@
                                     check_if_unit_is_per_solid_angle,
                                     create_angle_equivalencies_list,
                                     supported_sq_angle_units)
+from jdaviz.utils import _eqv_flux_to_sb_pixel, _eqv_pixar_sr
 
 __all__ = ['UnitConversion']
 
@@ -120,6 +123,8 @@ def __init__(self, *args, **kwargs):
 
         self.session.hub.subscribe(self, AddDataMessage,
                                    handler=self._on_add_data_to_viewer)
+        self.session.hub.subscribe(self, SliceValueUpdatedMessage,
+                                   handler=self._on_slice_changed)
 
         self.has_spectral = self.config in ('specviz', 'cubeviz', 'specviz2d', 'mosviz')
         self.spectral_unit = UnitSelectPluginComponent(self,
@@ -275,6 +280,11 @@ def _on_add_data_to_viewer(self, msg):
             self._handle_attribute_display_unit(self.sb_unit_selected, layers=layers)
             self._clear_cache('image_layers')
 
+    def _on_slice_changed(self, msg):
+        if self.config != "cubeviz":
+            return
+        self._cube_wave = u.Quantity(msg.value, msg.value_unit)
+
     @observe('spectral_unit_selected', 'flux_unit_selected',
              'angle_unit_selected', 'sb_unit_selected',
              'time_unit_selected')
@@ -376,6 +386,13 @@ def _handle_attribute_display_unit(self, attr_unit, layers=None):
                     or u.Unit(layer.layer.get_component("flux").units).physical_type != 'surface brightness'):  # noqa
                 continue
             if hasattr(layer.state, 'attribute_display_unit'):
-                layer.state.attribute_display_unit = _valid_glue_display_unit(attr_unit,
-                                                                              layer,
-                                                                              'attribute')
+                if self.config == "cubeviz":
+                    ctx = u.set_enabled_equivalencies(
+                        u.spectral() + u.spectral_density(self._cube_wave) +
+                        _eqv_flux_to_sb_pixel() +
+                        _eqv_pixar_sr(layer.layer.meta.get('_pixel_scale_factor', 1)))
+                else:
+                    ctx = nullcontext()
+                with ctx:
+                    layer.state.attribute_display_unit = _valid_glue_display_unit(
+                        attr_unit, layer, 'attribute')

jdaviz/utils.py

@@ -389,8 +389,11 @@ def flux_conversion(values, original_units, target_units, spec=None, eqv=None, s
         # the unit of the data collection item object, could be flux or surface brightness
         spec_unit = str(spec.flux.unit)
 
-        # Need this for setting the y-limits
-        eqv = u.spectral_density(spectral_values)
+        # Need this for setting the y-limits but values from viewer might be downscaled
+        if len(values) == spectral_values.size:
+            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