1# -*- coding: utf-8 -*-
2# Licensed under a 3-clause BSD style license - see LICENSE.rst
3
4"""
5Regression tests for coordinates-related bugs that don't have an obvious other
6place to live
7"""
8
9import io
10import copy
11import pytest
12import numpy as np
13from contextlib import nullcontext
14from erfa import ErfaWarning
15
16from astropy import units as u
17from astropy.coordinates import (
18    AltAz, EarthLocation, SkyCoord, get_sun, ICRS,
19    GeocentricMeanEcliptic, Longitude, Latitude, GCRS, HCRS, CIRS,
20    get_moon, FK4, FK4NoETerms, BaseCoordinateFrame, ITRS,
21    QuantityAttribute, UnitSphericalRepresentation,
22    SphericalRepresentation, CartesianRepresentation,
23    FunctionTransform, get_body,
24    CylindricalRepresentation, CylindricalDifferential,
25    CartesianDifferential)
26from astropy.coordinates.sites import get_builtin_sites
27from astropy.time import Time
28from astropy.utils import iers
29from astropy.table import Table
30
31from astropy.tests.helper import assert_quantity_allclose
32from astropy.utils.compat.optional_deps import HAS_SCIPY  # noqa
33from astropy.units import allclose as quantity_allclose
34
35
36def test_regression_5085():
37    """
38    PR #5085 was put in place to fix the following issue.
39
40    Issue: https://github.com/astropy/astropy/issues/5069
41    At root was the transformation of Ecliptic coordinates with
42    non-scalar times.
43    """
44    # Note: for regression test, we need to be sure that we use UTC for the
45    # epoch, even though more properly that should be TT; but the "expected"
46    # values were calculated using that.
47    j2000 = Time('J2000', scale='utc')
48    times = Time(["2015-08-28 03:30", "2015-09-05 10:30", "2015-09-15 18:35"])
49    latitudes = Latitude([3.9807075, -5.00733806, 1.69539491]*u.deg)
50    longitudes = Longitude([311.79678613, 72.86626741, 199.58698226]*u.deg)
51    distances = u.Quantity([0.00243266, 0.0025424, 0.00271296]*u.au)
52    coo = GeocentricMeanEcliptic(lat=latitudes,
53                                 lon=longitudes,
54                                 distance=distances, obstime=times, equinox=times)
55    # expected result
56    ras = Longitude([310.50095400, 314.67109920, 319.56507428]*u.deg)
57    decs = Latitude([-18.25190443, -17.1556676, -15.71616522]*u.deg)
58    distances = u.Quantity([1.78309901, 1.710874, 1.61326649]*u.au)
59    expected_result = GCRS(ra=ras, dec=decs,
60                           distance=distances, obstime=j2000).cartesian.xyz
61    actual_result = coo.transform_to(GCRS(obstime=j2000)).cartesian.xyz
62    assert_quantity_allclose(expected_result, actual_result)
63
64
65def test_regression_3920():
66    """
67    Issue: https://github.com/astropy/astropy/issues/3920
68    """
69    loc = EarthLocation.from_geodetic(0*u.deg, 0*u.deg, 0)
70    time = Time('2010-1-1')
71
72    aa = AltAz(location=loc, obstime=time)
73    sc = SkyCoord(10*u.deg, 3*u.deg)
74    assert sc.transform_to(aa).shape == tuple()
75    # That part makes sense: the input is a scalar so the output is too
76
77    sc2 = SkyCoord(10*u.deg, 3*u.deg, 1*u.AU)
78    assert sc2.transform_to(aa).shape == tuple()
79    # in 3920 that assert fails, because the shape is (1,)
80
81    # check that the same behavior occurs even if transform is from low-level classes
82    icoo = ICRS(sc.data)
83    icoo2 = ICRS(sc2.data)
84    assert icoo.transform_to(aa).shape == tuple()
85    assert icoo2.transform_to(aa).shape == tuple()
86
87
88def test_regression_3938():
89    """
90    Issue: https://github.com/astropy/astropy/issues/3938
91    """
92    # Set up list of targets - we don't use `from_name` here to avoid
93    # remote_data requirements, but it does the same thing
94    # vega = SkyCoord.from_name('Vega')
95    vega = SkyCoord(279.23473479*u.deg, 38.78368896*u.deg)
96    # capella = SkyCoord.from_name('Capella')
97    capella = SkyCoord(79.17232794*u.deg, 45.99799147*u.deg)
98    # sirius = SkyCoord.from_name('Sirius')
99    sirius = SkyCoord(101.28715533*u.deg, -16.71611586*u.deg)
100    targets = [vega, capella, sirius]
101
102    # Feed list of targets into SkyCoord
103    combined_coords = SkyCoord(targets)
104
105    # Set up AltAz frame
106    time = Time('2012-01-01 00:00:00')
107    location = EarthLocation('10d', '45d', 0)
108    aa = AltAz(location=location, obstime=time)
109
110    combined_coords.transform_to(aa)
111    # in 3938 the above yields ``UnitConversionError: '' (dimensionless) and 'pc' (length) are not convertible``
112
113
114def test_regression_3998():
115    """
116    Issue: https://github.com/astropy/astropy/issues/3998
117    """
118    time = Time('2012-01-01 00:00:00')
119    assert time.isscalar
120
121    sun = get_sun(time)
122    assert sun.isscalar
123    # in 3998, the above yields False - `sun` is a length-1 vector
124
125    assert sun.obstime is time
126
127
128def test_regression_4033():
129    """
130    Issue: https://github.com/astropy/astropy/issues/4033
131    """
132    # alb = SkyCoord.from_name('Albireo')
133    alb = SkyCoord(292.68033548*u.deg, 27.95968007*u.deg)
134    alb_wdist = SkyCoord(alb, distance=133*u.pc)
135
136    # de = SkyCoord.from_name('Deneb')
137    de = SkyCoord(310.35797975*u.deg, 45.28033881*u.deg)
138    de_wdist = SkyCoord(de, distance=802*u.pc)
139
140    aa = AltAz(location=EarthLocation(lat=45*u.deg, lon=0*u.deg), obstime='2010-1-1')
141    deaa = de.transform_to(aa)
142    albaa = alb.transform_to(aa)
143    alb_wdistaa = alb_wdist.transform_to(aa)
144    de_wdistaa = de_wdist.transform_to(aa)
145
146    # these work fine
147    sepnod = deaa.separation(albaa)
148    sepwd = deaa.separation(alb_wdistaa)
149    assert_quantity_allclose(sepnod, 22.2862*u.deg, rtol=1e-6)
150    assert_quantity_allclose(sepwd, 22.2862*u.deg, rtol=1e-6)
151    # parallax should be present when distance added
152    assert np.abs(sepnod - sepwd) > 1*u.marcsec
153
154    # in 4033, the following fail with a recursion error
155    assert_quantity_allclose(de_wdistaa.separation(alb_wdistaa), 22.2862*u.deg, rtol=1e-3)
156    assert_quantity_allclose(alb_wdistaa.separation(deaa), 22.2862*u.deg, rtol=1e-3)
157
158
159@pytest.mark.skipif(not HAS_SCIPY, reason='No Scipy')
160def test_regression_4082():
161    """
162    Issue: https://github.com/astropy/astropy/issues/4082
163    """
164    from astropy.coordinates import search_around_sky, search_around_3d
165    cat = SkyCoord([10.076, 10.00455], [18.54746, 18.54896], unit='deg')
166    search_around_sky(cat[0:1], cat, seplimit=u.arcsec * 60, storekdtree=False)
167    # in the issue, this raises a TypeError
168
169    # also check 3d for good measure, although it's not really affected by this bug directly
170    cat3d = SkyCoord([10.076, 10.00455]*u.deg, [18.54746, 18.54896]*u.deg, distance=[0.1, 1.5]*u.kpc)
171    search_around_3d(cat3d[0:1], cat3d, 1*u.kpc, storekdtree=False)
172
173
174def test_regression_4210():
175    """
176    Issue: https://github.com/astropy/astropy/issues/4210
177    Related PR with actual change: https://github.com/astropy/astropy/pull/4211
178    """
179    crd = SkyCoord(0*u.deg, 0*u.deg, distance=1*u.AU)
180    ecl = crd.geocentricmeanecliptic
181    # bug was that "lambda", which at the time was the name of the geocentric
182    # ecliptic longitude, is a reserved keyword. So this just makes sure the
183    # new name is are all valid
184    ecl.lon
185
186    # and for good measure, check the other ecliptic systems are all the same
187    # names for their attributes
188    from astropy.coordinates.builtin_frames import ecliptic
189    for frame_name in ecliptic.__all__:
190        eclcls = getattr(ecliptic, frame_name)
191        eclobj = eclcls(1*u.deg, 2*u.deg, 3*u.AU)
192
193        eclobj.lat
194        eclobj.lon
195        eclobj.distance
196
197
198def test_regression_futuretimes_4302():
199    """
200    Checks that an error is not raised for future times not covered by IERS
201    tables (at least in a simple transform like CIRS->ITRS that simply requires
202    the UTC<->UT1 conversion).
203
204    Relevant comment: https://github.com/astropy/astropy/pull/4302#discussion_r44836531
205    """
206    from astropy.utils.exceptions import AstropyWarning
207
208    # this is an ugly hack to get the warning to show up even if it has already
209    # appeared
210    from astropy.coordinates.builtin_frames import utils
211    if hasattr(utils, '__warningregistry__'):
212        utils.__warningregistry__.clear()
213
214    # check that out-of-range warning appears among any other warnings.  If
215    # tests are run with --remote-data then the IERS table will be an instance
216    # of IERS_Auto which is assured of being "fresh".  In this case getting
217    # times outside the range of the table does not raise an exception.  Only
218    # if using IERS_B (which happens without --remote-data, i.e. for all CI
219    # testing) do we expect another warning.
220    if isinstance(iers.earth_orientation_table.get(), iers.IERS_B):
221        ctx = pytest.warns(
222            AstropyWarning,
223            match=r'\(some\) times are outside of range covered by IERS table.*')
224    else:
225        ctx = nullcontext()
226    with ctx:
227        future_time = Time('2511-5-1')
228        c = CIRS(1*u.deg, 2*u.deg, obstime=future_time)
229        c.transform_to(ITRS(obstime=future_time))
230
231
232def test_regression_4996():
233    # this part is the actual regression test
234    deltat = np.linspace(-12, 12, 1000)*u.hour
235    times = Time('2012-7-13 00:00:00') + deltat
236    suncoo = get_sun(times)
237    assert suncoo.shape == (len(times),)
238
239    # and this is an additional test to make sure more complex arrays work
240    times2 = Time('2012-7-13 00:00:00') + deltat.reshape(10, 20, 5)
241    suncoo2 = get_sun(times2)
242    assert suncoo2.shape == times2.shape
243
244    # this is intentionally not allclose - they should be *exactly* the same
245    assert np.all(suncoo.ra.ravel() == suncoo2.ra.ravel())
246
247
248def test_regression_4293():
249    """Really just an extra test on FK4 no e, after finding that the units
250    were not always taken correctly.  This test is against explicitly doing
251    the transformations on pp170 of Explanatory Supplement to the Astronomical
252    Almanac (Seidelmann, 2005).
253
254    See https://github.com/astropy/astropy/pull/4293#issuecomment-234973086
255    """
256    # Check all over sky, but avoiding poles (note that FK4 did not ignore
257    # e terms within 10∘ of the poles...  see p170 of explan.supp.).
258    ra, dec = np.meshgrid(np.arange(0, 359, 45), np.arange(-80, 81, 40))
259    fk4 = FK4(ra.ravel() * u.deg, dec.ravel() * u.deg)
260
261    Dc = -0.065838*u.arcsec
262    Dd = +0.335299*u.arcsec
263    # Dc * tan(obliquity), as given on p.170
264    Dctano = -0.028553*u.arcsec
265
266    fk4noe_dec = (fk4.dec - (Dd*np.cos(fk4.ra) -
267                             Dc*np.sin(fk4.ra))*np.sin(fk4.dec) -
268                  Dctano*np.cos(fk4.dec))
269    fk4noe_ra = fk4.ra - (Dc*np.cos(fk4.ra) +
270                          Dd*np.sin(fk4.ra)) / np.cos(fk4.dec)
271
272    fk4noe = fk4.transform_to(FK4NoETerms())
273    # Tolerance here just set to how well the coordinates match, which is much
274    # better than the claimed accuracy of <1 mas for this first-order in
275    # v_earth/c approximation.
276    # Interestingly, if one divides by np.cos(fk4noe_dec) in the ra correction,
277    # the match becomes good to 2 μas.
278    assert_quantity_allclose(fk4noe.ra, fk4noe_ra, atol=11.*u.uas, rtol=0)
279    assert_quantity_allclose(fk4noe.dec, fk4noe_dec, atol=3.*u.uas, rtol=0)
280
281
282def test_regression_4926():
283    times = Time('2010-01-1') + np.arange(20)*u.day
284    green = get_builtin_sites()['greenwich']
285    # this is the regression test
286    moon = get_moon(times, green)
287
288    # this is an additional test to make sure the GCRS->ICRS transform works for complex shapes
289    moon.transform_to(ICRS())
290
291    # and some others to increase coverage of transforms
292    moon.transform_to(HCRS(obstime="J2000"))
293    moon.transform_to(HCRS(obstime=times))
294
295
296def test_regression_5209():
297    "check that distances are not lost on SkyCoord init"
298    time = Time('2015-01-01')
299    moon = get_moon(time)
300    new_coord = SkyCoord([moon])
301    assert_quantity_allclose(new_coord[0].distance, moon.distance)
302
303
304def test_regression_5133():
305    N = 1000
306    np.random.seed(12345)
307    lon = np.random.uniform(-10, 10, N) * u.deg
308    lat = np.random.uniform(50, 52, N) * u.deg
309    alt = np.random.uniform(0, 10., N) * u.km
310
311    time = Time('2010-1-1')
312
313    objects = EarthLocation.from_geodetic(lon, lat, height=alt)
314    itrs_coo = objects.get_itrs(time)
315
316    homes = [EarthLocation.from_geodetic(lon=-1 * u.deg, lat=52 * u.deg, height=h)
317             for h in (0, 1000, 10000)*u.km]
318
319    altaz_frames = [AltAz(obstime=time, location=h) for h in homes]
320    altaz_coos = [itrs_coo.transform_to(f) for f in altaz_frames]
321
322    # they should all be different
323    for coo in altaz_coos[1:]:
324        assert not quantity_allclose(coo.az, coo.az[0])
325        assert not quantity_allclose(coo.alt, coo.alt[0])
326
327
328def test_itrs_vals_5133():
329    """
330    Test to check if alt-az calculations respect height of observer
331
332    Because ITRS is geocentric and includes aberration, an object that
333    appears 'straight up' to a geocentric observer (ITRS) won't be
334    straight up to a topocentric observer - see
335
336    https://github.com/astropy/astropy/issues/10983
337
338    This is worse for small height above the Earth, which is why this test
339    uses large distances.
340    """
341    time = Time('2010-1-1')
342    height = 500000. * u.km
343    el = EarthLocation.from_geodetic(lon=20*u.deg, lat=45*u.deg, height=height)
344
345    lons = [20, 30, 20]*u.deg
346    lats = [44, 45, 45]*u.deg
347    alts = u.Quantity([height, height, 10*height])
348    coos = [EarthLocation.from_geodetic(lon, lat, height=alt).get_itrs(time)
349            for lon, lat, alt in zip(lons, lats, alts)]
350
351    aaf = AltAz(obstime=time, location=el)
352    aacs = [coo.transform_to(aaf) for coo in coos]
353
354    assert all([coo.isscalar for coo in aacs])
355
356    # the ~1 degree tolerance is b/c aberration makes it not exact
357    assert_quantity_allclose(aacs[0].az, 180*u.deg, atol=1*u.deg)
358    assert aacs[0].alt < 0*u.deg
359    assert aacs[0].distance > 5000*u.km
360
361    # it should *not* actually be 90 degrees, b/c constant latitude is not
362    # straight east anywhere except the equator... but should be close-ish
363    assert_quantity_allclose(aacs[1].az, 90*u.deg, atol=5*u.deg)
364    assert aacs[1].alt < 0*u.deg
365    assert aacs[1].distance > 5000*u.km
366
367    assert_quantity_allclose(aacs[2].alt, 90*u.deg, atol=1*u.arcminute)
368    assert_quantity_allclose(aacs[2].distance, 9*height)
369
370
371def test_regression_simple_5133():
372    """
373    Simple test to check if alt-az calculations respect height of observer
374
375    Because ITRS is geocentric and includes aberration, an object that
376    appears 'straight up' to a geocentric observer (ITRS) won't be
377    straight up to a topocentric observer - see
378
379    https://github.com/astropy/astropy/issues/10983
380
381    This is why we construct a topocentric GCRS SkyCoord before calculating AltAz
382    """
383    t = Time('J2010')
384    obj = EarthLocation(-1*u.deg, 52*u.deg, height=[10., 0.]*u.km)
385    home = EarthLocation(-1*u.deg, 52*u.deg, height=5.*u.km)
386
387    obsloc_gcrs, obsvel_gcrs = home.get_gcrs_posvel(t)
388    gcrs_geo = obj.get_itrs(t).transform_to(GCRS(obstime=t))
389    obsrepr = home.get_itrs(t).transform_to(GCRS(obstime=t)).cartesian
390    topo_gcrs_repr = gcrs_geo.cartesian - obsrepr
391    topocentric_gcrs_frame = GCRS(obstime=t, obsgeoloc=obsloc_gcrs, obsgeovel=obsvel_gcrs)
392    gcrs_topo = topocentric_gcrs_frame.realize_frame(topo_gcrs_repr)
393    aa = gcrs_topo.transform_to(AltAz(obstime=t, location=home))
394
395    # az is more-or-less undefined for straight up or down
396    assert_quantity_allclose(aa.alt, [90, -90]*u.deg, rtol=1e-7)
397    assert_quantity_allclose(aa.distance, 5*u.km)
398
399
400def test_regression_5743():
401    sc = SkyCoord([5, 10], [20, 30], unit=u.deg,
402                  obstime=['2017-01-01T00:00', '2017-01-01T00:10'])
403    assert sc[0].obstime.shape == tuple()
404
405
406def test_regression_5889_5890():
407    # ensure we can represent all Representations and transform to ND frames
408    greenwich = EarthLocation(
409        *u.Quantity([3980608.90246817, -102.47522911, 4966861.27310067],
410                    unit=u.m))
411    times = Time("2017-03-20T12:00:00") + np.linspace(-2, 2, 3)*u.hour
412    moon = get_moon(times, location=greenwich)
413    targets = SkyCoord([350.7*u.deg, 260.7*u.deg], [18.4*u.deg, 22.4*u.deg])
414    targs2d = targets[:, np.newaxis]
415    targs2d.transform_to(moon)
416
417
418def test_regression_6236():
419    # sunpy changes its representation upon initialisation of a frame,
420    # including via `realize_frame`. Ensure this works.
421    class MyFrame(BaseCoordinateFrame):
422        default_representation = CartesianRepresentation
423        my_attr = QuantityAttribute(default=0, unit=u.m)
424
425    class MySpecialFrame(MyFrame):
426        def __init__(self, *args, **kwargs):
427            _rep_kwarg = kwargs.get('representation_type', None)
428            super().__init__(*args, **kwargs)
429            if not _rep_kwarg:
430                self.representation_type = self.default_representation
431                self._data = self.data.represent_as(self.representation_type)
432
433    rep1 = UnitSphericalRepresentation([0., 1]*u.deg, [2., 3.]*u.deg)
434    rep2 = SphericalRepresentation([10., 11]*u.deg, [12., 13.]*u.deg,
435                                   [14., 15.]*u.kpc)
436    mf1 = MyFrame(rep1, my_attr=1.*u.km)
437    mf2 = mf1.realize_frame(rep2)
438    # Normally, data is stored as is, but the representation gets set to a
439    # default, even if a different representation instance was passed in.
440    # realize_frame should do the same. Just in case, check attrs are passed.
441    assert mf1.data is rep1
442    assert mf2.data is rep2
443    assert mf1.representation_type is CartesianRepresentation
444    assert mf2.representation_type is CartesianRepresentation
445    assert mf2.my_attr == mf1.my_attr
446    # It should be independent of whether I set the representation explicitly
447    mf3 = MyFrame(rep1, my_attr=1.*u.km, representation_type='unitspherical')
448    mf4 = mf3.realize_frame(rep2)
449    assert mf3.data is rep1
450    assert mf4.data is rep2
451    assert mf3.representation_type is UnitSphericalRepresentation
452    assert mf4.representation_type is CartesianRepresentation
453    assert mf4.my_attr == mf3.my_attr
454    # This should be enough to help sunpy, but just to be sure, a test
455    # even closer to what is done there, i.e., transform the representation.
456    msf1 = MySpecialFrame(rep1, my_attr=1.*u.km)
457    msf2 = msf1.realize_frame(rep2)
458    assert msf1.data is not rep1  # Gets transformed to Cartesian.
459    assert msf2.data is not rep2
460    assert type(msf1.data) is CartesianRepresentation
461    assert type(msf2.data) is CartesianRepresentation
462    assert msf1.representation_type is CartesianRepresentation
463    assert msf2.representation_type is CartesianRepresentation
464    assert msf2.my_attr == msf1.my_attr
465    # And finally a test where the input is not transformed.
466    msf3 = MySpecialFrame(rep1, my_attr=1.*u.km,
467                          representation_type='unitspherical')
468    msf4 = msf3.realize_frame(rep2)
469    assert msf3.data is rep1
470    assert msf4.data is not rep2
471    assert msf3.representation_type is UnitSphericalRepresentation
472    assert msf4.representation_type is CartesianRepresentation
473    assert msf4.my_attr == msf3.my_attr
474
475
476@pytest.mark.skipif(not HAS_SCIPY, reason='No Scipy')
477def test_regression_6347():
478    sc1 = SkyCoord([1, 2]*u.deg, [3, 4]*u.deg)
479    sc2 = SkyCoord([1.1, 2.1]*u.deg, [3.1, 4.1]*u.deg)
480    sc0 = sc1[:0]
481
482    idx1_10, idx2_10, d2d_10, d3d_10 = sc1.search_around_sky(sc2, 10*u.arcmin)
483    idx1_1, idx2_1, d2d_1, d3d_1 = sc1.search_around_sky(sc2, 1*u.arcmin)
484    idx1_0, idx2_0, d2d_0, d3d_0 = sc0.search_around_sky(sc2, 10*u.arcmin)
485
486    assert len(d2d_10) == 2
487
488    assert len(d2d_0) == 0
489    assert type(d2d_0) is type(d2d_10)
490
491    assert len(d2d_1) == 0
492    assert type(d2d_1) is type(d2d_10)
493
494
495@pytest.mark.skipif(not HAS_SCIPY, reason='No Scipy')
496def test_regression_6347_3d():
497    sc1 = SkyCoord([1, 2]*u.deg, [3, 4]*u.deg, [5, 6]*u.kpc)
498    sc2 = SkyCoord([1, 2]*u.deg, [3, 4]*u.deg, [5.1, 6.1]*u.kpc)
499    sc0 = sc1[:0]
500
501    idx1_10, idx2_10, d2d_10, d3d_10 = sc1.search_around_3d(sc2, 500*u.pc)
502    idx1_1, idx2_1, d2d_1, d3d_1 = sc1.search_around_3d(sc2, 50*u.pc)
503    idx1_0, idx2_0, d2d_0, d3d_0 = sc0.search_around_3d(sc2, 500*u.pc)
504
505    assert len(d2d_10) > 0
506
507    assert len(d2d_0) == 0
508    assert type(d2d_0) is type(d2d_10)
509
510    assert len(d2d_1) == 0
511    assert type(d2d_1) is type(d2d_10)
512
513
514def test_gcrs_itrs_cartesian_repr():
515    # issue 6436: transformation failed if coordinate representation was
516    # Cartesian
517    gcrs = GCRS(CartesianRepresentation((859.07256, -4137.20368,  5295.56871),
518                                        unit='km'), representation_type='cartesian')
519    gcrs.transform_to(ITRS())
520
521
522def test_regression_6446():
523    # this succeeds even before 6446:
524    sc1 = SkyCoord([1, 2], [3, 4], unit='deg')
525    t1 = Table([sc1])
526    sio1 = io.StringIO()
527    t1.write(sio1, format='ascii.ecsv')
528
529    # but this fails due to the 6446 bug
530    c1 = SkyCoord(1, 3, unit='deg')
531    c2 = SkyCoord(2, 4, unit='deg')
532    sc2 = SkyCoord([c1, c2])
533    t2 = Table([sc2])
534    sio2 = io.StringIO()
535    t2.write(sio2, format='ascii.ecsv')
536
537    assert sio1.getvalue() == sio2.getvalue()
538
539
540def test_regression_6597():
541    frame_name = 'galactic'
542    c1 = SkyCoord(1, 3, unit='deg', frame=frame_name)
543    c2 = SkyCoord(2, 4, unit='deg', frame=frame_name)
544    sc1 = SkyCoord([c1, c2])
545
546    assert sc1.frame.name == frame_name
547
548
549def test_regression_6597_2():
550    """
551    This tests the more subtle flaw that #6597 indirectly uncovered: that even
552    in the case that the frames are ra/dec, they still might be the wrong *kind*
553    """
554    frame = FK4(equinox='J1949')
555    c1 = SkyCoord(1, 3, unit='deg', frame=frame)
556    c2 = SkyCoord(2, 4, unit='deg', frame=frame)
557    sc1 = SkyCoord([c1, c2])
558
559    assert sc1.frame.name == frame.name
560
561
562def test_regression_6697():
563    """
564    Test for regression of a bug in get_gcrs_posvel that introduced errors at the 1m/s level.
565
566    Comparison data is derived from calculation in PINT
567    https://github.com/nanograv/PINT/blob/master/pint/erfautils.py
568    """
569    pint_vels = CartesianRepresentation(*(348.63632871, -212.31704928, -0.60154936), unit=u.m/u.s)
570    location = EarthLocation(*(5327448.9957829, -1718665.73869569,  3051566.90295403), unit=u.m)
571    t = Time(2458036.161966612, format='jd')
572    obsgeopos, obsgeovel = location.get_gcrs_posvel(t)
573    delta = (obsgeovel-pint_vels).norm()
574    assert delta < 1*u.cm/u.s
575
576
577def test_regression_8138():
578    sc = SkyCoord(1*u.deg, 2*u.deg)
579    newframe = GCRS()
580    sc2 = sc.transform_to(newframe)
581    assert newframe.is_equivalent_frame(sc2.frame)
582
583
584def test_regression_8276():
585    from astropy.coordinates import baseframe
586
587    class MyFrame(BaseCoordinateFrame):
588        a = QuantityAttribute(unit=u.m)
589
590    # we save the transform graph so that it doesn't acidentally mess with other tests
591    old_transform_graph = baseframe.frame_transform_graph
592    try:
593        baseframe.frame_transform_graph = copy.copy(baseframe.frame_transform_graph)
594
595        # as reported in 8276, this previously failed right here because
596        # registering the transform tries to create a frame attribute
597        @baseframe.frame_transform_graph.transform(FunctionTransform, MyFrame, AltAz)
598        def trans(my_frame_coord, altaz_frame):
599            pass
600
601        # should also be able to *create* the Frame at this point
602        MyFrame()
603    finally:
604        baseframe.frame_transform_graph = old_transform_graph
605
606
607def test_regression_8615():
608    # note this is a "higher-level" symptom of the problem that a test now moved
609    # to pyerfa (erfa/tests/test_erfa:test_float32_input) is testing for, but we keep
610    # it here as well due to being a more practical version of the issue.
611
612    crf = CartesianRepresentation(np.array([3, 0, 4], dtype=float) * u.pc)
613    srf = SphericalRepresentation.from_cartesian(crf)  # does not error in 8615
614
615    cr = CartesianRepresentation(np.array([3, 0, 4], dtype='f4') * u.pc)
616    sr = SphericalRepresentation.from_cartesian(cr)  # errors in 8615
617
618    assert_quantity_allclose(sr.distance, 5 * u.pc)
619    assert_quantity_allclose(srf.distance, 5 * u.pc)
620
621
622def test_regression_8924():
623    """This checks that the ValueError in
624    BaseRepresentation._re_represent_differentials is raised properly
625    """
626    # A case where the representation has a 's' differential, but we try to
627    # re-represent only with an 's2' differential
628    rep = CartesianRepresentation(1, 2, 3, unit=u.kpc)
629    dif = CartesianDifferential(4, 5, 6, u.km/u.s)
630    rep = rep.with_differentials(dif)
631
632    with pytest.raises(ValueError):
633        rep._re_represent_differentials(CylindricalRepresentation,
634                                        {'s2': CylindricalDifferential})
635
636
637def test_regression_10092():
638    """
639    Check that we still get a proper motion even for SkyCoords without distance
640    """
641    c = SkyCoord(l=10*u.degree, b=45*u.degree,
642                 pm_l_cosb=34*u.mas/u.yr, pm_b=-117*u.mas/u.yr,
643                 frame='galactic',
644                 obstime=Time('1988-12-18 05:11:23.5'))
645
646    with pytest.warns(ErfaWarning, match='ERFA function "pmsafe" yielded .*'):
647        # expect ErfaWarning here
648        newc = c.apply_space_motion(dt=10*u.year)
649    assert_quantity_allclose(newc.pm_l_cosb, 33.99980714*u.mas/u.yr,
650                             atol=1.0e-5*u.mas/u.yr)
651
652
653def test_regression_10226():
654    # Dictionary representation of SkyCoord should contain differentials.
655    sc = SkyCoord([270, 280]*u.deg, [30, 35]*u.deg, [10, 11]*u.pc,
656                  radial_velocity=[20, -20]*u.km/u.s)
657    sc_as_dict = sc.info._represent_as_dict()
658    assert 'radial_velocity' in sc_as_dict
659    # But only the components that have been specified.
660    assert 'pm_dec' not in sc_as_dict
661
662
663@pytest.mark.parametrize('mjd', (
664    52000, [52000], [[52000]], [52001, 52002], [[52001], [52002]]))
665def test_regression_10422(mjd):
666    """
667    Check that we can get a GCRS for a scalar EarthLocation and a
668    size=1 non-scalar Time.
669    """
670    # Avoid trying to download new IERS data.
671    with iers.earth_orientation_table.set(iers.IERS_B.open(iers.IERS_B_FILE)):
672        t = Time(mjd, format="mjd", scale="tai")
673        loc = EarthLocation(88258.0 * u.m, -4924882.2 * u.m, 3943729.0 * u.m)
674        p, v = loc.get_gcrs_posvel(obstime=t)
675        assert p.shape == v.shape == t.shape
676
677
678@pytest.mark.remote_data
679def test_regression_10291():
680    """
681    According to https://eclipse.gsfc.nasa.gov/OH/transit12.html,
682    the minimum separation between Venus and the Sun during the 2012
683    transit is 554 arcseconds for an observer at the Geocenter.
684
685    If light deflection from the Sun is incorrectly applied, this increases
686    to 557 arcseconds.
687    """
688    t = Time('2012-06-06 01:29:36')
689    sun = get_body('sun', t)
690    venus = get_body('venus', t)
691    assert_quantity_allclose(venus.separation(sun),
692                             554.427*u.arcsecond, atol=0.001*u.arcsecond)
693