基于MATLAB的MVDR波束形成的推导与仿真

写在前面

在本篇中将利用MATLAB进行MVDR波束形成的推导与仿真，其中最关键的内容为求解MVDR优化问题，不过如果有进行过矩阵论的相关学习应该不会很困难，对应的MATLAB实现也比较简单。
如果你对信通或电子科学方面感兴趣或是正在进行这方面的本科学习，可以去看一下我整理的这篇对信息与电子方向部分本科项目的总结（持续更新中），也许某些课程或是科研项目的总结能对你有所帮助。

1 推导过程

MVDR 的主要思想可以用以下优化问题进行描述，

其中${\mathbf{c}}_{MV}$为 MVDR 的驾驶向量，${\mathbf{R}}_x$为接收信号的自相关函数，$\mathbf{c}$为 CBF 的驾驶向量。由上式可知 MVDR 是通过最小化波束形成输出方差，同时约束观察方向不失真地全通所希望的信号来实现高分辨力且抑制旁瓣的效果。
为了实现MVDR，需要推导出最优的${\mathbf{c}}_{MV}$闭式解和对应的波束形成输出并尽可能化简。推导过程由三部分组成，第一部分为针对$\mathbf{c}$为向量的求解过程，第二部分为正则化方法求解过程，第三部分为针对$\mathbf{c}$为矩阵的求解过程，最终结果在第三部分给出。

1.1 驾驶向量$\mathbf{c}$为向量时的求解过程

对目标函数求Hessian矩阵，易得其为半正定矩阵，故目标函数为凸函数，可采用拉格朗日乘子法(Lagrange multiplier)求解最优解。由拉格朗日乘子法，可令函数

此时通过对$\mathbf{F}$求偏导，解$\mathbf{F}$的最小值，即可得到最优的${\mathbf{c}}_{MV}$闭式解，求解过程如下，

此时可求得最优解，

1.2 正则化方法求解过程

考虑到${\mathbf{R}}_x$未必是非奇异的，此时需要采取正则化方法，因此对正则化的结果进行了推导，推导过程和所求得解如下，

其中$\mathbf{\mu }$ 为很小的扰动，可以大大改善秩亏缺矩阵方程的数值稳定性。
不过从结果来看，通过正则化方法得到的效果并没有比直接求解得到的效果好，因此在本文中仅进行了推导，省略了结果及其分析。

1.3 驾驶向量$\mathbf{c}$为矩阵的求解过程

由于在具体实现过程中发现，我们在仿真过程中所面对的驾驶向量$\mathbf{c}$实际上是一个矩阵，即$\mathbf{C},{\mathbf{C}}_{MV}\in {\mathbf{C}}_{k\times n}$
因此对于推导过程，应该修改如下：
我们将${\mathbf{c}}_{MV}$和$\mathbf{c}$视为矩阵，令

此时，问题转换为

此处为方便起见，记$diag(\mathbf{A})$为一个由矩阵A的对角线元素所构成的对角阵，并且在下文中出现的$diag()$均为该含义，在MATLAB中由diag(diag(A))实现。
令

对上式求导，

此时即得最终解

2 MATLAB代码实现

根据最终推导得到的结果，可以编写相应的MATLAB仿真代码

%单目标MVDR和CBF波束形成，通过调整SNR改变信噪比，根据实验所需信号来选择性注释
clear;clc;fclose('all');
%声源、水平阵、环境参数设置
%环境参数
c = 1500;   %声速1500m/s
SNR = 20;   %信噪比为20dB
%声源（假定在很远处的点源，近场球面波,到达很远处的水平阵的时候可近似为平面波）
SL = 140;   %信号能量140dB
r = 7000;   %与水平阵首阵元相距7000m
f = 50;     %信号频率50Hz
lambda=c/f; %波长
angle=30;   %入射角
%水平阵
M = 9;             %水平阵阵元数9
d = c/(2*f);       %阵元间距15m
angles=-90:0.1:90; %检测角度范围
Nsnapshot=15;      %快拍数（结合多次发射信号）
%计算各阵元接收信号x=s+n,即信号加噪声
%阵列响应向量
v=sqrt(M)\exp(-1j*2*pi*(d*sin(angle*pi/180)/lambda)*(-(M-1)/2:(M-1)/2)'); 
%发射信号
s=sqrt(10^(SL/10))*exp(1j*2*pi*f*(1:Nsnapshot));
%高斯白噪声
n=sqrt(10^((SL-SNR)/10))*(randn(M,Nsnapshot)+1j*randn(M,Nsnapshot))/sqrt(2);
%接收信号
x=sqrt(M)*v*s+n;
%计算响应向量和波束形成响应
%Rx=(M*v*(s*s')*v'+10^((SL-SNR)/10)*eye(M))/Nsnapshot;%理想信号
Rx=(x*x')/Nsnapshot;%加入高斯白噪声的真实信号
%驾驶向量
c=sqrt(M)\exp(-1j*2*pi*(-(M-1)/2:(M-1)/2)'*(d*sin(angles*pi/180)/lambda));
%直接求解闭式解
Cmv1=(Rx\c)/(diag(diag((c'/Rx)*c)));
y1=diag(Cmv1'*Rx*Cmv1);     %MVDR响应
%y1=diag(c'*Rx*c);     %CBF响应
y1=abs(y1)/max(abs(y1));
%绘图
figure(1);
plot(angles,10*log10(y1),'k','linewidth',2);
xlabel('Angle(deg)','Fontsize',15);ylabel('Power Response(dB)','Fontsize',15) ;
title('单目标MVDR波束形成（实际信号）','Fontsize',20);

3 输出结果和分析

接下来在MATLAB内对理想信号和加入高斯白噪声的实际信号进行仿真并分析，同时为了体现MVDR的效果，将其与CBF的结果进行对比。

3.1 信噪比为20dB

从上图中可见虽然从CBF图像中能看出大致信号源方向，但主瓣较宽且旁瓣较大，与CBF相对的，理想波束MVDR结果的图像旁瓣基本消失，图像比较平滑，主瓣宽度也大大减小，分辨率很高，与理论预测相符合，对于增加了高斯白噪声的模拟实际信号，MVDR波束形成结果相对于理想波束虽然有些旁瓣存在，但仍可以有很好的分辨率，旁瓣数值低且主瓣窄，这也说明实验结果正确，起到了很好的优化效果。

3.2 不同信噪比结果比较

此处由于结果较多，不另外放图。从结果上看，当信噪比提升时，MVDR形成波束的主瓣宽度逐渐减小，30°处即信源方向变得越来越尖，分辨率越来越高，旁瓣数值也降低，曲线整体更加平滑。当信噪比降低时，从结果图来看MVDR形成波束稳定性逐渐降低，虽然仍可分辨出信号源方向，但主瓣宽度变大并且旁瓣数值升高，当信噪比取负数时已经有噪声能量大于信源发射信号的情况出现，旁瓣峰值接近主瓣，此时MVDR辨识将非常困难。

3.3 在20dB下增加不同角度目标结果比较

从结果上看，未加入白噪声时，如果两个目标相隔角度较大（例如-30°和30°），CBF和MVDR所形成的波束均能成功识别出两个目标的方向，其中MVDR分辨率高，两个主瓣宽度窄且旁瓣数值低，因此有很好的目标检测及定位的能力，与理论预测相符，但是两个目标角度小时，如20°和30°，MVDR仍能分辨两个主峰，但CBF无法进行辨识，这是由于该方法所产生旁瓣较大所导致的分辨率较低。
加入白噪声后，CBF的实验结果除了旁瓣增大以外没产生比较大的变化，但MVDR明显出现比较大的变动，不仅分辨率降低，在20°和30°的情况下无法辨识两个主峰，在目标相隔距离较大时（-30°和30°）也无法产生稳定的结果，很大概率导致旁瓣峰值接近甚至大于主瓣，导致其失去辨识的作用，因此可知在实际情况中MVDR并不适合多目标检测，还是使用CBF较为保险。
至此，仿真已全部完成。