Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

fix: SDSS-V SpectrumList format ambiguity, mwmVisit BOSS load fail #1185

Merged
merged 15 commits into from
Nov 1, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
8 changes: 8 additions & 0 deletions CHANGES.rst
Original file line number Diff line number Diff line change
Expand Up @@ -25,15 +25,23 @@ Bug Fixes
- Fixed ``Spectrum1D.with_flux_unit()`` not converting uncertainty along
with flux unit. [#1181]

- Fixed ``mwmVisit`` SDSS-V ``Spectrum1D`` and ``SpectrumList`` default loader
being unable to load files containing only BOSS instrument spectra. [#1185]

- Fixed automatic format detection for SDSS-V ``SpectrumList`` default loaders. [#1185]

- Fixed extracting a spectral region when one of spectrum/region is in wavelength
and the other is in frequency units. [#1187]


Other Changes and Additions
^^^^^^^^^^^^^^^^^^^^^^^^^^^

- Replaced ``LevMarLSQFitter`` with ``TRFLSQFitter`` as the former is no longer
recommended by ``astropy``. [#1180]

- "Multi" loaders have been removed from SDSS-V ``SpectrumList`` default loaders. [#1185]

1.17.0 (2024-10-04)
-------------------

Expand Down
71 changes: 40 additions & 31 deletions specutils/io/default_loaders/sdss_v.py
Original file line number Diff line number Diff line change
@@ -1,10 +1,12 @@
"""Register reader functions for various spectral formats."""
import warnings
from typing import Optional

import numpy as np
from astropy.units import Unit, Quantity, Angstrom
from astropy.nddata import StdDevUncertainty, InverseVariance
from astropy.io.fits import HDUList, BinTableHDU, ImageHDU
from astropy.utils.exceptions import AstropyUserWarning

from ...spectra import Spectrum1D, SpectrumList
from ..registers import data_loader
Expand Down Expand Up @@ -183,9 +185,10 @@ def load_sdss_apStar_1D(file_obj, idx: int = 0, **kwargs):


@data_loader(
"SDSS-V apStar multi",
"SDSS-V apStar",
identifier=apStar_identify,
dtype=SpectrumList,
force=True,
priority=10,
extensions=["fits"],
)
Expand Down Expand Up @@ -259,9 +262,10 @@ def load_sdss_apVisit_1D(file_obj, **kwargs):


@data_loader(
"SDSS-V apVisit multi",
"SDSS-V apVisit",
identifier=apVisit_identify,
dtype=SpectrumList,
force=True,
priority=10,
extensions=["fits"],
)
Expand Down Expand Up @@ -311,7 +315,6 @@ def load_sdss_apVisit_list(file_obj, **kwargs):


# BOSS REDUX products (specLite, specFull, custom coadd files, etc)

@data_loader(
"SDSS-V spec",
identifier=spec_sdss5_identify,
Expand All @@ -337,6 +340,8 @@ def load_sdss_spec_1D(file_obj, *args, hdu: Optional[int] = None, **kwargs):
"""
if hdu is None:
# TODO: how should we handle this -- multiple things in file, but the user cannot choose.
warnings.warn('HDU not specified. Loading coadd spectrum (HDU1)',
AstropyUserWarning)
hdu = 1 # defaulting to coadd
# raise ValueError("HDU not specified! Please specify a HDU to load.")
elif hdu in [2, 3, 4]:
Expand All @@ -347,9 +352,10 @@ def load_sdss_spec_1D(file_obj, *args, hdu: Optional[int] = None, **kwargs):


@data_loader(
"SDSS-V spec multi",
"SDSS-V spec",
identifier=spec_sdss5_identify,
dtype=SpectrumList,
force=True,
Comment on lines +355 to +358
Copy link
Contributor

@havok2063 havok2063 Oct 14, 2024

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

With the multi removed, what is the difference between this data loader and the above data loader at line 318? Same with the mwm data loaders?

With this removal, what would be the correct format to specify to load all data extensions? Is it now the default? Jdaviz only supports using format during data load to provide hints to the type of filepath.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

My goal with this was to simplify the interface. The difference previously was that multi loaded every valid HDU extension, whereas the single SpectrumList ones would load just a specified HDU. Since it's a SpectrumList, I felt the default loader for that type should just load every extension.

Across both the Spectrum1D and SpectrumList loaders there should be a single format (SDSS-V [DATATYPE]), which should be detected automagically and no longer needs to be manually specified on SpectrumList.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The issue might be mostly in Specviz. I'd like to retain the workflow of passing Specviz a single file and having it load all spectra from all extensions. Specviz uses the following https://github.com/spacetelescope/jdaviz/blob/main/jdaviz/configs/specviz/plugins/parsers.py#L78-L87 to load a specutils object from a file, which first attempts to load a Spectrum1D with a format option, and on fail, tries to load the SpectrumList with a format option.

With the multi option, I can explicitly tell Specviz to load using the SpectrumList loader, which will break will this change. Specviz only allows me to pass in a format keyword to specutils. Previously using SDSS-V mwm multi on the mwmStar files allowed me to load all spectra, since it triggered line 87. With this, it only loads the first spectrum from the first extension found.

In my mind, the easiest fix would be to reintroduce the multi here, but maybe it's more appropriate to fix Specviz instead. @rosteen since you're a maintainer/dev of both specutils and jdaviz, do you have thoughts on the best approach, and/or where the fix should go?

I am testing with the mwmStar and mwmVisit files for id 61731453, which has spectra in both the BOSS and APOGEE extensions. the mwmStar file has 1 spectrum per extension. The mwmVisit file has 1 spectrum in the BOSS extension, and 3 spectra in the APOGEE extension, loaded as a SpectrumList of 2 Spectrum1D objects, with flux.shapes of ((4648,), (3, 8575)).

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'm planning to catch up on this this afternoon, I might have thoughts then.

priority=5,
extensions=["fits"],
)
Expand Down Expand Up @@ -433,16 +439,21 @@ def _load_BOSS_HDU(hdulist: HDUList, hdu: int, **kwargs):
priority=20,
extensions=["fits"],
)
def load_sdss_mwm_1d(file_obj, hdu: Optional[int] = None, **kwargs):
def load_sdss_mwm_1d(file_obj,
hdu: Optional[int] = None,
visit: Optional[int] = None,
**kwargs):
"""
Load an unspecified spec file as a Spectrum1D.

Parameters
----------
file_obj : str, file-like, or HDUList
FITS file name, file object, or HDUList..
hdu : int
hdu : Optional[int]
Specified HDU to load.
visit : Optional[int]
Specified visit index to load.

Returns
-------
Expand All @@ -459,17 +470,22 @@ def load_sdss_mwm_1d(file_obj, hdu: Optional[int] = None, **kwargs):

# TODO: how should we handle this -- multiple things in file, but the user cannot choose.
if hdu is None:
for i in range(len(hdulist)):
for i in range(1, len(hdulist)):
if hdulist[i].header.get("DATASUM") != "0":
hdu = i
warnings.warn(
'HDU not specified. Loading spectrum at (HDU{})'.
format(i), AstropyUserWarning)
break

return _load_mwmVisit_or_mwmStar_hdu(hdulist, hdu, **kwargs)
# load spectra and return
return _load_mwmVisit_or_mwmStar_hdu(hdulist, hdu)


@data_loader(
"SDSS-V mwm multi",
"SDSS-V mwm",
identifier=mwm_identify,
force=True,
dtype=SpectrumList,
priority=20,
extensions=["fits"],
Expand All @@ -485,12 +501,8 @@ def load_sdss_mwm_list(file_obj, **kwargs):

Returns
-------
SpectrumList
The spectra contained in the file, where:
Spectrum1D
A given spectra of nD flux
None
If there are no spectra for that spectrograph/observatory
SpectrumList[Spectrum1D]
A list spectra from each visit with each instrument at each observatory (mwmVisit), or the coadd from each instrument/observatory (mwmStar).
"""
spectra = SpectrumList()
with read_fileobj_or_hdulist(file_obj, memmap=False, **kwargs) as hdulist:
Expand All @@ -505,10 +517,6 @@ def load_sdss_mwm_list(file_obj, **kwargs):
for hdu in range(1, len(hdulist)):
if hdulist[hdu].header.get("DATASUM") == "0":
# Skip zero data HDU's
# TODO: validate if we want this printed warning or not.
# it might get annoying & fill logs with useless alerts.
print("WARNING: HDU{} ({}) is empty.".format(
hdu, hdulist[hdu].name))
continue
spectra.append(_load_mwmVisit_or_mwmStar_hdu(hdulist, hdu))
return spectra
Expand All @@ -527,8 +535,8 @@ def _load_mwmVisit_or_mwmStar_hdu(hdulist: HDUList, hdu: int, **kwargs):

Returns
-------
Spectrum1D
The spectrum with nD flux contained in the HDU.
list[Spectrum1D]
List of spectrum with 1D flux contained in the HDU.

"""
if hdulist[hdu].header.get("DATASUM") == "0":
Expand Down Expand Up @@ -563,12 +571,11 @@ def _load_mwmVisit_or_mwmStar_hdu(hdulist: HDUList, hdu: int, **kwargs):
mask = mask != 0

# collapse shape if 1D spectra in 2D array
# NOTE: this fixes a jdaviz handling bug for 2D of shape 1,
# it could be that it's expected to be parsed this way.
if flux.shape[0] == 1:
flux = flux[0]
e_flux = e_flux[0]
mask = mask[0]
flux = np.ravel(flux)
e_flux = e_flux[0] # different class
mask = np.ravel(mask)

# Create metadata
meta = dict()
Expand All @@ -578,24 +585,26 @@ def _load_mwmVisit_or_mwmStar_hdu(hdulist: HDUList, hdu: int, **kwargs):
meta["snr"] = np.array(hdulist[hdu].data["snr"])

# Add identifiers (obj, telescope, mjd, datatype)
# TODO: need to see what metadata we're interested in for the MWM files.
meta["telescope"] = hdulist[hdu].data["telescope"]
meta["instrument"] = hdulist[hdu].header.get("INSTRMNT")
try:
try: # get obj if exists
meta["obj"] = hdulist[hdu].data["obj"]
except KeyError:
pass

# choose between mwmVisit/Star via KeyError except
try:
meta["date"] = hdulist[hdu].data["date_obs"]
meta["mjd"] = hdulist[hdu].data["mjd"]
meta['mjd'] = hdulist[hdu].data['mjd']
meta["datatype"] = "mwmVisit"
except KeyError:
meta["mjd"] = (str(hdulist[hdu].data["min_mjd"][0]) + "->" +
str(hdulist[hdu].data["max_mjd"][0]))
meta["min_mjd"] = str(hdulist[hdu].data["min_mjd"][0])
meta["max_mjd"] = str(hdulist[hdu].data["max_mjd"][0])
meta["datatype"] = "mwmStar"
finally:
meta["name"] = hdulist[hdu].name
meta["sdss_id"] = hdulist[hdu].data['sdss_id']

# drop back a list of Spectrum1Ds to unpack
return Spectrum1D(
spectral_axis=spectral_axis,
flux=flux,
Expand Down
Loading