CN101645701A - 一种基于滤波器组的时延估计方法及其系统 - Google Patents

Abstract

本发明涉及一种基于滤波器组的时延估计的方法及其系统，利用滤波器组分别估计两路目标信号的互谱和相应的加权谱，以实现时延估计，1)采集两路目标信号的采样序列x₁(n)，x₂(n)；2)预先生成L×N维傅里叶矩阵F＝[f₀，f₁，…，f_N－1]；3)分别计算协方差矩阵R₁₁，R₂₂和相应的各次逆矩阵R₁₁ ^－1，R₂₂ ^－1，以及互协方差矩阵R₁₂；4)计算得到N个长度为L相应数据子段的滤波器组；5)计算序列x₁(n)，x₂(n)的互谱

和加权函数谱

；6)计算采用广义互谱相位法的时延估计。本方案中还可以进行数据分段和加窗，目的是保证数据序列的平稳性，提高估计的稳健度。本发明原理明确，易于工程实现；不仅适合窄带和功率谱平坦的情况，也适合功率谱有较大起伏的情况，在信道起伏导致功率谱不平坦的情况下，仍可实现稳健的时延估计。

Description

一种基于滤波器组的时延估计方法及其系统

技术领域

本发明属于阵列信号处理领域，特别涉及一种基于滤波器组的时延估计方法及其系统。

背景技术

阵列信号处理领域中面临的一个基本问题就是时延估计，即确定两路传感器接收信号的时延差，该时延差与传感器间距，目标的方位、距离，传播环境密切相关。时延估计在声纳、雷达、通信、地震勘探和生物医学信号处理等领域中有着广泛的应用。

目前典型的时延估计方法有：广义互相关法；相位谱法；最小均方自适应滤波法等。

其中，广义互相关法(可参考文献：“Knapp C H等，The generalized correlationmethod for estimation of time delay，IEEE Transaction on Acoustics，Speech and SignalProcessing，24(4)：320～327，1976”)，即利用两路有时延差的信号的互相关函数在时延差位置处最大的特点来进行时延估计。该方法的特点是简单，易于实现，但受噪声影响较大，低信噪比下测量精度较差且不稳定。

相位谱法(可参考文献：“Piesol等，Time Delay Estimation Using Phase Data，IEEE Transaction on Acoustics，Speech and Signal Processing，，29：，471～477，1981”)，通过计算两路时延信号的互功率谱的相位谱的斜率来进行时延估计。该方法基于快速傅里叶变换(FFT)运算，处理较简单，但当功率谱有较大起伏而又没有考虑频谱加权时，性能将显著下降。

最小均方自适应滤波法(可参考文献：“吴慧娟等，低信噪比下的LMS自适应无偏时延估计，电子学报，37(3)：500～505，2009”)，通过设定迭代初值、参数和自适应学习来估计时延，不过该方法需要假设通道间的背景噪声为不相关的高斯白噪声，与实际环境未必一致。

除此以外，近年来，由于现代信号处理理论的迅速发展，高阶累积量、小波分析、循环平稳信号分析与处理等理论和方法在时延估计中也得到了广泛的研究。

总体而言，由于水声工作环境的恶劣性(噪声模型易失配，海洋传播信道的复杂，阵形易失配等)，以上算法在现阶段的实际应用中，效果不是很理想，实际应用中迫切需要一种稳健的时延估计方法。

发明内容

本发明的目的在于，为了克服现有技术的不足，解决现有技术中时延估计稳健性差、难以实用的问题，从而提供一种基于滤波器组的时延估计的方法及其系统。

为了实现上述目的，本发明提供的一种基于滤波器组的时延估计方法，该方法利用滤波器组分别估计两路目标信号的互谱和相应的加权谱，以实现时延估计，所述的基于滤波器组的时延估计方法包括以下步骤：

1)采集两路目标信号的采样序列x₁(n)，x₂(n)；所述的两路目标信号的到达接受传感器时存在时延差，即x₁(n)，x₂(n)之间存在时延τ，采样序列x₁(n)，x₂(n)的采样频率为F_s；

2)预先生成L×N维傅里叶矩阵F＝[f₀，f₁，…，f_N-1]，傅里叶矩阵的列向量 $f_n=\frac{1}{\sqrt{L}}{\left[\begin{array}{cccc}1& e^{i{\mathrm{\ω}}_n}& \·\·\·& e^{i{\mathrm{\ω}}_n(L-1)}\end{array}\right]}^T,$ ω_n＝2πn/N，其中，n＝0，1，…，N-1；Fs/N为频率分辨率；

3)分别计算两路目标信号的采样序列x₁(n)，x₂(n)的各子段相应的协方差矩阵R₁₁，R₂₂和相应的各次逆矩阵R₁₁ ^-1，R₂₂ ^-1，以及两路目标信号的采样序列x₁(n)和x₂(n)各子段相应的互协方差矩阵R₁₂；所述协方差矩阵R₁₁，R₂₂和互协方差矩阵R₁₂的计算公式为： $R_{11}=E\left\{x_1\left(n\right)x_1^H\left(n\right)\right\},$ $R_{22}=E\left\{x_2\left(n\right)x_2^H\left(n\right)\right\},$ $R_{12}=E\left\{x_1\left(n\right)x_2^H\left(n\right)\right\};$ 其中上标H表示共轭转置操作，E{·}表示求数学期望；

4)设定参数α，由设定参数α、协方差矩阵估计及其逆矩阵，得到相应的N个长度为L相应数据子段的滤波器组： $g_{1,i}^{\mathrm{\α}}=R_{11}^{-\mathrm{\α}}{f}_i/\left(f_i^H{R}_{11}^{-\mathrm{\α}}{f}_i\right),$ $g_{2,i}^{\mathrm{\α}}=R_{22}^{-\mathrm{\α}}{f}_i/\left(f_i^H{R}_{22}^{-\mathrm{\α}}{f}_i\right),$ 其中i＝0，1，…，N-1；f_i为频率点；所述参数α的取值范围为0≤α≤1；所述滤波器g_1，i ^α和g_2，i ^α分别表示第p个子段数据对应的滤波器系数；

5)计算序列x₁(n)，x₂(n)的互谱 $P_{12}^p\left(f_i\right)=g_{1,i}^H{R}_{12}{g}_{2,i}$ 和加权函数谱

6)计算采用广义互谱相位法的时延估计。

上述技术方案中，所述步骤1)中，所述的采样序列x₁(n)，x₂(n)可以是两个传感器的时间采样序列，也可以是由多个传感器采样数据进行波束形成所得到的波束序列。

上述技术方案中，所述步骤3)具体为：

31)首先，将采样序列x₁(n)，x₂(n)划分为P个相互重叠的长度为L的子段，重叠的长度一般在数据长度的一半以上，即重叠数据长度≥L/2；

32)然后，对每一个子段数据进行加窗处理后得到新的子段数据；所述的窗函数为汉宁窗或矩形窗；本方案中的数据分段和加窗是信号处理中常用的预处理手段，目的是保证数据序列的平稳性，提高估计的稳健度。

33)分别计算各子段的协方差矩阵和互协方差矩阵，然后平均得到协方差矩阵R₁₁，R₂₂和互协方差矩阵R₁₂。对应于该技术方案中，所述步骤6)具体为：

61)计算平均互谱的相位

其中Re[P₁₂(f)]和Im[P₁₂(f)]分别表示复数P₁₂(f)的实部和虚部；

62)由相位谱和加权谱计算时延估计

其中k_min，k_max分别对应处理带宽的起止频率。

上述技术方案中的另一种情况是，所述步骤3)具体为：

31)首先，将采样序列x₁(n)，x₂(n)划分为P个相互重叠的长度为L的子段，重叠的长度一般在数据长度的一半以上，即重叠数据长度≥L/2；

32)然后，对每一个子段数据进行加窗处理后得到新的子段数据；所述的窗函数为汉宁窗或矩形窗；

33)分别计算各子段的协方差矩阵和互协方差矩阵，得到协方差矩阵R₁₁，R₂₂和互协方差矩阵R₁₂。对应于该技术方案中，所述步骤6)具体为：

61′)计算互谱的相位其中Re[P₁₂(f)]和Im[P₁₂(f)]分别表示复数P₁₂(f)的实部和虚部；

62′)由相位谱和加权谱计算时延估计其中k_min，k_max分别对应处理带宽的起止频率；

63′)对各子段的时延估计平均后得到平均的时延估计。

上述技术方案中，所述步骤5)中，所述加权函数谱

为平均加权 ${\hat{W}}_a\left(f_i\right)=1,$ 或互谱幅度平方加权 ${\hat{W}}_{\mathrm{MS}}\left(f_i\right)={\left|P_{12}\left(f_i\right)\right|}^2$ 或最优加权 ${\hat{W}}_{\mathrm{MSC}}^p\left(f_i\right)=\frac{{\left|f_i^H{R}_{x_1{x}_1}^{-p}{R}_{x_1{x}_2}{R}_{x_2{x}_2}^{-p}{f}_i\right|}^2}{\left[f_i^H{R}_{x_1{x}_1}^{1-2p}{f}_i\right]\left[f_i^H{R}_{x_2{x}_2}^{1-2p}{f}_i\right]-{\left|f_i^H{R}_{x_1{x}_1}^{-p}{R}_{x_1{x}_2}{R}_{x_2{x}_2}^{-p}{f}_i\right|}^2}.$

为了实现本发明的另一目的，本发明还提供了一种基于滤波器组的时延估计系统，该系统利用滤波器组分别估计两路目标信号的互谱和相应的加权谱，以实现时延估计，所述的基于滤波器组的时延估计系统包括：

数据采集及划分模块，用于数据采集并划分P个相互重叠的长度为L的数据子段；

傅里叶矩阵产生模块，用于产生傅里叶矩阵，预先生成L×N维傅里叶矩阵F＝[f₀，f₁，…，f_N-1]，傅里叶矩阵的列向量 $f_n=\frac{1}{\sqrt{L}}{\left[\begin{array}{cccc}1& e^{i{\mathrm{\ω}}_n}& \·\·\·& e^{i{\mathrm{\ω}}_n(L-1)}\end{array}\right]}^T,$ ω_n＝2πn/N，其中，n＝0，1，…，N-1；Fs/N为频率分辨率；

滤波器组产生模块，用于依据计算得到的协方差和互协方差产生处理各子段数据的滤波器组；

互谱计算模块，用于计算序列x₁(n)，x₂(n)的互谱 $P_{12}^p\left(f_i\right)=g_{1,i}^H{R}_{12}{g}_{2,i}$ 和加权函数谱

以及

时延估计模块，用于根据获得的互谱和加权函数计算采用广义互谱相位法的估计时延

上述技术方案中，所述数据采集及划分模块包括：

数据采集单元，用于采集两路目标信号的采样序列x₁(n)，x₂(n)；所述的两路目标信号的到达存在时延差，即x₁(n)，x₂(n)之间存在时延τ，采样序列x₁(n)，x₂(n)的采样频率为F_s；

数据划分单元，用于将采样序列x₁(n)，x₂(n)划分为P个相互重叠的长度为L的子段，重叠的长度一般在数据长度的一半以上，即重叠数据长度≥L/2；以及

加窗处理单元，用于产生窗函数，对每一个子段数据进行加窗处理后得到新的子段数据；所述的窗函数为汉宁窗或矩形窗。

上述技术方案中，所述滤波器组产生模块包括：

协方差和互协方差计算单元，用于分别计算两路目标信号的采样序列x₁(n)，x₂(n)各自的协方差矩阵R₁₁，R₂₂和相应的各次逆矩阵R₁₁ ^-1，R₂₂ ^-1，以及两路目标信号的采样序列x₁(n)和x₂(n)的互协方差矩阵R₁₂；所述协方差矩阵R₁₁，R₂₂和互协方差矩阵R₁₂的计算公式为： $R_{11}=E\left\{x_1\left(n\right)x_1^H\left(n\right)\right\},$ $R_{22}=E\left\{x_2\left(n\right)x_2^H\left(n\right)\right\},$ $R_{12}=E\left\{x_1\left(n\right)x_2^H\left(n\right)\right\};$ 其中上标H表示共轭转置操作，E{·}表示求数学期望；以及

滤波器产生单元，用于根据设定参数α、协方差矩阵估计及其逆矩阵，得到相应的N个长度为L相应数据子段的滤波器组： $g_{1,i}^{\mathrm{\α}}=R_{11}^{-\mathrm{\α}}{f}_i/\left(f_i^H{R}_{11}^{-\mathrm{\α}}{f}_i\right),$ $g_{2,i}^{\mathrm{\α}}=R_{22}^{-\mathrm{\α}}{f}_i/\left(f_i^H{R}_{22}^{-\mathrm{\α}}{f}_i\right),$ 其中i＝0，1，…，N-1；f_i为频率点；所述参数α的取值范围为0≤α≤1；所述滤波器g_1，i ^α和g_2，i ^α分别表示第p个子段数据对应的滤波器系数。

上述技术方案中，所述互谱计算模块包括：

子段互谱计算单元，用于计算各数据子段的互谱 $P_{12}^p\left(f_i\right)=g_{1,i}^H{R}_{12}{g}_{2,i}$ 和互谱幅度平方函数；以及

平均互谱计算单元，用于平均各数据子段的互谱及互谱幅度平方函数。

与现有技术相比，本发明具有如下技术效果：

(1)原理明确，易于工程实现；

(2)不仅适合窄带和功率谱平坦的情况，也适合功率谱有较大起伏的情况，因此在信道起伏导致功率谱不平坦的情况下，仍可实现稳健的时延估计。

附图说明

以下，结合附图来详细说明本发明的实施例，其中：

图1是现有技术中的常规线列阵示意图；

图2是时延估计信号模型示意图；

图3是本发明的总体流程图；

图4是本发明的具体流程图；

图5是本发明对10kHz采样的时延0.5ms的信号进行100次蒙特卡洛仿真后得到的均方根误差示意图；

图6是本发明实施例的声纳数据处理和显控系统；

具体实施方式

本发明的基本构思为：利用滤波器组的方法，分别估计两路信号的互谱和相应的加权谱，进而实现稳健的时延估计。本发明的总体流程如图3所示。

下面结合附图和具体实施方式对本发明做进一步的详细描述：

实施例1

如图4所示，本实施例的2个水听器的基于滤波器组的时延估计方法及其系统的具体步骤如下：

步骤401：用2个水听器接收时空二维信号，从而获得2个数据序列x₁(n)和x₂(n)。

本步骤中的2个数据序列也可以是常规线列阵等声纳波束形成后得到的2个波束数据序列。假设可处理的数据长度K＝1000。本实施例中为两个传感器在Fs＝30kHz下采样的数据序列。

步骤402：预先产生傅里叶矩阵。

本步骤中傅里叶矩阵F＝[f₀，f₁，…，f_N-1]的维数为L×N维，其列向量 $f_n=\frac{1}{\sqrt{L}}{\left[\begin{array}{cccc}1& e^{i{\mathrm{\ω}}_n}& \·\·\·& e^{i{\mathrm{\ω}}_n(L-1)}\end{array}\right]}^T,$ ω_n＝2πn/N，其中，n＝0，1，…，N-1。本实施例中N＝100，即频率分辨率为100Hz。

步骤403：对所得数据进行分段，并加窗处理。

本步骤中的数据分段和加窗是信号处理中常用的预处理手段，目的是保证数据序列的平稳性，提高估计的稳健度。将序列x₁(n)，x₂(n)划分为P个相互重叠的长度为L的子段，重叠的长度一般在数据长度的一半以上，即重叠数据长度≥L/2。对每一个子段数据进行加窗处理后得到新的子段数据；常用的窗函数为汉宁窗，也可以为其他形式的窗函数如矩形窗等；本实施例中L＝100，P＝19，重叠一半数据，使用汉宁窗。

步骤404：计算平均协方差矩阵及其逆矩阵。

本步骤中序列x₁(n)协方差矩阵R₁₁和序列x₂(n)协方差矩阵R₂₂得定义分别为 $R_{11}=E\left\{x_1\left(n\right)x_1^H\left(n\right)\right\},$ $R_{22}=E\left\{x_2\left(n\right)x_2^H\left(n\right)\right\}.$ 对协方差矩阵的估计有多种方式，本实施例中协方差矩阵的估计为： $R_{11}=x_{1,p}\left(n\right)x_{1,p}^H\left(n\right)/L,$ $R_{22}=x_{2,p}\left(n\right)x_{2,p}^H\left(n\right)/L,$ 其中x_1，p(n)和x_2，p(n)分别表示2个序列的第p个子段数据。对求得的P个协方差矩阵进行平均就得到平均协方差矩阵，然后可得到相应的逆矩阵R₁₁ ^-1，R₂₂ ^-1。

步骤405：计算平均互协方差矩阵。

本步骤中序列x₁(n)和序列x₂(n)的互协方差矩阵R₁₂的定义为： $R_{12}=E\left\{x_1\left(n\right)x_2^H\left(n\right)\right\}.$ 对互协方差矩阵的估计有多种方式，本实施例中互协方差矩阵的估计为： $R_{12}=x_{1,p}\left(n\right)x_{2,p}^H\left(n\right)/L,$ 对求得的P个互协方差矩阵进行平均就得到平均协方差矩阵。

步骤406：设计滤波器组。

本步骤中首先确定参数α的值。本实施例中为了便于比较分别设定了α的值为：0，0.5和1。实际中只需设定一个值即可。由参数α、协方差矩阵估计及其逆矩阵，设计相应的N个长度为L的滤波器组为： $g_{1,i}^{\mathrm{\α}}=R_{11}^{-\mathrm{\α}}{f}_i/\left(f_i^H{R}_{11}^{-\mathrm{\α}}{f}_i\right),$ $g_{2,i}^{\mathrm{\α}}=R_{22}^{-\mathrm{\α}}{f}_i/\left(f_i^H{R}_{22}^{-\mathrm{\α}}{f}_i\right),$ 其中i＝0，1，…，N-1。

步骤407：计算互谱。

本步骤中互谱的计算公式为 $P_{12}^p\left(f_i\right)=g_{1,i}^H{R}_{12}{g}_{2,i}.$

步骤408：计算互谱的相位谱。

本步骤中相位谱其中和

分别表示P₁₂(f)的实部和虚部。

步骤409：计算加权谱

本步骤中加权谱

可以为平均加权 ${\hat{W}}_a\left(f_i\right)=1,$ 互谱幅度平方加权 ${\hat{W}}_{\mathrm{MS}}\left(f_i\right)={\left|P_{12}\left(f_i\right)\right|}^2,$ 最优加权 ${\hat{W}}_{\mathrm{MSC}}^p\left(f_i\right)=\frac{{\left|f_i^H{R}_{x_1{x}_1}^{-p}{R}_{x_1{x}_2}{R}_{x_2{x}_2}^{-p}{f}_i\right|}^2}{\left[f_i^H{R}_{x_1{x}_1}^{1-2p}{f}_i\right]\left[f_i^H{R}_{x_2{x}_2}^{1-2p}{f}_i\right]-{\left|f_i^H{R}_{x_1{x}_1}^{-p}{R}_{x_1{x}_2}{R}_{x_2{x}_2}^{-p}{f}_i\right|}^2};$ 本实施例中分别采用了平均加权和最优加权进行了比较。

步骤410：时延估计。

本步骤中有前述所得的相位谱和加权谱得到相应的时延估计为：

其中k_min，k_max分别对应处理带宽的起止频率。

图5是本发明对延时τ＝0.5ms，处理带宽为100Hz～1000Hz的情况下的时延估计均方根误差在各种信噪比下的示意图。可以看出在功率谱不平坦的情况下，最优加权相较于平均加权具有更稳健的性能。

下面对本实施例的硬件环境的一个示例进行描述。

本实施例可以直接使用现有技术中的常规线列阵，如图1所示，用于接收空间信号的线列阵声纳102的线阵是由32个市场所售常规的接收信号中心频率为5000KHz的无指向性水听器101组成，两两水听器间的间隔是半个波长，从而线阵声振段的总长度是4.65米。线列阵102安装在潜艇或无人潜航器上。

图6是本发明实施例的声纳数据处理和显控系统。水听器102和前置电路601从每个水听器接收到模拟信号。前置电路601包括前放、滤波和其它常规的电路。每个通道的模拟信号输入到A/D转换器602得到数字信号。从A/D转换器602出来的是多通道的数字数据流，每个通道的数据流对应于一个水听器接收到的模拟信号。将这些数据流输入到微型处理器603。经微型处理器603处理后的输出信息可以存储在数据存储器件605，比如磁盘存储设备中，或直接输出到显示设备606上显示。

微型处理器603首先将接收到的数据流存储到动态存取区604，在输入满足处理要求数量的数据流后就开始处理。图4是算法具体流程示意图，它包括数据采集401，傅里叶矩阵生成402，数据分段及加窗403，协方差矩阵及其逆矩阵计算404，互协方差矩阵计算405，滤波器组设计406，互谱估计407，互谱相位谱估计408，加权谱估计409和时延估计410。这些程序储存在动态存取区604中。

因为本发明中的A/D转换器602输出的是多通道的数据流，因此可以采用多片的微型处理器来并行处理。能够实现图6功能的其它一些硬件设备，比如专用硬件、基于应用的集成电路(ASIC)、DSP、ARM等都可以用来代替微型处理器603。

最后所应说明的是，以上实施例仅用以说明本发明的技术方案而非限制。尽管参照实施例对本发明进行了详细说明，本领域的普通技术人员应当理解，对本发明的技术方案进行修改或者等同替换，都不脱离本发明技术方案的精神和范围，其均应涵盖在本发明的权利要求范围当中。

Claims (11)

1、一种基于滤波器组的时延估计方法，该方法利用滤波器组分别估计两路目标信号的互谱和相应的加权谱，以实现时延估计，所述的基于滤波器组的时延估计方法包括以下步骤：

1)采集两路目标信号的采样序列x₁(n)，x₂(n)；所述的两路目标信号的到达接受传感器时存在时延差，即x₁(n)，x₂(n)之间存在时延τ，采样序列x₁(n)，x₂(n)的采样频率为F_s；

2)预先生成L×N维傅里叶矩阵F＝[f₀，f₁，…，f_N-1]，傅里叶矩阵的列向量 $f_n=\frac{1}{\sqrt{L}}{\left[\begin{array}{cccc}1& e^{i{\mathrm{\ω}}_n}& \·\·\·& e^{i{\mathrm{\ω}}_n(L-1)}\end{array}\right]}^T,$ ω_n＝2πn/N，其中，n＝0，1，…，N-1；Fs/N为频率分辨率；

3)分别计算两路目标信号的采样序列x₁(n)，x₂(n)的各子段相应的协方差矩阵R₁₁，R₂₂和相应的各次逆矩阵R₁₁ ^-1，R₂₂ ^-1，以及两路目标信号的采样序列x₁(n)和x₂(n)各子段相应的互协方差矩阵R₁₂；所述协方差矩阵R₁₁，R₂₂和互协方差矩阵R₁₂的计算公式为： $R_{11}=E\left\{x_1\left(n\right)x_1^H\left(n\right)\right\},$ $R_{22}=E\left\{x_2\left(n\right)x_2^H\left(n\right)\right\},$ $R_{12}=E\left\{x_1\left(n\right)x_2^H\left(n\right)\right\};$ 其中上标H表示共轭转置操作，E{·}表示求数学期望；

4)设定参数α，由设定参数α、协方差矩阵估计及其逆矩阵，得到相应的N个长度为L相应数据子段的滤波器组： $g_{1,i}^{\mathrm{\α}}=R_{11}^{-\mathrm{\α}}{f}_i/\left(f_i^H{R}_{11}^{-\mathrm{\α}}{f}_i\right),$ $g_{2,i}^{\mathrm{\α}}=R_{22}^{-\mathrm{\α}}{f}_i/\left(f_i^H{R}_{22}^{-\mathrm{\α}}{f}_i\right),$

其中i＝0，1，…，N-1；f_i为频率点；所述参数α的取值范围为0≤α≤1；所述滤波器g_1，i ^α和g_2，i ^α分别表示第p个子段数据对应的滤波器系数；

5)计算序列x₁(n)，x₂(n)的互谱 $P_{12}^p\left(f_i\right)=g_{1,i}^H{R}_{12}{g}_{2,i}$ 和加权函数谱

6)计算采用广义互谱相位法的时延估计。

2、根据权利要求1所述的基于滤波器组的时延估计方法，其特征在于，所述步骤1)中，所述的采样序列x₁(n)，x₂(n)可以是两个传感器的时间采样序列，也可以是由多个传感器采样数据进行波束形成所得到的波束序列。

3、根据权利要求1所述的基于滤波器组的时延估计方法，其特征在于，所述步骤3)具体为：

31)首先，将采样序列x₁(n)，x₂(n)划分为P个相互重叠的长度为L的子段，重叠的长度一般在数据长度的一半以上，即重叠数据长度≥L/2；

32)然后，对每一个子段数据进行加窗处理后得到新的子段数据；所述的窗函数为汉宁窗或矩形窗；

33)分别计算各子段的协方差矩阵和互协方差矩阵，然后平均得到协方差矩阵R₁₁，R₂₂和互协方差矩阵R₁₂。

4、根据权利要求3所述的基于滤波器组的时延估计方法，其特征在于，所述步骤6)具体为：

61)计算平均互谱的相位

其中Re[P₁₂(f)]和Im[P₁₂(f)]分别表示复数P₁₂(f)的实部和虚部；

62)由相位谱和加权谱计算时延估计

其中k_min，k_max分别对应处理带宽的起止频率。

5、根据权利要求1所述的基于滤波器组的时延估计方法，其特征在于，所述步骤3)具体为：

31)首先，将采样序列x₁(n)，x₂(n)划分为P个相互重叠的长度为L的子段，重叠的长度一般在数据长度的一半以上，即重叠数据长度≥L/2；

32)然后，对每一个子段数据进行加窗处理后得到新的子段数据；所述的窗函数为汉宁窗或矩形窗；

33)分别计算各子段的协方差矩阵和互协方差矩阵，得到协方差矩阵R₁₁，R₂₂和互协方差矩阵R₁₂。

6、根据权利要求5所述的基于滤波器组的时延估计方法，其特征在于，所述步骤6)具体为：

61′)计算互谱的相位

其中Re[P₁₂(f)]和Im[P₁₂(f)]分别表示复数P₁₂(f)的实部和虚部；

62′)由相位谱和加权谱计算时延估计其中k_min，k_max分别对应处理带宽的起止频率；

63′)对各子段的时延估计平均后得到平均的时延估计。

7、根据权利要求1所述的基于滤波器组的时延估计方法，其特征在于，所述步骤5)中，所述加权函数谱为平均加权 ${\hat{W}}_a\left(f_i\right)=1,$ 或互谱幅度平方加权 ${\hat{W}}_{\mathrm{MS}}\left(f_i\right)={\left|P_{12}\left(f_i\right)\right|}^2,$ 或最优加权 ${\hat{W}}_{\mathrm{MSC}}^p\left(f_i\right)=\frac{{\left|f_i^H{R}_{x_1{x}_1}^{-p}{R}_{x_1{x}_2}{R}_{x_2{x}_2}^{-p}{f}_i\right|}^2}{\left[f_i^H{R}_{x_1{x}_1}^{1-2p}{f}_i\right]\left[f_i^H{R}_{x_2{x}_2}^{1-2p}{f}_i\right]-{\left|f_i^H{R}_{x_1{x}_1}^{-p}{R}_{x_1{x}_2}{R}_{x_2{x}_2}^{-p}{f}_i\right|}^2}.$

8、一种基于滤波器组的时延估计系统，该系统利用滤波器组分别估计两路目标信号的互谱和相应的加权谱，以实现时延估计，所述的基于滤波器组的时延估计系统包括：

数据采集及划分模块，用于数据采集并划分P个相互重叠的长度为L的数据子段；

傅里叶矩阵产生模块，用于产生傅里叶矩阵，预先生成L×N维傅里叶矩阵F＝[f₀，f₁，…，f_N-1]，傅里叶矩阵的列向量 $f_n=\frac{1}{\sqrt{L}}{\left[\begin{array}{cccc}1& e^{i{\mathrm{\ω}}_n}& \·\·\·& e^{i{\mathrm{\ω}}_n(L-1)}\end{array}\right]}^T,$ ω_n＝2πn/N，其中，n＝0，1，…，N-1；Fs/N为频率分辨率；

滤波器组产生模块，用于依据计算得到的协方差和互协方差产生处理各子段数据的滤波器组；

互谱计算模块，用于计算序列x₁(n)，x₂(n)的互谱 $P_{12}^p\left(f_i\right)=g_{1,i}^H{R}_{12}{g}_{2,i}$ 和加权函数谱

；以及

时延估计模块，用于根据获得的互谱和加权函数计算采用广义互谱相位法的估计时延

9、根据权利要求

所述的基于滤波器组的时延估计系统，其特征在于，所述数据采集及划分模块包括：

数据采集单元，用于采集两路目标信号的采样序列x₁(n)，x₂(n)；所述的两路目标信号的到达存在时延差，即x₁(n)，x₂(n)之间存在时延τ，采样序列x₁(n)，x₂(n)的采样频率为F_s；

数据划分单元，用于将采样序列x₁(n)，x₂(n)划分为P个相互重叠的长度为L的子段，重叠的长度一般在数据长度的一半以上，即重叠数据长度≥L/2；以及

加窗处理单元，用于产生窗函数，对每一个子段数据进行加窗处理后得到新的子段数据；所述的窗函数为汉宁窗或矩形窗。

10、根据权利要求所述的基于滤波器组的时延估计系统，其特征在于，所述滤波器组产生模块包括：

协方差和互协方差计算单元，用于分别计算两路目标信号的采样序列x₁(n)，x₂(n)各自的协方差矩阵R₁₁，R₂₂和相应的各次逆矩阵R₁₁ ^-1，R₂₂ ^-1，以及两路目标信号的采样序列x₁(n)和x₂(n)的互协方差矩阵R₁₂；所述协方差矩阵R₁₁，R₂₂和互协方差矩阵R₁₂的计算公式为： $R_{11}=E\left\{x_1\left(n\right)x_1^H\left(n\right)\right\},$ $R_{22}=E\left\{x_2\left(n\right)x_2^H\left(n\right)\right\},$ $R_{12}=E\left\{x_1\left(n\right)x_2^H\left(n\right)\right\};$ 其中上标H表示共轭转置操作，E{·}表示求数学期望；以及滤波器产生单元，用于根据设定参数α、协方差矩阵估计及其逆矩阵，得到相应的N个长度为L相应数据子段的滤波器组： $g_{1,i}^{\mathrm{\α}}=R_{11}^{-\mathrm{\α}}{f}_i/\left(f_i^H{R}_{11}^{-\mathrm{\α}}{f}_i\right),$ $g_{2,i}^{\mathrm{\α}}=R_{22}^{-\mathrm{\α}}{f}_i/\left(f_i^H{R}_{22}^{-\mathrm{\α}}{f}_i\right),$ 其中i＝0，1，…，N-1；f_i为频率点；所述参数α的取值范围为0≤α≤1；所述滤波器g_1，i ^α和g_2，i ^α分别表示第p个子段数据对应的滤波器系数。

11、根据权利要求

所述的基于滤波器组的时延估计系统，其特征在于，所述互谱计算模块包括：

子段互谱计算单元，用于计算各数据子段的互谱 $P_{12}^p\left(f_i\right)=g_{1,i}^H{R}_{12}{g}_{2,i}$ 和互谱幅度平方函数；以及

平均互谱计算单元，用于平均各数据子段的互谱及互谱幅度平方函数。

Images