基于动态规划算法的移动声源到达时延估计方法与流程

本发明涉及阵列信号领域，特别是涉及麦克风阵列信号处理和到达时延估计。

背景技术：

本发明的背景基于实际需要产生。近年来，在无人机迅速成为研究热点时，也带来一系列问题，比如无人机黑飞，严重影响区域安全。因此无人机防御正在成为各国政府和军方重点关注的新领域。无人机的声音具有明显的特征，可以有效侦别出空中飞行的物体。通过架设多个麦克风阵列，能有效地收集到无人机的信号，但是实际环境存在的噪声、信号在电路设备传输过程中的噪声以及一些其他的混响、多径干扰等，收集到的声音信号非常嘈杂，导致对无人机进行定位时，估计信号到达麦克风时延(tdoa)会出现比较大的误差，为了提高tdoa估计的准确性，可以在估计tdoa时进行改进。

在目前的研究内容中，一般利用上一时刻的tdoa对本时刻tdoa估计进行修正，但这种方法仅依赖上一个时刻的tdoa，冗余性较差，容易出现误差。因此，急需一种新的方法结合前多个时刻的tdoa来估计本时刻tdoa，同时能够保证较快的运行速度和系统较低的成本。而本发明能有效的计算出tdoa，提升后续对信号分析的准确性。

技术实现要素：

为了实现对无人机移动时声音信号时延的跟踪，本发明采用麦克风阵列传感器对空中无人机声音信号进行处理，可以有效计算移动中的无人机声音信号到达麦克风tdoa。

本发明解决其技术问题所采用的技术方案是：一种基于动态规划算法的移动声源到达时延估计方法，包括如下步骤：

(1)根据广义互相关函数计算麦克风阵列t时刻采集的两组两组时域声信号x1(t)、x2(t)互功率谱密度分布函数

其中为x1(t)和x2(t)的傅立叶变换结果的乘积，为x1(t)和x2(t)的频域滤波器；

(2)提取出函数连续n个时刻中每个时刻最大的m个峰值对应的m个时延。

(3)根据移动声源时延估计的连贯性，利用如下公式找到n个时刻时延相连曲线中波动最小的曲线，并计算到达麦克风时延tdoa：

其中th代表在第h时刻的时延，th取自集合[th,1,th,2,...,th,m]，th,i代表第h时刻的第i个最大峰值对应的时延；ph代表th所取时延对应的互功率谱密度分布函数的峰值大小；a用来将ph与th的大小调至同一个数量级。

(4)通过动态规划算法求解公式(2)，得到第n个时刻的最终tdoa。

(5)利用滑动时间窗求得之后每个时刻的tdoa。

进一步地，步骤(1)中，采用相位变换(phat)加权函数

进一步地，步骤(4)中，利用如下动态规划算法求解公式(2)：

其中i,v＝1,…,m；h＝2,...,n；表示第h时刻的第i个计算中间变量，这里初值的取值设定为第一个时刻的互功率谱密度分布函数的第v个峰值的大小，迭代公式(3)得到最后结果

利用如下公式计算最终的tdoa：

本发明提出的基于动态规划算法的移动声源到达时延估计方法，可以计算移动中的无人机声音信号到达麦克风时延，具有成本低廉，计算迅速等特点。本发明具有以下优势：

(1)充分融合多时刻声源信号的信息，提高了tdoa估计的精度。

(2)在融合时刻信号时延时，引入了互功率谱密度分布函数的峰值，消除了多径效应对麦克风时延的影响。

(3)充分利用动态规划算法对整个计算过程进行优化，简化了计算步骤，降低了计算复杂度，加快了运算效率。

附图说明

图1为本发明方法流程图；

图2为多时刻多时延散点图；

图3为最小波动时延连线图；

图4为利用错误tdoa初值后后续tdoa跟踪结果；

图5为利用正确tdoa初值后后续tdoa跟踪结果。

具体实施方式

以下结合附图和具体实施例对本发明作进一步详细说明。

本发明提出的一种基于动态规划算法的移动声源到达时延估计方法，包括如下步骤：

(1)根据广义互相关函数计算麦克风阵列t时刻采集的两组两组时域声信号x1(t)、x2(t)互功率谱密度分布函数

其中为x1(t)和x2(t)的傅立叶变换结果的乘积，为x1(t)和x2(t)的频域滤波器；

可以采用相位变换(phat)加权函数

(2)提取出函数连续n个时刻中每个时刻最大的m个峰值对应的m个时延。

(3)根据移动声源时延估计的连贯性，利用如下公式找到n个时刻时延相连曲线中波动最小的曲线，并计算到达麦克风时延tdoa：

其中th代表在第h时刻的时延，th取自集合[th,1,th,2,...,th,m]，th,i代表第h时刻的第i个最大峰值对应的时延；ph代表th所取时延对应的互功率谱密度分布函数的峰值大小；a用来将ph与th的大小调至同一个数量级。

(4)通过动态规划算法求解公式(2)，得到第n个时刻的最终tdoa。

利用如下动态规划算法求解公式(2)：

其中i,v＝1,…,m；h＝2,...,n；表示第h时刻的第i个计算中间变量，这里初值的取值设定为第一个时刻的互功率谱密度分布函数的第v个峰值的大小，迭代公式(3)得到最后结果

利用如下公式计算最终的tdoa：

(5)利用滑动时间窗求得之后每个时刻的tdoa。

附图1是一个实施例的实施流程图，首先由麦克风阵列采集声音信号，先连续采集10秒左右的数据，将每一秒的数据都利用广义互相关算法进行处理，求出每一秒的广义互相关函数的前几个最大峰值对应的时延点，然后将多时刻这样的时延代入上述动态规划模型，求出最小波动时延连线以及tdoa，之后滑动时间窗口计算后续tdoa。

附图2是采集十秒数据，每秒从该秒数据的互功率密度分布函数上提取5个最大峰值对应的时延点，将这些时延点画成散点图。可以看到同一时刻的时延分布有多个值，一般的算法都是直接把最大峰值对应的时延作为tdoa，但这样做，对声源信号的信噪比要求较高，一旦有噪声或者多径效应导致tdoa不出现在最大峰值对应tdoa上，这样估计出的tdoa会有很大误差。

附图3是利用本发明在附图2中十秒的时延点中找到波动最小的时延点连线，黑色连线就是波动最小的时延点连线，从图中看出，在该连线下方几个点也可以连成一条比较小波动的时延点连线，但由于我们提出的公式：

中，加入了ph这个量，使我们选择连线时，优先选择对应互功率密度分布函数上峰值比较大的时延点所连成的线。

附图4和5展示了我们将此方法用于tdoa初值估计的结果。图4是我们利用将tdoa估计的初始值设置为0后，后续利用高斯函数结合gcc算法进行tdoa连续跟踪估计的效果，具体实施步骤如下：

根据上一时刻的时延估计结果以及互功率谱密度分布函数，计算高斯分布先验概率密度分布函数：

其中μ＝tpast，tpast表示上一时刻的时延估计的最终结果，b表示上一时刻互功率谱密度分布函数所对应的值左右两侧零点限定的区间范围；根据计算当前时刻的互功率谱密度分布函数，其中t满足t∈[-d/c,d/c]约束条件，表示当前时刻的互功率谱密度函数，d表示两麦克风之间的间距，c表示声速。可以看到，因为初值是错误的，后续利用高斯算法进行tdoa估计跟踪时，跟踪效果非常差，虽然后续由于跟踪算法不断迭代将tdoa重新跟回正确轨迹，但在这之前有很长一段时间tdoa估计是错误的。附图5展示了我们将tdoa初值设置为该发明提取的tdoa后，后续利用跟踪算法进行tdoa连续跟踪估计的效果，由于初值是正确的，所以这一次的后续tdoa估计都是正确的。

上述实施例用来解释说明本发明，而不是对本发明进行限制，在本发明的精神和权利要求的保护范围内，对本发明作出的任何修改和改变，都落入本发明的保护范围。