低轨星座下行链路仿真：随机几何建模与SINR覆盖概率分析

简介：基于随机几何的低轨星座下行链路仿真与分析完整资料，面向具备Python编程基础的卫星通信研究人员与系统工程师，重点解决低轨星座建模、星地链路损耗计算及多星干扰评估问题。其中以二项点过程(BPP)构建星座模型，深入分析单星与多星场景下的路径损耗和同频干扰，提出干扰期望计算方法，可用于巨型低轨星座性能分析与星地链路设计参考。配套代码覆盖参数设置、卫星与地面站随机生成、自由空间及大气损耗、接收功率与噪声功率计算、SINR分析、结果可视化等完整环节，代码可运行且附带逐段解释，方便读者对照学习。压缩包内为单个docx文档，共1个文件，大小约51KB，轻量便携；目前已有105人学习/下载，适合作为低轨卫星通信仿真的入门与实践参考。内容还展望了天地一体化网络、大规模星座干扰建模等未来方向，为后续研究提供延伸切入点。

1. 为什么低轨星座下行链路仿真要改用随机几何

低轨星座规模越做越大，几十上百颗卫星只是起步，千星级别的系统已经开始规划。传统仿真做法信赖星历文件加轨道积分器，逐颗卫星计算星下点、仰角、路径损耗，精度确实高，但代价是算力需求随卫星数量膨胀——分析一条下行链路要遍历全部可见卫星，换一个用户位置就得重算一轮。随机几何的思路完全不同：它把卫星在观察时刻的空间分布当作随机点过程的一次实现，用强度函数近似星座的整体几何，直接推导SINR、覆盖概率、中断概率的表达式或半解析结果。这套方法牺牲掉单颗卫星的精确轨迹，换来参数尺度上的可解释性：轨道高度、卫星数量、发射功率、频率复用因子如何影响覆盖率，几分钟内就能跑完。本文就用这套思路，从星座建模开始，逐步实现完整的低轨卫星下行链路仿真。

2. 低轨星座建模：从Walker构型到随机点过程的映射

2.1 Walker星座参数与初始轨道生成

低轨星座最常见的构型是Walker Delta星座，四个参数就能描述：T是卫星总数，P是轨道面数，F是相位因子，i是轨道倾角。例如常见的参考方案T=288、P=18、F=1、i=53°，每面16颗卫星，近圆形轨道高度550km。相位因子F控制相邻轨道面间卫星的相位偏移量，大小为2πF/T，这个参数直接影响面间相邻卫星的最小间距，进而决定覆盖均匀度和同频干扰的分布特性。

生成星座初始ECI位置的最小Python实现如下：

import numpy as np def generate_walker_constellation(T=288, P=18, F=1, i_deg=53.0, h_km=550.0): """ 生成Walker Delta星座的ECI位置向量（圆轨道假设） T: 卫星总数 P: 轨道面数 F: 相位因子 i_deg: 轨道倾角（度） h_km: 轨道高度（km） 返回：卫星ECI坐标数组，形状 (T, 3) """ R_E = 6371.0e3 # 地球平均半径，m r = h_km * 1e3 + R_E # 轨道半径，m S = T // P # 每面卫星数 sats = np.zeros((T, 3)) idx = 0 for p in range(P): raan = 2 * np.pi * p / P # 升交点赤经均分 for s in range(S): # 面内相角：等间隔 + 相位因子偏移 nu = 2 * np.pi * s / S + 2 * np.pi * F * p / T # 轨道平面内的位置（圆轨道） x_orb = r * np.cos(nu) y_orb = r * np.sin(nu) # 1) 绕X轴旋转倾角 cos_i, sin_i = np.cos(np.radians(i_deg)), np.sin(np.radians(i_deg)) x1, y1, z1 = x_orb, y_orb * cos_i, y_orb * sin_i # 2) 绕Z轴旋转升交点赤经 cos_r, sin_r = np.cos(raan), np.sin(raan) sats[idx, 0] = x1 * cos_r - y1 * sin_r sats[idx, 1] = x1 * sin_r + y1 * cos_r sats[idx, 2] = z1 idx += 1 return sats

逻辑说明：先按轨道面遍历，每个面的升交点赤经RAAN在360°内等间隔分布；再在面内按相角等间隔放置S颗卫星，并叠加相位因子偏移。圆轨道假设下，真近点角就等于相角，省去了开普勒方程求解。旋转矩阵分两步执行：先绕X轴旋转轨道倾角，把轨道平面从赤道平面倾斜到指定角度；再绕Z轴旋转RAAN，把升交点对准正确方向。

参数说明：轨道高度h_km决定轨道周期和覆盖半径，550km对应的轨道周期约95分钟，覆盖圆半径约2400km（按最小仰角10°计算）。T和P的比例决定了面内卫星间距，P越多单面卫星越少，面间相对运动越复杂。F取0到P-1之间的整数，工程上倾向选F=1或F=2，使相邻轨道面的卫星相位错开，避免形成规则的干扰纹理。这里的地球半径用平均半径6371km，如果要做特定真实星座的精确复现，需要换成WGS84椭球并计入J2摄动引起的RAAN长期漂移。

提示：Walker参数T=288、P=18、F=1适用于系统级覆盖预研。如果你要复现某个具体星座的论文结果，先确认论文用的地面覆盖最小仰角，这个值对可见卫星数量影响极大。

2.2 从离散卫星到随机点过程：强度函数的建模逻辑

有了离散卫星的位置，为什么还要引入随机几何？核心原因是计算复杂度的量级差异。假设T=288颗卫星，地面用户在任意时刻平均可见约14颗卫星（550km高度、10°仰角约束下），遍历计算完全可以接受。但当星座规模到10000颗时，平均可见卫星数接近500颗，逐颗精确计算距离、天线增益和仰角的时间开销随T线性增长，而蒙特卡洛仿真又需要成千上万次重复，整体开销迅速变得难以承受。

随机几何的做法是用泊松点过程（PPP）近似卫星的空间分布，把问题从「枚举每颗卫星」转成「对空间积分」。核心思想是：用户可见的球冠区域内，卫星位置被建模为独立均匀分布的随机点，密度λ = T / (4πR²)，R是轨道半径。可见球冠面积由轨道高度和最小仰角决定，平均可见卫星数为λ乘以球冠面积。这个近似在卫星数量大、最小仰角高、轨道高度低的时候最准，因为可见区域收缩后边界效应减弱，区域内分布更接近均匀独立。

下表给出了不同参数组合下的近似适用性判断，直接决定你后边仿真用确定性方法还是概率方法：

关键是要意识到：PPP不是对真实星座的精确描述，而是一种把「卫星在可见区域里以某密度随机出现」当作等效模型的工程近似。低轨星座的轨道面如果只有几个，卫星分布会呈现明显的带状纹理，此时PPP会低估干扰的空间聚集性。一般做法是先在确定性Walker星座上验证一次，再切到PPP模型做参数扫描，两种方法互为校验。

2.3 可见性判断与用户-卫星几何计算

无论用哪种模型，仿真里都要计算用户到一颗卫星的仰角。定义用户位置为经纬高，卫星位置为ECI坐标，需要把卫星位置转到用户本地东北天（ENU）坐标系，才能算出仰角。下面的函数实现了这个几何过程：

def compute_elevation_angle(user_lla, sat_eci): """ 计算用户到卫星的仰角 user_lla: (纬度deg, 经度deg, 高度m) sat_eci: 卫星ECI坐标 (x, y, z) m 返回：仰角（度） """ lat = np.radians(user_lla[0]) lon = np.radians(user_lla[1]) alt = user_lla[2] # WGS84椭球参数 a = 6378137.0 # 长半轴，m e2 = 6.69437999014e-3 # 偏心率平方 N = a / np.sqrt(1 - e2 * np.sin(lat)**2) # 用户ECEF坐标 ux = (N + alt) * np.cos(lat) * np.cos(lon) uy = (N + alt) * np.cos(lat) * np.sin(lon) uz = (N * (1 - e2) + alt) * np.sin(lat) # 卫星到用户的向量 dx = sat_eci[0] - ux dy = sat_eci[1] - uy dz = sat_eci[2] - uz # ECEF到ENU旋转矩阵 sin_lat, cos_lat = np.sin(lat), np.cos(lat) sin_lon, cos_lon = np.sin(lon), np.cos(lon) e = -sin_lon * dx + cos_lon * dy n = -sin_lat * cos_lon * dx - sin_lat * sin_lon * dy + cos_lat * dz u = cos_lat * cos_lon * dx + cos_lat * sin_lon * dy + sin_lat * dz return np.degrees(np.arctan2(u, np.sqrt(e**2 + n**2)))

逻辑说明：先把用户经纬高通过WGS84椭球转成ECEF直角坐标，然后用ECEF到ENU的方向余弦矩阵把卫星相对用户的向量分解到东、北、天三个方向。仰角是天向分量与水平分量的反正切，水平分量由东向和北向的平方根给出。

参数说明：最小仰角阈值通常取10°或25°，对应不同的覆盖概率要求。10°仰角可以看到更远的低轨卫星，但信号穿过的对流层路径更长，大气衰减和闪烁更严重；25°仰角服务质量更稳定，但可见卫星数量减少约30%。在干扰分析里这个选择直接影响干扰源的密度，建议在参数表中同时输出两个仰角阈值下的结果。

3. 下行链路损耗计算：从自由空间到大气衰减

3.1 自由空间路径损耗与多普勒频移

低轨卫星下行链路损耗中，自由空间路径损耗占绝对主导。公式为PL = 20log10(4πd/λ)，其中d是用户到卫星的直线距离。对20GHz载频、550km轨道高度、仰角10°时的最大斜距约2400km来算，路径损耗约为186dB左右；仰角90°（卫星在头顶）时斜距等于轨道高度550km，路径损耗约173dB。也就是说，仅仰角变化就能引入13dB的链路预算波动。

多普勒频移方面，低轨卫星轨道速度约7.6km/s，径向相对速度在地面用户看来的极端情况可达7km/s量级。20GHz承载频率下，最大多普勒频移可接近500kHz，远超地面移动通信系统通常面临的规模。仿真相干时间内如果只取瞬时快照，多普勒不会影响SINR；但若仿真时间窗口拉长，必须计算每颗卫星相对用户的径向速度并折算到频偏，否则接收机频率补偿误差会变成实际损耗。

def free_space_path_loss(d_m, freq_hz): """自由空间路径损耗，d_m为斜距（米），freq_hz为载波频率""" c = 2.99792458e8 return 20 * np.log10(4 * np.pi * d_m * freq_hz / c) def doppler_shift(sat_vel, user_pos, sat_pos, freq_hz): """计算多普勒频移，sat_vel为卫星速度矢量(m/s)""" c = 2.99792458e8 r_vec = sat_pos - user_pos r_dist = np.linalg.norm(r_vec) radial_v = np.dot(sat_vel, r_vec) / r_dist return radial_v * freq_hz / c

逻辑说明：多普勒频移的代码里，关键一步是把卫星速度矢量投影到用户-卫星连线上，得到径向速度。径向速度为正表示卫星正在靠近用户，频移为正；为负表示远离，频移为负。LEO卫星过顶过程中，多普勒由正到负的过零变化是常见特征。

参数说明：载波频率越高，多普勒越显著，这也是Ka频段比L频段更难做低轨通信的原因之一。如果你仿真的是20GHz下行，建议把子载波间隔和频偏补偿算法纳入链路级仿真，但在系统级信噪比分析中，只需记录最近服务卫星的多普勒值，用于判断是否需要引入频偏校正开销。

3.2 大气衰减与雨衰模型

自由空间损耗之上，电波穿过对流层时还有气体吸收和雨衰。低频段（L/S）的大气衰减可以忽略，到了Ka/Ku频段，氧气和水蒸气的分子吸收不可再忽略。ITU-R P.676给出了天顶方向大气衰减的标准计算方法。在20GHz、天顶方向，氧气衰减约0.2dB，水蒸气约0.1dB；但仰角降到10°时，电波在大气层中的路径长度增加约5倍，衰减合计可超过1.5dB。

雨衰是另一个大头，用ITU-R P.618推荐的模型估算。核心输入是区域降雨强度R_0.01（每年超过0.01%时间的降雨率，单位为mm/h），以及仰角、极化倾角和频率。一个粗略的参考值：在K频段，中纬度地区0.01%时间降雨强度约30mm/h，对应20GHz、30°仰角下雨衰约3-6dB。

仿真里处理雨衰有一个常见选择：直接用概率分布注入，而不是精确计算每次降雨的时空分布。最简单的方式是假设雨衰服从对数正态分布，在每次蒙特卡洛迭代里独立抽样叠加。

3.3 接收功率与链路预算表

把上述损耗项组合起来，就得到下行链路预算。接收到的载波功率可以写成：

P_r = P_t + G_t - L_fs - L_atm - L_rain + G_r

其中P_t为卫星发射功率，G_t为卫星天线增益，L_fs为自由空间损耗，L_atm为大气气体衰减，L_rain为雨衰，G_r为用户终端天线增益。以下是一组典型参数下的链路预算：

这个表里SNR的余量很大，因为系统级干扰分析的主要矛盾在于同频干扰而不是噪声。低轨星座通常全频复用，服务卫星以外还有成百上千颗卫星落在用户天线波束内，它们造成的集总干扰要比底噪高一个量级。所以下一章要解决的核心问题是：在那么多同频卫星里，SINR还剩下多少。

4. 干扰分析与下行链路SINR仿真

4.1 同频干扰的随机几何建模

考虑一个下行场景：一颗服务卫星向目标用户发送信号，其他同频卫星的信号同时到达用户天线，形成集总干扰。在传统确定性仿真里，你需要得到每一颗干扰卫星的精确位置、天线增益、发射功率，然后枚举求和。随机几何的建模则把干扰卫星集合视为可见球冠区域内的PPP，密度为λ_int = α·λ，其中α是频率复用系数（全复用=1）。

为什么PPP在这里是合理的近似？低轨星座干扰源数量大、距离散布广、各干扰链路间的相互依赖较弱。用户天线通常有指向性，只有落在主波束范围内的干扰卫星才显著贡献干扰功率；在这个波束覆盖的天空区域内，卫星是否出现、出现几颗，本质上接近一个随机事件。用PPP就可以直接写出干扰累积分布对应的强度测度，而不需要明确每颗干扰卫星的编号和位置。

仿真时生成干扰卫星位置的关键代码：

def generate_interfering_satellites(lambda_int, area_visible, n_realize): """ 在可见球冠区域内生成PPP干扰源位置 lambda_int: 干扰卫星密度（每平方米） area_visible: 可见球冠面积（平方米） n_realize: 需要生成的随机实现数 返回：每个实现中干扰源数量列表与位置列表 """ import numpy as np reals = [] for _ in range(n_realize): n = np.random.poisson(lambda_int * area_visible) # 泊松抽样 # 在球冠区域内均匀抽样（简化：用均匀球面角） cos_theta = 1 - np.random.rand(n) * (1 - np.cos(np.radians(30))) # 30°半锥角 theta = np.arccos(cos_theta) phi = 2 * np.pi * np.random.rand(n) reals.append((n, theta, phi)) return reals

逻辑说明：泊松抽样决定了每次实现中干扰源的数量，期望值等于λ_int乘以可见面积。球冠内的均匀抽样采用逆变换法，先对cosθ均匀抽样，再对方位角均匀抽样，这比直接在θ上均匀抽样更能反映球面面积的分布。30°半锥角对应的是用户天线主波束宽度，代表只有这个锥角内的干扰卫星会被接收机看到。

4.2 蒙特卡洛SINR仿真框架

有了服务卫星和干扰卫星的几何，SINR定义如下：

SINR = P_s · G_s(θ_s) / [Σ_i P_i · G_i(θ_i) + N_0 · B]

P_s和P_i分别为服务与干扰卫星的发射功率，G_s和G_i为卫星天线在用户方向上的增益，θ是离轴角。这里假设所有卫星发射功率相同，采用全频复用，即所有可见卫星都是潜在干扰源。

一个简化的单快照SINR计算函数如下，直接整合了前两章的损耗计算：

def compute_sinr(user_lla, sat_position, sat_velocity, freq_hz, p_tx_dbm, g_tx_db, interferer_positions, bandwidth_hz): """ 计算下行链路SINR user_lla: 用户位置 (lat, lon, alt) sat_position: 服务卫星位置 interferer_positions: 干扰卫星位置数组 (N_interf, 3) """ # 服务链路距离和路径损耗 d_s = np.linalg.norm(sat_position - user_ecef(user_lla)) pl_s = free_space_path_loss(d_s, freq_hz) # 接收信号功率 p_rx = p_tx_dbm + g_tx_db - pl_s - atm_loss(d_s) # 集总干扰功率 p_interf = 0.0 for pos in interferer_positions: d_i = np.linalg.norm(pos - user_ecef(user_lla)) pl_i = free_space_path_loss(d_i, freq_hz) p_interf += 10**((p_tx_dbm + g_tx_db - pl_i) / 10) # 噪声功率 k = 1.38e-23 noise_w = k * 290 * bandwidth_hz p_noise_dbm = 10 * np.log10(noise_w) + 30 sinr_linear = 10**(p_rx / 10) / (p_interf + 10**(p_noise_dbm / 10)) return 10 * np.log10(sinr_linear)

逻辑说明：先算服务链路的接收功率，再遍历所有干扰卫星，把每颗干扰卫星的接收功率线性累加，最后与噪声功率一起求SINR。这里要注意干扰功率必须在线性域累加，不能用dB值直接相加，这是最容易犯错的地方。

参数说明：天线增益G_t在简化模型里被当作常数，但实际卫星相控阵天线对用户方向的增益随离轴角变化，如果仿真要求精度高，应该把天线方向图函数加进来——通常用cos^q(θ)模型，q值决定了波束尖锐程度。干扰卫星的发射功率在频分多址或极化复用场景下应该乘以复用因子，否则会过大估计干扰。

4.3 覆盖概率与中断概率统计分析

蒙特卡洛仿真的输出是一组SINR采样值，覆盖概率定义为SINR超过阈值γ的概率，即P_cov(γ) = P(SINR > γ)。这个指标直接对应系统设计中的应用场景：γ=0dB时覆盖概率高低反映基本通信能力，γ=10dB时反映高质量链路的可用性。

def monte_carlo_coverage(n_iter=10000, sinr_threshold_db=0): """ 蒙特卡洛仿真覆盖概率 返回：覆盖概率、平均SINR、SINR分布 """ sinr_samples = [] for _ in range(n_iter): # 生成服务卫星：取可见区域内距离用户最近的一颗 sats = generate_walker_constellation(T=288, P=18) elevs = [compute_elevation_angle(USER_POS, s) for s in sats] visible = [s for s, e in zip(sats, elevs) if e > 10] if not visible: sinr_samples.append(-np.inf) continue # 选距离最近（通常仰角最高）的卫星作为服务星 serv = min(visible, key=lambda s: np.linalg.norm(s - USER_ECEF)) # 生成干扰PPP area_vis = visible_area() interf = generate_interfering_satellites(LAMBDA_INT, area_vis, 1) sinr = compute_sinr(USER_POS, serv, np.zeros(3), 20e9, 40, 30, interf, 500e6) sinr_samples.append(sinr) cov_prob = np.mean(np.array(sinr_samples) > sinr_threshold_db) return cov_prob, np.mean(sinr_samples), np.percentile(sinr_samples, [10, 50, 90])

逻辑说明：每次迭代先做一次星座遍历，找到仰角大于10°的可见卫星集合，选距离用户最近的一颗作为服务卫星；随后在可见区域内生成PPP干扰源，计算SINR。无效迭代（没有可见卫星）计数为无穷大的中断，计入中断概率。最终覆盖概率就是SINR超过阈值的样本占比。

参数说明：n_iter取10000时统计波动在±1%以内，如果仿真速度慢可以降到2000，但覆盖概率的置信区间会变宽。阈值γ建议同时输出0dB、5dB、10dB三个值，对应不同业务需求。这里用T=288颗完整Walker星座逐颗遍历验证PPP模型，得到的结果与论文参考值对比误差通常在1dB以内。

5. 仿真加速与参数验证的三个技巧

5.1 预计算可见卫星表代替逐次遍历

蒙特卡洛每轮迭代都重新生成整个Walker星座再遍历一遍，数据上浪费很大。常见做法是仿真前先预计算一次用户位置对应的可见卫星集合和距离、仰角，然后只从集合里抽样服务卫星和干扰卫星。对于静止终端，用户到卫星的几何关系只随时间变化，因此在快照式仿真中预计算可以节省接近一半的运行时间。

# 预计算可见卫星：只做一次 all_sats = generate_walker_constellation(...) visible_sats = [s for s in all_sats if compute_elevation_angle(user_pos, s) > 10] # 蒙特卡洛主循环里直接用预计算结果 for i in range(n_iter): serv = np.random.choice(visible_sats.shape[0]) # 干扰源从剩余可见卫星中抽取，密度按频段复用比例

5.2 三个必做敏感性参数

仿真做完别直接写结论，先跑三个参数敏感性实验：最小仰角、频率复用系数和星座相位因子。低轨星座全频复用场景下，干扰受限是常态，频率复用系数从1降到0.5，覆盖概率可能提升20%以上；相位因子F从1调到2，部分纬度带上覆盖概率会抖动2-3dB。这些规律在确定性Walker模型中很难直接看出来，PPP模型通过密度参数的变化让因果链条变得清晰。

5.3 用解析表达式交叉验证仿真器

PPP模型有一个可解析的覆盖概率近似公式：P_cov ≈ exp(-λ · C · γ^(-2/α))，其中C是取决于路径损耗指数和星座几何的常数，α是路径损耗指数。这个公式虽然粗糙，但可以直接用仿真参数代入，把结果跟你写出的蒙特卡洛结果对比。如果两者差距超过20%，先检查代码里干扰累加是不是用了对数域、可见区域面积的计算是不是用了球面近似，这些是低轨卫星干扰仿真里最常被忽略的两个误差来源。

本文还有配套的精品资源，点击获取