From 9373a9eb310cec97e2168ab8a4003db680033934 Mon Sep 17 00:00:00 2001 From: "P. L. Lim" <2090236+pllim@users.noreply.github.com> Date: Fri, 2 Feb 2024 13:18:33 -0500 Subject: [PATCH] Aperture photometry in Cubeviz (#2666) * Implement Cubeviz aperture photometry. Add tests and documentation. Attach WCS to some Cubeviz 2D products, where it is needed but missing. Minor clean-ups. * Add slice info to aper phot plugin for cubeviz * Clarify slice field in aper phot Co-authored-by: Kyle Conroy * Fix plugin data menu bug * Implement slice in output table and do not display slice for cubeviz 2D data, as requested. * Add doc for new slice column for Cubeviz phot table. * Fix table GUI for slice wavelength and display slice wavelength instead of index in read-only text field. * Display slice_wave in plot and text results --------- Co-authored-by: Kyle Conroy --- CHANGES.rst | 2 + docs/cubeviz/export_data.rst | 23 + docs/cubeviz/plugins.rst | 13 + jdaviz/components/plugin_dataset_select.vue | 8 +- jdaviz/configs/cubeviz/cubeviz.yaml | 1 + jdaviz/configs/cubeviz/helper.py | 14 + .../plugins/moment_maps/moment_maps.py | 10 +- .../moment_maps/tests/test_moment_maps.py | 7 +- jdaviz/configs/cubeviz/plugins/slice/slice.py | 1 + .../plugins/tests/test_cubeviz_aperphot.py | 184 +++++ .../cubeviz/plugins/tests/test_parsers.py | 2 +- .../default/plugins/collapse/collapse.py | 22 +- .../plugins/collapse/tests/test_collapse.py | 11 +- .../gaussian_smooth/gaussian_smooth.py | 17 +- .../aper_phot_simple/aper_phot_simple.py | 118 +++- .../aper_phot_simple/aper_phot_simple.vue | 14 + jdaviz/conftest.py | 15 +- jdaviz/core/events.py | 3 +- .../cubeviz_aperture_photometry.ipynb | 630 ++++++++++++++++++ 19 files changed, 1046 insertions(+), 49 deletions(-) create mode 100644 jdaviz/configs/cubeviz/plugins/tests/test_cubeviz_aperphot.py create mode 100644 notebooks/concepts/cubeviz_aperture_photometry.ipynb diff --git a/CHANGES.rst b/CHANGES.rst index 73590fb05b..223531c779 100644 --- a/CHANGES.rst +++ b/CHANGES.rst @@ -26,6 +26,8 @@ Cubeviz - Spectral extraction plugin re-organized into subsections to be more consistent with specviz2d. [#2676] +- New aperture photometry plugin that can perform aperture photometry on selected cube slice. [#2666] + Imviz ^^^^^ diff --git a/docs/cubeviz/export_data.rst b/docs/cubeviz/export_data.rst index 11c439e7ad..fb860fd915 100644 --- a/docs/cubeviz/export_data.rst +++ b/docs/cubeviz/export_data.rst @@ -133,3 +133,26 @@ Markers Table All mouseover information in the :ref:`markers plugin ` can be exported to an :ref:`astropy table ` by calling :meth:`~jdaviz.core.template_mixin.TableMixin.export_table` (see :ref:`plugin-apis`). + + +.. _cubeviz_export_photometry: + +Aperture Photometry +=================== + +Cubeviz can export photometry output table like Imviz: + +.. code-block:: python + + results = cubeviz.get_aperture_photometry_results() + +.. seealso:: + + :ref:`Imviz Aperture Photometry ` + Imviz documentation describing exporting of aperture photometry results in Jdaviz. + +In addition to the columns that :ref:`Imviz Aperture Photometry ` +would provide, the table from Cubeviz has this extra column after ``data_label``: + +* ``slice_wave``: Wavelength value at the selected slice of the cube used for computation. + If a 2D data (e.g., collapsed cube) is selected, the value would be NaN instead. diff --git a/docs/cubeviz/plugins.rst b/docs/cubeviz/plugins.rst index 46a57b2560..57f2f13db5 100644 --- a/docs/cubeviz/plugins.rst +++ b/docs/cubeviz/plugins.rst @@ -291,6 +291,19 @@ Click :guilabel:`EXTRACT` to produce a new 1D spectrum dataset from the spectral cube, which has uncertainties propagated by `astropy.nddata `_. +.. _cubeviz-aper-phot: + +Aperture Photometry +=================== + +Cubeviz allows aperture photometry on some 3D and 2D data, as long as they +have valid flux units. For 3D data, the current :ref:`slice` is used. + +.. seealso:: + + :ref:`Imviz Aperture Photometry ` + Imviz documentation describing the concept of aperture photometry in Jdaviz. + .. _cubeviz-export-plot: Export Plot diff --git a/jdaviz/components/plugin_dataset_select.vue b/jdaviz/components/plugin_dataset_select.vue index 64f3de685c..e0e73635dd 100644 --- a/jdaviz/components/plugin_dataset_select.vue +++ b/jdaviz/components/plugin_dataset_select.vue @@ -20,7 +20,7 @@
- + {{ data.item.label }} @@ -64,12 +64,6 @@ diff --git a/jdaviz/configs/cubeviz/cubeviz.yaml b/jdaviz/configs/cubeviz/cubeviz.yaml index 274a4abb75..ade379f001 100644 --- a/jdaviz/configs/cubeviz/cubeviz.yaml +++ b/jdaviz/configs/cubeviz/cubeviz.yaml @@ -30,6 +30,7 @@ tray: - specviz-line-analysis - cubeviz-moment-maps - cubeviz-spectral-extraction + - imviz-aper-phot-simple - g-export-plot viewer_area: - container: col diff --git a/jdaviz/configs/cubeviz/helper.py b/jdaviz/configs/cubeviz/helper.py index 2846cfe008..f059bb87e6 100644 --- a/jdaviz/configs/cubeviz/helper.py +++ b/jdaviz/configs/cubeviz/helper.py @@ -168,6 +168,20 @@ def get_data(self, data_label=None, spatial_subset=None, spectral_subset=None, f spectral_subset=spectral_subset, function=function, cls=cls, use_display_units=use_display_units) + # Need this method for Imviz Aperture Photometry plugin. + + def get_aperture_photometry_results(self): + """Return aperture photometry results, if any. + Results are calculated using :ref:`cubeviz-aper-phot` plugin. + + Returns + ------- + results : `~astropy.table.QTable` or `None` + Photometry results if available or `None` otherwise. + + """ + return self.plugins['Aperture Photometry']._obj.export_table() + def layer_is_cube_image_data(layer): return isinstance(layer, BaseData) and layer.ndim in (2, 3) diff --git a/jdaviz/configs/cubeviz/plugins/moment_maps/moment_maps.py b/jdaviz/configs/cubeviz/plugins/moment_maps/moment_maps.py index 3dd70e9a69..6db49c68f5 100644 --- a/jdaviz/configs/cubeviz/plugins/moment_maps/moment_maps.py +++ b/jdaviz/configs/cubeviz/plugins/moment_maps/moment_maps.py @@ -239,7 +239,15 @@ def calculate_moment(self, add_data=True): # Need transpose to align JWST mirror shape: This is because specutils # arrange the array shape to be (nx, ny, nz) but 2D visualization # assumes (ny, nx) as per row-major convention. - data_wcs = getattr(cube.wcs, 'celestial', None) + + # Extract 2D WCS from input cube. + data = self.dataset.selected_dc_item + # Similar to coords_info logic. + if '_orig_spec' in getattr(data, 'meta', {}): + w = data.meta['_orig_spec'].wcs + else: + w = data.coords + data_wcs = getattr(w, 'celestial', None) if data_wcs: data_wcs = data_wcs.swapaxes(0, 1) # We also transpose WCS to match. diff --git a/jdaviz/configs/cubeviz/plugins/moment_maps/tests/test_moment_maps.py b/jdaviz/configs/cubeviz/plugins/moment_maps/tests/test_moment_maps.py index 8c66bf7ec1..c9ac3599af 100644 --- a/jdaviz/configs/cubeviz/plugins/moment_maps/tests/test_moment_maps.py +++ b/jdaviz/configs/cubeviz/plugins/moment_maps/tests/test_moment_maps.py @@ -30,6 +30,7 @@ def test_user_api(cubeviz_helper, spectrum1d_cube): assert len(mm._obj.continuum_marks['center'].x) > 0 mom_sub = mm.calculate_moment() + assert isinstance(mom_sub.wcs, WCS) assert mom != mom_sub @@ -45,7 +46,7 @@ def test_user_api(cubeviz_helper, spectrum1d_cube): mm.calculate_moment() -def test_moment_calculation(cubeviz_helper, spectrum1d_cube, tmpdir): +def test_moment_calculation(cubeviz_helper, spectrum1d_cube, tmp_path): dc = cubeviz_helper.app.data_collection with warnings.catch_warnings(): warnings.filterwarnings("ignore", message="No observer defined on WCS.*") @@ -110,7 +111,7 @@ def test_moment_calculation(cubeviz_helper, spectrum1d_cube, tmpdir): "204.9998877673 27.0001000000 (deg)") assert mm._obj.filename == 'moment0_test_FLUX.fits' # Auto-populated on calculate. - mm._obj.filename = str(tmpdir.join(mm._obj.filename)) # But we want it in tmpdir for testing. + mm._obj.filename = str(tmp_path / mm._obj.filename) # But we want it in tmp_path for testing. mm._obj.vue_save_as_fits() assert os.path.isfile(mm._obj.filename) @@ -134,7 +135,7 @@ def test_moment_calculation(cubeviz_helper, spectrum1d_cube, tmpdir): "204.9998877673 27.0001000000 (deg)") -def test_moment_velocity_calculation(cubeviz_helper, spectrum1d_cube, tmpdir): +def test_moment_velocity_calculation(cubeviz_helper, spectrum1d_cube): with warnings.catch_warnings(): warnings.filterwarnings("ignore", message="No observer defined on WCS.*") cubeviz_helper.load_data(spectrum1d_cube, data_label='test') diff --git a/jdaviz/configs/cubeviz/plugins/slice/slice.py b/jdaviz/configs/cubeviz/plugins/slice/slice.py index 1b98e6458c..54b57246fe 100644 --- a/jdaviz/configs/cubeviz/plugins/slice/slice.py +++ b/jdaviz/configs/cubeviz/plugins/slice/slice.py @@ -220,6 +220,7 @@ def _on_slider_updated(self, event): self.hub.broadcast(SliceWavelengthUpdatedMessage(slice=value, wavelength=self.wavelength, + wavelength_unit=self.wavelength_unit, sender=self)) def vue_goto_first(self, *args): diff --git a/jdaviz/configs/cubeviz/plugins/tests/test_cubeviz_aperphot.py b/jdaviz/configs/cubeviz/plugins/tests/test_cubeviz_aperphot.py new file mode 100644 index 0000000000..891b069cab --- /dev/null +++ b/jdaviz/configs/cubeviz/plugins/tests/test_cubeviz_aperphot.py @@ -0,0 +1,184 @@ +import numpy as np +import pytest +from astropy import units as u +from astropy.tests.helper import assert_quantity_allclose +from astropy.utils.exceptions import AstropyUserWarning +from numpy.testing import assert_allclose +from regions import RectanglePixelRegion, PixCoord + + +def test_cubeviz_aperphot_cube_orig_flux(cubeviz_helper, image_cube_hdu_obj_microns): + cubeviz_helper.load_data(image_cube_hdu_obj_microns, data_label="test") + flux_unit = u.Unit("1E-17 erg*s^-1*cm^-2*Angstrom^-1") + + aper = RectanglePixelRegion(center=PixCoord(x=1, y=2), width=3, height=5) + cubeviz_helper.load_regions(aper) + + # Make sure MASK is not an option even when shown in viewer. + cubeviz_helper.app.add_data_to_viewer("flux-viewer", "test[MASK]", visible=True) + + plg = cubeviz_helper.plugins["Aperture Photometry"]._obj + assert plg.dataset.labels == ["test[FLUX]", "test[ERR]"] + assert plg.cube_slice == "4.894e+00 um" + + plg.dataset_selected = "test[FLUX]" + plg.aperture_selected = "Subset 1" + plg.vue_do_aper_phot() + row = cubeviz_helper.get_aperture_photometry_results()[0] + + # Basically, we should recover the input rectangle here. + assert_allclose(row["xcenter"], 1 * u.pix) + assert_allclose(row["ycenter"], 2 * u.pix) + sky = row["sky_center"] + assert_allclose(sky.ra.deg, 205.43985906934287) + assert_allclose(sky.dec.deg, 27.003490103642033) + assert_allclose(row["sum"], 75 * flux_unit) # 3 (w) x 5 (h) x 5 (v) + assert_allclose(row["sum_aper_area"], 15 * (u.pix * u.pix)) # 3 (w) x 5 (h) + assert_allclose(row["mean"], 5 * flux_unit) + assert_quantity_allclose(row["slice_wave"], 4.894499866699333 * u.um) + + # Move slider and make sure it recomputes for a new slice automatically. + cube_slice_plg = cubeviz_helper.plugins["Slice"]._obj + cube_slice_plg.slice = 0 + plg.vue_do_aper_phot() + row = cubeviz_helper.get_aperture_photometry_results()[1] + + # Same rectangle but different slice value. + assert_allclose(row["xcenter"], 1 * u.pix) + assert_allclose(row["ycenter"], 2 * u.pix) + sky = row["sky_center"] + assert_allclose(sky.ra.deg, 205.43985906934287) + assert_allclose(sky.dec.deg, 27.003490103642033) + assert_allclose(row["sum"], 15 * flux_unit) # 3 (w) x 5 (h) x 1 (v) + assert_allclose(row["sum_aper_area"], 15 * (u.pix * u.pix)) # 3 (w) x 5 (h) + assert_allclose(row["mean"], 1 * flux_unit) + assert_quantity_allclose(row["slice_wave"], 4.8904998665093435 * u.um) + + # We continue on with test_cubeviz_aperphot_generated_2d_collapse here + # because we want to make sure the result would append properly between 3D and 2D. + collapse_plg = cubeviz_helper.plugins["Collapse"]._obj + collapse_plg.vue_collapse() + + # Need this to make it available for photometry data drop-down. + cubeviz_helper.app.add_data_to_viewer("uncert-viewer", "test[FLUX] collapsed") + + plg = cubeviz_helper.plugins["Aperture Photometry"]._obj + plg.dataset_selected = "test[FLUX] collapsed" + plg.aperture_selected = "Subset 1" + plg.vue_do_aper_phot() + row = cubeviz_helper.get_aperture_photometry_results()[2] + + # Basically, we should recover the input rectangle here. + assert_allclose(row["xcenter"], 1 * u.pix) + assert_allclose(row["ycenter"], 2 * u.pix) + sky = row["sky_center"] + assert_allclose(sky.ra.deg, 205.43985906934287) + assert_allclose(sky.dec.deg, 27.003490103642033) + assert_allclose(row["sum"], 540 * flux_unit) # 3 (w) x 5 (h) x 36 (v) + assert_allclose(row["sum_aper_area"], 15 * (u.pix * u.pix)) # 3 (w) x 5 (h) + assert_allclose(row["mean"], 36 * flux_unit) + assert np.isnan(row["slice_wave"]) + + +def test_cubeviz_aperphot_generated_2d_moment(cubeviz_helper, image_cube_hdu_obj_microns): + cubeviz_helper.load_data(image_cube_hdu_obj_microns, data_label="test") + flux_unit = u.Unit("1E-17 erg*s^-1*cm^-2*Angstrom^-1") + + moment_plg = cubeviz_helper.plugins["Moment Maps"] + _ = moment_plg.calculate_moment() + + # Need this to make it available for photometry data drop-down. + cubeviz_helper.app.add_data_to_viewer("uncert-viewer", "test[FLUX] moment 0") + + aper = RectanglePixelRegion(center=PixCoord(x=1, y=2), width=3, height=5) + cubeviz_helper.load_regions(aper) + + plg = cubeviz_helper.plugins["Aperture Photometry"]._obj + plg.dataset_selected = "test[FLUX] moment 0" + plg.aperture_selected = "Subset 1" + plg.vue_do_aper_phot() + row = cubeviz_helper.get_aperture_photometry_results()[0] + + # Basically, we should recover the input rectangle here. + assert_allclose(row["xcenter"], 1 * u.pix) + assert_allclose(row["ycenter"], 2 * u.pix) + sky = row["sky_center"] + assert_allclose(sky.ra.deg, 205.43985906934287) + assert_allclose(sky.dec.deg, 27.003490103642033) + assert_allclose(row["sum"], 540 * flux_unit) # 3 (w) x 5 (h) x 36 (v) + assert_allclose(row["sum_aper_area"], 15 * (u.pix * u.pix)) # 3 (w) x 5 (h) + assert_allclose(row["mean"], 36 * flux_unit) + assert np.isnan(row["slice_wave"]) + + # Moment 1 has no compatible unit, so should not be available for photometry. + moment_plg.n_moment = 1 + moment_plg.reference_wavelength = 5 + _ = moment_plg.calculate_moment() + m1_lbl = "test[FLUX] moment 1" + cubeviz_helper.app.add_data_to_viewer("uncert-viewer", m1_lbl) + assert (m1_lbl in cubeviz_helper.app.data_collection.labels and + m1_lbl not in plg.dataset.choices) + + +def test_cubeviz_aperphot_generated_3d_gaussian_smooth(cubeviz_helper, image_cube_hdu_obj_microns): + cubeviz_helper.load_data(image_cube_hdu_obj_microns, data_label="test") + flux_unit = u.Unit("1E-17 erg*s^-1*cm^-2*Angstrom^-1") + + gauss_plg = cubeviz_helper.plugins["Gaussian Smooth"]._obj + gauss_plg.mode_selected = "Spatial" + with pytest.warns(AstropyUserWarning, match="The following attributes were set on the data"): + _ = gauss_plg.smooth() + + # Need this to make it available for photometry data drop-down. + cubeviz_helper.app.add_data_to_viewer("uncert-viewer", "test[FLUX] spatial-smooth stddev-1.0") + + aper = RectanglePixelRegion(center=PixCoord(x=1, y=2), width=3, height=5) + cubeviz_helper.load_regions(aper) + + plg = cubeviz_helper.plugins["Aperture Photometry"]._obj + plg.dataset_selected = "test[FLUX] spatial-smooth stddev-1.0" + plg.aperture_selected = "Subset 1" + plg.vue_do_aper_phot() + row = cubeviz_helper.get_aperture_photometry_results()[0] + + # Basically, we should recover the input rectangle here. + assert_allclose(row["xcenter"], 1 * u.pix) + assert_allclose(row["ycenter"], 2 * u.pix) + sky = row["sky_center"] + assert_allclose(sky.ra.deg, 205.43985906934287) + assert_allclose(sky.dec.deg, 27.003490103642033) + assert_allclose(row["sum"], 48.54973 * flux_unit) # 3 (w) x 5 (h) x <5 (v) + assert_allclose(row["sum_aper_area"], 15 * (u.pix * u.pix)) # 3 (w) x 5 (h) + assert_allclose(row["mean"], 3.236648941040039 * flux_unit) + assert_quantity_allclose(row["slice_wave"], 4.894499866699333 * u.um) + + +def test_cubeviz_aperphot_cube_orig_flux_mjysr(cubeviz_helper, spectrum1d_cube_custom_fluxunit): + cube = spectrum1d_cube_custom_fluxunit(fluxunit=u.MJy / u.sr) + cubeviz_helper.load_data(cube, data_label="test") + + aper = RectanglePixelRegion(center=PixCoord(x=1, y=3), width=1, height=1) + bg = RectanglePixelRegion(center=PixCoord(x=0, y=2), width=1, height=1) + cubeviz_helper.load_regions([aper, bg]) + + plg = cubeviz_helper.plugins["Aperture Photometry"]._obj + plg.dataset_selected = "test[FLUX]" + plg.aperture_selected = "Subset 1" + plg.background_selected = "Subset 2" + + # Make sure per steradian is handled properly. + assert_allclose(plg.pixel_area, 0.01) + assert_allclose(plg.flux_scaling, 0.003631) + + plg.vue_do_aper_phot() + row = cubeviz_helper.get_aperture_photometry_results()[0] + + # Basically, we should recover the input rectangle here, minus background. + assert_allclose(row["xcenter"], 1 * u.pix) + assert_allclose(row["ycenter"], 3 * u.pix) + assert_allclose(row["sum"], 1.1752215e-12 * u.MJy) # (15 - 10) MJy/sr x 2.3504431e-13 sr + assert_allclose(row["sum_aper_area"], 1 * (u.pix * u.pix)) + assert_allclose(row["pixarea_tot"], 2.350443053909789e-13 * u.sr) + assert_allclose(row["aperture_sum_mag"], 23.72476627732448 * u.mag) + assert_allclose(row["mean"], 5 * (u.MJy / u.sr)) + assert_quantity_allclose(row["slice_wave"], 0.46236 * u.um) diff --git a/jdaviz/configs/cubeviz/plugins/tests/test_parsers.py b/jdaviz/configs/cubeviz/plugins/tests/test_parsers.py index 12a10adebe..8a42cc422a 100644 --- a/jdaviz/configs/cubeviz/plugins/tests/test_parsers.py +++ b/jdaviz/configs/cubeviz/plugins/tests/test_parsers.py @@ -53,7 +53,7 @@ def test_fits_image_hdu_with_microns(image_cube_hdu_obj_microns, cubeviz_helper) flux_unit_str = "erg / (Angstrom cm2 s)" else: flux_unit_str = "erg / (Angstrom s cm2)" - assert label_mouseover.as_text() == (f'Pixel x=00.0 y=00.0 Value +1.00000e+00 1e-17 {flux_unit_str}', # noqa + assert label_mouseover.as_text() == (f'Pixel x=00.0 y=00.0 Value +5.00000e+00 1e-17 {flux_unit_str}', # noqa 'World 13h41m45.5759s +27d00m12.3044s (ICRS)', '205.4398995981 27.0034178810 (deg)') # noqa unc_viewer = cubeviz_helper.app.get_viewer('uncert-viewer') diff --git a/jdaviz/configs/default/plugins/collapse/collapse.py b/jdaviz/configs/default/plugins/collapse/collapse.py index b3cad6909d..671051a51b 100644 --- a/jdaviz/configs/default/plugins/collapse/collapse.py +++ b/jdaviz/configs/default/plugins/collapse/collapse.py @@ -4,7 +4,6 @@ from astropy.nddata import CCDData from glue.core import Data -from specutils import Spectrum1D from specutils.manipulation import spectral_slab from traitlets import Bool, List, Unicode, observe @@ -101,9 +100,20 @@ def collapse(self, add_data=True): """ # Collapsing over the spectral axis. Cut out the desired spectral # region. Defaults to the entire spectrum. - cube = self.dataset.get_object(cls=Spectrum1D, statistic=None) + cube = self.dataset.selected_obj spec_min, spec_max = self.spectral_subset.selected_min_max(cube) + # Extract 2D WCS from input cube. + data = self.dataset.selected_dc_item + # Similar to coords_info logic. + if '_orig_spec' in getattr(data, 'meta', {}): + w = data.meta['_orig_spec'].wcs + else: + w = data.coords + data_wcs = getattr(w, 'celestial', None) + if data_wcs: + data_wcs = data_wcs.swapaxes(0, 1) # We also transpose WCS to match. + with warnings.catch_warnings(): warnings.filterwarnings('ignore', message='No observer defined on WCS') spec = spectral_slab(cube, spec_min, spec_max) @@ -111,15 +121,15 @@ def collapse(self, add_data=True): collapsed_spec = spec.collapse(self.function_selected.lower(), axis=-1).T # Quantity # stuff for exporting to file - self.collapsed_spec = collapsed_spec + self.collapsed_spec = CCDData(collapsed_spec, wcs=data_wcs) self.collapsed_spec_available = True fname_label = self.dataset_selected.replace("[", "_").replace("]", "") self.filename = f"collapsed_{self.function_selected.lower()}_{fname_label}.fits" if add_data: - data = Data() + data = Data(coords=data_wcs) data['flux'] = collapsed_spec.value - data.get_component('flux').units = str(collapsed_spec.unit) + data.get_component('flux').units = collapsed_spec.unit.to_string() self.add_results.add_results_from_plugin(data) @@ -170,7 +180,7 @@ def _save_collapsed_spec_to_fits(self, overwrite=False, *args): return filename = str(filename) - CCDData(self.collapsed_spec).write(filename) + self.collapsed_spec.write(filename) # Let the user know where we saved the file. self.hub.broadcast(SnackbarMessage( diff --git a/jdaviz/configs/default/plugins/collapse/tests/test_collapse.py b/jdaviz/configs/default/plugins/collapse/tests/test_collapse.py index 021cc7c8b1..a6405d507f 100644 --- a/jdaviz/configs/default/plugins/collapse/tests/test_collapse.py +++ b/jdaviz/configs/default/plugins/collapse/tests/test_collapse.py @@ -1,6 +1,5 @@ import numpy as np import pytest -import os from astropy.nddata import CCDData from astropy import units as u from specutils import Spectrum1D @@ -40,7 +39,7 @@ def test_linking_after_collapse(cubeviz_helper, spectral_cube_wcs): assert dc.external_links[1].cids2[0] == dc[-1].pixel_component_ids[0] -def test_save_collapsed_to_fits(cubeviz_helper, spectral_cube_wcs, tmpdir): +def test_save_collapsed_to_fits(cubeviz_helper, spectral_cube_wcs, tmp_path): cubeviz_helper.load_data(Spectrum1D(flux=np.ones((3, 4, 5)) * u.nJy, wcs=spectral_cube_wcs)) @@ -58,15 +57,15 @@ def test_save_collapsed_to_fits(cubeviz_helper, spectral_cube_wcs, tmpdir): # check that default filename is correct, then change path fname = 'collapsed_sum_Unknown spectrum object_FLUX.fits' assert collapse_plugin._obj.filename == fname - collapse_plugin._obj.filename = os.path.join(tmpdir, fname) + collapse_plugin._obj.filename = str(tmp_path / fname) # save output file with default name, make sure it exists collapse_plugin._obj.vue_save_as_fits() - assert os.path.isfile(os.path.join(tmpdir, fname)) + assert (tmp_path / fname).is_file() # read file back in, make sure it matches - dat = CCDData.read(os.path.join(tmpdir, fname)) - assert np.all(dat.data == collapse_plugin._obj.collapsed_spec.value) + dat = CCDData.read(collapse_plugin._obj.filename) + np.testing.assert_array_equal(dat.data, collapse_plugin._obj.collapsed_spec.data) assert dat.unit == collapse_plugin._obj.collapsed_spec.unit # make sure correct error message is raised when export_enabled is False diff --git a/jdaviz/configs/default/plugins/gaussian_smooth/gaussian_smooth.py b/jdaviz/configs/default/plugins/gaussian_smooth/gaussian_smooth.py index c1ac6a6c45..e8cd99e3ae 100644 --- a/jdaviz/configs/default/plugins/gaussian_smooth/gaussian_smooth.py +++ b/jdaviz/configs/default/plugins/gaussian_smooth/gaussian_smooth.py @@ -214,12 +214,7 @@ def spatial_smooth(self): cube : `~specutils.Spectrum1D` The smoothed cube """ - # Get information from the flux component - attribute = self.dataset.selected_dc_item.main_components[0] - - cube = self.dataset.get_object(cls=Spectrum1D, - attribute=attribute, - statistic=None) + cube = self.dataset.selected_obj flux_unit = cube.flux.unit # Extend the 2D kernel to have a length 1 spectral dimension, so that @@ -228,8 +223,16 @@ def spatial_smooth(self): convolved_data = convolve(cube, kernel) + # Copy 3D WCS from input cube. + data = self.dataset.selected_dc_item + # Similar to coords_info logic. + if '_orig_spec' in getattr(data, 'meta', {}): + w = data.meta['_orig_spec'].wcs + else: + w = data.coords + # Create a new cube with the old metadata. Note that astropy # convolution generates values for masked (NaN) data. - newcube = Spectrum1D(flux=convolved_data * flux_unit, wcs=cube.wcs) + newcube = Spectrum1D(flux=convolved_data * flux_unit, wcs=w) return newcube diff --git a/jdaviz/configs/imviz/plugins/aper_phot_simple/aper_phot_simple.py b/jdaviz/configs/imviz/plugins/aper_phot_simple/aper_phot_simple.py index 19fa688357..c21667cd05 100644 --- a/jdaviz/configs/imviz/plugins/aper_phot_simple/aper_phot_simple.py +++ b/jdaviz/configs/imviz/plugins/aper_phot_simple/aper_phot_simple.py @@ -18,7 +18,7 @@ from traitlets import Any, Bool, Integer, List, Unicode, observe from jdaviz.core.custom_traitlets import FloatHandleEmpty -from jdaviz.core.events import SnackbarMessage, LinkUpdatedMessage +from jdaviz.core.events import SnackbarMessage, LinkUpdatedMessage, SliceWavelengthUpdatedMessage from jdaviz.core.region_translators import regions2aperture, _get_region_from_spatial_subset from jdaviz.core.registries import tray_registry from jdaviz.core.template_mixin import (PluginTemplateMixin, DatasetMultiSelectMixin, @@ -43,7 +43,7 @@ class SimpleAperturePhotometry(PluginTemplateMixin, ApertureSubsetSelectMixin, :ref:`public plugin API `: * :meth:`~jdaviz.core.template_mixin.PluginTemplateMixin.show` - * :meth:`~jdaviz.core.template_mixin.PluginTemplateMixin.open_in_tray` + * :meth:`~jdaviz.core.template_mixin.PluginTemplateMixin.open_in_tray` * :meth:`~jdaviz.core.template_mixin.PluginTemplateMixin.close_in_tray` """ template_file = __file__, "aper_phot_simple.vue" @@ -70,6 +70,10 @@ class SimpleAperturePhotometry(PluginTemplateMixin, ApertureSubsetSelectMixin, fit_radial_profile = Bool(False).tag(sync=True) fit_results = List().tag(sync=True) + # Cubeviz only + cube_slice = Unicode("").tag(sync=True) + is_cube = Bool(False).tag(sync=True) + icon_radialtocheck = Unicode(read_icon(os.path.join(ICON_DIR, 'radialtocheck.svg'), 'svg+xml')).tag(sync=True) # noqa icon_checktoradial = Unicode(read_icon(os.path.join(ICON_DIR, 'checktoradial.svg'), 'svg+xml')).tag(sync=True) # noqa @@ -107,6 +111,23 @@ def __init__(self, *args, **kwargs): self.session.hub.subscribe(self, SubsetUpdateMessage, handler=self._on_subset_update) self.session.hub.subscribe(self, LinkUpdatedMessage, handler=self._on_link_update) + # Custom dataset filters for Cubeviz + if self.config == "cubeviz": + def valid_cubeviz_datasets(data): + comp = data.get_component(data.main_components[0]) + img_unit = u.Unit(comp.units) if comp.units else u.dimensionless_unscaled + acceptable_types = ['spectral flux density wav', + 'photon flux density wav', + 'spectral flux density', + 'photon flux density'] + return ((data.ndim in (2, 3)) and + ((img_unit == (u.MJy / u.sr)) or + (img_unit.physical_type in acceptable_types))) + + self.dataset.add_filter(valid_cubeviz_datasets) + self.session.hub.subscribe(self, SliceWavelengthUpdatedMessage, + handler=self._on_slice_changed) + # TODO: expose public API once finalized # @property # def user_api(self): @@ -116,6 +137,23 @@ def __init__(self, *args, **kwargs): # 'calculate_photometry', # 'unpack_batch_options', 'calculate_batch_photometry')) + def _on_slice_changed(self, msg): + if self.config != "cubeviz": + return + self.cube_slice = f"{msg.wavelength:.3e} {msg.wavelength_unit}" + self._cube_wave = u.Quantity(msg.wavelength, msg.wavelength_unit) + self._cube_idx = int(msg.slice) + + @observe("dataset_selected") + def _on_dataset_selected_changed(self, event={}): + if self.config != "cubeviz": + return + # self.dataset might not exist when app is setting itself up. + if hasattr(self, "dataset") and self.dataset.selected_dc_item.ndim > 2: + self.is_cube = True + else: + self.is_cube = False + def _get_defaults_from_metadata(self, dataset=None): defaults = {} if dataset is None: @@ -132,11 +170,17 @@ def _get_defaults_from_metadata(self, dataset=None): else: telescope = meta.get('TELESCOP', '') if telescope == 'JWST': + # Hardcode the flux conversion factor from MJy to ABmag + mjy2abmag = 0.003631 if 'photometry' in meta and 'pixelarea_arcsecsq' in meta['photometry']: defaults['pixel_area'] = meta['photometry']['pixelarea_arcsecsq'] if 'bunit_data' in meta and meta['bunit_data'] == u.Unit("MJy/sr"): - # Hardcode the flux conversion factor from MJy to ABmag - defaults['flux_scaling'] = 0.003631 + defaults['flux_scaling'] = mjy2abmag + elif 'PIXAR_A2' in meta: + defaults['pixel_area'] = meta['PIXAR_A2'] + if 'BUNIT' in meta and meta['BUNIT'] == u.Unit("MJy/sr"): + defaults['flux_scaling'] = mjy2abmag + elif telescope == 'HST': # TODO: Add more HST support, as needed. # HST pixel scales are from instrument handbooks. @@ -269,11 +313,24 @@ def _calc_background_median(self, reg, data=None): raise ValueError("cannot calculate background median in multiselect") else: data = self.dataset.selected_dc_item + comp = data.get_component(data.main_components[0]) + + if self.config == "cubeviz" and data.ndim > 2: + comp_data = comp.data[:, :, self._cube_idx].T # nx, ny --> ny, nx + # Similar to coords_info logic. + if '_orig_spec' in getattr(data, 'meta', {}): + w = data.meta['_orig_spec'].wcs.celestial + else: + w = data.coords.celestial + else: # "imviz" + comp_data = comp.data # ny, nx + w = data.coords + if hasattr(reg, 'to_pixel'): - reg = reg.to_pixel(data.coords) + reg = reg.to_pixel(w) aper_mask_stat = reg.to_mask(mode='center') - img_stat = aper_mask_stat.get_values(comp.data, mask=None) + img_stat = aper_mask_stat.get_values(comp_data, mask=None) # photutils/background/_utils.py --> nanmedian() return np.nanmedian(img_stat) # Naturally in data unit @@ -380,14 +437,31 @@ def calculate_photometry(self, dataset=None, aperture=None, background=None, except ValueError: # Clearer error message raise ValueError('Missing or invalid background value') + if self.config == "cubeviz" and data.ndim > 2: + comp_data = comp.data[:, :, self._cube_idx].T # nx, ny --> ny, nx + # Similar to coords_info logic. + if '_orig_spec' in getattr(data, 'meta', {}): + w = data.meta['_orig_spec'].wcs + else: + w = data.coords + else: # "imviz" + comp_data = comp.data # ny, nx + w = data.coords + if hasattr(reg, 'to_pixel'): sky_center = reg.center - xcenter, ycenter = data.coords.world_to_pixel(sky_center) + if self.config == "cubeviz" and data.ndim > 2: + ycenter, xcenter = w.world_to_pixel(self._cube_wave, sky_center)[1] + else: # "imviz" + xcenter, ycenter = w.world_to_pixel(sky_center) else: xcenter = reg.center.x ycenter = reg.center.y if data.coords is not None: - sky_center = data.coords.pixel_to_world(xcenter, ycenter) + if self.config == "cubeviz" and data.ndim > 2: + sky_center = w.pixel_to_world(self._cube_idx, ycenter, xcenter)[1] + else: # "imviz" + sky_center = w.pixel_to_world(xcenter, ycenter) else: sky_center = None @@ -395,7 +469,6 @@ def calculate_photometry(self, dataset=None, aperture=None, background=None, include_pixarea_fac = False include_counts_fac = False include_flux_scale = False - comp_data = comp.data if comp.units: img_unit = u.Unit(comp.units) bg = bg * img_unit @@ -464,6 +537,13 @@ def calculate_photometry(self, dataset=None, aperture=None, background=None, 'data_label', 'subset_label', 'timestamp'], indexes=[1, 1, 1, 1, 3, 3, 3, 3, 3, 3, 18, 18, 18]) + if self.config == "cubeviz": + if data.ndim > 2: + slice_val = self._cube_wave + else: + slice_val = u.Quantity(np.nan, self._cube_wave.unit) + phot_table.add_column(slice_val, name="slice_wave", index=29) + if add_to_table: try: phot_table['id'][0] = self.table._qtable['id'].max() + 1 @@ -479,7 +559,10 @@ def calculate_photometry(self, dataset=None, aperture=None, background=None, # Plots. if update_plots: if self.current_plot_type == "Curve of Growth": - self.plot.figure.title = 'Curve of growth from aperture center' + if self.config == "cubeviz" and data.ndim > 2: + self.plot.figure.title = f'Curve of growth from aperture center at {slice_val:.4e}' # noqa: E501 + else: + self.plot.figure.title = 'Curve of growth from aperture center' x_arr, sum_arr, x_label, y_label = _curve_of_growth( comp_data, (xcenter, ycenter), aperture, phot_table['sum'][0], wcs=data.coords, background=bg, pixarea_fac=pixarea_fac) @@ -494,7 +577,10 @@ def calculate_photometry(self, dataset=None, aperture=None, background=None, self.plot.figure.axes[1].label = comp.units or 'Value' if self.current_plot_type == "Radial Profile": - self.plot.figure.title = 'Radial profile from aperture center' + if self.config == "cubeviz" and data.ndim > 2: + self.plot.figure.title = f'Radial profile from aperture center at {slice_val:.4e}' # noqa: E501 + else: + self.plot.figure.title = 'Radial profile from aperture center' x_data, y_data = _radial_profile( phot_aperstats.data_cutout, phot_aperstats.bbox, (xcenter, ycenter), raw=False) @@ -502,7 +588,10 @@ def calculate_photometry(self, dataset=None, aperture=None, background=None, self.plot.update_style('profile', line_visible=True, color='gray', size=32) else: # Radial Profile (Raw) - self.plot.figure.title = 'Raw radial profile from aperture center' + if self.config == "cubeviz" and data.ndim > 2: + self.plot.figure.title = f'Raw radial profile from aperture center at {slice_val:.4e}' # noqa: E501 + else: + self.plot.figure.title = 'Raw radial profile from aperture center' x_data, y_data = _radial_profile( phot_aperstats.data_cutout, phot_aperstats.bbox, (xcenter, ycenter), raw=True) @@ -549,7 +638,7 @@ def calculate_photometry(self, dataset=None, aperture=None, background=None, x = phot_table[key][0] if (isinstance(x, (int, float, u.Quantity)) and key not in ('xcenter', 'ycenter', 'sky_center', 'sum_aper_area', - 'aperture_sum_counts', 'aperture_sum_mag')): + 'aperture_sum_counts', 'aperture_sum_mag', 'slice_wave')): tmp.append({'function': key, 'result': f'{x:.4e}'}) elif key == 'sky_center' and x is not None: tmp.append({'function': 'RA center', 'result': f'{x.ra.deg:.6f} deg'}) @@ -561,6 +650,9 @@ def calculate_photometry(self, dataset=None, aperture=None, background=None, f'{x:.4e} ({phot_table["aperture_sum_counts_err"][0]:.4e})'}) elif key == 'aperture_sum_mag' and x is not None: tmp.append({'function': key, 'result': f'{x:.3f}'}) + elif key == 'slice_wave': + if data.ndim > 2: + tmp.append({'function': key, 'result': f'{slice_val:.4e}'}) else: tmp.append({'function': key, 'result': str(x)}) diff --git a/jdaviz/configs/imviz/plugins/aper_phot_simple/aper_phot_simple.vue b/jdaviz/configs/imviz/plugins/aper_phot_simple/aper_phot_simple.vue index 9d4236fe15..7d1185a63f 100644 --- a/jdaviz/configs/imviz/plugins/aper_phot_simple/aper_phot_simple.vue +++ b/jdaviz/configs/imviz/plugins/aper_phot_simple/aper_phot_simple.vue @@ -32,6 +32,20 @@ hint="Select the data for photometry." /> +
+ + + + + +
+