简介相控阵波束形成是雷达、通信与声纳领域的核心技术通过阵元相位控制实现波束电子扫描无需机械转动天线。这份Matlab仿真资源面向通信工程、信号处理方向的学生和工程师适合从理论走向实践的快速上手实验。压缩包内含3个文件包括1个Matlab仿真脚本.m与2张波束方向图未加窗与加窗对比整体仅78KB轻量实用脚本结构清晰便于直接运行和二次扩展。已有135人学习参考。通过运行脚本可直观观察不同窗函数对主瓣宽度与副瓣电平的影响理解时间延迟、相位延迟等波束形成算法的实现思路同时串起阵列模型设计、仿真参数设置与结果可视化流程帮助读者掌握从建模仿真到性能评估的完整链路为后续数字波束形成DBF和抗干扰研究打下基础。1. 相控阵波束形成Matlab仿真从阵列模型到自适应权值一次讲透接手相控阵波束形成Matlab仿真这个方向多数人第一反应是“不就是画个方向图吗”。实际跑起来才发现真正卡住你的不是方向图公式而是阵列流型怎么构建、权值怎么算、信号模型怎么设、以及为什么仿真结果和理论曲线对不上。这篇文章直接回答四个问题波束形成在解决什么、Matlab里最小可跑通的仿真怎么写、MVDR这类自适应算法怎么调参、以及工程上最常见的翻车点在哪。适合正在做阵列信号处理课程设计、雷达通信仿真预研、或者刚接触自适应波束形成的开发者和研究人员。2. 相控阵基础阵列模型、信号模型与波束形成的核心逻辑2.1 均匀线阵的几何模型与阵列流型向量相控阵波束形成的本质是通过控制每个阵元的相位补偿让某个方向的电磁波在接收端同相叠加。最常见的仿真对象是均匀线阵ULA也就是N个阵元等间距排成一条直线。设阵元间距为d信号入射方向与阵法线的夹角为θ那么相邻阵元之间的波程差是d·sinθ对应的相位差是2π·d·sinθ/λ。在Matlab里构建阵列流型向量时我一般会用这样的方式% 均匀线阵阵列流型向量生成 N 16; % 阵元数 d_lambda 0.5; % 阵元间距与波长比通常取0.5避免栅瓣 theta_deg 30; % 期望信号方向单位度 theta theta_deg * pi / 180; n (0:N-1); % 阵元索引列向量 a_theta exp(1j * 2 * pi * d_lambda * n * sin(theta));这段代码中最关键的是n采用列向量这样exp运算会得到N×1的复数向量对应N个阵元的相位响应。d_lambda 0.5意味着阵元间距是半波长这是相控阵仿真最常用的配置——既能保证不出现栅瓣又有足够的阵列孔径。当你把theta_deg从-90扫到90度逐点计算a_theta并画出来就得到了静态方向图。2.2 窄带信号模型与快拍数据的构造方法波束形成的仿真不能只算方向图还要模拟接收到的数据快拍。窄带假设下第k个快拍、N个阵元的接收向量可以写成期望信号、干扰和噪声的线性叠加。构造仿真数据时我一般采用以下方式% 构造接收数据快拍矩阵 M 500; % 快拍数 theta_s 30 * pi / 180; % 期望信号方向 theta_i -20 * pi / 180; % 干扰方向 a_s exp(1j * 2 * pi * 0.5 * (0:N-1) * sin(theta_s)); a_i exp(1j * 2 * pi * 0.5 * (0:N-1) * sin(theta_i)); s exp(1j * randn(M, 1)); % 期望信号复包络随相位 j exp(1j * randn(M, 1)); % 干扰复包络 noise sqrt(0.01) * (randn(N, M) 1j * randn(N, M)) / sqrt(2); X a_s * s a_i * j noise; % N×M 接收数据矩阵注意a_s * s中a_s是N×1s是1×M外积得到N×M数据矩阵每一列代表一个快拍。噪声功率设为0.01对应信噪比约20dB信号功率为1时。如果你期望调高干扰强度可以把1j * randn(M, 1)乘上大于1的系数。这个数据构造是整个仿真链路中最容易被忽视的部分——很多人直接从网上下载方向图脚本却不知道方向图和自适应波束之间的差距就在于有没有这个数据矩阵。2.3 相控阵波束形成的工作逻辑从相移补偿到权矢量波束形成的核心运算是一个加权求和y wᴴ·X。其中w是N×1的复权矢量上标H表示共轭转置。常规相控阵的权矢量就是阵列流型本身w a_θ它的作用是把各阵元接收信号做相位补偿后相加形成对θ方向的增益。这个操作在物理上对应移相器网络在仿真里就是一个点乘加求和。权矢量的选择决定了波束形成的性能。静态权值只在指定方向形成固定波束而自适应权值会根据接收数据的统计特性自动调整在期望方向保持增益的同时在干扰方向形成零陷。从仿真角度看静态波束只需一行w a_s就能画方向图自适应波束则要估计协方差矩阵并求解优化问题。理解了这一步你就知道为什么仿真里还要写协方差估计和矩阵求逆的代码。3. 自适应波束形成核心MVDR准则与Matlab实现3.1 MVDR的数学原理最小方差无失真响应MVDR最小方差无失真响应是自适应波束形成最经典的准则目标是在期望方向增益恒定的约束下使输出功率最小化。数学上写成min wᴴ·R·w约束条件是wᴴ·a_s 1。这个优化问题的闭式解是 w R⁻¹·a_s / (a_sᴴ·R⁻¹·a_s)。这里R通常是对接收数据协方差矩阵的估计值。仿真中用样本协方差矩阵R_hat X·Xᴴ / M替代理论协方差。MVDR之所以在工程仿真里被广泛使用是因为它只需要知道期望信号方向不需要知道干扰方向就能自动在干扰处形成零陷。如果要更稳健的版本可以在R的对角线上加一个加载量这就是后面要讲的对角加载技术。3.2 Matlab代码实现协方差估计、权值求解与波束方向图下面给出一个完整的MVDR波束形成仿真代码这段代码可以直接运行并输出方向图对比%% MVDR波束形成器仿真 clear; clc; N 16; % 阵元数 M 500; % 快拍数 d_lambda 0.5; % 半波长间距 theta_s 30 * pi / 180; % 信号方向 theta_i -20 * pi / 180; % 干扰方向 % 阵列流型 a_s exp(1j * 2 * pi * d_lambda * (0:N-1) * sin(theta_s)); a_i exp(1j * 2 * pi * d_lambda * (0:N-1) * sin(theta_i)); % 接收数据 s exp(1j * randn(M, 1)); j exp(1j * randn(M, 1)) * 100; % 强干扰功率是信号的40dB noise sqrt(0.01) * (randn(N, M) 1j * randn(N, M)) / sqrt(2); X a_s * s a_i * j noise; % 样本协方差矩阵 R_hat (X * X) / M; % MVDR权值 w_mvdr (R_hat \ a_s) / (a_s * (R_hat \ a_s)); % 注意使用左除 \ 而不是 inv 求逆数值稳定性更好 % 计算方向图 theta_scan (-90:0.5:90) * pi / 180; pattern_mvdr zeros(size(theta_scan)); pattern_cbf zeros(size(theta_scan)); for idx 1:length(theta_scan) a_scan exp(1j * 2 * pi * d_lambda * (0:N-1) * sin(theta_scan(idx))); pattern_mvdr(idx) abs(w_mvdr * a_scan); pattern_cbf(idx) abs(a_s * a_scan) / N; % 常规波束形成 end % 归一化方向图并绘图 pattern_mvdr pattern_mvdr / max(pattern_mvdr); pattern_cbf pattern_cbf / max(pattern_cbf); plot(theta_scan * 180 / pi, 20 * log10(pattern_mvdr), b-, ... theta_scan * 180 / pi, 20 * log10(pattern_cbf), r--); legend(MVDR, 常规波束); xlabel(角度(度)); ylabel(归一化增益(dB)); axis([-90 90 -60 0]); grid on;代码里的核心是w_mvdr (R_hat \ a_s) / (a_s * (R_hat \ a_s))这行做了一次矩阵左除和一次标量归一化。左除比inv(R_hat)更快也更稳定特别是当N变大、R_hat接近奇异的时候。方向图计算部分对整个扫描角范围逐点计算阵列流型并与权值做内积得到的是复增益的幅度。3.3 仿真结果解读零陷位置、主瓣宽度与旁瓣电平运行这段代码会发现三个关键特征。第一MVDR方向图在-20度干扰方向出现了一个深度超过-40dB的零陷而常规波束在那个方向只有大概-13dB的旁瓣。第二MVDR的主瓣方向和宽度与常规波束基本一致因为约束条件保证了信号方向的增益为1。第三MVDR方向图的旁瓣整体抬高了这是自适应波束的自由度被用来抑制干扰所付出的代价。这里有一个比较坑的现象当信噪比很高、快拍数不足时MVDR方向图会出现严重的畸变主瓣歪斜、旁瓣急剧升高。这个现象叫“信号自消”原因是期望信号分量被协方差矩阵中的高能量项误当作干扰而抑制掉了。后面第5章会详细讲这个问题的排查和处理。4. 参数扫描与性能评估阵元数、阵元间距、快拍数对波束的影响4.1 阵元数N对波束宽度和增益的影响阵元数直接决定阵列孔径孔径越大主瓣越窄。半波长间距下16阵元阵列的主瓣宽度大约是6.4度32阵元则压缩到约3.2度这个关系可以近似用公式BW ≈ 0.886·λ/(N·d)估算。仿真时可以用一个简单的循环来观察变化% 扫描阵元数对主瓣宽度的影响 N_list [8, 16, 32, 64]; beam_width zeros(size(N_list)); for idx 1:length(N_list) N N_list(idx); a exp(1j * 2 * pi * 0.5 * (0:N-1) * sin(30 * pi / 180)); theta_scan (-90:0.1:90) * pi / 180; pattern zeros(size(theta_scan)); for k 1:length(theta_scan) a_scan exp(1j * 2 * pi * 0.5 * (0:N-1) * sin(theta_scan(k))); pattern(k) abs(a * a_scan); end pattern pattern / max(pattern); % 找到3dB下降点对应的角度宽度 idx_3db find(pattern 0.707); beam_width(idx) (theta_scan(idx_3db(end)) - theta_scan(idx_3db(1))) * 180 / pi; end disp(beam_width);实际仿真中你会发现除了主瓣变窄阵元数增加还会让旁瓣数量变多第一旁瓣电平仍然维持在约-13dB对均匀加权来说。这解释了为什么实际系统里会用泰勒加权、切比雪夫加权来压低旁瓣——增加阵元数并不能压低相对旁瓣电平。4.2 阵元间距d的约束栅瓣问题与空间混叠阵元间距是相控阵仿真中决定成败的参数。当d/λ大于0.5时方向图会在某些角度出现和主瓣一样高的“栅瓣”这是因为相位差2π·d·sinθ/λ在某些角度又等于整数的2π倍。栅瓣会让波束形成产生模糊的方位判断必须避免。在Matlab里验证栅瓣非常直观把d_lambda 0.5改成1.0再画方向图你能看到在±90度扫描范围内出现第二个高度为1的瓣。实际工程中如果物理条件限制了阵元间距必须大于半波长通常会限制扫描范围来避免栅瓣进入可见区域。我在仿真中一般把d_lambda设为0.4到0.5之间留一点余量同时观察-90到90度全程是否有异常高旁瓣。4.3 快拍数M的影响协方差估计质量与收敛速度快拍数M直接影响协方差矩阵估计的准确性。理论分析表明MVDR的输出的信干噪比损失与(M-N2)/(M1)相关意味着快拍数越接近阵元数性能退化越明显。经验上M至少要是N的2到3倍最好在10倍以上方向图才会比较稳定。仿真时可以做一个有趣的实验固定N16分别跑M30、M100、M1000看零陷深度和主瓣形状的变化。M30时方向图往往出现剧烈抖动主瓣偏移好几度M1000时方向图已经和理论值几乎重合。这里顺带提一个技巧如果你用的信号源是随机相位如exp(1j*randn)每次运行结果都会略有差异这是正常的因为协方差估计有随机性。5. 相控阵波束形成仿真避坑信号自消、矩阵病态与角度失配排查5.1 信号自消导致主瓣塌陷现象、原因与对角加载处理现象加了期望信号的MVDR方向图在某些信噪比和快拍数组合下主瓣方向出现凹陷而不是峰值输出信噪比甚至低于常规波束形成。原因期望信号本身被包含在协方差矩阵R中当快拍数有限且信噪比高时样本协方差矩阵中期望信号分量与噪声分量的耦合使得优化过程把期望信号方向也当成了需要抑制的方向。这就是“信号自消”。解决最有效的工程手段是对角加载。把R替换为R γI其中γ是加载系数I是单位矩阵。加载量会降低自适应算法的自由度但能显著提高稳健性。一般γ取噪声功率的1到10倍例如噪声功率为0.01时γ0.1已经能稳定方向图。修改代码只需要一行gamma 0.1; % 对角加载量 R_loaded R_hat gamma * eye(N); w_mvdr (R_loaded \ a_s) / (a_s * (R_loaded \ a_s));加了加载后零陷会变浅一些但主瓣形状稳定得多。实际雷达系统中对角加载几乎是标配只不过加载量需要根据干扰场景和阵元数做标定。5.2 期望信号方向失配指向误差与波束畸变现象仿真时设定的信号实际方向是32度但代码里theta_s写的是30度结果MVDR在30度形成波束峰值同时在32度形成了一个零陷输出信号功率大幅下降。原因MVDR约束的是“设计方向”的增益为1如果实际信号方位与约束方向偏差超过波束宽度的几分之一信号就会被当作干扰信号抑制掉。解决有两种思路。一是尽可能精确估计信号方向通常用空间谱估计算法做到亚度级精度。二是将方向图展宽用多约束法在期望方向附近增加多个约束点。工程仿真中最简单的做法是加一个指向方向的搜索错位测试扫描从-0.5到0.5度的偏置观察输出信噪比随偏置的变化曲线从而确定系统允许的最大角度误差。我在仿真代码里经常会加这个测试循环每偏置0.1度跑一次MVDR画一条曲线比口头说“系统稳健”有说服力得多。5.3 协方差矩阵求逆的数值问题奇异与条件数过大现象当N16且M30时计算R_hat \ a_s偶尔报出警告“Matrix is singular or nearly singular”方向图变成一片噪声。原因快拍数小于或接近阵元数时样本协方差矩阵不满秩即使满秩条件数也可能大到让浮点求逆产生巨大误差。解决第一优先增加快拍数这是最根本的方案。第二优先使用对角加载加载量可以设到噪声功率的一半也就是0.005左右。如果不允许加载、快拍又少可以考虑降维处理——将全阵列投影到子空间只对主分量做自适应权值求解。但降维处理会改变主瓣宽度实际中需要权衡。5.4 宽带信号的误用把窄带模型直接套到宽带场景现象用这段仿真代码处理带宽200MHz的信号方向图出现严重的频率依赖偏移零陷位置在带宽内移动导致干扰抑制失效。原因窄带模型假设所有阵元的信号只差一个相位而宽带信号在阵元之间存在与频率相关的时延不再是简单的相移。对于宽带信号需要用到聚焦变换或频域子带处理。解决如果必须仿真宽带场景最简单的做法是分频段处理——把信号经FFT分成多个窄带子带每个子带独立做MVDR然后重新合成。Matlab代码中只需对频域数据矩阵重复执行第3章的MVDR流程即可但需要额外注意子带划分数量不能太少否则频率分辨率不够。5.5 方向图绘制时的坐标陷阱角度定义与dB转换现象画出来的方向图在0度位置出现“凹陷”但数组上0度对应的值明明不是最低的。原因很多Matlab方向图脚本里theta_scan用的是-90:90的向量但plot时直接用了向量下标作为横坐标相当于横轴单位变成了“扫描点数”不是角度。另外用20*log10(abs(pattern))时如果pattern里有0值会出现-Inf点导致整条曲线显示异常。解决绘制前统一做theta_deg theta_scan * 180 / pi作为横轴归一化后再取对数加一个eps防止log10(0)对于dB坐标限制纵轴范围在-60到0这样方向图的旁瓣结构才能看清。这个小问题看似低级实际上在很多人分享的仿真代码里反复出现。6. 波束扫描图与方向图动态展示用Matlab动画验证相控阵扫描过程做相控阵仿真不能只输出一张静态方向图最好能动态展示波束扫描过程这样既能验证权值计算的正确性也能直观看到主瓣移动和旁瓣变化。常见的做法是让期望信号方向从-60度逐步扫描到60度每个角度计算一次权值和方向图然后用Matlab的animatedline或者逐帧刷新plot来生成动画效果。% 波束扫描动态展示 theta_targets -60:2:60; % 扫描角度序列 figure; for t_idx 1:length(theta_targets) theta_s theta_targets(t_idx) * pi / 180; a_s exp(1j * 2 * pi * 0.5 * (0:N-1) * sin(theta_s)); % 计算MVDR权值这里复用R_hat或重新生成数据 w (R_hat \ a_s) / (a_s * (R_hat \ a_s)); theta_scan (-90:0.5:90) * pi / 180; pattern zeros(size(theta_scan)); for k 1:length(theta_scan) a_scan exp(1j * 2 * pi * 0.5 * (0:N-1) * sin(theta_scan(k))); pattern(k) abs(w * a_scan); end pattern 20 * log10(pattern / max(pattern) eps); plot(theta_scan * 180 / pi, pattern, b-); ylim([-60, 0]); xlabel(角度(度)); ylabel(增益(dB)); title(sprintf(波束指向 %.0f°, theta_targets(t_idx))); grid on; drawnow; pause(0.05); % 控制动画速度 end这段代码的核心思路是循环更新期望方向并重新计算权值通过drawnow强制刷新图形。如果你要保存成GIF可以在循环里用exportgraphics把每帧写入文件。但有一个性能隐患每次循环内嵌的扫描方向图计算要执行361次复数内积如果阵元数很大比如128阵元整体运行会明显变慢。这时候可以预先构造一个N×361的阵列流型矩阵一次性完成所有角度的响应计算把内层循环向量化% 向量化方向图计算预先构造扫描流型矩阵 theta_scan (-90:0.5:90) * pi / 180; A_scan exp(1j * 2 * pi * 0.5 * (0:N-1) * sin(theta_scan)); pattern abs(w * A_scan); % 一次矩阵运算得到所有角度响应这行代码把原来几百次的循环压缩成一次矩阵乘法在N64、扫描点数721的情况下有接近50倍的提速。这也是Matlab仿真中一个通用的习惯能用矩阵运算表达的逻辑绝不用循环。做这个动态仿真的过程中我发现一个值得注意的细节当波束从-60度扫描到0度再往正角度方向移动时主瓣的物理宽度在低角度时稍宽接近端射方向时明显变宽。这是因为方向图投影到sinθ域是等宽的映射到角度域就不均匀了。如果你在仿真报告里要讨论扫描时的波束展宽问题关注这个现象就对了。另一个经验是固定对角加载量在低角度扫描时表现尚可但扫描到接近±80度时方向图副瓣抬升明显必要时需要根据扫描角调整加载量。这些细节都是实际调试中一点点试出来的比背公式来得直接得多。最后说一个我自己的习惯每次写完段仿真代码第一件事就是固定随机种子也就是在文件开头加rng(2024)。否则每次运行方向图都不一样很难判断代码改动到底有没有生效。等你把算法调稳定了再删掉这行做多次蒙特卡洛验证。这个习惯帮我少走了很多弯路也希望帮到你。本文还有配套的精品资源点击获取