CN114781430A - 一种局部放电信号去噪方法 - Google Patents

Abstract

本发明涉及一种局部放电信号去噪方法，该方法首先采用量子粒子群寻优法确定VMD分解模态分量个数K和惩罚因子α，其次利用峭度阈值、相关系数阈值法确定保留的IMF分量并进行重构以去除窄带干扰，最后通过自适应阈值法的短时奇异值变换进行白噪声的去除。与现有技术相比，本发明具有在实现局部放电信号窄带干扰及白噪声干扰去除的同时，提高了算法自适应性和计算效率等优点。

Description

一种局部放电信号去噪方法

技术领域

本发明涉及一种局部放电信号去噪方法，尤其是涉及一种基于改进变分模态分解与自适应阈值短时奇异值分解结合的局部放电信号去噪方法。

背景技术

在电力系统中，由生产组装、运输、运行绝缘老化等原因产生的电力设备局部放电是导致电力设备故障甚至损坏的主要原因。在对局部放电的研究中，局部放电波形是判断故障种类、绝缘劣化程度的重要信息。然而，由于测量现场环境和测量仪器的影响，采集得到的局部放电波形会不可避免的受到噪声干扰。去噪方法一般集中在物理去噪方法和数字去噪方法，数字去噪因为其广泛的适用性和灵活性被许多技术人员应用。由于噪声干扰的存在，可能会导致在提取局部放电波形的特征时产生较大的偏差，给故障诊断带来负面影响。由此可见，研究局部放电的窄带干扰及白噪声去噪方法是必要的。

局部放电测量受到的干扰噪声种类可分为窄带干扰和白噪声干扰，窄带干扰分布于某一集中频段内，白噪声则存在于全频段内，由于各种数字去噪算法存在局限性，同一种去噪算法不能很好的兼顾窄带干扰和白噪声的去除。传统的去噪方法主要是基于快速傅里叶变换(IFFT)的去噪方法，通过观察局放信号和周期性窄带干扰在频域上幅值的区别，设定阈值来滤除窄带干扰后快速傅里叶逆变换来还原得到去噪后信号。快速傅里叶变换的缺点在于损失能量较大，且时域到频域的转换中丢失了时域中的信息。为了解决快速傅里叶变换在处理非平稳信号时的缺陷，又有多种去噪方法被提出并应用于局部放电去噪中：短时傅里叶变换(STFT)、小波去噪(WT)、形态学滤波、经验模态分解(EMD)等。

短时傅里叶变换在快速傅里叶变换的基础上增加时域窗，将时域分段进行傅里叶变换，可以很好的保留信号的时频域特征，然而选择时域窗的宽度会很大程度上影响变换后的分辨率。

小波变换以多种小波函数作为小波基，将原始信号看作多种小波和噪声的叠加，以小波系数作为各小波在真实信号中的含量。小波的能量在时域上非常集中，使得其适用于处理非平稳信号，原始信号中含有的窄带干扰在小波尺度增加后会趋于0，而非平稳的局放信号则能在小波基上保留大部分幅值。目前小波去噪的方法有模极大值法、相关系数法、阈值法和它们的改进形式等。然而，合适的小波基和小波系数阈值难以选取，限制了小波去噪的效果。因此，有学者提出了出了自适应阈值法的小波去噪法、小波多尺度变换去噪和小波包多阈值法去噪等方法来进行改进。

相对于小波去噪方法中小波基难以选择的问题，国外的学者提出了经验模态分解法(EMD)。EMD法认为任何信号都可以视为若干内涵模态分量(IMF)的叠加，将原始信号进行多次经验模态分解后舍弃一部分模态分量和残差分量，将内涵模态分量叠加重构便可得到去噪后的真实信号。EMD法的一大缺点是在分解时存在模态混叠和能量泄露问题，有学者提出了EEMD方法和CEEMD方法加以改进。

为了解决EMD方法中的模态混叠和端点效应问题，Dragomiretskiy于2014年提出了变分模态分解(VMD)方法。VMD方法认为原始信号由固有模态分量(BLIMF)组成，通过构造变分问题模型并进行搜索求解来确定每个BLIMF，实现信号从低频到高频的分离。虽然VMD方法弥补了EMD方法的不足，但是产生了选取固有模态分量数K和惩罚因子α的参数选取问题，使得该方法的自适应性较低；VMD对分离出的IMF分量的选择性保留也是一个问题。

另外，VMD对于白噪声的去除体现在最后对残差的舍弃上，这样仍然会导致复原后的波形含有噪声分量。采用短时奇异值分解(STSVD)可以有效的改善这种情况，但对于重构时奇异值的数量确定也是一个问题。

发明内容

本发明的目的就是为了克服上述现有技术存在的缺陷而提供一种基于改进变分模态分解与自适应阈值短时奇异值分解结合的局部放电信号去噪方法，在实现局部放电信号窄带干扰及白噪声干扰去除的同时，提高了算法自适应性和计算效率。

本发明的目的可以通过以下技术方案来实现：

根据本发明的一个方面，提供了一种局部放电信号去噪方法，该方法首先采用量子粒子群寻优法确定VMD分解模态分量个数K和惩罚因子α，其次利用峭度阈值、相关系数阈值法确定保留的IMF分量并进行重构以去除窄带干扰，最后通过自适应阈值法的短时奇异值变换进行白噪声的去除。

作为优选的技术方案，该方法具体以下步骤：

步骤1，获取待检测的局部放电信号f(t)；

步骤2，使用QPSO对VMD分解参数进行寻优，所述分解参数包括惩罚因子α和固有模态分量数K，以平均包络熵值最小化作为寻优函数，寻找合适的K、α值并计算其固有模态分量，其中平均包络熵值为分解得到IMF分量的算术平均值；

根据设定的VMD分解参数和获得的VMD模态分解个数，对含有噪声分量的局部放电信号f(t)进行VMD分解，得到K个模态分量u_k(t)及每个模态分量对应的中心频率ω_k；

步骤3，根据获得的固有模态分量，设定峭度阈值大于正峭度；窄带干扰的周期性强于局放信号分量，故计算各模态分量与它们中心频率所对应的正弦波形两者幅值归一化后的相关系数，相关系数绝对值小于0.5则为应保留的局放分量；对步骤2中得到的固有模态分量进行筛选并重构获得信号f₁(t)；

步骤4，将在步骤3得到的重构信号f₁(t)加信号窗进行截取；对截取的信号片段构建hankel矩阵并进行奇异值分解，计算每窗信号对应hankel矩阵的奇异值平均值，并进行降序排列，形成奇异值平均值序列；计算该序列的平均值，并取1.2倍裕度作为奇异值软阈值，在重构过程中对每窗信号构建hankel矩阵中高于阈值的奇异值对应的重构矩阵进行保留；对得到的重构矩阵进行叠加重构每窗信号对应hankel矩阵，并通过对角线叠加计算的方法还原每窗去噪后的信号片段；将信号序列相同位置上的元素进行算术平均，得到去噪后的局部放电信号。

作为优选的技术方案，所述的步骤1具体为：通过检测阻抗和高频电流互感器的耦合装置采集信号，并采用高压端并联耦合电容接检测阻抗或接地线接高频电流传感器来获得局部放电波形f(t)。

作为优选的技术方案，所述的局部放电波形f(t)为含有窄带干扰和白噪声干扰的电力信号f(t)，其表达式为：

其中等式的左半部分为含噪局部放电信号，右半部分依次为纯净局放信号，窄带干扰信号，白噪声干扰信号，A_i、f_i和Φ_i分别为构成窄带干扰信号的谐波的幅值、频率与相位，a为白噪声幅值。

作为优选的技术方案，所述的步骤2具体为：

步骤2-1：初始化QPSO寻优范围：K-[1,10]；α-[0,5000]、种群规模pop＝30及最大迭代次数maxgen＝30；

步骤2-2：初始化种群个体，在范围内随机初始化粒子的当前位置；以初始化粒子为输入参数进行VMD分解，得到分解后的固有模态分量u'_k(t)；

步骤2-3：以平均包络熵值最小化作为寻优条件进行寻优，寻优函数为

步骤2-4：记录适应度值，与最优适应度值进行比较；若优于全局最优适应度值，则更新最优适应度值f_gbest和全局最优粒子位置g_best，其中算法的迭代过程如下：

步骤2-5：计算收缩-扩张系数a：a＝(1-0.5)*(maxgen-t)/maxgen+0.5，maxgen为寻优算法的最大迭代次数，计算平均粒子历史最优位置M_best：

其中p_{best_i}为当前迭代中第i个粒子当前最优位置；

步骤2-6：粒子位置更新：

其中P_i为第i个粒子位置的更新值，

为(0，1)上均匀分布数值；

其中x_i为第i个粒子的位置，u为(0，1)上均匀分布数值，迭代次数+1；

迭代次数为最大值时查看最优适应度变化曲线，查看是否收敛；若收敛，则确定最终粒子为VMD最优参数K与α；

步骤2-7，初始化

和n，其中

表示第一次循环中分解的第k个模态分量，

表示第一次循环中

对应的中心频率，

表示第一次循环中

对应的拉格朗日乘法算子；

步骤2-8，将利用交替方向乘子法更新模态分量u_k及其对应的中心频率ω_k以及λ：

满足

时停止迭代

其中，

为第n+1次迭代循环中分解出的第k个模态分量、

代表原始信号f(t)的频域谱、

为除第k个模态分量之外的模态分量、

为拉格朗日算子、α为二次惩罚因子、ω_k为模态分量u_k对应的中心频率、

为

的简写、τ为双上升计算更新拉格朗日算子时的噪声容忍度(设置为0)、ω表示上述物理量均在频域范围内计算，ε为迭代停止精度约束条件(设置为1e-7)。

作为优选的技术方案，，所述的步骤3中的峭度的计算公式为：

其中Kurtosis为峭度、x为信号幅值、μ为信号平均幅值、n为信号长度点数。

作为优选的技术方案，所述的步骤3中的相关系数计算公式为：

其中r为相关系数值，f₁(t)、f₂(t)为计算相关系数的两信号。

作为优选的技术方案，所述的步骤4中的信号窗的长度T_w满足略小于波形持续时间。

作为优选的技术方案，所述的步骤4中的对含有白噪声的电力信号f₁(t)进行短时奇异值分解，具体为：

步骤4-1，设定参数值：设定短时奇异值分解参数信号窗长度T；

步骤4-2，设信号f₁(t)长度为N，以T为信号窗长度滑动截取信号片段y_i(t)，其中i＝1，2，……N；

步骤4-3，将y_i(t)变换为hankel矩阵Y_i，过程如下：

其中L的取值设置为T/2，需令L+1大于T-L；

步骤4-4，对矩阵Y_i进行奇异值分解：

Y_i＝U_iS_iV_i ^T

其中

为奇异值序列，K＝T-L，U_i的列向量为S_iS_i ^T对应特征列向量，即U_i为左特征向量矩阵、V_i ^T行向量为

S_i对应特征行向量，即V_i ^T为右特征向量矩阵；

步骤4-5，计算每个hankel矩阵的奇异值平均值

并对其进行降序排列形成奇异值平均值序列S'＝(s'₁,s'₂,…,s'_N)；由于奇异值的性质，信号的大部分能量集中在前几个奇异值上，故对奇异值平均值序列取平均值作为奇异值阈值s_r；

步骤4-6，对S_i进行阈值筛选得符合条件的奇异值矩阵

重构得到去噪后的hankel矩阵

步骤4-7，对去噪后的hankel矩阵Y_i ^r以取反对角线平均值的方法重构出去噪后信号片段

并对不同信号片段中在原信号中位置相同的元素进行叠加取平均值，计算得到去噪后信号s(t)。

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

1)本发明利用量子粒子群算法，结合平均包络熵值最小化作为寻优函数优化变分模态分解参数K与惩罚因子α，提高了本方法的自适应性，也减少了信号的信息泄露；

2)本发明设置峭度、相关度阈值对分解得到的固有模态分量进行筛选，以实现窄带干扰分量的去除；

3)本发明VMD靠去除残差来去除白噪声干扰，故采用短时奇异值变换，并计算平均奇异值序列作为自适应阈值对重构时的奇异值进行选取，提高了自适应性，解决了奇异值选取问题，也改善了去噪后信号在突变处的波形。

4)本发明仿真实例验证了该算法在局部放电信号去噪时的有效性和鲁棒性，对于局部放电去噪及其后处理有重要意义。

5)本发明提高了局部放电信号窄带干扰和白噪声干扰去除的效率和准确性，有利于对去噪结果对波形进行分析，提高局放波形在故障诊断中的应用，增进电力系统故障检测的准确性。

附图说明

图1是本发明局部放电去噪方法的整体流程图；

图2是本发明实施例的仿真含噪波形示意图；

图3是本发明实施例的仿真含噪波形去噪后波形示意图；

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明的一部分实施例，而不是全部实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例，都应属于本发明保护的范围。

本发明基于量子粒子群优化VMD结合自适应阈值选取的短时奇异值分解去噪方法，采用量子粒子群寻优法确定VMD分解模态分量个数K和惩罚因子α，利用峭度阈值、相关系数阈值法确定保留的IMF分量并进行重构以去除窄带干扰，最后通过自适应阈值法的短时奇异值变换进行白噪声的去除，在实现局部放电信号窄带干扰及白噪声干扰去除的同时，提高了算法自适应性和计算效率。

技术术语解释：

QPSO:量子粒子群优化算法，相对于PSO(粒子群算法)取消了粒子的移动方向属性，粒子位置的更新跟该粒子之前的运动没有任何关系，这样就增加了粒子位置的随机性。

VMD：变分模态分解，英文为Variational Mode Decomposition,简称为VMD。

SVD：奇异值分解(Singular Value Decomposition)是线性代数中一种重要的矩阵分解，是矩阵分析中正规矩阵酉对角化的推广。

STSVD：短时奇异值分解，在SVD(奇异值分解)的基础上加窗移动截取信号片段进行去噪，相对于SVD效果去噪效果更好。

如图1所示，一种基于改进变分模态分解与自适应阈值短时奇异值分解结合的局部放电信号去噪方法，包括如下步骤：

步骤1：获取待去噪的局部放电信号；

可选的，通过检测阻抗、高频电流互感器等耦合装置采集信号，可以采用高压端并联耦合电容接检测阻抗或接地线接高频电流传感器来获得局部放电波形f(t)；

步骤2：使用QPSO对VMD分解参数进行寻优，所述分解参数包括惩罚因子α和固有模态分量数K，平均包络熵值(分解得到IMF分量的算术平均值)最小化作为寻优函数，寻找合适的K、α值并计算其固有模态分量。根据设定的VMD分解参数和获得的VMD模态分解个数，对含有噪声分量的局部放电信号f(t)进行VMD分解，得到K个模态分量u_k(t)及每个模态分量对应的中心频率ω_k；每一个模态都是一个调幅调频函数，对应原始信号的某一分量。

步骤3：设定峭度阈值大于3(正峭度)；窄带干扰的周期性强于局放信号分量，故计算各模态分量与它们中心频率(在步骤2VMD分解时求得)所对应的正弦波形两者幅值归一化后的相关系数，相关系数小于0.5则为应保留的局放信号分量。

步骤4：将在步骤3得到的重构信号f₁(t)加信号窗进行截取，信号窗的长度T_w应满足略小于波形持续时间；对截取的信号片段构建hankel矩阵并进行奇异值分解，计算每窗信号对应hankel矩阵的奇异值平均值，并进行降序排列，形成奇异值平均值序列；计算该序列的平均值，并取1.2倍裕度作为奇异值软阈值，在重构过程中对每窗信号构建hankel矩阵中高于阈值的奇异值对应的重构矩阵进行保留。对得到的重构矩阵进行叠加重构每窗信号对应hankel矩阵，并通过对角线叠加计算的方法还原每窗去噪后的信号片段；将信号序列相同位置上的元素进行算术平均，得到去噪后的局部放电信号。

步骤1中，含有窄带干扰和白噪声干扰的电力信号f(t)的表达式可以为：

其中等式的左半部分为含噪局部放电信号，右半部分依次为纯净局放信号，窄带干扰信号，白噪声干扰信号，A_i、f_i和Φ_i为构成窄带干扰信号的谐波的幅值、频率与相位，a为白噪声幅值。

步骤2中，量子粒子群优化(Quantum Particle Swarm Optimization)，简称为QPSO；QPSO相较于PSO算法，只有一个参数影响寻优效果，即收缩-扩张因子a；并且a可以根据寻优的迭代数进行计算。

变分模态分解：Variational Mode Decomposition,简称为VMD；VMD参数惩罚因子α越大，分解得到的各分量信号带宽越小，反之，带宽越大，α的设置可能导致噪声去除效果弱或过度去噪损失信号能量；

可选的，所述步骤2中，模态分解个数K与惩罚因子α根据以平均包络熵值最小化为寻优条件的QPSO确定，可以如图2所示，步骤可以包括：

步骤2-1：初始化QPSO寻优范围：K-[1,10]；α-[0,5000]、种群规模pop＝30及最大迭代次数maxgen＝30；

步骤2-2：初始化种群个体，在范围内随机初始化粒子的当前位置(整数)；以初始化粒子为输入参数进行VMD分解，得到分解后的固有模态分量u'_k(t)；

步骤2-3：以平均包络熵值最小化作为寻优条件进行寻优，寻优函数为

步骤2-4：记录适应度值，与最优适应度值进行比较；若优于全局最优适应度值，则更新最优适应度值f_gbest和全局最优粒子位置gbset。算法的迭代过程如下：

步骤2-5：计算收缩-扩张系数a：a＝(1-0.5)*(maxgen-t)/maxgen+0.5，maxgen为寻优算法的最大迭代次数，计算平均粒子历史最优位置mbest：

其中p_{best_i}为当前迭代中第i个粒子当前最优位置pbest。

步骤2-6：粒子位置更新：

其中Pi用于第i个粒子位置的更新，

为(0，1)上均匀分布数值；

其中x_i为第i个粒子的位置，u为(0，1)上均匀分布数值。迭代次数+1。

迭代次数为最大值时查看最优适应度变化曲线，查看是否收敛；若收敛，则确定最终粒子为VMD最优参数K与α。

步骤2-7，初始化

和n，其中

表示第一次循环中分解的第k个模态分量，

表示第一次循环中

对应的中心频率，

表示第一次循环中

对应的拉格朗日乘法算子；

步骤2-8，将利用交替方向乘子法更新模态分量u_k及其对应的中心频率ω_k以及λ：

满足

时停止迭代

其中，

为第n+1次迭代循环中分解出的第k个模态分量、

代表原始信号f(t)的频域谱、

为除第k个模态分量之外的模态分量、

为拉格朗日算子、α为二次惩罚因子、ω_k为模态分量u_k对应的中心频率、

为

的简写、τ为双上升计算更新拉格朗日算子时的噪声容忍度(设置为0)、ω表示上述物理量均在频域范围内计算，ε为迭代停止精度约束条件(设置为1e-7)。

步骤3中峭度的计算公式为

其中Kurtosis为峭度、x为信号幅值、μ为信号平均幅值、n为信号长度点数；相关系数计算公式为

其中r为相关系数值、f₁(t)、f₂(t)为计算相关系数的两信号。对步骤2中得到的固有模态分量进行筛选并重构获得信号f₁(t)。

所述步骤4中对含有白噪声的电力信号f₁(t)进行短时奇异值分解，具体步骤为：

步骤4-1，设定参数值：设定短时奇异值分解参数信号窗长度T，其中T略小于波形长度；

步骤4-2，设信号f₁(t)长度为N，以T为信号窗长度滑动截取信号片段y_i(t)，其中i＝1，2，……N；

步骤4-3，将y_i(t)变换为hankel矩阵Y_i，过程如下：

其中L的取值设置为T/2左右，需令L+1大于T-L。

步骤4-4，对矩阵Y_i进行奇异值分解：

Y_i＝U_iS_iV_i ^T

其中

为奇异值序列，K＝T-L，U_i的列向量为S_i

对应特征列向量，即U_i为左特征向量矩阵、V_i ^T行向量为

S_i对应特征行向量，即V_i ^T为右特征向量矩阵；

步骤4-5，计算每个hankel矩阵的奇异值平均值

并对其进行降序排列形成奇异值平均值序列S'＝(s'₁,s'₂,…,s'_N)；由于奇异值的性质，信号的大部分能量集中在前几个奇异值上，故对奇异值平均值序列取平均值作为奇异值阈值s_r。

步骤4-6，对S_i进行阈值筛选得符合条件的奇异值矩阵

重构得到去噪后的hankel矩阵

步骤4-7，对去噪后的hankel矩阵Y_i ^r以取反对角线平均值的方法重构出去噪后信号片段

并对不同信号片段中在原信号中位置相同的元素进行叠加取平均值，计算得到去噪后信号s(t)。

本发明在实现窄带干扰和白噪声干扰去除的同时，提高了算法自适应性。本发明提高了局部放电信号去噪的有效性，根据去噪结果对局部放电进行更精确的特征提取，更有利于电力设备的故障诊断方面的准确度。

以上所述，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到各种等效的修改或替换，这些修改或替换都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应以权利要求的保护范围为准。

Claims (9)

1.一种局部放电信号去噪方法，其特征在于，该方法首先采用量子粒子群寻优法确定VMD分解模态分量个数K和惩罚因子α，其次利用峭度阈值、相关系数阈值法确定保留的IMF分量并进行重构以去除窄带干扰，最后通过自适应阈值法的短时奇异值变换进行白噪声的去除。

2.根据权利要求1所述的一种局部放电信号去噪方法，其特征在于，该方法具体以下步骤：

步骤1，获取待检测的局部放电信号f(t)；

步骤2，使用QPSO对VMD分解参数进行寻优，所述分解参数包括惩罚因子α和固有模态分量数K，以平均包络熵值最小化作为寻优函数，寻找合适的K、α值并计算其固有模态分量，其中平均包络熵值为分解得到IMF分量的算术平均值；

根据设定的VMD分解参数和获得的VMD模态分解个数，对含有噪声分量的局部放电信号f(t)进行VMD分解，得到K个模态分量u_k(t)及每个模态分量对应的中心频率ω_k；

步骤3，根据获得的固有模态分量，设定峭度阈值大于正峭度；窄带干扰的周期性强于局放信号分量，故计算各模态分量与它们中心频率所对应的正弦波形两者幅值归一化后的相关系数，相关系数绝对值小于0.5则为应保留的局放分量；对步骤2中得到的固有模态分量进行筛选并重构获得信号f₁(t)；

步骤4，将在步骤3得到的重构信号f₁(t)加信号窗进行截取；对截取的信号片段构建hankel矩阵并进行奇异值分解，计算每窗信号对应hankel矩阵的奇异值平均值，并进行降序排列，形成奇异值平均值序列；计算该序列的平均值，并取1.2倍裕度作为奇异值软阈值，在重构过程中对每窗信号构建hankel矩阵中高于阈值的奇异值对应的重构矩阵进行保留；对得到的重构矩阵进行叠加重构每窗信号对应hankel矩阵，并通过对角线叠加计算的方法还原每窗去噪后的信号片段；将信号序列相同位置上的元素进行算术平均，得到去噪后的局部放电信号。

3.根据权利要求2所述的一种局部放电信号去噪方法，其特征在于，所述的步骤1具体为：通过检测阻抗和高频电流互感器的耦合装置采集信号，并采用高压端并联耦合电容接检测阻抗或接地线接高频电流传感器来获得局部放电波形f(t)。

4.根据权利要求3所述的一种局部放电信号去噪方法，其特征在于，所述的局部放电波形f(t)为含有窄带干扰和白噪声干扰的电力信号f(t)，其表达式为：

其中等式的左半部分为含噪局部放电信号，右半部分依次为纯净局放信号，窄带干扰信号，白噪声干扰信号，A_i、f_i和Φ_i分别为构成窄带干扰信号的谐波的幅值、频率与相位，a为白噪声幅值。

5.根据权利要求2所述的一种局部放电信号去噪方法，其特征在于，所述的步骤2具体为：

步骤2-1：初始化QPSO寻优范围：K-[1,10]；α-[0,5000]、种群规模pop＝30及最大迭代次数maxgen＝30；

步骤2-2：初始化种群个体，在范围内随机初始化粒子的当前位置；以初始化粒子为输入参数进行VMD分解，得到分解后的固有模态分量u'_k(t)；

步骤2-3：以平均包络熵值最小化作为寻优条件进行寻优，寻优函数为

步骤2-4：记录适应度值，与最优适应度值进行比较；若优于全局最优适应度值，则更新最优适应度值f_gbest和全局最优粒子位置g_best，其中算法的迭代过程如下：

步骤2-5：计算收缩-扩张系数a：a＝(1-0.5)*(maxgen-t)/maxgen+0.5，maxgen为寻优算法的最大迭代次数，计算平均粒子历史最优位置M_best：

其中p_{best_i}为当前迭代中第i个粒子当前最优位置；

步骤2-6：粒子位置更新：

其中P_i为第i个粒子位置的更新值，

为(0，1)上均匀分布数值；

其中x_i为第i个粒子的位置，u为(0，1)上均匀分布数值，迭代次数+1；

迭代次数为最大值时查看最优适应度变化曲线，查看是否收敛；若收敛，则确定最终粒子为VMD最优参数K与α；

步骤2-7，初始化

和n，其中

表示第一次循环中分解的第k个模态分量，

表示第一次循环中

对应的中心频率，

表示第一次循环中

对应的拉格朗日乘法算子；

步骤2-8，将利用交替方向乘子法更新模态分量u_k及其对应的中心频率ω_k以及λ：

满足

时停止迭代

其中，

为第n+1次迭代循环中分解出的第k个模态分量、

代表原始信号f(t)的频域谱、

为除第k个模态分量之外的模态分量、

为拉格朗日算子、α为二次惩罚因子、ω_k为模态分量u_k对应的中心频率、

为

的简写、τ为双上升计算更新拉格朗日算子时的噪声容忍度、ω表示上述物理量均在频域范围内计算，ε为迭代停止精度约束条件。

6.根据权利要求2所述的一种局部放电信号去噪方法，其特征在于，所述的步骤3中的峭度的计算公式为：

其中Kurtosis为峭度、x为信号幅值、μ为信号平均幅值、n为信号长度点数。

7.根据权利要求2所述的一种局部放电信号去噪方法，其特征在于，所述的步骤3中的相关系数计算公式为：

其中r为相关系数值，f₁(t)、f₂(t)为计算相关系数的两信号。

8.根据权利要求2所述的一种局部放电信号去噪方法，其特征在于，所述的步骤4中的信号窗的长度T_w满足略小于波形持续时间。

9.根据权利要求2所述的一种局部放电信号去噪方法，其特征在于，所述的步骤4中的对含有白噪声的电力信号f₁(t)进行短时奇异值分解，具体为：

步骤4-1，设定参数值：设定短时奇异值分解参数信号窗长度T；

步骤4-2，设信号f₁(t)长度为N，以T为信号窗长度滑动截取信号片段y_i(t)，其中i＝1，2，……N；

步骤4-3，将y_i(t)变换为hankel矩阵Y_i，过程如下：

其中L的取值设置为T/2，需令L+1大于T-L；

步骤4-4，对矩阵Y_i进行奇异值分解：

Y_i＝U_iS_iV_i ^T

其中

为奇异值序列，K＝T-L，U_i的列向量为

对应特征列向量，即U_i为左特征向量矩阵、V_i ^T行向量为

对应特征行向量，即V_i ^T为右特征向量矩阵；

步骤4-5，计算每个hankel矩阵的奇异值平均值

并对其进行降序排列形成奇异值平均值序列S'＝(s′₁,s′₂,…,s′_N)；由于奇异值的性质，信号的大部分能量集中在前几个奇异值上，故对奇异值平均值序列取平均值作为奇异值阈值s_r；

步骤4-6，对S_i进行阈值筛选得符合条件的奇异值矩阵

重构得到去噪后的hankel矩阵

步骤4-7，对去噪后的hankel矩阵Y_i ^r以取反对角线平均值的方法重构出去噪后信号片段

并对不同信号片段中在原信号中位置相同的元素进行叠加取平均值，计算得到去噪后信号s(t)。

Images