抑制局部放电信号中周期性窄带干扰的子空间重构方法

段伟润，张宇辉，李天云

（1.国网天津市电力公司，天津 300010；2.东北电力大学 电气工程学院，吉林 吉林 132012）

0 引言

局部放电是造成变压器等电力设备绝缘劣化的重要原因［1］，局部放电信号的检测对保障电力设备正常运行具有重要意义。现场存在的各种干扰信号（如连续性周期窄带干扰、脉冲干扰、白噪声干扰等）增大了局部放电信号检测难度，如何削弱这3类干扰的影响，特别是连续性周期窄带干扰，是提高电力设备绝缘状态监测性能的一个关键问题。

已有学者针对窄带干扰抑制进行了大量研究工作。文献［2-4］采用快速傅里叶变换（FFT）处理局部放电信号中的窄带干扰，但FFT方法受其自身缺陷的影响（如频谱泄漏等）导致消噪效果不够理想。文献［5-8］利用小波、多小波变换良好的时频分析能力，在抑制窄带干扰方面显示出一定的优势，但选取合适的小波类型比较困难。经验模态分解完全根据信号自身进行自适应分解，文献［9-10］将其应用于局部放电信号窄带干扰抑制，取得了较好效果，但模态混叠现象使该算法的稳定性和准确性欠佳。集合经验模态分解［11］削弱了经验模态分解模态混叠现象，但当局部放电信号频率与叠加的窄带干扰频率接近时，分解效果不够理想。文献［12-13］采用混沌理论抑制周期性窄带干扰，失真较小，对白噪声和局部放电信号本身具有免疫性，但该方法需预知窄带干扰频率和预置系统周期策动力，缺乏灵活性。文献［14-15］从估计窄带干扰波形出发，采用傅里叶级数法抑制局部放电窄带干扰，其先通过选择参考段进行FFT频谱分析，再计算傅里叶系数，从而获得窄带干扰估计波形，该方法较好地保留了局部放电信号的波形及幅值，但受FFT的固有缺陷及随机干扰的影响，该类算法仍存在局限性。

本文提出一种基于子空间重构的周期性窄带干扰抑制方法。将奇异值分解［16］与 HANKEL 矩阵［16］相结合，对局部放电信号数据形成的HANKEL矩阵进行奇异值分解，实现窄带干扰子空间和局部放电信号子空间的划分，在窄带干扰子空间中重构干扰信号，通过与原始信号相减获得待检测局部放电信号。

1 HANKEL矩阵与奇异值分解

设局部放电信号的采样序列X为：

其中，s（k）为周期性窄带干扰信号；n（k）为局部放电信号；k=1，2，…，N，N 为采样点数。

由采样数据形成如下HANKEL矩阵：

其中，H（i，j）=x（i+j-1）；p+q-1=N，p≥q，q 一般在N／4～N／3 之间取值［17］。

式（3）可进一步表示为：

其中，Hs（i，j）=s（i+j-1）；Hn（i，j）=n（i+j-1）。

对构造的HANKEL矩阵H进行奇异值分解：

其中，正交矩阵 UϵRp×p、UsϵRp×rs和 UnϵRp×rn的列向量分别由的特征向量组成；正交矩阵 V ϵRq×q、VsϵRq×rs和 VnϵRq×rn的列向量分别由的特征向量组成；ΛϵRp×q、ΛsϵRrs×rs和为对角矩阵，对角元素 ε1≥ε2≥…≥εr≥0、σ1≥σ2≥…≥σrs≥0、γ1≥γ2≥…≥γrn≥0 分别为HANKEL 矩阵 H、Hs和 Hn的奇异值，r、rs和 rn分别为矩阵 H、Hs和 Hn的秩，且 r=rs+rn。

由式（5）知因此，通过获取 Hs（i，j）=s（i+j-1），可重构出周期性窄带干扰信号 s（k）。

式（5）可进一步写成：

其中，εi为H第i个奇异值；ui为HHT的第i个特征向量；vi为 HTH 的第 i个特征向量；s和n分别为Hs和Hn的最佳逼近矩阵。

由上述分析可知，通过奇异值分解并根据周期性窄带干扰和局部放电信号各自的特点，即窄带干扰信号s（k）与局部放电信号 n（k）之间的不相关性，以及窄带干扰信号能量比较集中而局部放电信号能量比较分散的特点，可以将由局部放电测量信号所构成的HANKEL矩阵H分成2个互不相关的子空间，即窄带干扰子空间s和局部放电信号子空间n。从HANKEL矩阵 H 中去除n得到s，进而得到窄带干扰信号 sˆ（k）。

2 窄带干扰子空间和局部放电子空间划分

2个信号子空间的正确划分对窄带干扰的重构精度至关重要。窄带干扰子空间和局部放电子空间的划分即是求H的有效秩l。本文根据奇异值大小判断l值，即满足式（7）的i的最大值记为有效秩。

其中，εi为第i个奇异值；c可根据奇异值变化曲线设置。

3 重构周期性窄带干扰信号

其中

子空间重构抑制周期性窄带干扰的步骤如下：

a.将采集到的局部放电信号数据形成HANKEL矩阵，对HANKEL矩阵进行奇异值分解；

b.确定参数l的值，划分窄带干扰子空间与局部放电信号子空间；

c.根据获得的窄带干扰子空间数据，由式（8）重构窄带干扰波形，通过与原始信号相减获得待检测局部放电信号。

4 仿真分析和实测处理

4.1 窄带干扰仿真分析

仿真中，周期性窄带干扰f（t）由4个不同频率和幅值的正弦波叠加而成，表达式为：

其中，ki、li（i=1，2，3，4）分别为窄带干扰幅值（mV）及频率（kHz）。考虑加入 l1、l2、l3、l4分别为 348 kHz、450 kHz、500kHz、800 kHz的窄带干扰信号，各频率成分幅值 k1、k2、k3、k4是随机的。采样频率 10 MHz，其时域波形如图1所示，计算数据窗长度为0.15 ms。

图1 窄带干扰仿真信号f（t）Fig.1 Simulative narrowband noises f（t）

采集窄带干扰信号数据形成HANKEL矩阵，根据大量仿真验证，文中HANKEL矩阵q取为N／4，p则取为3N／4+1，同时为保证计算的快速性，数据窗长度不宜过长。对HANKEL矩阵进行奇异值分解后的奇异值变化曲线如图2所示，从第9个奇异值开始无明显变化，可将HANKEL矩阵有效秩取为8。

图2 奇异值变化曲线1Fig.2 Curve of singular value（case 1）

采用本文方法重构窄带干扰信号，并将其与图1原始信号相减得到窄带干扰的重构误差，如图3所示。可见，窄带干扰重构误差在10-13数量级上，重构信号与原始信号基本一致，验证了子空间重构窄带干扰的可行性。

局部放电信号可用单指数衰减振荡模型和双指数衰减振荡模型描述：

图3 窄带干扰重构误差Fig.3 Reconstruction error of narrowband noises

其中，τ为衰减系数；fc为振荡频率；A为信号幅值。模拟4组放电脉冲，参数如表1所示。

表1 局部放电信号参数Table 1 Parameters of PD signals

仿真信号如图4所示，其中图4（a）为理想局部放电信号，幅值为0.9439 mV；图4（b）为叠加图1中窄带干扰后的信号，信噪比（SNR）为-15.0352 dB。

图4 仿真信号Fig.4 Simulative signals

采集局部放电信号数据形成HANKEL矩阵，对其进行奇异值分解后的奇异值变化曲线如图5所示，从第9个奇异值开始无明显变化，视为局部放电信号分量引起的奇异值。

图5 奇异值变化曲线2Fig.5 Curve of singular value（case 2）

此时，可将HANKEL矩阵有效秩取为8。采用本文方法抑制周期性窄带干扰后得到的局部放电信号如图6所示，信噪比为17.7336 dB。

图7为图6第3组脉冲信号的局部放大图，可以看出，本文方法能够较好地保留局放脉冲的高频分量、低频分量、幅值和极性，有效地抑制窄带干扰。

作为参考，分别用傅里叶级数法［14］、小波消噪方法［5］和多小波消噪方法［6］对上述数据进行分析。表2为采用傅里叶级数法时窄带干扰频率FFT估计值。

图6 抑制周期性窄带干扰后的局部放电信号Fig.6 PD signals after periodic narrowband noise suppression

图7 局部放电信号提取结果Fig.7 Extracted PD signals

表2 FFT频率估计值Table 2 Estimated frequencies by FFT kHz

图8为采用傅里叶级数法消除窄带干扰后的局部放电波形，信噪比为-3.0918 dB。由表2及图8的消噪结果可知，FFT算法本身的固有缺陷导致窄带干扰波形未能完全抵消，进而影响了结果的准确性。

图8 经傅里叶级数法提取的局部放电信号Fig.8 PD signals extracted by Fourier series

图9上、下波形分别对应小波消噪和多小波消噪处理结果，信噪比分别为0.7006 dB和2.5064 dB。

图9 小波和多小波处理结果Fig.9 PD signals extracted by wavelet and multi-wavelet

比较图6、图8和图9可知：小波和多小波方法将放电信号的部分信息和窄带干扰一起消除，导致放电信号的幅值和波形等发生严重变化；傅里叶级数方法在一定程度上抑制了窄带干扰，保留了放电信号的主要特征，但受FFT算法本身固有缺陷的影响，消噪效果不够理想；本文方法通过求解窄带干扰子空间数据来重构窄带干扰波形，避免了傅里叶级数方法频率估计引起的误差，自适应逼近能力较强。

通过不同的ki和li值组合，验证本文方法抑制窄带干扰信号的有效性，参数如表3所示。

表3 窄带干扰参数Table 3 Parameters of narrowband noises

本文方法与傅里叶级数法、小波和多小波方法的评价参数对比结果如表4所示。

表4 评价参数计算结果Table 4 Comparison of performances among different methods

从第4、5组样本评价参数对比结果中可以看出，相比于小波和多小波方法，本文方法和傅里叶级数法在抑制窄带干扰方面具有一定的优势；从这5组样本评价参数对比结果中可以看出，本文方法在抑制窄带干扰的稳定性和准确性上面，整体上优于傅里叶级数法、小波和多小波方法。

在图4（b）数据的基础上，考虑局部放电信号的衰减系数τ在 100 ns～2.5 μs（步长为 50 ns）范围内动态变化对本文方法抑制窄带干扰效果的影响。

图10 均方误差变化曲线Fig.10 Curve of mean square error

从图10可以看出，在整个衰减系数变化范围内，均方误差在10-4数量级上，本文方法较好地保留了局放波形的振荡特征。

4.2 随机干扰对抗干扰效果的影响

在图4（b）所示数据的基础上增加方差为0.2的随机噪声，利用本文方法进行抗干扰处理。对局部放电信号形成的HANKEL矩阵进行奇异值分解后的奇异值变化曲线如图11所示，从第9个奇异值开始无明显变化，视为局部放电信号和随机噪声分量引起的奇异值。

图11 奇异值变化曲线3Fig.11 Curve of singular value（case 3）

将HANKEL矩阵有效秩取为8，窄带干扰抑制结果如图12所示。可见，窄带干扰得到了很好的抑制，剩余干扰主要为白噪声，信噪比为-3.6617dB。

图12 抑制周期性窄带干扰后的局部放电信号Fig.12 PD signals after periodic narrowband noise suppression

图13给出了采用傅里叶级数法消除窄带干扰后的局部放电波形，信噪比为-6.4432 dB。

图13 经傅里叶级数法提取的局部放电信号Fig.13 PD signals extracted by Fourier series

由于随机干扰并不满足Dirichlet条件，傅里叶级数无法对其重构，因此引入了计算误差，再加上FFT频率估计偏差的影响导致窄带干扰波形不能完全抵消，信噪比较低。

4.3 现场实测数据分析

图14 实测信号Fig.14 Measured PD signals

图14为某变电站现场所采集到的一段局部放电信号数据，由于获取的局部放电信号窄带干扰不太明显，故在检测到的信号中加入式（9）形式的窄带干扰，初始相位随机设置为 60°、30°、-45°、70°，叠加干扰后的信号如图15所示。

图15 加入窄带干扰后的信号Fig.15 PD signals mixed with periodic narrowband noises

将图15的局部放电信号数据形成HANKEL矩阵，对其进行奇异值分解后的奇异值变化曲线如图16所示，HANKEL矩阵有效秩取为8。

图16 奇异值变化曲线4Fig.16 Curve of singular value（case 4）

采用本文方法抑制窄带干扰后的局部放电波形如图17所示。

图17 抑制周期性窄带干扰后的局部放电信号Fig.17 PD signals after periodic narrowband noise suppression

表5为采用傅里叶级数法时，窄带干扰频率FFT估计值，抑制窄带干扰后的局部放电波形见图18。

表5 FFT频率估计值Table 5 Estimated frequencies by FFT kHz

比较图14、图17和图18可知，本文提出的抑制局部放电周期性窄带干扰方法可以有效抑制窄带干扰，很好地保留了局部放电脉冲。从图18中可以看出，傅里叶级数法能在一定程度上抑制窄带干扰，但在估计347 kHz频率分量时，FFT发生了谱峰偏移，再加上随机干扰的影响，导致抗干扰效果不够理想。

图18 经傅里叶级数法提取的局部放电信号Fig.18 PD signals extracted by Fourier series

5 结论

局部放电信号的多态性及频谱具有几乎分布在整个频率区间上的特点，使基于波形匹配或频带划分的窄带干扰抑制方法存在适用局限性。仿真和实测信号处理结果表明：

a.从估计窄带干扰信号波形出发，逆向分析局部放电信号，能够有效抑制窄带干扰，更好地保留局部放电脉冲的高频分量、低频分量、幅值等信息；

b.本文提出的基于子空间重构的窄带干扰抑制方法，从矩阵的角度出发，利用窄带干扰信号能量比较集中、局部放电信号和随机干扰信号能量比较分散的特点，划分窄带干扰信号子空间和局部放电信号子空间（含随机干扰），通过子空间数据重构窄带干扰波形，相比于傅里叶级数法，方法简单、鲁棒性强，较适合实际使用，为局放信号窄带干扰抑制提供了一种新的选择。

参考文献：

［1］苑津莎，尚海昆.基于主成分分析和概率神经网络的变压器局部放电模式识别［J］.电力自动化设备，2013，33（6）：27-31.YUAN Jinsha，SHANG Haikun.Patternrecognitionbased on principal component analysis and probabilistic neural networks for partial discharge of power transformer［J］.Electric Power Automation Equipment，2013，33（6）：27-31.

［2］谢良聘，朱德恒.FFT频域分析算法抑制窄带干扰的研究［J］.高电压技术，2000，26（4）：6-8.XIE Liangpin，ZHU Deheng.Research ofspectrum analysis based on FFT for suppressing narrowband interference in PD signal［J］.High Voltage Engineering，2000，26（4）：6-8.

［3］卢毅，韩志锟，陈瑜.基于能量比预处理的FFT窄带滤波方法［J］.东南大学学报（自然科学版），2010，40（5）：948-951.LU Yi，HAN Zhikun，CHEN Yu.FFT narrow bandfiltering method based on energy ratio pretreatment［J］.Journal of Southeast University（Natural Science Edition），2010，40（5）：948-951.

［4］罗新，牛海清，胡日亮，等.一种改进的用于快速傅里叶变换功率谱中的窄带干扰抑制的方法［J］.中国电机工程学报，2013，33（12）：167-175.LUO Xin，NIU Haiqing，HU Riliang，et al.A modified method of suppressing narrow-band interference using fast Fourier transform power spectrum［J］.Proceedings of the CSEE，2013，33（12）：167-175.

［5］ZHENG Wei，ZHOU Xichao，XU Honglei，et al.Partial discharge signalobtain based on adaptive waveletwith mathematical morphology［J］.Physics Procedia，2012，24（Part B）：912-917.

［6］钱勇，黄成军，陈陈，等.多小波消噪算法在局部放电检测中的应用［J］.中国电机工程学报，2007，27（6）：89-95.QIAN Yong，HUANG Chengjun，CHEN Chen，et al.Application of multi-wavelet based on de-noising algorithm in partial discharge detection［J］.Proceedings of the CSEE，2007，27（6）：89-95.

［7］李化，杨新春，李剑，等.基于小波分解尺度系数能量最大原则的GIS局部放电超高频信号自适应小波去噪［J］.电工技术学报，2012，27（5）：84-91.LI Hua，YANG Xinchun，LI Jian，et al.The maximum energy of wavelet decomposition approximation-related adaptive wavelet denosing for partial discharge UHF pulse in GIS［J］.Transactions of China Electrotechnical Society，2012，27（5）：84-91.

［8］刘卫东，刘尚合.基于小波分解与滑动峰态的微弱放电信号联合检测方法［J］.电力自动化设备，2016，36（11）：152-156，164.LIU Weidong，LIU Shanghe. Weak partialdischarge signal detection based on wavelet decomposition and sliding kurtosis［J］.Electric Power Automation Equipment，2016，36 （11）：152-156，164.

［9］沈宏，张蒲，徐其惠，等.改进自适应噪声对消算法的窄带干扰抑制［J］.仪器仪表学报，2008，29（12）：2632-2636.SHEN Hong，ZHANG Pu，XU Qihui，et al.Narrow-band interference suppression based on improved adaptive noise cancellation algorithm［J］.Chinese Journal of Scientific Instrument，2008，29（12）：2632-2636.

［10］臧怀刚，李清志.改进的EMD方法在局部放电信号提取中的应用［J］.电力系统及其自动化学报，2014，26（11）：78-81.ZANG Huaigang，LIQingzhi.ApplicationofimprovedEMD method on extraction of partial discharge signal［J］.Proceedings of the CSU-EPSA，2014，26（11）：78-81.

［11］姚林朋，郑文栋，钱勇，等.基于集合经验模态分解的局部放电信号的窄带干扰抑制［J］.电力系统保护与控制，2011，39（22）：133-139.YAO Linpeng，ZHENG Wendong，QIAN Yong，et al.A narrowband interferencesuppression method based on EEMD for partial discharge［J］.Power System Protection and Control，2011，39（22）：133-139.

［12］李楠，廖瑞金，孙才新，等.一种用混沌振子去除局部放电信号中窄带干扰的新方法［J］.电工技术学报，2006，21（2）：88-92.LI Nan，LIAO Ruijin，SUN Caixin，et al.A new suppression method for narrow-band interference in partial discharge signals with chaotic oscillator［J］.Transactions of China Electrotechnical Society，2006，21（2）：88-92.

［13］舒娜，张晓星，孙才新，等.采用Van-der混沌振子抑制局部放电信号中周期性窄带干扰［J］.高电压技术，2012，38（1）：89-94.SHU Na，ZHANG Xiaoxing，SUN Caixin，et al.Van-der chaotic osillator to suppressing the periodic narrow-band interference from partial discharge pulse signal［J］.High Voltage Engineering，2012，38（1）：89-94.

［14］程养春，李成榕，王伟.用傅立叶级数法消除局部放电检测中窄带干扰的研究［J］.中国电机工程学报，2005，25（20）：106-111.CHENG Yangchun，LIChengrong，WANG Wei.Researchon eliminating narrow frequency noise by Fourier series in partial discharge detection［J］.Proceedings of the CSEE，2005，25（20）：106-111.

［15］刘卫东，刘尚合，魏明.抑制周期性窄带干扰的时域重构优化算法性能分析［J］.高电压技术，2010，36（5）：1211-1216.LIU Weidong，LIU Shanghe，WEI Ming.Analysis of the optimization reconstruction method in time domain forperiodic narrowband interference eliminating［J］.High Voltage Engineering，2010，36（5）：1211-1216.

［16］张贤达.矩阵分析与应用［M］.北京：清华大学出版社，2004.

［17］张宇辉，陈峰，李慧敏，等.基于小波变换和矩阵束算法的同步电机参数辨识［J］.电力系统保护与控制，2012，40（9）：87-92.ZHANG Yuhui，CHEN Feng，LI Huimin，et al.Parameter identification of synchronous machine based on wavelet transform and matrix pencil algorithm［J］.Power System Protection and Control，2012，40（9）：87-92.