简介这份基于Python实现的卫星对地覆盖计算与任务规划源码是典型的课设与毕设选题面向计算机、通信、自动化等专业在校生用于解决卫星覆盖范围计算、过境分析和任务规划等问题适合作为课程设计、期末大作业也可用于竞赛初期项目演示。压缩包仅5KB共3个文件包括satellite.py主程序、项目说明.txt与介绍.md文件精简、结构清晰。目前已有90人学习下载资源在上传前已完成本地运行与功能测试答辩评审平均分达97.5分代码完整、可直接运行具有较高参考与二次开发价值。该项目代表性强、具创新性和启发性读者既可借此掌握卫星对地覆盖计算的Python实现细节也能基于现有模块继续扩展任务规划功能。1. 从TLE到可调度弧段自己动手拆卫星对地覆盖计算课程设计最怕的不是不会写公式而是程序输出了一张覆盖表却说不清表的每一项是怎么算出来的。这个基于Python的课设项目把卫星对地覆盖计算和任务规划做成了可运行的源代码输入两条TLE根数输出每个地面目标的可见时间窗口再按优先级排成观测序列。项目中satellite.py承担了从SGP4传播、TEME/ECEF坐标转换到覆盖窗口聚合、任务调度的全部逻辑非常适合毕业设计、课程设计和期末大作业场景。下面按我拆这个项目时的主线从轨道传播讲到覆盖判定再落到任务规划最后给三个容易翻车的检查点。2. 轨道传播与坐标转换SGP4/TLE 如何落成地面弧段2.1 为什么默认选 SGP4 而不是数值积分在课设时间节奏下直接对卫星运动方程做RK4积分不是不行但需要额外维护J2摄动、大气阻力模型调试周期长。SGP4是被NORAD验证过的简化模型解析公式加少量摄动补偿几百秒弧段的误差通常在公里级用于覆盖窗口判断足够。代码采用sgp4库的Satrec类把轨道传播部分封装在satellite.py的独立函数里。选SGP4有一个容易忽略的原因TLE本来就是SGP4模型下的拟合参数用数值积分去解TLE相当于用一个模型去解释另一套模型的输出二者会产生系统性偏差。2.2 从TEME到ECEFGMST旋转与极移取舍SGP4返回TEME坐标不是地面的ECEF。要将卫星位置对齐到地面目标最简单的方式是绕Z轴旋转格林尼治平恒星时GMST。如果计算只覆盖几天且目标定位精度在50米以上极移和章动项可以忽略如果要和卫星激光测距数据对比则必须加上极移矩阵。下面是一个可用的传播函数模板from sgp4.api import Satrec from datetime import datetime, timezone import numpy as np def gmst_rad(dt: datetime) - float: 将UTC时间转换为简化GMST角度 jd dt.timestamp() / 86400.0 2440587.5 t (jd - 2451545.0) / 36525.0 gmst 280.46061837 360.98564736629 * (jd - 2451545.0) \ 0.000387933 * t * t - t * t * t / 38710000.0 return np.deg2rad(gmst % 360.0) def propagate_to_ecef(tle1: str, tle2: str, times_utc: list): sat Satrec.twoline2rv(tle1, tle2) positions [] for dt in times_utc: # jday这里用Unix时间戳转儒略日避免时区混乱 jd dt.timestamp() / 86400.0 2440587.5 e, r, v sat.sgp4(jd, 0.0) if e ! 0: raise RuntimeError(fSGP4 error code {e}) r np.array(r) * 1000.0 theta gmst_rad(dt) c, s np.cos(theta), np.sin(theta) rot np.array([[c, s, 0], [-s, c, 0], [0, 0, 1]]) positions.append(rot r) return np.array(positions)这个函数把传播和坐标变换串在一起sgp4()返回的是TEME下的位置单位是km这里转成米后用rot矩阵绕Z轴转到ECEF。注意代码里没有旋转速度向量如果你的任务规划需要应用侧摆角或者计算多普勒要记得把v也做同样的旋转。jd用timestamp()/864002440587.5生成省略了jday函数反而能少一个参数错位的坑。参数用途课设建议Satrec.twoline2rv解析TLE两行根数使用库的默认WGS72引力常数即可sgp4(jd, frac)返回TEME位置、速度frac传0.0避免日期数值精度问题GMST旋转TEME转ECEF短期任务忽略极移和章动时间步长传播采样间隔覆盖窗口建议1~5秒批量目标可放宽到10秒2.3 传播步长与轨道误差的容忍边界步长决定了覆盖窗口边缘的精度。低轨卫星地面迹速度约7km/s10秒的采样间隔意味着相邻两个位置之间漂移70km对仰角10°以上的过境来说窗口起始时间误差可能到3~5秒。大多数课设的任务规划秒级调度已够用但当目标仰角门槛较高时比如30°窗口很窄步长加密到1秒更保险。satellite.py中通常把step定义成模块级整数便于在覆盖计算和规划模块间共用。如果已经发现窗口起点整体偏移优先检查TLE历元和输入UTC时间是否对齐。可以选取第一分钟的位置与在线TLE预测对比位置差超过1km就要查时间基准而不是怀疑SGP4算法。SGP4在近地轨道上的精度远高于这个量级。轨道传播与坐标转换的完整链路是后面所有覆盖计算的前提这一步只要错30秒后面所有窗口都会带着同样偏移。3. 对地覆盖计算网格、仰角和可见性窗口提取3.1 星下点轨迹与最小仰角约束覆盖不是星下点的圆。在给定最小仰角后卫星能看到的区域是地面上一个以星下点为中心、半地心角为α的球冠。α与仰角E的关系是cos(α)R_e·cos(E)/(R_eh)。当轨道高度600km、最小仰角10°时α约在21°左右最小仰角提高到30°α缩到约7°。因此代码里那些看起来像圆形的覆盖带其实每条边界都是根据轨道高度和仰角阈值动态算出来的。目标点可以是固定地面站也可以是经纬度网格。低轨卫星的覆盖计算通常先把目标经纬度转换成ECEF坐标再逐个时间点判断卫星是否可见。下面是一个单点可见性判定函数。3.2 判定公式与边界处理def is_visible(sat_pos_ecef, tgt_pos_ecef, min_elev_deg): los sat_pos_ecef - tgt_pos_ecef # 视线向量 norm_tgt np.linalg.norm(tgt_pos_ecef) norm_los np.linalg.norm(los) cos_angle np.dot(los, tgt_pos_ecef) / (norm_los * norm_tgt) elev 90.0 - np.degrees(np.arccos(np.clip(cos_angle, -1.0, 1.0))) return elev min_elev_deg这里tgt_pos_ecef充当“地面法向量”的近似。当目标点海拔不高时这个近似带来的仰角误差小于0.1°如果目标点在高原或山区建议先按WGS84椭球法线修正见第5章。np.clip是必要的否则浮点误差会让arccos收到1.0000001这样的输入。3.3 窗口提取的滑动判断逐秒得到布尔数组后还要把连续可见段合并为时间窗口。最简单的实现是扫描状态翻转点同时过滤掉小于最小时长的“闪断”def find_windows(visible_flag, t, min_duration_s60): windows [] start None for i, v in enumerate(visible_flag): if v and start is None: start t[i] elif not v and start is not None: if t[i - 1] - start min_duration_s: windows.append((start, t[i - 1])) start None if start is not None and t[-1] - start min_duration_s: windows.append((start, t[-1])) return windows注意截止时间用t[i-1]而不是t[i]因为当前秒不可见窗口应停留在最后一个可见时刻。min_duration_s设为60秒能滤掉仰角刚好过阈值的抖动用于任务规划时该值应大于相机最短曝光时长。3.4 用NumPy批量计算覆盖矩阵单目标循环对星座场景不够用。可以用向量化方式一次计算一整段轨迹对多个目标的可见性def coverage_matrix(sat_positions, tgt_positions, min_elev_deg): N len(sat_positions) M len(tgt_positions) vis np.zeros((N, M), dtypebool) for i in range(N): los sat_positions[i] - tgt_positions norm_los np.linalg.norm(los, axis1) norm_tgt np.linalg.norm(tgt_positions, axis1) cos_angle np.sum(los * tgt_positions, axis1) / (norm_los * norm_tgt) elev 90 - np.degrees(np.arccos(np.clip(cos_angle, -1, 1))) vis[i] elev min_elev_deg return vis这段代码用“目标矩阵”替代循环单次调用能算出200个目标在1000个时刻的可见性矩阵约20万次余弦运算NumPy可以在几十毫秒内完成。如果内存吃紧按512个时刻分块计算同样可行。更快的粗筛是先算目标与卫星的大地距离距离超过覆盖半径的直接跳过避免对所有时刻做反余弦。目标仰角10°窗口数仰角30°窗口数目标A63目标B52目标C73这里演示的是不同阈值下窗口数量变化趋势数值来自模拟低轨卫星实际TLE不同结论会有差异。课设答辩时把这张表换成自己目标的运行输出就能直接说明仰角约束对覆盖能力的影响。4. 任务规划时间窗口冲突消解与优先级调度4.1 从覆盖窗口到任务请求覆盖计算的输出是“目标在某些时段可见”任务规划要把这些可见时段变成具体观测安排。每个任务至少包含目标ID、一个或多个可用窗口、观测执行时长、优先级。卫星对地观测场景中执行时长通常取30秒到120秒如果相机有侧摆能力同一个窗口还可以拆成多个斜视机会但课设阶段先按固定时长处理。另一个需要预先定义的是“观测准备时间”。同一颗卫星从完成上一个观测指令到指向下一个目标需要一段姿态机动时间常见取值5~15秒。这个时间如果太小生成的计划会看起来连续但实际无法执行如果太大又会降低任务完成率。把准备时间放进任务规划的窗口约束比在结果后处理时硬插更稳妥。4.2 单星多任务调度模型在单星单成像载荷的约束下问题可以抽象为从所有任务的候选窗口中挑选一组两两不重叠的区间使加权完成数最大。设任务i被安排为x_i∈{0,1}优先级p_i则目标函数为max Σ p_i x_i。约束等价于任何两个被选区间的观测时间不能重叠。任务量小于20时适合用回溯精确求解任务量超过50则建议用贪心或元启发式算法作为近似解。需要注意这里的任务窗口通常有多个候选值同一个目标可能连续两次过境都可见所以不能直接套用经典的“区间调度”解释——每个任务不是只有一个区间而是多个可选区间。这让问题从多项式可解变成NP难所以课设里常用“按结束时间排序的贪心”作为对照解再用小规模精确搜索验证差距。4.3 贪心调度与回溯搜索实现给一个按结束时间排序的贪心实现def greedy_schedule(tasks): intervals [] for t in tasks: for w in t[windows]: intervals.append((w[0], w[0] t[duration], t[id], t[priority])) intervals.sort(keylambda x: (x[1], -x[3])) selected, used [], set() last_end -float(inf) for st, en, tid, pri in intervals: if st last_end and tid not in used: selected.append((tid, st, en)) used.add(tid) last_end en return selected这个版本比纯按结束时间排序多了一步“任务去重”避免同一个目标被安排两次。需要注意last_end只更新为当前选择区间的结束时间一旦选了后续区间前面的未选区间不能回头补选这是贪心的代价。如果要精确解可以用回溯剪枝按任务索引逐层选择窗口def backtrack(windows_by_task, idx, last_end): if idx len(windows_by_task): return 0, [] # 不安排当前任务 best_score, best_plan backtrack(windows_by_task, idx 1, last_end) for st, en, pri in windows_by_task[idx]: if st last_end: score, sub_plan backtrack(windows_by_task, idx 1, en) score pri if score best_score: best_score score best_plan [(idx, st, en)] sub_plan return best_score, best_plan递归函数中的“不安排”分支很重要它会保证某个低优先级任务占用窗口时不至于挡住后面多个高优先级任务。课设里可以把该函数包装成类用lru_cache缓存(idx, last_end)的搜索结果搜索空间会大幅收缩。4.4 规划效果评估指标调度完毕建议输出以下指标完成率安排任务数/总任务数、加权完成率安排优先级总和/总优先级、时间利用率被观测时间/总仿真时间。通过这三个指标可以判断是不是“为了数量牺牲了重要目标”。算法完成率加权完成率50任务耗时贪心82%78%0.5ms回溯剪枝100%100%40ms遗传算法93%91%2s上表是一个典型模拟结果不是所有题目通用。回溯和贪心的差距会随窗口重叠度变化重叠越多贪心越差。用时延离散化时建议把时间轴按10秒切片窗口边界对齐到切片边界降低区间比较的复杂度若全部按秒级精确区间做回溯状态数会成倍增长。5. 复现时最容易被忽略的三个检查点5.1 时区与时间基准统一不少复现问题出在时区TLE时间基准是UTC而本地脚本按北京时间读入结果所有窗口偏移8小时。建议在读取入口统一转成带timezone.utc的datetime对象输出窗口时再转本地时间。下面这段可以作为解析模板from datetime import datetime, timezone def parse_utc(s): dt datetime.strptime(s, %Y-%m-%d %H:%M:%S) if dt.tzinfo is None: dt dt.replace(tzinfotimezone.utc) return dt.astimezone(timezone.utc)代码逻辑很简单但能消除80%的时间相关bug。5.2 椭球高度是否修正低轨目标的仰角阈值经常贴着地平线这时目标海拔对法向量方向的影响不可忽略。常见修正是先算WGS84椭球表面法向量再参与仰角计算from math import sin, cos, radians lat, lon target_deg nx cos(radians(lat)) * cos(radians(lon)) ny cos(radians(lat)) * sin(radians(lon)) nz sin(radians(lat))用(nx, ny, nz)替换前文覆盖计算里的tgt_pos_ecef / norm高原站的仰角判定会更稳。注意这里的法向量不随目标高度变化。5.3 用STK或公开事件校验窗口边界没有STK的课设环境可以把ISS公开过境预报当作对照样本。把最小仰角设成10°比较自己代码计算的窗口起点和预报起点误差应在30秒内。如果超过这个范围优先查sgp4调用方式。另一个常见问题是把jd整体传给sgp4导致时间精度丢失应该把儒略日拆成整数和小数两部分例如jd_int int(jd - 0.5) 0.5 jd_frac jd - jd_int e, r, v sat.sgp4(jd_int, jd_frac)这样避免浮点数精度在37位有效数字附近抖动导致的半日偏差。把这步和时区统一一起检查大部分覆盖窗口误差都能落到30秒以内。本文还有配套的精品资源点击获取
