一种基于稀疏采集模型的谐波检测方法

1.本发明涉及电子信息领域，特别是涉及一种基于稀疏采集模型的谐波检测方法。


背景技术：

2.随着智能电网和新能源发电的迅速发展，大量非线性装置应用于电力系统，导致电网中谐波间谐波及ddc分量增加。一方面，谐波和间谐波对关注频率存在频谱干扰，影响了电能的准确计量。另一方面，针对电力系统发生故障时电信号中衰减的直流分量，传统算法无法对其进行滤波，导致计算结果误差较大。因此，需要寻找一种有效的相量估计算法来准确测量宽频动态相量。
3.目前，相量测量方法主要分为两类：dft算法（离散傅里叶变换）和非dft算法。dft算法是对有限时域离散信号的频谱进行等间隔采样，并对频域函数进行离散化处理。当采用传统的dft变换算法测量相量时，实际频率发生偏移会出现信息冗余、互谐波干扰、信息泄露等问题，导致产生较大的测量误差。一些研究提出引入过零检测法、小波变换法、瞬时值法等来优化dft算法，但该类方法在非标称频率下仍存在较大的误差。相关文献介绍了改进的二次谐波滤波技术在非标称频率下单相相量测量中的应用。该方法可以在不改变传统相量处理方法的结构下集成均匀采样和固定窗口长度相量测量，具有良好性能。然而，静态模型下基于dft的相量估计算法存在动态性能差的问题。
4.对于破解dtf模型仅局限于静态信号分析这一劣点，一些探究提议将动态相量估算拓展到不同的数学框架，以得到基于该模型的新一代相量测量方法。这类方法运用特定的数学模型来定义信号在稳定采样时间间隔内的振幅和相位变化，摆脱了傅里叶分析中周期的局限。
5.目前，动态谐波分析方法，如泰勒-傅立叶变换（tft），在电力振荡情况下可以实现相量估计，但易受到间谐波和高次谐波的干扰。相关文献提出了一种基于sinc插值函数的动态谐波同步算法，本发明拓宽了tft的频率范围，克服了dft周期性的限制。然而，该方法系统参数多且较为复杂，多次拟合后会存在误差大的问题，且测量误差会随着阶数的增加而增大。为了大幅提升估量精准度，相关文献提出具备高分辨率的离散傅里叶变换插值法。本发明具有恢复阶段信息的能力，通过阻断迭代次数达到噪声抑制。然而，该方法中dft序列仅能以额定速率重采样，这毋庸置疑增添了效果的不确定度。
6.基于上述各类方法都是通过高阶矩阵求逆运算来获得矩阵变换系数。尽管有很多的高效减小误差算法，但仍不能解决高阶间谐波和噪声叠加带来的影响。由此，考虑引用机器学习算法来降低相量求解计算复杂度。相关文献通过识别信号中最相关的分量，有效跟踪谐波进行建模估计，从而限制间谐波的影响。然而，该方法仅适用于低频段谐波向量，由于高次谐波信号具有时变性，且求解高阶伪逆矩阵具有高复杂性，无法直接获得相量估计参数。进而考虑动态谐波分析领域，相关文献利用stwls算法对谐波参数进行迭代估计，并补偿谐波参数对估计的影响。该方法可以提高频率测量精度，为动态信息分配更多空间，但仍难以抵抗约束间谐波的干扰。相关文献提出了o样条动态谐波分析，该方法以闭合形式降
低dttft的计算复杂程度, 为振荡提供了最佳的数据压缩算法。然而，在非理想条件下，由于空间步长重建过程中存在较大的噪声，随着阶次的升高，间谐波干扰会越来越严重。
7.为解决同步相量估计求解问题，有公开相关文献提出基于快速傅里叶变换(fft)算法的复小波变换，对原始波形进行分割，分析信号重构。然而，在连续多尺度分析中，这可能会导致严重的相位失真，并丢失下一步信号处理的时序信息。也有公开相关文献通过改进taylor加权最小二乘算法来获得谐波中的频率、幅值和相位估计。由于该方法只适用于测量三阶连续波而不能用于高精度谐波分析，所以本发明仍不适用于高阶谐波测量。


技术实现要素：

8.基于此，针对上述技术问题，本发明提出了一种基于稀疏采集模型的谐波检测方法，首先利用sinc插值函数对衰减直流分量进行参数化建模，通过最小二乘法得到准确的ddc分量，然后针对宽频谐波分量估计，采用改进的泰勒-傅立叶模型，结合机器学习算法，以较少的数据点准确恢复特定信号，进而将谐波相量重构问题转化为稀疏采集模型求解问题。在基于正则化框架下，本发明将相量重构问题表述为最小化函数的形式，利用sbi体系下更迭极小化运算有效地求解上述欠定逆问题，由此高效准确地获得谐波相量估计值。为了证明该模型的合理性，利用交叉熵目标函数来度量估计值的误差范围。仿真结果表明，本发明对宽频动态相量的估计具有较高精度和速度，在对比算法中保持一定优越性，同时验证了本发明的动态性能和抗干扰能力。
9.本发明的技术方案是：一种基于稀疏采集模型的谐波检测方法，包括以下步骤：a、宽频动态信号中衰减直流分量估算步骤：建立宽频动态信号模型，计算衰减直流分量；b、宽频动态信号中谐波相量估算步骤：建立谐波相量估计模型，重构谐波相量估计模型，求解谐波估计模型；c、仿真分析步骤：建立测试场景，依次建立基本性能测试场景、频率斜坡测试场景、阶跃变换测试场景和干扰测试场景。
10.步骤a中，“建立宽频动态信号模型”的步骤如下：对于包含直流分量和基本分量的电力系统波形表达式如(1)所示:
ꢀꢀꢀꢀ
(1)其中，是含有基波、谐波以及间谐波分量的信号，是时变振幅，是时变相位角，参考频率是，指最高次谐波分量，是衰减直流分量信号，是衰减直流分量幅值，是时间常数；含有基波、谐波以及间谐波分量信号在时刻t的动态相量模型为:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(2) 将(2)式代入(1)式可得宽频动态信号的动态表达式:
ꢀꢀꢀꢀ
(3)其中，为共轭算子。
11.步骤a中，“计算衰减直流分量”的步骤如下：衰减直流分量由振幅和相位时变的余弦分量和表示:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(4)其中，，，分别表示第n个余弦分量频率、幅值、相位，t为观测间隔长度，将y设置为3；基于频域采样定理，利用sinc插值函数对衰减直流分量进行时域信号参数化建模，ddc分量的动态相量记为，具体表示为:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(5)其中，是采样频率，表示最大模型阶数，将设置为2，为(2)式中第n个余弦分量在时刻上的样本相量；对衰减直流分量进行采样，式(1)中拟合的ddc分量的离散形式用(6)表示；
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(6)其中，表示相量为余弦信号的列向量，为共轭算子，矩阵的元素为在采样窗中的样本点，表示矩阵，的定义如(7)所示：
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(7)通过最小二乘法求出相量的估计值，计算过程如(7)所示，在矩阵中，为该矩阵元素：
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(8)其中，h是共轭转置符号；利用式(8)构建每个余弦分量的，接着用重构矩阵，通过(5)求得每个余弦信号的幅值和相位，进而估计ddc分量。
12.步骤b中，“建立谐波相量估计模型”的步骤如下：利用阶泰勒展开模型来定义动态相量:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(9)
其中，为t=0时的阶导数，k为泰勒展开阶数，t为观测间隔长度；包含谐波和间谐波的正弦信号采用离散序列的形式，序列长度是n，且，通过多频相量分析，将动态相量方法扩展到实际序列模型，离散信号各相量分量的泰勒傅里叶系数表达式如下:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(10)其中，h是采样间隔，观测间隔长度，是平均谐波及间谐波相量，是相量频率在 处的k阶导数，是基波频率的整数倍代表谐波频率，非整数倍时代表间谐波频率；以大小为的采样率将谐波信号频率归一化，采样长度是，得到谐波分量频率，泰勒傅里叶基向量的归一化频率是， ，对应n个系数集合的频率分辨率为，每个i的泰勒开阶数都为k，则(10)写成矩阵形式:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(11)其中，是样本列相量，大小为的传感矩阵a是泰勒傅里叶基相量，e表示噪声或测量误差，x是长度为的列相量，描述了当谐波时，的集合；相量求解，利用最小二乘法转换为伪逆矩阵的求解，通过伪逆矩阵计算得到的是欧几里得范数最小时对应的解，此时是满足下列约束条件下的最优解：
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(12)其中，表示欧几里得范数，进一步通过引入拉格朗日算子并求导后，得动态相量系数估计矩阵：
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(13)步骤b中，“重构谐波相量估计模型”的步骤如下：设定修正因子，每组频率分辨率达到，采样序列长度为，采样长度h为5个基波周期，谐波分量归一化频率为，，得到改进的泰勒傅里叶系数表达式：
ꢀꢀꢀꢀꢀꢀꢀꢀ
(14)其中，，表示大小为的曲波变换矩阵，是的k阶导数；根据表达式(9)，得到近似二阶动态相量估计值、谐波频率和频率变化率:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(15)
ꢀꢀꢀꢀꢀꢀꢀ
(16)
ꢀꢀꢀꢀꢀꢀꢀ
(17)当修正因子时，若传感矩阵a与稀疏矩阵高度不相关，(14)转化为矩阵形式即为:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(18)其中，是测量矩阵，x表示长度为的每个数据块对应的样本序列，是对原始样本分割后提取的少量相量的叠加，为维数为的正交矩阵，该矩阵的每一列是离散傅里叶变换的基相量，e表示噪声信号相量；上述方法实现稀疏采样测量模型构建后，利用欧式搜索算法对x进行矩阵分块，在搜索范围内找到m个最佳匹配的数据点，并将所有数据点组合成一个矩阵, i为系数矩阵的坐标，对矩阵进行曲波变换，采用曲波阈值准则去除不显著的曲波系数，抽取的系数需要有效的去噪和压缩算法，进行曲波变换，采用曲波阈值准则去除不显著的曲波系数，抽取的系数需要有效的去噪和压缩算法；稀疏矩阵只包含一个非零有效项集，该有效项集按字典顺序排列后所得的变换系数列相量是，变换系数集中在零区附近且以细峰的形式分布，使用变换系数的分布特点来表征信号的可重复性，记为:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(19)其中，表示l1范数；每组谐波信号估计值可根据变换系数倒置得:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(20)
其中，是的逆算子。
13.步骤b中，“求解谐波估计模型”的步骤如下：求解正则化线性最小二乘代价函数，重构泰勒系数矩阵x，即:
ꢀꢀꢀꢀꢀ
(21)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(22)其中，表示l2范数，表示正则化项，是正则化系数，、分别表示局部和非局部稀疏项，是一个正则化参数，用于平衡两个稀疏项； (22)代入(21)得:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(23)其中，，引入两个辅助相量p和q，最终得到以下方案:
ꢀꢀꢀꢀ
(24)
ꢀꢀꢀꢀꢀꢀꢀ
(25)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(26)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(27)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(28)其中，和为固定值参数，作用是提高算法数值的稳定性，p，q是sbi算法辅助迭代相量；(24)是求解一个严格凸二次函数的极小化，将目标函数的梯度设为0，得到对应的封闭解:
ꢀꢀꢀꢀꢀꢀ
(29)其中，是大小为的单位矩阵；采用步长最优的最速下降法求解函数的最小化问题，表示为:，
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(30)其中，为目标函数的梯度，是最佳步长，，估计值为，运用最小均方误差定理将(25)的问题变成:
ꢀꢀꢀꢀꢀꢀꢀ
(31)同理，(26)的问题写成:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(32)这里（省略）。
14.得到一组封闭解，将看作有噪声的观测值的某种类型，用表示误差相量，每个元素的误差分别为（），然后假设e独立分布于，均值为0且方差为，由大数定律得，对于任意，均有:
ꢀꢀꢀꢀꢀ
(33)所以得到:
ꢀꢀꢀꢀꢀ
(34)变换后的误差相量为，，，表示每个元素的误差相量，m是最佳匹配的数据点个数，根据矩阵的正交性质，变换不改变每个组的方差，由此得出，每个组的都独立分布于，（零均值且方差为），大数定理得，对于任意，均有:
ꢀꢀꢀꢀꢀ
(35)由以上得: (36)得出:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(37)将其和(26)问题合并，得:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(38)未知量是上式中是可分离的，每一个组成部分化简为:
ꢀꢀ
(39)其中，是指对应相量的模；
在矢量收缩形式下，x的估计值表达式:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(40)上述得到泰勒系数矩阵x的估计值后，通过谐波采样方程，对谐波信号的重组，在logistic回归的框架下，引入交叉熵[27]目标函数来度量估计值和理论值x的概率分布差异:
ꢀꢀ
(41)假设误差为二元分布，视为预测概率分布与实际概率分布的关联性接近，证明假设与预期模型是一致的。
[0015]
步骤c中，“建立基本性能测试场景”的步骤如下：测试所提算法的基本性能，构造式(43)所示含ddc分量的宽频动态信号模型：
ꢀꢀꢀꢀꢀꢀꢀꢀ
(42)其中，为基波频率，此处设置为50hz，、分别表示基波和各次谐波相角，且在范围内任意取值，低频段谐波次数h取值为2-13，高频段谐波次数h分别为77、79、80、83，采样频率设置为10 khz，ddc分量幅值的取值为0.3，0.4，0.5，
…
1，ddc分量的时间常数从0.01s开始，按步长为0.01s变化到0.1s；运用bmp算法，迭代相量初始值为，正则化系数为，固定值参数，， 时，ζ值的范围为[0.05,0.3]，搜索窗口中最佳匹配的数据点数设为m；选择fft，prony，twls，sife作为对比算法。
[0016]
步骤c中，“建立频率斜坡测试场景
”ꢀ
的步骤如下：分析bmp在频率斜坡下的性能，提供的信号如下所示：
ꢀꢀꢀꢀ
(43)其中，为基波频率，且取值为50hz，为基频斜率，取值为1hz/s，和分别为基波相位和谐波相位，相位设置为范围内均匀分布的随机数。
[0017]
步骤c中，“建立阶跃变换测试场景”的步骤如下：
在测试开始时，将各分量的幅值设置为初始幅值的115%，并且相位变化到，提供的宽频动态信号如下所示：
ꢀꢀꢀꢀꢀꢀꢀ
(44)设定采样率为5khz，采样周期长度为的5个周期，在标准中，用响应时间的快慢来评估各算法在阶跃变化条件下的性能。
[0018]
步骤c中，“建立干扰测试场景”的步骤如下：信号引入信噪比为55db的高斯白噪声，具体的宽频动态信号如下：
ꢀꢀꢀꢀ
(45)其中，为间谐波频率，并且其取值依次为9652.5 hz、9751.5 hz、9850.5 hz、9949.5 hz、10048.5 hz、10147.5 hz，和为基波和谐波的相位，且取值为范围内的随机数，设定采样频率为10khz，采样窗长设为5个工频周期。
[0019]
本发明的有益效果是：a、本发明分为两部分,第一部分基于最小二乘法得到准确的ddc分量，第二部分建立宽频动态相量的泰勒-傅立叶模型，利用谐波向量估计方法解决稀疏采集模型的正则化优化问题，最后利用分段sbi迭代框架得到谐波相量测量的估计值，实现对原始信号的重构，通过仿真试验和性能测试，显著提高了相量测量和估计的精度，能够为pmu测量提供可靠的理论依据；b、本发明首先利用sinc插值函数对衰减直流分量进行参数化建模，通过最小二乘法得到准确的ddc分量，然后针对宽频谐波分量估计，采用改进的泰勒-傅立叶模型，结合机器学习算法，以较少的数据点准确恢复特定信号，进而将谐波相量重构问题转化为稀疏采集模型求解问题，在基于正则化框架下，本发明将相量重构问题表述为最小化函数的形式，利用sbi体系下更迭极小化运算有效地求解上述欠定逆问题，由此高效准确地获得谐波相量估计值，为了证明该模型的合理性，利用交叉熵目标函数来度量估计值的误差范围，仿真结果表明，对宽频动态相量的估计具有较高精度和速度，在对比算法中保持一定优越性，同时验证了算法的动态性能和抗干扰能力；c、本发明提出的一种新的宽频动态相量测量算法，解决了包含ddc分量的宽频动态相量的有效处理问题；基于正则化稀疏捕获矩阵的思想，通过求解稀疏正则化问题得到相量估计结果；仿真和实验结果显示，识别了宽频动态相量的关键信息，显著降低了计算复杂度；这代表着能够在短时间内获取更精确的结果，从而有效地检测宽频信号的瞬态特性；同时在噪声干扰、间谐波干扰、频率斜坡等静态和动态条件下，测量结果能够满足m级pmu的测试要求。
附图说明
图1为重构效果及运行时间范围随参数变化情况图；图2 [a]、图2 [b]、图2 [c]为在频率斜坡条件下tve、fe、rfe的最大误差对比图，图2 [a]为在频率斜坡条件下tve的最大误差对比图, 图2 [b]为在频率斜坡条件下fe的最大误差对比图, 图2 [c]为在频率斜坡条件下rfe的最大误差对比图；图3为阶跃下各方法的运行时间比例图；图4 [a]、图4 [b]、图4 [c]为间谐波干扰下tve、fe、rfe的最大误差图, 图4 [a]为间谐波干扰下tve、的最大误差图, 图4 [b]为间谐波干扰下fe的最大误差图, 图4 [c]为间谐波干扰下rfe的最大误差图。
具体实施方式
[0020]
下面结合附图对本发明的实施例进行详细说明。
[0021]
实施例：1、宽频动态信号中衰减直流分量估计1.1、宽频动态信号模型对于包含直流分量和基本分量的电力系统波形表达式如(1)所示:
ꢀꢀꢀ
(1)其中是含有基波、谐波以及间谐波分量的信号。是时变振幅，是时变相位角，参考频率是，指最高次谐波分量。是衰减直流分量信号，是衰减直流分量幅值，是时间常数。
[0022]
含有基波、谐波以及间谐波分量信号在时刻t的动态相量模型为:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(2) 将(2)式代入(1)式可得宽频动态信号的动态表达式:
ꢀꢀꢀ
(3)其中，为共轭算子。
[0023]
1.2、衰减直流分量估计衰减直流分量可以由振幅和相位时变的余弦分量和表示:
ꢀꢀꢀꢀꢀꢀꢀꢀ
(4)
其中，，，分别表示第n个余弦分量频率、幅值、相位。t为观测间隔长度。为了满足精度要求，本实施例将y设置为3。
[0024]
基于频域采样定理，利用sinc插值函数对衰减直流分量进行时域信号参数化建模，这里ddc分量的动态相量记为，具体表示为:
ꢀꢀꢀꢀꢀ
(5)其中，是采样频率，表示最大模型阶数，本实施例将设置为2。为(2)式中第n个余弦分量在时刻上的样本相量。
[0025]
对衰减直流分量进行采样，式(2)中拟合的ddc分量的离散形式可用(5)表示。
[0026]
ꢀꢀꢀꢀ
(6)其中表示相量为余弦信号的列向量。为共轭算子。矩阵的元素为在采样窗中的样本点。表示矩阵，的定义如(7)所示：
ꢀꢀꢀ
(7)通过最小二乘法求出相量的估计值。计算过程如(9)所示。在矩阵中，为该矩阵元素。
[0027]
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(8)其中h是共轭转置符号。
[0028]
利用式(8)构建每个余弦分量的。接着用重构矩阵。通过(5)求得每个余弦信号的幅值和相位，进而估计ddc分量。
[0029]
2、宽频动态信号中谐波相量估计2.1、谐波相量估计模型的建立考虑到泰勒傅里叶方法能够描述观测间隔内振幅和相位随时间的变换。由此，利用阶泰勒展开模型来定义动态相量:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(9)其中，为t=0时的阶导数，k为泰勒展开阶数，t为观测间隔长度。
[0030]
在实际应用中，包含谐波和间谐波的正弦信号一般采用离散序列的形式，
序列长度是n，且。通过多频相量分析，将动态相量方法扩展到实际序列模型，离散信号各相量分量的泰勒傅里叶系数表达式:
ꢀꢀꢀꢀꢀꢀ
(10)其中，h是采样间隔，观测间隔长度。是平均谐波及间谐波相量，是相量频率在 处的k阶导数。是基波频率的整数倍代表谐波频率，非整数倍时代表间谐波频率。
[0031]
以大小为的采样率将谐波信号频率归一化，采样长度是，得到谐波分量频率，泰勒傅里叶基向量的归一化频率是， 。对应n个系数集合的频率分辨率为。考虑到每个i的泰勒开阶数都为k，则(11)可以写成矩阵形式:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(11)其中，是样本列相量。大小为的传感矩阵a是泰勒傅里叶基相量。e表示噪声或测量误差。x是长度为的列相量，描述了当谐波时，的集合。
[0032]
通常相量求解可利用最小二乘法转换为伪逆矩阵的求解。通过伪逆矩阵计算得到的是欧几里得范数最小时对应的解，此时是满足下列约束条件下的最优解：
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(12)其中，表示欧几里得范数。进一步通过引入拉格朗日算子并求导后，可得动态相量系数估计矩阵：
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(13)然而，随着最大谐波阶数的增加，矩阵维数迅速增加，求解方程(14)对应的多元线性方程所需的计算量显著增加。对于高频谐波，基于最小二乘法的相量法的计算复杂度最低。因此，这种相量法在求解多维矩阵时存在效率低的问题。此外，由于观测时间和频率分辨率相互矛盾，会导致频率分析不够精准，时间分辨率低。
[0033]
2.2、重构谐波相量估计模型为了获得更准确的频域结果，本实施例在谐波相量分析中引入修正因子。每组频
率分辨率达到，采样序列长度为，采样长度h为5个基波周期。谐波分量归一化频率可表示为，。
[0034]
在更精准的修正条件下得到改进的泰勒傅里叶系数表达式：
ꢀꢀꢀ
(14)这里，表示大小为的曲波变换矩阵，是的k阶导数。
[0035]
根据表达式(9)，可得到近似二阶动态相量估计值、谐波频率和频率变化率:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(15)
ꢀꢀꢀꢀꢀꢀꢀ
(16)
ꢀꢀꢀꢀꢀꢀ
(17)由于谐波未知量的数量远大于观测值，因此，假设x在curvelet变换域中是可压缩的。当修正因子时，若传感矩阵a与稀疏矩阵高度不相关，(14)转化为矩阵形式即为:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(18)其中是测量矩阵。x表示长度为的每个数据块对应的样本序列，是对原始样本分割后提取的少量相量的叠加。为维数为的正交矩阵，该矩阵的每一列是离散傅里叶变换的基相量。e表示噪声信号相量。
[0036]
上述方法实现稀疏采样测量模型构建后，本实施例基于谐波频域分布的稀疏性，利用欧式搜索算法对x进行矩阵分块。在搜索范围内找到m个最佳匹配的数据点，并将所有数据点组合成一个矩阵, i为系数矩阵的坐标。对矩阵进行曲波变换，采用曲波阈值准则去除不显著的曲波系数，抽取的系数需要有效的去噪和压缩算法。进行曲波变换，采用曲波阈值准则去除不显著的曲波系数，抽取的系数需要有效的去噪和压缩算法。
[0037]
为保证谐波检测的局部平滑特征，稀疏矩阵只包含一个很小的非零有效项集。该有效项集按字典顺序排列后所得的变换系数列相量是。由于变换系数集中在零
区附近且以细峰的形式分布，因此使用变换系数的分布特点来表征信号的可重复性，记为:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(19)其中表示l1范数。
[0038]
每组谐波信号估计值可根据变换系数倒置得:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(20)其中是的逆算子。
[0039]
2.3、谐波估计模型的求解为了获得具有合理精度和抗噪声鲁棒性的估计，基于上述方法建立谐波相量恢复模型后，下文通过求解正则化线性最小二乘代价函数的优化问题，来重构泰勒系数矩阵x，即:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(21)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(22)其中，表示l2范数。表示正则化项，是正则化系数。、分别表示局部和非局部稀疏项。是一个正则化参数，用于平衡两个稀疏项。
[0040] (22)代入(21)得:
ꢀꢀꢀꢀꢀꢀꢀꢀ
(23)其中。
[0041]
为了解决上述极小化问题，本实施例采用一种交替的sbi算法。引入两个辅助相量p和q，最终得到以下方案:
ꢀꢀ
(24)
ꢀꢀꢀꢀꢀ
(25)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(26)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(27)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(28)其中，和为固定值参数，作用是提高算法数值的稳定性。p，q是sbi算法辅
助迭代相量。
[0042]
公式(24)是求解一个严格凸二次函数的极小化的问题。将目标函数的梯度设为0，得到对应的封闭解:
ꢀꢀꢀ
(29)其中是大小为的单位矩阵。
[0043]
采用步长最优的最速下降法求解函数的最小化问题，可表示为:，
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(30)其中为目标函数的梯度，是最佳步长。，考虑到估计值，本实施例运用最小均方误差定理可将公式(25)的问题变成:
ꢀꢀꢀꢀꢀꢀꢀ
(31)同理，(27)的问题可以写成:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(32)这里（省略）。
[0044]
由于的定义繁琐，很难直观的解出上式。本实施例通过提出合理的假设，由此得到一组封闭解。将看作有噪声的观测值的某种类型，用表示误差相量，每个元素的误差分别为（），然后假设e独立分布于，均值为0且方差为，由大数定律得，对于任意，均有:
ꢀꢀꢀꢀꢀ
(33)所以得到:
ꢀꢀꢀꢀꢀ
(34)变换后的误差相量为，，表示每个元素的误差相量，m是最佳匹配的数据点个数。根据矩阵的正交性质，变换不改变每个组的方差。由此得出，每个组的都独立分布于，（零均值且方差为），大数定理得，对于任意，均有:
ꢀꢀꢀꢀꢀꢀꢀ
(35)由以上得: (36)不难得出:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(37)将其和公式(26)问题合并，得:
ꢀꢀꢀꢀꢀꢀꢀ
(38)由于未知量是上式中是可分离的，每一个组成部分化简为:
ꢀꢀꢀ
(39)其中是指对应相量的模。
[0045]
在矢量收缩形式下，x的估计值表达式:
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(40)上述得到泰勒系数矩阵x的估计值后，通过谐波采样方程，实现了对谐波信号的重组，达到精确检测谐波相量变化的目的。为了保证模型的准确性，在logistic回归的框架下，引入交叉熵目标函数来度量估计值和理论值x的概率分布差异:
ꢀꢀ
(41)假设误差为二元分布，可认为预测概率分布与实际概率分布的关联性非常接近，这证明了假设与预期模型是一致的。算法流程见下表：
3、仿真分析本节在不同测试场景下对本实施例所提出的宽频动态相量测量算法(bmp)结果展开了讨论。以下4种算法作为对比算法：快速傅里叶方法(fft)，prony算法，泰勒加权最小二乘(twls)，基于sinc插值函数的动态同步相量测量算法（sife）。测试场景包括基本性能测试、干扰测试、频率斜坡测试以及阶跃变换测试。
[0046]
为比较不同方法的估计效果，引入总相量误差(tve)来描述理论相量与估计相量之间的相对偏差。tve与幅值误差、相角误差关系紧密，但不能单独反映一个方面的变化。因此，本实施例还引入了两个指标：频率误差(fe)和绝对频率变化率误差(rfe)来综合评价相量估计的效果。本实施例根据ieeec37.118.1标准，选择m类pmu测量要求作为不同工况下的标准限值。
[0047]
5种算法都使用相同的矩形观测窗口，基频带宽为1hz，采样窗长设为5个工频周期。为提高频域分辨率，本实施例综合考虑算法复杂度和算法精度，将的值设置为20。
[0048]
4.1、基本性能测试为测试所提算法的基本性能，构造式(42)所示含ddc分量的宽频动态信号模型：
ꢀꢀꢀ
(42)其中，为基波频率，此处设置为50hz，、分别表示基波和各次谐波相角，且在范围内任意取值。低频段谐波次数h取值为2-13，高频段谐波次数h分别为77、79、80、83，采样频率设置为10 khz。ddc分量幅值的取值为0.3，0.4，0.5，
…
1。ddc分量的时间常数从0.01s开始，按步长为0.01s变化到0.1s。
[0049]
在运用bmp算法时，迭代相量初始值，正则化系数。固定值
参数，， 时，它们的ζ值通常在一个范围[0.05,0.3]内变化。搜索窗口中最佳匹配的数据点数设为m。在进行矩阵分块过程中，若匹配块数目过多，块数组中必然存在噪声影响大且匹配度低的数据点；反之，则不可避免偶然性对构造矩阵的影响。本实施例设最大迭代次数为j，迭代次数越高，计算精度越高，但计算成本显著提高。本实施例结合在参数m、j变化情况下，对重构效果及算法运行时间进行分析，结果如图1所示。
[0050]
图1是重构效果及运行时间范围随参数变化情况图。由图1可知，随着参数m和j的增大，总相量误差呈现减小趋势，但是算法运行时间增加。当数据点数目m=240、迭代次数j小于100时，改变迭代次数对计算成本影响较小(运行时间变化小)且此时重构效果趋于稳定。时，迭代次数的递增导致计算复杂度大大增加(运行时间变化大)，tve也呈现不稳定的下降趋势。因此，适当减少m能够有效提高重构效果。综合重构性能和运行时间两方面考虑，本实施例应用算法中参数m=240，迭代次数j=100。
[0051]
选择fft，prony，twls，sife作为对比算法，本实施例与对比算法的总相量误差估计、频率和频率变化率误差估计结果如下表所示，各算法估计精度比较见下表：
ꢀꢀ
从上表中可以看出，bmp算法的tve、fe和rfe指数的最大值分别为2.79%、0.096hz和2.437hz/s。当包含谐波分量时，tve和fe的ieee标准限值分别为3%和0.1hz，rfe的限值为2.7hz/s。与其他比较算法相比，本发明的估计指标完全满足ieee标准的要求。结果表明，本发明在含ddc分量的宽频谐波条件下仍具有良好的检测效果，相量估计精度最高。bmp使用谐波频域分布的稀疏分布来识别信号中最相关的分量，这显著提高了测量结果的准确性。
[0052]
sife方法在低次谐波附近的估计误差不符合ieee的测量标准。这是因为sife方法是基于基带信号滤波的低通滤波器，所以很难得到零误差结果。但由于其具有宽通带和宽阻带，可以很好地估计谐波同步相量，所以本发明性能优于fft、twls。fft测量结果受频谱泄漏影响最大这将大大降低谐波参数辨识准确度。最大总相量误差超过8%，其fe和rfe测量结果的准确性也不理想。在动态条件下，傅里叶变换模型无法跟踪观测窗口中的相量变化，导致相量评估不正确。twls方法采用二阶泰勒阶拟合信号分量。然而，泰勒信号模型误差大，精度有限。增加泰勒模型阶数可以减小模型误差，但阶数越高，滤波器的通带性能越差。prony算法使用参数模型计算信号参数。然而，其估计阶数将限制估计的频率分量的数量，并且频率估计误差将逐渐增大。
[0053]
4.2、频率斜坡测试负载和发电机之间的功率不平衡会导致宽频信号的频率随着负载的增大而减小，随着输入功率的增大而增大。为了分析bmp(本实施例算法)在频率斜坡下的性能，提供的信号如下所示：
ꢀꢀꢀꢀ
(42)其中，为基波频率，且取值为50hz。为基频斜率，且本实施例取值为1hz/s。和分别为基波相位和谐波相位，相位设置为范围内均匀分布的随机数。本实施例和对比算法的谐波相量估计、频率估计和频率变化率估计结果如图1所示。假定采样频率为10khz，并且采样窗长设为5个工频周期。采用fft，prony，twls，sife作为对比算法对
(42)的信号进行分析，则图2 [a]、图2 [b]、图2 [c]为各方法在频率斜坡条件下产生的估计结果，图2 [a]、图2 [b]、图2 [c]为在频率斜坡条件下tve、fe、rfe的最大误差对比图，图2 [a]为在频率斜坡条件下tve的最大误差对比图, 图2 [b]为在频率斜坡条件下fe的最大误差对比图, 图2 [c]为在频率斜坡条件下rfe的最大误差对比图。
[0054]
从图2 [a]、图2 [b]、图2 [c]可知，本实施例的tve最大误差为1.06%，频率误差最大为0.0095hz，频率变换率最大误差为0.103hz/s，均满足ieee标准要求[14]。结果表明，本发明在基频变化较大、频率线性变化时仍能保持较高的精度，其估计效果优于其他四种算法。bmp算法采用泰勒阶数逼近动态信号模型，能够准确估计频率和相位角，且受线性频率变化的影响最小。该方法实现了最高的相量估计精度。其余算法中，fft方法不能在动态条件下实时跟踪频率变化，因此频率变化率大。因为twls是基于动态模型的估计算法，prony算法能够准确地提取主导模态的低频振荡特征值，所以prony和twls算法的误差计算结果均小于fft，且算法相位角和频率等参数的误差特性受频率变化的影响较小，但测量结果仍难以满足ieee的要求。若对信号的时变特性不够了解，prony算法很难获得精确的模式参数。而由于基波分量的影响，stwls也无法准确估计高阶谐波相量。
[0055]
4.3、阶跃变换测试为了模拟电压/电流信号幅值和相位突然变化的故障条件，有必要在这些条件下对所提出的算法进行仿真，以评估响应时间和延迟。在测试开始时，将各分量的幅值设置为初始幅值的115%，并且相位变化到。提供的宽频动态信号如下所示： (43)本节测试假设采样率为5khz，采样周期长度为的5个周期。在标准中，用响应时间的快慢来评估各算法在阶跃变化条件下的性能。它被定义为在第一次和最后一次大于给定阈值的瞬间之间的时间间隔。根据本节测试条件下的ieee标准，其最大tve，fe和rfe值的阈值为1.5%，0.13hz，和0.78hz/s。具体的仿真结果如图3所示。
[0056]
从图3可以看出，与其他算法相比，bmp方法在振幅和相位阶跃变换条件下达到ieee标准的tve、fe和rfe值所需的时间更少。可以认为，在阶跃变换条件下，本发明的估计精度最高且响应时间完全满足p类pmu的标准要求。对于分步光滑函数的整个傅里叶变换的求值，fft算法需要进行非常复杂的乘法运算。twls性能可以通过在每个报告帧中重新计算其系数来提高，然而这大大增加了计算负担，需要计算伪逆。为了实现sife算法达到精度要求，需要在内存中存储一个70行1400列的复杂矩阵，需要实时处理次实际乘法和次实际加法，这样的内存和处理能力是非常有限的。为了让prony算法得到较为准确的分析结果，需要加大本发明的模型阶数，从而加大了计算量。
[0057]
4.4、抗干扰测试通常的情况下，电力系统信号中都会含有一定的间谐波和噪声，这严重影响到谐波相量的估计。本节测试中，信号引入信噪比为55db的高斯白噪声。具体的宽频动态信号如下：
ꢀ
(43)其中，为间谐波频率，并且其取值依次为9652.5hz、9751.5hz、9850.5hz、9949.5hz、10048.5hz、10147.5hz。和为基波和谐波的相位，且取值为范围内的随机数。本节假定采样频率为10khz，并且采样窗长设为5个工频周期。具体的仿真结果如图4 [a]、图4 [b]、图4 [c]所示，图4 [a]、图4 [b]、图4 [c]为间谐波干扰下tve、fe、rfe的最大误差图, 图4 [a]为间谐波干扰下tve、的最大误差图, 图4 [b]为间谐波干扰下fe的最大误差图, 图4 [c]为间谐波干扰下rfe的最大误差图。
[0058]
从图4 [a]、图4 [b]、图4 [c]可知，本实施例所提算法的tvemax为2.43%，femax为0.071hz，rfemax为0.184hz/s。本实施例算法各指标计算结果皆优于其余对比算法。在存在间谐波干扰条件下，tve、fe与rfe的ieee标准限值分别为3.5%、0.2hz、3hz/s。同理sife算法tvemax、femax、rfemax参数结果分别为6.64%，0.115hz，7.39hz/s，prony算法为19.73%，0.136hz，9.207hz/s，twls算法为26.35%，0.171hz，19.634hz/s，fft算法为49.76%，17.294hz，24.270hz/s，本发明各指标的计算结果均优于其他比较算法。bmp算法的各项指标均满足ieee标准的要求。在该段存在噪声干扰的情况下，fft方法的相邻谐波之间会发生严重的干扰和频谱泄漏，影响分辨率和精度。在本段噪声干扰的情况下，twls和sife方法会增加其谐波滤波器的过渡带振幅，引起相邻谐波之间的干扰，导致较大的估计误差。此外，twls算法受间谐波干扰影响严重，难以准确估计每个高频分量。对于高阶分量，prony算法的模型参数会随着谐波阶数的增加而不断修改。因此，它在间谐波相量和频率估计方面具有良好的效果，并能在一定程度上抑制间谐波分量频谱泄漏的影响。然而，这种方法需要预先估计动态时变信号的阶数。
[0059]
bmp算法采用改进的taylor-fourier模型和机器学习算法，以较少的数据点准确地恢复特定信号，然后将谐波相量重构以稀疏采集模型来求解，有效地提高了算法的重构性能和抗噪声能力。此外，本发明重建的频谱分辨率为1hz，有利于间谐波分量的准确检测。
[0060]
4、结论本实施例提出的一种基于稀疏采集模型的谐波检测方法，用来解决包含ddc分量的宽频动态相量的有效处理问题。本发明基于正则化稀疏捕获矩阵的思想，通过求解稀疏正则化问题得到相量估计结果。仿真和实验结果显示，该方法识别了宽频动态相量的关键信息，显著降低了了计算复杂度。这代表着能够在短时间内获取更精确的结果，从而有效地检测宽频信号的瞬态特性。同时，在噪声干扰、间谐波干扰、频率斜坡等静态和动态条件下，测量结果能够满足m级pmu的测试要求。
[0061]
以上所述实施例仅表达了本发明的具体实施方式，其描述较为具体和详细，但并不能因此而理解为对本发明专利范围的限制。应当指出的是，对于本领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干变形和改进，这些都属于本发明的保护范围。