fix flux unit conversions in plugins (spacetelescope#3228)

fix flux unit conversions in plugins

CHANGES.rst

@@ -53,6 +53,10 @@ Specviz2d
 Bug Fixes
 ---------
 
+- Fixed broken flux unit conversions in all plugins that respond to changes in flux unit changes. These cases
+  occured when certain flux-to flux-conversions occured, as well as certain conversions between flux and surface
+  brightness. This PR also fixed an issue with unit string formatting in the aperture photometry plugin. [#3228]
+
 Cubeviz
 ^^^^^^^
 - Removed the deprecated ``save as fits`` option from the Collapse, Moment Maps, and Spectral Extraction plugins; use the Export plugin instead. [#3256]

docs/reference/api_nuts_bolts.rst

@@ -79,7 +79,11 @@ Utilities API
    :no-inheritance-diagram:
    :no-inherited-members:
 
-.. automodapi:: jdaviz.core.validunits
+.. automodapi:: jdaviz.core.unit_conversion_utils
+   :no-inheritance-diagram:
+   :no-inherited-members:
+
+.. automodapi:: jdaviz.core.custom_units_and_equivs
    :no-inheritance-diagram:
    :no-inherited-members:
 

jdaviz/app.py

@@ -52,10 +52,10 @@
 from jdaviz.utils import (SnackbarQueue, alpha_index, data_has_valid_wcs, layer_is_table_data,
                           MultiMaskSubsetState, _wcs_only_label, flux_conversion,
                           spectral_axis_conversion)
-from jdaviz.core.custom_units import SPEC_PHOTON_FLUX_DENSITY_UNITS
-from jdaviz.core.validunits import (check_if_unit_is_per_solid_angle,
-                                    combine_flux_and_angle_units,
-                                    supported_sq_angle_units)
+from jdaviz.core.custom_units_and_equivs import SPEC_PHOTON_FLUX_DENSITY_UNITS
+from jdaviz.core.unit_conversion_utils import (check_if_unit_is_per_solid_angle,
+                                               combine_flux_and_angle_units,
+                                               supported_sq_angle_units)
 
 __all__ = ['Application', 'ALL_JDAVIZ_CONFIGS', 'UnitConverterWithSpectral']
 

jdaviz/configs/cubeviz/plugins/moment_maps/moment_maps.py

@@ -17,6 +17,7 @@
                                         SpectralContinuumMixin,
                                         skip_if_no_updates_since_last_active,
                                         with_spinner)
+from jdaviz.core.unit_conversion_utils import convert_integrated_sb_unit
 from jdaviz.core.user_api import PluginUserApi
 
 __all__ = ['MomentMap']
@@ -324,7 +325,40 @@ def calculate_moment(self, add_data=True):
         # convert units for moment 0, which is the only currently supported
         # moment for using converted units.
         if n_moment == 0:
-            self.moment = self.moment.to(self.moment_zero_unit)
+            if self.moment_zero_unit != self.moment.unit:
+                spectral_axis_unit = u.Unit(self.spectrum_viewer.state.x_display_unit)
+
+                # if the flux unit is a per-frequency unit but the spectral axis unit
+                # is a wavelength, or vice versa, we need to convert the spectral axis
+                # unit that the flux was integrated over so they are compatible for
+                # unit conversion (e.g., Jy m / sr needs to become Jy Hz / sr, and
+                # (erg Hz)/(s * cm**2 * AA) needs to become (erg)/(s * cm**2)
+                desired_freq_unit = spectral_axis_unit if spectral_axis_unit.physical_type == 'frequency' else u.Hz  # noqa E501
+                desired_length_unit = spectral_axis_unit if spectral_axis_unit.physical_type == 'length' else u.AA  # noqa E501
+                moment_temp = convert_integrated_sb_unit(self.moment,
+                                                         spectral_axis_unit,
+                                                         desired_freq_unit,
+                                                         desired_length_unit)
+                moment_zero_unit_temp = convert_integrated_sb_unit(1 * u.Unit(self.moment_zero_unit),  # noqa E501
+                                                                   spectral_axis_unit,
+                                                                   desired_freq_unit,
+                                                                   desired_length_unit)
+
+                moment = moment_temp.to(moment_zero_unit_temp.unit, u.spectral())
+
+                # if flux and spectral axis units were incompatible in terms of freq/wav
+                # and needed to be converted to an intermediate unit for conversion, then
+                # re-instate the original chosen units (e.g Jy m /sr was converted to Jy Hz / sr
+                # for unit conversion, now back to Jy m / sr)
+                if spectral_axis_unit not in moment.unit.bases:
+                    if spectral_axis_unit.physical_type == 'frequency':
+                        moment *= (1*desired_length_unit).to(desired_freq_unit,
+                                                             u.spectral()) / desired_length_unit
+                    elif spectral_axis_unit.physical_type == 'length':
+                        moment *= (1*desired_freq_unit).to(desired_length_unit,
+                                                           u.spectral()) / desired_freq_unit
+
+                self.moment = moment
 
         # Reattach the WCS so we can load the result
         self.moment = CCDData(self.moment, wcs=data_wcs)

jdaviz/configs/cubeviz/plugins/moment_maps/tests/test_moment_maps.py

@@ -10,6 +10,8 @@
 from numpy.testing import assert_allclose
 from specutils import SpectralRegion
 
+from jdaviz.core.custom_units_and_equivs import PIX2, SPEC_PHOTON_FLUX_DENSITY_UNITS
+
 
 @pytest.mark.parametrize("cube_type", ["Surface Brightness", "Flux"])
 def test_user_api(cubeviz_helper, spectrum1d_cube, spectrum1d_cube_sb_unit, cube_type):
@@ -350,3 +352,80 @@ def test_correct_output_spectral_y_units(cubeviz_helper, spectrum1d_cube_custom_
 
     mm.calculate_moment()
     assert mm.moment.unit == moment_unit.replace('m', 'um')
+
+
+@pytest.mark.parametrize("flux_unit", [u.Unit(x) for x in SPEC_PHOTON_FLUX_DENSITY_UNITS][1:2])
+@pytest.mark.parametrize("angle_unit", [u.sr, PIX2])
+@pytest.mark.parametrize("new_flux_unit", [u.Unit(x) for x in SPEC_PHOTON_FLUX_DENSITY_UNITS][1:2])
+def test_moment_zero_unit_flux_conversions(cubeviz_helper,
+                                           spectrum1d_cube_custom_fluxunit,
+                                           flux_unit, angle_unit, new_flux_unit):
+    """
+    Test cubeviz moment maps with all possible unit conversions for
+    cubes in spectral/photon surface brightness units (e.g. Jy/sr, Jy/pix2).
+
+    The moment map plugin should respect the choice of flux and angle
+    unit selected in the Unit Conversion plugin, and inputs and results should
+    be converted based on selection. All conversions between units in the
+    flux dropdown menu in the unit conversion plugin should be supported
+    by moment maps.
+    """
+
+    if new_flux_unit == flux_unit:  # skip 'converting' to same unit
+        return
+    new_flux_unit_str = new_flux_unit.to_string()
+
+    cube_unit = flux_unit / angle_unit
+
+    sb_cube = spectrum1d_cube_custom_fluxunit(fluxunit=cube_unit)
+
+    # load surface brigtness cube
+    with warnings.catch_warnings():
+        warnings.filterwarnings("ignore", message="No observer defined on WCS.*")
+        cubeviz_helper.load_data(sb_cube, data_label='test')
+
+    # get plugins
+    uc = cubeviz_helper.plugins["Unit Conversion"]
+    mm = cubeviz_helper.plugins['Moment Maps']._obj
+
+    # and flux viewer for mouseover info
+    flux_viewer = cubeviz_helper.app.get_viewer(cubeviz_helper._default_flux_viewer_reference_name)
+    label_mouseover = cubeviz_helper.app.session.application._tools['g-coords-info']
+
+    # convert to new flux unit
+    uc.flux_unit.selected = new_flux_unit_str
+
+    new_mm_unit = (new_flux_unit * u.m / u.Unit(angle_unit)).to_string()
+    assert mm.output_unit_items[0]['label'] == 'Surface Brightness'
+    assert mm.output_unit_items[0]['unit_str'] == new_mm_unit
+
+    # calculate moment with new output label and plot in flux viewer
+    mm.add_results.label = new_flux_unit_str
+    mm.add_results.viewer.selected = cubeviz_helper._default_flux_viewer_reference_name
+    mm.calculate_moment()
+
+    assert mm.moment.unit == new_mm_unit
+
+    # make sure mouseover info in flux unit is new moment map unit
+    # which should be flux/sb unit times spectral axis unit (e.g. MJy m / sr)
+    label_mouseover._viewer_mouse_event(flux_viewer,
+                                        {'event': 'mousemove',
+                                         'domain': {'x': 0, 'y': 0}})
+    m_orig = label_mouseover.as_text()[0]
+    assert ((new_flux_unit / angle_unit) * u.m).to_string() in m_orig
+
+    # 'jiggle' mouse so we can move it back and compare original coordinate
+    label_mouseover._viewer_mouse_event(flux_viewer,
+                                        {'event': 'mousemove',
+                                         'domain': {'x': 1, 'y': 1}})
+
+    # when flux unit is changed, the mouseover unit conversion should be
+    # skipped so that the plotted moment map remains in its original
+    # unit. setting back to the original flux unit also ensures that
+    # each iteration begins on the same unit so that every comparison
+    # is tested
+    uc.flux_unit.selected = new_flux_unit_str
+    label_mouseover._viewer_mouse_event(flux_viewer,
+                                        {'event': 'mousemove',
+                                         'domain': {'x': 0, 'y': 0}})
+    assert m_orig == label_mouseover.as_text()[0]

jdaviz/configs/cubeviz/plugins/parsers.py

@@ -10,11 +10,10 @@
 from astropy.wcs import WCS
 from specutils import Spectrum1D
 
-from jdaviz.core.custom_units import PIX2
+from jdaviz.core.custom_units_and_equivs import PIX2, _eqv_flux_to_sb_pixel
 from jdaviz.core.registries import data_parser_registry
-from jdaviz.core.validunits import check_if_unit_is_per_solid_angle
-from jdaviz.utils import (standardize_metadata, PRIHDR_KEY, download_uri_to_path,
-                          _eqv_flux_to_sb_pixel)
+from jdaviz.core.unit_conversion_utils import check_if_unit_is_per_solid_angle
+from jdaviz.utils import standardize_metadata, PRIHDR_KEY, download_uri_to_path
 
 __all__ = ['parse_data']
 

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

@@ -20,10 +20,11 @@
                                         skip_if_no_updates_since_last_active,
                                         with_spinner, with_temp_disable)
 from jdaviz.core.user_api import PluginUserApi
-from jdaviz.core.validunits import check_if_unit_is_per_solid_angle
+from jdaviz.core.unit_conversion_utils import (all_flux_unit_conversion_equivs,
+                                               flux_conversion_general,
+                                               check_if_unit_is_per_solid_angle)
 from jdaviz.configs.cubeviz.plugins.parsers import _return_spectrum_with_correct_units
 from jdaviz.configs.cubeviz.plugins.viewers import WithSliceIndicator
-from jdaviz.utils import _eqv_pixar_sr
 
 
 __all__ = ['SpectralExtraction']
@@ -562,9 +563,20 @@ def _preview_x_from_extracted(self, extracted):
         return extracted.spectral_axis
 
     def _preview_y_from_extracted(self, extracted):
-        # TODO: use extracted.meta.get('PIXAR_SR') once populated
-        return extracted.flux.to(self.spectrum_y_units,
-                                 equivalencies=_eqv_pixar_sr(self.dataset.selected_obj.meta.get('PIXAR_SR', 1.0)))  # noqa:
+        """
+        Convert y-axis units of extraction preview to display units,
+        if necessary.
+        """
+
+        if extracted.flux.unit != self.spectrum_y_units:
+
+            eqv = all_flux_unit_conversion_equivs(self.dataset.selected_obj.meta.get('PIXAR_SR', 1.0),  # noqa
+                                                  self.dataset.selected_obj.spectral_axis)
+
+            return flux_conversion_general(extracted.flux.value, extracted.flux.unit,
+                                           self.spectrum_y_units, eqv)
+
+        return extracted.flux
 
     @with_spinner()
     def extract(self, return_bg=False, add_data=True, **kwargs):

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

@@ -8,14 +8,17 @@
 from astropy.table import QTable
 from astropy.tests.helper import assert_quantity_allclose
 from astropy.utils.exceptions import AstropyUserWarning
-
 from glue.core.roi import CircularROI, RectangularROI
 from numpy.testing import assert_allclose, assert_array_equal
 from regions import (CirclePixelRegion, CircleAnnulusPixelRegion, EllipsePixelRegion,
                      RectanglePixelRegion, PixCoord)
 from specutils import Spectrum1D
 from specutils.manipulation import FluxConservingResampler
 
+from jdaviz.core.custom_units_and_equivs import PIX2, SPEC_PHOTON_FLUX_DENSITY_UNITS
+from jdaviz.core.unit_conversion_utils import (all_flux_unit_conversion_equivs,
+                                               flux_conversion_general)
+
 calspec_url = "https://archive.stsci.edu/hlsps/reference-atlases/cdbs/current_calspec/"
 
 
@@ -650,3 +653,68 @@ def test_spectral_extraction_scientific_validation(
         ).to_value(u.dimensionless_unscaled) - 1
     ))
     assert median_abs_relative_dev < expected_rtol
+
+
+@pytest.mark.parametrize("flux_angle_unit", [(u.Unit(x), u.sr) for x in SPEC_PHOTON_FLUX_DENSITY_UNITS]  # noqa
+                                              + [(u.Unit(x), PIX2) for x in SPEC_PHOTON_FLUX_DENSITY_UNITS])  # noqa
+def test_spectral_extraction_flux_unit_conversions(cubeviz_helper,
+                                                   spectrum1d_cube_custom_fluxunit,
+                                                   flux_angle_unit):
+    """
+    Test that cubeviz spectral extraction plugin works with all possible
+    flux unit conversions for a cube loaded in units spectral/photon flux
+    density. The spectral extraction result will remain in the native
+    data unit, but the extraction preview should be converted to the
+    display unit.
+    """
+
+    flux_unit, angle_unit = flux_angle_unit
+
+    sb_cube = spectrum1d_cube_custom_fluxunit(fluxunit=flux_unit / angle_unit,
+                                              shape=(5, 4, 4),
+                                              with_uncerts=True)
+    cubeviz_helper.load_data(sb_cube)
+
+    spectrum_viewer = cubeviz_helper.app.get_viewer(
+        cubeviz_helper._default_spectrum_viewer_reference_name)
+
+    uc = cubeviz_helper.plugins["Unit Conversion"]
+    se = cubeviz_helper.plugins['Spectral Extraction']
+    se.keep_active = True  # keep active for access to preview markers
+
+    # equivalencies for unit conversion, for comparison of outputs
+    equivs = all_flux_unit_conversion_equivs(se.dataset.selected_obj.meta.get('PIXAR_SR', 1.0),
+                                             se.dataset.selected_obj.spectral_axis)
+
+    for new_flux_unit in SPEC_PHOTON_FLUX_DENSITY_UNITS:
+        if new_flux_unit != flux_unit:
+
+            uc.flux_unit.selected = flux_unit.to_string()
+            uc.spectral_y_type.selected = 'Flux'
+
+            # and set back to sum initially so units will always be flux not sb
+            se.function.selected = 'Sum'
+            se.extract()
+
+            original_sum_y_values = se._obj.marks['extract'].y
+
+            # set to new unit
+            uc.flux_unit.selected = new_flux_unit
+            assert spectrum_viewer.state.y_display_unit == new_flux_unit
+
+            # still using 'sum', results should be in flux
+            collapsed = se.extract()
+
+            # make sure extraction preview was translated to new display units
+            new_sum_y_values = se._obj.marks['extract'].y
+            new_converted_to_old_unit = flux_conversion_general(new_sum_y_values,
+                                                                u.Unit(new_flux_unit),
+                                                                u.Unit(flux_unit),
+                                                                equivs, with_unit=False)
+            np.testing.assert_allclose(original_sum_y_values, new_converted_to_old_unit)
+
+            # collapsed result will still have the native data flux unit
+            assert uc.spectral_y_type.selected == 'Flux'
+            assert collapsed.flux.unit == collapsed.uncertainty.unit == flux_unit
+            # but display units in spectrum viewer should reflect new flux unit selection
+            assert se._obj.spectrum_y_units == se._obj.results_units == new_flux_unit