Недействительные результаты со скриптом, который использует Skyfield

Я изучаю возможности великолепного программного обеспечения Skyfield от Брэндона Родса. Я сделал скрипт для вычисления союзов в прямом восхождении между случайными объектами. Я использую следующий скрипт:

from skyfield import almanac
from skyfield.searchlib import find_maxima, find_minima, find_discrete
from skyfield.api import Star, load
from datetime import datetime, date,timedelta
import pytz

planets = load('de430t.bsp')
earth = planets['earth']

x = [
['Aldebaran',[4, 35, 55.2],[16, 30, 33]],
['Regulus',[10, 8, 22.3],[11, 58, 2]],
['Pollux',[7, 45, 18.9],[28, 1, 34]],
['Antares',[16, 29, 24.4],[-26, 25, 55]],
]

ts = load.timescale(builtin=True)
t = ts.now()
tzn = 'Europe/Amsterdam'

tz = pytz.timezone(tzn)
now = datetime(2020, 1, 1, 0, 0, 0)
t0 = ts.utc(tz.localize(now))
t1 = ts.utc(tz.localize(now) + timedelta(days=+365))

def difference(t):
    e = earth.at(t)
    ra11, dec11, distance = e.observe(object).radec()
    ra12, dec12, distance2 = e.observe(planets[target1]).radec()
    diff = ra11.hours - ra12.hours
    return diff >= 0

difference.rough_period = 1.0

for count in range (len(x)):
    object = Star(ra_hours=(x[count][1][0], x[count][1][1], x[count][1][2]),dec_degrees=(x[count][2][0], x[count][2][1], x[count][2][2]))
    target1 = 'venus'

    t, b = find_discrete(t0, t1, difference)
    if len(t) > 0:
        print (f"{x[count][0]} and {target1}")

        for i, (a, b) in enumerate(zip(t, b)):
            e = earth.at(a)
            ra11, dec11, distance = e.observe(object).radec('date')
            ra12, dec12, distance2 = e.observe(planets[target1]).radec('date')
            print (f"Diff: {ra11.hours - ra12.hours}, ra_{x[count][0]}: {ra11.hours}, ra_{target1}: {ra12.hours}")
            print(f"{a.utc_iso()},{dec11._degrees - dec12._degrees}")
            print ("")

Я считаю, что это вычисляет момент времени, когда оба объекта имеют одинаковое RA.

К сожалению, я получаю такие результаты:

Aldebaran and venus
Diff: 4.600705365706519, ra_Aldebaran: 4.6176105612629375, ra_venus: 0.016905195556418524
2020-02-07T21:20:06Z,16.962942031863825

Diff: -0.0014205156698450239, ra_Aldebaran: 4.617748605529588, ra_venus: 4.619169121199433
2020-04-17T20:25:49Z,-10.136618155596008

Diff: -0.000670093218860579, ra_Aldebaran: 4.617892691238783, ra_venus: 4.618562784457644
2020-06-08T07:56:08Z,-4.921187478324768

Diff: -0.0001286749609095139, ra_Aldebaran: 4.618000948409605, ra_venus: 4.618129623370515
2020-07-12T06:44:16Z,-0.962286810883981

Regulus and venus
Diff: 10.140247344361857, ra_Regulus: 10.157152539918275, ra_venus: 0.016905195556418524
2020-02-07T21:20:06Z,12.28436748615726

Diff: 5.852858068422506e-06, ra_Regulus: 10.157702949500333, ra_venus: 10.157697096642265
2020-10-02T23:42:45Z,0.0903429562716429

Pollux and venus
Diff: 7.758742204719277, ra_Pollux: 7.7756474002756955, ra_venus: 0.016905195556418524
2020-02-07T21:20:06Z,28.391501492522927

Diff: 0.001226225278400328, ra_Pollux: 7.776229220287632, ra_venus: 7.775002995009232
2020-09-01T16:39:22Z,8.682000412217121

Antares and venus
Diff: 16.493491164600684, ra_Antares: 16.510396360157102, ra_venus: 0.016905195556418524
2020-02-07T21:20:06Z,-26.059118330110437

Diff: 0.000832014094154232, ra_Antares: 16.51126040187071, ra_venus: 16.510428387776557
2020-12-23T00:34:39Z,-5.652225571050259

Строка, начинающаяся с Diff, - это строка для контроля правильности вывода. Diff - рассчитанная разница в RA. Оно должно быть близко к нулю. Два других значения - это прямое восхождение обоих объектов. Они должны быть очень похожи. Вторая строка - это результат, который я хочу, это расчетное время и расстояние в градусах между объектами. Как вы можете видеть, по какой-то причине для каждого набора объектов я получаю недопустимый результат для экземпляра time: 2020-02-07T21: 20: 06Z, и значение разницы для этого экземпляра определенно не близко к нулю. Если я изменю объект Венеру на Луну, станет хуже, потому что каждый второй результат недействителен. Я проверил другие результаты по программе Skychart / Cartes du Ciel и по этим расчетам.

Я не могу понять, что здесь не так. Может кто-то помочь мне, пожалуйста?


person wie5Ooma    schedule 09.07.2020    source источник


Ответы (1)


Хороший вопрос! Я должен добавить новый раздел в https://rhodesmill.org/skyfield/searches.html объяснение этого обычного поведения, наблюдаемого при вычитании двух долгот или прямого восхождения. Ключ к разгадке тайны - посмотреть, что происходит с разницей углов в один из моментов, который проявляется в вашем выводе как фантомное соединение. Я прикрепил сценарий, который распечатывает это для самого первого события, которое вы распечатываете, между Венерой и Альдебараном:

2020-02-07 Difference: -19.33880215224481 Venus RA: 23.93746881891146
2020-02-08 Difference: 4.5908654129343995 Venus RA: 0.007801253732248779

Разница углов колеблется между -19,3 и 4,6, что сразу должно показаться нам подозрительным, поскольку это просто два разных названия для одного и того же угла! Вы можете подтвердить это, добавив 24,0 к -19,3, и вы получите угол, очень близкий к 4,6 (плюс-минус немного фактического движения, которое Венера совершила за один день).

Почему результат будет прыгать между двумя псевдонимами при точно такой же угловой разнице в небе?

Ответ заключается во втором факте, приведенном выше: прямое восхождение Венеры. Разрыв происходит точно в тот момент, когда Венера пересекает 0h прямого восхождения! Несмотря на то, что 23,93746881891146 и 0,007801253732248779 - это почти одинаковые углы, они отличаются на 24,0, потому что они охватывают то место в небе, где мы меняем название прямого восхождения.

Мой сценарий ниже также отображает сюжет, проясняющий ситуацию:

Участок прямых восхождений

На верхнем графике видно, что именно в тот момент, когда Венера сбрасывает свое RA обратно на ноль, разница RA совершает скачок на 24,0 часа от одного псевдонима к другому за те же четыре с половиной часа. разница в РА.

Решение?

Угловые различия должны быть ограничены диапазоном вроде [-12h, + 12h), чтобы принудительно выбрать один предпочтительный псевдоним для каждого возможного углового значения. В Python, как вы можете видеть в приведенном ниже скрипте, типичный маневр:

(ra1.hours - ra2.hours + 12.0) % 24.0 - 12.0

Это показано на втором графике выше как «Улучшенная разница» и не только правильно скрывает разрыв 7 февраля, больше не делая его похожим на событие, но и теперь правильно определяет противостояние между Венерой и Альдебараном ближе к концу. 2020 года (на правом краю графика), которая раньше казалась просто разницей, превышающей -12,0, но теперь сияет как решающий момент для угловой разницы.

Наконец, этот скрипт проверяет наличие возражений и отфильтровывает их из результатов поиска. Вы также найдете несколько возможных настроек вашего кода Python, которые вы можете рассмотреть, продолжая писать код на Python. Вот оно:

from skyfield.searchlib import find_maxima, find_minima, find_discrete
from skyfield.api import Star, load
from datetime import datetime, date,timedelta
import pytz

planets = load('de430t.bsp')
earth = planets['earth']

stars = [
    ['Aldebaran', (4, 35, 55.2), (16, 30, 33)],
    ['Regulus', (10, 8, 22.3), (11, 58, 2)],
    ['Pollux', (7, 45, 18.9), (28, 1, 34)],
    ['Antares', (16, 29, 24.4), (-26, 25, 55)],
]

ts = load.timescale(builtin=True)
t0 = ts.utc(2020, 1, 1)
t1 = ts.utc(2021, 1, 1)

# Exploring the first bad result for Venus and Aldebaran.

star = Star(ra_hours=stars[0][1], dec_degrees=stars[0][2])

for utc in (2020, 2, 7), (2020, 2, 8):
    t = ts.utc(*utc)
    ra1, _, _ = planets['earth'].at(t).observe(star).radec()
    ra2, _, _ = planets['earth'].at(t).observe(planets['venus']).radec()
    print(t.utc_strftime('%Y-%m-%d'),
          'Difference:', ra1.hours - ra2.hours,
          'Venus RA:', ra2.hours)

# Plot showing how to protect an angular difference against discontinuity.

import matplotlib.pyplot as plt

t = ts.utc(2020, 1, range(365))
e = planets['earth'].at(t)
star = Star(ra_hours=stars[0][1], dec_degrees=stars[0][2])
ra1, _, _, = e.observe(star).radec()
ra2, _, _, = e.observe(planets['venus']).radec()

fig, (ax, ax2) = plt.subplots(2, 1)
ax.plot(t.J, ra2.hours, label='Venus RA')
ax.plot(t.J, ra1.hours - ra2.hours, label='RA difference')
ax.set(xlabel='Year', ylabel='Hours')
ax.grid()
ax.legend()
ax2.plot(t.J, ra2.hours, label='Venus RA')
ax2.plot(t.J, (ra1.hours - ra2.hours + 12.0) % 24.0 - 12.0,
         label='Improved difference')
ax2.set(xlabel='Year', ylabel='Hours')
ax2.grid()
ax2.legend()
fig.savefig('tmp.png')

# So we need to force the difference into the range [-12 hours .. +12 hours]

def difference(t):
    e = earth.at(t)
    ra11, dec11, distance = e.observe(object).radec()
    ra12, dec12, distance2 = e.observe(planets[target1]).radec()
    diff = (ra11.hours - ra12.hours + 12.0) % 24.0 - 12.0
    return diff >= 0

difference.rough_period = 1.0

for name, ra_hms, dec_dms in stars:
    object = Star(ra_hours=ra_hms, dec_degrees=dec_dms)
    target1 = 'venus'

    t, b = find_discrete(t0, t1, difference)
    if len(t) == 0:
        break

    print (f"{name} and {target1}")

    for a, b in zip(t, b):
        e = earth.at(a)
        ra1, dec1, _ = e.observe(object).radec('date')
        ra2, dec2, _ = e.observe(planets[target1]).radec('date')
        if abs(ra1.hours - ra2.hours) > 6.0:
            continue  # ignore oppositions
        print(f"Diff: {ra1.hours - ra2.hours:.4f}, ra_{name}: {ra1.hours}, ra_{target1}: {ra2.hours}")
        print(f"{a.utc_iso()},{dec1._degrees - dec2._degrees}")
        print()

Теперь напечатаны следующие события:

Aldebaran and venus
Diff: -0.0014, ra_Aldebaran: 4.617748605529591, ra_venus: 4.619169121320681
2020-04-17T20:25:49Z,-10.136618155920797

Diff: -0.0007, ra_Aldebaran: 4.617892691238804, ra_venus: 4.618562784294269
2020-06-08T07:56:08Z,-4.921187477025711

Diff: -0.0001, ra_Aldebaran: 4.618000948409604, ra_venus: 4.618129623302464
2020-07-12T06:44:16Z,-0.9622868107603999

Regulus and venus
Diff: 0.0000, ra_Regulus: 10.157702949500333, ra_venus: 10.157697096607553
2020-10-02T23:42:45Z,0.09034295610918264

Pollux and venus
Diff: 0.0012, ra_Pollux: 7.776229220287636, ra_venus: 7.775002995218732
2020-09-01T16:39:22Z,8.68200041253418

Antares and venus
Diff: 0.0008, ra_Antares: 16.511260401870718, ra_venus: 16.51042838802181
2020-12-23T00:34:39Z,-5.652225570418828

Что, я считаю, решает и устраняет вашу проблему!

person Brandon Rhodes    schedule 10.07.2020
comment
Спасибо за уделенное время. Это отличное объяснение, которое действительно решает проблему. - person wie5Ooma; 11.07.2020