大阪市立科学館学芸員 江越航のホームページ  大阪市立科学館 学芸員

 江越 航のホームページ

Python で天文計算

日の出・日の入りの計算(下記)

Python で今日の月

Python で日食計算

Python で月食計算

Python で軌道計算


日の出・日の入りの計算 (Skyfield版)

(→ PyEphem を使った計算方法はこちら

Python の天文計算ライブラリ Skyfield を使用することで、簡単に日の出・日の入り時刻の計算ができます。

必要モジュール をインポートします。

from skyfield.api import load, Topos
from skyfield.almanac import find_discrete, sunrise_sunset
from datetime import timedelta
from pytz import timezone

Timescale オブジェクト、および天体暦のファイルを読み込みます。

ts = load.timescale()
eph = load('de421.bsp')

日の出・日の入りを調べる開始の時刻 t0 と終了の時刻 t1 を設定します。t0, t1 は世界時で設定します。
ここでは t0 に現在時刻、t1 に 1日後を設定しています。

t0 = ts.now()
t1 = ts.utc(t0.utc_datetime() + timedelta(days=1))

タイムゾーンを日本に設定します。

tz = timezone('Asia/Tokyo')

観測者の位置を指定する Topos オブジェクトを設定します。緯度、経度は文字列として設定します。

osaka = Topos('34.6914 N', '135.4917 E')

find_discrete(t0, t1, sunrise_sunset(eph, osaka)) で、t0, t1 の間に起こる日の出・日の入りの時刻を計算します。
結果は、t、updown に配列形式で代入されます。
t には日の出または日の入りの時刻が、updown には日の出の場合 True、日の入りの場合 False が、順番に配列形式で代入されます。

t, updown = find_discrete(t0, t1, sunrise_sunset(eph, osaka))

t、updown に代入された結果を、順番に表示します。

for yi, ti in zip(updown, t):
    print('次の日の出:' if yi else '次の日の入り:', ti.astimezone(tz).strftime('%Y-%m-%d %H:%M:%S'), 'JST')

全体では、次のようになります。

from skyfield.api import load, Topos
from skyfield.almanac import find_discrete, sunrise_sunset
from datetime import timedelta
from pytz import timezone

ts = load.timescale()
eph = load('de421.bsp')

t0 = ts.now()
t1 = ts.utc(t0.utc_datetime() + timedelta(days=1))
tz = timezone('Asia/Tokyo')

osaka = Topos('34.6914 N', '135.4917 E')

t, updown = find_discrete(t0, t1, sunrise_sunset(eph, osaka))

for yi, ti in zip(updown, t):
    print('次の日の出:' if yi else '次の日の入り:', ti.astimezone(tz).strftime('%Y-%m-%d %H:%M:%S'), 'JST')

表示結果

次の日の入り: 2020-01-22 17:16:14 JST
次の日の出: 2020-01-23 07:02:21 JST