自适应的时频支持的调频信号分解方法与系统

1.本发明涉及工程信号处理的技术领域，具体地，涉及自适应的时频支持的调频信号分解方法与系统。


背景技术：

2.在实际的工业领域，存在大量的瞬时频率和瞬时幅值随时间变化的调幅调频信号(am-fm信号)。比如旋转机械的变转速下的振动信号，切削装备的主轴振动信号，声音信号和雷达探测信号等。许多实际测量的调幅调频信号不可避免地受到噪声的干扰，特别是在恶劣的工况下，信号甚至受到强噪声的干扰。强噪声会对信号的瞬时频率(if)和瞬时幅值(ia)的估计产生很大的干扰，存在的信号分解方法不能准确地分解具有强噪声的复杂am-fm信号，在强噪声的干扰下容易产生模态混叠问题。对具有强噪声的复杂多成分am-fm信号进行分解能有效地估计每个子信号的瞬时频率并重构每个子信号，进而有效地滤除强噪声的干扰，揭示信号的时频变化规律，有效地理解信号所包含的关键信息。
3.针对上述中的相关技术，发明人认为存在信号分解方法的不足的问题。因此，需要提出一种新的技术方案以改善上述技术问题。


技术实现要素：

4.针对现有技术中的缺陷，本发明的目的是提供一种自适应的时频支持的调频信号分解方法与系统。
5.根据本发明提供的一种自适应的时频支持的调频信号分解方法，所述方法包括如下步骤：
6.步骤s1：计算信号y的参数化时频变换ptft；
7.步骤s2：从信号y的参数化时频变换的时频表示tfr中提取能量最大的子信号ym的粗糙的瞬时频率脊度曲线ifm；
8.步骤s3：将提取的该子信号ym的瞬时频率曲线ifm作为信号分解的初始瞬时频率fm；
9.步骤s4：根据该子信号ym的瞬时频率fm计算该子信号的瞬时频率核函数矩阵φm；
10.步骤s5：计算该子信号ym的解调信号um和vm；
11.步骤s6：计算该子信号ym的初始瞬时频率增量
12.步骤s7：计算该子信号ym的优化的初始瞬时频率增量δfm；
13.步骤s8：更新该子信号ym的瞬时频率fm；
14.步骤s9：重构该子信号ym；
15.步骤s10：重复步骤s3-步骤s9，直到重构的子信号ym满足终止迭代的条件；
16.步骤s11：采用k阶傅里叶级数拟合该子信号ym估计的瞬时频率fm；
17.步骤s12：更新信号y的参数化时频变换的核函数参数，并对信号y再次进行参数化时频变换；
18.步骤s13：重复步骤s1-步骤s12，直到该子信号ym的瞬时频率fm满足终止迭代的条件；
19.步骤s14：重构该子信号ym；
20.步骤s15：将该子信号ym从当前剩余信号rm中移除，然后继续分解下一个子信号y
m+1
，并重复步骤s1-步骤s14，直到原始信号中每个子信号被完全分解。
21.优选地，所述步骤s1中ptft的核函数参数为傅里叶级数，其计算公式为：
[0022][0023]
其中的计算公式为：
[0024][0025]
式中，k是傅里叶级数的阶数，f0为基频，f0＝1/(2t)，t为信号持续时间，＝1/(2t)，t为信号持续时间，为傅里叶级数的系数；wh的计算公式为：
[0026][0027]
所述步骤s2中通过求解e(ifm)的局部最大值得到子信号ym的瞬时频率脊度曲线ifm；e(ifm)的计算公式为：
[0028][0029]
式中s(t,ifm(t))表示信号的ptft，λ和β为惩罚系数；
[0030]
所述步骤s4中的瞬时频率核函数矩阵φm的计算公式为：
[0031]
φm＝[cm,sm]
[0032]
式中cm表示正弦函数矩阵，sm表示余弦函数矩阵；
[0033]cm
的计算公式为：
[0034]cm
＝diag[cos(θm(t1))
…
cos(θm(tn))]
[0035]
式中θm表示子信号ym的瞬时相位；
[0036]
瞬时相位θm的计算公式为：
[0037][0038]
式中表示子信号ym估计的瞬时频率；
[0039]
sm的计算公式为：
[0040]
sm＝diag[sin(θm(t1))
…
sin(θm(tn))]
[0041]
式中θm表示子信号ym的瞬时相位。
[0042]
优选地，所述步骤s5中子信号ym的解调信号um和vm的计算公式为：
[0043][0044]
式中n表示迭代次数，φm表示瞬时频率核函数矩阵，ρ表示惩罚系数，d表示二阶差分矩阵ω的对角矩阵rm表示当前的剩余信号；
[0045]
二阶差分矩阵ω的表达式为：
[0046][0047]
ω为n行n列矩阵，n表示信号的离散采样点数；
[0048]
所述步骤s6中子信号ym的初始瞬时频率增量的计算公式为：
[0049][0050]
式中和表示子信号ym在第n次迭代时的解调信号；
[0051]
所述步骤s7中子信号ym的优化的初始瞬时频率增量δfm的计算公式为：
[0052][0053]
式中ω表示二阶差分矩阵，i表示单位矩阵。
[0054]
优选地，所述步骤s8中子信号ym的瞬时频率fm的更新公式为：
[0055][0056]
式中n表示迭代次数；
[0057]
所述步骤s9中子信号ym的重构计算公式为：
[0058]
ym＝φmxm[0059]
式中φm表示子信号ym的瞬时频率核函数矩阵，xm表示子信号ym的解调信号；
[0060]
所述步骤s10中最开始设置迭代次数n为1，完成一次迭代后，将迭代次数n加1；重构的子信号ym的终止迭代条件为：
[0061][0062]
式中ε1表示终止迭代阈值；重构子信号ym的终止迭代阈值ε1设置为10-8
。
[0063]
优选地，所述步骤s11中用于拟合子信号ym的瞬时频率fm的k阶傅里叶级数的公式为：
[0064][0065]
式中为傅里叶级数的系数；
[0066]
所述步骤s13中，子信号ym的瞬时频率fm的终止迭代条件为：
[0067]
[0068]
式中h为迭代次数，ε2表示终止迭代阈值；终止迭代阈值ε2设置为10-8
；
[0069]
所述步骤s14中子信号ym的重构计算公式为：
[0070]
ym＝φmxm[0071]
式中φm表示子信号ym的瞬时频率核函数矩阵，xm表示子信号ym的解调信号；
[0072]
所述步骤s15中最开始时设置m＝1，当前剩余信号r1等于原始信号y；在完成对第m个子信号的分解后，移除子信号ym，此时当前剩余信号r
m+1
＝r
m-ym，并设置m＝m+1。
[0073]
本发明还提供一种自适应的时频支持的调频信号分解系统，所述系统包括如下模块：
[0074]
模块m1：计算信号y的参数化时频变换ptft；
[0075]
模块m2：从信号y的参数化时频变换的时频表示tfr中提取能量最大的子信号ym的粗糙的瞬时频率脊度曲线ifm；
[0076]
模块m3：将提取的该子信号ym的瞬时频率曲线ifm作为信号分解的初始瞬时频率fm；
[0077]
模块m4：根据该子信号ym的瞬时频率fm计算该子信号的瞬时频率核函数矩阵φm；
[0078]
模块m5：计算该子信号ym的解调信号um和vm；
[0079]
模块m6：计算该子信号ym的初始瞬时频率增量
[0080]
模块m7：计算该子信号ym的优化的初始瞬时频率增量
[0081]
模块m8：更新该子信号ym的瞬时频率fm；
[0082]
模块m9：重构该子信号ym；
[0083]
模块m10：重复调用模块m3-模块m9，直到重构的子信号ym满足终止迭代的条件；
[0084]
模块m11：采用k阶傅里叶级数拟合该子信号ym估计的瞬时频率fm；
[0085]
模块m12：更新信号y的参数化时频变换的核函数参数，并对信号y再次进行参数化时频变换；
[0086]
模块m13：重复调用模块m1-模块m12，直到该子信号ym的瞬时频率fm满足终止迭代的条件；
[0087]
模块m14：重构该子信号ym；
[0088]
模块m15：将该子信号ym从当前剩余信号rm中移除，然后继续分解下一个子信号y
m+1
，并重复调用模块m1-模块m14，直到原始信号中每个子信号被完全分解。
[0089]
优选地，所述模块m1中ptft的核函数参数为傅里叶级数，其计算公式为：
[0090][0091]
其中的计算公式为：
[0092]
[0093]
式中，k是傅里叶级数的阶数，f0为基频，f0＝1/(2t)，t为信号持续时间，＝1/(2t)，t为信号持续时间，为傅里叶级数的系数；wh的计算公式为：
[0094][0095]
所述模块m2中通过求解e(ifm)的局部最大值得到子信号ym的瞬时频率脊度曲线ifm；e(ifm)的计算公式为：
[0096][0097]
式中s(t,ifm(t))表示信号的ptft，λ和β为惩罚系数；
[0098]
所述模块m4中的瞬时频率核函数矩阵φm的计算公式为：
[0099]
φm＝[cm,sm]
[0100]
式中cm表示正弦函数矩阵，sm表示余弦函数矩阵；
[0101]cm
的计算公式为：
[0102]cm
＝diag[cos(θm(t1))
…
cos(θm(tn))]
[0103]
式中θm表示子信号ym的瞬时相位；
[0104]
瞬时相位θm的计算公式为：
[0105][0106]
式中表示子信号ym估计的瞬时频率；
[0107]
sm的计算公式为：
[0108]
sm＝diag[sin(θm(t1))
…
sin(θm(tn))]
[0109]
式中θm表示子信号ym的瞬时相位。
[0110]
优选地，所述模块m5中子信号ym的解调信号um和vm的计算公式为：
[0111][0112]
式中n表示迭代次数，φm表示瞬时频率核函数矩阵，ρ表示惩罚系数，d表示二阶差分矩阵ω的对角矩阵rm表示当前的剩余信号；
[0113]
二阶差分矩阵ω的表达式为：
[0114][0115]
ω为n行n列矩阵，n表示信号的离散采样点数；
[0116]
所述模块m6中子信号ym的初始瞬时频率增量的计算公式为：
[0117]
[0118]
式中和表示子信号ym在第n次迭代时的解调信号；
[0119]
所述模块m7中子信号ym的优化的初始瞬时频率增量δfm的计算公式为：
[0120][0121]
式中ω表示二阶差分矩阵，i表示单位矩阵。
[0122]
优选地，所述模块m8中子信号ym的瞬时频率fm的更新公式为：
[0123][0124]
式中n表示迭代次数；
[0125]
所述模块m9中子信号ym的重构计算公式为：
[0126]
ym＝φmxm[0127]
式中φm表示子信号ym的瞬时频率核函数矩阵，xm表示子信号ym的解调信号；
[0128]
所述模块m10中最开始设置迭代次数n为1，完成一次迭代后，将迭代次数n加1；重构的子信号ym的终止迭代条件为：
[0129][0130]
式中ε1表示终止迭代阈值；重构子信号ym的终止迭代阈值ε1设置为10-8
。
[0131]
优选地，所述模块m11中用于拟合子信号ym的瞬时频率fm的k阶傅里叶级数的公式为：
[0132][0133]
式中为傅里叶级数的系数；
[0134]
所述模块m13中，子信号ym的瞬时频率fm的终止迭代条件为：
[0135][0136]
式中h为迭代次数，ε2表示终止迭代阈值；终止迭代阈值ε2设置为10-8
；
[0137]
所述模块m14中子信号ym的重构计算公式为：
[0138]
ym＝φmxm[0139]
式中φm表示子信号ym的瞬时频率核函数矩阵，xm表示子信号ym的解调信号；
[0140]
所述模块m15中最开始时设置m＝1，当前剩余信号r1等于原始信号y；在完成对第m个子信号的分解后，移除子信号ym，此时当前剩余信号r
m+1
＝r
m-ym，并设置m＝m+1。
[0141]
与现有技术相比，本发明具有如下的有益效果：
[0142]
1、本发明的自适应的时频支持的调频信号分解方法与系统结合了if反正切解调优化和信号的时频表示迭代优化，采用信号的时频表示信息辅助每个子信号的if的迭代，能有效地抑制强噪声对信号分解的干扰，进而准确地估计信号的if；
[0143]
2、本发明的自适应的时频支持的调频信号分解方法与系统能有效地分解具有强噪声的am-fm信号，准确地估计信号的if和重构子信号。
附图说明
[0144]
通过阅读参照以下附图对非限制性实施例所作的详细描述，本发明的其它特征、
目的和优点将会变得更明显：
[0145]
图1为本发明的自适应的时频支持的调频信号分解方法的流程图；
[0146]
图2为本发明的具有强噪声的仿真信号的时域曲线图；
[0147]
图3为本发明的具有强噪声的仿真信号的短时傅里叶变换的时频表示图；
[0148]
图4为本发明的仿真信号的实际的瞬时频率图；
[0149]
图5为本发明的仿真信号的每个子信号的实际的时域曲线图；
[0150]
图6为本发明的自适应的时频支持的调频信号分解方法估计的仿真信号的if曲线图；
[0151]
图7为本发明的自适应的时频支持的调频信号分解方法重构的仿真信号2个子信号的时域曲线图；
[0152]
图8为本发明受到噪声干扰的外圈故障的轴承在变转速运行时的轴承座振动信号的时域曲线图；
[0153]
图9为本发明轴承座振动信号的短时傅里叶变换的时域曲线图；
[0154]
图10为本发明的自适应的时频支持的调频信号分解方法对这个轴承座振动信号进行信号分解得到的时频表示图。
具体实施方式
[0155]
下面结合具体实施例对本发明进行详细说明。以下实施例将有助于本领域的技术人员进一步理解本发明，但不以任何形式限制本发明。应当指出的是，对本领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干变化和改进。这些都属于本发明的保护范围。
[0156]
实施例1：
[0157]
根据本发明提供的一种自适应的时频支持的调频信号分解方法，所述方法包括如下步骤：
[0158]
步骤s1：计算信号y的参数化时频变换ptft；ptft的核函数参数为傅里叶级数，其计算公式为：
[0159][0160]
其中的计算公式为：
[0161][0162]
式中，k是傅里叶级数的阶数，f0为基频，f0＝1/(2t)，t为信号持续时间，＝1/(2t)，t为信号持续时间，为傅里叶级数的系数；wh的计算公式为：
[0163]
[0164]
步骤s2：从信号y的参数化时频变换的时频表示tfr中提取能量最大的子信号ym的粗糙的瞬时频率脊度曲线ifm；通过求解e(ifm)的局部最大值得到子信号ym的瞬时频率脊度曲线ifm；e(ifm)的计算公式为：
[0165][0166]
式中s(t,ifm(t))表示信号的ptft，λ和β为惩罚系数。
[0167]
步骤s3：将提取的该子信号ym的瞬时频率曲线ifm作为信号分解的初始瞬时频率fm；
[0168]
步骤s4：根据该子信号ym的瞬时频率fm计算该子信号的瞬时频率核函数矩阵φm；瞬时频率核函数矩阵φm的计算公式为：
[0169]
φm＝[cm,sm]
[0170]
式中cm表示正弦函数矩阵，sm表示余弦函数矩阵；
[0171]cm
的计算公式为：
[0172]cm
＝diag[cos(θm(t1))
…
cos(θm(tn))]
[0173]
式中θm表示子信号ym的瞬时相位；
[0174]
瞬时相位θm的计算公式为：
[0175][0176]
式中表示子信号ym估计的瞬时频率；
[0177]
sm的计算公式为：
[0178]
sm＝diag[sin(θm(t1))
…
sin(θm(tn))]
[0179]
式中θm表示子信号ym的瞬时相位。
[0180]
步骤s5：计算该子信号ym的解调信号um和vm；子信号ym的解调信号um和vm的计算公式为：
[0181][0182]
式中n表示迭代次数，φm表示瞬时频率核函数矩阵，ρ表示惩罚系数，d表示二阶差分矩阵ω的对角矩阵rm表示当前的剩余信号；
[0183]
二阶差分矩阵ω的表达式为：
[0184][0185]
ω为n行n列矩阵，n表示信号的离散采样点数。
[0186]
步骤s6：计算该子信号ym的初始瞬时频率增量子信号ym的初始瞬时频率增量的计算公式为：
[0187][0188]
式中和表示子信号ym在第n次迭代时的解调信号。
[0189]
步骤s7：计算该子信号ym的优化的初始瞬时频率增量δfm；子信号ym的优化的初始瞬时频率增量δfm的计算公式为：
[0190][0191]
式中ω表示二阶差分矩阵，i表示单位矩阵。
[0192]
步骤s8：更新该子信号ym的瞬时频率fm；子信号ym的瞬时频率fm的更新公式为：
[0193][0194]
式中n表示迭代次数。
[0195]
步骤s9：重构该子信号ym；子信号ym的重构计算公式为：
[0196]
ym＝φmxm[0197]
式中φm表示子信号ym的瞬时频率核函数矩阵，xm表示子信号ym的解调信号。
[0198]
步骤s10：重复步骤s3-步骤s9，直到重构的子信号ym满足终止迭代的条件；最开始设置迭代次数n为1，完成一次迭代后，将迭代次数n加1；重构的子信号ym的终止迭代条件为：
[0199][0200]
式中ε1表示终止迭代阈值；重构子信号ym的终止迭代阈值ε1设置为10-8
。
[0201]
步骤s11：采用k阶傅里叶级数拟合该子信号ym估计的瞬时频率fm；用于拟合子信号ym的瞬时频率fm的k阶傅里叶级数的公式为：
[0202][0203]
式中为傅里叶级数的系数。
[0204]
步骤s12：更新信号y的参数化时频变换的核函数参数，并对信号y再次进行参数化时频变换；
[0205]
步骤s13：重复步骤s1-步骤s12，直到该子信号ym的瞬时频率fm满足终止迭代的条件；子信号ym的瞬时频率fm的终止迭代条件为：
[0206][0207]
式中h为迭代次数，ε2表示终止迭代阈值；终止迭代阈值ε2设置为10-8
。
[0208]
步骤s14：重构该子信号ym；子信号ym的重构计算公式为：
[0209]
ym＝φmxm[0210]
式中φm表示子信号ym的瞬时频率核函数矩阵，xm表示子信号ym的解调信号。
[0211]
步骤s15：将该子信号ym从当前剩余信号rm中移除，然后继续分解下一个子信号y
m+1
，并重复步骤s1-步骤s14，直到原始信号中每个子信号被完全分解；所述步骤s15中最开始时设置m＝1，当前剩余信号r1等于原始信号y；在完成对第m个子信号的分解后，移除子信
号ym，此时当前剩余信号r
m+1
＝r
m-ym，并设置m＝m+1。
[0212]
本发明还提供一种自适应的时频支持的调频信号分解系统，所述自适应的时频支持的调频信号分解系统可以通过执行所述自适应的时频支持的调频信号分解方法的流程步骤予以实现，即本领域技术人员可以将所述自适应的时频支持的调频信号分解方法理解为所述自适应的时频支持的调频信号分解系统的优选实施方式。
[0213]
实施例2：
[0214]
本发明还提供一种自适应的时频支持的调频信号分解系统，所述系统包括如下模块：
[0215]
模块m1：计算信号y的参数化时频变换ptft；ptft的核函数参数为傅里叶级数，其计算公式为：
[0216][0217]
其中的计算公式为：
[0218][0219]
式中，k是傅里叶级数的阶数，f0为基频，f0＝1/(2t)，t为信号持续时间，＝1/(2t)，t为信号持续时间，为傅里叶级数的系数；wh的计算公式为：
[0220][0221]
模块m2：从信号y的参数化时频变换的时频表示tfr中提取能量最大的子信号ym的粗糙的瞬时频率脊度曲线ifm；通过求解e(ifm)的局部最大值得到子信号ym的瞬时频率脊度曲线ifm；e(ifm)的计算公式为：
[0222][0223]
式中s(t,ifm(t))表示信号的ptft，λ和β为惩罚系数。
[0224]
模块m3：将提取的该子信号ym的瞬时频率曲线ifm作为信号分解的初始瞬时频率fm；
[0225]
模块m4：根据该子信号ym的瞬时频率fm计算该子信号的瞬时频率核函数矩阵φm；瞬时频率核函数矩阵φm的计算公式为：
[0226]
φm＝[cm,sm]
[0227]
式中cm表示正弦函数矩阵，sm表示余弦函数矩阵；
[0228]cm
的计算公式为：
[0229]cm
＝diag[cos(θm(t1))
…
cos(θm(tn))]
[0230]
式中θm表示子信号ym的瞬时相位；
[0231]
瞬时相位θm的计算公式为：
[0232][0233]
式中表示子信号ym估计的瞬时频率；
[0234]
sm的计算公式为：
[0235]
sm＝diag[sin(θm(t1))
…
sin(θm(tn))]
[0236]
式中θm表示子信号ym的瞬时相位。
[0237]
模块m5：计算该子信号ym的解调信号um和vm；子信号ym的解调信号um和vm的计算公式为：
[0238][0239]
式中n表示迭代次数，φm表示瞬时频率核函数矩阵，ρ表示惩罚系数，d表示二阶差分矩阵ω的对角矩阵rm表示当前的剩余信号；
[0240]
二阶差分矩阵ω的表达式为：
[0241][0242]
ω为n行n列矩阵，n表示信号的离散采样点数。
[0243]
模块m6：计算该子信号ym的初始瞬时频率增量子信号ym的初始瞬时频率增量的计算公式为：
[0244][0245]
式中和表示子信号ym在第n次迭代时的解调信号。
[0246]
模块m7：计算该子信号ym的优化的初始瞬时频率增量δfm；子信号ym的优化的初始瞬时频率增量δfm的计算公式为：
[0247][0248]
式中ω表示二阶差分矩阵，i表示单位矩阵。
[0249]
模块m8：更新该子信号ym的瞬时频率fm；子信号ym的瞬时频率fm的更新公式为：
[0250][0251]
式中n表示迭代次数。
[0252]
模块m9：重构该子信号ym；子信号ym的重构计算公式为：
[0253]
ym＝φmxm[0254]
式中φm表示子信号ym的瞬时频率核函数矩阵，xm表示子信号ym的解调信号。
[0255]
模块m10：重复调用模块m3-模块m9，直到重构的子信号ym满足终止迭代的条件；最开始设置迭代次数n为1，完成一次迭代后，将迭代次数n加1；重构的子信号ym的终止迭代条
件为：
[0256][0257]
式中ε1表示终止迭代阈值；重构子信号ym的终止迭代阈值ε1设置为10-8
。
[0258]
模块m11：采用k阶傅里叶级数拟合该子信号ym估计的瞬时频率fm；用于拟合子信号ym的瞬时频率fm的k阶傅里叶级数的公式为：
[0259][0260]
式中为傅里叶级数的系数。
[0261]
模块m12：更新信号y的参数化时频变换的核函数参数，并对信号y再次进行参数化时频变换；
[0262]
模块m13：重复调用模块m1-模块m12，直到该子信号ym的瞬时频率fm满足终止迭代的条件；子信号ym的瞬时频率fm的终止迭代条件为：
[0263][0264]
式中h为迭代次数，ε2表示终止迭代阈值；终止迭代阈值ε2设置为10-8
。
[0265]
模块m14：重构该子信号ym；子信号ym的重构计算公式为：
[0266]
ym＝φmxm[0267]
式中φm表示子信号ym的瞬时频率核函数矩阵，xm表示子信号ym的解调信号。
[0268]
模块m15：将该子信号ym从当前剩余信号rm中移除，然后继续分解下一个子信号y
m+1
，并重复调用模块m1-模块m14，直到原始信号中每个子信号被完全分解；所述模块m15中最开始时设置m＝1，当前剩余信号r1等于原始信号y；在完成对第m个子信号的分解后，移除子信号ym，此时当前剩余信号r
m+1
＝r
m-ym，并设置m＝m+1。
[0269]
实施例3：
[0270]
为了克服存在的信号分解方法的不足，更有效地分解具有强噪声的复杂的am-fm信号，准确地估计每个子信号的瞬时频率与重构子信号，我们发明了一种自适应的时频支持的调频信号分解方法与系统，发明的信号分解方法优化了每个子信号的瞬时频率的迭代过程，将信号的时频表示信号用于辅助信号的瞬时频率的迭代，能有效地滤除强噪声的干扰。
[0271]
一种自适应的时频支持的调频信号分解方法，其步骤包括：
[0272]
步骤s1：计算信号y的参数化时频变换(ptft)。
[0273]
步骤s2：从信号y的参数化时频变换的时频表示(tfr)中提取能量最大的子信号ym的粗糙的瞬时频率脊度曲线ifm。
[0274]
步骤s3：将提取的该子信号ym的瞬时频率曲线ifm作为信号分解的初始瞬时频率fm。
[0275]
步骤s4：根据该子信号ym的瞬时频率fm计算该子信号的瞬时频率核函数矩阵φm。
[0276]
步骤s5：计算该子信号ym的解调信号um和vm。
[0277]
步骤s6：计算该子信号ym的初始瞬时频率增量
[0278]
步骤s7：计算该子信号ym的优化的初始瞬时频率增量δfm。
[0279]
步骤s8：更新该子信号ym的瞬时频率fm。
[0280]
步骤s9，重构该子信号ym。
[0281]
步骤s10，重复步骤s3至s9，直到重构的子信号ym满足终止迭代的条件。
[0282]
步骤s11，采用k阶傅里叶级数拟合该子信号ym估计的瞬时频率fm。
[0283]
步骤s12，更新信号y的参数化时频变换的核函数参数，并对信号y再次进行参数化时频变换。
[0284]
步骤s13，重复步骤s1至s12，直到该子信号ym的瞬时频率fm满足终止迭代的条件。
[0285]
步骤s14，重构该子信号ym。
[0286]
步骤s15，将该子信号ym从当前剩余信号rm中移除，然后继续分解下一个子信号y
m+1
，并重复步骤s1至s14，直到原始信号中每个子信号被完全分解。
[0287]
步骤s1中ptft的核函数参数为傅里叶级数，其计算公式为：其中的计算公式为：
[0288][0289]
式中，k是傅里叶级数的阶数，f0为基频，f0＝1/(2t)，t为信号持续时间，式中，k是傅里叶级数的阶数，f0为基频，f0＝1/(2t)，t为信号持续时间，为傅里叶级数的系数。wh的计算公式为：
[0290]
步骤s2中通过求解e(ifm)的局部最大值得到子信号ym的瞬时频率脊度曲线ifm。e(ifm)的计算公式为：
[0291][0292]
，式中s(t,ifm(t))表示信号的ptft，λ和β为惩罚系数。
[0293]
步骤s4中的瞬时频率核函数矩阵φm的计算公式为：φm＝[cm,sm]，式中cm表示正弦函数矩阵，sm表示余弦函数矩阵。
[0294]
步骤s4中的cm的计算公式为：cm＝diag[cos(θm(t1))
…
cos(θm(tn)))],式中θm表示子信号ym的瞬时相位。
[0295]
步骤s4中的瞬时相位θm的计算公式为：式中表示子信号ym估计的瞬时频率。
[0296]
步骤s4中的sm的计算公式为：sm＝diag[sin(θm(t1))
…
sin(θm(tn)))],式中θm表示子信号ym的瞬时相位。
[0297]
步骤s5中子信号ym的解调信号um和vm的计算公式为：式中n表示迭代次数，φm表示瞬时频率核函数矩阵，ρ
表示惩罚系数，d表示二阶差分矩阵ω的对角矩阵rm表示当前的剩余信号。
[0298]
步骤s5中二阶差分矩阵ω的表达式为：ω为n行n列矩阵，n表示信号的离散采样点数。
[0299]
步骤s6中子信号ym的初始瞬时频率增量的计算公式为：式中和表示子信号ym在第n次迭代时的解调信号。
[0300]
步骤s7中子信号ym的优化的初始瞬时频率增量δfm的计算公式为：式中ω表示二阶差分矩阵，i表示单位矩阵。
[0301]
步骤s8中子信号ym的瞬时频率fm的更新公式为：式中n表示迭代次数。
[0302]
步骤s9中子信号ym的重构计算公式为：ym＝φmxm，式中φm表示子信号ym的瞬时频率核函数矩阵，xm表示子信号ym的解调信号。
[0303]
步骤s10中最开始设置迭代次数n为1，完成一次迭代后，将迭代次数n加1。
[0304]
步骤s10中重构的子信号ym的终止迭代条件为：式中ε1表示终止迭代阈值。
[0305]
步骤s10中重构子信号ym的终止迭代阈值ε1一般设置为10-8
。
[0306]
步骤s11中用于拟合子信号ym的瞬时频率fm的k阶傅里叶级数的公式为：式中为傅里叶级数的系数。
[0307]
步骤s13中，子信号ym的瞬时频率fm的终止迭代条件为：式中h为迭代次数，ε2表示终止迭代阈值。
[0308]
步骤s13中,终止迭代阈值ε2一般设置为10-8
。
[0309]
步骤s14中子信号ym的重构计算公式为：ym＝φmxm，式中φm表示子信号ym的瞬时频率核函数矩阵，xm表示子信号ym的解调信号。
[0310]
步骤s15中最开始时设置m＝1，当前剩余信号r1等于原始信号y。在完成对第m个子信号的分解后，移除子信号ym，此时当前剩余信号r
m+1
＝r
m-ym，并设置m＝m+1。
[0311]
本发明提供的方法用于受到强噪声干扰的变化规律复杂的仿真信号和实际信号。图2和图3是仿真信号的时域曲线图和短时傅里叶变换的tfr。图4至图5是仿真信号的实际的if和和每个子信号的时域曲线。图6至图7是本发明的方法分解该仿真信号后的结果，结果表明本发明能精确地估计每个子信号的if和重构子信号。图8是实际的受到噪声干扰的
变转速轴承振动信号的时域曲线。图9是该轴承振动信号的短时傅里叶变换的tfr。图10是本发明的方法分解该轴承振动信号后得到的tfr，结果表明本发明能有效地滤除原始振动信号的噪声，精确地估计每个子信号的if和重构子信号，进而得到准确反映原始信号时频变化规律的能量集中性很强的tfr。
[0312]
本领域技术人员可以将本实施例理解为实施例1、实施例2的更为具体的说明。
[0313]
本领域技术人员知道，除了以纯计算机可读程序代码方式实现本发明提供的系统及其各个装置、模块、单元以外，完全可以通过将方法步骤进行逻辑编程来使得本发明提供的系统及其各个装置、模块、单元以逻辑门、开关、专用集成电路、可编程逻辑控制器以及嵌入式微控制器等的形式来实现相同功能。所以，本发明提供的系统及其各项装置、模块、单元可以被认为是一种硬件部件，而对其内包括的用于实现各种功能的装置、模块、单元也可以视为硬件部件内的结构；也可以将用于实现各种功能的装置、模块、单元视为既可以是实现方法的软件模块又可以是硬件部件内的结构。
[0314]
以上对本发明的具体实施例进行了描述。需要理解的是，本发明并不局限于上述特定实施方式，本领域技术人员可以在权利要求的范围内做出各种变化或修改，这并不影响本发明的实质内容。在不冲突的情况下，本技术的实施例和实施例中的特征可以任意相互组合。

技术特征：
1.一种自适应的时频支持的调频信号分解方法，其特征在于，所述方法包括如下步骤：步骤s1：计算信号y的参数化时频变换ptft；步骤s2：从信号y的参数化时频变换的时频表示tfr中提取能量最大的子信号y
m
的粗糙的瞬时频率脊度曲线if
m
；步骤s3：将提取的该子信号y
m
的瞬时频率曲线if
m
作为信号分解的初始瞬时频率f
m
；步骤s4：根据该子信号y
m
的瞬时频率f
m
计算该子信号的瞬时频率核函数矩阵φ
m
；步骤s5：计算该子信号y
m
的解调信号u
m
和v
m
；步骤s6：计算该子信号y
m
的初始瞬时频率增量步骤s7：计算该子信号y
m
的优化的初始瞬时频率增量δf
m
；步骤s8：更新该子信号y
m
的瞬时频率f
m
；步骤s9：重构该子信号y
m
；步骤s10：重复步骤s3-步骤s9，直到重构的子信号y
m
满足终止迭代的条件；步骤s11：采用k阶傅里叶级数拟合该子信号y
m
估计的瞬时频率f
m
；步骤s12：更新信号y的参数化时频变换的核函数参数，并对信号y再次进行参数化时频变换；步骤s13：重复步骤s1-步骤s12，直到该子信号y
m
的瞬时频率f
m
满足终止迭代的条件；步骤s14：重构该子信号y
m
；步骤s15：将该子信号y
m
从当前剩余信号r
m
中移除，然后继续分解下一个子信号y
m+1
，并重复步骤s1-步骤s14，直到原始信号中每个子信号被完全分解。2.根据权利要求1所述的自适应的时频支持的调频信号分解方法，其特征在于，所述步骤s1中ptft的核函数参数为傅里叶级数，其计算公式为：其中的计算公式为：式中，k是傅里叶级数的阶数，f0为基频，f0＝1/(2t)，t为信号持续时间，＝1/(2t)，t为信号持续时间，为傅里叶级数的系数；w
h
的计算公式为：所述步骤s2中通过求解e(if
m
)的局部最大值得到子信号y
m
的瞬时频率脊度曲线if
m
；e(if
m
)的计算公式为：式中s(t,if
m
(t))表示信号的ptft，λ和β为惩罚系数；
所述步骤s4中的瞬时频率核函数矩阵φ
m
的计算公式为：φ
m
＝[c
m
,s
m
]式中c
m
表示正弦函数矩阵，s
m
表示余弦函数矩阵；c
m
的计算公式为：c
m
＝diag[cos(θ
m
(t1))
…
cos(θ
m
(t
n
))]式中θ
m
表示子信号y
m
的瞬时相位；瞬时相位θ
m
的计算公式为：式中表示子信号y
m
估计的瞬时频率；s
m
的计算公式为：s
m
＝diag[sin(θ
m
(t1))
…
sin(θ
m
(t
n
))]式中θ
m
表示子信号y
m
的瞬时相位。3.根据权利要求1所述的自适应的时频支持的调频信号分解方法，其特征在于，所述步骤s5中子信号y
m
的解调信号u
m
和v
m
的计算公式为：式中n表示迭代次数，φ
m
表示瞬时频率核函数矩阵，ρ表示惩罚系数，d表示二阶差分矩阵ω的对角矩阵r
m
表示当前的剩余信号；二阶差分矩阵ω的表达式为：ω为n行n列矩阵，n表示信号的离散采样点数；所述步骤s6中子信号y
m
的初始瞬时频率增量的计算公式为：式中和表示子信号y
m
在第n次迭代时的解调信号；所述步骤s7中子信号y
m
的优化的初始瞬时频率增量δf
m
的计算公式为：式中ω表示二阶差分矩阵，i表示单位矩阵。4.根据权利要求1所述的自适应的时频支持的调频信号分解方法，其特征在于，所述步骤s8中子信号y
m
的瞬时频率f
m
的更新公式为：式中n表示迭代次数；
所述步骤s9中子信号y
m
的重构计算公式为：y
m
＝φ
m
x
m
式中φ
m
表示子信号y
m
的瞬时频率核函数矩阵，x
m
表示子信号y
m
的解调信号；所述步骤s10中最开始设置迭代次数n为1，完成一次迭代后，将迭代次数n加1；重构的子信号y
m
的终止迭代条件为：式中ε1表示终止迭代阈值；重构子信号y
m
的终止迭代阈值ε1设置为10-8
。5.根据权利要求1所述的自适应的时频支持的调频信号分解方法，其特征在于，所述步骤s11中用于拟合子信号y
m
的瞬时频率f
m
的k阶傅里叶级数的公式为：式中为傅里叶级数的系数；所述步骤s13中，子信号y
m
的瞬时频率f
m
的终止迭代条件为：式中h为迭代次数，ε2表示终止迭代阈值；终止迭代阈值ε2设置为10-8
；所述步骤s14中子信号y
m
的重构计算公式为：y
m
＝φ
m
x
m
式中φ
m
表示子信号y
m
的瞬时频率核函数矩阵，x
m
表示子信号y
m
的解调信号；所述步骤s15中最开始时设置m＝1，当前剩余信号r1等于原始信号y；在完成对第m个子信号的分解后，移除子信号y
m
，此时当前剩余信号r
m+1
＝r
m-y
m
，并设置m＝m+1。6.一种自适应的时频支持的调频信号分解系统，其特征在于，所述系统包括如下模块：模块m1：计算信号y的参数化时频变换ptft；模块m2：从信号y的参数化时频变换的时频表示tfr中提取能量最大的子信号y
m
的粗糙的瞬时频率脊度曲线if
m
；模块m3：将提取的该子信号y
m
的瞬时频率曲线if
m
作为信号分解的初始瞬时频率f
m
；模块m4：根据该子信号y
m
的瞬时频率f
m
计算该子信号的瞬时频率核函数矩阵φ
m
；模块m5：计算该子信号y
m
的解调信号u
m
和v
m
；模块m6：计算该子信号y
m
的初始瞬时频率增量模块m7：计算该子信号y
m
的优化的初始瞬时频率增量δf
m
；模块m8：更新该子信号y
m
的瞬时频率f
m
；模块m9：重构该子信号y
m
；模块m10：重复调用模块m3-模块m9，直到重构的子信号y
m
满足终止迭代的条件；模块m11：采用k阶傅里叶级数拟合该子信号y
m
估计的瞬时频率f
m
；模块m12：更新信号y的参数化时频变换的核函数参数，并对信号y再次进行参数化时频变换；模块m13：重复调用模块m1-模块m12，直到该子信号y
m
的瞬时频率f
m
满足终止迭代的条件；
模块m14：重构该子信号y
m
；模块m15：将该子信号y
m
从当前剩余信号r
m
中移除，然后继续分解下一个子信号y
m+1
，并重复调用模块m1-模块m14，直到原始信号中每个子信号被完全分解。7.根据权利要求6所述的自适应的时频支持的调频信号分解系统，其特征在于，所述模块m1中ptft的核函数参数为傅里叶级数，其计算公式为：其中的计算公式为：式中，k是傅里叶级数的阶数，f0为基频，f0＝1/(2t)，t为信号持续时间，＝1/(2t)，t为信号持续时间，为傅里叶级数的系数；w
h
的计算公式为：所述模块m2中通过求解e(if
m
)的局部最大值得到子信号y
m
的瞬时频率脊度曲线if
m
；e(if
m
)的计算公式为：式中s(t,if
m
(t))表示信号的ptft，λ和β为惩罚系数；所述模块m4中的瞬时频率核函数矩阵φ
m
的计算公式为：φ
m
＝[c
m
,s
m
]式中c
m
表示正弦函数矩阵，s
m
表示余弦函数矩阵；c
m
的计算公式为：c
m
＝diag[cos(θ
m
(t1))
…
cos(θ
m
(t
n
))]式中θ
m
表示子信号y
m
的瞬时相位；瞬时相位θ
m
的计算公式为：式中表示子信号y
m
估计的瞬时频率；s
m
的计算公式为：s
m
＝diag[sin(θ
m
(t1))
…
sin(θ
m
(t
n
))]式中θ
m
表示子信号y
m
的瞬时相位。8.根据权利要求6所述的自适应的时频支持的调频信号分解系统，其特征在于，所述模块m5中子信号y
m
的解调信号u
m
和v
m
的计算公式为：
式中n表示迭代次数，φ
m
表示瞬时频率核函数矩阵，ρ表示惩罚系数，d表示二阶差分矩阵ω的对角矩阵r
m
表示当前的剩余信号；二阶差分矩阵ω的表达式为：ω为n行n列矩阵，n表示信号的离散采样点数；所述模块m6中子信号y
m
的初始瞬时频率增量的计算公式为：式中和表示子信号y
m
在第n次迭代时的解调信号；所述模块m7中子信号y
m
的优化的初始瞬时频率增量δf
m
的计算公式为：式中ω表示二阶差分矩阵，i表示单位矩阵。9.根据权利要求6所述的自适应的时频支持的调频信号分解系统，其特征在于，所述模块m8中子信号y
m
的瞬时频率f
m
的更新公式为：式中n表示迭代次数；所述模块m9中子信号y
m
的重构计算公式为：y
m
＝φ
m
x
m
式中φ
m
表示子信号y
m
的瞬时频率核函数矩阵，x
m
表示子信号y
m
的解调信号；所述模块m10中最开始设置迭代次数n为1，完成一次迭代后，将迭代次数n加1；重构的子信号y
m
的终止迭代条件为：式中ε1表示终止迭代阈值；重构子信号y
m
的终止迭代阈值ε1设置为10-8
。10.根据权利要求6所述的自适应的时频支持的调频信号分解系统，其特征在于，所述模块m11中用于拟合子信号y
m
的瞬时频率f
m
的k阶傅里叶级数的公式为：式中为傅里叶级数的系数；所述模块m13中，子信号y
m
的瞬时频率f
m
的终止迭代条件为：式中h为迭代次数，ε2表示终止迭代阈值；终止迭代阈值ε2设置为10-8
；
所述模块m14中子信号y
m
的重构计算公式为：y
m
＝φ
m
x
m
式中φ
m
表示子信号y
m
的瞬时频率核函数矩阵，x
m
表示子信号y
m
的解调信号；所述模块m15中最开始时设置m＝1，当前剩余信号r1等于原始信号y；在完成对第m个子信号的分解后，移除子信号y
m
，此时当前剩余信号r
m+1
＝r
m-y
m
，并设置m＝m+1。

技术总结
本发明提供一种自适应的时频支持的调频信号分解方法与系统，包括：计算信号的参数化时频变换。从信号的参数化时频变换的时频表示中提取能量最大的子信号的粗糙的瞬时频率脊度曲线。将提取的该子信号的瞬时频率曲线作为信号分解的初始瞬时频率。根据该子信号的瞬时频率计算该子信号的瞬时频率核函数矩阵。计算该子信号的解调信号。计算该子信号的初始瞬时频率增量。计算该子信号的优化的初始瞬时频率增量。更新该子信号的瞬时频率。本发明实现了对具有强噪声的复杂多成分信号的分解，能准确地估计每个子信号的瞬时频率(IF)与重构每个子信号，有效地揭示了原始复杂信号的变化规律，能应用于许多领域的实际工程信号处理，有很强的工程实用价值。很强的工程实用价值。很强的工程实用价值。


技术研发人员：覃程锦 石岗 刘成良 陶建峰
受保护的技术使用者：上海交通大学
技术研发日：2023.07.14
技术公布日：2023/10/6