Method for measuring power frequency phasor of urban substation

Abstract

The invention discloses a method for measuring the power frequency phasor of an urban substation, and the method comprises the following steps: S1, building a mathematical model of a to-be-analyzed component in a power supply system of the urban substation, and inputting a to-be-measured signal component into the mathematical model; s2, fitting the power frequency phasor parameters by using linear Bessel interpolation to obtain a fundamental frequency, an amplitude and a phase; s3, compensating a frequency spectrum negative image through an iteration process by using a fitting result in the S2, and reducing the frequency spectrum leakage of the negative image; s4, constructing an energy threshold inequality, judging whether the ratio of the energy of the inter-harmonic component to the component of the signal to be detected meets a threshold or not according to the result of the S3, and if so, directly outputting the result; and if not, repeating the steps S2-S4 until the threshold value is met.

Description

Technical Field

[0001] The invention relates to the technical field of signal and data analysis, and in particular to a method for measuring power frequency phasors of urban substations. Background Art

[0002] At present, the access of a high proportion of renewable energy and power electronic equipment not only generates a large amount of harmonic current in the power system, but also generates a large amount of interharmonic interference, which seriously affects the measurement of power frequency phasors. The phasor measurement unit (PMU) is a basic tool for estimating synchronous phasors in modern power systems. The estimation results obtained by PMU are established in IEEE / IEC 60255-118-1.

[0003] In recent years, the accurate measurement of power frequency phasors in urban substations in multiple disturbance components has attracted much attention. This work proposes a power frequency phasor measurement technique with low response time and high accuracy, and meets the compliance of Class M PMU.

[0004] Traditional power frequency phasor measurement methods can be roughly divided into two categories: those based on Fourier transform and those not based on Fourier transform. Methods based on Fourier transform, such as interpolation discrete Fourier transform (IpDFT) and Taylor-Fourier transform (TFT), can obtain estimated values ​​of power frequency phasors to a certain extent, but when applied in the multi-disturbance environment of urban substations, the above methods cannot obtain satisfactory results. Methods not based on Fourier transform, such as Prony method and Kalman filter, usually have strict requirements on the model order and cannot dynamically track the power frequency phasors under multi-disturbance conditions. Summary of the invention

[0005] The object of the present invention is to provide a method for measuring power frequency phasors of urban substations, which is based on linear Bessel interpolation. Under the condition of considering multi-component interference such as harmonics and interharmonics, the power frequency phasors can be measured with a frequency resolution of 12.5 Hz in accordance with Class M PMU compliance through iterative process analysis.

[0006] To achieve the above object, the present invention provides the following technical solutions:

[0007] A method for measuring power frequency phasors in urban substations can measure power frequency phasors with a frequency resolution of 12.5 Hz in accordance with Class M PMU compliance, while taking into account multi-component interference such as harmonics and interharmonics. The method includes the following steps:

[0008] S1. Establish a mathematical model of the components to be analyzed in the power supply system of the urban substation, and input the signal components to be measured into the mathematical model;

[0009] S2, using linear Bessel interpolation to fit the power frequency phasor parameters to obtain the fundamental frequency, amplitude and phase;

[0010] S3, using the fitting result in S2, compensating the spectrum negative image through an iterative process to reduce the spectrum leakage of the negative image;

[0011] S4, construct an energy threshold inequality, and based on the result of S3, determine whether the ratio of the energy of the interharmonic component to the signal component to be measured meets the threshold. If so, directly output the result; if not, repeat S2-S4 until the threshold is met, and finally measure and obtain the power frequency vector.

[0012] Preferably, in S1, the mathematical model is:

[0013]

[0014] Where, h = 1, 2, ..., H is the highest order of harmonics; a h , are the amplitude and phase of the harmonic respectively; N is the data window length and is set to an even number; n is the sequence number of the discrete sampling point.

[0015] Preferably, in S1, usually, the signal component to be measured needs to be windowed to limit the calculation time, and the main lobe and side lobe characteristics of the window function are used to limit the spectrum leakage. The signal component to be measured is windowed using a window function, and after discrete Fourier transformation, the representation of the signal component to be measured in the frequency domain is obtained:

[0016]

[0017] Wherein, m = [0, N-1]; w(n) is the Kaiser window function; W(n) is the discrete Fourier transform form of the Kaiser window function; x(n) is the windowed signal.

[0018] Preferably, S2 comprises the following steps:

[0019] S21, establishing a relationship between the correction coefficient ξ, the phasor parameters (amplitude and phase) and the linear polynomial;

[0020] S22, solving the linear polynomial to obtain the correction coefficient ξ, amplitude and phase;

[0021] S23. Using the correction coefficient ξ, the fundamental frequency is obtained by solving the relationship between the correction coefficient ξ and the fundamental frequency.

[0022] Preferably, in S21, the relationship between the correction coefficient ξ, the phasor parameter and the linear polynomial is:

[0023]

[0024] a1≈P x (ξ)(|X w (m-σ)|+|Xw (m)|+|X w (m+σ)|)|m=m m ;

[0025]

[0026] The selection of “±” and the value of σ in the above relationship are determined by the position of the second largest spectral line. If the index number corresponding to the second largest spectral line position is m m +1, take "+", σ=1; if the index number is m m -1, take “-”, σ=-1.

[0027] Where P d , P x and P p are linear polynomials for calculating fundamental frequency, amplitude and phase respectively; m m M is the index number corresponding to the maximum spectrum line obtained by DFT; f is the frequency polynomial order; p mf is the polynomial coefficient; σ is a constant, taking the value of -1 or 1; m = [0, N-1], N is the data window length, and is set to an even number; a1 is the amplitude corresponding to the power frequency phasor; is the phase corresponding to the industrial frequency phasor.

[0028] Preferably, in S22, the linear polynomial needs to rely on the fundamental frequency dependency F fr , Amplitude dependence A f , phase dependence P f and the fitting factor F f The fitting solution is obtained;

[0029] Specifically, the frequency depends on the fundamental frequency F fr and the fitting factor F f Carry out M f The polynomial P is obtained by fitting the polynomial d , by solving the polynomial P d Get the correction coefficient ξ;

[0030] The amplitude a1 corresponding to the power frequency phase quantity is determined by the amplitude dependence A f and the fitting factor F f Carry out M a The polynomial P is obtained by fitting the polynomial x , substitute the correction coefficient ξ into the polynomial P x The measured value of the amplitude a1 corresponding to the power frequency phasor is obtained;

[0031] The phase corresponding to the power frequency phasor Through the phase dependence P f and the fitting factor F f Carry out Mp The polynomial P is obtained by fitting the polynomial p , substitute the correction coefficient ξ into the polynomial P p The phase corresponding to the power frequency phasor is obtained The measured value of

[0032] in,

[0033]

[0034] F f (n) = w(n)·N / 2π;

[0035]

[0036] p f (m) = arg{W(m)};

[0037]

[0038] Where n bin is the number of frequency bands required; W(m) is the spectrum of the Kaiser window; w(n) is the Kaiser window; I0(β) is the first kind of Bessel function; M a is the order of the amplitude polynomial; M p is the order of the phase polynomial.

[0039] Preferably, in S23, the relationship between the correction coefficient ξ and the fundamental frequency is:

[0040] f1=(m m ±ξ)Δf;

[0041] The selection of “±” and the value of σ in the relationship are determined by the position of the second largest spectral line. If the index number corresponding to the second largest spectral line position is m m +1, take "+", σ=1; if the index number is m m -1, take "-", σ = -1;

[0042] Where f1 is the fundamental frequency; Δf is the frequency resolution, Δf=1 / N.

[0043] Preferably, S3 comprises the following steps:

[0044] S31. The discrete Fourier transform (DFT) of the signal base frequency can be represented based on the positive and negative images of the spectrum:

[0045]

[0046] S32, compensating the spectrum negative image by iteration:

[0047]

[0048] In the formula, is the spectrum positive image; is the negative image of the spectrum; j is the number of iterations, j = 1, 2, ..., J; PID is the correction coefficient ξ, the amplitude a1 corresponding to the power frequency phasor, and the phase corresponding to the power frequency phasor A simplified representation of the calculation process.

[0049] Preferably, in S4, the result is made to satisfy the energy threshold inequality through an iterative process. The iterative process ends until the threshold is satisfied, and the calculation result of the power frequency phasor is returned. The energy threshold inequality is:

[0050]

[0051] in,

[0052] After eliminating the interference of the negative image, the interference of the intermediate harmonic of the measured component on the power frequency phasor measurement result is reduced through L iterations. The initial value of the iteration is:

[0053]

[0054] In the formula, is the interharmonic component; is the threshold; J is the maximum number of iterations of spectral negative image compensation.

[0055] Compared with the prior art, the present invention has the following beneficial effects:

[0056] The present invention is aimed at the multi-disturbance scenario of urban substations and aims to accurately measure the power frequency phasor parameters of urban substations. The present invention constructs a linear Bessel interpolation method based on the first-class Bessel function, and the proposed interpolation method can be used for power frequency phasor measurement in multi-disturbance scenarios of urban substations. Compared with the prior art, the main calculation steps of the present invention are less affected by frequency domain fluctuations through the fitting of linearized Bessel functions. Therefore, the method has good stability and is less affected by harmonics and interharmonics.

[0057] In summary, compared with the prior art, the present invention can significantly reduce the spectrum leakage effect, the interference of interharmonics on the measurement results and the mutual coupling between the two disturbances, and the measurement accuracy of the power frequency phasor amplitude and frequency is greatly improved. BRIEF DESCRIPTION OF THE DRAWINGS

[0058] Figure 1 It is a schematic diagram of the process of the present invention;

[0059] Figure 2 A schematic diagram of the relationship between positive and negative images of a spectrum according to an embodiment of the present invention;

[0060] Figure 3Schematic diagram of the results of the out-of-band interference test under static conditions of an embodiment of the present invention; wherein (a) is the total vector error TVE test result, (b) is the absolute frequency error FE test result, and (c) is the absolute frequency change rate error RFE test result. DETAILED DESCRIPTION

[0061] The technical solutions in the embodiments of the present invention will be clearly and completely described below in conjunction with the drawings in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, rather than all the embodiments; based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making creative work are within the scope of protection of the present invention.

[0062] See also Figure 1-Figure 3 , a method for measuring power frequency phasors of urban substations, comprising the following steps:

[0063] S1. Establish a mathematical model of the components to be analyzed in the power supply system of the urban substation, and input the signal components to be measured into the mathematical model, wherein the signal components to be measured include power frequency components and harmonics.

[0064]

[0065] Where, h = 1, 2, ..., H is the highest order of harmonics; a h , are the amplitude and phase of the harmonic respectively; N is the data window length and is set to an even number; n is the sequence number of the discrete sampling point.

[0066] Usually, the signal components to be measured need to be windowed to limit the calculation time of the algorithm, and the characteristics of the main lobe and side lobe of the window function are used to limit the spectrum leakage.

[0067] The windowed signal can be expressed as:

[0068] x w (n) = x(n)w(n);

[0069] The windowed signal is transformed through the discrete Fourier transform (DFT) process to obtain the representation of the signal component to be measured in the frequency domain:

[0070]

[0071] Wherein, m = [0, N-1]; w(n) is the Kaiser window function; W(n) is the discrete Fourier transform form of the Kaiser window function; x(n) is the windowed signal.

[0072] S2. Use linear Bessel interpolation to fit the power frequency phasor parameters to obtain the fundamental frequency, amplitude and phase.

[0073] Under the condition of multi-component aliasing, the choice of window function will affect the performance of the algorithm. The Kaiser window can change the proportion of the main lobe and side lobe by adjusting the window parameters, which is helpful for feature extraction in the process of power frequency phasor measurement and suppresses the spectrum leakage effect to a certain extent.

[0074] The calculation expression of the Kaiser window is:

[0075]

[0076]

[0077] Where I0(·) is the first kind Bessel function; β is the Kaiser window parameter; b is the sequence number of the series in the first kind Bessel function; G is the gamma function.

[0078] Due to the incoherent sampling in actual situations, the fundamental frequency of the signal component to be measured is located between the two consecutive highest spectral lines. What needs to be solved urgently is to determine a correction coefficient ξ to accurately locate the position of the fundamental frequency. The correction coefficient calculation method in the prior art is usually based on multi-spectral line interpolation based on a cosine family window, but it is not applicable to the first kind of Bessel function. The embodiment of the present invention uses a numerical method to process a class of Bessel functions to obtain the corresponding phasor parameters.

[0079] First, establish the correction coefficient ξ, the amplitude a1 corresponding to the power frequency phasor, and the phase corresponding to the power frequency phasor The relationship between the linear polynomial, the correction coefficient ξ, the amplitude a1 corresponding to the power frequency phasor, and the phase corresponding to the power frequency phasor The value of can be calculated by the number of the largest spectral line:

[0080]

[0081] a1≈P x (ξ)(|X w (m-σ)|+|X w (m)|+|X w (m+σ)|)|m=m m ;

[0082]

[0083] The selection of “±” and the value of σ in the above relationship are determined by the position of the second largest spectral line. If the index number corresponding to the second largest spectral line position is m m +1, take "+", σ=1; if the index number is m m -1, take “-”, σ=-1.

[0084] Where P d , P x and Pp are linear polynomials for calculating fundamental frequency, amplitude and phase respectively, obtained by spectrum fitting; m m M is the index number corresponding to the maximum spectrum line obtained by DFT; f is the frequency polynomial order; p mf are polynomial coefficients; σ is a constant, taking the value of -1 or 1; m = [0, N-1], N is the data window length, and is set to an even number.

[0085] Secondly, the linear polynomial is solved to obtain the correction coefficient ξ, amplitude and phase.

[0086] As shown in Table 1, the linear polynomial needs to rely on the fundamental frequency dependence F fr , Amplitude dependence A f , phase dependence P f and the fitting factor F f The fitting solution is obtained.

[0087] Table 1 Correspondence table of dependency, fitting factors and polynomials

[0088]

[0089] According to Table 1, the frequency depends on the fundamental frequency F fr and the fitting factor F f Carry out M f The polynomial P is obtained by fitting the polynomial of (frequency polynomial order) order. d , by solving the polynomial P d Get the correction coefficient ξ.

[0090] The amplitude a1 corresponding to the power frequency phase quantity is determined by the amplitude dependence A f and the fitting factor F f Carry out M a The polynomial P is obtained by fitting the polynomial of (amplitude polynomial order) order. x , substitute the correction coefficient ξ into the polynomial P x The measured value of the amplitude a1 corresponding to the power frequency phasor is obtained.

[0091] The phase corresponding to the power frequency phasor Through the phase dependence P f and the fitting factor F f Carry out M p The polynomial P is obtained by fitting the polynomial of (phase polynomial order) order. p , substitute the correction coefficient ξ into the polynomial P p The phase corresponding to the power frequency phasor is obtained The measured value.

[0092] in,

[0093]

[0094] F f (n) = w(n)·N / 2π;

[0095]

[0096] p f (m) = arg{W(m)};

[0097]

[0098] Where n bin is the required number of frequency bands; W(m) is the spectrum of the Kaiser window; w(n) is the Kaiser window; I0(β) is the first kind Bessel function.

[0099] Third, using the correction coefficient ξ, the fundamental frequency is obtained by solving the relationship between the correction coefficient ξ and the fundamental frequency;

[0100] The relationship between the correction coefficient ξ and the fundamental frequency is:

[0101] f1=(m m ±ξ)Δf;

[0102] Where f1 is the fundamental frequency; Δf is the frequency resolution, Δf=1 / N.

[0103] Substitute the correction coefficient ξ calculated by the above steps into the equation to solve for the fundamental frequency f1. The selection of "±" and the value of σ in the equation are determined according to the position of the second largest spectral line. If the index number corresponding to the second largest spectral line position is m m +1, take "+", σ=1; if the index number is m m -1, take “-”, σ=-1.

[0104] S3. Using the fitting result in S2, the spectrum negative image is compensated through an iterative process to reduce the spectrum leakage of the negative image.

[0105] like Figure 2 As shown, according to the convolution theorem, the characteristic of the windowed signal in the frequency domain is the transformation of the window function spectrum in the spatial position.

[0106] Therefore, the discrete Fourier transform (DFT) of the signal fundamental frequency can be represented based on the positive and negative images of the spectrum:

[0107]

[0108] Assuming that the maximum number of iterations used to compensate the negative image is J, the spectrum negative image can be expressed as follows:

[0109]

[0110] in,

[0111]

[0112] The iterative process of compensating the spectral negative image can be described as follows:

[0113]

[0114] In the formula, is the spectrum positive image; is the negative image of the spectrum; j is the number of iterations, j = 1, 2, ..., J; PID is the correction coefficient ξ, the amplitude a1 corresponding to the power frequency phasor, and the phase corresponding to the power frequency phasor A simplified representation of the calculation process.

[0115] The specific iteration process can be represented by the following table:

[0116] Table 2 Pseudo code table of spectrum negative image compensation based on PID

[0117]

[0118]

[0119] S4, construct an energy threshold inequality, and based on the result of S3, determine whether the ratio of the energy of the interharmonic component to the signal component to be measured meets the threshold. If so, there is no need to consider the interharmonic component and the result is output directly; if not, return to S2 to construct the interharmonic component and filter it out, and then perform S3-S4 until the threshold is met.

[0120] After eliminating the interference of the negative image, the interference of the intermediate harmonic of the measured signal component on the power frequency phasor measurement result is reduced through L iterations. The initial value of the iteration is:

[0121]

[0122] In the formula, is the constructed interharmonic component; is the threshold; J is the maximum number of iterations for compensating the negative image.

[0123] During the iteration process, when the constructed interharmonic components The energy of the measured signal component X w The ratio of (m) exceeds the threshold When , the accuracy of the result is not acceptable. Through L iterations, the result can satisfy the energy threshold inequality. When the result satisfies the energy threshold inequality, the iteration process ends, and the power frequency phasor of the substation is finally measured.

[0124] The energy threshold inequality is:

[0125]

[0126] The specific iterative calculation process is as follows:

[0127] Table 3 Energy threshold iteration pseudo code table

[0128]

[0129]

[0130] In the formula, E e is the interharmonic component energy; f OOB 、a OOB , are the frequency, amplitude and phase of the approximate interharmonics respectively; iPID is the iterative process of PID and the compensated negative image.

[0131] In a specific embodiment, in order to test the detection performance of the present invention, the basic estimation accuracy of the signal component by the evaluation algorithm is selected.

[0132] The single frequency signal is set to: frequency of 50Hz, interharmonic frequencies of 10~25Hz and 75~100Hz respectively. In addition, Gaussian white noise is added to the signal component to be measured, and the signal-to-noise ratio (SNR) is 60dB. Each method is calculated 10 times according to its nominal period, and the reporting rate is 50fps. The maximum values ​​of the total vector error (TVE), absolute frequency error (FE) and absolute frequency rate error (RFE) of the estimation results are shown in Figure 3 As shown in the figure, the red horizontal line is the measurement compliance of the M-type PMU in IEEE / IEC 60255-118-1. The method of this embodiment is represented by A, and the advanced algorithms in the prior art: TFT, MW-FIR, CS-TFM, iIpDFT and iIpDFT-dc are considered as comparison algorithms.

[0133] The maximum TVE, FE and RFE obtained by the present invention are 0.6667%, 0.0081 Hz and 0.0265 Hz / s respectively.

[0134] The maximum TVE, FE and RFE obtained by iIpDFT are 1.5757%, 0.0139Hz, 0.0544Hz / s respectively; the maximum TVE, FE and RFE obtained by TFT are 4.4162%, 0.0196Hz, 0.1628Hz / s respectively; the maximum TVE, FE and RFE obtained by CS-TFM are 1.2724%, 0.0095Hz, 0.0673Hz / s respectively; the maximum TVE, FE and RFE obtained by MW-FIR are 2.3331%, 0.0144Hz, 0.1011Hz / s respectively; the maximum TVE, FE and RFE obtained by iIpDFT-dc are 1.6267%, 0.0152Hz, 0.1192Hz / s respectively. Compared with the prior art, the errors of the present invention are all smaller than those of the existing algorithms, highlighting the superiority of the present invention in measurement accuracy.

[0135] The test signal is as follows:

[0136]

[0137] Where, f = 50HZ; f i is the frequency of the interharmonic; are the phases of the fundamental wave and interharmonics, respectively, and are set to random numbers ranging from 0 to 2π.

[0138] The calculation expression of the total vector error is:

[0139]

[0140] The calculation expression of the absolute frequency error is:

[0141] FE(n)=f measured (n)-f ref (n);

[0142] The calculation expression of the absolute frequency change rate error is:

[0143]

[0144] In the formula, are the real and imaginary PMU estimates at reporting time n, respectively; X r (n), X i (n) are the reference values ​​of the real and imaginary parts of the phasor at the reporting time n; n is the number corresponding to the reporting time; f measured (n), f ref (n) are the measured frequency and frequency reference value respectively.

Claims

1. A method for measuring power frequency phasors in urban substations, characterized in that: The following steps are involved: S1. Establish a mathematical model of the components to be analyzed in the power supply system of the urban substation, and input the signal components to be measured into the mathematical model; S2. Use linear Bessel interpolation to fit the power frequency phasor parameters to obtain the fundamental frequency, amplitude and phase; S3, using the fitting result in S2, compensating the spectrum negative image through an iterative process to reduce the spectrum leakage of the negative image; S4, construct an energy threshold inequality, and based on the result of S3, determine whether the ratio of the energy of the interharmonic component to the signal component to be measured meets the threshold. If so, directly output the result; if not, repeat S2-S4 until the threshold is met, and finally measure and obtain the power frequency vector.

2. A method for measuring power frequency phasors of urban substations according to claim 1, characterized in that: In S1, the mathematical model is: Where, h = 1, 2, ..., H is the highest order of harmonics; a h , are the amplitude and phase of the harmonic respectively; N is the data window length and is set to an even number; n is the sequence number of the discrete sampling point.

3. A method for measuring power frequency phasors of urban substations according to claim 2, characterized in that: In S1, the signal component to be measured is windowed using a window function, and after discrete Fourier transform, the representation of the component to be measured in the frequency domain is obtained: Wherein, m = [0, N-1]; w(n) is the Kaiser window function; W(n) is the discrete Fourier transform form of the Kaiser window function; x(n) is the windowed signal.

4. A method for measuring power frequency phasors of urban substations according to claim 1, characterized in that: S2 includes the following steps: S21, establishing a relationship between the correction coefficient ξ, the phasor parameter and the linear polynomial; S22, solving the linear polynomial to obtain the correction coefficient ξ, amplitude and phase; S23. Using the correction coefficient ξ, the fundamental frequency is obtained by solving the relationship between the correction coefficient ξ and the fundamental frequency.

5. A method for measuring power frequency phasors of urban substations according to claim 4, characterized in that: In S21, the relationship between the correction coefficient ξ, the phasor parameter and the linear polynomial is: a1≈P x (ξ)(|X) w (m-σ)|+|X W (m)|+|X w (m+σ)|)|m=m m ; The selection of "±" and the value of σ in the above relationship are determined by the position of the second largest spectral line. If the index number corresponding to the second largest spectral line position is m m +1, take "+", σ=1; if the index number is m m -1, take "-", σ = -1; Where P d , P x and P p are linear polynomials for calculating fundamental frequency, amplitude and phase respectively; m m N is the index number corresponding to the maximum spectrum line obtained by DFT; f is the frequency polynomial order; are polynomial coefficients; v is a constant, taking the value of -1 or 1; m = [0, N-1], N is the data window length, and is set to an even number; a1 is the amplitude corresponding to the power frequency phasor; is the phase corresponding to the industrial frequency phasor.

6. A method for measuring power frequency phasors of urban substations according to claim 4, characterized in that: In S22, the linear polynomial needs to rely on the fundamental frequency dependency F fr , Amplitude dependence A f , phase dependence P f and the fitting factor F f The fitting solution is obtained; Specifically, the frequency depends on the fundamental frequency F fr and the fitting factor F f Carry out M f The polynomial P is obtained by fitting the polynomial d , by solving the polynomial P d Get the correction coefficient ξ; The amplitude a1 corresponding to the power frequency phase quantity is determined by the amplitude dependence A f and the fitting factor F f Carry out M a The polynomial P is obtained by fitting the polynomial x , substitute the correction coefficient ξ into the polynomial P x The measured value of the amplitude a1 corresponding to the power frequency phasor is obtained; The phase corresponding to the power frequency phasor Through the phase dependence P f and the fitting factor F f Carry out M p The polynomial P is obtained by fitting the polynomial p , substitute the correction coefficient ξ into the polynomial P p The phase corresponding to the power frequency phasor is obtained The measured value of in, F f (n)=w(n)·N / 2π; P f (m)=arg{W(m)}; Where n bin is the number of frequency bands required; W(m) is the spectrum of the Kaiser window; w(n) is the Kaiser window; I0(β) is the first kind of Bessel function; M a is the order of the amplitude polynomial; M p is the order of the phase polynomial.

7. A method for measuring power frequency phasors of urban substations according to claim 4, characterized in that: In S23, the relationship between the correction coefficient ξ and the fundamental frequency is: f1=(m m ±ξ)Δf; The selection of "±" and the value of σ in the relationship are determined by the position of the second largest spectral line. If the index number corresponding to the second largest spectral line position is m m +1, take "+", σ=1; if the index number is m m -1, take "-", σ = -1; Where f1 is the fundamental frequency; Δf is the frequency resolution, Δf=1 / N.

8. A method for measuring power frequency phasors of urban substations according to claim 1, characterized in that: S3 includes the following steps: S31, the discrete Fourier transform DFT of the signal base frequency is represented by positive and negative images based on the spectrum: S32, compensating the spectrum negative image by iteration: In the formula, is the spectrum positive image; is the negative image of the spectrum; j is the number of iterations, j = 1, 2, ..., J; PID is the correction coefficient ξ, the amplitude a1 corresponding to the power frequency phasor, and the phase corresponding to the power frequency phasor A simplified representation of the calculation process.

9. A method for measuring power frequency phasors of urban substations according to claim 1, characterized in that: In S4, the energy threshold inequality is: in, The initial value of the iteration is: In the formula, is the interharmonic component; is the threshold; J is the maximum number of iterations of spectral negative image compensation.

Images