CN105425212A - 一种声源定位方法 - Google Patents

Abstract

本发明涉及一种声源定位方法，包括以下步骤：1.布放传感器节点，获得风场参数。2根据待测声源到各个传感器节点的相对位置，对不同传播方向的声速进行修正。3.利用修正的声速表示出待测声源到各传感器节点的传播时间。4.计算待测声源到达每两个传感器节点的第一时延差，利用信号检测技术测得声源到达对应两个传感器的第二时延差，第一时延差与第二时延差联立方程组，求解获得声源位置。本发明对由于风场造成的不同传播方向上声速不同的影响进行了根本性地修正，从而克服了传统方法在非均匀声速问题上的缺陷，并且由于本发明没有引入其他变量，没有增加对传感器节点个数的要求。

Description

一种声源定位方法

技术领域

本发明属于声学信号处理领域，涉及分布式传感器定位技术，尤其涉及一种声源定位方法。

背景技术

随着无线传感网络及电子技术的发展，声源定位技术因其成本低，安全性高，隐蔽性高等优点，在军事和民用领域得到了迅速发展，声源定位技术已经从对空中直升飞机、地面车辆和炮弹落点、水面舰船以及水下潜艇等目标定位逐渐延伸到了地震探测、医学成像等民用领域，呈现出越来越大的应用价值。

分布式声源定位技术是用多个传感器按照一定的规则布放，每个传感器对声音信号进行采集，根据采集到的信号利用信号处理算法确定声源的位置。

利用分布式声源定位技术进行定位时，其定位精度除了受到传感器本身的影响，例如传感器节点的不一致、节点布放规则等，还会受到周围环境的影响，例如温度、风等环境因素，其中，对于传感器节点的布放位置的研究已经比较广泛，不同的研究人员基于不同的系统要求提出了不同的布放形状。环境对定位精度的影响主要体现在不同的环境因素如温度、风速对声速的影响。

针对声速对定位性能的影响影响，主要是通过将声速视为未知量与声源位置进行联合估计，或者通过引入辅助设施或利用后处理技术进行修正。这些方法均假设声源到达各个传感器的声速是相同的。但声源以球面波的形式向四周传播，在风场条件下，声源向不同方向传播的速度就会由于风场存在而受到不同程度的影响，造成声源到各个传感器的声速不同。前期的研究之所以一般不考虑此问题，是因为节点布置一般比较集中，或声源位于测试系统的远场，而当传感器与声源均散布在较大区域内时(如数公里范围)，这个问题就凸显出来，成为影响定位性能的关键因素之一。

TDOA定位技术是一种常用的声源定位技术，该方法根据声源到达不同传感器的时间差(TDOA)，得到一组超定非线性方程。通过求解非线性方程组就可以得到估计的声源位置。理论上三个传感器节点就可以定出声源位置，但是由于各种测量误差的存在，仅使用三个传感器节点进行定位得到的结果一般不太理想，通常采用多个传感器按照一定的规则布放在测试区域进行声源的精确定位。

在气象学中，将距离地面50m-100m以下的气层称为近地层，在我们通常的应用中，声信号主要工作在这一层，风速一般为每秒几米至十几米，甚至更大些。风对声信号传播产生的影响主要有两方面：一是风速对实际声速的直接影响，二是风速梯度使声线产生弯曲；通过下面的分析可以看到，由于传感器的分布式部署以及风场的存在，声源到达各个传感器的传播途径上的声速将会不同。

发明内容

为解决现有技术存在的问题，本发明的目的是在传统的声源定位算法的基础上对声源到不同传感器的传播速度进行修正，进而进行精确的声源定位。

为达到上述目的，本发明提供的技术方案是：

一种声源定位方法，包括以下步骤：

步骤一：将多个传感器节点布放在测试区域内，获得各个传感器节点的位置坐标；利用风速风向传感器获得测试区域的风场参数，风场参数包括风向和风速；利用温度传感器获得测试区域的当前温度，计算当前温度下的无风时声速；

步骤二：根据待测声源到各个传感器节点的相对位置，计算风速在该方向上的分量，结合当前温度下无风时声速对不同传播方向的声速进行修正；

步骤三：利用步骤二修正后的声速计算待测声源到各传感器节点的传播时间；

步骤四：计算待测声源到达每两个传感器节点的第一时延差，利用信号检测技术测得声源到达对应两个传感器的第二时延差，第一时延差与第二时延差联立方程组，求解方程组，获得声源位置。

所述步骤二中根据待测声源到各个传感器节点的相对位置，计算风速在该方向上分量，对各个传感器节点的声速进行修正，其具体方法如下：

风向θ₀与参考方向x轴的夹角用α₀表示；风速用ν₀表示；当前温度下无风时的声速用c₀表示；待测声源的位置为S(x,y)，各个传感器节点的位置为A_i(x_i,y_i)，其中，i＝1,2,…,N；与参考方向x轴正向的角度为α_i∈[-π,π]；

风速在从待测声源到第i个传感器节点的传播方向上的分量为:v＝ν₀cosθ_i＝ν₀cos(α_i-α₀)；

其中， ${\mathrm{\&alpha;}}_i=\left\{\begin{array}{c}\arctan \left(\frac{y_i-y}{x_i-x}\right)-\mathrm{\&pi;},\left(y_i-y\&le;0,x_i-x<0\right)\\ \arctan \left(\frac{y_i-y}{x_i-x}\right),\left(x_i-x>0\right)\\ \arctan \left(\frac{y_i-y}{x_i-x}\right)+\mathrm{\&pi;},\left(y_i-y>0,x_i-x<0\right)\end{array}\right.;$

待测声源到不同传感器节点的修正声速为：c_i＝c₀+v₀cos(α_i-α₀)；

所述步骤三利用修正的风速计算待测声源到各传感器节点传播时间，其具体方法如下：

待测声源到传感器节点的时间为：其中，

所述步骤四联立第一时延差与第二时延差联立方程组，获得声源位置的具体方法如下：

以第一个传感器节点(i＝1)为参考节点，当声速进行修正时,声源实际到达第i个传感器节点与到达参考节点的第一时延差为

利用信号检测技术测得的声源到达对应两个传感器第二时延差为τ_i1，根据第一时延差和第二时延差可以得到方程：其中i＝2,3,...,N，e_i为测量误差；令则f_i(x,y)＝τ_i1-e_i是关于(x,y)的非线性方程；在大于2个传感器节点的情况下，可以得到N-1个关于S(x,y)的非线性方程，求解方程组得到声源位置。

所述利用信号检测技术测得的声源到达对应两个传感器第二时延差_τi1，其具体方法如下：

对传感器节点采集到的M个点的时域信号s(t)，求快速傅立叶变换(FFT)转换到频域s(f)，将M点分成n个频段进行处理，(M为2的幂函数)，其中n远远小于M；计算各个频段频谱能量的平均值，记为E＝[E₁E₂…E_n]。初始检测门限向量为T₀＝E+C，其中C为各频段允许超出的最大能量值；

信号检测门限向量为T＝ρT₀+(1-ρ)(E+C)，(ρ＞0.5)；

当检测到E中存在任意一个元素大于T中对应的门限值，即判断目标信号出现，否则以T作为新的T₀，进行迭代计算新的检测门限向量，继续检测目标信号。当检测到目标信号时，以该时段信号的时域最大值所对应的时刻作为信号到达传感器的时刻；将测量的信号到达两个传感器节点的时刻相减即可得到第二时延差为τ_i1。

本发明具有以下优点：

与现有定位技术相比，本方法有效的解决了传统TDOA声源定位算法由于不考虑测试场区内声速的差异，使用统一声速造成额外的定位误差。对该算法提出不同传播途径采用不同声速的改进措施，针对均匀风场，基于已知风场参数(风速和风向)，对声源到各传感器的不同传播方向的声速进行数学表述，然后利用修正的声速进行定位。本发明在不引入其他变量的情况下，对由于风场造成的不同传播方向上声速进行了根本性地修正，从而克服了传统算法非均匀声速问题上的缺陷。本发明的信号检测方法中检测门限根据背景信号的变化而缓慢地自适应变化，能够自始至终将门限与慢变的噪声背景相吻合，从而保证检测的灵敏度和可靠性。

附图说明

图1是本发明分布式声源定位方法的模型示意图；

具体实施方式

下面结合附图和具体实施方式对本发明进行进一步的说明。

下面以二维空间为例，例举本发明的详细实施步骤和计算方法：

如图1所示：

第一步：将传感器节点布放在测试区域内，获得各个传感器节点的位置坐标；利用风速风向传感器获得风场参数，包括风向和风速；同时利用温度传感器获得测试场区的温度T；根据温度传感器测量出的温度计算当前温度下的声速，温度(T)与声速的经验关系为:

c₀≈331.6+0.6T(1)

第二步：根据待测声源到各个传感器节点的相对位置，计算风速在该方向上的分量，结合当前温度下无风时声速对不同传播方向的声速进行修正；

1)风向与参考方向x轴的夹角为α₀，风速为ν₀；

2)目标声源位置为S(x,y)；各个节点的位置为A_i(x_i,y_i)，其中i＝1,2,…,N，参

考节点(i＝1)位于坐标原点，；

3)当前温度下无风时的声速为c₀；

4)与x轴(参考方向)正向的角度为α_i∈[-π,π]。

根据当前温度下无风时的声速c₀，利用最小二乘法计算泰勒级数展开法迭代的初始点(x⁽⁰⁾,y⁽⁰⁾)。其具体步骤如下:

计算声源到传感器节点的距离R_i，根据声源到达第i个传感器A_i(x_i,y_i)与到达参考节点A₁(x₁,y₁)的时间差τ_i,1，可得到方程

R_i-R₁＝c₀τ_i,1(2)

将式(2)两边平方并整理，可得到方程

$x_{i}x+y_{i}y+c_0{\mathrm{\&tau;}}_{i,1}{R}_1=x_i^2+y_i^2-{\left(c_0{\mathrm{\&tau;}}_{i,1}\right)}^2---\left(3\right)$

根据N‐1个方程，可以得到矩阵方程

ΦZ＝b(4)

$\mathrm{\&Phi;}=\left[\begin{array}{ccc}{x}_2& y_2& c_0{\mathrm{\&tau;}}_{2,1}\\ x_3& y_3& c_0{\mathrm{\&tau;}}_{3,1}\\ .& .& .\\ .& .& .\\ .& .& .\\ x_N& y_N& c_0{\mathrm{\&tau;}}_{N,1}\end{array}\right]---\left(5\right)$

Z＝[xyR₁]^Τ(6)

b＝[b₁b₂…b_N-1]^Τ(7)

$b_i=\frac{1}{2}(x_{i-1}^2+y_{i-1}^2-c_0^2{\mathrm{\&tau;}}_{i-1,1}^2)---\left(8\right)$

利用最小二乘法对式(8)求解可得

Z＝(Φ^ΤΦ)^-1Φ^Τb(9)

解出来的(x,y)的值作为迭代初始点(x⁽⁰⁾,y⁽⁰⁾)。若系数矩阵Φ不是列满秩矩阵，则将Φ^ΤΦ的元素加上适当的摄动，使其变为可逆矩阵。

根据(x⁽⁰⁾,y⁽⁰⁾)和各个传感器节点的相对位置(x_i,y_i)，根据式(10)对风速在该方向上的分量进行数学表述，利用风速的分量对各个传感器节点的声速进行修正。

声源到各个节点的传播方向α_i为

${\mathrm{\&alpha;}}_i=\left\{\begin{array}{c}\arctan \left(\frac{y_i-y}{x_i-x}\right)-\mathrm{\&pi;},\left(y_i-y\&le;0,x_i-x<0\right)\\ \arctan \left(\frac{y_i-y}{x_i-x}\right),\left(x_i-x>0\right)\\ \arctan \left(\frac{y_i-y}{x_i-x}\right)+\mathrm{\&pi;},\left(y_i-y>0,x_i-x<0\right)\end{array}\right.---\left(10\right)$

由式(10)可以看出，α_i为目标声源位置的函数，并且与待测声源和节点的相对位置有关。参考图1待测声源到各个传感器节点的实际声速为

c_i＝c₀+v₀cos(α_i-α₀)(11)

第三步：利用修正的风速c_i计算待测声源到各传感器节点的传播时间；待测声源到传感器节点的时间为：其中， $R_i=\sqrt{{(x-x_i)}^2+{(y-y_i)}^2}.$

第四步：计算待测声源到达每两个传感器节点的第一时延差，利用信号检测技术测得声源到达对应两个传感器的第二时延差，第一时延差与第二时延差联立方程组，求解方程组，获得声源位置。

当声速用式(11)进行修正时，声源实际到达第i个传感器节点A_i(x_i,y_i)与到达参考节点A₁(x₁,y₁)的第一时延差为

利用信号检测技术测得的声源到达对应两个传感器第二时延差为τ_i1，具体方法为：

所述利用信号检测技术测得的声源到达对应两个传感器第二时延差τ_i1，其具体方法如下：

对传感器节点采集到的M个点的时域信号s(t)，求快速傅立叶变换(FFT)转换到频域s(f)，将M点分成n个频段进行处理，(M为2的幂函数)，其中n远远小于M。计算各个频段频谱能量的平均值，记为E＝[E₁E₂…E_n]；初始检测门限向量为T₀＝E+C，其中C为各频段允许超出的最大能量值；

信号检测门限向量为T＝ρT₀+(1-ρ)(E+C)，(ρ＞0.5)；

当检测到E中存在任意一个元素大于T中对应的门限值，即判断目标信号出现，否则以T作为新的T₀，进行迭代计算新的检测门限向量，继续检测目标信号。当检测到目标信号时，以该时段信号的时域最大值所对应的时刻作为信号到达传感器的时刻；将测量的信号到达两个传感器节点的时刻相减即可得到第二时延差为τ_i1。

第一时延差与第二时延差联立方程可得

$\frac{R_i}{c_i}-\frac{R_1}{c_1}={\mathrm{\&tau;}}_{i1}-e_i---\left(12\right)$

其中i＝2,3,...,N，其中e_i为测量误差。令

$f_i(x,y)=\frac{R_i}{c_i}-\frac{R_1}{c_1}---\left(13\right)$

则

f_i(x,y)＝τ_i1-e_i(14)

此时，式(14)是关于(x,y)的非线性方程。假设有N(N>2)个传感器节点，则可以得到N-1个关于(x,y)的非线性方程，联立可得非线性方程组。

用泰勒级数展开法解非线性方程组；将式(14)在点(x⁽⁰⁾,y⁽⁰⁾)泰勒级数展开并忽略二阶以上分量，令δ_x＝x-x⁽⁰⁾,δ_y＝y-y⁽⁰⁾，则可得到

$f_i(x,y)\&ap;f_i(x^{\left(0\right)},y^{\left(0\right)})+\frac{\&part;f_i(x,y)}{\&part;x}{|}_{(x^{\left(0\right)},y^{\left(0\right)})}{\mathrm{\&delta;}}_x+\frac{\&part;f_i(x,y)}{\&part;y}{|}_{(x^{\left(0\right)},y^{\left(0\right)})}{\mathrm{\&delta;}}_y---\left(15\right)$

令

$a_{i1}=\frac{\&part;f_i(x,y)}{\&part;x}{|}_{(x^{\left(0\right)},y^{\left(0\right)})},a_{i2}=\frac{\&part;f_i(x,y)}{\&part;y}{|}_{(x^{\left(0\right)},y^{\left(0\right)})}---\left(16\right)$

$A=\left[\begin{array}{cc}{a}_{21}& a_{22}\\ a_{31}& a_{32}\\ .& .\\ .& .\\ .& .\\ a_{n1}& a_{n2}\end{array}\right]---\left(17\right)$

$\mathrm{\&delta;}=\left[\begin{array}{c}{\mathrm{\&delta;}}_x\\ {\mathrm{\&delta;}}_y\end{array}\right]---\left(18\right)$

$Z=\left[\begin{array}{c}{\mathrm{\&tau;}}_{21}-f_2(x^{\left(0\right)},y^{\left(0\right)})\\ {\mathrm{\&tau;}}_{31}-f_3(x^{\left(0\right)},y^{\left(0\right)})\\ .\\ .\\ .\\ {\mathrm{\&tau;}}_{n1}-f_n(x^{\left(0\right)},y^{\left(0\right)})\end{array}\right]---\left(19\right)$

$e=\left[\begin{array}{c}{e}_2\\ e_3\\ .\\ .\\ .\\ e_n\end{array}\right]---\left(20\right)$

上式可写为

$A\mathrm{\&delta;}\&cong;Z-e---\left(21\right)$

所以加权最小二乘估计为

δ＝(A^ΤQ^-1A)^-1A^ΤQ^-1Z(22)

其中Q＝cov(ee^Τ),cov(·)表示互相关运算。此时

x⁽¹⁾＝δ_x+x⁽⁰⁾,y⁽¹⁾＝δ_y+y⁽⁰⁾(23)

若ε为设定的门限值，则式(23)计算结果即为目标声源位置；反之，则以式(23)计算出来的值作为新的初始值，返回第二步进行重新计算。直到满足条件为止。

本发明的内容不限于实施例所列举，本领域普通技术人员通过阅读本发明说明书而对本发明技术方案采取的任何等效的变换，均为本发明的权利要求所涵盖。

Claims (5)

1.一种声源定位方法，其特征在于：

包括以下步骤：

步骤一：将多个传感器节点布放在测试区域内，获得各个传感器节点的位置坐标；利用风速风向传感器获得测试区域的风场参数，风场参数包括风向和风速；利用温度传感器获得测试区域的当前温度，计算当前温度下的无风时声速；

步骤二：根据待测声源到各个传感器节点的相对位置，计算风速在该方向上的分量，结合当前温度下无风时声速对不同传播方向的声速进行修正；

步骤三：利用步骤二修正后的声速计算待测声源到各传感器节点的传播时间；

步骤四：计算待测声源到达每两个传感器节点的第一时延差，利用信号检测技术测得声源到达对应两个传感器的第二时延差，第一时延差与第二时延差联立方程组，求解方程组，获得声源位置。

2.根据权利要求1所述的风场条件下的分布式声源定位方法，其特征在于：

所述步骤二中根据待测声源到各个传感器节点的相对位置，计算风速在该方向上分量，对各个传感器节点的声速进行修正，其具体方法如下：

风向θ₀与参考方向x轴的夹角用α₀表示；风速用ν₀表示；当前温度下无风时的声速用c₀表示；待测声源的位置为S(x,y)，各个传感器节点的位置为A_i(x_i,y_i)，其中，i＝1,2,…,N；与参考方向x轴正向的角度为α_i∈[-π,π]；

风速在从待测声源到第i个传感器节点的传播方向上的分量为:v＝ν₀cosθ_i＝ν₀cos(α_i-α₀)；其中，

${\mathrm{\&alpha;}}_i=\left\{\begin{array}{c}\arctan \left(\frac{y_i-y}{x_i-x}\right)-\mathrm{\&pi;},\left(y_i-y\&le;0,x_i-x<0\right)\\ \arctan \left(\frac{y_i-y}{x_i-x}\right),\left(x_i-x>0\right)\\ \arctan \left(\frac{y_i-y}{x_i-x}\right)+\mathrm{\&pi;},\left(y_i-y>0,x_i-x<0\right)\end{array}\right.;$

待测声源到不同传感器节点的修正声速为：c_i＝c₀+v₀cos(α_i-α₀)。

3.根据权利要求2所述的风场条件下的分布式声源定位方法，其特征在于：

所述步骤三利用修正的风速计算待测声源到各传感器节点传播时间，其具体方法如下：

待测声源到传感器节点的时间为：其中，

$R_i=\sqrt{{(x-x_i)}^2+{(y-y_i)}^2}.$

4.根据权利要求3所述的均匀风场条件下的分布式声源定位方法，其特征在于：

所述步骤四联立第一时延差与第二时延差联立方程组，获得声源位置的具体方法如下：

以第一个传感器节点(i＝1)为参考节点，当声速进行修正时,声源实际到达第i个传感器节点与到达参考节点的第一时延差为

利用信号检测技术测得的声源到达对应两个传感器第二时延差为τ_i1，根据第一时延差和第二时延差可以得到方程：其中i＝2,3,...,N，e_i为测量误差；令则f_i(x,y)＝τ_i1-e_i是关于(x,y)的非线性方程；在大于2个传感器节点的情况下，可以得到N-1个关于S(x,y)的非线性方程，求解方程组得到声源位置。

5.根据权利要求4所述的均匀风场条件下的分布式声源定位方法，其特征在于：

所述利用信号检测技术测得的声源到达对应两个传感器第二时延差τ_i1，其具体方法如下：

对传感器节点采集到的M个点的时域信号s(t)，求快速傅立叶变换(FFT)转换到频域s(f)，将M点分成n个频段进行处理，(M为2的幂函数)，其中n远远小于M；计算各个频段频谱能量的平均值，记为E＝[E₁E₂...E_n]。初始检测门限向量为T₀＝E+C，其中C为各频段允许超出的最大能量值；

信号检测门限向量为T＝ρT₀+(1-ρ)(E+C)，(ρ＞0.5)；

当检测到E中存在任意一个元素大于T中对应的门限值，即判断目标信号出现，否则以T作为新的T₀，进行迭代计算新的检测门限向量，继续检测目标信号。当检测到目标信号时，以该时段信号的时域最大值所对应的时刻作为信号到达传感器的时刻；将测量的信号到达两个传感器节点的时刻相减即可得到第二时延差为τ_i1。