Keywords: acoustic method, buried pipe, weakly heterogeneous medium, inverse problem

Acoustic waves are commonly used to locate buried polyethylene pipes. In this preliminary study we are particularly interested in pipes depth. To obtain depth information we are moving towards a multi-sensor solution. Several estimators are implemented and tested on real data. A depth estimator according to the relative delays between sensors is proposed. We compare two relative delays estimators : the method using the cross-correlation and the one using the coherence function. We will verify on real measurements that the second method is much more efficient than the first one. Before discussing the results we will present another approach which consists in adapting the MUSIC (MUltiple SIgnals Classification) algorithm to our problem.

INTRODUCTION

Detecting and locating gas pipelines is a major challenge for network operators. We are interested in buried polyethylene pipelines which is a non-conductive material. Methods have been developed and marketed for years to trace and locate pipes [START_REF] Liu | A Comprehensive Review of Acoustic Methods for Locating Underground Pipelines[END_REF]. Other methods exist, like seismic methods (interferometry [START_REF] Harmankaya | Locating near-surface scatterers using non-physical scattered waves resulting from seismic interferometry[END_REF], surface waves dispersion [START_REF] Filippi | New robust observables on Rayleigh waves affected by an underground cavity: from numerical to experimental modelling[END_REF] etc) for the investigation of near-surface structures, or invasive electromagnetic methods (insertion of metallic cable in the pipe). Here we are interested in acoustic methods. Different acoustic methods, including the pipe excitation and the vibration excitation of the ground surface have been developed [START_REF] Muggleton | Mapping the Uderworld : recent developments in vibro-acoustic techniques to locate buired infrastructure[END_REF]. We have chosen to focus on the pipe excitation method. The advantage of this one is to discern the interest pipe in a dense pipe network. Tools based on this method allow us to estimate the pipe plumb but do not provide any information on the pipe depth. The objective is to obtain some estimates on its depth with an accuracy of 0.1 m, with a non-destructive or noninvasive method and without 'a priori' information on the surrounding basement. Indeed, this method must be applicable regardless of the location. This preliminary study goes in this direction.

We know that the celerity of the sound in the pipe is equivalent to that of the sound in the air. Depending on the difference in wave propagation velocity between the pipe and the underground, the pressure waves contribution is more or less attenuated but the shear waves contribution stays high. Moreover, we consider that the pipe is analogous to a cylindrically secondary source which radiates through the heterogeneous ambiant medium to receivers (geophones) [START_REF] Liu | A Numerical Study of Axisymmetric Wave Propagation in Buried Fluid-Filled Pipes for Optimizing the Vibro-Acoustic Technique When Locating Gas Pipelines[END_REF][START_REF] Mcfadden | Radial Vibrations of Thick-Walled Hollow Cylinders[END_REF]. To support this assumption, we can observe the experimental measurements presented below (Sect. 5.ii).

We first present, section 2, the considerations about problem modeling. Relative time delays and depth estimators, will be developed section 3. We will present results of numerical simulations in section 4. Experimental set-up (electronic and acquisition system developed for our measurements in situ) will be described section 5, with experiments in semi-controlled conditions. Finally we will present results section 6, before concluding, section 7.

PROBLEM MODELING

We assume first, in accordance with our field experiments (Sect. 5.ii) that the cylindrical pipe is equivalent to a secondary source radiating in the surrounding heterogeneous medium when we inject an acoustic signal at one of the ends. The vibrating pipe is then a cylindrical source with resulting symmetries. Geophones located on the surface of the land can record the propagated signal from the pipe through the medium. Secondly, we consider that our medium is weakly heterogenous. We consider this source seen by receivers as a point source. The trench is not taken into account here. On figure 1, we can see the model scheme. Sensors are placed in a line following the perpendicular to the pipe plumb and sensors position is known. We consider a single weakly heterogeneous propagation medium with an unknown average velocity in the medium, noted V0. The coordinates of S are Sx, the Xcoordinate which represents the plumb deviation, and Sz, , the Z-coordinate which represents the pipe depth. For the receivers, S is seen as a source. We note the i th receiver Ri, its X-coordinate is riX and its Z-coordinate is riZ.

ESTIMATES i) Estimators of relative delay times

The source seen by sensors is the pipe (cylindrical secondary source), but we do not know the exact time when the signal is scattered from the pipe. For these reasons we are interested in the relative delay times between sensors. Estimating relative delay times is an important issue, because the estimate of the pipe depth depends on it.

We will use two different methods to estimate the relative delay times. First using the classical cross-correlation and then, using the smoothed coherence transform (SCOT) [START_REF] Carter | Coherence and Time Delay Estimation[END_REF][START_REF] Sun | The SCOT Weighted Adaptive Time Delay Estimation Algorithm Based on Minimum Dispersion Criterion[END_REF].

The SCOT weighted cross correlation method is based on coherence function. The coherence function GXY (f) between two stationary random processes x and y is equal to the crosscorrelation Fourier transform Γ XY (f ) divided by the square root of the product of the two autocorrelation Fourier transform Γ XX (f ) and Γ YY (f ) .

G XY (f ) = Γ XY (f ) √ Γ XX (f ) Γ YY ( f ) (1)
We will verify in the results section 6 that the estimator using the SCOT weighted improves the estimates on our real data.

ii) Depth estimator from relative delay times We want to estimate model unknown variables of models from the relative delay times. We will use generalized least squares to define an estimator [START_REF] Menke | Review of the Generalized Least Squares Method[END_REF][START_REF] Tarantola | Generalized Nonlinear Inverse Problems Solved Using the Least Squares Criterion[END_REF]. Here, sensor one is the reference, all relative delays are calculated from it. Let be θ the parameters vector to be estimated :

θ = [ S X S Z V 0 ] (2)
Sx is the pipe plumb, Sz the pipe depth and V0 the average propagation velocity in the medium. We note τ 1i the relative delay between receivers 1 and i :

τ 1i = | S R i |-| S R 1 | V 0 (3)
i = 2, … , K with K the number of receivers. The relative delays is then written as a function of the vector θ .

τ 1i (θ) = √[S X -r iX ] 2 +[S Z -r iZ ] 2 -√[ S X -r 1X ] 2 +[S Z -r 1Z ] 2 V 0 (4)
We define the following criterion :

C (θ) = [ τ m -τ(θ)] T [ τ m -τ(θ)] (5)
τ m is the relative delays vector estimated from the measurements and τ(θ) is the vector of τ 1i (θ) , the symbol [.] T means the transposed operator.

C(θ) is the sum of the squared errors. We minimize the criterion C(θ) in a nonlinear least squares sense. It will be interesting, to weight the criterion later, in order to take into account the correlation of the measurement noises between sensors since we consider relative delay times.

iii) Depth estimator based on MUSIC algorithm As we will see in section 4.i, with the model used we need a fine precision on the delay times between sensors (between 10 -5 and 10 -6 seconds). This is why we are interested in a so-called "high resolution" algorithm such as the MUSIC algorithm.

The MUSIC algorithm is classically used in the farfield to estimate sources number and its directions of arrival (DOA) [START_REF] Schmidt | Multiple Emitter Location and Signal Parameter Estimation[END_REF]. It can be adapted in near-field to estimate sources distance in addition to their DOA [START_REF] Wang | Mixed-Order MUSIC Algorithm for Localization of Far-Field and Near-Field Sources[END_REF][START_REF] Fossati | Chapter 3: Localization of Buried Objects in Sediment Using High Resolution Array Processing Methods[END_REF][START_REF] Han | Bearing and Range Estimation Algorithm for Buried Object in Underwater Acoustics[END_REF]. In our case, we are looking for a single source in the near-field, we want its 3D location, specially the depth, but we have no knowledge of the propagation velocity.

2

© 2021 by ASME In the far-field case, we can identify the following steps of the algorithm: a) Estimate the variance-covariance matrix of the system from signals received by the sensors. We note this matrix M. b) Decompose the matrix M into eigenvalues and eigenvectors.

c) Define the noise subspace with the eigenvectors corresponding to the smallest eigenvalues. We note Vb the matrix containing the eigenvectors which constitute the noise subspace.

d) From the information on the propagation medium and the knowledge of the antenna array (four geophones), we construct a family of vectors "a" parameterized by the theoretical relative phase state ϕ . In the case of a linear and regularly spaced antenna array we have :

a( ϕ) = [1 e -j ϕ e -j2 ϕ e -j 3 ϕ ... e -j( K-1) ϕ ] T (6)
e) Projection of the vector family a( ϕ) onto the noise subspace and definition of a criterion.

P MUSIC (ϕ) = 1 a (ϕ) H V b V b H a(ϕ) (7)
The symbol [.] H represents the conjugate transposed of [.].

From the ϕ estimate which maximizes the criterion PMUSIC, we can deduce the DOA. In our case, the step (d) is modified. We are in near field, and we do not only want the DOA of signal. We want to obtain the source position. Moreover, the propagation medium is unknown, and we integrate the variable V0 which represents the average propagation velocity of the signal. Delay times are expressed as a function of our model and are no longer proportional to each other as in the far field case. We construct the vector family "a" parameterized by S and V0 from equation

(3). a (S , V 0 ) = [ 1 e -j 2 π f0 | S R 2 |-| S R 1 | V 0 ... e -j 2 πf 0 | S R K |-| S R 1 | V 0 ] T (8)
Where f0 is the central frequency of the signal. We can define a criterion in the same way as equation [START_REF] Carter | Coherence and Time Delay Estimation[END_REF].

P MUSIC (S , V 0 ) = 1 a(S , V 0 ) H V b V b H a(S , V 0 ) (9)

NUMERICAL SIMULATIONS i) Cramer-Rao Bound

In this section we present results of numerical simulations. We calculate the Cramer-Rao Bound (CRB) to obtain information about the estimation accuracy [START_REF] Stoica | MUSIC, maximum likelihood, and Cramer-Rao bound[END_REF]. The CRB represents the smallest possible standard deviations of all unbiased estimates of the model parameters. From our model (equation 2 and 4), we calculate the BRC : [START_REF] Tarantola | Generalized Nonlinear Inverse Problems Solved Using the Least Squares Criterion[END_REF] with Var ( τ 1i ) the variance of a relative delay time.

CRB(θ) = Var ( τ 1i ) V 0 2 ( ∑ i=2 K [ S X -r iX √(SX-riX) 2 +(S Z -r iZ ) 2 - S X √(SX-r1X) 2 +(S Z -r 1Z ) 2 SZ-riX √(SX-riX) 2 +(SZ-riZ) 2 - SZ √(SX-r1X) 2 +(SZ-r1Z) 2 √(SX-r1X) 2 +(S Z -r 1Z ) 2 -√(S X -r iX ) 2 +(S Z -r iZ ) 2 V 0 ] ) -1
We observe, figure 2, the standard deviation evolution of the depth estimate as a function of the standard deviation of relative delay times estimates. This curve is obtained for a distance between the sensors of 0.2 m, an average propagation velocity (V0) in the basement of 500 m/s, a deviation on the ground from the plumb (SX) of 0 m and for the case of two depths (SZ) 0.4 m and 1 m. We note that to achieve 0.1 meter depth accuracy, relative delay times accuracy in the microsecond range is required. The standard deviation evolution of the depth estimate is presented figure 4, as a function of the standard deviation of the average velocity estimate. This curve is obtained for a distance between the sensors (geophones) of 0.2 m, a deviation of geophones from the plumb (SX) of 0 m, a depth (SZ) of 0.4 m and for different values of the average acoustic velocity (V0) through the basement. We note that to reach an accuracy of 0.1 meter of depth, it is necessary to have an accuracy on the average velocity which varies according to its real value. The lower the value of this velocity is, the more precise it must be estimated.

3

© 2021 by ASME We observe, figure 3, the standard deviation evolution of the depth estimate as a function of the standard deviation of relative delay times estimates. This curve is obtained for a depth (SZ) of 1 m, an average propagation velocity (V0) in the basement of 500 m/s, a deviation on the ground from the plumb (SX) of 0 m. We vary the distance between the sensors. We note that for the same accuracy on the relative delay times, as the distance between the sensors increases, the accuracy of the depth estimation also increases. This observation gives two issues. First, the signal-to-noise ratio degrades the further away the sensors are from the pipe. Second, when the distance between the sensors becomes too great, the signal propagation medium for each sensor is likely to be different. In the future, a compromise must be found.

ii) Statistics on estimators We present results of the estimators tests by Monte Carlo method. The principle is to run a large number of noise draws to characterize one estimator. Here, we evaluate the statistics of the estimators on 1,000 runs of the noise. To simulate the data on which the least squares estimator will run, we compute the theoretical delay times given by our model from the fixed simulation geometry. We noise these theoretical delay times with a centered Gaussian noise of standard deviation 10 -6. To simulate the data on which the MUSIC estimator will run, we construct the received theoretical signals from the previously calculated theoretical delay times. We noise these theoretical signals with a centered Gaussian noise of standard deviation 10 -1 . At the moment, the estimator constructed from the MUSIC algorithm does not give an estimate of the pipe plumb. We test this estimator for cases where the reference sensor 1 is placed at the pipe plumb. We can see, in Table 2, the behavior of the estimator on depth estimation as a function of the position error of the reference sensor. In the future, this estimator will have to be evolved so that it also estimates the pipe plumb deviation. Table 1 shows results obtained for the least squares estimator from the relative delay times. The case 1 considered for this simulation is a depth (SZ) of 0.42 m, a plumb deviation (SX) of 0.05 m, an average propagation velocity (V0) of 500 m/s and sensors are spaced at 0.2 m. We note that for an accuracy of microsecond range on the relative delay times we obtain the desired accuracy on the depth. Tables 3 and4 show, respectively, the results obtained for the least squares estimator and the MUSIC estimator. These results are obtained by considering the same case. The case 2 considered for these simulations is a depth (SZ) of 0.42 m, a plumb deviation (SX) of 0 m, an average propagation velocity (V0) of 500 m/s, and the sensors are spaced at 0.2 m.

EXPERIMENTAL SET-UP i) Presentation of the equipment

We carried out experimental measurements on a semicontrolled test area in order to test methods on real measurements. We have developed an acquisition system (emission (Fig. 5 and7) and multi-reception (Fig. 6 and8) chain) adapted to our experimental needs. The source signal is controlled by computer to adjust signal parameters (time evolution, duration, frequency band, sampling frequency, repetition of the sequence, period, ...). At the reception, acoustic signals are amplified and filtered analogically by the electronic cards. At the input of electronic cards, signals are filtered with a high-pass filter with a cut-off frequency of 100 Hz and at the output with a low-pass filter with a cut-off frequency of 100 kHz. This wide analog filtering allows us to digitally refine the filtering and to have an experimental flexibility. Moreover, the amplification chain is adjustable between 0 dB and 112 dB. Its also allows us to have experimental adaptability.

On figure 7, we can see the emission chain. We use the computer sound card to send the audio signal to an amplifier before it is transmitted to the loudspeaker. Controlling the emission with a computer allows sending different types of signals.

On figure 8, we can see the acquisition chain. Sensors use are geophones with a bandwidth up to 1 000Hz. Each sensor is connected to an electronic card that filters and amplifies the signal. The output electronic card signal is connected to an acquisition card which serves as an analog to digital converter. And that acquisition card is connected to a computer. A battery is required to power the electronic cards.

ii) Experiment to support the assumption of the cylindrical source

We conducted an experiment to support the cylindrical source assumption in our problem. A 0.1 second pulse was injected into the pipe and a geophone was placed at the pipe plumb. The transmitted signal is 100 milliseconds (ms) long, so its energy is spread over these 100 ms. If the energy of the signal transmitted by the pipe is not mainly radial then we expect to observe on the received signal an energy that is spread over a longer duration than the transmitted signal. Indeed, the received signal would be the sum of the signal emitted by the pipe before and after the passage of the signal under the sensor.

5

© 2021 by ASME

Geophones Electronic cards Filtering Amplification

Acquisition cards ADC Computer

On figure 9, we observe the signal received by the sensor. We note that the width of the received signal is approximately 110 ms, which is of the same order as the transmitted signal. On figure 10, we can see the power evolution of the received signal, and we also note that the peak width of the received energy is of the same order as the transmitted signal. We can conclude that the sensor observes the signal from a small portion of pipe closest to it. These results support the hypothesis that the direction of the scattered signal is mainly radial.

iii) Taking measurements to estimate the depth In a real situation, to focus on the depth estimation, we can use, first, estimates on the pipe plumb given by the Gas Tracker. This instrument has been developed several years ago, by MADE-SA company. The principle is to inject a low frequency monochromatic acoustic signal at one end of the pipeline. Information on propagated signals energy distribution can be then obtained from mobile sensor on ground surface. By locating the maximum energy measurement points, we get a pipe plumb estimate.

In this part, we are interested in an experiment on a semi-controlled test area. We qualified the depth and plumbness of the pipe using an invasive electromagnetic method. This method is not applicable to real pipe networks. However, it allows us to obtain an accuracy of the order of a cm on the depth and plumbness of the pipe in our test area. The semi-controlled test area is tarred, a pipe filled of air is buried at 0.42 meter depth and the pipe plumb is known. However, the propagation medium is unknown. To take experimental measurements a line of four geophones is placed perpendicular to the pipe plumb. The sensor 1 is at the pipe plumb and the others are at distances increasingly remote from the plumb. Distance between sensors is 0.2 m and pipe depth is 0.42 m. The source signal injected in the pipe is a sweep from 300 Hz to 1 kHz with a duration of 10 ms repeated every 0.5 s. Signals are synchronously recorded on the four sensors. The acquisition time is 10 seconds with a sampling frequency of 100 kHz. On figure 11, we observe the received signal. On figure 12, signals have been digitally filtered with a bandpass filter between 300 Hz and 1,000 Hz.

EXPERIMENTAL RESULTS & INTERPRETATION

We present the results obtained on real measurements. There are measures for which the algorithms do not work. For these cases where the algorithms do not converge, we can ask the question of signals consistency with the model used. In these cases the model used is not close enough to reality (strongly heterogeneous medium).

In Table 4, we observe the estimation of interest parameters with the least squares estimator (section 3.ii). Results of these estimates are given and compared for the two methods of relative times delays estimates (section 3.i). Using the cross-correlation to estimate the relative times delays, we obtain a depth estimate of 0.18 m. Using the coherence function, we obtain a depth estimate of 0.31 m. The method using the coherence function to estimate the relative times 6

© 2021 by ASME delays improves the depth estimate significantly. However, the depth estimate accuracy still needs to be improved. In table 5, we compare the estimated relative times delays with the theoretical ones according to the estimated velocity (V0 = 420 m/s). We can be seen when sensors are further away, the error of relative times delays estimation is higher. In Table 6, we present the result of the estimation using the MUSIC algorithm adapted to our problem. We obtain an estimate of 0.37 m for a reference value of 0.42 m, that is an error of 0.05 m. We obtained the best results with the high resolution algorithm MUSIC. In the future, we will be interested in evolving the model to cover a larger number of cases. In particular, we will be interested in the trench, which can represent a strong change in the propagation environment (Fig. 13).

CONCLUSION AND PERSPECTIVES

Results obtained are encouraging since we reach the expected accuracy. However, the algorithms do not work on all the measurements we were able to perform and this calls the model into question. It will be interesting to evolve the model, for example, by considering a vertically stratified heterogeneous medium (e.g. the trench). We have an equations number equal to unknowns number. To increase the equations number we need to increase the number of sensors or measurements. Relative delay times estimation with the coherence function shows that the error increases with the distance between sensors. The signal-to-noise ratio are strongly degraded. On the other hand, the depth estimation accuracy is improved by increasing the distance between sensors. We have to find the right compromise. Relative delay times estimate is also related to signal richness. We may increase the bandwidth of sensors or work on the source signal. We plan, at last, to include in MUSIC algorithm, the signal attenuation (e.g. for strongly heterogeneous medium). If it is interesting to evolve the model, we have to take into account the increase of parameters number (e.g. the trench). The trench

Figure 1 :

 1 Figure 1 : Model scheme, Ri geophones (mobile receivers), S the point source (pipe). A specific low frequency acoustic signal [100 -1,000 Hz] is injected at one end of the pipe.

Figure 2 :

 2 Figure 2 : Standard deviation evolution of the depth estimate function of the standard deviation of relative delay times estimates for different depths.

Figure 4 :

 4 Figure 4 : Standard deviation evolution of the depth estimate function of the standard deviation of the average velocity estimate.

Figure 3 :

 3 Figure 3 : Standard deviation evolution of the depth estimate function of the standard deviation of relative delay times estimates for different distances between sensors

Figure 9 :

 9 Figure 9 : Signal measured on a geophone placed at the pipe plumb, the emission is a 0.1 second pulse.

Figure 5 :Figure 8 :

 58 Figure 5 : Emission chain used to take measurements presented. Computer sound card Amplifier Loudspeaker (connected to the pipe)

Figure 7 :

 7 Figure 7 : Photo of the emission chain

Figure 6 :

 6 Figure 6 : Acquisition chain used to take measurements presented.

Figure 11 :

 11 Figure 11 : Signals received during the experiment.

Figure 12 :

 12 Figure 12 : Digitally filtered signals with a bandpass filter between 300 Hz and 1,000 Hz.

Figure 10 :

 10 Figure 10 : Signal power evolution over time, the emission is a 0.1 second pulse.

Figure 13 :

 13 Figure 13 : Scheme of an example of model evolution (with the trench).

  0) R 2 (r 2X ; r 2Z ) R 3 (r 3X ; r 3Z ) R 4 ( r 4X ; r 4Z )

of the ASME 2021 48th Annual Review of Progress in Quantitative Nondestructive Evaluation QNDE2021 July 28-30, 2021, Virtual, Online QNDE2021-74945 Preliminary Acoustic Study of 3D Localization of Buried Polyethylene Pipe. William XERRI 1,2 , Gineth SARACCO 1 , Ghislain GASSIER 1 , Laurent ZOMERO 2 , Philippe PICON 2

1

© 2021 by ASME Proceedings

Table 1 :

 1 Statistics on the least squares estimator by Monte Carlo method, case number 1.

		mean	Standard deviation	true value
	Depth SZ [m]	0,423	0,018	0,42
	Plumb deviation SX [m]	0,0497	0,0016	0,05
	Average velocity V0 [m/s]	498	10,5	500

Table 2 :

 2 Depth estimate by the estimator based on MUSIC algorithm for several plumb deviation.

	Plumb deviation [m]	0	0,01	0,02	0,03	0,04
	Depth estimate [m]	0,416	0,479	0,579	0,704	0,956
	Error [m]	0,004	0,059	0,159	0,284	0,536
	4				© 2021 by ASME

Table 3 :

 3 Statistics on the least squares estimator by Monte Carlo method, case number 2.

		mean	Standard deviation	true value
	Depth SZ [m]	0,4210	0,031	0,42
	Plumb deviation SX [m]	-5.10 -5	0,0048	0
	Average velocity V0 [m/s]	499,8	14,9	500

Table 4 :

 4 Statistics on the estimator based on MUSIC algorithm by Monte Carlo method, case number 2.

		mean	Standard deviation	true value
	Depth SZ [m]	0,414	0,013	0,42
	Average velocity V0 [m/s]	507	10,2	500

Table 4 :

 4 Interest parameters estimation from relative times delays.

		Estimated value	Estimated value
		using cross-	using coherence	Reference value
		correlation	function
	Depth (SZ)	0.18 m	0.31 m	0.42 m
	Distance to the pipe plumb (SX)	0.02 m	0.01 m	Approximately 0
	Average velocity (V0)	490 m/s	420 m/s	unknown

Table 5 :

 5 Comparison estimated relative delay times with the theoretical relative delay times calculated from the estimated propagation velocity.

		Estimated value using coherence function	Theoretical value (depth = 0,42 m and velocity = 420 m/s)
	Between sensors 1 and 2	0.150 ms	0.108 ms
	Between sensors 1 and 3	0.490 ms	0.381 ms
	Between sensors 1 and 4	0.890 ms	0.744 ms

Table 6 :

 6 Estimation with the MUSIC algorithm adapted to our problem.

		MUSIC algorithm adapted	Reference value
	Depth (SZ)	0.37 m	0.42 m
	Average velocity (V0)	450 m/s	unknown

CNRS-UMR7330 CEREGE, AMU, CdF, IRD, Europole de l'Arbois, BP 80, Aix-en-Provence, 13545, France e-mail : xerri@cerege.fr , ginet@cerege.fr , gassier@cerege.fr

MADE-SA, 167 Impasse de la Garrigue, La Farlède, 83210, France e-mail : l.zomero@made-sa.com , p.picon@made-sa.com

© 2021 by ASME

ACKNOWLEDGEMENTS

This project is supported by CIFRE contract between MADE-SA and CNRS-UMR7330, CEREGE, from the ANRT French program of CNRS (Agence National de la Recherche Technique). This project is partially supported by RICE-Research & Innovation Center of Energy, GRTgaz.