一种基于卡尔曼滤波的声源三维定位方法与流程

本发明涉及一种声源定位方法，尤其是涉及一种基于卡尔曼滤波的声源三维定位方法。

背景技术：

近年来，得益于现代计算机技术、现代数字信号处理技术的发展，声源定位技术有了长足的进步。该项技术的特点在于：系统本身仅依赖于目标声源的声音信号，通过对阵列信号的处理实现声源位置的检测和定位。传声器阵列在时域和频域的基础上增加了空间域，对声音信息的处理能力增强。基于麦克风阵列的声源定位研究一直是非常活跃的研究课题。许多国家如美国、法国、意大利、以色列和加拿大等相继开展了这方面的研究工作并卓有成效。研制的很多麦克风阵列声源定位系统已经广泛应用于社会生活的各种场合，发挥着极其重要的作用。

目前基于传声器阵列的声源定位算法主要包括三类：基于信号到达时差的定位算法(TDOA)、基于波束形成的定位算法(BF)和基于谱估计的定位算法。这三类算法仅能实现声源的方位角估计，无法实现声源的三维定位。近年来已经有研究将基于谱估计的定位算法拓展到三维定位领域，但定位效果受到测量过程中误差的干扰，精度有待进一步提高。TDOA和BF算法均无法应用于多元声源的定位。三类算法均建立在声源远场模型的基础上，故对于近场声源的定位更是无能为力。单一的余弦定理算法，由于应用过程中误差的存在，定位效果难以保证。所以必须提高相关算法的抗噪能力，提升算法的定位精度。首先，分析影响定位精度的影响因素。由计算过程可知，主要影响因素包括3个：1)全息面上的投影误差。由于重建、投影网格点与真实声源并不能完全重合，这就造成计算声源、真实声源在全息面上的投影距离存在误差。该误差无法通过加密重建网格消除；2)传声器采集误差。该误差存在于所有采集系统中，将使上述的投影距离进一步出现偏差；3)转动角度误差。该误差在采集过程中也是普遍存在的，并且随着转动角度的不同而变化。转动角度小时，该误差影响较大，转动角度大时，该误差影响较小。总体来说声源的三维定位方法目前还在探索之中。

技术实现要素：

本发明的目的就是为了克服上述现有技术存在的缺陷而提供一种基于卡尔曼滤波的声源三维定位方法。

本发明的目的可以通过以下技术方案来实现：

一种基于卡尔曼滤波的声源三维定位方法，包括下列步骤：

步骤一、利用传声器阵列采集声场中的声压信号。

采集声场中信号的传声器阵列为随机面阵列。阵元数为36，阵列直径0.5m，可以满足近场声源的定位需求。采集后的数据经过矩形窗处理，相当于对三维声场进行了局部二维平面采集。

步骤二、对采集的声压信号进行傅立叶变换，并选取所要分析的频率。

在进行傅立叶变换中，加窗处理的窗函数为汉宁窗(hanning window)，同时选取需要分析的频率数据。采集的数据长度为1s，那么频率分辨率为1Hz，满足大部分的分析需求。

步骤三、利用统计最优近场声全息计算声源映射在传声器阵列面上的二维坐标。

统计最优近场声全息是一种局部声全息算法，适用于大声源小阵列的情况。算法不强制要求阵列必须是规则阵列(如网格阵列)，并且算法对测量面以外的声源不敏感。算法实施过程中，声场被假定为单元波的叠加，单元波的幅值均假设等于1，单元波的个数应该大于阵元数。

步骤四、将传声器阵列沿纵轴旋转一个小角度，并利用统计最优近场声全息计算声源映射在旋转后阵列面上的二维坐标。具体方法如步骤三中所述。值得注意的是，由于两次计算的声源坐标有可能位于旋转轴的同一侧或者两侧，这两种情况的计算公式略有不同。

将起始测量面b与终止测量面a的夹角(α)进行n等分，等分角度α1＝α/n，将b全息面沿y轴旋转角度α1，定义为α1全息面(即：下次旋转的起始全息面)。利用统计最优近场声全息确定声源几何中心映射在b全息面、α1上的2D坐标。

优选地，α设定为60度，n定为60。

步骤五、重复步骤三和步骤四，获得两组二维坐标数据和一组旋转角度数据，利用两组二维坐标数据和对应的角度数据，计算出一组三维坐标数据，即声源相对b全息面的虚拟运动轨迹，将此三维坐标数据作为声源位置数据。过程中，需要将这组数据转换到初始测量面所对应的坐标系下，由运动的相对性原理可知，测量面绕垂直轴轴旋转相当于声源绕同一轴旋转，且总旋转角度相等。这样就形成一条声源绕轴旋转的虚拟运动轨迹。优选地，根据实际情况可以将总旋转角度定为60度。

声源到全息面距离的计算公式分两种情况：1)投影映射点在转轴的两侧，其计算公式如式(1)～式(5)；2)投影映射点在转轴的同侧，其计算公式如式(4)～式(9)。

γ＝180-α-β(3)

β＝90-α-θ (8)

γ＝90+θ (9)

式中，Xab为b全息面的映射点到a全息面的映射点之间的距离，Xa为a全息面的映射点到a全息面与b全息面的交叉点之间的距离，Xb为b全息面的映射点到b全息面与a全息面的交叉点之间的距离。Da为声源S到a全息面的距离，Db为声源S到b全息面的距离，β为Da与Xab之间的夹角，γ为Db与Xab之间的夹角，α为Da与Db之间的夹角。

将所有3D坐标转化为b全息面坐标系下，转换原理如下：声源在旋转后全息面x轴上投影坐标数组为(Xa1，Xa2，Xa3，…，Xan)、旋转角度数组为(α1，α2，α3，…，αn)、声源在b全息面x轴上的投影坐标Xb不变，声源在所有全息面上y轴的投影坐标Ya1＝Ya2＝Ya3＝…＝Yan＝Yb。基于上述数据及权利要求4中给出的公式可计算出：声源到的旋转后全息面的垂直距离(Da1，Da2，Da3，…，Dan)。在b全息面的坐标系下，坐标点(Xb，Yb，Db)，(Xa1，Y a1，Da1)，(Xa2，Ya2，Da2)，…，(Xan，Yan，Dan)即可组成声源相对于b全息面坐标系的虚拟运动轨迹线。

优选地，步骤五中重复步骤三和步骤四的次数为四十次。

步骤六、声源位置数据中含有测量误差和网格误差，利用卡尔曼滤波对上述声源位置数据进行滤波估计，以抑制上述误差，最终实现声源的三维定位。

与现有技术相比，本发明考虑了实践中的测量误差和网格误差，这些误差具有高斯分布特性，卡尔曼滤波作为高效率的递归滤波器，能够从一系列的不完全且包含高斯噪声的测量数据中，预测和估计出信号数据，本发明通过多次旋转全息面获得一组声源位置数据，并利用卡尔曼滤波抑制误差对定位结果的影响，提升了定位结果的精度。

附图说明

图1为本发明方法的流程示意图；

图2为本发明的声源重建网格示意图；

图3为投影映射点示意图，其中，图3(a)为投影映射点在转轴的同侧的示意图，图3(b)为投影映射点在转轴的两侧的示意图；

图4为声源位移推算示意图；

图5为本发明实施例中的定位效果图。

具体实施方式

下面结合附图和具体实施例对本发明进行详细说明。显然，所描述的实施例是本发明的一部分实施例，而不是全部实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例，都应属于本发明保护的范围。

如图1所示，本发明涉及一种基于卡尔曼滤波的声源三维定位方法，包括下列步骤：

(一)利用传声器阵列采集声场中的声压信号，测量阵列为直径0.5m的随机面阵列，阵元数为36，空间窗选定为空间矩形窗函数，采集仪器为LMS48通道采集前端。

(二)对采集的各声压信号进行傅立叶变换，每段信号的时间长度为1s。信号的窗函数为汉宁窗函数，采样频率为4096Hz。

(三)利用统计最优近场声全息计算声源映射在传声器阵列面上的二维坐标。

声源重建网格设定为规则网格阵列，示意图如图2所示。从图中可知，由于实际网格与声源真实位置并不能完全重合，导致实际投影坐标带有误差，且这种误差无法通过加密网格消除。

统计最优近场声全息的核心思想为：

假设p(rhn)为全息面上第n个测量点rhn＝(xn,yn,zh)处的复声压，p(r)为0≤z≤zh面上任意重建点r＝(x,y,z)处的复声压，则令：

式中，cn(r)为权重系数，N为测量点总数。

若上式成立，那么该假设对于单元平面波也成立，且约定在两种假设中的叠加系数不变，故权重系数cn(r)应满足：

式中，ΦK(r)为空间波数域的单元平面波，Km为波数矢量。

由上式确定的M(M≥N)个线性方程所构成的线性方程组表示成矩阵的形式为：

b＝Ac(r)

式中，c(r)为权重系数矩阵，b为单元平面波矩阵，A为空间传递函数矩阵。

上式的正则化解为：

c(r)＝(A⁺A+θ²I)^-1A⁺b

其中：A⁺为矩阵A的共轭转置矩阵；θ为正则化参数；I为单位对角矩阵。

将式(4)代入式(1)中可得：

(四)将传声器阵列沿着垂直轴旋转1度，形成新的测量面。利用新采集的数据及统计最优近场声全息计算声源相对于新测量面的二维坐标，具体方法如步骤(三)中所述。

将起始测量面b与终止测量面a的夹角(α)进行n等分，等分角度α1＝α/n，将b全息面沿y轴旋转角度α1，定义为α1全息面(即：下次旋转的起始全息面)。利用统计最优近场声全息确定声源几何中心映射在b全息面、α1上的2D坐标；该实例中α设定为60度，n定为60。

(五)重复步骤(三)和步骤(四)，获得一组二维坐标数据和旋转角度数据，利用余弦定理计算声源相对于传声器阵列面的虚拟运动轨迹坐标，此三维坐标数据即为声源位置数据，即：

声源到全息面距离的计算公式分两种情况：

1)投影映射点在转轴的两侧，参见图3(a)，其计算公式如下列式子所示：

γ＝180-α-β

2)投影映射点在转轴的同侧，参见图3(b)，其计算公式如下列式子所示：

β＝90-α-θ

γ＝90+θ

式中，Xab为b全息面的映射点到a全息面的映射点之间的距离，Xa为a全息面的映射点到a全息面与b全息面的交叉点之间的距离，Xb为b全息面的映射点到b全息面与a全息面的交叉点之间的距离。Da为声源S到a全息面的距离，Db为声源S到b全息面的距离，β为Da与Xab之间的夹角，γ为Db与Xab之间的夹角，α为Da与Db之间的夹角。

(六)利用卡尔曼滤波对该虚拟轨迹进行估计，由于测量误差和网格误差，上述步骤中获得的虚拟运动轨迹不够精确，误差很大。本步骤可以通过卡尔曼滤波的估计功能，抑制上述误差，得到更加精确的声源虚拟轨迹线，最终实现声源的三维定位。

基于旋转全息面的数据，采用数据平均算法初步估算出虚拟运动轨迹的半径R(参见图4)，则声源在xoz平面上的位移为：

其中：虚拟轨迹逆时针转动时，n＝1；虚拟轨迹顺时针转动时，n＝0。

本实施例针对真实路径进行了对比实验，图5为真实路径、原始定位数据及采用本发明方法滤波后的定位数据的对比图，由图中可以看出：1)真实轨迹(虚线)是一条圆弧线；2)由于数据存在误差，导致基于该数据的定位结果(实线)偏离真实轨迹，最大误差达到了300％；3)本发明的定位结果(星号线)明显好于已有的算法，最大误差不超过10％，实验结果验证的了本发明的有效性。

以上所述，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的工作人员在本发明揭露的技术范围内，可轻易想到各种等效的修改或替换，这些修改或替换都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应以权利要求的保护范围为准。