CN108595758A - 一种任意形式传感器阵列最优宽带波束图综合方法 - Google Patents

Abstract

本发明涉及一种任意形式传感器阵列最优宽带波束图综合方法，首先将传感器接收到的平面波信号展开成不同的级数形式，并利用不同频率处的权值向量可以表达成期望波束权值向量的函数，构建约束优化问题，获得满足约束条件的最优期望波束，最后通过解析的函数关系获得其它频率处的权值向量，合成得到宽带波束。本发明公开的方法可以在满足约束条件下获得最优的期望波束，从而得到更灵活的频率不变响应宽带波束，比现有技术公开的方法更优越。

Description

一种任意形式传感器阵列最优宽带波束图综合方法

技术领域

本发明属于一种波束图综合方法，涉及一种任意形式传感器阵列最优宽带波束图综合方法，适用于任意形式传感器阵列的宽带目标检测以及宽带目标方位估计，属于声学阵列信号处理、语音信号处理和声呐技术等领域。

背景技术

波束图综合在声呐、雷达、通信以及语音工程等领域有大量应用，其主要功能是设计具有特定形状的波束图，以达到抑制噪声和干扰、提高输出信噪比、改善方位估计精度等目的。对实际应用中常使用的宽带信号进行处理时，需要利用波束图综合方法设计频率不变响应波束形成器，以保证宽带信号的无畸变输出，从而有利于提取更精确的信息，如波形、频谱结构、量级以及方位等。由此可见，波束图综合在传感器阵列信号处理中发挥了重要作用，受到了很多研究者的关注，并提出了很多方法，主要有：文献1“任意结构形状传感器阵方向图的最佳化.中国造船,1984,87(4):78-85”公开的凹槽噪声场方法，该方法通过人为地在旁瓣区域增加干扰并利用迭代算法可以得到低旁瓣的宽带波束。文献2“Optimalarray pattern synthesis using semidefinite programming,IEEE Trans.SignalProcess.,2003,51(5):1172-1183”公开的波束图优化方法，该方法利用半正定规划可得到逼近于期望波束的宽带波束。文献3“Optimal array pattern synthesis for broadbandarrays.J.Acoust.Soc.Am.,2007,122(5):2686-2696”公开的宽带波束图综合方法，该方法在旁瓣满足一定约束条件下使得宽带波束主瓣之间误差最小化，结合二阶锥优化可以得到不同形式的宽带波束图。上述文献2公开的方法需要已知的期望波束，而文献1和3公开的方法则不需要期望波束，但三种方法均没有解析闭式解，需要数值计算，并且对于每个频率都需要重新计算，计算量较大。文献4“一种圆环形阵列恒定束宽波束形成器设计方法.中国:CN103903609B,2016”公开的期望波束拟合方法通过最小化合成波束和期望波束之间的均方误差，得到了圆环阵在不同频率上权值向量的精确闭式解，可方便地合成宽带频率不变响应波束。但该方法的期望波束只能事先给定，不能保证期望波束是满足需求的最佳选择，且仅适用于圆环阵。文献5“圆环阵最优频率不变响应超指向性波束形成方法.中国:CN107170441 A,2017”在文献3方法的基础上公开了一种最优频率不变响应宽带波束图综合方法，该方法可以获得满足要求的最优频率不变响应波束，但同样仅适用于圆环阵。

发明内容

要解决的技术问题

为了避免现有技术的不足之处，本发明提出一种任意形式传感器阵列最优宽带波束图综合方法，利用了下述性质克服了现有技术不够精确和适用范围有限的不足：依据阵列形式的不同，传感器接收到的平面波信号可以展开成不同的级数形式，且不同频率处的权值向量是一般化期望波束权值向量的函数。

技术方案

一种任意形式传感器阵列最优宽带波束图综合方法，其特征在于步骤如下：

步骤1、计算最优期望波束权值向量：

构建约束优化问题：

其中：Ω_SL是波束旁瓣区域，Ω₀＝(θ₀,φ₀)为期望波束指向角，|·|表示绝对值，||·||₂表示Euclidean范数，上标T表示转置；

所述B_d表示期望波束响应向量，其元素为：B_d(θ,φ)＝a^Tdiag{Ξ^*(θ₀,φ₀)}Ξ(θ,φ)

其中：a是实的权值向量，Ξ是与波束形状相关的向量，(θ₀,φ₀)是波束指向角，(θ,φ)∈Ω_SL，θ表示垂直俯仰角，φ表示水平方位角，上标“*”表示求共轭；

所述矩阵矩阵

其中：f为入射平面波的频率，f₁是感兴趣的频率范围内的最小值，f₂是感兴趣的频率范围内的最大值，L为f₁到f₂之间离散化点的个数，Σ₁＝Q-V^Hρ^-1V，Σ₂＝V^Hρ^-2V，其中ρ为任意形式传感器阵列的空域相关矩阵，上标H表示共轭转置；

对于三维波束：

矩阵Q＝diag{|Ξ(θ₀,φ₀)|²/(4π)}，

a＝[a₀,a_-1,1,a_0,1,a_1,1,...,a_-n,n,...,a_0,n,...,a_n,n,...,a_m,n,...,a_N,N]^T为(N+1)²×1维由实数构成的向量，同样为(N+1)²×1维的向量；

其中： 为缔合勒让德函数，N(N>0)为最高阶数，diag{·}表示对角阵，(·)！表示阶；

所述矩阵V＝[v₀ v_-1,1 v_0,1 v_1,1... v_m,n ... v_N,N]，

其中：b_n＝iⁿj_n(kr_s)，j_n(·)是第n阶球形Bessel函数，r_s是第s号阵元与坐标原点的距离；

所述矩阵ρ的第m行m′列的元素为ρ_mm′＝sinc(2πd_mm′/λ)，d_mm′为第m和m′号阵元间的距离，λ是信号波长；

对于二维波束、即与阵共面的波束：

所述矩阵Q＝I_{(2N+1)×(2N+1)}，I表示单位矩阵，

a＝[a_-N,a_-N+1,…,a₀,...,a_N-1,a_N]^T为(2N+1)×1维由实数构成的向量，Ξ(φ)＝[Y_-N(φ),Y_-N+1(φ),...,Y₀(φ)，...,Y_N-1(φ),Y_N(φ)]^T同样为(2N+1)×1维的向量；

其中Y_n(φ)＝e^imφ，N(N>0)为最高阶数，θ₀＝π/2；

所述矩阵V＝[v_-N v_-N+1 … v₀ … v_N-1 v_N]，

v_n＝[b_nY_n(φ₀),b_nY_n(φ₁),...,b_nY_n(φ_M-1)]^T，

其中：b_n＝iⁿJ_n(kr_s)，J_n(·)是第n阶柱形Bessel函数，r_s是第s号阵元与坐标原点的距离；

所述矩阵ρ的第m行m′列的元素为ρ_mm′＝J₀(2πd_mm′/λ)，d_mm′为第m和m′号阵元间的距离，λ是信号波长，J₀(·)是第0阶柱形Bessel函数；

步骤2、计算在所需频率范围内合成波束的权值向量：

合成波束的权值向量的计算：w＝ρ^-1Tdiag{Ξ(θ₀,φ₀)}a

所述w＝[w₀,w₁,...,w_M-1]^T；将权值向量a代入上式即得到所需要的合成波束权值向量w；

步骤3、合成最终波束：

将步骤2得到的权值向量w代入下式，即得到合成波束：

所述P(θ,φ)＝[p₀(θ,φ),p₁(θ,φ),...,p_M-1(θ,φ)]^T为阵列流形向量；p_s(θ,φ)为第s个阵元接收到的从方向(θ,φ)入射的单位幅度平面波信号，表达式为：

k＝ku为波数向量，k＝2π/λ，

u＝-[sinθcosφ,sinθsinφ,cosθ]^T，r_s＝[r_ssinθ_scosφ_s,r_ssinθ_ssinφ_s,r_scosθ_s]^T表示第s号阵元的位置向量，(r_s,θ_s,φ_s)为第s号阵元的球坐标。

有益效果

本发明提出的一种任意形式传感器阵列最优宽带波束图综合方法，首先将传感器接收到的平面波信号展开成不同的级数形式，并利用不同频率处的权值向量可以表达成期望波束权值向量的函数，构建约束优化问题，获得满足约束条件的最优期望波束，最后通过解析的函数关系获得其它频率处的权值向量，合成得到宽带波束。

有益效果体现在：

本发明提出一种适用于任意形式传感器阵列的最优宽带波束图综合方法，有益效果体现在：

1.本发明公开的方法对于所有频率仅需要一次数值计算，比文献1、2和3公开的方法所需计算量更小。

2.本发明公开的方法不受阵列形式限制，适用于任意形式传感器阵列，相比于文献4和5公开的方法，其使用范围更广。

本发明公开的方法可以在满足约束条件下获得最优的期望波束，从而得到更灵活的频率不变响应宽带波束，比文献1、2、3和4公开的方法更优越。

附图说明

图1是直线形传感器阵列示意图。

图2是直线形阵列最优期望波束。

图3是直线形阵列不同方法的指向性指数(三维噪声场中)。

图4是直线形阵列不同方法的白噪声增益。

图5是直线形阵列合成的频率不变响应宽带波束。图5(a)是波束随频率变化的三维显示，图5(b)是波束随频率变化的叠加显示。

图6是“V”形传感器阵列示意图。

图7是“V”形阵列最优期望波束。

图8是“V”形阵列不同方法的指向性指数(二维噪声场中)。

图9是“V”形阵列不同方法的白噪声增益。

图10是“V”形阵列合成的频率不变响应宽带波束。图10(a)是水平波束(θ＝90°)随频率变化的三维显示，图10(b)是水平波束(θ＝90°)随频率变化的叠加显示。

图11是半球形传感器阵列示意图。

图12是半球形阵列最优期望波束。图12(a)是三维显示，图12(b)是水平波束(θ＝60°)，图12(c)是垂直波束(φ＝180°)。

图13是半球形阵列不同方法的指向性指数(三维噪声场中)。

图14是半球形阵列不同方法的白噪声增益。

图15是半球形阵列合成的频率不变响应宽带波束。图15(a)是水平波束(θ＝60°)随频率变化的三维显示，图15(b)是水平波束(θ＝60°)随频率变化的叠加显示，图15(c)是垂直波束(φ＝180°)随频率变化的三维显示，图15(d)是垂直波束(φ＝180°)随频率变化的叠加显示。

具体实施方式

现结合实施例、附图对本发明作进一步描述：

本发明涉及的一种适用于任意形式传感器阵列的最优宽带波束图综合方法，可在最小均方误差准则下给出不同频率处的权值向量与期望波束权值向量的函数，并可获得满足约束条件的最优期望波束，进而合成得到最优宽带波束。其过程为：

1.计算最优期望波束权值向量。

构建如下约束优化问题：

其中Ω_SL是波束旁瓣区域，Ω₀＝(θ₀,φ₀)为期望波束指向角，|·|表示绝对值，||·||₂表示Euclidean范数，上标T表示转置。所述B_d表示期望波束响应向量，其元素为

B_d(θ,φ)＝a^Tdiag{Ξ^*(θ₀,φ₀)}Ξ(θ,φ) (2)

其中a是实的权值向量，Ξ是与波束形状相关的向量，(θ₀,φ₀)是波束指向角，(θ,φ)∈Ω_SL，θ表示垂直俯仰角，φ表示水平方位角，上标“*”表示求共轭。

所述矩阵矩阵其中f为入射平面波的频率，f₁是感兴趣的频率范围内的最小值，f₂是感兴趣的频率范围内的最大值，L为f₁到f₂之间离散化点的个数，Σ₁＝Q-V^Hρ^-1V，Σ₂＝V^Hρ^-2V，其中ρ为任意形式传感器阵列的空域相关矩阵，上标H表示共轭转置。

对于三维波束：

所述矩阵Q＝diag{|Ξ(θ₀,φ₀)|²/(4π)}，

a＝[a₀,a_-1,1,a_0,1,a_1,1,...,a_-n,n,...,a_0,n,...,a_n,n,...,a_m,n,...,a_N,N]^T为(N+1)²×1维由实数构成的向量，同样为(N+1)²×1维的向量，其中 为缔合勒让德函数，N(N>0)为最高阶数，diag{·}表示对角阵，(·)！表示阶乘。

所述矩阵V＝[v₀ v_-1,1 v_0,1 v_1,1 ... v_m,n ... v_N,N]，

其中b_n＝iⁿj_n(kr_s)，j_n(·)是第n阶球形Bessel函数，r_s是第s号阵元与坐标原点的距离。

所述矩阵ρ的第m行m′列的元素为ρ_mm′＝sinc(2πd_mm′/λ)，d_mm′为第m和m′号阵元间的距离，λ是信号波长。

对于二维波束，即与阵共面的波束：

所述矩阵Q＝I_{(2N+1)×(2N+1)}，I表示单位矩阵，a＝[a_-N,a_-N+1,...,a₀,...，,a_N-1,a_N]^T为(2N+1)×1维由实数构成的向量，Ξ(φ)＝[Y_-N(φ),Y_-N+1(φ),...,Y₀(φ),...,Y_N-1(φ),Y_N(φ)^T]同样为(2N+1)×1维的向量，其中Y_n(φ)＝e^imφ，N(N>0)为最高阶数，θ₀＝π/2。

所述矩阵V＝[v_-N v_-N+1 ... v₀ ... v_N-1 v_N]，

v_n＝[b_nY_n(φ₀),b_nY_n(φ₁),...,b_nY_n(φ_M-1)]^T，其中b_n＝iⁿJ_n(kr_s)，J_n(·)是第n阶柱形Bessel函数，r_s是第s号阵元与坐标原点的距离。

所述矩阵ρ的第m行m′列的元素为ρ_mm′＝J₀(2πd_mm′/λ)，d_mm′为第m和m′号阵元间的距离，λ是信号波长，J₀(·)是第0阶柱形Bessel函数。

对于式1的使用，有以下说明：约束目标是使得合成波束与期望波束之间的均方误差最小，参数δ用于控制波束旁瓣响应，参数μ用于控制稳健性，参数σ用于控制波束指向性。上式所示的约束优化问题可由二阶锥规划进行求解。

参照图1和2。考虑一个均匀直线形阵列，如图1所示，相关参数为：M＝4，N＝3，(θ₀,φ₀)＝(0°,180°)，声速c＝1500m/s，阵元间距0.1m，f₁＝100Hz，f₂＝1600Hz，(相当于无旁瓣约束)，δ＝∞(相当于无旁瓣约束)，μ＝∞(相当于无白噪声增益约束)，σ＝10^-9/10，L＝31。

由式(1)计算得到的直线形阵列最优期望权值向量a见表1，未在表中列出的均为0。由该权值向量得到的最优期望波束如图2所示，其指向性指数约等于9dB，满足约束条件。

表1直线形阵列最优期望权值系数

a0	a0,1	a0,2	a0,3
				2.0108	2.0442	0.8101	0.0532

参照图6和7。考虑一个“V”形阵列，阵元所在位置的角度坐标列于表2，阵列形式如图6所示。

表2“V”形阵列各阵元所在位置的角度坐标

阵元序号	1	2	3	4	5	6	7
								X坐标(m)	0	-0.0050	-0.0100	-0.0150	0.0050	0.0100	0.0150
Y坐标(m)	0	-0.0087	-0.0173	-0.0260	-0.0087	-0.0173	-0.0260
								Z坐标(m)	0	0	0	0	0	0	0

其它相关参数为：M＝7，N＝2，(θ₀,φ₀)＝(90°,130°)，声速c＝340m/s，f₁＝200Hz，f₂＝6400Hz，Ω_SL＝{(θ,φ)|θ＝90°,φ∈[0°,φ₀-Δ]∪[φ₀+Δ,360°]}，Δ＝102°，δ＝10^-30/20，μ＝∞(相当于无白噪声增益约束)，σ＝10^-6/10，L＝32。

由式(1)计算得到的最优期望权值系数向量a列于表3，由该权值向量得到的最优期望波束如图7所示，其旁瓣级均在-30dB以下，指向性指数大于6dB，满足约束条件。

表3“V”形阵列最优期望权值系数

a-2	a-1	a0	a1	a2
					0.1000	0.2424	0.3152	0.2424	0.1000

参照图11和12。考虑一个半球形阵列，半径a＝0.042m，阵元所在位置的角度坐标列于表4，阵列形式如图11所示。

表4半球形阵列各阵元所在位置的角度坐标

阵元序号	1	2	3	4	5	6	7	8	9	10	11
												方位角(°)	180	0	72	144	216	-72	-36	36	108	180	252
俯仰角(°)	0	63.4	63.4	63.4	63.4	63.4	37.4	37.4	37.4	37.4	37.4
												阵元序号	12	13	14	15	16	17	18	19	20	21
方位角(°)	36	72	108	144	180	216	252	-72	-36	0
												俯仰角(°)	79.2	100.8	79.2	100.8	79.2	100.8	79.2	100.8	79.2	100.8

其它相关参数为：M＝21，N＝4，(θ₀,φ₀)＝(60°,180°)，声速c＝340m/s，f₁＝400Hz，f₂＝3200Hz，Ω_SL＝{(0°,0°),(60°,0°),(60°,70°),(60°,120°),(120°,180°),(140°,0°)}，δ＝10^-30/²⁰，μ＝10⁶⁰/¹⁰(相当于白噪声增益下限为-60dB)，σ＝10^-12/¹⁰，L＝29。

由式(1)计算得到的最优期望权值系数向量a列于表5，由该权值向量得到的最优期望波束如图12所示，其所约束的旁瓣级均在-30dB以下，指向性指数约等于12dB，最小白噪声增益等于-60dB，满足约束条件。

表5半球形阵列最优期望权值系数

a0	a-1,1	a0,1	a1,1	a-2,2	a-1,2	a0,2	a1,2	a2,2	a-3,3	a-2,3	a-1,3	a0,3
													1.0681	0.6886	0.6456	0.6886	0.5786	0.9298	1.0861	0.9298	0.5786	1.0912	0.6919	0.4089	0.3328
a1,3	a2,3	a3,3	a-4,4	a-3,4	a-2,4	a-1,4	a0,4	a1,4	a2,4	a3,4	a4,4
													0.4089	0.6919	1.0912	0.1123	0.0722	0.1445	-0.1014	0.0575	-0.1014	0.1445	0.0722	0.1123

2.计算在所需频率范围内合成波束的权值向量。

合成波束的权值向量的由下式计算：

w＝ρ^-1Tdiag{Ξ(θ₀,φ₀)}a (3)

所述w＝[w₀,w₁,...,w_M-1]^T。将式1计算得到的权值向量a代入式3，即可得到所需要的合成波束权值向量w。

3.合成最终波束。

将式3得到的权值向量代入下式，即可得到合成波束：

所述P(θ,φ)＝[p₀(θ,φ),p₁(θ,φ),...,p_M-1(θ,φ)]^T为阵列流形向量，p_s(θ,φ)为第s个阵元接收到的从方向(θ,φ)入射的单位幅度平面波信号，表达式为：k＝ku为波数向量，k＝2π/λ，u＝-[sinθcosφ,sinθsinφ,cosθ]^T，r_s＝[r_ssinθ_scosφ_s,r_ssinθ_ssinφ_s,r_scosθ_s]^T表示第s号阵元的位置向量，(r_s,θ_s,φ_s)为第s号阵元的球坐标。

参照图3、4和5。针对如图1所示的直线形阵列合成得到的频率不变响应宽带波束如图5所示，显示了较好的频率不变特性。对应的指向性指数和白噪声增益分别如图3和图4所示，其中指向性指数控制在9dB左右，稳健性较MVDR方法也有所改善。

参照图8、9和10。针对如图6所示的“V”形阵列合成得到的频率不变响应宽带波束如图10所示，显示了较好的频率不变特性。对应的指向性指数和白噪声增益分别如图8和图9所示，其中指向性指数控制在6dB以上，部分波束的旁瓣级略高于-30dB，稳健性较MVDR方法也有所改善。

参照图13、14和15。针对如图11所示的半球形阵列合成得到的频率不变响应宽带波束如图15所示，显示了较好的频率不变特性。对应的指向性指数和白噪声增益分别如图13和图14所示，其中所约束的旁瓣级基本都在-30dB左右，指向性指数控制在12dB左右，白噪声增益大于等于-60dB。

Claims (1)

1.一种任意形式传感器阵列最优宽带波束图综合方法，其特征在于步骤如下：

步骤1、计算最优期望波束权值向量：

构建约束优化问题：

其中：Ω_SL是波束旁瓣区域，Ω₀＝(θ₀,φ₀)为期望波束指向角，|·|表示绝对值，||·||₂表示Euclidean范数，上标T表示转置；

所述B_d表示期望波束响应向量，其元素为：B_d(θ,φ)＝a^Tdiag{Ξ*(θ₀,φ₀)}Ξ(θ,φ)

其中：a是实的权值向量，Ξ是与波束形状相关的向量，(θ₀,φ₀)是波束指向角，(θ,φ)∈Ω_SL，θ表示垂直俯仰角，φ表示水平方位角，上标“*”表示求共轭；

所述矩阵矩阵

其中：f为入射平面波的频率，f₁是感兴趣的频率范围内的最小值，f₂是感兴趣的频率范围内的最大值，L为f₁到f₂之间离散化点的个数，Σ₁＝Q-V^Hρ^-1V，Σ₂＝V^Hρ^-2V，其中ρ为任意形式传感器阵列的空域相关矩阵，上标H表示共轭转置；

对于三维波束：

矩阵Q＝diag{|Ξ(θ₀,φ₀)|²/(4π)}，

a＝[a₀,a_-1,1,a_0,1,a_1,1,...,a_-n,n,...,a_0,n,...,a_n,n,...,a_m,n,...,a_N,N]^T为(N+1)²×1维由实数构成的向量，同样为(N+1)²×1维的向量；

其中： 为缔合勒让德函数，N(N>0)为最高阶数，diag{·}表示对角阵，(·)！表示阶；

所述矩阵V＝[v₀ v_-1,1 v_0,1 v_1,1 ... v_m,n ... v_N,N]，

其中：b_n＝iⁿj_n(kr_s)，j_n(·)是第n阶球形Bessel函数，r_s是第s号阵元与坐标原点的距离；

所述矩阵ρ的第m行m′列的元素为ρ_mm′＝sinc(2πd_mm′/λ)，d_mm′为第m和m′号阵元间的距离，λ是信号波长；

对于二维波束、即与阵共面的波束：

所述矩阵Q＝I_{(2N+1)×(2N+1)}，I表示单位矩阵，

a＝[a_-N,a_-N+1,...,a₀,...,a_N-1,a_N]^T为(2N+1)×1维由实数构成的向量，Ξ(φ)＝[Y_-N(φ),Y_-N+1(φ),...,Y₀(φ),...,Y_N-1(φ),Y_N(φ)]^T同样为(2N+1)×1维的向量；

其中Y_n(φ)＝e^imφ，N(N>0)为最高阶数，θ₀＝π/2；

所述矩阵V＝[v_-N v_-N+1 ... v₀ ... v_N-1 v_N]，

v_n＝[b_nY_n(φ₀),b_nY_n(φ₁),...,b_nY_n(φ_M-1)]^T，

其中：b_n＝iⁿJ_n(kr_s)，J_n(·)是第n阶柱形Bessel函数，r_s是第s号阵元与坐标原点的距离；

所述矩阵ρ的第m行m′列的元素为ρ_mm′＝J₀(2πd_mm′/λ)，d_mm′为第m和m′号阵元间的距离，λ是信号波长，J₀(·)是第0阶柱形Bessel函数；

步骤2、计算在所需频率范围内合成波束的权值向量：

合成波束的权值向量的计算：w＝ρ^-1Tdiag{Ξ(θ₀,φ₀)}a

所述w＝[w₀,w₁,...,w_M-1]^T；将权值向量a代入上式即得到所需要的合成波束权值向量w；

步骤3、合成最终波束：

将步骤2得到的权值向量w代入下式，即得到合成波束：

所述P(θ,φ)＝[p₀(θ,φ),p₁(θ,φ),...,p_M-1(θ,φ)]^T为阵列流形向量；p_s(θ,φ)为第s个阵元接收到的从方向(θ,φ)入射的单位幅度平面波信号，表达式为：

k＝ku为波数向量，k＝2π/λ，

u＝-[sinθcosφ,sinθsinφ,cosθ]^T，r_s＝[r_ssinθ_scosφ_s,r_ssinθ_ssinφ_s,r_scosθ_s]^T表示第s号阵元的位置向量，(r_s,θ_s,φ_s)为第s号阵元的球坐标。