稀疏双时空非凸罚自适应Chirp模态交叉混叠分解方法

本发明属于信息处理科学与，尤其涉及工业控制系统中的性能评估与故障诊断领域，具体是稀疏双时空非凸罚自适应chirp模态交叉混叠分解方法。

背景技术：

1、目前，在工程实践中，非线性、非平稳、非高斯振动信号承载着机电装备（如大型减速机齿轮箱、核发电机组、汽轮机组）的服役状态与性能退化信息，振动信号处理（如滤波、降噪、重构等）与信号特征提取（如时域特征因子提取，频域特征因子提取等）是机电装备故障诊断、寿命预测以及智能运维的核心任务。

2、然而，由于设备内部结构复杂，导致测量振动信号各个成分相互叠加复杂，多源信号传递路径衰减与延迟，使故障信息来源辨识困难；外边载荷波动与工况随机复杂多变（如变载荷、变速）导致动态测量振动信号被严重污染，增加了故障分离的难度；内部故障激励与外部工况相互耦合，弱化了故障响应信号的动态特性，严重影响了故障振动信号的特征提取、故障分离及分类识别，使得装备故障诊断、寿命预测以及智能运维问题非常棘手。

3、在复杂多源异构故障信号耦合机理上，叠加了外界干扰的振动信号往往与纯净故障振动信号分布差异极大，原始纯净故障振动信号完全淹没在外界干扰分量中，原始纯净故障振动信号频带与干扰分量频带存在部分或全部重叠，导致严重的模态频率交叉混叠问题，而传统的频谱与时频谱信号处理方法难以实现纯净故障振动分量的剥离与提取。因此，如何从模态频率交叉混叠的测量振动信号中剥离出有用的低频故障特征信息是目前机电装备健康智能运维亟待解决的挑战。

4、目前，传统信号分解方法主要包括：奇异值分解方法、经验模态分解法、局部均值分解方法、小波/小波包分解方法、变分模态分解方法、自适应chirp模态分解方法等。

5、然而，上述传统信号分解方法均存在不同程度的算法与应用的局限性，例如，小波/小波包分解方法存在小波基选择问题，小波基与故障响应物理波形匹配问题等；经验模态分解方法与局部均值分解方法存在严重的端点效应、干扰项、模态混叠固有病态问题；变分模态分解方法与变分非线性chirp模态分解方法需要人为设置固有模态分量的个数等；对于多通道信号分解问题，自适应chirp模态分解方法没有考虑多通道信号的时间-空间耦合特性，导致分解降噪效果不明显。

技术实现思路

1、为克服上述瓶颈技术问题，本发明提出了稀疏双时空非凸罚自适应chirp模态交叉混叠分解方法，可自适应提取出模态频率交叉混叠信号分量，实现固有模态分量与外界干扰分量的剥离，解决了频率交叉混叠分解问题。

2、为实现上述目的，本发明采用了以下技术方案：

3、稀疏双时空非凸罚自适应chirp模态交叉混叠分解方法，所述的分解方法包括以下步骤：

4、获取待分解的快变混叠调频信号；

5、计算所述快变混叠调频信号时间序列的分数阶数，构建目标成本函数模型；

6、利用分裂布雷格曼迭代算法求解构建的目标成本函数模型，并依次解调得到快变混叠调频信号中的各个隐含固有模态分量和时频图。

7、进一步的，所述快变混叠调频信号的模型表达式为：

8、     (1)

9、式（1）中， s( t)为快变混叠调频信号，为快变混叠调频信号中的第 i个chirp信号分量， k为chirp信号分量的个数。

10、进一步的，所述计算所述快变混叠调频信号时间序列的分数阶数，构建目标成本函数模型的步骤包括：

11、根据频率调制原理，调整快变混叠调频信号表达为，

12、   (2)

13、式（2）中，第个chirp信号分量的瞬时幅值，瞬时频率，初始相位为，为0均值高斯白噪声；为与的频率函数；

14、幅值 ai( t) 表示为：

15、；

16、幅值 bi( t) 表示为：

17、；

18、从信号观测分量中估计，构建所述目标成本函数模型表示为：

19、 (3)

20、式（3）中，为从观测分量 s( t)去除估计分量 si( t)后的剩余能量， α>0为权重系数， p为分数阶数， λ1与 λ2为正则化参数；

21、双时间域矩阵包括矩阵 dt与矩阵 mt；

22、双空间域矩阵包括矩阵 ds与 ms；

23、矩阵 dt， mt， ds与 ms大小均为 mn× mn，其中 m表示为通道数， n表示为每个通道的采样点数。

24、进一步的，所述利用分裂布雷格曼迭代算法求解构建的目标成本函数模型，并依次解调得到快变混叠调频信号中的各个隐含固有模态分量和时频图的步骤包括；

25、利用分裂布雷格曼迭代算法求解构建的目标成本函数模型，解调得到快变混叠调频信号中所有隐含分量的瞬时频率，得到原始快变混叠调频信号的瞬时时频图。

26、进一步的，所述构建所述目标成本函数模型的具体步骤包括：

27、利用kronecker积方法构造双时间域矩阵 dt与 mt，空间域矩阵 ds与 ms，其中：

28、；

29、；

30、；

31、；

32、其中，符号为kronecker积；

33、将式（3）构建的目标成本函数模型用约束矩阵形式可改写为，

34、(4)

35、式（4）中，对角矩阵， ω为二阶微分算子；

36、；

37、；

38、；

39、  ；

40、；

41、利用重标极差法计算快变混叠调频信号的分数阶 p。

42、进一步的，所述利用重标极差法计算快变混叠调频信号的分数阶 p的具体方法如下：

43、快变混叠调频信号序列的hurst指数值 h可计算为：

44、    (5)

45、式（5）中， r( n)为数据重整化范围， s( n)为标准差，， h指数可通过在对数图上绘制 与 的曲线得到，则分数阶数 p=h-0.5。

46、进一步的，所述利用分裂布雷格曼迭代算法求解构建的目标成本函数模型，解调得到快变混叠调频信号中所有隐含分量的瞬时频率的步骤包括：

47、利用分裂布雷格曼迭代算法求解目标成本函数模型，式(4)可改写为：

48、 (6)

49、式（6）中，参数 μ1， μ2， μ3， μ4， μ5与 μ6为拉格朗日乘子， b1, b2, b3, b4, b5与 b6为布雷格曼变量；

50、其中，布雷格曼变量 b1， b2， b3， b4， b5， b6的更新规则如下：

51、 (7)

52、将式（6）分解为以下7个子问题：

53、  (8a)

54、        (8b)

55、     (8c)

56、     (8d)

57、   (8e)

58、   (8f)

59、   (8g)

60、根据lp-norm最小优化问题，信号 x可利用广义软阈值算法计算得到，

61、；

62、利用广义软阈值算法，子问题 (8a) 的计算过程如下：

63、(9a)

64、(9b)

65、  (9c)

66、其中，子问题 (8a) 的求解结果为：

67、(9d)

68、同理，子问题 (8b), (8c), (8d), (8e), (8f) 与 (8g) 可利用广义软阈值算法得到：

69、      (10a)

70、    (10b)

71、     (10c)

72、    (10d)

73、     (10e)

74、   (10f)

75、chirp信号分量可解调为：

76、；      (11)

77、在解调得到第1个chirp信号分量后，从快变混叠调频信号中去除第1个chirp信号分量，可得到：

78、                (12)

79、式（12）中， r1( t) 为从快变混叠调频信号中去除第1个chirp估计信号后的剩余分量；

80、将剩余分量 r1( t)看作新的快变混叠调频信号，重复上述步骤，依次得到，第2个chirp估计分量，第3个chirp估计分量，… ，第 k个chirp估计分量，直到余量信号 r k( t) 满足预先设定阈值；

81、原始快变混叠调频信号 s( t)可表达为：     (13)

82、式（13）中， r k( t) 为余量信号。

83、进一步的，所述解调得到快变混叠调频信号中所有隐含分量的瞬时频率，得到原始快变混叠调频信号的瞬时时频图的步骤包括：

84、各个chirp估计分量的瞬时频率计算为：；

85、各个chirp估计分量的增量瞬时频率计算为：

86、    (14)

87、最终原始快变混叠调频信号的瞬时频率计算为：

88、      (15)

89、式（15）中，， i为单位矩阵， ω为二阶导数算子， ζ为一常数。

90、与现有技术相比，本发明提供的稀疏双时空非凸罚自适应chirp模态交叉混叠分解方法具有以下优点：

91、1.相比传统信号分解方法，本发明在模态频率交叉混叠分量提取过程中不需要人为选择参数，不需要设置固有模态分量的个数，弥补了变分模态分解方法中先验参数选择的瓶颈问题；

92、2.本发明提出的稀疏双时空非凸罚自适应chirp模态交叉混叠分解方法，考虑原始模态频率交叉混叠信号的时间-空间耦合特性与长相关特性，与现有的技术相比如自适应chirp模态分解方法，可有效去除噪声干扰分量，实现良好的去噪效果；

93、3.本发明提出的分解方法，与现有的技术相比如经验模态分解方法与局部均值分解方法，可有效缓解传统分解方法的断点效应问题，且该方法计算复杂度低，响应速度快，可操作性强，具有良好的工业应用价值。