一种基于频率预估计的双阵元正弦信号来波方向估计方法与流程

本发明属于信号处理领域，具体来说，涉及一种基于频率预估计的双阵元正弦信号来波方向估计方法。

背景技术：

在非合作条件下，对被噪声污染的双阵元正弦信号的来波方向进行精确估计是是信号处理中的研究热点之一，它在语音信号处理、雷达、声纳以及电子战等领域具有广泛的应用，尤其在水声和电子侦察脉内信号处理中扮演了极其重要的角色。

双阵元正弦信号来波方向估计的实质是将来波方向估计问题转化为目标正弦信号到达两个接收基元的时延差估计问题，其关键是要获得高精度的时延估计。目前，国内外学者提出了很多时延差估计算法，主要分为两大类：互相关法和互谱法。

互相关法的实质是是比较双阵元正弦信号的波形，求出两个接收正弦信号的相对时延差。互相关法实现简单，但受正弦信号波形自身周期性的影响，互相关波形结果不可避免的存在多值性，导致方向估计误差大，工程实用性差。

互谱法的实质是通过求取两阵元正弦信号在信号频点处的相位差来估计信号来波方向，不存在信号波形周期性重复所带来的相位模糊问题，在实际工程中得到了广泛的应用。但常规互谱方法均是基于双阵元正弦信号的离散频谱直接测向，当正弦信号的相对频率偏差较大时，信号频谱泄露严重，对测向精度影响较大。本发明提出的方法是在常规互谱法基础上的优化，其利用双阵元正弦信号互功率谱模值预估信号频率相对偏差和信号频率，依据相对偏差预估计值与门限值的大小，自动选择在预估频点处傅里叶变换共轭相乘的实部和虚部，从而进一步利用双阵元正弦信号的时延模型估计出双阵元正弦信号的来波方向，有效的克服了常规互谱方法的不足，运算量小，具有很强的工程实用性。

技术实现要素：

技术问题：本发明提供一种利用双阵元正弦信号互功率谱模值预估计信号频率的，在仅增加两次单点傅里叶变换运算的基础上，来波方向估计性能较常规互谱测向法有明显改进的基于频率预估计的双阵元正弦信号来波方向估计方法。

技术方案：本发明的基于频率预估计的双阵元正弦信号来波方向估计方法，包括以下步骤：

第一步：获取待处理的双阵元正弦信号采样数据序列x1(n1),n1＝0,1,…,N-1和x2(n2),n2＝0,1,…,N-1，其中x1(n1)为双阵元基元1的正弦信号采样数据序列，x2(n2)为双阵元基元2的正弦信号采样数据序列：从所述双阵元传感器接收N个采样点的实时采集数据作为待处理的数据序列x1(n1),n1＝0,1,…,N-1和x2(n2),n2＝0,1,…,N-1，或从存储器中提取从检测到信号时刻起始的N个采样点的数据作为待处理的数据序列x1(n1),n1＝0,1,…,N-1和x2(n2),n2＝0,1,…,N-1，所述的N为检测到的正弦信号脉宽长度所对应的采样点个数，取值为2的整数次幂；

第二步：分别计算所述采样数据序列x1(n1)和x2(n2)的离散傅里叶变换X1(l1)和X2(l2)以及两者的互功率谱模值P(k)，计算过程如下：

首先，分别计算x1(n1)和x2(n2)的离散傅里叶变换X1(l1)和X2(l2)，即

其中l1和l2分别为X1(l1)和X2(l2)的离散频率索引，j表示虚数单位，即

然后，根据下式计算所述离散傅里叶变换X1(l1)和X2(l2)的互功率谱模值：

其中k为P(k)的离散频率索引，||代表取模值运算，代表开根号运算；

第三步：计算两阵元正弦信号在离散频点处共轭相乘的实部ξd与虚部ηd、正弦信号频率相对偏差δ与预估频率计算过程如下：

首先，根据下式搜索互功率谱模值P(k)最大值所对应的索引k0：

其中表示在1≤k≤N/2-1范围内搜索P(k)的最大值所对应的离散频率索引；

然后，计算两阵元信号在离散频点k0处共轭相乘的实部ξd与虚部ηd，即：

ξd＝Re[X1(k0)]Re[X2(k0)]+Im[X1(k0)]Im[X2(k0)] 式(5)

ηd＝Re[X1(k0)]Im[X2(k0)]-Im[X1(k0)]Re[X2(k0)] 式(6)

其中Re[]代表取实部运算，Im[]代表取虚部运算；

最后，估计正弦信号频率相对偏差δ和预估频率计算过程如下：

取P(k)在离散频率索引k0-1，k0和k0+1处的值，分别记为Al，Am和Ar，即

Al＝P(k0-1),Am＝P(k0),Ar＝P(k0+1) 式(7)

利用所述值Al，Am和Ar计算正弦信号频率相对偏差δ和预估频率即

其中△f为长度为N的离散傅里叶变换的频率分辨率，Δf＝fs/N，fs为采样频率；δ0为设定的相对频率偏差阈值，δ0为任一大于0且小于0.5的数；

第四步：比较|δ|与设定相对偏差阈值δ0的大小，如果|δ|<δ0，则令两阵元信号在预估频点处傅里叶变换共轭相乘的实部ξ＝ξd，虚部η＝ηd，并转到第六步；否则转到第五步；

第五步：计算两阵元信号在预估频点处傅里叶变换共轭相乘的实部ξ和虚部η，过程如下：

首先，分别计算所述采样数据序列x1(n)和x2(n)在预估频点处的单点傅里叶变换Z1和Z2，即：

然后，计算两阵元信号在预估频点处傅里叶变换共轭相乘的实部ξ和虚部η，即：

ξ＝Re[Z1]Re[Z2]+Im[Z1]Im[Z2] 式(12)

η＝Re[Z1]Im[Z2]-Im[Z1]Re[Z2] 式(13)

第六步：估计双阵元正弦信号来波方向过程如下：

首先，计算双阵元正弦信号的相位差α，即：

其中atan()代表反正切运算，π代表圆周率；

然后，估计双阵元正弦信号的时延差τ，并依据所述双阵元信号的最大时延差对τ进行修正，得到修正时延差τR，即：

其中，d为两阵元的阵间距，c为正弦信号在介质中的传播速度；

最后，估计双阵元正弦信号来波方向即：

其中，acos()为反余弦运算。

进一步的，本发明方法的第二步中，x1(n1)和x2(n2)的离散傅里叶变换X1(l1)和X2(l2)是采用快速傅里叶变换。

进一步的，本发明方法的第三步中，通过搜索互功率谱模值P(k)最大值所对应的离散频率索引k0，然后计算两阵元正弦信号在离散频点k0处共轭相乘的实部ξd与虚部ηd和频率相对偏差δ与预估频率

进一步的，本发明方法的优选方案中，第三步和第四步中，当相对偏差阈值δ0＝0.1时，估计效果较好。

本发明方法首先利用双阵元正弦信号互功率谱模值预估信号频率相对偏差和信号频率，然后依据相对偏差预估计值与门限值的大小关系，自动选择在预估频点处傅里叶变换共轭相乘的实部和虚部，从而进一步利用双阵元正弦信号的时延模型估计出双阵元正弦信号的来波方向。

有益效果：本发明与现有技术相比，具有以下优点：

1.本发明的估计方法利用双阵元正弦信号互功率谱模值对正弦信号的频率进行预估计，与常规方法仅利用双阵元正弦信号中的一路来预估信号频率相比，本发明方法充分利用了双阵元采样数据正弦信号的相关性和噪声的非相关性，同一信噪比条件下，本发明方法估计精度更高；不同信噪比条件下，本发明方法可达到的信噪比下限更低。

2.本发明的估计方法通过估计频率相对偏差δ的大小，在不同的频段采用不同的估计方法：当δ较小时，本发明方法直接使用常规互谱测向方法；当δ较大时，利用两阵元单点傅里叶变换的互谱进行测向。与常规互谱测向方法相比，常规方法不分估计频率相对偏差δ的大小，均直接使用两阵元正弦信号的离散频谱进行互谱测向，当δ较大时，信号频谱泄漏严重，信号的离散频谱与信号的实际频谱相差较大，来波方向估计误差较大；而本发明方法利用单点傅里叶变换对正弦信号离散傅里叶变换的频谱泄露进行补偿，充分利用了正弦信号的频谱特征，提高了双阵元正弦信号的来波方向估计精度。

附图说明

图1所示为本发明的流程图。

图2所示为双阵元接收正弦信号示意图。

图3所示为实施例1中仿真双阵元正弦信号的互功率谱模值。

图4所示为实施例2中仿真双阵元正弦信号的互功率谱模值。

具体实施方式

下面结合实施例和说明书附图对本发明作进一步的说明。

如图1所示，一种基于频率预估计的双阵元正弦信号来波方向估计方法，包括以下步骤：

第一步：获取待处理的双阵元正弦信号采样数据序列x1(n1),n1＝0,1,…,N-1和x2(n2),n2＝0,1,…,N-1，其中x1(n1)为双阵元基元1的正弦信号采样数据序列，x2(n2)为双阵元基元2的正弦信号采样数据序列：从所述双阵元传感器接收N个采样点的实时采集数据作为待处理的数据序列x1(n1),n1＝0,1,…,N-1和x2(n2),n2＝0,1,…,N-1，或从存储器中提取从检测到信号时刻起始的N个采样点的数据作为待处理的数据序列x1(n1),n1＝0,1,…,N-1和x2(n2),n2＝0,1,…,N-1，所述的N为检测到的正弦信号脉宽长度所对应的采样点个数，取值为2的整数次幂；

第二步：分别计算所述采样数据序列x1(n1)和x2(n2)的离散傅里叶变换X1(l1)和X2(l2)以及两者的互功率谱模值P(k)，计算过程如下：

首先，分别计算x1(n1)和x2(n2)的离散傅里叶变换X1(l1)和X2(l2)，即

其中l1和l2分别为X1(l1)和X2(l2)的离散频率索引，j表示虚数单位，即

然后，根据下式计算所述离散傅里叶变换X1(l1)和X2(l2)的互功率谱模值：

其中k为P(k)的离散频率索引，||代表取模值运算，代表开根号运算；

在第二步中，列x1(n1)和x2(n2)的的离散傅里叶变换即式(1)和式(2)，是通过快速傅里叶变换实现的，利用快速傅里叶变换可以降低算法的运算量，提高算法的计算效率；式(3)中k＝0,1,2…,N/2是因为实数据序列的离散傅里叶变换关于中心对称，因此k可以只取前N/2+1个点。

第三步：计算两阵元正弦信号在离散频点处共轭相乘的实部ξd与虚部ηd、正弦信号频率相对偏差δ与预估频率计算过程如下：

首先，根据下式搜索互功率谱模值P(k)最大值所对应的索引k0：

其中表示在1≤k≤N/2-1范围内搜索P(k)的最大值所对应的离散频率索引；

然后，计算两阵元信号在离散频点k0处共轭相乘的实部ξd与虚部ηd，即：

ξd＝Re[X1(k0)]Re[X2(k0)]+Im[X1(k0)]Im[X2(k0)] 式(5)

ηd＝Re[X1(k0)]Im[X2(k0)]-Im[X1(k0)]Re[X2(k0)] 式(6)

其中Re[]代表取实部运算，Im[]代表取虚部运算；

最后，估计正弦信号频率相对偏差δ和预估频率计算过程如下：

取P(k)在离散频率索引k0-1，k0和k0+1处的值，分别记为Al，Am和Ar，即

Al＝P(k0-1),Am＝P(k0),Ar＝P(k0+1) 式(7)

利用所述值Al，Am和Ar计算正弦信号频率相对偏差δ和预估频率即

其中△f为长度为N的离散傅里叶变换的频率分辨率，Δf＝fs/N，fs为采样频率；δ0为设定相对频率偏差阈值，δ0为任一大于0且小于0.5的数；

在第三步中，计算频率相对偏差δ，例如《电子学报》2004年，32(4)的第625页至第628页中公开的内容，本文发明的方法将其作为双阵元正弦信号频率相对偏差的预估值；

在第三步中，分两步实现的：第一步，搜索互功率谱模值P(k)最大值所对应的离散频率索引k0；第二步，计算两阵元正弦信号在离散频点k0处共轭相乘的实部ξd与虚部ηd和频率相对偏差δ与预估频率

第四步：比较|δ|与设定相对偏差阈值δ0的大小，如果|δ|<δ0，则令两阵元信号在预估频点处傅里叶变换共轭相乘的实部ξ＝ξd，虚部η＝ηd，并转到第六步；否则转到第五步；

在第四步中，频率相对偏差阈值δ0可以取大于0且小于0.5之间的任一值，但是通过蒙特卡洛仿真实验得到的优选值为δ0＝0.1；

第五步：计算两阵元信号在预估频点处傅里叶变换共轭相乘的实部ξ和虚部η，过程如下：

首先，分别计算所述采样数据序列x1(n)和x2(n)在预估频点处的单点傅里叶变换Z1和Z2，即：

然后，计算两阵元信号在预估频点处傅里叶变换共轭相乘的实部ξ和虚部η，即：

ξ＝Re[Z1]Re[Z2]+Im[Z1]Im[Z2] 式(12)

η＝Re[Z1]Im[Z2]-Im[Z1]Re[Z2] 式(13)

在第五步中，分两步实现的：第一步，计算x1(n)和x2(n)在预估频点处的单点傅里叶变换Z1和Z2；第二步，计算两阵元信号在预估频点处傅里叶变换共轭相乘的实部ξ和虚部η。

第六步：估计双阵元正弦信号来波方向过程如下：

首先，计算双阵元正弦信号的相位差α，即：

其中atan()代表反正切运算，π代表圆周率；

然后，估计双阵元正弦信号的时延差τ，并依据所述双阵元信号的最大时延差对τ进行修正，得到修正时延差τR，即：

其中，d为两阵元的阵间距，c为正弦信号在介质中的传播速度；

最后，估计双阵元正弦信号来波方向即：

其中，acos()为反余弦运算。

在第六步中，为防止双阵元正弦信号存在相位模糊从而导致估计方向存在多值性的问题，要求预估频率若预估频率本发明方法不再适用。

在第六步中，式(16)和式(17)中的c为正弦信号在介质中的传播速度，如正弦信号表现为声信信号，则在水介质中的传播速度为c＝1500m/s，在空气中的传播速度为c＝340m/s；

在第六步中，分三步实现的：第一步，计算双阵元正弦信号的相位差α；第二步，估计双阵元正弦信号的时延差τ，并依据双阵元信号的最大时延差对τ进行修正得修正时延差τR；第三步，估计双阵元正弦信号来波方向

本发明的一种基于频率预估计的双阵元正弦信号来波方向估计方法：首先利用双阵元正弦信号互功率谱预估信号频率相对偏差和信号频率，然后依据相对偏差预估计值与门限值的大小，自动选择在预估频点处傅里叶变换共轭相乘的实部和虚部，从而进一步利用双阵元正弦信号的时延模型估计出双阵元正弦信号的来波方向。

本发明的实施例中，双阵元接收正弦信号示意图如图2所示，仿真双阵元接收正弦信号模型为：

基元1接收正弦信号模型为

基元2接收正弦信号模型为

其中A为双阵元接收正弦信号的幅度，为双阵元接收正弦信号的初始相位，N为正弦信号采样点数，f0为正弦信号频率，fs为采样频率，d为两阵元的间距，c为正弦信号在介质中的传播速度，θ0为信号来波方向，即待估计值，w1(n1)和w2(n2)分别为基元1和基元2接收的相互独立的高斯白噪声，两者的均值为0，方差为σ²且，方差σ²的大小由信噪比SNR决定：SNR＝10log10[A²/(2σ²)]。

实施例1

仿真信号参数分别设置为：信号幅度A＝1，初始相位采样点数N＝1024，采样频率fs＝4000Hz，则离散傅里叶变换的频率分辨率Δf＝fs/N＝3.9063Hz，信号频率f0＝291Hz，阵元间距d＝0.5m，正弦信号在介质中的传播速度c＝1500m/s(声音在水中的传播速度)，信号来波方向θ＝π/4rad，信噪比SNR＝0dB，设定相对偏差阈值δ0＝0.1。

图3所示是频率为291Hz的仿真双阵元正弦信号的互功率谱模值，从图3可以看出，互功率谱模值P(k)最大值所对应的离散频率索引为75，次大值索引为74。

搜索互功率谱模值P(k)最大值所对应的索引得k0＝75，取P(k)在离散频率索引k0-1，k0和k0+1处的模值Al，Am和Ar，分别为

Al＝|P(k0-1)|＝0.2829，Am＝|P(k0)|＝0.369，Ar＝|P(k0+1)|＝0.1303

由此可知Al>Ar，因此依据式(8)可以计算出频率相对偏差

由此可知|δ|>δ0，进一步依据式(9)计算出预估频率

然后，比较|δ|与设定相对偏差阈值δ0的大小，δ0的设定阈值为0.1，由此可得|δ|>δ0，依据第四步的判决条件转入第五步。利用式(10)和式(11)求得：

Z1＝0.539-j0.1456 Z2＝0.4856+j0.088

进一步利用式(12)和式(13)求得

ξ＝0.2489 η＝0.1181

然后，转到第六步，因为ξ＝0.2489>0且η＝0.1181>0，依据式(14)可知双阵元正弦信号的相位差

然后，利用式(15)和式(16)得修正时延差

τR＝2.4208*10^-4

最后，依据第六步的式(17)双阵元正弦信号来波方向

估计来波方向相对误差为

实施例2

仿真信号参数分别设置为：信号幅度A＝1，初始相位采样点数N＝1024，采样频率fs＝4000Hz，则离散傅里叶变换的频率分辨率Δf＝fs/N＝3.9063Hz，信号频率f0＝390Hz，阵元间距d＝0.1m，正弦信号在介质中的传播速度c＝340m/s(声音在空气中的传播速度)，信号来波方向θ＝π/3rad，信噪比SNR＝3dB，设定相对偏差阈值δ0＝0.1。

图4所示是频率为390Hz的仿真双阵元正弦信号的互功率谱模值，从图4可以看出，互功率谱模值P(k)最大值所对应的离散频率索引为100，次大值索引为99。

搜索互功率谱模值P(k)最大值所对应的索引得k0＝100，取P(k)在离散频率索引k0-1，k0和k0+1处的模值Al，Am和Ar，分别为

Al＝|P(k0-1)|＝0.0973，Am＝|P(k0)|＝0.4375，Ar＝|P(k0+1)|＝0.0526

由此可知Al>Ar，因此依据式(8)可以计算出频率相对偏差

由此可知|δ|>δ0，进一步依据式(9)计算出预估频率

然后，比较|δ|与设定相对偏差阈值δ0的大小，δ0的设定阈值为0.1，由此可得δ|>δ0，依据第四步的判决条件转入第五步。利用式(10)和式(11)求得：

Z1＝0.4409+j0.0205 Z2＝-0.1715+j0.4505

进一步利用式(12)和式(13)求得

ξ＝-0.0664 η＝0.2021

然后，转到第六步，因为ξ＝-0.0664<0且η＝0.2021>0，依据式(14)可知双阵元正弦信号的相位差

然后，利用式(15)和式(16)得修正时延差

τR＝7.7068*10^-4

最后，依据第六步的式(17)双阵元正弦信号来波方向

估计来波方向相对误差为

从实施例1和例2的结果可以看出，本发明估计方法可以获得良好的估计精度，而且计算简单，计算量小，适用于高精度快速估计双阵元正弦信号来波方向的场合。

上述实施例仅是本发明的优选实施方式，应当指出：对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和等同替换，这些对本发明权利要求进行改进和等同替换后的技术方案，均落入本发明的保护范围。