无人机运动适应波束成形与人工噪声反窃听仿真解析

无人机通信、运动适应光束控制、人工噪声反窃听这三个词叠在一起基本就能判断这不是那种跑通就算完事的波束成形 demo。这期源码14912 期我完整跑过一遍核心思路很明确无人机在高速移动时波束方向按运动状态动态修正始终对准合法接收者同时在发射信号里注入按零空间投影的人工噪声把窃听者压在低信噪比区间。整套方案在 Matlab 里用蒙特卡洛仿真验证它的用法不只是复现一条安全速率曲线更是一套可以自由扩展参数的物理层安全研究框架。适合正在做无人机通信、物理层安全、波束成形相关课题的研究生也适合需要在应急通信、野外作业或活动保障场景中验证链路安全性的工程师。下面我把建模、代码、调参和几个坑一次说清楚。1. 问题拆解移动节点的安全传输到底难在哪1.1 空-地信道为什么“藏不住”信号地面蜂窝通信里多径、遮挡、小区规划都能帮你把信号控制在某个空间范围。无人机通信则完全是另一套逻辑。无人机飞得高和地面用户之间绝大多数时间存在视距LOS路径而且空地信道往往由一根主径主导信道能量非常集中。对合法接收者是好事增益高但对窃听者也同样是好事只要它出现在波束覆盖范围内或者干脆就在旁瓣方向截获信号的难度比地面通信低得多。物理层安全在这个场景里比加密更值得优先考虑因为加密依赖上层协议和密钥分发无人机节点算力有限、链路时延敏感把安全做进物理层信号设计里代价更小、响应更快。1.2 移动性毁掉的不只是链路增益你肯定会想波束成形不就能把能量集中在合法用户方向吗对但别忘了无人机是动的。无人机一边飞自己和地面合法用户之间的方位角、俯仰角都在持续变化。如果波束权向量还是用起飞时刻的角度算出来的过几秒钟方向就跑偏了。波束方向图的主瓣增益会明显下降甚至原来指向合法用户的方向掉到旁瓣或零陷里。这时候不只是合法用户接收变差窃听者如果凑巧在原来主瓣方向附近反而会拿到更高增益。这个坑比单纯信道衰减更隐蔽系统看起来在正常发射功率也没降但安全速率已经悄悄崩了。1.3 人工噪声为什么是性价比最高的反窃听手段反制窃听者可以走两条路一是让合法链路的信号质量尽可能高二是让窃听链路的信干噪比SINR尽可能低。人工噪声就是第二条路的经典实现发射端在数据信号之外额外发射一路已知的伪随机噪声信号通过波束成形设计让这路噪声落到合法用户信道的正交方向上合法用户几乎不受影响但窃听者只要不在同一个方向就会收到一波额外的干扰。它的最大优势是不需要知道窃听者的信道状态信息CSI也不知道它的位置只要合法信道和窃听信道有差异就能有效果。在无人机这种窃听者可能随时移动、位置不确定的场景里这个特性非常实用。2. 系统建模与关键数学定义2.1 场景布局与坐标约定仿真场景我建议按最常见的“发射无人机—地面用户—地面/空中窃听者”来搭。发射端无人机标记为 T合法用户为 Bob窃听者为 Eve。坐标上无人机初始位置设为 (0, 0, H)H 取 100 m 左右Bob 固定在地面比如 (200, 0, 0)Eve 可以设成地面窃听终端比如 (180, 60, 0)也可以放在空中模拟另一架无人机。做蒙特卡洛仿真时Bob 位置可以固定Eve 每次随机撒点最后统计平均安全速率。这样既能看出人工噪声的稳定压制效果又能模拟窃听者位置不确定的真实场景。2.2 运动模型与方向角更新无人机按匀速直线运动建模t 时刻位置为p_t p0 v_uav * t;速度矢量 v_uav 是水平面上的二维向量。设无人机朝向 x 正方向飞行则从无人机到 Bob 的方位角和俯仰角分别用 atan2 和 atan 计算% Bob_pos 固定在地面 dx Bob_pos(1) - p_t(1); dy Bob_pos(2) - p_t(2); dz Bob_pos(3) - p_t(3); % 近似 -H theta_az atan2(dy, dx); % 方位角 theta_el atan(dz / sqrt(dx^2 dy^2)); % 俯仰角负值表示向下看问题里有个很关键的点角度不是线性变化的。无人机匀速直线飞过用户上空时方向和距离都在变特别是飞得比较近时角度变化率很大飞得很远时角度变化相对平缓。这就是“运动适应光束控制”要解决的时变特性。2.3 信道模型与阵列导向矢量空地链路用 LOS 主导信道建模。接收端收到的信号相位由传播距离决定阵列上每个阵元因为位置不同相位也不同。对水平放置的均匀平面阵UPA阵元坐标取 (mxd, myd, 0)其中 mx 0...Nx-1my 0...Ny-1d λ/2。以阵列中心为参考点某方向 (θ_az, θ_el) 的单位方向向量分量为u_x cos(theta_el) * cos(theta_az); u_y cos(theta_el) * sin(theta_az);第 (mx, my) 个阵元相对参考点的波程差为 mxdu_x mydu_y相位为 2π/λ 乘以这个波程差。导向矢量元素写成phase exp(-1j * 2 * pi / lambda * (mx * d * u_x my * d * u_y));这里的负号表示入射波到达不同阵元的相位延迟方向图翻转问题多半就出在这个符号上。合法信道向量写为 h β * a(θ_az_b, θ_el_b)窃听信道向量写为 g β_e * a(θ_az_e, θ_el_e)。β 为复用路径损耗的小尺度系数在做归一化仿真时先取 β 1重点观察波束增益和人工噪声的影响。2.4 安全速率的完整表达式发射端总功率为 P功率分配因子 ρ 表示数据信号占比则数据功率 Ps ρP人工噪声功率 Pj (1-ρ)P。合法用户接收 SINR 为γ_b Ps * |h^H w|^2 / (Pj * |h^H z|^2 σ_b^2)因为 z 设计在 h 的零空间内|h^H z|^2 理论上为零所以 γ_b ≈ ρP * |h^H w|^2 / σ_b^2。窃听者接收 SINR 为γ_e Ps * |g^H w|^2 / (Pj * |g^H z|^2 σ_e^2)这里 |g^H z|^2 通常不为零人工噪声就体现在这个分母里。安全速率定义为Rs max(0, log2(1 γ_b) - log2(1 γ_e))当窃听信道太好时Rs 可能为负物理上表示此时无法保证安全传输取 0。这个公式是整个仿真代码的输出指标后面调参全部围绕提高 Rs 展开。3. 运动适应光束控制Matlab 实现详解3.1 先写一个可复用的 UPA 导向矢量函数不要在主循环里现算相位建议单独写一个函数方便后面频繁调用function a_vec upa_steering(theta_az, theta_el, Nx, Ny, d, lambda) mx (0:Nx-1).; my (0:Ny-1).; u_x cos(theta_el) * cos(theta_az); u_y cos(theta_el) * sin(theta_az); phase exp(-1j * 2 * pi / lambda * (mx * d * u_x my * d * u_y)); a_vec phase(:); a_vec a_vec / norm(a_vec); % 归一化 end注意这里用 phase(:) 把 Nx x Ny 矩阵拉成 N x 1 列向量。波束成形权向量通常取 w a(θ_az, θ_el)这样主瓣方向就是合法用户方向阵列增益达到 N。上面这段把 u_x 和 u_y 分开算是因为 UPA 同时控制方位角和俯仰角比 ULA 更适合空地三维场景。3.2 运动状态更新与角度预测“运动适应”的核心是实现每个更新周期内的方向预估。一个常用的做法是已知当前时刻无人机位置、速度矢量和更新周期 T_slot直接预测下一时刻位置再用预测位置计算角度。这样可以抵消一部分因处理、计算带来的延迟p_next p_uav v_uav * T_slot; dx Bob_pos(1) - p_next(1); dy Bob_pos(2) - p_next(2); dz Bob_pos(3) - p_next(3); theta_az atan2(dy, dx); theta_el atan(dz / sqrt(dx^2 dy^2));如果无人机有姿态变化或速度变化可以在状态方程里加加速度项甚至用卡尔曼滤波做轨迹预测。对于匀速直线飞行的基准场景用线性外推就够了。3.3 波束权向量更新与固定波束对照每次更新周期内执行a_b upa_steering(theta_az_b, theta_el_b, Nx, Ny, d, lambda); w a_b; % 对准合法用户的运动适应波束为了做对比实验固定波束方案要保存初始时刻角度算出的 w0w_fixed w0; % 整个仿真过程中不更新这里有个细节对比时必须保证两种方案的总发射功率相同否则结果没有意义。我习惯的做法是先归一化 w 和 z再按 ρ 做功率分配。运动适应方案不是每时每刻都重新算波束而是每隔 T_slot 更新一次。更新周期本身就是个可调参数完全可以和速度、距离一起做敏感性分析。3.4 如何判断波束是否真的“对准”了跑仿真的时候除了安全速率我建议额外输出一个诊断指标合法方向上的归一化波束增益 |h^H w|^2 / N。这个值在理想对准时为 1角度偏了之后会掉到 0.5 甚至更低。实测下来10 m/s 速度、200 m 距离的典型场景里如果 1 秒才更新一次波束增益损失可能超过 3 dB安全速率损失更大。加了运动适应后增益能稳定在 0.95 以上。这个指标比单纯看安全速率更直观能帮你快速定位问题是出在波束对准还是人工噪声分配。4. 人工噪声反窃听核心代码与功率分配4.1 零空间投影复数共轭转置是最大的坑人工噪声要设计在合法信道 h 的零空间里也就是要找一个 z使得 h^H z 0。在 Matlab 里很多人图省事写 null(h.)结果仿真里人工噪声没起到“只干扰窃听者”的作用合法用户的 SINR 也被拉低还找不到原因。问题出在 h. 是普通转置不是共轭转置。复数导向矢量下必须用 h也就是 ctranspose。正确的是Z_null null(h); % 返回 N x (N-1) 矩阵满足 h * Z_null 0原因是 null(A) 返回满足 A*z 0 的基而我们要的是 h^H z h * z 0所以 A 必须取 h。这个小细节直接决定人工噪声到底能不能“精准绕过”合法用户。如果你用 h.相位反了噪声会漏进合法信道系统性能整体恶化而且不仔细看根本发现不了。4.2 生成人工噪声向量并合成发射信号拿到零空间基后生成随机噪声向量并投影即可v_noise (randn(size(Z_null, 2), 1) 1j * randn(size(Z_null, 2), 1)) / sqrt(2); z_an Z_null * v_noise; z_an z_an / norm(z_an); % 归一化 % 发射信号 x sqrt(rho * P_total) * w * s sqrt((1 - rho) * P_total) * z_an;这里 s 是归一化数据符号。合法的接收 SINR 由 |h^H w|^2 决定噪声项因为 h * z_an 0 而不出现窃听者方向 g 和 h 不一致时|g^H z_an|^2 给出了干扰功率。功率分配因子 ρ 不要乱取我实测在 0.5 到 0.8 之间有较好表现后面第 5 节细说。4.3 这里有个前提人工噪声不是万能的人工噪声的压制效果建立在合法信道和窃听信道不完全一致这个前提上。如果窃听者恰好站在合法用户旁边两个信道的导向矢量高度相关|g^H z_an|^2 同样趋近于零人工噪声对窃听者几乎没有影响。这时只能靠波束本身的空间分辨率区分两者阵元数越大波束越窄能区分的最小角度差越小。所以人工噪声不是替代阵列设计而是和波束成形配合使用。4.4 主仿真循环框架整套蒙特卡洛循环可以按下面骨架搭MC_count 500; num_slots round(sim_T / T_slot); Rs_matrix zeros(MC_count, num_slots); for mc 1:MC_count % 重新随机化 Eve 位置视实验目的可选 for k 1:num_slots p_uav p0 v_uav * (k * T_slot); [theta_az_b, theta_el_b] calc_direction(p_uav, Bob_pos); [theta_az_e, theta_el_e] calc_direction(p_uav, Eve_pos); h beta_b * upa_steering(theta_az_b, theta_el_b, Nx, Ny, d, lambda); g beta_e * upa_steering(theta_az_e, theta_el_e, Nx, Ny, d, lambda); w upa_steering(theta_az_b, theta_el_b, Nx, Ny, d, lambda); % 运动适应 Z_null null(h); v_noise (randn(size(Z_null,2),1) 1j*randn(size(Z_null,2),1)) / sqrt(2); z_an Z_null * v_noise; z_an z_an / norm(z_an); gamma_b rho * P_total * abs(h * w)^2 / sigma2_b; gamma_e rho * P_total * abs(g * w)^2 / ... ((1 - rho) * P_total * abs(g * z_an)^2 sigma2_e); Rs max(0, log2(1 gamma_b) - log2(1 gamma_e)); Rs_matrix(mc, k) Rs; end end Rs_avg mean(Rs_matrix, 1);注意每个蒙特卡洛循环里窃听者位置要不要重新撒点取决于你研究的是“平均窃听位置下的安全速率”还是“固定窃听者下的安全速率”。两个问题答案不同别混在一起出图。5. 结果分析与参数敏感性5.1 固定波束 vs 运动适应差距比想象中大固定波束方案在无人机起飞初期还能维持不错的安全速率随着无人机向前飞波束方向逐渐偏离 Bob合法用户方向的增益从主瓣掉向旁瓣。在 4x4 UPA、速度 15 m/s、初始距离 200 m 的设定下5 秒后固定波束的安全速率可能下降 30% 以上而运动适应方案基本能维持平稳。这背后的原因不只是增益损失还有人工噪声的泄漏特性变化波束指向偏了之后合法方向的等效 SINR 降低而窃听者方向可能恰好处在波束旁瓣抬升区两个因素叠加安全速率掉得比单个指标还快。5.2 功率分配因子 ρ 存在平台区把 ρ 从 0.1 扫到 0.9安全速率曲线大致呈“先升后降”的单峰形态。ρ 太小数据功率不足合法用户 SINR 上不去安全速率自然低ρ 太大人工噪声功率不足窃听者的 SINR 不会被压住安全速率也低。平台区通常出现在 0.5 到 0.7 之间具体位置受路径损耗、噪声基底和窃听者位置影响。窃听者离合法用户越近人工噪声效果越弱最优 ρ 会右移因为你更需要把功率给数据信号。跑参数扫描时建议对每个 ρ 都做蒙特卡洛平均单次信道实现的噪声太大看不出规律。5.3 速度、距离和更新周期怎么匹配无人机运动速度直接影响波束方向角的变化率。角速度近似为 v / RR 是无人机与用户水平距离。比如速度 20 m/s、距离 200 m方向角变化率约 0.1 rad/s也就是每秒约 5.7°。如果阵元数较少、波束本来就很宽更新周期可以长一点阵元数多、波束窄就必须加密更新。经验公式上建议保证每个更新周期内角度变化不超过 3 dB 波束宽度的三分之一到四分之一。4x4 阵列的 3 dB 波束宽度大约在十几度量级所以 T_slot 取 0.05 s 到 0.1 s 是充分安全的。表不同参数下的稳定性和更新周期建议场景速度 (m/s)近距离 R (m)建议 T_slot (s)固定波束损失情况慢速巡检51500.2轻微常规巡航152000.1明显高速跟踪251000.02~0.05严重高空盘旋103000.2轻微5.4 窃听者位置远近和角度都要看窃听者如果离合法用户很远但恰好在波束主瓣方向人工噪声仍会因为它“和合法信道方向接近”而失效。真正决定安全速率的不是距离而是窃听信道的导向矢量与合法信道导向矢量的相关性。仿真时要画一张二维热力图以 Bob 位置为原点x 轴是窃听者的水平 x 偏移y 轴是 y 偏移颜色是安全速率。你会发现安全速率在 Bob 周围存在一个“暗区”暗区大小由波束宽度决定越靠近 Bob人工噪声和波束都分不开安全速率趋近于 0。这张图对理解系统的物理极限很有帮助。6. 踩坑记录与调试速查表6.1 高频问题速查跑这套源码的过程中我梳理了几个出现频率最高的问题直接列成表供你排查现象可能原因排查/解决方法人工噪声毫无作用用了 h. 而不是 h 做零空间改用 h检查 h * z_an 是否为 0安全速率持续为 0功率分配 ρ 太小或窃听者接近 Bob增大 ρ重新设置 Eve 位置曲线抖动剧烈蒙特卡洛次数不足至少 500 次取平均固定随机种子 rng(2024)波束增益始终上不去导向矢量相位符号写反检查 exp 的负号用诊断指标固定波束收益反而高飞行距离太近、波束太宽拉大初始距离或增加阵元数缩小波束角度计算结果异常atan2 参数顺序写错atan2(dy, dx)不是 atan2(dx, dy)6.2 几个值得养成的调试习惯第一个习惯把合法方向的归一化波束增益打出来别的先不管。如果增益没有保持在 0.95 以上波束控制这块肯定有问题人工噪声调得再好也没有意义。第二个习惯设置 rng 固定种子先跑一条单次信道实现用 debug 模式逐行看 h * z_an 是不是接近机器精度级别的小数再上蒙特卡洛。第三个习惯把噪声功率从线性单位换算成 dB脚本里统一线性计算只在出图时转 dB别混着算否则很容易出现 10^7 量级的数字溢出问题。6.3 说实话这套方案最值钱的地方在哪我最深的体会是运动适应光束控制和人工噪声不是两个孤立模块它们的配合逻辑才是整套方案的核心。波束控制负责“让合法用户始终在照亮区域里”人工噪声负责“让合法用户之外的方向尽量黑”。单独用波束控制窃听者只要在旁瓣附近还是能截获单独用人工噪声合法用户也可能被误伤。只有把两者放在同一个发射信号结构里才能实现安全速率的稳定提升。这个思路你后续完全可以扩展到多无人机协作、IRS 反射面辅助等场景只需要在现有框架上改信道建模和波束优化部分。最后再补一个我实测的建议如果你想把这个仿真作为论文或毕设的基础不要只贴一张安全速率曲线。把第 5 节的敏感性分析做全特别是 ρ 和窃听者位置的热力图加上与固定波束的对比审稿人看到这种系统性对比基本不会再问“为什么不用更简单的固定波束”。仿真代码本身不难难的是把每个参数背后的物理意义讲清楚。这一点做到了这套源码的价值就真正体现出来了。