CN102435845A - Harmonic energy metering method of Nuttall self-convolution weighted Fourier transform - Google Patents

Abstract

The invention relates to a harmonic energy metering method of Nuttall self-convolution weighted Fourier transform. The method comprises the following steps that: a discrete Nuttall window sequence with a length of N is constructed; convolution operation is carried out on the discrete Nuttall window so as to obtain a convolution sequence; zero filling is carried out on the head or the tail of the convolution sequence so as to obtain a discrete Nuttall self-convolution window; the Nuttall self-convolution window sequence is utilized to carry out weighted operation on a discrete power signal; Fourier transform is employed to obtain a discrete frequency spectrum of the power signal; spectral lines with a greatest amplitude and a second great amplitude are searched nearby a harmonic frequency corresponded to the discrete frequency spectrum; fitting is carried out on a peak value spectral line by utilizing a least square method so as to obtain an amplitude, a frequency and an initial phase angle of a power harmonic wave, so that integral harmonic energy metering in a complex power signal is realized.

Description

A kind of harmonic electric energy metering method of Nuttall self-convolution window weighting Fourier transform

Technical field

The present invention relates to a kind of harmonic electric energy metering method, specifically is a kind of electric harmonic electric energy gauging method of Nuttall self-convolution window weighting Fourier transform.Belong to the signal Processing field.

Background technology

The high-precision signal analysis be treated to harmonic trend calculating, Equipment Inspection, Harmonious Waves in Power Systems compensation and inhibition, analysis of vibration signal and random fault processing etc. scientific basis be provided.When adopting the Fourier transform theory to carry out signal Processing, because there is fluctuation in signal frequency or disturbs, strict synchronized sampling can't realize that spectrum leakage and fence effect can be introduced than mistake.

Find that through the document retrieval adopt window function can effectively reduce spectrum leakage to signal weighting, the main lobe of window function frequency spectrum is directly relevant with frequency resolution, main lobe is wide, frequency resolution is low; Secondary lobe is directly relevant with leakage, and secondary lobe is big, and leakage is many; The speed of side lobe attenuation slope reflection side lobe attenuation, side lobe attenuation is fast more, and is strong more to leak suppressing.Chinese scholars has proposed a series of window functions, like Hanning window, Blackman-Harris window, Rife-Vincent (I) window, Nuttall window and rectangle convolution window etc., and they is applied in the signal weighting Processing Algorithm, has improved the signal Processing precision.

Signal weighting analytical algorithm based on classical window function (Hanning window, Blackman-Harris window, Rife-Vincent (I) window, Nuttall window and rectangle convolution window etc.) utilizes the frequency spectrum secondary lobe dropping characteristic of window function can reduce spectrum leakage to a certain extent.But because the sidelobe performance of classical window is still not ideal enough, still limited to the inhibiting effect of spectrum leakage, the phase mutual interference between harmonic wave can not ignore, and the signal analysis and processing precision is restricted.

Summary of the invention

The object of the invention is to overcome the defective of existing windowing Fourier transform harmonic electric energy metering method; A kind of easy Nuttall self-convolution window weighting Fourier transform harmonic electric energy metering method is provided; Solving signal weighting processing procedure intermediate frequency spectrum leaks excessive; The problem of phase mutual interference between harmonic wave, thus accuracy and the practicality that signal harmonic is analyzed improved, for signal parameter identification and the harmonic electric energy metering of further carrying out provides reliable basis.

The present invention realizes above-mentioned purpose through following technical scheme: a kind of harmonic electric energy metering method of Nuttall self-convolution window weighting Fourier transform comprises the steps:

At first make up discrete Nuttall self-convolution window; Adopting length is that the Nuttall self-convolution window sequence of N is that the discrete electrical force signal of N computes weighted to length; Adopt Fourier transform to obtain the discrete spectrum of electric power signal again; The utilization least square method is carried out the spectrum peak match, obtains amplitude, frequency and the initial phase angle of electric harmonic, and then calculates harmonic electric energy.

The said structure Nuttall self-convolution window that disperses; Be that structure length is the discrete Nuttall window sequence of M earlier; Then p discrete Nuttall window sequence made convolution algorithm p-1 time; Obtain convolution sequence, zero in the head or p-1 of the tail benefit of convolution sequence, promptly obtaining length is the discrete Nuttall self-convolution window in p rank of N=pM.

Said utilization least square method is carried out the spectrum peak match; It is near the maximum and inferior big spectral line of the searching amplitude corresponding harmonic frequency of discrete spectrum; The utilization least square method is carried out match to the peak value spectral line, obtains amplitude, frequency and the initial phase angle of electric harmonic, and then calculates harmonic electric energy.

Said Nuttall self-convolution window weighting Fourier transform is the increase along with the convolution exponent number, and the sidelobe level of Nuttall self-convolution window is reduced rapidly.

Principle of the present invention is:

The Nuttall window is the cosine composite window, and its time domain expression formula does

$w\left(n\right)=\underset{g=0}{\overset{G-1}{\mathrm{\Σ}}}{(-1)}^g{b}_g\cos (2\mathrm{\πn}\·g/N)$

In the formula, M is the item number of window function; N=0,1, Λ, N-1; Bg satisfies constraint condition

$\left\{\begin{array}{c}\underset{g=0}{\overset{G-1}{\mathrm{\Σ}}}{b}_g=1\\ \underset{g=0}{\overset{G-1}{\mathrm{\Σ}}}{(-1)}^g{b}_g=0\end{array}\right.$

Table 1 has provided the coefficient of typical Nuttall window function

The coefficient of table 1 Nuttall window function

The Nuttall self-convolution window of time domain discrete does

Wherein, Subscript p representes to participate in the number of the Nuttall self-convolution window of convolution algorithm, and subscript N representes to participate in the Nuttall window of convolution, for ease of the realization of Fast Fourier Transform (FFT); Convolution results is carried out the zero padding operation, and then the length of discrete Nuttall self-convolution window is N=pM.

If only contain single-frequency signal components x (t), be the analog to digital conversion of fs through SF after, obtain discrete series:

$x\left(m\right)=A_0\sin (2\mathrm{\π}\frac{f_0}{f_s}m+{\mathrm{\φ}}_0),m=\mathrm{0,1,2},L,+\∞$

In the formula, A ₀, f ₀, φ ₀Be respectively amplitude, frequency and the initial phase angle of signal.

Signal after the discretize is added the p rank Nuttall self-convolution window w that length is N _N-p(n) (n=0,1, L, N-1), the sequence that obtains after the brachymemma does

x(n)＝x(m)w _N-p(n)，n＝0，1，L，N-1

After sequence x (n) carried out discrete Fourier transformation, obtain its discrete spectrum and do

$X\left(k\right)=A_0{e}^{{\mathrm{j\φ}}_0}\frac{2^p{e}^{-j\left({k-k}_0\right)\mathrm{\π}}}{M^p}{\left\{\frac{\sin [\mathrm{\π}(k-k_0)/2p]}{\sin [\mathrm{\π}(k-k_0)/N]}\right\}}^{2p}$

In the formula, k ₀=f ₀N/f _s

If in the discrete spectrum, frequency f ₀Near local amplitude is maximum to be respectively k with time maximum spectral line ₁And k ₂Root satisfies k ₁≤k ₀≤k ₂=k ₁+ 1.In frequency f ₀Near the local peaking's search strategy that adopts finds this two spectral lines, can confirm k ₁And k ₂If these two spectral line amplitudes are respectively y ₁And y ₂, promptly

${{\mathrm{<figure-callout\; id="1"\; label="y"\; filenames="HSA00000596614000012.png"\; state="\{\{state\}\}">y</figure-callout>}}_1=\left|X\left(k_1\right)\right|=\frac{A_0{2}^p}{M^p}\left\{\frac{\sin [\mathrm{\&pi;}(k_1-k_0)/2p]}{\sin [\mathrm{\&pi;}(k_1-k_0)/N]}\right\}}^{2p}$

${{\mathrm{<figure-callout\; id="2"\; label="y"\; filenames="HSA00000596614000012.png"\; state="\{\{state\}\}">y</figure-callout>}}_2=\left|X\left(k_2\right)\right|=\frac{A_0{2}^p}{M^p}\left\{\frac{\sin [\mathrm{\&pi;}(k_2-k_0)/2p]}{\sin [\mathrm{\&pi;}(k_2-k_0)/N]}\right\}}^{2p}$

Consider 0≤k ₀-k ₁≤1, the definition alpha does

α＝k ₀-k ₁-0.5，-0.5≤α≤0.5

Y then ₁And y ₂Can be rewritten as function about α:

${{\mathrm{<figure-callout\; id="1"\; label="y"\; filenames="HSA00000596614000012.png"\; state="\{\{state\}\}">y</figure-callout>}}_1=\frac{A_0{2}^p}{M^p}\left\{\frac{\sin [\mathrm{\&pi;}(-\mathrm{\&alpha;}-0.5)/2p]}{\sin [\mathrm{\&pi;}(-\mathrm{\&alpha;}-0.5)/N]}\right\}}^{2p}$

${{\mathrm{<figure-callout\; id="2"\; label="y"\; filenames="HSA00000596614000012.png"\; state="\{\{state\}\}">y</figure-callout>}}_2=\frac{A_0{2}^p}{M^p}\left\{\frac{\sin [\mathrm{\&pi;}(-\mathrm{\&alpha;}+0.5)/2p]}{\sin [\mathrm{\&pi;}(-\mathrm{\&alpha;}+0.5)/N]}\right\}}^{2p}$

The definition factor beta does

$\mathrm{\&beta;}=\frac{y_2-y_1}{y_2+y_1}$

With the y in the formula following formula ₁And y ₂Replace, then β can be written as the function about α:

$\mathrm{\&beta;}=\frac{\left|W_{N-p}\right[2\mathrm{\&pi;}(-\mathrm{\&alpha;}+0.5)/N\left]\right|-\left|W_{N-p}\right[2\mathrm{\&pi;}(-\mathrm{\&alpha;}-0.5)/N\left]\right|}{\left|W_{N-p}\right[2\mathrm{\&pi;}(-\mathrm{\&alpha;}-0.5)/N\left]\right|+\left|W_{N-p}\right[2\mathrm{\&pi;}(-\mathrm{\&alpha;}+0.5)/N\left]\right|}=g\left(\mathrm{\&alpha;}\right)$

In α ∈ [0.5,0.5] scope,, utilize each α with the periodic sampling some spots _iValue calculates corresponding β _i, confirm to treat fitting of a polynomial number of times K, utilize and find the solution coefficient a on the formula _k(k=0,1, L K-1), thereby sets up frequency spectrum interpolation polynomial expression S _K(x), can be designated as

α＝g ^-1(β)≈a ₁β+a ₂β ²+a ₃β ³+L+a _Kβ ^K

In the formula, a ₁, a ₂..., a _KCoefficient for polynomial fitting.

Frequency f ₀Calculating formula do

$f_0=k_0\frac{f_s}{N}=(\mathrm{\&alpha;}+k_1+0.5)\frac{f_s}{N}$

Amplitude does

$A_0=\frac{{2\mathrm{<figure-callout\; id="1"\; label="y"\; filenames="HSA00000596614000012.png"\; state="\{\{state\}\}">y</figure-callout>}}_1}{\left|W_{N-p}\right[2\mathrm{\&pi;}(-\mathrm{\&alpha;}-0.5)/N\left]\right|}$

Or

$A_0=\frac{{2\mathrm{<figure-callout\; id="2"\; label="y"\; filenames="HSA00000596614000012.png"\; state="\{\{state\}\}">y</figure-callout>}}_2}{\left|W_{N-p}\right[2\mathrm{\&pi;}(-\mathrm{\&alpha;}+0.5)/N\left]\right|}$

Initial phase angle does

${\mathrm{\&phi;}}_0=\arg \left[X\left(k_1\right)\right]+\frac{\mathrm{\&pi;}}{2}-\arg \left\{W_{N-p}\right[\frac{2\mathrm{\&pi;}(-\mathrm{\&alpha;}-0.5)}{N}\left]\right\}$

For the each harmonic component, all can find the solution parameter as stated above.

Description of drawings

Fig. 1 is the process flow diagram of the harmonic electric energy metering method of Nuttall self-convolution window weighting Fourier transform of the present invention.

Fig. 2 is a discrete spectrum interpolation processing synoptic diagram of the present invention.

Specific embodiments

Below through embodiment technical scheme of the present invention is further specified.

As shown in Figure 1, behind the signal sampling (discretize), it is added the computing of Nuttall self-convolution window; Promptly the signal of discretize is carried out brachymemma, carry out Fast Fourier Transform (FFT) (FFT) then, obtain discrete spectrum with the Nuttall self-convolution window; Near its local spectrum peak of the search h subfrequency finds maximum spectral line of the local amplitude of discrete spectrum and time big spectral line, and confirms corresponding spectral line amplitude; Utilizing the least square interpolation algorithm to carry out frequency spectrum interpolation calculates; Obtain frequency offset, obtain the each harmonic parameter, thus output harmonic wave energy value.The harmonic electric energy metering method of Nuttall self-convolution window weighting Fourier transform of the present invention comprises the steps:

At first set up discrete 4 minimum secondary lobe Nuttall windows of M=64, G=3, promptly select b0=0.3635819, b1=0.4891775, b2=0.1365995, b3=0.0106411, through the sequence of 4 minimum secondary lobe Nuttall windows of computes:

$w\left(n\right)=\underset{g=0}{\overset{G-1}{\mathrm{\&Sigma;}}}{(-1)}^g{b}_g\cos (2\mathrm{\&pi;n}\&CenterDot;g/M)$

Secondly, w (n) is carried out 2～4 rank from convolution algorithm, promptly p respectively value be 2,3 and 4, carry out corresponding zero padding operation again, make the length of 4 the minimum secondary lobe Nuttall self-convolution windows in 2～4 rank be respectively 128,192 and 256.

The voltage signal that is provided with a multi-frequency composition does

$x\left(n\right)=A_0+A_1\sin (2\mathrm{\&pi;}\frac{f_{1}n}{N}+{\mathrm{\&phi;}}_1)+A_3\sin (2\mathrm{\&pi;}\frac{f_{3}n}{N}+{\mathrm{\&phi;}}_3)$

In the formula, A ₀=0.2V; A ₁=6V; f ₁=20.2Hz; φ ₁=0.1rad; A ₃=1; f ₃=60.6Hz; φ ₃=0rad.

Adopt 4 the minimum secondary lobe Nuttall self-convolution windows in 4 rank set up that discrete series x (n) is carried out weighting, carry out discrete Fourier transformation again, obtain the discrete spectrum after the voltage signal weighting.Discrete spectrum is carried out normalization handle, near the spectrum amplitude the frequency f 1 distributes as shown in Figure 2.Signal frequency f ₁Do not overlap, and be positioned at the maximum spectral line k of discrete spectrum amplitude with discrete spectral line ₁=4 with time big spectral line k ₂Between=5.After confirming that k1 and k2 are respectively maximum and second largest peak value spectral line; Obtain corresponding range value y1 and y2; According to the least square peaks spectrum fit procedure in the aforementioned summary of the invention; Calculate the high reps of polynomial expression when being set to K=5, add 4 the minimum secondary lobe Nuttall self-convolution windows in 4 rank after, the discrete spectrum interpolation polynomial does

α _4th＝0.4026β ⁵+0.7644β ³+11.9048β

Will

Bring following formula into, promptly calculating f1 is 20.1999999996Hz, and A0 is 0.20000000007V, and A1 is 5.9999999999V, φ ₁Be 0.10000008rad.And the like, can calculate A3 is 0.99999999997V, f3 is 60.60000000001Hz, φ ₃For-0.00000002rad.

As shown in Figure 2, discrete spectrum interpolation processing principle is:

Be located in the discrete spectrum frequency f ₀Near local amplitude is maximum to be respectively k with time maximum spectral line ₁And k ₂Root, establishing these two spectral line amplitudes is respectively y ₁And y ₂Consider 0≤k ₀-k ₁≤1, the definition alpha does

α＝k ₀-k ₁-0.5，-0.5≤α≤0.5

Y then ₁And y ₂Can be rewritten as the function about α, the definition factor beta does

$\mathrm{\&beta;}=\frac{y_2-y_1}{{\mathrm{<figure-callout\; id="2"\; label="y"\; filenames="HSA00000596614000012.png"\; state="\{\{state\}\}">y</figure-callout>}}_2+{\mathrm{<figure-callout\; id="1"\; label="y"\; filenames="HSA00000596614000012.png"\; state="\{\{state\}\}">y</figure-callout>}}_1}$

With the y in the following formula ₁And y ₂Replace, then β can be written as the function about α:

$\mathrm{\&beta;}=\frac{\left|W_{N-p}\right[2\mathrm{\&pi;}(-\mathrm{\&alpha;}+0.5)/N\left]\right|-\left|W_{N-p}\right[2\mathrm{\&pi;}(-\mathrm{\&alpha;}-0.5)/N\left]\right|}{\left|W_{N-p}\right[2\mathrm{\&pi;}(-\mathrm{\&alpha;}-0.5)/N\left]\right|+\left|W_{N-p}\right[2\mathrm{\&pi;}(-\mathrm{\&alpha;}+0.5)/N\left]\right|}=g\left(\mathrm{\&alpha;}\right)$

In α ∈ [0.5,0.5] scope,, utilize each α with the periodic sampling some spots _iValue calculates corresponding β _i, confirm to treat fitting of a polynomial number of times K, utilize and find the solution coefficient a on the formula _k(k=0,1, L K-1), thereby sets up frequency spectrum interpolation polynomial expression S _K(x), can be designated as

α＝g ^-1(β)≈a ₁β+a ₂β ²+a ₃β ³+L+a _Kβ ^K

In the formula, a ₁, a ₂..., a _KCoefficient for polynomial fitting.Thereby realize the discrete spectrum interpolation processing.

Claims (4)

1. the harmonic electric energy metering method of a Nuttall self-convolution window weighting Fourier transform is characterized in that, this method comprises the steps:

At first make up discrete Nuttall self-convolution window; Adopting length is that the Nuttall self-convolution window sequence of N is that the discrete electrical force signal of N computes weighted to length; Adopt Fourier transform to obtain the discrete spectrum of electric power signal again; The utilization least square method is carried out the spectrum peak match, obtains amplitude, frequency and the initial phase angle of electric harmonic, and then calculates harmonic electric energy.

2. the harmonic electric energy metering method of Nuttall self-convolution window weighting Fourier transform according to claim 1; It is characterized in that the said structure Nuttall self-convolution window that disperses is that to make up length earlier be the discrete Nuttall window sequence of M; Then p discrete Nuttall window sequence made convolution algorithm p-1 time; Obtain convolution sequence, zero in the head or p-1 of the tail benefit of convolution sequence, promptly obtaining length is the discrete Nuttall self-convolution window in p rank of N=pM.

3. the harmonic electric energy metering method of Nuttall self-convolution window weighting Fourier transform according to claim 1; It is characterized in that; Said utilization least square method is carried out the spectrum peak match, is near the maximum and inferior big spectral line of the searching amplitude corresponding harmonic frequency of discrete spectrum, and the utilization least square method is carried out match to the peak value spectral line; Obtain amplitude, frequency and the initial phase angle of electric harmonic, and then calculate harmonic electric energy.

4. the harmonic electric energy metering method of Nuttall self-convolution window weighting Fourier transform according to claim 1; It is characterized in that; Said Nuttall self-convolution window weighting Fourier transform is the increase along with the convolution exponent number, and the sidelobe level of Nuttall self-convolution window is reduced rapidly.

Images