Journal of Modern Power Systems and Clean Energy

ISSN 2196-5625 CN 32-1884/TK

网刊加载中。。。

使用Chrome浏览器效果最佳,继续浏览,你可能不会看到最佳的展示效果,

确定继续浏览么?

复制成功,请在其他浏览器进行阅读

Parameter Estimation of Sub-/Super-synchronous Oscillation Based on Interpolated All-phase Fast Fourier Transform with Optimized Window Function  PDF

  • Bo Sun 1
  • Xi Wu 1 (Senior Member, IEEE)
  • Xi Chen 1
  • Zixiao Zou 1
  • Qiang Li 2
  • Bixing Ren 2
1. School of Electrical Engineering, Southeast University, Nanjing, China; 2. State Grid Jiangsu Electric Power Company Limited Research InstituteNanjing, China

Updated:2024-07-26

DOI:10.35833/MPCE.2023.000179

  • Full Text
  • Figs & Tabs
  • References
  • Authors
  • About
CITE
OUTLINE

Abstract

In recent years, with increasing amounts of renewable energy sources connecting to power networks, sub-/super-synchronous oscillations (SSOs) have occurred more frequently. Due to the time-variant nature of SSO magnitudes and frequencies, as well as the mutual interferences among SSO modes with close frequencies, the accurate parameter estimation of SSO has become a particularly challenging topic. To solve this issue, this paper proposes an improved spectrum analysis method by improving the window function and a spectrum correction method to achieve higher precision. First, by aiming at the sidelobe characteristics of the window function as evaluation criteria, a combined cosine function is optimized using a genetic algorithm (GA). Furthermore, the obtained window function is self-convolved to extend its excellent characteristics, which have better performance in reducing mutual interference from other SSO modes. Subsequently, a new form of interpolated all-phase fast Fourier transform (IpApFFT) using the optimized window function is proposed to estimate the parameters of SSO. This method allows for phase-unbiased estimation while maintaining algorithmic simplicity and expedience. The performance of the proposed method is demonstrated under various conditions, compared with other estimation methods. Simulation results validate the effectiveness and superiority of the proposed method.

I. Introduction

IN recent years, with the large-scale grid connection of renewable energy resources such as wind power and photovoltaic, the proportion of power electronics in power systems is increasing. Sub-/super-synchronous oscillations (SSOs) have occurred frequently. According to the involved elements, SSO can be classified as sub-/super-synchronous resonance (SSR), sub-synchronous torsional interaction (SSTI), and sub-synchronous control interaction (SSCI) [

1]. SSR is a phenomenon that a series compensated system exchanges significant energy with a turbine generator at a frequency below the synchronous frequency [2]. SSCI is a pure electrical oscillation phenomenon caused by interactions between control systems and networks. The adverse effect of high-voltage direct current (HVDC) rectifiers on the sub-synchronous torsional oscillations, referred to as SSTI [3]. Recent SSO events were caused by the interaction between wind turbine clusters and series-compensated transmission lines in southern Texas of the USA and in Guyuan of China [4], [5]. During these incidents, SSO became so strong that it caused serious damage to the rotor shaft system of turbine generators or even weakened the safety and reliability of power systems. These accidents have shown that SSOs have characteristics of time-variant and mode superposition, posing great challenges for detecting SSOs [6]. Therefore, it is imperative and significant to identify the parameters of SSO such as frequency, amplitude, and phase in order to establish control strategies and mitigate SSOs.

Phasor measurement units (PMUs) are widely used in power systems for synchronous phasor estimation. However, traditional PMUs are primarily designed to measure fundamental phasors [

7]. To accurately obtain the fundamental phasors, PMU usually filters out sub-/super-synchronous components as interferences or noise and only retains the fundamental wave. Consequently, the sub-/super-synchronous information can hardly be obtained through traditional PMUs. Therefore, it is necessary to modify the filtering setting to capture sub-/super-synchronous components. Similar modifications have been done in [8]. What’s more, compared to purpose-built protective relaying intelligent electronic devices (IEDs), PMUs have been widely adopted in power systems, and modifying them for sub-/super-synchronous component detection can yield cost savings. The algorithm of modified PMU has a wider detection range. It can not only detect the fundamental frequency but also the SSO frequency.

Discrete Fourier transform (DFT) based methods are the most commonly employed in the identification and monitoring of SSOs, owing to their straightforward algorithms, rapid computation speeds, and ease of implementation in embedded devices [

9]-[11]. However, DFT can be subjected to fence effect and spectral leakage, leading to decreased detection accuracy. To overcome these issues, some researchers have improved the traditional algorithms. Iterative [12], window, and interpolated DFT (IpDFT) algorithms have been proposed [13], [14]. In [15], an iterative Taylor-Fourier multi-frequency parameter estimation model of SSO is proposed. However, its computation burden increases dramatically with the number of oscillation modes. Window techniques are also used to suppress spectral leakage [16]. Numerous window types have been thoroughly categorized and explained [17]. A well-designed window function with a small peak sidelobe level and a large asymptotic sidelobe attenuation rate can be effective in suppressing spectral leakage. In [18], it only considered the generalized maximum asymptotic sidelobe attenuation window, which means that the effect of optimization performance is limited, and the spectrum of SSO mode will be significantly interfered by the spectral leakage from other close SSO modes. In addition to improving the windows, increasing the number of interpolated spectral lines can also enhance the accuracy of DFT. Currently, two-spectral line, three-spectral line, and four-spectral line IpDFT algorithms are used to reduce the spectral interference [19]-[21]. In [22], a three-point IpDFT algorithm was proposed to capture the sub-/sup-synchronous dynamics. However, the accuracy of phase measurements is limited by both white noise and spectral interference. In the domain of information processing, a novel signal analysis algorithm named all-phase FFT (ApFFT) has been introduced in the past few years [23]-[25]. Due to its favorable phase characteristics and enhanced leakage suppression, it offers the potential to further improve the precise parameter estimation of SSO.

Currently, DFT-based methods require enhancements in window techniques to suppress spectral leakage, thereby improving the accuracy of SSO detection. Concerning phase measurement, it is more susceptible to noise and spectral interference, resulting in less precise data for monitoring. Therefore, it is imperative to develop a more accurate method for SSO detection.

This paper proposes an improved interpolated ApFFT (IpApFFT) based on the optimized combined cosine self-convolution window (OCCSCW) function to estimate the parameters of SSO. Compared with traditional IpDFT, the proposed IpApFFT with OCCSCW function in this paper has the following advantages.

1) Our significant contribution lies in the optimization of window parameters using sidelobe peak and sidelobe attenuation rate as objective functions. Superior to traditional windows, the optimized OCCSCW function greatly improves the ability to suppress spectral leakage. Furthermore, it reduces the mutual interference of close SSO modes on the parameter estimation of SSO. Thus, it significantly improves the precision of parameter estimation for SSO.

2) An ApFFT correction algorithm based on spectral line interpolation is proposed. This algorithm provides improved accuracy in suppressing spectral leakage compared with traditional DFT, and effectively reduces errors induced by the fence effect. Additionally, the ApFFT has the phase-invariance property, resulting in zero error in phase estimation.

3) Simulation tests with close frequency of SSO modes have been conducted, and the results demonstrate that the proposed method outperforms IpDFT under interharmonic amplitude modulation (AM), phase modulation (PM), frequency ramp (FR), and harmonic interference (HI) conditions. It exhibits good reliability and stability in the presence of interference from other close SSO modes.

This paper is structured as follows. First, Section II describes an optimized combined cosine (OCC) window function based on the genetic algorithm (GA) in detail. Next, the principles of the ApFFT spectrum analysis technique are introduced in Section III. Then, the IpApFFT based on the OCCSCW is proposed to estimate the frequency, amplitude, and phase of SSO in Section IV. Section V presents verification studies that demonstrate the advantages of the proposed method under conditions of sub-/super-synchronous parameter variations, noise, dynamic conditions, and simulation test. Finally, the main conclusions of the proposed method are summarized in Section VI.

II. OCC Window Function Based on GA

In this section, the combined cosine window function is introduced and the parameter of the window function is optimized using a GA to obtain the optimized window function. Subsequently, the self-convolution window function can be obtained. It is noteworthy that the parameter optimization of window function is an offline preparatory work that only needs to be optimized once, and does not occupy any computation time.

A. Combined Cosine Window Function

Applying a window function to sampled data is an effective method to suppress spectral leakage. The performance of the window function in suppressing spectral leakage is heavily influenced by its sidelobe characteristics. Different window functions have different suppression effects. It is crucial to find a new window function, which is more effective than most window functions.

The combined cosine window consists of various terms and coefficients that lead to different window functions with different performances. By changing the terms and coefficients, it is possible to find a better target window function. The combined cosine window comprises several terms and coefficients, which can generate different window functions with different performances. Altering these terms and coefficients enables identification of a more optimal target window function. The combined cosine window function in its general time-domain form can be rewritten as:

ω(n)=m=0M-1(-1mamcos(2πnm/N) (1)

where n=0,1,,N-1, and N is the length of window function; m=0,1,,M-1, and M is the number of window function terms; and am  is the coefficient which needs to meet the following constraints.

m=0M-1(-1)mam=0m=0M-1am=1 (2)

The spectrum function of the combined cosine window function can be written as:

W(ω)=m=0M-1(-1)mam2WRω-2πmN+WRω+2πmN (3)

where WR(ω)=(sinωN/2)/(sinω/2)e-jω(N-1)/2, which is a spectrum function of a rectangular window.

B. Parameter Optimization of Combined Cosine Window Function

The peak sidelobe level and asymptotic sidelobe attenuation rate serve as crucial metrics for evaluating the performance of window functions. Window functions characterized by low peak sidelobe levels and high asymptotic sidelobe attenuation rates contribute to enhance the accuracy in signal parameter estimation [

26]. To find such a window function, the GA can be used.

The first step in using GA is to determine the optimization constraints. The input variables of GA are the coefficients am of the combined cosine window function. The constraints of optimization are shown in (4).

m=1M-1am=1m=0M-1-am=00<am<1m=0,1,,M-1 (4)

The peak sidelobe level and the asymptotic sidelobe attenuation rate are taken as the parameter optimization targets of the window function. The objective optimization function can be constructed as:

Fobj=min[-θ(|ασ|+|βλ|)] (5)

where σ and λ are the peak sidelobe level and asymptotic sidelobe attenuation rate, respectively; α and β are the weighting coefficients for the peak sidelobe level and asymptotic sidelobe attenuation rate, respectively, with the constraint α+β=1; and θ is penalty factor.

The values of σ and λ are determined as follows. The time-domain function of the combined cosine window is derived by calculating its coefficients am. The spectrum function is obtained by applying a Fourier transform to the time-domain function, resulting in the corresponding amplitude frequency response curve. The maximum value of each sidelobe is determined by analyzing the curve. Finally, the value of σ is obtained by comparing the sidelobe peak closest to the main lobe, and λ is calculated by determining the slope between each adjacent sidelobe peak points.

The values of θ, α, and β are determined as follows. θ is generally set to be 1, while it is set to be 0 if both λ<0 and σ>0. For α=0.6 and β=0.4, the testing of the algorithm has shown better optimization results.

The parameters for the GA are set as follows. The population size is 100; the number of iterations is 1000; the stall generation is 300; and the function tolerance is 1×10-15.

The width of the main lobe is determined by the number of terms in the combined cosine window function, and the performance of the sidelobes is affected by the coefficients. To achieve optimal results, it is essential to strike a balance and select an appropriate number of terms, neither excessive nor insufficient. In this paper, we achieve this balance by using six terms, and the time-domain expression is calculated as:

ω(n)=0.2625-0.4266cos(2πn/M)+0.2250cos(4πn/M)-0.0726cos(6πn/M)+0.0125cos(8πn/M)-0.0008cos(10πn/M) (6)

C. Self-convolution Window Function and its Spectrum Characteristic

The window function can be self-convoluted to improve the sidelobe performance, and this improvement is further amplified with increasing convolution orders. However, the higher the order, the longer the width of the window function. Considering the sidelobe performance of the window function, the second-order self-convolution window function is adopted in this paper. The new window function can be written as:

ω2(n')=ω(n)*ω(n) (7)

where n'=2n-1 is the length of the convolution window. To ensure that the new window has a length of n'=2n, it is necessary to apply zero to the window.

The spectrum function of the OCCSCW can be derived as:

W2=m=0M-1(-1)mam2WRω-2πmN+WRω+2πmN2 (8)

where ω=2πk/N,k=0,1,...,N-1.

Equation (8) can be simplified as:

W2=e-jπkm=0M-1am2sinkπ2-mπsinkπ2-2mπN+sinkπ2+mπsinkπ2+2mπN2 (9)

Equation (9) can be simplified as:

W2=e-jπksin2kπ2·m=0M-1am2sin2kπNsinkπ2-2mπNsinkπ2+2mπN2 (10)

The cosine self-convolution window function has improved the sidelobe performance with a smaller peak sidelobe peak and a higher sidelobe attenuation rate. This sidelobe performance can better suppress the spectral leakage that occurs in FFT algorithms, thus improving the parameter estimation accuracy of SSO. To demonstrate the advantages of the cosine self-convolution window function compared to traditional window functions, the normalized spectrum of traditional window functions, OCCSCW function, and OCC window function is depicted in Fig. 1. From Fig. 1, it is evident that the OCCSCW function has the smallest peak sidelobe and the highest sidelobe attenuation rate among all the window functions, indicating the excellent characteristics of the OCCSCW function. These characteristics significantly suppress spectral leakage, resulting in a more precise parameter estimation of SSO.

Fig. 1  Frequency response of different window functions.

Note that the window function parameters depend on the peak sidelobe and the sidelobe attenuation rate of spectrum. The parameter optimization process is an offline preparatory work and it only needs to be optimized once. The optimized window function can ensure excellent parameter estimation performance within the sub-/super-synchronous frequency range. The tests are provided to demonstrate the excellent parameter estimation performance within the sub-synchronous frequency range in Section V-A.

III. ApFFT Spectrum Analysis Technique

In this section, an ApFFT spectrum analysis technique is introduced. The superiority of ApFFT in suppressing leakage and enhancing phase estimation can be demonstrated by comparing the spectral functions of FFT and ApFFT.

The signal model of power grid can be modeled as:

x(t)=A0cos(2πf0t+φ0) (11)

where A0, f0, and φ0 are the amplitude, frequency, and phase of the signal, respectively.

Its discrete expression can be written as:

x(n)=A0cos(2πnf0/fs+φ0) (12)

where fs is the sampling rate.

The complex exponential series of (12) can be expressed as:

x(n)=ej(ω0n+φ0)=ej2πk0nN+φ0 (13)

where ω0=2πk0/N; and k0 is the theoretical frequency index, which may not necessarily be an integer.

The traditional N-point FFT spectrum x(k) of x(n) is:

X(k)=1Nn=0N-1ejφ0ej2πk0n/Ne-j2πkn/N=1Nejφ01-ej2π(k0-k)1-ej2π(k0-k)/N=1Nsin(π(k0-k))sin(π(k0-k)/N)ejφ0+N-1Nπ(k0-k) (14)

where k=0,1,,N-1.

ApFFT is a recent signal processing technique. It is well-suited for analyzing close frequency distributions and effectively suppressing spectral leakage. The all-phase data sampling sequence can be written as a data vector.

x=[x(-N+1)        x(0)        x(N-1)] (15)

All N-length shifted data vectors that contain x(0) are constructed as:

x0=[x(0)    x(1)        x(N-1)]x1=[x(-1)    x(0)        x(N-2)]xN-1=[x(-N+1)    x(-N+2)        x(0)] (16)

With the central sampling point x(0) as the reference point, the N vectors obtained from (16) are cyclically shifted.

x0'=[x(0)    x(1)        x(N-1)]x1'=[x(0)    x(1)        x(-1)]xN-1'=[x(0)    x(-N+1)        x(-1)] (17)

From (17), a new data vector with a length of N can be created by calculating the average of x0',x1',,xN-1'.

xap'=1N[Nx(0)    (N-1)x(1)+x(-N+1)        x(N-1)+(N-1)x(-1)] (18)

Based on the shift property of FFT, the spectrum of (16) and (17) has the following relationship.

Xi'(k)=Xi(k)ej2πNiki,k=0,1,,N-1 (19)

As a result, the spectrum of xap' can be derived as:

Y(k)=1Ni=0N-1Xi'(k)=1Ni=0N-1Xi(k)ej2πNik=1N2i=0N-1n=0N-1x(n-i)e-j2πNknej2πNik=1N2i=0N-1n=0N-1ejφ0ej2π(n-i)k0/Ne-j2πNknej2πNik=ejφ0N2i=0N-1e-j2π(k0-k)iNn=0N-1ej2π(k0-k)nN=ejφ0N21-e-j2π(k0-k)1-e-j2π(k0-k)/N1-ej2π(k0-k)1-ej2π(k0-k)/N=1N2sin2(π(k0-k))sin2(π(k0-k)/N)ejφ0 (20)

By comparing (14) and (20), |Y(k)| is equal to the square of |X(k)|. The square relationship results in a faster asymptotic attenuation rate of the sidelobes in the ApFFT spectrum compared to the FFT spectrum. This means that the main spectral line of ApFFT is more prominent, resulting in less spectral leakage and thereby offering superior spectrum analysis capabilities compared with traditional FFT with respect to spectral leakage.

According to (20), each phase value of the traditional FFT is closely related to the frequency deviation value k0-k. The phase value obtained through the traditional FFT is closely associated with the frequency deviation value k0-k. The phase value obtained through the ApFFT method is determined by the value of φ0 presented in (11). This value corresponds to the theoretical phase value of the central sampling point x(0) and remains unchanged by variations in the frequency deviation value k0-k. This characteristic is commonly referred to as the phase-invariance property of ApFFT.

IV. Proposed IpApFFT Based on OCCSCW Function

This section proposes the IpApFFT based on the OCCSCW function to estimate the frequency, amplitude, and phase of SSO.

The operation of applying the OCCSCW function to (12) is performed, and then DFT is used to obtain the spectrum function.

X(k)=A02ejφ0W22π(k-k0)N+e-jφ0W22π(k+k0)N (21)

In synchronous sampling, the frequency f0 falls on the k0th frequency bin of the spectrum in synchronous sampling. However, in asynchronous sampling, the frequency f0 cannot be precisely aligned with a specific bin in the spectrum. As a result, the DFT during asynchronous sampling can cause the fence effect and spectral leakage, which may lead to an inaccurate value. In other words, the actual spectral line k0 is located between km-1 and km or between km and km+1, as shown in Fig. 2, where km represents the bin with the highest magnitude. To minimize errors arising from asynchronous sampling, a double-point spectral interpolation method is employed.

Fig. 2  Discrete spectrum of signal under noncoherent sampling.

In this paper, the two largest bins are defined as km and km+1, respectively, satisfying km<k0<km+1. The amplitudes of the two bins are y1=Y(km) and y2=Y(km+1), respectively. The expressions of y1 and y2 are obtained as:

y1=A12W2(2π(km-k0)/N) (22)
y2=A22W2(2π(km-k0+1)/N) (23)

where A1 and A2 are the amplitudes corresponding to km and km+1, respectively.

A variable parameter a is introduced to represent the deviation between k0 and km. Set a=k0-km-0.5, so a[-0.5, 0.5], and the ratio β is defined as:

β=y2-y1y2+y1 (24)

By substituting (22) and (23) into (24), the ratio β is calculated as:

β=W2(2π(-a+0.5)/N)2-W2(2π(-a-0.5)/N)2W2(2π(-a+0.5)/N)2+W2(2π(-a-0.5)/N)2 (25)

Thus, the relationship between a and β in (25) can be expressed as:

β=g(a) (26)

The specific method for determining the relationship between a and β is to substitute a range of values for a into (25) to obtain a corresponding set of β. Then, a curve fitting function (26) can be obtained by fitting a polynomial, allowing the inverse function of (26) to be calculated as:

α=g-1(β)=0.2181β7+0.3206β5+0.6250β3+2.4269β (27)

Finally, the variable parameter a can be calculated by substituting β, which can be obtained from (24). This enables the derivation of correction formulas for frequency, amplitude, and phase. The frequency correction formula can be derived as:

fest=k0fsN=(a+k1+0.5)fsN (28)

The amplitude correction formula can be derived from (22) and (23) as:

A=A1W2(2π(k1-k0)/N)2+A2W2(2π(k2-k0)/N)2W2(2π(k1-k0)/N)2+W2(2π(k2-k0)/N)2=2(y1+y2)W2(2π(k1-k0)/N)2+W2(2π(k1-k0)/N)2 (29)

The phase value of the maximum spectral bin k1 can be used as the initial phase of the signal due to the phase-invariance property of ApFFT.

The implementation steps of the proposed method can be summarized into three parts.

1) The optimization of window functions is achieved using GAs. It is worth noting that once the window function is obtained, this process is no longer needed, making it an offline preparation work.

2) The sampled signal is windowed with the optimized window function and then processed with ApFFT to obtain the spectral function.

3) The positions of the two most prominent bins in the spectral function are determined. Subsequently, the frequency, amplitude, and phase of SSO are corrected using the double-point spectral interpolation method.

Note that the window function parameters depend on the peak sidelobe and the sidelobe attenuation rate of its spectrum. The window function parameters are optimized by considering the objective function with sidelobe peak and sidelobe attenuation rate. Additionally, the optimized window function parameters do not change with variations in the system parameters. Variations in the system parameters like total connected loads, generators, and lines may lead to changes in oscillation frequencies, but the oscillation frequencies are still within the sub-/super-synchronous frequency range. The optimized window function can ensure excellent parameter estimation performance within this range. The tests are provided to demonstrate the excellent parameter estimation performance within the sub-synchronous frequency range in Section V-A.

V. Verification Studies

To demonstrate the effectiveness and accuracy of the proposed method, verification studies are carried out using synthetic and simulated data under various conditions.

A. Ideal Signal Test

The ideal voltage signal that consists of the fundamental and two sets of complementary sub-/super-synchronous components is represented as:

x(t)=A0cos(2πf0t+φ0)+Asubcos(2πfsubt+φsub)+           Asupcos(2πfsupt+φsup) (30)

where Asub, fsub, φsub and Asup, fsup, φsup are the amplitudes, frequencies, and initial phases of the sub- and super-synchronous components, respectively.

The synthetic signal is modeled as (30), and its parameters are specified as: A0 is set to be 100 V, and fsup is set to be 75.5 Hz. φ0, φsub, and φsup are set to be 90°, 60°, and 30°, respectively. Both nominal and off-nominal frequencies within the range [49.

5, 50]Hz (i.e., f0=49.5, 49.7, 49.9, 50.0, 50.1, 50.3, and 50.5 Hz) are considered. The range of fsub is [10, 45]Hz with a step of 1 Hz. The ranges of Asub and Asup are [10, 100]V with a step of 5 V.

The sampling rate is 3200 Hz, and the sampling data point N is 3200. To assess the accuracy of the results, the estimation errors EE are determined as the relative errors, which are defined as:

EE=p^-pp×100% (31)

where p and p^ are the true and estimated parameter values, respectively.

1) The impact of fsub: with the fixed Asub=10 V and Asup=14 V, the estimation errors of frequencies and amplitudes with various fsub in the range of [

10,45]Hz are shown in Fig. 3.

Fig. 3  Estimation errors of frequency and amplitude with various fsub in range of [

10, 45]Hz. (a) |f0-f^0|/f0. (b) |A0-A^0|/A0. (c) |fsub-f^sub|/fsub. (d) |Asub-A^sub|/Asub. (e) |fsup-f^sup|/fsup. (f) |Asup-A^sup|/Asup.

The relative errors for the frequencies of the fundamental, sub-synchronous, and super-synchronous components are less than 1×10-6, 5×10-4, and 1×10-15, respectively. The relative errors for the amplitudes of fundamental, sub-synchronous, and super-synchronous components are less than 2×10-5, 9×10-4, and 2×10-15, respectively. These results demonstrate that the proposed method accurately identifies the frequencies and amplitudes of each component under different values of fsub.

2) The impact of Asub: with the fixed Asup=14 V and fsub=24.5 Hz, the estimation errors of frequencies and amplitudes for various values of Asub are depicted in Fig. 4. The relative errors for the frequencies of fundamental, sub-synchronous, and super-synchronous components are less than 6×10-7, 1×10-15, and 1×10-16, respectively. The relative errors for the amplitudes of fundamental, sub-synchronous, and super-synchronous components are less than 2×10-5, 4×10-15, and 2×10-15, respectively. These results demonstrate that the proposed method accurately identifies the frequencies and amplitudes of each component under different values of Asub.

Fig. 4  Estimation errors of frequency and amplitude with various Asub in range of [

10, 100]V. (a) |f0-f^0|/f0. (b) |A0-A^0|/A0. (c) |fsub-f^sub|/fsub. (d) |Asub-A^sub|/Asub. (e) |fsup-f^sup|/fsup. (f) |Asup-A^sup|/Asup.

3) The impact of Asup: with the fixed Asub=10 V and fsub=24.5 Hz, the estimation errors of frequencies and amplitudes with various Asup in the range of [

10, 100]V are shown in Fig. 5. The relative errors for the frequencies of the fundamental, sub-synchronous, and super-synchronous components are less than 5×10-7, 2×10-6, and 1×10-16, respectively. The relative errors for the amplitudes of fundamental, sub-synchronous, and super-synchronous components are less than 2×10-5, 4×10-15, and 2×10-15, respectively. These results demonstrate that the proposed method accurately identifies the frequencies and amplitudes of each component under different values of Asup.

Fig. 5  Estimation errors of frequency and amplitude with various Asup in range of [

10, 100]V. (a) |f0-f^0|/f0. (b) |A0-A^0|/A0. (c) |fsub-f^sub|/fsub. (d) |Asub-A^sub|/Asub. (e) |fsup-f^sup|/fsup. (f) |Asup-A^sup|/Asup.

The parameter estimation results under ideal conditions demonstrate that the proposed method achieves high accuracy in estimating parameters of SSO. This indicates that the optimized window function based on the IpApFFT exhibits excellent detection performance within the sub-/super-synchronous frequency range.

To further compare and analyze the estimation accuracy of the proposed method, several sub-/super-synchronous components with close frequency distributions are set up. The Hanning-based IpDFT, a widely-used technique in measurement devices [

27], is used for comparison with the proposed method. This comparison is conducted under conditions where the signal contains non-noise, noise, and non-nominal fundamental frequency variations. The voltage signal model, as shown in (33), contains a fundamental component and four sub-/super-synchronous components, with frequencies of 11.5 Hz, 41.5 Hz, 58.5 Hz, and 88.5 Hz, corresponding to oscillation modes 1-4, respectively.

x(t)=A0cos(2πf0t+φ0)+        Asub1cos(2πfsub1t+φsub1) +Asub2cos(2πfsub2t+φsub2)+        Asup1cos(2πfsup1t+φsup1) +Asup2cos(2πfsup2t+φsup2) (33)

where A0=100 V, f0=49.7 Hz, φ0=15°; Asub1=10 V, fsub1=11.5 Hz, φsub1=60°; Asub2=11 V, fsub2=41.5 Hz, φsub2=30°; Asup1=14 V, fsup1=58.5 Hz, φsup1=45°; and Asup2=17 V, fsup2=88.5 Hz, φsup2=90°.

The comparison results of the proposed method and IpDFT in terms of estimation accuracy are presented in Table I under ideal conditions. The results show that the proposed method has higher accuracy than IpDFT. The relative errors of the frequencies estimated by the proposed method are at least 7 orders lower than those of IpDFT, and the relative errors of the amplitudes estimated by the proposed method are at least 2 orders lower than those of IpDFT. Additionally, the proposed method can achieve zero error in estimating phase, while the maximum error of IpDFT is approximately 4°. These results indicate that the IpDFT has a higher error than the proposed method in estimating the parameters of sub-/super-synchronous components under close frequency conditions. Besides, the proposed method demonstrates better performance in suppressing leakage and significantly improving the estimation accuracy of sub-/super-synchronous components.

TABLE I  Comparison Results of Proposed Method and IpDFT Under Ideal Conditions
ModeParameterProposed methodIpDFT

Estimated

value

Relative error (%)

Estimated

value

Relative error(%)
1 fsub1 11.5000 Hz 2.03×10-12 11.4987 Hz 1.56×10-2
Asub1 9.9998 V 1.58×10-3 10.0139 V 1.40×10-1
φsub1 60.0000° 0 60.4062° 6.80×10-1
2 fsub2 41.5000 Hz 5.64×10-12 41.4771 Hz 5.53×10-2
Asub2 10.9998 V 1.58×10-3 11.2462 V 2.24×100
φsub2 30.0000° 0 34.0914° 13.64×100
3 fsup1 58.5000 Hz 2.64×10-12 58.4975 Hz 8.58×10-3
Asup1 13.9998 V 1.58×10-3 14.1059 V 7.60×10-1
φsup1 45.0000° 0 45.9803° 2.19×100
4 fsup2 88.5000 Hz 4.00×10-12 88.4999 Hz 4.89×10-5
Asup2 16.9997 V 1.58×10-3 17.0120 V 7.00×10-2
φsup2 90.0000° 0 90.0665° 7.38×10-2

The reason for the improved detection accuracy of proposed method compared with traditional ones can be attributed to the following two factors. Firstly, the optimized window exhibits superior sidelobe characteristics compared with traditional ones, which not only greatly improves the ability to suppress spectral leakage, but also reduces the mutual interference with close SSO modes on the parameter estimation of SSO. Secondly, the ApFFT has better spectral leakage suppression capabilities than traditional FFT, and it has the phase-invariance property, resulting in zero error in phase estimation. By combining these two advantages, the proposed method achieves better parameter estimation performance than traditional IpDFT.

B. Time Latency Analysis

Due to the extremely short sampling time of approximately 0.0003 ms, the algorithm of the device does not need to consider the effect of sampling latency [

28], [29]. The latency of the algorithm is mainly determined by the computation time. To determine the computation time, 100 sets of tests are carried out in MATLAB platform (MATLAB R2021(b), PC with 16 GB RAM and 1.6 GHz Intel Core CPU). The results show that the computation time of the proposed method is approximately 7 ms. If this method is applied to the protective relaying system, the latency of this method is still very short compared with the communication latency of 20-50 ms [29].

C. Noise Interference Test

To investigate the performance of the proposed method under noisy conditions, (33) is redone by adding zero-mean Gaussian noise to the test signals [

30]. The corresponding results are given in Tables II-IV under noisy conditions. Noise is the main source of error when the signal-to-noise ratio (SNR) is 20 dB. The results described in Tables II-IV indicate that the proposed method performs better than the traditional IpDFT in the presence of noise. When SNR is 20 dB, the relative errors of frequencies estimated by the proposed method are at least 2 orders lower than those of IpDFT. Similarly, the relative errors of the amplitudes are at least 1 order lower. Additionally, the proposed method for phase estimation has a maximum error of only 4°, which is significantly lower compared with 26.3039° of IpDFT. When the SNR is 40 dB, the noise and mutual interference from other SSO modes are the main sources of errors. The proposed method has relative errors in the frequencies and amplitudes that are at least 2 and 1 orders lower, respectively, compared with those of IpDFT. Furthermore, the proposed method can achieve zero error in phase estimation, whereas the maximum error of IpDFT is 9.2765°.

TABLE II  Frequency Estimation Results of Proposed Method and IpDFT Under Noisy Conditions

SNR

(dB)

Methodfsub1fsub2fsup1fsup2
Estimated value (Hz)Relative error (%)Estimated value (Hz)Relative error (%)Estimated value (Hz)Relative error (%)Estimated value (Hz)Relative error (%)
40 Proposed method 11.5019 1.65×10-2 41.4981 4.58×10-3 58.5011 1.88×10-3 88.5011 1.24×10-3
IpDFT [27] 11.4964 3.13×10-2 41.1827 7.60×10-1 58.2822 3.70×10-1 88.4974 2.94×10-3
20 Proposed method 11.5072 6.26×10-2 41.4764 5.69×10-2 58.5034 5.81×10-3 88.4820 2.03×10-2
IpDFT [27] 11.4175 7.20×10-1 40.9733 1.27×100 58.0832 7.10×10-1 88.3026 2.20×10-1
TABLE III  Amplitude Estimation Results of Proposed Method and IpDFT Under Noisy Conditions

SNR

(dB)

MethodAsub1Asub2Asup1Asup2
Estimated value (V)Relative error (%)Estimated value (V)Relative error (%)Estimated value (V)Relative error (%)Estimated value (V)Relative error (%)
40 Proposed method 10.0024 0.02 11.0329 0.30 13.9375 0.45 16.9969 0.02
IpDFT [27] 10.0129 0.13 11.9020 8.20 14.7447 5.32 17.0312 0.18
20 Proposed method 9.7362 2.64 11.3013 2.74 14.6093 4.35 16.6353 2.15
IpDFT [27] 10.4633 4.63 13.7220 24.75 16.1630 15.45 17.5841 3.44
TABLE IV  Phase Estimation Results of Proposed Method And IpDFT Under Noisy Conditions

SNR

(dB)

Methodφsub1φsub2φsup1φsup2
Estimated value (°)Relative error (%)Estimated value (°)Relative error (%)Estimated value (°)Relative error (%)Estimated value (°)Relative error (%)
40 Proposed method 60.0000 0 30.0000 0 45.0000 0 90.0000 0
IpDFT [27] 60.7883 1.31 39.2765 30.92 47.5567 5.68 91.0379 1.15
20 Proposed method 59.0000 1.67 29.0000 3.33 46.0000 2.22 88.0000 2.22
IpDFT [27] 65.3214 8.87 56.3039 87.68 57.5667 27.93 87.1782 3.14

The analysis conducted under different SNR conditions demonstrates excellent noise resistance in estimating SSO. Even under very low SNR, the proposed method consistently exhibits significantly higher phase estimation accuracy compared with the IpDFT.

D. Dynamic Condition Test

Dynamic features in SSO will affect parameter estimation. To assess the performance of the proposed method, several test cases are provided that consider factors such as PM, AM, FR, and HI. A set of SSO modes from (33) is chosen, specifically frequencies of 41.5 Hz and 58.5 Hz, to simulate dynamic conditions. Other parameters remain consistent with those in (33). The modulation signal model is represented as (18). It comprises the fundamental and two sets of complementary sub-/super-synchronous components, which include AM, PM, FR, and HI.

x(t)=Asub(1+kxcos(2πfmodt))cos(2π(fsub+Rt)t+kacos(2πfmodt))+Asup(1+kxcos(2πfmodt))cos(2π(fsup+Rt)t+kacos(2πfmodt))+A0cos(2πf0t)+Ahcos(2πfht) (34)

where kx is the AM coefficient; ka is the PM coefficient; fmod is the frequency modulation; R is the ramp rate; Ah is the amplitude of the harmonic; and fh is the frequency of the harmonic.

The specific test conditions are as follows.

1) When an SSO occurs, the amplitudes of SSO may exhibit modulations. For AM, we set kx=0.15, ka=0, and fmod=0.15 Hz.

2) PM may also occur in an SSO. For PM, we set kx=0, ka=0.15, and fmod=0.15 Hz.

3) The frequencies of a voltage/current signal can exhibit linear changes in an SSO. The sub-synchronous frequency changes at a ramp rate of R=0.02 Hz/s, and the super-synchronous frequency changes at a ramp rate of R=-0.02 Hz/s.

4) Harmonics can also be present in the voltage signal. In this case, the voltage signal is superimposed with a 5% third harmonic and a 3% fifth harmonic.

The total vector error (TVE) and frequency error (FE) as defined in the IEEE Synchrophasor Standard [

31] are used to evaluate estimation results. The results in Table V demonstrate that the proposed method improves the estimation accuracy in dynamic conditions compared with IpDFT, with lower errors in all tested cases. The maximum TVE and FE of the proposed method are 0.0374% and 0.9 mHz, while those of IpDFT are 3.21% and 25.4 mHz. In particular, the proposed method performs better under the condition of FR, with lower errors compared with IpDFT. In summary, the proposed method demonstrates higher precision in estimation.

TABLE V  Results Under Different Test Conditions
Frequency (Hz)MethodAMPMFRHI
TVE (%)FE (mHz)TVE (%)FE (mHz)TVE (%)FE (mHz)TVE (%)FE (mHz)
41.5 Proposed method 0.0009 0.9 0.09 0 10.60 64.1 0.0018 0
IpDFT [27] 3.2100 17.1 2.56 25.4 66.91 641.8 1.9300 4.8
58.5 Proposed method 0.0374 0.9 0.09 0 10.22 1.2 0.0018 0
IpDFT [27] 0.1600 10.8 0.86 7.2 61.86 629.1 0.1500 0.1

E. Application to Simulated Data

The proposed method is further validated using simulated data of SSO. As shown in Fig. 6 [

32], a grid-connected directed-drive wind farm (DDWF) via a voltage source converter-based high-voltage direct current (VSC-HVDC) transmission system is modeled in MATLAB/Simulink platform. The rating of permanent magnet synchronous generator (PMSG) is 50 MW. The length of the AC transmission line is 10 km, and the unit equivalent resistance, inductance, and capacitance parameters of the line are 0.01273 Ω/km, 0.934 mH/km, and 0.01274 μF/km, respectively. The parameters of the PMSG and VSC-HVDC in Fig. 6 are shown in [32].

Fig. 6  System diagram of DDWFs via VSC-HVDC transmission systems.

Figure 7(a) shows the waveform of the phase-A voltage obtained through simulation. The simulated data between 16 s and 20 s are adopted for analysis, with a sampling frequency of 1200 Hz. Figure 7(b) presents the corresponding spectrum. The spectrum in Fig. 7(b) reveals the presence of three sub-/super-synchronous components in the voltage signal. A sub-synchronous frequency (fsub29 Hz) and a super-synchronous frequency (fsup71 Hz) are selected as a set of frequency components in the close frequency distributions.

Fig. 7  Simulation results of phase-A voltage. (a) Voltage waveform. (b) Spectrum.

Since the actual measurement results of the signal are unknown, the Hanning-based two-point IpDFT [

26] is used to compare its performance with the proposed method. The length of the data window for both methods is set to be 1 s. The frequency and amplitude of the sub-/super-synchronous components are continuously measured in Figs. 8 and 9, respectively. The proposed method exhibits a smoother frequency and amplitude response curve, and it yields voltage signal estimates that are closer to the actual values than those obtained using IpDFT. The primary source of significant fluctuations in IpDFT estimates is interference from other mode components. To further reinforce the superiority of the proposed method, the standard deviations of estimates are given in Table VI. The suggested approach has standard deviations for frequency and amplitude that are 9 and 32 times smaller, respectively, than IpDFT. In conclusion, the proposed method has higher precision and stability.

Fig. 8  Frequency response curves of proposed method and IpDFT. (a) Sub-synchronous frequency. (b) Super-synchronous frequency.

Fig. 9  Amplitude response curves of proposed method and IpDFT. (a) Sub-synchronous amplitude. (b) Super-synchronous amplitude.

TABLE VI  Standard Deviation Statistics
ComponentMethodStandard deviation
Frequency (Hz)Amplitude (V)
Sub-synchronous Proposed method 0.0059 0.27
IpDFT [27] 0.1900 2.46
Super-synchronous Proposed method 0.0060 0.28
IpDFT [27] 0.2000 2.51

VI. Conclusion

In this paper, a novel method, i.e., IpApFFT with an OCCSCW, is proposed for estimating parameters of SSO. The proposed method has better performance than the traditional one in suppressing spectral leakage and the fence effect. Furthermore, it provides accurate phase estimation of SSO due to its phase-invariance property. Simulation results demonstrate that the proposed method provides high precision and noise resistance in estimating amplitude and frequency of SSO. Under ideal conditions, the proposed method offers 107 times higher precision in estimating frequency and 200 times higher precision in estimating amplitude compared to IpDFT. Under dynamic conditions, the TVEs and FEs of the proposed method are 86 times and 26 times higher, respectively, compared with IpDFT.

References

1

J. Shair, X. Xie, L. Wang et al., “Overview of emerging subsynchronous oscillations in practical wind power systems,” Renewable and Sustainable Energy Reviews, vol. 99, pp. 159-168, Jan. 2019. [Baidu Scholar] 

2

X. Wu, M. Wang, M. Shahidehpour et al., “Model-free adaptive control of STATCOM for SSO mitigation in DFIG-based wind farm,” IEEE Transactions on Power Systems, vol. 36, no. 6, pp. 5282-5293, Nov. 2021. [Baidu Scholar] 

3

S. Kovacevic, D. Jovcic, P. Rault et al., “Three-way subsynchronous torsional interactions between LCC-HVDC, MMC-HVDC and a thermal generator,” Journal of Modern Power Systems and Clean Energy, vol. 11, no. 4, pp. 1331-1340, Jul. 2023. [Baidu Scholar] 

4

T. Rauhala, A. M. Gole, and P. Järventausta, “Detection of subsynchronous torsional oscillation frequencies using phasor measurement,” IEEE Transactions on Power Delivery, vol. 31, no. 1, pp. 11-19, Feb. 2016. [Baidu Scholar] 

5

Y. Li, L. Fan, and Z. Miao, “Replicating real-world wind farm SSR events,” IEEE Transactions on Power Delivery, vol. 35, no. 1, pp. 339-348, Feb. 2020. [Baidu Scholar] 

6

X. Xie, X. Zhang, H. Liu et al., “Characteristic analysis of subsynchronous resonance in practical wind farms connected to series-compensated transmissions,” IEEE Transactions on Energy Conversion, vol. 32, no. 3, pp. 1117-1126, Sept. 2017. [Baidu Scholar] 

7

S. Xu, H. Liu, and T. Bi, “Field PMU test and calibration method– part I: general framework and algorithms for PMU calibrator,” Journal of Modern Power Systems and Clean Energy, vol. 10, no. 6, pp. 1507-1518, Nov. 2022. [Baidu Scholar] 

8

N. Ma, X. Xie, P. Kang et al., “Wide-area monitoring and analysis of subsynchronous oscillation in power systems with high-penetration of wind power,” Proceedings of the CSEE, vol. 41, no. 1, pp. 65-74, Jan. 2021. [Baidu Scholar] 

9

P. Gopakumar, B. Mallikajuna, M. Jaya et al., “Remote monitoring system for real time detection and classification of transmission line faults in a power grid using PMU measurements,” Protection and Control of Modern Power Systems, vol. 3, no. 1, pp. 159-168, Jun. 2018. [Baidu Scholar] 

10

A. Carta, N. Locci, and C. Muscas, “A PMU for the measurement of synchronized harmonic phasors in three-phase distribution networks,” IEEE Transaction on Instrumentation and Measurement, vol. 58, no. 10, pp. 3723-3730, Oct. 2009. [Baidu Scholar] 

11

M. Chakir, I. Kamwa, and H. L. Huy, “Extended C37.118.1 PMU algorithms for joint tracking of fundamental and harmonic phasors in stressed power systems and microgrids,” IEEE Transaction on Power Delivery, vol. 29, no. 3, pp. 1465-1480, Jun. 2014. [Baidu Scholar] 

12

H. C. Lin, “Power harmonics and interharmonics measurement using recursive group-harmonic power minimizing algorithm,” IEEE Transactions on Industrial Electronics, vol. 59, no. 2, pp. 1184-1193, Feb. 2012. [Baidu Scholar] 

13

D. Macii, D. Petri, and A. Zorat, “Accuracy analysis and enhancement of DFT-based synchrophasor estimators in off-nominal conditions,” IEEE Transactions on Instrumentation and Measurement, vol. 61, no. 10, pp. 2653-2664, Oct. 2012. [Baidu Scholar] 

14

D. Agrez, “Weighted multipoint interpolated DFT to improve amplitude estimation of multifrequency signal,” IEEE Transactions on Instrumentation Measurement, vol. 51, no. 2, pp. 287-292, Apr. 2002. [Baidu Scholar] 

15

L. Chen, W. Zhao, F. Wang et al., “An interharmonic phasor and frequency estimator for subsynchronous oscillation identification and monitoring,” IEEE Transactions on Instrumentation Measurement, vol. 68, no. 6, pp. 1714-1723, Jun. 2019. [Baidu Scholar] 

16

H. Wen, C. Li, and L. Tang, “Novel three-point interpolation DFT method for frequency measurement of sine-wave,” IEEE Transaction on Industrial Informatics, vol. 13, no. 5, pp. 2333-2338, Oct. 2017. [Baidu Scholar] 

17

T. Jin and W. Zhang, “A novel interpolated DFT synchrophasor estimation algorithm with an optimized combined cosine self-convolution window,” IEEE Transaction on Instrumentation Measurement, vol. 70, pp. 1-10, Jan. 2021. [Baidu Scholar] 

18

J. Borkowski, D. Kania, and J. Mroczka, “Interpolated-DFT-based fast and accurate frequency estimation for the control of power,” IEEE Transaction on Industrial Electronics, vol. 61, no. 12, pp. 7026-7034, Dec. 2014. [Baidu Scholar] 

19

H. Wen, C. Li, and L. Tang, “Novel three-point interpolation DFT method for frequency measurement of sine-wave,” IEEE Transaction on Industrial Informatics, vol. 13, no. 5, pp. 2333-2338, Oct. 2017. [Baidu Scholar] 

20

H. Wen, C. Li, and W. Yao, “Power system frequency estimation of sine-wave corrupted with noise by windowed three-point interpolated DFT,” IEEE Transactions on Smart Grid, vol. 9, no. 5, pp. 5163-5172, Sept. 2018. [Baidu Scholar] 

21

D. Belega and D. Petri, “Accuracy analysis of the multicycle synchrophasor estimator provided by the interpolated DFT algorithm,” IEEE Transaction on Instrumentation Measurement, vol. 62, no. 5, pp. 942-953, May 2013. [Baidu Scholar] 

22

X. Xie, Y. Zhan, H. Liu et al., “Improved synchrophasor measurement to capture sub/super-synchronous dynamics in power systems with renewable generation,” IET Renewable Power Generation, vol. 13, no. 1, pp. 49-56, Jan. 2019. [Baidu Scholar] 

23

T. Su, M. Yang, T. Jin et al., “Power harmonic and interharmonic detection method in renewable power based on Nuttall double-window all-phase FFT algorithm,” IET Renewable Power Generation, vol. 12, no. 8, pp. 953-961, May 2018. [Baidu Scholar] 

24

S. Liu, N. Lyu, J. Cui et al., “Improved blind timing skew estimation based on spectrum sparsity and ApFFT in time-interleaved ADCs,” IEEE Transactions on Instrumentation Measurement, vol. 68, no. 1, pp. 73-86, Jan. 2019. [Baidu Scholar] 

25

J. Zhao, Y. Zhou, J. Zhao et al., “Precision position measurement of PMSLM based on ApFFT and temporal sinusoidal fringe pattern phase retrieval,” IEEE Transactions on Industrial Informatics, vol. 16, no. 12, pp. 7591-7601, Dec. 2020. [Baidu Scholar] 

26

D. Belega, D. Dallet, and D. Petri, “Statistical description of the sine-wave frequency estimator provided by the interpolated DFT method,” Measurement, vol. 45, no. 1, pp. 109-117, Jan. 2012. [Baidu Scholar] 

27

P. Romano and M. Paolone, “Enhanced interpolated-DFT for synchro phasor estimation in FPGAs: theory, implementation, and validation of a PMU prototype,” IEEE Transactions on Instrumentation Measurement, vol. 63, no. 12, pp. 2824-2836, Dec. 2014. [Baidu Scholar] 

28

J. Li, S. Xu, H. Liu et al., “High-accuracy and low-complexity phasor estimation method for PMU calibration,” CSEE Journal of Power and Energy Systems, vol. 7, no. 6, pp. 1202-1212, Nov. 2021. [Baidu Scholar] 

29

X. Xie, Y. Zhan, H. Liu et al., “Wide-area monitoring and early-warning of subsynchronous oscillation in power systems with high-penetration of renewables,” International Journal of Electrical Power & Energy Systems, vol. 108, pp. 31-39, Jun. 2019. [Baidu Scholar] 

30

M. Brown, M. Biswal, S. Brahma et al., “Characterizing and quantifying noise in PMU data,” in Proceeding of 2016 IEEE PES General Meeting (PESGM), Boston, USA, Jul. 2016, pp. 1-5. [Baidu Scholar] 

31

IEEE Standard for Synchrophasor Measurements for Power Systems-Amendment 1: Modification of Selected Performance Requirements, IEEE Standard C37.118.1a-2014, 2014. [Baidu Scholar] 

32

X. Guo, Y. Li, X. Xie et al., “Sub-synchronous oscillation characteristics caused by PMSG-based wind plant farm integrated via flexible HVDC system,” Proceedings of the CSEE, vol. 40, no. 4, pp. 1149-1160, Feb. 2020. [Baidu Scholar]