Implement spectral coordinate object

[nmearl]

Resolves #422. See the astropy PR for further discussion (which should take place here now).

[keflavich]

This looks good, but it's also way too much for me to digest right now. Is this the object I should be using as a mixin in spectral-cube or as the spectral_axis object? It's the latter. This is great.

[nmearl]

@eteq @astrofrog Any further opinions on this approach?

[astrofrog]

I ran into an issue with the following example:

In [1]: from astropy import units as u 
   ...: from astropy.coordinates import SkyCoord                                                                                                       

In [2]: LSRD = SkyCoord(u=0 * u.km, v=0 * u.km, w=0 * u.km, 
   ...:                 U=9 * u.km / u.s, V=12 * u.km / u.s, W=7 * u.km / u.s, 
   ...:                 frame='galactic', representation_type='cartesian', differential_type='cartesian')                                              

In [3]: from specutils.spectra.spectral_coordinate import SpectralCoord                                                                                

In [4]: spectral_axis = SpectralCoord([3, 4, 5] * u.micron, observer=LSRD)                                                                             
---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
<ipython-input-4-947d16e15fcb> in <module>
----> 1 spectral_axis = SpectralCoord([3, 4, 5] * u.micron, observer=LSRD)

~/Dropbox/Code/Astropy/specutils/specutils/spectra/spectral_coordinate.py in __new__(cls, value, unit, radial_velocity, redshift, rest, velocity_convention, observer, target, **kwargs)
     58         obj.velocity_convention = velocity_convention
     59 
---> 60         obj.observer = observer
     61         obj.target = target
     62 

~/Dropbox/Code/Astropy/specutils/specutils/spectra/spectral_coordinate.py in observer(self, value)
    145     @observer.setter
    146     def observer(self, value):
--> 147         value = self._validate_frame(value, True)
    148         self._observer = value
    149 

~/Dropbox/Code/Astropy/specutils/specutils/spectra/spectral_coordinate.py in _validate_frame(frame, is_observer)
    110             # If no distance value is defined on the frame, assume a default
    111             # distance of either zero (for an observer), or 1000 kpc for target
--> 112             if not isinstance(frame.distance, Distance):
    113                 auto_dist = 0 * u.AU if is_observer else 1000 * u.kpc
    114 

~/python/dev/lib/python3.7/site-packages/astropy/coordinates/baseframe.py in __getattr__(self, attr)
   1543             return val
   1544 
-> 1545         return self.__getattribute__(attr)  # Raise AttributeError.
   1546 
   1547     def __setattr__(self, attr, value):

AttributeError: 'Galactic' object has no attribute 'distance'

I haven't had a chance to investigate - am I doing something that shouldn't work?

[astrofrog]

Another thing which isn't an issue as such but would be nice to have work: if we pass an EarthLocation to observer it would be nice to translate on-the-fly to a coordinate object:

spectral_axis = SpectralCoord([3, 4, 5] * u.micron, observer=EarthLocation.of_site('alma'))

At the moment this gives:

ValueError: `(2225015.30883296, -5440016.41799762, -2481631.27428014) m` is not a subclass of `BaseCoordinateFrame` or `SkyCoord`.

[astrofrog]

If I modify the example with EarthLocation I still get an error related to distance:

In [5]: from astropy.coordinates import EarthLocation                                                                                                  

In [6]: alma = EarthLocation.of_site('ALMA').itrs                                                                                                      

In [7]: from specutils.spectra.spectral_coordinate import SpectralCoord                                                                                

In [8]: spectral_axis = SpectralCoord([3, 4, 5] * u.micron, observer=alma)                                                                             
---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
<ipython-input-8-9017f7b95060> in <module>
----> 1 spectral_axis = SpectralCoord([3, 4, 5] * u.micron, observer=alma)

~/Dropbox/Code/Astropy/specutils/specutils/spectra/spectral_coordinate.py in __new__(cls, value, unit, radial_velocity, redshift, rest, velocity_convention, observer, target, **kwargs)
     58         obj.velocity_convention = velocity_convention
     59 
---> 60         obj.observer = observer
     61         obj.target = target
     62 

~/Dropbox/Code/Astropy/specutils/specutils/spectra/spectral_coordinate.py in observer(self, value)
    145     @observer.setter
    146     def observer(self, value):
--> 147         value = self._validate_frame(value, True)
    148         self._observer = value
    149 

~/Dropbox/Code/Astropy/specutils/specutils/spectra/spectral_coordinate.py in _validate_frame(frame, is_observer)
    110             # If no distance value is defined on the frame, assume a default
    111             # distance of either zero (for an observer), or 1000 kpc for target
--> 112             if not isinstance(frame.distance, Distance):
    113                 auto_dist = 0 * u.AU if is_observer else 1000 * u.kpc
    114 

~/python/dev/lib/python3.7/site-packages/astropy/coordinates/baseframe.py in __getattr__(self, attr)
   1543             return val
   1544 
-> 1545         return self.__getattribute__(attr)  # Raise AttributeError.
   1546 
   1547     def __setattr__(self, attr, value):

AttributeError: 'ITRS' object has no attribute 'distance'

[nmearl]

Another thing which isn't an issue as such but would be nice to have work: if we pass an EarthLocation to observer it would be nice to translate on-the-fly to a coordinate object:

@astrofrog I don't think we can support this, as currently to convert an EarthLocation to a SkyCoord requires that the user provide an observation time and frame:

alma = EarthLocation.of_site('alma')
obstime = time.Time('2018-12-13 9:00')

observer_gcrs = alma.get_gcrs(obstime)

[nmearl]

Ah, well, I stand corrected: if the frame doesn't require it, like in your example (ITRS), there's no need for the observation time.

[astrofrog]

@nmearl - thanks for the fixes! I'm running into an issue related to (I guess?) the repr for SpectralCoord:

This is using the same notebook as before, just with an extra spectral_axis cell. I think this should probably work since I didn't actually ask for any velocities? The second example in my notebook has a slightly different error:

And another repr-related failure:

[astrofrog]

@nmearl - I got further this time! See

https://gist.github.com/astrofrog/b11eb8e0213b682658789b984f63cd9f

for the remaining issues. In short:

• repr often errors even if the info mentioned by the error shouldn't be needed
• I can't figure out how to set the observation time when creating an observer frame on Earth
• Adding zero velocities to observatories on the Earth's surface is tedious at the moment, and maybe EarthLocations, when converted to ITRS, should have zero velocities by default?
• I was having issues transforming SpectralCoord from an Earth-based frame to another frame

Thanks for continuing to work on this!

[nmearl]

@astrofrog Thanks for looking this over! A few restrictions come to mind: it sounded like we shouldn't be in the business of transforming frames (@eteq?). In the case where we want to change observer, we should really only allow SkyCoordinates to the with_observer method, or another observer-like object?

For your last line: spectral_axis.with_observer(LSRD), I'm not sure what behavior you're trying to get? I can't even transform the alma_with_vel to the LSRD frame without astropy giving me an error:

alma_lsrd = alma_with_vel.transform_to(LSRD)

---------------------------------------------------------------------------
NotImplementedError                       Traceback (most recent call last)
<ipython-input-30-d9e9fba3f3c1> in <module>
      1 print(LSRD.velocity)
----> 2 alma_lsrd = alma_with_vel.transform_to(LSRD)
      3 spectral_axis.with_observer(alma_lsrd)

~/anaconda3/envs/specutils_dev/lib/python3.7/site-packages/astropy/coordinates/baseframe.py in transform_to(self, new_frame)
   1161            hasattr(self, 'obstime') and hasattr(new_frame, 'obstime') and
   1162            np.any(self.obstime != new_frame.obstime)):
-> 1163             raise NotImplementedError('You cannot transform a frame that has '
   1164                                       'velocities to another frame at a '
   1165                                       'different obstime. If you think this '

NotImplementedError: You cannot transform a frame that has velocities to another frame at a different obstime. If you think this should (or should not) be possible, please comment at https://github.com/astropy/astropy/issues/6280

[astrofrog]

Just to be clear I didn't mean changing the frame but indeed changing the observer via SkyCoord.

Regarding the last line:

alma_lsrd = alma_with_vel.transform_to(LSRD)

what I'm really trying to accomplish is to transform the velocities to an observer that is stationary compared to LSRD. I'm surprised by the error you showed when directly transforming the coordinates

In general, I think what I've realized by thinking about this some more is that I don't think ever really want to change the actual position of the observer, but it's more that we want to change to an observer at the same position that is stationary relative to the specified SkyCoord?

I think it might be productive to have a short voice chat about all this between the three of us to iron out some of these remaining issues?

[astrofrog]

I think the transformation to LSRD should really work, the LSRD obstime is not set: astropy/astropy#6280 (comment)

[astrofrog]

I've now played around with this more again since we last spoke. Latest notebook:

https://gist.github.com/astrofrog/11d925ab02282440709f63114221e650

The bottom line is that I still had issues with the case where the observer location changes and we take into account geometrical issues, but I managed to get the case where the observer location stays the same and the velocities change, although there were issues with the units on the result.

[eteq]

I've left several inline comments, some of which are relatively minor but a couple of which are rather significant implementation items.

I'm also pushing up a few tests that test out the new behavior as discussed in #422 (comment) . Sadly, some of those are currently failing. Some of that is fixed by the two suggestions I made (or a similar fix should do it), but some of it looks like a different set of failures. It's not clear to me exactly what the origin of those failures are, but hopefully the tests provide enough of a guidance on the direction to finish it out. So I think this will be good once the review comments are addressed and the new tests are passing.

[eteq]

(oh and note that the tests I pushed up also include 8b512a0 which I think fixes the mysterious failure due to not finding keck by switching to a builtin site that requires no download)