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

Add ability to load 2D flux arrays in Specviz #3229

Merged
merged 16 commits into from
Nov 1, 2024
Merged
1 change: 1 addition & 0 deletions CHANGES.rst
Original file line number Diff line number Diff line change
Expand Up @@ -15,6 +15,7 @@ Mosviz

Specviz
^^^^^^^
- Specviz parser will now split a spectrum with a 2D flux array into multiple spectra on load. [#3229]
Copy link
Contributor

Choose a reason for hiding this comment

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

If this is added to support a specific telescope/instrument, please state also.


Specviz2d
^^^^^^^^^
Expand Down
6 changes: 4 additions & 2 deletions jdaviz/configs/specviz/helper.py
Original file line number Diff line number Diff line change
Expand Up @@ -41,7 +41,8 @@ def __init__(self, *args, **kwargs):
handler=self._redshift_listener)

def load_data(self, data, data_label=None, format=None, show_in_viewer=True,
concat_by_file=False, cache=None, local_path=None, timeout=None):
concat_by_file=False, cache=None, local_path=None, timeout=None,
load_as_list=False):
"""
Load data into Specviz.
Expand Down Expand Up @@ -80,7 +81,8 @@ def load_data(self, data, data_label=None, format=None, show_in_viewer=True,
concat_by_file=concat_by_file,
cache=cache,
local_path=local_path,
timeout=timeout)
timeout=timeout,
load_as_list=load_as_list)

def get_spectra(self, data_label=None, spectral_subset=None, apply_slider_redshift="Warn"):
"""Returns the current data loaded into the main viewer
Expand Down
90 changes: 70 additions & 20 deletions jdaviz/configs/specviz/plugins/parsers.py
Original file line number Diff line number Diff line change
Expand Up @@ -18,7 +18,8 @@

@data_parser_registry("specviz-spectrum1d-parser")
def specviz_spectrum1d_parser(app, data, data_label=None, format=None, show_in_viewer=True,
concat_by_file=False, cache=None, local_path=os.curdir, timeout=None):
concat_by_file=False, cache=None, local_path=os.curdir, timeout=None,
load_as_list=False):
"""
Loads a data file or `~specutils.Spectrum1D` object into Specviz.

Expand Down Expand Up @@ -46,59 +47,81 @@
remote requests in seconds (passed to
`~astropy.utils.data.download_file` or
`~astroquery.mast.Conf.timeout`).
load_as_list : bool, optional
Force the parser to load the input file with the `~specutils.SpectrumList` read function
instead of `~specutils.Spectrum1D`.
"""

spectrum_viewer_reference_name = app._jdaviz_helper._default_spectrum_viewer_reference_name
# If no data label is assigned, give it a unique name
if not data_label:
data_label = app.return_data_label(data, alt_name="specviz_data")
# Still need to standardize the label
elif isinstance(data_label, (list, tuple)):
data_label = [app.return_data_label(label, alt_name="specviz_data") for label in data_label] # noqa
else:
data_label = app.return_data_label(data_label, alt_name="specviz_data")

if isinstance(data, SpectrumCollection):
raise TypeError("SpectrumCollection detected."
" Please provide a Spectrum1D or SpectrumList")
elif isinstance(data, Spectrum1D):
data_label = [app.return_data_label(data_label, alt_name="specviz_data")]
data = [data]
# Handle the possibility of 2D spectra by splitting into separate spectra
if data.flux.ndim == 1:
data_label = [data_label]
data = [data]
elif data.flux.ndim == 2:
data = split_spectrum_with_2D_flux_array(data)
data_label = [f"{data_label} [{i}]" for i in range(len(data))]
# No special processing is needed in this case, but we include it for completeness
elif isinstance(data, SpectrumList):
pass
elif isinstance(data, list):
# special processing for HDUList
if isinstance(data, fits.HDUList):
data = [Spectrum1D.read(data)]
data_label = [app.return_data_label(data_label, alt_name="specviz_data")]
data_label = [data_label]
else:
# list treated as SpectrumList if not an HDUList
data = SpectrumList.read(data, format=format)
else:
# try parsing file_obj as a URI/URL:
data = download_uri_to_path(
data, cache=cache, local_path=local_path, timeout=timeout
)
data = download_uri_to_path(data, cache=cache, local_path=local_path, timeout=timeout)

path = pathlib.Path(data)

if path.is_file():
if path.is_dir() or load_as_list:
data = SpectrumList.read(str(path), format=format)
if data == []:
raise ValueError(f"`specutils.SpectrumList.read('{str(path)}')` "

Check warning on line 96 in jdaviz/configs/specviz/plugins/parsers.py

View check run for this annotation

Codecov / codecov/patch

jdaviz/configs/specviz/plugins/parsers.py#L96

Added line #L96 was not covered by tests
"returned an empty list")
elif path.is_file():
try:
data = [Spectrum1D.read(str(path), format=format)]
data_label = [app.return_data_label(data_label, alt_name="specviz_data")]
data = Spectrum1D.read(str(path), format=format)
if data.flux.ndim == 2:
data = split_spectrum_with_2D_flux_array(data)

Check warning on line 102 in jdaviz/configs/specviz/plugins/parsers.py

View check run for this annotation

Codecov / codecov/patch

jdaviz/configs/specviz/plugins/parsers.py#L102

Added line #L102 was not covered by tests
else:
data = [data]
data_label = [app.return_data_label(data_label, alt_name="specviz_data")]

except IORegistryError:
# Multi-extension files may throw a registry error
data = SpectrumList.read(str(path), format=format)
elif path.is_dir():
data = SpectrumList.read(str(path), format=format)
if data == []:
raise ValueError(f"`specutils.SpectrumList.read('{str(path)}')` "
"returned an empty list")
else:
raise FileNotFoundError("No such file: " + str(path))

# step through SpectrumList and convert any 2D spectra to 1D spectra
if isinstance(data, SpectrumList):
new_data = SpectrumList()
for spec in data:
if spec.flux.ndim == 2:
new_data.extend(split_spectrum_with_2D_flux_array(spec))
else:
new_data.append(spec)
data = SpectrumList(new_data)
rosteen marked this conversation as resolved.
Show resolved Hide resolved

if not isinstance(data_label, (list, tuple)):
temp_labels = []
for i in range(len(data)):
temp_labels.append(f"{data_label} {i}")
data_label = temp_labels
data_label = [f"{data_label} [{i}]" for i in range(len(data))]
elif len(data_label) != len(data):
raise ValueError(f"Length of data labels list ({len(data_label)}) is different"
f" than length of list of data ({len(data)})")
Expand All @@ -121,7 +144,7 @@
# if the concatenated list of uncertanties is all nan (meaning
# they were all nan to begin with, or all None), it will be set
# to None on the final Spectrum1D
if spec.uncertainty[wlind] is not None:
if spec.uncertainty is not None and spec.uncertainty[wlind] is not None:
dfnuallorig.append(spec.uncertainty[wlind].array)
else:
dfnuallorig.append(np.nan)
Expand Down Expand Up @@ -246,3 +269,30 @@
spec = Spectrum1D(flux=fnuall * flux_units, spectral_axis=wlall * wave_units,
uncertainty=unc)
return spec


def split_spectrum_with_2D_flux_array(data):
"""
Helper function to split Spectrum1D of 2D flux to a SpectrumList of nD objects.
rosteen marked this conversation as resolved.
Show resolved Hide resolved
Parameters
----------
data : `~specutils.Spectrum1D`
Spectrum with 2D flux array

Returns
-------
new_data : `~specutils.SpectrumList`
list of unpacked spectra
rosteen marked this conversation as resolved.
Show resolved Hide resolved
"""
new_data = []
for i in range(data.flux.shape[0]):
unc = None
mask = None
if data.uncertainty is not None:
unc = data.uncertainty[i, :]

Check warning on line 292 in jdaviz/configs/specviz/plugins/parsers.py

View check run for this annotation

Codecov / codecov/patch

jdaviz/configs/specviz/plugins/parsers.py#L292

Added line #L292 was not covered by tests
if data.mask is not None:
mask = data.mask[i, :]

Check warning on line 294 in jdaviz/configs/specviz/plugins/parsers.py

View check run for this annotation

Codecov / codecov/patch

jdaviz/configs/specviz/plugins/parsers.py#L294

Added line #L294 was not covered by tests
new_data.append(Spectrum1D(flux=data.flux[i, :], spectral_axis=data.spectral_axis,
Copy link
Contributor

Choose a reason for hiding this comment

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

What about the other properties like redshift, metadata, etc? Do they not have to be propagated?

Copy link
Contributor

Choose a reason for hiding this comment

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

Also... is an utility function like this better suited for specutils upstream?

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

I just added the meta, I had discussed with Brian and we agreed that it's ok to attach the meta from the original to each split spectrum at least for now. Redshift comes along with the spectral_axis.

I think I'd prefer to keep this here for now, but we could potentially upstream it at some point. I don't know exactly where I'd want to put it in specutils.

Copy link
Contributor

Choose a reason for hiding this comment

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

OK, I opened astropy/specutils#1188

uncertainty=unc, mask=mask))

return new_data
18 changes: 18 additions & 0 deletions jdaviz/configs/specviz/tests/test_helper.py
Original file line number Diff line number Diff line change
Expand Up @@ -405,6 +405,24 @@ def test_load_spectrum_list_directory_concat(tmpdir, specviz_helper):
assert len(specviz_helper.app.data_collection) == 41


def test_load_2d_flux(specviz_helper):
spec = Spectrum1D(spectral_axis=np.linspace(4000, 6000, 10)*u.Angstrom,
Copy link
Contributor

Choose a reason for hiding this comment

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

Maybe add comments to explain what is going on here for future devs. Hard to tell from reading the code without more context.

flux=np.ones((4, 10))*u.Unit("1e-17 erg / (Angstrom cm2 s)"))
specviz_helper.load_data(spec, data_label="test")

assert len(specviz_helper.app.data_collection) == 4
assert specviz_helper.app.data_collection[0].label == "test [0]"

spec2 = Spectrum1D(spectral_axis=np.linspace(4000, 6000, 10)*u.Angstrom,
flux=np.ones((2, 10))*u.Unit("1e-17 erg / (Angstrom cm2 s)"))

spec_list = SpectrumList([spec, spec2])
specviz_helper.load_data(spec_list, data_label="second test")

assert len(specviz_helper.app.data_collection) == 10
assert specviz_helper.app.data_collection[-1].label == "second test [5]"


def test_plot_uncertainties(specviz_helper, spectrum1d):
specviz_helper.load_data(spectrum1d)

Expand Down