CN104218973A - 基于Myriad滤波的跳频信号参数估计方法 - Google Patents

Abstract

本发明公开了一种基于Myriad滤波的跳频信号参数估计方法。本发明首先利用Myriad滤波器对脉冲噪声进行抑制，然后结合线性时频分析方法STFT，完成α稳定分布噪声下跳频信号的跳频周期、跳变时刻和跳频频率参数的估计。具体步骤包括：1、采集信号，2、加权麦瑞德Myriad滤波，3、时频分析，4、搜索最大值，5、提取参数信息。本发明克服了已有技术无法抑制强脉冲噪声的缺陷，提高了低广义信噪比下跳频信号参数估计的精度。

Description

基于Myriad滤波的跳频信号参数估计方法

技术领域

本发明属于通信技术领域，更进一步涉及无线通信技术、信号处理技术领域中基于麦瑞德Myriad滤波的跳频信号参数估计方法。本发明利用麦瑞德Myriad滤波器对脉冲噪声进行抑制，用短时傅里叶变换STFT估计跳频信号的跳频周期、跳变时刻和跳频频率参数，实现非高斯脉冲噪声环境中跳频信号的参数提取。

背景技术

跳频是一种重要的扩频通信方式，具有抗干扰能力强、频谱利用率高，优秀的组网能力等特点，广泛应用于非协作通信领域中。因此，在非协作通信系统的接收端，对跳频信号的跳频周期、跳变时刻和跳频频率参数进行精确估计，具有非常重要的意义。目前，跳频信号的参数估计方法主要有短时傅里叶变换(STFT)、维格纳-威尔分布(WVD)及其改进方法等线性时频分析方法。

以短时傅里叶变换(STFT)为代表的线性时频分析方法，将背景噪声建立为高斯模型，直观地给出了跳频信号的参数信息，可以在不知道任何先验知识的情况下估计跳频信号的参数。但是大量的研究发现，实际环境中的干扰和噪声，例如通信多通道干扰、低频空气噪声、水声和雷达杂波等均服从非高斯分布，并且具有一定的脉冲特性，更适合用α稳定分布来描述。在这类具有脉冲特性的背景噪声下，基于高斯噪声模型的线性时频分析方法对跳频信号的参数估计性能严重下降。因此，研究脉冲噪声下跳频信号的参数估计成为非协作通信中亟待解决的问题。

中国工程物理研究院电子工程研究所拥有的专利技术“跳频通信系统跳频间隔的一种盲识别方法”(申请号201110158600.7，申请日2011.06.14，授权号CN102223155B，授权日2013.10.09)中提出了一种基于误差函数的跳频间隔估计方法。该方法构造了一个误差函数，通过搜索频率集合中的频率值，将使得误差函数最小的频率值作为跳频间隔的估计值，无需利用传统的线性时频分析方法提取参数信息。该专利技术存在的不足是，脉冲噪声对误差函数的影响很大，尤其在当噪声的脉冲性较强时，已经无法通过该误差函数得到跳频间隔的准确估计值。

刘杰等人发表的文章“Parameter Estimation of Frequency Hopping Signals inAlpha Stable Noise Environment”(IEEE 11th International Conference on SignalProcessing,2012,1:250-253)中提出了一种基于分数低阶时频分析的方法FLOSTFT。该方法首先对接收信号做分数低阶运算，然后采用短时傅里叶变换STFT方法提取跳频信号的参数信息。该方法在一定程度上能够抑制脉冲噪声对跳频参数估计的影响，提高脉冲噪声下跳频信号参数估计的精度。但是，该方法存在的不足是，由于分数低阶运算对跳频信号本身也有一定的抑制作用，因而在低广义信噪比下，该方法的参数估计性能会明显下降。

发明内容

本发明的目的在于克服上述已有跳频信号参数估计方法的不足，提供了一种基于麦瑞德Myriad滤波的跳频信号参数估计方法，以抑制脉冲噪声对跳频信号参数估计的影响，提高低广义信噪比下跳频参数的估计性能。

实现本发明目的的具体思路是：首先对采集信号进行加权麦瑞德Myriad滤波处理，然后采用基于短时傅里叶变换STFT的时频分析方法，完成跳频信号的参数估计，提高非高斯脉冲噪声下跳频信号参数估计的性能。

实现本发明目的的具体步骤如下：

(1)采集信号：

信号采集系统通过接收天线，选取任意一段含有原始跳频信号和脉冲噪声的接收信号，将所选取的接收信号作为采集信号；

(2)加权麦瑞德Myriad滤波：

(2a)采用加权麦瑞德Myriad滤波器，对采集信号进行滤波处理，得到加权麦瑞德Myriad滤波器的初始输出值；

(2b)采用基于平均误差最小准则的自适应方法，通过迭代，获得加权麦瑞德Myriad滤波器的最优权值；

(2c)采用最优权值的麦瑞德Myriad滤波器，对采集信号进行滤波处理，得到加权麦瑞德Myriad滤波器的最终输出值；

(3)时频分析：

采用基于短时傅里叶变换STFT的时频分析方法，对加权麦瑞德Myriad滤波器的最终输出值做短时傅里叶变换，得到滤波处理后采集信号的时频分布二维矩阵；

(4)搜索最大值：

从滤波处理后采集信号的时频分布二维矩阵中，搜索每一时刻沿频率轴的最大值，将搜索得到的所有最大值构成一个最值向量

(5)提取参数信息：

(5a)对最值向量做快速傅里叶变换FFT，得到跳频间隔的估计值和跳频周期的估计值；

(5b)按照下式，搜索最值向量的最大值，得到与该最大值所对应的最值坐标值：

$p_i=\arg \underset{n}{\max }\left\{\hat{y}\right\}$

其中，p_i表示最值向量的最大值所对应的最值坐标值；表示取最大值所对应的下标值操作；表示最值向量；n表示最值向量的下标值，n∈[H+1,N-H]；H表示跳频间隔的估计值；N表示加权麦瑞德Myriad滤波器的窗口长度；

(5c)按照下式，计算跳频信号跳变时刻的估计值：

$n_i=p_i-\frac{1}{2}H$

其中，n_i表示跳频信号跳变时刻的估计值；p_i表示最值向量的最大值所对应的最值坐标值；i表示最值向量的频率，i＝1,2…P；P表示最值向量的长度；H表示跳频间隔的估计值；

(5d)按照下式，计算跳频信号跳频频率的估计值：

${\hat{f}}_i=\arg \underset{l}{\max }\left\{\underset{n_i}{\overset{n_{i+1}}{\mathrm{\Σ}}}{\mathrm{MYRSTFT}}_x(k,l)\right\}*f_s/\left(2N\right)$

其中，表示跳频信号跳频频率的估计值；表示取最大值所对应的下标值操作；n_i表示跳频信号跳变时刻的估计值；i表示最值向量的频率，i＝1,2…P；P表示最值向量的长度；MYRSTFT_x(k,l)表示滤波处理后采集信号的时频分布二维矩阵；k表示采集信号的时刻；l表示采集信号的频率；f_s表示采集信号的采样频率；N表示采集信号的长度。

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

第一，本发明采用了麦瑞德Myriad滤波器，其过程可行性高，可有效抑制具有脉冲特性的非高斯噪声，尤其是α稳定分布噪声，克服了现有技术中脉冲噪声会对结果造成很大影响的不足，使得本发明对脉冲性很强的噪声具有很好的抑制效果。

第二，本发明将麦瑞德Myriad滤波方法与线性时频分析方法STFT相结合，由于麦瑞德Myriad滤波方法对跳频信号本身的影响较小，因而克服了现有技术在低广义信噪比下跳频参数估计性能明显下降的缺点，使得本发明在低广义信噪比下仍然保持较高的估计精度。

附图说明

图1是本发明的流程图；

图2是对称稳定SαS分布噪声的时域波形示意图；

图3是本发明最值向量的示意图；

图4是本发明与现有方法所得跳频周期的估计性能比较图。

具体实施方式

下面结合附图对本发明做进一步的描述。

参照附图1，本发明的具体步骤如下。

步骤1，采集信号。

信号采集系统通过接收天线，选取任意一段含有原始跳频信号和脉冲噪声的接收信号，将所选取的接收信号作为采集信号。

本发明采集的接收信号中所含有的脉冲噪声为对称稳定SαS分布噪声，对称稳定SαS分布噪声的时域波形如图2所示。

图2(a)所示的时域波形为特征指数α＝1.5的对称稳定SαS分布噪声，图2(a)中以实线段连起的曲线表示对称稳定SαS分布噪声的时域波形，横坐标表示对称稳定SαS分布噪声的时域采样点，纵坐标表示对称稳定SαS分布噪声的振幅。

图2(b)所示的时域波形为特征指数α＝0.8的对称稳定SαS分布噪声，图2(b)中以实线段连起的曲线表示对称稳定SαS分布噪声的时域波形，横坐标表示对称稳定SαS分布噪声的时域采样点，纵坐标表示对称稳定SαS分布噪声的振幅。

图2(a)和图2(b)中的对称稳定SαS分布噪声在某些时域采样点的振幅很大，图2(b)比图2(a)中具有大振幅的时域采样点的个数更多。因此，可以用对称稳定SαS分布噪声表示具有不同强度的脉冲噪声。

步骤2，加权麦瑞德Myriad滤波。

采用加权麦瑞德Myriad滤波器，对采集信号进行滤波处理，具体的实施步骤如下：

选取加权麦瑞德Myriad滤波器的窗口长度N，N的取值范围为5～10；选取加权麦瑞德Myriad滤波器的权向量w，w＝[w₁,w₂,...,w_N]，初始权向量w(0)＝[0,...,2N,...,0]，其中N表示在5～10范围内所选定的加权麦瑞德Myriad滤波器的窗口长度；选取加权麦瑞德Myriad滤波器的尺度参数K，其中α表示稳定分布噪声的特征指数，α的取值区间为(0,2]。

用柯西分布公式对采集信号构造似然函数，在跳频信号幅度的取值范围内等间隔选取100个样值，分别代入似然函数，将似然函数取得最大值时所对应的样值作为加权麦瑞德Myriad滤波器的输出值。

用加权麦瑞德Myriad滤波器的初始输出值减去原始跳频信号的幅度值，作为加权麦瑞德Myriad滤波器的滤波误差。采用基于平均绝对误差最小准则，对加权麦瑞德Myriad滤波器的初始权向量进行L次迭代，L的取值范围为50～100，权值的迭代公式如下：

$w_i(n+1)=P[w_i\left(n\right)+\mathrm{\μsgn}\left(e\right)\frac{(y-x_i)}{{[1+\frac{w_i}{K^2}{(y-x_i)}^2]}^2}]$

其中，w_i(n)和w_i(n+1)分别表示加权麦瑞德Myriad滤波器第i个加权值的第n次与第n+1次迭代值；μ表示加权麦瑞德Myriad滤波器的步长因子，μ的取值范围为[0.001,0.1]；e表示加权麦瑞德Myriad滤波器的滤波误差；y表示加权麦瑞德Myriad滤波器的输出值；x表示采集信号的幅度；K表示加权麦瑞德Myriad滤波器的尺度参数；sgn(·)表示符号函数；P[·]表示矩阵函数。

采用最优权值的麦瑞德Myriad滤波器，对采集信号进行滤波处理，得到加权麦瑞德Myriad滤波器的最终输出值。

步骤3，时频分析。

采用基于短时傅里叶变换STFT的时频分析方法，对加权麦瑞德Myriad滤波器的最终输出值做短时傅里叶变换，得到滤波处理后采集信号的时频分布二维矩阵。

步骤4，搜索最大值。

从滤波处理后采集信号的时频分布二维矩阵中，搜索每一时刻沿频率轴的最大值，将搜索得到的所有最大值构成一个最值向量可表示如下：

$\hat{y}=\underset{l}{\max }\left\{{\mathrm{MYRSTFT}}_x(k,l)\right\}$

其中，max{·}表示取最大值操作；MYRSTFT_x(k,l)表示滤波处理后采集信号的时频分布二维矩阵；k表示采集信号的时刻；l表示采集信号的频率。

图3所示为本发明最值向量的示意图，图3中以实线段连起的曲线表示最值向量的幅度，横坐标表示最值向量的采样点，纵坐标表示最值向量的幅度值。图3中的曲线有8个具有周期性的峰值，每个峰值所对应的采样点代表跳频信号的频率，每个峰值之间的距离代表跳频周期，根据峰值的周期性可以得到跳频间隔的估计值。

步骤5，提取参数信息。

采用加权麦瑞德Myriad滤波器，对采集信号进行滤波处理，具体的实施步骤如下：

用快速傅里叶变换FFT对最值向量做谱分析，得到一组跳频间隔的估计值；对相邻的跳频间隔估计值作差，将所有相邻跳频间隔估计值差值的平均值，作为跳频周期的估计值。

按照下式，搜索最值向量的最大值，得到与该最大值所对应的最值坐标值：

$p_i=\arg \underset{n}{\max }\left\{\hat{y}\right\}$

其中，pi表示最值向量的最大值所对应的最值坐标值；表示取最大值所对应的下标值操作；表示最值向量；n表示最值向量的下标值，n∈[H+1,N-H]；H表示跳频间隔的估计值；N表示加权麦瑞德Myriad滤波器的窗口长度。

按照下式，计算跳频信号跳变时刻的估计值：

$n_i=p_i-\frac{1}{2}H$

其中，n_i表示跳频信号跳变时刻的估计值；p_i表示最值向量的最大值所对应的最值坐标值；i表示最值向量的频率，i＝1,2…P；P表示最值向量的长度；H表示跳频间隔的估计值。

按照下式，计算跳频信号跳频频率的估计值：

${\hat{f}}_i=\arg \underset{l}{\max }\left\{\underset{n_i}{\overset{n_{i+1}}{\mathrm{\Σ}}}{\mathrm{MYRSTFT}}_x(k,l)\right\}*f_s/\left(2N\right)$

其中，表示跳频信号跳频频率的估计值；表示取最大值所对应的下标值操作；n_i表示跳频信号跳变时刻的估计值；i表示最值向量的频率，i＝1,2…P；P表示最值向量的长度；MYRSTFT_x(k,l)表示滤波处理后采集信号的时频分布二维矩阵；k表示采集信号的时刻；l表示采集信号的频率；f_s表示采集信号的采样频率；N表示采集信号的长度。

下面结合仿真图对本发明做进一步的描述。

1.仿真条件：

本发明仿真实验的运行系统为Intel(R)Core(TM)i5 CPU 6503.20GHz，32位Windows操作系统，仿真软件采用MATLAB R(2010b)，仿真参数设置如下所示。

选取跳频信号的跳频频率f_k＝{1.1,1.3,1.6,1.0,1.2,1.5,1.7,1.4}(kHz)，跳频周期T＝0.05s；采样频率f_s＝4kHz，采样点数M＝1600。选取加权麦瑞德Myriad滤波器的参数为：窗口长度N＝5，初始权向量w(0)＝[0,0,10,0,0]，迭代次数L＝75，步长因子μ＝0.01。

2.仿真内容与结果分析：

仿真1：

在对称稳定SαS分布噪声下采用本发明与现有方法估计跳频周期，得到的的估计性能比较图如图4所示。

图4(a)为特征指数α＝1.5的对称稳定SαS分布噪声下采用不同方法估计跳频周期的性能比较图。图4(a)中的横坐标表示广义信噪比，纵坐标表示标准差。图4(a)中以三角形连起的曲线表示本发明所得跳频周期估计性能曲线，以正方形连起的曲线表示基于分数低阶时频分析的方法所得跳频周期估计性能曲线，以圆圈连起的曲线表示基于时频分析的方法所得跳频周期估计性能曲线。

图4(b)为特征指数α＝0.8的对称稳定SαS分布噪声下采用不同方法估计跳频周期的性能比较图。图4(b)中的横坐标表示广义信噪比，纵坐标表示标准差。图4(a)中以三角形连起的曲线表示本发明所得跳频周期估计性能曲线，以正方形连起的曲线表示基于分数低阶时频分析的方法所得跳频周期估计性能曲线，以圆圈连起的曲线表示基于时频分析的方法所得跳频周期估计性能曲线。

由图4(a)可见，当广义信噪比GSNR≥0dB时，本发明可以准确估计出跳频周期。当广义信噪比GSNR≥2dB时，基于分数低阶时频分析的方法可以准确估计出跳频周期，基于时频分析的方法估计性能严重下降。

由图4(b)可见，当广义信噪比GSNR≥6dB时，本发明可以准确估计出跳频周期。当广义信噪比GSNR≥9dB时，基于分数低阶时频分析的方法可以准确估计出跳频周期，基于时频分析的方法估计性能严重下降。

比较图4(a)和图4(b)可见，本发明在低广义信噪比下的跳频周期估计性能优于基于分数低阶时频分析的方法和基于时频分析的方法。

仿真2：

在特征指数α＝1.5的对称稳定SαS分布噪声下，分别采用基于分数低阶时频分析的方法和本发明估计跳变时刻，得到如表1所示的3dB对称稳定SαS分布噪声下跳变时刻估计值的相对误差对比表。

表1 跳变时刻估计值的相对误差表

由表1可见，采用本发明估计跳变时刻的最大相对误差仅为2.8％，采用基于分数低阶时频分析的方法估计跳变时刻的最大相对误差为3.0％。可见，本发明比基于分数低阶时频分析的方法的跳变时刻估计性能更好。

在特征指数α＝1.5的对称稳定SαS分布噪声下，分别采用基于分数低阶时频分析的方法和本发明估计跳变时刻，得到如表2所示的3dB对称稳定SαS分布噪声下跳频频率估计值的标准误差对比。

表2 跳频频率估计值的标准误差

由表2可见，本发明的跳频频率估计值比基于分数低阶时频分析的方法标准误差小，说明本发明的跳频频率估计精度比基于分数低阶时频分析的方法更高。

综上所述，由两个仿真实验所获得的四个结果表明，采用本发明能够有效抑制脉冲噪声对跳频参数估计的影响，提高低广义信噪比下跳频参数的估计性能。

Claims (5)

1.一种基于Myriad滤波的跳频信号参数估计方法，包括如下步骤：

(1)采集信号：

信号采集系统通过接收天线，选取任意一段含有原始跳频信号和脉冲噪声的接收信号，将所选取的接收信号作为采集信号；

(2)加权麦瑞德Myriad滤波：

(2a)采用加权麦瑞德Myriad滤波器，对采集信号进行滤波处理，得到加权麦瑞德Myriad滤波器的初始输出值；

(2b)采用基于平均误差最小准则的自适应方法，通过迭代，获得加权麦瑞德Myriad滤波器的最优权值；

(2c)采用最优权值的麦瑞德Myriad滤波器，对采集信号进行滤波处理，得到加权麦瑞德Myriad滤波器的最终输出值；

(3)时频分析：

采用基于短时傅里叶变换STFT的时频分析方法，对加权麦瑞德Myriad滤波器的最终输出值做短时傅里叶变换，得到滤波处理后采集信号的时频分布二维矩阵；

(4)搜索最大值：

从滤波处理后采集信号的时频分布二维矩阵中，搜索每一时刻沿频率轴的最大值，将搜索得到的所有最大值构成一个最值向量

(5)提取参数信息：

(5a)对最值向量做快速傅里叶变换FFT，得到跳频间隔的估计值和跳频周期的估计值；

(5b)按照下式，搜索最值向量的最大值，得到与该最大值所对应的最值坐标值：

$p_i=\arg \underset{n}{\max }\left\{\hat{y}\right\}$

其中，p_i表示最值向量的最大值所对应的最值坐标值；表示取最大值所对应的下标值操作；表示最值向量；n表示最值向量的下标值，n∈[H+1,N-H]；H表示跳频间隔的估计值；N表示加权麦瑞德Myriad滤波器的窗口长度；

(5c)按照下式，计算跳频信号跳变时刻的估计值：

$n_i=p_i-\frac{1}{2}H$

其中，n_i表示跳频信号跳变时刻的估计值；p_i表示最值向量的最大值所对应的最值坐标值；i表示最值向量的频率，i＝1,2…P；P表示最值向量的长度；H表示跳频间隔的估计值；

(5d)按照下式，计算跳频信号跳频频率的估计值：

${\hat{f}}_i=\arg \underset{l}{\max }\left\{\underset{n_i}{\overset{n_{i+1}}{\mathrm{\Σ}}}{\mathrm{MYRSTFT}}_x(k,l)\right\}*f_s/\left(2N\right)$

其中，表示跳频信号跳频频率的估计值；表示取最大值所对应的下标值操作；n_i表示跳频信号跳变时刻的估计值；i表示最值向量的频率，i＝1,2…P；P表示最值向量的长度；MYRSTFT_x(k,l)表示滤波处理后采集信号的时频分布二维矩阵；k表示采集信号的时刻；l表示采集信号的频率；f_s表示采集信号的采样频率；N表示采集信号的长度。

2.根据权利要求1所述的基于麦瑞德Myriad滤波的跳频信号参数估计方法，其特征在于：步骤(2a)、步骤(2c)所述采集信号滤波的具体步骤如下：

第一步，选取加权麦瑞德Myriad滤波器的窗口长度N，N的取值范围为5～10；选取加权麦瑞德Myriad滤波器的权向量w，w＝[w₁,w₂,...,w_N]，其中N表示在5～10范围内所选定的加权麦瑞德Myriad滤波器的窗口长度；选取加权麦瑞德Myriad滤波器的尺度参数K，其中α表示稳定分布噪声的特征指数，α的取值范围为(0,2]；

第二步，用柯西分布公式对采集信号构造似然函数，在跳频信号幅度的取值范围内等间隔选取100个样值，分别代入似然函数，将似然函数取得最大值时所对应的样值作为加权麦瑞德Myriad滤波器的输出值。

3.根据权利要求1所述的基于麦瑞德Myriad滤波的跳频信号参数估计方法，其特征在于：步骤(2b)所述基于平均误差最小准则的自适应方法的具体步骤如下：

第一步，用加权麦瑞德Myriad滤波器的初始输出值减去原始跳频信号的幅度值，作为加权麦瑞德Myriad滤波器的滤波误差；

第二步，按照下式，选取加权麦瑞德Myriad滤波器的初始权向量：

w(0)＝[0,...,2N,...,0]

其中，w(0)表示加权麦瑞德Myriad滤波器的初始权向量；N表示加权麦瑞德Myriad滤波器的窗口长度；

第三步，采用基于平均绝对误差最小准则，对加权麦瑞德Myriad滤波器的初始权向量进行L次迭代，L的取值范围为50～100，权值的迭代公式如下：

$w_i(n+1)=P[w_i\left(n\right)+\mathrm{\μsgn}\left(e\right)\frac{(y-x_i)}{{[1+\frac{w_i}{K^2}{(y-x_i)}^2]}^2}]$

其中，w_i(n)和w_i(n+1)分别表示加权麦瑞德Myriad滤波器第i个加权值的第n次与第n+1次迭代值；μ表示加权麦瑞德Myriad滤波器的步长因子，μ的取值范围为[0.001,0.1]；e表示加权麦瑞德Myriad滤波器的滤波误差；y表示加权麦瑞德Myriad滤波器的输出值；x表示采集信号的幅度；K表示加权麦瑞德Myriad滤波器的尺度参数；sgn(·)表示符号函数；P[·]表示矩阵函数。

4.根据权利要求1所述的基于麦瑞德Myriad滤波的跳频信号参数估计方法，其特征在于：步骤(4)所述搜索最大值按照下式进行：

$\hat{y}=\underset{l}{\max }\left\{{\mathrm{MYRSTFT}}_x(k,l)\right\}$

其中，表示由滤波处理后采集信号的时频分布二维矩阵得到的最值向量；MYRSTFT_x(k,l)表示滤波处理后采集信号的时频分布二维矩阵；max{·}表示取最大值；x表示采集信号；k表示时刻；l表示频率。

5.根据权利要求1所述的基于麦瑞德Myriad滤波的跳频信号参数估计方法，其特征在于：步骤(5a)所述快速傅里叶变换FFT的具体步骤如下：

第一步，用快速傅里叶变换FFT对最值向量做谱分析，得到一组跳频间隔的估计值；

第二步，对相邻的跳频间隔估计值作差，将所有相邻跳频间隔估计值差值的平均值，作为跳频周期的估计值。