CN105181121A - 采用加权迭代等效源法的高精度近场声全息算法 - Google Patents

Abstract

本发明公开了一种采用加权迭代等效源法的高精度近场声全息算法，其特征是在声源近场辐射区域内布置全息面H并测量H上的声压P_H；在目标重建面T远离H的一侧布置等效源面S_e，并在S_e上布置等效源；利用等效源与H之间的声压传递矩阵建立P_H与各等效源间的关系；采用带有后验加权范数约束惩罚项的新型迭代正则化算法求解等效源源强Q，再用求出的Q以及等效源与目标重建面T之间的传递矩阵计算T上的声场数据。本发明方法通过带有后验加权范数约束惩罚项的新型迭代正则化算法实现等效源源强的精确求解，避免了Tikhonov正则化过程中2范数惩罚项导致的源强能量泄漏，与基于常规等效源的近场声全息相比，本发明方法获得的计算结果更加精确。

Description

采用加权迭代等效源法的高精度近场声全息算法

技术领域

本发明涉及噪声源识别与分析领域中的高精度、高分辨率近场声全息算法。

背景技术

近场声全息是一种具有强大的噪声源识别、定位及声场可视化功能的声学前沿技术，在机械故障诊断、噪声治理工程中有着广泛的应用前景。经过多年研究，目前已经发展出适用于平面、柱面、球面等可分离坐标曲面的空间Fourier变换法和统计最优算法，以及适用于任意形状声源的边界元法和等效源法(又被称为波叠加法)等多种近场声全息算法。其中等效源法由于其实现过程简单、几何适应性强且克服了Fourier方法相关的窗效应、卷绕误差的影响，具有较高的计算精度等，在工程应用中具有更大的潜力。

但基于等效源法的近场声全息在计算精度和分辨率方面还能进一步提高。在等效源法中，等效源源强的求解是方法的关键所在，只有求出的等效源源强与实际情况相符才能保证对声场的高精度重建与分析。然而等效源源强的求解过程属于数学意义上的逆问题，具有求解的不适定性，微小的输入误差可能会导致解的巨大变化，为此在等效源源强的求解过程中通常采用正则化方法稳定求解过程。目前常规等效源法的源强求解过程中，使用最多的是在最小2范数准则下的Tikhonov正则化方法，然而由于最小2范数约束是一个光滑性约束条件，它倾向于将源强能量分散到所有的等效源中，从而获得一个光滑解。这样往往导致等效源源强能量的泄漏，降低源强峰值，造成近场声全息计算结果分辨率与计算精度的损失。

发明内容

本发明是为避免上述现有技术所存在的不足之处，提供一种采用加权迭代等效源法的高精度近场声全息算法，通过带有后验加权范数约束惩罚项的迭代正则化算法过程，替代现有等效源法近场声全息中采用的Tikhonov正则化过程，实现等效源源强的精确求解，避免Tikhonov正则化过程中2范数惩罚项导致的源强能量泄漏，使获得的等效源源强峰值更加突出，从而获得具有更高分辨率和更高重建精度的近场声全息计算结果。

本发明解决技术问题采用的技术方案是：

本发明采用加权迭代等效源法的高精度近场声全息算法的特点是按如下步骤进行：

步骤a、在声源近场辐射区域内布置全息面并测量全息面上的声压

在由声源产生的声场中，于近场区域内布置全息面H，所述全息面H为任意形状曲面，在所述全息面H上分布各测量网格点，利用传感器测量各测量网格点处的声压P_H；

步骤b、目标重建面T位于声源与全息面H之间，在目标重建面T远离全息面H的一侧布置等效源面S_e，所述等效源面S_e与目标重建面T之间的垂直距离为d_h，布置在所述等效源面S_e上的等效源的个数不大于全息面H上测量网格点的个数；所述等效源为单极子、偶极子或四极子，所述等效源为点源、面源或体源；

步骤c、建立各等效源与测量获得的全息面H上的声压P_H之间的关系为：

P_H＝G_HpQ

其中，G_Hp为等效源面S_e上各等效源与全息面H上各测量网格点之间的声压传递矩阵，Q为等效源源强向量；

步骤d、按如下过程求解等效源源强向量Q：

第一步：采用常规等效源法求解获得等效源源强向量作为等效源源强向量初始解Q⁰，令初始迭代次数为1，设定最大迭代次数为J；

第二步：利用向量Q^j-1作为后验加权系数构建加权矩阵W，W＝diag(Q^j-1)，diag(Q^j-1)是利用向量Q^j-1作为主对角元素构建的对角阵，所述向量Q^j-1为第j-1次迭代计算获得的等效源源强向量，对于第1次迭代计算过程，向量Q^j-1即为Q⁰；

第三步：利用所述加权矩阵W的逆矩阵W^-1构造向量Q^j的加权范数为||W^-1Q^j||；则在加权范数最小化准则意义下，通过如下极小化过程获得等效源源强向量Q^j的正则解：

min{||G_HpQ^j-P_H||²+λ²||W^-1Q^j||²}

通过预条件处理令W^-1Q^j≡X^j，将求解Q^j的加权范数正则化过程化为求解中间变量X^j标准的Tikhonov正则化过程：

min{||G_HpWX^j-P_H||²+λ²||X^j||²}

采用奇异值分解法求解中间变量X^j的正则解，正则化参数λ的选择采用广义交叉验证法或L曲线法来实现；经过第j次迭代获得等效源源强向量Q^j为：Q^j＝WX^j；

第五步：若j≤J则将迭代次数赋值为j+1，并转入第二步继续进行迭代；若j＞J则终止迭代，并令Q^j为等效源源强向量Q_f；

步骤e、利用P_T＝G_TpQ_f计算获得目标重建面T上的声压P_T；利用V_T＝G_TvQ_f计算获得目标重建面T上的法向振速V_T，其中，G_Tp为各等效源与目标重建面之间的声压传递矩阵，G_Tv为各等效源与目标重建面T之间的振速传递矩阵。

本发明采用加权迭代等效源法的高精度近场声全息算法的特点也在于：所述d_h的取值为0.02m-0.1m。

本发明采用加权迭代等效源法的高精度近场声全息算法的特点也在于：所述最大迭代次数J不超过10次。

与已有技术相比，本发明有益效果体现在：

1、本发明方法中通过带有后验加权范数约束惩罚项的新型迭代正则化算法实现等效源源强的精确求解，避免了Tikhonov正则化过程中2范数惩罚项导致的源强能量泄漏，因此与基于常规等效源的近场声全息相比，本发明方法获得的等效源源强更加精确。

2、本发明方法在整个适宜于近场声全息运用的频段内，其计算精度均明显优于基于常规等效源法的近场声全息。

3、本发明方法显著提高了近场声全息重建结果的分辨率，增强了噪声源识别与定位能力。

4、本发明方法具备优良的鲁棒性和抗干扰能力，在低信噪比条件下计算精度明显优于常规等效源法，性能稳定。

附图说明

图1为本发明采用的等效源计算模型；

图2为单极子声源频率为500Hz时，本发明方法和常规等效源法求解所得单极子声源源强与实际源强比较；

图3为单极子声源频率为1000Hz时，本发明方法和常规等效源法求解所得单极子声源源强与实际源强比较；

图4为单极子声源频率为1500Hz时，本发明方法和常规等效源法求解所得单极子声源源强与实际源强比较；

图5为单极子声源频率为2000Hz时，采用加权迭代等效源法和常规等效源法求解所得单极子声源源强与实际源强的比较；

图6a为脉动球声源频率为100Hz时，本发明方法计算所得的声压分布；

图6b为脉动球声源频率为100Hz时，常规等效源法计算所得的声压分布；

图6c为脉动球声源频率为100Hz时，目标重建面上的理论声压分布；

图7a为脉动球声源频率为400Hz时，本发明方法计算所得的声压分布；

图7b为脉动球声源频率为400Hz时，常规等效源法计算所得的声压分布；

图7c为脉动球声源频率为400Hz时，目标重建面上的理论声压分布；

图8a为脉动球声源频率为700Hz时，本发明方法计算所得的声压分布；

图8b为脉动球声源频率为700Hz时，常规等效源法计算所得的声压分布；

图8c为脉动球声源频率为700Hz时，目标重建面上的理论声压分布；

图9为不同频率下本发明方法重建误差与常规等效源法重建误差的比较；

图10为信噪比10dB的低信噪比条件下本发明方法与常规等效源法重建误差的比较。

具体实施方式

本实施例中采用加权迭代等效源法的高精度近场声全息算法是按如下步骤进行：

步骤a、在声源近场辐射区域内布置全息面并测量全息面上的声压

在由声源产生的声场中，于近场区域内布置全息面H，所述全息面H为任意形状曲面，在所述全息面H上分布各测量网格点，利用传感器测量各测量网格点处的声压P_H。

步骤b、目标重建面T位于声源与全息面H之间，在目标重面T远离全息面H的一侧布置等效源面S_e，等效源面S_e与目标重建面之间的垂直距离为d_h，d_h的取值既不能太大也不能太小，d_h的取值太大将导致声压传递矩阵列向量间的相关性增强，造成矩阵条件数迅速增加，从而影响求解精度；而d_h的取值太小又会因为等效源距离目标重建面T太近而产生奇异性，本实施例中d_h取值范围设定为0.02m-0.1m；布置在等效源面S_e上的等效源的个数不大于全息面H上测量网格点的个数，这样可保证建立的方程组不是欠定的，从而可以增加求解的稳定性；等效源为单极子、偶极子或四极子，等效源为点源、面源或体源。

步骤c、建立各等效源与测量获得的全息面H上的声压P_H之间的关系为：

P_H＝G_HpQ

其中，G_Hp为等效源面S_e上各等效源与全息面H上各测量网格点之间的声压传递矩阵，Q为等效源源强向量。

步骤d、按如下过程求解等效源源强向量Q：

第一步：采用常规等效源法求解获得等效源源强向量作为等效源源强向量初始解Q⁰，令初始迭代次数为1，设定最大迭代次数为J。

由于本算法收敛迅速，通常只需3-5次迭代计算就能稳定地收敛到等效源源强向量Q的高精度解，此时如果继续增加迭代次数求解精度将基本保持不变，因此本算法在实施过程中J通常取为3-5次，最大不超过10次。

常规等效源法是指由学者A.Sarkissian等提出的，目前广泛用于近场声全息计算的等效源法，是采用以等效源源强向量Q的2范数||Q||²作为惩罚项的标准Tikhonov正则化方法来求解等效源源强向量Q。

第二步：利用向量Q^j-1作为后验加权系数构建加权矩阵W，W＝diag(Q^j-1)，diag(Q^j-1)是利用向量Q^j-1作为主对角元素构建的对角阵，所述向量Q^j-1为第j-1次迭代计算获得的等效源源强向量，对于第1次迭代计算过程，向量Q^j-1即为Q⁰。

第三步：利用所述加权矩阵W的逆矩阵W^-1构造向量Q^j的加权范数为||W^-1Q^j||；则在加权范数最小化准则意义下，通过如下极小化过程获得等效源源强向量Q^j的正则解：

min{||G_HpQ^j-P_H||²+λ²||W^-1Q^j||²}

通过预条件处理令W^-1Q^j≡X^j，将求解Q^j的加权范数正则化过程化为求解中间变量X^j标准的Tikhonov正则化过程：

min{||G_HpWX^j-P_H||²+λ²||X^j||²}

采用奇异值分解法求解中间变量X^j的正则解，正则化参数λ的选择采用广义交叉验证法或L曲线法来实现；经过第j次迭代获得等效源源强向量Q^j为：Q^j＝WX^j。

第五步：若j≤J则将迭代次数赋值为j+1，并转入第二步继续进行迭代；若j＞J则终止迭代，并令Q^j为等效源源强向量Q_f。

本实施例中，通过前一步迭代中获得的等效源源强向量的估计值作为后验加权系数，通过对正则化过程中惩罚项进行加权来进一步限制其解空间的大小，从而避免等效源源强能量的泄漏，使获得的等效源源强峰值更加突出，从而提高计算结果的分辨率和精度。

步骤e、利用P_T＝G_TpQ_f计算获得目标重建面T上的声压P_T；利用V_T＝G_TvQ_f计算获得目标重建面T上的法向振速V_T，其中，G_Tp为各等效源与目标重建面之间的声压传递矩阵，G_Tv为各等效源与目标重建面T之间的振速传递矩阵。

理论模型

本发明方法在声场计算过程中仍采用等效源模型，参见图1，其基本思想是采用在目标重建面后部背离分析域的区域内布置的一系列等效源，通过这些等效源产生的声场的线性叠加来近似实际声场；但与常规等效源法不同的是，在等效源源强求解过程中，本发明中提出了一种带有后验加权范数约束惩罚项的新型迭代正则化模型，并以此来代替传统的Tikhonov正则化模型，实现等效源源强的精确求解。其具体细节如下：

根据等效源法基本原理，图1中声源产生的声场可由等效源面上的一系列不同源强的等效源辐射声场的线性叠加来近似。假设等效源面S_e上分布有N个等效源，则场点r处的声压p(r)可由等效源表示为：

$p\left(r\right)=\underset{n=1}{\overset{N}{\mathrm{\Σ}}}i{\mathrm{\ρ}}_0\mathrm{ckg}(r,r_{\mathrm{on}})q\left(r_{\mathrm{on}}\right)---\left(1\right)$

其中，i为虚数单位，ρ₀为媒质密度，c为声速，k为声波波数，r_on为等效源面S_e上第n个等效源所在的位置，q(r_on)为第n个等效源对应的源强；g(r,r_on)为场点r与第n个等效源所在位置r_on之间的格林函数，并有：

$g(r,r_{\mathrm{on}})=\frac{1}{4\mathrm{\π}|r-r_{\mathrm{on}}|}{e}^{\mathrm{ik}|r-r_{\mathrm{on}}|}---\left(2\right)$

其中，|r-r_on|为场点r与第n个等效源所在位置r_on之间的距离。

由于全息面H位于声源的声场中，则由式(1)可知全息面H上M个测量网格点处的声压可表示为

$\begin{array}{cc}{p}_H\left(r_m\right)=\underset{n=1}{\overset{N}{\mathrm{\Σ}}}i{\mathrm{\ρ}}_0\mathrm{ckg}(r_m,r_{\mathrm{on}})q\left(r_{\mathrm{on}}\right)& m=\mathrm{1,2},...,M\end{array}---\left(3\right)$

其中，p_H(r_m)为全息面H上测量网格点r_m处的声压，g(r_m,r_on)为全息面H上第m个测量网格点与第n个等效源所在位置r_on之间的格林函数，将式(3)表达为矩阵的形式如式(4)：

P_H＝G_HpQ(4)

其中，P_H为全息面H上的M个测量网格点处的声压组成的列向量，Q为N个等效源源强组成的列向量，G_Hp为N个等效源到全息面上M个测量网格点处的声压传递矩阵，并有：

G_Hp(m,n)＝iρ₀ckg(r_m,r_on)(5)

Q＝[q(r_o1),q(r_o2),…,q(r_on),…,q(r_oN)](6)

式中，G_Hp(m,n)为各等效源与全息面H之间的声压传递矩阵G_Hp中位于第m行、第n列的矩阵元素。

对于式(4)逆向求解获得等效源源强向量Q为：

$Q=G_{\mathrm{Hp}}^{+}{P}_H---\left(7\right)$

其中，为G_Hp的广义逆矩阵。

等效源源强向量Q的求解精度直接决定了近场声全息计算结果的精度和分辨率。由于等效源源强向量Q的求解过程对应于数学上的逆问题，具有不适定性，微小的输入误差可能导致解的巨大变化；而实际中噪声误差是不可避免的，因此实际中并不采用式(5)进行计算，而是采用正则化方法求解源强向量Q。

常规等效源法采用Tikhonov正则化方法求解Q，其正则化过程可归结为求解如下极小化问题：

min{||G_HpQ-P_H||²+λ²||Q||²}(8)

其中：||||表示2范数，λ为正则化参数；惩罚项||Q||²是一个2范数约束，由于2范数约束是一个强光滑性约束条件，它倾向于将信号能量分散到所有的等效源中，从而获得一个光滑解。这样虽然可以较好地抑制噪声误差的影响，但往往导致等效源源强能量的泄漏，降低源强峰值，造成近场声全息计算结果分辨率与计算精度的损失。

因此，本发明方法构建了一种新型的加权迭代正则化方法，并以此代替Tikhonov正则化来求解等效源源强向量。该方法采用迭代的算法流程分为如下五个步骤：

第一步：采用常规等效源法求解获得等效源源强向量作为等效源源强向量初始解Q⁰，令初始迭代次数为1，设定最大迭代次数为J；

第二步：利用向量Q^j-1作为后验加权系数构建加权矩阵W，

W＝diag(Q^j-1)(9)

式(9)中diag(Q^j-1)是利用向量Q^j-1作为主对角元素构建的对角阵，向量Q^j-1为第j-1次迭代计算获得的等效源源强向量，对于第1次迭代计算过程，向量Q^j-1即为Q⁰；

第三步：利用加权矩阵W的逆矩阵W^-1构造向量Q^j的加权范数为||W^-1Q^j||；则在加权范数最小化准则意义下，通过如下极小化过程获得等效源源强向量Q^j的正则解：

min{||G_HpQ^j-P_H||²+λ²||W^-1Q^j||²}(10)

由于加权范数惩罚项能保持主要等效源的源强能量，同时能抑制由于噪声误差因素造成的源强很小的那些等效源的能量，因此不但能避免2范数惩罚项导致的源强能量泄漏，使获得的等效源源强峰值更加突出，还能在迭代过程中不断削弱对应于噪声误差因素的等效源能量，抑制噪声误差，提高算法精度和算法在低信噪比条件下的性能；同时由于每次迭代中均通过加权范数引入了前次迭代获得的有关等效源源强向量的后验信息来限制本次正则化过程中解空间的大小，从而能获得更加精确的、与实际情况更加符合的等效源源强向量正则解。

由于直接求解式(10)所示的正则化模型并不方便，实际计算中为了简化计算，通过预条件处理令W^-1Q^j≡X^j，将求解Q^j的加权范数正则化过程化为求解中间变量X^j标准的Tikhonov正则化过程：

min{||G_HpWX^j-P_H||²+λ²||X^j||²}(11)

采用奇异值分解法求解式(11)获得中间变量X^j的正则解，正则化参数λ的选择采用广义交叉验证法或L曲线法来实现；经过第j次迭代获得等效源源强向量Q^j为：

Q^j＝WX^j(12)

第五步：若j≤J则将迭代次数赋值为j+1，并转入第二步继续进行迭代；若j＞J则终止迭代，并令Q^j为等效源源强向量Q_f；

在获得Q_f后，通过Q_f结合等效源与目标重建面之间的声压、振速传递矩阵计算目标重建面T上的声压、质点振速等声场数据，从而实现声场的重建和预测。例如目标重建面T上的声压P_T和法向振速V_T可分别按式(13)和式(14)进行计算

P_T＝G_TpQ_f(13)

V_T＝G_TvQ_f(14)

其中，G_Tp和G_Tv分别为各等效源与目标重建面T之间的声压传递矩阵和振速传递矩阵，并有：

G_Tp(m,n)＝iρ₀ckg(r′_m,r_on)(15)

$G_{\mathrm{Tv}}(m,n)=\frac{\∂g(r_m^{\′},r_{\mathrm{on}})}{\∂n}---\left(16\right)$

其中，G_Tp(m,n)为各等效源与目标重建面T之间的声压传递矩阵G_Tp中位于第m行、第n列的矩阵元素，G_Tv(m,n)为各等效源与目标重建面T之间的声压传递矩阵G_Tv中位于第m行、第n列的矩阵元素，r′_m为目标重建面T上第m个网格点的位置，g(r′_m,r_on)为目标重建面T上第m个网格点与第n个等效源所在位置r_on之间的格林函数，n为目标重建面T的外法线方向，表示求法向导数。

方法的检验：

仿真1

验证本发明方法相对于常规等效源法具有更高的等效源源强求解精度。

仿真中声源为两个位于Z＝0m的平面内的单极子声源，声源位置分别为(0.1m,0,0)、(-0.1m,0,0)，强度为给定值0.1。

全息面H位于声源上方Z＝0.06m的平面内，H的尺寸为0.3m×0.3m,H上均匀分布7×7个测量网格点，因此全息面H上的测量间隔为0.05m。

为了比较等效源源强求解精度，将等效源面S_e同样布置在单极子声源所在的Z＝0m的平面内，并且S_e的尺寸大小以及其上等效源的分布与全息面H上测量网格点的分布一致，因此刚好有两个等效源与实际单极子声源位置重合。

为使仿真与实际实施中存在测量噪声的情况更加一致，全息面声压添加了高斯白噪声，信噪比为30dB。

分别采用本发明方法和常规等效源法求解S_e上各等效源源强，并与给定的理论值进行比较，以验证本发明方法的优势。采用本发明方法计算时，最大迭代次数取为3次。

根据给定的计算模型，等效源源强的理论分布应为如下形式：与实际单极子声源重合的等效源源强应为给定值0.1，而其他等效源强度应为0。图2、图3、图4以及图5分别给出了500、1000、1500、2000四个不同频率下，本发明方法与常规等效源法求得的等效源源强分布与理论分布间的对比，图2、图3、图4以及图5中曲线a为等效源源强的理论分布，曲线b为采用本发明方法获得的等效源源强分布，曲线c为采用常规等效源法获得的等效源源强分布。通过图2、图3、图4以及图5中曲线b、曲线c与曲线a的对比可知，采用本发明方法获得的等效源源强分布显然与理论分布更为一致,仅在等效源与实际单极子声源重合处有强度接近0.1的峰值，其他等效源源强接近于0；而常规等效源法求得的源强分布与理论值并不一致，其主峰值远小于0.1，并在其它本应为0的位置出现了多个小峰，这说明常规等效源法在求解过程中产生了源强能量泄漏。

上述仿真表明本发明方法可以克服常规等效源法源强求解过程中的能量泄漏问题，获得更高的源强求解精度，从而获得具有更高分辨率和更高重建精度的近场声全息计算结果。

仿真2

验证本发明方法相对于常规等效源法具有更高的分辨率以及更高的声场重建精度。

声源为两个位于Z＝0m的平面内的脉动球声源，声源位置分别为(0.075m,0,0)、(-0.075m,0,0),脉动球半径为0.05m，振速为0.1m/s。

全息面H位于声源上方Z＝0.09m的平面内，H的尺寸为0.6m×0.6m,H上均匀分布25×25个测量网格点，因此全息面H上的测量间隔为0.025m。

目标重建面T位于Z＝0.05m的平面内，并且目标重建面T的尺寸大小以及其上网格点的分布与全息面H上测量网格点的分布一致。

等效源面S_e布置在目标重建面T远离全息面H的一侧，与目标重建面T之间的垂直距离d_h为0.03m，因此等效源面S_e位于Z＝0.02m的平面内，并且S_e的尺寸大小以及其上等效源的分布与全息面H上测量网格点的分布一致。

为进一步增加计算难度，并检验方法的鲁棒性，对全息面声压添加了高斯白噪声，信噪比为20dB。

仿真中分别采用本发明方法以及常规等效源法对100-2000Hz范围内20个频率(频率间隔100Hz)的脉动球声源进行了声压重建计算，所有计算中最大迭代次数J均取为3次。

图6a、图6b、图6c、图7a、图7b、图7c以及图8a、图8b和图8c分别给出了脉动球声源频率为100Hz、400Hz、700Hz时本发明方法和常规等效源法重建获得的目标重建面上声压分布图，以及相应的理论声压分布图，其中，图6a、图7a和图8a为本发明方法获得的目标重建面上声压分布，图6b、图7b和图8b为常规等效源法获得的目标重建面上声压分布，图6c、图7c和图8c为按照脉动球声压解析表达式计算得到的目标重建面理论声压。

通过分别对比图6a、图6b、图6c、图7a、图7b、图7c以及图8a、图8b和图8c可知，采用本发明方法获得的目标重建面声压分布与理论值非常一致，重建结果中清晰地显示出了两个声压峰值，并且峰值位置与模型中脉动球声源的布置位置完全一致；而常规等效源法获得的重建结果却由于等效源源强能量的泄漏，丢失了细节信息，重建结果中只有一个光滑的峰值，根本无法看出两个脉动球声源的存在，与理论值相差很大。

上述仿真结果表明本发明方法相对于常规等效源法能够更准确地识别声源，更精确地重建出声源的细节信息，也就是说本发明方法具有更高的分辨率和更好的精度。为了定量描述重建结果的精确程度，定义声压重建误差E_p为

E_p＝||P_True-P_{reconstructed}||/||P_True||(17)

其中，P_True表示目标重建面上的理论声压，P_{reconstructed}表示计算获得的重建声压。

按式(17)对本发明方法在100Hz-2000Hz频率范围内的重建误差进行了计算，并与常规等效源法在相同频率的重建误差进行了比较，参见图9。

图9表明在所有计算频率处，本发明方法的重建误差均小于常规等效源法，并且在绝大多数频率处，本发明方法的重建误差不足常规等效源法误差的二分之一。

上述仿真表明在本发明方法在分辨率与重建精度方面均明显优于基于常规等效源法的近场声全息。

仿真3

验证本发明方法在低信噪比条件下的优越性

仿真条件基本与仿真2中相同，仅将全息面声压的信噪比降低为10dB，分别采用本发明方法以及常规等效源法计算目标重建面声压分布，并与理论值进行比较。计算中最大迭代次数J取为3次，所得结果按式(17)计算两种方法相应的重建误差，如图10所示，图中曲线a为常规等效源法在各频率处的重建误差，曲线b为本发明方法在各频率处的重建误差。对比曲线a和曲线b可知在100Hz-2000Hz的频带内所有分析频率处，本发明方法的重建误差都明显低于常规等效源法的重建误差，也就是说采用本发明方法使近场声全息重建精度获得了显著提高。

上述仿真表明在低信噪比条件下，本发明方法具备更加良好的抗噪声干扰能力，其重建误差在整个分析频带内均明显低于常规等效源法的重建误差。

Claims (3)

1.采用加权迭代等效源法的高精度近场声全息算法，其特征是按如下步骤进行：

步骤a、在声源近场辐射区域内布置全息面并测量全息面上的声压

在由声源产生的声场中，于近场区域内布置全息面H，所述全息面H为任意形状曲面，在所述全息面H上分布各测量网格点，利用传感器测量各测量网格点处的声压P_H；

步骤b、目标重建面T位于声源与全息面H之间，在目标重建面T远离全息面H的一侧布置等效源面S_e，所述等效源面S_e与目标重建面T之间的垂直距离为d_h，布置在所述等效源面S_e上的等效源的个数不大于全息面H上测量网格点的个数；所述等效源为单极子、偶极子或四极子，所述等效源为点源、面源或体源；

步骤c、建立各等效源与测量获得的全息面H上的声压P_H之间的关系为：

P_H＝G_HpQ

其中，G_Hp为等效源面S_e上各等效源与全息面H上各测量网格点之间的声压传递矩阵，Q为等效源源强向量；

步骤d、按如下过程求解等效源源强向量Q：

第一步：采用常规等效源法求解获得等效源源强向量作为等效源源强向量初始解Q⁰，令初始迭代次数为1，设定最大迭代次数为J；

第二步：利用向量Q^j-1作为后验加权系数构建加权矩阵W，W＝diag(Q^j-1)，diag(Q^j-1)是利用向量Q^j-1作为主对角元素构建的对角阵，所述向量Q^j-1为第j-1次迭代计算获得的等效源源强向量，对于第1次迭代计算过程，向量Q^j-1即为Q⁰；

第三步：利用所述加权矩阵W的逆矩阵W^-1构造向量Q^j的加权范数为‖W^-1Q^j‖；则在加权范数最小化准则意义下，通过如下极小化过程获得等效源源强向量Q^j的正则解：

min{‖G_HpQ^j-P_H‖²+λ²‖W^-1Q^j‖²}

通过预条件处理令W^-1Q^j≡X^j，将求解Q^j的加权范数正则化过程化为求解中间变量X^j标准的Tikhonov正则化过程：

min{‖G_HpWX^j-P_H‖²+λ²‖X^j‖²}

采用奇异值分解法求解中间变量X^j的正则解，正则化参数λ的选择采用广义交叉验证法或L曲线法来实现；经过第j次迭代获得等效源源强向量Q^j为：Q^j＝WX^j；

第五步：若j≤J则将迭代次数赋值为j+1，并转入第二步继续进行迭代；若j＞J则终止迭代，并令Q^j为等效源源强向量Q_f；

步骤e、利用P_T＝G_TpQ_f计算获得目标重建面T上的声压P_T；利用V_T＝G_TvQ_f计算获得目标重建面T上的法向振速V_T，其中，G_Tp为各等效源与目标重建面之间的声压传递矩阵，G_Tv为各等效源与目标重建面T之间的振速传递矩阵。

2.根据权利要求1所述的采用加权迭代等效源法的高精度近场声全息算法，其特征是：所述d_h的取值为0.02m-0.1m。

3.根据权利要求1所述的采用加权迭代等效源法的高精度近场声全息算法，其特征是：所述最大迭代次数J不超过10次。