mzechmeister / serval

calculate radial velocities from stellar spectra
MIT License
36 stars 9 forks source link

astropy WARNING: ErfaWarning: ERFA function "dtf2d" yielded 1 of "dubious year (Note 6)" [astropy._erfa.core] #48

Open mzechmeister opened 2 years ago

mzechmeister commented 2 years ago

A new astropy warning (only python2)

WARNING: ErfaWarning: ERFA function "dtf2d" yielded 1 of "dubious year (Note 6)" [astropy._erfa.core]

https://github.com/nu-radio/NuRadioReco/issues/99 https://github.com/astropy/astropy/issues/5809 https://github.com/astropy/astropy/issues/9603

Arising in https://github.com/mzechmeister/serval/blob/f26706f77b2505a94029c8bad3f994b9aed03392/src/brv_we14py.py#L75

Minimal script to create the warning:

from astropy.utils.iers import conf as iers_conf
iers_conf.iers_auto_url = 'https://datacenter.iers.org/data/9/finals2000A.all'

from astropy.time import Time
from astropy import coordinates as coord
targ = coord.SkyCoord(0, 0, unit=('deg', 'deg'))
loc = coord.EarthLocation.from_geodetic(0, 0)
Time(2457395.2, format='jd', scale='utc').light_travel_time(targ, location=loc)
mzechmeister commented 2 years ago

The stack

.light_travel_time(targ, location=loc)
  ~/.local/lib/python2.7/site-packages/astropy/time/core.py(682)light_travel_time()
-> gcrs_coo = itrs.transform_to(GCRS(obstime=self))
  ~/.local/lib/python2.7/site-packages/astropy/coordinates/baseframe.py(934)transform_to()
-> return trans(self, new_frame)
  ~/.local/lib/python2.7/site-packages/astropy/coordinates/transformations.py(1314)__call__()
-> curr_coord = t(curr_coord, curr_toframe)
  ~/.local/lib/python2.7/site-packages/astropy/coordinates/transformations.py(914)__call__()
-> return supcall(fromcoord, toframe)
  ~/.local/lib/python2.7/site-packages/astropy/coordinates/builtin_frames/intermediate_rotation_transforms.py(89)itrs_to_cirs()
-> pmat = cirs_to_itrs_mat(itrs_coo.obstime)
  ~/.local/lib/python2.7/site-packages/astropy/coordinates/builtin_frames/intermediate_rotation_transforms.py(33)cirs_to_itrs_mat()
-> xp, yp = get_polar_motion(time)
> ~/.local/lib/python2.7/site-packages/astropy/coordinates/builtin_frames/utils.py(43)get_polar_motion()
-> xp, yp, status = iers.IERS_Auto.open().pm_xy(time, return_status=True)

More minimal:

from astropy.utils import iers
iers.conf.iers_auto_url = 'https://datacenter.iers.org/data/9/finals2000A.all'
from astropy.time import Time

t = Time(2457395.2, format='jd', scale='utc')
iers.IERS_Auto.open().pm_xy(t)

Now we are here:

from astropy.time import Time
Time.now()

Further, in core.py (line 19603, https://github.com/astropy/astropy/blob/4ea556e14340760bd58e53ea2cbccde7ce5567c3/astropy/_erfa/core.py.templ#L78)

    stat_ok = _core._dtf2d(it)

calls a shared library.

(Pdb) _core._dtf2d(it)
False
(Pdb) it
<numpy.nditer object at 0x7f54d0201170>
(Pdb) arrs
[array('UTC', dtype='|S16'), array(2022, dtype=int32), array(1, dtype=int32), array(13, dtype=int32), array(14, dtype=int32), array(38, dtype=int32), array(24.352779), array(2459592.5), array(0.61000408), array(1, dtype=int32)]
_core._dtf2d(it)

I found the core.c in https://files.pythonhosted.org/packages/2b/5d/d38bcc7dfa6d84b8cd720f944a4cfdb893cda10a926786f55350ab86c137/astropy-2.0.3.tar.gz

l4700: static PyObject *Py_dtf2d(PyObject *self, PyObject *args, PyObject *kwds)
l4730:    _c_retval = eraDtf2d(_scale, *_iy, *_im, *_id, *_ihr, *_imn, *_sec, _d1, _d2);

The function may look like: https://github.com/liberfa/erfa/blob/5233e79cd97c2d3fdae4b7a8f802d46be46867a5/src/dtf2d.c#L5

mzechmeister commented 2 years ago

As noted in https://github.com/astropy/astropy/issues/9603#issuecomment-1014861407 it might be difficult to get barycentric correction working for astropy 2.