AltAz transformation fails for vector SkyCoord

[bmorris3]

I'm calling a coord.transform_to(AltAz(...)) transformation on a vector of coordinates, like this:

from astropy.coordinates import SkyCoord, EarthLocation, AltAz
from astropy.time import Time

# Set up list of targets
vega = SkyCoord.from_name('Vega')
capella = SkyCoord.from_name('Capella')
sirius = SkyCoord.from_name('Sirius')
targets = [vega, capella, sirius]

# Feed list of targets into SkyCoord
combined_coords = SkyCoord(targets)
print(combined_coords)
# Set up AltAz frame
time = Time('2015-01-01 00:00:00')
location = EarthLocation('10d', '45d', 0)
aa = AltAz(location=location, obstime=time)
print(combined_coords.transform_to(aa))

and getting the following error:

In [12]: run exp
<SkyCoord (ICRS): (ra, dec, distance) in (deg, deg, )
    [(279.23473479, 38.78368896, 1.0), (79.17232794, 45.99799147, 1.0),
     (101.28715533, -16.71611586, 1.0)]>
---------------------------------------------------------------------------
UnitConversionError                       Traceback (most recent call last)
/Users/bmorris/git/astroplan/exp.py in <module>()
     22 location = EarthLocation('10d', '45d', 0)
     23 aa = AltAz(location=location, obstime=time)
---> 24 print(combined_coords.transform_to(aa))
     25 
     26 

/Users/bmorris/git/astropy/astropy/coordinates/sky_coordinate.pyc in transform_to(self, frame)
    385         # Do the transformation, returning a coordinate frame of the desired
    386         # final type (not generic).
--> 387         new_coord = trans(self.frame, generic_frame)
    388 
    389         # Finally make the new SkyCoord object from the `new_coord` and

/Users/bmorris/git/astropy/astropy/coordinates/transformations.pyc in __call__(self, fromcoord, toframe)
    918 
    919             curr_toframe = t.tosys(**frattrs)
--> 920             curr_coord = t(curr_coord, curr_toframe)
    921 
    922         # this is safe even in the case where self.transforms is empty, because

/Users/bmorris/git/astropy/astropy/coordinates/transformations.pyc in __call__(self, fromcoord, toframe)
    706 
    707     def __call__(self, fromcoord, toframe):
--> 708         res = self.func(fromcoord, toframe)
    709         if not isinstance(res, self.tosys):
    710             raise TypeError('the transformation function yielded {0} but '

/Users/bmorris/git/astropy/astropy/coordinates/builtin_frames/icrs_cirs_transforms.pyc in icrs_to_cirs(icrs_coo, cirs_frame)
     29         px = 0
     30     else:
---> 31         px = 1 / icrs_coo.distance.to(u.parsec).value
     32     srepr = icrs_coo.represent_as(UnitSphericalRepresentation)
     33     i_ra = srepr.lon.to(u.radian).value

/Users/bmorris/git/astropy/astropy/units/quantity.pyc in to(self, unit, equivalencies)
    619         unit = Unit(unit)
    620         new_val = np.asarray(
--> 621             self.unit.to(unit, self.value, equivalencies=equivalencies))
    622         return self._new_view(new_val, unit)
    623 

/Users/bmorris/git/astropy/astropy/units/core.pyc in to(self, other, value, equivalencies)
    970             If units are inconsistent
    971         """
--> 972         return self._get_converter(other, equivalencies=equivalencies)(value)
    973 
    974     def in_units(self, other, value=1.0, equivalencies=[]):

/Users/bmorris/git/astropy/astropy/units/core.pyc in _get_converter(self, other, equivalencies)
    871         except UnitsError:
    872             return self._apply_equivalencies(
--> 873                 self, other, self._normalize_equivalencies(equivalencies))
    874         return lambda val: scale * _condition_arg(val)
    875 

/Users/bmorris/git/astropy/astropy/units/core.pyc in _apply_equivalencies(self, unit, other, equivalencies)
    862         raise UnitConversionError(
    863             "{0} and {1} are not convertible".format(
--> 864                 unit_str, other_str))
    865 
    866     def _get_converter(self, other, equivalencies=[]):

UnitConversionError: '' (dimensionless) and 'pc' (length) are not convertible

As you can see by the print statement that does not fail, the SkyCoord constructor unpacks the coordinates (which only have RA, Dec) into a vector with distances. What's that about? @eteq @cdeil @adrn

[adrn]

Looks like it's because your EarthLocation height kwarg has no units (dimensionless). I'd say this is a bug since 0 can be interpreted unambiguously even without units. @bmorris3 best way to get around this is to always pass in quantities:

location = EarthLocation('10d', '45d', 0*u.meter)

or, while we're at it:

location = EarthLocation(10*u.degree, 45*u.degree, 0*u.meter)

you crazy people with your strings...:)

[adrn]

BTW I don't know why this only happens for non-scalar SkyCoord's -- that's weird...

[bmorris3]

Hmm, every other time I've used EarthLocations with a unitless elevation I haven't had issues. The EarthLocation docs say:

height : Quantity or float, optional
Height above reference ellipsoid (if float, in meters; default: 0).

So I think what I typed should have worked.

[bmorris3]

@adrn 's suggested work-around of using a unit in the height kwarg to EarthLocation actually doesn't help. Using location = EarthLocation(10*u.deg, 45*u.deg, 0*u.m) I still get the same error.

[adrn]

@bmorris3 oops, you're right! when I tested it on my machine I used a scalar coordinate...doh

[bmorris3]

I think the error comes from /coordinates/builtin_frames/icrs_cirs_transforms.py with the expectation that there is a distance associated with the SkyCoord. The code should be going into the px=0 case since there is no distance, but it is instead trying to calculate the parallax.

[bmorris3]

It appears that SkyCoords initialized without a distance get a default dimensionless distance of 1 assigned to them:

In [36]: vega = SkyCoord.from_name('Vega')

In [37]: vega.distance
Out[37]: <Quantity 1.0>

Could it be possible that the default distance of <Quantity 1.0> for UnitSphericalRepresentation is being misinterpreted in the transformation, @eteq ?

[adrn]

Ah, no, it's actually that combined_coords for some reason doesn't have a UnitSphericalRepresentation:

>>> combined_coords = SkyCoord(targets)
>>> print(combined_coords.representation)
astropy.coordinates.representation.SphericalRepresentation

[adrn]

Hm, well, that's because the default representation is spherical not unitspherical. @eteq I don't remember why this is, and I don't understand that...but in any case, the frames get initialized with dimensionless distances which causes trouble for the CIRS transformation because the way it checks for an uninitialized distance is to do:

if isinstance(icrs_coo.data, UnitSphericalRepresentation):  # no distance

[eteq]

Investigated this a bit. Turns out the confusion is that representation is not necessarily the internal value of .data - it's the "how it appears" representation. So indeed the underlying problem is that the UnitSphericalRepresentation of the three inputs gets silently converted to SphericalRepresentation, which means a distance that doesn't make sense.

I think the solution is to just correctly convince the new combined_coords to be a UnitSphericalRepresentation. Self-assigned and I'll work on a PR for this.