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

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.


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 multidata 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.

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 1 as the reference, the propagation distance difference ∆r i1 is: (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].

Optimized Energy Accumulation Method
The energy accumulation curve E is calculated based on the UHF partial discharge signal: 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 first and the i-th sensors. The coordinates of Si (xi, yi, zi) and c are known. When m than three Δti1 are achieved, the coordinate of P (x, y, z) can be obtained by solving equation groups.
To suppress the periodic interference or white noises, the wavelet denoising and tering are performed on the collected UHF signals before the TDOA extraction [19,20]

Optimized Energy Accumulation Method
The energy accumulation curve E is calculated based on the UHF partial discha signal: where uk denotes the voltage value of the k-th point in the PD signal segment, and n notes the sampling point number. Two channels of the UHF signals and corresponding energy accumulation curves shown in Figure 1. The PD pulse period in Figure 1a corresponds to the rising segmen the energy accumulation curve in Figure 1b. The traditional energy accumulation method takes the starting point of the cu slope as the reference point to calculate the TDOA [12]. This is easily affected by inter ence signals before the first wave peak.
In this paper, the wavelet modulus maximum calculation is employed to detect gularities so that the inflection point of the slope of the curve can be determined. inflection point of the curve slope is then used as a reference point to extract the TDO Wavelet transform is applied firstly to smooth the energy accumulation curve f(x reduce the impact of interference signals and noises. The cubic B-spline function θ(x 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 Define ϕ (x) = dθ(x)/dx, and then perform the convolution calculation with f (x) and ϕ (x). The calculated wavelet transform function is W f s (x), which can be expressed as: where ⊗ denotes convolution calculation, and s denotes the transform scale. W f s (x) 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 f s (x). 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.
where ⊗ denotes convolution calculation, and s denotes the transform scale.
Then the inflection point of the slope of the energy accumulation curve can be d mined by calculating the partial modulus maximum of the absolute value of ' s W f Figure 2 shows the wavelet modulus maximum lines corresponding to the energy a mulation curves in Figure 1. It can be seen that the modulus maximum line for the wave peak of the PD pulse is significantly higher than the lines for the subsequent w peaks and noises. Therefore, the time point corresponding to the modulus maximum with the highest amplitude can be used as a reference point to calculate the TDOA. optimized energy accumulation method overcomes the effect of the interference si before the wave peak and achieves high accuracy when the PD pulse is prominent.

Secondary Correlation Method
The UHF partial discharge signals of two sensors from the same PD source are orded as s1(t) and s2(t), respectively, which are expressed as:

t a s t t n t s t a s t t n t
where s(t) denotes signal of the PD source; a1 and a2 denote attenuation constants, res tively; t1 and t2 denote time differences, Δt = t2 − t1 denotes the TDOA of the two sig n1(t) and n2(t) are signal distortions and noises, respectively.
Employ the self-correlation and the cross-correlation calculation to signals s1(t) s2(t). The self-correlation function R s 1 s 1 (τ) and cross-correlation function R s 1 s 2 (τ) are

Secondary Correlation Method
The UHF partial discharge signals of two sensors from the same PD source are recorded as s 1 (t) and s 2 (t), respectively, which are expressed as: where s(t) denotes signal of the PD source; a 1 and a 2 denote attenuation constants, respectively; t 1 and t 2 denote time differences, ∆t = t 2 − t 1 denotes the TDOA of the two signals; n 1 (t) and n 2 (t) are signal distortions and noises, respectively. Employ the self-correlation and the cross-correlation calculation to signals s 1 (t) and s 2 (t). The self-correlation function R s 1 s 1 (τ) and cross-correlation function R s 1 s 2 (τ) are: Entropy 2023, 25, 572 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 crosscorrelation function R s 1 s 2 (τ) is calculated as: where R n 1 n 2 (τ) ≈ 0 . If s(t) is unrelated to n 1 (t) and n 2 (t), then: In Equation (8), when τ = ∆t, R ss reaches its maximum, so does R s 1 s 2 (τ). Thus, ∆t can be determined by: 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, 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 1 (t) and perform self-correlation calculation to s 1 (t). The self-correlation function R s 1 s 1 (τ) is expressed as: Then employ cross-correlation calculation with reference signals s 1 (t) and s 2 (t), as in Equation (8). Since R s 1 s 1 (τ) and R s 1 s 2 (τ) are functions with τ as the independent variable, τ can be substituted by t.
Then, perform the cross-correlation calculation with R s 1 s 1 (t) and R s 1 s 2 (t) and obtain the secondary cross-correlation function R R 1 R 2 (τ), in which the secondary correlation calculation of the first-round correlation functions (such as R ss (t) and R sn 1 (t)) will appear. In R s 1 s 1 (t) and R s 1 s 2 (t), according to Equations (8) and (11), R n 1 n 1 (t) are unrelated to other correlation functions. If R ss (t), R sn 1 (t) and R sn 2 (t) are unrelated, then the cross-correlation Thus, the secondary cross-correlation function R R 1 R 2 (τ) can be simplified as: where R ss,ss (τ) denotes the self-correlation function of R ss (t). Then, ∆t can be expressed as: In the secondary correlation method, the weight of original signal s(t) is strengthened, while the influence of noises n 1 (t) and n 2 (t) are weakened in the first round of correlation calculations. The effect of R sn 1 (τ) and R sn 2 (τ) are weakened in the second round of the correlation calculation.

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 microstrip 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 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.

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.

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 l1, l2 and l3 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 = 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.
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 1 (0.

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 1 , l 2 and l 3 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 = 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.

Application of Optimized Energy Accumulation Method
The application of the optimized energy accumulation method to the MSPD ca based on time window scanning. In this paper, the width of the main PD pulse (def by 30% of the peak value) is about 13 ns, as shown in Figure 5, where l1, l2 and l3 de the length, width and height of the tank, respectively. The theoretical maximum TDO the signals from the same PD source is ∆t max = l 1 2 +l 2 2 +l 3 2 /c ≈ 20 ns, where c denotes propagation velocity of the electromagnetic wave in air and its approximate value i cm/ns.

Time (ns)
Amplitude (V) Figure 5. Single UHF partial discharge waveform of PD source.
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. 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. 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. 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. First, for the UHF signal in a time window of 0-40 ns, the energy accumulatio are conversed. The wavelet transform and the modulus maximization are perfo sufficiently high modulus maximum lines do not appear simultaneously in the ti dow of the signals, then the time window is advanced by the scanning step of 10 above procedure is repeated for the time window of 10-50 ns. When sufficien modulus maximum lines appear simultaneously on all signals, record the corres time points.
Take a set of multi-channel synchronized signals as an example. For time w shown in Figure 6, sufficiently high modulus maximum lines of all signals appea taneously in the 40 ns time window, as shown in Figure 7b. Since interference noise is weak and the main PD pulse is well-defined, each energy accumulation one inflection point, and the TDOAs can be computed from the corresponding tim 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. However, for time window 1 in Figure 6, when the interferences tiple peaks occur in a time window, several energy accumulation curv one inflection point, as shown in Figure 8, leading to hard determinat maximum value lines. Since the interference is accidental, the TDOA data can still be multiple sets of repeated experimental data and summarizing the calcu ever, the TDOA extraction is less efficient. At the same time, this ap cope with the extreme case that multiple PD main pulses are exactly w window.

Application of Secondary Correlation Method
The application of the secondary correlation method in the MSP tract the feature signal segment from the reference PD signal and the ondary correlation with the signals of other channels.
Take the PD signals of channels 1 and 4 shown in Figure 9 as exa steps of the signal feature segment and the TDOA extraction are descr 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.

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: 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.

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 1 (t). First, the highest value of s 1 (t) is searched for and marked as s 1 (t m1 ). Then, search for the nearest points whose amplitude equals 30%·s 1 (t m1 ), and the corresponding time coordinates are marked as t s1 , t e1 . The fragment of s 1 (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 4 (t). From s 4 (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 (τ) of f 1 (t) and T 1 (t) is deduced. The function curve of R R 1 R 2 (τ) 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 1 (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 4 (t), as shown in Figure 10. Then, the second TDOA is extracted. The TDOA extracting process above will be repeated until s 1 (t mi ) is lower than the threshold s t = 30%·s 1 (t m1 ).  (1) Set the signal of channel 1 as the reference, which is marked as s1(t). First, the highest value of s1(t) is searched for and marked as s1(tm1). Then, search for the nearest points whose amplitude equals 30%·s1(tm1), and the corresponding time coordinates are marked as ts1, te1. The fragment of s1(t) in the time interval [ts1, te1] is defined as the first reference fragment s1r1(t).
(2) Signal of channel 4 is the contrast signal, which is marked as s4(t). From s4(t), the contrast fragment s4c1(t) corresponding to the reference fragment s1r1(t) is selected. Considering the theoretical maximum time difference Δtmax, the time interval of s4c1(t) should be [ts1 − Δtmax, te1 + Δtmax]. Then, based on Equations (11)-(13), the secondary function R R 1 R 2 (τ) of f1(t) and T1(t) is deduced. The function curve of R R 1 R 2 (τ) is shown in Figure 9b. The first TDOA is then extracted. (3) After completing the first TDOA extraction, set the function value of [ts1, te1] in s1 (t) to 0 and obtain a new reference signal s1′(t). Then, repeat the above steps to extract the new reference fragment s1r2(t) from s1′(t) and the contrast fragment s4c2(t) from s4(t), as shown in Figure 10. Then, the second TDOA is extracted. The TDOA extracting process above will be repeated until s1(tmi) is lower than the threshold st = 30%·s1(tm1). 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.

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 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.

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 Tables 2 and 3.  According to Tables 2 and 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 1 (0.2 m, 0, 0.2 m) and PD source P 1 (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 1 has a large amplitude when it arrives at the reference sensor S 1 and can be clearly distinguished from the interference signal. Therefore, the arrival time of the signal at the reference node T 1 is relatively accurate, and the TDOA extraction accuracy is relatively high.
When the distance between the reference node S 1 (0.2 m, 0, 0.2 m) and the PD source P 2 (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 Tables 2 and 3. Therefore, the PD sources localized by only one set of TDOAs may not be accurate.

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.

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: 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 /3 denotes the volume of the sphere. The density threshold value λ 0 will be defined according to the number of initial source points. The error initial source points (λ(p) ≤ λ 0 ) will be removed first, such as circuited point k in Figure 11a.

Linear Classification
In addition to the error points, the number of initial source points N0 is counted. Then, the point with the highest λ(p) is searched for and defined as p1 (x1, y1, z1). p (x, y, z) denotes one of the remaining initial source points. Define the point set S1. For any point p, if a sphere with p as its center and r as its radius contains p1, then p belongs to S1. d(p1, p) denotes the distance between p1 and p. To facilitate linear classification, the distance con-

Linear Classification
In addition to the error points, the number of initial source points N 0 is counted. Then, the point with the highest λ(p) is searched for and defined as p 1 (x 1 , y 1 , z 1 ). p (x, y, z) denotes one of the remaining initial source points. Define the point set S 1 . For any point p, if a sphere with p as its center and r as its radius contains p 1 , then p belongs to S 1 . d(p 1 , p) denotes the distance between p 1 and p. To facilitate linear classification, the distance constant D is defined, and the expression of D is: Then the relative magnitude of d(p 1 , p) and radius r can be judged by (e D − 1). If (e D − 1) ≥ 0, then p belongs to set S 1 .
Then, the points belonging to set S 1 are removed from all initial source points. In the remaining points, search for the point with the highest λ(p) and define it as p 2 , and then find out the corresponding point set S 2 . 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.

Dynamic Weighting
After classification, n point sets are obtained (S 1~Sn ). The number of points in S i is recorded as N i . n high density initial source points (p 1~pn ) 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: 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:

Application Results and Analysis
The obtained 20 sets of TDOAs (some of which are shown in Tables 2 and 3) are used to verify the dynamic weighting algorithm. r is set as 5 cm, and λ 0 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. 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 3 (220, 135, 55) and P 4 (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.  Figure 12. UHF partial discharge signals in the second set of experiments.
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.

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. 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.

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.

Data Availability Statement:
No new data were created or analyzed in this study. Data sharing is not applicable to this article.

Conflicts of Interest:
The authors declare that they have no known competing financial interest or personal relationships that could have appeared to influence the work reported in this paper.