Multi-Source Partial Discharge Fault Location with Comprehensive Arrival Time Difference Extraction Method and Multi-Data Dynamic Weighting Algorithm

Abstract

The location of the partial discharge source is an important part of fault diagnosis inside power equipment. As a key step of the ultra-high frequency location method, the extraction of the time difference of arrival can generate large errors due to interference. To achieve accurate time difference extraction and further multi-source partial discharge location, a location method with comprehensive time difference extraction and a multi-data dynamic weighting algorithm is proposed. For time difference extraction, the optimized energy accumulation curve method applies wavelet transform and mode maximization calculations such that it overcomes the effect of interference signals before the wave peak. The secondary correlation method improves the interference capability by performing two rounds of correlation calculations. Both extraction methods are combined to reduce the error in time difference extraction. Then, the dynamic weighting algorithm effectively utilizes multiple data and improves the location accuracy. Experimental results on multi-source partial discharge locations performed in a transformer tank validate the accuracy of the proposed method.

1. Introduction

Partial discharge (PD) is a common fault in the internal insulation of power equipment. The PD source location is an important part of PD diagnosis. Conventional PD source location methods include acoustic, optical, electrical and ultra-high frequency (UHF) methods [1,2,3,4,5]. Because of its high sensitivity and interference suppression, the UHF method has been widely studied and applied. The accuracy of the time difference of arrival (TDOA) extraction is an important factor in the location effect of the UHF method [6]. TDOA extraction methods for UHF signals mainly include the first peak method, the wavelet transform method, the cross-correlation method and the energy accumulation method [7,8,9].

The first peak method uses the starting point or some threshold value of the first pulse peak as a reference point to calculate the TDOA [10,11]. When the first peak protrudes, this method is convenient with high accuracy [7]. Unfortunately, the interference signal before the first peak can affect the determination of the reference point. The wavelet transform method is performed by applying a wavelet transform to the UHF signal and determining the mode maximum point as the reference point [12,13]. However, interference-generated pulse peaks can easily lead to a shift of the mode maximum point, and a large TDOA error may still exist. The cross-correlation method establishes the cross-correlation function of two UHF partial discharge signals [14,15]. Then, when the cross-correlation function reaches its peak value, the time-independent variable corresponding to the peak point is equal to the TDOA. If there is an interference signal similar to the discharge pulse in the UHF signal, then the cross-correlation function may have multiple peaks, resulting in difficulty in determining the TDOA. The energy accumulation method reduces the effect of interference signals by converting UHF partial discharge signals into energy accumulation curves. Existing studies take the starting point of the slope of the energy accumulation curve as the reference point to calculate TDOA [8,16]. Although the energy accumulation method is the most robust to noise from single PD sources [8], it is still susceptible to large interference signals before the first peak.

In addition, the multi-source partial discharge (MSPD) location is more difficult than that of a single PD source. When there are multiple PD sources in space-limited power equipment, the UHF signals generated by different PD sources interfere with each other. At the same time, signals propagating in other directions can also begin to interfere after reflecting on structures such as metal tank walls. To distinguish UHF signals from different PD sources, an additional acoustic measurement is adopted [10,17], which increases measurement complexity. In [18], the MSPD location has been realized based on the simulation model and the TDOA database. However, it is not convenient to use because it only applies to specific equipment.

To address the issues described, this paper proposes a PD location method for MSPD inside power equipment by integrating two TDOA extraction methods and using multi-data dynamic weighting. The proposed PD fault location method has the following salient characteristics:

(1)

The optimized energy accumulation method and the secondary correlation method are proposed and combined to reduce the error in TDOA extraction. In addition, dynamic weighting algorithms are used to effectively utilize multiple groups of TDOA data in PD source calculations and improve the location accuracy. Compared with the other methods, the localization method has better interference and accuracy.

(2)

Compared to the conventional method, the optimized energy accumulation curve method obtains the inflection point of the energy curve as a reference point through wavelet transform and mode maximization calculations, which overcomes the effect of the interference signal before the wave peak.

(3)

In the secondary correlation method, the effect of interfering signals is mitigated by the two rounds of correlation calculations. The interference capability of TDOA extraction is further improved compared to the cross-correlation method.

The remainder of this paper is organized as follows. Two TDOA extraction methods are described in Section 2. Practical applications of these two methods to TDOA extraction of the MSPD are described, and their comprehensive application strategy is provided in Section 3. Following, the dynamic weighting algorithm for multiple TDOA data is proposed and the experimental results of the MSPD location are analyzed in Section 4. Finally, conclusions are drawn in Section 5.

2. TDOA Extraction Methods

In practice, the TDOA is the time difference between different detection sensors receiving the same UHF partial discharge signal. The PD source position can be calculated using three TDOAs. In the calculation, the UHF signals are considered to be transmitted along a straight line. Define P (x, y, z) as the PD source coordinates and S_i (x_i, y_i, z_i) as the i-th sensor coordinates. Taking S₁ as the reference, the propagation distance difference Δr_i1 is:

$$\Delta r_{i1}=\sqrt{{\left(x-x_i\right)}^2+{\left(y-y_i\right)}^2+{\left(z-z_i\right)}^2}-\sqrt{{\left(x-x_1\right)}^2+{\left(y-y_1\right)}^2+{\left(z-z_1\right)}^2}=c\Delta t_{i1}$$

(1)

where c denotes the velocity of the electromagnetic wave; Δt_i1 denotes the TDOA between the first and the i-th sensors. The coordinates of S_i (x_i, y_i, z_i) and c are known. When more than three Δt_i1 are achieved, the coordinate of P (x, y, z) can be obtained by solving the equation groups.

To suppress the periodic interference or white noises, the wavelet denoising and filtering are performed on the collected UHF signals before the TDOA extraction [19,20].

2.1. Optimized Energy Accumulation Method

The energy accumulation curve E is calculated based on the UHF partial discharge signal:

$$E={\displaystyle \sum _{k=0}^n{u}_k{}^2}$$

(2)

where u_k denotes the voltage value of the k-th point in the PD signal segment, and n denotes the sampling point number.

Two channels of the UHF signals and corresponding energy accumulation curves are shown in Figure 1. The PD pulse period in Figure 1a corresponds to the rising segment of the energy accumulation curve in Figure 1b.

The traditional energy accumulation method takes the starting point of the curve slope as the reference point to calculate the TDOA [12]. This is easily affected by interference signals before the first wave peak.

In this paper, the wavelet modulus maximum calculation is employed to detect singularities so that the inflection point of the slope of the curve can be determined. The inflection point of the curve slope is then used as a reference point to extract the TDOA.

Wavelet transform is applied firstly to smooth the energy accumulation curve f(x) to reduce the impact of interference signals and noises. The cubic B-spline function θ(x) is applied in this paper, which can be regarded as the response function of the low-pass filter and can be approximated as:

$$\widehat{\theta }\left(x\right)={\left(\frac{sin\left(\frac{x}{4}\right)}{\frac{x}{4}}\right)}^4$$

(3)

Define ${\phi }^{\prime }(x)=d\theta (x)/dx$, and then perform the convolution calculation with f(x) and ${\phi }^{\prime }(x)$. The calculated wavelet transform function is $W^{\prime }{f}_s\left(x\right)$, which can be expressed as:

$$W^{\prime }{f}_s\left(x\right)=f\left(x\right)\otimes {\phi }^{\prime }\left(x\right)=f\left(x\right)\otimes \left(\frac{d\theta \left(x\right)}{dx}\right)=\frac{d}{dx}\left(f\otimes \theta \right)\left(x\right),x\in s$$

(4)

where $\otimes$ denotes convolution calculation, and s denotes the transform scale. $W^{\prime }{f}_s\left(x\right)$ is proportional to the first derivative of smoothed f(x).

Then the inflection point of the slope of the energy accumulation curve can be determined by calculating the partial modulus maximum of the absolute value of $W^{\prime }{f}_s\left(x\right)$. Figure 2 shows the wavelet modulus maximum lines corresponding to the energy accumulation curves in Figure 1. It can be seen that the modulus maximum line for the first wave peak of the PD pulse is significantly higher than the lines for the subsequent wave peaks and noises. Therefore, the time point corresponding to the modulus maximum line with the highest amplitude can be used as a reference point to calculate the TDOA. The optimized energy accumulation method overcomes the effect of the interference signal before the wave peak and achieves high accuracy when the PD pulse is prominent.

2.2. Secondary Correlation Method

The UHF partial discharge signals of two sensors from the same PD source are recorded as s₁(t) and s₂(t), respectively, which are expressed as:

$$\{\begin{array}{l}{s}_1\left(t\right)=a_1\cdot s\left(t-t_1\right)+n_1\left(t\right)\\ s_2\left(t\right)=a_2\cdot s\left(t-t_2\right)+n_2\left(t\right)\end{array}$$

(5)

where s(t) denotes signal of the PD source; a₁ and a₂ denote attenuation constants, respectively; t₁ and t₂ denote time differences, Δt = t₂ − t₁ denotes the TDOA of the two signals; n₁(t) and n₂(t) are signal distortions and noises, respectively.

Employ the self-correlation and the cross-correlation calculation to signals s₁(t) and s₂(t). The self-correlation function $R_{s_1{s}_1}\left(\tau \right)$ and cross-correlation function $R_{s_1{s}_2}\left(\tau \right)$ are:

$$R_{s_1{s}_1}\left(\tau \right)={\displaystyle {\int }_{-\infty }^{+\infty }{s}_1\left(t\right)}\cdot s_1\left(t+\tau \right)dt=E\left[s_1\left(t\right)s_1\left(t+\tau \right)\right]$$

(6)

$$R_{s_1{s}_2}\left(\tau \right)={\displaystyle {\int }_{-\infty }^{+\infty }{s}_1\left(t\right)}\cdot s_2\left(t+\tau \right)dt=E\left[s_1\left(t\right)s_2\left(t+\tau \right)\right]$$

(7)

The traditional cross-correlation method extracts the TDOA by calculating the peak value of the cross-correlation function. According to Equations (5) and (7), the cross-correlation function $R_{s_1{s}_2}\left(\tau \right)$ is calculated as:

$$\begin{array}{ll}\hfill R_{s_1{s}_2}\left(\tau \right)& =E\left[s_1\left(t\right)s_2\left(t+\tau \right)\right]\hfill \\ & =E\left\{\left[a_{1}s\left(t-t_1\right)+n_1\left(t\right)\right]\cdot \left[a_{2}s\left(t-t_2+\tau \right)+n_2\left(t+\tau \right)\right]\right\}\hfill \\ & =a_1{a}_2{R}_{ss}\left(\tau -\mathsf{\Delta }t\right)+a_2{R}_{n_{1}s}\left(\tau -t_2\right)+a_1{R}_{s{n}_2}\left(\tau +t_1\right)+R_{n_1{n}_2}\left(\tau \right)\hfill \end{array}$$

(8)

where $R_{n_1{n}_2}(\tau ) \approx 0$. If s(t) is unrelated to n₁(t) and n₂(t), then:

$$R_{s_1{s}_2}\left(\tau \right)=a_1{a}_2{R}_{ss}\left(\tau -\mathsf{\Delta }t\right)$$

(9)

In Equation (8), when τ = Δt, R_ss reaches its maximum, so does $R_{s_1{s}_2}\left(\tau \right)$. Thus, Δt can be determined by:

$$\mathsf{\Delta }t=\mathit{arg}\left\{max\left[R_{ss}\left(\tau -\mathsf{\Delta }t\right)\right]\right\}=\mathit{arg}\left\{max\left[R_{s_1{s}_2}\left(\tau \right)\right]\right\}$$

(10)

where arg() denotes taking independent variables.

In practice, the cross-correlation function may have multiple peaks due to attenuation, oscillation of the UHF signal and the effect of interference. That is, $a_2{R}_{n_{1}s}(\tau -t_2)+a_1{R}_{{\mathit{sn}}_2}(\tau +t_1)$ is large. In this case, the first peak of the PD pulse is not prominent, making the reference point difficult to determine.

To improve the resolution ratio of the first peak and reduce the effect of interference, a secondary correlation method is proposed.

Define the reference signal s₁(t) and perform self-correlation calculation to s₁(t). The self-correlation function $R_{s_1{s}_1}\left(\tau \right)$ is expressed as:

$$\begin{array}{ll}\hfill R_{s_1{s}_1}\left(\tau \right)& =E\left[x\left(t\right)x\left(t+\tau \right)\right]\hfill \\ & =E\left\{\left[a_{1}s\left(t-t_1\right)+n_1\left(t\right)\right]\cdot \left[a_{1}s\left(t-t_1+\tau \right)+n_1\left(t+\tau \right)\right]\right\}\hfill \\ & =a_1{}^2{R}_{ss}\left(\tau \right)+a_1\left[R_{s{n}_1}\left(\tau +t_1\right)+R_{n_{1}s}\left(\tau -t_1\right)\right]+R_{n_1}{}_{n_1}\left(\tau \right)\hfill \end{array}$$

(11)

Then employ cross-correlation calculation with reference signals s₁(t) and s₂(t), as in Equation (8). Since $R_{s_1{s}_1}\left(\tau \right)$ and $R_{s_1{s}_2}(\tau )$ are functions with τ as the independent variable, τ can be substituted by t.

Then, perform the cross-correlation calculation with $R_{s_1{s}_1}\left(t\right)$ and $R_{s_1{s}_2}\left(t\right)$ and obtain the secondary cross-correlation function $R_{R_1{R}_2}\left(\tau \right)$, in which the secondary correlation calculation of the first-round correlation functions (such as $R_{\mathit{ss}}\left(t\right)$ and $R_{s{n}_1}\left(t\right)$) will appear. In $R_{s_1{s}_1}\left(t\right)$ and $R_{s_1{s}_2}\left(t\right)$, according to Equations (8) and (11), $R_{n_1{n}_1}(t)$ are unrelated to other correlation functions. If $R_{\mathit{ss}}\left(t\right)$, $R_{s{n}_1}\left(t\right)$ and $R_{s{n}_2}\left(t\right)$ are unrelated, then the cross-correlation functions (such as $E\left[R_{\mathit{ss}}\right(t-∆t\left)R_{s{n}_1}\right(t+t_1\left)\right]$ and $E\left[R_{\mathit{ss}}\right(t\left)R_{s{n}_2}\right(t+t_1\left)\right]$) are equal to zero. Thus, the secondary cross-correlation function $R_{R_1{R}_2}\left(\tau \right)$ can be simplified as:

$$\begin{array}{ll}\hfill R_{R_1{R}_2}\left(\tau \right)& =E\left[R_{s_1{s}_1}\left(t\right)R_{s_1{s}_2}\left(t\right)\right]\hfill \\ & =a_1{}^3{a}_{2}E\left[R_{ss}\left(t\right)R_{ss}\left(t-\mathsf{\Delta }t+\tau \right)\right]+a_1{}^2{a}_{2}E\left[R_{ss}\left(t\right)R_{n_{1}s}\left(t-t_2+\tau \right)\right]+\cdots \cdots \hfill \\ & =a_1{}^3{a}_{2}E\left[R_{ss}\left(t\right)R_{ss}\left(t-\mathsf{\Delta }t+\tau \right)\right]=a_1{}^3{a}_2{R}_{ss,ss}\left(\tau -\mathsf{\Delta }t\right)\hfill \end{array}$$

(12)

where R_ss,ss(τ) denotes the self-correlation function of R_ss(t). Then, Δt can be expressed as:

$$\mathsf{\Delta }t=\mathit{arg}\left\{max\left[R_{ss,ss}\left(\tau -\mathsf{\Delta }t\right)\right]\right\}=\mathit{arg}\left\{max\left[R_{R_1{R}_2}\left(\tau \right)\right]\right\}$$

(13)

In the secondary correlation method, the weight of original signal s(t) is strengthened, while the influence of noises n₁(t) and n₂(t) are weakened in the first round of correlation calculations. The effect of $R_{s{n}_1}\left(\tau \right)$ and $R_{s{n}_2}\left(\tau \right)$ are weakened in the second round of the correlation calculation.

2.3. Preliminary Comparison of TDOA Extraction Method

A preliminary comparison of the extraction effects of different extraction methods is presented. Experiments are carried out in laboratory conditions and in air, respectively. The PD sources are modeled by a needle-plate PD model. The UHF sensors adopt micro-strip antennas, whose absolute bandwidth is 340~440 MHz when the standing wave ratio (SWR) is lower than 2, and the comparative bandwidth reaches 25.6% [21]. The applied oscilloscope is the LeCroy7100, which features a simulation bandwidth of 1 GHz and a maximum sampling rate of 20 GS/s.

In the experiments, the distance between the reference sensor and the PD source is set as 1 m. Different values are used for the distance difference between the sensors. Multiple repeated experiments are performed under each set of distance differences. The TDOAs are extracted separately using the first peak method, the secondary correlation method and optimized energy accumulation method, respectively. The mean values of the TDOAs and their errors compared to theoretical values are calculated. As the results in

Table 1. TDOA extraction methods comparison.

	Error (ns)	First Peak Method	Secondary Correlation Method	Optimized Energy Accumulation Method
Theoretical TDOA (ns)
1.6		0.16	0.47	0.19
7.3		0.68	0.48	0.45
11.0		1.35	0.30	0.17
16.8		1.26	0.52	0.43

show, in the presence of air and without significant interference, the error of the optimized energy accumulation method and the secondary correlation method are, in general, smaller than that of the first peak method. Moreover, the optimized energy accumulation method has the highest accuracy.

3. MSPD TDOA Extraction

The PD locations are usually used to diagnose PD faults in power equipment such as transformers. Within a piece of power equipment, electromagnetic waves reflect off tank walls and other structures and form interference sources. Meanwhile, in the situation of multiple localized discharge sources, the signals of the discharge pulses interfere with each other. To further study the practical application of the proposed TDOA extraction method and put forward the optimal application strategy, an experimental platform is built on the basis of a real transformer tank (5.2 × 2.6 × 1.8 m). The experimental platform diagram is shown in Figure 3, and the transformer tank is shown in Figure 4a.

The needle-plate PD models are used for the PD sources to simulate PD caused by metal protrusions in a transformer, as shown in Figure 4b. The sensors use the UHF Archimedean spiral antennas [22], as shown in Figure 4c. Sensors are mounted on different surfaces of the tank. S_i denotes the i-th sensor. Taking the corner of the tank as the origin, a virtual three-dimensional coordinate system is established. Consequently, the sensor coordinates are S₁ (0.2 m, 0, 0.2 m), S₂ (5.2 m, 1.8 m, 1 m), S₃ (0, 1.3 m, 0.1 m) and S₄ (5 m, 2.6 m, 1.8 m), respectively. The PD source coordinates are P₁ (0.2 m, 0.3 m, 0.2 m) and P₂ (5 m, 2.4 m, 1.8 m), respectively.

3.1. Application of Optimized Energy Accumulation Method

The application of the optimized energy accumulation method to the MSPD case is based on time window scanning. In this paper, the width of the main PD pulse (defined by 30% of the peak value) is about 13 ns, as shown in Figure 5, where l₁, l₂ and l₃ denote the length, width and height of the tank, respectively. The theoretical maximum TDOA of the signals from the same PD source is $∆t_{\max }=\sqrt{{l_1}^2+{l_2}^2+{l_3}^2}/c$ ≈ 20 ns, where c denotes the propagation velocity of the electromagnetic wave in air and its approximate value is 30 cm/ns.

Based on the main pulse width and the theoretical maximum TDOA, the time window and scanning time step are set as 40 and 10 ns, respectively. After a discharge, the synchronization signals collected by the four sensors are shown in Figure 6.

First, for the UHF signal in a time window of 0–40 ns, the energy accumulation curves are conversed. The wavelet transform and the modulus maximization are performed. If sufficiently high modulus maximum lines do not appear simultaneously in the time window of the signals, then the time window is advanced by the scanning step of 10 ns. The above procedure is repeated for the time window of 10–50 ns. When sufficiently high modulus maximum lines appear simultaneously on all signals, record the corresponding time points.

Take a set of multi-channel synchronized signals as an example. For time window 2 shown in Figure 6, sufficiently high modulus maximum lines of all signals appear simultaneously in the 40 ns time window, as shown in Figure 7b. Since interference such as noise is weak and the main PD pulse is well-defined, each energy accumulation has only one inflection point, and the TDOAs can be computed from the corresponding time points of the modulus maximum lines.

However, for time window 1 in Figure 6, when the interferences are heavy and multiple peaks occur in a time window, several energy accumulation curves have more than one inflection point, as shown in Figure 8, leading to hard determination of the modulus maximum value lines.

Since the interference is accidental, the TDOA data can still be obtained by taking multiple sets of repeated experimental data and summarizing the calculated results. However, the TDOA extraction is less efficient. At the same time, this approach still fails to cope with the extreme case that multiple PD main pulses are exactly within the same time window.

3.2. Application of Secondary Correlation Method

The application of the secondary correlation method in the MSPD problem is to extract the feature signal segment from the reference PD signal and then conduct the secondary correlation with the signals of other channels.

Take the PD signals of channels 1 and 4 shown in Figure 9 as examples. The specific steps of the signal feature segment and the TDOA extraction are described as follows:

(1)

Set the signal of channel 1 as the reference, which is marked as s₁(t). First, the highest value of s₁(t) is searched for and marked as s₁(t_m1). Then, search for the nearest points whose amplitude equals 30%·s₁(t_m1), and the corresponding time coordinates are marked as t_s1, t_e1. The fragment of s₁(t) in the time interval [t_s1, t_e1] is defined as the first reference fragment s_1r1(t).

(2)

Signal of channel 4 is the contrast signal, which is marked as s₄(t). From s₄(t), the contrast fragment s_4c1(t) corresponding to the reference fragment s_1r1(t) is selected. Considering the theoretical maximum time difference Δt_max, the time interval of s_4c1(t) should be [t_s1 − Δt_max, t_e1 + Δt_max]. Then, based on Equations (11)–(13), the secondary function $R_{R_1{R}_2}\left(\tau \right)$ of f₁(t) and T₁(t) is deduced. The function curve of $R_{R_1{R}_2}\left(\tau \right)$ is shown in Figure 9b. The first TDOA is then extracted.

(3)

After completing the first TDOA extraction, set the function value of [t_s1, t_e1] in s₁(t) to 0 and obtain a new reference signal s_1′(t). Then, repeat the above steps to extract the new reference fragment s_1r2(t) from s_1′(t) and the contrast fragment s_4c2(t) from s₄(t), as shown in Figure 10. Then, the second TDOA is extracted. The TDOA extracting process above will be repeated until s₁(t_mi) is lower than the threshold s_t = 30%·s₁(t_m1).

Compared to the optimized energy accumulation method, the TDOA extracted by the secondary correlation method is less affected by the interference signal. The TDOAs can also be extracted when the amplitude of the interference pulse is close to the first peak of the PD pulse.

3.3. Method Comparison and Application Strategy

According to the above studies, both the optimized energy accumulation method and the secondary correlation method still have weaknesses. The accuracy of the optimized energy accumulation method decreases when the amplitude of the interference pulse approaches that of the first peak. The accuracy of the secondary correlation method will still be affected by the waveform distortion of the discharge pulse.

To compare the differences between the above two extraction methods more intuitively and to make full use of them to form a complement, they are both used to extract TDOAs. Several TDOA calculation data are listed, as shown in

Table 2. TDOAs obtained by the optimized energy accumulation method.

	TDOA (ns)	P1 (0.2 m, 0.3 m, 0.2 m)			P2 (5 m, 2.4 m, 1.8 m)
Group		Δt12	Δt13	Δt14	Δt12	Δt13	Δt14
Theoretical value		15.3	0.7	18.0	−16.6	−2.4	−17.5
1		15.3	0.6	18.1	−18.0	−3.8	−19.5
2		15.3	0.7	18.1	−17.8	−3.2	−18.7
3		-	-	-	−15.7	−2.1	−16.6
4		-	-	-	−17.6	−3.4	−18.6
5		15.7	0.7	18.5	−16.9	−2.7	−17.8
6		15.3	0.5	18.1	-	-	-
7		14.8	0.6	17.0	−16.8	−2.0	−18.7
8		-	-	-	−17.4	−3.5	−18.8
9		15.4	0.7	18.2	−17.1	−2.9	−18
10		14.8	0.4	16.7	−18.0	−3.2	−19.2
11		14.6	0.5	16.5	−18.0	−3.7	−17.2
12		-	-	-	−17.2	−2.2	−18.8
Average absolute error		0.28	0.11	0.60	0.88	0.74	1.07

and

Table 3. TDOAs obtained by the secondary correlation method.

	TDOA (ns)	P1 (0.2 m, 0.3 m, 0.2 m)			P2 (5 m, 2.4 m, 1.8 m)
Group		Δt12	Δt13	Δt14	Δt12	Δt13	Δt14
Theoretical value		15.3	0.7	18.0	−16.6	−2.4	−17.5
1		15.0	0.6	19.5	−16.2	−2.4	−17.0
2		15.3	0.7	18.0	−15.7	−1.5	−16.6
3		15.7	0.5	19.0	−17.8	−3.2	−18.8
4		15.0	0.1	17.2	−17.0	−2.1	−18.0
5		14.3	0.5	17.1	−18.0	−2.4	−18.9
6		14.8	0.3	17.0	−15.1	−1.7	−18.3
7		15.4	0.7	18.5	−16.0	−2.1	−18.9
8		14.7	0.5	16.9	−17.5	−3.0	−18.9
9		14.6	0.5	16.5	−17.0	−2.9	−17.0
10		15.0	0.6	18.5	−18.0	−3.6	−18.0
Average absolute error		0.42	0.21	0.88	0.91	0.53	0.99

.

According to Table 2 and Table 3, both the optimized energy accumulation method and the secondary correlation method are able to perform relatively accurate TDOA extraction under the conditions of internal equipment and multiple PD sources. For the optimized energy accumulation method, a part of the TDOA is not exactly solved, as shown in the blank column of Table 2.

By comparing the TDOA extraction errors of the two methods, it can be seen that when the distance between reference sensor S₁ (0.2 m, 0, 0.2 m) and PD source P₁ (0.2 m, 0.3 m, 0.2 m) is close, the accuracy of the optimized energy accumulation method is relatively higher. The possible reason is that, in this case, the signal from the PD source P₁ has a large amplitude when it arrives at the reference sensor S₁ and can be clearly distinguished from the interference signal. Therefore, the arrival time of the signal at the reference node T₁ is relatively accurate, and the TDOA extraction accuracy is relatively high.

When the distance between the reference node S₁ (0.2 m, 0, 0.2 m) and the PD source P₂ (5 m, 2.4 m, 1.8 m) is far away, propagation of the pulse signal is relatively more disturbed. In this case, the secondary correlation method is relatively more robust against interference, and hence, the extracted TDOA is more accurate.

Based on the above analysis and the statistics of the experimental results, the strategy to apply these two TDOA extraction methods is described as follows.

(a)

If the signal–noise ratio (SNR) of the PD signal is low, the secondary correlation method is used.

(b)

If the distance between the reference sensor and PD source is estimated to be shorter than 1 m, and SNR of the PD signal is not very low, the optimized energy accumulation method is used.

(c)

Otherwise, the results of the two methods should be compared.

In addition, the TDOA data with large errors can be seen in Table 2 and Table 3. Therefore, the PD sources localized by only one set of TDOAs may not be accurate.

4. Multi-Data Dynamic Weighting Algorithm

Multi-data statistics help to reduce the effect of accidental TDOA data and improve location accuracy [23]. In this paper, a dynamic weighting algorithm is proposed to improve the fault location accuracy for the MSPD by fully exploiting the multiple sets of the TDOAs obtained. First, the point density estimation is used to further eliminate errors due to interference. The initial points are then linearly classified into multisets. Finally, the points in each set are dynamically weighted to compute the coordinates of multiple PD sources.

4.1. Point Density Estimation

The PD source is defined as the primary source point. The initial source points are distributed in a 3D coordinate system based on the actual spatial size. Their densities are higher in regions close to the true PD sources. The point density estimation method is used to count the number of points in a sphere centered at a point. The density λ(p) of the initial source point p(x, y, z) is given by:

$$\lambda \left(p\right)=\frac{N(p,r)}{4\pi r^3/3}$$

(14)

where N(p, r) denotes the number of points involved inside the sphere, whose center is the point p and the radius is r; 4πr³/3 denotes the volume of the sphere. The density threshold value λ₀ will be defined according to the number of initial source points. The error initial source points (λ(p) ≤ λ₀) will be removed first, such as circuited point k in Figure 11a.

4.2. Linear Classification

In addition to the error points, the number of initial source points N₀ is counted. Then, the point with the highest λ(p) is searched for and defined as p₁ (x₁, y₁, z₁). p (x, y, z) denotes one of the remaining initial source points. Define the point set S₁. For any point p, if a sphere with p as its center and r as its radius contains p₁, then p belongs to S₁. d(p₁, p) denotes the distance between p₁ and p. To facilitate linear classification, the distance constant D is defined, and the expression of D is:

$$\begin{array}{l}D=\frac{r^2-d^2\left(p_1,p\right)}{r^2},\\ d^2\left(p_1,p\right)={\left(x_1-x\right)}^2+{\left(y_1-y\right)}^2+{\left(z_1-z\right)}^2\end{array}$$

(15)

Then the relative magnitude of d(p₁, p) and radius r can be judged by (e^D − 1). If (e^D − 1) ≥ 0, then p belongs to set S₁.

Then, the points belonging to set S₁ are removed from all initial source points. In the remaining points, search for the point with the highest λ(p) and define it as p₂, and then find out the corresponding point set S₂. The above linear classification procedure will be repeated until the number of unclassified points is less than the threshold N_t. An example is shown in Figure 11b.

4.3. Dynamic Weighting

After classification, n point sets are obtained (S₁~S_n). The number of points in S_i is recorded as N_i. n high density initial source points (p₁~p_n) are obtained. The coordinates of P_i are (x_i, y_i, z_i). The remaining points in S_i are denoted as p_in (x_in, y_in, z_in), n ∈ [1, N_i−1]. Then the dynamic weight coefficient α_in of p_in is defined as:

$$\begin{array}{l}{\alpha }_{in}=\frac{1-g^2\left(p_i,p_{in}\right)/G}{N_i-1},G={\displaystyle \sum _{n=1}^{N_i-1}{g}^2\left(p_i,p_{in}\right)},\\ g\left(p_i,p_{in}\right)=\sqrt{r^2-{\left(x_i-x_{in}\right)}^2-{\left(y_i-y_{in}\right)}^2-{\left(z_i-z_{in}\right)}^2}\end{array}$$

(16)

The final PD source P_i (X_i, Y_i, Z_i) is determined by the initial source point in S_i. Finally, the calculated coordinate of the PD source is:

$$X_i={\displaystyle \sum _{n=1}^{N_i-1}{\alpha }_{in}{x}_{in}},Y_i={\displaystyle \sum _{n=1}^{N_i-1}{\alpha }_{in}{y}_{in}},Z_i={\displaystyle \sum _{n=1}^{N_i-1}{\alpha }_{in}{z}_{in}}$$

(17)

4.4. Application Results and Analysis

The obtained 20 sets of TDOAs (some of which are shown in Table 2 and Table 3) are used to verify the dynamic weighting algorithm. r is set as 5 cm, and λ₀ is set as 2. The results of the secondary correlation TDOA extraction and dynamic weighting algorithms are taken as examples. The actual PD source (marked with red solid dots) and the calculated PD source (marked with black five-pointed stars) are shown in Figure 11c.

The detailed PD source coordinates calculated from the TDOAs obtained with both methods are shown in

Table 4. Location results and error.

Method	PD Source (cm)	Location Result (cm)	Error (cm)
Optimized energy accumulation method	P1 (20, 30, 20)	(25.5, 36.7, 24.9)	10.0
	P2 (500, 240, 180)	(494.7, 233.5, 178.7)	8.5
Secondary correlation method	P1 (20, 30, 20)	(29.8, 28.1, 26.9)	12.1
	P2 (500, 240, 180)	(503.7, 235.5, 185.2)	8.4
	P3 (220, 135, 55)	(233.7, 141.0, 52.4)	15.2
	P4 (285, 155, 50)	(272.1, 148.0, 52.0)	14.8

. The location error is within the acceptable limits considering that the propagation speed of the electromagnetic wave is 30 cm/ns.

The second set of experiments used two other PD sources P₃ (220, 135, 55) and P₄ (285, 155, 50), which are deeper in the tank and are close to each other. One of the measured UHF partial discharge signals is shown in Figure 12. Since the PD signal fluctuates and oscillates while the interference signal is evident, the secondary correlation TDOA extraction method is used. The location results for the second set of experiments are also presented in Table 4.

Comparing the results of the two sets of experiments with the secondary correlation method, the location error is relatively larger in the second set. In addition to the effect of PD pulse oscillations, the main causes of the errors could be the proximity of the two PD sources and the cross-influence of the UHF partial discharge signals.

Experimental results show that the proposed PD location method can maintain high accuracy when applied to the MSPD location inside power equipment with simple internal structures.

5. Conclusions

In this paper, a PD fault location method based on the UHF signals is proposed. To obtain accurate TDOA, the optimized energy accumulation method and the secondary correlation method are proposed, and a comprehensive application strategy is provided to fully exploit the advantages of both methods. Moreover, the dynamic weighting algorithm is used to further improve the location accuracy, and finally, the high-precision MSPD location inside the power equipment is achieved.

The main conclusions of this paper are as follows:

(1)

The optimized energy accumulation curve method overcomes the effect of the interference signal before the wave peak. The secondary correlation method further improves the interference capability of the TDOA extraction. Both methods have smaller errors than the first peak method for a single PD source.

(2)

For the MSPD location inside a piece of power equipment, the optimized energy accumulation method should be applied when the interference signal is weak and the distance between the reference sensor and the PD source is estimated to be small. The secondary correlation method should be applied when the interference signal is strong and the distance between the reference sensor and the PD source is estimated to be large.

(3)

The proposed dynamic weighting algorithm can fully utilize multiple TDOA data to reduce the effect of accidental data and improve location accuracy.

The UHF partial discharge signal can still be measured in the presence of internal structures, such as iron cores and windings. However, due to the obstruction of the structures, the propagation path of the electromagnetic wave is more complex, and the attenuation of the signal wave head is serious [22]. For practical application, the proposed PD location method can be further investigated by taking into account transmission velocity variations, electromagnetic wave refraction and other factors.