An adaptive method for oscillations monitoring in power systems with high penetration of renewable energy

This paper proposes an adaptive monitoring method for wideband, multimodal, and time-varying oscillations in power systems with high penetration of renewable energy. The method is based on the sliding windowed interpolation FFT (SWIFFT) technique and adaptive threshold spectral filtering method to accurately identify the amplitude, frequency, phase, and damping coefficients of the oscillation modes in each segment. Additionally, the time-domain relationships of these parameters are established, enabling real-time monitoring of the oscillations. Finally, the proposed method is tested using a 5-mode example signal and a simulation system measurement signal. The results demonstrate that this method can accurately identify the parameters of complex oscillation signals, enabling real-time monitoring of oscillations. It exhibits flexibility in adapting to various oscillation scenarios and meets the signal characteristic requirements in modern power systems.


Introduction
With the development of modern power systems with high penetration of renewable energy [1], the interactions and mutual influences between power electronic devices and the power grid may give rise to a risk of complex oscillations [2] [3], which consist of low-frequency (0.2-2.5 Hz), sub/supersynchronous, and medium-to-high-frequency oscillatory components.These phenomena are increasingly exhibiting strong nonlinearity, multimodality, and time-varying characteristics [4] [5] [6].
Accurately identifying each oscillatory mode and capturing their temporal variations in the time domain would contribute to the development of targeted suppression measures based on different frequency bands.However, existing methods face challenges in achieving precise real-time monitoring of complex signals that exhibit wide frequency spectra, multimodality, and time-varying characteristics.Consequently, their direct applicability in practical environments of advanced highpenetration renewable energy power systems is limited [7] [8].
This paper introduces a novel adaptive analysis algorithm called sliding window interpolation FFT (SWIFFT) by incorporating a sliding window approach into the windowed interpolation FFT algorithm and integrating an adaptive filtering method for signal spectrum processing.Section II and III presented the structure of the oscillatory signal and the main methods used in this paper.The combination of the principles presented before was proposed in Section III to give a complete set of steps for monitoring methods.In Sections IV and V, a example signal and a measurement signal of the simulated system model were used to evaluate the proposed method.Upon examination, it was determined that the method presented in this paper exhibits feasibility, applicability, and flexibility.Section VI provides a comprehensive summary of the paper.

Windowed signal model
Let the oscillatory signal model with M frequency components be represented after discrete sampling at a sampling frequency fs as follows: The amplitude, frequency, and phase of each component are represented as Am, fm, and θm.By setting ωm=2πfm/fs, the following expression can be obtained through the Fourier transform: The cosine window function is represented as a weighted sum of R cosine terms, and expressed as: The weighting coefficients in this expression satisfy two constraints:


In the selection of window functions, existing high-order window function expressions are often complex, requiring the solution of high-order equations.The classical Blackman-Harris window, a cosine window function, is widely used due to its ability to minimize sidelobe levels, achieve ideal main lobe width, and exhibit excellent sidelobe suppression.It finds applications in scenarios with high dynamic range, such as acoustic sensors and radar signal processing, where accurate localization of peak frequency points is desired.
The amplitude function of Blackman-Harris window can be expressed as follows.
According to the convolution theorem, , the spectrum function of the discrete signal after applying the cosine window can be obtained.
Considering the symmetry of the spectrum function with respect to positive and negative frequencies, we can neglect the sidelobe effects in the negative frequency range.Let the frequency spacing of the spectrum be denoted as Δf=fs/N, the positive frequency point fm corresponds to the positive integer sampling point km, m m f k f    , then:

Three-spectral line interpolation correction method
To suppress the aliasing effects caused by non-integer sample rates (km being a non-integer), a spectrum interpolation algorithm is used for correction.The choice of the number of spectrum lines determines the efficiency of the algorithm.To strike a balance between recognition accuracy, algorithm complexity, and computational time, selecting the three spectral lines closest to the peak frequency point for interpolation correction can yield more optimal results.
Let km1 represent the line number which is closest to the peak frequency point of the spectrum, with km0 being the adjacent line on the left and km2 being the adjacent line on the right.Introducing the correction coefficient , =km-km1, and -0.5<<0.5.When disregarding the influence of other modes on the m-th component, l0=|YW(km0)|, l1=|YW(km1)|, and l2=|YW(km2)| are magnitudes of the three spectral lines.
Then set another parameter  , , which can be expressed as: From this expression, it is evident that β is a function of α, α can be considered as the inverse function of β.In other words, ,and α can be expressed as 1 ( ) g     .After being fitted, α can be expressed as a polynomial weighted sum, given by: Where, 0 = 0.938919, 1 = 0.082038, 2 = 0.015411, 3 = 0.0031699.
Weighted averaging is performed on the amplitudes of the three spectral lines.The spectral line with the central index, which exhibits the highest amplitude, is assigned a higher weight.The expression for calculating the weight is as follows: ) The expression for calculating the amplitudes Am of each mode is expressed as follows: By utilizing the polynomial fitting method as well as the aforementioned relationship, we can obtain that The correction formula for the frequencies fm and phase m of each mode can be expressed as follows: 1 ( ) For signals with time-varying amplitudes of each mode, there exists a functional relationship between Am and time t.This relationship describes how the amplitudes of the different modes change over time.By analyzing this relationship, we can gain insights into the temporal dynamics of the signal and accurately capture the variations in amplitudes for each mode.
Each damping coefficient can be determined by fitting the curve of the amplitudes in the time domain.After fitting, the expression of the amplitude can be represented as Am(t), ( ) . Where the damping coefficient of each mode is σm.Taking the logarithm of the equation above: So the expression of σm can be obtained.

Sliding windowed interpolation FFT
Considering the increasing time-varying characteristics of power systems, an improved approach called SWIFFT has been developed by using sliding windows into the WIFFT methodology to capture localized information effectively.The main principle is to use a sliding window to segment the signal along the time axis.Each window has a relatively short length, assuming that the signal's waveform characteristics remain constant within each window.This allows us to treat the truncated signal as stationary.Consequently, the original signal can be seen as a collection of short-time stationary signals.By applying interpolated FFT, we obtain local spectra.As the window slides along the time axis, the identified parameter information is accumulated, resulting in the temporal variation of the parameters.
Firstly, the signal extracted from the first time window is analyzed by applying a windowed interpolation WIFFT.Subsequently, shift the time window.All the parameters are correspondingly matched with the time coordinates, using the moment window's endpoint as the time coordinate.Finally, the variation of each parameter with time is obtained, and the damping coefficient is determined based on the amplitude Am(t).
For a discrete signal y(n) with a sampling frequency of fs and a total of N samples, the set window length is L seconds in the time domain, and each window will contain N0 = Lfs samples.A step size of d samples, resulting in a total of (N-N0)/d window shifts.The number of windows is (N-N0)/d +1.Furthermore, the width of the window function determines the time resolution and frequency resolution.Therefore, it is necessary to determine the width of the window function based on the characteristics of the signal, striking a balance between time resolution and frequency resolution.

Adaptive threshold spectral filtering method
In power systems with high penetration of renewable energy, the spectrum of signals may contain numerous false spectral peaks with very small amplitudes due to the presence of noise and the inherent asynchronous sampling of FFT.Based on the difference in quantity between the two types of peaks, a threshold-based method can be used to distinguish the actual components in the signal.If the threshold is set to a fixed value, it lacks flexibility in distinguishing.When the SNR is low, a low threshold may classify some noise as actual components or consider actual components as noise.The threshold should be adaptively adjusted based on the noise level.Here introduces a spectral filtering method with adaptive thresholding, which follows the specific procedure outlined below:  Calculate the mean value μ and standard deviation value σ of the amplitude values at each frequency point across the entire spectrum. Calculate the preliminary threshold value Pth , Pth = μ + K × σ, K is a adjustment coefficient. Compare each amplitude value in the spectrum with the threshold Pth.If there are frequency components with larger amplitudes than the threshold Pth, they are considered actual components.Remove these actual components from the spectrum and return to the previous step.If no frequency components remain after the previous step, it indicates that all the remaining components are false. By subtracting the spectrum containing only false components from the original signal spectrum, we obtain the spectrum containing only the dominant components.

Specific steps of the method
 The continuous-time signal is globally discretely sampled, and the sliding window ought to be set with a specific width and step size. The position of the sliding window is initialized, and the signal segment obtained from the beginning is windowed using the Blackman-Harris window. The Fourier transform is performed on the windowed signal, resulting in the complete frequency spectrum of the signal. Adaptive threshold filtering is applied to remove noise and false spectral peaks caused by asynchronous sampling from the frequency spectrum, thereby identifying the dominant modes. Three-spectral line interpolation is applied at the peak frequency points to correct the amplitude, frequency, and initial phase parameters of each mode in the signal.
 Move the sliding window, and the same processing steps are applied to the signal segment captured by the window, continuing until the entire signal is traversed in the time domain. Consequently, the relationships between the parameters and time can be obtained, and the damping coefficient can be determined based on the relationship between amplitude and time.

Test of example signals
To validate the feasibility and applicability of the proposed method, we first construct a signal as follows, assuming a total of 5 modes: The white noise component is ξ, which is introduced into the signal model to better simulate the real environment of power systems with high penetration of renewable energy.The white noise is generated based on a signal-to-noise ratio (SNR) of 50.The signal parameters are set as follows： Assuming a signal length of 10s, fs = 1024 Hz.The sliding window has a width of 1s in the time domain, containing 1024 discrete points.The sliding window moves with a fixed step size of d=5 points.Since the analysis is conducted on the signal captured by the first window from 0s to 1s, the parameter analysis results of this window represent the characteristics at 1s.The relationship between each parameter and time will start from 1s.For a stationary signal, the parameter identification results of the initial window are presented in Table 2 and the relative errors are illustrated in Figure 1.The method based on the Blackman-Harris window and three-spectrum line interpolation is found to be reliable, with modal parameter identification exhibiting very small relative errors in the range of 10 -4 to 10 -8 .These results indicate that the proposed method in this study achieves high precision.
After applying the sliding window, the amplitude and frequency curves analysis shows that the amplitude decreases based on the preset coefficient, while the frequency remains constant.The phase of the sinusoidal signal only considers the initial phase, this study focuses on the changes in amplitude and frequency over time.

Test of simulated system model signal
A simulation model was constructed in Figure 3 in the DIgSILENT/PowerFactory software.The model represents a single-machine infinite bus power system operating at a frequency of 60Hz.
Figure 3. Simulation model of a small system.This system consists of a series of compensating lines between buses A and B, and one generator.During abnormal operating conditions, the capacitance of the series compensating line and the synchronous generator's shaft system may form an oscillatory loop.In this scenario, one or multiple oscillatory modes typically exist within the loop.Analysis of the voltage signal at bus A reveals the relationship between amplitude and frequency over time, as shown in Figure 4.   shows the time-varying frequency.This system exhibits two dominant modes, they are the fundamental mode at 60Hz and the sub-synchronous oscillation mode at 39.7Hz.Obviously, in the presence of sub-synchronous oscillation mode with negative damping characteristics, small disturbances trigger unstable sub-synchronous oscillations.

Conclusion
This paper proposes an oscillation monitoring method for high-proportion renewable energy power systems.The method utilizes a windowed interpolated FFT algorithm and improved sliding window technique to estimate parameter values and capture their time-varying patterns, facilitating effective oscillation detection.Extensive testing confirms the theoretical feasibility of the proposed algorithm.In the future, our method will be enhanced to accommodate a higher number of modes and improve precision, while also focusing on improving noise resistance capabilities.These improvements aim to better adapt to the requirements of new energy power systems and ensure accurate and reliable analysis for their effective monitoring and integration into the grid.

Figure 1 .Figure 2 .
Figure 1.Relative errors of each Parameter.The estimated parameters vary with time by sliding the window as presented in Figure 2.

Figure 4 .
Figure 4.Each amplitude and frequency of the model signal vary with time.Figure 4 (a) shows the time-varying amplitude of each mode obtained through sliding window after recognizing the simulated model signals.Figure 4 (b)shows the time-varying frequency.This system exhibits two dominant modes, they are the fundamental mode at 60Hz and the sub-synchronous oscillation mode at 39.7Hz.Obviously, in the presence of sub-synchronous oscillation mode with negative damping characteristics, small disturbances trigger unstable sub-synchronous oscillations.

Figure 4 (
Figure 4.Each amplitude and frequency of the model signal vary with time.Figure 4 (a) shows the time-varying amplitude of each mode obtained through sliding window after recognizing the simulated model signals.Figure 4 (b)shows the time-varying frequency.This system exhibits two dominant modes, they are the fundamental mode at 60Hz and the sub-synchronous oscillation mode at 39.7Hz.Obviously, in the presence of sub-synchronous oscillation mode with negative damping characteristics, small disturbances trigger unstable sub-synchronous oscillations.

Table 1 .
Parameters of each mode of the signal The oscillatory signal is used to simulate the sustained and stable oscillation stage.

Table 2 .
Parameter identification results