Power Data

Location:

College Station, TX

Posted:

November 12, 2012

Contact this candidate

Resume:

Elimination of DC Offset in Accurate Phasor

Estimation Using Recursive Wavelet Transform

J. Ren, Student Member, IEEE, M. Kezunovic, Fellow, IEEE

This method requires both voltage and current inputs. As a

result, it is not applicable to the current-based protection

Abstract DC offset has significant effect on extracting

schemes.

fundamental frequency components. It directly impacts the

Recursive wavelet approach has been introduced in

accuracy of the fundamental frequency component based

protective relaying algorithms. This paper proposes a novel protective relaying for a long time [11]-[13]. The improved

method for estimating the fundamental phasor while eliminating

model with single-direction recursive equations is more

the dc offset using improved recursive wavelet transform. The

suitable for the application to real-time signal processing [14].

proposed approach converges to the actual value in less than one

The band energy of any center frequency can be extracted

cycle, and the computation burden is fairly low because of the

through improved recursive wavelet transform (IRWT) with

recursive formula. Studies indicate that to achieve a certain level

moderately low computation burden. Recursive wavelet

of accuracy, the higher sampling frequency one uses, the shorter

data window the computation needs, and vice versa. Comparing features band-pass filter achieves good performance in

with conventional DFT filter, simulation results demonstrate that

suppressing the sub-harmonic components. A novel filtering

the proposed phasor estimation method achieves better

approach using IRWT is proposed for estimating the

performance.

fundamental frequency phasor while eliminating the effect of

dc offset.

Index Terms protective relaying, fundamental frequency

This paper first introduces the recursive wavelet transform

component, phasor estimation, dc offset, improved recursive

and its characteristics in the time and frequency domain. The

wavelet transform, data window, conventional DFT filter

phasor estimation and DC offset removal algorithm are

I. INTRODUCTION presented next. The relationship between the convergence and

I the sampling rate is studied as well. It indicates that to achieve

N protective relaying application, Discrete Fourier

a certain level of accuracy, the higher sampling rate one uses,

Transform (DFT) is widely used as a filtering algorithm for

the shorter data window it needs, and vice versa. Simulation

extracting fundamental phasors [1], [2]. Conventional DFT

results demonstrate the effectiveness of the proposed phasor

algorithm achieves excellent performance when the signals

estimation method. It converges to the actual value in less than

contain only fundamental frequency and integer harmonic

one cycle, and the computation burden is fairly low because of

frequency components. Since in most cases the currents

the recursive formula, thus it can satisfy the time response

contain DC offsets this may introduce fairly large errors in the

requirement of the high speed relaying schemes.

estimation of fundamental frequency phasor [3], [4].

Many techniques have been proposed to eliminate the DC

offset in waveforms. A digital mimic filter based method was II. FILTERING ALGORITHM

proposed in [5]. This filter features high-pass frequency

A. Characteristics of Conventional DFT Filters

response which results in bringing high frequency noise to the

outcome. It performs well when its time constant matches the In terms of the length of the data window used for the

time constant of the exponentially decaying component. filtering calculation, conventional Discrete Fourier Transform

Theoretically, the decaying component can be completely (DFT) can be classified into two categories: Full Cycle DFT

removed from the original waveform once its parameters can (FCDFT) and Half Cycle DFT (HCDFT). Frequency

be obtained. Based on this idea, [6], [7] utilize additional responses of FCDFT and HCDFT shown in Fig. 1 and Fig. 2

samples to calculate the parameters of the decaying respectively indicate the performance in suppressing integer

component. Reference [8] uses the simultaneous equations frequency harmonics. The sub-harmonic components and DC

derived from the harmonics. A new Fourier algorithm and offset can not be easily eliminated with conventional DFT

three simplified algorithms based on Taylor expansion were filters. This can be seen from Fig. 3, which presents the

proposed to eliminate the decaying component in [9]. In [10], frequency spectrum of a set of exponentially decaying signals

author estimates the parameters of the decaying component by with a broad range of time constants (0.5 to 5 cycles). In

using the phase angle difference between voltage and current. power system, DC offset widely exists in voltage and current

signals when various disturbances occur, such as fault or

oscillation, and it may cause the calculated amplitude deviate

J. Ren and M. Kezunovic are with the Department of Electrical and from the real value 15% in worst case [9].

Computer Engineering, Texas A&M University, College Station, TX 77843-

3128, USA (e-mails: *.*.***@***.****.***, *******@***.****.***).

imaginary part in time domain. The real part, imaginary part

and magnitude of the frequency characteristics are given in

Fig. 5, where a = 1/ f0. As we can see in Fig. 5, IRW features a

band-pass filter with the center frequency f0.

Fig. 1. Frequency response of the FCDFT

Fig. 4. Waveforms of IRW in time domain

Fig. 2. Frequency response of the FCDFT

Fig. 5. Waveforms of IRW in frequency domain

Assume the input signal x(k) and the sampling period T,

the formula of improved recursive wavelet transform is given

as follows:

f { x(k 1) 2 x(k 2)

W xIR ) ( f, k ) T

3 x(k 3) 4 x(k 4) 5 x(k 5)}

Fig. 3. Amplitude-frequency characteristic of exponentially decaying

signals

1W xIR ) ( f, k 1) 2W xIR ) ( f, k 2) (2)

(k (k

B. Improved Recursive Wavelet Transform

3W xIR ) ( f, k 3) 4W xIR ) ( f, k 4)

(k (k

The mother wavelet is given as:

5W xIR ) ( f, k 5) 6W xIR ) ( f, k 6)

3t 3 4t 4 5t 5 ( j 0 ) t (k (k

(t ) ( u ( t )

where,

3 6 15

e f T j

A set of wavelet functions can be created by dilating and 0

shifting the mother wavelet, as given below:

1 [( f T ) 3 / 3 ( f T ) 4 / 6 ( f T ) 5 / 15]

t b

a,b (t ) a 1 / 2 2 2 [2( f T ) 3 / 3 5( f T ) 4 / 3 26( f T ) 5 / 15]

3 3 [ 6( f T ) 3 / 3 22( f T ) 5 / 5]

The wavelet coefficient for a given signal x(t) can be express

as below:

4 4 2( f T ) 3 / 3 5( f T ) 4 / 3 26( f T ) 5 / 15

t b

W x (t ) (a, b) a 1 / 2 x(t ) dt 5 5 [( f T ) 3 / 3 ( f T ) 4 / 6 ( f T ) 5 / 15]

(1)

1 6, 2 15 2, 3 20 3

Assign 2 / 3, 0 2 thus (t ) is admissible.

4 15 4, 5 6 5, 6 6

Besides, we have a 1 / f, that is the scale coefficient a is

reciprocal to the frequency f used in the analysis.

Improved Recursive Wavelet (IRW) exhibits good time-

frequency characteristics. Fig. 4 shows its real part and

g D e T ( k 1) / h(k 1)

III. PHASOR ESTIMATION SCHEME (6)

where

A. Estimating Phasor Using IRWT

h(k 1) W yIRk 1) (a, k 1) I (a, k ) y (k 1)

Consider a sinusoidal signal expressed in complex form: (

x(t ) Am e j ( t ) t 0 From (5) and (6), we obtain,

T h( k )

where Am is the amplitude, is the phase angle., D

g e T k /

log h(k ) / h(k 1)

Apply IRWT to signal x(t) using (1). As derived in

Appendix, we obtain,

It should be noticed that g is known after is calculated.

W x (t ) (a, b) Am e j ( b ) I (a, b) (3) Then, the decaying component can be removed completely by

the following formulae:

Assume the sampling frequency f equals N times the

W yIRk ) (a, k ) I (a, k ) D e T k /

fundamental frequency f0. That is f = N f0, and the sample

e j 1

j ( 2 k / N ) (

Am e

period is T = 1/ f0 N, t = T k. The discrete expression of

I ( a, k )

equation (3) is:

W x ( k ) (a, k ) Am e j ( k T ) I (a, k ) Therefore, we have

W yIRk ) (a, k ) I (a, k ) D e T k /

To extract the fundamental component, simply let a = 1/f0. (

From (2) we have the coefficient W x ( k ) ( a, k ) . Select proper I ( a, k ) (7)

sampling rate so that the error resulting from the discrete data

1 2 k / N

computation is within the limit of tolerance, we

have W x ( k ) ( a, k ) W

(a, k ) . That is: C. Analysis of the Impact of Data Window Length

x(k )

j ( 2 k / N )

W IR Theoretically, the amplitude and phase angle of the input

Am e ( a, k ) / I ( a, k )

x(k )

signal can be accurately estimated using (7) by two samples,

W xIR ) (a, k ) / I (a, k ) e j which means the length of the data window would be 2 T. In

reality two factors should be considered when selecting the

Thus, data window length. One is that formula (4) and (7) are

Am W xIR ) (a, k ) / I (a, k ), 2 k / N (4) derived based on the assumption that the error resulting from

the discrete computation is negligible. Another is the error

In above formula, W x ( k ) ( a, k ) is the wavelet transform introduced by the recursive calculation.

In practical application, the length of the data window

coefficient calculated by the recursive equation (2). I(a,k) is

should be selected in terms of the sampling frequency. Fig. 6

constant given in Appendix which can be calculated in

presents the convergence characteristics of the proposed

advance.

filtering algorithm vs. the sampling frequency. In Fig. 6, the

B. Eliminating Exponentially Decaying Component window length is one cycle of the fundamental frequency; f is

Consider a signal consisting of exponentially decaying the sampling frequency while f0 is the fundamental frequency.

component,

y (t ) D e t / x(t ) t 0

where D is the amplitude, is the time constant.

Apply IRWT to signal y(t) using (1). The IRWT coefficient

of signal y(t) derived in Appendix is given as

W y (t ) (a, b) I (a, b) D e b / W x (t ) (a, b)

As we can see from above equation, D and can be calculated

by subtracting I(a,b) x(t) from Wy(t)(a,b). Introducing the

wavelet transform coefficient of signal y(t) by applying

formula (2) with a = 1/f0, in discrete form, we have

W yIRk ) (a, k ) I (a, k ) y (k )

(

[ I (a, k ) I (a, k )] D e T k /

Designate Fig. 6. Convergence characteristics vs. sampling frequency

g (, k ) I (a, k ) I (a, k )

To obtain certain accuracy of 1% of the real value, the

h(k ) W yIRk ) (a, k ) I (a, k ) y (k ), relationship between sampling frequency f and data window

(

length Ts is studied considering various phase angles of input

We obtain

signals, as shown in Fig. 7, where is the phase angle of the

g D e T k / h(k ) (5) input signal. We can conclude that higher sampling rate

Introduce additional sample y(k+1)

shortens the data window and expedites the convergence

process for the phasor estimation.

Fig. 9 Phase angle of the estimated phasor

Table I summarizes the test results, in which the estimated

values are obtained at half cycle, full cycle and 0.75 cycle for

Fig. 7. Sampling frequency vs. data window length HCDFT, FCDFT and proposed filter respectively. The phase

error is calculated in full scale (360 ).

IV. PERFORMANCE TEST

TABLE I

In this section, an input signal comprising fundamental SUMMARY OF TEST RESULTS

frequency component and exponentially decaying component

is used to test the performance of the proposed phasor Estimate Value Error estimate and DC offset elimination schemes. The time Filter Type

(cycle) Am (p.u.) (deg) ErrAm Err

constant is studied over a broad range 0.5 5 cycles. The

HCDFT 0.7623 5.6721 23.7734 15.0911

results are compared with the conventional DFT (full cycle

0.5 FCDFT 0.8473 46.6454 15.2655 3.7096

and half cycle) method. Select N = 400, Ts = 0.75 cycle, f0 =

Proposed 1.0034 59.2125 0.3387 0.2187

60 Hz, Am = 1 p.u., = 60 . Consider the severe condition, that HCDFT 0.6796 -11.1511 32.0379 19.7642

is the fundamental frequency component and decaying 1 FCDFT 0.8559 51.4980 14.4133 2.3617

component have the same amplitude. The input signal is Proposed 1.0027 59.2304 0.2662 0.2138

described in complex form: HCDFT 0.6522 -23.4222 34.7817 23.1728

y (k ) Am e k T / Am e j ( k / 200 ) 2 FCDFT 0.9005 55.4360 9.9481 1.2678

Proposed 1.0025 59.2037 0.2521 0.2212

Apply formula (7) and obtain the amplitude and phase angle HCDFT 0.6483 -28.1817 35.1729 24.4949

of the input signal y(k). Fig. 8 and Fig. 9 give the calculated 3 FCDFT 0.9262 56.9217 7.3844 0.8551

Proposed 1.0024 59.1922 0.2391 0.2244

amplitudes and phase angles in four cycles over = 1.0 cycle

HCDFT 0.6477 -30.6813 35.2314 25.1893

respectively using the full cycle DFT, half cycle DFT and the

4 FCDFT 0.9416 57.6840 5.8434 0.6433

proposed algorithm. The convergence time of the proposed Proposed 1.0024 59.1859 0.2411 0.2261

algorithm is 0.75 cycle, and the calculation errors are 0.2662% HCDFT 0.6478 -32.2173 35.2169 25.6159

and 0.2138% for amplitude and phase angle respectively. 5 FCDFT 0.9517 58.1454 4.8275 0.5152

Proposed 1.0024 59.1788 0.2442 0.2281

V. CONCLUSIONS

DC offset has significant effect on extracting the

fundamental frequency components. Conventional DFT filters

can not easily eliminate it when estimating the components of

interest. This directly impacts the accuracy of the fundamental

frequency component based protective relaying algorithms.

This paper proposes a novel method for estimating the

fundamental phasor and eliminating the DC offset using

improved recursive wavelet transform. Studies indicate that

the convergence of proposed algorithm is related to the

sampling frequency. To achieve a certain level of accuracy,

the higher sampling rate one uses, the shorter data window it

Fig. 8 Amplitude of the estimated phasor

needs, and vice versa. The proposed method converges to

correct results within one cycle, and the computation burden is

fairly low because of the recursive formula. Comparing with

conventional DFT filter, performance tests demonstrate that

the proposed method achieves very good results.

I (a, b) D e t / Wx ( t ) (a, b)

VI. APPENDIX

The IRWT coefficient of the input signal x(t) is: where I (a,b) has the same with I(a,b) by replacing 1 with 2.

t b

W x (t ) (a, b) a x(t ) (

1 / 2

)dt b 0

VII. REFERENCES

3 t b 3 4 t b 4 [1] A. G. Phadke and J. S. Thorp, Computer Relaying for Power Systems.

a 1 / 2 Am e j ( t ) [ New York: John Wiley and Sons, 1988.

3 a 6 a [2] M. V. V. S. Yalla, A Digital Multifunction Protective Relays, IEEE

Trans. On Power Delivery, vol. 7, no. 1, pp. 193-201, 1992.

t b

5 t b 5 ( j 0 ) ( a ) [3] A. G. Phadke, T. Hlibka, M. Ibrahim, A Digital Computer System for ]e dt EHV Substation: Analysis and Field Tests, IEEE Trans on Power

15 a Apparatus and Systems, Vol. PAS-95, pp. 291-301, Jan. /Feb. 1976.

[4] N. T. Stringer, The Effect of DC Offset on Current-operated Relays,

Denote IEEE Trans on Industry Applications, vol. 34, no. 1, pp. 30-34, Jan/Feb

t k a b, k [ b / a,0] 1998.

[5] G. Benmouyal, Removal of DC-offset in Current Waveforms Using

1 j (a 0 ) Digital Mimic Filtering, IEEE Trans on Power Delivery, vol. 10, no. 2,

pp. 621-630, April 1995.

We have [6] Jyh-Cherng Gu, Sun-Li Yu, Removal of DC Offset in Current and

b a e [ j ( a 0 )] k

j ( b ) Voltage Signals Using a Novel Fourier Filter Algorithm, IEEE Trans

W x (t ) (a, b) Am e on Power Delivery, vol. 15, no. 1, pp. 73-79, Jan. 2000.

a [7] Jun-Zhe Yang, Chih-Wen Liu, Complete Elimination of DC Offset in

Current Signal for Relaying Applications, in Proc. 2000 IEEE Power

3 4 5

( k3 k4 k 5 )dk Engineering Society Winter Meeting, vol. 3, pp. 1933-1938, Jan 2000.

[8] T. S. Sidhu, X. Zhang, F. Albasri, M. S. Sachdev, Discrete-Fourier-

3 6 15 Transform-Based Technique for Removal of Decaying DC Offset from

j ( b )

Am e I ( a, b) Phasor Estimates, IEE Proc. Generation, Transmission and

Distribuation, vol. 150, no. 6, pp. 745-752, Nov. 2003.

where [9] Y. Guo, M. Kezunovic, Simplified Algorithms for Removal of the

b b Effect of Exponentially Decaying DC-Offset on the Fourier 3 ( b ) 3 2 Algorithms, IEEE Trans on Power Delivery, vol. 18, no. 3, pp. 711-

[ a a 717, July 2003.

I ( a, b) a { e a [10] Chi-Shan Yu, A Discrete Fourier Transform-Based Adaptive Mimic

1 12

3 Phasor Estimator for Distance Relaying Applications, IEEE Trans on

Power Delivery, vol. 21, no. 4, pp. 1836-1846, Oct. 2006.

b [11] O. Chaari, M. Meunier, F. Brouaye, Wavelets: A New Tool For the

6 ( b )

b b

a e 1 a 6 (1 e 1 ( a ) )]

1 Resonant Grounded Power Distribution Systems Relaying, IEEE Trans

on Power Delivery, vol. 11, no. 3, pp. 1301-1308, July 1996.

1 1

3 4

[12] Chuan-li Zhang, Yi-zhuang Huang, et al., A New Approach to Detect

Transformer Inrush Current by Applying Wavelet Transform, in Proc.

b 4 b 1998 POWERCON, vol. 2, pp. 1040-1044, Aug. 1998.

4 3 ( b ) b

4 [13] Xiang-ning Lin, Hai-feng Liu, A Fast Recursive Wavelet Based

1 [ a a 1

e a e a Boundary Protection Scheme, in Proc. 2005 IEEE Power Engineering

1 12

6 Society General Meeting, vol. 1, pp. 722-727, June 2005.

[14] C. Zhang, Y. Huang, X. Ma, W. Lu, G. Wang, A New Approach to

b b Detect Transformer Inrush Current by Applying Wavelet Transform, in

12 2 ( b ) 24 ( b ) Proc. 1998 POWERCON, vol. 2, pp. 1040-1044, Aug. 1998.

a a e 1 a

e a

13 14 VIII. BIOGRAPHIES

b 5 Jinfeng Ren (S 07) received his B.S. degree from Xi an Jiaotong University, Xi an, China, in electrical engineering in 2004. He has been with

b b

24 1 ( 1 [ a

) Texas A&M University pursuing his Ph.D degree since Jan. 2007. His

5 (1 e )] e a

research interests are digital simulator and automated simulation method for

1 15 protective relay and phasor measurement unit testing and applications in

power system control and protection.

b b

5 4 ( b ) 20 3 ( b ) Mladen Kezunovic (S 77, M 80, SM 85, F 99) received the Dipl. Ing.,

a a M.S. and Ph.D. degrees in electrical engineering in 1974, 1977 and 1980,

1 1

e a e a respectively. Currently, he is the Eugene E. Webb Professor and Site Director

12 13 of Power Engineering Research Center (PSerc), an NSF I/UCRC.at Texas

A&M University He worked for Westinghouse Electric Corp., Pittsburgh, PA,

b 2 b 1979-1980 and the Energoinvest Company, in Europe 1980-1986, and spent a

60 120 ( b )

b b

a e 1 a 120 (1 e 1 ( a ) )]}

1 a sabbatical at EdF in Clamart 1999-2000. He was also a Visiting Professor at

e a Washington State University, Pullman, 1986-1987 and The University of

1 1 1

4 5 6

Hong Kong, fall of 2007. His main research interests are digital simulators and

Similarly as deriving the wavelet coefficient of signal x(t), simulation methods for relay testing as well as application of intelligent

methods to power system monitoring, control, and protection. Dr. Kezunovic

for signal y(t), denote 2 a / j 0, we have is a Fellow of the IEEE, member of CIGRE and Registered Professional

t b Engineer in Texas.

Wy ( t ) (a, b) a 1 / 2 y (t ) dt

Contact this candidate