一种基于自适应窗长的改进同步提取时频分析方法与流程

1.本发明属于信号处理技术，具体涉及一种基于自适应窗长的改进同步提取时频分析方法。


背景技术：

2.时频分析(tfa)方法是描述非平稳信号时变特征的有效工具，在通信、雷达和声纳等工程领域得到了广泛的应用。经典的线性方法，如短时傅里叶变换(stft)和小波变换(wt)，可以将一维时间序列信号扩展到二维时频(tf)平面，从这个时频平面中可以观察信号的时变特征，并实现多分量信号的分解。然而，受海森伯格不确定原理的限制，传统方法产生的时频表示通常是模糊的，不能为时变信号提供精确的时频描述。
3.传统tfa方法的缺点严重限制了它们在实际数据处理中的应用，为解决这一问题，近年来发展了许多先进的后处理方法，如重分配法(rm)、同步压缩变换(sst)、和同步提取变换(set)。rm首先基于时频相位信息计算每个时频点重新分配的位置，再沿二维时频方向将时频系数重新分配到瞬时轨迹中，因此可以大大提高时频分辨率，但rm缺乏重构信号的能力。sst仅在频率方向上将时频系数挤压到瞬时轨迹中，因此在提高经典tfa方法性能的同时，可以精确地重构信号。和sst相比，set只保留与信号时变特征最相关的时频系数，具有更高的能量聚集性和噪声鲁棒性。以上后处理方法都是基于固定窗长的stft，对于强调频信号和多分量信号的分析性能难以令人满意。
4.公开号为cn 107576943a的中国专利，公开了一种基于瑞利熵的自适应时频同步压缩方法，该方法首先利用瑞利熵衡量信号局部能量的聚集性,估计出适合信号时频分析的最佳窗口参数序列,通过时频同步压缩提高聚集性,实现多分量信号的分离。该方法具有较高的时频分辨率，自适应时窗也能适应更多形式的信号，但对于频率变化差异大的多分量信号，难以同时适应各个分量的变化。
5.公开号为cn108694392a的中国专利，公开了一种高精度同步提取广义s变换时频分析方法，该方法将三参数广义s变换和同步提取相结合，可以根据信号频率调节窗长，具有较高的时频分辨率，但其性能受参数影响很大，难以确定最优参数，不适用于强调频信号，且丢失了信号强度信息。


技术实现要素：

6.本发明的目的在于提供一种基于自适应窗长的改进同步提取时频分析方法。
7.实现本发明目的的技术解决方案为：一种基于自适应窗长的改进同步提取时频分析方法，包括以下步骤：
8.步骤1、对输入信号做短时傅里叶变换；
9.步骤2、根据二阶局部调制算子估计每个时频点的频率变化率；
10.步骤3、通过剥离算法迭代求解优化问题，搜索局部能量峰值，提取时频脊线，估计信号分量个数，得到各分量的瞬时频率变化率；
11.步骤4、对于单分量信号，根据信号的瞬时频率变化率设计自适应时窗；
12.步骤5、对于多分量信号，根据提取出的脊线将时频域划分成不同区域，针对不同分量的变化趋势设计自适应时频窗；
13.步骤6、以自适应窗函数构造新的短时傅里叶变换，得到聚集性更高的时频系数；
14.步骤7、通过检测时频系数在频率方向的局部极大值，构造新的频率重分配算子，得到信号瞬时频率的估计值；
15.步骤8、运用改进的同步提取算子提取时频系数，去除模糊的时频能量，并且保留信号强度信息。
16.一种电子设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述程序时实现上述的基于自适应窗长的改进同步提取时频分析方法。
17.一种计算机可读存储介质，其上存储有计算机程序，该程序被处理器执行时实现上述的基于自适应窗长的改进同步提取时频分析方法。
18.本发明与现有技术相比，其显著优点为：1)本发明具有自适应时频窗，可以根据各个分量的调频强度自适应改变窗长，使得各分量时频表示清晰，时频聚集性高，适应多种形式的信号；2)本发明可以更好地提取弱幅度信号，并保留各分量的强度信息。
19.下面结合附图对本发明作进一步详细描述。
附图说明
20.图1为本发明基于自适应窗长的改进同步提取时频分析方法的流程图。
21.图2是本发明实施例1的时频分析结果图。
22.图3是同步压缩实施例1的时频分析结果图。
23.图4是本发明和同步压缩施例1的瞬时频率估计绝对误差对比图。
24.图5是本发明实施例2的时频分析结果图。
25.图6是同步压缩实施例2的时频分析结果图。
具体实施方式
26.结合图1所示，本发明为一种基于自适应窗长的改进同步提取时频分析方法，包括以下步骤：
27.步骤1、对输入信号x(t)做初步短时傅里叶变换，得到短时傅里叶系数具体为：
[0028][0029]
其中，t是时间变量；ω是频率变量；τ是积分变量；i是复数单位；g(t)为高斯窗函数，取为这里σ是高斯窗的标准方差，代表高斯窗窗长；
[0030]
步骤2、根据二阶局部调制算了估计每个时频点的频率变化率，具体为：
[0031][0032]
其中，分别表示窗函数为g
′
(t)，g
′
(t)，tg(t)，tg
′
(t)的短时傅里叶变换；
[0033]
步骤3、通过剥离算法迭代求解下述优化问题，搜索局部能量峰值，提取信号k条时频脊线：
[0034][0035][0036][0037]
其中，t
n
是第n个时刻，n∈{1，...，n}，n是总时刻数；l
k
是第k条时频脊线，k∈{1，...，k}，k是信号总分量数；l
k
(t
n
)是第k条时频脊线在t
n
时刻的频点；表示剥离k
‑
1条脊线及其邻域后的短时傅里叶系数，t为时间，ω是频率变量；是第k条时频脊线在t
n
时刻的短时傅里叶系数；λ是放大系数；是二阶局部调制算子；γ
n
是t
n
时刻的自定义阈值，可以设置为具体步骤为：
[0038]
步骤3
‑
1、在t1时刻，沿频率方向搜索全局能量峰值确定阈值γ1；
[0039]
步骤3
‑
2、从t1到t
n
时刻重复步骤3
‑
1，得到各个时刻的阈值γ
n
；
[0040]
步骤3
‑
3、在t1时刻，沿频率方向搜索局部峰值，且幅度需大于阈值γ1，得到满足条件的频点(ω1、ω2…
)；
[0041]
步骤3
‑
4、确定脊线起始位置(t1，ω1)，在t2时刻及约束频域区间时刻及约束频域区间搜索能量最大处，且幅度需大于阈值γ2，得到满足条件的时频点(t2，ω2)；
[0042]
步骤3
‑
5、从t1到t
n
‑1时刻重复步骤3
‑
4，得到第1条时频脊线l1；
[0043]
步骤3
‑
6、将脊线l1及其领域系数从初始短时傅里叶系数中去除，记为
[0044]
步骤3
‑
7、重复步骤3
‑
3至步骤3
‑
6，得到k条时频脊线；
[0045]
步骤4、对于单分量信号，根据信号的瞬时频率变化率设计自适应时窗，具体步骤为：
[0046]
步骤4
‑
1、构造时变窗长σ(t)
[0047][0048]
其中，l(t)是信号的时频脊线；是二阶局部调制算子；θ是窗长调节因子；
[0049]
步骤4
‑
2、将σ(t)代入高斯窗函数构造自适应窗函数g
σ
(t)为：
[0050][0051]
步骤5、对于多分量信号，根据提取出的脊线将时频域划分成不同区域，针对不同分量的变化趋势设计自适应时频窗，具体步骤为：
[0052]
步骤5
‑
1、构造随时间、频率变化的窗长σ
k
(t，ω)
[0053][0054]
其中，是二阶局部调制算子；θ是窗长调节因子；每个分量在其脊点附近的频率区间[b
k
，b
k 1
]内窗长是固定，l
k
(t)是提取的第k条脊线，k∈{1，
…
，k}，k是信号总分量数，假设l1(t)≤l2(t)≤
…
l
k
(t)，那么可取
[0055][0056]
其中m为总频点数；
[0057]
步骤5
‑
2、将σ
k
(t，ω)代入高斯窗函数构造自适应窗函数g
σ
(t，ω)为：
[0058][0059]
步骤6、以自适应窗函数构造新的短时傅里叶变换，得到聚集性更强的时频系数f
a
(t，ω)，具体为：
[0060][0061]
其中，t是时间变量；ω是频率变量；τ是积分变量；i是复数单位，g
σ
(t，ω)是自适应窗函数；
[0062]
步骤7、通过检测时频系数在频率方向的局部极大值，构造新的频率重分配算子，得到信号瞬时频率的估计值具体为；
[0063][0064]
其中，f
a
(t，ω)是自适应时频系数，δ是定义的频谱间隔；
[0065]
步骤8、运用改进的同步提取算子aseo提取精确的时频系数，去除模糊的时频能量，并且保留信号强度信息，具体步骤为：
[0066]
步骤8
‑
1、根据得到的瞬时频率估计，计算改进同步提取算子
[0067][0068]
[0069]
其中，η是频率变量；是信号瞬时频率估计值；f
a
(t，ω)是自适应时频系数；m是时频系数的均值。
[0070]
步骤8
‑
2、与自适应短时傅里叶系数f
a
(t，ω)级联得到自适应改进同步提取变换：
[0071][0072]
其中，t是时间变量；η是频率变量；ω是积分变量；为信号瞬时频率估计值；f
a
(t，ω)是自适应时频系数。
[0073]
本发明还提供一种电子设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述程序时实现上述的基于自适应窗长的改进同步提取时频分析方法。
[0074]
以及一种计算机可读存储介质，其上存储有计算机程序，该程序被处理器执行时实现上述的基于自适应窗长的改进同步提取时频分析方法。
[0075]
下面结合两个实施例对本发明作进一步详细描述。
[0076]
实施例1
[0077]
仿真信号为两个非线调频信号x1(t)、x2(t)的叠加，解析式为：
[0078]
x1(t)＝(1
‑
0.6e
‑
0.3t
)sin(2π(38t 3sin(1.5t)))
[0079]
x2(t)＝e
‑
0.4t
sin(2π(15t 2sin(4.5t)))
[0080]
采样频率为100hz。图2为采用自适应窗长的改进同步提取变换得到的时频谱，可以看出本发明可以同时清晰地刻画两个分量的频率变化，并且能够反映信号的强度变化。图3为原同步压缩变换的分析结果，对比图2、图3可以看出同步压缩变换对于弱调频分量具有较清晰的时频表示，但对于强调频分量，时频谱在频率快速变换处能量发散，而本发明的时频结果聚集性更强。图4是本发明和同步压缩变换估计强调频分量瞬时频率的绝对误差曲线，可以看出本发明的瞬时频率估计误差更小，更精确。
[0081]
实施例2
[0082]
仿真信号为三个具有明显不同幅度的非线调频信号x1(t)、x2(t)、x3(t)的叠加，解析式为：
[0083]
x1(t)＝sin(2π(38t 3sin(1.5t)))
[0084]
x2(t)＝0.1sin(2π(25t 3sin(1.5t)))
[0085]
x3(t)＝0.01sin(2π(13t 3sin(1.5t)))
[0086]
采样频率为100hz。图5为采用自适应窗长的改进同步提取变换得到的时频谱，可以看出本发明可以清晰地表示弱分量信号。图6为原同步压缩变换的分析结果，对比图5、图6可以看出同步压缩变换对于弱调频分量已经几乎观察不到，而本发明可以清楚地表征每个分量的时变轨迹，时频聚集性高，并且可以反映各分量的强度大小。
[0087]
本发明可以适应不同分量信号的调频强度，在时间、频率方向上自适应调整窗长，可以清晰表征弱信号分量，能够反映信号强度大小，适应多种信号形式，时变特征刻画精细，具有很高的时频分辨率。
[0088]
以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为
本发明的保护范围。