[发明专利]一种声矢量圆阵稳健宽带MVDR方位估计方法

权利要求书

1.一种声矢量圆阵稳健宽带MVDR方位估计方法，其特征在于：

(1)将二维矢量圆阵接收的声压数据x_p以及振速数据x_vx和x_vy，进行子频带分解，分别生成声压、v_x振速和v_y振速的宽带频域信号矩阵和得到矢量阵的宽带频域快拍数据矩阵X_j(f_k)，

$\left\{\begin{array}{c}{X}_j^p\left(f_k\right)=A^p\left(f_k\right)S_j\left(f_k\right)+N_j^p\left(f_k\right)\\ X_j^{vx}\left(f_k\right)=A^{vx}\left(f_k\right)S_j\left(f_k\right)+N_j^{vx}\left(f_k\right)\\ X_j^{vy}\left(f_k\right)=A^{vy}\left(f_k\right)S_j\left(f_k\right)+N_j^{vy}\left(f_k\right)\end{array}\right.$

其中，S_j(f_k)为源信号的频域快拍数据矩阵；和分别为声压、v_x振速及v_y振速信号的频域快拍数据矩阵；和分别为声压、v_x振速及v_y振速通道噪声信号的频域快拍数据矩阵；A^p(f_k)、A^vx(f_k)和A^vy(fk)分别为声压、v_x振速及v_y振速通道的阵列流型矢量；j＝1,2,…,J，k＝l,l+1,…,h-1,h，f_l为下限频率，f_h为上限频率，下标“_l、_h和_k”代表频率编号，

得到矢量阵第j个频域快拍数据矩阵X_j(f_k)为：

$X_j\left(f_k\right)=\left[\begin{array}{c}{X}_j^p\left(f_k\right)\\ X_j^{vx}\left(f_k\right)\\ X_j^{vy}\left(f_k\right)\end{array}\right];$

(2)生成各个子频带上的矢量阵互谱矩阵R(f_k)，

$R\left(f_k\right)=\frac{1}{J}\underset{j=1}{\overset{J}{\Σ}}{X}_j\left(f_k\right)X_j^H\left(f_k\right);$

(3)设置步长，实施方位角扫描，对任意一个方位角在子频带f_k上构造矢量阵导向矢量

$\hat{a}\left(f_k\right|\widehat{\mathrm{\θ}})=\left[\begin{array}{c}{\hat{a}}_p\left(f_k\right|\widehat{\mathrm{\θ}})\\ {\hat{a}}_{vx}\left(f_k\right|\widehat{\mathrm{\θ}})\\ {\hat{a}}_{vy}\left(f_k\right|\widehat{\mathrm{\θ}})\end{array}\right]$

其中，

${\hat{a}}_p\left(f_k\right|\widehat{\mathrm{\θ}})=\left[\begin{array}{c}{\hat{a}}_1\left(f_k\right|\widehat{\mathrm{\θ}})\\ {\hat{a}}_2\left(f_k\right|\widehat{\mathrm{\θ}})\\ \begin{array}{c}.\\ .\\ .\end{array}\\ {\hat{a}}_N\left(f_k\right|\widehat{\mathrm{\θ}})\end{array}\right]=\left[\begin{array}{c}\exp \left(-i2{\mathrm{\πf}}_k\left(x_1\cos \widehat{\mathrm{\θ}}+y_1\sin \widehat{\mathrm{\θ}}\right)/c\right)\\ \exp \left(-i2{\mathrm{\πf}}_k\left(x_2\cos \widehat{\mathrm{\θ}}+y_2\sin \widehat{\mathrm{\θ}}\right)/c\right)\\ \begin{array}{c}.\\ .\\ .\end{array}\\ \exp \left(-i2{\mathrm{\πf}}_k\left(x_N\cos \widehat{\mathrm{\θ}}+y_N\sin \widehat{\mathrm{\θ}}\right)/c\right)\end{array}\right];$

c为声速，i表示复数，i²＝-1；

(4)选取f₀为聚焦参考频率点，在相同的方位角上构造矢量阵聚焦导向矢量

$\hat{a}\left(f_0\right|\widehat{\mathrm{\θ}})=\left[\begin{array}{c}{\hat{a}}_p\left(f_0\right|\widehat{\mathrm{\θ}})\\ {\hat{a}}_{vx}\left(f_0\right|\widehat{\mathrm{\θ}})\\ {\hat{a}}_{vy}\left(f_0\right|\widehat{\mathrm{\θ}})\end{array}\right],$

其中，

${\hat{a}}_{vx}\left(f_0\right|\widehat{\mathrm{\θ}})=cos\widehat{\mathrm{\θ}}\·{\hat{a}}_p\left(f_0\right|\widehat{\mathrm{\θ}});{\hat{a}}_{vy}\left(f_0\right|\widehat{\mathrm{\θ}})=\sin \widehat{\mathrm{\θ}}\·{\hat{a}}_p\left(f_0\right|\widehat{\mathrm{\θ}});$

(5)采用相干信号子空间CSS聚焦变换方法，得到f_k上的矢量阵聚焦变换矩阵

$T\left(f_k|\widehat{\mathrm{\θ}}\right)=\left[\begin{array}{ccc}{T}_p\left(f_k|\widehat{\mathrm{\θ}}\right)& 0_N& 0_N\\ 0_N& T_{vx}\left(f_k|\widehat{\mathrm{\θ}}\right)& 0_N\\ 0_N& 0_N& T_{vy}\left(f_k|\widehat{\mathrm{\θ}}\right)\end{array}\right],$

其中，

0_N表示N×N维零矩阵；

(6)对(2)和(5)进行导向矢量聚焦变换后，得到f_k上的协方差矩阵进一步将总共K个频带的进行累加后得到宽带聚焦协方差矩阵

$\hat{R}\left(f_k|\widehat{\mathrm{\θ}}\right)=T\left(f_k|\widehat{\mathrm{\θ}}\right)R\left(f_k\right)T^H\left(f_k|\widehat{\mathrm{\θ}}\right),$

$\hat{R}\left(\widehat{\mathrm{\θ}}\right)=\frac{1}{K}\underset{k=l}{\overset{h}{\Σ}}\hat{R}\left(f_k\right|\widehat{\mathrm{\θ}});$

(7)对实际存在失配误差的矢量阵导向矢量施加稳健性约束优化条件，得到最优权矢量

至分别为矢量中的1至2N号元素；

理想情况下的矢量阵导向矢量为与实际的矢量阵导向矢量之间存在误差量e，则有针对存在的失配误差e进行约束||e||≤ε，ε为约束参数，“||·||”表示取矢量的2范数，则属于下面的集合：

$a\left(\mathrm{\ϵ}\right)=\left\{a\left(f_0\right|\widehat{\mathrm{\θ}})\right|a\left(f_0\right|\widehat{\mathrm{\θ}})=\hat{a}\left(f_0\right|\widehat{\mathrm{\θ}})+e,\left|\right|e\left|\right|\≤\mathrm{\ϵ}\}$

稳健性约束优化问题可表示为

$\left\{\begin{array}{cc}\underset{w\left(\widehat{\mathrm{\θ}}\right)}{\min }& w^H\left(\widehat{\mathrm{\θ}}\right)\hat{R}\left(\widehat{\mathrm{\θ}}\right)w\left(\widehat{\mathrm{\θ}}\right)\\ s.t.& \left|w^H\left(\widehat{\mathrm{\θ}}\right)a\left(f_0\right|\widehat{\mathrm{\θ}})\right|\≥1foralla\left(f_0|\widehat{\mathrm{\θ}}\right)\∈a\left(\mathrm{\ϵ}\right)\end{array}\right.$

其中，“|·|”表示取模，“∈”表示属于，为待优化得到的最优权矢量，“s.t.”为Subject to的缩写，表示约束条件，式中大括号的第一行表示目标函数，第二行表示约束条件；

通过化简，可得上式的等价形式为：

$\left\{\begin{array}{cc}\underset{w\left(\widehat{\mathrm{\θ}}\right)}{\min }& w^H\left(\widehat{\mathrm{\θ}}\right)\hat{R}\left(\widehat{\mathrm{\θ}}\right)w\left(\widehat{\mathrm{\θ}}\right)\\ s.t.& w^H\left(\widehat{\mathrm{\θ}}\right)\hat{a}\left(f_0\right|\widehat{\mathrm{\θ}})\≥\mathrm{\ϵ}\left|\right|w\left(\widehat{\mathrm{\θ}}\right)\left|\right|+1,\mathrm{Im}\left\{w^H\left(\widehat{\mathrm{\θ}}\right)\hat{a}\left(f_0|\widehat{\mathrm{\θ}}\right)\right\}=0\end{array}\right.$

其中，“Im”表示取虚部；

对进行Cholesky分解得到分解因子引进一个非负标量并构造一个新的约束形式则有：

$\left\{\begin{array}{cc}\underset{\mathrm{\τ},w\left(\widehat{\mathrm{\θ}}\right)}{\min }& \mathrm{\τ}\\ s.t.& \mathrm{\ϵ}\left|\right|w\left(\widehat{\mathrm{\θ}}\right)\left|\right|\≤w^H\left(\widehat{\mathrm{\θ}}\right)\hat{a}\left(f_0|\widehat{\mathrm{\θ}}\right)-1,\mathrm{Im}\left\{w^H\left(\widehat{\mathrm{\θ}}\right)\hat{a}\left(f_0|\widehat{\mathrm{\θ}}\right)\right\}=0,\left|\right|Uw\left(\widehat{\mathrm{\θ}}\right)\left|\right|\≤\mathrm{\τ}\end{array}\right.$

其实值形式为：

其中，

“Re”表示取实部，

最终写成二阶锥规划问题的求解形式：

其中，

0表示元素全为0的矢量，I表示元素全为1的矢量，R表示实数域，R右上角表示维数，SOC₁^2N+1表示第一个2N+1维的二阶锥，SOC₂^2N+1表示第二个2N+1维的二阶锥，{0}表示零锥；

在扫描角度上的最优化权矢量表示为：

(8)将代入目标函数中，得在优化后的阵列平均输出功率

(9)重复(3)至(8)的步骤直至完成全部的空间角度扫描，绘制声矢量圆阵稳健宽带空间谱图，通过空间谱的谱峰位置确定声源来波方向。