简介：这是一份面向通信、雷达与信号处理学习者的 MATLAB 波束赋形（Beamforming）技术文档，以单个 doc 格式文件封装，大小约 1.23MB。文档围绕均匀线阵方向图展开，完整提供 8 阵元、16/128/1024 阵元等多种配置下的仿真代码，系统比较不同阵元数对波束宽度、分辨率的影响，并进一步演示栅瓣产生条件、最优权傅里叶变换等效性，以及基于最大信噪比准则的方向图和功率谱生成流程，内容由浅入深，覆盖波束赋形核心原理。文件内不仅包含可直接运行的 MATLAB 源码，也配有仿真结果分析和文字讲解，读者可据此理解相位加权、阵列响应、自适应权值优化等关键概念，并自行修改角度、阵元数、干扰方向等参数进行扩展实验。目前已有 335 人学习使用，特别适合需要快速掌握多种波束成形算法比较、天线数量影响分析以及干扰抑制技巧的高校学生和工程技术人员。
1. 从 MATLAB 仿真实例拆解波束成形：这份代码能让你少走一半弯路
做阵列信号处理的人八成都有过这种经历：理论推导头头是道，一写 MATLAB 仿真就翻车——画出来的方向图跟教科书长得不一样，明明设了零陷却对不上干扰方向，改个阵元数波束宽度变化跟预期不符。这份 matlab beamforming 文档正好是这么一套能直接跑的波束成形仿真实例集，覆盖了均匀线阵方向图、最大信噪比准则、MSE 准则、LCMV 准则和 Capon 波束形成，每种算法都附了完整代码和仿真图。它不是教科书式的原理重讲，而是把波束成形从阵元配置到权矢量求解的仿真链路打通了，适合正在做课程设计、毕业设计或者刚接手阵列信号处理项目的人拿来做基线参考。我拆完这份代码后发现，它的价值不止在能跑通，更在于参数怎么调、阵元数怎么影响波束宽度、栅瓣什么时候出现、不同准则下的协方差矩阵估计有什么区别，这些在代码注释和对比仿真里都能直接看到。
2. 均匀线阵方向图与波束宽度：从阵元配置看懂仿真背后的物理规律
2.1 阵元数与半波长间距的选择依据
打开文档里的第一段代码，8 阵元均匀线阵方向图，核心变量就三个：element_num=8定义阵元数，d_lamda=1/2定义阵元间距与波长的比值，theta0=0定义来波方向。代码用扫描向量theta=linspace(-pi/2,pi/2,200)遍历空域角度，每个角度下计算导向矢量a，再与最优权w做内积得到方向图响应p(j)=w'*a。
clc; clear all; close all; imag=sqrt(-1); element_num=8; % 阵元数为 8 d_lamda=1/2; % 阵元间距 d 与波长 lamda 的关系 theta=linspace(-pi/2,pi/2,200); theta0=0; % 来波方向 w=exp(imag*2*pi*d_lamda*sin(theta0)*[0:element_num-1]'); for j=1:length(theta) a=exp(imag*2*pi*d_lamda*sin(theta(j))*[0:element_num-1]'); p(j)=w'*a; end figure; plot(theta,abs(p)), grid on xlabel('theta/radian') ylabel('amplitude') title('8 阵元均匀线阵方向图')这里d_lamda=1/2是最常见的配置，物理意义是阵元间距取半波长。选择半波长的理由有两个：一是避免栅瓣，二是保证空间采样不模糊。如果间距大于半波长，方向图会出现与主瓣等高的栅瓣，这是个必须记住的边界条件。w的构造是把来波方向的相位补偿预先算好，相当于把阵列指向theta0，然后扫描整个空域看响应。注意w'*a的结果是复数，取abs(p)才是幅度方向图，很多新手直接plot(theta, p)会看到震荡的实部虚部混叠图，这不是代码错了，是复数绘图没取模。
2.2 波束宽度如何随阵元数收窄
文档第二节专门做了 16、128、1024 三个阵元数量的对比，代码用asin(sin(theta*pi/180) - lamda/(element_num*d))计算主瓣零点位置，从而得到波束宽度。这段代码的核心是推导波束零点的解析表达式，不是靠仿真曲线去量，而是用公式直接算：
clc; clear all; close all; ima=sqrt(-1); element_num1=16; % 阵元数 element_num2=128; element_num3=1024; lamda=0.03; % 波长为 0.03 米 d=1/2*lamda; % 阵元间距与波长的关系 theta=0:0.5:90; for j=1:length(theta) fai(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_num1*d)); psi(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_num2*d)); beta(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_num3*d)); end figure; plot(theta, fai, 'r', theta, psi, 'b', theta, beta, 'g'), grid on xlabel('theta'); ylabel('Width in radians') title('波束宽度与波达方向及阵元数的关系')从仿真结果能直观看到：阵元数从 16 增加到 1024，波束宽度从零点几弧度收窄到接近零。这条曲线是衡量阵列分辨力的关键依据。实际做系统设计时，阵元数不是越多越好，因为阵元数翻倍意味着通道数翻倍，硬件成本和计算量同步上涨。折中做法是先根据目标分辨力反推波束宽度要求，再确定最小阵元数。另一个细节是波束宽度还跟扫描角有关——波束指向阵列法线方向时最窄，越往端射方向扫描波束越宽，这个规律在代码里通过theta从 0 扫到 90 度的曲线斜率变化也能看出来。
2.3 栅瓣现象与空间模糊的边界条件
文档提到当阵元间距大于半波长时会出现栅瓣，导致空间模糊。栅瓣的本质是空间混叠：阵元间距太大，不同方向的来波在阵列上产生的相位差可能相差 2π 的整数倍，阵列无法区分这两个方向。最典型的反面例子就是d_lamda=1时，在端射方向附近会出现与主瓣等高的栅瓣。仿真时想复现这个现象，把d_lamda改成 1 或者更大，重新跑第一节的方向图代码就能看到多个等高的波峰。
实际工程里栅瓣不总是坏事——在某些稀疏阵列设计中，栅瓣位置是已知的，可以通过适当设计避开。但在常规均匀线阵里，栅瓣意味着方位估计出现模糊，必须避免。判断标准很简单：阵元间距必须满足d <= lamda/2。如果由于物理条件限制不得不大于半波长，就需要用非均匀阵列或者超分辨算法来弥补，这就超出这份代码的范畴了。
3. 自适应波束形成算法实现：从最大信噪比到 Capon 波束形成的完整仿真链路
3.1 最大信噪比准则：广义特征分解求最优权
第五段代码是最大信噪比准则的完整实现。场景设置为 8 阵元、来波方向 0 度、干扰方向 20 度，采样单元数 L=512。代码用随机幅度模拟信号和干扰，构造信号自相关矩阵Rs和干扰加噪声自相关矩阵Rnj，然后通过广义特征分解求最优权：
for i=1:L amp0=10*randn(1); amp1=200*randn(1); ampn=1; s(:,i)=amp0*exp(ima*2*pi*1/2*sin(theta0*pi/180)*[0:element_num-1]'); j(:,i)=amp1*exp(ima*2*pi*1/2*sin(theta1*pi/180)*[0:element_num-1]'); n(:,i)=ampn*(randn(element_num,1)+ima*randn(element_num,1)); end Rs=1/L*s*s'; % 信号自相关矩阵 Rnj=1/L*(j*j'+n*n'); % 干扰+噪声的自相关矩阵 [V,D]=eig(Rs,Rnj); % (Rs,Rnj) 的广义特征值和特征向量 [D,I]=sort(diag(D)); % 排序 Wopt=V(:,I(8)); % 最优权矢量这段代码的关键是eig(Rs, Rnj)的用法。MATLAB 的eig(A, B)返回的是广义特征值问题的解，即A*v = lambda*B*v。最大信噪比准则下，最大广义特征值对应的特征向量就是最优权。排序后取最后一个特征向量V(:, I(8))，因为sort默认升序，I(8)对应最大特征值的位置。这里有个容易踩的坑：如果矩阵是奇异的，eig可能返回 NaN 或 Inf，所以采样数 L 必须大于阵元数，否则自相关矩阵不满秩。
得到最优权后，方向图的绘制方式和第一节不同——这里用了subplot(121)和subplot(122)分别画方向图和功率谱。方向图用20*log10(abs(f)/max(abs(f)))归一化到 dB 单位，方便观察旁瓣电平和零陷深度。功率谱用p(j)=a'*Rs*a+a'*Rnj*a计算，实际上是在每个扫描方向计算信号加干扰噪声的总功率。仿真结果里能看到在 20 度干扰方向形成较深的零陷，这是最大信噪比准则的典型表现。
3.2 Capon 波束形成：最小方差无失真响应的实现细节
第八段代码是标准的 Capon 波束形成，也叫最小方差无失真响应（MVDR）波束形成。核心公式是w = R_inv * steer / (steer' * R_inv * steer)，其中steer是期望方向的导向矢量，R_inv是接收数据自相关矩阵的逆：
Rx=1/L*x*x'; R=inv(Rx); steer=exp(ima*2*pi*1/2*sin(theta0*pi/180)*[0:element_num-1]'); w=R*steer/(steer'*R*steer); % Capon 最优权矢量 for j=1:length(theta) a=exp(ima*2*pi*d_lamda*sin(theta(j)*pi/180)*[0:element_num-1]'); f(j)=w'*a; p(j)=1/(a'*R*a); end与最大信噪比准则相比，Capon 的权矢量求解不需要知道干扰方向，只需要知道期望方向。它通过最小化输出功率来抑制干扰，同时保证期望方向增益为 1。这份代码里有个细节值得注意：p(j)=1/(a'*R*a)是 Capon 功率谱估计的公式，它用扫描导向矢量遍历整个空域，在每个方向计算对应的输出功率，峰值位置就是信号到达角估计值。
需要提醒的是，inv(Rx)在实际系统中不是首选做法，因为矩阵求逆的数值稳定性不如线性方程求解。特别是当阵元数较多、快拍数有限时，Rx可能是病态的，直接求逆会放大数值误差。更稳妥的做法是用Rx\steer代替inv(Rx)*steer，或者在Rx的对角线上加一个小的对角加载因子。文档里的代码为了教学清晰用了inv，实际做项目时可以改进。
3.3 LCMV 准则：多点约束下的波束形成
第七段代码实现了线性约束最小方差（LCMV）准则，它把 Capon 扩展到多个约束方向。场景是 8 阵元、0 度、30 度、60 度三个信号方向，约束矩阵C由三个方向的导向矢量组成，约束响应向量F=[1 0 1]'表示 0 度和 60 度方向增益为 1，30 度方向增益为 0：
C=[steer1 steer2 steer3]; F=[1 0 1]'; % 把三个方向都作为来波方向 w=inv(Rx)*C*(inv(C'*inv(Rx)*C))*F;LCMV 的解是w = R_inv * C * (C' * R_inv * C)^(-1) * F。这段代码的约束设置很有意思——把 30 度方向约束为零增益，相当于在 30 度方向形成零陷，但代价是 30 度方向的信号完全被抑制。实际应用中，约束向量的设计需要权衡：约束方向越多，自由度消耗越大，旁瓣电平可能会升高。这份代码的仿真图对比了 LCMV 与 Capon 的方向图和功率谱，LCMV 在三个约束点都精确满足了增益要求。
我一般做 LCMV 仿真时会把约束矩阵C的列做归一化，因为不同方向的导向矢量模长相同，但数值范围可能随阵元数变化。归一化后约束条件的数值稳定性更好，避免C'*inv(Rx)*C条件数过大导致求逆失败。
4. 波束成形仿真避坑指南：五个最容易翻车的细节与排查方法
4.1 方向图出现多根等高峰值：先查阵元间距
现象：仿真出的方向图在多个角度出现幅度相同的波峰，看起来像有几个主瓣。
原因：阵元间距d_lamda设置大于 0.5，导致栅瓣出现。这是空间混叠的直接表现，阵列无法区分栅瓣方向和真实来波方向。
解决：检查代码里的d_lamda是否为1/2。如果确实需要大于半波长的间距，考虑改用稀疏阵列并配合解模糊算法，但常规仿真务必保持d_lamda <= 1/2。另外注意theta=linspace(-pi/2, pi/2, 200)的采样间隔，如果栅瓣恰好落在采样点之间，方向图显示可能不完整，可以把采样点数增加到 400 再看看。
4.2 零陷位置偏了几个度：协方差矩阵估计方式要背锅
现象：Capon 或 LCMV 方向图中干扰方向的零陷不在预设角度，偏移了 2-5 度。
原因：协方差矩阵Rx是用有限快拍数估计的，Rx=1/L*x*x'是真实协方差矩阵的近似。快拍数 L 越小，估计误差越大，零陷位置偏移越明显。另外，如果信号和干扰的幅度不是随机独立产生的，可能出现信号与干扰的相干性，导致协方差矩阵秩亏。
解决：把 L 从 512 增大到 1000 或 2048，观察零陷是否落在正确位置。文档后面专门有一节对比了不同协方差矩阵估计方法对 Capon 波束形成的影响——直接用接收信号估计Rx和只用干扰加噪声估计Rnj的零陷深度有明显差异。如果必须用小快拍数，可以在Rx上做对角加载：Rx_loaded = Rx + epsilon*eye(element_num)，epsilon 取trace(Rx)/element_num * 0.01量级。
4.3 矩阵奇异导致 NaN 输出：采样数不够或信源相干
现象：运行最大信噪比或 Capon 代码时，w或方向图输出出现 NaN 或 Inf，程序不报错但结果不可用。
原因：eig(Rs, Rnj)或inv(Rx)时矩阵奇异。常见原因有两个：一是采样数 L 小于阵元数，自相关矩阵秩不足；二是多个信号完全相干（比如幅度不成随机变化），导致协方差矩阵的秩进一步下降。
解决：确保L > element_num，最好L > 2*element_num。信号幅度用randn随机产生，避免相干。检查amp0、amp1是否设置成固定常数——如果固定不变，每次采样都一样，信号就完全相干了，这时可以用空间平滑技术解相干。
4.4 方向图不对称：阵元权矢量相位没对齐
现象：方向图主瓣不在 0 度或预设的theta0方向，整体向左或向右偏。
原因：构造w时用的导向矢量与扫描时用的导向矢量表达式不一致。比如w构造时用了sin(theta0)而扫描时用了sin(theta(j)*pi/180)，角度单位混用（弧度与度混用）导致相位计算错误。
解决：统一角度单位。文档代码中theta0=0时sin(theta0)和sin(0*pi/180)结果相同所以没出问题，但改成theta0=45就暴露了：一个用弧度一个用度，计算结果完全不同。建议全篇统一用弧度，或者统一在构造导向矢量的函数里先转弧度。我的习惯是写一个独立的导向矢量函数，输入角度自动转弧度，避免每次手写*pi/180。
4.5 旁瓣电平太高：权矢量没有归一化或阵元数太少
现象：方向图主瓣和旁瓣幅度接近，看起来没有明显的波束形状。
原因：归一化方式不对。20*log10(abs(f)/max(max(abs(f))))是常见的归一化方式，但如果max(abs(f))不是主瓣峰值而是旁瓣峰值，整个方向图的比例就会失真。另外 8 阵元的均匀线阵本身旁瓣电平就在 -13dB 左右（相对于主瓣），这是理论极限，不是代码问题。
解决：先画出线性幅度方向图确认主瓣位置，再用 dB 归一化。如果旁瓣实在太高影响视觉效果，可以用窗函数加权（比如汉明窗）压低旁瓣，但主瓣会相应变宽。文档里的代码没有做窗函数处理，保持的是均匀加权的最原始状态，这也是理解波束形成原理的最佳起点——先看均匀加权的效果，再做改进。
5. 自适应波束形成大作业变体：从固定参数到全场景参数扫描的调参技巧
文档最后一部分是自适应波束形成的大作业，和前面的示例不同，它用的是实际载频 f=500e6、波长lamda=c/f，阵元间距d=1/2*lamda。这部分的价值在于它做了多组参数对比：L=512 和 L=2048 的采样数对比，不同信号幅度amp0、amp1、amp2组合下的方向图对比。把这些变体组合起来，就能系统性地看出三个关键参数对自适应波束形成的影响。
c=3e8; f=500e6; lamda=c/f; d=1/2*lamda; theta0=0; theta1=45; theta2=60; L=512; for i=1:L amp0=10*randn(1); amp1=100*randn(1); amp2=100*randn(1); ampn=10; x(:,i)=amp0*exp(ima*2*pi*1/2*sin(theta0*pi/180)*[0:element_num-1]')+... amp1*exp(ima*2*pi*1/2*sin(theta1*pi/180)*[0:element_num-1]')+... amp2*exp(ima*2*pi*1/2*sin(theta2*pi/180)*[0:element_num-1]')+... ampn*(randn(element_num,1)+ima*randn(element_num,1)); end Rx=1/L*x*x'; steer1=exp(ima*2*pi*1/2*sin(theta0*pi/180)*[0:element_num-1]'); steer2=exp(ima*2*pi*1/2*sin(theta1*pi/180)*[0:element_num-1]'); steer3=exp(ima*2*pi*1/2*sin(theta2*pi/180)*[0:element_num-1]'); C=[steer1 steer2 steer3]; F=[1 0 0]'+0.001*randn(3,1); % 约束向量加微扰，模拟实际系统不完美约束 w=inv(Rx)*C*(inv(C'*inv(Rx)*C))*F;我的做法是保留文档原始代码作为基线，然后在此基础上做三个方向的扩展：一是把L从 512 改到 2048，对比协方差矩阵估计精度对零陷深度的影响；二是把amp0从 10 改到 100，观察信号功率增强后约束方向增益的保持情况；三是把干扰幅度从 100 降到 10，观察自适应算法在低干扰功率下是否还能形成有效零陷。每组参数跑完都保存方向图数据，方便对比。
我还习惯在主循环外面加一层参数遍历，把阵元数作为外层循环变量，自动生成不同阵元数下的方向图叠加对比图。这样能一次性看清阵元数从 8 到 32 变化时主瓣宽度、旁瓣电平和零陷深度的趋势。用一个小技巧可以提升效率——预先分配矩阵，把x的维度设为[element_num, L]而非在循环里动态增长，L 较大时速度能快好几倍。
element_nums = [8, 16, 32]; figure; hold on; for e_idx = 1:length(element_nums) element_num = element_nums(e_idx); x = zeros(element_num, L); for i=1:L % 填充 x 的代码，与文档一致 end Rx = 1/L * x * x'; % 计算权矢量 w % 扫描方向图，保存到 f_e plot(theta, 20*log10(abs(f_e)/max(abs(f_e)))); end legend('8 阵元', '16 阵元', '32 阵元');这类参数扫描式的跑法，是我拿到一份新仿真代码后必做的第一步复现验证。如果只跑通原代码就收工，你根本不知道哪些参数是临界值。把那几个关键参数——阵元数、采样快拍数、信号幅度比值——分别扫一遍，这份代码才算真正吃透了。从那以后我每拿到一份波束成形代码，都强制自己走一遍这个流程，先复现基线，再扫关键参数，最后再考虑改进算法。希望这些拆解和调参经验对你有帮助。
本文还有配套的精品资源，点击获取
