Skyfieldを使った地上駅パススケジュールの作成

ローカルパススケジューラが必要な理由

SATnogs ステーション、アマチュア無線シャックトラッキング AO-91/ISS/NOAA、または初期段階のキューブサット地上セグメントを運営している場合は、自動パス予測が必要です。商用ツール (STK、NOVA) はやり過ぎです。無料のウェブプレディクター (AMSAT、Heavens-Above) はうまくスクリプト化されません。STKクラウドが2026年3月に日没する中、最もクリーンで無料の代替手段は スカイフィールド です。これは Pyephem の後継となる最新の Python です。

この記事では、約80行のPythonで完全な地上局スケジューラについて説明します。新しいTLEを取得し、衛星リストの次の24時間のパスを計算し、最小高度でフィルタリングし、ダウンリンク頻度のドップラーを計算し、パスがカレンダーアプリに表示されるようにiCalendarを送信します。

インストール

pip install skyfield requests icalendar

それでおしまいです。スカイフィールドは初回実行時に独自のエフェメライド (DE421、IERS) を取得します。約 17 MB がローカルにキャッシュされます。

ステップ 1: 現在の TLE を取得する

TLE ドリフト — エポックから 7 日以内に使用してください。セレストラックから現在のエレメントセットを取り出してください:

from skyfield.api import load

stations_url = 'https://celestrak.org/NORAD/elements/gp.php?GROUP=amateur&FORMAT=tle'
sats = load.tle_file(stations_url)
by_name = {s.name: s for s in sats}
print(f'Loaded {len(sats)} amateur-radio satellites')

ターゲットを名前で選んでください (Celestrak ファイルでは大文字と小文字が区別されます):

targets = [
    by_name['AO-91 (FOX-1B)'],
    by_name['ISS (ZARYA)'],
    by_name['NOAA 19'],
    by_name['GreenCube (IO-117)'],
]

## ステップ 2: 地上局の定義

from skyfield.api import wgs84, load

ts = load.timescale()
my_station = wgs84.latlon(40.0150, -105.2705, elevation_m=1624)  # Boulder, CO
min_elevation_deg = 10  # horizon mask to skip low passes

elevation_mパラメータを使用してください。低い仰角では傾斜範囲が大きく影響します。

ステップ 3: 上昇/最高潮/セットイベントを検索するSatellite.find_events()はイベントのフラットリスト (0 = 上昇、1 = クライミング、2 = セット) を返します。これらを 3 つに分けて処理します。

from datetime import timedelta

t0 = ts.now()
t1 = ts.from_datetime(t0.utc_datetime() + timedelta(hours=24))

for sat in targets:
    t, events = sat.find_events(my_station, t0, t1,
                                altitude_degrees=min_elevation_deg)
    # Group into rise/peak/set triples
    triples = zip(t[0::3], t[1::3], t[2::3])
    for rise, culm, setting in triples:
        # Peak elevation
        topo = (sat - my_station).at(culm)
        alt, az, dist = topo.altaz()
        print(f'{sat.name:30s}  {rise.utc_iso():20s}  peak {alt.degrees:4.1f}°  slant {dist.km:4.0f} km')

t0で既に地平線上にあるサテライトは無視してください (イベントリストは最高潮とセットのペアで始まります)。

ステップ 4: パス上でのドップラーシフト

アマチュアのUHF (435 MHz) では、パスによってキャリアが±10 kHzずれる可能性があります。レシーバーはトラッキングする必要があります。レンジレートから瞬時ドップラーを計算:

import numpy as np

FREQ_HZ = 435_800_000  # AO-91 downlink
C = 299_792_458.0      # m/s

def doppler_shift(sat, observer, t):
    """Positive when satellite approaching (received frequency is higher)."""
    # Numerical range rate from 1s apart
    dt = 1.0  # seconds
    t_next = ts.from_datetime(t.utc_datetime() + timedelta(seconds=dt))
    r0 = (sat - observer).at(t).distance().m
    r1 = (sat - observer).at(t_next).distance().m
    range_rate = (r1 - r0) / dt  # m/s, positive = moving away
    return -FREQ_HZ * range_rate / C  # Hz, positive = approaching

# Sample across a pass
for pct in [0, 25, 50, 75, 100]:
    t_sample = ts.from_datetime(
        rise.utc_datetime() + (setting.utc_datetime() - rise.utc_datetime()) * pct / 100
    )
    d = doppler_shift(sat, my_station, t_sample)
    print(f'  {pct:3d}% pass:  doppler = {d:+8.0f} Hz')

一般的な 6 分間の AO-91 オーバーヘッドパスでは、ドップラーが約 +10 kHz (上昇) から 0 (天頂)、-10 kHz (セット) まで流れるのがわかります。これをラジオのメイン VFO でトラッキングすれば、信号が途切れることはありません。

ステップ 5: 各パスのリンクバジェット

傾斜範囲がわかったら、それを rftools 自由空間パス損失計算ツール に入力するか、インラインで計算します。

def fspl_db(distance_m, freq_hz):
    return 20 * np.log10(4 * np.pi * distance_m * freq_hz / C)

# At zenith of an AO-91 pass (~800 km slant at 40° elevation from Boulder)
fspl = fspl_db(800_000, FREQ_HZ)
print(f'FSPL = {fspl:.1f} dB')  # ~143 dB

ITU-R 伝搬モデルを使用してモンテカルロリンクバジェットを完成させるには、衛星リンクバジェットアナライザー を使用し、シナリオ URL をステーションログブックと共有してください。

ステップ 6: iCalendar を発行して、パスがカレンダーに表示されるようにする

これがスケジューラが人間にとって役立つところです。

from icalendar import Calendar, Event
from datetime import datetime
import pytz

cal = Calendar()
cal.add('prodid', '-//rftools ground station//')
cal.add('version', '2.0')

for sat in targets:
    t, events = sat.find_events(my_station, t0, t1, altitude_degrees=10)
    triples = zip(t[0::3], t[1::3], t[2::3])
    for rise, culm, setting in triples:
        topo = (sat - my_station).at(culm)
        alt, az, _ = topo.altaz()
        ev = Event()
        ev.add('summary', f'{sat.name} (peak {alt.degrees:.0f}°)')
        ev.add('dtstart', rise.utc_datetime())
        ev.add('dtend', setting.utc_datetime())
        ev.add('description', f'Culminate at {culm.utc_iso()}, AZ {az.degrees:.0f}°')
        cal.add_component(ev)

with open('passes.ics', 'wb') as f:
    f.write(cal.to_ical())
print('Wrote passes.ics — import into Google Calendar / Apple Calendar / Outlook')

passes.icsを Apple カレンダー/Google カレンダー/Outlook にインポートすると、各パスの 5 分前にプッシュ通知が届きます。

ステップ 7: cron で実行してください~/bin/station-passes.pyにスクリプトをドロップして、1 日に 1 回 cron してください:

# Refresh pass schedule nightly at 03:00 local
0 3 * * * /usr/bin/python3 ~/bin/station-passes.py > ~/passes.ics 2>> ~/station.log

WebDAV カレンダーにcurl --upload-fileとペアリングすると、すべてのパスが自動的にスマートフォンに表示されます。

STKクラウドを1つずつ置き換える

アマチュア/キューブサット/小規模の商業事業では、統合ツールキットがSTK Cloudの日常的な機能に取って代わります。

ミッション・デザインの貿易研究や機関プログラムの場合は、やはりデスクトップSTKまたはNASA GMATが必要です。それ以外は、AMSATの設計レビューを含め、すべてSkyfield + rftoolsスタックが揃っています。

さらに読む

-スカイフィールド・ドキュメンテーション -セレストラック TLE カタログ -rftools サテライトパスプレディクター — 最も一般的なアマチュア衛星用のこのスケジューラーのホストバージョン -STK クラウドからの移行:無料の代替手段 — STK ワークフロー全体の置き換えに関する関連記事 -UHF CubeSat リンクバジェットのサイズ設定 — エンドツーエンドの作業例