The Half-Sine Method: A New Accurate Location Method Based on Wavelet Transform for Transmission-Line Protection from Single-Ended Measurements

In this work, a new and accurate method based on the wavelet transform is proposed for fault location in transmission-line systems. The proposed wavelet method consists of the analysis of the transient signal measured at a single end of the transmission line. Aerial current modes are used, and zero modes are included in the fault-detection scheme for low fault-inception angles. The fault distance is evaluated using the wavelet modulus maxima technique and a method based on the response to a half-sine voltage is proposed to overcome drawbacks arising from the limited sampling frequency and low fault-inception angle. The fault distance is calculated using the difference between the time when a 100 kHz half-sine signal is sent and the time when the derivative signal is received. The proposed algorithm is tested considering harmonic distortion and varying fault resistance, ground resistivity, location and inception angle. The high accuracy of the proposed algorithm is obtained even for faults close to the bus and low inception angles.


Introduction
Transmission-line protection is critical to power system security. High speed and selective isolation of faults in transmission lines help to improve the reliability of the power supply. Accuracy in locating permanent faults is also important.
The most widely used fault-detection algorithms calculate the fundamental frequency of current and voltage waveforms and then impedance is determined. Although this solution features a low computational cost and high reliability, it is time-consuming (at least one cycle). As lines are repaired by crews, a precise fault location minimizes repair time and benefits the security and quality of the supply.
The accuracy of the impedance method for fault location is affected by fault resistance, the mutual coupling effect and vague line-parameter knowledge [1]. Parameter-free fault-location algorithms provide more accurate results [2,3].
When a fault happens, transients are generated and propagate through the transmission line. Fault information is included in the waves traveling toward both ends. Furthermore, traveling wave-based methods can provide accurate fault location in shorter time frames [4,5], which also make them suitable for DC systems [6].
Acquiring the arrival time of the traveling wave is difficult, but it can be done accurately using wavelet transform (WT) [7,8]. These methods can be single-ended [9] or double-ended [10]. In the first case, data are independent from the communication equipment and also make it possible to avoid the complexities of synchronization between terminals [11].

WT with Multiresolution Analysis
The discrete wavelet transform (DWT) decomposes a sampled signal f (t) into a set of wavelets ψ: where the integers k and j are the translation and dilation factors, respectively, and DW j k the DWT coefficients in the wavelet basis considered.
Alternatively, multiresolution analysis (MRA) provides a different time-and-frequency resolution at each level of analysis. MRA employs two sets of functions (scaling functions ϕ(t) and wavelet functions ψ(t)), which allow the sampled signal f (t) to expand thus: where A j k is the approximation coefficients and D j k the detail coefficients, which can be calculated in terms of the scaling and wavelet basis: MRA builds a wavelet decomposition tree as shown in Figure 1. At each level, the signal is decomposed by applying successive pairs of low-pass (h(n)) and high-pass (g(n)) filters. After each filtering, half of the samples are eliminated and a level of decomposition is achieved at each end row; Energies 2019, 12, 3293 3 of 15 detail coefficients D1, D2, D3 and approximation coefficients A3 are obtained as a result. It allows both time and frequency information to be extracted simultaneously, even though the resolution differs; consequently, high time resolution is achieved for lower frequencies. Once the current fault transients are decomposed into their modal components (Section 5), the MRA of these signals allow to obtain WT coefficients of the aerial mode and zero mode current which are used to detect faults. For unsymmetrical faults to earth, zero-mode coefficients (i 0 ) are considered. In the other cases, aerial-mode currents (i α , i β ) are considered. WT is also used to locate the fault, except for single-line-ground (SLG) faults with low inception angle and when faults occur close to the relay; the proposed half-sine method is applied in these cases. transients are decomposed into their modal components (Section 5), the MRA of these signals allow to obtain WT coefficients of the aerial mode and zero mode current which are used to detect faults. For unsymmetrical faults to earth, zero-mode coefficients (i0) are considered. In the other cases, aerialmode currents (iα , iβ) are considered. WT is also used to locate the fault, except for single-line-ground (SLG) faults with low inception angle and when faults occur close to the relay; the proposed half-sine method is applied in these cases.

Half-Sine Method
The transient response to a disturbance in the transmission system is obtained by solving the differential equations governing signal propagation: where the column vectors [V(t)] and [i(t)] represent the voltages and currents along the conductors, and [R], [L], [G] and [C] represent the per-unit-length resistance, inductance, conductance, and capacitance, respectively. The solution to these differential equations is the sum of an incident and a reflected wave. In time domain analysis, the response wave R(t) is the convolution of the incident wave Vin(t) with the impulse response h(t), which can be expressed in matrix notation as follows where [R(t)] and [Vin(t)] are two sequences of n samples, and [h(t)] is a cyclic convolution matrix (nxn). In the proposed method, a half-sine (per phase) is injected into the transmission line under normal conditions, before the fault event. The half-sine Hs is defined as where u(t) is the unit step function and f the frequency. After fault detection, a new half-sine (Hs) is injected, which propagates along the line until it reaches the fault point, where it is reflected back.
The difference (∆Hs) between the signals measured when half-sines are injected in normal conditions (Hsn) and in fault conditions (HsF) is given by: The travel time to the fault point (Δt), and the propagation speed (v), taking the instants when half-sines are injected as a time reference, allow us to calculate the distance to the fault (d) easily.
To inject the half-sine, the coupling capacitors of a PLC (power-line communication) connection can be used as shown in Figure 2, where the injection is done by means of the fault locator device (Hsrelay). The coupling filter and the capacitance of the coupling capacitors are a broad band-pass filter matching the transmission-line surge impedance with the device impedance. To keep the bandwidth as wide as possible and to reduce the loss of energy on the transmission, the capacitance in use usually varies from 1 to 50 nF. PLCs are in common use, which makes the approach more applicable.

Half-Sine Method
The transient response to a disturbance in the transmission system is obtained by solving the differential equations governing signal propagation: where the column vectors [V(t)] and [i(t)] represent the voltages and currents along the conductors, and [R], [L], [G] and [C] represent the per-unit-length resistance, inductance, conductance, and capacitance, respectively. The solution to these differential equations is the sum of an incident and a reflected wave.
In time domain analysis, the response wave R(t) is the convolution of the incident wave V in (t) with the impulse response h(t), which can be expressed in matrix notation as follows where [R(t)] and [V in (t)] are two sequences of n samples, and [h(t)] is a cyclic convolution matrix (nxn). In the proposed method, a half-sine (per phase) is injected into the transmission line under normal conditions, before the fault event. The half-sine H s is defined as where u(t) is the unit step function and f the frequency. After fault detection, a new half-sine (H s ) is injected, which propagates along the line until it reaches the fault point, where it is reflected back.
The difference (∆H s ) between the signals measured when half-sines are injected in normal conditions (H sn ) and in fault conditions (H sF ) is given by: The travel time to the fault point (∆t), and the propagation speed (v), taking the instants when half-sines are injected as a time reference, allow us to calculate the distance to the fault (d) easily.
To inject the half-sine, the coupling capacitors of a PLC (power-line communication) connection can be used as shown in Figure 2, where the injection is done by means of the fault locator device (H s -relay). The coupling filter and the capacitance of the coupling capacitors are a broad band-pass filter matching the transmission-line surge impedance with the device impedance. To keep the bandwidth as wide as possible and to reduce the loss of energy on the transmission, the capacitance in use usually varies from 1 to 50 nF. PLCs are in common use, which makes the approach more applicable.

Modeled Power System
The power system considered here is the 132-kV system shown in Figure 3, which has been modeled in PSCAD/EMTDC [23]. The Hs-relay analysed in this study is installed at the line position B1L at G1. The transmission-line lengths are 100 km for the transmission line TL1 and 120 km for the transmission line TL2 and, in both cases, the frequency-dependent phase model is applied. Harmonic distortion in voltage waveform is also considered; the expected THD in voltage waveform (THDV) must be lower than 2.5% for normal operations, and for conditions lasting less than one hour, this limit may be exceeded by 50%. Therefore, the THDV is set to a value of 3.75%.  The propagation speed for aerial modes and for the earth mode is calculated by:

Modeled Power System
The power system considered here is the 132-kV system shown in Figure 3, which has been modeled in PSCAD/EMTDC [23]. The H s -relay analysed in this study is installed at the line position B1L at G1. The transmission-line lengths are 100 km for the transmission line TL1 and 120 km for the transmission line TL2 and, in both cases, the frequency-dependent phase model is applied. Harmonic distortion in voltage waveform is also considered; the expected THD in voltage waveform (THDV) must be lower than 2.5% for normal operations, and for conditions lasting less than one hour, this limit may be exceeded by 50%. Therefore, the THDV is set to a value of 3.75%.

Modeled Power System
The power system considered here is the 132-kV system shown in Figure 3, which has been modeled in PSCAD/EMTDC [23]. The Hs-relay analysed in this study is installed at the line position B1L at G1. The transmission-line lengths are 100 km for the transmission line TL1 and 120 km for the transmission line TL2 and, in both cases, the frequency-dependent phase model is applied. Harmonic distortion in voltage waveform is also considered; the expected THD in voltage waveform (THDV) must be lower than 2.5% for normal operations, and for conditions lasting less than one hour, this limit may be exceeded by 50%. Therefore, the THDV is set to a value of 3.75%.  The propagation speed for aerial modes and for the earth mode is calculated by: The electrical networks connected to the buses G1 and G2 are represented by their respective Thévenin equivalent voltage and impedance ( Figure 3). Figure 4 shows the transmission tower dimensions.
The propagation speed for aerial modes and for the earth mode is calculated by: where L 1 and C 1 are the positive sequence inductance and capacitance, and L 0 and C 0 are the zero sequence inductance and capacitance whose values are shown in Table 1.  The electrical networks connected to the buses G1 and G2 are represented by their respective Thévenin equivalent voltage and impedance ( Figure 3). Figure 4 shows the transmission tower dimensions. The propagation speed for aerial modes and for the earth mode is calculated by: These values of L 0 and C 0 , and the zero sequence speed v 0 given by Equation (11), have been calculated considering the ground as a perfect conductor, which is true for seawater; however, the value of ground resistivity will vary for different soil types and, while C 0 is constant, L 0 and, consequently, v 0 are affected by this parameter.

Fault Detection and Location Algorithm
A flowchart of the proposed fault detection and location algorithm is shown in Figure 5. where L1 and C1 are the positive sequence inductance and capacitance, and L0 and C0 are the zero sequence inductance and capacitance whose values are shown in Table 1.
These values of L0 and C0, and the zero sequence speed v0 given by Equation (11), have been calculated considering the ground as a perfect conductor, which is true for seawater; however, the value of ground resistivity will vary for different soil types and, while C0 is constant, L0 and, consequently, v0 are affected by this parameter.

Fault Detection and Location Algorithm
A flowchart of the proposed fault detection and location algorithm is shown in Figure 5.

Fault Detection
Fault inception gives rise to a variation of current modes, which allows detection schemes based on the WT to determine whether faults are present. Detection schemes based on WT use the aerialmode current. However, when the fault-inception angle is close to zero, the fault signal detected is weak and can be masked by noise. To avoid this drawback, fault detection can be improved by Data

Fault Detection
Fault inception gives rise to a variation of current modes, which allows detection schemes based on the WT to determine whether faults are present. Detection schemes based on WT use the aerial-mode current. However, when the fault-inception angle is close to zero, the fault signal detected is weak and can be masked by noise. To avoid this drawback, fault detection can be improved by considering zero-mode coefficients, which provide more relevant information and appreciable magnitudes with fault-inception angles close to the zero-crossing region. In a high-impedance grounded distribution network, the accuracy of this detection method would drop due to problems such as a weak fault signal, noise or unbalanced operation, and other techniques should, therefore, be considered [21]. However, this is not the case of (usually solidly) grounded transmission systems.
The current fault transients are decomposed into their modal components and, to improve fault detection, both aerial-and zero-mode coefficients are employed for fault detection. For unsymmetrical faults to earth, zero-mode coefficients are considered. In the other cases, aerial-mode currents are considered. Table 2 shows the current modes considered to detect the fault type (single-line-to-ground (SLG), line-to-line (LL), line-to-line-to-ground (LLG) or three-phase faults (3L/3LG)).

Wavelet Analysis
Modal transformation allows a three-phase coupled line to be decomposed into three independent propagation modes. Out of the various modal transformation methods, Clarke transformation is employed to calculate the modal components of currents: where i 0 is the zero or earth mode, and i α , i β are the aerial modes of the current. The MRA of a signal provides WT coefficients; the absolute local maximum values of these coefficients and the instants of occurrence (i.e., the arrival times of the traveling waves from the fault point to the substation) are obtained using wavelet modulus maxima (WMM). Figure 6 shows the result of applying WMM to the detail coefficients. where i0 is the zero or earth mode, and iα, iβ are the aerial modes of the current.
The MRA of a signal provides WT coefficients; the absolute local maximum values of these coefficients and the instants of occurrence (i.e., the arrival times of the traveling waves from the fault point to the substation) are obtained using wavelet modulus maxima (WMM). Figure 6 shows the result of applying WMM to the detail coefficients. The mother wavelet type, as well as the number of decomposition scales, may differ from one application or condition to another and, therefore, cannot be generalized to all cases. Several wavelet families with a variety of vanishing moments, such as Haar, Daubechies (db), Symlets, Coiflets, Biorthogonal, Reverse Biorthogonal and Discrete Meyer, were tested to select the most suitable mother wavelet. According to the results, and considering Daubechies family properties in faultlocation applications and their low computational cost, db3 wavelet was selected as the mother wavelet. The mother wavelet type, as well as the number of decomposition scales, may differ from one application or condition to another and, therefore, cannot be generalized to all cases. Several wavelet families with a variety of vanishing moments, such as Haar, Daubechies (db), Symlets, Coiflets, Biorthogonal, Reverse Biorthogonal and Discrete Meyer, were tested to select the most suitable mother wavelet. According to the results, and considering Daubechies family properties in fault-location applications and their low computational cost, db3 wavelet was selected as the mother wavelet.
Applying WMM allows us to clearly identify the wave peaks, and, therefore, fault distance is determined as follows: where d_wt is the estimated distance to fault, t 2 − t 1 is the time difference between two consecutive maximum values of WT coefficients and v is the propagation speed through the line or ground. The propagation speed in aerial modes is calculated from positive sequence parameters by means of Equation (10). For SLG faults with a low inception angle, aerial-mode coefficients are not effective. The zero mode is not used in the proposed algorithm because fault location using zero-mode coefficients leads to errors. Measuring zero-mode velocity with enough accuracy is difficult as ground resistivity varies along the line.
When SLG faults occur with a low inception angle, the algorithm applies the half-sine H' s method, as explained below.

H s and H' s Methods
The fault-locator device injects H s signals and registers the response of the network under normal conditions, allowing the fault locator to detect any irregularity in the transmission line.

H' s Method
When a fault is detected, H s are injected in each phase of the transmission line and recorded. Therefore, the algorithm calculates the signal ∆H s , given by Equation (9). The maximum of ∆H s can be used to determine the time ∆t HsF elapsed between the instant a H s signal is injected and the reception of the reflected signal when the fault has happened, taking the instants when H s are injected as a time reference. Finally, the distance to the fault (d HsF ) is obtained as follows For faults near to a substation, both high and low frequency H s allow for fault location. For long-distance faults and high-frequency H s , the signal is affected by attenuation. However, for low frequency H s there is also a clear signal deformation and a deviation of the maximum (Figure 7), which finally results in greater inaccuracies in the fault location.
For faults near to a substation, both high and low frequency Hs allow for fault location. For longdistance faults and high-frequency Hs, the signal is affected by attenuation. However, for low frequency Hs there is also a clear signal deformation and a deviation of the maximum (Figure 7), which finally results in greater inaccuracies in the fault location.
Signals are sampled with a frequency of 1.7 MHz. Hence, a 5 kHz Hs, which has a duration of 0.1 ms, is described by 170 samples. However, in a 100 kHz Hs (duration of 5 μs) only 8 points are sampled. As can be seen in Figure 8, compared with the 5 kHz Hs, the 100 kHz, Hs is of shorter duration and, consequently, the attenuation affects the 100 kHz Hs amplitude without causing a deviation of the signal maximum, as happens for a low frequency Hs. Signals are sampled with a frequency of 1.7 MHz. Hence, a 5 kHz H s , which has a duration of 0.1 ms, is described by 170 samples. However, in a 100 kHz H s (duration of 5 µs) only 8 points are sampled. As can be seen in Figure 8, compared with the 5 kHz H s , the 100 kHz, H s is of shorter duration and, consequently, the attenuation affects the 100 kHz H s amplitude without causing a deviation of the signal maximum, as happens for a low frequency H s .

H's Method
A preferable alternative to the detection of the maximum in the Hs method is to consider the time derivative of ∆Hs which can be expressed as the convolution of H's and the variation of the system impulse response (ΔnF) under fault compared with normal conditions

H' s Method
A preferable alternative to the detection of the maximum in the H s method is to consider the time derivative of ∆H s which can be expressed as the convolution of H' s and the variation of the system impulse response (∆ nF ) under fault compared with normal conditions where H' s starts at a maximum value making it easier to locate the distance to the fault (d H'sF ). Likewise, it is easier to inject a sine than a cosine, whose difficulty is equivalent to generating an ideal step. In Figure 9, H' s is compared with H s .  Table 3 shows a comparison of results obtained using the aforementioned methods for a longdistance fault (100 km, at the end of the line) and considering a constant soil resistivity of 200 Ω•m: It can be noticed that the H's method shows good results, especially for high-frequency halfsines. Therefore, to improve accuracy, the fault distance is evaluated using a 100-kHz half-sine, which is a frequency compatible with PLC equipment that can operate over a frequency band up to 500 kHz.   Table 3 shows a comparison of results obtained using the aforementioned methods for a long-distance fault (100 km, at the end of the line) and considering a constant soil resistivity of 200 Ω·m:  It can be noticed that the H' s method shows good results, especially for high-frequency half-sines. Therefore, to improve accuracy, the fault distance is evaluated using a 100-kHz half-sine, which is a frequency compatible with PLC equipment that can operate over a frequency band up to 500 kHz.

Discussion
Although WT methods offer high accuracy in fault location, a drawback associated with this method is that faults close to the substation cannot be located due to the limitation of sampling frequency.
To ensure an accurate fault location, the time difference between two consecutive peaks of the WT coefficients is considered to be, at least, twice the period of the applied WT. Therefore, in the case of this fault locator (current signals are sampled at 1.7 MHz), 2 km is the threshold distance (Blind_d) from which the WT method gives accurate fault-location estimation, and below this distance there is a blind zone. To overcome this shortcoming, the sampling frequency could be increased, but that would make it too expensive. The proposed algorithm tackles this drawback by applying the H' s method for calculating the fault distance with single-phase-to-ground faults and with low inception-angle faults, achieving good accuracy in both cases.
In overhead lines, signals lose energy as they are transmitted. Attenuation per unit length increases with frequency due to the skin effect, but it is also affected by ground resistivity (ρ) [4]. In propagation in overhead lines with ground return, the series impedance per unit length is the sum of the series impedances considering ground as a perfect conductor and of the ground impedance integral [24], which is evaluated in PSCAD through direct numerical integration. Alternatively, the ground impedance integral can be estimated by an analytical approximation [25,26]. If ground is a perfect conductor (ρ ≈ 0 Ω·m), the attenuation is negligible. For low values of ρ (<100 Ω·m), the attenuation rises rapidly with increasing ρ, and, from this point (ρ = 100 Ω·m), it remains almost constant as shown in Figure 10. The ground effect return affects the injected half-sines as the length of the line increases; therefore, it will determine the half-sine amplitude.
Energies 2019, 12, x FOR PEER REVIEW 10 of 16 a blind zone. To overcome this shortcoming, the sampling frequency could be increased, but that would make it too expensive. The proposed algorithm tackles this drawback by applying the H's method for calculating the fault distance with single-phase-to-ground faults and with low inceptionangle faults, achieving good accuracy in both cases. In overhead lines, signals lose energy as they are transmitted. Attenuation per unit length increases with frequency due to the skin effect, but it is also affected by ground resistivity (ρ) [4]. In propagation in overhead lines with ground return, the series impedance per unit length is the sum of the series impedances considering ground as a perfect conductor and of the ground impedance integral [24], which is evaluated in PSCAD through direct numerical integration. Alternatively, the ground impedance integral can be estimated by an analytical approximation [25,26]. If ground is a perfect conductor (ρ ≈ 0 Ω•m), the attenuation is negligible. For low values of ρ (<100 Ω•m), the attenuation rises rapidly with increasing ρ, and, from this point (ρ = 100 Ω•m), it remains almost constant as shown in Figure 10. The ground effect return affects the injected half-sines as the length of the line increases; therefore, it will determine the half-sine amplitude.

Simulation Results
The algorithm performance is evaluated by varying fault resistance, location, and fault type on the power system of Figure 3. Figure 11 and 12 show the error in fault distance estimation for single-

Simulation Results
The algorithm performance is evaluated by varying fault resistance, location, and fault type on the power system of Figure 3. Figures 11 and 12 show the error in fault distance estimation for single-line-ground (SLG), double-line-to-ground (LLG), line-to-line (LL) and three-phase (3L) faults, varying fault location (1,8,20,30,40,50, 60 and 100 km). It should be noted that these errors are almost insensitive to different fault types, being slightly higher for SLG faults. Figure 13 depicts the error in fault distance estimation for different fault locations (1,8,20,40,50, 60 and 100 km) and fault resistance values (between 0 Ω and 150 Ω), in case of SLG faults.          In Table 5 the proposed method is compared with several methods: a traditional distance algorithm (DFFT distance algorithm) [27], a distance relaying scheme based on the current phase For faults close to the bus, below the threshold distance Blind_d and for the sampling frequency considered, the error (Table 4) is lower than 133 m (0.13%) and independent of the angle of inception, fault resistance or ground resistivity. For longer distances, the maximum relative error detected is 0.14%, which is obtained for a fault resistance value of 150 Ω when an SLG fault is applied. These results were obtained considering the worst harmonic distortion and inception angle (0 • ) conditions.  Figure 14 shows how the error for different location faults is not affected by the inception angle which is varied between 0 • and 270 • . Additionally, these results are not altered when ground resistivity values vary (between 0 and 10,000 Ω·m), as can be noticed in Figure 15. Ground resistivity just affects waveform attenuation and, for sufficient high-frequency half-sine, the impact is negligible.     In Table 5 the proposed method is compared with several methods: a traditional distance algorithm (DFFT distance algorithm) [27], a distance relaying scheme based on the current phase In Table 5 the proposed method is compared with several methods: a traditional distance algorithm (DFFT distance algorithm) [27], a distance relaying scheme based on the current phase jump behaviour during fault conditions to improve the apparent impedance estimated (DFFT distance algorithm + JCF) [27], a fault location method based on time difference of arrival of ground-mode and aerial-mode traveling waves [4], a wavelet-based method [20] and the proposed method. The comparison shows that the proposed method has better accuracy with respect to the traditional methods described in [20] and [27]. The proposed algorithm accuracy is obtained even for faults close to the bus (<2%), low inception angles and it is independent of the soil resistivity value. when a close fault happens in the adjacent line TL2 (Figure 3) and measuring errors are considered, there could be a small probability of incorrect operation. A time-delayed zone could be defined to solve this problem, but the protection should, nevertheless, be fast, and any intentional time delay should be avoided. There is a better alternative, the reception of the H s signal from the remote end can be used to block the trip. The proposed scheme will provide selective and fast clearance of faults. The response time of the proposed algorithm to detect faults is lower than 0.5 ms.

Conclusions
An accurate fault-location scheme for transmission lines is presented based on single-ended measurements.
Faults are detected by a wavelet multi-resolution analysis of the current transients. Aerial current modes are used, and zero modes are considered to improve fault detection when fault-inception angles are close to zero.
Fault location is also based on the wavelet analysis, except for faults close to the bus and fault-inception angles near to zero for SLG faults. In such cases, the fault locator employs the H' s method to improve accuracy. In the proposed H' s method, the fault distance is calculated by measuring the difference between the time when a 100-kHz half-sine signal is sent and the time when the derivative response is received.
The proposed algorithm was checked considering harmonic distortion and varying fault resistance, ground resistivity, location and inception angle. According to these results, the distance location can be determined with an accuracy lower than 150 m in all the cases. Funding: This research was funded by Spanish Economy, Industry and Competitiveness Ministry and by European Regional Development Fund (ERDF), grant number RTC-2016-5234-3, in the frame of the Project MHiRED "New technologies for wind-photovoltaic hybrid mini networks managed with storage in connection to the grid and with synchronous support in off-grid operation".