【2025年认证杯数学建模挑战赛】A题
该题为典型的空间几何建模+轨道动力学建模+预测问题。
⚙ 问题一：利用多个天文台的同步观测，确定小行星与地球的相对距离
问题分析
- 已知若干地面天文台的观测数据：方位角 (Azimuth) 和 高度角 (Altitude)；
- 所有天文台在同一时间观测同一颗小行星；
- 目标是根据这些观测结果三角测量出小行星的空间坐标，进而得出相对于地球的位置（特别是距离）。
解题思路
- 获取每个天文台的地理坐标（经纬度、海拔）；
- 将每个天文台的观测方位角、高度角转换为单位观测方向向量（以地心为参考）；
- 将地面站的位置用地心坐标表示；
- 利用两台或多台天文台的观测方向和坐标，建立**多射线最近点法（最小距离点）**来确定目标的小行星的位置。
可能用到的数学模型与方法
- 地理坐标系 → 地心坐标系转换；
- 球面坐标系 → 笛卡尔单位向量；
- 多射线几何定位（三维空间中线与线之间的最短距离问题）；
- 最小二乘法求解最优交点。
初步代码框架（Python）
import numpy as npdef geo_to_ecef(lat, lon, alt):# WGS84参数a = 6378137.0 # 地球长半轴e = 8.1819190842622e-2 # 偏心率lat, lon = np.radians(lat), np.radians(lon)N = a / np.sqrt(1 - e**2 * np.sin(lat)**2)x = (N + alt) * np.cos(lat) * np.cos(lon)y = (N + alt) * np.cos(lat) * np.sin(lon)z = (N * (1 - e**2) + alt) * np.sin(lat)return np.array([x, y, z])def azel_to_direction(az, el, lat, lon):az, el = np.radians(az), np.radians(el)x = np.cos(el) * np.sin(az)y = np.cos(el) * np.cos(az)z = np.sin(el)# 本地坐标转为地心坐标（简略版）return np.array([x, y, z]) # 需进一步旋转变换def triangulate(positions, directions):# 最小二乘法解多条射线的最近点A, b = [], []for p, d in zip(positions, directions):d = d / np.linalg.norm(d)I = np.eye(3)A.append(I - np.outer(d, d))b.append((I - np.outer(d, d)) @ p)A, b = np.sum(A, axis=0), np.sum(b, axis=0)return np.linalg.lstsq(A, b, rcond=None)[0]# 示例：三台天文台
stations = [{'lat': 34.0, 'lon': -118.0, 'alt': 1000, 'az': 130, 'el': 45},{'lat': 40.0, 'lon': -75.0, 'alt': 500, 'az': 120, 'el': 40},{'lat': 35.0, 'lon': 139.0, 'alt': 20, 'az': 150, 'el': 50},
]positions, directions = [], []
for s in stations:pos = geo_to_ecef(s['lat'], s['lon'], s['alt'])dir = azel_to_direction(s['az'], s['el'], s['lat'], s['lon']) # 简化positions.append(pos)directions.append(dir)asteroid_pos = triangulate(positions, directions)
distance = np.linalg.norm(asteroid_pos) / 1.496e+11 # 转为AU单位
print("小行星距离地球中心：", distance, "AU")
⚙ 问题二：利用多次观测数据，预测小行星未来15或30天的轨道及天文台观测角
问题分析
- 输入：多台天文台对某颗小行星的多时刻观测数据（时刻、方位角、高度角）；
- 输出：未来15或30天内，每天某固定时刻的方位角和高度角预测值。
解题思路
- 将观测数据反演为轨道初始状态（即轨道六要素或位置-速度）；
- 根据轨道力学模型（两体问题/摄动模型），进行轨道积分，得到未来轨道；
- 对于预测日的每个时间点，将轨道位置转换成各天文台的视角坐标，计算方位角和高度角。
可能用到的数学模型与方法
- Gauss法或Laplace法进行轨道确定（初轨道估计）；
- 开普勒轨道传播（Kepler Propagation）；
- JPL HORIZONS、NASA SPICE工具可模拟（也可调用简化库如 Skyfield）；
- 天球坐标计算公式（向量到方位角、高度角转换）；
初步代码框架（Skyfield库）
from skyfield.api import load, Toposts = load.timescale()
ephemeris = load('de421.bsp')
earth = ephemeris['earth']# 假设已知小行星轨道初值（可从MPC或观测数据拟合）
from skyfield.keplerlib import KeplerOrbit
# 伪代码: 构造轨道
# orbit = KeplerOrbit(elements, epoch=ts.utc(2024, 1, 1))# 假设某台站
station = earth + Topos(latitude_degrees=34.0, longitude_degrees=-118.0, elevation_m=1000)# 预测未来30天
t = ts.utc(2024, 1, range(1, 31))
# astrometric = station.at(t).observe(orbit) # 实际需从轨道构造天体
# alt, az, d = astrometric.apparent().altaz()# print("每天12:00的高度角：", alt.degrees)
# print("每天12:00的方位角：", az.degrees)
✅ 总结建议
| 问题 | 本质 | 方法建议 | 数据需求 |
|---|---|---|---|
| 问题一 | 三维空间定位 | 三角测量 + 多射线最短距离 | 多台站的同时观测数据 |
| 问题二 | 轨道传播与预测 | 轨道反演 + 两体运动传播 + 坐标转换 | 多次观测 + 天文台坐标 |
后续会更新详细的解题思路、论文与代码，请看下方~