# Performance Comparison of Windowed Interpolation FFT and Quasisynchronous Sampling Algorithm for Frequency Estimation

Authors: He Wen; Huifang Dai; Zhaosheng Teng; Yuxiang Yang; Fuhai Li
Journal: Mathematical Problems in Engineering
 (2014)

---

## Abstract

The DFT-based frequency estimations have inherent limitations such as spectral leakage and picket-fence effect due to asynchronous sampling. This paper focuses on comparing the windowed interpolation FFT (WIFFT) and quasisynchronous sampling algorithm (QSSA) for frequency estimation. The WIFFT uses windows to reduce spectral leakage and employs interpolation algorithm to eliminate picket-fence effect. And the QSSA utilizes quasisynchronous weighted iterations for frequency estimation. The accuracy and time complexity of WIFFT and QSSA are theoretically studied. Computer simulations of frequency estimations with noise and fluctuations by using WIFFT and QSSA are performed. Simulations results show that the wideband noise sensitivity of QSSA is lower than that of WIFFT. However, the WIFFT exhibits less time complexity than QSSA.


---

## Body

## 1. Introduction

Frequency, as we all know, is an important parameter in power system monitoring, control, and protection [1–3]. In power systems, typical applications of frequency estimation are for protection against loss of synchronism, underfrequency relaying, power-quality monitoring, and power system stabilization [4–6]. However, fluctuation of fundamental frequency has become a curse to the safe operation of power system which deteriorates the accuracy of frequency estimation. Thus, accurate frequency estimation of distorted and noisy signals in industrial power systems is a challenge that has attracted much attention.Recently, various frequency estimation algorithms have been proposed in the literature, which can be classified as time-domain (parametric) and frequency-domain (nonparametric) methods [7–10]. Time-domain methods, including iterative method, Kalman filtering, and neural network, are model-based and have very good selectivity and statistical efficiency [11, 12]. However, time-domain methods employ computationally intensive algorithms and require very good model agreement with a real multicomponent signal.Frequency-domain approaches are based on the well-known discrete Fourier transform (DFT), which exhibit a lower computational effort due to the availability of fast Fourier transform (FFT). To obtain the accurate frequency, a direct application of the DFT has the requirement of the periodicity of signal and the synchronization of sampling process [13, 14]. However, fluctuations of signals in power system, including fundamental frequency instability, fundamental variations, time varying harmonics, and harmonic interference, make synchronous sampling unattainable. The performances of DFT-based frequency estimations have inherent limitations such as spectral leakage and picket-fence effect due to asynchronous sampling. It is well known that a considerable reduction of spectral leakage can be achieved by weighting samples with a suitable window, and the picket-fence effect can be eliminated by adopting interpolation FFT algorithms with dual-spectrum line or multispectrum line. Although the drawback of DFT-based frequency estimations can be improved by the windowed interpolation FFT (WIFFT) [15–17], additional computational costs will be introduced. Also, traditional DFT-based approaches require more than one cycle of the signal samples for computation. Therefore, the performance of real-time monitoring may be deteriorated.Quasisynchronous sampling algorithm (QSSA) is a popular approach for reducing spectral leakage [18–20]. It employs quasisynchronous weighted iterations for the spectrum analysis and estimates frequency by the phase angle difference between two periodic signals. The QSSA has the advantages of requiring simple hardware and high accuracy. However, the time complexity of the QSSA has not been analyzed. Also, the noise sensitivity of the QSSA is still not clear.The worldwide increasing applications of frequency measurement have made the choosing of frequency estimation methods a vital issue than ever before. Unfortunately, in the scientific literature the accuracy achieved by the WIFFT has not yet been compared with that of the QSSA. This is the goal of this paper. In Section2, the WIFFT and QSSA for frequency estimation are discussed in brief. Comparisons of the WIFFT and QSSA, computational burden, are presented in Section 3. Simulations results are provided in Section 4. The conclusions are drawn in Section 5.


## 2. Frequency Estimation Algorithm

In order to better compare the performance of the QSSA and WIFFT, brief introductions are provided in this section.


### 2.1. Quasisynchronous Sampling Algorithm

Consider a sine-wavex
(
m
) of amplitude A, frequency f
=
1
/
T, phase φ, and sampled at frequency f
s
=
1
/
T
s, which can be expressed as
(1)
x
(
m
)
=
A
sin
⁡
(
2
π
f
m
f
s
+
φ
)
,
m
∈
[
0
,
+
∞
)
.Let the sequencex
1
(
n
) be selected from x
(
m
) that start at m
=
b
1 and consist of the L
N
+
1 samples {
x
1
(
1
)
=
x
(
b
1
)
,
x
1
(
2
)
=
x
(
b
1
+
1
)
,
…
,
x
1
(
L
N
+
1
)
=
x
(
b
1
+
L
N
)
}. According to the definition of Fourier transform of the x
1
(
n
) it can be calculated as
(2)
X
(
k
)
|
b
1
=
2
L
N
+
1
∑
n
=
1
L
N
+
1
x
1
(
n
)
exp
⁡
(
-
j
2
π
n
L
N
+
1
k
)
,
d
d
d
d
d
d
d
d
d
d
d
d
d
d
d
d
i
i
d
d
d
d
k
∈
[
0
,
L
N
]
.
When k
=
1, we can obtain the real part and image part of the second spectral line of the DFT of x
1
(
n
), which can, respectively, be calculated as
(3)
R
|
b
1
=
2
L
N
+
1
∑
n
=
1
L
N
+
1
x
1
(
n
)
cos
⁡
(
2
π
n
L
N
+
1
)
,
I
|
b
1
=
2
L
N
+
1
∑
n
=
1
L
N
+
1
x
1
(
n
)
sin
(
2
π
n
L
N
+
1
)
.If the sine-wave is sampled synchronously, the frequency can be accurately obtained from the from the spectral lines, that is,R
|
b
1 and I
|
b
1. Otherwise, the spectral leakage caused by the asynchronous sampling will deteriorate the accuracy of frequency estimation. To reduce the effect of the asynchronous sampling, a quasisynchronous iteration process for calculating R
|
b
1 and I
|
b
1 can be applied as shown in Figure 1.Figure 1
TheLth iteration process based on the quasisynchronous sampling algorithm with L
N
+
1 samples.In Figure1, the upper index l represents the number of recursions (l
=
1,2
,
3
,
…
,
L) and q
(
n
) represents the original sequence of length L
N
+
1. The weighted recursion can be expressed as [20]
(4)
Q
i
1
=
1
∑
r
=
i
N
+
i
ρ
r
∑
r
=
i
N
+
i
ρ
r
q
(
r
)
Q
i
l
=
1
∑
r
=
i
N
+
i
ρ
r
∑
r
=
i
N
+
i
ρ
r
Q
r
l
-
1
,
i
=
1,2
,
…
,
(
L
-
l
)
N
+
1
,
where ρ
r
(
r
=
i
,
i
+
1
,
…
,
i
+
N
) are weights according to the window function.The real partR
|
b
1 and image part I
|
b
1 of the DFT of the x
1
(
n
) can be obtained from Figure 1 when the x
1
(
n
)
sin
⁡
(
2
π
n
/
(
L
N
+
1
)
) and x
1
(
n
)
cos
⁡
(
2
π
n
/
(
L
N
+
1
)
) are, respectively, submitted into q
(
n
). The weighted recursion results of R
|
b
1 and I
|
b
1 can be expressed as
(5)
R
|
b
1
=
Q
i
=
b
1
L
|
x
1
cos
⁡
(
)
≈
(
γ
m
m
)
L
A
sin
⁡
(
2
π
b
1
f
f
s
+
φ
)
I
|
b
1
=
Q
i
=
b
1
L
|
x
1
sin
(
)
≈
(
γ
m
m
)
L
A
cos
⁡
⁡
(
2
π
b
1
f
f
s
+
φ
)
,
where γ
m
m is the iteration factor.Another sequencex
2
(
n
) can also be selected from x
(
m
) that start at m
=
b
2 and consist of the L
N
+
1 sample {
x
2
(
1
)
=
x
(
b
2
)
,
x
2
(
2
)
=
x
(
b
2
+
1
)
,
…
,
x
2
(
L
N
+
1
)
=
x
(
b
2
+
L
N
)
}. Using the same procedure as mentioned above, the real part and image part of the second spectral line of the DFT of x
2
(
n
), that is, R
|
b
2 and I
|
b
2, can be obtained similarly.The equivalent initial phase angles ofx
1(n) and x
2(n) can thus be calculated as
(6)
ϑ
|
b
1
≈
arctan
(
R
|
b
1
I
|
b
1
)
=
2
π
b
1
f
f
s
+
φ
,
ϑ
|
b
2
≈
arctan
(
R
|
b
2
I
|
b
2
)
=
2
π
b
2
f
f
s
+
φ
.
By subtracting ϑ
|
b
2 from ϑ
|
b
1, we can calculate the frequency of x
(
m
) as
(7)
f
=
(
ϑ
|
b
1
-
ϑ
|
b
2
)
f
s
2
π
(
b
1
-
b
2
)
.


### 2.2. Windowed Interpolation FFT Algorithm

The windowed interpolated FFT is an effective approach to eliminate the spectral leakage and picket-fence effect [21, 22]. The signal as shown in formula (1) is multiplied by a window function of length N; that is, x
w
(
n
)
=
x
(
m
)
w
(
n
) with n
=
1,2
,
…
,
N. The DFT of the windowed signal can be expressed as
(8)
X
(
k
)
=
A
2
j
[
e
j
φ
W
N
(
k
-
λ
)
-
e
-
j
φ
W
N
(
k
+
λ
)
]
,
d
d
d
d
d
d
d
d
d
d
d
d
d
d
d
d
d
d
d
k
∈
[
0
,
N
-
1
]
,
where W
(
) is the DFT of the adopted window w
(
n
) and λ
=
f
N
/
f
s is the signal frequency divided by the frequency resolution Δ
f
=
f
s
/
N. The λ would be an integer number if synchronous sampling is achieved. Otherwise the λ can be written as the sum of an integer part η and a fractional part δ under asynchronous sampling; that is, λ
=
η
+
δ. Without loss of generality, assume that the minimum distance between spectral lines is larger than the main lobe width of the adopted window, which is sufficient for avoiding spectral interferences from the imaginary part of the spectrum. So, the imaginary part of the spectrum can be ignored; (8) can thus be simplified as
(9)
X
(
k
)
=
A
2
j
e
j
φ
W
(
k
-
λ
)
.The WIFFT applies a peak-search procedure toX
(
k
) and thus obtains two spectral lines with the maximum and the second maximum magnitudes, that is, the k
1th and k
2th spectral lines with k
2
=
k
1
+
1. As shown in Figure 2, the λ should be larger than k
1 and smaller than k
2
(
k
1
≤
λ
≤
k
2
=
k
1
+
1
). By defining a coefficient α
=
λ
-
k
1
-
0
.
5 with -
0.5
<
α
<
0
.
5, we can get k
1
-
λ
=
-
0.5
-
α and k
2
-
λ
=
0.5
-
α. Thus the magnitudes of the k
1th and k
2th spectral lines are |
X
(
k
1
)
|
=
A
|
W
(
-
0.5
-
α
)
|
/
2 and |
X
(
k
2
)
|
=
A
|
W
(
0.5
-
α
)
|
/
2, respectively. Since the inverse of an odd and symmetric function does not contain any even order items, an odd and symmetric coefficient β can be defined as a function of α to simplify the calculation procedure, which can be written as
(10)
β
=
|
X
(
k
2
)
|
-
|
X
(
k
1
)
|
|
X
(
k
2
)
|
+
|
X
(
k
1
)
|
=
|
X
(
0.5
-
α
)
|
-
|
X
(
-
0.5
-
α
)
|
|
X
(
0.5
-
α
)
|
+
|
X
(
-
0.5
-
α
)
|
=
|
W
(
0.5
-
α
)
|
-
|
W
(
-
0.5
-
α
)
|
|
W
(
0.5
-
α
)
|
+
|
W
(
-
0.5
-
α
)
|
=
g
(
α
)
.Figure 2
Discrete spectrum of nonsynchronized sampling.As the expression of (10) is complex, it is hard to get the inverse of (10) directly. The least square method is a good choice to obtain an approximate polynomial. Firstly, a set of data points (β
i, α
i) can be obtained by giving a series of values α
i by using (10). Then, the β
i and α
i are referred to as the independent variable and dependent variable, respectively. Finally, the least square method is applied to obtain the approximate polynomial
(11)
α
=
g
-
1
(
β
)
=
∑
j
=
0
J
-
1
v
j
β
j
,
where J and v
j are the order and coefficients of the approximate polynomial, respectively. Once the coefficients v
j are determined, the α can be calculated by submitting β into (11). Consequently, the signal frequency f can be estimated by
(12)
f
=
λ
f
s
N
=
(
k
1
+
0.5
+
α
)
f
s
N
.


## 2.1. Quasisynchronous Sampling Algorithm

Consider a sine-wavex
(
m
) of amplitude A, frequency f
=
1
/
T, phase φ, and sampled at frequency f
s
=
1
/
T
s, which can be expressed as
(1)
x
(
m
)
=
A
sin
⁡
(
2
π
f
m
f
s
+
φ
)
,
m
∈
[
0
,
+
∞
)
.Let the sequencex
1
(
n
) be selected from x
(
m
) that start at m
=
b
1 and consist of the L
N
+
1 samples {
x
1
(
1
)
=
x
(
b
1
)
,
x
1
(
2
)
=
x
(
b
1
+
1
)
,
…
,
x
1
(
L
N
+
1
)
=
x
(
b
1
+
L
N
)
}. According to the definition of Fourier transform of the x
1
(
n
) it can be calculated as
(2)
X
(
k
)
|
b
1
=
2
L
N
+
1
∑
n
=
1
L
N
+
1
x
1
(
n
)
exp
⁡
(
-
j
2
π
n
L
N
+
1
k
)
,
d
d
d
d
d
d
d
d
d
d
d
d
d
d
d
d
i
i
d
d
d
d
k
∈
[
0
,
L
N
]
.
When k
=
1, we can obtain the real part and image part of the second spectral line of the DFT of x
1
(
n
), which can, respectively, be calculated as
(3)
R
|
b
1
=
2
L
N
+
1
∑
n
=
1
L
N
+
1
x
1
(
n
)
cos
⁡
(
2
π
n
L
N
+
1
)
,
I
|
b
1
=
2
L
N
+
1
∑
n
=
1
L
N
+
1
x
1
(
n
)
sin
(
2
π
n
L
N
+
1
)
.If the sine-wave is sampled synchronously, the frequency can be accurately obtained from the from the spectral lines, that is,R
|
b
1 and I
|
b
1. Otherwise, the spectral leakage caused by the asynchronous sampling will deteriorate the accuracy of frequency estimation. To reduce the effect of the asynchronous sampling, a quasisynchronous iteration process for calculating R
|
b
1 and I
|
b
1 can be applied as shown in Figure 1.Figure 1
TheLth iteration process based on the quasisynchronous sampling algorithm with L
N
+
1 samples.In Figure1, the upper index l represents the number of recursions (l
=
1,2
,
3
,
…
,
L) and q
(
n
) represents the original sequence of length L
N
+
1. The weighted recursion can be expressed as [20]
(4)
Q
i
1
=
1
∑
r
=
i
N
+
i
ρ
r
∑
r
=
i
N
+
i
ρ
r
q
(
r
)
Q
i
l
=
1
∑
r
=
i
N
+
i
ρ
r
∑
r
=
i
N
+
i
ρ
r
Q
r
l
-
1
,
i
=
1,2
,
…
,
(
L
-
l
)
N
+
1
,
where ρ
r
(
r
=
i
,
i
+
1
,
…
,
i
+
N
) are weights according to the window function.The real partR
|
b
1 and image part I
|
b
1 of the DFT of the x
1
(
n
) can be obtained from Figure 1 when the x
1
(
n
)
sin
⁡
(
2
π
n
/
(
L
N
+
1
)
) and x
1
(
n
)
cos
⁡
(
2
π
n
/
(
L
N
+
1
)
) are, respectively, submitted into q
(
n
). The weighted recursion results of R
|
b
1 and I
|
b
1 can be expressed as
(5)
R
|
b
1
=
Q
i
=
b
1
L
|
x
1
cos
⁡
(
)
≈
(
γ
m
m
)
L
A
sin
⁡
(
2
π
b
1
f
f
s
+
φ
)
I
|
b
1
=
Q
i
=
b
1
L
|
x
1
sin
(
)
≈
(
γ
m
m
)
L
A
cos
⁡
⁡
(
2
π
b
1
f
f
s
+
φ
)
,
where γ
m
m is the iteration factor.Another sequencex
2
(
n
) can also be selected from x
(
m
) that start at m
=
b
2 and consist of the L
N
+
1 sample {
x
2
(
1
)
=
x
(
b
2
)
,
x
2
(
2
)
=
x
(
b
2
+
1
)
,
…
,
x
2
(
L
N
+
1
)
=
x
(
b
2
+
L
N
)
}. Using the same procedure as mentioned above, the real part and image part of the second spectral line of the DFT of x
2
(
n
), that is, R
|
b
2 and I
|
b
2, can be obtained similarly.The equivalent initial phase angles ofx
1(n) and x
2(n) can thus be calculated as
(6)
ϑ
|
b
1
≈
arctan
(
R
|
b
1
I
|
b
1
)
=
2
π
b
1
f
f
s
+
φ
,
ϑ
|
b
2
≈
arctan
(
R
|
b
2
I
|
b
2
)
=
2
π
b
2
f
f
s
+
φ
.
By subtracting ϑ
|
b
2 from ϑ
|
b
1, we can calculate the frequency of x
(
m
) as
(7)
f
=
(
ϑ
|
b
1
-
ϑ
|
b
2
)
f
s
2
π
(
b
1
-
b
2
)
.


## 2.2. Windowed Interpolation FFT Algorithm

The windowed interpolated FFT is an effective approach to eliminate the spectral leakage and picket-fence effect [21, 22]. The signal as shown in formula (1) is multiplied by a window function of length N; that is, x
w
(
n
)
=
x
(
m
)
w
(
n
) with n
=
1,2
,
…
,
N. The DFT of the windowed signal can be expressed as
(8)
X
(
k
)
=
A
2
j
[
e
j
φ
W
N
(
k
-
λ
)
-
e
-
j
φ
W
N
(
k
+
λ
)
]
,
d
d
d
d
d
d
d
d
d
d
d
d
d
d
d
d
d
d
d
k
∈
[
0
,
N
-
1
]
,
where W
(
) is the DFT of the adopted window w
(
n
) and λ
=
f
N
/
f
s is the signal frequency divided by the frequency resolution Δ
f
=
f
s
/
N. The λ would be an integer number if synchronous sampling is achieved. Otherwise the λ can be written as the sum of an integer part η and a fractional part δ under asynchronous sampling; that is, λ
=
η
+
δ. Without loss of generality, assume that the minimum distance between spectral lines is larger than the main lobe width of the adopted window, which is sufficient for avoiding spectral interferences from the imaginary part of the spectrum. So, the imaginary part of the spectrum can be ignored; (8) can thus be simplified as
(9)
X
(
k
)
=
A
2
j
e
j
φ
W
(
k
-
λ
)
.The WIFFT applies a peak-search procedure toX
(
k
) and thus obtains two spectral lines with the maximum and the second maximum magnitudes, that is, the k
1th and k
2th spectral lines with k
2
=
k
1
+
1. As shown in Figure 2, the λ should be larger than k
1 and smaller than k
2
(
k
1
≤
λ
≤
k
2
=
k
1
+
1
). By defining a coefficient α
=
λ
-
k
1
-
0
.
5 with -
0.5
<
α
<
0
.
5, we can get k
1
-
λ
=
-
0.5
-
α and k
2
-
λ
=
0.5
-
α. Thus the magnitudes of the k
1th and k
2th spectral lines are |
X
(
k
1
)
|
=
A
|
W
(
-
0.5
-
α
)
|
/
2 and |
X
(
k
2
)
|
=
A
|
W
(
0.5
-
α
)
|
/
2, respectively. Since the inverse of an odd and symmetric function does not contain any even order items, an odd and symmetric coefficient β can be defined as a function of α to simplify the calculation procedure, which can be written as
(10)
β
=
|
X
(
k
2
)
|
-
|
X
(
k
1
)
|
|
X
(
k
2
)
|
+
|
X
(
k
1
)
|
=
|
X
(
0.5
-
α
)
|
-
|
X
(
-
0.5
-
α
)
|
|
X
(
0.5
-
α
)
|
+
|
X
(
-
0.5
-
α
)
|
=
|
W
(
0.5
-
α
)
|
-
|
W
(
-
0.5
-
α
)
|
|
W
(
0.5
-
α
)
|
+
|
W
(
-
0.5
-
α
)
|
=
g
(
α
)
.Figure 2
Discrete spectrum of nonsynchronized sampling.As the expression of (10) is complex, it is hard to get the inverse of (10) directly. The least square method is a good choice to obtain an approximate polynomial. Firstly, a set of data points (β
i, α
i) can be obtained by giving a series of values α
i by using (10). Then, the β
i and α
i are referred to as the independent variable and dependent variable, respectively. Finally, the least square method is applied to obtain the approximate polynomial
(11)
α
=
g
-
1
(
β
)
=
∑
j
=
0
J
-
1
v
j
β
j
,
where J and v
j are the order and coefficients of the approximate polynomial, respectively. Once the coefficients v
j are determined, the α can be calculated by submitting β into (11). Consequently, the signal frequency f can be estimated by
(12)
f
=
λ
f
s
N
=
(
k
1
+
0.5
+
α
)
f
s
N
.


## 3. Comparisons of Time Complexity

As the real-time frequency estimation is of great importance in power systems, it is of interest to know the time complexity of different algorithms of frequency estimations.The time complexity of an algorithm quantifies the amount of time taken by an algorithm to run as a function of the length of the string representing the input and is commonly expressed usingO
(
) notation, which excludes coefficients and lower order terms. Time complexity is commonly estimated by counting the number of elementary operations performed by the algorithm, where an elementary operation takes a fixed amount of time to perform. Thus the amount of time taken and the number of elementary operations performed by the algorithm differ by at most a constant factor. The elementary operations include addition, subtraction, multiplication, and division. Normally, divisions and subtractions are, respectively, counted into multiplications and additions. Comparisons of time complexities of WIFFT and QSSA are provided in Table 1.Table 1
Comparisons of time complexities.
QSSA
WIFFT
Addition
L
2
N
3
(
L
N
+
1
)
log
⁡
2
(
L
N
+
1
)
Multiplication
L
N
/
2
+
L
N
2
/
2
+
L
+
1
2
(
L
N
+
1
)
log
⁡
2
(
L
N
+
1
)
Time complexity
O
(
N
2
)
O
(
N
log
⁡
2
N
)The time complexity of the WIFFT is caused by the FFT operations. Thus, according to the definition of FFT, the WIFFT withL
N
+
1 samples requires 3
(
L
N
+
1
)
log
⁡
2
⁡
(
L
N
+
1
) additions and 2
(
L
N
+
1
)
log
⁡
2
⁡
(
L
N
+
1
) multiplications. So, the time complexity of the WIFFT is O
(
N
log
⁡
⁡
N
).The time complexity of the QSSA is caused by the iteration procedure, which is listed in Pseudocode1. As shown in Figure 1 and (4), for the QSSA with L
N
+
1 samples, there are L
N
/
2
+
L
N
2
/
2
+
L
+
1 multiplications and L
2
N additions. Thus the time complexity of the QSSA is O
(
N
2
).Pseudocode 1:Pseudocode description of QSSA.
Algorithm: iteration procedure of QSSA
Input:q
(
t
)
=
x
(
t
)
sin
⁡
⁡
(
2
π
t
/
(
L
N
-
L
+
1
)
) for image part
q
(
t
)
=
x
(
t
)
cos
⁡
⁡
(
2
π
t
/
(
L
N
-
L
+
1
)
) for real part
Output:Q—calculated image or real part
(1) forl
=
1 : L
(2)for  i
=
1 : (
L
-
l
)
N
+
1
(3)if  i
=
=
1
(4)    forn
=
1 : N
+
1
(5)  Q
(
i
)
=
Q
(
i
)
+
q
(
n
)
(6)    end for
(7)else Q
(
i
)
=
Q
(
i
-
1
)
-
q
(
i
-
1
)
+
q
(
N
+
i
)
(8)end if
(9) end for
(10)for  i
=
1 : (
L
-
l
)
N
+
1
(11)q
(
i
)
=
Q
(
i
)
(12) end for
(13) end for
(14) return


## 4. Simulation and Analysis

The aim of this section is to compare frequency estimations provided by the QSSA and WIFFT by using MATLAB. First, frequency estimations of a sine-wave signal are simulated under asynchronous sampling. Second, the influence of frequency fluctuation on frequency estimation is considered. Third, the sine-wave signal with white noise is simulated. Last, the joint influence of frequency fluctuation and white noise on frequency estimation is analyzed. Comparison summary is also presented at the end of this section. Notice that, in the following simulation,b
2
-
b
1
=
6 is adopted in QSSA.


### 4.1. With Asynchronous Sampling

Suppose a sine-wave signal as shown in (1) with A
=
22
0 V, φ
=
π
/
3, and f
=
50
.
2 Hz is simulated. The sampling rate is set as f
s = 3.2 kHz and f
s = 6.4 kHz to make comparisons. The QSSA with 513 samples and 641 samples and the WIFFT based on the 4-term maximum decay windows (MDW 4-III) with 512 and 1024 samples are adopted. The absolute error of frequency estimation and computational burden of the two adopted algorithms are shown in Table 2.Table 2
Comparisons of absolute errors of frequency estimation without frequency fluctuations.
QSSA
WIFFT (MDW 4-III)
N
=
128; L
=
4
N
=
128; L
=
5
512
1024
f
s
=
3.2 kHz
1
.
45
E
-
09
2
.
93
E
-
1
2
5
.
88
E
-
0
8
5
.
83
E
-
1
0
f
s
=
6.4 kHz
1
.
21
E
-
09
2
.
25
E
-
1
2
2
.
01
E
-
0
5
5
.
88
E
-
0
8
Additions
2048
3200
13824
30720
Multiplications
33029
41286
9216
20480When the sampling frequency isf
s = 3.2 kHz, 513 samples, 641 samples, and 1024 samples mean about 8, 10, and 16 periods of the sine-wave signal are used, respectively. When the sampling frequency is f
s = 6.4 kHz, 513 samples, 641 samples, and 1024 samples mean about 4, 5, and 8 periods of the sine-wave signal are used, respectively. As shown in Table 2, with the increasing of number of signal periods, that is, the decreasing of sampling frequency, the absolute errors of frequency estimation provided by the QSSA and WIFFT both decrease. Notice that the absolute errors of frequency estimation provided by WIFFT are more affected by the sampling frequency or the number of signal periods than that of the QSSA.As shown in Table2, the absolute errors of frequency estimation provided by the QSSA with 513 samples are much lower than that of the WIFFT with 512 samples. That is to say, by using almost the same numbers of samples, the QSSA outperforms the WIFFT in absolute errors of frequency estimation. However, Table 2 shows that the computational burden of the WIFFT is much less than that of the QSSA.


### 4.2. With Frequency Fluctuations

To investigate the influence of frequency fluctuations on frequency estimation, simulations are done with the frequency changing from 49.5 Hz to 50.5 Hz with a step of 0.1 Hz. The sampling frequency is set asf
s = 3.2 kHz. The QSSA with 513 samples and 641 samples and the WIFFT based on the 4-term maximum decay windows (MDW 4-III) with 512 and 1024 samples are adopted. Results of frequency estimations with frequency fluctuations are shown in Figure 3.Figure 3
Results of frequency estimations with frequency fluctuations when the sampling frequency is set asf
s = 3.2 kHz.Obviously, Figure3 shows that both QSSA and WIFFT can overcome the influence of frequency fluctuations on frequency estimation. Notice that, in Figure 3, the WIFFT (MDW 4-III) shows the lowest errors when the frequency is f
=
5
0 Hz with f
s = 3.2 kHz. However, it is well known that the f
=
5
0 Hz means synchronous sampling which is hard to be achieved in power system.Moreover, Figure3 also provides comparisons of two adopted algorithms when the asynchronous sampling procedure is considered; that is, f
≠
5
0 Hz with f
s = 3.2 kHz. The simulation results show that the absolute errors provided by the QSSA are less than that of the WIFFT (MDW 4-III) when the two adopted algorithms use almost the same numbers of samples.


### 4.3. With White Noise

Frequency estimation system may be influenced by external noise (electromagnetic interference), and the detection process cannot completely eliminate noise but is subject to the statistical law (normal or uniform distribution). We can regard noise in frequency estimation as white noise and study the impact of noise on the stability of frequency estimation by simulations with different signal-to-noise ratio (SNR).The frequency and sampling frequency are set asf
=
50
.
2 Hz and f
s = 3.2 kHz, respectively. The SNR changes from 10 dB to 100 dB with a step of 10 dB. For each SNR, 500 times of simulation are performed. The biases and variances of frequency estimations provided by the QSSA with 513 samples and 641 samples and the WIFFT based on the 4-term maximum decay windows (MDW 4-III) with 512 and 1024 samples are shown in Figures 4 and 5, respectively.Figure 4
Biases of frequency estimations with white noise when the SNR changes from 20 dB to 100 dB.Figure 5
Variances of frequency estimations with white noise when the SNR changes from 20 dB to 100 dB.As shown in Figure4, when the SNR ≤ 60 dB, the WIFFT (MDW 4-III) with 1024 samples provides the minimum biases among four adopted methods. However, when the SNR > 60 dB, the WIFFT (MDW 4-III) with 1024 samples and QSSA with 513 samples and 641 samples have almost the same biases of frequency estimation. This shows that the power of white noise has significant influence on the accuracy of frequency estimation. In Figure 5, the WIFFT (MDW 4-III) with 1024 samples provides the minimum variances among four adopted methods. Also, compared with the QSSA with 513 samples, the QSSA with 641 samples provides less variance of frequency estimation. This is because the increasing of samples with a fixed sampling frequency means the increasing of frequency resolution and thus can suppress the influence of white noise. Figure 5 shows that the variances provided by the QSSA with 513 samples are less than that of the WIFFT (MDW 4-III) with 512 samples. That is to say, by using almost the same numbers of samples, the QSSA outperforms the WIFFT in variances of frequency estimation with white noise.


### 4.4. With Frequency Fluctuations and White Noise

In this part, simulations are performed with SNR = 40 dB when the frequency changes from 49.5 Hz to 50.5 Hz. The sampling frequency is set asf
s = 3.2 kHz. The biases and variances of frequency estimations provided by the QSSA with 513 samples and 641 samples and the WIFFT based on the 4-term maximum decay windows (MDW 4-III) with 512 and 1024 samples are shown in Figures 6 and 7, respectively.Figure 6
Biases of frequency estimations with SNR = 40 dB when the frequency changes from 49.5 Hz to 50.5 Hz.Figure 7
Variances of frequency estimations with SNR = 40 dB when the frequency changes from 49.5 Hz to 50.5 Hz.In Figure6, the biases of frequency estimation obtained by the four adopted methods are less than 1
E
-
5 with SNR = 40 dB when the frequency changes from 49.5 Hz to 50.5 Hz. Figure 7 shows that the variances of frequency estimation obtained by the four adopted methods are less than 1
E
-
8 with SNR = 40 dB when the frequency changes from 49.5 Hz to 50.5 Hz. And the WIFFT (MDW 4-III) with 1024 samples has almost the minimum biases and variances among four adopted methods. Notice that the increasing of numbers of samples will lead to the significant increment of computational burden. Figures 6 and 7 show that the number of samples plays an important role in frequency estimation. As shown in Figure 7, when using almost the same numbers of samples, the variances of frequency estimation provided by the QSSA are less than that of the WIFFT (MDW 4-III).


## 4.1. With Asynchronous Sampling

Suppose a sine-wave signal as shown in (1) with A
=
22
0 V, φ
=
π
/
3, and f
=
50
.
2 Hz is simulated. The sampling rate is set as f
s = 3.2 kHz and f
s = 6.4 kHz to make comparisons. The QSSA with 513 samples and 641 samples and the WIFFT based on the 4-term maximum decay windows (MDW 4-III) with 512 and 1024 samples are adopted. The absolute error of frequency estimation and computational burden of the two adopted algorithms are shown in Table 2.Table 2
Comparisons of absolute errors of frequency estimation without frequency fluctuations.
QSSA
WIFFT (MDW 4-III)
N
=
128; L
=
4
N
=
128; L
=
5
512
1024
f
s
=
3.2 kHz
1
.
45
E
-
09
2
.
93
E
-
1
2
5
.
88
E
-
0
8
5
.
83
E
-
1
0
f
s
=
6.4 kHz
1
.
21
E
-
09
2
.
25
E
-
1
2
2
.
01
E
-
0
5
5
.
88
E
-
0
8
Additions
2048
3200
13824
30720
Multiplications
33029
41286
9216
20480When the sampling frequency isf
s = 3.2 kHz, 513 samples, 641 samples, and 1024 samples mean about 8, 10, and 16 periods of the sine-wave signal are used, respectively. When the sampling frequency is f
s = 6.4 kHz, 513 samples, 641 samples, and 1024 samples mean about 4, 5, and 8 periods of the sine-wave signal are used, respectively. As shown in Table 2, with the increasing of number of signal periods, that is, the decreasing of sampling frequency, the absolute errors of frequency estimation provided by the QSSA and WIFFT both decrease. Notice that the absolute errors of frequency estimation provided by WIFFT are more affected by the sampling frequency or the number of signal periods than that of the QSSA.As shown in Table2, the absolute errors of frequency estimation provided by the QSSA with 513 samples are much lower than that of the WIFFT with 512 samples. That is to say, by using almost the same numbers of samples, the QSSA outperforms the WIFFT in absolute errors of frequency estimation. However, Table 2 shows that the computational burden of the WIFFT is much less than that of the QSSA.


## 4.2. With Frequency Fluctuations

To investigate the influence of frequency fluctuations on frequency estimation, simulations are done with the frequency changing from 49.5 Hz to 50.5 Hz with a step of 0.1 Hz. The sampling frequency is set asf
s = 3.2 kHz. The QSSA with 513 samples and 641 samples and the WIFFT based on the 4-term maximum decay windows (MDW 4-III) with 512 and 1024 samples are adopted. Results of frequency estimations with frequency fluctuations are shown in Figure 3.Figure 3
Results of frequency estimations with frequency fluctuations when the sampling frequency is set asf
s = 3.2 kHz.Obviously, Figure3 shows that both QSSA and WIFFT can overcome the influence of frequency fluctuations on frequency estimation. Notice that, in Figure 3, the WIFFT (MDW 4-III) shows the lowest errors when the frequency is f
=
5
0 Hz with f
s = 3.2 kHz. However, it is well known that the f
=
5
0 Hz means synchronous sampling which is hard to be achieved in power system.Moreover, Figure3 also provides comparisons of two adopted algorithms when the asynchronous sampling procedure is considered; that is, f
≠
5
0 Hz with f
s = 3.2 kHz. The simulation results show that the absolute errors provided by the QSSA are less than that of the WIFFT (MDW 4-III) when the two adopted algorithms use almost the same numbers of samples.


## 4.3. With White Noise

Frequency estimation system may be influenced by external noise (electromagnetic interference), and the detection process cannot completely eliminate noise but is subject to the statistical law (normal or uniform distribution). We can regard noise in frequency estimation as white noise and study the impact of noise on the stability of frequency estimation by simulations with different signal-to-noise ratio (SNR).The frequency and sampling frequency are set asf
=
50
.
2 Hz and f
s = 3.2 kHz, respectively. The SNR changes from 10 dB to 100 dB with a step of 10 dB. For each SNR, 500 times of simulation are performed. The biases and variances of frequency estimations provided by the QSSA with 513 samples and 641 samples and the WIFFT based on the 4-term maximum decay windows (MDW 4-III) with 512 and 1024 samples are shown in Figures 4 and 5, respectively.Figure 4
Biases of frequency estimations with white noise when the SNR changes from 20 dB to 100 dB.Figure 5
Variances of frequency estimations with white noise when the SNR changes from 20 dB to 100 dB.As shown in Figure4, when the SNR ≤ 60 dB, the WIFFT (MDW 4-III) with 1024 samples provides the minimum biases among four adopted methods. However, when the SNR > 60 dB, the WIFFT (MDW 4-III) with 1024 samples and QSSA with 513 samples and 641 samples have almost the same biases of frequency estimation. This shows that the power of white noise has significant influence on the accuracy of frequency estimation. In Figure 5, the WIFFT (MDW 4-III) with 1024 samples provides the minimum variances among four adopted methods. Also, compared with the QSSA with 513 samples, the QSSA with 641 samples provides less variance of frequency estimation. This is because the increasing of samples with a fixed sampling frequency means the increasing of frequency resolution and thus can suppress the influence of white noise. Figure 5 shows that the variances provided by the QSSA with 513 samples are less than that of the WIFFT (MDW 4-III) with 512 samples. That is to say, by using almost the same numbers of samples, the QSSA outperforms the WIFFT in variances of frequency estimation with white noise.


## 4.4. With Frequency Fluctuations and White Noise

In this part, simulations are performed with SNR = 40 dB when the frequency changes from 49.5 Hz to 50.5 Hz. The sampling frequency is set asf
s = 3.2 kHz. The biases and variances of frequency estimations provided by the QSSA with 513 samples and 641 samples and the WIFFT based on the 4-term maximum decay windows (MDW 4-III) with 512 and 1024 samples are shown in Figures 6 and 7, respectively.Figure 6
Biases of frequency estimations with SNR = 40 dB when the frequency changes from 49.5 Hz to 50.5 Hz.Figure 7
Variances of frequency estimations with SNR = 40 dB when the frequency changes from 49.5 Hz to 50.5 Hz.In Figure6, the biases of frequency estimation obtained by the four adopted methods are less than 1
E
-
5 with SNR = 40 dB when the frequency changes from 49.5 Hz to 50.5 Hz. Figure 7 shows that the variances of frequency estimation obtained by the four adopted methods are less than 1
E
-
8 with SNR = 40 dB when the frequency changes from 49.5 Hz to 50.5 Hz. And the WIFFT (MDW 4-III) with 1024 samples has almost the minimum biases and variances among four adopted methods. Notice that the increasing of numbers of samples will lead to the significant increment of computational burden. Figures 6 and 7 show that the number of samples plays an important role in frequency estimation. As shown in Figure 7, when using almost the same numbers of samples, the variances of frequency estimation provided by the QSSA are less than that of the WIFFT (MDW 4-III).


## 5. Conclusion

This paper compares the performances of the QSSA and the WIFFT for frequency estimation. It has been shown that the absolute errors provided by the QSSA are less than that of the WIFFT (MDW 4-III) when the two adopted algorithms use about 512 samples when the fundamental frequency has fluctuations. Moreover, the variances of frequency estimation provided by the QSSA are less than that of the WIFFT (MDW 4-III) when the number of acquired samples is about 512. That is to say, the wideband noise sensitivity of QSSA is lower than that of WIFFT. However, we should also notice that the computational burden required by the QSSA is higher than that of the WIFFT. Thus, for frequency estimation with high accuracy and less required samples, the QSSA could be preferred. For fast frequency estimation with large numbers of samples, WIFFT would be a good choice.



---

*Source: 745830-2014-08-13.xml*