目前，常见的辐射噪声测量方法以单水听器测量为主^[1-3]。但随着减振降噪技术和新型材料的发展应用，舰船辐射声级逐年降低^[4-5]，而随着海洋背景噪声增加^[6-7]，导致接收信号的信噪比(signal to noise ratio, SNR)越来越小，最终影响辐射噪声测量结果的准确性。

为了解决低信噪比时辐射噪声测量效果差的问题，学者们提出通过阵处理的方法提高测量性能^[8-11]。小孔径阵具有成本低、布放方便等优点，在实际中应用较广，但是其阵增益小，低信噪比下噪声抑制效果差，最终辐射噪声测量精度较低。

为了使用小孔径阵测量低信噪比下的辐射噪声，可以对数据协方差矩阵进行预处理，先提高接收信噪比，再进行阵处理以测量辐射声源级。对于非相关噪声，当快拍数较多时，噪声分量主要位于协方差矩阵的对角线上，则可以采用对角减载类方法进行降噪处理。常见的对角减载方法为均匀对角减载(homogeneous diagonal unloading, HDU)^[12]，即各对角线元素减去相同的减载量，但是由于各阵元接收噪声的功率不相等，均匀对角减载仅能消除部分非相关噪声的影响。针对这种情况，有学者研究了非均匀对角减载(inhomogeneous diagonal unloading, IDU)方法^[13-14]，即各对角线元素的减载量不同。当噪声协方差矩阵为对角阵时，IDU方法可以完全消除噪声分量，但是当快拍数较少时，噪声协方差矩阵中非对角线元素不再为零，IDU方法的降噪性能急剧下降。稳健主成分分析(robust principal component analysis, RPCA)方法可以将由低秩矩阵和稀疏矩阵组成的矩阵分离开^[15-16]，当声源数小于阵元数且噪声相关性较小时，采样数据协方差矩阵就是由低秩的信号协方差矩阵和稀疏的噪声协方差矩阵组成，因此可以利用RPCA方法直接提取出信号协方差矩阵，达到降噪的目的。

本文针对低信噪比下小孔径垂直阵辐射噪声测量问题，利用RPCA方法提取信号协方差矩阵，提高接收信噪比，再利用阵处理方法估计辐射信号的声源级。相比于对角减载类方法，RPCA方法在快拍数较少时也有较好的降噪效果。

1 信号模型与测量参数

测量舰船辐射噪声时，常选择远离航道的水域，则可忽略其他干扰源，仅考虑单目标的情况。假设目标为点声源，一个M元垂直线列阵用来接收辐射声信号，则第i个阵元上的接收信号可以表示为：

$ {x_i}(t) = {h_i}(t) * s(t) + {n_i}(t),i = 1,2, \cdots ,M $

(1)

式中：s(t)为距声源等效声中心1 m处的辐射声信号；h_i(t)为声源到第i个阵元的信道冲激响应函数，以自由场中距声源等效声中心1 m处的信道响应函数归一化；n_i(t)为第i个阵元接收到的环境噪声。对接收信号进行短时傅里叶变换，并用矩阵表示：

$ \mathit{\boldsymbol{X}}(n,f) = \mathit{\boldsymbol{H}}(f)S(n,f) + \mathit{\boldsymbol{N}}(n,f) $

(2)

式中：n表示不同的时域信号段；f表示不同的信号频率；H(f)为频域信道传输函数，以自由场中距声源中心1 m处信道传输函数的模值归一化；X(n, f)、S(n, f)和N(n, f)分别为第n段接收信号、距声中心1 m处辐射信号和环境噪声的离散时间傅里叶变化。为表述简洁，下文所提声源辐射信号即为距声源等效声中心1 m处的辐射声信号。

辐射噪声测量的主要目的是测量舰船辐射噪声的总声源级、声源频带级和声源频谱级。其中，总声源级(source level, SL)定义为距等效声源中心1 m处声强与参考声强的比值，取分贝数，即：

$ {\rm{SL}} = 10\lg \frac{{{{\bar I}_1}}}{{{I_{{\rm{ref}}}}}} = 10\lg \frac{{\bar p_1^2}}{{p_{{\rm{ref}}}^2}} $

式中：I₁和p₁分别表示距等效声中心1 m处的有效声强和有效声压，参考声强为平面波下参考声压p_ref=1 μPa对应的声强，即I_ref=0.67×10^-18 W/m²。声源频带级(source band level, SBL)为声源辐射信号在一定频带内的声强级，通常频带宽度按照1/3倍频程变化。声源频谱级(source spectral level, SSL)为声源辐射信号在1 Hz带宽内的声强级。

针对式(1)的信号模型，假设经过处理已经得到声源辐射信号s(t)，快拍数为L，其有效声压为$\bar{p}=\sqrt{\sum\limits_{i=0}^{L-1}|s(t)|^{2} / L}$，则对应的声源级表示为：

$ {\rm{SL}} = 10\lg \frac{{{{\bar p}^2}}}{{p_{{\rm{ref}}}^2}} = 10\lg \frac{{\frac{1}{L}\sum\limits_{t = 0}^{L - 1} | s(t){|^2}}}{{p_{{\rm{ref}}}^2}} $

(3)

实际中，通常在频域处理宽带信号，利用帕塞瓦尔定理可以得到时域信号和频域信号之间的能量关系为：

$ \sum\limits_{t = 0}^{L - 1} | s(t){|^2} = \frac{1}{L}\sum\limits_{k = 0}^{L - 1} | S(k){|^2} $

(4)

将式(4)代入式(3)，得：

$ {\rm{SL}} = 10\lg \frac{{\sum\limits_{k = 0}^{L - 1} | S(k){|^2}}}{{{L^2}p_{{\rm{ref}}}^2}} $

(5)

信号采样频率为f_s时，傅里叶变换的频域分辨率为f_s/L，则在B=f_s/L带宽内的声源频带级为：

$ \begin{array}{l} {\rm{SBL}} (k) = 10\lg \frac{{|S(k){|^2}}}{{{L^2}p_{{\rm{ ref }}}^2}} = 10\lg \left( {\frac{{|S(k){|^2}}}{{L{f_s}p_{{\rm{ ref }}}^2}}\frac{{{f_s}}}{L}} \right) = \\ {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} 10\lg \frac{{{P_S}(k)}}{{p_{{\rm{ ref }}}^2}} + 10\lg \frac{{{f_s}}}{L} \end{array} $

式中：P_S(k)=|S(k)|²/(Lf_s)表示声源辐射信号的功率谱密度(power spectral density, PSD)。声源频谱级则表示为：

$ {\rm{SSL}}(k) = 10\lg \frac{{{P_S}(k)}}{{p_{{\rm{ref}}}^2}} $

(6)

显然，声源频谱级即为声源辐射信号的PSD在(1 μPa)²参考下的分贝表示，因此当测量的辐射信号单位是μPa时，声源频谱级和声源辐射信号的PSD等价。1/3倍频程内的声源频带级为：

$ {\rm{SB}}{{\rm{L}}_{\frac{{\rm{1}}}{{\rm{3}}}{\rm{OCT}}}}({B_k}) = 10\lg (\sum\limits_{k = {k_l}}^{{k_h}} 1 {0^{{\rm{SSL}}(k)/10}}) + 10\lg \frac{{{f_s}}}{L} $

(7)

式中k_l和k_h分别表示第B_k个1/3倍频程的上、下限频率点。

通过上面推导可知总声源级、声源频带级、声源频谱级和功率谱密度之间的关系，利用辐射信号的功率谱密度即可求出其他量，因此本文后续以估计辐射信号功率谱密度(即声源频谱级)为主。

2 基于稳健主成分分析的辐射信号测量 2.1 辐射信号功率谱密度估计

对于式(2)表示的频域信号模型，假设环境噪声和舰船辐射信号不相关，并且各阵元接收的环境噪声互不相关。恒定束宽波束形成器的归一化加权向量用W(f)表示，则波束输出信号表示为：

$ \begin{array}{*{20}{c}} {Y(n,f) = \mathit{\boldsymbol{W}}{{(f)}^{\rm{H}}}\mathit{\boldsymbol{X}}(n,f) = }\\ {\mathit{\boldsymbol{W}}{{(f)}^{\rm{H}}}\mathit{\boldsymbol{H}}S(n,f) + \mathit{\boldsymbol{W}}{{(f)}^{\rm{H}}}\mathit{\boldsymbol{N}}(n,f)} \end{array} $

为简化，后文公式中省略n和f。波束输出信号在频率f上的PSD表示为：

$ \begin{array}{*{20}{l}} {{P_Y} = 10\lg \left( {\frac{{|Y{|^2}}}{{L{f_s}}}} \right) = 10\lg \left( {\frac{{{\mathit{\boldsymbol{W}}^{\rm{H}}}\mathit{\boldsymbol{RW}}}}{{L{f_s}}}} \right) \approx }\\ {10\lg \left( {\frac{{|{\mathit{\boldsymbol{W}}^{\rm{H}}}\mathit{\boldsymbol{H}}S{|^2}}}{{L{f_s}}}} \right) = {P_S} + 10\lg |{\mathit{\boldsymbol{W}}^{\rm{H}}}\mathit{\boldsymbol{H}}{|^2}} \end{array} $

(8)

式中：R=(XX^H)/N为接收数据的频域协方差矩阵，这里N表示频域快拍数，由时域数据段的数量确定；P_S=10lg(|S|²/Lf_s)表示辐射声信号的功率谱密度。因此，辐射声信号的PSD可以表示为：

$ \begin{array}{*{20}{c}} {{P_S} \approx {P_Y} - 10\lg |{\mathit{\boldsymbol{W}}^{\rm{H}}}\mathit{\boldsymbol{H}}{|^2} = }\\ {10\lg \left( {\frac{{{\mathit{\boldsymbol{W}}^{\rm{H}}}\mathit{\boldsymbol{RW}}}}{{L{f_s}}}} \right) - 10\lg (A{G_s}) + {\rm{TL}}} \end{array} $

(9)

式中：AG_s=|W^HH|²M/(H^HH)，表示阵处理时的信号增益，用以修正阵处理对信号功率的影响；TL=-10lg(H^HH/M)为M个阵元的平均声传播损失，用以补偿声源到整个基阵的声传播损失。

实际中，真实的信道传输函数未知，所以需要估计声场传递函数。常规辐射噪声测量方法将直达波以外的多途信号都当作干扰，通过对直达波方向进行波束形成，以滤除其他方向的多途信号和环境噪声，最后再按照球面扩展补偿接收信号的声传播损失，估计出直达波方向辐射信号的功率谱密度。这种处理方法仅将直达波当作信号，因此声场传输函数可以简化为直达波的阵列流行向量。阵处理时按照恒定束宽波束形成进行加权，声传播损失按照球面扩展近似表示为TL≈20lg r，则辐射信号PSD的估计表示为：

$ {\hat P_S} \buildrel \Delta \over = 10\lg \left( {\frac{{{\mathit{\boldsymbol{W}}^{\rm{H}}}\mathit{\boldsymbol{RW}}}}{{L{f_s}}}} \right) - 10\lg ({\widehat {AG}_s}) + 20\lg r $

(10)

式中：$\hat{P}_{S}$表示辐射噪声信号真实PSD的估计；r表示声源和基阵的距离。远场测量时，归一化恒定束宽加权向量的信号增益$\hat{AG}$_s≈1。

根据式(8)，波束输出信号的PSD既包含接收信号的PSD，也包含了环境噪声的PSD，所以直接通过信号增益和传播损失修正波束输出信号的PSD来计算辐射信号功率谱密度的方法，受环境噪声影响较大，尤其是当接收信噪比较低时，会导致测量的辐射信号PSD大于实际值。

2.2 稳健主成分分析提取接收信号协方差矩阵

对于式(2)表示的频域信号模型，当信号和噪声之间互不相关，理想的数据协方差矩阵表示为：

$ \mathit{\boldsymbol{R}} = {\rm{E}}(\mathit{\boldsymbol{X}}{\mathit{\boldsymbol{X}}^{\rm{H}}}) = {\mathit{\boldsymbol{R}}_{rs}} + {\mathit{\boldsymbol{R}}_n} $

式中：R_rs为基阵接收信号协方差矩阵；R_n为基阵接收环境噪声协方差矩阵。对于远场单目标辐射信号测量，接收信号协方差矩阵的秩为1，并且当接收的环境噪声不相关时，噪声的协方差矩阵为对角阵。

稳健主成分分析是通过优化方法将一个矩阵分解为一个低秩矩阵和一个稀疏矩阵。根据上面对协方差矩阵的分析可知，接收信号协方差矩阵R_rs的秩等于1，且远小于阵元数，则R_rs为低秩矩阵，而噪声协方差矩阵R_n为对角阵，满足稀疏矩阵的要求，即矩阵中大部分元素为0。因此，通过RPCA方法可以将数据协方差矩阵分解为接收信号协方差矩阵和噪声协方差矩阵，其数学表示为：

$ \begin{array}{*{20}{c}} {\mathop {\min }\limits_{{\mathit{\boldsymbol{R}}_{rs}},{\mathit{\boldsymbol{R}}_n}} {\rm{rank}} ({\mathit{\boldsymbol{R}}_{rs}}) + \beta {{\left\| {{\mathit{\boldsymbol{R}}_n}} \right\|}_0}}\\ {{\rm{ s}}{\rm{.t}}{\rm{. }}{\mathit{\boldsymbol{R}}_{rs}} + {\mathit{\boldsymbol{R}}_n} = \mathit{\boldsymbol{R}}} \end{array} $

(11)

式中：rank(R_rs)表示矩阵的秩；‖R_n‖₀表示矩阵的L₀范数；β为正则化参数，用来平衡信号协方差矩阵低秩性和噪声协方差矩阵的稀疏性。

通过求解式(11)可以提取出接收信号协方差矩阵和噪声协方差矩阵，但是由于式(11)是非凸优化问题，求解过程很难实现。可以利用L₁范数作为L₀范数的松弛表示，将非凸优化问题转化为凸优化问题。对于任意矩阵R，矩阵的秩可以用矩阵的核范数表示为‖R‖_*=$\sum\limits_{k} \lambda_{k}$，其中λ_k表示矩阵R的第k个特征值。矩阵的L₀范数可以松弛表示为矩阵中所有元素绝对值的和‖R‖₁=$\sum\limits_{i} \sum\limits_{j}$|R_ij|。另外，由于真实的信号协方差矩阵和噪声协方差矩阵是共轭对称的半正定矩阵，则可以对分解后的矩阵作相应的约束，得到最终的表达式为：

$ \begin{array}{l} \begin{array}{*{20}{l}} {\{ {{\mathit{\boldsymbol{\hat R}}}_{rs}},{{\mathit{\boldsymbol{\hat R}}}_n}\} = {\rm{arg}}\mathop {{\rm{min}}}\limits_{{\mathit{\boldsymbol{R}}_{rs}},{\mathit{\boldsymbol{R}}_n}} {{\left\| {{\mathit{\boldsymbol{R}}_{rs}}} \right\|}_*} + \beta {{\left\| {{\mathit{\boldsymbol{R}}_n}} \right\|}_1}}\\ {{\rm{ }}{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\rm{s}}{\rm{.t}}{\rm{. }}{\mathit{\boldsymbol{R}}_{rs}} + {\mathit{\boldsymbol{R}}_n} = \mathit{\boldsymbol{R}}} \end{array}\\ {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \begin{array}{*{20}{l}} {{\mathit{\boldsymbol{R}}_{rs}} = \mathit{\boldsymbol{R}}_{rs}^{\rm{H}} \ge 0}\\ {{\mathit{\boldsymbol{R}}_n} = \mathit{\boldsymbol{R}}_n^{\rm{H}} \ge 0} \end{array} \end{array} $

(12)

根据文献[17]，取β=1/$\sqrt{M}$，其中M为阵元数。

接收信号协方差矩阵R_rs被提取出后，将其代入式(10)以取代采样数据协方差矩阵R，则可以降低环境噪声对测量结果的影响，提高低信噪比下辐射噪声测量的性能。

3 数值仿真与分析

考虑各阵元上接收信号的声传播损失和接收的环境噪声功率互不相同，这里将接收信号的谱级信噪比定义为：

$ {\rm{SNR}} = 10\lg \frac{{\sigma _s^2\frac{{{\mathit{\boldsymbol{H}}^{\rm{H}}}\mathit{\boldsymbol{H}}}}{M}}}{{\frac{{\sigma _n^2}}{M}\sum\limits_{i = 1}^M {{\rho _{n,ii}}} }} = 10\lg \frac{{\sigma _s^2{\mathit{\boldsymbol{H}}^{\rm{H}}}\mathit{\boldsymbol{H}}}}{{\sigma _n^2\sum\limits_{i = 1}^M {{\rho _{n,ii}}} }} $

式中：σ_s²为辐射声信号的功率谱；σ_n²为环境噪声的功率谱；ρ_{n, ii}表示矩阵ρ_n的第i个对角线元素。

3.1 数据协方差矩阵分解

本节使用11元均匀线列阵，阵元间距为半波长。假设信号为远场平面波，3个非相关信号的来波方向分别为[75°, 95°, 130°]，对应的信号功率分别为[100, 100, 100]，所有阵元接收的平均噪声功率为100，各阵元接收环境噪声的功率是平均功率±3 dB内的均匀分布，则基阵不同方向入射信号的接收信噪比为0 dB。首先使用理想数据协方差矩阵进行仿真，它由图 1所示的信号协方差矩阵和噪声协方差矩阵组成。利用RPCA方法对协方差矩阵进行分解，得到的结果如图 2所示。

比较图 1和图 2发现，RPCA方法分离的信号协方差矩阵和噪声协方差矩阵与真实的基本一致。为了定量分析RPCA方法分离的信号协方差矩阵和真实接收信号协方差矩阵的误差，定义信号协方差矩阵的估计误差为：

$ {E_R} = \frac{{{{\left\| {{{\mathit{\boldsymbol{\hat R}}}_{rs}} - {\mathit{\boldsymbol{R}}_{rs}}} \right\|}_{\rm{F}}}}}{{{{\left\| {{\mathit{\boldsymbol{R}}_{rs}}} \right\|}_{\rm{F}}}}} $

式中：$\hat{\boldsymbol{R}}_{r s}$表示RPCA提取的接收信号协方差矩阵；‖·‖ _F表示矩阵的F范数。图 2中，RPCA提取的信号协方差矩阵的误差为1.87×10^-9，基本可以忽略。

首先分析噪声相关性对数据协方差矩阵分离结果的影响。假设第i和第j个阵元接收环境噪声的相关系数为γ^|i-j|，其中γ为相邻阵元接收的环境噪声相关系数。令γ在0~1变化，快拍数为20，其余仿真条件不变，采用100次独立仿真结果计算平均误差，则信号协方差矩阵的估计误差如图 3所示。作为比较，仿真中还包括IDU方法的估计误差以及未降噪的数据协方差矩阵与信号协方差矩阵的误差。由仿真结果看出，随着噪声相关性增加，信号协方差矩阵的估计误差逐渐增大，并且当环境噪声完全相关时，IDU方法和RPCA方法失效。而对于部分相关的环境噪声(特别是当噪声相关系数小于0.6时)，虽然RPCA方法对信号协方差矩阵的估计误差增大，但是仍能消除部分噪声分量，仍具有一定的降噪效果。

下面分析不同信噪比和快拍数时RPCA方法对采样数据协方差矩阵的分离效果。信号和环境噪声都采用正态分布的随机信号，信噪比在-20~10 dB内变化，采样数据协方差矩阵通过时域采样数据生成，快拍数分别取5、10、15、20和无穷(即理想协方差矩阵)，采用多次独立仿真的结果估计平均性能。不同方法对信号协方差矩阵的估计误差随快拍数和信噪比变化情况如图 4所示，每个结果都是采用100次独立仿真取平均所得。

根据仿真结果可知，RPCA方法和IDU方法对信号协方差矩阵的估计误差随信噪比和快拍数增加而减小。对于理想协方差矩阵，无论信噪比多大，2种方法都可以准确估计出信号协方差矩阵，完全消除噪声的影响。而当数据快拍数有限时，由于噪声协方差矩阵不再是严格的对角阵，使用IDU方法仅能消除对角线上的噪声成分，因此估计的信号协方差矩阵误差较大，与图 4(a)相比，IDU方法对协方差矩阵的降噪效果并不明显。但是，对于RPCA方法，由于不限制噪声协方差矩阵为对角阵，所以不仅能消除对角线上的噪声分量，还能消除部分非对角线上的噪声分量，因此RPCA方法提取的信号协方差矩阵误差更小。

根据以上分析，当信噪比较低时，为了最大化降低环境噪声的影响，减小RPCA方法的估计误差，快拍数越大越好，环境噪声的相关性越小越好。

3.2 RPCA方法对辐射噪声测量影响分析

仿真使用的浅海波导环境和基阵声源布放方式如图 5所示。接收阵为11元均匀线列阵，阵元间距1 m。信道传输函数由AT (acoustic toolbox)工具箱中的KrakenC模型计算，用于生成基阵的接收信号。

为保证不同频率时各阵元接收信号的谱级信噪比基本一致，假设辐射噪声信号为宽带随机信号，如图 6所示，其PSD为125 dB，信号采样频率10 kHz，时域信号长2 s，FFT变换时采用矩形窗，窗长度0.1 s，时域信号窗重叠50%，共计产生39个频域快拍数。环境噪声为高斯白噪声，对应的PSD为90 dB，由于平均声传播损失大约为40 dB，则基阵上各阵元接收信号的谱级信噪比大约为-5 dB。

在400~800 Hz带宽内，使用恒定束宽波束形成，分别采用常规测量方法、基于IDU的测量方法和基于RPCA的测量方法计算辐射信号的PSD。为了量化分析数值仿真中PSD估计结果与真实值之间的误差，用均方根误差(root mean square error, RMSE)表示为：

$ {\rm{RMSE}} = \sqrt {{\rm{E}}{{({{\hat P}_S} - {P_S})}^2}} $

期望可以用多次独立的仿真结果求平均近似计算。使用100次独立仿真计算测量结果的均方根误差，3种测量方法的均方根误差如图 7所示。

根据仿真结果，在400~800 Hz，基于IDU的测量方法也能改善辐射信号测量性能，但是对性能的提升较小，而基于RPCA的测量方法测量误差最小，并且不同频率的测量误差波动也较小。

辐射信号测量受信噪比影响较大，而RPCA方法提取信号协方差矩阵时对频域快拍数也有一定要求，因此分别针对信噪比和频域快拍数进行仿真，对比不同测量方法对声源辐射信号PSD估计的RMSE。由于频域处理是针对单独的频点计算其PSD，因此仿真中只计算750 Hz处的估计结果。

假设环境噪声级不变，可以通过改变辐射信号PSD来调整基阵的平均接收信噪比。令辐射信号PSD在110~140 dB变化，对应的信噪比变化范围约为-20~10 dB。令信号时长在0.55~5.05 s变化，其余条件不变，对应的频域快拍数在10~100变化。分别使用3种测量方法计算辐射信号的PSD，再通过100次仿真计算不同方法测量结果的均方根误差，不同信噪比和快拍数时的测量结果如图 8所示。

根据图 8(a)的仿真结果可知，随着接收信号信噪比增加，辐射信号PSD的估计误差逐渐减小。通常要求测量的辐射信号PSD与真实值的误差不超过3 dB，对于常规测量方法，在此仿真条件下可以有效测量的最低接收信噪比为-7 dB，即辐射信号PSD应大于123 dB，基于IDU的测量方法性能提升不明显，要求接收信噪比大于-9 dB。但是基于RPCA的测量方法，可以测量PSD更小的辐射信号，最低接收信噪比可达到-13 dB。因此，基于RPCA的测量方法可以有效测量低信噪比下辐射信号的功率谱密度。

根据图 8(b)的仿真结果可知，基于RPCA的测量方法性能最好，其测量误差随快拍数增加而减小。常规测量方法在快拍数大于40后，辐射信号PSD的测量误差基本不随信号时间增加而减小，而基于IDU的测量方法，在快拍数为10时，其测量误差和常规方法相当，随着信号时长增加，其测量误差逐渐降低，但仍高于RPCA方法的测量误差。

根据以上仿真结果，基于RPCA的辐射信号PSD估计方法性能提升显著，可以测量更低信噪比下的辐射信号，并且随着快拍数增加，能够准确测量的信噪比门限会进一步降低。

4 结论

1) 理想条件下，使用RPCA方法分离的信号协方差矩阵可以完全消除非相关噪声的影响，而即使快拍数较少时，使用RPCA方法也能消除部分非相关噪声。

2) 相比于直接使用数据协方差矩阵的辐射噪声测量方法，基于RPCA的测量方法受环境噪声影响更小，相同接收信噪比时测量精度更高，相同测量精度要求下，所要求的接收信噪比更低。

3) 本文研究了小孔径垂直阵辐射噪声测量方法，但RPCA方法对于基阵孔径并无具体要求，只要基阵接收的环境噪声相关性较小，即可应用RPCA方法进行降噪处理。

未来工作中还需进一步研究噪声相关系数较高时，如何利用RPCA方法进行噪声分离。