A DC offset removal algorithm based on AR model

ABSTRACT

This work proposes an auto-regressive (AR) model-based algorithm that can eliminate the adverse influence of the exponentially decaying DC offset in the phasor estimation process using discrete Fourier transform (DFT). The proposed algorithm comprises three stages: DFT, AR model implementation, and compensation. First, DFT can eliminate all harmonic components and estimate the phasor of the fundamental frequency component. Second, the AR model is used to calculate the error caused by the exponentially decaying DC offset in the phasor value estimated by DFT. A second-order AR model that uses three successively estimated phasor values is utilized. Finally, the phasor of the fundamental frequency component can be accurately estimated by compensating the error caused by the exponentially decaying DC offset. The performance of the proposed algorithm is evaluated for a-phase to ground fault on a 154 kV 25 km overhead transmission line. Electromagnetic Transients Program is used to generate fault signals. The evaluation result suggests that the proposed algorithm can effectively suppress the adverse influence of the DC offset.

KEYWORDS:

• Dc offset
• phasor estimation
• auto regressive model
• Fourier transformation

1. Introduction

Discrete Fourier transform (DFT) is generally used to calculate the phasor of the fundamental frequency component in digital protective relays. However, the phasor value estimated by DFT contains an error caused by exponentially decaying DC offsets. For this reason, the DC offset component in a fault relaying signal adversely affects the accuracy and the operating speed of a digital protection relay. Therefore, the decaying DC component should be carefully considered in calculating the phasor of the fundamental frequency component of a relaying signal.

Several numerical techniques have been proposed to reduce or remove the adverse effect of the decaying DC offset component. In [1], a DFT algorithm with a preconditioning digital mimic filter was proposed to suppress the decaying DC offset in a current waveform. The digital mimic filter can completely remove the decaying DC offset only when the actual time constant of a DC offset in a fault relaying current is the same as that of the mimic filter. A modified Fourier algorithm [2] was proposed to remove the influence of the decaying DC offset on a fault relaying signal. This method can reduce the influence of but not completely remove the decaying DC offset. An AR model-based method was proposed in [3] to remove the effect of the DC offset on a relaying signal. However, this method can remove the decaying DC offset only when a relaying signal consists of the fundamental frequency component and a decaying DC offset component. A DFT-based modified phasor estimation method was proposed in [4]. The error due to the exponentially decaying DC offset on the DFT is calculated and eliminated using the outputs of the even- and odd-sample-set DFTs. However, this method is disadvantageous because the high-frequency noise components are amplified by the subtraction process of calculating the error caused by the DC offset.

We propose an AR model-based algorithm to eliminate the adverse influence of the exponentially decaying DC offset. The proposed AR model is implemented with phasor values estimated by the DFT to improve the disadvantage of the existing AR model in [3], which is vulnerable to harmonics. The error due to the DC offset in the phasor value estimated by the DFT is calculated and eliminated using three successively estimated phasor values to accurately estimate the phasor of the fundamental frequency component. The performance of the proposed algorithm is evaluated for a-phase to ground (a–g) faults on a 154 kV 25 km overhead transmission line. Electromagnetic Transients Program (EMTP) was used to generate the fault relaying signals. The results of the proposed algorithm are compared to those of the conventional DFT and Modified DFT in [2] to demonstrate its performance. Evaluation results indicate that the proposed algorithm can effectively eliminate the adverse influence of the DC offset.

2. Propsed algorithm

2.1. Auto-regressive model

In an AR model, series data can be described by past values and expressed as follows:

(1) $y\left(k\right)=a_{1}y\left(k-1\right)+a_{2}y\left(k-2\right)$(1)

This model is a second-order AR model because the value at time (k) is predicted from the values at times (k-1) and (k-2). A p-th-order AR model is generally a multiple linear regression in which the value of the series at any time (k) is a weighted sum of the values until the time (k-p) passes; it can be expressed as follows:

(2) $y\left(k\right)=a_{1}y\left(k-1\right)+a_{2}y\left(k-2\right)+\dots a_{p}y\left(k-p\right)$(2)

2.2. Error caused by exponentially decaying DC offset in phasor value estimated by DFT

When a fault occurs in a power system, the current signal can be generally considered the combination of an exponentially decaying DC offset and sinusoidal components. If sinusoidal components with frequencies higher than the (N/2–1) order are assumed to be eliminated by an anti-aliasing low-pass filter, the current signal can be expressed as follows:

(3) $i\left[k\right]=A_0{e}^{-\frac{k\mathrm{\Delta }t}{\tau }}+\sum _{n=1}^{N/2-1}{A}_n\mathrm{c}\mathrm{o}\mathrm{s}({\theta }_n+\frac{2\pi n}{N}\mathrm{k})$(3)

where $A_0$ and $\tau$ are the magnitude and the time constant of an exponential component, respectively; $A_n$and ${\theta }_n$ are the amplitude and the phase angle of the n-th harmonic component, respectively; $\mathrm{\Delta }t$ is the sampling interval; and $N$ is the number of samples per cycle.

The phasor of the fundamental frequency component is calculated by the DFT as follows:

(4) $I_I\left[k\right]=\frac{2}{N}\sum _{i=k-N+1}^{k}i\left[j\right]\cdot sin\left(\frac{2\pi }{N}j\right)$(4)

$I_R\left[k\right]=\frac{2}{N}\sum _{i=k-N+1}^{k}i\left[j\right]\cdot sin\left(\frac{2\pi }{N}j\right)$

where $I_I$ and $I_R$ are the imaginary and real parts of the fundamental frequency phasor, respectively.

The estimated phasor value by Equation (4) is as follows:

(5) $I_I\left[k\right]=A_{1}sin\left({\theta }_1+\frac{2\pi }{N}k\right)+\frac{2}{N}{A}_0\frac{1-E^N}{1-2Ecos\left(\frac{2\pi }{N}\right)+E^2}Esin\left(\frac{2\pi }{N}\right)E^k$(5)

(6) $I_R\left[k\right]=A_{1}cos\left({\theta }_1+\frac{2\pi }{N}k\right)+\frac{2}{N}{A}_0\frac{1-E^N}{1-2Ecos\left(\frac{2\pi }{N}\right)+E^2}\left\{1-Ecos\left(\frac{2\pi }{N}\right)\right\}{E}^k$(6)

where $E=e^{-\mathrm{\Delta }t/\tau }$

As shown in Equations (5) and (6), the estimated phasor value contains an error that is caused by the exponentially decaying DC offset in addition to the fundamental frequency component. This equation can be rewritten as follows:

(7) $I_I\left[k\right]=A_1\cdot sin\left({\theta }_1+\frac{2\pi }{N}k\right)+D_I\cdot e^{-\frac{\mathrm{\Delta }t}{\tau }}=I_I^{1st}\left[k\right]+I_I^{dc}\left[k\right]$(7)

(8) $I_R\left[k\right]=A_1\cdot cos\left({\theta }_1+\frac{2\pi }{N}k\right)+D_R\cdot e^{-\frac{\mathrm{\Delta }t}{\tau }}=I_R^{1st}\left[k\right]+I_R^{dc}\left[k\right]$(8)

where $D_I=\frac{2}{N}{A}_0\frac{1-E^N}{1-2Ecos\left(\frac{2\pi }{N}\right)+E^2}Esin\left(\frac{2\pi }{N}\right)$

$D_R=\frac{2}{N}{A}_0\frac{1-E^N}{1-2Ecos\left(\frac{2\pi }{N}\right)+E^2}\{1-Ecos\left(\frac{2\pi }{N}\right).$

Therefore, to obtain the accurate fundamental frequency component, the error caused by the DC offset included in the estimated phasor value must be eliminated.

2.3. AR model for fundamental frequency component

The second derivative of the imaginary part of the fundamental frequency component can be calculated using a central difference approximation formula of Equation (9).

(9) $I_I^{1st}{ }^{\prime \prime }\left[k\right]=\frac{\left(I_I^{1\mathrm{s}\mathrm{t}}\left[k+1\right]-I_I^{1\mathrm{s}\mathrm{t}}\left[k\right]\right)-\left(I_I^{1\mathrm{s}\mathrm{t}}\left[k\right]-I_I^{1\mathrm{s}\mathrm{t}}\left[k-1\right]\right)}{\mathrm{\Delta }{t}^2}$(9)

The second derivative term can be expressed as follows:

(10) $I_I^{1st{ }^{\prime \prime }}\left[k\right]=-{\left\{\frac{2sin\left(\frac{\pi }{N}\right)}{\mathrm{\Delta }t}\right\}}^2{I}_I^{1st}\left[k\right]$(10)

Substituting Equation (10) into the left side of Equation (9) obtains the following:

(11) $-{\left\{\frac{2sin\left(\frac{\pi }{N}\right)}{\mathrm{\Delta }t}\right\}}^2{I}_I^{1st}\left[k\right]=\frac{\left(I_I^{1st}\left[k+1\right]-I_I^{1st}\left[k\right]\right)-\left(I_I^{1st}\left[k\right]-I_I^{1st}\left[k-1\right]\right)}{\mathrm{\Delta }{t}^2}$(11)

From Equation (11), the imaginary part of the fundamental frequency component can be represented by the following second AR model:

(12) $I_I^{1st}\left[k+1\right]=\left[2-{\left\{2sin\left(\frac{\pi }{N}\right)\right\}}^2\right]I_I^{1st}\left[k\right]-I_I^{1st}\left[k-1\right]$(12)

Similarly, the real part of the fundamental frequency component can be represented by the following second AR model;

(13) $I_R^{1st}\left[k+1\right]=\left[2-{\left\{2sin\left(\frac{\pi }{N}\right)\right\}}^2\right]I_R^{1st}\left[k\right]-I_R^{1st}\left[k-1\right]$(13)

2.4. Calculation of error caused by exponentially decaying DC offset in estimated phasor value

The imaginary part of the estimated phasor value can be expressed with the error caused by the exponentially decaying DC offset as shown in Equation (14).

(14) $I_I\left[k-1\right]=I_I^{1st}\left[k-1\right]+I_I^{dc}\left[k-1\right]$(14)

$I_I\left[k\right]=I_I^{1st}\left[k\right]+I_I^{dc}\left[k\right]$

$I_I\left[k+1\right]=I_I^{st}\left[k+1\right]+I_I^{dc}\left[k+1\right]$

The pure fundamental frequency component from the estimated phasor value should be eliminated to estimate the error caused by the exponentially decaying DC offset in the estimated phasor value.

From Equations (12) and (14), the following equation can be obtained:

(15) $I_I\left[k+1\right]-\left[2-{\left\{2sin\left(\frac{\pi }{N}\right)\right\}}^2\right]I_I\left[k\right]+I_I\left[k-1\right]$(15)

$=I_I^{dc}\left[k+1\right]-\left[2-{\left\{2sin\left(\frac{\pi }{N}\right)\right\}}^2\right]I_I^{dc}\left[k\right]+I_I^{dc}\left[k-1\right]$

$\begin{array}{rl}& =D_I\cdot e^{-\frac{\mathrm{\Delta }t}{\tau }\left(k+1\right)}-\left[2+{\left\{2sin\left(\frac{\pi }{N}\right)\right\}}^2\right]D_I\cdot e^{-\frac{\mathrm{\Delta }t}{\tau }\left(k\right)}\\ & +D_I\cdot e^{-\frac{\mathrm{\Delta }t}{\tau }\left(k-1\right)}\end{array}$

Equation (15) can be rewritten with the following equation consisting only of the error caused by the exponentially decaying DC offset.

(16) $I_I\left[k+1\right]-\left[2-{\left\{2sin\left(\frac{\pi }{N}\right)\right\}}^2\right]I_I\left[k\right]+I_I\left[k-1\right]$(16)

$=D_I\cdot e^{-\frac{\mathrm{\Delta }t}{\tau }k}\left[e^{-\frac{\mathrm{\Delta }t}{\tau }}-2+{\left\{2sin\left(\frac{\pi }{N}\right)\right\}}^2+e^{\frac{\mathrm{\Delta }t}{\tau }}\right]=K_I\left[k\right]$

From Equation (16), the error due to the exponentially decaying DC offset in the estimated phasor value can be expressed as follows:

(17) $D_I\cdot e^{-\frac{\mathrm{\Delta }t}{\tau }k}=\frac{1}{e^{-\frac{\mathrm{\Delta }t}{\tau }}-2+{\left\{2sin\left(\frac{\pi }{N}\right)\right\}}^2+e^{\frac{\mathrm{\Delta }t}{\tau }}}{K}_I\left[k\right]=C_I\cdot K_I\left[k\right]$(17)

From Equation (17), the error due to the exponentially decaying DC offset in the calculated phasor value can be calculated by $C_I$ and $K_I$, where $K_I$ can be calculated by three imaginary parts of DFT, as shown in Equation (16), and $C_I$ can be calculated by the constants ($\mathrm{\Delta }t$, $N$) and unknown variable ($\tau$). The time constant of the DC offset is required to calculate $C_I$ given that the time constant ($\tau$) is an unknown variable. From Equation (16), the following equation can be obtained:

(18) $\frac{K_I\left[k\right]}{K_I\left[k-1\right]}=\frac{D_I\cdot e^{-\frac{\mathrm{\Delta }t}{\tau }k}\left[e^{-\frac{\mathrm{\Delta }t}{\tau }}-2+{\left\{2sin\left(\frac{\pi }{N}\right)\right\}}^2+e^{\frac{\mathrm{\Delta }t}{\tau }}\right]}{D_I\cdot e^{-\frac{\mathrm{\Delta }t}{\tau }\left(k-1\right)}\left[e^{-\frac{\mathrm{\Delta }t}{\tau }}-2+{\left\{2sin\left(\frac{\pi }{N}\right)\right\}}^2+e^{\frac{\mathrm{\Delta }t}{\tau }}\right]}=e^{-\frac{\mathrm{\Delta }t}{\tau }}$(18)

The error due to the exponentially decaying DC offset in the imaginary part of the estimated phasor value can be calculated in Equation (17) using the calculated time constant.

Similarly, the error due to the exponentially decaying DC offset in the real part of the estimated phasor value can be calculated by Equation (19).

(19) $D_R\cdot e^{-\frac{\mathrm{\Delta }t}{\tau }k}=\frac{1}{e^{-\frac{\mathrm{\Delta }t}{\tau }}-2+{\left\{2sin\left(\frac{\pi }{N}\right)\right\}}^2+e^{\frac{\mathrm{\Delta }t}{\tau }}}{K}_R\left[k\right]=C_R\cdot K_R\left[k\right]$(19)

$\mathrm{w}\mathrm{h}\mathrm{e}\mathrm{r}\mathrm{e}{K}_R\left[k\right]=I_R\left[k+1\right]-\left[2-{\left\{2sin\left(\frac{\pi }{N}\right)\right\}}^2\right]I_R\left[k\right]+I_R\left[k-1\right]$

$K_R$ can be calculated by three real parts of the estimated phasor value, and $C_R$ can be calculated from the constants ($\mathrm{\Delta }t$, $N$) and unknown variable (τ). Then, $C_R$ can be calculated because the time constant of the DC offset is already obtained by Equation (18). Therefore, the error due to the exponentially decaying DC offset in the real part of the estimated phasor value can be calculated by Equation (19).

2.5. Compensation of error caused by exponentially decaying DC offset in estimated phasor value

The phasor of the fundamental frequency component can be accurately estimated by subtracting the error caused by the exponentially decaying DC offset calculated from Equations (17) and (19) from the estimated phasor value.

(20) $I_I^{1st}\left[k\right]=I_I\left[k\right]-D_I\cdot e^{-\frac{\mathrm{\Delta }t}{\tau }k}$(20)

(21) $I_R^{1st}\left[k\right]=I_R\left[k\right]-D_R\cdot e^{-\frac{\mathrm{\Delta }t}{\tau }k}$(21)

Figure 1 shows the flowchart of the proposed algorithm. A low-pass filter is used for anti-aliasing, and the DFT can eliminate all harmonic components and estimate the phasor of the fundamental frequency component. The time constant of the DC offset is calculated using K_I, and the error due to the exponentially decaying DC offset in the estimated phasor value is calculated using the AR model that reflects the calculated time constant. Finally, the phasor of the fundamental frequency component can be accurately estimated by compensating the error caused by the exponentially decaying DC offset.

Figure 1. Flowchart of proposed algorithm.

3. Performance evaluation

3.1. System configuration

In this section, the performance evaluation of the proposed algorithm is described. The current signals are generated by EMTP. The model system for the simulations is a 154 kV 25 km overhead transmission with sources at both ends, as shown in Figure 2. The transmission line parameters used in the simulation are presented in Table 1.

Table 1. Overhead transmission line parameters.

Sequence	Parameter	Value	Units
Positive & Negative	R1, R2	0.0419	Ω/km
	L1, L2	0.8799	mH/km
	C1, C2	0.0128	μF/km
Zero	R0	0.2294	Ω/km
	L0	2.6666	mH/km
	C0	0.0043	μF/km

Figure 2. Single-line diagram of model system.

The time constant and magnitude of the decaying DC offset vary depending on the system configuration and fault conditions, such as the fault location, resistance, and inception angle. The time constant of the DC offset is generally determined by the fault resistance and the impedance ratio (X/R) of the system. As shown in Table 2, different source impedances are used for Cases 1 and 2 to consider the effect of the source impedance ratio on the time constant. A fault inception angle of 0° is considered for the maximum magnitude of the DC component.

Table 2. Equivalent source parameters.

`	Case 1		Case 2
Sequence	S [GVA]	X/R	S [GVA]	X/R
ZS0	5.0	10.7	5.0	15.0
ZS1	11.5	5.1	11.5	7.5
Sequence	R [GVA]	X/R	R [GVA]	X/R
ZR0	2.1	18.1	2.1	18.1
ZR1	6.4	4.9	6.4	4.9

3.2. Case studies

Figures 3-5 show the estimated RMS values of the fault currents. At this time, the effect of the fault distance is considered. The time constant of the DC offset is changed by the fault distance from the relaying point. For comparison, the results obtained with the conventional DFT and Modified DFT are also shown in each figure.

Figure 3. Test results according to change in fault distance (fault distance: 10 [%], source impedance ratio: case 1).

Figure 4. Test results according to change in fault distance (fault distance: 50 [%], source impedance ratio: case 1).

Figure 5. Test results according to change in fault distance (fault distance: 90 [%], source impedance ratio: case 1).

Figure 3 shows the test results for a solid a-phase to ground fault at 10[%] of the total line length from the relaying point. For a solid fault, the fault current flowing through a power line contains the exponentially decaying DC offset which depends on the local source impedance, the line impedance, and the impedance of the ground. Compared to the other algorithms, the result of the conventional DFT shows not only the largest transient overshoot but also the largest oscillatory response. Figure 3(b) shows only the proposed algorithm and the Modified DFT for accurate comparison. As shown in Figure 3(b), the proposed algorithm can more accurately remove the adverse influence of the DC offset than the Modified DFT.

Figures 4 and 5 show the test results for a solid a-phase to ground fault at 50[%] and 90[%] of the total line length from the relaying point, respectively. As shown in Figures 4(a) and 5(a), the result of the conventional DFT shows not only the largest transient overshoot but also the largest oscillatory response. Figures 4(b) and 5(b) shows only the proposed algorithm and the Modified DFT for accurate comparison. As shown in Figures 4(b) and 5(b), the proposed algorithm can more accurately remove the adverse influence of the DC offset than the Modified DFT.

Figures 6 and 7 show the test results for a solid a-phase to ground fault at 0.1 km from the relaying point. At this time, the effect of the source impedance is considered. The time constant of the DC offset is changed by the X/R of the local source impedance. As shown in Figures 6(a) and 7(a), the result of the conventional DFT shows not only the largest transient overshoot but also the largest oscillatory response. Figures 6(b) and 7(b) shows only the proposed algorithm and the Modified DFT. As shown in Figures 6(b) and 7(b), the proposed algorithm can more accurately remove the adverse influence of the DC offset than the Modified DFT.

Figure 6. Test results according to change of the source impedance ratio (fault distance: 0.1 [km], source impedance ratio: case 1).

Figure 7. Test results according to change of the source impedance ratio (fault distance: 0.1 [km], source impedance ratio: case 2).

4. Conclusions

An algorithm based on a proposed AR model is proposed to eliminate the adverse influence of exponentially decaying DC offset. The error due to the DC offset in the phasor value obtained from DFT is calculated and eliminated with the AR model. The performance results show that the proposed algorithm is not affected by the system and the fault conditions. The output of the proposed algorithm is stable and does not have oscillatory response, unlike that of the other algorithms. Thus, the proposed algorithm can effectively suppress the adverse influence of the DC offset in a relaying signal.

Disclosure statement

No potential conflict of interest was reported by the authors.

Additional information

Funding

This research was supported by the Korea Electric Power Corporation [R17XA05-2].

Notes on contributors

Woo-Joong Kim

Woo-Joong Kim: He received his B.S. and M.S. degrees from Myongji University, Korea, in 2012 and 2014, respectively. He is now studying for his Ph.D. degree at Myongji University. He has been with Next-generation Power Technology Center, Korea, since 2012. His main research interest is power system protection.

Sang-Hee Kang

Sang-Hee Kang: He is a professor at Myongji University, Korea. He received the B.S., M.S. and Ph.D. degrees from Seoul National University, Korea in 1985, 1987 and 1993, respectively. He was a visiting fellow and a visiting scholar at the University of Bath, UK in 1991 and 1999. His research interest is to develop digital protection systems for power systems using digital signal processing techniques.