CN103944174B - 基于互相关函数滤噪算法的低频振荡在线辨识方法 - Google Patents

Abstract

本发明公开了一种基于互相关函数滤噪算法的低频振荡在线辨识方法，其特点是以互相关函数代替实测信号，抑制色噪声，并结合总体最小二乘-旋转不变技术参数估计(TLS-ESPRIT)算法进行低频振荡的模态辨识。该算法能够快速、准确辨识出色噪声环境下电力系统低频振荡的模态，为系统的安全稳定分析及抑制措施提供有效依据。

Description

基于互相关函数滤噪算法的低频振荡在线辨识方法

技术领域

本发明涉及一种基于互相关函数滤噪算法的低频振荡在线辨识方法，具体地说，本发明是以互相关函数代替实测信号，从而抑制色噪声，并结合总体最小二乘-旋转不变技术参数估计(TLS-ESPRIT)算法进行低频振荡的模态辨识，属于电气信息领域。

背景技术

随着全国电网的互联，区域间低频振荡会直接影响区域间的功率交换，制约联络线传输能力。传统的低频振荡模态分析方法计算量大，只适合离线分析。随着广域测量系统(WAMS)，尤其是相量测量单元(PMU)在电力系统中的广泛应用，使低频振荡在线辨识得到实现。

近年来，快速傅里叶变换(FFT)、卡尔曼滤波以及希尔伯特—黄变换(HHT)等量测方法应用很广。FFT与卡尔曼滤波不能有效提取振荡的衰减特征。而HHT能提取低频振荡频率和衰减系数并且可自适应处理非线性、非平稳信号，但此算法采用经验模态分解(EMD)，复杂度高，不适合在线辨识；且EMD分解存在模态混叠现象，对于多模态存在的实际系统不太实用。

主流的模态辨识算法可分为两类：线性预测方法和子空间旋转不变方法。Prony以及改进的Prony算法属于线性预测方法。矩阵束(MP)算法以及各种ESPRIT算法属于子空间旋转不变方法。线性预测方法在低信噪比的情况下，性能会有所下降，而子空间旋转不变方法在运算效率和抗噪声性能方面更具优势。

电力系统的测量、传输环节通常都会引入噪声干扰，主要来源于系统负荷的随机变化，考虑采样误差及计算的截断误差，将这些误差看成高斯白噪声。PMU采集信号中的高斯白噪声经低通滤波处理后，会形成高度相关的高斯有色噪声。高斯有色噪声的存在，会产生估计偏差，从而影响低频振荡模态辨识的精度。目前，电力系统模态辨识一般都是通过自带的奇异值分解进行噪声的消除，这种情况对有色噪声的估计并不足，尤其对低信噪比情况下的噪声处理能力有限。

发明内容

本发明的目的是针对现有技术的不足而提供一种基于互相关函数滤噪算法的低频振荡在线辨识方法，其特点是采用互相关函数(CCF)来抑制高斯色噪声的影响，在此基础上通过总体最小二乘-旋转不变技术参数估计(TLS-ESPRIT)算法来在线辨识出低频振荡模态信息，以下统称为CCF-TLS-ESPRIT。

本发明的目的由以下技术措施实现

基于互相关函数滤噪算法的低频振荡在线辨识方法包括以下步骤：

1)去除实测信号的均值，保留振荡分量，电力系统低频振荡信号模型为：

式中，Q为假设的信号模式数，B_k为振幅，α_k为阻尼因子，f_k为频率，为初始相位。进一步表示为指数模型：

$X\left(t\right)=\underset{k=1}{\overset{k=2Q}{\Σ}}{R}_k{e}^{p_{k}t}---\left(2\right)$

式中，p_k＝2πf_kj±α_k， 为低频振荡信号的相位。

由此得两个不同时刻的有限长信号指数模型为：

$x\left(n\right)=\underset{i=1}{\overset{n}{\Σ}}{a}_i{e}^{p_{i}n},(n=0,1,2,...N_1-1)---\left(3\right)$

$y\left(m\right)=\underset{j=1}{\overset{m}{\Σ}}{b}_j{e}^{p_{j}m},(m=0,1,2,...N_2-1)---\left(4\right)$

式中，a_i、b_j分别是相应指数模型的幅值。对信号的不同采样序列，求取其互相关函数R_xy(τ)；

平稳采样序列{x(n)}和{y(m)}受噪声污染后的非平稳混合过程表示为：

$\tilde{x}\left(n\right)=x\left(n\right)+w\left(n\right)=\underset{i=1}{\overset{m}{\Σ}}{a}_i{e}^{p_{i}n}+w\left(n\right)---\left(5\right)$

$\tilde{y}\left(m\right)=y\left(m\right)+v\left(m\right)=\underset{i=1}{\overset{m}{\Σ}}{b}_j{e}^{p_{j}m}+v\left(m\right)---\left(6\right)$

其中，w(n)和v(m)是相互独立的噪声，高斯白噪声或者高斯有色噪声；对电力系统，其统计分布为零均值。

2)针对有限长序列{x(n)}和{y(n)}，得其互相关函数表达式为：

$R_{xy}\left(\mathrm{\τ}\right)=\frac{1}{N}\underset{n=S_1}{\overset{S_2}{\Σ}}x\left(n\right)y^{*}(n+\mathrm{\τ})---\left(7\right)$

类似于四阶自相关函数的协方差估计子，取S₁＝max(0，-τ)，

S₂＝min(N-1-τ_max，N-1+τ_max-τ)，N＝min(N₁，N₂)，0＜τ≤τ_max；

将式(5)(6)带入式(7)，得：

${\tilde{R}}_{xy}\left(\mathrm{\τ}\right)=R_{xy}\left(\mathrm{\τ}\right)+\mathrm{\ρ}\left(\mathrm{\τ}\right)+R_{wv}\left(t\right)---\left(8\right)$

其中，为含噪声的互相关函数，R_wv(t)＝0，

$\begin{array}{c}\mathrm{\ρ}\left(\mathrm{\τ}\right)=\frac{1}{N}\underset{n=S_1}{\overset{S_2}{\Σ}}x\left(n\right)w^{*}(m+\mathrm{\τ})+\frac{1}{N}\underset{n=S_1}{\overset{S_2}{\Σ}}v\left(n\right)y^{*}(m+\mathrm{\τ})\\ =E\[x\left(n\right)w^{*}(m+\mathrm{\τ})\]+E\[v\left(n\right)y^{*}(m+\mathrm{\τ})\]\\ \≈E\[x\left(n\right)\]E\[w^{*}(m+\mathrm{\τ})\]+E\[v\left(n\right)\]E\[y^{*}(m+\mathrm{\τ})\]\\ =0\end{array}---\left(9\right)$

故 ${\tilde{R}}_{xy}\left(\mathrm{\τ}\right)\≈R_{xy}\left(\mathrm{\τ}\right)$

振荡信号中的高斯噪声能够被有效抑制，互相关序列保留了原信号的极点信息；

3)利用互相关函数R_xy(τ)形成Hankel矩阵为：

τ_max、N取值为2)中所示，L＝N/3～N/2；

4)对矩阵Y进行奇异值分解，取出其前n个奇异值：σ₁、σ₂、…σ_n，按照奇异值范数得到有效秩，生成信号子空间V_S和噪声子空间V_N；

5)V_S删除第一行和第二行剩下的矩阵分别为V₁、V₂，由V₂＝V₁Ψ得到旋转算子Ψ；对[V₁,V₂]进行奇异值分解得到右特征向量P， $P=\left[\begin{array}{cc}{P}_{11}& P_{12}\\ P_{21}& P_{22}\end{array}\right];$

6)计算的特征根得到信号的频率、衰减系数；

7)利用总体最小二乘法求出信号的幅值和主导模态，通过ESPRIT算法使实测数据组成的向量经旋转后得到新的向量，保持两种向量对应的信号子空间的不变性，再通过求旋转算子的广义特征值来求取信号极点，应用总体最小二乘法求出信号幅值，得到主导模态。

本发明具有如下优点：

本发明采用互相关函数代替电力系统实测有限长信号，使其对高斯色噪声不敏感，能够有效滤除信号中的高斯色噪声，该算法收敛性较好，运算精度高，运算速度快，具有很高的工程实用价值。

附图说明

图1为联络线AC7功率信号图

图2为三种算法运算时间比较图

图3为发电机主振模态图

图4为配置PSS前后信号图

具体实施方式

下面通过实施例对本发明进行具体的描述，有必要在此指出的是本实施例只用于对本发明进行进一步说明，不能理解为对本发明包括范围的限制，该领域的技术熟练人员可以根据上述本发明的内容做出一些非本质的改进和调整。

实施例

选取四机两区域系统作为电力系统仿真算例，扰动设为交流线路AC5三相短路，持续0.05s。取本系统两区域间的交流联络线(AC7)的功率振荡信号作为辨识信号。

首先对交流联络线AC7功率信号进行无噪声辨识，为便于比较CCF与FOMMC算法的性能，分别采用TLS-ESPRIT、FOMMC-TLS-ESPRIT、CCF-TLS-ESPRIT算法进行辨识，取幅值最大的模态，如表1所示。

表1不加噪声辨识结果

由表1得知，此系统有一个弱阻尼的区间振荡主振模态，除此以外，还有一个0.96Hz的弱阻尼非主振模态。在无噪声干扰下，三种算法辨识的精度差别不大。

对联络线AC7的功率信号叠加15dB的高斯白噪声，并经过第3节所示的低通滤波器进行滤波，辨识时间窗为10s。原始信号、含噪信号以及CCF-TLS-ESPRIT算法的拟合信号如图1所示。

由图1，CCF-TLS-ESPRIT能够在高斯色噪声环境下的采集信号中准确拟合出原始信号的模态信息。拟合信号相对原始信号，除幅值减小了若干倍以外，其余信息均未发生改变。

采用TLS-ESPRIT、FOMMC-TLS-ESPRIT、CCF-TLS-ESPRIT算法分别对含有高斯色噪声的功率信号进行辨识，结果如表2所示。三种算法的运算时间如图2所示。

表2含高斯色噪声辨识结果

由表2得知，TLS-ESPRIT算法受色噪声影响较大，色噪声下辨识结果与无噪声时相差较多，尤其是非主振模态的辨识存在较大误差。FOMMC-TLS-ESPRIT与CCF-TLS-ESPRIT算法受色噪声影响很小，辨识结果与无噪声时接近。结合图2可看出，TLS-ESPRIT虽然运算时间最短，但其在色噪声环境下精度不高，不如其余两种算法。CCF-TLS-ESPRIT算法运算时间比FOMMC-TLS-ESPRIT算法短，具有运算效率的优势。

要提高系统阻尼，首先需抑制系统的主振模态振荡，取母线BUS8为相角参考点，取各发电机母线与其相对相角为CCF-TLS-ESPRIT的辨识信号，生成各发电机主振模态图，如图3所示。

由图3得知，两区域内的发电机G1、G2与G3、G4振形相反、振幅接近。考虑在每台发电机上配置电力系统稳定器(PSS)来抑制系统的弱阻尼区间振荡，PSS的设计方法采用传统的留数法，效果如图4所示。

结果表明：在有效在线辨识出系统低频振荡模态的基础上，合理配置PSS能有效提高系统阻尼。

Claims (1)

1.基于互相关函数滤噪算法的低频振荡在线辨识方法，其特征在于该方法包括以下步骤：

1)去除实测信号的均值，保留振荡分量，电力系统低频振荡信号模型为：

式中，Q为假设的信号模式数，B_k为振幅，α_k为阻尼因子，f_k为频率，为初始相位，进一步表示为指数模型：

$X\left(t\right)=\underset{k=1}{\overset{k=2Q}{\Σ}}{R}_k{e}^{p_{k}t}---\left(2\right)$

式中，p_k＝2πf_kj±α_k，为低频振荡信号的相位，由此得两个不同时刻的有限长信号指数模型为：

$x\left(n\right)=\underset{i=1}{\overset{n}{\Σ}}{a}_i{e}^{p_{i}n},(n=0,1,2,...N_1-1)---\left(3\right)$

$y\left(m\right)=\underset{j=1}{\overset{m}{\Σ}}{b}_j{e}^{p_{j}m},(m=0,1,2,...N_2-1)---\left(4\right)$

式中，a_i、b_j分别是相应指数模型的幅值，对信号的不同采样序列，求取其互相关函数R_xy(τ)；

平稳采样序列{x(n)}和{y(m)}受噪声污染后的非平稳混合过程表示为：

$\tilde{x}\left(n\right)=x\left(n\right)+w\left(n\right)=\underset{i=1}{\overset{m}{\Σ}}{a}_i{e}^{p_{i}n}+w\left(n\right)---\left(5\right)$

$\tilde{y}\left(m\right)=y\left(m\right)+v\left(m\right)=\underset{i=1}{\overset{m}{\Σ}}{b}_j{e}^{p_{j}m}+v\left(m\right)---\left(6\right)$

其中，w(n)和v(m)是相互独立的噪声，高斯白噪声或者高斯有色噪声；

2)针对有限长序列{x(n)}和{y(n)}，得其互相关函数表达式为：

$R_{xy}\left(\mathrm{\τ}\right)=\frac{1}{N}\underset{n=S_1}{\overset{S_2}{\Σ}}x\left(n\right)y^{*}(n+\mathrm{\τ})---\left(7\right)$

类似于四阶自相关函数的协方差估计子，取S₁＝max(0，-τ)，

S₂＝min(N-1-τ_max，N-1+τ_max-τ)，N＝min(N₁，N₂)，0＜τ≤τ_max；

将式(5)(6)带入式(7)，得：

${\tilde{R}}_{xy}\left(\mathrm{\τ}\right)=R_{xy}\left(\mathrm{\τ}\right)+\mathrm{\ρ}\left(\mathrm{\τ}\right)+R_{wv}\left(t\right)---\left(8\right)$

其中，为含噪声的互相关函数，R_wv(t)＝0，

$\begin{array}{c}\mathrm{\ρ}\left(\mathrm{\τ}\right)=\frac{1}{N}\underset{n=S_1}{\overset{S_2}{\Σ}}x\left(n\right)w^{*}(m+\mathrm{\τ})+\frac{1}{N}\underset{n=S_1}{\overset{S_2}{\Σ}}v\left(n\right)y^{*}(m+\mathrm{\τ})\\ =E\[x\left(n\right)w^{*}\left(m+\mathrm{\τ}\right)\]+E\[v\left(n\right)y^{*}\left(m+\mathrm{\τ}\right)\]\\ \≈E\[x\left(n\right)\]E\[w^{*}\left(m+\mathrm{\τ}\right)\]+E\[v\left(n\right)\]E\[y^{*}\left(m+\mathrm{\τ}\right)\]\\ =0\end{array}$

故 ${\tilde{R}}_{xy}\left(\mathrm{\τ}\right)\≈R_{xy}\left(\mathrm{\τ}\right)---\left(9\right)$

振荡信号中的高斯噪声能够被有效抑制，互相关序列保留了原信号的极点信息；

3)利用互相关函数R_xy(τ)形成Hankel矩阵为：

τ_max、N取值为2)中所示，L＝N/3～N/2；

4)对矩阵Y进行奇异值分解，取出其前n个奇异值：σ₁、σ₂、…σ_n，按照奇异值范数得到有效秩，生成信号子空间V_S和噪声子空间V_N；

5)V_S删除第一行和第二行剩下的矩阵分别为V₁、V₂，由V₂＝V₁Ψ得到旋转算子Ψ；对[V₁,V₂]进行奇异值分解得到右特征向量P， $P=\left[\begin{array}{cc}{P}_{11}& P_{12}\\ P_{21}& P_{22}\end{array}\right];$

6)计算的特征根得到信号的频率、衰减系数；

7)利用总体最小二乘法求出信号的幅值和主导模态，通过ESPRIT算法使实测数据组成的向量经旋转后得到新的向量，保持两种向量对应的信号子空间的不变性，再通过求旋转算子的广义特征值来求取信号极点，应用总体最小二乘法求出信号幅值，得到主导模态。