# An Analytical Solution for Acoustic Emission Source Location for Known P Wave Velocity System

Authors: Longjun Dong; Xibing Li; Gongnan Xie
Journal: Mathematical Problems in Engineering
 (2014)

---

## Abstract

This paper presents a three-dimensional analytical solution for acoustic emission source location using time difference of arrival (TDOA) measurements fromN receivers, N⩾5. The nonlinear location equations for TDOA are simplified to linear equations, and the direct analytical solution is obtained by solving the linear equations. There are not calculations of square roots in solution equations. The method solved the problems of the existence and multiplicity of solutions induced by the calculations of square roots in existed close-form methods. Simulations are included to study the algorithms' performance and compare with the existing technique.


---

## Body

## 1. Introduction

The solution of the problem of locating a signal source using time difference of arrival (TDOA) measurements has numerous applications in aerospace, surveillance, structural health, nondestructive testing, navigation, industrial process, speaker location, machine condition, the monitoring of nuclear explosions, and mining induced areal seismology [1–17]. Many authors have discussed and faced numerous problems connected with the location of acoustic emission. The time difference of arrival TDOA method, based on estimates of time delay for a correlated signal as detected by spatially distributed sensor elements in an array, remains a commonly used technique for source location [5].The TDOAs are proportional to the differences in sensor-source range, called range differences (RDs). Conventionally, the source location is estimated from the intersection of a set of hyperboloids defined by the RD measurements and the known sensor locations [18]. The inverse problem for TDOA source location is usually solved by an iterative technique such as nonlinear least squares, minimum error, or an optimization method in recognizing that the equations are nonlinear with respect to source location [19–23]. Kalman filtering has also been used to iteratively solve the TDOA source location problem for microphone speaker location [24, 25]. Although these iterative algorithms are resilient to varying extents to errors in arrival time data, they may be computationally expensive. This is a key consideration in some real-time applications [21, 22].Closed-form solutions are usually less computationally burdensome than iterative, nonlinear minimization, or the ML method and achieve good accuracy [18]. Several closed-form analytical solutions to the TDOA source location problem have been developed [4, 21, 26, 27]. Closed-form solutions have been found in terms of intersecting spheres of distance from each sensor to an arbitrarily located source spherical intersection method, for a monitoring array of four or more sensors, in some cases admitting dual source location solutions [28, 29]. A spherical interpolation (SI) method based on linear least-squares equation error minimization with respect to source range has also been developed [30]. The two-dimensional (2D) TDOA source location problem has also been solved [21, 31]. A technique has been successfully developed by Kundu et al. [32] for localizing acoustic source in anisotropic plates that avoids the need to solve a system of nonlinear equations. The advantage of the technique is that the knowledge of the wave velocity in isotropic or anisotropic plate is not required in two-dimensional structure, but the three-dimensional conditions are not considered. A family of exact solutions for 2D TDOA based on arrival times at four sensors in a plane is also known [33]. A linear intersection (LI) estimator for 3D source location has also been proposed, based on the closest intersection of source bearing lines calculated for four sensor element subarrays under a far-field assumption [22, 34]. This method is claimed to closely model the behavior of the maximum likelihood estimator and to be consistently less sensitive to noise in the time differences than the SI method. A similar but less heuristic closed-form solution based on the far-field assumption has also been proposed [4]. A simultaneous analytical solution for the spatial coordinates of a source and the redundant radial distance of the source/time of origin of the signal has also been reported in the literature [29]. A direct and short derivation of an algorithm based on the closed-form solution of the nonlinear equations for emitter location using time difference of arrival (TDOA) is given [18], and the influences of imaginary roots on solutions were discussed. Spencer developed closed-form solutions for both two-dimensional and three-dimensional source locations by formulating the TDOA equations in, respectively, polar and spherical coordinate systems, with the radial direction coincident with the assumed geodesic path of signal propagation to a reference sensor [5].Nevertheless, the existence and multiplicity of the closed-form solutions are affected by square root calculation in solution equation. In this paper, we simplify the nonlinear location equations for TDOA to linear equations and develop a unique three-dimensional analytical solution without having square root calculation.


## 2. Analytical Solutions

An acoustic emission source is located atsource(x,y,z) and the 5 geosensors are located at S1(x1,y1,z1), S2(x2,y2,z2), S3(x3,y3,z3), S4(x4,y4,z4), and S5(x5,y5,z5). The equations governing the position of acoustic emission source and sensors are
(1)(x1-x)2+(y1-y)2+(z1-z)2=v2t02,(2)(x2-x)2+(y2-y)2+(z2-z)2=v2(t0+t12)2,(3)(x3-x)2+(y3-y)2+(z3-z)2=v2(t0+t13)2,(4)(x4-x)2+(y4-y)2+(z4-z)2=v2(t0+t14)2,(5)(x5-x)2+(y5-y)2+(z5-z)2=v2(t0+t15)2,
where t0 is acoustic wave propagation time from source to the nearest sensor. The signal arrives first to the nearest sensor. After identification of the nearest sensor and calculation of propagation time delay between nearest sensor and other sensors, the sensors are numbered. The nearest sensor from source is numbered as sensor 1 (S1). t12, t13, t14, and t15 are the propagation time delay between sensor 1 and sensors 2, 3, 4, and 5, respectively. v is the P wave velocity. Equations (1)–(5) represent spheres, with the center at respective sensor position, passing through the source. Any two of the spheres (1)–(5) intersect and the source is located on the intersecting circle. The equation of the intersecting plane for the spheres (1) and (2) can be obtained by taking difference of (1) and (2). This equation of intersecting plane is given in (6). In the following equations, li is a constant:
(6)2(x2-x1)x+2(y2-y1)y+2(z2-z1)z+2t12v2t0=l1,
where l1=x22-x12+y22-y12+z22-z12-v2t122.Similarly, the equations for the intersecting planes for spheres (1) and (3), for spheres (1) and (4), and for spheres (1) and (5) are obtained by taking difference of (1) and (3) and (1) and (4), as well as (1) and (5), respectively. These equations are given in (7), (8), and (9), respectively:
(7)2(x3-x1)x+2(y3-y1)y+2(z3-z1)z+2v2t13t0=l2,(8)2(x4-x1)x+2(y4-y1)y+2(z4-z1)z+2v2t14t0=l3,(9)2(x5-x1)x+2(y5-y1)y+2(z5-z1)z+2v2t15t0=l4,
where l2=x32-x12+y32-y12+z32-z12-v2t132, l3=x42-x12+y42-y12+z42-z12-v2t142, and l4=x52-x12+y52-y12+z52-z12-v2t152.Equations (6), (7), (8), and (9) constitute a set of linear equations, which were rewritten as
(10)l5x+l6y+l7z+l8t0=l1,l9x+l10y+l11z+l12t0=l2,l13x+l14y+l15z+l16t0=l3,l17x+l18y+l19z+l20t0=l4,
where l5=2(x2-x1), l6=2(y2-y1), l7=2(z2-z1), l8=2t12v2, l9=2(x3-x1), l10=2(y3-y1), l11=2(z3-z1), l12=2t13v2, l13=2(x4-x1), l14=2(y4-y1), l15=2(z4-z1), l16=2t14v2, l17=2(x5-x1), l18=2(y3-y1), l19=2(z5-z1), and l20=2t15v2.Therefore, thex,y,z, which can be obtained by solving the linear equations, can be expressed and simplified as
(11)x=[l1(l10l16l19-l11l16l18-l12l14l19+l12l15l18-l10l15l20+l11l14l20)+l6(l15l2l20-l16l19l2+l12l19l3-l11l20l3+l11l16l4-l12l15l4)+l7(l16l18l2-l14l2l20-l12l18l3-l10l20l3-l10l16l4+l12l14l4)+l8(l14l19l2-l15l18l2-l10l19ll3+l11l18l3+l10l15l4-l11l14l4)]×[l5(l10l16l19-l11l16l18-l12l14l19+l12l15l18-l10l15l20+l11l14l20)+l6(l11l16l17+l12l13l19-l12l15l17-l11l13l20)+l7(l10l13l20-l10l16l17-l12l13l18+l12l14l17)+l8(l10l15l17-l10l13l19+l11l13l18-l11l14l17)+l9(l15l6l20-l19l16l6+l7l16l18-l20l14l17+l19l14l8-l15l18l8)]-1,y=[l1(l11l16l17-l12l13l19-l12l15l17+l11l13l20-l16l19l9+l15l19l20)+l5(l16l2l19-l5l20l2+l12l19l3-l11l20l3+l11l16l4-l12l15l4)+l7(l16l17l2-l13l2l20-l12l17l3-l12l13l4)+l9(l20l3l7-l16l4l7-l19l8l3-l15l14l8)+l8(l13l19l2-l15l17l2-l11l17l3+l11l13l4)]×[l5(l10l16l19-l11l16l18-l12l14l19+l12l15l18-l10l15l20+l11l14l20)+l6(l11l16l17+l12l13l19-l12l15l17-l11l13l20)+l7(l10l13l20-l10l16l17-l12l13l18+l12l14l17)+l8(l10l15l17-l10l13l19+l11l13l18-l11l14l17)+l9(l15l6l20-l19l16l6+l7l16l18-l20l14l17+l19l14l8-l15l18l8)]-1,z=[l5(l18l16l2-l14l2l20-l12l18l3+l10l3l20-l16l10l4+l12l14l4)+l6(l13l20l2-l16l2l17+l12l17l3-l12l13l4)+l8(l14l17l2-l13l2l18-l10l17l3+l10l13l4)+l9(l1l14l20-l1l16l18-l20l6l3+l16l6l4+l18l3l8-l14l8l4)]×[l5(l10l16l19-l11l16l18-l12l14l19+l12l15l18-l10l15l20+l11l14l20)+l6(l11l16l17+l12l13l19-l12l15l17-l11l13l20)+l7(l10l13l20-l10l16l17-l12l13l18+l12l14l17)+l8(l10l15l17-l10l13l19+l11l13l18-l11l14l17)+l9(l15l6l20-l19l16l6+l7l16l18-l20l14l17+l19l14l8-l15l18l8)]-1.
Equation (10) can also rewrite as
(12)AS=B,
where
(13)A=[l5l6l7l8l9l10l11l12l13l14l15l16l17l18l19l20],S=[xyzt0],B=[l1l2l3l4].
The x,y, and z can be also obtained by solving the linear systems. The proposed method is suitable for known wave velocity system. The closed-form solution for unknown wave velocity system is reported in [35].


## 3. Validated Examples

### 3.1. Example 1: Simulations for Sound Sources

For example, consider five receivers and a speed of sound of 330 m/s. The Cartesian coordinates of five sensors are S1 (250.25, 450.17, and 150.28), S2 (290.11, 180.13, and 140.18), S3 (150.35, 190.21, and 30.33), S4 (380.38, 250.28, and 98.15), and S5 (500.55, 350.25, and 210.45) m. The two source locations with Cartesian coordinates of P1 (310.75, 200.33, and 125.11) and P2 (210.88, 290.99, and 175.89) m and the original times of two sources are 0 s. The trigger times are obtained according to the distance formula between spatial two points as well as distance and velocity formula. The trigger times of five receivers (S1–S5) of P1 are 0.782697, 0.098713, 0.565408, 0.272225, and 0.777218 s, respectively. The trigger times of five receivers (S1–S5) of P2 are 0.502922, 0.426859, 0.566984, 0.578392, and 0.902068 s, respectively.To validate the proposed method, the coordinates of receivers and their trigger time were used to calculateli(i=1,2,…,20) according to (7)–(10). Then, the calculated li(i=1,2,…,20) were substituted in (11) to resolve the source location. The calculated results of sources P1 and P2 are (310.7499, 200.3302, and 125.1101) and (210.8797, 290.9899, and 175.8904) m, respectively. Results show that the calculated coordinates of P1 and P2 are consistent with authentic results.


### 3.2.  Example 2: Simulations of a Mine Acoustic Emission Location System and Location Error Analysis

A numerical test for an acoustic emission/microseismic monitoring system of a mine was used to validate the proposed method and analyze location errors. The locations of five sensors are shown in Figure1, and the coordinates are (200, 500, and 300), (300, 100, and 100), (400, 350, and 280), (600, 250, and 250), and (700, 200, and 180). The locations of acoustic emission/microseismic sources E1, E2, E3, E4, and E5 are 1 (250, 300, and 120), 2 (310, 450, and 130), 3 (450, 260, and 110), 4 (650, 180, and 135), and 5 (780, 250, and 125), respectively, and the unit is m; P-wave velocity is 5000 m/s. Errors of wave velocity and arrivals were discussed.Sketch map of acoustic emission/microseismic source location: A indicates analytical solution; N indicates numerical solution.
(a)
Locations of sensors and sources
(b)
Spatial locations of sensors and sourcesFirstly, the acoustic emission sources are located by proposed method using true wave velocity. The results are shown in Figure2. It can be seen that the calculated results are consistent with the true locations. Therefore, it can be concluded that proposed method is reasonable. To investigate the location errors induced by errors of velocities, the method STT was used to compare and analyze. The locations of acoustic emissions by STT method without considering the errors of velocities and arrivals are consistent with true locations (Figure 2). The errors of absolute distant are shown in Figure 2. It also can be seen that the results of STT are consistent with true locations and results of the proposed analytical solution method. It is proved that the STT method is stable and reliable. The intervals of velocity errors are from 4500 m/s to 5500 m/s. True velocity is 5000 m/s. The variable interval is 200 m/s. Location results with analytical solutions and numerical methods STT are shown in Figure 1. The absolute distance errors between the real coordinates and calculated coordinates are shown in Figure 2 under different levels of velocity error. The left vertical axis is the absolute distance error. The original time errors are also plotted in the right vertical axis of Figure 2. It can be seen from positioning results of Figures 1 and 2 that the positioning accuracy is affected by the velocity value. In view of the spatial distribution of positions, symmetrical changes of velocity will result in a positioning results symmetry spatial location using analytical solution positioning, since analytical solutions are based on the results of each positioning operation parameters of the true value. The location errors are completely dependent on the operator error parameters which are offset from the true value. Symmetrical changes of velocity will induce irregular spatial position using numerical iterative method, STT. If velocity error is large, it may cause serious positioning errors, which will locate a serious deviation from the true value, such as the E4-N. It can be seen from Figure 2 that higher positioning analytical location method requires measurement data with high quality, but when measurement data is reasonable, it can obtain accurate positioning results. Throughout all of located microseismic sources, it can be easily found that the greatest absolute distance error is from the numerical solution of STT. A velocity error of 500 m/s can result in a big absolute distance error of 700 m using STT, while the positioning error of the analytical solution is 450 m.Figure 2
Errors of distance error from located source to real source and original time induced by insufficient know velocity error: D indicates absolute distant error; TS indicates original time error.Secondly, to investigate the location errors induced by errors of arrivals, the method STT was also used to compare and analyze. The error of arrivals is limited between −1000μs and +1000 μs. The 1000 errors are generated by the random function of Microsoft Excel; then arrivals with errors were obtained by true arrivals plus generated errors. In this study, the minimum, average, and maximum of arrivals with errors were used to investigate the location errors induced by arrival errors. The calculated results are shown in Figure 3. The absolute distance error is expressed in the left vertical axis. Original time error is expressed in the right vertical axis. It can be clearly seen in Figure 3 that arrival errors also significantly affect the location errors. The absolute distant errors and original errors of STT are always larger than those of analytical method. The maximums of absolute distant errors of STT and analytical method are 12.5 m and 6.5 m, respectively. The maximums of original time errors of STT and analytical method are 0.001 s and 0.00023 s, respectively.Figure 3
Errors of Distance and original time induced by arrival time reading error: the lefty axil expressing the distance error from located sources to real sources and right y axil expressing the error of original time, E indicates absolute distant error; TS indicates original time errors; A indicates analytical solutions; and N indicates STT results.


## 3.1. Example 1: Simulations for Sound Sources

For example, consider five receivers and a speed of sound of 330 m/s. The Cartesian coordinates of five sensors are S1 (250.25, 450.17, and 150.28), S2 (290.11, 180.13, and 140.18), S3 (150.35, 190.21, and 30.33), S4 (380.38, 250.28, and 98.15), and S5 (500.55, 350.25, and 210.45) m. The two source locations with Cartesian coordinates of P1 (310.75, 200.33, and 125.11) and P2 (210.88, 290.99, and 175.89) m and the original times of two sources are 0 s. The trigger times are obtained according to the distance formula between spatial two points as well as distance and velocity formula. The trigger times of five receivers (S1–S5) of P1 are 0.782697, 0.098713, 0.565408, 0.272225, and 0.777218 s, respectively. The trigger times of five receivers (S1–S5) of P2 are 0.502922, 0.426859, 0.566984, 0.578392, and 0.902068 s, respectively.To validate the proposed method, the coordinates of receivers and their trigger time were used to calculateli(i=1,2,…,20) according to (7)–(10). Then, the calculated li(i=1,2,…,20) were substituted in (11) to resolve the source location. The calculated results of sources P1 and P2 are (310.7499, 200.3302, and 125.1101) and (210.8797, 290.9899, and 175.8904) m, respectively. Results show that the calculated coordinates of P1 and P2 are consistent with authentic results.


## 3.2.  Example 2: Simulations of a Mine Acoustic Emission Location System and Location Error Analysis

A numerical test for an acoustic emission/microseismic monitoring system of a mine was used to validate the proposed method and analyze location errors. The locations of five sensors are shown in Figure1, and the coordinates are (200, 500, and 300), (300, 100, and 100), (400, 350, and 280), (600, 250, and 250), and (700, 200, and 180). The locations of acoustic emission/microseismic sources E1, E2, E3, E4, and E5 are 1 (250, 300, and 120), 2 (310, 450, and 130), 3 (450, 260, and 110), 4 (650, 180, and 135), and 5 (780, 250, and 125), respectively, and the unit is m; P-wave velocity is 5000 m/s. Errors of wave velocity and arrivals were discussed.Sketch map of acoustic emission/microseismic source location: A indicates analytical solution; N indicates numerical solution.
(a)
Locations of sensors and sources
(b)
Spatial locations of sensors and sourcesFirstly, the acoustic emission sources are located by proposed method using true wave velocity. The results are shown in Figure2. It can be seen that the calculated results are consistent with the true locations. Therefore, it can be concluded that proposed method is reasonable. To investigate the location errors induced by errors of velocities, the method STT was used to compare and analyze. The locations of acoustic emissions by STT method without considering the errors of velocities and arrivals are consistent with true locations (Figure 2). The errors of absolute distant are shown in Figure 2. It also can be seen that the results of STT are consistent with true locations and results of the proposed analytical solution method. It is proved that the STT method is stable and reliable. The intervals of velocity errors are from 4500 m/s to 5500 m/s. True velocity is 5000 m/s. The variable interval is 200 m/s. Location results with analytical solutions and numerical methods STT are shown in Figure 1. The absolute distance errors between the real coordinates and calculated coordinates are shown in Figure 2 under different levels of velocity error. The left vertical axis is the absolute distance error. The original time errors are also plotted in the right vertical axis of Figure 2. It can be seen from positioning results of Figures 1 and 2 that the positioning accuracy is affected by the velocity value. In view of the spatial distribution of positions, symmetrical changes of velocity will result in a positioning results symmetry spatial location using analytical solution positioning, since analytical solutions are based on the results of each positioning operation parameters of the true value. The location errors are completely dependent on the operator error parameters which are offset from the true value. Symmetrical changes of velocity will induce irregular spatial position using numerical iterative method, STT. If velocity error is large, it may cause serious positioning errors, which will locate a serious deviation from the true value, such as the E4-N. It can be seen from Figure 2 that higher positioning analytical location method requires measurement data with high quality, but when measurement data is reasonable, it can obtain accurate positioning results. Throughout all of located microseismic sources, it can be easily found that the greatest absolute distance error is from the numerical solution of STT. A velocity error of 500 m/s can result in a big absolute distance error of 700 m using STT, while the positioning error of the analytical solution is 450 m.Figure 2
Errors of distance error from located source to real source and original time induced by insufficient know velocity error: D indicates absolute distant error; TS indicates original time error.Secondly, to investigate the location errors induced by errors of arrivals, the method STT was also used to compare and analyze. The error of arrivals is limited between −1000μs and +1000 μs. The 1000 errors are generated by the random function of Microsoft Excel; then arrivals with errors were obtained by true arrivals plus generated errors. In this study, the minimum, average, and maximum of arrivals with errors were used to investigate the location errors induced by arrival errors. The calculated results are shown in Figure 3. The absolute distance error is expressed in the left vertical axis. Original time error is expressed in the right vertical axis. It can be clearly seen in Figure 3 that arrival errors also significantly affect the location errors. The absolute distant errors and original errors of STT are always larger than those of analytical method. The maximums of absolute distant errors of STT and analytical method are 12.5 m and 6.5 m, respectively. The maximums of original time errors of STT and analytical method are 0.001 s and 0.00023 s, respectively.Figure 3
Errors of Distance and original time induced by arrival time reading error: the lefty axil expressing the distance error from located sources to real sources and right y axil expressing the error of original time, E indicates absolute distant error; TS indicates original time errors; A indicates analytical solutions; and N indicates STT results.


## 4. Conclusions

We simplify the nonlinear location equations for TDOA to linear equations. Based on solving simplified linear equations, we obtain a unique three-dimensional analytical solution for acoustic emission source location using TDOA measurements fromN receivers, N⩾5. There are no calculations of square roots in the solution equations of the proposed method. The problems of the existence and multiplicity induced by calculations of square roots in existed close-form methods were solved successfully. Simulations are included to study the algorithms' performance and compare with existing STT method.



---

*Source: 290686-2014-03-27.xml*