ESPRIT associated with filter bank for power-line harmonics, sub-harmonics and inter-harmonics parameters estimation

Electrical Power and Energy Systems 118 (2020) 105731

Contents lists available at ScienceDirect

Electrical Power and Energy Systems journal homepage: www.elsevier.com/locate/ijepes

ESPRIT associated with filter bank for power-line harmonics, sub-harmonics and inter-harmonics parameters estimation

T

⁎

Elaine Santosa, Mahdi Khosravya,b, , Marcelo A.A. Limaa, Augusto S. Cerqueiraa, Carlos A. Duquea a b

Electrical Engineering Department, Federal University of Juiz de Fora, Brazil Electrical Engineering Department, University of the Ryukyus, Okinawa, Japan

A R T I C LE I N FO

A B S T R A C T

Keywords: Power quality Harmonics estimation Inter-harmonics estimation Sub-harmonic estimation ESPRIT Spectrum spread Filter bank

Detection of power-line frequency components and estimating their parameters are crucial tasks for maintenance of applied power quality, especially due to uncertainties originated from growing utilization of renewable power sources, hybrid energy storages and increasing demand for electrical vehicles. This manuscript introduces the association of filter bank to ESPRIT (FB-ESPRIT) for higher efficiency in detection and estimation of the harmonics, inter-harmonics, and sub-harmonics parameters in power quality measurement. FB-ESPRIT initially applies an analysis filter bank of parameter L to the signal, then individually implements ESPRIT to each of L resultant signal sub-band while it is spread over the full bandwidth of the signal. Due to individual ESPRIT implementation on L times spectrum spread of subbands, FB-ESPRIT is expected to have L times higher accuracy and robustness in parameter estimation compared to ESPRIT. Its efficiency analysis has been carried out over three case studies composed of multiple harmonics, inter-harmonics, sub-harmonic and fundamental components; a synthesized signal with under 500 Hz bandwidth, a synthesized signal with a bandwidth up to 64th harmonic in 50, 000 Monte Carlo simulations (MCS) at different noise levels, and a real power signal measured at the output of a photovoltaic solar power plant in a partly cloudy windy summer morning. The comparative evaluations approve the multiple times accuracy and robustness dominance of efficiency of ESPRIT with versus without filter bank association. It has been theoretically illustrated that a FB-ESPRIT calibrated to give the same efficiency as ESPRIT has the advantage of L2 times less complexity order, as well approved via the mean result of extensive simulations runs.

1. Introduction The estimation of frequency harmonic components in the Electric Power System (EPS) is a crucial task not only for the evaluation of the power quality [1–3] but also to assist engineers in power system planning and designing like measuring the power quality disturbances [4], wind energy quality assessment [5,6], improving the losses in hybrid energy storage [7], electric-vehicle charger connection to the power grid [8], etc. In one side, the recent growing use of distributed, renewable and intermittent generations with their high uncertainty, and on the other side the proliferation of power electronic equipment and nonlinear loads are responsible for a significant amount of harmonic distortion in the power systems, and thus justifying the increasing interest of the researchers on harmonic analysis [9,10]. In the point of signal processing view, the harmonic analysis consists of two main challenges: (i) the detection of behavior of the power system

⁎

signal such as frequency deviations, variable harmonic content, presence of disturbances, etc. as faster as possible and (ii) providing an efficient estimation of harmonics content through a short length of signal with few samples for some application such as protection. Several algorithms have been developed over the years to estimate harmonics in the power system signal as Refs. [11,12] present a sequential evolution of them. IEC 61000-4-7 standard [13] recommends the Fast Fourier Transform (FFT) as it is widely in use, however, it has resolution limitation in its application to harmonic analysis, especially over short length signal measurements. Pisarenko was the first subspace eigendecomposition parametric method of spectrum estimation which could overcome FFT in giving higher resolution [14]. The drawback of Pisarenko was in its need for exact a priori information about the model order, also it estimates the autocorrelation lags in a statistical way which results in model incompatibility and inaccuracy. Because of that, as presented by a comparative analysis [15], the Prony technique was

Corresponding author. E-mail addresses: [email protected] (E. Santos), [email protected] (M. Khosravy), [email protected] (M.A.A. Lima), [email protected] (A.S. Cerqueira), [email protected] (C.A. Duque). https://doi.org/10.1016/j.ijepes.2019.105731 Received 30 July 2019; Accepted 21 November 2019 0142-0615/ © 2019 Elsevier Ltd. All rights reserved.

Electrical Power and Energy Systems 118 (2020) 105731

E. Santos, et al.

accuracy and robustness comparative evaluation via the mean and standard deviation of estimations in 50, 000 Monte Carlo runs, (vi) discussion over the down-sampling effect of aliasing on the signal and the resultant interference in sub-bands, (vii) illustration of the required filter bank characteristics to minimize the sub-bands interference, (viii) statistical evaluation over four different categories of sub-harmonics, fundamental frequency, middle frequencies, and edge frequencies close to sampling theorem limitation, (ix) performance evaluation on the real data of a photo-voltaic power plant output with multiple sub-harmonics due to the time of measurement in a partly cloudy windy morning, and (x) noise sensitivity comparative evaluation by statistical analysis of 50,000 Monte Carlo simulation at different levels of SNR values of 5 dB, 15 dB, 25 dB, and 35 dB. The organization of the manuscript is as follows. Section 2, after a brief explanation to ESPRIT, illustrates the structure of the proposed FBESPRIT. Section 3 illustrates the complexity advantage of the proposed technique and theoretically obtains the optimal number of filters in the complexity aspect. Section 4 analyses the efficiency of FB-ESPRIT in terms of accuracy, robustness, and computational load in comparison to the conventional ESPRIT over synthetic signal and real power data. Section 5 presents the future scope, and finally, Section 6 concludes the manuscript.

found more preferable which was also due to its lower computational load. Also, it was less suffering from spurious components because of using squared error of a recursive process, and a simple checking approach for model order estimation. However, Prony was prone to noise and computationally inefficient [16]. The subspace technique was generalized in MUSIC algorithm [17], which has high resolution over a short length of data without any side lobes effect [18]. However, despite the attractive advantages of MUSIC, it has two main drawbacks of high computation load and the need for a large storage [19]. At this stage, the well-known parametric method of harmonic estimation based on rotational invariance technique (ESPRIT) was introduced and drew the researchers attention [19–23]. ESPRIT possesses several main advantages [16] as (i) estimation efficiency, (ii) high accuracy , (iii) no need of memory storage, (iv) high resolution and capability of detecting inter-harmonics, etc. Spectral estimation by ESPRIT has been used in a variety of application in power quality as Ref. [24] employs ESPRIT for monitoring the frequency of generator as a feature for islanding detection, and Ref. [25] applies a modification of ESPRIT on synchrophasor measurements for identification of low-frequency mode of power systems. ESPRIT has also been applied via a virtual instrument [26] in power quality indices defined by IEEE std 14592010 [27]. The required di0mension of autocorrelation matrix in ESPRIT depends on the model order, which is according to the signal attributes [28]. While Ref. [29] uses ESPRIT for power frequency and parameter estimation, its model order is accurately estimated based on relative difference of consecutive eigenvalues, and the adaptive application of the same measure for model order estimation through an iterative formulation has been used to modify ESPRIT for time-efficient harmonic estimation over non-stationary signals [30]. Since ESPRIT has the ability of precisely estimating the signal parameters by a short length of measurements, it has been deployed for effective estimation of dynamic phasor of alternating current signal [31]. Also, ESPRIT has been implemented for harmonics parameters estimation when the signal is contaminated by multiplicative noise [32], efficient assessment of powerline waveform distortions [33], and online monitoring and measurement of and synchronized phasor for reliable power quality identification [34,35], harmonic source identification in distribution system [36]. This manuscript suggests the implementation of ESPRIT with the association of filter bank for power systems harmonics, inter-harmonics, and sub-harmonics detection and parameters estimation in a structure which results in a multi-times improvement in performance accuracy, robustness, and speed. In the suggested association of filter bank with ESPRIT (FB-ESPRIT), a uniform filter bank provides sub-band signals of equal bandwidth, then in a parallel structure, each sub-band after going through a down-sampling spread spectrum is fed to an individual ESPRIT parameter estimation. FB-ESPRIT is evaluated by detection and estimation of parameters of fundamental, harmonics, interharmonics, and sub-harmonic components of (i) a sinusoidal of 500 Hz low bandwidth signal with multiple sub-harmonics, harmonics and inter-harmonics, (ii) a synthesized multi-component sinusoidal signal of 3840 Hz bandwidth under different levels of noise, and (iii) a real data of a photovoltaic solar power plant with multiple sub-harmonic interferences. The accuracy and robustness of FB-ESPRIT are numerically evaluated through mean bias and mean standard deviation of estimations in 50, 000 Monte Carlo simulation (MCS) runs for each SNR level at the range from 5 to 35 dB. The manuscript extends the initial idea presented in Ref. [37] in the following aspects: (i) theoretical justification of higher efficiency of ESPRIT with the association of filter bank, (ii) theoretical illustration and calculation of the extent of improvement in time complexity as given in Lemma 1, (iii) a theoretical discussion on the optimum number of filters L constructing the associated uniform filter bank as presented in Lemma 2, (iv) statistical analysis in the comparative evaluation of ESPRIT performance with versus without filter bank association through 50, 000 Monte Carlo runs in different noisy conditions, (v)

2. Filter Bank associated ESPRIT Filter banked ESPRIT (FB-ESPRIT) is the algorithm of parallel individual implementations of ESPRIT to the spread spectrum of each sub-band of the signal. The proposed method associates a filter bank to the conventional ESPRIT to provide equally spaced sub-bands, then it spreads each sub-band by a down-sampling process to cover all the signal estimation bandwidth. Indeed, FB-ESPRIT estimates the parameters of each sub-band of the signal individually in the absence of the other ones. Dedication of all signal bandwidth to each multiple times narrower sub-band results in higher accuracy, and robustness of parameter estimation. After a brief explanation to ESPRIT technique, we go in technical details of FB-ESPRIT. 2.1. Estimation of Spectral Parameters via Rotational Invariance Technique (ESPRIT) ESPRIT is a parametric spectral estimation which explores the deterministic relationship between subspaces via rotational invariance, and it was originally developed by Roy et al. [19–21]. The rotational invariance technique is the heart of ESPRIT technique. The signal frequency components are made of complex exponentials, and according to time-shifting theory, each single complex exponential s0 [n] = e j2πfn fulfills the following property

s0 [n + 1] = s0 [n] e j2πf

(1)

which indicates that the next sample value is a phase-shifted version of the current value, or in other words, its rotation on the unit circle e j2πf . Considering a signal in length-M time-window vectorial form as its M−1 current and future samples as x [n] = [x [n] x [n + 1] …x [n + M − 1]]T to be indicated in a model of complex exponentials s [n] in addition to noise vector w [n] P

x [n] =

∑ p=1

αp v (fp ) e j2πnfp + w [n] = V Φnα + w [n]

(2)

where V is M × P matrix made up of length-M time-window frequency vectors v (fp ) = [1 e j2πfp … e j2π (M − 1) fp]T corresponding to P frequencies as V = [v (f1 ) v (f2 ) … v (fp )], the α = [α1 α2 … αP ]T is the vector consisted of amplitudes of each exponential αp , and the Φ is the diagonal matrix of phase shifts between adjacent time samples of each complex exponential component of s [n] 2

Electrical Power and Energy Systems 118 (2020) 105731

E. Santos, et al.

j2πf 0 ⎡e 1 j2πf2 ⎢ 0 e =⎢ ⋮ ⋮ ⎢ 0 ⎣ 0

0 ⎤ … 0 ⎥ … ⋱ ⋮ ⎥ ⎥ … e j2πfP ⎦

where Us is the matrix made of right-hand singular vectors corresponding to the P biggest singular values in magnitude. Please notice that according to the model of Equ. Eq. (2), all the frequency vectors v (f ) for frequencies f = f1 , f2 , …, fP locate in signal subspace. Thus, the V and Us span the same subspace, and there exists an invertible transformation T mapping Us into V :

(3)

where ϕi = e j2πfi . As it is observable, Φ is a rotational invariant characteristic of the signal known as rotation matrix which fully describes the frequencies of the signal complex exponentials. Using this fact, ESPRIT estimates the frequency components by finding the rotation matrix Φ. The approach is by consideration of two overlapping M − 1 length sub-windows sM − 1 [n] and sM − 1 [n + 1] within the length-M timewindow vector as follows

sM − 1 [n] ⎤ ⎡ s (n) ⎤ s [n] = ⎡ = ⎢ ⎥ ⎢ [n + 1]⎥ s [ 1] s n M + − M − 1 ⎦ ⎣ ⎦ ⎣

sM − 1 [n] = VM − 1 Φnα

As Equ. Eq. (8) partitions the matrix V Φ into two smaller –dimensional (M − 1) –dimensional subspaces, we partition the Us in the same way … 1 = ⎡ Us = ⎡U (16) ⎣…⎤ ⎦ ⎣U2 ⎤ ⎦ Note that U1 and V1 correspond to the unstaggered subspace, and U2 and V2 correspond to staggered subspace. The relation in (15) is hold between corresponding subspaces as follows

(4)

V1 = U1 T

V2 = U2 T

(17)

Similar to subspace rotation relation between V1 and V2 in Eq. (7), there should be a similar, though different subspace rotation relation between U1 and U2 as

(5)

The matrix VM − 1 is made in the same way as V with this difference that its time-window vectors vM − 1 (f ) are of M − 1 length, VM − 1 = [vM − 1 (f1 ) vM − 1 (f2 ) … vM − 1 (fP )]. From Eq. (5), the following definitions for V1 and V2 , the subspaces corresponding to unstaggered and staggered windows and their relations are obtained;

V2 = VM − 1 Φn + 1

(15)

V = Us T.

sM − 1 [n] is a length-(M − 1) subwindow of s (n) and it can be expressed as follows

V1 = VM − 1 Φn

(14)

U = [Us |Uw]

Φ=diag {ϕ1 ϕ2 …ϕP }

U2 = U1 Ψ.

(18)

SVD on data matrix X gives the subspaces U1 and U2 . Then, using a least square technique, we obtain Ψ from Eq. (18)

(6)

Ψ = (U1H U1)−1U1H U2

V2 = V1 Φ

(7)

… 1 = ⎡ ⎤ V Φ = ⎡V ⎣…⎤ ⎦ ⎣V2 ⎦

(8)

Recall that ESPRIT frequency estimation is via searching solution for subspace rotation matrix Φ. The Eqs. (7), (17) and (18) together lead to estimation of Φ . Solving for V2 in two ways (i) substituting (18) in (17) and (ii) substituting (17) in (7) gives V2 = U1 ΨT and V2 = U1 TΦ, where equality of right-hand sides results in

In the other side, there is a relation between V and the autocorrelation matrix of x [n] in Equ. Eq. (2) as illustrated in the following equations.

Rx = E {x [n] x H [n]} = Rs + Rw

ΨT = T Φ

(9)

∑

|αp |2 v (fp ) v H (fp ) + σw2 I = VAV H + σw2 I

p=1

VAV H ,

(10)

fp ̂ =

σw2 I

and Rw = are respectively the autocorrelation Rx , Rs = matrices of x [n], s [n] and the white noise w [n]. The matrix A is a diagonal matrix of the powers of each of the respective complex exponentials, A = diag{|α1 |2 |α2 |2 … |αP |2 } . In practice, the correlation matrix is not known and it is estimated from the available data samples x [n] at time n. Practically, Rx is estimated by a data matrix X

1 Rx̂ = X H X N

XN × M

2π

(21)

2.2. FB-ESPRIT In FB-ESPRIT, first, a filter bank decomposes the signal to its frequency sub-bands with equal bandwidth. Then each sub-band in the absence of the other signal sub-bands goes individually through the ESPRIT algorithm. Since the other sub-bands are not present, ESPRIT can dedicate whole the bandwidth deployed for parameter estimation to the sub-band under process, thereof a multiple times efficiency in parameter estimation is expected. The sub-band signal is frequency spread over the whole bandwidth of the signal. This aim is achieved by down sampling of the signal by factor L that is the same as the number of acquired sub-bands. The down-sampled sub-band covers the empty spectrum of the other sub-bands in the frequency band; thus it is expected that the accuracy and robustness of ESPRIT in the detection of the corresponding sub-band parameters L times increase. The FBESPRIT algorithm is divided into the following steps:

(11)

(12)

Having the data matrix X , and applying the singular value decomposition (SVD)

X = LΣU H

∡ϕp

where ϕp for p = 1, 2, …, P are the eigenvalues of Ψ , and ∡ϕp is the phase of ϕp . Fig. 1 depicts an illustrative block diagram of theory and algorithm of ESPRIT.

where the data matrix X is made up of stacking the rows of N measured length-M time window data vectors x [n] samples as T x [1] … x [M − 1] ⎤ ⎡ X [0] ⎤ ⎡ x [0] ⎢ X T [1] ⎥ ⎢ x [1] ⎥ x [2] … x [M ] =⎢ ⎥=⎢ ⎥ ⋮ ⋮ ⋱ ⋮ ⋮ ⎥ ⎢ ⎢ ⎥ T [N − 1]⎥ − … + − [ 1] [ ] [ 2] x N x N x N M ⎢ X ⎦ ⎦ ⎣ ⎣

(20)

where indicates the relationship between eigenvectors and eigenvalues of the matrix Ψ . Thus, the diagonal elements of Φ are the eigenvalues of Ψ , thus the estimate of frequencies are as follows

P

=

(19)

(13)

where L and U are N × N and M × M unitary matrices of left and right singular vectors, respectively. From Equ. Eq. (11), the squared magnitudes of the singular values equal to the eigenvalues of Rx̂ scaled by a factor of N, and the corresponding eigenvectors are the columns of U . Therefore, U shapes an orthonormal basis for the underlying M dimensional vector space. This subspace can be subdivided into signal and noise subspaces

Step 1: L bandpass filters of the filter bank separately filter the signal s (t ) , and give the L sub-banded signals s0 (t ), s1 (t ), …, sL (t ) . Step 2: Each of the i-th sub-banded signals si are down-sampled by parameter L equivalent to the number of band-pass filters used

3

Electrical Power and Energy Systems 118 (2020) 105731

E. Santos, et al.

Fig. 1. Demonstrating block diagram of ESPRIT theory and algorithm. The upper part with dashed boxes and ovals is with just theoretical relations without calculation, but the bottom part with solid boxes and ovals is with practical calculation on measured data.

in the filter bank, and resulting in siL ↓. Step 3: ESPRIT individually estimates the parameters of each downsampled sub-band signal siL ↓. Step 4: The estimated parameters are relocated to the corresponding sub-band. The union of the parameters obtained for all subbands makes the set of the signal parameters.

Our preference for filter bank implementation is by using Infinite Impulse Response (IIR) filters due to their considerable lower order than FIR filter with equivalent magnitude response. Especially, filters with sharp transition width are required, and as the transition width of an FIR filter decreases, its order inversely increases. Another important issue about the filter bank is decomposing the signal s (t ) to the subband si (t ) signals of equal bandwidth.

Fig. 2 illustrates the FB-ESPRIT block-diagram wherein L parallel ESPRITs estimate the parameters of L spread spectrum sub-bands.

2.2.2. Down-sampling considerations In a down-sampler operation, only the samples of x (n) where n is a multiple of L are used, and the rest of samples are decimated;

2.2.1. Filter bank considerations The filters used in the filter bank affect the efficiency of the proposed FB-ESPRIT. To clarify the role of each bandpass filter characteristics in FB-ESPRIT, notice that each separated signal sub-band by a filter is composed of the filter passband output besides the stop-band output.The stop-band of the filter attenuates the other sub-bands effect but does not eliminate it completely, and due to that as down-sampling by the parameter L spreads the spectrum of each sub-band to cover the F full range of sampling limit [0 − 2s ] Hz, the attenuated spectrum of other sub-bands alias the corresponding sub-band spectrum by superposing its frequency content. Therefore, to maximally avoid the interference of the other sub-bands and the adjacent sub-band contents, the filters should have (i) an efficient stop-band attenuation, and (ii) a sharp transition width, respectively.

x L ↓ (n) = x (Ln).

(22)

In the proposed FB-ESPRIT method, the filter bank output signals are decimated by a factor of L, since the filter bank is composed of L filters where the output of each is a sub-band occupying 1 of the signal L spectrum. Due to down-sampling by the factor of L, the spectrum of each sub-band signals is spread L times. When it is the first sub-band B0 = [0Hz − fcHz ], the bandwidth of the sub-band as spread by the factor of L, it covers from zero to L × fc Hertz and covers full bandwidth of the signal. The other sub-bands cover full bandwidth of the signal after down-sampling too, but in a different way. It is due to the aliasing phenomenon. As an instance, while already B0 is cut out, by downsampling by the factor of L, the corresponding sub-band spectrum

Fig. 2. The filter bank associated ESPRIT block diagram. 4

Electrical Power and Energy Systems 118 (2020) 105731

E. Santos, et al.

B1 = [fc − 2fc ] is spread over full spectrum of the signal but in flipped way wherein fc locates on Lfc, and 2fc locates on 0 Hz, and the middle frequencies in a mirror pattern locate from L × fc to 0Hz . Considering the effect mentioned above of down-sampling in mapping the location of the frequencies in the spread spectrum, after their estimation, the effect of spread spectrum is removed from them by the following relocation of the estimated frequencies: L↓

fi ̂ =

⎧i × fc + 1 f i ̂ , L

=

O (FB-ESPRIT) ≈

3.1. On the choice of L Although the Lemma 1 indicates that a FB-ESPRIT with signal spectrum division to more sub-bands has lower complexity, there are issues limiting the factor L. A practical limitation for increment of L is the filters efficiency. Even a high order filter does not have a stop-band with hundred percent attenuation. After L times spreading the spectrum of a sub-band, the attenuated content of the other L − 1 sub-bands superpose the sub-band due to aliasing effect. Therefore in noisy condition, the power of the noise will be added to L − 1 times attenuated noise of the other sub-bands. In the case of using a big number of filters L, such a noise increment can affect SNR condition and thereof the FBESPRIT efficiency. Another limitation for L comes from the complexity of the filters. Besides the effect of L on reducing complexity order via the first term in Eq. (27), it increases the computation load of the second term in the equation proportionally to L. Having filters with an order comparable to the signal snapshot length, the optimum L in the context of computational load is given by the following lemma.

3. Computational complexity Most of the parametric methods have high computation burden, mainly due to involving either the eigenvalue decomposition or matrix inversion, each having the computation load requirement of O (m2n) and O (n3) , respectively for n × m and n × n matrices. For example, the parametric method of Prony requires singular value decomposition of a Toeplitz matrix, finding the roots of polynomial and then matrix inversion by least-squares estimation of Vandermonde matrix [38]. Amongst the parametric methods, MUSIC [18] and ESPRIT are known for their high resolution and comparative computational efficiency, and as illustrated in Ref. [20], ESPRIT has 105 times complexity advantage compared to MUSIC. It is very important for a spectrum estimation technique to have low computational complexity, as possessing a faster performance makes the technique more suitable for tracking the spectrum deviations in a non-stationary condition where the harmonics and inter-harmonics parameters changes by time. Although the focus of the manuscript is on the improvement of accuracy and robustness by the association of filter bank, this section theoretically compares the FB-ESPRIT complexity with ESPRIT while they are set to give equivalent performance. The computational load of each ESPRIT implementation equivalent to the number of multiplications needed in a parametric technique is as follows [39].

Lemma 2. Assume the FB-ESPRIT algorithm with a filter bank of L filters with order p comparable to n the signal snapshot length, p ≈ n , then the number of filters in a FB-ESPRIT with optimal complexity is

Lopt = 2⌊log2

(29)

Proof. Please see the Appendix.

□

From Lemma 2, we have obtained the optimum L values for different lengths of signal snapshot n and the autocorrelation vector M, and shown in Table 1. Please notice that all L values given in the Table 1 are powers of two, it is for the sake of using the polyphase technique in filter bank implementation [40] which requires the number of filters to be a power of two.

(24)

Lemma 1. Supposing the distribution of frequency components along the signal spectrum is uniform, then a FB-ESPRIT with performance close to ESPRIT has the advantage of L2 times lower computational order than ESPRIT

4. Results analysis and discussion This section analyzes the performance of ESPRIT associated with filter bank and compares its efficiency with respect to the ESPRIT without filter bank in detection and quantification of the parameters of harmonics, inter-harmonics and sub-harmonics in power line system. We have used the multirate implementation for the filter banks composed of IIR filters of order 24 with 80 dB stopband attenuation and

(25)

Proof. Due to the assumption of uniform distribution of frequency components, since FB-ESPRIT divides the signal spectrum to L subbands which each goes separately through an ESPRIT estimation, the M can be divided by L for each individual implementation of ESPRIT. On the other hand, since the spectrum is spread L times for each ESPRIT block in FB-ESPRIT, subsequently as expected, the frequency detection resolution increases L times. Therefore, acquiring an L times shorter snapshot length n of the signal should results in the same resolution as L ESPRIT without spectrum spread. Thereof, from Eq. (24), to have similar performance as ESPRIT, the computational order of FB-ESPRIT implementation is expected to be as follows

n M n O (FB-ESPRIT) = L ⎛9( )3 + 12 ( )2⎞ + Lnp L L ⎠ ⎝ L

3 6(3n + 4M ) ⌉

where ⌊. ⌉ indicates the closest integer.

where n is the snapshot length, and M is the autocorrelation vector length. Lemma 1 demonstrates the complexity load of a FB-ESPRIT algorithm compared to ESPRIT, while both are set to give equivalent efficiency.

1 O (ESPRIT) L2

(28)

(23)

where f i ̂ is as estimated frequency for i-th spread sub-band, and fi is the corrected frequency by removing the effect of down-sampling.

O (FB-ESPRIT) ≈

1 (9n3 + 12Mn2) L2

□

L↓

O (ESPRIT) = 9n3 + 12Mn2.

(27)

The second term is due to computational load of implementation of L filters of order p. The second term is ignorable compared to the first term of O (n3) , thus

i = 0, 2, 4, ⋯

⎨ (i + 1) × f − 1 f L̂ ↓ , i = 1, 3, 5, ⋯ c L i ⎩

1 (9n3 + 12Mn2) + Lnp L2

Table 1 Optimum L from Lemma 2 for different snapshot length n, and autocorrelation vector length M. n

(26) 5

M

2

4

8

16

32

64

128

256

512

1024

4 16 64 256 1024

4 8 16 16 32

4 8 16 16 32

4 8 16 16 32

8 8 16 16 32

8 8 16 16 32

8 16 16 16 32

16 16 16 16 32

16 16 16 16 32

16 16 16 32 32

32 32 32 32 32

Electrical Power and Energy Systems 118 (2020) 105731

E. Santos, et al.

sharp transition width of 20 Hz.The performance analysis is over three different case studies; (i) synthetic power signal with under 500 Hz content of multiple harmonics, inter-harmonics located up to 2 Hz close to each other, and sub-harmonics with 1 Hz frequency difference, (ii) synthetic power signal with multiple harmonics, inter-harmonics, and sub-harmonics with frequency content up to 3840 Hz with different SNR levels of 5 dB, 15 dB, 25 dB, and 35 dB, and (iii) real data from a photovoltaic power plant in a windy partly cloudy morning. The filter bank factor L is taken as L = 4 for all the case studies mentioned above. Please note the optimum values of L in Table 1 are without consideration of filters stopband attenuation gain limitation and resulting interference of the sub-bands. IIR filters deployed here are close to ideal to have a clear analysis of FB-ESPRIT advantages in dividing and spreading the signal spectrum.

C0 = (384 − 512) Hz. The estimated parameters of s1 (t ) by ESPRIT without and with association of filter bank have been shown in Table 2. The comparison is over the signal lengths of 35, 37, 40, 48, 56, and 64 1.0938, samples which are respectively equivalent to 1.1563, 1.25, 1.5, 1.75, and 2 cycles of the signal. As it is observable from the Table 2, while in all cases of detection and estimation of frequencies, amplitudes and phases, FB-ESPRIT performs with zero error for the signal lengths N ⩾ 37 samples (i.e. N ⩾ 1.1563 cycles) ESPRIT without association of filter bank for the signal length N < 48 samples (i.e. N < 1.5 cycles) is not applicable (NA) in most of the cases, and for longer signal length the estimation is with a considerable error. Despite perfect detection of 30 Hz and 31 Hz sub-harmonics frequencies from each other by FB-ESPRIT as well 345 Hz and 347 Hz when N ⩾ 1.1563 cycles, ESPRIT is not able to distinguish these closely located sub-harmonics and inter-harmonics when the signal length is less than 2 cycles and 1.75 cycles, respectively. The detailed comparison of estimation relative error is given in Fig. 3. The relative estimation error of 31 Hz parameters at the presence of 30 Hz (ΔF = 1 Hz) by two techniques has been compared in top row curves of Fig. 3. As it is seen, FB-ESPRIT performs with zero error from N ⩾ 37 samples (i.e. N ⩾ 1.1563 cycles), ESPRIT from N = 64 samples (i.e. N = 2 cycles) can detect the 31 Hz with an acceptable relative error. Middle row curves of Fig. 3 compares the relative error of estimation 347 Hz parameters at the presence of similar power of 345 Hz component (ΔF = 2 Hz). Similarly, in this case from N ⩾ 37 samples (i.e. N ⩾ 1.1563 cycles) FB-ESPRIT performs with zero error, and ESPRIT start to perform with an acceptable relative error for N ⩾ 56 samples (i.e. N ⩾ 1.75 cycles). The bottom curves compares the two techniques in general case by the average of relative error of the estimated parameters corresponding to the other frequencies in Table 2. In general cases that frequencies are located with more than 2 Hz difference from each other, the ESPRIT performs with an acceptable error for N ⩾ 48 samples (i.e. N ⩾ 1.5 cycles). As mentioned in introduction and as the above comparison results

4.1. Synthetic signal with frequency contents under 500 Hz The first case study signal is as follows:

s1 (t )= cos(2π × 60t + 15°) + 0.4cos(2π × 30t + 25°)+ 0.5cos(2π × 31t + 35°) + 0.4cos(2π × 180t + 45°)+ 0.3cos(2π × 300t + 60°) + 0.25cos(2π × 345t + 75°)+ 0.3cos(2π × 347t + 25°) + 0.2cos(2π × 420t + 60°)+ 0.1cos(2π × 455t + 68°).

(30)

It is made of a fundamental, third, fifth, and seventh harmonics, two sub-harmonics of 30 Hz and 31 Hz with only 1 Hz difference from each other, and three inter-harmonics which two of them are 345 Hz and 347 Hz with only 2 Hz difference from each other. Since the highest frequency component is less than 500 Hz, the sampling frequency is taken as 1024 Hz. The number of sample points per cycle NPPC is taken as 32 points. The filter bank is made by multirate implementation of IIR filters, and the filter bank parameter L = 4 , thereof the sub-bands are C0 = (0 − 128) Hz, C0 = (128 − 256) Hz, C0 = (256 − 384) Hz, and

Table 2 Case study 1: Detection and parameters estimations for harmonics, inter-harmonics, and sub-harmonics with frequencies under 500 Hz. Signal length:

35 = 1.0938 cycles

37 = 1.1563 cycles

40 = 1.25 cycles

48 = 1.5 cycles

56 = 1.75 cycles

64 = 2 cycles

Technique:

ESP.

FBE.

ESP.

FBE.

ESP.

FBE.

ESP.

FBE.

ESP.

FBE.

ESP.

FBE.

Frequencies

Sub. 30 Hz Sub. 31 Hz Fnd. 60 Hz Har. 180 Hz Har. 300 Hz Int. 345 Hz Int. 347 Hz Har. 420 Hz Int. 455 Hz

NA NA NA NA NA NA NA NA NA

30.5541 30.5541 60.0046 180.0000 300.0000 346.0705 346.0705 420.0009 455.0000

NA NA 59.4584 NA NA NA NA 402.0637 442.1448

30.0000 31.0000 60.0000 180.0000 300.0000 345.0000 347.0000 420.0000 455.0000

NA NA 49.6768 191.4954 NA NA 353.2364 425.5114 459.0704

30.0000 31.0000 60.0000 180.0000 300.0000 345.0000 347.0000 420.0000 455.0000

30.5653 30.5653 59.9302 179.9934 299.9530 346.1332 346.1332 419.9826 455.0085

30.0000 31.0000 60.0000 180.0000 300.0000 345.0000 347.0000 420.0000 455.0000

30.5624 30.5624 60.0009 180.0000 300.0000 345.0201 347.0293 420.0000 455.0000

30.0000 31.0000 60.0000 180.0000 300.0000 345.0000 347.0000 420.0000 455.0000

29.9505 30.9210 60.0000 180.0000 300.0000 345.0001 347.0001 420.0000 455.0000

30.0000 31.0000 60.0000 180.0000 300.0000 345.0000 347.0000 420.0000 455.0000

Amplitudes

30 Hz, 0.4 31 Hz, 0.5 60 Hz, 1 180 Hz, 0.4 300 Hz, 0.3 345 Hz,.25 347 Hz, 0.3 420 Hz, 0.2 455 Hz, 0.1

NA NA NA NA NA NA NA NA NA

0.9124 0.9124 0.9981 0.4000 0.2984 0.5423 0.5423 0.2000 0.1000

NA NA NA NA NA NA NA NA NA

0.4000 0.5000 1.0000 0.4000 0.3000 0.2500 0.3000 0.2000 0.1000

NA NA NA NA NA 0.7939 0.7939 NA NA

0.4000 0.5000 1.0000 0.4000 0.3000 0.2500 0.3000 0.2000 0.1000

0.8970 0.8970 1.0076 0.3999 0.2990 0.5001 0.5001 0.1999 0.1000

0.4000 0.5000 1.0000 0.4000 0.3000 0.2500 0.3000 0.2000 0.1000

0.8841 0.8841 0.9999 0.4000 0.3000 0.2542 0.2910 0.2000 0.1000

0.4000 0.5000 1.0000 0.4000 0.3000 0.2500 0.3000 0.2000 0.1000

0.3531 0.5850 1.0000 0.4000 0.3000 0.2500 0.3000 0.2000 0.1000

0.4000 0.5000 1.0000 0.4000 0.3000 0.2500 0.3000 0.2000 0.1000

Phases

30 Hz, 25° 31 Hz, 35°

NA

30.6763

NA

25.0000

NA

25.0000

30.2488

25.0000

30.5231

25.0000

8.3836

25.0000

NA

30.6763

NA

35.0000

NA

35.0000

30.2488

35.0000

30.5231

35.0000

33.7256

35.0000

60 Hz, 15° 180 Hz, 45°

NA

14.8674

NA

15.0000

NA

15.0000

15.1537

15.0000

14.9892

15.0000

14.9998

15.0000

NA

45.0000

NA

45.0000

NA

45.0000

45.0246

45.0000

45.0000

45.0000

45.0000

45.0000

300 Hz, 60° 345 Hz, 75°

NA

59.8941

NA

60.0000

NA

60.0000

60.1699

60.0000

60.0000

60.0000

60.0000

60.0000

NA

48.2235

NA

75.0000

NA

75.0000

47.4459

75.0000

73.1246

75.0000

74.9930

75.0000

347 Hz, 25° 420 Hz, 60°

NA

48.2235

NA

25.0000

25.6714

25.0000

47.4459

25.0000

25.4369

25.0000

25.0000

25.0000

NA

60.0000

NA

60.0000

NA

60.0000

60.0761

60.0000

59.9999

60.0000

60.0000

60.0000

455 Hz, 68°

NA

68.0000

NA

68.0000

NA

68.0000

67.9345

68.0000

67.9999

68.0000

68.0000

68.0000

6

Electrical Power and Energy Systems 118 (2020) 105731

E. Santos, et al.

Fig. 3. Relative error of estimation of Under 500 Hz components parameters. Upper row: 31 Hz at the presence of 30 frequency (left), amplitude (middle), and phase (right). Middle row: 347 Hz at the presence of 345 Hz frequency (left), amplitude (middle), and phase (right). Bottom row: general case: frequency (left), amplitude (middle), and (right) phase.

show, a strong point of the parametric techniques is their ability in detection of frequency components with high resolution over a short length of the signal. As observed in the above case study, FB-ESPRIT and ESPRIT have respectively the ability of frequency detection with resolution of 1 Hz while the signal length is N ⩾ 2 cycles and N ⩾ 1.1563 cycles, respectively for ESPRIT and FB-ESPRIT. But in the case of using discrete Fourier transform (DFT), the number of cycles of 60 Hz fundamental frequency required for having resolution Δf is given by Ncycles = 60/Δf . Therefore, to detect the frequency contents with resolution less than 1 Hz by DFT, the signal should have the length of N > 60 cycles. Similarly, the required signal length for frequency detection with resolutions of 2 Hz, 3 Hz, 4 Hz, and 5 Hz by DFT is respectively N = 30, 20, 15, and 12 cycles. Apart from the considerable longer required signal length, the detectable frequencies by DFT just can be located on the multiples of fs / N . If real frequency is different from the bin frequency it leaks for the adjacent slots and not accurate detection is possible.

for both FB-ESPRIT and ESPRIT is acquired as 32 points. The goal of this section is comparison of accuracy and robustness of ESPRIT with and without association of filter bank under different noise conditions while the acquired signal length is taken equal to 1.25 cycle. The acquired sampling frequency is 7680 Hz, which guarantees the ability of measuring frequencies up to 3840 Hz. The filter bank divides the signal spectrum to four sub-bands with equal bandwidth by four filters as (i) C0 (z ) a low pass filter with cutoff frequency of 960 Hz, (ii) C1 (z ) a band pass filter with low and high cutoff frequencies at 960 Hz and 1920, (iii) another band pass filter with low and high cutoff frequencies at 1920 Hz and 2880 Hz, and (iv) a high pass filter with cutoff frequency at 2880 Hz. Fig. 4 illustrates frequency response magnitude of the filter bank. The equal bandwidth sub-bands at the output of the filter bank B0 = (0 − fc ), B1 = (fc − 2fc ), B2 = (2fc − 3fc ) , are as and B3 = (3fc − 4fc ) where fc = 960 Hz. Tables 3–5 respectively indicate the mean value, standard deviation (std.) and mean bias of estimated parameters of s2 (t ) in the abovementioned 50, 000 MCS runs at each signal to noise ratio (SNR) level of 5 dB, 15 dB, 25 dB and 35 dB with different random noise for each iteration. Figs. 5 and 6 visually compare the accuracy and robustness of FB-ESPRIT and ESPRIT via density distribution of the estimated frequencies of s2 (t ) in 50, 000 MCS runs for the SNRs of 5 dB and 25 dB,

4.2. Wide bandwidth signal at different noisy conditions To have an exact evaluation of association of filter bank with ESPRIT on accuracy and robustness in parameter estimation, we have acquire a synthesized sinusoidal signal with multiple harmonics, interharmonics and sub-harmonic components with wide bandwidth of 3300 Hz under different noise conditions. This case study signal is as follows: The synthesized signal s2 (t ) is as follows

s2 (t )= cos(2π × 60t + 15°) + 0.2cos(2π × 30t + 25°)+ 0.3cos(2π × 27.5 × 60t + 45°) + 0.2cos(2π × 43.5 × 60t + 60°) + 0.1cos(2π × 55 × 60t + 68°)

(31)

Signal s2 (t ) in Eq. (31) comprises a fundamental component, two interharmonics, one sub-harmonic and a harmonic, and the parameters which are targeted to be estimated are the frequencies, amplitudes and the phases of sinusoidal components of s2 (t ) . The evaluative analysis of s2 (t ) parameters estimation is carried out via 50, 000 Monte Carlo simulation (MCS) runs for different level of additive white Gaussian noise while the number of points per cycle Nppc

Fig. 4. Frequency response magnitude of the Filter Bank used for the analysis the proposed FB-ESPRIT. It comprises of (a) C0 (z ) low pass filter with cutoff frequency fc = 960 Hz, (b) C1 (z ) band pass filter with cutoff frequencies fc1 = 960 and fc 2 = 1920 , (c) band pass filter with fc1 = 1920 Hz and fc 2 = 2880 Hz, and (c) high pass filter with fc = 2880 Hz. 7

Electrical Power and Energy Systems 118 (2020) 105731

E. Santos, et al.

Table 3 The mean of the estimated parameters of the signal s2 (t ) (Eq. (31)) in 50, 000 Monte Carlo runs of the proposed FB-ESPRIT and ESPRIT at different SNR levels: 5, 15, 25, and 35 dB, while the snapshot length for both is 1.25 cycles. SNR Estimator Fund. F.: 60 Hz Sub. F.: 30 Hz Int. F.: 1650 Int. F.: 2610 Har. F.: 3300

5 dB

15 dB

25 dB

35 dB

FB-ESP.

ESPRIT

FB-ESP.

ESPRIT

FB-ESP.

ESPRIT

FB-ESP.

ESPRIT

59. 9596 30. 0835 1650. 0022 2609. 9731 3300. 0410

59.4771 NA 1650.0281 2609.9455 3298.7704

59. 9964 30. 0060 1650. 0014 2610. 0006 3300. 0063

60.3138 NA 1650.0761 2610.0462 3299.9222

60. 0018 30. 0005 1650. 0009 2610. 00016 3299. 9997

59.7658 NA 1650.0131 2610.0181 3300.0373

60. 0018 29. 9999 1650. 0006 2610. 00019 3300. 00002

59.7657 NA 1650.0123 2610.012 3300.0343

Fund. 60 A.: 1 Sub. 30 A.: 0. 2 Int. 1650 A.: 0. 3 Int. 2610 A.: 0. 2 Har. 3300 A.: 0. 1

1. 0. 0. 0. 0.

0217 3291 3655 3421 2407

1.2177 NA 0.3728 0.3352 0.2687

0. 0. 0. 0. 0.

9972 2287 3067 2541 1581

0.2160 NA 0.3378 0.3162 0.2359

0. 0. 0. 0. 0.

9989 2075 3006 2010 1023

0.1101 NA 0.3006 0.2076 0.1876

0. 0. 0. 0. 0.

9996 2022 3001 2001 1002

0.1037 NA 0.2999 0.2022 0.1983

Fund. 60 ϕ : 15° Sub. 30 ϕ : 25° Int. 1650 ϕ : 45° Int. 2610 ϕ : 60° Har. 3300 ϕ : 68°

15. 33. 45. 60. 68.

5655° 0339° 0800° 3699° 1831°

11.7016° NA 45.0776° 60.5723° 68.7135°

15. 25. 45. 60. 68.

0526° 1045° 0454° 0037° 2296°

7.9896° NA 44.9145° 60.1144° 68.3361°

15. 25. 45. 60. 68.

0058° 0040° 0132° 0014° 0773°

9.3676° NA 44.9442° 59.9460° 67.8625°

15. 25. 45. 60. 68.

0008° 0012° 0025° 0014° 0081°

9.5429° NA 44.9669° 59.9514° 67.8703°

In our analysis, the amplitude of 30 Hz component of s2 (t ) has been acquired 0.2 which is five times less than the amplitude of the fundamental component. While ESPRIT can not detect this component, FBESPRIT in addition to frequency, estimates the amplitude and the phase of 30 Hz component, as it can be observed with details in Tables 3–5. The accuracy of 30 Hz amplitude estimation is respectively 0.1291, 0.0282, 0.0075, and 0.0022 for the SNRs 5 dB, 15 dB, 25 dB, and 35 dB, and robustness indicated by mean Std of estimation for the corresponding SNRs is respectively 0.1364, 0.0275, 0.0091, and 0.0028. In the case of phase estimation, apart from very low SNR of 5 dB, that the mean estimation bias is 8.0339 degree, for the SNRs 15 dB, 25 dB, and 35 dB, FB-ESPRIT has a high accuracy as it is respectively 0.1045, 0.0045, and 0.0012 . The robustness indicated by mean Std of MCS runs for the corresponding SNRs is respectively 18.8042, 8.9455, 2.8323, and 0.8952 degree.

respectively. The data vector used for the autocorrelation in ESPRIT and each of internal ESPRITs of FB-ESPRIT has the length of 96 and 24 samples, respectively. The analysis of the results of MCS runs in Tables 3–5 and Figs. 5 and 6 is presented in the following subsections. 4.2.1. 30 Hz sub-harmonic parameters estimation Over the same signal lenght of 1.25 fundamental cycles, while FBESPRIT detects the sub-harmonic frequency of 30 Hz even in very low SNR of 5 dB, the conventional ESPRIT is unable of its detection even in high SNR of 35 dB. On the other hand, the detection of frequency of 30 Hz by FB-ESPRIT in 50, 000 MCS runs for the SNR levels of 5 dB, 15 dB, 25 dB and 35 dB is respectively with mean standard deviation (Std) values of 1.9858, 0.5587, 0.1769 and 0.0650 which indicate high robustness of the estimation. Also, in term of estimation accuracy, the mean estimation bias values in MCS runs (Table 5) are respectively 0.0835, 0.0060, 0.0005 and 9.88 × 10−5 for SNRs 5 dB, 15 dB, 25 dB, and 35 dB which indicate high accuracy of FB-ESPRIT in detection of 30 Hz. The density distribution of the estimated frequencies by both ESPRIT and FB-ESPRIT depicted in Figs. 5.up-left and 6.up-left show respectively for the SNRs of 5 dB and 25 dB while FB-ESPRIT detects the frequencies in a normal distribution around the mean of 30 Hz, ESPRIT does not detect any frequency around 30 Hz.

4.2.2. 60 Hz fundamental frequency parameters estimation In estimation of 60 Hz fundamental frequency, the FB-ESPRIT mean Std values (Table 4) of MCS runs for different SNR values show 2.5 times less in the case of low SNR of 5, and around 10 times less in the case of middle SNRs of 15 dB and 25 dB and up to 25 times less mean std in the case of SNR of 35 dB compared to ESPRIT. Apart from this considerable

Table 4 The standard deviation of the estimated parameters of the signal s2 (t ) (Eq. (31)) in 50, 000 Monte Carlo runs of the FB-ESPRIT and ESPRIT at different SNR levels: 5, 15, 25, and 35 dB, while the snapshot length for both methods is 1.25 cycles. SNR

5 dB

15 dB

25 dB

35 dB

Estimator

FB-ESP.

ESPRIT

FB-ESP.

ESPRIT

FB-ESP.

ESPRIT

FB-ESP.

ESPRIT

Fund. F.: 60 Hz Sub. F.: 30 Hz Int. F.: 1650 Hz Int. F.: 2610 Hz Har. F.: 3300 Hz

0. 7001 1. 9858 1. 3686 2. 8718 12. 3301

1.6736 NA 5.4396 5.3712 24.8458

0. 0. 0. 0. 1.

1200 5587 3641 5520 2799

1.0603 NA 1.3941 2.1181 5.1071

0. 0. 0. 0. 0.

0232 1769 1349 1731 3653

0.2934 NA 0.4447 0.6586 1.3410

0. 0. 0. 0. 0.

0113 0650 0839 0935 1137

0.2934 NA 0.1584 0.2318 0.4327

Fund. 60 A.: 1 Sub. 30 A.: 0. 2 Int. 1650 A.: 0. 3 Int. 2610 A.: 0. 2 Har. 3300 A.: 0. 1

0. 0. 0. 0. 0.

0.1357 NA 0.1235 0.1155 0.0918

0. 0. 0. 0. 0.

0343 0275 0340 1228 1316

0.0970 NA 0.1293 0.1234 0.1055

0. 0. 0. 0. 0.

0108 0091 0103 0108 0116

0.0441 NA 0.0090 0.0080 0.0287

0. 0. 0. 0. 0.

0034 0028 0033 0034 0040

0.0300 NA 0.0033 0.0032 0.0040

3.5881° NA 5.8780° 8.9052° 20.7514°

0. 2. 1. 2. 5.

6096° 8323° 9570° 9173° 9310°

1.4197° NA 1.8344° 2.7577° 5.6795°

0. 0. 0. 0. 1.

1927° 8952° 7252° 9161° 8673°

0.7812° NA 0.5804° 0.8719° 1.7831°

Fund. 60 ϕ : 15° Sub. 30 ϕ : 25° Int. 1650 ϕ : 45° Int. 2610 ϕ : 60° Har. 3300 ϕ : 68°

1256 1364 1335 1377 1013

6. 1562° 18. 8042° 19. 9218° 36. 2477° 53. 2489°

5.1832° NA 19.5359° 23.1362° 43.0474°

1. 9323° 8. 9455° 6. 2611° 9. 4082° 23. 1052°

8

Electrical Power and Energy Systems 118 (2020) 105731

E. Santos, et al.

Table 5 The mean bias of the estimated parameters of the signal s2 (t ) (Eq. (31)) in 50, 000 Monte Carlo runs of the proposed FB-ESPRIT and ESPRIT at different SNR levels: 5, 15, 25, and 35 dB, while the snapshot length for both is 1.25 cycles. SNR

5 dB

15 dB

25 dB

35 dB

Estimator

FB-ESP.

ESPRIT

FB-ESP.

ESPRIT

FB-ESP.

ESPRIT

FB-ESP.

ESPRIT

Fund. F.: 60 Hz Sub. F.: 30 Hz Int. F.: 1650 Hz Int. F.: 2610 Hz Har. F.: 3300 Hz

0. 0. 0. 0. 0.

0404 0835 0022 0269 0410

0.5229 NA 0.0281 0.0545 1.2296

0. 0. 0. 0. 0.

0036 0060 0014 0006 0063

0.3138 NA 0.0761 0.0462 0.0778

0. 0018 0. 0005 0. 0009 1. 57e − 4 0. 0003

0.2342 NA 0.0131 0.0181 0.0373

0. 0001 9. 88e − 5 0. 0006 1. 92e − 4 1. 7e − 5

0.2343 NA 0.0123 0.0127 .0343

Fund. 60 A.: 1 Sub. 30 A.: 0. 2 Int. 1650 A.: 0. 3 Int. 2610 A.: 0. 2 Har. 3300 A.: 0. 1

0. 0. 0. 0. 0.

0217 1291 0655 1421 1407

0.2177 NA 0.0728 0.1352 0.1687

0. 0. 0. 0. 0.

0028 0287 0067 0541 0581

0.7840 NA 0.0378 0.1162 0.1359

0. 0. 0. 0. 0.

0011 0075 0006 0010 0023

0.8899 NA 0.0006 0.0076 0.0876

0. 0. 0. 0. 0.

0004 0022 0001 0001 0002

0.8963 NA 0.0001 0.0022 0.0983

Fund. 60 ϕ : 15° Sub. 30 ϕ : 25° Int. 1650 ϕ : 45° Int. 2610 ϕ : 60° Har. 3300 ϕ : 68°

0. 8. 0. 0. 0.

5655° 0339° 0800° 3699° 1831°

3.2984° NA 0.0776° 0.5723° 0.7135°

0. 0. 0. 0. 0.

0526° 1045° 0454° 0037° 2296°

7.0104° NA 0.0855° 0.1144° 0.3361°

0. 0. 0. 0. 0.

0058° 0040° 0132° 0014° 0773°

5.6324° NA 0.0558° 0.0540° 0.1375°

0. 0. 0. 0. 0.

0008° 0012° 0025° 0014° 0081°

5.4571° NA 0.0331° 0.0486° 0.1297°

4.3. Harmonics close to limit of sampling rate

robustness dominance in detection of fundamental frequency, the bias values in Table 5 respectively show 12.94, 87.17, 130.11, and 2343 times higher accuracy for SNRs of 5 dB, 15 dB, 25 dB and 35 dB. In amplitude estimation of the fundamental frequency 60 Hz, there is a big gap between accuracy of FB-ESPRIT and the conventional ESPRIT; respectively for SNRs of 5 dB, 15 dB, 25 dB, and 35 dB, the mean bias of estimations is 10, 280, 809, and 2241 times less for FBESPRIT. Also, in term of robustness of estimations as observed, the mean Std of estimations in MCS runs under the same above SNRs are respectively 1.0804, 2.8280, 4.0833, and 8.8235 times less for FB-ESPRIT. Regarding the robustness in estimation of the phase of the 60 Hz fundamental, under the SNRs of 15 dB, 25 dB, and 35 dB the mean Std of estimations by FB-ESPRIT is respectively 1.8569, 2.3289, and 4.0540 times less compared to ESPRIT, however under very low SNR of 5 dB ESPRIT has 1.1877 better robustness. Interestingly, regarding the accuracy in estimation of the phase of the 60 Hz fundamental, there is a long distance between two methods; under the SNRs of 5 dB, 15 dB, 25 dB, and 35 dB, the mean bias of MCS runs by FB-ESPRIT is 5.8, 133.3, 971.1, and 6821.4 times less than conventional ESPRIT.

In estimation of 3300 Hz frequency which is close to the sampling frequency limitation of 3840 Hz, the mean Std of estimations by FBESPRIT for different SNRs is 3.8056 times less compared to ESPRIT which indicates higher robustness. More interesting, there is a considerable gap between the mean estimation bias by two methods. Respectively, for the SNRs 5 dB, 10 dB, 25 dB, and 35 dB, the mean estimation bias is 30, 12.3, 124.3, and 2017.6 times less for the FBESPRIT. It is due to spectrum spread which not only brings higher robustness and accuracy of estimations but also makes the edge frequency farther than usual from the sampling limitation of 3840 Hz. While in ESPRIT the frequency of 3300 Hz is 540 Hz far from the sampling frequency limitation, in FB-ESPRIT, 3300 maps to 2160 after filter bank and down-sampling which is 1680 Hz farther from the sampling frequency limit, and considerably more than 540 Hz in ESPRIT. Interestingly, the accuracy of FB-ESPRIT in estimation of the amplitude of 3300 Hz component is considerably higher than ESPRIT as indicated by 1.1990, 2.3391, 38.0870 , and 491.5 times less mean bias of estimations in MCS runs at the SNRs of 5 dB, 15 dB, 25 dB, and 35 dB, respectively. In estimation of the phase of 3300 Hz, it behaves as explained for the middle frequencies phases.

4.2.3. Middle frequencies parameters estimation In estimation of the middle frequencies of the spectrum, 1650 Hz and 2610 Hz, in higher noisy conditions including SNRs 5 dB, 15 dB, and 25 dB, in average the mean Std of estimations by FB-ESPRIT is 4.1033 times less compared to conventional ESPRIT. In close to ideal condition of SNR equal to 35 dB, this dominance ratio reduces to 2.1836 times. Also, the accuracy of FB-ESPRIT in estimation of middle harmonics frequencies has been observed through the mean estimation bias in MCS runs between 14 to 115 times higher than conventional ESPRIT for the SNRs of 15 dB, 25 dB, and 35 dB. For SNR of 5 dB this ratio decreases to 2.026 and 12.7727 times more accuracy for the 1650 Hz and 2610 Hz, respectively. In amplitude estimation of the middle frequencies, the robustness of FB-ESPRIT is partially better than ESPRIT as the Std values of estimations in MCS runs are in average 1.2977 times less under all SNR conditions. In term of accuracy in amplitude estimation of the middle frequencies 1650 Hz and 2610 Hz, FB-ESPRIT has in average 5.1816 times less mean bias. In the case of phase estimation of the middle frequencies, while the mean Std of the estimations show similar robustness for different SNRs for both methods, in term of accuracy of the estimation the mean of bias estimation shows 2.1380, 11.4220, 14.8592, and 21.3222 times higher accuracy for FB-ESPRIT under the SNRs of 5 dB, 15 dB, 25 dB, and 35 dB.

4.4. Evaluation on real data of a photovoltaic power plant To have more rigorous evaluation, a performance comparison has been carried out on the power signal of photovoltaic (PV) power plant located in Federal University of Juiz de Fora (UFJF) shown in Fig. 7(a). Fig. 7(b) illustrates the PV signal measuring set-up where the signal is collected after inverter and before the transformer and distribution lines. The measurement has been carried out in the morning of a summer partly cloudy windy day while the sunlight exposure to the solar panels was subjected to fast fluctuations due to the wind and cloud pieces. One second of the measured PV signal current and its spectrum by DFT taken from 60 cycles signal length has been demonstrated in Fig. 8. As it is observable from the spectrum of the PV signal, it contains considerable sub-harmonic disturbances of under 25 Hz frequencies. To compare the efficiency of ESPRIT with and without association of filter bank, both techniques have been applied to 3.5 cycle length of the signal to detect and estimate the sub-harmonics and their parameters.Since we do not have any a priori knowledge about the subharmonic contents of the signal, we have used DFT applied to 60 cycles as a reference. Although DFT over 60 cycles length of the signal gives the frequency content with resolution of 1 Hz, it can not be used as a 9

Electrical Power and Energy Systems 118 (2020) 105731

E. Santos, et al.

Fig. 5. Density distributions of estimated frequencies 30 Hz (left-up), 60 Hz (right-up), 1650 Hz (left-down), 2610 Hz (middle-down), and 3300 Hz (right-down) by FB-EPRIT and ESPRIT in low SNR of 5 dB. The snapshot length for both methods is 1.25 cycles.

for the estimated phase the differences from DFT estimation are respectively 41.1606° and 11.7294°. In the case of frequencies where ESPRIT is not applicable (NA), we can see a considerable correlation between the estimations of DFT (over 60 cycles) and FB-ESPRIT (over 3.5 cycles) for the frequencies of 1, 10 , and 17 Hertz. For the other frequencies detected by FB-ESPRIT we can not judge by DFT values due to DFT limitations, as an example in the case of 4 Hz, FB-ESPRIT detects 3.8457 Hz with almost twice amplitude of 4 Hz detected by DFT. It can be due to DFT limitations where this frequency is estimated in a scattered way between 3, 4 and 5 Hertz by DFT.

reference for evaluating the estimation error. It is because of the phenomenon of frequency leakage in DFT where the frequencies located between the slots are partially leaked and detected at the adjacent slots.Table 6 shows the estimated frequencies and their parameters by DFT, and ESPRIT without and with filter bank association. Just for having a comparison, the DFT measures corresponding to the frequencies detected and estimated by two other methods are given in the table. As it can be seen, while the FB-ESPRIT detects the frequencies in vicinity of 1, 4, 7, 10, 14, 17, 21, 49, 54 , and 60 Hertz, ESPRIT which uses the same length of the signal just can detect the frequencies around 7 and 60 Hertz. About the amplitude and phase of the detected frequencies, due to DFT limitations, just being close the estimated parameters by DFT has been used as a criterion. As it is observable in the case of frequencies 7 Hz and 60 Hz detected by both ESPRIT with and without FB, the parameters detected by FB-ESPRIT are closer to DFT estimations. This issue is especially firmly observable in amplitude and phase of 7 Hz where the differences from DFT estimation in amplitude case are respectively 13.6743 and 3.2945 by ESPRIT and FB-ESPRIT, and

4.5. Complexity evaluation Over three case study signals, both FB-ESPRIT and ESPRIT are implemented while the signal snapshot length n and the autocorrelation vector length M are individually set for each of them in a way both methodologies give equivalent efficiencies in term of accuracy and robustness in detection of the frequencies. Table 7 shows the 10

Electrical Power and Energy Systems 118 (2020) 105731

E. Santos, et al.

Fig. 6. Density distributions of estimated frequencies 30 Hz (left-up), 60 Hz (right-up), 1650 Hz (left-down), 2610 Hz (middle-down), and 3300 Hz (right-down) by FB-EPRIT and ESPRIT in low SNR of 25 dB. The snapshot length for both methods is 1.25 cycles.

Fig. 7. (a) Solar power plant of Federal University of Juiz de Fora (UFJF), (b) Power measuring set-up at UFJF solar power plant. 11

Electrical Power and Energy Systems 118 (2020) 105731

E. Santos, et al.

Fig. 8. (a) Current signal of UFJF PV power plant in a partly cloudy windy morning, (b) UFJF PV power plant signal spectrum by FFT of 60 cycles. Table 6 Evaluation on current output of UFJF photovoltaic power plant. Frequency

Amplitude

Phase

DFT

ESPRIT

FB-ESP.

DFT

ESPRIT

FB-ESP.

DFT

ESPRIT

FB-ESP.

1 4 7 10 14 17 21 49 54 60

NA NA 7.3850 NA NA NA NA NA NA 60.0511

1.0317 3.8457 7.0823 10.1324 13.7883 17.5203 21.4276 49.3723 54.4575 59.9769

21.7738 0.9707 15.0595 4.1054 1.3002 0.1866 0.7273 0.6555 0.5544 34.9640

NA NA 28.7338 NA A NA NA NA NA 33.6139

20.1455 2.0252 18.3540 5.1732 0.4235 0.4139 0.3148 1.0914 0.6087 33.7309

−21.1484 −32.4593 −118.0406 121.1077 −114.9467 −105.2372 84.6625 −114.1721 162.7442 29.7229

NA NA −159.2012 NA NA NA NA NA NA 25.9615

−26.5236 97.0703 −129.7700 117.0557 −78.7982 −73.0094 −13.4289 −160.1991 62.2054 34.0749

optimal solution via trade-off between efficiency and complexity.

computational load comparison of FB-ESPRIT and ESPRIT over the condition mentioned above. As it is observable from the table, the FB-ESPRIT gains equivalent efficiency with much lower computational load than ESPRIT. Interestingly, in SNR condition of 25 dB, it approves the theoretical L2 times complexity improvement by FB-ESPRIT as in Eq. (26). However, as Table 7 indicates, in lower SNRs of 5 dB and 15 dB, this ratio is less than theoretical expectation of 16, and in practice 7.7500 times lower complexity is observed.

6. Conclusion ESPRIT with the association of a filter bank (FB-ESPRIT) has been suggested for multiple times improvement of the efficiency of ESPRIT in the detection and parameter estimation of power line harmonics, interharmonics, and sub-harmonic components. FB-ESPRIT is indeed a parallel structure of the individual implementation of ESPRIT to the sub-bands of the signal after spreading their spectrum over the full bandwidth of the signal. The equal bandwidth sub-bands are provided by a uniform filter bank composed of L filters, and the spectrum spread of the sub-band is via down-sampling by the same parameter L. When each small section of the spectrum is widened L times and goes through parameter estimation, theoretically multiple times efficiency is expected. Furthermore, since the components are divided between the sub-bands, each sub-band individual ESPRIT implementation requires L times shorter autocorrelation vector. The manuscript theoretically illustrates the L2 times lower complexity order of FB-ESPRIT while its parameters set to give equivalent efficiency to ESPRIT. The mean result of Monte Carlo simulation runs approves the average L times higher efficiency of FB-ESPRIT in both terms of accuracy and robustness compared to ESPRIT, as well as the theoretical L2 lower complexity in

5. Future scope The proposed association of filter bank to the ESPRIT improves its robustness, accuracy, and speed in parameter estimation of signal power frequency components. Besides, the manuscript theoretically and practically illustrates that in the case that FB-ESPRIT and ESPRIT are set for equivalent performance, FB-ESPRIT has L2 times lower complexity burden. According to this advantage of faster performance over a shorter length of the signal, FB-ESPRIT has the potential for efficient tracking of the signal spectrum variations under non-stationary condition. The future scope of this research work is the analysis of filter bank association with ESPRIT and other parametric techniques for spectrum estimation of non-stationary signals where the goal is finding a Pareto

Table 7 Computation cost for estimation of frequency of the 1650 Hz Inter-harmonic component of s2 (t ) , while ESPRIT and FB-ESPRIT give their closest efficiency in term of accuracy and robustness of estimation in 50, 000 MCS runs. (One cycle of s (t ) is 32 samples.)

12

Electrical Power and Energy Systems 118 (2020) 105731

E. Santos, et al.

Declaration of Competing Interest

practice.

The authors declared that there is no conflict of interest. Appendix A. Appendix: Proof of Lemma 2 From Eq. (27), the complexity order of FB-ESPRIT with a filter-bank composed of L filters of order p comparable to O (n) the order of signal snapshot length, O (p) ≈ O (n) , can be expressed as

OFB-ESPRIT (L) =

1 (9n3 + 12Mn2) + Ln2. L2

(32)

The L which brings the minimum complexity order for FB-ESPRIT is the solution of the following optimization problem

Minimize OFB-ESPRIT (L),

L ∈ {2i}, i ∈ ,

(33)

where  is the set of natural numbers. Acquiring the optimum L as a power two, L ∈ {2i} , is for the sake of being beneficiary of technical advantages of polyphase technique in filter-bank implementation [40]. To solve the above problem we search for a L min ∈ + minimizing OFB-ESPRIT (L), where + is the set of positive real numbers, then we acquire the closest power of two as the optimum solution, Lopt . To find L min ∈ + minimizing

OFB-ESPRIT (L) , we solve the stationary condition of d d (OFB-ESPRIT (L))=dL dL

(

9n3 + 12Mn2 + L3n2 L2

n2 (L3 − 2(9n + 12M )) = L3

Thus

L min =

3

dOFB-ESPRIT (L) dL

= 0 or

)

= 0.

(34)

6(3n + 4M ).

(35)

The second derivative of OFB-ESPRIT (L) at the above-obtained L min is d2 (OFB-ESPRIT (L)) dL2 Lmin

=

d2 dL2

(

L3n2 − 2n2 (09n + 12M ) L3

=

6(9n3 + 12Mn2) 4 Lmin

)

Lmin

> 0,

(36)

which indicates the L min truly gives the minimum of the function. The lemma is proved as the closest power of two to L min is the optimum solution as follows

Lopt = 2⌊log2Lmin⌉ = 2⌊log2

3 6(3n + 4M ) ⌉

.

(37)

[15] Marple L.. Spectral line analysis by pisarenko and prony methods. In: ICASSP’79. IEEE international conference on acoustics, speech, and signal processing. vol. 4. IEEE, 1979. p. 159–61. [16] Jain SK, Singh S. Harmonics estimation in emerging power system: Key issues and challenges. Electric Power Syst Res 2011;81(9):1754–66. [17] Schmidt RO. A signal subspace approach to multiple emitter location and spectral estimation Ph.D. dissertation Stanford Univ.; 1981. [18] Schmidt R. Multiple emitter location and signal parameter estimation. IEEE Trans Antennas Propagat 1986;34(3):276–80. [19] Roy R, Paulraj A, Kailath T. Esprit–a subspace rotation approach to estimation of parameters of cisoids in noise. IEEE Trans Acoust Speech Signal Process 1986;34(5):1340–2. [20] Roy RH, Kailath T. Esprit–estimation of signal parameters via rotational invariance techniques. Opt Eng 1990;29(4):296–314. [21] Ottersten B, Viberg M, Kailath T. Performance analysis of the total least squares esprit algorithm. IEEE Trans Signal Process 1991;39(5):1122–35. [22] Swindlehurst AL, Ottersten B, Roy R, Kailath T. Multiple invariance esprit. IEEE Trans Signal Process 1992;40(4):867–81. [23] Zhang X-D, Liang Y-C. Prefiltering-based esprit for estimating sinusoidal parameters in non-gaussian arma noise. IEEE Trans Signal Process 1995;43(1):349–53. [24] Zeineldin H, Abdel-Galil T, El-Saadany E, Salama M. Islanding detection of grid connected distributed generators using tls-esprit. Electric Power Syst Res 2007;77(2):155–62. [25] Tripathy P, Srivastava S, Singh S. A modified tls-esprit-based method for low-frequency mode identification in power systems utilizing synchrophasor measurements. IEEE Trans Power Syst 2010;26(2):719–27. [26] Zolfaghari R, Shrivastava Y, Agelidis VG. Evaluation of windowed esprit virtual instrument for estimating power quality indices. Electric Power Syst Res 2012;83(1):58–65. [27] IEEE. Ieee standard definitions for the measurement of electric power qualities under sinusoidal, nonsinusoidal, balanced, or unbalanced conditions, 2010. [28] Jain SK, Singh S. Impact of signal attributes on autocorrelation matrix dimension for smart grid solutions. ISGT2011-India. IEEE; 2011. p. 43–8. [29] Hai Y, Chen J-Y. Exact model order esprit technique for frequency and power estimation. 2012 IEEE 14th international conference on communication technology. IEEE; 2012. p. 1281–5. [30] Jain SK, Singh S, Singh J. An adaptive time-efficient technique for harmonic estimation of nonstationary signals. IEEE Trans Industr Electron 2012;60(8):3295–303.

References [1] Ribeiro PF, Duque CA, Ribeiro PM, Cerqueira AS. Power systems signal processing for smart grids. John Wiley & Sons; 2013. [2] Terzija VV, Djuric MB, Kovacevic BD. Voltage phasor and local system frequency estimation using newton type algorithm. IEEE Trans Power Deliv 1994;9(3):1368–74. [3] Saini MK, Kapoor R. Classification of power quality events–a review. Int J Electrical Power Energy Syst 2012;43(1):11–9. [4] Terzija VV. Adaptive underfrequency load shedding based on the magnitude of the disturbance estimation. IEEE Trans Power Syst 2006;21(3):1260–6. [5] Liu G, Zhou J, Jia B, He F, Yang Y, Sun N. Advance short-term wind energy quality assessment based on instantaneous standard deviation and variogram of wind speed by a hybrid method. Appl Energy 2019;238:643–67. [6] Novanda H, Regulski P, Stanojevic V, Terzija V. Assessment of frequency and harmonic distortions during wind farm rejection test. IEEE Trans Sustain Energy 2013;4(3):698–705. [7] Parwal A, Fregelius M, Temiz I, Göteman M, de Oliveira JG, Boström C, Leijon M. Energy management for a grid-connected wave energy park through a hybrid energy storage system. Appl Energy 2018;231:399–411. [8] Taghizadeh S, Hossain M, Lu J, Water W. A unified multi-functional on-board ev charger for power-quality control in household networks. Appl Energy 2018;215:186–201. [9] Ray PK, Puhan PS, Panda G. Real time harmonics estimation of distorted power system signal. Int J Electr Power Energy Syst 2016;75:91–8. [10] Chen C-I, Chang GW. Virtual instrumentation and educational platform for timevarying harmonic and interharmonic detection. IEEE Trans Industr Electron 2010;57(10):3334–42. [11] Lin HC. Development of interharmonics identification using enhanced-fft algorithm. J Eng 2017;2017(7):333–42. [12] Jin Z, Zhang H, Shi F, Shi X. A robust and adaptive interharmonics detection method. 2017 IEEE conference on energy internet and energy system integration (EI2). IEEE; 2017. p. 1–5. [13] Chang GW, Chen C-I. An accurate time-domain procedure for harmonics and interharmonics detection. IEEE Trans Power Deliv 2010;25(3):1787–95. [14] Pisarenko VF. The retrieval of harmonics from a covariance function. Geophys J Int 1973;33(3):347–66.

13

Electrical Power and Energy Systems 118 (2020) 105731

E. Santos, et al.

[36] Jain P, Tiwari AK, Jain SK. Harmonic source identification in distribution system using estimation of signal parameters via rotational invariance technique-total harmonic power method. Trans Inst Measur Control 2018;40(12):3415–23. [37] Santos E, Duque CA, Cerqueira AS. Harmonic and inter-harmonic estimation using esprit associated with filter bank. In: 2018 18th International conference on harmonics and quality of power (ICHQP). IEEE, 2018. p. 1–6. [38] Zygarlicki J, Zygarlicka M, Mroczka J, Latawiec KJ. A reduced prony’s method in power-quality analysis—parameters selection. IEEE Trans Power Deliv 2010;25(2):979–86. [39] Moon TK, Stirling WC. Mathematical methods and algorithms for signal processing vol. 1. Upper Saddle River, NJ: Prentice hall; 2000. [40] Mitra SK, Kuo Y. Digital signal processing: a computer-based approach vol. 2. New York: McGraw-Hill; 2006.

[31] Banerjee P, Srivastava S. An effective dynamic current phasor estimator for synchrophasor measurements. IEEE Trans Instrum Meas 2015;64(3):625–37. [32] Yang S. Frequency estimation of harmonics with multiplicative noise using esprit. In: 2015 8th International congress on image and signal processing (CISP). IEEE, 2015. p. 1251–55. [33] Alfieri L, Carpinelli G, Bracale A. New esprit-based method for an efficient assessment of waveform distortions in power systems. Electric Power Syst Res 2015;122:130–9. [34] Jain SK. Memo-esprit based synchronized harmonic measurement for online applications. In: 2016 International conference on emerging trends in electrical electronics & sustainable energy systems (ICETEESES). IEEE, 2016. p. 152–8. [35] Jain SK, Jain P, Singh SN. A fast harmonic phasor measurement method for smart grid applications. IEEE Trans Smart Grid 2016;8(1):493–502.

14