阅读原文还请移步个人知乎专栏:
https://www.zhihu.com/column/c_1287066237843951616
首先力荐鄢社锋老师的书籍《优化阵列信号处理》,对于麦克风阵列信号处理的讲解很是全面,对阵列的基本概念与优化方法讲的非常透彻,并且最重要事情固然是全中文的啦,并且鄢社锋是本领域的顶尖专家,比不少翻译过来的外文书籍要好不少,值得阅读。我的感受比陈老师的《Microphone Array Signal Processing》要更容易理解、更容易上手不少,由于这本书里面包含了不少的案例用以对比演示。可是很遗憾,我并无找到鄢老师的相关代码,因此干脆就本身动手,看书的同时顺带实现了几个。优化
因此,本文是后续系列文章的第一篇,记录案例实现的过程与代码,对于案例的解释主要为原书内容,本身进行了简化,便于后续理解。本文讲述第二章的案例2.二、2.3与2.4,对比试验了常规波束造成与MVDR最优波束造成的波束图、方位扫描图与阵增益,以及本身实现的代码。spa
例2.2 常规波束图、常规波束扫描方位谱
考虑一个10元标准线列阵,假设一功率为1的单频信号从 方向入射到基阵,不考虑噪声,用以下公式计算数据协方差矩阵:翻译
以基阵在 方向的响应向量做为导向向量,即
,采用式设计
计算常规波束加权向量,用式3d
波束图code
计算波束响应。将获得的波束图显示于图 2(a)中。图中能够看见波束主瓣指向方位。这意味着该波束造成器对
方向到达的信号响应最大。附代码实现,在图片上方。orm
%% 2(a) 常规波束造成波束图 代码实现部分 clear; close all; clc; c=340; %声速 f=1000; %频率 lambda = c / f; k = 2*pi*f / c; space = lambda/2; %麦克风间距 M = 10; %麦克风数量 theta_d = 80*pi/180; %入射角度 theta_angle = 0:0.1:180; theta = theta_angle*pi/180; Wc = exp(-1i*k*space*[0:M-1]*cos(theta_d)) / M; B = zeros(size(theta)); %% 计算波束图 for i = 1:length(theta) p = exp(1i*k*space*[0:M-1]*cos(theta(i))); B(i) = conj(Wc)*p'; end B_db = 20*log10(B); limit_dB = -50; index = B_db < limit_dB; B_db(index) = limit_dB; figure; plot(theta_angle, B_db, 'linewidth', 1.5); grid on;xlabel('方位(^o)'); ylabel('20lg(B)/dB'); title('波束图');
2(a) 常规波束造成波束图blog
波束扫描方位谱
进行波束扫描,即让 在
内变化,取导向向量为
。用下式
图片
计算出波束扫描方位谱,显示于图 2(b)中。从图中能够看出,在 位存在一个主峰值,该峰值表示在该方位存在一个信号。峰值大小为 0dB,表示该方位波束输出功率为1,即该方向到达信号功率估计值为1。同时注意到,虽然只有
方向存在信号,其余方位没有信号,但从方位谱上观察,其余方位的波束输出功率并不为0,好像是出现了“能量泄漏”。这是因为波束旁瓣引发的,对于远离信号方向的波束,因为常规波束的旁瓣比较高,即便信号位于这些波束的旁瓣,旁瓣对该信号有必定的抑制,可是波束输出并不为0,所以表现为在方位谱上非信号方向存在输出功率。能够预见,因为常规波束较高的旁瓣,可能会由于强信号的“能量泄漏”而淹没了其余方位的弱信号,这一点在后文的仿真中将会看到。ci
比较图2(a)与图2(b)发现,指向信号方向的常规波束图与常规波束扫描方位谱曲线彻底相同。这是在仅存在单位功率平面波信号时常规波束造成特有的现象,其余波束造成方法不存在此现象。假设基阵接收到功率为 0dB的空间白噪声,信号功率增长到 30dB,即输入信噪比 。采用常规波束扫描获得的方位谱如图 2(c)所示,图中指示在
方向方位谱大约为 30dB,表示该方位信号功率估计值约为 30dB。代码实现参考 2(b),再也不赘述。
%% 2(b) 常规波束扫描方位谱,代码实现部分 p = exp(1i*k*space*[0:M-1]*cos(theta_d)); Rx = (p')*conj(p); B = zeros(size(theta)); for i = 1:length(theta) Wcc = exp(-1i*k*space*[0:M-1]*cos(theta(i))) / M; B(i) = conj(Wcc)*Rx*Wcc'; end B_db = 10*log10(B); limit_dB = -50; index = B_db < limit_dB; B_db(index) = limit_dB; figure; plot(theta_angle, B_db, 'linewidth', 1.5); grid on; xlabel('方位(^o)'); ylabel('10lg(P)/dB'); title('波束扫描方位谱');
2(b) 常规波束扫描方位谱
2(c) 常规波束扫描方位图
波束阵增益
应用公式:
计算扫描到各方位时波束阵增益,显示于图 2(d)中。图中可见,该曲线形状与图2(a)中波束形状相同,但高出10dB ,即为最高阵增益。从该图能够知,当波束对准信号方向时,波束造成器的阵增益最大;波束观察方位与信号方位存在误差时,波束阵增益减少。只要方位误差不太大(在半功率束宽之内),阵增益比最大值降低得并很少。
%% 2(d) 常规波束阵增益 代码实现部分 p_ = exp(1i*k*space*[0:M-1]*cos(theta_d)); B = zeros(size(theta)); for i = 1:length(theta) p = exp(1i*k*space*[0:M-1]*cos(theta(i))); B(i) = (conj(p_) * p').^2 / (conj(p_)*p_'); end B_db = 10*log10(B); limit_dB = -50; index = B_db < limit_dB; B_db(index) = limit_dB; figure; plot(theta_angle, B_db, 'linewidth', 1.5); grid on; xlabel('方位(^o)'); ylabel('10lg(G)/dB'); title('不一样方向波束阵增益');
2(d) 常规波束阵增益
例2.3 MVDR波束造成
考虑一个10元标准线列阵,假设基阵接收噪声是功率为0dB 高斯白噪声,即 。一信噪比为
的单频信号从
方向入射到基阵,用下式计算数据协方差矩阵
。
MVDR波束造成器波束图
采用下式
设计MVDR 波束造成器加权向量。图 3(a)显示了观察方向分别为与
的两个波束图。对于
方向波束造成器,它至关于
的状况,
方向平面波被视做干扰;对于
方向波束,
方向平面波为指望信号,对应于
的状况。
对于方向的波束图,它在
方向的响应为
。因为MVDR波束造成器的特性,为了使波束输出功率最小,该波束在
方向造成零点以最大限度抑制来自该方向的干扰。对于
方向的波束图,假想导向向量与真实基阵响应向量相等,又因为
,该MVDR波束造成器蜕化为常规波束造成器,因此该波束图与图 2(a)中的常规波束图相同。
%% 3(a) MVDR波束造成器波束图,代码实现 c=340; %声速 f=1000; %频率 lambda = c / f; k = 2*pi*f / c; space = lambda/2; %麦克风间距 M = 10; %麦克风数量 theta_angle = 0:0.1:180; theta = theta_angle*pi/180; theta_d1 = 80*pi/180; %入射角度 p_1 = exp(1i*k*space*[0:M-1]*cos(theta_d1))'; Rx = conj(p_1)*p_1'*10^(30/10) + eye(M)*10^(0/10); Wc = Rx\p_1 / (conj(p_1')*inv(Rx)*p_1); B = zeros(size(theta)); %% 计算波束图 80度 for i = 1:length(theta) p = exp(1i*k*space*[0:M-1]*cos(theta(i)))'; B(i) = conj(Wc')*p; end B_db = 20*log10(B); limit_dB = -50; index = B_db < limit_dB; B_db(index) = limit_dB; figure; plot(theta_angle, B_db, 'linewidth', 1.5); hold on;
3(a) MVDR波束造成器波束图
MVDR波束扫描方位谱
采用MVDR波束扫描估计方位谱,其中扫描方位间隔为 ,获得的方位谱显示于图 3(b)中。与图 2(b)相比较可见,MVDR波束扫描获得的方位谱峰值更尖锐。可见,MVDR波束扫描能提供比常规波束扫描较高的方位分辨能力。
%% 3(b) MVDR波束造成器波束扫描方位谱图,代码实现 theta_angle = 0:1:180; theta = theta_angle*pi/180; theta_d1 = 80*pi/180; %入射角度1 p_1 = exp(1i*k*space*[0:M-1]*cos(theta_d1))'; Rx = (p_1)*conj(p_1')*10^(30/10); B = zeros(size(theta)); for i = 1:length(theta) p = exp(1i*k*space*[0:M-1]*cos(theta(i)))'; Wcc = Rx\p / (conj(p')*(Rx\p)); B(i) = conj(Wcc')*Rx*Wcc; end B_db = 10*log10(B); limit_dB = -10; index = B_db < limit_dB; B_db(index) = limit_dB; figure; plot(theta_angle, B_db, 'linewidth', 1.5); hold on; plot(80, 30, 'ro', 'linewidth', 1.5); grid on; legend('方位谱', '真实信号'); xlabel('方位(^o)'); ylabel('10lg(P)/dB'); title('波束扫描方位谱');
3(b) MVDR波束造成器波束扫描方位谱图
MVDR波束造成器加权向量范数
图 3(v)显示了各扫描方位波束的加权向量范数 。从图中能够看出,在波束观察方向距离信号方向较远时(例如本例中与信号间隔
以上),波束造成器加权向量范数较小,接近于最小值。间隔较近时,随着间隔减少,加权向量范数逐渐增长;而当波束方向刚好等于信号方向时,加权向量范数忽然降低到与间隔较远的方位同一大小(由于此时波束造成器蜕化为常规波束造成器)。
%% 3(c) MVDR波束造成器加权向量范数,代码实现 theta_angle = 0:0.1:180; theta = theta_angle*pi/180; theta_d1 = 90*pi/180; %入射角度1 p_1 = exp(1i*k*space*[0:M-1]*cos(theta_d1))'; Rx = (p_1)*conj(p_1')*10^(30/10); B = zeros(size(theta)); for i = 1:length(theta) p = exp(1i*k*space*[0:M-1]*cos(theta(i)))'; Wcc = Rx\p / (conj(p')*(Rx\p)); B(i) = norm(Wcc)*norm(Wcc); end B_db = 10*log10(B); limit_dB = -50; index = B_db < limit_dB; B_db(index) = limit_dB; figure; plot(theta_angle, B_db, 'linewidth', 1.5); grid on; xlabel('方位(^o)'); ylabel('10lg(w^2)/dB'); title('加权向量范数');
3(c) MVDR波束造成器加权向量范数
例2.4 两信号源状况下MVDR与常规波束造成器的比较
考虑一个 10元标准线列阵,基阵接收噪声是功率为 0dB高斯白噪声,两个信噪比分别为 15dB与 30dB的单频信号分别从 与
方向入射到基阵。
分别采用常规波束造成与MVDR波束造成扫描,获得的方位谱分别如图4(a) 与图4(b) 所示。
4(a) 两信号源常规波束造成器
4(b) 两信号源MVDR波束造成器
参考资料
一、《优化阵列信号处理》,鄢社锋