A Robust Steered Response Power Localization Method for Wireless Acoustic Sensor Networks in an Outdoor Environment

Abstract

The localization of outdoor acoustic sources has attracted attention in wireless sensor networks. In this paper, the steered response power (SRP) localization of band-pass signal associated with steering time delay uncertainty and coarser spatial grids is considered. We propose a modified SRP-based source localization method for enhancing the localization robustness in outdoor scenarios. In particular, we derive a sufficient condition dependent on the generalized cross-correlation (GCC) waveform function for robust on-grid source localization and show that the SRP function with GCCs satisfying this condition can suppress the disturbances induced by the grid distance and the uncertain steering time delays. Then a GCC refinement procedure for band-pass GCCs is designed, which uses complex wavelet functions in multiple sub-bands to filter the GCCs and averages the envelopes of the filtered GCCs as the equivalent GCC to match the sufficient condition. Simulation results and field experiments demonstrate the excellent performance of the proposed method against the existing SRP-based methods.

1. Introduction

With the rapid development of communication technology and mobile computing devices, applications of wireless acoustic sensor networks (WASNs) are becoming popular in acoustic signal processing. Particularly, WASN-based sound source localization has captured researchers’ attention in the last two decades [1,2,3,4,5]. The existing methods available for passive source localization in WASNs include (1) the received energy-based approaches [6,7,8,9]; (2) the direction of arrival (DOA)-based approaches [10,11]; (3) the time of arrival (TOA)-based approaches [12]; (4) the time difference of arrival (TDOA)-based approaches [13,14,15] and (5) the steered response power (SRP)-based approaches [16,17,18,19,20,21,22].

Most methods require a pre-processing stage in which specific modalities are measured from sensor signals before the location-estimating stage. In contrast, the SRP-based approaches locate the source position or direction by maximizing the power of spatially steered filter and sum beamformer of a group of sensors and contain only one decision step in processing sensor signals to estimate location. Without information compression and disturbances resulting from partial mistakes in the front-end stage, the SRP-based solutions can usually yield more robust performance in noisy and reverberant acoustic environments. Practical implementations commonly use the generalized cross-correlation [23]-based form of the SRP function [16] to reduce computation. The methods similar to the GCC-expression of SRP function are also called a “global coherence field (GCF)” in several references [24,25].

In practice, the primary constraint of the SRP-based approaches is the time-consuming on-grid searching procedure for finding their global maximums. Hence, it has been a hot issue to reduce the computational cost for the SRP-based approaches. In [17], a stochastic region construction (SRC) method is proposed to avoid global grid searching. However, this strategy also causes information loss. In [26], a geometrically sampled grid set based on the TDOA gradient is proposed to improve the SRP performances. An alternative strategy to solve the high-cost searching problem is adopting some adaptive SRP functions regarding the grid resolution to apply a coarse or a hierarchical searching. In [27], the authors use the low-frequency component of GCC for coarse grid resolution and the high-frequency component for fine grids in the SRP-based DOA estimation. In [28], the authors adopt a Gaussian low-pass filter to the GCC for coarse grids. For full-band signals, a similar kind of modification is proposed both in microphone arrays [29] and WASNs [18,19], respectively, in which the spatial spectrum of a given grid is calculated from the sum of the phase-transform weighted GCCs (GCC-Phase Transform (PHAT)s) within a time window containing the TDOA values in the volume surrounding the grid, instead of the original GCC-PHAT in the SRP function.

The SRP-based approaches can provide a robust solution in DOA estimation and source localization tasks in confined spaces. However, they could lose their robustness in an outdoor WASN scenario due to the synthetic effect of the following factors. (1) Grid size, since the monitoring area in outdoor cases may become much more extensive than the area of indoor applications, and the proper searching grids would be much coarser (e.g., meter-level grids outdoors compared with centimeter-level grids indoors). (2) Steering time delay uncertainty; in the classical SRP-based localization frame, the steering time delay at a given position is generated from an ideal propagation model and is always assumed to be entirely right. However, the steering time delay to the source position is different from the actual propagation time. Such a difference becomes no more negligible in the outdoor environment and causes a defocus effect, even though the WASN system is well synchronized. (3) Signal passband; when processing the acoustic data collected in outdoor environments, high-pass or band-pass filtering is indispensable because the environmental noise is intense in the low-frequency range, and the source signals in the real world often possess the band-pass characteristic. The synthetic effect of these three factors would make it difficult to achieve stable localization results. The Modified-SRP functional (MSRP) method introduced in [18,19] provides an elegant solution for scalable grids but it is not suitable for band-pass signals. In [21], the authors elaborate on the SRP in band-pass situations and use the GCC-PHAT envelope or frequency-shifted GCC-PHAT to enhance the robustness in such situations. Nevertheless, the above two methods hardly consider the other two factors (the grid and the steering time uncertainty). In [30], the authors propose a Frequency-Sliding GCC (FSGCC) method, which uses singular value decomposition (SVD) or weighted SVD (WSVD) on the FSGCC matrix and can intelligently extract time delay information of the source signal from multiple sub-band GCCs. The authors adopt the WSVD-FSGCC to the MSRP functional for source localization. This solution can provide excellent localization performance in the band-pass situation with scalable grids. However, in outdoor applications, the high computation cost of the SVD of giant matrices is inevitable due to the long GCC range.

Previously, several common acoustic source placements have been proposed in outdoor scenarios. They mostly focus on localizing the source from TDOA [31] and DOA [32,33] measurements. Some uncertainties are then introduced by the estimation error of TDOA or DOA estimating algorithms. Moreover, some useful information is also compressed, which results in unstable performance. In this direction, in this paper, a robust SRP-based outdoor source localization problem is discussed.

In this paper, a modified SRP-based method is proposed, in which the systematic influence of the above inevitable factors in outdoor WASNs scenarios is considered. The localization performance is analyzed using the normalized contribution of the signal components in the SRP function. A sufficient condition dependent on the GCC waveform function for robust on-grid SRP-based source localization is derived by geometrical analysis. The SRP function using GCCs satisfying this condition can suppress the disturbances induced by the grid distance and the uncertain steering time delay. A GCC refinement procedure for band-pass GCCs is then designed, which uses the complex wavelet functions in multiple sub-bands to filter the GCC and averages the envelopes of the filtered GCCs as the equivalent GCC to match the sufficient condition. Simulation results and field experiments demonstrate the excellent performance of the proposed method against the existing SRP-based methods.

The rest of this paper is organized as follows. In Section 2, the outdoor SRP-based source localization problem is formulated. Section 3 gives the sufficient condition in brief and introduces the GCC refinement procedure. The results of the simulation and the field experiment are presented in Section 4. Conclusions are given in Section 5.

2. SRP-Based Localization in Outdoor Acoustic Sensor Network

2.1. System Models

We discuss the acoustic source localization problem in an N-dimensional Euclidean space with M distributed microphones ($M>N$). Let $\mathbf{x}\in {\mathbb{R}}^N$ be a spatial coordinate vector. Specifically, define ${\mathbf{x}}_s$ as the source location and ${\mathbf{z}}_m$ as the position of the $m_{th}$ sensor ($m=1,2,\cdots ,M$). Let $s\left(t\right)$ be the source signal in the time domain, and the received signal of the microphone at ${\mathbf{z}}_m$ can be modeled as

$$y_m\left[n\right]=\left[h_m\left(t\right)\ast s\left(t\right)+w_m\left(t\right)\right]\delta (t-n/F_s),$$

(1)

where $h_m\left(t\right)$ is the impulse response function representing the propagation of sound from ${\mathbf{x}}_s$ to ${\mathbf{z}}_m$, the operator “∗” represents the convolution operation, $w_m\left(t\right)$ stands for the additive noise signal, and $\delta (t-n/F_s)$ denotes the sampling process at rate $F_s$. When the multi-path delay and non-linear distortion are neglected, the propagation function in the frequency domain can be simplified as

$$H_m\left(\omega \right)=A_m{e}^{-j\omega t_m},$$

(2)

where $A_m\in \mathbb{R}$ is the amplitude-attenuation factor and $t_m$ is the time delay factor. In the frequency domain Equation (1) can be denoted as

$$Y_m\left(\mathsf{\Omega }\right)=A_{m}S\left(\mathsf{\Omega }\right)e^{-j\mathsf{\Omega }{F}_s{t}_m}+W_m\left(\mathsf{\Omega }\right),$$

(3)

where $\mathsf{\Omega }=\omega /F_s\in [-\pi ,\pi ]$ is the normalized angular frequency, $Y_m\left(\mathsf{\Omega }\right)$ is the discrete-time Fourier transform (DTFT) of $y_m\left[n\right]$, $S\left(\mathsf{\Omega }\right)$ and $W_m\left(\mathsf{\Omega }\right)$ are the Fourier transforms of $s\left(t\right)$ and $w_m\left(t\right)$, respectively.

Let ${\eta }_m\left(\mathbf{x}\right)\in \mathbb{R}$ be the steering time delay function describing the time delay associated with sound propagation from a given location $\mathbf{x}$ to ${\mathbf{z}}_m$. In practice, it is commonly modeled as the sound traveling time going through the line-of-sight (LOS) path with a constant sound speed $v_s$; i.e.,

$${\eta }_m\left(\mathbf{x}\right)=\left|\right|\mathbf{x}-{\mathbf{z}}_m\left|\right|/v_s,$$

(4)

where “$∥.∥$” denotes the Euclidean distance. Note that ${\eta }_m\left(\mathbf{x}\right)$ is not exactly the sound propagation in reality. Then the SRP function, which is defined as the output power of the filtered-and-sum beam-former, is given by:

$$P\left(\mathbf{x}\right)={\int }_{-\pi }^{\pi }{\left|\sum _{m=1}^M{G}_m\left(\mathsf{\Omega }\right)Y_m\left(\mathsf{\Omega }\right)e^{j{\mathsf{\Omega }}_m{F}_s{\eta }_m\left(\mathbf{x}\right)}\right|}^{2}d\mathsf{\Omega },$$

(5)

where $G_m\left(\mathsf{\Omega }\right)e^{j{\mathsf{\Omega }}_m{F}_s{\eta }_m\left(\mathbf{x}\right)}$ is the filter associated with the $m_{th}$ sensor. It can be equivalently expressed in term of GCCs [16]:

$$P\left(\mathbf{x}\right)=2\pi \sum _{l=1}^M\sum _{m=1}^M{R}_{l,m}\left({\eta }_m\left(\mathbf{x}\right)-{\eta }_l\left(\mathbf{x}\right)\right),$$

(6)

where

$$R_{l,m}\left(\tau \right)=\frac{1}{2\pi }{\int }_{-\pi }^{\pi }{\mathsf{\Psi }}_{l,m}\left(\mathsf{\Omega }\right)Y_l\left(\mathsf{\Omega }\right)Y_m^{*}\left(\mathsf{\Omega }\right)e^{j\mathsf{\Omega }{F}_s\tau }d\mathsf{\Omega }$$

(7)

denotes the GCC of the sensor pair $\{l,m\}$, $\tau$ is the time lag, superscript “${(.)}^{*}$” represents the conjugate operation, ${\mathsf{\Psi }}_{l,m}\left(\mathsf{\Omega }\right)=G_l\left(\mathsf{\Omega }\right)G_m^{*}\left(\mathsf{\Omega }\right)$ and denotes the weight function of the associated GCC. Ideally, each $R_{l,m}\left(\tau \right)$ achieves its peak at $\tau =t_m-t_l$ so that the SRP function is supposed to achieve its maximum value at the source position ${\mathbf{x}}_s$, as shown in Figure 1a,b. The Phase Transform (PHAT) weight function

$${\mathsf{\Psi }}_{l,m}^{PHAT}\left(\mathsf{\Omega }\right)=1/\left|Y_l\left(\mathsf{\Omega }\right)Y_m^{*}\left(\mathsf{\Omega }\right)\right|$$

(8)

is widely used in the TDOA- and SRP-based localization applications. The PHAT-weighted GCC is generally referred to as the GCC-PHAT, and the SRP using the GCC-PHAT is generally referred to as the SRP-PHAT.

Removing those irrelevant and repetitive terms in Equation (6), the effective component for source localization can be simplified as

$$P_E\left(\mathbf{x}\right)=\sum _{l=1}^{M-1}\sum _{m=l+1}^M{R}_{l,m}\left({\eta }_m\left(\mathbf{x}\right)-{\eta }_l\left(\mathbf{x}\right)\right)=\sum _{p=1}^{C_M^2}{R}_p\left({\tau }_p\left(\mathbf{x}\right)\right),$$

(9)

where p is the sequence number of the valid sensor pair ${\mathbf{c}}_p=\{l,m\}(l<m)$ and is deduced to be $p=(2M-l)(l-1)/2+m-l$, varying from one to a combinatorial number $C_M^2$; ${\tau }_p\left(\mathbf{x}\right)={\eta }_m\left(\mathbf{x}\right)-{\eta }_l\left(\mathbf{x}\right)$ and can be referred to as the steering TDOA function.

2.2. Problem Formulation

The classical SRP-based localization method often lacks robustness in outdoor scenarios. The steering time delay function ${\eta }_m\left(\mathbf{x}\right)$ in the SRP function is different from the sound propagation in reality denoted as ${\eta }_m^0\left(\mathbf{x}\right)$, and $\mathsf{\Delta }{\eta }_m\left(\mathbf{x}\right)={\eta }_m\left(\mathbf{x}\right)-{\eta }_m^0\left(\mathbf{x}\right)$ is denoted as the steering time-uncertainty function. Similarly, the steering TDOA-uncertainty functions in a pair of sensors can be expressed as

$$\mathsf{\Delta }{\tau }_p\left(\mathbf{x}\right) = \mathsf{\Delta }{\eta }_m\left(\mathbf{x}\right) - \mathsf{\Delta }{\eta }_l\left(\mathbf{x}\right) = {\tau }_p\left(\mathbf{x}\right) - {\tau }_p^0\left(\mathbf{x}\right),$$

(10)

where ${\tau }_p^0\left(\mathbf{x}\right) = {\eta }_m^0\left(\mathbf{x}\right) - {\eta }_l^0\left(\mathbf{x}\right)$, representing the real steering TDOA function for a given sensor pair ${\mathbf{c}}_p$. This term is usually negligible within a confined space, so it has been rarely discussed in classical SRP models. However, in outdoor applications, the sound propagation is much more unpredictable, resulting in enlarged uncertainty with the increase in distances. The steering time uncertainty can easily be influenced by the geography, temperature, wind, and self-localization error among sensors, and then yields a noticeable defocus effect on the SRP map, as shown in Figure 1c. The GCCs would intersect with each other dispersedly around ${\mathbf{x}}_s$.

Since the spatial spectrum generated by the SRP function contains many local extrema and ridged areas, the maximal value of $P\left(\mathbf{x}\right)$ is usually found through a grid-searching process. Consider a uniform sampling grid (USG) case in ${\mathbb{R}}^N$. Define ${\mathbf{X}}_g$ as the set of grid points in the candidate searching region ($\mathcal{V}\in {\mathbb{R}}^N$), and $d_g\in \mathbb{R}$, $N_g\in \mathbb{R}$ as the grid distance and the total number of the grids in ${\mathbf{X}}_g$, respectively, then the estimated on-grid location is formulated as

$${\widehat{\mathbf{x}}}_s=arg\underset{\mathbf{x}\in {\mathbf{X}}_g}{max}P\left(\mathbf{x}\right)=arg\underset{\mathbf{x}\in {\mathbf{X}}_g}{max}{P}_E\left(\mathbf{x}\right).$$

(11)

Note that the localization precision depends on the gird resolution. A more accurate estimation usually requires a smaller $d_g$. This will leads to a larger $N_g$ and significantly increased calculation burden because the number of grids is inversely proportional to the $N_{th}$ power of $d_g$ (i.e., $N_g\propto {\left(d_g\right)}^{-N}$). Hence, the accuracy and feasibility can hardly be balanced in an outdoor WASN system confronting a large search region, for which the minimal grid resolution limited by computing power is much coarser than that in indoor applications. However, most SRP approaches usually work well at subtle grid resolutions, and coarser grid resolution has an undersampled effect, as shown in Figure 1d. The searching process probably would miss the source peak.

It is known that the background noise always dominates at low frequencies in the field environment, and real sound sources often show band-pass characteristics. Thus a band-pass GCC is indeed required. However, the SRP-PHAT with a band-pass source would cause a rippling effect [21], as shown in Figure 1e. The rippling effect does not alter the location of the maximal value of the SRP function. However, it may lead to local extrema and even fake peaks such that the SRP spectrum is susceptible to the two other factors and shows a lack of robustness.

Under the influence of the synthetic effect of the above inevitable factors, the real-world SRP output is illustrated in Figure 1f. It shows that classical SRP implementations hardly deal with all these factors outdoors and yield a divergent localization result.

3. A Robust Outdoor SRP-Based Source Localization Method

3.1. On-Grid SRP-Based Localization Error Bound Condition

It is known that the SRP-based spatial spectra mainly depend on the phase information of the source components. It is always reasonable to assume that the additive noise of sensors is independent of each other and the source signal, and then it has no spatial preference (which means that they have zero mean in the phase domain). Their contributions to the SRP spectrum can be neglected and not related to the grid resolution and the steering time uncertainty. Therefore, only the contribution of the source signal is considered in analyzing the SRP function. With the terms of additive noise $w_m\left(\tau \right)$ neglected, the weight functions ${\mathsf{\Psi }}_p\left(\mathsf{\Omega }\right)$ of the sensor pair ${\mathbf{c}}_p$ usually can be expressed as

$${\mathsf{\Psi }}_p\left(\mathsf{\Omega }\right)=B_p{\mathsf{\Psi }}_0\left(\mathsf{\Omega }\right),$$

(12)

where $B_p\in \mathbb{R}$ is an amplitude-scaling factor irrelevant to the frequency, and ${\mathsf{\Psi }}_0\left(\mathsf{\Omega }\right) = {\mathsf{\Psi }}_0\left(-\mathsf{\Omega }\right)\in \mathbb{R}$ is a real function irrelevant to sensors. Substituting Equation (12) into Equation (7), the GCC $R_p\left(\tau \right)$ can be rewritten as

$$\begin{array}{cc}\hfill R_p\left(\tau \right)& =\frac{B_p{A}_l{A}_m}{2\pi }{\int }_{-\pi }^{\pi }{\mathsf{\Psi }}_0\left(\mathsf{\Omega }\right)S\left(\mathsf{\Omega }\right)S^{*}\left(\mathsf{\Omega }\right)e^{j\mathsf{\Omega }{F}_s(\tau -{\tau }_p^0\left({\mathbf{x}}_s\right))}d\mathsf{\Omega }\hfill \\ \hfill & =\frac{B_p{A}_l{A}_m{C}_0}{2\pi }{R}_0\left(\tau -{\tau }_p^0\left({\mathbf{x}}_s\right)\right),\hfill \end{array}$$

(13)

where $C_0=max\left|{\int }_{-\pi }^{\pi }{\mathsf{\Psi }}_0\left(\mathsf{\Omega }\right)S\left(\mathsf{\Omega }\right)S^{*}{\left(\mathsf{\Omega }\right)}^{j\mathsf{\Omega }{F}_s\tau }d\mathsf{\Omega }\right|$, and

$$R_0\left(\tau \right)=\frac{1}{C_0}{\int }_{-\pi }^{\pi }{\mathsf{\Psi }}_0\left(\mathsf{\Omega }\right)S\left(\mathsf{\Omega }\right)S^{*}\left(\mathsf{\Omega }\right)e^{j\mathsf{\Omega }{F}_s\tau }d\mathsf{\Omega },$$

(14)

is the amplitude-normalized version of the weighted self-correlation function of the source signal $s\left(t\right)$. Hence, each GCC contains the same waveform function $R_0\left(\tau \right)$ with different time-shifting factors ${\tau }_p^0\left({\mathbf{x}}_s\right)$ and amplitude factors $B_p{A}_l{A}_m/C_0$. In practice, the range information in amplitude is usually less stable or accurate than in time delay. Thus, a normalized mapping function representing the contribution of the source component in the SRP function can be constructed as

$$F_E(\mathbf{x},{\mathbf{x}}_s)=\frac{1}{C_M^2}\sum _{p=1}^{C_M^2}{R}_0({\tau }_p\left(\mathbf{x}\right)-{\tau }_p^0\left({\mathbf{x}}_s\right)).$$

(15)

In the above equation, the amplitude factors $B_p{A}_l{A}_m/C_0$ between different sensor pairs are removed. Thus, each pair yields an equal contribution to the SRP function. Note that $F_E\left(\mathbf{x}\right)\in \left[-1,1\right]$ has a definite value range regardless of the sensor number M.

For a given grid distance $d_g\in {\mathbb{R}}_{>0}$, an arbitrary uniform sampling grid set in ${\mathbb{R}}^N$ can be expressed as

$$\mathcal{X}(d_g,{\mathbf{x}}_g^o)= \left\{\mathbf{x}+{\mathbf{x}}_g^o:\mathbf{x}={[n_1{d}_g,\cdots ,n_N{d}_g]}^T;n_1,\cdots ,n_N\in \mathbb{Z}\right\},$$

(16)

where ${\mathbf{x}}_g^o\in {\mathbb{R}}^N$ is the position of the origin of the set. Then the on-grid location estimation is given by

$$\begin{array}{cc}\hfill {\widehat{\mathbf{x}}}_s^g& = arg\underset{\mathbf{x}\in \mathcal{X}(d_g,{\mathbf{x}}_o)}{max}{F}_E\left(\mathbf{x},{\mathbf{x}}_s\right)\hfill \\ \hfill & = arg\underset{\mathbf{x}\in \mathcal{X}(d_g,{\mathbf{x}}_g^o)}{max}\frac{1}{C_M^2}\sum _{p=1}^{C_M^2}{R}_0\left({\tau }_p\left(\mathbf{x}\right)-{\tau }_p\left({\mathbf{x}}_{\mathbf{s}}\right)+\mathsf{\Delta }{\tau }_p\left(\mathbf{x}\right)\right).\hfill \end{array}$$

(17)

It is worth pointing out that the grid resolution, the steering time uncertainty, and band-pass issues are comprehensively considered in the above-simplified SRP function.

The grid issue should be unrelated to the origin position ${\mathbf{x}}_g^o$. In the real world, the uncertainty functions $\mathsf{\Delta }{\tau }_p\left(\mathbf{x}\right)$ are hard to closely describe due to many interference factors, and it is reasonable to assume that they have an upper bound $\mathsf{\Delta }{\tau }_{max}$ (i.e., $\left|\mathsf{\Delta }{\tau }_p\left(\mathbf{x}\right)\right|\le \mathsf{\Delta }{\tau }_{max}$). $\mathsf{\Delta }{\tau }_{max}$ indicates the steering time delay uncertainty level and can be estimated from the environmental and devices’ conditions. Thus, the robustness of the on-grid localization problem can be described as: given a $d_g$ and a $\mathsf{\Delta }{\tau }_{max}$, there exists a $\epsilon \in (0,\infty )$ such that

$$\parallel {\widehat{\mathbf{x}}}_s^g-{\mathbf{x}}_s\parallel \le \epsilon .$$

(18)

Define a level-passed area based on $F_E\left(\mathbf{x},{\mathbf{x}}_s\right)$:

$$\mathcal{M}(\alpha ,{\mathbf{x}}_s)\triangleq \{\mathbf{x}:F_E(\mathbf{x},{\mathbf{x}}_s)\ge \alpha \}\subseteq {\mathbb{R}}^N,$$

(19)

where $\alpha \in \mathbb{R}$ is the level-pass threshold. Then a sufficient condition can be obtained in the following Proposition:

Proposition 1.

if $\mathcal{M}(\alpha ,{\mathbf{x}}_s)\cap \mathcal{X}(d_g,{\mathbf{x}}_g^o)\ne \varnothing$ and $\mathcal{M}(\alpha ,{\mathbf{x}}_s)$ is a bounded set (i.e., there exists a ${\epsilon }_M\in (0,\infty )$ such that $\parallel {\mathbf{x}}_1-{\mathbf{x}}_2\parallel \le {\epsilon }_M$ for all ${\mathbf{x}}_1,{\mathbf{x}}_2\in \mathcal{M}(\alpha ,{\mathbf{x}}_s)$), then Inequality (18) is satisfied.

The proof is given in Appendix A.1. Thus, the robustness of the on-grid source localization problem can be analyzed in terms of $\mathcal{M}\left(\alpha ,{\mathbf{x}}_s\right)$.

A practical example of $\mathcal{M}\left(\alpha ,{\mathbf{x}}_s\right)$ is depicted in Figure 2, and its area shrinks inwards when $\alpha$ increases. The first sub-condition ($\mathcal{M}(\alpha ,{\mathbf{x}}_s)\cap \mathcal{X}(d_g,{\mathbf{x}}_g^o)\ne \varnothing$) can be satisfied when $\mathcal{M}\left(\alpha ,{\mathbf{x}}_s\right)$ covers enough areas. The shape of $\mathcal{M}(\alpha ,{\mathbf{x}}_s)$ relates to $\alpha$, $R_0\left(\tau \right)$, $\mathsf{\Delta }{\tau }_p\left(\mathbf{x}\right)$, and sensor distribution, and it is generally irregular. Consider a closed ball ${\mathcal{B}}^N({\mathbf{x}}_0,r)\triangleq \left\{\mathbf{x}:|\mathbf{x}-{\mathbf{x}}_0|\le r;{\mathbf{x}}_0,\mathbf{x}\in {\mathbb{R}}^N\right\}$ with center ${\mathbf{x}}_0$ and radius r. If

$$r\ge d_g\sqrt{N}/2,$$

(20)

then ${\mathcal{B}}^N({\mathbf{x}}_0,r)\cap \mathcal{X}(d_g,{\mathbf{x}}_g^o)\ne \varnothing$ is satisfied. The proof can be seen in Appendix A.2. Consequently, if ${\mathcal{B}}^N({\mathbf{x}}_{\mathbf{s}},d_g\sqrt{N}/2)\subseteq \mathcal{M}(\alpha ,{\mathbf{x}}_s)$, then the first sub-condition is satisfied.

Figure 3 illustrates a typical waveform of $R_0\left(\tau \right)$, the GCC-PHAT of the passband $\left[{{\mathsf{\Omega }}_C-\mathsf{\Omega }}_B,{{\mathsf{\Omega }}_C+\mathsf{\Omega }}_B\right] \subset \left(0,\pi \right]$, which can be expressed by

$$R_0^{PHAT-BP}\left(\tau \right)=\mathrm{sinc}\left(\frac{{\mathsf{\Omega }}_B{F}_s}{\pi }\tau \right)cos\left({\mathsf{\Omega }}_C{F}_s\tau \right).$$

(21)

A valid $R_0\left(\tau \right)$ is an even and bounded function (i.e., $R_0\left(\tau \right) = R_0\left(-\tau \right)$ and $R_0\left(\tau \right) \in \left[-1,1\right]$) and contains a main-lobe around $\tau =0$, where its maximum $a_m$ lies. The maximum side-lobe height (or the maximum value outside the main-lobe area if $R_0\left(\tau \right)$ has no side-lobes) can be denoted as $a_s$, where $a_s<a_m$.

Let us define a function based on $R_0\left(\tau \right)$ by

$$T_R\left(a_T\right)\triangleq inf\left\{\right|\tau |:R_0\left(\tau \right)<a_T\},$$

(22)

where $a_T\in \left[a_S,a_M\right]$ is the level-pass threshold of GCC, “$inf\{.\}$” represents the infimum. $T_R\left(a_T\right)$ represents the half-width of the level-passed section of $R_0\left(\tau \right)$ within its main-lobe. It follows that $R_0\left(\tau \right) \ge a_T$ if and only if $\tau \in \left(-T_R\left(a_T\right),T_R\left(a_T\right)\right)$.

Based on a geometrical analysis in Appendix A.3, if $R_0\left(\tau \right)$ possesses the following property:

$$T_R\left(\alpha \right)\ge d_g\sqrt{N}/v_s+\mathsf{\Delta }{\tau }_{max},$$

(23)

then $\mathcal{M}(\alpha ,{\mathbf{x}}_s)\supset {\mathcal{B}}^N({\mathbf{x}}_{\mathbf{s}},d_g\sqrt{N}/2)$. Therefore, the first sub-condition can be satisfied.

For all $\alpha$ such that $\alpha >{max}_{\parallel \mathbf{x}\parallel \to +\infty }\left\{F_E(\mathbf{x},{\mathbf{x}}_s)\right\}$, the second sub-condition ($\mathcal{M}(\alpha ,{\mathbf{x}}_s)$ is a bounded set) is satisfied. The area of $\mathcal{M}(\alpha ,{\mathbf{x}}_s)$ is mainly the superposition of the projection area of the main-lobe sections of GCCs belonging to individual sensor pairs. Denote

$${\mathsf{\Lambda }}_p({\tau }^c,T)=\{\mathbf{x}:|{\tau }_p\left(\mathbf{x}\right)-{\tau }^c|\le T\}$$

to be the projection area of the TDOA section $\left[{\tau }^c-T,{\tau }^c+T\right]$ of sensor pair ${\mathbf{c}}_p$, where $T\in [0,\infty )$ and ${\tau }^c\in \left[-{\tau }_p^{max},{\tau }_p^{max}\right]$ are the half-width and the central TDOA of the section, respectively, and ${\tau }_p^{max}=\parallel {\mathbf{z}}_l-{\mathbf{z}}_m\parallel /v_s$ is the maximal TDOA value that this sensor pair can produce.

For each sensor pair ${\mathbf{c}}_p$, the solution set of the half hyperbolic equation ${\tau }_p\left(\mathbf{x}\right)={\tau }^c$ can be denoted as ${\mathsf{\Lambda }}_p\left({\tau }^c,0\right)$ and extends to infinity (i.e., there exists an $\mathbf{x}$ such that $\parallel \mathbf{x}\parallel =\infty$ and $\mathbf{x}\in {\mathsf{\Lambda }}_p\left({\tau }^c,0\right)$ ). For two different sensor pairs ${\mathbf{c}}_i$ and ${\mathbf{c}}_j$, if there exist a ${\tau }_i^c\in \left[-{\tau }_i^{max},{\tau }_i^{max}\right]$ and a ${\tau }_j^c\in \left[-{\tau }_j^{max},{\tau }_j^{max}\right]$ such that ${\mathsf{\Lambda }}_i\left({\tau }_i^c,0\right)\subseteq {\mathsf{\Lambda }}_j\left({\tau }_j^c,0\right)$ or ${\mathsf{\Lambda }}_i\left({\tau }_i^c,0\right)⊉{\mathsf{\Lambda }}_j\left({\tau }_j^c,0\right)$, then the half hyperbolic functions ${\tau }_i\left(\mathbf{x}\right)={\tau }_i^c$ and ${\tau }_j\left(\mathbf{x}\right)={\tau }_j^c$ are not independent. The sense might occur when the sensors of these two pairs are co-linear or have the same axis of symmetry; in the meantime, both ${\tau }_i^c$ and ${\tau }_j^c$ reach their extremum or become zero. In WASNs, this case rarely happens because the sensor distributions are often irregular. Despite this sense for all sensor pairs, the maximal value of $F_E\left(\mathbf{x},{\mathbf{x}}_s\right)$ at infinity does not exceed a linear combination of $a_m$ and $a_s$, which is given as

$${\alpha }_{inf}=\frac{C_N^2{a}_m+\left(C_M^2-C_N^2\right)a_s}{C_M^2}.$$

(24)

The detailed derivation can be found in Appendix A.4. If $\alpha >{\alpha }_{inf}$, then $\mathcal{M}(\alpha ,{\mathbf{x}}_s)$ is bounded.

Combining Inequality (23) and Equation (24) together, a sufficient condition for robust on-grid source localization is given by

$$T_R\left({\alpha }_{inf}\right) >d_g\sqrt{N}/v_s+\mathsf{\Delta }{\tau }_{max}.$$

(25)

It means that for a given grid distance $d_g$ and steering TDOA uncertainties within $\mathsf{\Delta }{\tau }_{max}$, if the GCC waveform function $R_0\left(\tau \right)$ has a wide main-lobe satisfying this condition, then the divergent on-grid location estimation can be avoided.

The SRP-PHAT generates a sharp GCC to increase the TDOA resolution for cases with reverberation or multiple sources. However, as shown in Figure 3, the band-pass effect would bring a narrow main-lobe section and strong side-lobes to the GCC waveform function. It can hardly satisfy the requirement Inequality(25), which is also shown by the poor performance of SRP-PHAT in Figure 1f. Next, we will introduce a GCC waveform refinement procedure for the band-pass SRP.

3.2. Robust SRP-Based Source Localization with Refined GCC Waveform

The condition in Inequality (25) is too strict for band-pass GCC situations with coarse grid resolution and perceptible steering TDOA uncertainties. Some classical GCC methods utilized low-pass filtering to meet a broader main-lobe requirement, but they are not applicable for band-pass signals. In this section, the GCC is refined to obtain a suitable waveform to modify the SRP function.

Consider a complex wavelet function ${\psi }_e(\tau ,{\mathsf{\Omega }}_C)=u_e\left(\tau \right)e^{-j{\mathsf{\Omega }}_C{F}_s\tau }$, where $u_e\left(\tau \right)\in {\mathbb{L}}^2\left(\mathbb{R}\right)$ is an even symmetrical function. Applying ${\psi }_e(\tau ,{\mathsf{\Omega }}_C)$ as the filtering function on the GCC-PHAT, the filtered output of ${\mathbf{c}}_p$ can be denoted as

$$R_p^{CF}(\tau ,{\mathsf{\Omega }}_C)=R_p^{PHAT}\left(\tau \right)\ast {\psi }_e(\tau ,{\mathsf{\Omega }}_C),$$

(26)

where $R_p^{PHAT}\left(\tau \right)$ is the GCC-PHAT of ${\mathbf{c}}_p$.

When the real function $u_e\left(\tau \right)$ has an effective support $[-{\mathsf{\Omega }}_B,{\mathsf{\Omega }}_B]\subset [-\pi ,\pi ]$ in the frequency domain, i.e.,

$${\int }_{-\infty }^{\infty }|U_e{\left(\mathsf{\Omega }\right)|}^{2}d\mathsf{\Omega }-{\int }_{-{\mathsf{\Omega }}_B}^{{\mathsf{\Omega }}_B}|U_e{\left(\mathsf{\Omega }\right)|}^{2}d\mathsf{\Omega }\ll {\int }_{-{\mathsf{\Omega }}_B}^{{\mathsf{\Omega }}_B}{\left|U_e\left(\mathsf{\Omega }\right)\right|}^{2}d\mathsf{\Omega },$$

(27)

where $U_e\left(\mathsf{\Omega }\right)$ is the Fourier Transform of $u_e\left(\tau \right)$, and if the source is dominant in the frequency band $[{\mathsf{\Omega }}_C-{\mathsf{\Omega }}_B,{\mathsf{\Omega }}_C+{\mathsf{\Omega }}_B]\subseteq (0,\pi ]$, then the approximation

$$\begin{array}{cc}\hfill R_p^{CF}(\tau ,{\mathsf{\Omega }}_C)& =\frac{1}{2\pi }{\int }_{-\infty }^{\infty }\frac{Y_{p_l}\left(\mathsf{\Omega }\right)Y_{p_m}^{*}\left(\mathsf{\Omega }\right)}{|Y_{p_l}\left(\mathsf{\Omega }\right)Y_{p_m}^{*}\left(\mathsf{\Omega }\right)|}{U}_e(\mathsf{\Omega }-{\mathsf{\Omega }}_C)e^{j\mathsf{\Omega }{F}_s\tau }d\mathsf{\Omega }\hfill \\ \hfill & \approx \frac{1}{2\pi }{\int }_{{\mathsf{\Omega }}_C-{\mathsf{\Omega }}_B}^{{\mathsf{\Omega }}_C+{\mathsf{\Omega }}_B}\frac{Y_{p_l}\left(\mathsf{\Omega }\right)Y_{p_m}^{*}\left(\mathsf{\Omega }\right)}{|Y_{p_l}\left(\mathsf{\Omega }\right)Y_{p_m}^{*}\left(\mathsf{\Omega }\right)|}{U}_e(\mathsf{\Omega }-{\mathsf{\Omega }}_C)e^{j\mathsf{\Omega }{F}_s\tau }d\mathsf{\Omega }\hfill \\ \hfill & \approx \frac{1}{2\pi }{\int }_{{\mathsf{\Omega }}_C-{\mathsf{\Omega }}_B}^{{\mathsf{\Omega }}_C+{\mathsf{\Omega }}_B}{e}^{-j\mathsf{\Omega }{F}_s{\tau }_p^0\left({\mathbf{x}}_s\right)}{U}_e(\mathsf{\Omega }-{\mathsf{\Omega }}_C)e^{j\mathsf{\Omega }{F}_s\tau }d\mathsf{\Omega }\hfill \\ \hfill & \approx \frac{1}{2\pi }{\int }_{-\infty }^{\infty }{e}^{-j\mathsf{\Omega }{F}_s{\tau }_p^0\left({\mathbf{x}}_s\right)}{U}_e(\mathsf{\Omega }-{\mathsf{\Omega }}_C)e^{j\mathsf{\Omega }{F}_s\tau }d\mathsf{\Omega }\hfill \\ \hfill & =u_e(\tau -{\tau }_p^0\left({\mathbf{x}}_s\right))e^{j{\mathsf{\Omega }}_C{F}_s(\tau -{\tau }_p^0\left({\mathbf{x}}_s\right))}\hfill \end{array}$$

(28)

exists. It can be observed that the approximate function carries the same envelope as $u_e\left(\tau \right)$ and extracts the TDOA information in $[{\mathsf{\Omega }}_C-{\mathsf{\Omega }}_B,{\mathsf{\Omega }}_C+{\mathsf{\Omega }}_B]$.

Note that the $R_p^{CF}(\tau ,{\mathsf{\Omega }}_C)$ is equal to the time domain approach of the sub-band GCC defined in [30]. Since the main goal is to obtain an equivalent GCC to match the sufficient condition in Inequality (25), a lightweight approach is to average the envelope of those filtered GCCs of multiple sub-bands in high SNR conditions. According to the power spectral density (PSD) of source signal or other prior knowledge, $N_q$ valid sub-bands can be selected with individual central frequency ${\mathsf{\Omega }}_q$. The final refined GCC is given by

$$R_p^{WR}\left(\tau \right)=\frac{1}{N_q}\sum _q\left|R_p^{CF}(\tau ,{\mathsf{\Omega }}_q)\right|\approx \left|u_e(\tau -{\tau }_p^0\left({\mathbf{x}}_s\right))\right|,$$

which has a specific waveform function $R_0\left(\tau \right)\approx \left|u_e\left(\tau \right)\right|$. Furthermore, the improved spatial function is calculated as

$$P^{WR}\left(\mathbf{x}\right)=\frac{1}{C_M^2}\sum _{p=1}^{C_M^2}{R}_p^{WR}\left({\tau }_p\left(\mathbf{x}\right)\right)=\frac{1}{C_M^2{N}_q}\sum _{p=1}^{C_M^2}\sum _{q=1}^{N_q}|R_p^{CF}({\tau }_p\left(\mathbf{x}\right),{\mathsf{\Omega }}_q)|.$$

(29)

The selection $u_e\left(\tau \right)$ has a significant influence on the refinement of GCC. Its envelope $|u_e\left(\tau \right)|$ provides the waveform function of refined GCCs. The suitable envelope of a suitable $u_e\left(\tau \right)$ should have no side-lobes, i.e., $\left|u_e\left({\tau }_1\right)\right|>\left|u_e\left({\tau }_2\right)\right|\ge 0$ for all $\left|{\tau }_1\right|<\left|{\tau }_2\right|$. Meanwhile, each $U_e\left(\mathsf{\Omega }-{\mathsf{\Omega }}_q\right)$ in the frequency domain serves as a band-pass filter, thus the spectral distribution of $U_e\left(\mathsf{\Omega }\right)$ should be concentrated to satisfy Inequality (27). Gaussian function given by

$$u_e\left(\tau \right)=e^{-{\left({\mathsf{\Omega }}_d{F}_s\tau \right)}^2}$$

(30)

which possesses the required properties both in the time domain and in the frequency domain. Then the corresponding complex filtering function ${\psi }_e\left(\tau ,{\mathsf{\Omega }}_C\right)$ can be regarded as a complex Morlet wavelet. According to (25), for a given grid distance $d_g$ and steering TDOA uncertainty level $\mathsf{\Delta }{\tau }_{max}$, the parameter ${\mathsf{\Omega }}_d$ can be given by

$${\mathsf{\Omega }}_d=v_s\sqrt{-ln\alpha }/\left(F_s{d}_g\sqrt{N}+v_s\mathsf{\Delta }{\tau }_{max}\right),$$

(31)

where N is the space dimension, $\alpha$ is the threshold value, which usually can be set as $\alpha =0.5$. Taking Equation (31) into Inequality (27) and dividing (27) by its right side term, it yields

$$\left[{\int }_{-\infty }^{\infty }{e}^{-{\left(\frac{\mathsf{\Omega }}{2{\mathsf{\Omega }}_d}\right)}^4}d\mathsf{\Omega }-{\int }_{-{\mathsf{\Omega }}_B}^{{\mathsf{\Omega }}_B}{e}^{-{\left(\frac{\mathsf{\Omega }}{2{\mathsf{\Omega }}_d}\right)}^4}d\mathsf{\Omega }\right]/\left({\int }_{-{\mathsf{\Omega }}_B}^{{\mathsf{\Omega }}_B}{e}^{-{\left(\frac{\mathsf{\Omega }}{2{\mathsf{\Omega }}_d}\right)}^4}d\mathsf{\Omega }\right)\ll 1.$$

Thus, the relation of ${\mathsf{\Omega }}_d$ and ${\mathsf{\Omega }}_B$ can be obtained by the following equivalent equation:

$$2{\mathsf{\Omega }}_d\left[{\int }_0^{\infty }{e}^{-{\mathsf{\Omega }}^4}d\mathsf{\Omega }-{\int }_0^{{\mathsf{\Omega }}_B}{e}^{-{\mathsf{\Omega }}^4}d\mathsf{\Omega }\right]/\left({\int }_0^{{\mathsf{\Omega }}_B}{e}^{-{\mathsf{\Omega }}^4}d\mathsf{\Omega }\right)=c,$$

where c is an extremely small number. Then, it can be obtained that

$${\mathsf{\Omega }}_B=2c_e{\mathsf{\Omega }}_d,$$

(32)

where $c_e$ is the positive solution of the following equation:

$$x{E}_{\frac{3}{4}}\left(x^4\right)=\frac{4c}{1+c}\mathsf{\Gamma }\left(\frac{5}{4}\right),$$

where $E_n\left(x\right)={\int }_1^{+\infty }\frac{e^{-xt}}{t^n}dt,(x>0)$ and $\mathsf{\Gamma }\left(x\right)={\int }_0^{+\infty }{t}^{x-1}{e}^{-1}dt,(x>0)$. When c is set as 0.001(−30 dB), $c_e$ in Equation (32) can be obtained as 2.89.

A simulation is performed to illustrate the effect of the GCC waveform refinement procedure on on-grid SRP-based source localization. As shown in Figure 4, the dot-dashed box shows the range of TDOA within the volume of the nearest gird ${\mathbf{x}}_g$, the dashed line with “$\mathsf{\Delta }$” shows the real TDOA, which should coincide with the peak of the GCC; the dotted line with “∇” marks $R_p\left({\tau }_p\left({\mathbf{x}}_g\right)\right)$, corresponding to the nearest gird ${\mathbf{x}}_g$. The $R_p\left({\tau }_p\left({\mathbf{x}}_g\right)\right)$ of the traditional GCC-PHAT is small, thus leading to poor performance in grid searching. In contrast, the proposed refining method generates a smooth waveform and high values throughout the TDOA region indicated by the box in the figure.

The modified algorithm with the GCC refinement procedure is shown in Algorithm 1, in which $u_e\left(\tau \right)=e^{-{\left({\mathsf{\Omega }}_d{F}_s\tau \right)}^2}$ is taken as the target waveform function.

Algorithm 1: SRP with the waveform refinement procedure

Parameter Setting (1) Set the maximum steering TDOA error $\mathsf{\Delta }{\tau }_{max}=\mathsf{\Delta }{\tau }_{max}^C+\mathsf{\Delta }{\tau }_{max}^S$, where the sub-items $\mathsf{\Delta }{\tau }_{max}^C$ and $\mathsf{\Delta }{\tau }_{max}^S$ are determined by the wind and the synchronization error of sensors, respectively. (2) Set the grid distance $d_g$ and searching region $\mathcal{V}$ that meet the system requirement. Then the searching grid set ${\mathbf{X}}_g$ is generated. (3) Set the waveform function $u_e\left(\tau \right)=e^{-{\left({\mathsf{\Omega }}_d{F}_s\tau \right)}^2}$ and $\alpha =$0.5. (4) Set $c=$0.001 and compute the bandwidth ${\mathsf{\Omega }}_B$ using Equation (32). Band selecting (1) Set up the passband $\left[{\mathsf{\Omega }}_L,{\mathsf{\Omega }}_U\right]$ (2) Pick up $N_q$ highest PSD bands of the source or divide the passband uniformly. Source Localization (1) Calculate the refinement waveform (WR)-SRP function $P^{WR}\left(\mathbf{x}\right)$ by Equation (29) at all $\mathbf{x}\in {\mathbf{X}}_g$. (2) Estimate the source location ${\widehat{\mathbf{x}}}_s$ by Equation (11).

4. Experiment Results and Discussion

4.1. Numerical Simulations

In this section, we use Monte Carlo simulations to analyze the efficiency of the proposed SRP-based localization method (the SRP functional with the refinement waveform, referred to as WR), compared with the traditional SRP functional with GCC-PHAT (PS), the SRP functional—the envelope of GCC-PHAT (PES) that is designed for acoustic band-pass signals [21], the modified-SRP (M-SRP) functional with GCC-PHAT (PM) [18] in which grid resolution is considered, and the M-SRP functional with the envelope of GCC-PHAT (PEM) in which both band-pass and grid resolution are considered.

In this setup, M = 8 sensors and one source are randomly deployed in a monitored area of 200 m by 200 m. The propagation model is set to be the line-of-sight path with a constant sound speed of 345 m/s. The input GCCs are generated by the waveform function in Equation (21) with passband of $\left[0.15\pi ,0.4\pi \right]$. The steering TDOA uncertainty $\mathsf{\Delta }{\tau }_p\left(\mathbf{x}\right)$ uniformly distributes over $\left[-\mathsf{\Delta }{\tau }_{max},\mathsf{\Delta }{\tau }_{max}\right]$, where $\mathsf{\Delta }{\tau }_{max}$ is the maximal time uncertainty dependent on the sound-propagation model error and the synchronization error.

We consider four different conditions in WASNs to test the algorithms: (a) a small steering TDOA uncertainty and small grid distance (STSG) condition with $\mathsf{\Delta }{\tau }_{max}=0.1$ ms, $d_g=0.1$ m, (b) a large steering TDOA uncertainty and small grid distance (LTSG) condition with $\mathsf{\Delta }{\tau }_{max}=100$ ms, $d_g=0.1$ m, (c) a small steering TDOA uncertainty and large grid distance (STLG) condition with $\mathsf{\Delta }{\tau }_{max}=0.1$ ms, $d_g=10$ m, (d) a large steering TDOA uncertainty and large grid distance (LTLG) condition with $\mathsf{\Delta }{\tau }_{max}=100$ ms or $d_g=10$ m.

The mean absolute error (MAE) $E\left\{\parallel {\widehat{\mathbf{x}}}_s-{\mathbf{x}}_s\parallel \right\}$ of distance and the cumulative distribution function (CDF) of estimation errors of relative distance are calculated to evaluate the accuracy and robustness of these algorithms, where the relative distance in the cumulative distribution function (CDF) is normalized by the grid distance, i.e.,

$$F\left(e_u\right)=P\left\{\parallel {\widehat{\mathbf{x}}}_s-{\mathbf{x}}_s\parallel /d_g\le e_u\right\},$$

(33)

where $e_u$ is the relative positioning error that is determined as the system requirement. Specifically, the 95th percentile of the localization error in meters is computed as $F^{-1}\left(0.95\right)·d_g$.

The MAE and 95th percentile results are listed in

Table 1. Mean absolute error (MAE) and 95th percentile under different conditions in the simulation.

MAE (m)
Condition	PS	PES	PM	PEM	WR
STSG	0.81	0.07	1.01	0.07	0.06
LTSG	44.53	29.27	52.04	36.37	13.16
STLG	51.90	15.39	42.97	4.07	4.46
LTLG	77.64	50.74	70.37	22.88	13.65
95th percentile (m)
Condition	PS	PES	PM	PEM	WR
STSG	2.83	0.17	2.99	0.18	0.17
LTSG	123.13	82.61	128.10	118.61	33.43
STLG	147.04	58.81	124.39	7.11	9.24
LTLG	172.37	139.73	163.95	74.07	34.68

. All the localization algorithms can obtain the best estimation accuracy in the STSD condition in which the defocus effect and undersampled effect are slight. When the steering TDOA uncertainty or the grid distance increases, the MAE would increase. However, compared with the PS, PES, PM, and PEM methods, the MAE in the WR has almost the smallest estimate error because all these factors have been considered. The 95th percentile has similar results with the MAE, which indicates that the proposed WR method has a stable localization performance in outdoor conditions.

Figure 5a–d depict the CDF of each algorithm in the range $e_u\in [0.5,100\mathrm{m}/d_g]$ under the four conditions. Specifically, the CDF curves will increase rapidly with the location error in the fine condition, and then the estimate errors are the smallest for all the algorithms in the STSG. The CDF curve will move down as the grid distance $d_g$ and steering TDOA uncertainty $\mathsf{\Delta }{\tau }_{max}$ increase, such as in the LTSG, STLG, and LTLG. Since the steering TDOA uncertainty is not considered in PES and PEM, their descent range of CDF in the SDLG is lower than that in the LDSG. Among these localization algorithms, the CDF of the WR is the highest or very close to the highest (STLG), and the PEM method is better than the PS, PES, and PM. The proposed WR method is very robust even though the condition becomes abominable.

Furthermore, Figure 6 presents the MAE in four situations: (a) fixed small steering TDOA uncertainty (ST) with $\mathsf{\Delta }{\tau }_{max}$ = 0.1 ms, $d_g$ ranges from 0.1 m to 50 m; (b) fixed large steering TDOA uncertainty level (LT) with $\mathsf{\Delta }{\tau }_{max}$ = 100 ms, $d_g$ ranges from 0.1 m to 50 m; (c) fixed small grid distance (SG) with $d_g$ = 0.1 m, $\mathsf{\Delta }{\tau }_{max}$ range from 0.1 ms to 100 ms; (d) fixed large grid distance (LG) with $d_g$ = 10 m, $\mathsf{\Delta }{\tau }_{max}$ range from 0.1 ms to 100 ms. The MAE increases with $d_g$ or $\mathsf{\Delta }{\tau }_{max}$ significantly, and this indicates that the steering TDOA uncertainty and grid distance have a severe influence on the performance of source localization. In each situation, the PS and PM produce larger MAE than the other algorithms when $d_g$ and $\mathsf{\Delta }{\tau }_{max}$ are small because they are not applied to band-pass signals. Since the scalable grid sampling and steering TDOA uncertainty are not considered in the PES, it shows reliable performance only when $d_g\le 1$ m and $\mathsf{\Delta }{\tau }_{max}\le 1$ ms. The PEM considered both grid size and band-pass effect; thus, it achieves the best performance in the small $\mathsf{\Delta }{\tau }_{max}$ case. However, the MAE becomes worse when the influence caused by the steering TDOA uncertainties is more significant than by the grid size. The WR obtains the MAE close to the PEM when $\mathsf{\Delta }{\tau }_{max}$ is small. Moreover, it is the smallest in all the other situations. These results abundantly demonstrate its excellent robust performance.

4.2. Field Experiment

In this experiment, seven nodes are distributed in a park, as shown in Figure 7a,b. Each node consists of a microphone sensor, a Wi-Fi module, and a GPS module for self-localization and time calibration. The monitoring area has the same 200 m × 200 m in addition with a hillock. A portable speaker generates the sound signals at 12 positions inside the area, such as the Gaussian signal (S-G), the whistle of vehicles (S-V) representing an urban source, and birdsong (S-B) representing a field source. The temperature was approximately 30 °C, and the wind speed is slower than 3 m/s. Therefore, in the proposed method $\mathsf{\Delta }{\tau }_{max}$ can be set to be 10 ms fully considering the self-localization error of the sensors and the effect of wind.

The sampling frequency is 10,000 Hz and Figure 7c shows the PSDs of both the background noise and received source signals, which are obtained with the Burg method of 50 order number and 2048 FFT length. The PSDs of the source signals are collected at about 30 m away from the speaker. Because the environmental noise is mainly distributed in the frequency bands below 1500 Hz, the passband is set to be (1500 Hz, 3500 Hz) for all sources. The estimated SNRs are shown in Figure 7d, and the SNRs of the full band (0, 5000 Hz) and of the passband (1500 Hz, 3500 Hz) are plotted in solid lines and dashed lines, respectively. For the three source types, the SNR is improved by 20 dB∼30 dB.

The recorded data are divided into 1242 two-second audio frames. SRP algorithms with full-band and band-pass cross-correlation (referred to as CSF and CSB) are added to analyze the necessity of band-pass signals. The PS and PM are not included since they have been proven unreliable in the simulation. Then the candidate SRP-based locators compared in this sub-section include: (1) SRP with full-band GCC (CSF), (2) SRP with band-pass GCC (CSB), (3) SRP with the envelope of band-pass GCC-PHAT (PES), (4) MSRP with the envelope of band-pass GCC-PHAT (PEM) and (5) WR-SRP with band-pass GCC (WR). A well known TDOA-based localization method [13] (referred to as TC) is also compared as a reference in which the TDOAs are obtained by band-pass GCC-PHATs.

The MAE and the 95th percentile of the localization errors of the TC method and the SRP-based methods with different grid distances ($d_g\in \{0.1,1,10\}$ m) are listed in

Table 2. Mean absolute error (MAE) and 95th percentile under different conditions in the field experiment.

MAE (m)
Condition	TC	CSF	CSB	PES	PEM	WR
no grid	102.2	-	-	-	-	-
$d_g$ = 0.1 m	-	79.2	23.5	7.1	18.7	1.4
$d_g$ = 1 m	-	83.0	33.0	12.6	27.4	2.0
$d_g$ = 10 m	-	93.3	66.0	42.6	46.1	7.2
95th percentile (m)
Condition	TC	CSF	CSB	PES	PEM	WR
no grid	322.8	-	-	-	-	-
$d_g$ = 0.1 m	-	146.5	100.8	53.7	105.0	5.4
$d_g$ = 1 m	-	150.4	113.1	91.6	105.1	6.0
$d_g$ = 10 m	-	171.8	149.0	138.5	104.6	21.0

. Moreover the MAEs with grid distance $d_g$ ranging from 0.1 m to 50 m are presented in Figure 8a. Figure 8b–d give the CDF curves at the three grid distances ($d_g\in \{0.1,1,10\}$ m).

Like the simulation, the MAEs increase and the CDF curves move down as the grid distance increases. The MAE of the TC method is the highest because some sensor pairs might produce very severe TDOA measurements in noisy acoustic environments. Its CDF curve also shows that the solution is not stable. By comparing the result of CSF and CSB, the band-pass GCC can significantly enhance the SNR and the localization performance. The PES and PEM obtain more significant localization errors and lack robustness, which indicates the influence of the steering TDOA uncertainty is very remarkable. The proposed WR method achieves the best estimation for all the grid distances, which thoroughly verifies its effectiveness.

5. Conclusions

In this work, a novel and robust Steered Response Power (SRP)-based source localization approach is proposed to localize the band-pass source in outdoor WASNs with steering time delay uncertainty and coarser spatial grids. The robustness of on-grid source localization is analyzed by a sufficient condition, in which the relation between GCC signal waveform and on-grid localization error is demonstrated. A band-pass GCC refinement procedure is designed to meet the sufficient condition for enhancing the on-grid source localization performance. The Monte Carlo simulation and field experiment show that the proposed method has a robust performance in outdoor WASNs scenarios, compared with some state-of-the-art SRP-based methods.