CN105137185B - 一种基于离散傅里叶变换的频域插值电力谐波分析方法 - Google Patents
一种基于离散傅里叶变换的频域插值电力谐波分析方法 Download PDFInfo
- Publication number
- CN105137185B CN105137185B CN201510439330.5A CN201510439330A CN105137185B CN 105137185 B CN105137185 B CN 105137185B CN 201510439330 A CN201510439330 A CN 201510439330A CN 105137185 B CN105137185 B CN 105137185B
- Authority
- CN
- China
- Prior art keywords
- mrow
- msub
- frequency
- mfrac
- lambda
- 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.)
- Expired - Fee Related
Links
Abstract
本发明公开了一种基于离散傅里叶变换的频域插值电力谐波分析方法，包括以下步骤，(1) 时域连续信号离散化，(2) 根据测量要求，选取合适的组合余弦窗对离散采样序列加权，(3) 对序列进行快速傅里叶变换（4）以谱峰搜索的形式寻找所测谐波附近最大的两条谱线，并计算他们的谱值比，（5）根据牛顿插值方法获得频率偏差，（6）根据求得的频率偏差，进一步求解该次谐波的幅值、频率和初相位。本发明的基于查表的快速插值与窗函数无关，具有普适性，可与任何窗函数结合进行谐波分析，同时，该方法无复杂运算，计算量较小，易于嵌入式系统中实现。
Description
技术领域
本发明涉及一种基于离散傅里叶变换的频域插值电力谐波分析方法，属于电力系统谐波测量领域。
背景技术
随着电网智能化的发展，由于越来越多的分布式电源并网运行，使得谐波污染越来越严重，也使谐波问题更加复杂。准确的谐波测量是实施管理和控制电网谐波的前提，因此，谐波测量技术受到了广泛关注。
目前，研究者已提出较多谐波参数估计的方法。在这些方法中，基于DFT的频谱分析方法，是当前应用最广泛的电力系统谐波参数估计方法。DFT方法物理意义明确，并有FFT提高信号估计的实时性，故被IEC(International Electrotechnical Commission,IEC)标准所推荐。
DFT方法假定信号是周期性的。但由于实际电网负荷变动导致基波频率波动和间谐波的影响，使得整周期采样不可实现。即使采用数字锁相环技术等硬件同步，也会有一个或几个采样周期的延时，不能完全实现同步采样。因此，常采用加窗插值方法对FFT分析结果进行频谱校正。旁瓣性能较好的窗函数可以抑制信号由于非整周期截断造成的频谱泄漏，插值方法则能够消除栅栏效应误差。
加窗插值FFT方法在一定程度上提高了谐波分析的精度。目前，已有较多的窗函数和插值方法被提出。但是，在这些加窗插值方法中，只有基于Rife-Vincent(I)窗的插值方法，在频率偏差与谱值比之间有解析解，其他形式的加窗插值方法只能通过求解高阶方程、迭代等方式，求得频率偏差的近似解，求解高阶方程或迭代运算，大大增加了谐波分析过程中的计算量，影响了相关方法在实际中的应用。
发明内容
本发明提出一种基于离散傅里叶变换的频域插值电力谐波分析方法，该方法独立于窗函数项数和类型，根据查表和牛顿插值方法来获得非整周期采样条件下的频率偏差，进而求得谐波幅值、频率和初始相位，该插值方法具有较高的谐波分析准确度，同时降低谐波分析过程的计算量，易于在嵌入式系统中实现。
为达到上述目的，本发明采用的技术方案如下：
一种基于离散傅里叶变换的频域插值电力谐波分析方法，包括以下步骤：
(1)时域连续信号离散化：将电力系统中含有谐波成分的电压/电流时域连续信号x(t)通过低通滤波器后，以固定的采样频率fs进行采样，得到离散采样序列x(n)，其中，n表示第n个采样点；
(2)根据测量要求，选取合适的组合余弦窗对离散采样序列x(n)加权，得到N点长序列hN(n)，其中，N表示采样点数；
(3)对序列hN(n)进行快速傅里叶变换，然后对其计算结果取模求绝对值，获得包括基波在内的各次谐波幅值；
(4)以谱峰搜索的形式寻找所测谐波附近最大的两条谱线，并计算谱值比γ：
γ＝次最大谱线值/最大谱线值；
(5)根据牛顿插值方法获得频率偏差；
(6)根据所述步骤(5)求得的频率偏差，进一步求解该次谐波的幅值、频率和初相位。前述的步骤(3)中，各次谐波幅值的求解方法如下：
3-1)对hN(n)进行快速傅里叶变换，根据傅里叶变换频域卷积定理，第r次谐波的离散频谱H(λ)为：
式中：G(λ)为余弦窗函数g(n)的频域函数；Ar和φr分别为第r次谐波的幅值和初始相位，R为谐波最高次数，λr为第r次谐波频率被频率分辨率归一化后的值，λ为被频率分辨率归一化的频率；
3-2)在不计谐波间相互干扰及负频次谐波影响情况下，式(2)可整理为：
式中，σ＝λ-λr为频率偏差，根据式(3)，第r次谐波幅值Ar表示为：
前述的步骤(5)根据牛顿插值方法获得频率偏差，包括以下步骤：
5-1)在组合余弦窗从中心位置到半个归一化频率[0,0.5]区间内，将频率均分为P等份，则均分点处对应的频率值λi等于i＝0,1...P，计算均分点处频率λi和λi-1处的主瓣幅值，计算方法如下：
组合余弦窗的离散时域表达式g(n)为：
其中，M为窗函数的项数，am为窗函数系数，
式(5)经离散傅里叶变换，得到余弦窗频域函数G(λ)：
式中，λ为被频率分辨率归一化的频率，
把λi和λi-1的值作为自变量λ代入式(6)求绝对值，即可获得对应的主瓣幅值；
5-2)计算λi-1处主瓣幅值与λi处主瓣幅值的比值，该比值为一组离散数据，用ζi表示：
5-3)以ζi为自变量，λi为因变量，构建牛顿插值多项式f(x)如下：
f(x)＝f[x0]+f[x0,x1](x-x0)
+f[x0,x1,x2](x-x0)(x-x1)+··· (8)
+f[x0,x1,···,xi](x-x0)···(x-xi-1)
式中，
其中，x0,x1,···,xi分别对应ζ0,ζ1,···,ζi的计算值，f(x0),f(x1),···,f(xi)分别对应λ0,λ1,···,λi的数值；
5-4)利用构造的牛顿插值多项式，计算谱值比为γ时对应的频率偏差σ。
前述的步骤5-1)中，P≤10。
前述的步骤(6)中，谐波的幅值、频率和初相位的求解方法为：
将所述步骤(5)求得的σ作为自变量λ代入式(6)求绝对值，即可获得G(σ)的值，H(λ)为最接近第r次谐波离散频谱H(λ)的幅值，也为已知，故根据式求得第r次谐波幅值Ar；
设频域轴上的第p根谱线最接近第r次谐波谱线，则第r次谐波频率fr可表示为：
其中，f0为基波频率，Ts为数据采集系统的采样周期；
最后，第r次谐波的初相位φr为：
φr＝φp+πσ (10)
式中，φp为频域轴上最接近第r次谐波的第p根谱线的初相角。
本发明的优点是：1、具有普遍适用性，在本发明的公式推导过程中并没有指定谐波次数，故该方法可用于任意次数谐波的参数求解；2、方法较为简单，无复杂运算，计算量小，易于在嵌入式系统中实现。
附图说明
图1本发明的实施例中余弦窗主瓣[0,0.5]区间内，i取不同值时的λi和λi-1位置关系；
图2本发明的实施例中四项Blackman-Harris窗归一化频率[0,0.5]区间内，均分为10份条件下，由本发明方法计算得到频率偏差与理论值之间的误差。
具体实施方式
下面结合附图和具体实施方式对本发明作进一步的说明。
本发明的基于离散傅里叶变换的频域插值电力谐波分析方法包括以下步骤：
(1)时域连续信号离散化：将电力系统中含有谐波成分的电压/电流时域连续信号x(t)通过低通滤波器后，以固定的采样频率fs进行采样，得到离散采样序列x(n)，其中，n表示第n个采样点。
非整周期采样下，一段电网采样信号，其时域离散采样序列x(n)可表示为：
式中：R为谐波最高次数，f0为基波频率，Ar和φr分别为第r次谐波的幅值和初始相位，Ts为数据采集系统的采样周期，采样频率1/Ts满足奈奎斯特抽样定理。
(2)根据测量要求，选取合适的组合余弦窗对离散采样序列x(n)加权，得到N点长序列hN(n)，其中，N表示采样点数。
用长度为N的余弦窗函数g(n)对离散采样序列x(n)加权，得到序列hN(n)，hN(n)＝x(n)·g(n)。
(3)对序列hN(n)进行快速傅里叶变换(fast Fourier transform,FFT)，然后对其计算结果取模求绝对值，获得包括基波在内的各次谐波幅值。
对hN(n)进行FFT变换，根据傅里叶变换频域卷积定理，第r次谐波的离散频谱H(λ)为：
式中：G(λ)为余弦窗函数g(n)的频域函数；λr为第r次谐波频率被频率分辨率归一化后的值，λ为被频率分辨率归一化的频率，
记频率分辨率为F0，F0＝1/NTs，
λr＝rf0/F0。
在不计谐波间相互干扰及负频次谐波影响情况下，式(2)可整理为：
式中，σ＝λ-λr为频率偏差。根据式(3)，第r次谐波幅值Ar表示为：
(4)以谱峰搜索的形式寻找所测谐波附近最大的两条谱线，并计算他们的谱值比γ，
γ＝次最大谱线值/最大谱线值。
(5)根据牛顿插值方法获得频率偏差σ，
由于在数据采样时频率分辨率满足各次谐波频率处于不同的频率间隔中，故频率偏差小于半个频率分辨率，将频率用频率分辨率归一化。在组合余弦窗从中心位置到半个归一化频率[0,0.5]区间内，将频率均分为P等份(P≤10)，则均分点处对应的频率值λi等于i＝0,1...P。
频域采样时，如果一个采样点落在λi处，由于相邻采样点之间的差值为1，则在λi减去1处，也会有一个采样点，记为ρi，则ρi＝λi-1。λi、ρi及余弦窗频域函数G(λ)(λ为频率分辨率归一化频率)之间的关系如图1所示。
计算ρi和λi处的主瓣幅值G(ρi)、G(λi)，求两者的比值并将其存储在内存中。
这里，计算主瓣幅值的方法如下：
组合余弦窗的离散时域表达式g(n)为：
其中，M为窗函数的项数，am为窗函数系数。
式(5)经离散傅里叶变换，得到余弦窗频域函数G(λ)：
式中，λ为被频率分辨率归一化的频率，
把λi和ρi的值作为自变量λ代入式(6)求绝对值，即可获得对应频率处的主瓣幅值。
计算ρi对应的主瓣幅值与λi处主瓣幅值的比值，该比值为一组离散数据，用ζi表示：
以ζi为自变量，λi为因变量，构建牛顿插值多项式f(x)如下：
f(x)＝f[x0]+f[x0,x1](x-x0)
+f[x0,x1,x2](x-x0)(x-x1)+··· (8)
+f[x0,x1,···,xi](x-x0)···(x-xi-1)
式中，
其中，x0,x1,···,xi分别对应ζ0,ζ1,···,ζi的计算值，f(x0),f(x1),···,f(xi)分别对应λ0,λ1,···,λi的数值。
假设被测信号是第r次谐波，非整周期采样下，频域内接近第r次谐波的次最大谱线与最大谱线之比为γ。显然，γ的值介于数组ζi中最小值和最大值之间，故可以用上述构造的牛顿插值多项式，计算谱值比为γ时对应的频率偏差σ。
(6)根据求得的频率偏差，进一步求解该次谐波的幅值、频率和初相位。
在步骤(3)，已经得到第r次谐波幅值Ar表示为式(4)：
式中的σ为频率偏差，由步骤(5)求得，则将σ作为自变量λ代入式(6)求绝对值，即可获得|G(σ)|的值，|H(λ)|为最接近第r次谐波离散频谱H(λ)的幅值，也为已知，故可以根据(4)求得第r次谐波幅值Ar。
设频域轴上的第p根谱线最接近第r次谐波谱线。则第r次谐波频率fr可表示为：
最后，由式(3)可知，第r次谐波的初相位φr为：
φr＝φp+πσ (10)
式中，φp为频域轴上最接近第r次谐波的第p根谱线的初相角。
以上即为该频域插值方法求解电力谐波参数的过程。
为了进一步说明本发明的具体实施方式，以四项Blackman-Harris窗为例，选用该窗对采样数据加权，然后用本发明所提频域插值方法进行谐波分析。
将四项Blackman-Harris窗在归一化频率[0,0.5]区间内，均分为10份，则各均分点处频率为λi＝0.05·i(i＝0,1,···,10)，分别计算λi和λi-1处的主瓣幅值；然后计算i取不同值时，λi-1与λi对应主瓣幅值的比值ζi。其主瓣幅值和比值如表1所示。将λi和ζi以内存形式存入内存。以ζi为自变量，λi为因变量，编写牛顿插值程序。
表1为四项Blackman-Harris窗归一化频率[0,0.5]区间内，均分为10份时，各均分点处频率为λi＝0.05·i(i＝0,1,···,10)与(λi-1)处的主瓣幅值，及两者主瓣幅值之比。
表1
从表1可知，λi和ζi都是单调递增的，故可以用插值方法求解自变量的某个内插值点对应的函数值。以电力谐波为例，频域内某次谐波附近，次最大谱线与最大谱线的幅值比为0.8849，把该值代入以ζi为自变量λi为因变量的牛顿插值多项式程序，求得频率偏差为0.3400。然后，根据该频率偏差就可以求得该次谐波的幅值、频率和初相位。
实际中，任意次数谐波附近的次最大谱线与最大谱线之比在[0.6805,1]之间任意取值时，不考虑随机噪声和间谐波等影响，及四项Blackman-Harris窗在归一化频率[0,0.5]区间内均分为10等份条件下，由本发明方法计算得到频率偏差与理论值之间的误差如图2所示，稳定在10-13数量级，测量误差非常小。
Claims (3)
1.一种基于离散傅里叶变换的频域插值电力谐波分析方法，其特征在于，包括以下步骤：
(1)时域连续信号离散化：将电力系统中含有谐波成分的电压/电流时域连续信号x(t)通过低通滤波器后，以固定的采样频率fs进行采样，得到离散采样序列x(n)，其中，n表示第n个采样点；
(2)根据测量要求，选取合适的组合余弦窗对离散采样序列x(n)加权，得到N点长序列hN(n)，其中，N表示采样点数；
(3)对序列hN(n)进行快速傅里叶变换，然后对其计算结果取模求绝对值，获得包括基波在内的各次谐波幅值；
(4)以谱峰搜索的形式寻找所测谐波附近最大的两条谱线，并计算谱值比γ：
γ＝次最大谱线值/最大谱线值；
(5)根据牛顿插值方法获得频率偏差，包括以下步骤：
5-1)在组合余弦窗从中心位置到半个归一化频率[0,0.5]区间内，将频率均分为P等份，则均分点处对应的频率值λi等于计算均分点处频率λi和λi-1处的主瓣幅值，计算方法如下：
组合余弦窗的离散时域表达式g(n)为：
<mrow>
<mi>g</mi>
<mrow>
<mo>(</mo>
<mi>n</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>m</mi>
<mo>=</mo>
<mn>0</mn>
</mrow>
<mrow>
<mi>M</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</munderover>
<msup>
<mrow>
<mo>(</mo>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
<mi>m</mi>
</msup>
<msub>
<mi>a</mi>
<mi>m</mi>
</msub>
<mi>c</mi>
<mi>o</mi>
<mi>s</mi>
<mrow>
<mo>(</mo>
<mfrac>
<mrow>
<mn>2</mn>
<mi>&pi;</mi>
</mrow>
<mi>N</mi>
</mfrac>
<mi>m</mi>
<mi>n</mi>
<mo>)</mo>
</mrow>
<mo>,</mo>
<mi>n</mi>
<mo>=</mo>
<mn>0</mn>
<mo>,</mo>
<mn>1</mn>
<mo>,</mo>
<mn>2</mn>
<mo>...</mo>
<mi>N</mi>
<mo>-</mo>
<mn>1</mn>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>5</mn>
<mo>)</mo>
</mrow>
</mrow>
其中，M为窗函数的项数，am为窗函数系数，
式(5)经离散傅里叶变换，得到余弦窗频域函数G(λ)：
<mrow>
<mi>G</mi>
<mrow>
<mo>(</mo>
<mi>&lambda;</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>m</mi>
<mo>=</mo>
<mn>0</mn>
</mrow>
<mrow>
<mi>M</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</munderover>
<msup>
<mrow>
<mo>(</mo>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
<mi>m</mi>
</msup>
<mfrac>
<msub>
<mi>a</mi>
<mi>m</mi>
</msub>
<mn>2</mn>
</mfrac>
<mo>&lsqb;</mo>
<mi>D</mi>
<mrow>
<mo>(</mo>
<mi>&lambda;</mi>
<mo>-</mo>
<mi>m</mi>
<mo>)</mo>
</mrow>
<mo>+</mo>
<mi>D</mi>
<mrow>
<mo>(</mo>
<mi>&lambda;</mi>
<mo>+</mo>
<mi>m</mi>
<mo>)</mo>
</mrow>
<mo>&rsqb;</mo>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>6</mn>
<mo>)</mo>
</mrow>
</mrow>
式中，λ为被频率分辨率归一化的频率，
把λi和λi-1的值作为自变量λ代入式(6)求绝对值，即可获得对应的主瓣幅值；
5-2)计算λi-1处主瓣幅值与λi处主瓣幅值的比值，该比值为一组离散数据，用ζi表示：
<mrow>
<msub>
<mi>&zeta;</mi>
<mi>i</mi>
</msub>
<mo>=</mo>
<mfrac>
<mrow>
<mo>|</mo>
<mi>G</mi>
<mrow>
<mo>(</mo>
<msub>
<mi>&lambda;</mi>
<mi>i</mi>
</msub>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
<mo>|</mo>
</mrow>
<mrow>
<mo>|</mo>
<mi>G</mi>
<mrow>
<mo>(</mo>
<msub>
<mi>&lambda;</mi>
<mi>i</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>|</mo>
</mrow>
</mfrac>
<mo>=</mo>
<mfrac>
<mrow>
<mo>|</mo>
<mi>G</mi>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>&CenterDot;</mo>
<mfrac>
<mn>0.5</mn>
<mi>P</mi>
</mfrac>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
<mo>|</mo>
</mrow>
<mrow>
<mo>|</mo>
<mi>G</mi>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>&CenterDot;</mo>
<mfrac>
<mn>0.5</mn>
<mi>P</mi>
</mfrac>
<mo>)</mo>
</mrow>
<mo>|</mo>
</mrow>
</mfrac>
<mo>,</mo>
<mi>i</mi>
<mo>=</mo>
<mn>0</mn>
<mo>,</mo>
<mn>1.</mn>
<mi>P</mi>
<mo>.</mo>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>7</mn>
<mo>)</mo>
</mrow>
<mo>;</mo>
</mrow>
5-3)以ζi为自变量，λi为因变量，构建牛顿插值多项式f(x)如下：
<mrow>
<mtable>
<mtr>
<mtd>
<mrow>
<mi>f</mi>
<mrow>
<mo>(</mo>
<mi>x</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mi>f</mi>
<mo>&lsqb;</mo>
<msub>
<mi>x</mi>
<mn>0</mn>
</msub>
<mo>&rsqb;</mo>
<mo>+</mo>
<mi>f</mi>
<mo>&lsqb;</mo>
<msub>
<mi>x</mi>
<mn>0</mn>
</msub>
<mo>,</mo>
<msub>
<mi>x</mi>
<mn>1</mn>
</msub>
<mo>&rsqb;</mo>
<mrow>
<mo>(</mo>
<mi>x</mi>
<mo>-</mo>
<msub>
<mi>x</mi>
<mn>0</mn>
</msub>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>+</mo>
<mi>f</mi>
<mo>&lsqb;</mo>
<msub>
<mi>x</mi>
<mn>0</mn>
</msub>
<mo>,</mo>
<msub>
<mi>x</mi>
<mn>1</mn>
</msub>
<mo>,</mo>
<msub>
<mi>x</mi>
<mn>2</mn>
</msub>
<mo>&rsqb;</mo>
<mrow>
<mo>(</mo>
<mi>x</mi>
<mo>-</mo>
<msub>
<mi>x</mi>
<mn>0</mn>
</msub>
<mo>)</mo>
</mrow>
<mrow>
<mo>(</mo>
<mi>x</mi>
<mo>-</mo>
<msub>
<mi>x</mi>
<mn>1</mn>
</msub>
<mo>)</mo>
</mrow>
<mo>+</mo>
<mn>...</mn>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>+</mo>
<mi>f</mi>
<mo>&lsqb;</mo>
<msub>
<mi>x</mi>
<mn>0</mn>
</msub>
<mo>,</mo>
<msub>
<mi>x</mi>
<mn>1</mn>
</msub>
<mo>,</mo>
<mn>...</mn>
<mo>,</mo>
<msub>
<mi>x</mi>
<mi>i</mi>
</msub>
<mo>&rsqb;</mo>
<mrow>
<mo>(</mo>
<mi>x</mi>
<mo>-</mo>
<msub>
<mi>x</mi>
<mn>0</mn>
</msub>
<mo>)</mo>
</mrow>
<mn>...</mn>
<mrow>
<mo>(</mo>
<mi>x</mi>
<mo>-</mo>
<msub>
<mi>x</mi>
<mrow>
<mi>i</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
</mtable>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>8</mn>
<mo>)</mo>
</mrow>
</mrow>
式中，
<mrow>
<mi>f</mi>
<mo>&lsqb;</mo>
<msub>
<mi>x</mi>
<mn>0</mn>
</msub>
<mo>,</mo>
<msub>
<mi>x</mi>
<mn>1</mn>
</msub>
<mo>,</mo>
<mo>...</mo>
<mo>,</mo>
<msub>
<mi>x</mi>
<mi>i</mi>
</msub>
<mo>&rsqb;</mo>
<mo>=</mo>
<mfrac>
<mrow>
<mi>f</mi>
<mo>&lsqb;</mo>
<msub>
<mi>x</mi>
<mn>0</mn>
</msub>
<mo>,</mo>
<msub>
<mi>x</mi>
<mn>1</mn>
</msub>
<mo>...</mo>
<mo>,</mo>
<msub>
<mi>x</mi>
<mrow>
<mi>i</mi>
<mo>-</mo>
<mn>2</mn>
</mrow>
</msub>
<mo>,</mo>
<msub>
<mi>x</mi>
<mi>i</mi>
</msub>
<mo>&rsqb;</mo>
<mo>-</mo>
<mi>f</mi>
<mo>&lsqb;</mo>
<msub>
<mi>x</mi>
<mn>0</mn>
</msub>
<mo>,</mo>
<msub>
<mi>x</mi>
<mn>1</mn>
</msub>
<mo>...</mo>
<mo>,</mo>
<msub>
<mi>x</mi>
<mrow>
<mi>i</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mo>&rsqb;</mo>
</mrow>
<mrow>
<msub>
<mi>x</mi>
<mi>i</mi>
</msub>
<mo>-</mo>
<msub>
<mi>x</mi>
<mrow>
<mi>i</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
</mrow>
</mfrac>
</mrow>
1
其中，x0,x1,…,xi分别对应ζ0,ζ1,…,ζi的计算值，f(x0),f(x1),…,f(xi)分别对应λ0,λ1,…,λi的数值；
5-4)利用构造的牛顿插值多项式，计算谱值比为γ时对应的频率偏差σ；
(6)根据所述步骤(5)求得的频率偏差，进一步求解该次谐波的幅值、频率和初相位，求解方法为：
将所述步骤(5)求得的σ作为自变量λ代入式(6)求绝对值，即可获得|G(σ)|的值，|H(λ)|为最接近第r次谐波离散频谱H(λ)的幅值，也为已知，故根据式求得第r次谐波幅值Ar；
设频域轴上的第p根谱线最接近第r次谐波谱线，则第r次谐波频率fr可表示为：
<mrow>
<msub>
<mi>f</mi>
<mi>r</mi>
</msub>
<mo>=</mo>
<msub>
<mi>rf</mi>
<mn>0</mn>
</msub>
<mo>=</mo>
<mfrac>
<mrow>
<mo>(</mo>
<mi>p</mi>
<mo>-</mo>
<mn>1</mn>
<mo>-</mo>
<mi>&sigma;</mi>
<mo>)</mo>
</mrow>
<mrow>
<msub>
<mi>T</mi>
<mi>s</mi>
</msub>
<mi>N</mi>
</mrow>
</mfrac>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>9</mn>
<mo>)</mo>
</mrow>
</mrow>
其中，f0为基波频率，Ts为数据采集系统的采样周期；
最后，第r次谐波的初相位φr为：
φr＝φp+πσ (10)
式中，φp为频域轴上最接近第r次谐波的第p根谱线的初相角。
2.根据权利要求1所述的一种基于离散傅里叶变换的频域插值电力谐波分析方法，其特征在于，所述步骤(3)中，各次谐波幅值的求解方法如下：
3-1)对hN(n)进行快速傅里叶变换，根据傅里叶变换频域卷积定理，第r次谐波的离散频谱H(λ)为：
<mrow>
<mi>H</mi>
<mrow>
<mo>(</mo>
<mi>&lambda;</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<msub>
<mi>A</mi>
<mi>r</mi>
</msub>
<mn>2</mn>
</mfrac>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>r</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>R</mi>
</munderover>
<mo>&lsqb;</mo>
<mi>G</mi>
<mrow>
<mo>(</mo>
<mi>&lambda;</mi>
<mo>-</mo>
<msub>
<mi>&lambda;</mi>
<mi>r</mi>
</msub>
<mo>)</mo>
</mrow>
<msup>
<mi>e</mi>
<mrow>
<mo>-</mo>
<mi>j</mi>
<mo>&lsqb;</mo>
<mi>&pi;</mi>
<mrow>
<mo>(</mo>
<mi>&lambda;</mi>
<mo>-</mo>
<msub>
<mi>&lambda;</mi>
<mi>r</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>-</mo>
<msub>
<mi>&phi;</mi>
<mi>r</mi>
</msub>
<mo>&rsqb;</mo>
</mrow>
</msup>
<mo>+</mo>
<mi>G</mi>
<mrow>
<mo>(</mo>
<mi>&lambda;</mi>
<mo>+</mo>
<msub>
<mi>&lambda;</mi>
<mi>r</mi>
</msub>
<mo>)</mo>
</mrow>
<msup>
<mi>e</mi>
<mrow>
<mo>-</mo>
<mi>j</mi>
<mo>&lsqb;</mo>
<mi>&pi;</mi>
<mrow>
<mo>(</mo>
<mi>&lambda;</mi>
<mo>+</mo>
<msub>
<mi>&lambda;</mi>
<mi>r</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>+</mo>
<msub>
<mi>&phi;</mi>
<mi>r</mi>
</msub>
<mo>&rsqb;</mo>
</mrow>
</msup>
<mo>&rsqb;</mo>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>2</mn>
<mo>)</mo>
</mrow>
</mrow>
式中：G(λ)为余弦窗函数g(n)的频域函数；Ar和φr分别为第r次谐波的幅值和初始相位，R为谐波最高次数，λr为第r次谐波频率被频率分辨率归一化后的值，λ为被频率分辨率归一化的频率；
3-2)在不计谐波间相互干扰及负频次谐波影响情况下，式(2)可整理为：
<mrow>
<mi>H</mi>
<mrow>
<mo>(</mo>
<mi>&lambda;</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<msub>
<mi>A</mi>
<mi>r</mi>
</msub>
<mn>2</mn>
</mfrac>
<mi>G</mi>
<mrow>
<mo>(</mo>
<mi>&sigma;</mi>
<mo>)</mo>
</mrow>
<msup>
<mi>e</mi>
<mrow>
<mo>-</mo>
<mi>j</mi>
<mrow>
<mo>(</mo>
<mi>&pi;</mi>
<mi>&sigma;</mi>
<mo>-</mo>
<msub>
<mi>&phi;</mi>
<mi>r</mi>
</msub>
<mo>)</mo>
</mrow>
</mrow>
</msup>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>3</mn>
<mo>)</mo>
</mrow>
</mrow>
式中，σ＝λ-λr为频率偏差，根据式(3)，第r次谐波幅值Ar表示为：
<mrow>
<msub>
<mi>A</mi>
<mi>r</mi>
</msub>
<mo>=</mo>
<mfrac>
<mrow>
<mn>2</mn>
<mo>|</mo>
<mi>H</mi>
<mrow>
<mo>(</mo>
<mi>&lambda;</mi>
<mo>)</mo>
</mrow>
<mo>|</mo>
</mrow>
<mrow>
<mo>|</mo>
<mi>G</mi>
<mrow>
<mo>(</mo>
<mi>&sigma;</mi>
<mo>)</mo>
</mrow>
<mo>|</mo>
</mrow>
</mfrac>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>4</mn>
<mo>)</mo>
</mrow>
<mo>.</mo>
</mrow>
3.根据权利要求1所述的一种基于离散傅里叶变换的频域插值电力谐波分析方法，其特征在于，所述步骤5-1)中，P≤10。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201510439330.5A CN105137185B (zh) | 2015-07-23 | 2015-07-23 | 一种基于离散傅里叶变换的频域插值电力谐波分析方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201510439330.5A CN105137185B (zh) | 2015-07-23 | 2015-07-23 | 一种基于离散傅里叶变换的频域插值电力谐波分析方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN105137185A CN105137185A (zh) | 2015-12-09 |
| CN105137185B true CN105137185B (zh) | 2017-09-22 |
Family
ID=54722601
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201510439330.5A Expired - Fee Related CN105137185B (zh) | 2015-07-23 | 2015-07-23 | 一种基于离散傅里叶变换的频域插值电力谐波分析方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN105137185B (zh) |
Families Citing this family (18)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN105675962A (zh) * | 2016-03-02 | 2016-06-15 | 国网上海市电力公司 | 一种用于变压器绕组状态在线监测的电压电流采样方法 |
| CN105822289B (zh) * | 2016-03-25 | 2019-03-19 | 重庆科技学院 | 用于油井动液面检测的频率估算方法 |
| CN106053940B (zh) * | 2016-08-09 | 2018-04-10 | 重庆大学 | 一种基于方波傅里叶级数分解的谐波分析方法 |
| CN106324342B (zh) * | 2016-08-22 | 2023-06-02 | 江苏南自通华智慧能源股份有限公司 | 一种基于查表的谐波检测方法 |
| CN106405229B (zh) * | 2016-08-30 | 2019-04-30 | 威胜集团有限公司 | 一种基波和谐波电能计量方法 |
| CN107064634B (zh) * | 2017-04-18 | 2019-07-02 | 辽宁大学 | 电力系统谐波的检测方法 |
| CN108490257B (zh) * | 2018-03-26 | 2020-08-21 | 华北电力大学 | 一种基于频谱曲线拟合的短时窗间谐波测量方法 |
| CN108776263B (zh) * | 2018-05-02 | 2020-07-28 | 三峡大学 | 基于高阶汉宁自卷积窗及改进插值算法的谐波检测方法 |
| CN112908345B (zh) * | 2019-01-29 | 2022-05-31 | 桂林理工大学南宁分校 | 一种物联网语音压缩与解压方法 |
| CN110007129B (zh) * | 2019-04-17 | 2019-12-03 | 贵州电网有限责任公司 | 一种应用于动态电能计量的三相电压实时估计方法 |
| CN110333389B (zh) * | 2019-05-07 | 2020-12-29 | 东南大学 | 基于插值dft的正弦信号频率估计方法 |
| CN110455223B (zh) * | 2019-07-31 | 2020-12-25 | 中国科学院西安光学精密机械研究所 | 一种基于条纹相位的干涉图倾斜角测量方法 |
| CN112129983A (zh) * | 2020-09-25 | 2020-12-25 | 中北大学 | 一种基于等时间间隔等效取样的波形恢复数据处理方法 |
| CN112327234B (zh) * | 2020-10-29 | 2021-09-14 | 厦门大学 | 可转位刀具切削动态信号的工频干扰高精度补偿方法 |
| CN114252699B (zh) * | 2021-12-24 | 2023-12-15 | 中国电子科技集团公司第十四研究所 | 一种离散傅里叶变换比幅测频方法 |
| CN115389816B (zh) * | 2022-10-27 | 2023-04-07 | 苏州大学 | 一种电网谐波与间谐波的检测方法 |
| CN117169590B (zh) * | 2023-08-15 | 2024-03-29 | 嘉兴市科讯电子有限公司 | 一种基于软件变采样率的电力谐波分析的方法和装置 |
| CN117311300B (zh) * | 2023-11-29 | 2024-02-13 | 西安热工研究院有限公司 | 分散控制系统动态调整采样频率的方法、装置及电子设备 |
Citations (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JPH08248070A (ja) * | 1995-03-08 | 1996-09-27 | Anritsu Corp | 周波数スペクトル分析装置 |
| CN101261292A (zh) * | 2008-04-14 | 2008-09-10 | 湖南大学 | 基于5项Rife-Vincent（I）窗双谱线插值FFT的基波与谐波检测方法 |
| CN101334431A (zh) * | 2007-12-20 | 2008-12-31 | 复旦大学 | 电网谐波的频谱插值校正分析方法 |
| CN101441233A (zh) * | 2008-12-02 | 2009-05-27 | 湖南海兴电器有限责任公司 | 基于Kaiser窗双谱线插值FFT的基波与谐波检测方法 |
| CN202339381U (zh) * | 2011-10-19 | 2012-07-18 | 广西电网公司电力科学研究院 | 一种纳托尔自卷积窗加权傅里叶变换的谐波电能计量系统 |
| CN103207319A (zh) * | 2013-03-12 | 2013-07-17 | 广东电网公司电力科学研究院 | 数字化变电站电力信号非同步采样条件下的谐波测量方法 |
| CN103575984A (zh) * | 2012-08-02 | 2014-02-12 | 西安元朔科技有限公司 | 基于凯塞窗双谱线插值fft的谐波分析方法 |
-
2015
- 2015-07-23 CN CN201510439330.5A patent/CN105137185B/zh not_active Expired - Fee Related
Patent Citations (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JPH08248070A (ja) * | 1995-03-08 | 1996-09-27 | Anritsu Corp | 周波数スペクトル分析装置 |
| US5706202A (en) * | 1995-03-08 | 1998-01-06 | Anritsu Corporation | Frequency spectrum analyzing apparatus and transmitter characteristics measuring apparatus using the same |
| CN101334431A (zh) * | 2007-12-20 | 2008-12-31 | 复旦大学 | 电网谐波的频谱插值校正分析方法 |
| CN101261292A (zh) * | 2008-04-14 | 2008-09-10 | 湖南大学 | 基于5项Rife-Vincent（I）窗双谱线插值FFT的基波与谐波检测方法 |
| CN101441233A (zh) * | 2008-12-02 | 2009-05-27 | 湖南海兴电器有限责任公司 | 基于Kaiser窗双谱线插值FFT的基波与谐波检测方法 |
| CN202339381U (zh) * | 2011-10-19 | 2012-07-18 | 广西电网公司电力科学研究院 | 一种纳托尔自卷积窗加权傅里叶变换的谐波电能计量系统 |
| CN103575984A (zh) * | 2012-08-02 | 2014-02-12 | 西安元朔科技有限公司 | 基于凯塞窗双谱线插值fft的谐波分析方法 |
| CN103207319A (zh) * | 2013-03-12 | 2013-07-17 | 广东电网公司电力科学研究院 | 数字化变电站电力信号非同步采样条件下的谐波测量方法 |
Non-Patent Citations (1)
| Title |
|---|
| 基于Blackman窗函数的插值FFT在电网谐波信号分析中的应用;周俊 等;《浙江大学学报（理学版）》;20061130;第33卷(第6期);第650-653页 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN105137185A (zh) | 2015-12-09 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN105137185B (zh) | 一种基于离散傅里叶变换的频域插值电力谐波分析方法 | |
| CN102435844B (zh) | 一种频率无关的正弦信号相量计算方法 | |
| CN103454497B (zh) | 基于改进加窗离散傅立叶变换的相位差测量方法 | |
| CN103576002B (zh) | 一种容性绝缘设备介质损耗角的计算方法 | |
| CN102539915B (zh) | 时延傅立叶变换测频法精确计算电力谐波参数方法 | |
| CN104391178B (zh) | 一种基于Nuttall窗的时移相位差稳态谐波信号校正方法 | |
| Wen et al. | Triangular self-convolution window with desirable sidelobe behaviors for harmonic analysis of power system | |
| Wen et al. | Hanning self-convolution window and its application to harmonic analysis | |
| CN107102255B (zh) | 单一adc采集通道动态特性测试方法 | |
| CN109521275B (zh) | 一种同步相量确定方法、系统、装置及可读存储介质 | |
| CN103308804A (zh) | 基于快速k-s变换电能质量扰动信号时频参数提取方法 | |
| CN108776263B (zh) | 基于高阶汉宁自卷积窗及改进插值算法的谐波检测方法 | |
| CN111222088B (zh) | 一种改进的平顶自卷积窗加权电力谐波幅值估计方法 | |
| CN107643446B (zh) | 一种基于主瓣宽度的多谱线插值谐波分析方法及系统 | |
| CN105785124A (zh) | 一种采用谱估计和互相关的电力系统谐波和间谐波测量方法 | |
| CN109374966A (zh) | 一种电网频率估计方法 | |
| CN103018555A (zh) | 一种高精度的电力参数软件同步采样方法 | |
| CN110837003B (zh) | 一种基于三角窗的双窗全相位dft同步相量测量方法及系统 | |
| CN109541304B (zh) | 基于六项最小旁瓣窗插值的电网高次弱幅值谐波检测方法 | |
| CN103543331A (zh) | 一种计算电信号谐波和间谐波的方法 | |
| CN109254202B (zh) | 一种应用于配电网的同步相量测量装置 | |
| CN102072987B (zh) | 短区间正弦信号的相位估计法及其实验装置 | |
| CN105372492A (zh) | 基于三条dft复数谱线的信号频率测量方法 | |
| Levashov et al. | The Ways of Reducing of the Active Power Measurement Error for the Method of Averaging of the Instantaneous Power | |
| Chen et al. | A phase difference measurement method based on the extended kalman filter for Coriolis mass flowmeters |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant | ||
| CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20170922 Termination date: 20210723 |
|
| CF01 | Termination of patent right due to non-payment of annual fee |