CN104749432B - 基于聚焦s变换的多分量非平稳信号瞬时频率估计方法 - Google Patents

Abstract

本发明公开了一种基于聚焦S变换的多分量非平稳信号瞬时频率估计方法，主要解决现有技术的能量聚集性不高、时频分辨率无法灵活调整且无法实现高精度的多分量信号瞬时频率估计问题，其实现步骤：1.输入待分析的非平稳信号；2.根据聚集性准则优化聚焦S变换中窗函数的控制参数；3.计算优化后的聚焦S变换并获得信号时频分布能量图；4.对非平稳信号的时频分布能量图进行二值化处理；5.根据二值化处理后的二值化时频图得到信号的瞬时频率估计。本发明通过改进S变换的窗函数，实现了窗函数控制参数的自适应调节，提高了时频分布能量聚集性和瞬时频率估计精准性，可用于雷达截获、通信对抗、语音识别、医学脑电波信号分析。

Description

基于聚焦S变换的多分量非平稳信号瞬时频率估计方法

技术领域

本发明属于信号与信息处理技术领域，特别涉及一种瞬时频率估计方法，可用于雷达截获、通信对抗、语音识别、医学脑电波信号分析。

背景技术

实际工程应用领域普遍存在非平稳信号，如雷达信号、通信信号、语音信号等，其共性是拥有时变的频率特征。传统的傅里叶变换FT针对周期性平稳信号，不能反映非平稳信号的局部特征。因此，在分析这类信号结构特点时，需要采用瞬时频率IF估计方法。当截获到的非平稳信号进入信号分析系统中，瞬时频率IF估计可以表征这类信号频率随时间的变化情况，对于分析非平稳信号具有瞬时有效性。

时频分析是一种信号处理方法，利用信号的时间—频率联合分布来表征信号能量的分布特征。将时频分析方法应用于瞬时频率估计中可以更精确地表征信号特征，并获得更高的能量聚集性和更强的抑噪能力。许多基于时频域的分析方法在文献中已被提出，其中常用的时频分析方法有短时傅里叶变换STFT，连续小波变换CWT，Wigner-Ville分布WVD，S变换ST，广义S变换GST等。短时傅里叶变换是一种线性时频分析方法且易于实现，但其能量聚集性很低。连续小波有着多分辨率特性但是其尺度参量缺乏明确物理意义。Wigner-Ville分布WVD是一种双线性变换，针对单分量信号具有很强的能量聚集性但是在多分量信号分析中会受到交叉项的影响，使得信号难以分辨。S变换是一种将小波与短时傅里叶变换结合的分析方法，它保留了信号的相位信息，但是由于基本小波形态固定，限制了S变换在工程应用中的灵活性。随着S变换在工程中的推广应用，各种改进的广义S变换被提出，并广泛应用于雷达信号、脑电图信号EEG，地震信号和电压信号等的分析中。但广义S变换并没有改变时频聚集性的变化趋势，不能灵活的控制特定时频区域的能量聚集性。

综上所述，现有瞬时频率估计技术中普遍存在能量聚集性不高、时频分辨率无法灵活调整且无法实现高精度的多分量信号瞬时频率估计。

发明内容

本发明的目的在于提供一种基于聚焦S变换的多分量非平稳信号的瞬时频率估计方法，以解决现有技术对该类信号的瞬时频率估计中能量聚集性较低、时频分辨率调整性较差且对多分量信号估计精度较低的缺陷。

本发明是基于聚焦S变换的时频分析方法，其技术方案通过改进广义S变换中的窗函数，加入中心频率控制参数shiftN并进行平移变换，同时加入窗长控制参数r和窗长变化趋势控制参数p并进行尺度变换，从而控制不同时频区域的能量聚集性，实现对目标时频区域的自适应时频分辨率调整。即使在低信噪比环境下，也能获得高聚集性和高分辨率的时频分布能量图。其实现步骤包括如下：

(1)输入非平稳信号的时域形式h(t)，其中t代表时间，非平稳信号包括线性调频信号、正弦调频信号、线性调频与二次调频的多分量混合信号；

(2)对非平稳信号的时域形式h(t)进行傅里叶变换，得到非平稳信号的频域形式H(f)，其中f代表频率；

(3)定义聚焦S变换窗函数的频域形式为：

W(α,f)＝k·exp((f+shiftN)^2pα^2r/2) <1>

其中，k代表归一化参数，shiftN代表中心频率控制参数，r代表窗长控制参数，p代表窗长变化趋势控制参数，α代表频移；

(4)根据式<1>求出聚焦S变换的时频表达式：

其中，H(α+f)代表加入频移α后非平稳信号的频域形式，W(α,f)代表聚焦S变换窗函数的频域形式，i代表虚数单位；

(5)根据式<2>求出聚集性准则判别式CM如下：

(6)依照式<3>优化中心频率控制参数shiftN，窗长控制参数r和窗长变化趋势控制参数p，即在聚集性准则判别式CM取最大值时，得到最优化中心频率控制参数shiftN'，最优化窗长控制参数r'和最优化窗长变化趋势控制参数p'；

(7)将最优化中心频率控制参数shiftN'，最优化窗长控制参数r'和最优化窗长变化趋势控制参数p'代入式<1>，得到聚焦S变换窗函数的频域聚焦形式W'(α,f)；

(8)根据聚焦S变换窗函数的频域聚焦形式W'(α,f)，按照式<2>得到非平稳信号最优化聚焦S变换的时频表达式S'(t,f)，根据S'(t,f)获得非平稳信号h(t)的时频分布能量图；

(9)根据时频域支撑区间置信度准则设定门限γ，以门限γ为二值化处理门限对非平稳信号h(t)的时频分布能量图进行二值化处理，获得表征非平稳信号h(t)时频分布的时频曲线，将该时频曲线作为非平稳信号h(t)的瞬时频率估计结果。

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

1.本发明由于采用基于聚焦S变换的时频分析方法，可通过自适应调节聚焦S变换窗函数的时频特性，使其与待分析信号的时频特性相匹配，从而有效提高信号时频分析的能量聚集性。

2.本发明由于采用三种不同的窗函数优化控制参数，可灵活调整目标信号所在时频域的时频分辨率。

3.本发明由于采用聚焦S变换，其线性变换特性可有效抑制交叉项的干扰，能够实现高精度的多分量非平稳信号的瞬时频率估计。

附图说明

图1为本发明的实现流程图；

图2为对线性调频信号进行聚焦S变换后的时频分布能量图；

图3为对线性调频信号进行聚焦S变换后的瞬时频率估计图；

图4为对正弦调频信号进行聚焦S变换后的时频分布能量图；

图5为对正弦调频信号进行聚焦S变换后的瞬时频率估计图；

图6为对线性调频和二次调频的多分量混合信号进行聚焦S变换后的时频分布能量图；

图7为对线性调频和二次调频的多分量混合信号进行聚焦S变换后的瞬时频率估计图。

具体实施方式

以下结合附图对本发明的技术方案和效果作进一步详细描述。

参照图1，本发明的实施步骤如下：

步骤1，输入非平稳信号的时域形式。

非平稳信号，包括线性调频信号、正弦调频信号、线性调频与二次调频的多分量混合信号，其时域形式表示为h(t)，其中t代表时间。

步骤2，根据聚集性准则优化聚焦S变换中窗函数的控制参数。

2a)对非平稳信号的时域形式h(t)进行傅里叶变换，得到非平稳信号的频域形式H(f)，其中f代表频率；

2b)根据聚焦S变换窗函数的频域形式设定窗函数初始值，包括：初始中心频率控制参数shiftN＝0，初始窗长控制参数r＝2，初始窗长变化趋势控制参数p＝2；设定聚焦S变换窗函数的频域形式为：

W(α,f)＝k·exp((f+shiftN)^2pα^2r/2)<1>

其中，k代表归一化参数，shiftN代表中心频率控制参数，r代表窗长控制参数，p代表窗长变化趋势控制参数，α代表频移；

2c)根据式<1>求出聚焦S变换的时频表达式：

其中，H(α+f)代表加入频移α后非平稳信号的频域形式，i代表虚数单位；

2d)根据式<2>计算聚集性准则判别式CM如下：

聚集性准则是时频分布能量聚集性的判定标准，其判别式CM越大代表与之对应的时频分布能量聚集性越高，由该时频分布所获得的信号瞬时频率估计越精确；

2e)在聚集性准则判别式CM最大化的约束条件下，在各自的优化区间内，搜索得到最优中心频率控制参数shiftN'，最优窗长控制参数r'和最优窗长变化趋势控制参数p'，其中，参数优化区间设置为：中心频率控制参数shiftN优化区间为[-2fs，2fs]，窗长控制参数r和窗长变化趋势控制参数p的优化区间均为[0.6，1.3]，其中fs为采样频率。

步骤3，计算优化后的聚焦S变换并获得信号时频分布能量图。

3a)将最优化中心频率控制参数shiftN'，最优化窗长控制参数r'和最优化窗长变化趋势控制参数p'代入式<1>，得到聚焦S变换窗函数的频域聚焦形式W'(α,f)；

3b)根据聚焦S变换窗函数的频域聚焦形式W'(α,f)，按照式<2>得到非平稳信号最优化聚焦S变换的二维时频表达式S'(t,f)，根据该二维时频表达式S'(t,f)获得非平稳信号h(t)的时频分布能量图。

步骤4，对非平稳信号的时频分布能量图进行二值化处理。

4a)根据非平稳信号h(t)的时频分布能量图得到二维矩阵S_x；

4b)将该二维矩阵S_x中所有M个元素按降序排列获得一维向量A，并对一维向量A中所有M个元素求和，获得4c)顺序累加一维向量A中前m个元素，直至取第m个元素值A(m)作为门限γ的取值，即γ＝A(m)，其中η代表残差能量控制参数，取值为η＝0.7；

4d)用非平稳信号h(t)的时频分布能量图中每个点与门限γ进行比较，将所有小于门限γ的点取值为0，其余点取值为1，得到非平稳信号h(t)的二值化时频图。

步骤5，确定信号的瞬时频率估计。

根据二值化处理后的非平稳信号h(t)的二值化时频图，获得表征非平稳信号h(t)时频分布的时频曲线，将该时频曲线作为非平稳信号h(t)的瞬时频率估计结果。

本发明的效果可以通过以下仿真结果进一步说明：

1、仿真条件

下述仿真中均取非平稳信号的时间t为0＜t≤2s，采样频率fs＝512Hz，附加噪声w(t)为高斯白噪声。

2、仿真内容

仿真1，输入信噪比SNR＝0dB的线性调频信号h(t)＝sin(2π(40t²+50t))+w(t)，根据输入信号，在最优化中心频率控制参数shiftN'＝450，最优化窗长控制参数r'＝0.9及最优化窗长变化趋势控制参数p'＝0.775的条件下，仿真信号的聚焦S变换时频分布能量，结果如图2所示；再对该时频分布能量图进行二值化处理，获得信号的瞬时频率估计，如图3所示。

仿真2，输入信噪比SNR＝5dB的正弦调频信号h(t)＝sin(2π·(5sin(10t)+120t))，根据输入信号，在最优化中心频率控制参数shiftN'＝500，最优化窗长控制参数r'＝1及最优化窗长变化趋势控制参数p'＝0.62的条件下，仿真信号的聚焦S变换时频分布能量，结果如图4所示；再对该时频分布能量图进行二值化处理，获得信号的瞬时频率估计，如图5所示。

仿真3，输入多分量信号h(t)＝sin(2π·(20t²+30t))+sin(2π·(50(t-1)³+200t))，其信噪比为SNR＝5dB；根据输入信号，在最优化中心频率控制参数shiftN'＝650，最优化窗长控制参数r'＝1及最优化窗长变化趋势控制参数p'＝0.7的条件下，仿真信号的聚焦S变换时频分布能量，结果如图6所示；再对该时频分布能量图进行二值化处理，获得信号的瞬时频率估计，如图7所示。

图2-图5表明，对于复杂调频模式的单分量信号，本发明能够在低信噪比下获得对该类信号的精确瞬时频率估计。

图6-图7表明，对于混合调频模式的多分量信号，本发明能够在精确瞬时频率估计的同时，去除交叉项的干扰。

以上描述仅是本发明的一个具体实例，并未构成对本发明的任何限制。显然对于本领域的专业人员来说，在了解了本发明内容和原理后，都可能在不背离本发明原理、结构的情况下，进行形式和细节上的各种修正和改变，但是这些基于本发明思想的修正和改变仍在本发明的权利要求保护范围之内。

Claims (4)

1.一种基于聚焦S变换的多分量非平稳信号瞬时频率估计方法，其特征在于：包括以下步骤：

(1)输入非平稳信号的时域形式h(t)，其中t代表时间，非平稳信号包括线性调频信号、正弦调频信号、线性调频与二次调频的多分量混合信号；

(2)对非平稳信号的时域形式h(t)进行傅里叶变换，得到非平稳信号的频域形式H(f)，其中f代表频率；

(3)定义聚焦S变换窗函数的频域形式为：

W(α,f)＝k·exp((f+shiftN)^2pα^2r/2) <1>

其中，k代表归一化参数，shiftN代表中心频率控制参数，r代表窗长控制参数，p代表窗长变化趋势控制参数，α代表频移；

(4)根据式<1>求出聚焦S变换的时频表达式：

$S(t,f)={\&Integral;}_{-\mathrm{\&infin;}}^{\mathrm{\&infin;}}H(\mathrm{\&alpha;}+f)\&CenterDot;W(\mathrm{\&alpha;},f)e^{2\mathrm{\&pi;}i\mathrm{\&alpha;}t}d\mathrm{\&alpha;}---<2>$

其中，H(α+f)代表加入频移α后非平稳信号的频域形式，W(α,f)代表聚焦S变换窗函数的频域形式，i代表虚数单位；

(5)根据式<2>求出聚集性准则判别式CM如下：

$CM=\frac{1}{{\&Integral;}_0^{+\mathrm{\&infin;}}{\&Integral;}_{-\mathrm{\&infin;}}^{+\mathrm{\&infin;}}|S(t,f){|}^{0.2}dtdf};---<3>$

(6)依照式<3>优化中心频率控制参数shiftN，窗长控制参数r和窗长变化趋势控制参数p，即在聚集性准则判别式CM取最大值时，得到最优化中心频率控制参数shiftN'，最优化窗长控制参数r'和最优化窗长变化趋势控制参数p'；

(7)将最优化中心频率控制参数shiftN'，最优化窗长控制参数r'和最优化窗长变化趋势控制参数p'代入式<1>，得到聚焦S变换窗函数的频域聚焦形式W'(α,f)；

(8)根据聚焦S变换窗函数的频域聚焦形式W'(α,f)，按照式<2>得到非平稳信号最优化聚焦S变换的时频表达式S'(t,f)，根据S'(t,f)获得非平稳信号h(t)的时频分布能量图；

(9)根据时频域支撑区间置信度准则设定门限γ，以门限γ为二值化处理门限对非平稳信号h(t)的时频分布能量图进行二值化处理，获得表征非平稳信号h(t)时频分布的时频曲线，将该时频曲线作为非平稳信号h(t)的瞬时频率估计结果。

2.根据权利要求1所述的基于聚焦S变换的多分量非平稳信号瞬时频率估计方法，其特征在于，步骤(6)所述的依照聚集性准则判别式CM优化中心频率控制参数shiftN，窗长控制参数r和窗长变化趋势控制参数p，按如下步骤进行：

(6a)设定初始值，包括：初始中心频率控制参数shiftN＝0，初始窗长控制参数r＝2，初始窗长变化趋势控制参数p＝2，初始聚集性准则判别式CM＝inf，其中inf为无限大；

(6b)设定参数优化区间，包括：中心频率控制参数shiftN优化区间为[-2fs，2fs]，窗长控制参数r和窗长变化趋势控制参数p的优化区间均为[0.6，1.3]，其中fs为采样频率；

(6c)在聚集性准则判别式CM最大化的约束条件下，在各自的优化区间内，搜索得到最优中心频率控制参数shiftN'，最优窗长控制参数r'和最优窗长变化趋势控制参数p'。

3.根据权利要求1所述的基于聚焦S变换的多分量非平稳信号瞬时频率估计方法，其特征在于，步骤(9)所述的根据时频域支撑区间置信度准则设定门限γ，按以下步骤设定：

(9a)根据非平稳信号h(t)的时频分布能量图得到二维矩阵S_x；

(9b)将该二维矩阵S_x中所有M个元素按降序排列获得一维向量A，并对一维向量A中所有M个元素求和，获得

(9c)顺序累加一维向量A中前m个元素，直至，取第m个元素值A(m)作为门限γ的取值，即γ＝A(m)，其中η代表残差能量控制参数，取值为η＝0.7。

4.根据权利要求1所述的基于聚焦S变换的多分量非平稳信号瞬时频率估计方法，其特征在于，步骤(9)所述的以门限γ为二值化处理门限对非平稳信号h(t)的时频分布能量图进行二值化处理，是用非平稳信号h(t)的时频分布能量图中每个点与门限γ进行比较，将所有小于门限γ的点取值为0，其余点取值为1。