CN111855192B - 一种用于编码器信号去噪的奇异值分解方法 - Google Patents

Abstract

本发明提供一种用于编码器信号去噪的奇异值分解方法，用于机械设备的故障诊断和信号处理，其包括以下步骤：S1、读取输出编码器信号并进行高频采样，得到输出轴的角位置测试信号θ(t_n)；S2、将角位置信号θ(t_n)进行二阶差分得到瞬时角加速度信号a(t_n)；S3、对瞬时角加速度信号a(t_n)进行时域同步平均得到信号x(n)，去除无关的噪声与干扰分量；S4、利用一维信号x(n)构造轨迹矩阵H并进行奇异值分解，再利用对角线平均得到瞬时角加速度信号m个的一维信号分量X_i[n]；S5、设计构建谐波的归一化比例指标NPH，计算每个信号分量X_i[n]对应的谐波归一化比例NPH，选择具有最大NPH的信号成分作为瞬时角加速度信号中微弱故障成分，根据故障成分中冲击的间隔距离得到故障类型分析结果。

Description

一种用于编码器信号去噪的奇异值分解方法

技术领域

本发明涉及机械设备故障诊断和信号处理分析技术领域，特别涉及一种用于编码器信号去噪的奇异值分解方法。

背景技术

齿轮箱作为机械设备传动链系统关键设备，其故障是导致设备停机甚至发生严重事故的主要原因。因此，监测齿轮箱健康状态，诊断齿轮箱早期故障已经成为机械故障诊断领域热点之一。现阶段，振动分析为机械设备故障诊断最为有效的途径之一，源于机械设备的状态劣化往往表现为振动信息的变化或异常，但由于现代机械设备更加复杂、更加集成的发展趋势，振动分析方法在齿轮箱故障诊断中存在一定局限性，例如测试信号传递路径长、多成分耦合等难题。因此，寻找新的诊断信息源成为当前机械故障诊断亟需解决的问题。

随着现代机械设备的智能化发展趋势，编码器已经被广泛安装于机械系统中用于设备运转速度和位置控制，同时，编码器安装位置通常距离齿轮箱系统较近，因此内置编码器信号中包含丰富的设备运行动态信息，当齿轮箱的齿轮发生故障时，故障齿轮啮合刚度的周期性变化在编码器信号中往往表现为周期瞬变特征。相比于振动分析，基于内置编码器信号的分析方法具有故障敏感度高、传递路径短、适用范围广、测试成本低等特点。

然而，编码器原始输出信号不仅仅包含由于早期故障引起的瞬态冲击，同时包含齿轮啮合等周期成分，且易于受到测试环境噪声的干扰，给编码器信号中的微弱故障特征提取带来困难。奇异值分解(SVD,Singular value decomposition)方法是近年来在故障诊断行业兴起的测试信号分析方法，其根据一维测试信号构造轨迹矩阵，并通过矩阵分解将信号分解为一组特征向量，并通过奇异值排序、设定阈值、截断重构等操作对信号进行去噪。然而，传统SVD方法需要人工预先设定阈值，且这一阈值的选择对于信号重构至关重要，同时，针对不同的应用场景，阈值往往不同且无统一的设定规则。因此，需要一种在复杂背景噪声干扰下能够提取出编码器信号中微弱故障成分的方法，并且能够尽可能减少人为设置参数的影响。

发明内容

为了克服现有技术的缺陷，本发明的目的在于提供一种在复杂背景噪声干扰下能够有效提取编码器信号中的微弱故障特征成分的方法，能精确提取设备运转情况下的故障信息，同时很大程度上简化了数据采集程序、降低了测试费用，能够实现机械设备故障特征提取和状态监测的自动化。

本发明提供一种用于编码器信号去噪的奇异值分解方法，用于机械设备的故障诊断和信号处理，其包括以下步骤：

S1：利用编码器数据采集卡读取机械设备传动装置上的输出编码器信号，并对测试得到的信号进行高频采样，得到输出轴的角位置测试信号θ(t_n)；

S2：将角位置信号θ(t_n)进行二阶差分得到瞬时角加速度信号a(t_n)；

S3：对瞬时角加速度信号a(t_n)进行时域同步平均得到信号x(n)；

S4：利用S3中一维信号x(n)构造轨迹矩阵H，并进行奇异值分解，得到m个矩阵H_i,再对矩阵H_i进行对角线平均得到瞬时角加速度信号m个的一维信号分量X_i[n]，具体步骤如下：

S41、根据表达式(3)构造S3中一维信号x(n)的轨迹矩阵H，所述表达式(3)为：

其中，N表示信号长度，m表示分解成分数量；

S42、根据表达式(4)对轨迹矩阵H进行奇异值分解，得到m个矩阵H_i，所述表达式(4)为：

其中，u_i，

分别表示正交矩阵U,V^T的列向量与行向量，D表示由奇异值σ_i构成的对角矩阵；

S43、对矩阵H_i进行对角线平均得到瞬时角加速度信号m个的一维信号分量X_i[n]，所述X_i[n]为：

X_i[n]＝[x_i(1),x_i(2),…,x_i(N)]；

S5：设计构建谐波的归一化比例指标NPH，计算每个信号分量X_i[n]对应的谐波归一化比例NPH，选择具有最大NPH的信号成分作为瞬时角加速度信号中微弱故障成分，根据故障成分中冲击的间隔距离给出机械设备的故障类型分析结果；

谐波的归一化比例指标NPH为：

其中，F是快速傅立叶变换运算，abs表示取绝对值运算,Hilbert表示希尔伯特变换，X_i[n]表示第i个信号成分，

表示第i个信号成分X_i[n]的希尔伯特包络谱，T表示机械设备故障周期，

是奈奎斯特频率，

表示取整操作,K表示包络谱中关注频率个数，j表示虚数单位,在NPH的计算中考虑了k/T±1处的点，以补偿由输入的周期序列的步长引起的误差；并通过表达式(7)找到信号的最佳候选周期；

其中，arg max()表示函数取最大值时对应的自变量取值；

根据表达式(8)选择最大的谐波的归一化比例指标NPH对应的奇异值分量，作为降噪后的信号，并以此作为瞬时角加速度信号中的微弱故障成分，根据故障成分的冲击间隔进行机械设备的故障分析，并得出机械设备的故障类型分析结果。

优选地，所述S2中瞬时角加速度信号a(t_n)为如下表达式(1)：

其中，t_n表示第n个数据采样时间，Δt表示时间间隔。

优选地，所述S3的具体方法为采用表达式(2)对瞬时角加速度信号a(t_n)进行时域同步平均得到x(n)：

其中，L表示平均次数，T是机械设备故障周期。

本发明相比于现有技术，具有以下有益效果：

a)本发明能够有效提取复杂干扰下编码器信号中的微弱故障特征；

b)本发明源于传统SVD方法，采用新指标谐波归一化比例代替传统指标评估编码器信号分解成分的故障信息，能够选择更可靠的分解成分；

c)本发明不需要任何先验周期知识，也无需传统SVD的奇异值排序、设定阈值、截断重构等过程，方法鲁棒性更强；

d)本发明只利用机械设备内置编码器信息，能精确提取设备运转情况下的故障信息，同时很大程度上简化了数据采集程序、降低了测试费用，能够实现故障特征提取和状态监测的自动化。

附图说明

图1为本发明实施例中试验台结构示意图；

图2为本发明方法流程图；

图3为本发明实施例中内置编码器角位置测试信号θ(t_n)示意图；

图4为本发明实施例中二次差分后瞬时角加速度信号a(t_n)示意图；

图5为本发明实施例中经过时域同步平均的瞬时角加速度信号x(n)示意图；

图6为本发明实施例中最优分解成分；

图7为本发明实施例中经最小熵解卷积滤波后信号；

图8为本发明实施例中传统奇异值分解前12个奇异值分量；以及

图9为本发明实施例中本方法的前12个奇异值分量。

图中：

1-磁粉制动器；2-第二联轴器；3-输出编码器；4-第二轴承；5-行星齿轮箱；6-第一轴承；7-输入编码器；8-第一联轴器；9-驱动电机。

具体实施方式

为更好的理解本发明的技术方案，下面结合附图1-9和具体实施例，对本发明的具体实施方式作进一步详细描述。

以某行星齿轮箱故障检测试验台为例进行说明，如图1所示，该试验台至少包括驱动电机9、第一联轴器8、第二联轴器2、磁粉制动器1、第一轴承6、第二轴承4，以及行星齿轮箱5。驱动电机9通过第一轴承6和第一联轴器8与行星齿轮箱5的输入轴相连接，磁粉制动器1通过第二轴承4和第二联轴器2与行星齿轮箱5的输出轴相连接。行星齿轮箱5的输入轴上安装有输入编码器7，行星齿轮箱5的输出轴上安装有输出编码器3。整个行星齿轮箱5试验台由驱动电机9驱动，将扭矩通过行星齿轮箱5传递到磁粉制动器1，磁粉制动器1完成加载过程，通过行星齿轮箱5的输入编码器7和输出编码器3进行信号采集。预制行星轮齿面剥落故障置于行星齿轮箱5中用于测试实验。现有技术中的行星齿轮箱5结构和原理都基本相同。在本实施例中的行星齿轮箱5行星轮故障检测试验台具体参数如下：

驱动电机9额定功率：1.2kW，额定转速：40Hz；

行星齿轮箱5传动比：5.1:1，内齿圈齿数：82，模数：1，行星轮齿数：31，模数：1，太阳轮齿数：20，模数：1；

磁粉制动器1额定功率下的扭矩0.06N*m。

如图2所示，一种用于编码器信号去噪的奇异值分解方法，用于机械设备的故障诊断和信号处理，其包括以下步骤：

S1：利用编码器数据采集卡读取行星齿轮箱5输出轴上的输出编码器3信号，并对测试得到的信号进行高频采样，得到输出轴的角位置测试信号θ(t_n)；如图3所示；

S2：瞬时角加速度信号受机器速度变化和噪声的影响较小，因此它是齿轮特征提取的信号源是一个合适的选择。将角位置信号θ(t_n)利用表达式(1)进行二阶差分得到瞬时角加速度信号a(t_n)，表达式(1)为：

其中，t_n表示第n个数据采样时间，Δt表示采样时间间隔。

S3：根据微分过程原理，经过S2中的二阶差分过程，虽然故障信息得到了增强，但随着采样时间间隔的减小，误差和噪声也同步增大。时域同步平均能够通过使用片段信号的整体平均值消除异步干扰和噪声来增强指定特征。由于齿轮故障周期是固定的并且在转速恒定时容易获得，因此将时域同步平均用于齿轮故障信息的增强。具体为采用表达式(2)对瞬时角加速度信号a(t_n)进行时域同步平均得到x(n)；表达式(2)为：

其中，L表示平均次数，本实施例中L设置为6，T是齿轮故障周期。

S4：经过时域同步平均后，与其他信号成分相比，故障信息仍然较弱。因此，进一步的增强操作对于特征提取至关重要。根据表达式(3)利用S3中一维信号x(n)构造轨迹矩阵H，根据表达式(4)对轨迹矩阵H进行奇异值分解，得到m个矩阵H_i,H_i与轨迹矩阵矩阵H结构相同，进一步对矩阵H_i进行对角线平均得到瞬时角加速度信号m个的一维信号分量X_i[n]，X_i[n]＝[x_i(1),x_i(2),…,x_i(N)]。表达式(3)为：

表达式(4)为：

其中，N表示信号长度，m表示分解成分数量。其中，u_i，

分别表示正交矩阵U,V^T的列向量与行向量，D表示由奇异值σ_i构成的对角矩阵。m是奇异值分解中控制奇异分量数量的重要参数，同时奇异值分解的分解性能与该参数的选择直接相关。通常，较大的m倾向于产生更好的传统奇异值分解性能。但是，另一方面，它也带来了更大的计算成本。因此，需要在计算效率和分解性能之间进行权衡。然而，与以往研究不同的是，本文着重于编码器信号而不是振动信号。在时域同步平均之后，减少了干扰分量。考虑到这些因素，本实施例中选择m为40。

S5：传统奇异值分解方法通过选择较大奇异值对应的奇异分量的方式进行信号的去噪与重构，然而根据奇异值大小的方式选择的奇异分量往往不一定具有最多的故障信息，因此，为了克服传统奇异值分解方法缺点并选择具有最多故障信息的最佳奇异分量，设计了一种新的指标NPH，即谐波的归一化比例NPH，用于识别编码器信号中的齿轮故障信号。通过计算每个信号分量X_i[n]对应的谐波归一化比例NPH，选择具有最大NPH的信号成分作为瞬时角加速度信号中微弱故障成分，根据故障成分中冲击的间隔距离给出齿轮箱的故障分析结果。

谐波的归一化比例NPH定义如下：

其中，F是快速傅立叶变换运算，abs表示取绝对值运算,Hilbert表示希尔伯特变换，X_i[n]表示第i个信号成分，

表示第i个信号成分X_i[n]的希尔伯特包络谱。T表示周期。

是奈奎斯特频率。

表示取整操作,K表示包络谱中关注频率个数,j表示虚数单位。应该强调的是，在NPH的计算中考虑了k/T±1处的点，以补偿由输入的周期序列的步长引起的误差。

并通过表达式(7)找到信号的最佳候选周期，表达式(7)为：

其中，arg max()表示函数取最大值时对应的自变量取值。

根据表达式(8)选择最大的谐波的归一化比例指标NPH对应的奇异值分量，作为降噪后的信号，并以此作为瞬时角加速度信号中的微弱故障成分，根据故障成分的冲击间隔进行机械设备的故障分析，并得出机械设备的故障类型分析结果。

表达式(8)为：

参照图4和图5可知，图4和图5分别是编码器测试信号经二次差分后的瞬时角加速度信号a(t_n)和瞬时角加速度经时域同步平均后的信号x(n)，这两幅图中均无法识别出周期性瞬态冲击故障成分。参照图6，图6显示了本方法应用于信号x(n)的结果，即根据NPH得到的最优信号分量，能够清晰辨别出间隔为31齿的周期冲击，根据行星齿轮箱5结构，只有行星轮的齿数为31，即分析结果证明行星齿轮箱5每间隔31齿会出现可反映在编码器信号中的刚度变化，据此，说明齿轮箱中行星轮发生故障。为凸显本方法在处理内置编码器信号的优势，利用最小熵解卷积办法增强信号x(n)中的瞬态冲击成分，参照图7，图7为基于最小熵解卷积得到的滤波后信号，对比两种方法的处理结果，发现最小熵解卷积专注于单个瞬变并且无法识别齿轮故障，而本发明提出的方法具有明显优势，能够清晰判断出周期故障冲击。参照图8和图9，图8和9分别示出了传统奇异值分解和改进奇异值分解的前12个奇异分量。从图8中可以看出，周期性干扰在信号能量中占主导地位。显然，在这种情况下，基于传统奇异值分解的传统方法无法为机器诊断提供任何帮助。从图9可以看出，最佳奇异分量是使用传统奇异值分解的第38个奇异分量。这意味着与其他组件相比，最佳奇异分量的能量非常低。因此进一步突出了所提出的NPH在所提出的方法中的重要作用。

最后应说明的是：以上所述的各实施例仅用于说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述实施例所记载的技术方案进行修改，或者对其中部分或全部技术特征进行等同替换；而这些修改或替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的范围。

Claims (3)

1.一种用于编码器信号去噪的奇异值分解方法，用于机械设备的故障诊断和信号处理，其特征在于：其包括以下步骤：

S1：利用编码器数据采集卡读取机械设备传动装置上的输出编码器信号，并对测试得到的信号进行高频采样，得到输出轴的角位置测试信号θ(t_n)，其中t_n表示第n个数据采样时间；

S2：将角位置测试信号θ(t_n)利用进行二阶差分得到瞬时角加速度信号a(t_n)；

S3：对瞬时角加速度信号a(t_n)进行时域同步平均得到信号x(n)；

S4：利用S3中一维信号x(n)构造轨迹矩阵H，并进行奇异值分解，得到m个矩阵H_i,再对矩阵H_i进行对角线平均得到瞬时角加速度信号m个的一维信号分量X_i[n]，具体步骤如下：

S41、根据表达式(3)构造S3中一维信号x(n)的轨迹矩阵H，所述表达式(3)为：

其中，N表示信号长度，m表示分解成分数量；

S42、根据表达式(4)对轨迹矩阵H进行奇异值分解，得到m个矩阵H_i，所述表达式(4)为：

其中，u_i，

分别表示正交矩阵U,V^T的列向量与行向量，D表示由奇异值σ_i构成的对角矩阵；

S43、对矩阵H_i进行对角线平均得到瞬时角加速度信号m个的一维信号分量X_i[n]，所述X_i[n]为：

X_i[n]＝[x_i(1),x_i(2),…,x_i(N)]；

S5：设计构建谐波的归一化比例指标NPH，计算每个信号分量X_i[n]对应的谐波归一化比例NPH，选择具有最大NPH的信号成分作为瞬时角加速度信号中微弱故障成分，根据故障成分中冲击的间隔距离给出机械设备的故障类型分析结果；

谐波的归一化比例指标NPH为：

其中，F是快速傅立叶变换运算，abs表示取绝对值运算,Hilbert表示希尔伯特变换，X_i[n]表示第i个信号成分，

表示第i个信号成分X_i[n]的希尔伯特包络谱，T表示机械设备故障周期，

是奈奎斯特频率，

表示取整操作,K表示包络谱中关注频率个数，j表示虚数单位,在NPH的计算中考虑了(k/T)±1处的点，以补偿由输入的周期序列的步长引起的误差；并通过表达式(7)找到信号的最佳候选周期T_op；

其中，arg max()表示函数取最大值时对应的自变量取值；

根据表达式(8)选择最大的谐波的归一化比例指标NPH对应的奇异值分量，作为降噪后的信号，并以此作为瞬时角加速度信号中的微弱故障成分，根据故障成分的冲击间隔进行机械设备的故障分析，并得出机械设备的故障类型分析结果。

2.根据权利要求1所述的用于编码器信号去噪的奇异值分解方法，其特征在于：所述S2中瞬时角加速度信号a(t_n)为如下表达式(1)：

其中，t_n表示第n个数据采样时间，Δt表示时间间隔。

3.根据权利要求1所述的用于编码器信号去噪的奇异值分解方法，其特征在于：所述S3的具体方法为采用表达式(2)对瞬时角加速度信号a(t_n)进行时域同步平均得到x(n)：

其中，L表示平均次数，T是机械设备故障周期。

Images