skyfielders / python-skyfield

Elegant astronomy for Python
MIT License
1.38k stars 208 forks source link

Oppositions/Conjuctions precision issue #876

Closed orchard0 closed 12 months ago

orchard0 commented 12 months ago

Consider the following code:

from datetime import datetime, timedelta, timezone
from skyfield.api import load
from skyfield import almanac

ts = load.timescale()

def findConjuctions(t0, t1):
    f = almanac.oppositions_conjunctions(eph, eph['moon'])
    ta, y = almanac.find_discrete(t0, t1, f)
    print(*ta)

eph = load('de441_part-2.bsp')

t = datetime(2023, 6, 1, tzinfo=timezone.utc)

for delta in [31, 60, 100]:
    t0 = ts.from_datetime(t)
    t1 = ts.from_datetime(t + timedelta(days=delta))
    findConjuctions(t0, t1) 

Results:

<Time tt=2460099.6547848135> <Time tt=2460113.6932674637>
<Time tt=2460099.654784815> <Time tt=2460113.6932674656> <Time tt=2460128.986006175> <Time tt=2460143.272909314>
<Time tt=2460099.654784815> <Time tt=2460113.6932674646> <Time tt=2460128.986006178> <Time tt=2460143.2729093144> <Time tt=2460158.2727949396> <Time tt=2460172.902319756> <Time tt=2460187.56721168>

The expected outcome for the opposition and conjunctions found is to have exactly the same value. But tt=2460099.6547848135 obtained when delta is 31 is not exactly the same as tt=2460099.654784815 when delta is 60 or 100. Or tt=2460113.6932674637 is not exactly the same as tt=2460113.6932674656 or tt=2460113.6932674646.

This becomes an issue if the results are store in a database and need to be retrieved; they can't be found because the values aren't exactly the same. Why does this happen and how to fix it?

brandon-rhodes commented 12 months ago

@orchard0 — Unfortunately, from any given starting point, all a computer can do to find a conjunction is to jump closer and closer to the moment of conjunction, until it gets ‘close enough’ where ‘close’ is defined by the user and what they are trying to accomplish. So two attempts to find the same conjunction, if they use even slightly different starting points, are going to land two different places within the realm of ‘close enough’ answers. If we compare your first two times, and convert to seconds:

(2460099.6547848135 - 2460099.654784815) * 3600 * 24

—then we see that they are within between one and two tenths of a millisecond of each other:

-0.00012069941

You can look at the EPSILON value and the epsilon keyword arguments in skyfield/searchlib.py if you want to see where find_discrete() gets its default idea of ‘close enough’. You can supply your own value for epsilon= when calling it if you want to use one of your own.

So the answer to your question is that, alas, you can't reliably give names to conjunctions that always work. But you can get close:

  1. What if you round each conjunction to the nearest second, and call it the conjunction of 2023-01-02-03:04:05? Then most of the time, two separate searches for it will still land in the same minute, and get the same name. Only once out of every thousand cases — because the search routine goes to within a millisecond of the answer — will a conjunction's millisecond-wide ‘close enough’ boundary happen to straddle the border between one second and the next, and cause it to get two different names from two different searches.
  2. If you name the conjunction by its minute, like 2023-01-02-03:04, then only one out of every 60,000 conjunctions will happen to have its millisecond fall across a minute boundary and give it two names.
  3. What if you just name each conjunction by its date, like 2023-01-02? Then only one out of every 86,400,000 conjunctions would have its millisecond-of-close-enough happen to fall across the boundary between two names, and alias the conjunction to two different dates.

Maybe try out the last option and see how it goes. Better yet, actively have your code check whether the conjunction falls within a millisecond of midnight, and if it does, then go ahead and proactively enroll it under two different names in the database — the previous date and the following date — and you'll be guaranteed to find it when you look it up later!

I'm going to close this issue as Skyfield is behaving as expected, but feel free to ask more questions by making further comments here even though it's closed!