Add cube writer to properly write out cube generated by model fitting

[pllim]

Description

This pull request is to add a custom specutils writer so one can write out spectral cube (as Spectrum1D object) to a FITS file, and then have it read back in as a cube (not a table).

This is necessary until the following issue/PR are resolved upstream:

• Writer needed for cube data astropy/specutils#592
• Fix and extend wcs1d-fits loader for multi-D WCS astropy/specutils#1009

With this patch, you can do cube fitting, put the fitted cube in uncert-viewer, and then also do this:

from specutils import Spectrum1D
sp = cubeviz.app.get_data_from_viewer("uncert-viewer", data_label="cube-fit model", cls=Spectrum1D)
sp.write("mydata.fits", format="jdaviz-cube-writer")

Then you can read it back into Cubeviz (a new instance) like this:

cubeviz2 = Cubeviz()
cubeviz2.load_data("mydata.fits")
cubeviz2.show()

TODO

• Fix test failure.
• Adam G: Is it safe to assume EXTNAME=MASK is a standard for HDU with mask info? Safe enough for Jdaviz.
• Adam G: Is the mask being written has 0=good or 0=bad? We want 0=good.

Change log entry

• Is a change log needed? If yes, is it added to CHANGES.rst? If you want to avoid merge conflicts,
  list the proposed change log here for review and add to CHANGES.rst before merge. If no, maintainer
  should add a no-changelog-entry-needed label.

Checklist for package maintainer(s)

This checklist is meant to remind the package maintainer(s) who will review this pull request of some common things to look for. This list is not exhaustive.

• Are two approvals required? Branch protection rule does not check for the second approval. If a second approval is not necessary, please apply the trivial label.
• Do the proposed changes actually accomplish desired goals? Also manually run the affected example notebooks, if necessary.
• Do the proposed changes follow the STScI Style Guides?
• Are tests added/updated as required? If so, do they follow the STScI Style Guides?
• Are docs added/updated as required? If so, do they follow the STScI Style Guides?
• Did the CI pass? If not, are the failures related?
• Is a milestone set? Set this to bugfix milestone if this is a bug fix and needs to be released ASAP; otherwise, set this to the next major release milestone.
• After merge, any internal documentations need updating (e.g., JIRA, Innerspace)? 🐱

[keflavich]

Coming from the other thread: this lgtm, but I'm not giving a formal review approval since I'm not part of the organization. I just don't see any obvious problems with this approach.

I'll note that this effectively establishes a standard in which the 'MASK' HDU stores the mask, which isn't necessarily widely adopted. I like it and can immediately add it to the spectral-cube FITS reader. Could you add to the docs, though, where the rules for the MASK are? I think it should be listed in the docstring of the writer. I'm pretty sure 0 means 'good', and anything nonzero is some form of 'bad', but past that I think it is in the JWST docs somewhere?

[pllim]

Could you add to the docs, though, where the rules for the MASK are?

Good question! I didn't even consider that the mask might be flipped because it is using the numpy.ma convention. I'll ask the other devs tomorrow. 😱 I think @bmorris3 is the most recent dev to do the mask stuff but he is out currently.

[pllim]

a standard in which the 'MASK' HDU stores the mask

This is also a good point. I picked MASK from this existing code but that happened before I joined the project, so I would also need to consult with other devs.

jdaviz/jdaviz/configs/cubeviz/plugins/parsers.py

Line 19 in 2d767af

mask=['mask', 'dq', 'quality'])

Masking was discussed in APE 11 that never got merged and also that discussion was not really for spectra, though you would think masking would be general enough to also cover it.

In real HST and JWST data, they actually contained bit-plane mask, not just simple 0s and 1s. AFAIK, that is currently not supported by Jdaviz parsers.

Regardless, naively, I would assume 0 means "good" because it is defined as such for the bit-plane stuff.

[dhomeier]

If this is the file structure you need, I don't think wcs1d-fits can do it in its present form – with or without the new PR. You can write to the first extension instead of the primary with write(..., hdu=1, format='wcs1d-fits'), but currently there is no provision to write the mask. It could certainly be added, but last I recall specutils had rather poor support to write to a new HDU in an existing HDUList.
If this loader is working as it should, might be best to use it as is and perhaps then contribute to the default_loaders upstream. Would this make sense to provide as a writer for jwst_s3d_loader, or is that yet another format?

[pllim]

This has nothing to do with JWST, it is a "generic" format for a cube output from model fitting plugin that could or could not have mask populated. So putting this in jwst_s3d_loader would be very misleading.

[dhomeier]

In that case it would probably belong into wcs1d-fits. I have not found much on mask treatment apart from SDSS/MaNGA (which is storing it in hdulist['MASK']), so otherwise don't know what would be a standard for storing it.

[pllim]

I dunno either. If this is not generic enough, personally, I am okay with this writer just living in Jdaviz forever.

[codecov]

Codecov Report

Base: 92.06% // Head: 92.09% // Increases project coverage by +0.02% 🎉

Coverage data is based on head (e05119d) compared to base (65b7779).
Patch coverage: 98.94% of modified lines in pull request are covered.

Additional details and impacted files

@@            Coverage Diff             @@
##             main    #2012      +/-   ##
==========================================
+ Coverage   92.06%   92.09%   +0.02%     
==========================================
  Files         140      140              
  Lines       15305    15373      +68     
==========================================
+ Hits        14091    14158      +67     
- Misses       1214     1215       +1     

Impacted Files	Coverage Δ
...igs/default/plugins/model_fitting/model_fitting.py	83.69% <66.66%> (-0.21%)	⬇️
jdaviz/configs/cubeviz/helper.py	98.82% <100.00%> (+0.41%)	⬆️
jdaviz/configs/cubeviz/plugins/parsers.py	68.60% <100.00%> (+1.32%)	⬆️
...efault/plugins/model_fitting/tests/test_fitting.py	99.35% <100.00%> (+0.19%)	⬆️

Help us with your feedback. Take ten seconds to tell us how you rate us. Have a feature suggestion? Share it here.

☔ View full report at Codecov.
📢 Do you have feedback about the report comment? Let us know in this issue.

[pllim]

The ASDF warnings in dev job are unrelated. ASDF is changing upstream.

[kecnry]

In general this seems like a reasonable approach to me for now. What are other valid format names (from upstream)? Would it make more sense to be just "jdaviz-cube" since that describes the format of the data (rather than "-writer" which to me more describes the method to export the data)?

[kecnry]

should there be an example on how to re-load the resulting fits file back into a Spectrum1D or is that just considered obvious? Does that require providing the format again?

[pllim]

I think you can just use vanilla Spectrum1D.read. Existing Cubeviz parser can read it back.

[kecnry]

Writer and reader cannot share the same identifier.

Is that true? It looks like other cases use the same format string for both (astropy/specutils#1009, for example... although that is implemented differently).

I think you can just use vanilla Spectrum1D.read. Existing Cubeviz parser can read it back.

Should this be covered in the test? The test currently just loads with fits.open.

As a user, I would just expect the same string to be required (or at least to work) when writing and then reading, since the kwargs is format and not writer. So even if it isn't required, does Spectrum1D.read(..., format='jdaviz-cube-writer') raise an error now?

[pllim]

Is that true?

Actually now that I think about it, if this really follows the astropy.io.registry stuff, then the name can indeed be shared between reader and writer. I think the specutils custom writer example confused me. I'll change the name as you suggested. Thanks!

[pllim]

Should this be covered in the test?

Sure, doesn't hurt.

The reader would raise error because that format isn't register on the reader side. I can also add a test for that.

[dhomeier]

Actually now that I think about it, if this really follows the astropy.io.registry stuff, then the name can indeed be shared between reader and writer. I think the specutils custom writer example confused me. I'll change the name as you suggested.

Yes, the typical pattern in the default loaders is to provide both a @data_loader and @custom_writer for the format name, and have the identifier check for origin ('read' or 'write') as required. Though I think some of those things are indeed handled differently in custom loaders.

[bmorris3]

Looks good overall!

I agree that we should be careful to write down somewhere that we're sticking with these mask conventions, as mentioned in #2012 (comment). I agree that zero/False should correspond to "unmasked," and non-zero/True to "masked."

[pllim]

to write down somewhere

Is the additional text linking to specutils doc in the diff here not enough, @bmorris3 ?

[bmorris3]

Is the additional text linking to specutils doc in the diff here not enough, @bmorris3 ?

I expect users will tell us if this is confusing. I doubt most of our users are even aware that glue masks are the inverse of this convention, so it's probably not a top priority to inject that confusion into our docs.

In summary – I don't think we should worry about it. 💯

[pllim]

I think I addressed all the comments? Can I merge when CI pass since there are 2 approvals or is there something else? Thanks!

[pllim]

Thanks, all!