WO2011089421A1 - Apparatus and method for measuring a phasor of an electrical power system signal - Google Patents

Abstract

An apparatus and a method for measuring a phasor of an electrical power system signal (20) is described. The power system signal (20) has a cycle having a predetermined period and is sampled to produce a sample stream (22) of the electrical power system signal (22) The apparatus comprises a processor ( 14) configured to measure a phasor in phase-space by using the value of a first sample from the sample stream (22) to represent the x-coordinate of the phasor and by using the value of a second sample from the sample stream (22) to represent the y-coordinate of the phasor The second sample is delayed by a quarter of the predetermined period from the first sample. In further implementations, the apparatus and method may take into account a DC offset of the power system signal (20) The apparatus may then further comprise a DC offset calculation unit(24) for calculating a DC offset (26) of the sample stream (22) The processor (14) is further configured to subtract the DC offset (26) from the sample stream (22), thereby generating a DC-corrected sample stream The DC-corrected sample stream is then used to determine the x and y coordinates of the phasor

Description

APPARATUS AND METHOD FOR MEASURING A PHASOR OF AN ELECTRICAL

POWER SYSTEM SIGNAL

The present invention relates to an apparatus and method that can measure a phasor of an electrical power system signal. The apparatus and method may be used in phasor measurement units, relays and many other electrical power system devices.

Phasor measurement is used in various devices that form part of electrical power systems, whether for transmission or distribution. Electrical power systems transmit and distribute electrical power using AC signals with a defined frequency. For example, in the UK the electricity system has a nominal frequency of 50 Hz. Phasor measurement expresses the current and/or voltage of the electrical power system at the point of measurement as a phasor with an amplitude and phase angle representing the fundamental component of the signal. Phasor measurement has a number of different applications in electrical power systems. It can be used to determine whether a fault has occurred, for example in a protection relay. Phasor measurement can also be used to monitor the overall quality and stability of power system voltage measured by phasor measurement units at various points on the network, in particular the power system dynamic process. Measuring a phasor by analysing an electrical power system signal (voltage or current) using a Fourier Transform is known. However, this is computationally complex and susceptible to errors and transients, and subject to a significant time delay. It is desirable to measure a phasor of an electrical power system signal with reduced complexity and without a large time delay.

According to the present invention, an apparatus and method make use of the surprising result that if a power system signal is sampled the signal can be represented in phase-space by an x-coordinate represented by a first sample and a y-coordinate represented by a second sample delayed by one quarter of a cycle of the power system signal.

This result arises when Delay Coordinate Embedding techniques are applied to power system signals. Although Delay Coordinate Embedding is known, it has not been applied to electrical power system signals. Examples of Delay Coordinate Embedding are discussed in F. Takens, "Detecting Strange Attractors in turbulence", Lecture notes in mathematics, 898:366-381 , 1981 and in J.D. Fanner and J.J. Sidorowich "Predicting chaotic time series", Physical Review Letters, 59(8):845-848, August 1987. Apply delay coordinate embedding to a power system signal is discussed in more detail in the following detailed description. It leads to the surprising result that a sampled power system signal can be represented in phase-space by simply taking a first sample as the x coordinate and a second sample, delayed by one quarter of a cycle from the first sample, as the y coordinate. According to a first aspect of the present invention, there is provided an apparatus for estimating a phasor of an electrical power system signal with a cycle having a

predetermined period from a sample stream of the electrical power system signal, wherein the apparatus comprises:

a processor configured to estimate a phasor in phase-space by using the value of a first sample from the sample stream to represent the x-coordinate of the phasor and by using the value of a second sample from the sample stream to represent the y-coordinate of the phasor, wherein the second sample is delayed by a quarter of the predetermined period from the first sample. The electrical power system signal may be current or voltage. Within an electrical power system a nominal frequency is set, for example 50 Hz in the UK, which enables the period to be determined.

The sample stream is preferably sampled at a rate at least four times the nominal frequency, so that a sample delayed by a quarter of the predetermined period is available. The first and second sample need not be contiguous within the sample stream. Higher sampling rates may be advantageous to improve accuracy and speed up response time. Although it is preferred that the sampling rate is chosen so that samples may be taken directly from the sample stream this is not essential. For example the sample stream may be subsequently oversampled or interpolated so a sample delayed by a quarter of the predetermined period is available. The processor may be a microprocessor, FPGA, ASIC, system-on-chip or any other suitable integrated circuit.

The apparatus of the present invention can measure the phasor of a power system signal in phase-space in an efficient and quick manner. It requires a delay of only a quarter cycle for the second sample to be available. Far less computational resources are required than prior art methods based on Fourier Transfonns. The representation in phase space can be used directly in some embodiments, for example by comparing a Euclidean norm against a threshold to determine the presence of a fault.

In a preferred embodiment, the processor may be further configured to measure the phase, φ, of the phasor by evaluating φ = arctan — , wherein x is the first sample and y is the x J

second sample. This allows the angle to be calculated in a simple, accurate and efficient manner. The simple and accurate measurement of phasor angle in not easily available in prior art techniques. For example prior art analysis using Fourier Transfonns can only get accurate results in the case where the power system signal is static in a fixed frequency process without any transients involved. If required, in some embodiments the angle may also be expressed relative to an arbitrary sine wave calculated by: arctan (A mod(iy/₎₍ ,2^) where τηοά(ωί_η ,2π) is the residue of cot,, divided by 2π.

[ x j

In one embodiment, the apparatus may further comprise a median filter for post-processing the measured phasor. This can improve the accuracy of the phasor measurement. Use of a median filter has been found to give a better rise time because it is better at preserving edges and removing spikes than a mean filter. Preferably the median filter is a half-period median filter to give a suitable trade-off between the speed reduction introduced by the filter and the improved accuracy. In a further advantageous embodiment the median filter is cleared after a fault to reduce rise time and speed up response time when a fault signal is present. The processor ma also be further configured to measure the amplitude, A, of the phasor by evaluating A =

, wherein x is the first sample and j is the second sample. This allows the magnitude of the fundamental component to measured quickly and simply. An electrical power system is affected by transient events such as faults, switching operations, change in rotor angle of generators, lightening strikes, high frequency resonances, harmonics and reflections. These transient events alter the signal and present problems for measuring the phasor. High frequency transients may be removed by anti-aliasing and other filters. Low- frequency transients present more difficulty. A particular problem is presented by step changes caused by faults and switching operations, as well as the low frequency oscillations following such a disturbance due to rotor angle shift. The step change can produce a DC offset to the signal and causes the result of phasor measurement to be inaccurate for several cycles. Although acceptable for some applications, such as phasor measurement units, this is undesirable for many applications, particularly for protection relays which need to detect a fault as soon as possible.

It is known to remove the effects of a DC offset from a sampled power system signal. Generally, known techniques fall into two categories. The first category includes methods which try to reject the DC offset based on a model. Examples include digital mimic filters and LES filter models which model the DC offset by a Taylor series expansion. In general methods in the first category have limitations on accuracy, are sensitive to noise and computationally complex. The second category includes methods which try to estimate the DC offset so that it can then be removed from the signal. Examples include methods based on the Discrete Fourier Transform (DFT) and Half-cycle Discrete Fourier Transform

(HCDFT).

In an embodiment of the present invention, errors introduced by a DC offset in a fault component may be compensated for by the apparatus further comprising a DC offset calculation unit for calculating a DC offset of the sample stream; and wherein the processor is further configured to subtract the DC offset from the sample stream, thereby generating a DC-corrected sample stream and wherein the processor is configured to use the DC-corrected sample stream to determine the x and y coordinates of the phasor. This allows improved accuracy. The calculation of a DC-corrected sample stream can add complexity and increase the delay in the measurement of the phasor. It is desirable to reduce the complexity and delay for calculation of the DC offset. In an advantageous embodiment the DC-offset is calculated using Mathematical Morphology techniques. This allows estimate of the DC offset with less complexity than method based on Fourier transforms.

The paper "Morphological transform for removal of exponentially decaying DC-offset", Electronics Letters, 44(9):595-U104, April 2008 by Z. Lu, T.Y. Ji and Q.H. Wu proposes using a mathematical morphological transform to extract decaying DC offset components which are contained in fault currents in electrical power system signals. It discloses a mathematical morphological transform that can be used to extract a DC offset from a sampled signal. Methods are described which require as little as a quarter of a cycle to compute the DC offset. The arithmetical calculations involved in the Morphological Transform are simple, reducing the computational burden. The paper describes using the Morphological Transform to calculate the DC offset which is then subtracted from the signal. However, the resulting DC-Corrected signal is then processed using a Fourier Transform to determine the amplitude of the fundamental frequency.

It has been found that the combination of DC-offset calculation by mathematical morphology techniques has synergistic advantages when combined with phase-space representation to measure the phasor. Both techniques can be carried out with simple operations, enabling them to be implemented in simple and inexpensive hardware. They are well suited to implementation in an FPGA, for example.

In some embodiments, the apparatus may further comprise a fault detection unit for detecting a fault on the electrical power system by monitoring the phasor. Various techniques may be used. In one embodiment the fault detection unit is configured to monitor the Euclidean norm of the phasor and compare it to a threshold. Other embodiments may monitor the angle, φ, of the phasor. For example, the fault detection unit may be configured to monitor the rate of change of the phase shift between the measured phase and a nominal signal with the predetermined period to detect a fault. This may be compared to a threshold and give an accurate detection of the start of fault.

The performance may be improved in some embodiments by compensating the DC offset depending on when the fault is detected. If the fault detection unit detects a fault after the peak or trough in the half cycle of the electrical power system signal then the DC offset calculated by the DC offset calculation unit is held constant for the duration of the half cycle. If the fault detection unit detects a fault before the peak in the half cycle of the electrical power system signal, then the DC offset calculated by the DC offset calculation unit is held constant until the point at which the fault is detected, following that point the DC offset calculated from the second half of the cycle is used. These measures enable DC offset calculation to use mathematical morphology techniques based on a quarter cycle.

The apparatus of the invention may be used in an electrical relay, which may be a protection relay.

The apparatus of the invention can also be used in a phasor measurement unit in combination with a transmitter for transmitting the estimated phasor to a control centre. The transmitter may be wired or wireless.

In another embodiment of the invention, an electrical power network comprises a plurality of the phasor measurement units described above.

According to another aspect of the present invention, there is provided a method of measuring a phasor of an electrical power system signal with a cycle having a

predetermined period from a sample stream of the electrical power system signal, the method comprising:

using a first sample from the sample stream to represent the x-coordinate of the phasor in phase-space; and using a second sample from the sample stream to represent the y-coordinate o phasor in phase-space, wherein the second sample is delayed by a quarter of the predetermined period from the first sample.

The method may further comprise measuring the phase, φ, of the phasor by evaluating φ = arctan wherein x is the first sample and y is the second sample.

The method may further comprise post-processing the measured phasor with a median filter.

In some embodiments the method may further comprise monitoring the rate of change of the phase shift between the measured phase and a nominal signal with the predetemiined period for detection of a fault on the electrical power system. The amplitude, A, of phasor by may be measured by evaluating A = ^x + y , wherein x is the first sample and y is the second sample.

The method may further comprise determining a DC-corrected sample stream by:

calculating a DC offset of the sample stream; and

subtracting the DC offset from the sample stream thereby generating the DC-corrected sample stream; and

wherein the steps of using a first sample and using a second sample use samples from the DC-corrected sample stream. The DC offset may be calculated using Mathematical Morphology techniques.

Embodiments of the invention will now be described by way of example only and not limitation with reference to the accompanying drawings, in which: Figure 1 is a graph of fault occurrence detection using the Euclidean norm of a coordinate embedded signal; Figure 2 depicts the embedding of a fault signal in a two-dimensional phase space;

Figure 3 is a graph which depicts the result of removing the DC offset from the fault current;

Figure 4 is a graph which depicts an estimate of the amplitude of the fault current;

Figure 5 is a graph which depicts an estimation of the phase of the fault current;

Figure 6 is a graph which depicts the estimation of the phase difference;

Figure 7 is a graph which depicts the calculation of p;

Figure 8 is a graph depicting the estimate of the phase of the fault current in one embodiment;

Figure 9 is a graph depicting the result of a mathematical morphological DC offset operator;

Figure 10 is a graph depicting the mirror technique in mathematical morphology of power system signals;

Figure 1 lA-1 1C depicts different fault onset points in the cycle of a power system signal; Figure 12 depicts a transmission line model used for simulation;

Figure 13 depicts an equivalent circuit of a transmission line with a fault;

Figure 14 depicts a comparison of FCDFT and method of an embodiment of the present invention with the DC offset removed against use of FCDFT of the original signal;

Figures 15A and 15B depict a comparison of FCDFT and method of an embodiment of the present invention with the DC offset removed against use of FCDFT of the original signal at a frequency of 53Hz and 47Hz, respectively;

Figure 16 depicts a comparison of FCDFT and method of an embodiment of the present invention with the DC offset removed against use of FCDFT of the original signal with a 64 Samples Per Cycle sampling rate;

Figures 17A-17C depict a comparison of FCDFT and method of an embodiment of the present invention with the DC offset removed against use of FCDFT of the original signal for different time constants;

Figures 18A and 18B depict a comparison of FCDFT and method of an embodiment of the present invention with the DC offset removed against use of FCDFT of the original signal under different phase shift condition; Figures 19A and 19B depict a comparison of FCDFT and method of an embodiment of the present invention with the DC offset removed against use of FCDFT of the original signal using real data collected from a power system;

Figure 20 is a diagrammatic representation of a protection relay according to an embodiment of the present invention.

The theory behind phasor measurement in electrical power system signals through embedding used embodiments of the present invention will first be described. According to delay coordinate embedding theory, a ίϊΈ-dimensional phase space can be constructed given the delay coordinate, τ, and the signal forms a matrix described

by:

J^' l + r ^{■ χ} 1 + (ά_Ε - 1 )τ

x2 X-2 X + T ^{■ ■}i^;2 + ( d_E - 1 )r

X = =

_ XM _ X- M XM +T ^{■ "} ·¾/+ {<½ - 1 ) 7 The columns construct the coordinates of the phase space and each row vector Xj of the matrix represents a point in the phase space.

In electrical power systems, a source current or voltage signal has the standard sinusoidal waveform:

where ω = 27r/Twith T the period. On the other hand, a fault current or voltage signal is considered to comprise a fundamental sinusoidal component, a series of harmonics and an exponential decaying DC offset, which can be expressed as:

hit)

<k(i) (6) where 7Ί is the fundamental component and it has the same frequency as the source signal but has a phase shift of an angle of φ - Φ, Ι_\λ is the combination of the harmonics, and _c denotes the DC offset. Therefore, the signal including a source part and a fault part can be expressed by: h(t), o<t<t_s

/₅ (/ /,). / /.

(7) where t_s denotes the time instant that the fault occurs. An example of such a signal is given in Figure 1 in dotted line. In this case, up to 15 harmonics are included in the signal, and a Gaussian noise of a signal-to-noise ratio of 15 dB is added.

Embedded to a 2-dimensional phase space, fault signal / forms the following matrix:

I(t_n) I{t_n +r)

(8) where t_n = t + nAt with t the beginning and Δί the sampling interval. Matrix I can be considered as a 2-dimensional signal in the phase space with the left column being its x- values and the right column its corresponding v-values. The detection of a fault on the signal will now be described. Due to the features of circular function, the source signal, has the following embedding when τ = 774:

Io - [ /... (*) h (t + r) }

= [ x y ] ₍₉₎

where

./· = A{) cos (ά /· + φ)

2 2 2

As x +y = AQ , it shows that the source signal forms a circle whose radium is A₀ and whose centre is at (0, 0) in the phase space. In other words, a pair of samples, I (t„) and Ι (ί_η+τ), form a point in the phase space and the Euclidean norm of the point is AQ.

However, when the fault occurs, the fault part forms some other shape and the Euclidean norm of the point (I(t_s - τ), I(t_s)) suddenly increases. Therefore, a threshold can be set to check if the fault occurs. The embedding of the signal shown Figure 1 in dotted line is given in Figure 2. Dots depict the embedding of the source part and crosses depict the embedding of the fault part. The Euclidean norms of the points are plotted in Figure 1 in solid line, and the estimated fault occurrence point is highlighted by a dot. In this case, it is simulated that the fault occurs at t_s = t₈₇ and the detection result is also i_s = t_g7 .

In an alternative embodiment a mathematical morphological filter may be used to remove the DC offset. This is described in more detail below. In the following description of more detailed embodiments measuring particular characteristics of a signal, the estimation of amplitude and phase angle is based on the signal whose DC offset is removed by the morphological filter. An example is given in Figure 3. In Figure 3 the dotted line is the fault current, the dashed line is the extracted DC offset and the solid line is the result of removing the DC offset from the fault current. Estimation of the fundamental amplitude in one embodiment of the invention will now be described. For the fundamental component, 7i(t), its embedding in the phase space can be expressed by:

+ τ) -f φ)

= ..4₁ cos Θ cos(ijt + φ) — Αχ sin Θ sin(u,'i + φ) where θ = ωτ = 2πτ/Τ. From the equation (12) and (13) for x and y, the amplitude of the fundamental component can be calculated from:

where A_x denotes the estimation οϊΑ_\. Equations (12)-(14) show that any two samples from the fault current are enough to estimate its amplitude. Usually τ is selected to make Θ a common angle. For example, when τ = 778, θ = π/4. In this case, as the signal is digitised at N= T/At samples per cycle, it uses N/8 samples to calculate A i . When τ = 774, we have θ = 2πτ/Τ= I2 and the embedding of Ij becomes

y = f] {t + r) = -A ^'n i L / + φ) . The estimation therefore becomes:

Since

{ A_\ os( t + ψ) )² + { -Ai sm{ t + φ) )^'2 — A

(16) the embedded signal

fonns a circle in the phase space, and the centre and the radius of the circle are the origin of the phase space and A _\, respectively. In practice, due to the presence of noise and haraionics, the estimated amplitude calculated from different pairs of samples varies slightly:

1-1 I + vl

where A_hl denotes the estimation of_^4 | calculated from x_n = I(t_n) and j^„ = I(t„ + τ). Assume that the onset of the fault is detected at time f_s . From the time instant of t = i_s + At , the estimated amplitude is calculated using (14), where Θ = a>{t_n - i_s), until t_n = t_s + τ - At . In this manner, the estimation of the fault current amplitude is not affected by the source current amplitude, except the amplitude at time instant t_s . In order to reduce the estimation error caused by noise, the average of A over a user-defined window can be used. An example is given in Fig. 4. A fault signal that has its DC offset removed is considered as the input signal, as plotted in dotted line in Figure 4. The signal-to-noise ratio of the fault signal is 10.67 dB. Using the method described above to estimate the amplitude of its fundamental component, the result is given in solid line in Figure 4. As it can be seen, the sudden increase in amplitude has been successfully estimated as the fault occurs. The vibration caused by the noise is avoided by taking the average of the estimated amplitude. As a comparison, the estimation result by the halfcycle Fourier Transform (FT) is also included, which is shown in dashed line in Figure 4. The half-cycle FT uses the samples over a half of a cycle to calculate the amplitude. Hence, it causes a half-cycle delay to accurately estimate the amplitude of the fault signal. On the other hand, the embedding- based method does not introduce such an error in the estimation.

Estimation of phase in another embodiment of the invention will now be described. As stated previously, when τ = Γ/4, the embedded fundamental component can be defined by: COS(O'/O + φ) — A y sm(u io

Ai cosf ,, - ) — A [ sin (^f„ + )

(1 ί

When to = 0, the phase of the fundamental component can be calculated from the first point of l, :

Ii (1. 2)

ψ = arctan

i i f i , i i (19)

The phase can also be calculated from an arbitrary point of Ii using:

arctan I— ^ ^' ~ } 1 — mod(i ,'/:„ . 2π) .

1 ₁ ( n . 1 (20) where mod( oi,„ 2π) is the residue of cot_n divided by 2π. This is to make sure that φ_η falls in the range of [0, 2π]. Ideally, for different n, φ_η should be the same. In practice, it is more accurate to use the average of φ_η over a certain window to eliminate the influence of any possible noise:

The test is also carried out using the above method and the half-cycle FT, respectively, and the results are given in Figure 5. In this case, the signal-to-noise ratio of the fault signal is also 10.67 dB. The actual phase is 78.7500°, and the estimated results by embedding and the half-cycle FT are 79.4820° and 79.7806°, respectively. The method of this embodiment is slightly better with an estimation error of 0.93%, while the error of the half-cycle FT is 1.31 %. However, as it does not need to involve samples from a whole cycle in calculation, the method of this embodiment is more computationally efficient. Phase difference detection in another embodiment of the invention will now be described. In some cases, there is a phase difference between the cun-ent signal and the voltage signal. The difference can be detected through the embedding of the two signals to the phase space with τ = T/4. Define the current and voltage signals as follows:

f_y = A_v cos(iJt + φ_ν )

(23) where A_c and A_v are their amplitudes and φ₀ and φ_ν are their phases, respectively.

When τ = T/4, an arbitrary point from the embedded current signal has the coordinates of (x_cn, y_en) = (A_c cosfft)/,, + <p_c), -A_c η(ωί„ + <p_c)) and has a phase angle of

φ_αι = arctan(- y_cn I x_ai ) . At the same sampling time, the point from the embedded voltage signal has the coordinates of (x_vn, y_vn) = (A_r cos(cot_n + φ^ -Α_γ sin(<¾i„ + < >,,)) and has a phase angle of φ_νη = arctan(- y_m I x_m ) . Therefore, the phase difference between the two signals, denoted by Δφ, can be calculated from Αφ = φ_αι - φ_νιι . Considering the interference of noise, it is more accurate to use the average of the phase differences over a certain window:

Figure 6 shows the estimation result of a test, where the input voltage and current signals both have a signal-to-noise ratio of 10.67 dB. The simulated phase difference is 307.7440° and the estimated phase difference at each sampling instant varies between 306.9407° and 310.1625°, as the figure shows. The average estimation error is 0.12%.

Fundamental frequency shift estimation in another embodiment of the invention will now be described. In some cases, power system signals are influenced by fundamental frequency shift, which can be up to ±5%. A considerable eiTor would be introduced to the result of phasor measurement if the signal is embedded according to its nominal period. To estimate the actual fundamental frequency, a method is proposed as follows. The signal is first embedded with τ = 772, where T is the nominal period. Hence, the embedding in the phase space is expressed by: x = {t) = Λ _\ ι 'ί + ψ)

y = I₁ (t + T) = A₁ CS>s{ J^, {t + T/2) + )

= A_} cos(u f + ψ +— )

(26) where ω' denotes the actual fundamental frequency. As the function of

F(t) = x + y

= 2,4_! cos(—— J cos{ t + φ +—— )

"^ω (27) fonris a sinusoidal signal, ω' can be calculated from ω'— 2ω arccos(±p/2 4i) /7T where p is the amplitude of F and is calculated from p = {F @ g - F Q g) /2 (29) where + and Θ denote the operators of dilation and erosion used in mathematical morphological techniques, respectively, and g is a flat structuring element of a cycle 1 An example is given in Fig. 7 to show the calculation of p. In equation (28), if ω' > ω, used; otherwise, p is used.

To determine if ω' is larger or smaller than co, the following strategy is used. As stated above, if there is no fundamental frequency shift, the estimated phase, φ_η , is a flat function. However, if the fundamental frequency has a positive variation, i.e. ω' > co, ψ_η increases gradually. On the other hand, if co' < co, <p_n decreases. For the fault current shown in Fig. 7, the actual fundamental frequency is f = 51 .2542 Hz. Estimating the phase when the delay coordinate is still set at τ = 774, where T is the nominal period, the value of φ_η is depicts in Fig. 8. As φ_η is basically an increasing function, it is considered that co' > co and -p is used in equation (28).

According to equation (28), the estimation of co' also requires the value oiA _\ . The estimation process is the same as described above, and the final result οϊΑ^" _\ is calculated as the average ϊΑ _η over a cycle. In this example, A

= 3.0124. The estimated fundamental frequency is ^"/ = 51.2794 Hz and the estimation error is 4.92 ^χ 10 ⁴. On the estimation of the fundamental frequency, the fault current is re-embedded according to the estimated actual fundamental frequency and phasor measurement can be carried out afterwards. As samples from at least a cycle of the fault current are used to estimate the fundamental frequency, the method causes a delay of at least a cycle at the very beginning.

Thus, an embedding-based method has been described in the above embodiments for phasor measurement in electrical power systems. The method takes advantage of the mathematical properties of a power system signal to transfer it to a 2-dimensional phase space through delay coordinate embedding. In this manner, the amplitude and phase angle of a current/voltage signal and the phase difference of current and voltage signals can be calculated sample by sample. The calculation involves two or four samples only, unlike traditional FT-based methods that use samples of half a cycle or an entire cycle. Moreover, the method can also be used to estimate the actual fundamental frequency when it deviates from its nominal value.

As mentioned above, some of the embodiments of the invention use mathematical morphology (MM) to extract the DC offset from the power system signal, so that a DC- corrected power system signal is produced. The principle is that there are two sets: a larger set representing the image or signal which will be denoted S and a smaller set representing the structuring element (SE) which will be denoted B. The structuring element is then moved through the image/signal pixel by pixel performing the required operator. The basic operators are erosion and dilation, the greyscale operators are defined as:

(S Θ B) ( x) =

( S θ B ) ( x) = max { Six - v) + B( v ) \

·'·— "e.¾ '-e J ^■ (31 ) where x corresponds to a pixel in the image/signal and v a pixel in the structuring element.

To remove the DC offset a one dimensional structuring element is used, as a current or voltage signal is one dimensional. The following operator is developed to extract the DC offset component, where the SE is a flat SE with length equal to half a cycle of the fundamental frequency.

DC Orfset = :

(32) Where ° is an opening, an erosion followed by a dilation and defined by equation (33) below. · is a closing, a dilation followed by an erosion and defined by equation (34) below. The opening extracts the DC offset from the peak region, and the closing from the trough region. This is illustrated in Figure 9. The delays between the different operators have not been shown, so as to clearly display the result of each operator. It can be seen that the DC offset operator removes the DC offset, and the components opening and closing, remove the DC offset in the peak part and the trough part respectively. The delay for this DC offset operator is 3/4 of a cycle, the delay for just a closing or opening is a 1/2 cycle.

B = (S O B)

B = (S Φ B) (34)

To reduce the delay of the MM DC offset removal operator given in equation (32) the following steps may be performed in some embodiments. Firstly the delay increases in relation to the number of operators performed in series, as the window needs to be larger. Therefore if we use the opening operator in the peak half cycle and the closing operator in the trough half cycle the delay can be reduced to a 1/2 cycle as mentioned above. Secondly if the window is set to half the structuring element size, and a mirroring is performed at the peak or trough point the delay can be reduced to 1/4 cycle as proposed in the paper "Morphological Transform for removal of exponentially decaying dc-offset" by Lu et al acknowledged above. This is illustrated in Fig. 10. Though to improve the performance the result of the two halves of the mirror are averaged resulting in a delay of a ½ cycle but better performance. Impulse noise, as that caused by ADC misread, interference, etc, can cause the peak detection to mislocate the peak. To solve this problem the input signal may be firstly processed by a morphological filter to remove this noise. The filter used is the average of an opening and a closing given by:

_ S o B + S » B

^{V ~~} 2 (35) where y is the output of the filter, and the length of the structuring element and the window are 3 samples. This causes an additional delay of 1 sample. Therefore the delay is, 1 sample for the filter, and one half cycle for the morphology DC offset detection, giving a total delay of a half cycle and 1 sample. Some embodiments may handling the phase shift at fault onset. When a fault occurs, there is a change of impedance, which will introduce a phase shift. The phase shift will change the width of the period during the fault onset. The effects of the phase shift can be divided into three categories as illustrated in Fig 1 1 A-C. In case 1 , depicted in Figure 1 1 A, the fault instant is following the peak in the current half cycle. In case 2, depicted in Figure 1 I B, the fault instant is proceeding the peak in the current half cycle. In case 3, depicted in Figure 1 1 C, the fault instant is at the same instant as the peak in the current half cycle.

In case 1 the DC offset value should be held, for the duration of the half cycle. For case 2 the DC offset value should be held until the fault time, after the fault time the value from the second quarter of the half period should be used. In case 3 since the value of the fault instant is the same as the peak instant, the remaining part of the half cycle must be less than or equal to a quarter cycle. Hence the structuring element is reduced to a quarter cycle. Case 1 and 2 accurate represent the signal, in case 3 the error is minimized to an acceptable level.

As discussed above, following the removal of the DC offset the signal can be embedded to a 2 dimensional phase space as given by the following matrix:

I ( to ) I(to + τ)

(36) where t_n = t₀ + nAt with to the beginning and At the sampling interval. The left column represents the x value and the right column the y value in the phase space. Each row represents a different sampling point.

If we assume the signal only consists of the fundamental frequency the coordinates in the phase space are given by:

A i cosf 'i + φ)

Ai cos(c / + φ + θ) (38) where θ = ωτ = 2πτ/Τ , T is the samples per cycle and A] is the fundamental amplitude. From (37) and (38) the amplitude of the sinusoidal signal can be calculated as:

Ai ² = (x²— 2xy cos Θ + ij²)/ sur Θ

If τ = T/4 then the x and y coordinates are given by: x = A\ cos(tfctf + ψ)

y — Ai sm( t + φ) Since these are orthogonal components the magnitude is given by:

A_x = v^² + ir ( ) The angle is given by:

Ψ— arctan (—

■ xJ (43) The phase can also be given relative to an arbitrary sine wave: vrctan (— j— ΐΆ ά(ωί_η, 2π)

(44)

Hence the Magnitude and phase can be calculated with a delay of ¼ cycle instead of one whole cycle if using the FCDFT, or a 1/2 cycle using the HCDFT. To improve the accuracy a median filter of length ½ cycle may be used to post-process the measurement. The median filter gives a better rise time as it is better at preserving edges, and is better at removing spikes than a mean filter. To give a quick rise time the filter is preferably cleared when a fault occurs as given by the following equation:

A{x) = median, ... . .. a{S(.r + v)} ₍₄₅₎ where

B = 6.1.6₂ · ·,/>/- J -fy : / .." or // · .r

B = ^h _.f- bf₊i...bi-i. hi : .r ^~ I < ft < ,r. (₄₆₎ where / = 772 and f=x -fi.fi is the fault time, and/ the position of the fault in the window. In one embodiment, to detect the fault start point, the rate of change of phase shift between the signal and a nominal 50Hz signal is used. This is derived from equation 44.

When the rate of change of phase shift is greater than a threshold, the fault start is detected. To prevent the fault start being detected multiple times, the rate of change of phase shift preferably has to settle below a lower threshold before a new fault start can be detected. In another embodiment, frequency compensation is used. When the frequency changes the length of a period becomes longer or shorter, and therefore the length of the structuring element is not quite the right size. If the frequency increases the length of the period becomes shorter, hence the structuring element is longer than half a cycle. This means the DC offset value is lower in the peak half cycle and higher in the trough half cycle, causing the magnitude of the signal to be increased. If the frequency is decreased the converse is true.

To compensate for the change in frequency, firstly the frequency deviation is measured. The frequency deviation is the average rate of change of phase shift between the signal and a nominal 50Hz signal derived from equation (44).

¾H J- Δ ₂ (*)

t=x-T ( 9)

Secondly the magnitude calculated by the Phasor Measurement is multiplied by twice the proportion change in frequency. This is due to the increase in magnitude in the peak half cycle added to the increase in magnitude in the trough half cycle.

1 + (&ψ₂ (χ)/Μ' )) * 2 (50) Thirdly the Phasor Measurement is "tuned" to 50Hz, which causes oscillations in the output. To prevent this the Phasor Measurement can be "retuned" by altering the Θ in Eqn. 39. When using the assumption τ = T/4 such as to calculate the phase the value of τ can be updated accordingly assuming a significantly high sampling frequency is used.

The result of using embodiments of the invention will now be discussed with reference to simulation results. The simulation is preformed in Matlab, a commercially available computer program. Simulations were carried out using data generated from PSCAD (a computer commercially available from Manitoba HVDC Research Centre Inc. of

Manitoba, Canada), data based on the equations, and real-data gathered from

measurements of a real power system. The model used in the PSCAD software is shown in Figure 12. The fault is located equidistant between the two buses. The transmission line can also be modelled by the lumped parameter model shown in Figure 13 which leads to the following simple equations for current before and after the fault respectively.

1₂ = \ , si u + θ + [)) + Be^' (52)

Tests are designed to check the operation of the above described embodiment, they include testing over a range of frequencies from 47Hz to 53 Hz, testing with different values of the time constant τ, testing with real noisy data, testing with different phase shift conditions and testing with different sampling frequencies. Details and results of the simulation testing are given below.

First a general comparison of the methods of the present invention against existing techniques will be given. Figure 14 shows a comparison of the Full Cycle Discrete Fourier Transform (FCDFT) and the phasor measurement by embedding in phase space of the present invention (MEPS) with DC removed compared to the FCDFT of the original signal. The data is 50Hz, 32 samples per cycle test data generated from PSCAD. It shows that both the FCDFT and MEPS with DC removed, do not have the large over and under shoot, in the FCDFT of the original signal. It also shows that the MEPS is significantly quicker than the FCDFT. Next the methods of the present invention were tested in the presence of frequency deviation. The frequency is tested from 47Hz to 53Hz, 53Hz is shown in Figure 15A, 47Hz is shown in Figure 15B and 50Hz is shown in Figure 14. Comparison of these figures show the method of the invention works correctly over this frequency range. In Figures 15 A and 15B the offset caused by the basic morphology has been compensated for and hence does not appear in the output. It also shows the oscillations, as in the FCDFT of the original signal, due to wrong tuning have been removed in the output by correct tuning. Next the method of the invention was tested at different sampling rates. Three different sampling rates are tested: 32 samples per cycle(SPC) as shown in Figure 14, a higher 64 SPC sampling rate as shown in Figure 16, finally the real data shown in Figure 19 uses a sampling rate of 100 SPC. In all these cases the large overshoots and undershoots are removed by the method of the invention. In principle a higher sampling rate will produce a higher degree of accuracy inline with the more accurate representation of the signal.

Next, the method of the invention was tested with different time constants. Using equation (51 ) and (52) the time constant is varied between 0.5 cycles and 5 cycles and the effect on the DC offset removal compared. The results are shown in Figure 17A-17C. It shows that the large oscillations are successfully removed. For τ = 0.5 the value for the method of the invention in the first half period is a little lower, though for the FCDFT with DC removed the output is still raising, so the performance is good.

Next, the method of the invention was tested with different cases of phase shift. Three cases of phase shift have been tested and are depicted in Figure 1 8A and 18B. Case 1 and case 2 shown in Figure 18B can be illustrated using the same figure because a case 1 phase shift is usually followed by a case 2 phase shift. Figures 1 8A and 18B illustrate what is happening to the signals as well as what happens to the magnitude measurement. It can be seen that the correct value has been found for the DC offset, to successfully remove it from the signal. In the magnitude measurements, it can be seen that most of the large oscillations are removed by the method of the invention. Though for easel 12 the first oscillation is still significant, this may be due to the steep increase at the beginning of the fault. Finally, the method of the invention was tested with real data collected from a power system. Figure 19A and 19B depict an example with real data (provided by Siemens). Looking at Figure 19B the noise especially in the peak and trough regions can be seen. This is largely reduced in the output by the impulse noise preprocessor, and hence the method of the invention correctly identifies and removes the DC offset. Thus in the magnitude measurement as seen in Figure 19A the large oscillations are correctly removed.

To summarise the results of testing the method of the invention, Tables 1 and 2 below show the rise times and normalized maximum and minimum oscillation values

respectively.

Table 1 - Rise Times

Max oscillation Min oscillation

1 C 1 )1 1 , , M I PS ,. 1 ( 1 )1 1 , FCDFT, _;, MEPS,. FCDFT,

50Hz 32SPC 1. 1 14 1.01 1 1 .006 0.915 0.998 0.999

53Hz 1.066 1.013 1 .009 0.905 0.992 0.993

47Hz 1 . 126 1.031 1.013 0.947 0.994 0.990

64SPC 1.1 15 1.008 1.005 0.916 0.993 0.995 r = 0.5 1 . 130 1.021 1.01 1 0.953 0.997 0.998 r = 2 1 .070 1 .009 1 .005 0.945 0.991 0.998 r = 5 1 .033 1 .006 1 .005 0.970 0.993 0.998 casel/2 1.060 1.048 1.019 0.953 0.991 0.998 case3 1.051 1.008 1.003 0.934 0.992 0.998 real 1.066 1.019 1.012 0.941 0.979 0.986

Table 2 - Normalised Maximum and Minimum Oscillation Values

The rise time is defined as the sample time when the signal reaches 90% of the settled value. The maximum and minimum oscillation values, are the size of the largest peak oscillation and largest trough oscillation respectively, normalized to the settled value. From Table 1 it can be seen that the MEPS (the method of the invention) with DC removed (MEPSr) is 0.5 to 0.625 cycles faster than the FCDFTwith DC removed(FCDFTr). The delay between the MEPSr and the FCDFT of original signal(FCDFTo) is 0.24 cycles on average ranging between -0.03 cycles and 0.57cycles. From Table II it can be seen that by removing the DC offset improved the performance: the overshoot is significantly reduced from an average of 8.3% to an average of 1.7% for MEPSr. It can also be seen that the undershoot is significantly reduced from an average of 6.2% to an average of 0.8% for MEPSr. The results also show the benefits of the method of the present invention on rise times, with the MEPSr results showing significantly quicker rise times which will give improved performance when fast response is required, such as protection relays.

Thus, as discussed above, the method of the invention gives advantages when combined with mathematical Morphological removal of the DC offset. The DC offset can be calculated in one half cycle plus a few samples as a pre-processor to the phasor measurement of the invention. The results show that the MM preprocessor can successfully remove the DC offset. It has been shown that DC offset removal improves the results of the magnitude calculations, with significantly less overshoot. It has also been shown that by embedding the signal in phase space according to the method of the invention to perform the phasor measurement the delay is much shorter than using the conventional FCDFT. It will be appreciated that the method of the phasor measurement of the invention can also be used without DC offset removal, but performance is improved when DC offset removal is combined with the method of the invention.

Having discussed the theory and method of the present invention, its application to an example apparatus will now be explained by way of example with reference to Figure 20. Figure 20 depicts a high level diagram of a protection relay 10. The relay includes an analogue-to-digital converter (ADC) 12, a processor 14, a fault condition calculation unit 16 and an output stage 18.

The ADC 12 is connected to inputs of the voltage and current signals 20 measured from an electrical power system (not shown). The inputs may include one or more of these signals, for example measured at different points on the power system. In the simplest

embodiment a single voltage or current signal is received. The input signals are measured using known techniques, for example a current transformer (CT) to measure the current. The ADC 12 converts the input analogue power system signal(s) 20 into a digital sample stream 22. The digital sample stream 22 is provided as an input to the processor 14. Within the processor 14, the sample stream 22 is provided as an input to both a DC Offset calculation unit 24 and a phasor measurement unit 28. The DC Offset calculation unit 24 calculates the DC offset of the signal from the sample stream 22 and outputs the DC offset 26. In a preferred embodiment the DC offset is calculated using Mathematical Morphological techniques as described above, although this is not essential.

The phasor measurement unit 28 receives an input of the sample stream 22 and the DC offset 26. It subtracts the DC offset 26 from the sample stream 22 to generate a DC- corrected sample stream. The phasor measurement unit 28 then generates a phasor measurement from the DC-corrected sample stream by taking a first sample as the x coordinate and a second sample, delayed by one quarter of a cycle of the power system signal as the y coordinate. A detailed explanation of the theory of this is given above. Optionally the processor 14 may apply further correction to the signal as discussed above, for example expressing phase relative to an arbitrary sine wave or compensating for the point in the cycle at which a fault occurs. The phasor measurement unit outputs a phasor measurement 30 which forms the output from the processor 14. The phasor measurement 30 is provided as an input to the Fault Condition calculation unit 16. The fault condition calculation unit monitors the phasor measurement 30 for characteristics indicating a fault. For example the rate of change of phase relative to a nominal signal can be calculated and compared to a threshold; a fault is present if the rate of change exceeds the threshold. Methods for detecting faults were described in detail above and any or all of these may be used by the fault condition calculation unit 16, for example using the Euclidean norm of the phasor or a specific fault detection algorithm subject to the relaying functions. The fault condition calculation unit 16 outputs its determination of the fault condition to the output stage 18. In the event that a fault is detected the output stage 18 sends a signal 32 to trip a circuit breaker (not shown) to disconnect the affected circuit.

This embodiment can monitor one or more signals. For example it may be implemented to control a single set of circuit breakers for a single circuit. In alternate embodiments it can monitor several signals, for example several circuits and their associated circuit breakers present at a substation.

The present embodiment is well suited to production on an integrated circuit because the techniques of the present invention can be implemented with simple mathematical operators requiring less processing power than prior art methods. The various functional blocks described above may be provided as discrete components or integrated into a single integrated circuit. In an alternative embodiment, the DC offset calculation unit 24 is omitted. This simplifies construction but is less preferred because it reduces accuracy. In further alternate embodiments the phasor measurement unit may be provided separately, for example to function as a phasor measure unit without a protection function. In these embodiments the ADC and DC offset calculation unit may also optionally be included. When used a phasor measurement unit, a further embodiment (not shown) may incorporate a transmitter to transmit the phasor measurement to a controller, for example a central control centre of the power system. The unit may also include a time stamp unit to add an accurate timestamp (for example derived from GPS data) to the phasor measurement to assist comparison of measurements taken at different parts of a power system. A power system provided with such phasor measurement units distributed throughout the system can be operated more efficiently using the phasor measurements.

Claims

1 . An apparatus for measuring a phasor of an electrical power system signal with a cycle having a predetermined period from a sample stream of the electrical power system signal, wherein the apparatus comprises:

a processor configured to measure a phasor in phase-space by using the value of a first sample from the sample stream to represent the x-coordinate of the phasor and by using the value of a second sample from the sample stream to represent the y-coordinate of the phasor, wherein the second sample is delayed by a quarter of the predetermined period from the first sample.

2. An apparatus according to claim 1 , wherein the processor is further configured to measure the phase, φ, of the phasor by evaluating φ (A , wherein x is the first

[ x j

sample and j> is the second sample.

3. An apparatus according to claim 1 or 2, further comprising a median filter for post-processing the measured phasor.

4. An apparatus according to any one of the preceding claims, wherein the processor is further configured to measure the amplitude, A, of the phasor by evaluating

A -

wherein x is the first sample and v is the second sample.

5. An apparatus according to any one of the preceding claims, further comprising: a DC offset calculation unit for calculating a DC offset of the sample stream; and wherein the processor is further configured to subtract the DC offset from the sample stream, thereby generating a DC-corrected sample stream and wherein the processor is configured to use the DC-corrected sample stream to determine the x and y coordinates of the phasor.

6. An apparatus according to claim 5, wherein the DC offset calculation unit is configured to calculate the DC offset using Mathematical Morphology techniques.

7. An apparatus according to any one of the preceding claims, further comprising a fault detection unit for analysing the phasor and detecting a fault on the electrical power system.

8. An electrical relay comprising an apparatus according to any one of the preceding claims.

9. A phasor measurement unit comprising an apparatus according to any one of the preceding claims and a transmitter for transmitting the measured phasor to a control centre.

10. An electrical power network comprising a plurality of phasor measurement units according to claim 9.

1 1. A method of measuring a phasor of an electrical power system signal with a cycle having a predetermined period from a sample stream of the electrical power system signal, the method comprising:

using a first sample from the sample stream to represent the x-coordinate of the phasor in phase-space; and

using a second sample from the sample stream to represent the y-coordinate of the phasor in phase-space, wherein the second sample is delayed by a quarter of the predetermined period from the first sample.

12. A method according to claim 1 1 , further comprising measuring the phase, φ, of the phasor by evaluating φ = arctani— , wherein x is the first sample and _y is the second sample.

13. A method according to claim 1 1 or 12, further comprising post-processing the measured phasor with a median filter.

14. A method according to any one of claims 1 1 tol 3, further comprising monitoring the rate of change of the phase shift between the estimated phase and a nominal signal with the predetermined period for detection of a fault on the electrical power system.

15. A method according to any one of claims 1 1 to 14, further comprising measuring the amplitude, A, of the phasor by evaluating A = ^x² + y² , wherein x is the first sample and y is the second sample.

16. A method according to any one of claims 1 1 to 15, wherein the method further comprises determining a DC-corrected sample stream by:

calculating a DC offset of the sample stream; and

subtracting the DC offset from the sample stream thereby generating the DC-corrected sample stream; and

wherein the steps of using a first sample and using a second sample use samples from the DC-corrected sample stream.

17. A method according to claim 16, wherein the DC offset is calculated using Mathematical Morphology techniques.