转载：Capon波束形成器

首发于 麦克风阵列学习

写文章

MVDR波束形成器

longRookie3

划水

​ 关注他

​

目录

收起

一、背景

二、公式推导

三、代码尝试

一、背景

常规波束形成器的阵增益有限，为了最大限度地提高阵增益，Capon于1969年提出MVDR波束形成器，又称Capon波束形成器。

二、公式推导

MVDR波束形成器的设计原理是对感兴趣方位的信号无失真地输出，而且波束输出噪声（可以是包含干扰的广义噪声）方差最小。

这里无失真地输出，需要满足无失真约束（ ）,波束输出噪声方差最小（对于高斯白噪声的平均功率会等于它的方差）。

波束输出噪声方差为

于是，MVDR波束形成加权向量设计问题表述为：

采用Lagrange算子，定义函数：

将该函数对 求导，令该导数为 ,得到：

代入 ,得到：

将公式（2.5）代入公式（2.4）得到：

若噪声场为空间白噪声，则退化为常规波束形成器加权向量。

但是，这样有个问题，我们无法预先估计噪声协方差矩阵 ，这时的解决方法往往是直接用数据协方差矩阵 来代替噪声协方差矩阵 。那么，该波束设计问题为：

有些文献将该方法称为最小功率无失真响应（mininum power distortionless response,MPDB）波束形成器，也有很多将该方法仍称为MVDR波束形成法，公式（2.2）可以看做公式（2.8）的特例，后面会提到。

---

考虑如下阵列接收数据模型：

若存在平面波干扰，可以将干扰纳入噪声 中。

数据协方差矩阵表示为：

式中，

式中， 是高斯白噪声； 表示非白噪声（有色噪声）分量，例如干扰是一种非白噪声。

当 时，公式（2.8）和公式（2.2）等价。

当 时，此时的MVDR波束形成器加权向量为：

该MVDR波束响应式为

波束输出功率谱为

阵增益为

三、代码尝试

10元标准线列阵，假设基阵接收噪声是功率为0dB的高斯白噪声，即 。一个信噪比 的单频信号从 方向入射到基阵。

数据协方差计算公式

加权向量计算采用公式（2.12）

环境设置：

clear; close all; clc;

%噪声：空间白噪声
c=343; %声速
f=16000; %频率
lambda = c / f;
k = 2pif / c;
space = lambda/2; %麦克风间距
M = 10; %麦克风数量
P=[1:M];

观察方向为 和 的波束图：

%% 计算波束图

theta_angle = -90:1:90;

theta = theta_angle*pi/180;

theta_d1 = 10pi/180; %入射角度1
p_1 = exp(1jkspacePsin(theta_d1))';
Rx = (p_1)p_1'10^(30/10)+eye(M);
B = zeros(size(theta));
Wcc = Rx</span>p_1 / (p_1'/Rxp_1);

for i = 1:length(theta)
p = exp(1ikspacePsin(theta(i)))‘;
B(i) = Wcc’p;
end
B_db = 20log10(B);
figure;
plot(theta_angle, B_db,‘-.’, ‘linewidth’, 1.5);
hold on;

theta_d2 = -10pi/180; %观察方向
p_2 = exp(1jkspacePsin(theta_d2))‘;
Wcc = Rx</span>p_2 / (p_2’/Rxp_2);
for i = 1:length(theta)
p = exp(1ikspacePsin(theta(i)))‘;
B(i) = Wcc’p;
end
B_db = 20log10(B);
plot(theta_angle, B_db, ‘linewidth’, 1.5);

%plot(10, 30, ‘ro’, ‘linewidth’, 1.5);
grid on;
xlim([-90,90]);
legend(‘\theta_0=10^\circ’, ‘\theta_0=-10^\circ’);
xlabel(‘\theta /(\circ)’); ylabel(‘20lg B(\theta)/dB’);
title(‘波束图’);

MVDR波束扫描估计方位谱：

%% 计算波束扫描方位谱

theta_angle = -90:1:90;

theta = theta_angle*pi/180;

theta_d1 = 10pi/180; %入射角度1
p_1 = exp(1jkspacePsin(theta_d1))';
Rx = (p_1)p_1'*10^(30/10)+eye(M);
B = zeros(size(theta));

for i = 1:length(theta)
p = exp(1ikspacePsin(theta(i)))‘;
Wcc = Rx</span>p / (p’/Rxp);
B(i) = Wcc'RxWcc;
end
B_db = 10log10(B);
figure;
plot(theta_angle, B_db, ‘linewidth’, 1.5);
hold on;
plot(10, 30, ‘ro’, ‘linewidth’, 1.5);
grid on;
ylim([-15,40]);
legend(‘方位谱’, ‘真实信号’);
xlabel(‘方位（^o）’); ylabel(‘10lg§/dB’);
title(‘波束扫描方位谱’);

加权向量范数：

%% 计算加权向量范数

theta_angle = -90:0.1:90;

theta = theta_angle*pi/180;

theta_d1 = 10pi/180; %入射角度1
p_1 = exp(1jkspacePsin(theta_d1))';
Rx = (p_1)p_1'*10^(30/10)+eye(M);
B = zeros(size(theta));

for i = 1:length(theta)
p = exp(1ikspacePsin(theta(i)))‘;
Wcc = Rx</span>p / (p’/Rxp);
B(i) = norm(Wcc)norm(Wcc);
end
B_db = 10*log10(B);
figure;
plot(theta_angle, B_db, ‘linewidth’, 1.5);

grid on;
ylim([-15,30]);
xlabel(‘方位（^o）’); ylabel(‘10lg§/dB’);
title(‘加权向量范数’);