一种基于数学形态学的电力信号频率测量方法与流程

1.本发明涉及电网交流信号频率测量的技术领域，尤其是指一种基于数学形态学的电力信号频率测量方法。


背景技术：

2.在电力信号分析中经常需要测量电力系统的频率及其偏差。现有的基波频率测量方法中，常采用整周期采样点计数法。该方法首先确定基波的周期，在确定基波周期后，再根据基波周期内采样点之间的时间间隔确定基波的频率。该方法在采样周期较少时，所测量的频率精度很差，而在采样周期较长时，又不能满足实时性的要求。此外，在电力系统发生扰动时，采用该方法测量频率时，过零点会存在畸变，从而影响整个基波周期的确定，最终导致所测量的频率不准确。基于傅里叶变换的方法需要足够长的数据窗口，才能足够小的频率分辨率，以满足对频率测量精度的要求。
3.在信号处理领域，数学形态学对于滤波、信号分解、梯度提取、峰谷提取等方面的操作都极其有效。与傅里叶变换和小波变换等其它信号处理方法比较，数学形态学由于只有加减和比较运算的参与，计算量较小。形态学滤波器已广泛利用于心电图信号的基线校正及噪音抑制，滚动轴承缺陷特征提取的噪声抑制等。在形态学的发展中，有学者提出了软形态学的概念，将结构元素分为内核和软边缘两部分，增强了形态学滤波器在噪声环境下的表现。其它形态算子，如形态学梯度，形态小波和多分辨率形态学也被用于心电图，变压器励磁涌流，输电线路故障等场景的奇异性检测。


技术实现要素：

4.本发明的目的在于克服现有技术的缺点与不足，提出了一种基于数学形态学的电力信号频率测量方法，能够在稳态和非稳态情况下，提高信号频率测量的精度，且该方法的计算量小、实时性好。
5.为实现上述目的，本发明所提供的技术方案为：一种基于数学形态学的电力信号频率测量方法，包括以下步骤：
6.1)根据数学形态学运算准则，采用软结构元素，构建数学形态学滤波器对信号进行滤波降噪；
7.2)利用数学形态学的高帽变换和低帽变换对滤波后的信号进行处理，并提取处理结果的零点；
8.3)在每一个周期中，搜索经高帽变换结果的最左侧零点和低帽变换的最右侧零点，求取两者平均值，作为该周期的频率计算参考点；
9.4)根据相邻周期的频率计算参考点之间的时间差继而获取信号的频率测量值。
10.进一步，在步骤1)中，数学形态学的基本运算包括腐蚀运算和膨胀运算，设f(x)是一维待处理的信号，g(s)是用于提取信号特征的结构元素，则利用结构元素g(s)对信号f(x)进行腐蚀运算和膨胀运算，分别表示为：
[0011][0012][0013]
式中，表示腐蚀运算，表示膨胀运算，表示利用结构元素g(s)对待处理信号f(x)作腐蚀运算，表示利用结构元素g(s)对待处理信号f(x)作膨胀运算，df和dg分别为待处理信号f(x)和结构元素g(s)的定义域，x和s分别表示f(x)和g(s)的自变量，都是离散点，x∈df表示x在f(x)的定义域内，s∈dg表示s在g(s)的定义域内，f(x s)表示自变量由x向右移s时待处理信号的值，同理f(x-s)表示自变量由x向左移s时待处理信号的值，移动的范围由结构元素的定义域决定；
[0014]
所述软结构元素，即将结构元素分为两个部分，其中一部分是核b1，其对应权重大于或等于1，另一个部分是软边缘b2，其对应权重等于1；根据数学形态学的腐蚀运算和膨胀运算准则，软形态学腐蚀运算和膨胀运算分别表示为：
[0015][0016][0017]
式中，表示利用软结构元素b1和b2对待处理信号f(x)作软形态学腐蚀运算，表示利用软结构元素b1和b2对待处理信号f(x)作软形态学膨胀运算，k为核b1所对应的权重，kth为迭代次数，核b1所对应的权重与迭代次数相等，y和z分别表示核b1和软边缘b2的自变量，都是离散点，y∈b1表示y在b1的定义域内，z∈b2表示z在b2的定义域内,∪是集合并运算，k
◇
f(x)表示由k个f(x)形成的集合，即k
◇
f(x)＝{f(x),f(x),
…
,f(x)}；
[0018]
当结构元素为扁平结构时，即结构元素的值域为0，数学形态学的基本运算实际上是在长度为结构元素长度的数据窗口内求取极值；对于腐蚀运算是取最小值，对于膨胀运算则是取最大值，因此，数学形态学是一种非线性运算，膨胀运算与腐蚀运算不是互为逆运算，将这两种基本运算级联使用而得到开运算与闭运算；开运算是先腐蚀后膨胀的运算，对信号有削减波峰的作用；闭运算则是先膨胀后腐蚀的运算，对信号有填充波谷的作用；所述数学形态学滤波器采用软结构元素，待处理信号进行开、闭运算，并取均值，表示为：
[0019][0020][0021][0022]
式中，表示开运算，
·
表示闭运算，表示利用软结构元素b1和b2对待处理信号f(x)作开运算，f(x)
·
[b1,b2,k]表示利用软结构元素b1和b2对待处理信号f(x)作闭运算，soft_filter[f(x)]表示利用具有软结构元素的数学形态学滤波器对信号f(x)进行滤波降噪处理。
[0023]
进一步，在步骤2)中，利用数学形态学的高帽top-hat变换和低帽bottom-hat变换对滤波后的信号进行处理，分别提取信号的波峰和波谷部分；top-hat算子用于提取信号的波峰部分，表示为：
[0024]
th[f(x)]＝f(x)-f(x)οg(s)
[0025]
而bottom-hat算子则是用于提取信号的波谷部分，表示为：
[0026]
bh[f(x)]＝f(x)-f(x)
·
g(s)
[0027]
式中，th表示高帽变换，bh表示低帽变换，f(x)是一维待处理的信号，th[f(x)]表示对待处理信号f(x)作高帽变换，bh[f(x)]表示对待处理信号f(x)作低帽变换，g(s)是用于提取信号特征的结构元素，长度取信号额定周期的二分之一；
[0028]
在每一个周期中，高帽变换和低帽变换的结果均有一个零点区间，提取处理结果的零点，提取的零点与信号的频率存在以下关系：在同一个周期中，当信号的真实频率等于额定频率时，高帽变换处理结果的最左侧的零点和低帽变换最右侧的零点重合；当真实频率小于额定频率时，则高帽变换的零点被低帽变换的零点所包含；当真实频率大于额定频率时，低帽变换的零点被高帽变换的零点所包含。
[0029]
进一步，在步骤3)中，所述频率计算参考点的计算方法根据信号的极性，既能选择同属于上升沿的参考点，也能选择同属于下降沿的参考点；当选择同属于上升沿的参考点时，取每一个周期内高帽变换处理结果的最左侧的零点和低帽变换最右侧的零点的平均值，表示为：
[0030][0031]
式中，n
up,left,i
表示以上升沿为参考时，在第i个周期的高帽变换处理结果的最左侧的零点所对应的采样序列编号；n
up,right,i
表示以上升沿为参考时，在第i个周期的低帽变换处理结果的最右侧的零点所对应的编号；n
up,ref,i
表示以上升沿为参考时，在第i个周期的频率计算参考点所对应的编号；
[0032]
同理，当选择同属于下降的参考点时，取每一个周期内高帽变换处理结果的最右侧的零点和低帽变换最左侧的零点的平均值来求取频率计算参考点，表示为：
[0033][0034]
式中，n
down,left,i
表示以下降沿为参考时，在第i个周期的低帽变换处理结果的最左侧的零点所对应的采样序列编号；n
down,right,i
表示以下降沿为参考时，在第i个周期的高帽变换处理结果的最右侧的零点所对应的编号；n
down,ref,i
表示以下降沿为参考时，在第i个周期的频率计算参考点所对应的编号；
[0035]
因此，第i个周期的频率计算参考点表示为：
[0036][0037]
式中，n
ref,i
表示在第i个周期的频率计算参考点所对应的编号。
[0038]
进一步，在步骤4)中，选择两个相邻周期的频率计算参考点，并计算其时间差从而获得频率测量值，表示为：
[0039][0040]
式中，f
mea,i
为第i个周期的频率测量值，fs为信号的采样频率，n
ref,i
和n
ref,i-1
分别表示第i个周期和第i-1个周期的频率计算参考点所对应的编号。
[0041]
本发明与现有技术相比，具有如下优点与有益效果：
[0042]
1、本发明方法首次提出了基于数学形态学的频率测量方法，采用软结构元素的数学形态学滤波器滤除噪声的干扰，响应速度快，无相位偏。
[0043]
2、本发明方法首次提出了基于数学形态学的频率测量方法，采用形态学高帽变换和低帽变换对波形进行处理，根据输出结果的零点的关系来确定频率计算参考点，克服了传统过零点检测方法因受到干扰产生过零点漂移导致频率测量不准确的难题。
[0044]
3、本发明方法所需计算量小，不受直流分量的影响，适用于非稳态情况下的信号频率测量。
附图说明
[0045]
图1为形态滤波前后的信号波形对比图。
[0046]
图2为当稳态信号频率等于额定频率(50hz)时的形态学处理结果图。
[0047]
图3为当稳态信号频率(40hz)小于额定频率(50hz)时的形态学处理结果图。
[0048]
图4为当稳态信号频率(60hz)小于额定频率(50hz)时的形态学处理结果图。
[0049]
图5为当非稳态情况下信号的形态学处理结果图。
具体实施方式
[0050]
下面结合实施例及附图对本发明作进一步详细的描述，但本发明的实施方式不限于此。
[0051]
本实施例公开了一种基于数学形态学的电力信号频率测量方法，使用了软结构元素构建数学形态学滤波器，并对滤波后的信号进行高帽变换和低帽变换处理，其包括以下步骤：
[0052]
1)输入待处理信号，根据数学形态学运算准则，采用软结构元素，构建数学形态学滤波器对信号进行滤波降噪。
[0053]
数学形态学的基本运算包括腐蚀运算和膨胀运算，设f(x)是一维待处理的信号，g(s)是用于提取信号特征的结构元素，则利用结构元素g(s)对信号f(x)进行腐蚀运算和膨胀运算，可分别表示为：
[0054][0055][0056]
式中，表示腐蚀运算，表示膨胀运算，表示利用结构元素g(s)对待处理信号f(x)作腐蚀运算，表示利用结构元素g(s)对待处理信号f(x)作膨胀运算，df和dg分别为待处理信号f(x)和结构元素g(s)的定义域，x和s分别表示f(x)和g(s)的自变量，都是离散点，x∈df表示x在f(x)的定义域内，s∈dg表示s在g(s)的定义域内，f(x s)表
示自变量由x向右移s时待处理信号的值，同理f(x-s)表示自变量由x向左移s时待处理信号的值，移动的范围由结构元素的定义域决定；
[0057]
所述软结构元素，即将结构元素分为两个部分，其中一部分是核b1，其对应权重大于或等于1，另一个部分是软边缘b2，其对应权重等于1。根据数学形态学的腐蚀运算和膨胀运算准则，软形态学腐蚀运算和膨胀运算分别表示为：
[0058][0059][0060]
式中，表示利用软结构元素b1和b2对待处理信号f(x)作软形态学腐蚀运算，表示利用软结构元素b1和b2对待处理信号f(x)作软形态学膨胀运算，k为核b1所对应的权重，kth为迭代次数，核b1所对应的权重与迭代次数相等，y和z分别表示核b1和软边缘b2的自变量，都是离散点，y∈b1表示y在b1的定义域内，z∈b2表示z在b2的定义域内,∪是集合并运算，k
◇
f(x)表示由k个f(x)形成的集合，即k
◇
f(x)＝{f(x),f(x),
…
,f(x)}；
[0061]
当结构元素为扁平结构时，即结构元素的值域为0，数学形态学的基本运算实际上是在长度为结构元素长度的数据窗口内求取极值。对于腐蚀运算是取最小值，对于膨胀运算则是取最大值，因此，数学形态学是一种非线性运算，膨胀运算与腐蚀运算不是互为逆运算。将这两种基本运算级联使用而得到开运算与闭运算。开运算是先腐蚀后膨胀的运算，对信号有削减波峰的作用。闭运算则是先膨胀后腐蚀的运算，对信号有填充波谷的作用。所述数学形态学滤波器采用软结构元素，待处理信号进行开、闭运算，并取均值，表示为：
[0062][0063][0064][0065]
式中，表示开运算，
·
表示闭运算，表示利用软结构元素b1和b2对待处理信号f(x)作开运算，f(x)
·
[b1,b2,k]表示利用软结构元素b1和b2对待处理信号f(x)作闭运算，soft_filter[f(x)]表示利用具有软结构元素的数学形态学滤波器对信号f(x)进行滤波降噪处理。
[0066]
案例1中待测量信号为添加了20db白噪声的正弦信号，信号的频率为50hz，采样频率为10khz。进行数学形态学滤波时，所用的软结构元素为扁平结构元素，b1长度为5，b2长度为2。由图1可见，滤波前信号有许多毛刺，而滤波后波形变得很光滑，说明具有软结构元素的数学形态学滤波器对含噪的信号滤波效果良好。
[0067]
2)利用数学形态学的高帽变换和低帽变换对滤波后的信号进行处理，并提取处理结果的零点。
[0068]
利用数学形态学的高帽(top-hat)变换和低帽(bottom-hat)变换对滤波后的信号进行处理，分别提取信号的波峰和波谷部分。top-hat算子用于提取信号的波峰部分，表示
为：
[0069][0070]
而bottom-hat算子则是用于提取信号的波谷部分，表示为：
[0071]
bh[f(x)]＝f(x)-f(x)
·
g(s)
[0072]
式中，th表示高帽变换，bh表示低帽变换，f(x)是一维待处理的信号，th[f(x)]表示对待处理信号f(x)作高帽变换，bh[f(x)]表示对待处理信号f(x)作低帽变换，g(s)是用于提取信号特征的结构元素，长度取信号额定周期的二分之一；
[0073]
在每一个周期中，高帽变换和低帽变换的结果均有一个零点区间，提取处理结果的零点。提取的零点与信号的频率存在以下关系：在同一个周期中，当信号的真实频率等于额定频率时，高帽变换处理结果的最左侧的零点和低帽变换最右侧的零点重合；当真实频率小于额定频率时，则高帽变换的零点被低帽变换的零点所包含；当真实频率大于额定频率时，低帽变换的零点被高帽变换的零点所包含。
[0074]
3)在每一个周期中，搜索经高帽变换结果的最左侧零点和低帽变换的最右侧零点，求取两者平均值，作为该周期的频率计算参考点。
[0075]
频率计算参考点的计算方法根据信号的极性，既能选择同属于上升沿的参考点，也能选择同属于下降沿的参考点。当选择同属于上升沿的参考点时，取每一个周期内高帽变换处理结果的最左侧的零点和低帽变换最右侧的零点的平均值，表示为：
[0076][0077]
式中，n
up,left,i
表示以上升沿为参考时，在第i个周期的高帽变换处理结果的最左侧的零点所对应的采样序列编号；n
up,right,i
表示以上升沿为参考时，在第i个周期的低帽变换处理结果的最右侧的零点所对应的编号；n
up,ref,i
表示以上升沿为参考时，在第i个周期的频率计算参考点所对应的编号；
[0078]
同理，当选择同属于下降的参考点时，取每一个周期内高帽变换处理结果的最右侧的零点和低帽变换最左侧的零点的平均值来计算频率计算参考点，表示为：
[0079]
式中，n
down,left,i
表示以下降沿为参考时，在第i个周期的低帽变换处理结果的最左侧的零点所对应的采样序列编号；n
down,right,i
表示以下降沿为参考时，在第i个周期的高帽变换处理结果的最右侧的零点所对应的编号；n
down,ref,i
表示以下降沿为参考时，在第i个周期的频率计算参考点所对应的编号；
[0080]
因此，第i个周期的频率计算参考点可以表示为：
[0081][0082]
式中，n
ref,i
表示在第i个周期的频率计算参考点所对应的编号。
[0083]
4)根据相邻周期的频率计算参考点之间的时间差继而获得信号的频率测量值，表示为：
[0084][0085]
式中，f
mea,i
为第i个周期的频率测量值，fs为信号的采样频率，n
ref,i
和n
ref,i-1
分别表示第i个周期和第i-1个周期的频率计算参考点所对应的编号。
[0086]
下面我们以具体的案例对上述频率测量算法进行详细说明：
[0087]
案例2为对案例1中滤波后输出的信号，进行高帽变换和低帽变换，得到的结果如图2所示，其中矩形框代表一个周期。在该周期内，可以看到高帽变换的零点区间的最左侧零点和低帽变换的零点区间的最右侧点重合，取二者平均值后得到的参考点与真实过零点一致，与本发明的分析吻合，即在同一个周期中，当信号的真实频率等于额定频率时，高帽变换处理结果的最左侧的零点和低帽变换最右侧的零点重合。其它周期的结果相类似，最终得到频率测量值为50.07hz，相对误差为0.140％。
[0088]
案例3中待测量信号为45hz的正弦信号，信号采样频率为10khz，添加了20db白噪声，进行数学形态学滤波时，所用的软结构元素为扁平结构元素，b1长度为5，b2长度为2。对滤波后的信号进行高帽变换和低帽变换，如图3所示，其中矩形框代表一个周期。由于信号的真实频率小于额定频率(50hz)，使得高帽变换的零点区间的最左侧零点滞后于滤波后信号的过零点，而低帽变换的零点区间的最右侧零点超前于滤波后信号的过零点，与本发明的分析吻合，即在同一个周期中，当真实频率小于额定频率时，则高帽变换的零点被低帽变换的零点所包含。取高帽变换和低帽变换二者平均值得到的频率计算参考点与真实过零点仍一致。在其它周期中也可得到类似的结果，最终计算得到为44.96hz，相对误差0.089％。
[0089]
案例4中待测量信号为55hz的正弦信号，信号采样频率为10khz，添加了20db白噪声。进行数学形态学滤波时，所用的软结构元素为扁平结构元素，b1长度为5，b2长度为2。对滤波后的信号进行高帽变换和低帽变换，如图4所示，其中矩形框代表一个周期。由于信号的真实频率大于额定频率(50hz)，使得高帽变换的零点区间的最左侧零点超前于滤波后信号的过零点，而低帽变换的零点区间的最右侧零点滞后于滤波后信号的过零点，与本发明的分析吻合，即在同一个周期中，当真实频率大于额定频率时，低帽变换的零点被高帽变换的零点所包含。取高帽变换和低帽变换二者平均值得到的频率计算参考点与真实过零点仍一致。在其它周期中也可得到类似的结果，最终计算得到频率为55.03hz，相对误差0.054％。
[0090]
案例5中待测量信号为基波频率为50hz的非稳态信号，信号采样频率为10khz，添加了20db白噪声。进行数学形态学滤波时，所用的软结构元素为扁平结构元素，b1长度为5，b2长度为2。对滤波后的信号进行高帽变换和低帽变换，如图5所示，其中矩形框代表一个周期。由于叠加了衰减直流分量产生了正向偏移，使得滤波后信号的过零点发生漂移。高帽变换的零点区间的最左侧零点和低帽变换的零点区间的最右侧零点均滞后于滤波后信号的过零点，取高帽变换和低帽变换二者平均值得到的频率计算参考点与基波频率信号的过零点更接近，说明本发明方法不受过零点漂移的影响，适用于非稳态情况下的信号频率测量。在其它周期中也可得到类似的结果，最终计算得到频率为49.98hz，0.040％。
[0091]
上述实施例为本发明较佳的实施方式，但本发明的实施方式并不受上述实施例的限制，其他的任何未背离本发明的精神实质与原理下所作的改变、修饰、替代、组合、简化，均应为等效的置换方式，都包含在本发明的保护范围之内。