An acoustic method for simultaneous condition detection, localization, and classification in air-filled pipes is proposed. The contribution of this work is threefold: (1) a microphone array is used to extend the usable acoustic frequency range to estimate the reflection coefficient from blockages and lateral connections; (2) a robust regularization method of sparse representation based on a wavelet basis function is adapted to reduce the background noise in acoustical data; and (3) the wavelet components are used to localize and classify the condition of the pipe. The microphone array and sparse representation method enhance the acoustical signal reflected from blockages and lateral connections and suppress unwanted higher-order modes. Based on the sparse representation results, higher-level wavelet functions representing the impulse response are used to localize the position of the sensor corresponding to a blockage or lateral connection with higher spatial resolution. It is shown that the wavelet components can be used to train and to test a support vector machine (SVM) classifier for the condition identification more accurately than with a time domain SVM classifier. This work paves the way for the development of simultaneous condition classification and localization methods to be deployed on autonomous robots working in buried pipes.

Buried pipe infrastructure is important to urban life and forms a vital part of many engineering structures for transporting fluids and gases. In the United Kingdom alone, there are over 600 000 km of sewer pipes.1 The US Environmental Protection Agency estimates that water collection systems in the United States have a total replacement value between $1 and $2 trillion. In Europe, buried water pipe networks are much longer and have a much higher replacement value. These networks are aging rapidly and becoming more heavily used due to population growth, increasing demand for water, and climate change, which leads to an increased rate and severity of faults in these pipes. Therefore, reliable techniques for condition monitoring and fault detection are required for the inspection and targeted maintenance of pipe infrastructure.

Autonomous robotic sensing systems working in buried pipes for condition monitoring and fault detection offer the opportunity to capitalise on recent advances in acoustic and ultrasonic sensing techniques.1 Acoustic methods have been investigated for blockage detection and condition assessment in sewage pipes in recent decades.2 These methods are a very attractive alternative to traditional visual closed-circuit television (CCTV) inspection methods because they are rapid and highly efficient computationally. Acoustically reflective artefacts, including blockages, can be localized remotely with respect to the robot position using the time delay of acoustic echoes measured with a microphone.3 In sewer pipes, the power reflection ratio and signal phase measured with the microphone can be used to discriminate between various in-pipe conditions, e.g., blockage, lateral connection, or pipe end.4 

Although acoustic methods are well suited for use on an autonomous robotic platform, they are complicated by the multi-modal sound wave propagation in a partially filled sewer pipe.5 As a result, in this class of applications, it is common to limit the frequency range to the so-called plane wave regime only, i.e., to the range below the first eigen-frequency of the round sewer pipe, f 10 = 0.59 c / 2 R ,6 where c is the sound speed in air and R is the radius of the pipe. In the case of a typical 300 mm sewer pipe, this frequency is 669 Hz for c = 340 m/s. To radiate a sufficient acoustic power in this frequency range, a large powerful speaker is generally required, which is difficult to deploy on a small robot that would operate in a typical sewer. Furthermore, such a low-frequency range limits the condition localization and classification accuracy due to a relatively long wavelength and restricted frequency band. The main contribution of this paper is to overcome this limitation by proposing a new microphone array processing and machine learning method that extends the frequency range well above the first eigen-frequency to achieve much higher spatial resolution for condition detection and classification in sewer pipes.

A microphone array consists of a set of microphones positioned in a specific way to capture the spatial information about the sound field that can be used for various purposes, e.g., spatial filtering, noise reduction, and dereverberation problems for audio processing.7 This paper uses the microphone array located at the same cross section to capture the acoustic signal containing the spatial information about the first four modes. This information is then processed to extract the fundamental mode (plane wave) to enable defect localization and classification. This method makes use of wavelets, which are well suited to reconstruct a transient signal in the presence of background noise. The idea of using wavelets to deal with transient signals is motivated by prior work, such as that of Ferrante et al.,8 who used the wavelet transform to analyze the transient pressure signal for leakage detection, and Owowo and Oyadiji,9 who used wavelets and the soft threshold method to cancel background noise from acoustic signals for leakage detection in an air-filled pipe. However, little or no work has been performed on the use of wavelets to identify and localize conditions in sewer pipes. This paper proposes a sparse representation method that uses a wavelet basis to cancel the background noise, improve the resolution for condition localization, and increase the accuracy of classification between blockages and lateral connections through post-processing of enhanced acoustical data with support vector machine (SVM).

A health monitoring system should seek to answer a number of key questions, including the presence of a fault (i.e., a blockage or lateral connection) and the location of a fault,10 that are needed to target repair or cleanup operation. An advantage of the microphone sensing array we propose here is that it has a dual use for both condition detection and localization. The method we propose for localization uses the Kalman filter, which operates on acoustic features extracted from the signals measured by the microphone array. This makes the approach highly efficient as a single sensing method is used for both tasks.

The structure of this paper is organised as follows. Section II discusses the theory of acoustic wave propagation in a cylindrical pipe and signal processing theory, including wavelets, sparse representation, and SVM. Section III presents the simulation results of the microphone array processing and acoustic reflection from blockages and lateral connections in the pipe. The experimental setup is discussed in Sec. IV. Experimental results for blockage localization and identification are discussed in Sec. V.

The acoustic field in a rigid cylindrical pipe filled with air is the solution of the wave equation written in cylindrical coordinates ( r , θ , z) as illustrated in Fig. 1. A convenient representation of the total acoustic field in the frequency domain is the normal mode decomposition as suggested by Morse and Ingard,6,
p r , θ , z , ω = A m n Ψ m n r , θ e i γ m n z ,
(1)
where ω is the angular frequency, m and n are the mode indices, and Ψ m n is the mode shape function of a duct cross section. The acoustic field in the pipe is multi-modal, i.e., it is a superposition of an infinite number of waves with the modal amplitude A m n and coordinate-dependent shape functions Ψ m n ( r , θ ). For a cylindrical pipe with the radius R, the mode shape function is given by6,
Ψ m n r , θ = cos m θ J m k m n r ,
(2)
where J m · denotes the mth Bessel function. The eigen-number γ m n can be obtained from the equation for zero velocity on a rigid wall pipe,6 yielding
J m k m n r | r = R = 0.
(3)
FIG. 1.

The system of coordinates in a cylindrical pipe.

FIG. 1.

The system of coordinates in a cylindrical pipe.

Close modal
In the above equation, denotes partial derivative with respect to r. z axis wavenumber is given by6 
γ m n = k 0 2 k m n 2 ,
(4)
where k0 is the wavenumber in a free space ( k 0 = ω / c 0, where c0 is sound velocity in air).

Equation (4) predicts the wavenumber for different modes at different frequencies, which means that the sound velocity in each mode (except in the case of plane wave when k 00 = 0 ) is frequency- and mode-dependent. When the free field wavenumber k0 is larger than the eigen-number kmn, or the frequency is above the corresponding eigen-frequency fmn, the particular acoustic mode can propagate along the pipe with relatively little attenuation. Figure 2 shows schematically the angular and radial dependence of the first four mode shapes in the cylindrical pipe. In this figure, the plus or minus corresponds to the sign taken by the modal shape ( Ψ m n) in Eq. (1).

FIG. 2.

An illustration of the behavior of the first four mode shapes in the cylindrical pipe.

FIG. 2.

An illustration of the behavior of the first four mode shapes in the cylindrical pipe.

Close modal
The frequency response function (FRF) between the acoustic pressure p at r , θ , z and input point source at r s , 0 , 0 with volume velocity Q can be calculated from5,6,11
p ( ω ) Q ω = ω ρ 0 π R e 2 m = 0 n = 0 J m k m n r s cos  m θ J m k m n r e i z γ m n δ m 0 + 1 γ m n J m 2 k m n R 1 m k m n R 2 .
(5)
In particular, for the point source located at the centre of the pipe cross section, only the axisymmetric mode can be excited when m =0 in Eq. (5). Similarly, only the axisymmetric modes can be separated when the receiver is located at the centre.
The sound pressure in an infinite pipe with a blockage or lateral connection is equal to the summation between the incident ( p i) and reflected ( p r) components, given by
p r , θ , z , ω = p i + p r = ( A m n e i k z z + R m n A m n e i k z z ) Ψ m n r , θ ,
(6)
where R m n is the reflection coefficient from the artefact. Using the modal orthogonality, the acoustic pressure of mode ( m , n ) can be estimated from6,
P m n z = A m n e i k z z + R m n A m n e i k z z = p r , θ , z , ω Ψ m n r , θ r drd θ / S ,
(7)
where S is the cross-sectional area of the pipe. The modal amplitude P m n can be used to predict the absolute value of the acoustic reflection coefficient R m n,6 
R m n = max P m n ( z ) min P m n max P m n ( z ) + min P m n .
(8)
Using Eq. (8), the reflection coefficient of a blockage and lateral junction can be predicted via numerical simulations, e.g., finite element method (FEM), to compare against the experimental result as shown in Sec. V.
Note that in this paper, only the plane wave [mode (0,0)] is analyzed for the identification and classification of blockages and lateral connections. This is because the acoustic plane wave is not dispersive, i.e., its speed does not depend on the frequency of sound.6 As a result, the behavior of the reflection coefficient of the plane wave is easier to interpret and to use for condition classification and localization. For this mode, the integral in Eq. (7) can be simplified as
P 00 ( z ) = p r , θ , z , ω r drd θ / S .
(9)
In Eq. (9), the accuracy of the integration depends on the number and distribution of the sensing points and on the choice of the frequency range. In practical situations, there are usually not enough sensing points to cover the whole cross section over which the integral in Eq. (9) is taken. This is particularly important at high frequencies, leading to unwanted contribution from higher modes in P 00 ( z ). To remove the unwanted residue of higher modes as well as some background noise from the measured or predicted acoustic pressure, a sparse representation method using wavelet basis function is proposed and studied in this paper.

Acoustic wave in the pipe can be complicated with multiple modes that travel at different velocities. Using the plane wave mode reconstruction method introduced in Sec. II A [i.e., Eq. (9)] a single mode can be extracted across a broad frequency range and used for in-pipe condition detection, localization, and classification. This paper proposes a sparse representation with wavelets for the simplified impulse response to clean up the higher mode residue and to cancel background noise. Different from the wavelet decomposition and soft shrinkage for noise cancelation proposed in Ref. 6, this paper uses sparse wavelet representation to cancel background noise and to clear up some higher mode residue after the plane wave reconstruction with Eq. (9).

There are two main reasons for using the sparse representation. First, it is possible to assume that the acoustic echoes from the pipe artefacts (e.g., blockages/junctions) have relatively short duration time. This means that the impulse response actively measured on the robot is considered to contain a large number of zero components apart from the initial pulse and reflected echo wave packs, which leads to the sufficient sparsity relative to its dimension in time domain. The sparse nature of the impulse responses has been illustrated in Ref. 3 as an example and is explained further in our paper. Second, the impulse response can be written in terms of appropriate basis vectors where only a few vectors are active, hence, reducing the number of the time domain signals required to store for an accurate signal representation in accordance with the Nyquist sampling theorem.12 The basis vectors used in this paper are wavelet functions.

Wavelets is an extended Fourier analysis that makes use of more general orthogonal bases rather than classical sinusoidal functions. This approach is helpful to overcome the uncertainty principle by exploiting a multi-resolution decomposition12 to balance between different time and frequency fidelities in different frequency bands. It is particularly useful for decomposing complex signals that arise from multi-scale processes such as images and audio signals.12 For a given mother wavelet function φ t, the subspace of scale a is generated by the functions (sometimes called child wavelets),12 
φ j , q t = a j φ a j t q b ,
(10)
where j , q Z; a >1 defines the scale, and b is any real number that defines the shift. j is denoted as the level number of wavelet. In this paper, the typical Daubechies wavelet family was used with a = 2 and b = 1.
Acoustic impulse response measured in pipes can be represented accurately using a limited number of wavelet basis functions. An acoustic impulse response of the pipe, x may be written as a sparse vector s (containing mostly zeros) in the form of a transform matrix W,12,
x = W s .
(11)
As shown in Fig. 3, the wavelet matrix W is generated using Symlets wavelet functions, sym4 (i.e., generated by the matlab function @wmpdictionary). Symlets wavelets are a modified version of Daubechies wavelets with increased symmetry.12 The sparsest solution s ̂ (i.e., the solution represented with a minimum number of wavelet basis functions) satisfies the following optimization problem:12,
s ̂ = arg min s 0 s subject to W s x = 0 ,
(12)
where · 0 denotes the 0 pseudo-norm, which is the number of the non-zero components of the vector. This is also referred to as the cardinality of s .12 Effectively, s ̂ is the vector that is composed of the amplitudes of the wavelet basis functions at different levels.
FIG. 3.

(Color online) Fourth order basis functions of sym4.

FIG. 3.

(Color online) Fourth order basis functions of sym4.

Close modal
The optimization problem in Eq. (12) is non-convex, and its solution is usually found by using a brute-force search,12 which can be computationally expensive. Fortunately, it is possible to relax the optimization in Eq. (12) to a convex 1-minimization,12 
s ̂ = arg min s 1 s subject to W s x = 0 ,
(13)
where · 1 denotes the l1-norm, which describes the sum of absolute values of the vector.
A related convex optimization problem is the following:
s ̂ = arg min s W s x 2 2 + λ s 1 ,
(14)
where λ > 0 is a parameter that weights the importance of sparsity. The estimated impulse response can then be represented by
x ̂ = W s ̂ .
(15)
In this study, the SpaRSA algorithm11 was used to solve Eq. (14). The SpaRSA algorithm provides a general framework for solving the sparse representation problem.13 It can be viewed as an accelerated version of the iterative shrinkage thresholding algorithm.14 The regularization parameter λ is usually set as a specific small constant. The SpaRSA algorithm uses an adaptive continuation technique to optimize and update the value of λ for a more efficient convergence.11 Key steps in this algorithm are presented in the  Appendix.

Different levels of wavelets have different frequency components. Wavelets with higher frequency components can provide a higher spatial resolution to the problem of condition and/or robot localization. For example, the fourth order Symlets function, sym4, has five different levels, where s1 corresponds to the lower frequency and s5 to the higher frequency components. Therefore, it is convenient to use the higher level of wavelet domain vector s ̂ to predict the location of the robot with higher precision. Furthermore, wavelet domain vector s ̂ can also be used to fit into the SVM trainer to identify the blockage from junctions.

After the sparse representation of the impulse response, the higher-level wavelet components can be used to localize the robot position along the pipe. Robot localization is the means by which a robot estimates its position with respect to the surrounding environment. Localization is required for robot control and autonomous navigation, reporting the location of conditions detected in a pipe network, and mapping unknown parts of the pipe network. Normally, information from sensing of the robot's motion and surroundings is input into a localization function. In typical robotic applications, vision and range-finding sensors such as scanning lidar are popular means of making perceptions, as they are able to acquire a large amount of information from the arbitrary environment. In the pipe environment, however, these sensors are limited in scope and only able to observe nearby artefacts due to the confined space within the cross section that is very limited compared to a relatively long length of the pipe and scale of the overall pipe network. This limit in scope means that a localization estimate will accumulate uncertainty over time and that the estimate will likely drift from the true robot position. Acoustic echo sensing is able to perceive more distant artefacts in the pipe environment and offers a means of perception that will not cause an accumulation of uncertainty. In a previous study,3 the robot localization has been validated with a speaker and single microphone sensing system using the plane wave below the cut-off frequency. This paper uses the microphone array to extend the frequency range of the signal to localize the robot and artefacts more precisely.

Robot localization typically takes a probabilistic approach,15 where the uncertainty in each measurement is acknowledged, and a localization estimate is the most likely value of the robot's state in the probability distribution computed over all possible states. Many robot localization approaches are derived from a Bayes filter, a mathematical tool that facilitates the incorporation of prior knowledge and measurements to produce a posterior estimate. A practical implementation of this is the Kalman filter described below.

The robot is modelled as moving in one dimension along the z axis of a pipe with the position at the time t given by
z t = z t 1 + u t + w t ,
(16)
where ut is the command robot motion, and wt is additive error in the motion, drawn from a normal distribution with standard deviation Σu,t. After moving, the robot makes an acoustic echo measurement of the distance to N nearby acoustically reflective artefacts, ξ t n, described by
ξ t = ξ t 1 , ξ t 2 , , ξ t n , , ξ t N .
(17)
These measurements correspond to artefacts in the pipe environment, the positions of which are to be estimated along with the robot position. Therefore, the state is given by
z t = z t , z 1 , z 2 , , z m , , z M T .
(18)
The state space model is, therefore, given by
z t = A z t 1 + B ( u t + w t ) ,
(19)
where A is the identity matrix, and B = 1 , 0 , 0 , , 0 T.
The Kalman filter estimate is parameterised as a multivariable Gaussian distribution described by the mean vector μt and covariance matrix Σt. The prediction step of the filter incorporates the motion information. The predicted mean and covariance are given by
μ ¯ t = A μ t 1 + B u t ,
(20)
Σ ¯ t = A Σ t 1 A T + Σ u , t ,
(21)
where Σu,t is the motion uncertainty.
The acoustic echoes are only used to measure the distance to some artefacts in the pipe environment, so data association needs to be computed between echo measurements ξ t n and artefacts in the environment zm. The negative log-likelihood πn,m,t of a match is found for each combination of measurement and artefact, by computing
Ψ n , m , t = C m Σ ¯ t C m T + Σ ξ , t ,
(22)
π n , m , t = ( ξ t n ξ ¯ t m ) 2 Ψ n , m , t 1 ,
(23)
where Cm is the state space output vector that describes the relative direction of the robot and artefact m, Σξ,t is the measurement uncertainty, and ξ t m is the expected measurement for the artefact m. The smallest value of πn,m,t is found, which corresponds to the most likely match between measurement and artefact. If the value is larger than a threshold for negative log-likelihood ∏, then it is likely that the measurement corresponds to a previously unobserved artefact, so a new artefact is added to the estimate ξ t. If the value is smaller than the threshold, it is likely that the measurement n corresponds to artefact m.
Once the data association has been computed, the update step of the Kalman filter provides an estimate and uncertainty using the following equations:
K n , m , t = Σ ¯ t C m T Ψ n , m , t 1 ,
(24)
μ t = μ ¯ t + K n , m , t ξ t n ξ ¯ t m ,
(25)
Σ t = I K n , m , t C m Σ ¯ t .
(26)
These two steps, prediction and update, are computed at each time point t to recursively estimate the robot position and position of surrounding artefacts using the previous estimate and new information.

This estimation process can be improved through improved acoustic sensing and processing described in this paper. When high level wavelet components are used, as described in Sec. II B, the measurement uncertainty, Σξ,t, and subsequent estimate uncertainty, Σt, will be reduced, and the likelihood of correct data association will be higher, improving robustness. If classification of each artefact in the environment is possible, as described in Sec. II D, it can be incorporated into this estimation of data association. Again, this can improve the robustness of the estimation process.

SVM is a special type of feed-forward neural network used in machine learning. Given a set of training samples consisting of pairs of a co-occurrence feature vector and class of pipe defect labels as ( X m, Y m), in which the subscript m indexes the lengths of the vectors ( Y m ∈ {1,−1}, X mRn), the training approach of SVM is to adjust the weights α and biases b (b ∈ R) to search for an optimal hyperplane and maximum margin. The latter is defined as the distance of the closest vectors in both classes to the hyperplane. In this paper, two different data sets were used to train and test the binary SVM classifier: (1) the time domain impulse response from the echo [i.e., x ̂ in Eq. (15)] is used straightforwardly as Xm; (2) the wavelet domain data s ̂ [in Eq. (15)] associated with the pipe artefacts (e.g., blockage/junctions) is applied as the training and testing input. To find the hyperplane, the following quadratic programming problem was solved:10 
min α , b ϕ α = 1 2 α 2 such that Y m α · X m + b 1 for   all   m = 1 , , N .
(27)
However, it is difficult to straightforwardly solve Eq. (27). Therefore, the Lagrange function was introduced,
max α min α , b L α , b , λ = 1 2 α 2 λ m Y m ( X m · α + b ) 1 ,
(28)
where λ m is Lagrange multipliers. A maximum of the Lagrange function L must satisfy the following conditions by applying a first-order derivation to Eq. (28) respective to b and α as
L α , b , λ b = λ m Y m = 0 , λ m > 0 ,
(29)
L α , b , λ α = α Y m λ m X m = 0 , λ m > 0 .
(30)
The patterns X m for which λ m> 0 are also called support vectors (SVs), which lie exactly on the margins, and all the remaining training samples satisfy automatically their constraints [Eq. (27)]. Substituting the conditions for Eqs. (29) and (30) into the Lagrange function [Eq. (28)], the following dual form is obtained:
max  L p α , b , λ = λ m 1 2 Y m Y n λ m λ n X m · X n ,
(31)
such that
λ m Y n = 0 ; λ n 0 .
(32)
From Eqs. (31) and (32), the optimal solution for λ m can be determined. In pattern recognition, a decision function that correctly classifies the labelled samples ( X m, Y m) is defined as
f α , b X = sgn α · X + b ,
(33)
where
α = λ m Y m X m ; b = Y m λ m Y m X m · X m .
(34)
To make SVM appropriate for real-world decision, radial basis function (RBF) and polynomial kernels are used to nonlinearly map the input data, Xm,
RBF : K X , X m = exp X X m 2 / 2 σ 2 ,
(35)
Polynomial : K X , X m = X T · X m + β d .
(36)
This paper focuses on the binary classification of blockage in sewer pipes, namely, to identify the acoustic echo from a blockage and to discriminate it from other pipe artefacts. For a complete binary classification, a class can be labelled as either +1 or −1 for Ym. On the other hand, the wavelet components of the artefact wave pack can be used as the training and testing input Xm for the SVM classifier. To illustrate the advantage of a wavelet, the time domain artefact wave pack, Xm, was obtained through the following process:
  1. Acquisition of impulse response: Speaker sent a chirp signal and simultaneously recorded the response using six microphones. After deconvolution and bandpass filtering (200–3000 Hz), the six-channel impulse response was obtained, xm, where m = 1:6.

  2. Plane wave construction: Averaging the six-channel impulse response provided the pre-processed plane wave impulse response x.

  3. Denoising and feature extraction using wavelets and sparse representation: After generating the wavelet matrix W using sysm4 level-5 wavelets, plane wave impulse response x was constructed using sparse representation algorithm (SpaRSA) to obtain the denoised signal x ̂ = W s ̂ and wavelet components s ̂ = s ̂ 1 s ̂ 2 s ̂ 3 s ̂ 4 s ̂ 5 T.

  4. Localization: The higher-level wavelet components s ̂ 3 s ̂ 4 s ̂ 5 T were used to represent the higher frequency signal x ̂ h = W s ̂ 3 s ̂ 4 s ̂ 5 T and to apply the Hilbert transform to x ̂ h, where the coordinates of the peaks of the envelope were associated with the location of pipe artefacts referred to the robot. For sequential robotic localization, the measured coordinates were then imported to the ξ t for Kalman filter to predict the locations.

  5. Classification: The wavelet components ( s ̂) associated with the artefacts as the input (Xm) and the label Ym (blockage uses 1 and non-blockage uses −1) were used for SVM training and testing.

Figure 4 shows the processing steps for denoising, localization, and classification.

FIG. 4.

Diagram of the process for condition detection and localization.

FIG. 4.

Diagram of the process for condition detection and localization.

Close modal

This section discusses the analytical and numerical simulations of the microphone array processing used to extract the plane wave from the overall acoustic pressure [Eq. (1)] measured on the microphone array and to estimate the reflection coefficient for an artefact in the pipe. The sensor placement and position uncertainties due to the robot movement in the pipe are discussed in Sec. III A. This can provide evidence in support of the adopted sensor placement strategy for the plane wave reconstruction. The reflection coefficient from blockage and lateral connections will be obtained from numerical simulations to validate the plane wave reconstruction method proposed in the paper via comparison with the experimental results in Sec. V.

A numerical simulation was implemented based on the transfer function described by Eq. (5). The excitation point source was located close to the pipe wall so that all the acoustic modes were excited. The setup used in this numerical simulation is shown in Fig. 5(a). The six virtual microphones were positioned circumferentially and equidistantly spaced at 0.628 R. At this radial position, J 0 k 01 r = 0, so that the amplitude of the first axisymmetric mode is equal to zero as illustrated in Fig. 2. However, the microphones may not be at the ideal positions in a practical situation, e.g., when the robot platform cannot be perfectly located. Therefore, a simulation using slightly shifted microphones (the centres of the six microphones were shifted at a distance of 0.02 R) was also implemented [see Fig. 5(b)].

FIG. 5.

(Color online) (a) A diagram of the simulation setup with the point source excitation and six-microphone array; (b) a diagram explaining the simulated microphone array shift; (c) a comparison of the frequency response function between a single microphone and averaged from the six microphones estimated at the ideal or shifted positions.

FIG. 5.

(Color online) (a) A diagram of the simulation setup with the point source excitation and six-microphone array; (b) a diagram explaining the simulated microphone array shift; (c) a comparison of the frequency response function between a single microphone and averaged from the six microphones estimated at the ideal or shifted positions.

Close modal

As shown in Fig. 5(c), averaging the acoustic sound pressures predicted for the six-microphone array removes the higher modes over the frequency range of 0–5 kHz (within 1 dB fluctuation for plane wave mode) if the microphones were ideally positioned circumferentially at 0.628 R. When the sensors shifted slightly at the distance of 0.02 R, the first four modes were cancelled significantly over the frequency range 0–3.7 kHz (within 1.5 dB fluctuation) apart from the first eigen-frequency around 1.3 kHz. At higher frequencies, the spatial information collected by six microphones tends to be more sensitive to the shifted distance, resulting in a more significant error in the plane wave reconstruction. In this case, the fluctuation in the mode (4, 0) and mode (5, 0) is larger than 3 dB as shown in Fig. 5. Furthermore, it is also observed from additional simulations that the plane wave reconstruction error tends to increase with the microphone shift distance, although the dependence of the error on the microphone shift distance it is not discussed in detail in this paper.

In this paper, the acoustic wave reflection from a blockage or lateral junction was studied using the FEM available in commercial software comsol. Figure 6 shows the simulation setup for sound propagation in the presence of a blockage and lateral connection. The pipe diameter was 0.15 m, which is consistent with that used in the experiments. The height of the blockage was set as 0.6 times the pipe diameter, i.e., h/R =1.2, and it was a diameter long [see Fig. 5(a)]. The maximum mesh size in this numerical study was below 9.5 mm, which corresponded to 1/12 of the acoustic wavelength at 3 kHz. Plane wave excitation was used in the simulations. A perfectly matched layer (PML) was set up at the ends of the pipe to absorb sound to simulate an infinite pipe length. The surface of the blockage in this study was assumed solid, i.e., its acoustic characteristic impedance was much larger than that of air ( Z blockage Z air). The pipe wall was also assumed as rigid.

FIG. 6.

(Color online) Illustration of the simulation using FE modeling for a pipe with (a) a blockage (h/R = 1.2) and (b) a 90° lateral connection.

FIG. 6.

(Color online) Illustration of the simulation using FE modeling for a pipe with (a) a blockage (h/R = 1.2) and (b) a 90° lateral connection.

Close modal

As discussed in Sec. II A, the sound pressure in the plane wave mode was predicted by Eq. (9) using two-dimensional (2D) integration over the cross section. The incident and reflected plane wave interfere with each other, resulting in the fluctuation of sound amplitude in z axis. Using the peak and trough values of the fluctuating acoustic pressure at different axial coordinates, the reflection coefficient can be estimated from Eq. (8). The 2D integration over the cross section was implemented repeatedly with 0.005 m intervals and over 1 m range.

The sound pressure in the plane wave mode, P00, for three frequencies obtained through the simulation and integration in Eq. (9) is shown in Fig. 7. These frequencies were chosen to be between the eigen-frequencies and to illustrate the dependence of the sound pressure as a function of the axial direction. Note that the curves shown in Fig. 7 can be used to determine the amplitude and phase of the complex reflection coefficient for the plane wave mode at these particular frequencies, although only the amplitude of the reflection coefficient is discussed in this paper [see Eqs. (6)–(9)].

FIG. 7.

(Color online) The absolute value of the acoustic pressure of the plane wave mode at three different frequencies as a function of the axial distance in the pipe with a blockage (h/R = 1.2).

FIG. 7.

(Color online) The absolute value of the acoustic pressure of the plane wave mode at three different frequencies as a function of the axial distance in the pipe with a blockage (h/R = 1.2).

Close modal

It is also worth noting that, although the incident wave was a plane wave, the reflection contains higher modes due to the wave scattering at the artefacts (see Fig. 6). There is also a complex relation between the mode number and modal excitation coefficients, depending on the nature of an artefact. Using the integral from Eq. (9), the higher-order modes can then be cancelled as discussed in Sec. III A. Therefore, the amplitude of P00 is a combination of the direct and reflected plane waves. It is frequency-dependent because of the complexity of Eq. (5) and integral (9). The modal pressure oscillates as a function of z with the period determined by the wavelength. As shown in Fig. 7, the amplitude of the plane wave mode at 1 kHz (below the first eigen-frequency for the 0.15 m diameter pipe, f < f10) over the axial direction is significantly higher than the amplitude of the 2 kHz wave (between the first and second eigen-frequency: f10 < f < f20). This can be understood intuitively as a part of the plane wave excitation energy being converted into the higher modes after scattering from the blockage.

The acoustic reflection coefficient can be calculated from Eq. (8) using the predicted spectra of P 00. The reflection coefficient simulation results are discussed in Sec. V and compared against the experimental data.

The acoustic sensing system used in this work consisted of a loudspeaker, six-microphone array, and processor [including power amplifier for loudspeaker, analog-to-digital converter (ADC), digital-to-analog converter (DAC), and Raspberry Pi 4 (Raspberry Pi, Cambridge, UK) for data acquisition] as shown in Fig. 8. This system was installed on a remotely controlled robot [iRobot Looj 330 by iRobot (Bedford, MA)]. The sampling rate was 16 kHz. A bandpass filter with the frequency response of 200–3000 Hz was used to reduce noise. A 100–4000 Hz sweep sine with 10 s duration was used as the excitation signal. The speaker and microphone array were located at the centre of the pipe within 5 mm positional error initially, although this could change due to the robot movement inside the pipe. The radial coordinates of the microphone array were around 60 mm from the pipe centre [see also Fig. 5(a)]. The microphone type used in this test was MSM321A3729H9CP by MEMSensing Microsystems Co., Ltd. (Suzhou, China), and the speaker (Visaton 2242) size was 32 mm diameter driven with a 3 W power supply.

FIG. 8.

(Color online) Robotic platform and acoustic sensing system.

FIG. 8.

(Color online) Robotic platform and acoustic sensing system.

Close modal

In this work, different sizes of blockages were used in a 150 mm diameter polyvinyl chloride (PVC) pipe laid in the iCAIR laboratory at the University of Sheffield. These blockages are described by the ratio of the height of the blockage to the pipe radius h / R = 0, 0.4, 0.8, 1.2, 1.6, and 2, as shown in Fig. 9(a). Figure 9(b) shows an impression of a sandbag blockage in the 150 mm pipe. Other kinds of blockages were also used in the experiment, e.g., acoustic absorbent foam and plastic block, as shown in Figs. 9(c) and 9(d), respectively. To simulate a full 100% blockage, a heavy wooden board was put at the end of the pipe. Efforts were made to seal the circumferential gap between the pipe and board. The straight pipe was constructed from several pipe sections connected with joints at different angles as illustrated in Figs. 9(e) and 9(f). The pipes were not perfectly joined, and joints were not perfectly sealed, so some energy in the acoustic wave was able to reflect and leak out due to the discontinuity at a joint.

FIG. 9.

(Color online) Blockage and lateral junctions simulated in the iCAIR laboratory: (a) different size of concrete blockages; (b) sandbag blockage in the pipe; (c) foam blockage; (d) plastic blockages three-dimensionally (3D) printed; (e) 90° lateral connection; (e) 45°/135° lateral connection.

FIG. 9.

(Color online) Blockage and lateral junctions simulated in the iCAIR laboratory: (a) different size of concrete blockages; (b) sandbag blockage in the pipe; (c) foam blockage; (d) plastic blockages three-dimensionally (3D) printed; (e) 90° lateral connection; (e) 45°/135° lateral connection.

Close modal

The impulse response measured using the acoustic system in the pipe with a blockage and lateral connection is shown in Figs. 10 and 11, respectively. Note that time domain impulse response was converted into the distance domain response by multiplying the time by sound velocity (343 m/s). For a single microphone, the wave dispersion into the higher-order modes significantly complicates the impulse response, causing high frequency noise in the data. This noise can cause difficulties in identifying the condition, particularly when the condition is a small blockage, e.g., h / R = 0.2, as shown in Fig. 10. Averaging the six-microphone data removes the higher-order modes and provides a cleaner signal than that obtained on a single microphone. This is consistent with the theoretical study presented in Sec. II A. Sparse representation cleans up the data further, making it more convenient to apply the localization and classification algorithms detailed in Sec. II. This is even more evident in the case of the data obtained for the pipe with a lateral connection as illustrated in Fig. 11.

FIG. 10.

Impulse response measured in a straight pipe with blockage (h/R = 0.2) located at 4.1 m away from the robot.

FIG. 10.

Impulse response measured in a straight pipe with blockage (h/R = 0.2) located at 4.1 m away from the robot.

Close modal
FIG. 11.

Impulse response measured in a straight pipe with 90° lateral connection located at 4.1 m away from the robot.

FIG. 11.

Impulse response measured in a straight pipe with 90° lateral connection located at 4.1 m away from the robot.

Close modal

As discussed in Sec. II A, the amplitude of the reflection coefficient from a blockage and lateral connection for the plane wave mode can be estimated from Eq. (8). The simulation results for the reflection coefficient were obtained based on the FE modeling in comsol (see Sec. III). For experimental measurements, the impulse response from the blockage and lateral connection (e.g., see the first echo pulse in Figs. 10 and 11, respectively) were extracted using time windowing, zero-padded, and transferred into the frequency domain. A comparison between the simulation results and measurements of the reflection coefficient spectra from blockage and lateral connection with different setups is shown in Figs. 12–14. The measured data were obtained from the denoised plane wave echo using the sparse representation algorithm described in Sec. V A.

FIG. 12.

(Color online) A comparison of the amplitude of the predicted and measured reflection coefficient from blockages of different size (h/2R = 20%, 40%, 60%, 80%, 100%); dashed vertical lines represent the cut-off frequencies.

FIG. 12.

(Color online) A comparison of the amplitude of the predicted and measured reflection coefficient from blockages of different size (h/2R = 20%, 40%, 60%, 80%, 100%); dashed vertical lines represent the cut-off frequencies.

Close modal
FIG. 13.

(Color online) A comparison of the amplitude of reflection coefficient from 90° lateral with different branch sizes between the simulation and measurement.

FIG. 13.

(Color online) A comparison of the amplitude of reflection coefficient from 90° lateral with different branch sizes between the simulation and measurement.

Close modal
FIG. 14.

(Color online) Comparison of the amplitude of reflection coefficient from a 45°/135° lateral (100 mm branch) with different angle between the simulation and measurement.

FIG. 14.

(Color online) Comparison of the amplitude of reflection coefficient from a 45°/135° lateral (100 mm branch) with different angle between the simulation and measurement.

Close modal

Figure 12 presents the reflection coefficient from blockages with different sizes [ 0.2 h / 2 R 0.8; see Fig. 9(a)] obtained from the simulation and experiments. The full blockage was used as a reference, where all the other reflection coefficient were normalized by the full blockage echo obtained experimentally. The simulated results had less than 0.06 average discrepancy with the measurement over the frequency range 200–3000 Hz. As shown in Fig. 12, the acoustic reflection becomes stronger when the blockage size increases. This information can be used to estimate the size of the blockage. Furthermore, the acoustic reflection coefficient becomes larger when the frequency approaches the first cut-off frequency. These properties of the frequency domain spectra can be used to classify blockages from other artefacts such as junctions.

Figure 13 shows the predicted and measured acoustic reflection coefficient spectra from the lateral connection of different diameters attached to the main pipe perpendicularly. The 150 mm diameter branch lateral connection results in a higher reflection coefficient than the 100 mm diameter branch. The lateral connection works like a high-pass filter that allows for the propagation of higher frequency sound wave. The reflection coefficient from a lateral connection drops significantly as the frequency of sound approaches the first cut-off frequency. This highlights the importance of extending the frequency range in the proposed analysis to enable the localization and classification of conditions beyond a lateral connection.

Figure 14 shows the predicted and measured acoustic reflection coefficient from the lateral connection installed at different angles. An increase in the angle of the lateral connection results reduces the sound pressure in the reflected plane wave mode over the frequency range 200–3000 Hz. For a lateral connection, the reflection coefficient reduces significantly as the frequency of sound approaches the first cut-off frequency. This slope in the frequency-dependent reflection coefficient becomes steeper when the lateral angle gets smaller. At higher frequencies beyond the first cut-off frequency, the reflection coefficient for the lateral connection with an increased angle tends to develop a local peak, beyond which it decreases gradually. For example, in the case of a 100 mm lateral connected at 90°, the measured reflection coefficient increases until reaching a peak around 2 kHz and then decreases gradually to almost zero at 2.7 kHz (see Fig. 13).

This section provides the knowledgebase for further identification of blockages and lateral connections according to their acoustic reflection properties. The close agreement (less than 0.06 error on average) between the simulated and measured reflection coefficient for blockages and lateral connections in the frequency domain (300–3000 Hz) also validates the signal processing methods of sparse representation. The denoised blockage/lateral signal is useful for the localization and classification algorithms.

After the reconstruction of plane wave mode using multiple microphone processing and sparse representation method, a robot can process the denoised signal to localize its position with respect to a blockage or lateral connection. In the previous work,3 the Hilbert transform was used to obtain the envelope of the time domain impulse response, where the coordinates of the peaks of the envelope correspond to the relative distance between the pipe artefacts and robot. In this paper, high frequency components of the plane wave impulse response were used for a more precise acoustic localization achieved with microphone array processing. This was accomplished by using the envelope of the higher-level wavelet representation of the impulse response.

As shown in Fig. 15(a), the original plane wave impulse response x after the averaging from six-microphone data can be sparse represented and denoised using wavelets to obtain a clearer signal x ̂ for post-processing. Wavelet components ( s ̂ = s ̂ 1 s ̂ 2 s ̂ 3 s ̂ 4 s ̂ 5 T) and their representation ( x ̂ = W s ̂) for the impulse response are shown in Fig. 15(b). Highest-level wavelet components s ̂ 5 are zero in this signal after the shrinkage from the sparse representation algorithm. The third and fourth levels of wavelet representation show higher resolution in the time domain data than the lower levels (see Fig. 15). This is because the wavelet components s ̂ 3 , s ̂ 4 correspond to the amplitude of the higher frequency signal. After the representation, W s ̂ 3 and W s ̂ 4 contain the acoustic feature wave packs with shorter duration in the time domain, whereas the lower frequency representations, i.e., W s ̂ 1 and W s ̂ 2, provide wider wave packs. Therefore, in this work, the higher frequency representations ( W s ̂ 3 + W s ̂ 4) were used to estimate the location of artefacts in the pipe. Specifically, the locations of artefacts were determined using the coordinates of the peaks of the envelope of the signal at high frequencies W s ̂ 3 + W s ̂ 4. The envelope of the signal W s ̂ 3 + W s ̂ 4 was calculated using the magnitude of its analytic signal, which was computed by filtering W s ̂ 3 + W s ̂ 4 with a Hilbert FIR filter of five points length16 (implemented using the function @envelope in matlab). The envelope results are shown in Fig. 16.

FIG. 15.

(Color online) The wavelet components for the impulse response used for the localization in a pipe with a blockage (h/2R = 0.6) at 4 m.

FIG. 15.

(Color online) The wavelet components for the impulse response used for the localization in a pipe with a blockage (h/2R = 0.6) at 4 m.

Close modal
FIG. 16.

The envelope of (a) the impulse response and (b) the wavelet representation at higher-level (Ws3 + Ws4).

FIG. 16.

The envelope of (a) the impulse response and (b) the wavelet representation at higher-level (Ws3 + Ws4).

Close modal

The blockage was located at 4 m away from the robot with 0.5% prediction error using wavelet representation, whereas the localization error was 4.2% when the impulse response was used. More experiments were carried out using different sizes of blockages and lateral connections. The prediction error of higher-level wavelet representation was below 0.7%. This demonstrates an advantage of using wavelets for the robotic localization with higher accuracy and precision.

Although the pipe artefacts can be localized with respect to the position of the robot, the directions of echo pulses are unknown. Sound intensity measurement can be a solution to determine whether the echo comes from front or back of the robot.17 This paper takes the robot position uncertainty into consideration to use the sequential measurement for the localization of artefacts and pipe mapping. The localization measurement using wavelet can be applied to this sequential robotic localization using the Kalman filter as discussed in Sec. II C, where the robot moves in the pipe and takes the acoustic measurement sequentially every several meters (i.e., every 2 m in this paper).

Photographs of the robotic localization test rig are shown in Fig. 9, and Fig. 17(a) illustrates this rig schematically. In this experiment, a heavy wooden board was installed at the far end of the pipe to represent a full blockage. The robot was moved toward the lateral connection and stopped every 2 m to measure the impulse response. The measured impulse responses are shown in Figs. 17(b) and 17(c) with and without the sparse representation method, respectively. The sparse representation method removes a significant amount of background noise, including dispersive higher modes and some unwanted reflections from the pipe joints.

FIG. 17.

The robotic sensing system in a pipe network for the localization of blockage/lateral connection: (a) an illustration of the experimental system; (b) six-microphone average impulse responses when the robot was at different positions, without sparse representation processing; (c) the impulse responses after the averaging of six-microphone data and with sparse representation process when the robot was at different positions in the pipe.

FIG. 17.

The robotic sensing system in a pipe network for the localization of blockage/lateral connection: (a) an illustration of the experimental system; (b) six-microphone average impulse responses when the robot was at different positions, without sparse representation processing; (c) the impulse responses after the averaging of six-microphone data and with sparse representation process when the robot was at different positions in the pipe.

Close modal

Using the acoustic echoes reflecting from different features in the environment, the robot can estimate its position while it moves along the pipe. Using the process described in Sec. II D, an estimate of the robot's position can be made by combining the prior estimate of the robot's position with traditional odometry and with new acoustic information obtained from the reported measurements. The uncertainty in new information is incorporated into the estimate, so that the process is designed to have robustness to noise in measurements. However, the measurement noise does have an effect on the precision of the estimate, which is investigated here.

Figure 18 shows the results from the simulation of the variation in the mean estimate error for the robot's trajectory (over 100 trajectories) along a pipe obtained for a range of measurement noise levels. This measurement noise was the standard deviation of the Gaussian noise added to each continuous value of distance measurement found from an acoustic echo. For comparison, the estimate made with traditional odometry without using acoustics is also shown. This result illustrates the impact of the uncertainty in robot motion along the pipe. As the acoustic measurement precision increases, the measurement uncertainty decreases, and the median estimate error is seen to decrease from close to the benchmark estimate error of 0.6 m at 1 m of measurement uncertainty to 0.25 m at 0.1 m of measurement uncertainty. This illustrates the strong impact of improved acoustic echo measurement precision on robot localization that can be achieved using the approaches described in this paper.

FIG. 18.

(Color online) Localization error per time step with varying measurement noise.

FIG. 18.

(Color online) Localization error per time step with varying measurement noise.

Close modal

In this work, the two classes of pipe artefacts were identified: (1) a blockage and (2) a lateral connection. The blockages (35 different types in total) shown in Fig. 9 were used for the training and testing of the SVM model. Junctions (25 different types in total) included pipe joints, lateral connections, T-junction, and corner junction. As discussed in Sec. II D, the time domain wave packs x ̂ were used directly as input (Xi) for the training and testing. Wavelet components ( s ̂) associated with the artefacts were also used as the input (Xi), where the classifier was expressed as a wavelet-SVM classifier. In this study, cross-validation was used with five folds, i.e., groups that data samples are split into for the evaluation of SVM modeling to protect against overfitting via data partitioning.

For the wavelet-based SVM classifier, the accuracy of the blockage detection in front of the robot at the first echo was 88% (53/60) based on the provided cases as shown with the confusion matrix in Table I. The time domain SVM classifier enabled us to achieve 78% (47/60) accuracy based on the provided cases. It is worth noting that a linear SVM was also implemented in this study with 65% (39/60) and 53% (32/60) accuracy using the wavelet components and time domain data, respectively. The detailed accuracy, precision, recall, and F1 score12 for these four classifiers are shown in Table II. These results provide the evidence that using wavelet components and non-linear kernel (RBF) improves the classification accuracy relative to using linear kernel and raw time domain data.

TABLE I.

Confusion matrix of the wavelet-SVM classifier from the testing data. TP, true positive; FP, false positive; FN, false negative; TN, true negative.

Predicted blockage Predicted junction
Actual blockage  TP = 31  FP = 4  35 
Actual junction  FN = 3  TN = 22  25 
Total  34  26   
Predicted blockage Predicted junction
Actual blockage  TP = 31  FP = 4  35 
Actual junction  FN = 3  TN = 22  25 
Total  34  26   
TABLE II.

Comparison of time domain linear SVM, wavelet linear SVM, time domain RBF SVM, and wavelet RBF SVM.

Metric Time domain linear SVM Wavelet linear SVM Time domain RBF SVM Wavelet RBF SVM
Accuracy  53%  65%  78%  88% 
Precision  0.571  0.686  0.829  0.886 
Recall  0.625  0.615  0.806  0.912 
F1 score  0.597  0.649  0.817  0.899 
Metric Time domain linear SVM Wavelet linear SVM Time domain RBF SVM Wavelet RBF SVM
Accuracy  53%  65%  78%  88% 
Precision  0.571  0.686  0.829  0.886 
Recall  0.625  0.615  0.806  0.912 
F1 score  0.597  0.649  0.817  0.899 

Eight different settings of pipe network have been used as the demonstration examples, which are shown in Fig. 19. In Fig. 19, time domain SVM classifier can estimate the first artefacts close to the robot accurately, apart from a small blockage (blockage 1), which presents smaller reflection energy, whereas the SVM model using the wavelet components as the training and testing data shows a more accurate (10% accuracy improvement, particularly for small blockages) classification result than the time domain SVM classifier. Furthermore, the reflection from joints/lateral connections 7 m behind the robot can also be predicted by wavelet-SVM classifier with around 92% accuracy. However, the time domain SVM enabled us to achieve only 50% prediction accuracy. This provides evidence that the wavelet method, which takes advantage of the sparsity of the impulse response, can be used to improve the prediction accuracy in comparison to the time domain SVM classification method.

FIG. 19.

(Color online) A comparison of the SVM classifier test between time domain impulse response and wavelet components.

FIG. 19.

(Color online) A comparison of the SVM classifier test between time domain impulse response and wavelet components.

Close modal

Although the prediction using wavelet-SVM classifier tends to be accurate in the testing examples (Fig. 19), there are some cases in the experiment resulting in error classifications:

  1. Small blockages or sound absorbent blockage materials [e.g., acoustic absorption foam in this study; see Fig. 9(c)]. This is because these kinds of blockages do not reflect enough acoustic energy, which leads to a negligibly small amplitude of their impulse response and distortion in the assumed reflection coefficient spectra. The smallest successful blockage identified in this work using wavelet-SVM classifier was 20% blockage (blockage 1 in Fig. 19).

  2. Blockages located close to a junction (<1 m). The acoustic echo from blockage overlaps with the reflection from the junction or another artefact, making it difficult to separate multiple reflections and to classify each of them.

  3. A robot located too close to the blockage or junction (<1 m). This is similar to case 2, where the multiple reflections occur and overlap.

  4. Blockages located behind any artefacts. For example, if the blockage is behind a lateral connection, then the reflected signal can be colored by the presence of this lateral connection. As shown in Fig. 19, the lateral connection after blockage 3 is mistakenly classified.

The classification method provided in this paper uses a limited number of blockage and lateral connection cases simulated in the laboratory. More experimental data and realistic environmental testing will be needed to extend this method for multiple classification of different types of blockages and junctions.

This paper proposed a new acoustic method to simultaneously detect, localize, and identify the conditions in an air-filled pipe. Compared with previous studies, the main novel contributions of this paper are (1) the use of a microphone array to extend the usable acoustic frequency range to estimate the reflection coefficient from blockages and lateral connections, (2) a robust regularization method of sparse representation based on a wavelet basis function adapted to reduce the background noise in the acoustical data, and (3) the use of wavelet components to localize and classify the blockages.

In particular, multiple microphones have been used to reconstruct the plane wave mode beyond the first three eigen-frequencies to support more accurate condition detection, localization, and classification. Numerical and experimental results for the modal reflection coefficient from a blockage and lateral connection have been predicted and compared with measurements. This information has been used to support condition detection and classification.

Wavelet basis functions have been used to sparsely represent the plane wave mode impulse response for the condition detection and classification using the l 1-norm regularization method. The higher-level wavelet functions referring to the higher frequency components of the impulse response have been used to localize the robot and blockage/lateral connection with a higher resolution and accuracy. It has been shown that the wavelet components can also be used to train and to test the SVM classifier for the blockage identification with higher accuracy than using the time domain SVM classifier.

This work is supported by the UK Engineering and Physical Sciences Research Council (EPSRC) Programme Grant No. EP/S016813/1. The authors would like to gratefully thank Gavin Sailor for kindly helping with the design of the robotic platform to support the acoustic sensing system. The authors would also like to gratefully thank Dr. Will Shepherd and Paul Osbourne for kindly helping with the design of the blockages and providing other experimental facilities.

The detailed key steps of the SpaRSA algorithm are presented in Table III.

TABLE III.

The algorithm of sparse reconstruction by separable approximation [SpaRSA (Ref. 18)] for l1-norm regularization.

Task: To solve the problem s ̂ = arg min 1 2 W s x 2 2 + λ s 1 
Input: Response signal e , wavelet dictionary W, parameter λ = 0.001 
Initialization: k =1, A = W, x 1 = x, τ 1 I = A T A, tolerance ε = 10 5 
Iteration: 1. λ k = max 0.2 A T x k , λ
 2. Exploit soft shrinkage: s k + 1 = shrink s k A T A s k x / τ k , λ k / τ k 
(where shrink s i , λ = sign s i max s i λ , 0
 3. Update the step size: τ k = s k + 1 s k T ϑ s k + 1 ϑ s k s k + 1 s k T s k + 1 s k 
 4. If s k + 1 s k s k ε, go to step 5. Otherwise, return to step 2 
 5. s k + 1 = x A s k + 1 
 6. If λ k = λ, stop; Otherwise k = k +1, and return to step 1. 
Output: s ̂ = s k, x ̂ = W s ̂ 
Task: To solve the problem s ̂ = arg min 1 2 W s x 2 2 + λ s 1 
Input: Response signal e , wavelet dictionary W, parameter λ = 0.001 
Initialization: k =1, A = W, x 1 = x, τ 1 I = A T A, tolerance ε = 10 5 
Iteration: 1. λ k = max 0.2 A T x k , λ
 2. Exploit soft shrinkage: s k + 1 = shrink s k A T A s k x / τ k , λ k / τ k 
(where shrink s i , λ = sign s i max s i λ , 0
 3. Update the step size: τ k = s k + 1 s k T ϑ s k + 1 ϑ s k s k + 1 s k T s k + 1 s k 
 4. If s k + 1 s k s k ε, go to step 5. Otherwise, return to step 2 
 5. s k + 1 = x A s k + 1 
 6. If λ k = λ, stop; Otherwise k = k +1, and return to step 1. 
Output: s ̂ = s k, x ̂ = W s ̂ 
1.
Discover water: Treating sewage
,” https://discoverwater.co.uk/treating-sewage (Last accessed 1/17/2023).
2.
Y.
Yu
,
A.
Safari
,
X.
Niu
,
B.
Drinkwater
, and
K. V.
Horoshenkov
, “
Acoustic and ultrasonic techniques for defect detection and condition monitoring in water and sewerage pipes: A review
,”
Appl. Acoust.
183
,
108282
(
2021
).
3.
R.
Worley
,
Y.
Yu
, and
S.
Anderson
, “
Acoustic echo-localization for pipe inspection robots
,” in
Proceedings of the 2020 IEEE International Conference on Multisensor Fusion and Integration for Intelligent Systems (MFI)
, Karlsruhe, Germany (September 14–16,
2020
).
4.
W.
Duan
,
R.
Kirby
,
J.
Prisutova
, and
K. V.
Horoshenkov
, “
On the use of power reflection ratio and phase change to determine the geometry of a blockage in a pipe
,”
Appl. Acoust.
87
,
190
197
(
2015
).
5.
Y.
Yu
,
A.
Krynkin
,
Z.
Li
, and
K. V.
Horoshenkov
, “
Analytical and empirical models for the acoustic dispersion relations in partially filled water pipes
,”
Appl. Acoust.
179
,
108076
(
2021
).
6.
P. M.
Morse
and
K. U.
Ingard
,
Theoretical Acoustics
(
Princeton University
,
Princeton, NJ
,
1986
).
7.
J.
Benesty
,
J.
Chen
, and
Y.
Huang
,
Microphone Array Signal Processing
(
Springer
,
New York
,
2008
), Vol. 1.
8.
M.
Ferrante
,
B.
Brunone
, and
S.
Meniconi
, “
Wavelets for the analysis of transient pressure signals for leak detection
,”
J. Hydraul. Eng.
133
(
11
),
1274
1282
(
2007
).
9.
J.
Owowo
and
S. O.
Oyadiji
, “
Finite element analysis and experimental measurement of acoustic wave propagation for leakage detection in an air-filled pipe
,”
Int. J. Struct. Integr.
8
,
452
467
(
2017
).
10.
C. R.
Farrar
and
K.
Worden
, “
An introduction to structural health monitoring
,”
Philos. Trans. R. Soc. A Math. Phys. Eng. Sci.
365
(
1851
),
303
315
(
2007
).
11.
G.
Kokossalakis
, “
Acoustic data communication system for in-pipe wireless sensor networks
,” Ph.D. dissertation,
Massachusetts Institute of Technology
,
Cambridge, MA
,
2006
.
12.
S. L.
Brunton
and
J. N.
Kutz
,
Data-Driven Science and Engineering: Machine Learning, Dynamical Systems, and Control
(
Cambridge University
,
Cambridge, UK
,
2019
).
13.
Z.
Zhang
,
Y.
Xu
,
J.
Yang
,
X.
Li
, and
D.
Zhang
, “
A survey of sparse representation: Algorithms and applications
,”
IEEE Access
3
,
490
530
(
2015
).
14.
M. A. T.
Figueiredo
and
R. D.
Nowak
, “
A bound optimization approach to wavelet-based image deconvolution
,” in
Proceedings of the EEE International Conference on Image Processing 2005
, Genova, Italy (September 14,
2005
), Vol. 2.
15.
S.
Thrun
,
W.
Burgard
, and
F.
Dieter
,
Probabilistic Robotics
(
MIT
,
Cambridge, MA
,
2006
).
16.
A. V.
Oppenheim
,
R. W.
Schafer
, and
J. R.
Buck
,
Discrete-Time Signal Processing
, 2nd ed. (
Pearson Education India
,
Chennai, India
,
1999
).
17.
M. T. B.
Ali
,
K.
Horoshenkov
, and
S. J.
Tait
, “
Rapid detection of sewer defects and blockages using acoustic-based instrumentation
,”
Water Sci. Technol.
64
(
8
),
1700
1707
(
2011
).
18.
S. J.
Wright
,
R. D.
Nowak
, and
M. A.
Figueiredo
, “
Sparse reconstruction by separable approximation
,”
IEEE Trans. Signal Process.
57
(
7
),
2479
2493
(
2009
).
Published open access through an agreement withThe University of Sheffield