02/04/2026
GoogleのAI検索モードで、SOHO/LASCO C3の視野内に入るC/2026 A1をプロットするpythonコードを生成してもらいました。座標の反転やC3の視野円が倍になってたりの間違いはありましたが、ちゃんとJPL Horizonsのデータを引用して、まずますの結果が得られました。
彗星の名前と計算開始、終了時刻を入力するだけのシンプルなコードです。
表示が気に入らない点があるので手直しして、今後使えるようにしておこうかと思います。
更新1：近日点を追記し、ちょっと見やすくしました。
更新2：LASCOの光学系とイメージセンサーの仕様から計算した画角に変更、SOHOと太陽の距離から視直径を計算するように変更しました。コードが長くなるので、下記コードには反映させていません。
https://lasco-www.nrl.navy.mil/index.php?p=content/handbook/hndbk_6
更新3：細かい修正を AIにお願いすると、前の修正を忘れたり余計な修正をしてきたので、手作業でコードを変更して、好みのプロットにしました。下記コードには反映させていません。
---
import numpy as np
import matplotlib.pyplot as plt
import astropy.units as u
from astropy.coordinates import SkyCoord
from astropy.time import Time, TimeDelta
from astroquery.jplhorizons import Horizons
from sunpy.coordinates import frames
# 1. 設定
start_time = Time("2026-04-02T00:00:00")
end_time = Time("2026-04-07T00:00:00")
# 1時間ステップに変更
time_step_hours = 0.25
delta_t = TimeDelta(time_step_hours * u.h)
# 時間配列の作成
num_steps = int((end_time - start_time).jd * 24 / time_step_hours) + 1
times = start_time + np.arange(num_steps) * delta_t
comet_name = "C/2026 A1"
tx_list = []
ty_list = []
r_sun_list = []
print(f"{time_step_hours}時間おきに計算中（全{len(times)}点）: {comet_name}...")
for t in times:
# SOHOの座標取得
soho_query = Horizons(id='SOHO', location='', epochs=t.jd)
soho_table = soho_query.vectors()
soho_coord = SkyCoord(x=soho_table['x'], y=soho_table['y'], z=soho_table['z'],
unit='au', representation_type='cartesian',
frame='heliocentrictrueecliptic', obstime=t)
# 彗星の座標取得
comet_query = Horizons(id=comet_name, location='', epochs=t.jd)
comet_table = comet_query.vectors()
comet_coord = SkyCoord(x=comet_table['x'], y=comet_table['y'], z=comet_table['z'],
unit='au', representation_type='cartesian',
frame='heliocentrictrueecliptic', obstime=t)
r_au = np.sqrt(comet_table['x']**2 + comet_table['y']**2 + comet_table['z']**2).item()
r_sun_list.append(r_au)
# SOHOから見た太陽投影座標(HPC)に変換
hpc_frame = frames.Helioprojective(observer=soho_coord, obstime=t)
comet_hpc = comet_coord.transform_to(hpc_frame)
# 配列から数値を取り出してリストに追加 (.item()で確実にスカラー化)
tx_list.append(comet_hpc.Tx.to('deg').value.item())
ty_list.append(comet_hpc.Ty.to('deg').value.item())
# 2. グラフの作成
fig, ax = plt.subplots(figsize=(10, 10))
# LASCO C3 視野円 (約15.8度)
c3_fov = plt.Circle((0, 0), 15.8/2, color='blue', fill=False, linestyle='--', alpha=0.5, label='LASCO C3 FOV (15.8°)')
ax.add_patch(c3_fov)
# 太陽の視直径 (約0.25度半径)
sun_circle = plt.Circle((0, 0), 0.25, color='orange', fill=True, label='Sun')
ax.add_patch(sun_circle)
# 彗星の軌跡 (1時間ごとの細かな線、マーカーは小さく)
ax.plot(tx_list, ty_list, 'r-', linewidth=1.5, alpha=0.8, label=f'Path ({time_step_hours}h step)')
# 12時間おきに小さなドットを表示して時間経過を分かりやすくする
ax.scatter(tx_list[::int(12/time_step_hours)], ty_list[::int(12/time_step_hours)], color='red', s=10, zorder=3)
# 日付ラベルの追加 (24時間 = 24ステップおきに表示)
for i, t in enumerate(times):
if i % int(24/time_step_hours) == 0:
ax.annotate(t.strftime('%m-%d %H:%M'), (tx_list[i], ty_list[i]),
textcoords="offset points", xytext=(5/2,5/2), fontsize=8)
# グラフ設定
ax.set_aspect('equal')
ax.set_xlim(-20/2, 20/2) # East Left
ax.set_ylim(-20/2, 20/2)
ax.set_xlabel('Solar X (HPC Tx) [deg]')
ax.set_ylabel('Solar Y (HPC Ty) [deg]')
ax.set_title(f'SOHO/LASCO C3 View: {comet_name}')
ax.grid(True, linestyle=':', alpha=0.4)
ax.legend(loc='upper right')
# 日心距離が最小のインデックスを取得
peri_idx = np.argmin(r_sun_list)
peri_time = times[peri_idx]
peri_dist = r_sun_list[peri_idx]
peri_tx = tx_list[peri_idx]
peri_ty = ty_list[peri_idx]
# 結果の表示
print("\n=== 近日点（Perihelion）の結果 ===")
print(f"時刻 (UTC): {peri_time.iso}")
print(f"近日点距離: {peri_dist:.6f} AU")
print(f"SOHOから見た位置 (HPC): Tx={peri_tx:.4f}, Ty={peri_ty:.4f}")
# グラフ上に近日点地点を強調表示
ax.scatter(peri_tx, peri_ty, color='magenta', s=120, edgecolors='black',
marker='P', zorder=6, label='Perihelion')
ax.annotate(f'Perihelion: {peri_dist:.3f} AU\n{peri_time.strftime("%y-%m-%d %H:%M UT")}',
(peri_tx, peri_ty), xytext=(40, 20), textcoords="offset points",
arrowprops=dict(arrowstyle="->", color='magenta'),
fontsize=9, color='darkmagenta', fontweight='bold')
# 凡例を更新
ax.legend(loc='upper right')
plt.show()
