CN109145476A - 用于电力系统信号处理的时域自适应分段复指数级数法 - Google Patents

Abstract

本发明公开了用于电力系统信号处理的时域自适应分段复指数级数法，涉及电力系统和现代信号处理技术领域，解决了现有电力系统信号处理方法不能直接计算出电力系统信号分量的幅值、频率、相位和衰减因子4个特性参数及不能处理带不连续非平稳分量的电力系统信号的问题，其技术方案要点是：首先，采用复指数级数对电力系统信号的时域采样序列进行拟合；然后，以方均根误差最小为目标自适应调整所截取采样序列的长度，实现采样序列的分段；最后，计算各子段上的4个特性参数。本发明不仅具有能同时获取电力系统信号分量的4个特性参数及处理带不连续非平稳分量的电力系统信号的优点，还可用于电力系统信号重构，在电力系统领域具有广阔的应用前景。

Description

用于电力系统信号处理的时域自适应分段复指数级数法

技术领域

本发明属于电力系统和现代信号处理技术领域，更具体地说，其涉及用于电力系统信号处理的时域自适应分段复指数级数法。

背景技术

电力系统的信号包括电压、电流、频率、有功功率、无功功率和视在功率等。电力系统信号处理是电力系统分析的基础，高质量的电力系统信号处理是实现电力系统安全、稳定、经济、优质运行的前提条件。电力系统信号处理的本质是用数学方法将电力系统信号所含的各种分量分离开来，并采用特定的数学工具予以描述，从而提取或揭示重要的信息用以指导决策，包括时域和频域两类方法。

目前广泛采用的电力系统信号处理方法有离散傅里叶变换(discrete Fouriertransform，DFT)、小波变换(wavelet transform，WT)和自回归滑动平均建模(autoregressive moving average，ARMA)等。DFT及其各种变化形式(如傅里叶级数和快速傅里叶变换)仅适于处理平稳信号，当采样序列的时间长度不能准确地为各信号分量周期的整数倍或当信号含有随时间衰减或增长的直流分量、周期分量及非平稳噪声分量时，计算结果将受到极大影响。WT是时频处理方法，适于处理非平稳信号，但小波基函数的选取带有经验性，且分量的时频表征方式不能直观反映电力系统的特性，实际应用不方便。ARMA模型适宜描述平稳随机过程，故对具有非平稳性的电力系统信号进行建模缺乏严格的理论支持，且ARMA模型仅能给出信号分量的频率和幅值，在应用中具有局限性。

综上，为满足电力系统各分支专业的应用需求，要求信号处理方法应能同时获取信号分量的幅值、频率、相位和衰减因子4个特性参数，且须在时域内进行表征。因此，如何设计出能同时获取信号分量的幅值、频率、相位和衰减因子4个特性参数并能处理带不连续非平稳分量的电力系统信号的信号处理方法是目前迫切需要解决的问题。

发明内容

针对现有信号处理方法存在的问题，本发明提供了用于电力系统信号处理的时域自适应分段复指数级数法，不仅具有能同时获取电力系统信号分量的幅值、频率、相位和衰减因子4个特性参数及处理带不连续非平稳分量的电力系统信号的优点，还可用于电力系统信号的重构，即根据需要对已有电力系统信号进行插值或延拓，在电力系统领域具有广阔的应用前景。

本发明的上述技术目的是通过以下技术方案实现的：首先，采用复指数级数对电力系统信号的时域采样序列进行拟合；然后，以方均根误差(root mean square error，RMSE)最小为目标自适应地调整所截取的采样序列的长度实现采样序列的分段；最后，计算给出各子段上电力系统信号分量的幅值、频率、相位和衰减因子4个特性参数。

进一步，所述用于电力系统信号处理的时域自适应分段复指数级数法包括：

对电力系统的时域连续信号x(t)按采样频率f_s等距采样后，得到N个时域采样数据x(nT_s)，简记为时域采样序列x_n，n＝0,1,···,N–1；采用p个复指数组成级数对各x_n进行线性表示，即拟合：

式中，0<p≤N/2；T_s＝1/f_s为采样间距/周期(单位s)；j为虚数单位；A_i、θ_i、α_i和f_i分别为第i个复指数的幅值(单位视电力系统信号的具体物理意义而定)、初相角(单位°)、衰减因子(单位s^–1)和频率(单位Hz)；设：

b_i为复留数，s_i为复频率。

问题归结为已知x_n，求解未知参数b_i和z_i。方法是视x_n为p阶齐次前向差分方程：

的解。该齐次差分方程的通解形式为x_n＝czⁿ(c为待定常数)，将其代入式(3)可得一特征方程：

约去常数c，两边同乘z^–n，式(4)变为：

求解式(5)可得p个特征根z_i，并得到式(3)的齐次解：

式中，c_i为待定常数，由边界条件确定。令c_i＝b_i，则式(6)就与式(1)相同。

综上，复指数级数法的基本步骤如下：先由式(3)解出a_k，k＝0,1,···,p–1；再由式(5)解出z_i；最后由式(1)或式(6)解出b_i。

由式(3)求解a_k时，为使用全部N个采样数据以最大限度地利用x_n所携带的信息，采用x_n的后向排序序列进行后向预测，即将x_n依次向后延时(平移)一个时间单位(序号)后代入式(3)，得一方程组，其矩阵形式为：

简记为：

Xa＝x_A (8)

其中，X为(N–p)×p矩阵，a和x_A均为(N–p)×1列向量。式(8)为超定方程，须采用最小二乘法(least squares，LS)求解，得到LS解为：

a_LS＝(X^TX)^–1X^Tx_A (9)

式中，上标T表示矩阵转置运算，X^TX常严重秩亏和病态，无法正常求逆，需用Moore-Penrose的广义逆(伪逆)替代计算(方法1)，或依据LS与奇异值分解(singularvaluedecomposition，SVD)的关系进行数值估计(方法2)。方法1以Frobenius范数||(X^TX)(X^TX)⁺–I||_F最小化为目标，其中I为p×p单位矩阵，(X^TX)⁺为X^TX的伪逆，用以替代(X^TX)^–1；方法2以2范数||Xa–x_A||₂最小化为目标，与LS的目标等价。方法2比方法1在理论上更严谨，而实验也表明方法2的精度比方法1高3～4个数量级，故本发明采用方法2计算a_LS，步骤为：先求X的SVD，得到X＝UΣV^T，其中U和V分别为(N–p)×(N–p)和p×p正交矩阵，Σ为(N–p)×p矩阵(其对角元素为X的p个奇异值σ_i，其余元素为0)；再将U和V改写为列向量的形式U＝[u₁,u₂,···,u_N–p]和V＝[v₁,v₂,···,v_p]；最后由下式计算a_LS：

将求得的a_k代入式(5)求解z_i。为克服传统迭代法求解的低效性及其无法求取重根或复根的缺点，本发明提出一快速求解方法。设p×p友矩阵：

则A的p阶特征多项式为：

式中，I为p×p单位矩阵；符号|·|表示行列式计算；z为A的特征值。将B_p按第1行展开，得：

B_p＝zB_p–1+(–1)¹⁺p(–1)^p–1a₀＝zB_p–1+a₀ (13)

用相同的方法求取B_p–1，并对式(13)进行递推可得：

B_p＝zB_p–1+a₀＝z(zB_p–2+a₁)+a₀＝z²(zB_p–3+a₂)+a₁z+a₀＝···

＝zp+a_p–1z^p–1+···+a2z²+a₁z+a₀ (14)

式(14)的系数取负后即为A第p列的元素，故可由式(14)的系数a_k按式(11)写出A；令式(14)等于0，即得到A的特征方程，且该方程与式(5)相同，故式(5)的根就是A的特征值。因此，对A进行特征值分解(eigenvalue decomposition，EVD)，所得特征值就是方程(5)的根z_i。该方法经一次EVD即可求出方程(5)的所有根，故求解效率和速度比传统迭代法高很多。另外，SVD可由两步EVD来实现，故信号处理程序仅需编制一个EVD子程序，即可实现SVD和方程(5)的求解。

求得z_i后，将其与x_n一起代入式(1)得一方程组，其矩阵形式为：

简记为：

Zb＝x_B (16)

Z为Vandermonde矩阵。可用求解式(8)的方法求解式(16)，故不再赘述。

最后，由求得的b_i和z_i根据式(2)求出各信号分量的4个特性参数：

A_i＝|b_i|，θ_i＝argb_i，α_i＝ln|z_i|/T_s，f_i＝argz_i/(2πT_s) (17)

式中arg[·]表示求复数的幅角主值，有argb_i、arg z_i∈(–π,π]。求幅角主值时须判断幅角所处的象限，以argb_i为例，设其实部和虚部为Re(b_i)和Im(b_i)，有：

式中，arctan[·]∈(–π/2,π/2)。注意，按式(17)计算得到的是复指数级数的参数，并非电力系统常用的实数型正、余弦信号的参数，还需通过Euler公式将2个复指数相加合成为1个实余弦信号，故除衰减因子α_i外，实余弦信号的幅值、相角和频率分别为2A_i、|θ_i|和|f_i|(对于实余弦信号分量，θ_i和f_i是正负对称出现的)。

经实验发现，当复指数的个数p取最大值N/2时，可完全保证复指数级数对采样序列x_n的拟合性，不会出现算法失稳的问题。但考虑到N可能为奇数，应用中应取p＝int_[N/2]，符号int_表示向负向取整，比如int_[5.95]＝5。然而，经实验还发现，p＝int_[N/2]虽能保证拟合性，但可能出现过拟合现象，即拟合数据在采样序列的尾部出现轻微振荡，使拟合误差有所增加。而改取：

p＝int_[N/2]–1 (19)

则可在保证拟合性的同时避免出现过拟合现象。

进一步，所述用于电力系统信号处理的时域自适应分段复指数级数法的自适应分段策略具体包括：

因为复指数是关于时间t的连续函数，故要求被拟合的电力系统信号在采样序列占据的整个时间段内各个分量持续存在且特性参数不会发生突变。但部分电力系统信号，比如由数字故障录波器记录的输电线路自动重合闸的断续信号，可能不满足该要求。因此，本发明提出自适应分段策略，用以自动确定不连续信号分量的起点和终点，在起点和终点间的各个子段上分别采用复指数级数法实现拟合，从而保证复指数级数法的适用性和精确性。

设n_s和n_e分别为子段起点和终点的序号，n_Δ为搜索步长(子段每次增加的点数)，e_f为子段的RMSE，L_min、N_Δm和E_m分别为预先指定的子段的最小长度、最大搜索步长和最大RMSE。RMSE的计算公式如下：

式中，L＝n_e–n_s+1为子段长度；||·||₂为2-范数运算；为原时域采样序列子段的向量；为用复指数级数拟合后的子段的向量；的元素按式(1)计算：

式中，以符号∧标识拟合值；子段复指数的个数按式(19)确定，p_s＝int_[(n_e–n_s+1)/2]–1；注意，n的起始值不是n_s，而是与式(1)中相同的0，相当于将x_n的0起点右移到了当前子段的起点n_s，也就是将每个子段均视为一由0起始的新序列。

进一步，所述用于电力系统信号处理的时域自适应分段复指数级数法的算法流程如图1，具体步骤为：

开始时，取子段长度为L_min，可能因数据点少使e_f≥E_m；此时，保持n_Δ＝1不变并执行分支④，避免n_Δ>1使当前子段的n_e跨越分界点延伸至下一子段，继续使e_f≥E_m，陷入死循环。在循环执行若干次分支④后，找到使e_f首次小于E_m的n_e，接着在满足e_f<E_m的条件下，取n_Δ＝N_Δm，经循环执行分支②搜索子段的最大长度(即终点)。当再次出现e_f≥E_m时，转而执行分支③，以确定出精确的分界点，方法是将n_Δ从N_Δm逐步减半，采用图中的函数int₄₅[n_Δ/2]实现(符号int₄₅表示4舍5入取整，比如int₄₅[3.4]＝3、int₄₅[3.5]＝4)。当e_f<E_m再次满足时，程序又转回执行分支②，因此时n_Δ已改变故不再置其为N_Δm。分支②和③交替执行，直到n_Δ被减为0，表明精确分界点已找到，则计算并存储当前子段的特性参数，完成后将n_Δ重置为1并调整n_s和n_e经分支①进行下一子段的搜索。

附图说明

图1为本发明提供的用于电力系统信号处理的时域自适应分段复指数级数法的算法流程；

图2为本发明实施例1提供的电力系统仿真信号及其拟合信号的波形；

图3为本发明实施例1提供的5种信号插值方法的插值结果及理想插值结果；

图4为本发明实施例1提供的采用时域自适应分段复指数级数法对电力系统仿真信号同时进行延拓和插值后的结果；

图5为本发明实施例2提供的真实电力系统故障信号及其拟合信号的波形；

图6为本发明实施例2提供的采用时域自适应分段复指数级数法对真实电力系统信号同时进行延拓和插值后的结果。

具体实施方式

以下结合附图1～6和实施例，对本发明作进一步说明。实施例1将时域自适应分段复指数级数法应用于一个由仿真产生的已知各分量特性参数且带有不连续非平稳分量的电力系统信号，以验证其有效性；实施例2将时域自适应分段复指数级数法应用于一个各分量特性参数未知的真实电力系统信号，以验证其实用性。

实施例1：包括电力系统仿真信号的拟合及其特性参数求取、电力系统仿真信号的插值及延拓3个子例。

1、电力系统仿真信号的拟合及其特性参数求取。

电力系统仿真信号，以下简称仿真信号，其离散表达式为：

式中，N＝600，T_s＝1ms，f_s＝1kHz；各信号分量x_i(nT_s)的表达式为：

具体参数如表1，其包含4个间谐波、1个3次谐波和1个衰减直流分量。

表1 6个信号分量的参数

序号i	A<sub>i</sub>	α<sub>i</sub>/s<sup>–1</sup>	f<sub>i</sub>/Hz	θ<sub>i</sub>/(°)
					1	3.0	–2.0	12.0	30.0
2	2.0	–1.0	10.0	20.0
					3	1.0	–0.5	5.0	10.0
4	0.7	–3.0	8.0	15.0
					5	3.5	–5.0	250.0	160.0
6	1.5	–1.5	0.0	0.0

为模拟突然接入的及非持续的信号分量，设置：

式中，序号乘以T_s即得到相应的时间值；x₁(t)从0s接入，持续到0.122s消失；x₅(t)在0～0.122s内不存在，从0.123s接入，持续到0.299s消失；x₆(t)在0～0.299s内不存在，从0.3s接入，一直持续到仿真信号的终止时刻；其余3个信号分量均由0s接入，一直持续到仿真信号的终止时刻。6个信号分量按式(22)合成的仿真信号波形如图2(a)。

首先，暂不加入自适应分段策略，直接将复指数级数法用于对仿真信号进行整体的拟合。在式(21)中取n_e–n_s＝N–1，将计算得到的拟合波形叠加绘制在仿真信号上，如图2(b)。可见，虽然拟合数据点在整体上很接近仿真信号，但在某些局部(如图中的虚线圈部分)，可看到拟合数据点偏离了仿真信号，此时p＝299，已达到复指数个数允许的最大值300，但e_f≈1.731×10³，表明拟合精度仍非常低；若减小p，则e_f更大，甚至出现无法拟合的情况。此时计算出的特性参数误差过大，故不再列出。因此，如果不加入自适应分段策略，复指数级数法在处理不连续信号分量时会存在缺陷。

接下来，用时域自适应分段复指数级数法对仿真信号进行拟合，选取实验参数L_min＝20、N_Δm＝20、E_m＝1×10^–6。同样，将由式(21)计算得到的各子段的拟合信号依次叠加绘制在仿真信号上，如图2(c)。图中的虚竖线为各子段的分界线，可见，时域自适应分段复指数级数法将仿真信号自动分为了4个子段，拟合数据点全部都在仿真信号上，从视觉上已无法看出差别。为评判分段效果和拟合精度，在表2中列出了各子段的详细信息及特性参数的计算值。可见，时域自适应分段复指数级数法准确地寻找到了式(24)设定的序号122、299和300，且各子段的e_f数量级均在10^–20左右，说明拟合精度非常高。

最后，通过对比表2和表1，进一步确认计算得出的特性参数的精度。先将表2子段1的4组特性参数与表1的前4行参数对比(根据式(24)的设定子段1仅包含x₁(t)～x₄(t))，可见两者相互对应，仅排列顺序相反，表明时域自适应分段复指数级数法准确地找到了子段1应包含的4个信号分量，且定量分析显示，计算得出的特性参数与其真实值间的相对误差在10^–7量级。再将表2子段2的4组特性参数与表1的2～5行参数对比(根据式(24)的设定子段2仅包含x₂(t)～x₅(t))，可见两者的α_i、f_i相差很小，而A_i、θ_i却相差甚远。

表2各子段的详细信息及特性参数计算值

注意，这并不是计算有误，而是由时间偏移导致的“正确结果”：子段2各信号分量的A_i、θ_i实际上是x₂(t)～x₅(t)的A_i、θ_i经0.122s后的结果，比如，x₄(t)的幅值初值为0.7，初相角为15.0°，在α_i、f_i的作用下经0.122s后，分别变为0.7e^{–3.0×0.122}和arg[15.0°+2π×8.0×0.122×180°/π]，具体数值见表3子段2数据的第1行；同样处理后，将x₃(t)、x₂(t)和x₅(t)的A_i、θ_i经0.122s后的数值列于表3子段2数据的2～4行，这时再次将表3子段2的数据与表1的2～5行参数对比，可见两者的数值一致了，定量分析显示两者的相对误差也在10^–7量级。后续的第3、第4子段也有类似情况，可仿照上述方法处理，其各信号分量经相应时间偏移后的A_i、θ_i的数值也列于表3中。总之，若已知一信号分量的接入时

表3考虑时间偏移量后各子段A_i、θ_i的计算值

间，欲使其在各子段上的分析结果一致，在计算A_i、θ_i时应计及相应的时间偏移量，这是在具体实施时要特别注意的。

2、电力系统仿真信号的插值。

信号插值是多采样率系统统一采样频率的主要手段，用于将低采样频率序列向高采样频率序列转换；此外，对一采样周期较大的采样序列，有时需获知其采样点间某些特定时刻的采样点，此时也需对原采样序列进行插值，以补出那些未被记录下来的采样点。

此处算例的目标是通过插值将仿真信号的采样频率由f_s＝1kHz提高为f_s1＝10kHz，采样周期由T_s＝1ms缩小为T_s1＝0.1ms。这需要在仿真信号的每2个点间插入9个点，即将1个T_s等分为10份，每1份为1个T_s1。因仿真信号各分量对f_s、f_s1均满足采样定理，故插值不会造成频谱混叠，无需增加滤波环节。

采用时域自适应分段复指数级数法实现信号插值的方法为：对于仿真信号的数据点，将其直接拷贝至相应的时间点上；而对于仿真信号数据点间的点，则由式(21)的如下变形形式产生：

上式的实质是在仿真信号的2个整数序号间插入9个浮点数序号，比如0,0.1,···,0.9,1和1,1.1,···,1.9,2等，而与这些浮点数序号对应的数据点和即为0～1和1～2间应插入的数据点。由式(1)可知，指数n+0.1k在计算时会与T_s相乘形成对应的时间(n+0.1k)T_s，这正是各整数序号所对应时间之间更细的时间，故即为这些更细时间上的数据。

以下将时域自适应分段复指数级数法的信号插值性能与多项式插值、填0滤波插值这两类目前使用最广的传统信号插值方法进行对比，以验证其有效性。图3依次示出了以上信号插值方法的插值结果和理想的插值结果，图中以实线表示f_s时的仿真信号，以点线表示f_s1时的插值信号。说明：①因理论证明高次插值不可取，且实验也表明4次及以上多项式的插值结果并不比3次多项式好，故图中仅给出3例低次多项式的插值结果；②默认以f_s1时的仿真波形作为理想的插值结果即插入数据的真实值，由此计算各信号插值方法插值结果的RMSE(e_f)，e_f也于各子图中示出；③因插值后数据总量增加10倍，插入数据点过于密集，从图上无法看出差异，故仅截取了0.115～0.135s的一段波形(因其包含0.122s这一分界点，很具代表性)进行横向展宽，以较清晰地揭示插值的细节。

图3(a)为1次Lagrange多项式的插值结果。因两插值结点(仿真信号的数据点)间用直线连接，故插入数据点在仿真信号两点间的直线上排列齐整，看似很完美，然而e_f却很大，对比图3(f)可知，这是因为除插值结点外，多数插入的数据点均非f_s1时波形上的数据点。

图3(b)为3次Hermite多项式的插值结果。Hermite方法在Lagrange方法的基础上增加了插值结点处的导数约束，比Lagrange方法多使用了一些数据信息，故插值结果与理想插值结果稍接近了一些，e_f也比图3(a)小，但仍很大。

图3(c)为3次样条函数的插值结果。样条插值的波形不仅光滑，且其e_f较之前两种方法进一步减小，波形已比较接近图3(f)的理想波形，但主要缺陷在于各子段的末端会因样条插值固有的末端效应而出现摆动(比如第1子段的末端)。

图3(d)为填0滤波的插值结果。该方法是在拟插入数据的位置上填上0，用0替代其数值，然后经抗混叠滤波器滤除折叠频率以上因填0作用复制产生的频率分量，插值效果与滤波器密切相关。经实验发现，长度10、截止数字域频率0.5(对应模拟域频率0.5f_s1/(2π)≈795.77Hz)的数字滤波器可使e_f基本为最小，且使插值精度又较图3(c)提高了一些。然而，这是在理想插值结果已知的情况下通过多次尝试设计的数字滤波器，在实际应用中不具通用性，难以保证数字滤波器参数对所有信号都是最优的；其次，数字滤波器也存在末端效应。

图3(e)为时域自适应分段复指数级数法的插值结果。其插值精度较之图3(d)有大幅提高，e_f的量级降为10^–7。

小结：造成4种传统插值方法插值误差大的主要原因在于其均不具备信号信息的提取和利用能力；而时域自适应分段复指数级数法能充分提取和利用信号信息，自适应地、目标明确地实现插值。

3、电力系统仿真信号的延拓。

由一段已知信号推知其终止时刻后或起始时刻前的波形，称为信号的正向延拓(也称为预测)或负向延拓(也称为反演)。若电力系统信号经时域自适应分段复指数级数法实现了分段，则正向延拓波形与最后一个子段的相关性最强，负向延拓波形与第1子段的相关性最强，故分别按最后一个子段和第1个子段的复指数级数拟合式实现正向和负向延拓。因延拓和插值同时进行的实用价值更大，仅延拓、不插值的情况可视为其简化形式，故以仿真信号延拓、插值同时进行的情况为例。

此处算例的目标是对仿真信号分别进行正、负向延拓，各延拓0.1s时长，并同时进行插值。设N_F和N_B分别为正、负向延拓的点数，因0.1s/T_s1＝1000，故正、负向需各延拓N_F＝1000点和N_B＝1000点。

先进行负向延拓。因仿真信号的起始时刻为0，延拓信号的时间值为负值，故各数据点的序号不仅为浮点数，且为负数。负向延拓的公式为：

式中，参数为计算得到的第1子段的参数；注意，计算式(26)前，应先提取第1子段的主导信号分量(p_s取主导信号分量的个数)，去除虚假分量，用主导分量的参数进行计算，因的指数为负数，相当于做除法，若某些虚假分量的过小，将使计算出的有很大误差。

接下来进行正向延拓。延拓数据点的序号为正浮点数，无负向延拓的问题，故无需提取主导分量，但为保证波形光滑，须从最后一个子段的最后一个数据点开始添加延拓数据。设最后一个数据点的序号为n_last(相当于最后一个子段的n_s，由第4子段的信息知n_last＝146)，正向延拓的公式为：

经正、负向延拓并插值后的仿真信号如图4的点线。图中，实线为由式(22)～(24)按T_s1在–0.1～0.7s上计算得到的精确波形，可见，–0.1～0s的负向延拓波形和0.6～0.7s的正向延拓波形与精确波形重合得很好，定量分析表明其e_f均在10^–7量级。

实施例2：包括真实电力系统信号的拟合及其特性参数求取、插值和延拓2个子例。

1、真实电力系统信号的拟合及其特性参数求取。

图5(a)为一由数字故障录波器记录的电力系统中某条输电线路发生短路故障前后的实际故障电流信号，为便于处理对数值进行了归一化处理即取标么值。

首先，暂不加入自适应分段策略，直接将复指数级数法用于对实际故障电流信号进行整体的拟合，结果如图5(b)。可见，在0～0.25s拟合信号与实际故障电流信号重合较好，但到末端时拟合信号出现了发散，相应的e_f很大，求取特性参数已无意义。这再次说明，如果不加入自适应分段策略，复指数级数法在处理含有不连续信号分量的电力系统信号时会存在缺陷。

接下来，用时域自适应分段复指数级数法对实际故障电流信号进行拟合，选取实验参数L_min＝20、N_Δm＝20、E_m＝5×10^–4，适当加大E_m是因为实际故障电流信号成分复杂，寻求过高的拟合精度是不现实的，若E_m过小，则判据e_f<E_m可能永远无法得到满足。拟合结果如图5(c)，图中用虚竖线示出了各子段的分界线。因实际故障电流信号所含噪声分量较多，故仅在表4中列出各子段的e_f及主导分量的特性参数。时域自适应分段复指数级数法将实际故障电流信号分为了6个子段：1～3子段的基波的幅值基本相同且较小，据此可判定此3子段为故障前的正常运行电流；到第4子段，基波的幅值迅速增大为前3段的约5倍，据此可判定此子段为短路暂态电流；到第5子段，基波的幅值较之第4子段减小很多，可判定继电保护在动作切除故障输电线路；到第6子段，基波缺失了，仅剩下1个伴随整个故障过程的3次谐波，这是由电力系统的非线性元件产生的1个环绕于电力变压器三角形绕组的固有环流。可见，正是因为实际故障电流信号中基波的非持续性导致了不加入自适应分段策略时算法的失效。注意，表中前5个子段的RMSE均满足e_f<E_m，而第6子段的e_f>E_m，这是因为其为最后一段，无再多的数据使e_f继续下降，只能将其独立为一段，但其e_f还是可接受的。本子例还可视为是时域自适应分段复指数级数法在电力系统故障诊断方面的一个应用案例。

表4故障信号各子段的e_f及若干主导分量的特性参数

2、真实电力系统信号的插值和延拓。

本子例的目标是将实际故障电流信号的采样频率从原先的1kHz经插值后提高到10kHz，且在正、负向延拓各0.1s(即5个工频周波)。将式(26)和(27)用于实际故障电流信号，经插值、延拓后的实际故障电流信号如图6，注意，因实际故障电流信号在采样前会经过抗混叠滤波，而时域自适应分段复指数级数法不会无端地产生频率分量，故插值、延拓不会造成频谱混叠。图中，0s虚竖线左侧为负向延拓信号，0.3s虚竖线右侧为正向延拓信号。可见，正向延拓信号有逐渐增大且出现振荡的趋势，这表明与经典预测方法一样，随着预测时长的增加，预测的精度和可信度都将下降，因为当前的已知信息绝不能完全取代未来的信息；同样，对负向延拓(反演)也如此。不过，若对延拓时长加以限制，比如仅延拓0.1s左右，则精度和可信度均可得到保证。

以上所述仅为本发明的实施例，并不用以限制本发明，凡在本发明精神和原则内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围内。

Claims (7)

1.用于电力系统信号处理的时域自适应分段复指数级数法，其特征在于，首先，采用复指数级数对电力系统信号的时域采样序列进行拟合；然后，以方均根误差最小为目标自适应地调整所截取的采样序列的长度，实现采样序列的分段；最后，计算给出各子段上电力系统信号分量的幅值、频率、相位和衰减因子4个特性参数。

2.如权利要求1所述的用于电力系统信号处理的时域自适应分段复指数级数法，其特征在于，所述复指数级数法具体包括：

对电力系统的时域连续信号x(t)按采样频率f_s等距采样后，得到N个时域采样数据x(nT_s)，简记为时域采样序列x_n，n＝0,1,···,N–1；采用p个复指数组成级数对各x_n进行线性表示，即拟合：

式中，0<p≤N/2；T_s＝1/f_s为采样间距/周期(单位s)；j为虚数单位；A_i、θ_i、α_i和f_i分别为第i个复指数的幅值(单位视电力系统信号的具体物理意义而定)、初相角(单位°)、衰减因子(单位s^–1)和频率(单位Hz)；设：

已知x_n求解未知参数b_i和z_i，设x_n为p阶齐次前向差分方程的解：

该齐次差分方程的特征方程为：

复指数级数法的基本步骤如下：先由式(3)解出a_k，k＝0,1,···,p–1；再由式(4)解出z_i；最后由式(1)解出b_i；

由式(3)求解a_k时，采用x_n的后向排序序列进行后向预测，即将x_n依次向后延时、平移一个时间单位、序号后代入式(3)，得一方程组，其矩阵形式为：

简记为：

Xa＝x_A (6)

其中，X为(N–p)×p矩阵，a和x_A均为(N–p)×1列向量；式(6)为超定方程，采用最小二乘法(least squares，LS)求解，得到LS解为：

a_LS＝(X^TX)^–1X^Tx_A (7)

式中，上标T表示矩阵转置运算；X^TX常严重秩亏和病态，无法正常求逆；依据LS与奇异值分解(singularvalue decomposition，SVD)的关系对a_LS进行数值估计，该方法以2范数||Xa–x_A||₂最小化为目标，与LS的目标等价，步骤为：

先求X的SVD，得到X＝UΣV^T，其中U和V分别为(N–p)×(N–p)和p×p正交矩阵，Σ为(N–p)×p矩阵，其对角元素为X的p个奇异值σ_i，其余元素为0；再将U和V改写为列向量的形式U＝[u₁,u₂,···,u_N–p]和V＝[v₁,v₂,···,v_p]；最后由下式计算a_LS：

将求得的a_k代入式(4)求解z_i，所采用的快速求解方法为：首先，用a_k组成p×p友矩阵：

然后，对A进行特征值分解(eigenvalue decomposition，EVD)，所得特征值就是方程(4)的根z_i；

求得z_i后，将其与x_n一起代入式(1)得一方程组，其矩阵形式为：

简记为：

Zb＝x_B (11)

用求解式(6)的方法求解式(11)；

最后，由求得的b_i和z_i根据式(2)求出各信号分量的4个特性参数：

A_i＝|b_i|，θ_i＝arg b_i，α_i＝ln|z_i|/T_s，f_i＝arg z_i/(2πT_s) (12)

式中arg[·]表示求复数的幅角主值，有arg b_i、arg z_i∈(–π,π]；求幅角主值时须判断幅角所处的象限，以arg b_i为例，设其实部和虚部为Re(b_i)和Im(b_i)，有：

式中，arctan[·]∈(–π/2,π/2)；按式(12)计算得到的是复指数级数的参数，通过Euler公式将2个复指数相加合成为1个实余弦信号，除衰减因子α_i外，实余弦信号的幅值、相角和频率分别为2A_i、|θ_i|和|f_i|；

复指数的个数p按下式确定：

p＝int-[N/2]–1 (14)

式中，符号int-表示向负向取整。

3.如权利要求2所述的用于电力系统信号处理的时域自适应分段复指数级数法，其特征在于，所述自适应分段策略具体包括：

设n_s和n_e分别为子段起点和终点的序号，n_Δ为搜索步长，即子段每次增加的点数，e_f为子段的RMSE，L_min、N_Δm和E_m分别为预先指定的子段的最小长度、最大搜索步长和最大RMSE，RMSE的计算公式如下：

式中，L＝n_e–n_s+1为子段长度；||·||₂为2-范数运算；为原时域采样序列子段的向量；为用复指数级数拟合后的子段的向量；的元素按式(1)计算：

式中，以符号^标识拟合值；子段复指数的个数按式(14)确定，p_s＝int–[(n_e–n_s+1)/2]–1；其中n的起始值不是n_s，而是与式(1)中相同的0。

4.如权利要求3所述的用于电力系统信号处理的时域自适应分段复指数级数法，其特征在于，所述时域自适应分段复指数级数法的算法流程如下：

开始时，取子段长度为L_min，当数据点少e_f≥E_m时，保持n_Δ＝1不变并执行分支④；

在循环执行若干次分支④后，找到使e_f首次小于E_m的n_e，在满足e_f<E_m的条件下，取n_Δ＝N_Δm，经循环执行分支②搜索子段的最大长度，即子段的终点；

当再次出现e_f≥E_m时，转而执行分支③，采用函数int₄₅[n_Δ/2]将n_Δ从N_Δm逐步减半，确定出精确的分界点，此处符号int₄₅表示4舍5入取整；

当e_f<E_m再次满足时，又转回执行分支②，但n_Δ不再置为N_Δm；

分支②和③交替执行，当n_Δ被减为0，表明精确分界点已找到，则计算并存储当前子段的特性参数，完成后将n_Δ重置为1并调整n_s和n_e经分支①进行下一子段的搜索。

5.如权利要求1～4任意一项所述的用于电力系统信号处理的时域自适应分段复指数级数法，其特征在于，所述电力系统信号的插值方法为：

对于已知信号的数据点，将其直接拷贝至相应的时间点上；而对于已知信号数据点间的点，则由式(16)的如下变形形式产生：

上式在已知信号的2个整数序号间插入K个浮点数序号，与浮点数序号对应的数据点即为应插入的数据点。

6.如权利要求1～4任意一项所述的用于电力系统信号处理的时域自适应分段复指数级数法，其特征在于，所述电力系统信号的延拓方法为：

设N_F和N_B分别为正、负向延拓的点数；

负向延拓的公式为：

式中，参数为计算得到的第1子段的参数；计算前，应先提取第1子段的主导信号分量，p_s取主导信号分量的个数，去除虚假分量，用主导分量的参数进行计算；

正向延拓的公式为：

式中，n_last为最后一个数据点的序号，即最后一个子段的n_s。

7.实现权利要求书1～6任意一项所述的用于电力系统信号处理的时域自适应分段复指数级数法的计算机程序。