shbhuk / barycorrpy

Python version of Barycorr
GNU General Public License v3.0
37 stars 6 forks source link

USNO IERS is offline - barycorrpy frozen. #27

Closed njcuk9999 closed 2 years ago

njcuk9999 commented 4 years ago

Took me several hours to figure out what was wrong.

This issue here explain more of the problem one would think is initially an astropy.utils.iers problem.

Basically barycorrpy is unable to fetch anything from this server (as of Oct 24th 2019 and will continue to be this way until at least April 2020!!), for the finals2000A.all (iers.IERS_A) this is solved by updating astropy (they now provide a mirror and different link for IERS_A) i.e.:

IERS_A_URL = 'https://datacenter.iers.org/data/9/finals2000A.all'
IERS_A_URL_MIRROR = 'ftp://cddis.gsfc.nasa.gov/pub/products/iers/finals2000A.all'

and/or providing an offline file in a hacky way:

        from astropy.utils import iers
        iers_a_file = 'my/path/to/finals2000A.all`
        iers.IERS.iers_table = iers.IERS_A.open(iers_a_file)
        import barycorrpy

which is not the best solution (should be settable via barycorrpy)

however in barycorrpy.utc_tdb.py (Line 59) you have a hard coded link to http://maia.usno.navy.mil/ser7/tai-utc.dat which is and will be offline for a long time and hence on my machine just causes a code freeze and does not allow any calculation, and I cannot find anywhere else to download this file. I know I can turn off updating for leap seconds but this seems against the point.

shbhuk commented 4 years ago

Hi @njcuk9999 , Thanks for raising this issue. I did come across this notice and was waiting to see how Astropy deals with it, and what the USNO recommendations are.

As far as leap seconds area concerned it will be safe to turn off the update until the server comes back up since there will be no new leap seconds at least till June 2020 - https://datacenter.iers.org/data/latestVersion/16_BULLETIN_C16.txt. The next leap second would have been for December 2019 onwards, but there won't be one as this notice mentions.

I'll look into the IERS table update provided by Astropy..

shbhuk commented 4 years ago

I'm surprised barycorrpy freezes when it tries to update the leap second, The way I have coded it is that it should just default to the Astropy Time object converter when the leap second update from the USNO site fails.

Could you please send me a minimum working example that I can run and troubleshoot.

njcuk9999 commented 4 years ago

Sorry for not getting back to you sooner. So it doesn't freeze it just takes way longer than before to finish (i.e. ~20 second for just 3 times), this is on a large linux machine with 24 cores and 256 GB ram so resources isn't a problem.

Minimum working example:


import barycorrpy
import numpy as np
import time

bkwargs = {'ra': 269.4486143580427,
 'dec': 4.737980765585507,
 'epoch': 2457205.9999942128,
 'px': 547.4506196448664,
 'pmra': -802.802647878488,
 'pmdec': 10362.542163273869,
 'lat': 19.825252,
 'longi': -155.468876,
 'alt': 4204.0,
 'rv': 0.0}

times = np.array([2458787.70556851, 2458789.37223517, 2458791.03890184])

print('barycorrpy.get_BC_vel')
start = time.time()
out1 = barycorrpy.get_BC_vel(JDUTC=times, zmeas=0.0, **bkwargs)
end = time.time()
print('\tTime taken = {0:.3f} s'.format(end - start))

print('barycorrpy.utc_tdb.JDUTC_to_BJDTDB')
start = time.time()
out2 = barycorrpy.utc_tdb.JDUTC_to_BJDTDB(times, **bkwargs)
end = time.time()
print('\tTime taken = {0:.3f} s'.format(end - start))

output is this:

barycorrpy.get_BC_vel
        Time taken = 21.170 s
barycorrpy.utc_tdb.JDUTC_to_BJDTDB
        Time taken = 12.676 s

Cancelling it (with Ctrl+C) brings up this dialog:

^C---------------------------------------------------------------------------
KeyboardInterrupt                         Traceback (most recent call last)
<ipython-input-1-96541a78d561> in <module>
     17 
     18 print('barycorrpy.get_BC_vel')
---> 19 out1 = barycorrpy.get_BC_vel(JDUTC=times, zmeas=0.0, **bkwargs)
     20 print('barycorrpy.utc_tdb.JDUTC_to_BJDTDB')
     21 out2 = barycorrpy.utc_tdb.JDUTC_to_BJDTDB(times, **bkwargs)

~/anaconda3/lib/python3.7/site-packages/barycorrpy/barycorrpy.py in get_BC_vel(JDUTC, starname, hip_id, ra, dec, epoch, pmra, pmdec, px, rv, obsname, lat, longi, alt, zmeas, ephemeris, leap_dir, leap_update)               
    157                  zmeas=zm,
    158                  loc=loc,
--> 159                  ephemeris=ephemeris, leap_dir=leap_dir, leap_update=leap_update,**star_output)
    160         vel.append(a[0])
    161         warning.append(a[1])

~/anaconda3/lib/python3.7/site-packages/barycorrpy/barycorrpy.py in BCPy(JDUTC, ra, dec, epoch, pmra, pmdec, px, rv, zmeas, loc, ephemeris, leap_dir, leap_update)                                                            
    185 
    186     # Convert times to obtain TDB and TT
--> 187     JDTDB, JDTT, warning, error = utc_tdb.JDUTC_to_JDTDB(JDUTC, fpath=leap_dir, leap_update=leap_update)                                                                                                              
    188 
    189     ##### NUTATION, PRECESSION, ETC. #####

~/anaconda3/lib/python3.7/site-packages/barycorrpy/utc_tdb.py in JDUTC_to_JDTDB(utctime, leap_update, fpath)
    199 
    200     # Call function to check leap second file and find offset between UTC and TAI.
--> 201     tai_utc,warning,error=leap_manage(utctime=utctime,fpath=fpath,leap_update=leap_update)
    202 
    203     if utctime.scale != 'utc':

~/anaconda3/lib/python3.7/site-packages/barycorrpy/utc_tdb.py in leap_manage(utctime, fpath, leap_update)
    117     if (os.path.isfile(log_fpath)==False or os.path.isfile(ls_fpath)==False):
    118         if leap_update==True:
--> 119             flag=leap_download(ls_fpath=ls_fpath,log_fpath=log_fpath)
    120             if flag==0:
    121                 warning+=['WARNING JD = '+str(utctime.jd)+' : Downloaded leap second file from http://maia.usno.navy.mil/ser7/tai-utc.dat ']                                                                          

~/anaconda3/lib/python3.7/site-packages/barycorrpy/utc_tdb.py in leap_download(ls_fpath, log_fpath)
     61     if sys.version_info.major==3:
     62         try:
---> 63             urllib.request.urlretrieve( url,ls_fpath)
     64             flag=0
     65             with open(log_fpath,'w') as f:

~/anaconda3/lib/python3.7/urllib/request.py in urlretrieve(url, filename, reporthook, data)
    245     url_type, path = splittype(url)
    246 
--> 247     with contextlib.closing(urlopen(url, data)) as fp:
    248         headers = fp.info()
    249 

~/anaconda3/lib/python3.7/urllib/request.py in urlopen(url, data, timeout, cafile, capath, cadefault, context)
    220     else:
    221         opener = _opener
--> 222     return opener.open(url, data, timeout)
    223 
    224 def install_opener(opener):

~/anaconda3/lib/python3.7/urllib/request.py in open(self, fullurl, data, timeout)
    523             req = meth(req)
    524 
--> 525         response = self._open(req, data)
    526 
    527         # post-process response

~/anaconda3/lib/python3.7/urllib/request.py in _open(self, req, data)
    541         protocol = req.type
    542         result = self._call_chain(self.handle_open, protocol, protocol +
--> 543                                   '_open', req)
    544         if result:
    545             return result

~/anaconda3/lib/python3.7/urllib/request.py in _call_chain(self, chain, kind, meth_name, *args)
    501         for handler in handlers:
    502             func = getattr(handler, meth_name)
--> 503             result = func(*args)
    504             if result is not None:
    505                 return result

~/anaconda3/lib/python3.7/urllib/request.py in http_open(self, req)
   1343 
   1344     def http_open(self, req):
-> 1345         return self.do_open(http.client.HTTPConnection, req)
   1346 
   1347     http_request = AbstractHTTPHandler.do_request_

~/anaconda3/lib/python3.7/urllib/request.py in do_open(self, http_class, req, **http_conn_args)
   1318             except OSError as err: # timeout error
   1319                 raise URLError(err)
-> 1320             r = h.getresponse()
   1321         except:
   1322             h.close()

~/anaconda3/lib/python3.7/http/client.py in getresponse(self)
   1334         try:
   1335             try:
-> 1336                 response.begin()
   1337             except ConnectionError:
   1338                 self.close()

~/anaconda3/lib/python3.7/http/client.py in begin(self)
    304         # read until we get a non-100 response
    305         while True:
--> 306             version, status, reason = self._read_status()
    307             if status != CONTINUE:
    308                 break

~/anaconda3/lib/python3.7/http/client.py in _read_status(self)
    265 
    266     def _read_status(self):
--> 267         line = str(self.fp.readline(_MAXLINE + 1), "iso-8859-1")
    268         if len(line) > _MAXLINE:
    269             raise LineTooLong("status line")

~/anaconda3/lib/python3.7/socket.py in readinto(self, b)
    587         while True:
    588             try:
--> 589                 return self._sock.recv_into(b)
    590             except timeout:
    591                 self._timeout_occurred = True
njcuk9999 commented 4 years ago

Note that I thought it froze completely because normal the length of "times" is 200 to 300 so something is going on inside a loop around individual "time" calculations, and that setting leap_update=False means the same code is much faster:

barycorrpy.get_BC_vel
        Time taken = 0.026 s
barycorrpy.utc_tdb.JDUTC_to_BJDTDB
        Time taken = 0.019 s
shbhuk commented 4 years ago

What is your numpy and astropy version numbers?

njcuk9999 commented 4 years ago

astropy 3.2.3 (just updated for new mirror for iers) numpy 1.17.2

platform: Linux platform-release: 3.10.0-514.26.2.el7.x86_64 platform-version: 1 SMP Tue Jul 4 15:04:05 UTC 2017 architecture: x86_64 ram: 252 GB

shbhuk commented 4 years ago

Hi, That is interesting -

I first updated my astropy from 3.2.1 to 3.2.3, and then numpy from 1.15.0 to 1.17.3, then I ran your code.

See below -

import barycorrpy
import numpy as np
import time

import astropy
print(astropy.__version__)
print(np.__version__)

bkwargs = {'ra': 269.4486143580427,
 'dec': 4.737980765585507,
 'epoch': 2457205.9999942128,
 'px': 547.4506196448664,
 'pmra': -802.802647878488,
 'pmdec': 10362.542163273869,
 'lat': 19.825252,
 'longi': -155.468876,
 'alt': 4204.0,
 'rv': 0.0}

times = np.array([2458787.70556851, 2458789.37223517, 2458791.03890184])
# times = np.array([2458000])

print('barycorrpy.get_BC_vel')
start = time.time()
out1 = barycorrpy.get_BC_vel(JDUTC=times, zmeas=0.0, **bkwargs)
end = time.time()
print('\tTime taken = {0:.3f} s'.format(end - start))

print('barycorrpy.utc_tdb.JDUTC_to_BJDTDB')
start = time.time()
out2 = barycorrpy.utc_tdb.JDUTC_to_BJDTDB(times, **bkwargs)
end = time.time()
print('\tTime taken = {0:.3f} s'.format(end - start))

The output is -

3.2.3
1.17.3
barycorrpy.get_BC_vel
        Time taken = 0.119 s
barycorrpy.utc_tdb.JDUTC_to_BJDTDB
        Time taken = 0.030 

I'll look into this more and issue an update soon.

shbhuk commented 4 years ago

The issue is that if I require Astropy>=3.2.3 in the setup.py file then I am basically excluding Python 2.x and closing support for it.

njcuk9999 commented 4 years ago

It might be that you already have a leap second file stored somewhere so it doesn't try querying the server (which is hardcoded and is down until April 2020) - this could be where the massive amount of time lag comes from for each loop operation

cusher commented 4 years ago

The issue is that if I require Astropy>=3.2.3 in the setup.py file then I am basically excluding Python 2.x and closing support for it.

Note that astropy 2.0.16 also includes the fixed IERS URLs. I don't think that would do much to solve the problem of issues caused by the offline leap second file though.

shbhuk commented 4 years ago

Exactly! The best way to make the code go faster might just be too turn off leap_update by setting it to False

leap_update=False when you call get_BC_vel()

I have summarized this in the wiki, let me know if there are any lingering doubts or issues that I can help with.

@cusher you're right about the Astropy update 2.0.16. If my next release is before the servers come back up, I may add this requirement of Astropy >=3.2.3 for Python 3.x and >=2.0.16 for Python 2.x, which would necessitate releasing separate versions for the two.

shbhuk commented 4 years ago

I'll just leave this issue open for others to see when they get this error, in case they do not go to the wiki immediately.

shbhuk commented 3 years ago

As per PR #41 , since leap seconds are now being handled by astropy from barycorrpy v0.4.0 onwards, this is only a concern for IERS tables, and not leap updates through barycorrpy.