CN103902832B - 一种圆环形阵列实数加权超指向性波束形成方法 - Google Patents

Abstract

本发明涉及一种圆环形阵列实数加权超指向性波束形成方法，利用具有偶数个阵元的均匀圆环形阵列的对称性质，将波束图、指向性因子和误差敏感度函数均变换为实数域的计算形式，且权值向量的维数减小，然后将超指向性波束形成构造为多约束优化问题，在满足旁瓣级和稳健性约束的条件下使得指向性因子最大，最后利用成熟的二阶锥规划方法求解实数权值向量。本发明利用了下述性质克服了计算量偏大的不足：具有偶数个阵元的均匀圆环形阵列的波束图、指向性因子和误差敏感度函数都可以在实数域计算，并且计算维度减小；二阶锥规划在实数域的计算效率高于在复数域的效率。

Description

一种圆环形阵列实数加权超指向性波束形成方法

技术领域

本发明属于一种波束形成方法，涉及一种圆环形阵列实数加权超指向性波束形成方法，适用于圆环形阵列的低信噪比目标检测以及目标方位的高分辨率估计，属于水声学、阵列信号处理和声纳技术等领域。

背景技术

圆环形阵列在声纳、雷达、通信以及语音工程等领域都有十分广泛的应用，关于该阵列的波束形成方法也是层出不穷，其中超指向性波束形成方法更是受到了极大的关注。超指向性方法是相对于常规的延迟求和方法而言的，后者对各通道信号进行延时求和处理以提取所需信息，虽然操作简单，性能稳健，但能提供的阵增益和空间指向性十分有限，尤其在较低频段更是如此。另外估计目标方位时存在“瑞利限”，难以较好地分辨低频目标。相比之下，超指向性方法可以在不增加阵列尺寸的前提下，获得超过延迟求和方法所能得到的指向性，或者在同样的指向性条件下大大减小阵列孔径，因此在探测低信噪比目标以及提高目标方位分辨率方面都具有很大的应用潜力。

已有的关于圆环形阵列的超指向性波束形成方法主要有：文献1“任意几何形状和阵元指向性的传感器阵列优化波束形成方法,声学学报,2005,vol.30(3),p.264-270”公开的稳健超指向性方法，此方法虽适用于任意形状的阵列，但对圆环形阵列而言，并不最优。这是因为其未充分利用圆环形阵列的性质，由二阶锥规划求解复数权值向量的计算量偏大，效率较低。文献2“Theoretical and practical solutions for high-ordersuperdirectivity of circular sensor arrays,IEEE Trans.Ind.Electron.,2013,Vol.60(1),p.203-209”公开的特征波束分解与综合方法，该方法给出了圆环形阵列高阶超指向性的精确闭式解，通过一种降秩处理技术获得了不错的超指向性结果，其计算简单，但未考虑超指向性波束旁瓣过高的问题。

发明内容

要解决的技术问题

为了避免现有技术的不足之处，本发明提出一种圆环形阵列实数加权超指向性波束形成方法，针对均匀圆环形阵列提出一种实数加权超指向性波束形成方法，可有效减小计算量，提高效率。

技术方案

一种圆环形阵列实数加权超指向性波束形成方法，其特征在于步骤如下：

步骤1：均匀圆环形阵列第s个阵元接收到的从方向(θ,φ)入射的单位幅度平面波信号为p_s(θ,φ)=e^{-ikasinθcos(φ-φs)}，其中为圆环形阵列半径，k=2π/λ，λ表示入射平面波的波长，φ_s=sβ，β=2π/M；

所述均匀圆环形阵列包含M个阵元，且M为偶数；

步骤2：将p_s(θ,φ)构成阵列流形向量P(θ,φ)=[p₀(θ,φ),p₁(θ,φ),...p_M-1(θ，φ)]^T，然后变换为新的阵列流形向量 $\hat{E}(\mathrm{\θ},\mathrm{\φ})=[E_0(\mathrm{\θ},\mathrm{\φ}),E_1(\mathrm{\θ},\mathrm{\φ}),...,E_{M-1}(\mathrm{\θ},\mathrm{\φ}){]}^T,$ 相应元素的变换关系为：

$E_m(\mathrm{\θ},\mathrm{\φ})=v_m^{H}P(\mathrm{\θ},\mathrm{\φ})=\frac{1}{\sqrt{M}}\underset{s=0}{\overset{M-1}{\mathrm{\Σ}}}{e}^{-\mathrm{ism\β}}\·p_s(\mathrm{\θ},\mathrm{\φ})$

其中v_m=M^-1/2[1 e^imβ...e^i(M-1)mβ]^T，上标T表示转置，上标H表示共轭转置；新阵列流形向量的元素具有对称性：其中“*”表示复数共轭；

步骤3：将波束图变换为：B(θ,φ)=η^TF(θ,φ)，其中η=[η₀,η₁,...,η_M/2]^T为所要求解的实数权值向量，表示Hadamard积，(θ₀,φ₀)为预先设定的波束指向角，E=[E₀,E₁,...,E_M/2]^T，ε=[ε₀,ε₁,...,ε_M/2]^T；

所述向量ε元素的取值为： ${\mathrm{\ϵ}}_m=\left\{\begin{array}{cc}1,& m=0,M/2\\ 2,& m=\mathrm{1,2},...,M/2-1\end{array}\right.;$

步骤4：将指向性因子变换为其中

Σ=diag{ε₀|E₀|²,ε₁|E₁|²,...,ε_M/2|E_M/2|²}，Λ_n=diag{λ₀,λ₁,...,λ_M/2}；

所述矩阵Λ_n元素为实数，其值为：

式中ρ_s=sinc(k·Δr_s),Δr_s=2asin(sβ/2)为第m和m'号阵元间的距离，s=|m-m′|；

步骤5：将误差敏感度函数变换为

步骤6：将步骤2～步骤5变换的参量代入下式

s.t.F(Ω₀)^Tη=1，η^TΣη≤σ，

F(Ω_j)^Tη≤δ_j，Ω_j∈Ω_SL,j=I,...,N_sL

其中Ω=(θ,φ)，Ω_SL为选定的旁瓣区域，σ是误差敏感度函数的上界，即稳健性约束值，δ_j是期望旁瓣高度，δ_j是以主瓣为中心左右对称取值的；

步骤7：采用两种方式合成最终波束图

1、将计算得到的实数权值向量η直接代入B(θ,φ)=η^TF(θ,φ)得到最终需要的波束图；

2、由以下步骤：

a）由η得到 $\widehat{\mathrm{\η}}={[{\mathrm{\η}}_0,{\mathrm{\η}}_1,...,{\mathrm{\η}}_M]}^T,$ 其中η_M-m=η_m；

b）由式得到

c）由式求出一般化的权值向量w，其中V=[v₀,v₁,...,v_M-1]；

d）将w代入式B(θ,φ)=w^HP(θ,φ)得到最终所需要的波束图。

有益效果

本发明提出的一种圆环形阵列实数加权超指向性波束形成方法，利用具有偶数个阵元的均匀圆环形阵列的对称性质，将波束图、指向性因子和误差敏感度函数均变换为实数域的计算形式，且权值向量的维数减小，然后将超指向性波束形成构造为多约束优化问题，在满足旁瓣级和稳健性约束的条件下使得指向性因子最大，最后利用成熟的二阶锥规划方法求解实数权值向量。本发明利用了下述性质克服了计算量偏大的不足：具有偶数个阵元的均匀圆环形阵列的波束图、指向性因子和误差敏感度函数都可以在实数域计算，并且计算维度减小；二阶锥规划在实数域的计算效率高于在复数域的效率。

有益效果体现在：

1.本发明利用具有偶数阵元的均匀圆环形阵列的对称性质，将波束图、指向性因子和误差敏感度函数均变换为实数域的计算形式，并且权值向量的维数得以减小，相比于文献1的方法，利用二阶锥规划计算的效率提高了。

2.本发明将超指向性方法构造为一个多约束优化问题，可以在指向性、旁瓣级和稳健性之间获得最好的折中。相比于文献2，本发明虽然没有闭式解，但计算更灵活，可以得到更符合实际需要的结果。

附图说明

图1是本发明方法的计算过程示意图。

图2是本发明方法所用的均匀圆环形阵列示意图。

图3是本发明方法和文献1方法得到的等旁瓣超指向性波束图与常规方法得到的ka=2时波束图；

图4是本发明方法和文献1方法得到的非等旁瓣超指向性波束图与常规方法得到的ka=2时波束图。

具体实施方式

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

本发明的圆环形阵列实数加权超指向性波束形成方法，圆环形阵列的阵元数为偶数，利用对称性质将超指向性波束形成转换到实数域，并将其构造为简单的多约束优化问题，通过二阶锥规划计算权值向量，而需要求解的权值向量均为实数，且维数减少。其过程为：

（1）均匀圆环形阵列包含M个阵元，且M为偶数。记p_s(θ,φ)为第s个阵元接收到的从方向(θ,φ)入射的单位幅度平面波信号，表达式为：

p_s(θ,φ)=e^{-ikasinθcos(φ-φs)}, （1）其中a为圆环形阵列半径，k=2π/λ，λ表示入射平面波的波长，φ_s=sβ，β=2π/M。转换的具体步骤如下：

a）将阵列流形向量P(θ,φ)=[p₀(θ,φ),p₁(θ,φ)，...，p_M-1(θ,φ)^T]变换为新的阵列流形向量 $\hat{E}(\mathrm{\θ},\mathrm{\φ})=[E_0(\mathrm{\θ},\mathrm{\φ}),E_1(\mathrm{\θ},\mathrm{\φ}),...,E_{M-1}(\mathrm{\θ},\mathrm{\φ}){]}^T,$ 相应元素的变换关系为：

$E_m(\mathrm{\θ},\mathrm{\φ})=v_m^{H}P(\mathrm{\θ},\mathrm{\φ})=\frac{1}{\sqrt{M}}\underset{s=0}{\overset{M-1}{\mathrm{\Σ}}}{e}^{-\mathrm{ism\β}}\·p_s(\mathrm{\θ},\mathrm{\φ}),---\left(2\right)$

其中v_m=M^-1/2[1 e^imβ...e^i(M-1)mβ]^T，上标T表示转置，上标H表示共轭转置。新阵列流形向量的元素具有对称性：其中“*”表示复数共轭。

b）将波束图变换为如下形式：

B(θ,φ)=η^TF(θ,φ), （3）

其中η=[η₀,η₁,...,η_M/2]^T为所要求解的实数权值向量，(θ₀,φ₀)为预先设定的波束指向角，E=[E₀,E₁,...,E_M/2]^T，ε=[ε₀,ε₁,...,ε_M/2]^T。向量ε元素的取值为：

${\mathrm{\ϵ}}_m=\left\{\begin{array}{cc}1& (m=0,M/2)\\ 2& (m=\mathrm{1,2},...,M/2-1).\end{array}\right.----\left(4\right)$

c）将指向性因子变换为如下形式：

其中Σ=diag{ε₀|E₀|²,ε₁|E₁|²,...,ε_M/2|E_M/2|²}，Λ_n=diag{λ₀,λ₁,...,λ_M/2}。矩阵Λ_n元素的取值为：

${\mathrm{\λ}}_m=\underset{s=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ρ}}_s{e}^{\mathrm{ism\β}},---\left(6\right)$

式中ρ_s=sinc(k·Δr_s),Δr_s=2asin(sβ/2)为第m和m'号阵元间的距离，s=|m-m′|。

d）将误差敏感度函数变换为如下形式：

$T={\mathrm{\η}}^T\mathrm{\Σ\η}=\underset{m=0}{\overset{M/2}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_m{\mathrm{\η}}_m^2{\left|E_m({\mathrm{\θ}}_0,{\mathrm{\φ}}_0)\right|}^2.---\left(7\right)$

（2）基于以上变换，将超指向性波束形成方法表述为一个简单的多约束优化问题：

s.t.F(Ω₀)^Tη=1,η^TΣη≤σ, (8)

F(Ω_j)^Tη≤δ_j(Ω_j∈Ω_SL,j=1,...，N_SL).

其中Ω=(θ,φ)，表示Hadamard积，Ω_SL为选定的旁瓣区域，σ是误差敏感度函数的上界，即稳健性约束值，δ_j是期望旁瓣高度。注意，对于主瓣左右对称的旁瓣区域，δ_j的取值也应是对称的。

（3）上述优化问题是典型的二阶锥规划问题，由成熟的内点算法进行求解，实际中直接利用软件工具包SeDuMi进行计算。

（4）有两种方式合成最终波束图：

第一种将计算得到的实数权值向量η直接代入式（3）即得到最终需要的波束图。

第二种分为以下步骤：

a）由η得到 $\widehat{\mathrm{\η}}={[{\mathrm{\η}}_0,{\mathrm{\η}}_1,...,{\mathrm{\η}}_M]}^T,$ 其中η_M-m=η_m；

b）由式得到

c）由式求出一般化的权值向量w，其中V=[v₀,v₁,...,v_M-1]；

d）将w代入式B(θ,φ)=w^HP(θ,φ)得到最终所需要的波束图。

参照图1。本发明先将波束形成变换到实数域，然后将其构造为简单的多约束优化问题，利用软件工具包进行求解，最后利用求得的权值向量根据具体情况选用不同方式合成超指向性波束图。各个过程的具体实施步骤如下：

（1）参照图2。均匀圆环形阵列包含16个阵元，半径为0.25m。利用该圆环形阵列具有的对称性质，将其波束图、指向性因子和误差敏感度函数变换为实数域的计算形式。

第s个阵元接收到的从方向(θ,φ)入射的单位幅度平面波信号的表达式为：

p_s(θ,φ)=e^{-ikasinθcos(φ-φs)}, （9）其中a为圆环形阵列半径，k=2π/λ，λ表示入射平面波的波长，φ_s=sβ，β=2π/M。转换的具体步骤如下：

a）将阵列流形向量P(θ,φ)=[p₀(θ,φ),p₁(θ,φ),...，p_M-1(θ,φ)^T]变换为新的阵列流形向量 $\hat{E}(\mathrm{\θ},\mathrm{\φ})=[E_0(\mathrm{\θ},\mathrm{\φ}),E_1(\mathrm{\θ},\mathrm{\φ}),...,E_{M-1}(\mathrm{\θ},\mathrm{\φ}){]}^T,$ 相应元素的变换关系为：

$E_m(\mathrm{\θ},\mathrm{\φ})=v_m^{H}P(\mathrm{\θ},\mathrm{\φ})=\frac{1}{\sqrt{M}}\underset{s=0}{\overset{M-1}{\mathrm{\Σ}}}{e}^{-\mathrm{ism\β}}\·p_s(\mathrm{\θ},\mathrm{\φ}),---\left(10\right)$

其中v_m=M^-1/2[1 e^imβ...e^i(M-1)mβ]^T，上标T表示转置，上标H表示共轭转置。新阵列流形向量的元素具有对称性：其中“*”表示复数共轭。

b）将波束图变换为如下形式：

B(θ,φ)=η^TF(θ,φ), （11）

其中η=[η₀,η₁,...,η_M/2]^T为所要求解的实数权值向量，(θ₀,φ₀)为预先设定的波束指向角，设定为(θ₀,φ₀)=(90°,180°)，E=[E₀,E₁,...,E_M/2]^T，ε=[ε₀,ε₁,...,ε_M/2]^T。向量ε元素的取值为：

${\mathrm{\ϵ}}_m=\left\{\begin{array}{cc}1& (m=0,M/2)\\ 2& (m=\mathrm{1,2},...,M/2-1).\end{array}\right.---\left(12\right)$

c）将指向性因子变换为如下形式：

其中Σ=diag{ε₀|E₀|²,ε₁|E₁|²,...,ε_M/2|E_M/2|²}，Λ_n=diag{λ₀,λ₁,...,λ_M/2}。矩阵Λ_n元素的取值为：

${\mathrm{\λ}}_m=\underset{s=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ρ}}_s{e}^{\mathrm{ism\β}},---\left(14\right)$

式中ρ_s=sinc(k·Δr_s),Δr_s=2asin(sβ/2)为第m和m'号阵元间的距离，s=|m-m′|。

d）将误差敏感度函数变换为如下形式：

$T={\mathrm{\η}}^T\mathrm{\Σ\η}=\underset{m=0}{\overset{M/2}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_m{\mathrm{\η}}_m^2{\left|E_m({\mathrm{\θ}}_0,{\mathrm{\φ}}_0)\right|}^2.---\left(15\right)$

（2）基于以上变换，将超指向性波束形成方法表述为一个简单的多约束优化问题：

s.t.F(Ω₀)^Tη=1,η^TΣη≤σ, (16)

F(Ω_j)^Tη≤δ_j(Ω_j∈Ω_SL，j=1，...，N_SL).

其中Ω=(θ,φ)，表示Hadamard积，Ω_SL为选定的旁瓣区域，σ是误差敏感度函数的上界，即稳健性约束值，δ_j是期望旁瓣高度。注意，对于主瓣左右对称的旁瓣区域，δ_j的取值也应是对称的。

参照图3。相关参数为：Ω_SL={(θ,φ)|θ=θ₀,φ∈[0°,φ₀-50°]∪[φ₀+50°,360°]}，σ=10^0/10，δ_j=10^-15/20。

参照图4。相关参数为：Ω_SL={(θ,φ)|θ=θ₀,φ∈[0°,φ₀-24°]∪[φ₀+24°,360°]}，σ=∞，即没有稳健性约束。δ_j的取值在Ω₁={(θ,φ)|θ=θ₀,0°≤φ≤φ₀-24°}中随方位角φ的增大从10^-25/20线性递增到10^-15/20，在Ω₂={(θ,φ)|θ=θ₀,φ₀+24°≤φ≤360°}中随方位角φ的增大从10^-15/20线性递减到10^-25/20。

（3）上述优化问题是典型的二阶锥规划问题，由成熟的内点算法进行求解，实际中直接利用软件工具包SeDuMi进行计算。

（4）有两种方式合成最终波束图：

第一种将计算得到的实数权值向量η直接代入式（3）即得到最终需要的波束图。

第二种分为以下步骤：

a）由η得到 $\widehat{\mathrm{\η}}={[{\mathrm{\η}}_0,{\mathrm{\η}}_1,...,{\mathrm{\η}}_M]}^T,$ 其中η_M-m=η_m；

b）由式得到

c）由式求出一般化的权值向量w，其中V=[v₀,v₁,...,v_M-1]；

d）将w代入式B(θ,φ)=w^HP(θ,φ)得到最终所需要的波束图。

注意，两种方法得到的结果一样，而第一种方法较直观简单，第二种方法更利于实际操作。

参照图3。由本发明采用第二种方式得到的ka=2时的等旁瓣超指向性波束与文献1方法得到的波束是重合的，说明两种方法在这个实施例中得到的结果是一样的。然而本发明方法计算所需要的时间减少了，效率得以提高。具体而言，在此实施例中，本发明方法的运行时间为202.7毫秒，文献1方法的运行时间为326.2毫秒，前者是后者的62.14%。

参照图4。由本发明采用第二种方式得到的ka=2时的非等旁瓣超指向性波束与文献1方法得到的波束也是重合的，说明两种方法在这个实施例中得到的结果也是一样的。在此实施例中，本发明方法的运行时间为149.0毫秒，文献1方法的运行时间为358.7毫秒，前者是后者的41.54%，本发明方法同样比文献1方法的效率高。

以上实施例中，所用计算机的CPU为Intel Core2，主频为2.33GHz，计算时间都由50次运行后平均得到。

参照表1。与常规方法比较，本发明方法在以上两个实施例中得到的超指向性波束具有更大的指向性指数，更窄的半功率波束以及更低的旁瓣级。

表1本发明方法得到超指向性波束和常规波束的性能指标

ka=2，等旁瓣	本发明方法	常规方法
			指向性指数（dB）	10.4	6.5
半功率波束宽度（°）	37.8	65.4
			旁瓣级（dB）	-15.0	-7.9
ka=2，非等旁瓣	本发明方法	常规方法
			指向性指数（dB）	14.3	同上
半功率波束宽度（°）	21.0	同上
			旁瓣级（dB）	-15.0	同上

相比于文献2的方法，本发明方法虽然没有闭式解，但可以在指向性、稳健性和旁瓣级之间进行更好折中，得到更符合实际需要的结果，使用更灵活方便。

Claims (1)

1.一种圆环形阵列实数加权超指向性波束形成方法，其特征在于步骤如下：

步骤1：均匀圆环形阵列第s个阵元接收到的从方向(θ,φ)入射的单位幅度平面波信号为其中，a为圆环形阵列半径，k＝2π/λ，λ表示入射平面波的波长，φ_s＝sβ，β＝2π/M；

所述均匀圆环形阵列包含M个阵元，且M为偶数；

步骤2：将p_s(θ,φ)构成阵列流形向量P(θ,φ)＝[p₀(θ,φ),p₁(θ,φ),…,p_M-1(θ,φ)]^T，然后变换为新的阵列流形向量相应元素的变换关系为：

$E_m(\mathrm{\θ},\mathrm{\φ})=v_m^{H}P(\mathrm{\θ},\mathrm{\φ})=\frac{1}{\sqrt{M}}\underset{s=0}{\overset{M-1}{\mathrm{\Σ}}}{e}^{-ism\mathrm{\β}}\·p_s(\mathrm{\θ},\mathrm{\φ})$

其中v_m＝M^-1/2[1 e^imβ … e^i(M-1)mβ]^T，上标T表示转置，上标H表示共轭转置；新阵列流形向量的元素具有对称性：其中“*”表示复数共轭；

步骤3：将波束图变换为：B(θ,φ)＝η^TF(θ,φ)，其中η＝[η₀,η₁,…,η_M/2]^T为所要求解的实数权值向量，F(θ,φ)＝Re{εοE^*(θ₀,φ₀)οE(θ,φ)}，ο表示Hadamard积，(θ₀,φ₀)为预先设定的波束指向角，E＝[E₀,E₁,…,E_M/2]^T，ε＝[ε₀,ε₁,…,ε_M/2]^T；

所述向量ε元素的取值为：

步骤4：将指向性因子变换为其中

Σ＝diag{ε₀|E₀|²,ε₁|E₁|²,…,ε_M/2|E_M/2|²}，Λ_n＝diag{λ₀,λ₁,…,λ_M/2}；

所述矩阵Λ_n元素为实数，其值为：

式中ρ_h＝sinc(k·Δr_h),Δr_h＝2asin(hβ/2)为第m和m'号阵元间的距离，h＝|m-m′|；

步骤5：将误差敏感度函数变换为

步骤6：将步骤2～步骤5变换的参量代入下式

s.t.F(Ω₀)^Tη＝1，η^TΣη≤σ，

F(Ω_j)^Tη≤δ_j，Ω_j∈Ω_SL,j＝1，…,N_SL

其中Ω＝(θ,φ)，Ω_SL为选定的旁瓣区域，σ是误差敏感度函数的上界，即稳健性约束值，δ_j是期望旁瓣高度，δ_j是以主瓣为中心左右对称取值的；

步骤7：采用两种方式合成最终波束图

1、将计算得到的实数权值向量η直接代入B(θ,φ)＝η^TF(θ,φ)得到最终需要的波束图；

2、由以下步骤：

a)由η得到其中η_M-m＝η_m；

b)由式得到

c)由式求出一般化的权值向量w，其中V＝[v₀,v₁,…,v_M-1]；

d)将w代入式B(θ,φ)＝w^HP(θ,φ)得到最终所需要的波束图。