CN107832777A - 一种采用时域压缩多分辨率快速s变换特征提取的电能质量扰动识别方法 - Google Patents

一种采用时域压缩多分辨率快速s变换特征提取的电能质量扰动识别方法 Download PDF

Info

Publication number
CN107832777A
CN107832777A CN201710948335.XA CN201710948335A CN107832777A CN 107832777 A CN107832777 A CN 107832777A CN 201710948335 A CN201710948335 A CN 201710948335A CN 107832777 A CN107832777 A CN 107832777A
Authority
CN
China
Prior art keywords
mrow
frequency
msub
point
mtd
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201710948335.XA
Other languages
English (en)
Other versions
CN107832777B (zh
Inventor
林琳
高兴泉
韩光信
孙明革
陈玲玲
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Dragon Totem Technology Hefei Co ltd
Luoyi Technology Shanghai Co ltd
Original Assignee
Jilin Institute of Chemical Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Jilin Institute of Chemical Technology filed Critical Jilin Institute of Chemical Technology
Priority to CN201710948335.XA priority Critical patent/CN107832777B/zh
Publication of CN107832777A publication Critical patent/CN107832777A/zh
Application granted granted Critical
Publication of CN107832777B publication Critical patent/CN107832777B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/24—Classification techniques
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
    • G06F18/213—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
    • G06F18/214—Generating training patterns; Bootstrap methods, e.g. bagging or boosting
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
    • G06Q50/00—Information and communication technology [ICT] specially adapted for implementation of business processes of specific business sectors, e.g. utilities or tourism
    • G06Q50/06—Energy or water supply
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
    • G06F2218/08—Feature extraction
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
    • G06F2218/12—Classification; Matching

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Evolutionary Computation (AREA)
  • General Engineering & Computer Science (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Business, Economics & Management (AREA)
  • Evolutionary Biology (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Artificial Intelligence (AREA)
  • Economics (AREA)
  • Health & Medical Sciences (AREA)
  • Human Resources & Organizations (AREA)
  • Strategic Management (AREA)
  • Tourism & Hospitality (AREA)
  • General Business, Economics & Management (AREA)
  • Primary Health Care (AREA)
  • Geometry (AREA)
  • Marketing (AREA)
  • General Health & Medical Sciences (AREA)
  • Water Supply & Treatment (AREA)
  • Computer Hardware Design (AREA)
  • Public Health (AREA)
  • Other Investigation Or Analysis Of Materials By Electrical Means (AREA)
  • Complex Calculations (AREA)

Abstract

本发明提供一种采用时域压缩多分辨率快速S变换特征提取的电能质量扰动识别方法，应用于电能质量扰动信号分析与研究，其特点是：包括电能质量扰动信号仿真，对扰动信号进行多分辨率快速S变换处理，在确定特征的前提下，对快速S变换过程中每个主要频率点的加窗傅里叶变换逆变换结果计算所需保留信息，构造中间矩阵，对中间矩阵提取有效特征构建特征向量，用于分类器构建。具有科学合理，适用性强，效果佳，能够完成对复杂电能质量扰动信号进行识别等优点。

Description

一种采用时域压缩多分辨率快速S变换特征提取的电能质量 扰动识别方法
技术领域
本发明是一种采用时域压缩多分辨率快速S变换(Time-domain CompressionMultiresolution Fast S-transform,TCMFST)特征提取的电能质量扰动识别方法，应用于电能质量扰动信号分析与研究。
背景技术
随着可持续发展及低碳理念的提出，新能源发电规模不断扩大，新能源发电的并网以及大量电力电子设备的应用将影响电网的电能质量(Power Quality,PQ)，在科学技术及工业生产快速发展的同时，对供电质量的关注日益提升，因此，需要对电能质量开展深入的监控与分析，而大量监控点收集的海量电能质量数据，对电能质量信号分类系统提出了更高的要求。目前，电能质量扰动信号已从传统的单一扰动过渡为复合扰动，扰动种类更多，且实测信号持续时间长，实时监测需要的数据存储空间大。为了提高电能质量分析效果，同时改善扰动识别过程中嵌入式系统内存的限制，应致力于算法的研究与开发，以降低硬件内存对扰动识别的影响。
电能质量扰动信号分类一般是在电能质量扰动信号的时-频分析(Time-Frequency Analysis，TFA)结果中提取特征，之后将特征向量提供给分类器进行电能质量扰动类型识别。
现有研究对电能质量信号的时-频分析进行了大量的深入分析，如Hilbert-Huang变换 (Hilbert-Huang transform，HHT)、小波变换(Wavelet Transform，WT)和S变换(S-transform， ST)等。其中，S变换被证明具有良好的抗噪性，采用S变换及其改进形式的电能质量信号特征提取受到了广泛的重视。快速S变换(Fast S-transform，FST)可以实现频域上的压缩，降低运算复杂度，相较于S变换，运算量有明显降低。但是在实际应用中，当信号采样率过高，信号过长时，经快速S变换得到的模时频矩阵规模仍然较大，空间复杂度较高，需要数据存储空间大，对硬件要求高。分类器方面，传统的支持向量机、BP神经网络以及决策树等具有较强的鲁棒性，但是不能满足高效、准确的分类需求，而极限学习机作为单隐含层神经网络的优化算法克服了这个缺点。但是，由于极限学习机的输入权值和隐藏层节点是随机选择的，必然存在一系列非最优或者不必要的取值，一方面会增加极限学习机算法对测试数据的反应时间，另一方面也会增加网络结构的复杂度。
发明内容
本发明的目的在于，克服现有技术的不足，提供一种科学合理，适用性强，效果佳，能够完成对复杂电能质量扰动信号进行识别的采用时域压缩多分辨率快速S变换特征提取的电能质量扰动识别方法。
一种采用时域压缩多分辨率快速S变换特征提取的电能质量扰动识别方法，其特征是：它包括电能质量扰动信号仿真，对扰动信号进行多分辨率快速S变换处理，在确定特征的前提下，对快速S变换过程中每个主要频率点的加窗傅里叶变换逆变换结果计算所需保留信息，构造中间矩阵，对中间矩阵提取有效特征构建特征向量，用于分类器构建，具体包括的步骤有：
1)电能质量扰动信号仿真：
使用Matlab按照数学模型仿真生成信噪比为20dB～50dB间随机值的12种电能质量扰动信号，包括暂降、暂升、中断、尖峰、闪变、谐波、暂态振荡、切痕、谐波含暂降、谐波含暂升、谐波含闪变以及谐波含振荡，信号采样率为3.2kHz，采样时间为1秒，基频为50Hz；
2)将电能质量扰动信号利用多分辨率快速S变换方法进行处理：
在多分辨率快速S变换的基础上，首先，对仿真扰动信号进行快速傅里叶变换，根据傅立叶谱在频域范围内确定主要频率点；之后，只针对主要频率点对应的部分开展傅里叶变换逆变换，得到时频矩阵；
设t为时间，f为频率，τ为时移参数，当f→nx/NT,τ→jT时，多分辨率快速S变换的离散表达式为：
式中，N为信号长度(本发明中N＝3200)，T为采样时间间隔，j和m均为0至N-1范围内的整数，λx为窗宽调整因子，nx为根据傅里叶谱值确定的扰动信号主要频率点，x为0 至M-1范围内的整数，M为主要频率点的个数。由原始信号经快速傅立叶变换后位移获得。各频段频率范围划分为：频率点在1Hz至100Hz范围内为低频频段，在101Hz至 700Hz范围内为中频频段，在701Hz至(N/2)Hz范围内为高频频段。在此基础上，将λx在不同频率范围内分别定义。低频频段内λx＝2；中频频段内高频频段内，通过基频 (50Hz)傅里叶谱分析判断信号中是否含有扰动基频扰动，若无扰动，取否则，取
3)确定要提取的特征种类：
经过对各特征分类效果分析，确定对扰动信号提取4类特征用于构建特征向量，包括： F1、原始信号1/4周期能量跌落幅度；F2、基频幅值标准差；F3、700Hz以上高频频段能量； F4、100Hz以上各频率对应最大幅值的最大值，
各特征的具体计算方法分别为：
(a)原始信号1/4周期能量跌落幅度：
式中，h(t)为原始信号，t为信号采样点，R(z)为原始信号1/4周期的均方根值(Root Mean Square，RMS)，z＝1,2……200，R0为1/4周期无噪声标准电能质量信号的均方根值，取值为0.7071；
(b)基频50Hz幅值标准差：
式中，t为信号采样点，f0为信号基频频率，50Hz，N为信号长度，N＝3200，S(t,f0)为信号经快速傅里叶变换后，对基频50Hz频率点进行加窗傅里叶变换逆变换后得到的行向量；
(c)700Hz以上高频频段保留能量：
式中，S(t,f)为模时频矩阵，t为采样点，f为频率，N为信号长度，N＝3200；
(d)100Hz以上各频率对应最大幅值的最大值：
F4＝max{max(S(t,f))} 101Hz≤f≤(N/2)Hz (6)
式中，S(t,f)为模时频矩阵，t为采样点，f为频率，N为信号长度，N＝3200；
对以上4类特征分别以三种方式进行提取，
方式一：以式(2)、式(3)为依据，对原始信号进行数学运算后提取特征F1；
方式二：以式(4)为依据，对原始信号进行快速傅立叶变换后，对基频(50Hz)频率点对应部分进行加窗傅里叶变换逆变换后得到的行相量直接进行数学运算，获得特征F2，避免了对模时频矩阵整体运算和储存的复杂过程；
方式三：在对信号进行快速S变换处理的过程中，频率点在700Hz以上的，以式(5)为依据，对此频率范围内主要频率点在加窗傅里叶变换逆变换得到行相量同时提取各行向量所需保留的特征计算中间值，逐个保留频率对应行计算完毕后获得计算特征F3所需中间特征相量；频率点在100Hz以上的，以式(6)为依据，对此频率范围内主要频率点在加窗傅里叶变换逆变换得到行相量同时提取各行向量所需保留的特征计算中间值，逐个保留频率对应行计算完毕后获得计算特征F4所需中间特征相量，随即对以上各中间特征相量进行储存，由此形成的特征矩阵定义为中间矩阵，最后，对中间矩阵进行数学运算，分别提取特征F3和 F4；
4)计及时域压缩与特征提取需要的中间向量与中间矩阵构造过程：
为进一步减小时频矩阵的时域维数，降低矩阵储存空间，采用基于时域压缩多分辨率快速S变换构造中间矩阵；在快速S变换处理信号的过程中，对各主要频率点在加窗傅里叶变换逆变换得到行相量同时提取各行向量所需保留的特征计算中间值，用于构建中间矩阵；经压缩后的中间矩阵时域维数减少，所需的存储空间减小，降低了相关分析的空间存储要求与存储成本；
设原始信号为h(t)，t为信号采样点，构造中间矩阵的运算流程为：
①对原始信号h(t)进行快速傅里叶变换，得到快速傅里叶变换频谱H(m)，即 H(m)＝FFT(h(t))，其中m为0至N-1范围内的整数，N为信号长度，N＝3200；
②对傅里叶谱取模平方得|H(m)|2，设|H(m)|2共有L个极大点，所有极大值序列记为 |H[l]|2，l＝0,1,2,...,L-1，若|H(m)|2中相邻的极大点之间分别有X1,...,Xl,...,XL-1个频率点(Xl为|H(m)|2中第l个极大点与第l+1个极大值之间的频率点数目)，则|H(m)|2的极大值包络计算方式如下：
式中，x＝0,1,...,Xl-1；l＝0,1,...,L-1，为|H(m)|2中第l个极大点对应的数值，由此得到|H(m)|2的极大值包络谱之后，利用动态测度法求中每一个极大点的动态测度。动态测度的求解方法如下：
设M为的一个极大点，如果M两侧各有一个或多个比点M更高的极大点时，点M两侧也一定存在两个与点M同高度的点N1、N2，由点M通向点N1的路径记为P1，M通向 N2的路径记为P2，则极大点M的动态测度等于路径P1的动态测度与路径P2的动态测度中的较小者。设M的动态测度为Dyn(M)，则有：
式中，halt(M)和halt(N)分别表示点M、N的高度，P(M,N)为点M、N之间的路径，表示这两点之间的部分；Dyn(P(M,N))为路径P(M,N)的动态测度，表示此路径上最高点和最低点的高度差，inf表示下确界，表示一个集合中的最小值，通过以上计算，得到极大|H(m)|2值包络的动态测度，通过设置阈值，将大于该阈值的动态测度所对应的频率点定为主要频率点nx，其中，x为0至M-1范围内的整数，M为主要频率点的个数；
③根据公式确定主要频率点nx对应的高斯窗函数W(m,nx)；
④将H(m)移位到H(m+nx)，对H(m+nx)进行加窗处理，得到H(m+nx)·W(m,nx)；
⑤对加窗的H(m+nx)进行快速傅里叶变换逆变换，得到频率点nx对应的加窗傅立叶变换逆变换的结果S(j,nx)，即S(j,nx)＝IFFT(H(m+nx)·W(m,nx))，其中j为采样点，随即对该行提取特征并存储，形成中间矩阵，中间矩阵构造过程具体步骤如下：
根据700Hz以上高频频段能量特征F3和100Hz以上各频率对应最大幅值的最大值F4的原理，在700Hz以上频率范围内，对主要频率点对应加窗傅立叶变换逆变换的结果按照式(9) 进行数学运算，提取特征相量；在100Hz以上频率范围内，对主要频率点对应加窗傅立叶变换逆变换的结果按照式(10)进行数学运算，提取特征相量，对相量进行储存，构建中间矩阵；
设S(t,f)为信号经快速傅里叶变换后，对频率f进行加窗傅里叶变换逆变换后得到的行向量，t为采样点，N为扰动信号采样点数目，N＝3200，中间矩阵构建规则为：
xmax＝max{S(t,f)} 101Hz≤f≤(N/2)Hz (10)
式中，xmax为频率f下对应幅值中的最大值，E为频率f下的能量，由此分别求得中间矩阵中的特征向量m1和m2，m1＝[m11,m12,…,m1d]T为100Hz以上特定频率各采样点对应幅值中的最大值，m2＝[m21,m22,…,m2d]T为700Hz以上特定频率对应的能量，d为此频段内主要频率点数目，由以上特征相量分别构建各频段下的中间矩阵，中间矩阵结构为：
中频频段，101Hz≤f≤700Hz的中间矩阵结构如式(11)所示，式中：Vm为中频频段的中间矩阵，L2为中频频段主要频率点数目，
高频频段(701Hz≤f≤(N/2)Hz，N＝3200)的中间矩阵结构如式(12)所示，式中：Vh为高频频频段的中间矩阵，L3为高频频段主要频率点数目，
经过以上计算，最终得到规模为L2×1的中频中间矩阵Vm，和规模为L3×2的高频中间矩阵Vh，经时域压缩快速S变换后所得中频频段和高频频段中间矩阵列数分别为1列和2列，与未时域压缩的3200列时频矩阵相比，所需储存空间明显减小；
5)基于时域压缩快速S变换中间矩阵的特征向量构建：
在步骤4)中获得各频段中间矩阵的基础上，对中间矩阵按照式(13)、式(14)进行数学运算，最终提取特征F3和F4：
式中，Vh为高频频段中间矩阵，i为中间矩阵中的各行，k为中间矩阵中能量特征对应的列(k＝2)，L3为高频频段中主要频率点的数目，
F4＝max[max(Vm(i,k1)),max(Vh(j,k2))] (14)
式中，Vm和Vh分别为中频频段和高频频段的中间矩阵，i和j分别为两个矩阵中的行， i∈[1,L2]，J∈[1,L3]，L2和L3分别为中频频段和高频频段主要频率点的数目，k1＝1和k2＝1分别为计算特征F4所需中间向量在两个中间矩阵中的对应列，
经过以上对中间矩阵的计算，得到特征F3和F4，与之前提取的特征F1和F2一起，构建扰动信号的特征相量F＝[F1,F2,F3,F4]；
6)基于粒子群算法的极限学习机分类器构建：
极限学习机输入权值和隐藏层节点的任意选取，一方面会增加算法对测试数据的反应时间，另一方面也增加了网络结构的复杂度，因此，应用粒子群算法优化的极限学习机以改善以往极限学习机存在的弊端，粒子群算法通过优化极限学习机的输入权重和隐藏层结点偏置实现对神经网络的改进，从而获得最优模型；对扰动信号使用特征提取方法得到4维特征向量，训练分类器用于识别电能质量扰动信号类型，当需要识别未知电能质量扰动信号时，将采集的信号通过特征提取方法得到输入向量，输入到已训练好的分类器中，判断该信号属于哪一类电能质量扰动。
本发明提供一种采用时域压缩多分辨率快速S变换特征提取的电能质量扰动识别方法, 首先，针对电能质量扰动信号进行快速傅里叶变换(Fast Fourier Transform,FFT)，利用得到的傅里叶谱值确定扰动信号的主要频率点；之后，对每个主要频率点进行加窗傅里叶变换逆变换(Inverse Fast Fourier Transform,IFFT)，同时在确定特征的前提下，对每个主要频率点的加窗傅里叶变换逆变换结果计算所需保留信息，构造中间矩阵，通过对中间矩阵的数学运算实现有效特征的提取；最后，设计基于粒子群算法优化的极限学习机分类器，完成对复杂电能质量扰动信号的识别。具有方法科学合理，适用性强，效果佳等优点。
附图说明
图1为本发明的总体流程图；
图2为基于时域压缩快速S变换获取中间矩阵流程图。
具体实施方式
本发明采用基于时域压缩快速S变换特征提取的电能质量扰动识别方法对扰动信号进行特征提取与识别，其具体实施包括以下步骤：
1)使用Matlab仿真生成12类电能质量扰动信号；
2)使用多分辨率快速S变换处理扰动信号；
3)确定所提取的特征种类；
4)基于时域压缩快速S变换构造中间向量与中间矩阵；
5)对中间矩阵进一步提取特征，构建特征向量；
6)构建基于粒子群算法的极限学习机分类器，对扰动信号进行分类。
下面结合附图和具体实施例对本发明做进一步描述。
图1为本发明的总体流程图。
A、Matlab仿真生成电能质量扰动信号
使用Matlab按照数学模型仿真生成信噪比为20dB～50dB间随机值的12种电能质量扰动信号，包括暂降、暂升、中断、尖峰、闪变、谐波、暂态振荡、切痕、谐波含暂降、谐波含暂升、谐波含闪变以及谐波含振荡，信号采样率为3.2kHz，采样时间为1秒，基频为50Hz；
B、使用多分辨率快速S变换处理扰动信号
在多分辨率快速S变换的基础上，首先，对仿真扰动信号进行快速傅里叶变换，根据傅立叶谱在频域范围内确定主要频率点；之后，只针对主要频率点对应的部分开展傅里叶变换逆变换，得到时频矩阵；
设t为时间，f为频率，τ为时移参数，当f→nx/NT,τ→jT时，多分辨率快速S变换的离散表达式为：
式中，N为信号长度(本发明中N＝3200)，T为采样时间间隔，j和m均为0至N-1范围内的整数，λx为窗宽调整因子，nx为根据傅里叶谱值确定的扰动信号主要频率点，x为0 至M-1范围内的整数，M为主要频率点的个数。由原始信号经快速傅立叶变换后位移获得。各频段频率范围划分为：频率点在1Hz至100Hz范围内为低频频段，在101Hz至 700Hz范围内为中频频段，在701Hz至(N/2)Hz范围内为高频频段。在此基础上，将λx在不同频率范围内分别定义。低频频段内λx＝2；中频频段内高频频段内，通过基频(50Hz)傅里叶谱分析判断信号中是否含有扰动基频扰动，若无扰动，取否则，取
C、确定需要提取的扰动特征
经过对各特征分类效果分析，确定对扰动信号提取4类特征用于构建特征向量，包括： F1、原始信号1/4周期能量跌落幅度；F2、基频幅值标准差；F3、700Hz以上高频频段能量；F4、100Hz以上各频率对应最大幅值的最大值，
各特征的具体计算方法分别为：
(a)原始信号1/4周期能量跌落幅度：
式中，h(t)为原始信号，t为信号采样点，R(z)为原始信号1/4周期的均方根值(Root Mean Square，RMS)，z＝1,2……200，R0为1/4周期无噪声标准电能质量信号的均方根值，取值为0.7071；
(b)基频50Hz幅值标准差：
式中，t为信号采样点，f0为信号基频频率，50Hz，N为信号长度，N＝3200，S(t,f0)为信号经快速傅里叶变换后，对基频50Hz频率点进行加窗傅里叶变换逆变换后得到的行向量；
(c)700Hz以上高频频段保留能量：
式中，S(t,f)为模时频矩阵，t为采样点，f为频率，N为信号长度，N＝3200；
(d)100Hz以上各频率对应最大幅值的最大值：
F4＝max{max(S(t,f))} 101Hz≤f≤(N/2)Hz (6)
式中，S(t,f)为模时频矩阵，t为采样点，f为频率，N为信号长度，N＝3200；
对以上4类特征分别以三种方式进行提取，
方式一：以式(2)、式(3)为依据，对原始信号进行数学运算后提取特征F1；
方式二：以式(4)为依据，对原始信号进行快速傅立叶变换后，对基频(50Hz)频率点对应部分进行加窗傅里叶变换逆变换后得到的行相量直接进行数学运算，获得特征F2，避免了对模时频矩阵整体运算和储存的复杂过程；
方式三：在对信号进行快速S变换处理的过程中，频率点在700Hz以上的，以式(5)为依据，对此频率范围内主要频率点在加窗傅里叶变换逆变换得到行相量同时提取各行向量所需保留的特征计算中间值，逐个保留频率对应行计算完毕后获得计算特征F3所需中间特征相量；频率点在100Hz以上的，以式(6)为依据，对此频率范围内主要频率点在加窗傅里叶变换逆变换得到行相量同时提取各行向量所需保留的特征计算中间值，逐个保留频率对应行计算完毕后获得计算特征F4所需中间特征相量，随即对以上各中间特征相量进行储存，由此形成的特征矩阵定义为中间矩阵，最后，对中间矩阵进行数学运算，分别提取特征F3和 F4；
D、计及时域压缩与特征提取需要的快速S变换构造中间向量与中间矩阵；
为进一步减小时频矩阵的时域维数，降低矩阵储存空间，采用基于时域压缩多分辨率快速S变换构造中间矩阵；在快速S变换处理信号的过程中，对各主要频率点在加窗傅里叶变换逆变换得到行相量同时提取各行向量所需保留的特征计算中间值，用于构建中间矩阵；经压缩后的中间矩阵时域维数减少，所需的存储空间减小，降低了相关分析的空间存储要求与存储成本；
设原始信号为h(t)，t为信号采样点，构造中间矩阵的运算流程为：
①对原始信号h(t)进行快速傅里叶变换，得到快速傅里叶变换频谱H(m)，即 H(m)＝FFT(h(t))，其中m为0至N-1范围内的整数，N为信号长度，N＝3200；
②对傅里叶谱取模平方得|H(m)|2，设|H(m)|2共有L个极大点，所有极大值序列记为 |H[l]|2，l＝0,1,2,...,L-1，若|H(m)|2中相邻的极大点之间分别有X1,...,Xl,...,XL-1个频率点(Xl为|H(m)|2中第l个极大点与第l+1个极大值之间的频率点数目)，则|H(m)|2的极大值包络计算方式如下：
式中，x＝0,1,...,Xl-1；l＝0,1,...,L-1，为|H(m)|2中第l个极大点对应的数值，由此得到|H(m)|2的极大值包络谱之后，利用动态测度法求中每一个极大点的动态测度。动态测度的求解方法如下：
设M为的一个极大点，如果M两侧各有一个或多个比点M更高的极大点时，点M两侧也一定存在两个与点M同高度的点N1、N2，由点M通向点N1的路径记为P1，M通向 N2的路径记为P2，则极大点M的动态测度等于路径P1的动态测度与路径P2的动态测度中的较小者，设M的动态测度为Dyn(M)，则有：
式中，halt(M)和halt(N)分别表示点M、N的高度，P(M,N)为点M、N之间的路径，表示这两点之间的部分；Dyn(P(M,N))为路径P(M,N)的动态测度，表示此路径上最高点和最低点的高度差，inf表示下确界，表示一个集合中的最小值，通过以上计算，得到极大|H(m)|2值包络的动态测度，通过设置阈值，将大于该阈值的动态测度所对应的频率点定为主要频率点nx，其中，x为0至M-1范围内的整数，M为主要频率点的个数；
③根据公式确定主要频率点nx对应的高斯窗函数W(m,nx)；
④将H(m)移位到H(m+nx)，对H(m+nx)进行加窗处理，得到H(m+nx)·W(m,nx)；
⑤对加窗的H(m+nx)进行快速傅里叶变换逆变换，得到频率点nx对应的加窗傅立叶变换逆变换的结果S(j,nx)，即S(j,nx)＝IFFT(H(m+nx)·W(m,nx))，其中j为采样点，随即对该行提取特征并存储，形成中间矩阵，中间矩阵构造过程具体步骤如下：
根据700Hz以上高频频段能量特征F3和100Hz以上各频率对应最大幅值的最大值F4的原理，在700Hz以上频率范围内，对主要频率点对应加窗傅立叶变换逆变换的结果按照式(9) 进行数学运算，提取特征相量；在100Hz以上频率范围内，对主要频率点对应加窗傅立叶变换逆变换的结果按照式(10)进行数学运算，提取特征相量，对相量进行储存，构建中间矩阵；
设S(t,f)为信号经快速傅里叶变换后，对频率f进行加窗傅里叶变换逆变换后得到的行向量，t为采样点，N为扰动信号采样点数目，N＝3200，中间矩阵构建规则为：
xmax＝max{S(t,f)} 101Hz≤f≤(N/2)Hz (10)
式中，xmax为频率f下对应幅值中的最大值，E为频率f下的能量，由此分别求得中间矩阵中的特征向量m1和m2，m1＝[m11,m12,…,m1d]T为100Hz以上特定频率各采样点对应幅值中的最大值，m2＝[m21,m22,…,m2d]T为700Hz以上特定频率对应的能量，d为此频段内主要频率点数目，由以上特征相量分别构建各频段下的中间矩阵，中间矩阵结构为：
中频频段，101Hz≤f≤700Hz的中间矩阵结构如式(11)所示，式中：Vm为中频频段的中间矩阵，L2为中频频段主要频率点数目，
高频频段(701Hz≤f≤(N/2)Hz，N＝3200)的中间矩阵结构如式(12)所示，式中：Vh为高频频频段的中间矩阵，L3为高频频段主要频率点数目，
经过以上计算，最终得到规模为L2×1的中频中间矩阵Vm，和规模为L3×2的高频中间矩阵Vh，经时域压缩快速S变换后所得中频频段和高频频段中间矩阵列数分别为1列和2列，与未时域压缩的3200列时频矩阵相比，所需储存空间明显减小；
E、基于时域压缩快速S变换中间矩阵的特征向量构建
在D中获得各频段中间矩阵的基础上，对中间矩阵按照式(13)、式(14)进行数学运算，最终提取特征F3和F4：
式中，Vh为高频频段中间矩阵，i为中间矩阵中的各行，k为中间矩阵中能量特征对应的列(k＝2)，L3为高频频段中主要频率点的数目，
F4＝max[max(Vm(i,k1)),max(Vh(j,k2))] (14)
式中，Vm和Vh分别为中频频段和高频频段的中间矩阵，i和j分别为两个矩阵中的行， i∈[1,L2]，J∈[1,L3]，L2和L3分别为中频频段和高频频段主要频率点的数目，k1＝1和k2＝1分别为计算特征F4所需中间向量在两个中间矩阵中的对应列，
经过以上对中间矩阵的计算，得到特征F3和F4，与之前提取的特征F1和F2一起，构建扰动信号的特征相量F＝[F1,F2,F3,F4]；
F、基于粒子群算法的极限学习机分类器构建
极限学习机输入权值和隐藏层节点的任意选取，一方面会增加算法对测试数据的反应时间，另一方面也增加了网络结构的复杂度，因此，应用粒子群算法优化的极限学习机以改善以往极限学习机存在的弊端，粒子群算法通过优化极限学习机的输入权重和隐藏层结点偏置实现对神经网络的改进，从而获得最优模型；对扰动信号使用特征提取方法得到4维特征向量，训练分类器用于识别电能质量扰动信号类型，当需要识别未知电能质量扰动信号时，将采集的信号通过特征提取方法得到输入向量，输入到已训练好的分类器中，判断该信号属于哪一类电能质量扰动。
为了体现本发明特征提取方法与分类器结合时的分类效果，使用Matlab仿真生成信噪比为20dB～50dB间随机值的12类各扰动信号各300组，包括暂降(C1)、暂升(C2)、中断(C3)、尖峰(C4)、闪变(C5)、谐波(C6)、暂态振荡(C7)、切痕(C8)、谐波含暂降(C9)、谐波含暂升(C10)、谐波含闪变(C11)以及谐波含振荡(C12)。使用本发明特征提取方法对各扰动信号提取特征，构建特征集合，用于训练基于粒子群算法优化的极限学习机分类器；此外，生成信噪比分别为50dB、40dB、30dB、20dB的扰动信号每类各100组，使用本发明特征提取方法对各扰动信号提取特征，构建各测试集合。测试结果如上表所示：
经以上分析可知，计及时域压缩快速S变换所得中频频段和高频频段中间矩阵列数分别为1列和2列，相较于时域压缩前的3200列时频矩阵而言，所需数据存储空间大大减小，存储成本大大降低，另外，本发明特征提取方法与分类器结合可以得到良好的分类效果，所以，本发明具有实际工业应用价值。

Claims (1)

1.一种采用时域压缩多分辨率快速S变换特征提取的电能质量扰动识别方法，其特征是：它包括电能质量扰动信号仿真，对扰动信号进行多分辨率快速S变换处理，在确定特征的前提下，对快速S变换过程中每个主要频率点的加窗傅里叶变换逆变换结果计算所需保留信息，构造中间矩阵，对中间矩阵提取有效特征构建特征向量，用于分类器构建，具体包括的步骤有：
1)电能质量扰动信号仿真：
使用Matlab按照数学模型仿真生成信噪比为20dB～50dB间随机值的12种电能质量扰动信号，包括暂降、暂升、中断、尖峰、闪变、谐波、暂态振荡、切痕、谐波含暂降、谐波含暂升、谐波含闪变以及谐波含振荡，信号采样率为3.2kHz，采样时间为1秒，基频为50Hz；
2)将电能质量扰动信号利用多分辨率快速S变换方法进行处理：
在多分辨率快速S变换的基础上，首先，对仿真扰动信号进行快速傅里叶变换，根据傅立叶谱在频域范围内确定主要频率点；之后，只针对主要频率点对应的部分开展傅里叶变换逆变换，得到时频矩阵；
设t为时间，f为频率，τ为时移参数，当f→nx/NT,τ→jT时，多分辨率快速S变换的离散表达式为：
<mrow> <mi>S</mi> <mrow> <mo>(</mo> <mi>j</mi> <mi>T</mi> <mo>,</mo> <mfrac> <msub> <mi>n</mi> <mi>x</mi> </msub> <mrow> <mi>N</mi> <mi>T</mi> </mrow> </mfrac> <mo>)</mo> </mrow> <mo>=</mo> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>m</mi> <mo>=</mo> <mn>0</mn> </mrow> <mrow> <mi>N</mi> <mo>-</mo> <mn>1</mn> </mrow> </munderover> <mi>H</mi> <mrow> <mo>(</mo> <mfrac> <mrow> <mi>m</mi> <mo>+</mo> <msub> <mi>n</mi> <mi>x</mi> </msub> </mrow> <mrow> <mi>N</mi> <mi>T</mi> </mrow> </mfrac> <mo>)</mo> </mrow> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mn>2</mn> <msup> <mi>&amp;pi;</mi> <mn>2</mn> </msup> <msup> <mi>m</mi> <mn>2</mn> </msup> <mo>/</mo> <msup> <msub> <mi>&amp;lambda;</mi> <mi>x</mi> </msub> <mn>2</mn> </msup> <msup> <msub> <mi>n</mi> <mi>x</mi> </msub> <mn>2</mn> </msup> </mrow> </msup> <msup> <mi>e</mi> <mrow> <mi>i</mi> <mn>2</mn> <mi>&amp;pi;</mi> <mi>m</mi> <mi>j</mi> <mo>/</mo> <mi>N</mi> </mrow> </msup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </mrow>
式中，N为信号长度(本发明中N＝3200)，T为采样时间间隔，j和m均为0至N-1范围内的整数，λx为窗宽调整因子，nx为根据傅里叶谱值确定的扰动信号主要频率点，x为0至M-1范围内的整数，M为主要频率点的个数。由原始信号经快速傅立叶变换后位移获得。各频段频率范围划分为：频率点在1Hz至100Hz范围内为低频频段，在101Hz至700Hz范围内为中频频段，在701Hz至(N/2)Hz范围内为高频频段。在此基础上，将λx在不同频率范围内分别定义。低频频段内λx＝2；中频频段内高频频段内，通过基频(50Hz)傅里叶谱分析判断信号中是否含有扰动基频扰动，若无扰动，取否则，取
3)确定要提取的特征种类：
经过对各特征分类效果分析，确定对扰动信号提取4类特征用于构建特征向量，包括：F1、原始信号1/4周期能量跌落幅度；F2、基频幅值标准差；F3、700Hz以上高频频段能量；F4、100Hz以上各频率对应最大幅值的最大值，
各特征的具体计算方法分别为：
(a)原始信号1/4周期能量跌落幅度：
<mrow> <mi>F</mi> <mn>1</mn> <mo>=</mo> <mfrac> <mrow> <mi>m</mi> <mi>i</mi> <mi>n</mi> <mo>&amp;lsqb;</mo> <mi>R</mi> <mrow> <mo>(</mo> <mi>z</mi> <mo>)</mo> </mrow> <mo>&amp;rsqb;</mo> </mrow> <msub> <mi>R</mi> <mn>0</mn> </msub> </mfrac> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <mi>R</mi> <mrow> <mo>(</mo> <mi>z</mi> <mo>)</mo> </mrow> <mo>=</mo> <msqrt> <mrow> <mfrac> <mn>1</mn> <mn>16</mn> </mfrac> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>t</mi> <mo>=</mo> <mn>16</mn> <mi>z</mi> <mo>-</mo> <mn>15</mn> </mrow> <mrow> <mi>t</mi> <mo>=</mo> <mn>16</mn> <mi>z</mi> </mrow> </munderover> <msup> <mi>h</mi> <mn>2</mn> </msup> <mrow> <mo>(</mo> <mi>t</mi> <mo>)</mo> </mrow> </mrow> </msqrt> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo> </mrow> </mrow>
式中，h(t)为原始信号，t为信号采样点，R(z)为原始信号1/4周期的均方根值(RootMean Square，RMS)，z＝1,2……200，R0为1/4周期无噪声标准电能质量信号的均方根值，取值为0.7071；
(b)基频50Hz幅值标准差：
<mrow> <mi>F</mi> <mn>2</mn> <mo>=</mo> <msqrt> <mrow> <mfrac> <mn>1</mn> <mi>N</mi> </mfrac> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>t</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>N</mi> </munderover> <msup> <mrow> <mo>(</mo> <mi>S</mi> <mo>(</mo> <mi>t</mi> <mo>,</mo> <msub> <mi>f</mi> <mn>0</mn> </msub> <mo>)</mo> <mo>-</mo> <mfrac> <mn>1</mn> <mi>n</mi> </mfrac> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>t</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>N</mi> </munderover> <mi>S</mi> <mo>(</mo> <mi>t</mi> <mo>,</mo> <msub> <mi>f</mi> <mn>0</mn> </msub> <mo>)</mo> <mo>)</mo> </mrow> <mn>2</mn> </msup> </mrow> </msqrt> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </mrow>
式中，t为信号采样点，f0为信号基频频率，50Hz，N为信号长度，N＝3200，S(t,f0)为信号经快速傅里叶变换后，对基频50Hz频率点进行加窗傅里叶变换逆变换后得到的行向量；
(c)700Hz以上高频频段保留能量：
<mrow> <mtable> <mtr> <mtd> <mrow> <mi>F</mi> <mn>3</mn> <mo>=</mo> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>f</mi> <mo>=</mo> <mn>701</mn> </mrow> <mrow> <mi>N</mi> <mo>/</mo> <mn>2</mn> </mrow> </munderover> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>t</mi> <mo>=</mo> <mn>1</mn> </mrow> <mn>3200</mn> </munderover> <msup> <mrow> <mo>&amp;lsqb;</mo> <mi>S</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>,</mo> <mi>f</mi> <mo>)</mo> </mrow> <mo>&amp;rsqb;</mo> </mrow> <mn>2</mn> </msup> </mrow> </mtd> <mtd> <mrow> <mn>701</mn> <mi>H</mi> <mi>z</mi> <mo>&amp;le;</mo> <mi>f</mi> <mo>&amp;le;</mo> <mrow> <mo>(</mo> <mi>N</mi> <mo>/</mo> <mn>2</mn> <mo>)</mo> </mrow> <mi>H</mi> <mi>z</mi> </mrow> </mtd> </mtr> </mtable> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>5</mn> <mo>)</mo> </mrow> </mrow>
式中，S(t,f)为模时频矩阵，t为采样点，f为频率，N为信号长度，N＝3200；
(d)100Hz以上各频率对应最大幅值的最大值：
F4＝max{max(S(t,f))}101Hz≤f≤(N/2)Hz (6)
式中，S(t,f)为模时频矩阵，t为采样点，f为频率，N为信号长度，N＝3200；
对以上4类特征分别以三种方式进行提取，
方式一：以式(2)、式(3)为依据，对原始信号进行数学运算后提取特征F1；
方式二：以式(4)为依据，对原始信号进行快速傅立叶变换后，对基频(50Hz)频率点对应部分进行加窗傅里叶变换逆变换后得到的行相量直接进行数学运算，获得特征F2，避免了对模时频矩阵整体运算和储存的复杂过程；
方式三：在对信号进行快速S变换处理的过程中，频率点在700Hz以上的，以式(5)为依据，对此频率范围内主要频率点在加窗傅里叶变换逆变换得到行相量同时提取各行向量所需保留的特征计算中间值，逐个保留频率对应行计算完毕后获得计算特征F3所需中间特征相量；频率点在100Hz以上的，以式(6)为依据，对此频率范围内主要频率点在加窗傅里叶变换逆变换得到行相量同时提取各行向量所需保留的特征计算中间值，逐个保留频率对应行计算完毕后获得计算特征F4所需中间特征相量，随即对以上各中间特征相量进行储存，由此形成的特征矩阵定义为中间矩阵，最后，对中间矩阵进行数学运算，分别提取特征F3和F4；
4)计及时域压缩与特征提取需要的中间向量与中间矩阵构造过程：
为进一步减小时频矩阵的时域维数，降低矩阵储存空间，采用基于时域压缩多分辨率快速S变换构造中间矩阵；在快速S变换处理信号的过程中，对各主要频率点在加窗傅里叶变换逆变换得到行相量同时提取各行向量所需保留的特征计算中间值，用于构建中间矩阵；经压缩后的中间矩阵时域维数减少，所需的存储空间减小，降低了相关分析的空间存储要求与存储成本；
设原始信号为h(t)，t为信号采样点，构造中间矩阵的运算流程为：
①对原始信号h(t)进行快速傅里叶变换，得到快速傅里叶变换频谱H(m)，即H(m)＝FFT(h(t))，其中m为0至N-1范围内的整数，N为信号长度，N＝3200；
②对傅里叶谱取模平方得|H(m)|2，设|H(m)|2共有L个极大点，所有极大值序列记为|H[l]|2，l＝0,1,2,...,L-1，若|H(m)|2中相邻的极大点之间分别有X1,...,Xl,...,XL-1个频率点(Xl为|H(m)|2中第l个极大点与第l+1个极大值之间的频率点数目)，则|H(m)|2的极大值包络计算方式如下：
<mrow> <mo>|</mo> <mover> <mi>H</mi> <mo>^</mo> </mover> <mrow> <mo>(</mo> <mi>m</mi> <mo>)</mo> </mrow> <msup> <mo>|</mo> <mn>2</mn> </msup> <mo>=</mo> <mo>|</mo> <mover> <mi>H</mi> <mo>~</mo> </mover> <mrow> <mo>(</mo> <mi>l</mi> <mo>-</mo> <mn>1</mn> <mo>)</mo> </mrow> <msup> <mo>|</mo> <mn>2</mn> </msup> <mo>+</mo> <mi>x</mi> <mfrac> <mrow> <mo>|</mo> <mover> <mi>H</mi> <mo>~</mo> </mover> <mrow> <mo>(</mo> <mi>l</mi> <mo>)</mo> </mrow> <msup> <mo>|</mo> <mn>2</mn> </msup> <mo>-</mo> <mo>|</mo> <mover> <mi>H</mi> <mo>~</mo> </mover> <mrow> <mo>(</mo> <mi>l</mi> <mo>-</mo> <mn>1</mn> <mo>)</mo> </mrow> <msup> <mo>|</mo> <mn>2</mn> </msup> </mrow> <msub> <mi>X</mi> <mi>l</mi> </msub> </mfrac> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>7</mn> <mo>)</mo> </mrow> </mrow>
式中，x＝0,1,...,Xl-1；l＝0,1,...,L-1，为|H(m)|2中第l个极大点对应的数值，由此得到|H(m)|2的极大值包络谱之后，利用动态测度法求中每一个极大点的动态测度。动态测度的求解方法如下：
设M为的一个极大点，如果M两侧各有一个或多个比点M更高的极大点时，点M两侧也一定存在两个与点M同高度的点N1、N2，由点M通向点N1的路径记为P1，M通向N2的路径记为P2，则极大点M的动态测度等于路径P1的动态测度与路径P2的动态测度中的较小者。设M的动态测度为Dyn(M)，则有：
<mrow> <mtable> <mtr> <mtd> <mrow> <msub> <mi>D</mi> <mrow> <mi>y</mi> <mi>n</mi> </mrow> </msub> <mrow> <mo>(</mo> <mi>M</mi> <mo>)</mo> </mrow> <mo>=</mo> <mrow> <mo>(</mo> <msub> <mi>i</mi> <mrow> <mi>n</mi> <mi>f</mi> </mrow> </msub> <mo>(</mo> <msub> <mi>D</mi> <mrow> <mi>y</mi> <mi>n</mi> </mrow> </msub> <mo>(</mo> <mi>P</mi> <mo>(</mo> <mi>M</mi> <mo>,</mo> <mi>N</mi> <mo>)</mo> </mrow> <mo>)</mo> <mo>)</mo> <mo>;</mo> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mi>h</mi> <mrow> <mi>a</mi> <mi>l</mi> <mi>t</mi> </mrow> </msub> <mrow> <mo>(</mo> <mi>M</mi> <mo>)</mo> </mrow> <mo>=</mo> <msub> <mi>h</mi> <mrow> <mi>a</mi> <mi>l</mi> <mi>t</mi> </mrow> </msub> <mrow> <mo>(</mo> <mi>N</mi> <mo>)</mo> </mrow> </mrow> </mtd> </mtr> </mtable> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>8</mn> <mo>)</mo> </mrow> </mrow>
式中，halt(M)和halt(N)分别表示点M、N的高度，P(M,N)为点M、N之间的路径，表示这两点之间的部分；Dyn(P(M,N))为路径P(M,N)的动态测度，表示此路径上最高点和最低点的高度差，inf表示下确界，表示一个集合中的最小值，通过以上计算，得到极大|H(m)|2值包络的动态测度，通过设置阈值，将大于该阈值的动态测度所对应的频率点定为主要频率点nx，其中，x为0至M-1范围内的整数，M为主要频率点的个数；
③根据公式确定主要频率点nx对应的高斯窗函数W(m,nx)；
④将H(m)移位到H(m+nx)，对H(m+nx)进行加窗处理，得到H(m+nx)·W(m,nx)；
⑤对加窗的H(m+nx)进行快速傅里叶变换逆变换，得到频率点nx对应的加窗傅立叶变换逆变换的结果S(j,nx)，即S(j,nx)＝IFFT(H(m+nx)·W(m,nx))，其中j为采样点，随即对该行提取特征并存储，形成中间矩阵，中间矩阵构造过程具体步骤如下：
根据700Hz以上高频频段能量特征F3和100Hz以上各频率对应最大幅值的最大值F4的原理，在700Hz以上频率范围内，对主要频率点对应加窗傅立叶变换逆变换的结果按照式(9)进行数学运算，提取特征相量；在100Hz以上频率范围内，对主要频率点对应加窗傅立叶变换逆变换的结果按照式(10)进行数学运算，提取特征相量，对相量进行储存，构建中间矩阵；
设S(t,f)为信号经快速傅里叶变换后，对频率f进行加窗傅里叶变换逆变换后得到的行向量，t为采样点，N为扰动信号采样点数目，N＝3200，中间矩阵构建规则为：
<mrow> <mtable> <mtr> <mtd> <mrow> <mi>E</mi> <mo>=</mo> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>t</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>N</mi> </munderover> <msup> <mrow> <mo>{</mo> <mi>S</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>,</mo> <mi>f</mi> <mo>)</mo> </mrow> <mo>}</mo> </mrow> <mn>2</mn> </msup> </mrow> </mtd> <mtd> <mrow> <mn>701</mn> <mi>H</mi> <mi>z</mi> <mo>&amp;le;</mo> <mi>f</mi> <mo>&amp;le;</mo> <mrow> <mo>(</mo> <mi>N</mi> <mo>/</mo> <mn>2</mn> <mo>)</mo> </mrow> <mi>H</mi> <mi>z</mi> </mrow> </mtd> </mtr> </mtable> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>9</mn> <mo>)</mo> </mrow> </mrow>
xmax＝max{S(t,f)}101Hz≤f≤(N/2)Hz (10)
式中，xmax为频率f下对应幅值中的最大值，E为频率f下的能量，由此分别求得中间矩阵中的特征向量m1和m2，m1＝[m11,m12,…,m1d]T为100Hz以上特定频率各采样点对应幅值中的最大值，m2＝[m21,m22,…,m2d]T为700Hz以上特定频率对应的能量，d为此频段内主要频率点数目，由以上特征相量分别构建各频段下的中间矩阵，中间矩阵结构为：
中频频段，101Hz≤f≤700Hz的中间矩阵结构如式(11)所示，式中：Vm为中频频段的中间矩阵，L2为中频频段主要频率点数目，
<mrow> <msub> <mi>V</mi> <mi>m</mi> </msub> <mo>=</mo> <msub> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <msub> <mi>m</mi> <mn>11</mn> </msub> </mtd> </mtr> <mtr> <mtd> <mo>.</mo> </mtd> </mtr> <mtr> <mtd> <mo>.</mo> </mtd> </mtr> <mtr> <mtd> <mo>.</mo> </mtd> </mtr> <mtr> <mtd> <msub> <mi>m</mi> <mrow> <mn>1</mn> <msub> <mi>L</mi> <mn>2</mn> </msub> </mrow> </msub> </mtd> </mtr> </mtable> </mfenced> <mrow> <msub> <mi>L</mi> <mn>2</mn> </msub> <mo>&amp;times;</mo> <mn>1</mn> </mrow> </msub> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>11</mn> <mo>)</mo> </mrow> </mrow>
高频频段(701Hz≤f≤(N/2)Hz，N＝3200)的中间矩阵结构如式(12)所示，式中：Vh为高频频频段的中间矩阵，L3为高频频段主要频率点数目，
<mrow> <msub> <mi>V</mi> <mi>h</mi> </msub> <mo>=</mo> <msub> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <msub> <mi>m</mi> <mn>11</mn> </msub> </mtd> <mtd> <msub> <mi>m</mi> <mn>21</mn> </msub> </mtd> </mtr> <mtr> <mtd> <mo>.</mo> </mtd> <mtd> <mo>.</mo> </mtd> </mtr> <mtr> <mtd> <mo>.</mo> </mtd> <mtd> <mo>.</mo> </mtd> </mtr> <mtr> <mtd> <mo>.</mo> </mtd> <mtd> <mo>.</mo> </mtd> </mtr> <mtr> <mtd> <msub> <mi>m</mi> <mrow> <mn>1</mn> <msub> <mi>L</mi> <mn>3</mn> </msub> </mrow> </msub> </mtd> <mtd> <msub> <mi>m</mi> <mrow> <mn>2</mn> <msub> <mi>L</mi> <mn>3</mn> </msub> </mrow> </msub> </mtd> </mtr> </mtable> </mfenced> <mrow> <msub> <mi>L</mi> <mn>3</mn> </msub> <mo>&amp;times;</mo> <mn>2</mn> </mrow> </msub> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>12</mn> <mo>)</mo> </mrow> </mrow>
经过以上计算，最终得到规模为L2×1的中频中间矩阵Vm，和规模为L3×2的高频中间矩阵Vh，经时域压缩快速S变换后所得中频频段和高频频段中间矩阵列数分别为1列和2列，与未时域压缩的3200列时频矩阵相比，所需储存空间明显减小；
5)基于时域压缩快速S变换中间矩阵的特征向量构建：
在步骤4)中获得各频段中间矩阵的基础上，对中间矩阵按照式(13)、式(14)进行数学运算，最终提取特征F3和F4：
<mrow> <mi>F</mi> <mn>3</mn> <mo>=</mo> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mrow> <mi>L</mi> <mn>3</mn> </mrow> </munderover> <msup> <mrow> <mo>&amp;lsqb;</mo> <msub> <mi>V</mi> <mi>h</mi> </msub> <mrow> <mo>(</mo> <mi>i</mi> <mo>,</mo> <mi>k</mi> <mo>)</mo> </mrow> <mo>&amp;rsqb;</mo> </mrow> <mn>2</mn> </msup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>13</mn> <mo>)</mo> </mrow> </mrow>
式中，Vh为高频频段中间矩阵，i为中间矩阵中的各行，k为中间矩阵中能量特征对应的列(k＝2)，L3为高频频段中主要频率点的数目，
F4＝max[max(Vm(i,k1)),max(Vh(j,k2))] (14)
式中，Vm和Vh分别为中频频段和高频频段的中间矩阵，i和j分别为两个矩阵中的行，i∈[1,L2]，J∈[1,L3]，L2和L3分别为中频频段和高频频段主要频率点的数目，k1＝1和k2＝1分别为计算特征F4所需中间向量在两个中间矩阵中的对应列，
经过以上对中间矩阵的计算，得到特征F3和F4，与之前提取的特征F1和F2一起，构建扰动信号的特征相量F＝[F1,F2,F3,F4]；
6)基于粒子群算法的极限学习机分类器构建：
极限学习机输入权值和隐藏层节点的任意选取，一方面会增加算法对测试数据的反应时间，另一方面也增加了网络结构的复杂度，因此，应用粒子群算法优化的极限学习机以改善以往极限学习机存在的弊端，粒子群算法通过优化极限学习机的输入权重和隐藏层结点偏置实现对神经网络的改进，从而获得最优模型；对扰动信号使用特征提取方法得到4维特征向量，训练分类器用于识别电能质量扰动信号类型，当需要识别未知电能质量扰动信号时，将采集的信号通过特征提取方法得到输入向量，输入到已训练好的分类器中，判断该信号属于哪一类电能质量扰动。
CN201710948335.XA 2017-10-12 2017-10-12 一种采用时域压缩多分辨率快速s变换特征提取的电能质量扰动识别方法 Active CN107832777B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201710948335.XA CN107832777B (zh) 2017-10-12 2017-10-12 一种采用时域压缩多分辨率快速s变换特征提取的电能质量扰动识别方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201710948335.XA CN107832777B (zh) 2017-10-12 2017-10-12 一种采用时域压缩多分辨率快速s变换特征提取的电能质量扰动识别方法

Publications (2)

Publication Number Publication Date
CN107832777A true CN107832777A (zh) 2018-03-23
CN107832777B CN107832777B (zh) 2021-01-26

Family

ID=61647913

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201710948335.XA Active CN107832777B (zh) 2017-10-12 2017-10-12 一种采用时域压缩多分辨率快速s变换特征提取的电能质量扰动识别方法

Country Status (1)

Country Link
CN (1) CN107832777B (zh)

Cited By (16)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109239458A (zh) * 2018-09-06 2019-01-18 华北电力大学（保定） 一种高噪声背景下电能质量扰动信号降噪方法
CN109271975A (zh) * 2018-11-19 2019-01-25 燕山大学 一种基于大数据多特征提取协同分类的电能质量扰动识别方法
CN109856503A (zh) * 2018-12-27 2019-06-07 国网江苏省电力有限公司检修分公司 一种基于s变换及同步相量测量的输电线路故障定位方法
CN110046593A (zh) * 2019-04-22 2019-07-23 三峡大学 基于分段改进s变换和随机森林的复合电能质量扰动识别方法
CN110147760A (zh) * 2019-05-20 2019-08-20 吉林化工学院 一种高效电能质量扰动图像特征提取与识别新方法
CN110222953A (zh) * 2018-12-29 2019-09-10 北京理工大学 一种基于深度学习的电能质量混合扰动分析方法
CN110231117A (zh) * 2019-06-10 2019-09-13 中设设计集团股份有限公司 一种基于s变换的拉索基频特征辨识方法
CN110728195A (zh) * 2019-09-18 2020-01-24 武汉大学 基于yolo算法的电能质量扰动检测方法
CN111241902A (zh) * 2019-07-24 2020-06-05 成都理工大学 一种高精度多重同步压缩广义s变换时频分析方法
CN111325485A (zh) * 2020-03-22 2020-06-23 东北电力大学 计及物联带宽约束的轻量级梯度提升机电能质量扰动识别方法
CN111999607A (zh) * 2020-07-20 2020-11-27 中国南方电网有限责任公司超高压输电公司广州局 一种单通道信号下局部放电窄带干扰盲源分离方法及装置
CN113238110A (zh) * 2021-05-10 2021-08-10 合肥工业大学 一种电能质量扰动诊断方法及系统
CN113361427A (zh) * 2021-06-11 2021-09-07 温州大学 基于wt与ga-pnn的复合电能质量扰动分类方法
CN113392813A (zh) * 2021-07-30 2021-09-14 湖北工业大学 精确识别振动信号主频的方法及系统
CN114638270A (zh) * 2022-03-31 2022-06-17 合肥工业大学 一种电能质量扰动特征精确快速提取方法
CN114799610A (zh) * 2022-06-24 2022-07-29 苏芯物联技术(南京)有限公司 一种基于傅里叶逆变换及自编码器的焊接质量实时检测方法及系统

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101924354A (zh) * 2010-04-19 2010-12-22 昆明理工大学 利用s变换能量样本属性的配电网故障的人工神经网络选线方法
CN102944773A (zh) * 2012-10-29 2013-02-27 华南理工大学 基于空间变换的电能扰动检测及分类方法
CN103308804A (zh) * 2013-06-17 2013-09-18 湖南大学 基于快速k-s变换电能质量扰动信号时频参数提取方法
CN104034974A (zh) * 2014-05-09 2014-09-10 西南交通大学 一种复合电能质量扰动信号的识别方法
CN104459398A (zh) * 2014-12-08 2015-03-25 东北电力大学 一种采用二维形态学降噪的电能质量复合扰动识别方法
CN104730384A (zh) * 2015-03-16 2015-06-24 华南理工大学 基于不完全s变换的电能扰动识别与定位方法

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101924354A (zh) * 2010-04-19 2010-12-22 昆明理工大学 利用s变换能量样本属性的配电网故障的人工神经网络选线方法
CN102944773A (zh) * 2012-10-29 2013-02-27 华南理工大学 基于空间变换的电能扰动检测及分类方法
CN103308804A (zh) * 2013-06-17 2013-09-18 湖南大学 基于快速k-s变换电能质量扰动信号时频参数提取方法
CN104034974A (zh) * 2014-05-09 2014-09-10 西南交通大学 一种复合电能质量扰动信号的识别方法
CN104459398A (zh) * 2014-12-08 2015-03-25 东北电力大学 一种采用二维形态学降噪的电能质量复合扰动识别方法
CN104730384A (zh) * 2015-03-16 2015-06-24 华南理工大学 基于不完全s变换的电能扰动识别与定位方法

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
XIANYONG XIAO: "Short duration disturbance classifying based on S-transform maximum similarity", 《ELECTRICAL POWER AND ENERGY SYSTEMS》 *
黄南天 等: "基于粒子群与极限学习机的电能质量信号特征选择与识别", 《电工电能新技术》 *
黄南天 等: "采用改进多分辨率快速S 变换的电能质量扰动识别", 《电网技术》 *
黄南天 等: "采用最优多分辨率快速S 变换的电能质量分析", 《仪器仪表学报》 *

Cited By (23)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109239458A (zh) * 2018-09-06 2019-01-18 华北电力大学（保定） 一种高噪声背景下电能质量扰动信号降噪方法
CN109239458B (zh) * 2018-09-06 2021-03-19 华北电力大学（保定） 一种高噪声背景下电能质量扰动信号降噪方法
CN109271975A (zh) * 2018-11-19 2019-01-25 燕山大学 一种基于大数据多特征提取协同分类的电能质量扰动识别方法
CN109271975B (zh) * 2018-11-19 2020-08-21 燕山大学 一种基于大数据多特征提取协同分类的电能质量扰动识别方法
CN109856503A (zh) * 2018-12-27 2019-06-07 国网江苏省电力有限公司检修分公司 一种基于s变换及同步相量测量的输电线路故障定位方法
CN110222953A (zh) * 2018-12-29 2019-09-10 北京理工大学 一种基于深度学习的电能质量混合扰动分析方法
CN110046593B (zh) * 2019-04-22 2023-08-01 三峡大学 基于分段改进s变换和随机森林的复合电能质量扰动识别方法
CN110046593A (zh) * 2019-04-22 2019-07-23 三峡大学 基于分段改进s变换和随机森林的复合电能质量扰动识别方法
CN110147760A (zh) * 2019-05-20 2019-08-20 吉林化工学院 一种高效电能质量扰动图像特征提取与识别新方法
CN110147760B (zh) * 2019-05-20 2022-08-12 吉林化工学院 一种高效电能质量扰动图像特征提取与识别新方法
CN110231117A (zh) * 2019-06-10 2019-09-13 中设设计集团股份有限公司 一种基于s变换的拉索基频特征辨识方法
CN111241902A (zh) * 2019-07-24 2020-06-05 成都理工大学 一种高精度多重同步压缩广义s变换时频分析方法
CN110728195A (zh) * 2019-09-18 2020-01-24 武汉大学 基于yolo算法的电能质量扰动检测方法
CN111325485A (zh) * 2020-03-22 2020-06-23 东北电力大学 计及物联带宽约束的轻量级梯度提升机电能质量扰动识别方法
CN111999607B (zh) * 2020-07-20 2021-08-20 中国南方电网有限责任公司超高压输电公司广州局 一种单通道信号下局部放电窄带干扰盲源分离方法及装置
CN111999607A (zh) * 2020-07-20 2020-11-27 中国南方电网有限责任公司超高压输电公司广州局 一种单通道信号下局部放电窄带干扰盲源分离方法及装置
CN113238110A (zh) * 2021-05-10 2021-08-10 合肥工业大学 一种电能质量扰动诊断方法及系统
CN113361427A (zh) * 2021-06-11 2021-09-07 温州大学 基于wt与ga-pnn的复合电能质量扰动分类方法
CN113361427B (zh) * 2021-06-11 2024-03-15 温州大学 基于wt与ga-pnn的复合电能质量扰动分类方法
CN113392813A (zh) * 2021-07-30 2021-09-14 湖北工业大学 精确识别振动信号主频的方法及系统
CN114638270A (zh) * 2022-03-31 2022-06-17 合肥工业大学 一种电能质量扰动特征精确快速提取方法
CN114799610A (zh) * 2022-06-24 2022-07-29 苏芯物联技术(南京)有限公司 一种基于傅里叶逆变换及自编码器的焊接质量实时检测方法及系统
CN114799610B (zh) * 2022-06-24 2022-10-04 苏芯物联技术(南京)有限公司 一种基于傅里叶逆变换及自编码器的焊接质量实时检测方法及系统

Also Published As

Publication number Publication date
CN107832777B (zh) 2021-01-26

Similar Documents

Publication Publication Date Title
CN107832777A (zh) 一种采用时域压缩多分辨率快速s变换特征提取的电能质量扰动识别方法
Wang et al. Wind speed forecasting based on the hybrid ensemble empirical mode decomposition and GA-BP neural network method
CN106846173A (zh) 基于ewt‑esn的短期风电功率预测方法
CN101604356B (zh) 一种不确定性中长期水文预报模型的建立方法
CN101975910B (zh) 一种特高压直流输电线路故障智能分类与测距方法
CN106875041A (zh) 一种短期风速预测方法
CN106933778A (zh) 一种基于爬坡事件特征识别的风电功率组合预测方法
CN104459398B (zh) 一种采用二维形态学降噪的电能质量复合扰动识别方法
Bigdeli et al. A comparative study of optimal hybrid methods for wind power prediction in wind farm of Alberta, Canada
CN104792522A (zh) 一种基于分数阶小波变换和bp神经网络的齿轮缺陷智能分析方法
CN105488466A (zh) 一种深层神经网络和水声目标声纹特征提取方法
CN102609766B (zh) 一种风电场风速智能预测方法
CN113238110B (zh) 一种电能质量扰动诊断方法
CN103876736B (zh) 一种基于功率谱划分的复杂度谱脑电信号分类识别方法
CN111325485B (zh) 计及物联带宽约束的轻量级梯度提升机电能质量扰动识别方法
CN105160204A (zh) 一种碳排放价格组合预测方法
CN104459397A (zh) 采用自适应多分辨率广义s变换的电能质量扰动识别方法
CN110399796A (zh) 一种基于改进深度学习算法的电能质量扰动识别方法
CN103413134A (zh) 基于稀疏分解的地面移动目标微动信号特征提取
CN106068017B (zh) 基于小波神经元网络的无线链路质量预测方法
CN102315642B (zh) 一种电网规划方案的环境影响评价方法
CN106548136A (zh) 一种无线信道场景分类方法
CN104376402A (zh) 一种基于频域指标的负荷特性分类与综合方法
CN103577877B (zh) 一种基于时频分析和bp神经网络的船舶运动预报方法
Zhang et al. Second-order coupled tristable stochastic resonance and its application in bearing fault detection under different noises

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant
TR01 Transfer of patent right

Effective date of registration: 20231204

Address after: 201100 2nd floor, building 11, 988 Zhongchun Road, Minhang District, Shanghai

Patentee after: Luoyi Technology (Shanghai) Co.,Ltd.

Address before: 230000 floor 1, building 2, phase I, e-commerce Park, Jinggang Road, Shushan Economic Development Zone, Hefei City, Anhui Province

Patentee before: Dragon totem Technology (Hefei) Co.,Ltd.

Effective date of registration: 20231204

Address after: 230000 floor 1, building 2, phase I, e-commerce Park, Jinggang Road, Shushan Economic Development Zone, Hefei City, Anhui Province

Patentee after: Dragon totem Technology (Hefei) Co.,Ltd.

Address before: 132022 No. 45 Chengde street, Longtan District, Jilin, Jilin

Patentee before: JILIN INSTITUTE OF CHEMICAL TECHNOLOGY

TR01 Transfer of patent right