0 引　言

波达方向（DOA）估计为阵列信号处理领域的一个研究热点，其在水声通信中也有着广泛应用。矢量水听器^[1]与传统的标量水听器不同，前者可额外获取到振速通道信息，因此使用矢量水听器进行DOA估计正在逐步成为一种新兴水下声呐技术。

传统的DOA估计方法，如快速求根多重信号分类（ROOTMUSIC）算法^[2]等，在水下低信噪比环境中，鲁棒性较差。而压缩感知技术（Com-pressive Sesing，CS）^[3]很好解决了这个问题，但压缩感知方法可能会出现离网问题（off-grid），这严重影响了DOA估计算法的性能。基于原子范数最小化的无网格压缩感知^[4]进一步推进了DOA估计的发展，它利用原子范数具有半正定规划性质（SDP），可将原子范数最小化问题（A-NM）转化为SDP问题，并用凸优化工具箱求解。但SDP问题属于非光滑凸优化问题，利用CVX求解该问题所要消耗的时间会随着矢量水听器阵列规模的增加，逐渐增大。

为解决上述问题，本文提出基于交替方向乘子法（ADMM）^[5]的矢量水听器DOA估计方案，该方案将交替方向乘子法应用到矢量水听器阵列的DOA估计中，考虑到海洋环境噪声，本文使用原子范数去噪方法（AST）来估计线谱参数^[6]，将原子范数最小化问题转化为SDP问题，使用ADMM对SDP问题进行求解，最后使用对偶多项式估计角度。通过仿真对比实验，结果表明，ADMM既可保证算法收敛性又可提高算法的效率。

1 系统模型

如图1所示，在远场条件下，假设三维空间内有$ k (k = 1,2, \ldots K) $个窄带水声信号入射到$ M $个矢量水听器组成的均匀线阵上。矢量水听器可以同时接收一个声压通道信息$ p $和3个振速通道信息$ v $，$ v $在三维空间中的投影分量可以表示为$ ({v_x},{v_y},{v_z}) $。

在远场条件下，参考阵元位于坐标原点，$ {\Theta _k} = \left( {{\theta _k},{\phi _k}} \right) $为第$ k $个水声信号的二维到达角，$ {\theta _k} $为水声信号的水平方位角，$ {\phi _k} $为俯仰角。假设$ \phi = \displaystyle\frac{\text{π} }{2} $，矢量水听器阵列接收到的声压$ p $和2个振速分量$ {v_x} $和$ {v_y} $，记为:

$ \left\{ \begin{gathered} p = s ，\\ {v_x} = \cos \theta \cdot s，\\ {v_y} = \sin \theta \cdot s 。\\ \end{gathered} \right. $

(1)

式中，$ s $为水声信号。

矢量水听器的相位延迟向量表示为：

$ {a_v}(\theta ) ={\left[{e^{ j2\text{π} \frac{{0d\cos \theta }}{\lambda }}},{e^{ j2\text{π} \frac{{1d\cos \theta }}{\lambda }}}, \ldots ,{e^{ j2\text{π} \frac{{\left( {M - 1} \right)d\cos \theta }}{\lambda }}}\right]^{\text T}}。$

(2)

式中：$ \lambda $为水声信号波长，$ d $为阵元间距，$ {\mathrm{T}} $为转置。

$ M $个矢量阵元的接收水声信号为：

$ Y = AS + N 。$

(3)

式中：$ Y $为$ 3 M \times 1 $维矢量阵输出向量；$ {\boldsymbol A} = [a\left( {{\theta _1}} \right), a\left( {{\theta _2}} \right), \ldots , a\left( {{\theta _K}} \right)] $为$ 3 M\times K $维阵列流型矩阵；$ a(\theta ) = {a_v}(\theta ) \otimes h $为矢量水声器的导向向量；$ h = {[1,\cos \theta , \sin \theta ]^{\text T}} $为单个矢量水听器的方向向量；$ \otimes $为Kronecker 积；$ S\left( t \right) = [{s_1}, {s_2}, \ldots ,{s_K}]^{\mathrm{T}} $为$ K \times 1 $维水声信号；$ N(t) = [{n_0}, {n_1}, \ldots , {n_{M - 1}}]^{\mathrm{T}} $为$ 3 M\times1 $维水声噪声向量。

2 本文方案 2.1 原子范数

在理想无噪声的单快拍模型下，矢量阵列接收水声信号的模型如下：

$ X = \left[ \begin{gathered} \sum\limits_{k = 1}^K {{s_k} \cdot {e^{ \displaystyle j2\text{π} \frac{{0d\cos {\theta _k}}}{\lambda }}} \cdot {h_k}} \\ \sum\limits_{k = 1}^K {{s_k} \cdot {e^{ \displaystyle j2\text{π} \frac{{1d\cos {\theta _k}}}{\lambda }}} \cdot {h_k}} \\ \cdots \\ \sum\limits_{k = 1}^K {{s_k} \cdot {e^{ \displaystyle j2\text{π} \frac{{\left( {M - 1} \right)d\cos {\theta _k}}}{\lambda }}} \cdot {h_k}} \\ \end{gathered} \right] 。$

(4)

由文献[7]，类似的可将$ {s_k} $表示为 $ \left| {{s_k}} \right|{e^{j{\varphi _k}}} $的形式，$ \left| {{s_k}} \right| $为水声信号的模值，$ {\varphi _{{k}}} = \left[ {0,2{\text π} } \right) $为水声信号的初始相位。

$ X(t) $中的第$ m $个矢量水听器的单快拍接收数据模型可表示：

$ \begin{split} {X_m}=\ &\sum\limits_{k = 1}^K {\left| {{s_k}} \right| \cdot {e^{j\left[ \displaystyle {2\text{π} \frac{{(m - 1)d\cos {\theta _k}}}{\lambda } - {\varphi _k}} \right]}} \cdot {h_k}} = \\ &\sum\limits_{k = 1}^K {\left| {{s_k}} \right|\cdot a({\theta _k},{\varphi_{{k}}})} 。\end{split} $

(5)

这里，可将$ X(t) $看作原子集元素的稀疏组合，由此定义原子集合$ { A} $，它为一个无限字典^[8]:

$ { A} = \left\{ {a(\theta ,\varphi ):\theta \in \left[ {0,2{\text π} } \right),\varphi \in \left[ {0,2{\text π} } \right)} \right\} 。$

(6)

给出该原子集合上原子范数的定义：

$ {\left\| X \right\|_{ A}} = \mathop {\inf }\limits_{a({\theta _k},{\varphi _k}) \in { A}} \left\{ {\sum\limits_{k = 1}^K {\left| {{s_k}} \right|:X = \sum\limits_{k = 1}^K {\left| {{s_k}} \right|a({\theta _k},{\varphi _k})} } } \right\}。$

(7)

任何范数都有与之对应的对偶函数，原子范数的对偶范数等价于原子集合$ \rm{\mathit{A}} $的支撑集函数，$ {\left\| X \right\|_A} $的对偶范数为：

$ \left\| Q \right\|_A^ * = \mathop {\sup }\limits_{a(\theta ,\varphi ) \in {\rm A}} {Re} \left( {\left\langle {Q,a({\theta _k},{\varphi _k})} \right\rangle } \right) 。$

(8)

式中，$ {Re(}\langle ·\rangle \text{)} $为实内积。

2.2 AST

考虑到实际噪声，本文使用AST^[6]来估计线谱参数。原子范数软阈值的定义：

$ \mathop {\min imize}\limits_X \left\| {X - Y} \right\|_2^2 + 2\tau {\left\| X \right\|_{\rm A}}。$

(9)

式中：$ {\Vert ·\Vert }_{2} $为2范数；$ \tau $为一个与噪声模型有关，用于正则化的参数，若$ n\sim N(0,{\sigma ^2}) $，则：

$ \tau = \sigma (1 + \frac{1}{{\log M}})\sqrt {M\log M + M\log (4\text{π} \log M)} 。$

(10)

从本质上讲，求解原子范数最小化问题与求解其对偶范数最大化问题是等价的，AST的对偶问题：

$ \begin{gathered} \mathop {\max imize}\limits_X = \left\| Y \right\|_2^2 - \left\| {Y - Q} \right\|_2^2，\\ {\mathrm{subject}}\;{\mathrm{to}}\quad \left\| Q \right\|_{\rm A}^ * \leqslant 1 。\\ \end{gathered} $

(11)

2.3 SDP与ADMM

原子范数$ {\left\| X \right\|_A} $近似等于下面问题的最优解，即AST的SDP^[9]表达式为:

$ \begin{split} & \mathop {\min imize}\limits_{X,t,u} \left\| {X - Y} \right\|_2^2 + \tau \left( {t + {\omega ^{\rm T}}u} \right) ，\\ & {\mathrm{subject}}\;{\mathrm{to}}\quad \left[ {\begin{array}{*{20}{c}} {T(u)}&X \\ {{X^H}}&t \end{array}} \right] \geqslant 0 。\end{split}$

(12)

AST是通过在式（12）中选择$ w = 2{e_0} $来获得的，其中$ {e_0} $为第一个元素为1，其他都为0的向量。另外，$ {\boldsymbol u} \in {{C}^{3M \times 1}} $，$ {\boldsymbol {T(u)}} \in {{C}^{3M \times 3M}} $为Toe-plitz矩阵。

${\boldsymbol{T(u)}} = \left[ {\begin{array}{*{20}{c}} {{u_1}}&{{u_2}}& \ldots &{{u_{3M}}} \\ {u_2^{\mathrm{H}}}&{{u_1}}& \ldots &{{u_{3M - 1}}} \\ \vdots & \vdots & \ddots & \vdots \\ {u_{3M}^{\mathrm{H}}}&{u_{3M - 1}^{\mathrm{H}}}& \ldots &{{u_1}} \end{array}} \right] 。$

(13)

式中：$ {u_j} $为$ u $的第$ j $个值；$ {\left| \cdot \right|^{\mathrm{H}}} $为共轭转置。

原子范数最小化问题是个NP（Non-Deter-ministic Polynormial）问题，很难直接求解，但原子范数具有SDP，可将式（9）转化为半正定规化问题（12）。SDP问题属于非光滑凸优化问题，随着矢量水听器阵列规模的增大，SDP问题所涉及的变量个数也随之增多，使得要优化的正定矩阵也随之增大，计算复杂度上升。而且求解SDP问题通常需使用凸优化工具箱（如CVX），CVX一般采用内点法进行迭代求解。当矢量水听器阵列规模的增大时，变量规模增大，内点法求解的迭代次数也会增加，导致计算时间增长。ADMM就没有上述问题，ADMM通过将一个大问题转化为多个小问题交替求解，这种方法即可保证算法收敛性，又通过分解子问题降低了计算量，提高算法效率。

将SDP优化问题为增广拉格朗日函数^[10]形式：

$ \begin{split} &{L_\rho }(t,u,X,Z,\Lambda ) = \left\| {X - Y} \right\|_2^2 + \tau \left( {t + {\omega ^{\rm T}}u} \right)+\\ & 2\left\langle {\Lambda ,Z - \underbrace {\left[ {\begin{array}{*{20}{c}} {T(u)}&X \\ {{X^H}}&t \end{array}} \right]}_T} \right\rangle + \rho \left\| {Z - \left[ {\begin{array}{*{20}{c}} {T(u)}&X \\ {{X^H}}&t \end{array}} \right]} \right\|_F^2 。\end{split} $

(14)

式中：$ \langle ·\rangle $为内积，惩罚项参数$ \rho > 0 $，惩罚项参数影响迭代效率；$ \Lambda \geqslant 0 $为拉格朗日乘子，$ \Lambda $和$ Z $中块的划分与T中的块划分相匹配。

第$ l $次迭代，对各个变量的更新如下：

$ ({t}^{l+1},{u}^{l+1},{X}^{l+1})\leftarrow \mathrm{arg}\underset{t,u,x}{\mathrm{min}}{L}_{\rho }（t,u,X,{Z}^{l},{\Lambda }^{l}），$

(15)

$ {Z}^{l+1}\leftarrow \mathrm{arg}\underset{Z\geqslant 0}{\mathrm{min}}{L}_{\rho }（{t}^{l+1},{u}^{l+1},{X}^{l+1},Z,{\Lambda }^{l}），$

(16)

$ {\Lambda ^{l + 1}} \leftarrow {\Lambda ^l} + \rho \left( {{Z^{l + 1}} - {T^{l + 1}}} \right) 。$

(17)

在更新变量的时候，AD-MM采用交替迭代的方式，通过不断递归调整参数，直到达到收敛状态。在文献[10]中，可找到以上模型中主要参数的闭式解法。

2.4 对偶多项式恢复角度

如果强对偶性成立，则原始解$ \widehat X $和对偶解$ \widehat Q $满足：

$ Y = \widehat X + \tau \widehat Q ，$

(18)

对偶解$ \widehat Q $的对偶多项式^[8]为：

$ \widehat q(\theta ,\varphi ) = \left\langle {\widehat Q,a(\theta ,\varphi )} \right\rangle ，$

(19)

通过寻找$ \left| {\widehat q(\theta ,\varphi )} \right| $的峰值可得到目标角度：

$ \left| {\widehat q(\theta ,\varphi )} \right| = 1 。$

(20)

图2为使用对偶多项式进行DOA估计的结果：水声信号中$ \cos (\theta ) $的实际位置显示为实线，虚线为求解式（20）得到的对偶多项式模值，如图2可知，当对偶多项式模值=1的时候，可有效估计目标角度。其中，矢量水听器阵元数为32，信噪比为10 dB。

3 仿真实验

在仿真对比实验中，本文使用ADMM在均方根误差（Root Mean Square Error，RMSE）、运行时间等方面与ROOTMUSIC和CVX算法进行比较。

实验设置：在远场条件下，假设有$ K = 2 $个窄带水声信号入射到$ M = [8,16,32,64,128] $个矢量水听器组成的均匀线阵。其中，阵元分布均满足阵元间距为$ d = {\lambda }/{2} $，即水声信号的半波长，水声信号信噪比为$ SNR = [0,5,10,15,20] $，快拍数$ L = 1 $，所有仿真实验重复100次，并且对结果取平均值。

RMSE的计算方式为：

$ {RMSE}=\sqrt{\frac{1}{m}{\displaystyle \sum _{i=1}^{m}({\widehat{\theta }}_{i}-{\theta }_{i}{)}^{2}}} 。$

(21)

式中：$ \theta $为真实波达角度；$ \hat \theta $为估计的波达角度；$ m $为需估计的数据个数。均方根误差可很好反映出算法的精度和稳定性。

图3为 ROOTMUSIC、CVX和ADMM这3种算法在不同信噪比条件下的均方根误差。可看出，由于存在水下噪声，3 种方法均受到不同程度的影响。其中，传统子空间分解方法所受影响最大，ROOTMUSIC算法的DOA估计效果较差。而CVX和ADMM算法均方根误差总体优于ROOTMUSIC 算法，且随着信噪比降低，CVX和ADMM算法性能一直保持稳定。在高信噪比情况下，3 种算法均具备较好的分辨能力。总体来说，CVX、ADMM在广泛的信噪比范围下具有比其他算法更加稳定的算法精度与分辨能力。

图4为ROOTMUSIC、CVX和ADMM 等3 种算法在不同矢量阵元数条件下的运行时间。可以看出，随着矢量水听器阵列规模逐渐增加，CVX求解问题的速度越来越慢，而ROOTMUS-IC、ADMM几乎不受影响。

综上，ADMM在保证算法估计性能的同时，算法效率也有着明显优势。

4 结　语

为了提高无网格压缩感知的矢量水听器DOA估计的精度和速度，本文将ADMM应用到矢量水听器阵列的DOA估计中，考虑到海洋环境噪声，本文使用AST来估计线谱参数，将原子范数最小化问题转化为SDP问题，使用ADMM对SDP问题进行求解，最后使用对偶多项式估计角度。通过仿真对比实验，结果表明，ADMM在保证DOA估计模型收敛性的同时，提高了效率。