SpectralCube to specutils transition in cubeviz

[javerbukh]

Resolves #544

Requires Specutils 811

Also requires a fix to glue-astronomy 28

Edit: Requires glue-astronomy 29

2nd Edit: Requires Specutils 822 and glue-astronomy 36

• Also fix [FEAT] Cubeviz: Update cubeviz plugins to use Spectrum1D instead of SpectralCube. #587
• Drop spectral-cube as a dependency
• Check if [BUG] Cubeviz: Cannot Open Cubes Produced in JWST ETC #935 is resolved by this

[javerbukh]

Latest changes can be tested by running the CubevizExample notebook to the 4th cell. Then, go to Viewer and change the x_axis to Declination and y_axis to Right Ascension. This change needs to be done because the latest specutils updates have a different shape from SpectralCube. I'm not sure yet how to do this programmatically and may need some assistance from @rosteen and @ibusko on how to implement this. Also, the slider seems to snap back to 0 after dragging, not sure why this is happening either.

[ibusko]

I am seeing this output, but no data gets loaded:

FLUX
Wrong number of dimensions in input array.  Expected 2.
IVAR
Wrong number of dimensions in input array.  Expected 2.
MASK
Wrong number of dimensions in input array.  Expected 2.
WAVE
SPECRES
SPECRESD

But isn't that expected, since _new_parse_hdu is just a no-op pass?

[javerbukh]

That is correct about _new_parse_hdu but the 4th cell should manually load data into data_collection and it should be available in the data drop down for a viewer.

[javerbukh]

Latest changes have a manga cube load into cubeviz using Spectrum1D. Some hacky things had to be done for this to work, like using np.moveaxis() on some axes in hdu.data and using Spectrum1D to load individual hdu's rather than the entire data cube at once. Also, if you check the viewer tab, you can see that the x axis is WAVE and the y axis is Declination. Not sure why this is happening but using swapaxes on the wcs does not change this either. I am looking for some feedback and advice on the previously mentioned issues.

[duytnguyendtn]

Here's what I get when I change the axes:

which is probably not what it's supposed to look like. I spent a few hours trying different things to no avail. My only possible theory is whether us moving the axes around (to accommodate for the spectral axis) is either confusing us, or the glue backend) and causing the axes to be mislabeled. I was going to try and see if I could plot it in matplotlib to see if I understood the axes correctly, but promptly got confused...

[javerbukh]

TL;DR: App-wide conversion to Spectrum1D will take longer if wavelength axis is forced to be last, and the time it takes for data to load is significantly longer. Do we still want to continue development in this direction with that in mind? Is there any chance we can change the decision on the wavelength axis being last in Spectrum1D? @eteq

The latest commit has working axes (check _fix_axes()) and uses u.count if no flux units are found, so that solves two of the concerns with this PR. The next concern would be getting the cubeviz plugins to work (i.e. use Spectrum1D instead of SpectralCube) with these changes to the parser, which is tracked here. Since Spectrum1D forces a different order of axes from SpectralCube, we need to make sure that every place that uses a particular axis for wavelength, RA, or Dec now instead either figures out programmatically which is correct or now uses the Spectrum1D location for axes, which make take some time.

Another thought is that there may be work in glue-astronomy that can make this conversion simpler, but that development work is made more complicated by the "wavelength as last axis" decision, in my opinion. That said, I will proceed in this direction for now.

[pllim]

I can't keep track anymore. Can we completely remove spectral-cube yet, or do we have to wait for "Spiderman: Far from Home"?

Suggested change

spectral-cube>=0.5

[javerbukh]

Waiting for Far from Home (Spectrum1D collapse, for posterity) unfortunately!

[pllim]

p.s. I am basically a co-author by now, so I'll abstain from approving.

[PatrickOgle]

There appears to be a problem with spatial Gaussian smoothing, which only smooths in the y direction. Spectral smoothing works just fine.

[javerbukh]

@PatrickOgle I'm not sure I understand, does smoothing need to happen in multiple axes, or in a different axis? This is what I see after smoothing

:

[rosteen]

The dev test failure looks unrelated, so I'm approving at this point. Nice job! I would still like to see @PatrickOgle comment that the moment maps output looks reasonable before merging (especially since he caught a problem with the smoothing!).

[PatrickOgle]

The Gaussian spatial smoothing now appears to be working properly.

I have verified the following plugin functionality:

1. Gaussian smooth spectral
2. Gaussian smooth spatial
3. Collapse (Sum) over spectral region
4. Moment 0 Map
5. Moment 1 Map

However, Moment 2 Crashes, even though it works fine with jdaviz2.0.0...

[PatrickOgle]

For this PR, Moment 2 map on either the whole spectrum or a spectral region crashes on NIRSpec IFU spectral cube processed with JWST pipeline. Moment 0 and Moment 1 work properly. Same error does not occur with jdaviz2.0.0.
Here is the output:

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
~/miniconda3/envs/jdaviz2.0.0/lib/python3.8/site-packages/ipyvue/VueTemplateWidget.py in _handle_event(self, _, content, buffers)
     58                 getattr(self, "vue_" + event)(data, buffers)
     59             else:
---> 60                 getattr(self, "vue_" + event)(data)
     61 
     62 

~/Desktop/jdat/jdaviz/jdaviz/configs/cubeviz/plugins/moment_maps/moment_maps.py in vue_calculate_moment(self, *args)
    134         except ValueError:
    135             raise ValueError("Moment must be a positive integer")
--> 136         self.moment = analysis.moment(slab, order=n_moment)
    137 
    138         moment_ccd = CCDData(self.moment, unit=self.moment.unit)

~/miniconda3/envs/jdaviz2.0.0/lib/python3.8/site-packages/specutils/analysis/moment.py in moment(spectrum, regions, order, axis)
     40 
     41     """
---> 42     return computation_wrapper(_compute_moment, spectrum, regions,
     43                                order=order, axis=axis)
     44 

~/miniconda3/envs/jdaviz2.0.0/lib/python3.8/site-packages/specutils/analysis/utils.py in computation_wrapper(func, spectrum, region, **kwargs)
     17     # No region, therefore whole spectrum.
     18     if region is None:
---> 19         return func(spectrum, **kwargs)
     20 
     21     # Single region

~/miniconda3/envs/jdaviz2.0.0/lib/python3.8/site-packages/specutils/analysis/moment.py in _compute_moment(spectrum, regions, order, axis)
     80             m1 = np.tile(m1, _shape).T
     81 
---> 82         return np.sum(flux * (dispersion - m1) ** order, axis=axis) / m0

~/miniconda3/envs/jdaviz2.0.0/lib/python3.8/site-packages/astropy/coordinates/spectral_quantity.py in __array_ufunc__(self, function, method, *inputs, **kwargs)
     86     def __array_ufunc__(self, function, method, *inputs, **kwargs):
     87         # We always return Quantity except in a few specific cases
---> 88         result = super().__array_ufunc__(function, method, *inputs, **kwargs)
     89         if ((function is np.multiply
     90             or function is np.true_divide and inputs[0] is self)

~/miniconda3/envs/jdaviz2.0.0/lib/python3.8/site-packages/astropy/units/quantity.py in __array_ufunc__(self, function, method, *inputs, **kwargs)
    484 
    485         # Call our superclass's __array_ufunc__
--> 486         result = super().__array_ufunc__(function, method, *arrays, **kwargs)
    487         # If unit is None, a plain array is expected (e.g., comparisons), which
    488         # means we're done.

ValueError: operands could not be broadcast together with shapes (33,39,2059) (39,33,2059) 

[javerbukh]

@rosteen Is the bug @PatrickOgle found in specutils? It looks like the axes of the slab are getting reversed at some point.

[PatrickOgle]

Filed this jdaviz issue report for Moment = 2 failure:
#947

[PatrickOgle]

As I previously commented, the Gaussian Smooth, Collapse and Moment Maps Plugins appear to work properly, except for an issue with Moment 2 for a NIRSpec IFU dataset.
This lien notwithstanding, the rest looks good.