1
Maximum Likelihood Wireless Sensor Network Source Localization Using Acoustic Signal Energy Measurements
Xiaohong Sheng and Yu-Hen Hu Fellow Department of Electrical and Computer Engineering University of Wisconsin - Madison, WI 53706, USA, Email: sheng@ece.wisc.edu, hu@engr.wisc.edu
Abstract A maximum likelihood (ML) acoustic source location estimation method is presented. This method uses acoustic signal energy measurements taken at individual sensors of an ad hoc wireless sensor network to estimate the locations of multiple acoustic sources. Compared to existing acoustic energy based source localization methods, this proposed ML method delivers more accurate results and offers the enhanced capability of multiple source localization. A multi-resolution search algorithm and an expectation- maximization (EM ) like iterative algorithm are proposed to expediate the computation of source locations. The Cramer-Rao Bound (CRB) of the ML source location estimate has been derived. When there is only a single source in the sensor ?eld, the corresponding CRB formulation can be used to analyze the impacts of sensor placement to the accuracy of location estimates. Extensive simulations have been conducted. Empirically, it is observed that this proposed ML method consistently outperforms existing acoustic energy based source localization methods. An example applying this method to track military vehicles using real world experiment data also demonstrates the performance advantage of this proposed method over a previously proposed acoustic energy source localization method. Index Terms Source Localization, Sensor Network, Acoustic Energy, Acoustic Intensity, Maximum Likelihood, Nonlinear Least Square
This project is supported by DARPA under grant no. F 30602-00-2-0555
2
I. I NTRODUCTION The emergence of miniature low-power devices that integrate micro-sensing and actuation with onboard processing and wireless communication capabilities has stimulated great interests in developing wireless ad hoc sensor network [1], [2]. In an ad hoc network, the sensors do not form a structured formation but often are deployed randomly. Here we assume that the senosr locations are known in advance. A wireless ad hoc sensor network often performs monitoring tasks such as detection, classi?cation, localization and tracking of one or more targets in the sensor ?eld. The sensors are typically batterypowered and have limited wireless communication bandwidth. Therefore, ef?cient collaborative signal processing algorithms that consume less energy for computation and less bandwidth for communication are needed [3]. In this paper, we focus on the task of source location estimation using passive and stationary acoustic sensors (microphones) in a wireless ad hoc sensor network. Source localization using acoustic sensors has found numerous applications. In sonar signal processing, the focus is on locating under-water acoustic sources using an array of hydrophones [4]. In video conference and multimedia human computer interface applications, microphone arrays have been developed to locate and track speakers head position in a room environment [6], [7], [8], [9], [10],[11], [12]. Acoustic signatures have also been used to estimate vehicle locations in an open-?eld sensor network [13], [14]. Existing acoustic source localization methods make use of three types of physical measurements: time delay of arrival (TDOA), direction of arrival (DOA) and source signal strength or energy. DOA can be estimated by exploiting the phase difference measured at receiving sensors [15], [16], [17], [18],[19] and is applicable when the acoustic source emits a coherent, narrow band signal. TDOA is suitable for broadband acoustic source localization and has been extensively investigated [5], [6] [13], [14], [20], [21]. It requires accurate measurements of the relative time delay between sensor nodes. It is known that the intensity or equivalently the energy of acoustic signal attenuates as a function of distance from the source. Using this property, recently, an energy-based acoustic source localization method has been reported [22] for locating single target in an open-terrain acoustic sensor ?eld. For wireless sensor network applications, energy (intensity) based acoustic features is an appropriate choice since the acoustic power emitted by targets such as moving vehicles usually varies slowly with respect to time. As such, the acoustic energy time series can be sampled at a much lower rate compared to the raw acoustic time series. Furthermore, the positions of ground moving targets need not be up-
3
dated too often. Therefore, little data will need to be transmitted to a data fusion center via the often congested wireless communication channels. This, in terms, will reduce the energy consumption for data transmission on individual sensor nodes, and will reduce the demand of communication bandwidth over shared wireless channels. Moreover, acoustic signal intensity measurements can also be used to detect the presence of a target. Hence there is no need to compute additional features for the task of source localization. The method to be presented in this paper is based on a maximum likelihood estimation of both the source locations and corresponding acoustic energy readings. Compared to a previously reported method [22], this proposed method promises two signi?cant advantages:
?
The ML method can handle more than one targets within a sensor ?eld. The previously reported method can handle only single target localization.
?
The ML method yields higher accuracy in terms of source location estimates compared to the earlier methods.
In this paper, an acoustic signal energy attenuation model as a function of source-to-sensor distance is ?rst established. Based on this model, the acoustic source localization problem is formulated as a maximum likelihood estimation problem. Two complementary methods then are proposed to solve this nonlinear optimization problem. The ?rst method is based on a projection formulation and uses multiresolution search to expediate computation. The second method is based on an Expectation-Maximization (EM ) like iterative algorithm. This method has less computation complexity but may converge to a local minimum. Next, a formulation of the Cramer-Rao Bound (CRB) of the variance of the location estimates is derived. When there is only a single target, the corresponding CRB can be used to analyze the impact of sensor deployment strategy on the localization accuracy. Speci?cally, the formula clearly veri?ed the conventional belief that smaller estimation error is achieved if the sensors are densely and uniformly laid out over the sensor ?eld. The remaining of this paper is organized as follows: In section II, the acoustic energy attenuation model is developed. In section III, the ML source location estimation method is derived. Both the projection and the EM solution methods are also presented. In section IV, the CRB is derived. In particular, in the single source case, explicit formulation of CRB is derived. This expression allows one to explore the impacts of different sensor deployment strategies to the localization accuracy. Experiments and simulations are provided in section V. Extensive simulation have been conducted to compare the performance of the proposed ML method to three existing energy based source localization methods. Then, time series data
4
samples obtained from a ?eld experiment conducted in California during November 2001 are used to demonstrate the feasibility of using this proposed method to solve real world sensor network source localization problem. II. ACOUSTIC E NERGY ATTENUATION M ODEL Let N be the number of acoustic sensors, and K be the number of acoustic sources in the sensor ?eld. In theory [23], the intensity of an acoustic signal emitted omni-directionally from a point sound source and propagating through ground surface will attenuate at a rate that is inversely proportional to the distance from the source. In this paper, we will further assume that the acoustic intensities of the K sources will be linearly superimposed without any interaction among them. The acoustic signal received at the ith sensor, i = 1, 2, ...N , during time interval n then can be expressed as:
xi (n) = si (n) + νi (n)
(1)
where
si (n) = γi
K k=1
ak (n ? tki ) ρk (n ? tki ) ? ri
(2)
is the acoustic intensity measured at the ith sensor due to all K acoustic sources, and νi (n) is the background noise. In this paper, νi (n) is modeled as a zero-mean additive white Gaussian (AWGN)
2 . a (n ? t ) is the intensity of the k th acoustic source measured noise random variable with variance σn k ki
at one meter from that source, and tki is the propagation delay of the acoustic signal from the k th source to the ith sensor. ak (n ? tki ) will be modeled as a zero-mean random variable uncorrelated with each other. That is, E {ak (n ? tki )} = 0 and E {ak (n ? tki )aj (n ? tki )} = E {ak (n ? tki )} · E {aj (n ? tji )} = 0 if k = j . ρk is an unknown p × 1 vector denoting the position vector of the k th target; and ri is a given
p × 1 vector denoting the position vector of the ith stationary sensor. γi is sensor gain factor of the ith
acoustic sensor. We remark that this model does not account for the potential echoes and reverberation [24] due to the presence of obstacles such as man-made walls or rocky hills. It also does not account for the potential impacts of wind direction [25] [26] and dense vegetation [27] in the sensor ?eld. Moreover, in practice, some of the assumptions made in deriving this model may not be always true and hence may affect its accuracy. For example, the engine sound of a vehicle may not be omni-directional, and will be biased toward the side closer to the engine. The physical size of the acoustic source may be too large compared to the sensor ?eld to be adequately modelled as a point source. In an outdoor environment, strong background noise, including wind-gust may be encountered during operation. In addition, the gain of
5
individual microphones will need to be calibrated to yield consistent acoustic energy readings. Some of these issues will be addressed in a future work. Assume that si (n) and νi (n) are uncorrelated such that E {si (n)νi (n)} = E {si (n)}E {νi (n)} = 0.
2 and S (n ? t ) = E a2 (n ? t ) ): Then, one may represent the acoustic energy as (setting gi = γi k ki ki k K k=1 K
E
s2 i (n)
=
2 γi
E a2 k (n ? tki ) ρk (n ? tki ) ? ri Sk (n ? tki ) ρk (n ? tki ) ? ri
2
= gi
k=1
2
In practice, the expectation is realized using the ensemble average over a time window T = M/fs , where M is the number of sample points used for estimating the acoustic energy received by the sensor during this time interval, and fs is the sampling frequency. Denote this average energy measurements over the time window [t ? T /2, t + T /2] as yi (t), it leads to
yi (t) = 1 T 1 T
t +T / 2 fs
x2 i (n)
t+T /2 fs
T /2 n= t? f s t+T /2 fs
=
s2 i (n)
T /2 n= t ? f s
1 + T
2 νi (n)
(3)
T /2 n= t ? f s
signal energy, ys (t)
noise energy, yn =εi (t)
Since the sources are assumed to be within a sensor network, it is safe to assume that the propagation delay tki is small enough such that a(n ? tki ) ≈ a(n) and ρ(n ? tki ) ≈ ρ(n). This leads to a more concise acoustic energy decay model
yi (t) = ys (t) + εi (t) = gi Sj (t) + εi (t) d2 (t) j =1 ij
K
(4)
where dij (t) = ρj (t) ? ri is the Euclidean distance between the ith sensor and the j th source. The
2 (n), will have a χ2 distribution with mean equals to E [ν 2 (n)] = σ 2 square of the background noise, νi i i 4 /M . If M is suf?ciently large (M >> 30), according to the central limit and variance equals to 2σi 2 , 2σ 4 /M ). For theorem, εi can be approximated well with a normal distribution, namely, εi ? N (σi i 2 , and σ 2 = 2σ 4 /M in later derivations. The validity of this energy convenience, we shall denote ?i = σi i i
attenuation model has been veri?ed with a simple experiment. Details can be found in [22]. III. M AXIMUM L IKELIHOOD L OCATION E STIMATION In this section, we will introduce the ML estimation with different solutions to estimate the source location. Note that the estimation is based on single frame of energy readings from different individual
6
sensors. Let us now de?ne the following matrix notations. The time index t is omitted in the interests of brevity.
T
Z =
(y1 ? ?1 )/σ1 . . . (yN ? ?N )/σN g1 /σ1 g2 /σ2 . . . gN /σN 1/d2 12 1/d2 22 ... ... .. . 1/d2 1K 1/d2 2K
? ? ? ? ? ? ? ? ?
G = diag
?
D =
? ? ? 1/d2 21 ? ? . ? . . ? ?
1/d2 11
. . .
. . .
(5)
2 2 1/d2 N 1 1/dN 2 . . . 1/dN K T
S =
S1 S2 · · · SK
T
H = GD ξ = ξ1 ξ2 · · · ξN
where ξi = ( i ? ?i )/σi ? N (0, 1) are independent Gaussian random variables. Using these notations, equation (4) can be represented as:
Z = GDS + ξ = HS + ξ
(6)
The joint probability density function of Z then can be expressed as:
N 1 f (Z|θ ) = (2π )? 2 exp ? (Z ? HS)T (Z ? HS) 2
(7)
T
where
θ= ρT 1 ρT 2 · · · ρT K S1 S2 · · · SK
is the vector of unknown parameters. The negative log-likelihood function is proportional to a nonlinear quadratic form:
(θ ) = Z ? GDS
2
(8)
Thus the maximum likelihood parameter estimation of θ can be obtained by minimizing (θ ). Equation (8) is a nonlinear least square cost function because the D matrix that contains N K elements is a nonlinear function of the Kp unknown source location coordinates {ρj ; 1 ≤ j ≤ K }. The S vector also contains K unknown parameters. Since there are K (p + 1) unknown parameters, there must be at least K (p + 1) or more sensors reporting acoustic energy readings to yield an unique solution.
7
To minimize (θ ), the solution must lie on a stationary point where
? (θ ) = 0, ?Sk
and
?ρk (θ ) = 0.
for k = 1, 2, · · · K . These conditions lead to the following set of relations
S = H? Z BT G(Z ? HS) = 0
(9) (10)
B1 B2 · · · BN
T
where H? is the pseudo-inverse of the matrix H, and B = with its k th column de?ned as:
Bk =
def
is a K × N matrix
? (DS)T = ?2Sk ? ρk
3 b1k /d3 1k . . . bN k /dN k
(11)
and bij = ?dij /? ρj = ρj ? ri /dij is a unit vector from the j th source to the ith sensor. Equation (10) may be expressed alternatively as:
N K
?ρk (θ ) = 2Sk
i=1
gi σi
ρk ? rk d4 ik
zi ?
gi σi
Sm d2 m=1 im
= 0
(12)
Solving above equation, ρk can be represented as a linearly weighted sum of all sensor locations ri :
N N
ρk =
i=1
αi ri /
i=1
αi
K
(13)
where
αi =
N i=1
gi σi
1 d4 ik
gi zi ? σi
Sm d2 m=1 im
We caution readers that equation (13) is not an explicit expression for ρk since the distance dim is a function of ρm on the right hand side of this equation. Note that equation (9) gives explicit expression of S as a function of H and Z. However, equation (10) is an implicit relation where the unknown {ρk } appear in both sides of the equation. These two equations motivated the development of two different approaches to solve for the maximum likelihood source localization problem as described below:
8
A. Multi-Resolution Projection Solution Substitute equation (9) into equation (8), the variables {Sk }K k=1 can be elimited, giving a modi?ed negative log-likelihood function
L (ρ1 , ρ2 , . . . , ρK ) = ZT (I ? PH )T (I ? PH )Z = ZT (I ? PH )Z = ZT Z ? ZT PH Z
where
PH = H(HT H)?1 HT = UH UT H.
is a projection matrix, and UH is the left singular vectors of the H matrix. The properties PH = PT H and PH · PH = PH have been used during derivation. Since ZT Z are known (normalized) energy measurements, minimizing L is equivalent to maximizing
T L (ρ1 , ρ2 , . . . , ρK ) = ZT UH UT H Z = UH Z 2
(14)
that depends only on source coordinates {ρK k=1 }. In a way, equation (14) is the log-likelihood function with the constraint of equation (9) imposed. As such, L contains K fewer unknown variables than (θ ). We caution the reader that although L is expressed as a quadratic form, maximizing it still requires the solution of a nonlinear least square problem as UH is still a nonlinear function of the source locations
{ρk }K k=1 .
1) Multi-resolution (MR) search: A straight-forward method to ?nd a solution that maximizes L is exhaustive search. However, the computation cost is extremely high especially when there are multiple sources. For example, let there be K sources and q grid points to be searched in each dimension. Then the total number of search points with a p-dimensional sensor ?eld will be equal to q pK . While the computation complexity may be feasible for a desktop computer, it is likely to be overly excessive for low power sensor nodes with limited computing capabilities. This exponentially growing computation complexity can be mitigated with the use of multi-resolution (MR) search method. Among several choices, a logrithmic MR search strategy will examine only w points in each dimension per iteration where q = wm , m being a positive integer. Hence, in each iteration, only wpK grid points needs to be searched. Then, another iteration of search will be con?ned in the neighborhood of the current best solution by sub-dividing the coarser mesh around the current solution into w sub-divisions, and perform search. After m iterations, the MR method will search at a grid size equal to that of exhaustive search. However, the total search points will be m · wpK rather than
q pK = wmpK . To appreciate the amount of saving, let q = 128 = 27 = wm , and p = 2. Using exhaustive
9
search, the total search points will be 214K . Using MR search, the total number of search points will be 7 · 22K . For K = 1, the difference is 4096 points versus 28 points. For K = 2, the difference is 268435456 points versus 112 points. Obviously, due to the coarser search grid at earlier iterations, the MR method may be trapped in a local minimum and yields an inferior solution. B. Expectation-Maximization Iterative Solution Equations (9) and (13) can be used together to yield an iterative solution similar to the expectationmaximization algorithm. With this solution, we assume the source intensity vector S is the missed data rather than unknown parameter. We initiate the K unknown source location ρk at the beginning. During the iteration, we expect the missed data according to equation (9). And then we maximize the loglikelihood function using equation (13) to get the updated estimate ρk . The iteration keeps on going until it convergence. This EM algorithm has much less computation complexity compared to the projection solution. However, it will trap into a local minimum. Hence, they may be applied to re?ne the MR search results. C. Comparison with Other Energy Based Source Localization Methods There are other energy based acoustic source localization methods besides the maximum likelihood method presented in this paper. We now summarize these algorithms and exploit their relationship to the proposed ML method. 1) Closest Point Approach (CPA) Methods: The closest point approach (CPA) is a nevigation term that describes the closest position of two objects moving along non-intersecting straight lines. Here, we borrow the term to refer to a nearest neighbor localization method: Identify a sensor with largest acoustic energy measurement:
i? = arg max yi
i
Assign the source location to be the sensor location
ρ = ri?
(15)
In general, the CPA method is suitable for single source situation. When there are multiple sources, the algorithm must ?rst identify all local maxima among all sensor acoustic energy readings. As a matter of fact, the CPA algorithm can be deduced from equation (13) by setting K = 1. Speci?cally, let di = dik when K = 1, and set di? → 0. Then equation (13) becomes equation (15). Another implementation detail is that the actual source location can not be overlap with the sensor location. Hence, it is often chosen as a location that is nearest to the sensor with largest energy reading.
10
D. Energy Ratios Source Localization, Nonlinear Least Square (ER-NLS) and Least Square (ER-LS) Formulations When there is only a single source (K = 1) within the sensor ?eld, the H matrix become a vector
H= g1 g2 gn , ,···, 2 2 σn d2 σ1 d1 σ2 d2 n
T
,
As such, UH = H/ H will be a unit vector. If each entry of the UH vector were an independent variable, then an obviously solution to maximize the modi?ed cost function L in equation (14) would be
UH = c · Z where c is a proportional constant. Equivalently, for i = 1, · · · , N , yi ? ?i gi gi =c· =c· 2 σi σi ρ ? ri σi di
2
(16)
where ρ denotes the source location. Note that although there are N equalities, there are actually p+1 < N unknowns, including the proportional constant c. One way to solve these set of nonlinear equalities is to ?rst eliminate the unknown constant c by computing the energy ratio ?ij of the ith and the j th sensors as follows:
(yi ? ?i )/(yj ? ?j ) ?ij = gi /gj
?1/2
=
ρ ? ri ρ ? rj
(17)
By sorting the calibrated energy readings (yi ? ?i )/gi , for 0 < ?ij = 1, all the possible source coordinates
ρ that satisfy equation (17) reside on a p-dimensional hyper-sphere described by the equation: ρ ? cij
2 2 = ζij
(18)
where the center cij and the radius ζij of this hyper-sphere associated with sensor i and j are given by:
cij = ri ? ?2 ?ij ri ? rj ij rj , ζij = 2 1 ? ?ij 1 ? ?2 ij
2
(19)
If ?ij = 1, the solution of equation (17) form a hyper-plane between ri and rj , i.e.:
ρ(t)ιij = τij
(20)
where ιij = ri ? rj , τij = |ri |2 ? |rj |2 /2. With all the energy ratios computed, the target location can be solved by minimizing a nonlienar least square cost function:
L1 L2 l=12
J (ρ) =
l1 =1
(|ρ ? cl1 | ? ζl1 )2 +
ιT l2 ρ ? τl2
2
(21)
where L1 and L2 are the number of hyper-spheres and the number of hyper-planes respectively. We call it the energy-ratio nonlinear least square (ER-NLS) method.
11
Since every pair of hyper-spheres (with double indices ij replaced by a single index m for the brevity of notations) ρ ? cm common terms:
(cm ? cn )T ρ = ( cm
2 2 2 ? cn 2 ) ? (ζm ? ζn ) /2 2 2 and ρ ? c = ζm n 2 2 , a hyper plane can be determined by eliminating the = ζn
(22)
Combining equations (22) with (20), the source location can be solved using a least square solution without lengthy nonlinear optimization search. We call it the energy-ratio least square (ER-LS) method. The ER-NLS and ER-LS methods have been reported previously in [22]. They are summarized here for the convenience of readers. IV. C RAMER -R AO B OUNDS (CRB) CRB is a theoretical lower bound of the variance of an unbiased parameter estimate. It is de?ned as the inverse of the Fisher information Matrix:
J = ?E = E ? ?θ ? ln f (Z|θ ) ?θ
T
? ln f (Z|θ ) ?θ
? ln f (Z|θ ) ?θ
T
(23)
Substitute equation (7) into equation (23), one has
J= ? (DS)T T ? (DS) G G ?θ ? θT
(24)
where ? (DS)/? θ T = equation (11). Hence
B D
since ? (DS)/? ST = D, and ? (DS)T /? ρj = Bj as de?ned in
? ?
J=?
?
BT DT
? T ?G G
B
D
(25)
The lower bound of the variance of the source location estimates var ρij can be expressed as
var ρij ≥ J ?1 1 ≤ i ≤ K , and 1 ≤ j ≤ p.
(i?1)p+j,(i?1)p+j
A. CRB for the Single Source Case For convenience, we will focus on the single source (K = 1) and two-dimensional sensor ?eld (p = 2) situation. Furthermore, we set gi = g and σi = σ for i = 1, ..., n. Then, the Fisher Information matrix
12
becomes
?
?
J = ?
? J11 ?
J12 ?
? ?
n T 6 i=1 bi bi /di n 5 i=1 bi /di n 5 i=1 bi /di n 4 i=1 1/di
J21 J22
? ?
=
g2 σ2
4s2
·
?2s
?2s
? ? ?
? ?
De?ne E = J11 ? J12 JT 12 /J22 , it can be factored as:
? ? ?
E = κ?
LT X LT Y
LX LY
?
= κ?
?
|LX
|2
|LX | · |LY | cos β ? |LY |2
?
|LX | · |LY | cos β
(26)
where κ = 4s2 g 2 /σ 2 , LX and LY are both N × 1 vectors, and β = cos?1 LT X LY /|LX | |LY | . The
j th column of the 2 × N matrix LX LY
T
LX LY
T
can be expressed as:
1/d2 j
N i=1 N
=
j
bj ? d3 j
2 1/d2 i=1 i
1 d2 i
bi d3 i
2 = bj /d3 j ? b/dj
(27)
where b is a weighted average of the individual (weighted) source-to-sensor direction vector. Using the matrix inversion formula, the lower bound of the variance of the ML source location estimates can be found in terms of the CRB as:
?1 V ar(ρ ? x ) ≥ {E }11 =
1 1 2 κ |LX | (1 ? cos2 β ) 1 1 κ |LY |2 (1 ? cos2 β )
(28)
?1 V ar(ρ ? y ) ≥ {E }22 =
(29)
Due to the weight factors of 1/d3 i , sensors that are close to the source will have a much higher impact on the CRB than those far away. Speci?cally, the fact that the CRBs are the diagonal elements of E?1 implies that the E matrix must be strongly non-singular. That in turn implies that at least two of the sensors must be close enough to the source so that the [LX LY ]T matrix will have full column rank
(= p).
To reduce the CRB and hence to improve the accuracy of the source localization results, three approaches may be taken: (i) To keep κ large, (ii) To increase both |LX |2 and |LY |2 , and (iii) to let β → 90o . Since source energy emission level s and background noise level σ can not be controlled,
13
to increase κ, the sensor gain g needs to be increased. To increase the magnitude of both |LX |2 and
|LY |2 , the set of sensor target distance {di ; 1 ≤ i ≤ N } must be reduced. This implies that there must be
suf?cient number of sensors close to the acoustic source whenevern the source is within the sensor ?eld. Since the target may be moving within the sensor ?eld, this requirement translates to the requirement of dense and uniform sensor deployment within the sensor ?eld. As for (iii), to make β → 90o , it would require the sensors within the sensor ?eld be laid out in a particular formation. Nonetheless, this requirement is not suitable for an ad hoc sensor network application where sensor deployment is assumed to have no speci?c formation or structure. V. P ERFORMANCE E VALUATION A. Performance Comparison using Simulation – Single Source Case Extensive simulation runs have been conducted to compare the performance of the proposed ML energy based source localization algorithm to other energy based source localization algorithms. We use equation (4) to generate the acoustic energy readings of a two-dimensional (p = 2) sensor ?eld of size 100 meters by 100 meters. Since existing localization algorithms are all for single target situation, we set K = 1. The source location and the sensor locations are randomly chosen from within the sensor ?eld in each run. The source energy is set at S = 5000, and the background noise level is set at σi = 1 for all sensors in the sensor ?eld. Note that although the SNR is 37dB at the source location, the actual SNR at different sensors depends very much on the sensor to source distance. For example, for a sensor that is 50 meters away from the source, its SNR is merely 10 × log10 (5000/502 ) = 3dB . The energy variation
2 , σ 4 /M ) with M = 100. modeled as a Gaussian random variable N (σi i i (t)
is
2000 repeated trials were conducted. In each trial, all four energy based acoustic source localization methods, ML, CPA, ER-NLS, and ER-LS will be used to locate the source location, and the error is recorded. For the CPA method, the target location is the location of the maximum energy sensor added by a Gaussian random variable with zero mean and unit variance. For the ML and ER-NLS methods, exhaustive search is used at a resolution of 5 meters by 5 meters per search grid. Three different sensor densities have been used: 4, 10, 25. In the last case, there is approximately one sensor in every 20 by 20 meter cells in the sensor ?eld. From simulation results, we observe that the mean values of these methods do not show any statistically signi?cant bias, and hence yield unbiased estimates. Furthermore, the estimation error in different dimensions are also uncorrelated. In addition,
14
4 sensors 1 CPA mean = 27.8825 stdev= 16.059 0.5 0.5 1 mean = 17.3775 stdev= 10.3028 0.5 10 sensors 1 mean = 10.75 stdev= 6.3295 25 sensors
0 1 ML
0
20
40
60
80
0 1
0
20
40
60
80
0 1
0
20
40
60
80
mean = 9.2625 stdev= 11.5869 0.5 0.5
mean = 4.1725 stdev= 3.7869 0.5
mean = 3.425 stdev= 2.4687
0 1 ER?NLS
0
20
40
60
80
0 1
0
20
40
60
80
0 1
0
20
40
60
80
mean = 15.4925 stdev= 19.6636 0.5 0.5
mean = 8.44 stdev= 13.9371 0.5
mean = 7.6875 stdev= 12.8784
0 1 ER?LS
0
20
40
60
80
0 1
0
20
40
60
80
0 1
0
20
40
60
80
mean = 22.6775 stdev= 25.529 0.5 0.5
mean = 7.5875 stdev= 11.414 0.5
mean = 7.4325 stdev= 11.7562
0
0
20
40
60
80
0
0
20
40
60
80
0
0
20
40
60
80
Fig. 1.
Comparison of four acoustic energy-based source localization algorithms
?
The ML method consistently outperforms all other three existing methods in terms of minimum estimation error variances and minimum mean estimation error.
?
The CPA method bene?ts most when the number of sensors increases. Its rankings jump from the last the the second best as a result.
?
With both the 10-sensor or the 25-sensor cases, the least square formulation (ER-LS) out-performs the computation-wise more expensive nonlinear least square formulation (ER-NLS).
Another way to analyze the simulation results is to examine the distribution of the magnitude of location estimation error. The results are summarized in Figure 1. In this ?gure, each column represent results obtained from a particular method. Each row represents results from a particular sensor density. The histograms of the magnitudes of the localization error are plotted using a bin of 5 meters increment. Since the histogram can be regarded as an approximation of the probability density function, the mean and standard deviation of the magnitude of the localization error are calculated and listed in each ?gure. B. Performance Analysis: Two Sources Case We have also performed simulation to examine the performance of the proposed algorithm when there are two targets present in the sensor ?eld. Again, we assume a sensor ?eld of 100 by 100 meter squares. Three different source separations are used: 50 meters, 20 meters, and 10 meters. These pre-de?ned source locations are further perturbed randomly during simulation over a ±5 meter square area. The background noise parameters are ? = 1 and σ = 0.1. Both source intensity at 1 meter distance are de?ned as 5000 units. Three sensor network con?gurations are used: 6, 12, and 30 sensors. For each source location and
15
sources 50m apart 1 mean = 20.295 6 sensors stdev= 17.1205 0.5 0.5 1 mean = 18.235 stdev= 14.9661 0.5 sources 25m apart 1 mean = 16.075 stdev= 15.5369 sources 10m apart
0 1
5 15 25 35 45 55 65 75 85 95
0 1
5 15 25 35 45 55 65 75 85 95
0 1
5 15 25 35 45 55 65 75 85 95
mean = 13.305 12 sensors stdev= 14.0313 0.5 0.5
mean = 11.255 stdev= 9.3126 0.5
mean = 9.505 stdev= 8.2858
0 1
5 15 25 35 45 55 65 75 85 95
0 1
5 15 25 35 45 55 65 75 85 95
0 1
5 15 25 35 45 55 65 75 85 95
mean = 11.535 30 sensors stdev= 13.9335 0.5 0.5
mean = 9.265 stdev= 7.3932 0.5
mean = 7.84 stdev= 6.0112
0
5 15 25 35 45 55 65 75 85 95
0
5 15 25 35 45 55 65 75 85 95
0
5 15 25 35 45 55 65 75 85 95
Fig. 2.
Localization error of ?rst target of the two targets
sources 50m apart 1 mean = 21.005 6 sensors stdev= 16.8372 0.5 0.5 1 mean = 18.18 stdev= 14.7644 0.5 sources 25m apart 1 mean = 15.035 stdev= 14.4931 sources 10m apart
0 1
5 15 25 35 45 55 65 75 85 95
0 1
5 15 25 35 45 55 65 75 85 95
0 1
5 15 25 35 45 55 65 75 85 95
mean = 13.005 12 sensors stdev= 13.3892 0.5 0.5
mean = 11.19 stdev= 8.876 0.5
mean = 9.4 stdev= 8.3391
0 1
5 15 25 35 45 55 65 75 85 95
0 1
5 15 25 35 45 55 65 75 85 95
0 1
5 15 25 35 45 55 65 75 85 95
mean = 11.565 30 sensors stdev= 13.9015 0.5 0.5
mean = 9.1 stdev= 6.8257 0.5
mean = 8.005 stdev= 5.7114
0
5 15 25 35 45 55 65 75 85 95
0
5 15 25 35 45 55 65 75 85 95
0
5 15 25 35 45 55 65 75 85 95
Fig. 3.
Localization error of second target of the two targets
each sensor con?guration, 2000 simulations are performed. The mean and variance are computed. Since there are two targets, the localized results will return two positions. We assume the target association is done properly so that the location estimation error is the minimum of the two possible assignments. The histogram of the location estimation error of each source are listed in ?gures 2 and 3. C. Experiments:Application to Moving Vehicle Localization An application of the proposed EBSL methods to locate moving vehicle using distributed microphone sensor nodes will be reported in this section. First, we describe the overall system. Then the experiment results will be presented.
16
sensor deployment, road location and region specification for real experiment 100 sensor location road location
50 region 1 42 0 Y coordinates 1 region 2 58 52 53 ?50 49 ?100 47 48 50 51 54 56 55 59 41
Node 1 is the manager node for region 1 ?150 the detection node for region 1 are 41 42 51 54 55 56 59 Node 58 is the manager node for region 2 the detection node for region 2 are 47 48 49 50 52 53 ?200 ?200 ?150 ?100 ?50 0 X coordinates 50 100 150
Fig. 4.
sensor deployment, road coordinate and region speci?cation for experiments
1) Sensor Network System: In November 2001, a ?eld experiment sponsored by the DARPA ITO SensIT project has been carried out at 29 Palms California, USA. Custom-made prototype sensor nodes are laid out along side a road. Each sensor node is equipped with acoustic, seismic, and polorized infrared sensors, a 16-bit micro-processor, and radio transceiver and modem. It is powered by external car battery. During experiment, military vehicles such as AAV (Amphibious Assault Vehicle), DW (dragon wagon) were driven through the road, and sensors sampled the corresponding multi-modality data. The acoustic signal is sampled at 4960 Hz at 16-bit resolution. The set of data segments reported below are taken from the acoustic signatures of a single AAV travelling from east to west along the road during a time period of approximately 2 minutes. The sensor ?eld is partitioned into two regions with some details given in Figure 4. 2) Acoustic Energy-based Localization Experiment: The energy reading collected from all sensor nodes within the region at the same 0.75 second time interval were used for acoustic source localization at the manager node. Fig. 5 shows the AAV ground truth and the localization results based on the ML algorithm with MR projection solution and ER-NLS algorithm. The ground truth is obtained by interpolating an on-board GPS recordings that recorded a position ?x every 15 seconds. To use multi-resolution search, 3 grid sizes of 4, 2 and 1 The localization results are summarized in Figures 5 and 6. From ?gure 6, it is clear that the ML method out-perform the ER-NLS method. The data sampled during this experiment is very noisy. During the experiments, Gusty wind often blew directly into the microphone, creating isolated energy spikes in some of the sensors. Many microphones are not properly calibrated, and the gain factors gi estimated from the time series manually, and can be grossly inaccurate. The ground truth is not necessarily correct either due to the lack of differential GPS
17
experiment for aav3 single target localization 60 ground truth ebl ml road coordniates
40
miss detection, so localization is not active here
20
0 Y coordinate
?20
?40
?60
?80
?100 ?100
?50
0 X coordinate
50
100
150
Fig. 5.
AAV ground truth and localization estimation results based on ML algorithm with projection solution and NLS
algorithm (MR search is used, grid size is 4*4, 2*2 and 1*1. Estimation results look bias from the ground-truth, see discussion for reasoning)
histogram of localization err for ebl and ml estimation for aav3 experiment 35 ebl ml number of estimation points within the coresponding error range 30
25
20
15
10
5
0
10
20 30 40 estimation error range (m), 0~10,10~20,20~30,30~40,40~50
50
Fig. 6.
Estimation error histogram for AAV experiment data
settings and the lack of long term averaging at each position ?x. This is evident from ?gure 5 that the GPS marking is consistently off the road. Hence, the superior performance of the proposed ML source localization method in this experiment clearly demonstrate the feasibility of applying this method to handle real world data. VI. C ONCLUSION In this paper, a novel ML source localization method is presented. This method promises superior performance, multiple source localization, and is easy to compute. Besides algorithm derivation, the CRB of this algorithm has also been reported. Extensive simulations show that the proposed algorithm consistently out-performs other existing energy-based localization methods.
18
Future works include parameter sensitivity analysis, and sequential Bayesian estimation. R EFERENCES
[1] D. Estrin, Culler, D., Pister, K., and Sukhatme, G.,: Connecting the Physical World with Pervasive Networks, IEEE Pervasive Computing, vol. 1, Issue 1, 2002, pp. 59-69. [2] C. Savarese, Rabaey, J. M., and Reutel, J.: Localization in distributed Ad-hoc wireless sensor networks, Proc. ICASSP’2001, Salt Lake City, UT, (2001), 2037-2040 [3] D. Li, Wong, K.D., Hu, Y. H., Sayeed, A. M.,: Detection, classi?cation, and tracking of targets. IEEE Signal Processing Magazine, 19, (2002), 17–29 [4] Tolstory, A.,: Matched-Field Processing for underwater acoustics, Sigapore: World Scienti?c, (1993) [5] Brandstein, M., and Silverman, H.,: A localization-error-based method for microphone-array design, Proc. ICASSP’96, Atlanta, GA, (1996), 901-904 [6] Brandstein, M. S., Adcock, J. E., and Silverman,H. F.,: A closed form location estimator for use with room environment microphone arrays, IEEE Trans. Speech and Audio Processing, vol. 5, (1997), 45-50 [7] Huang, J., Ohnishi, N., and Sugie, N.,: Sound localization in reverberant enviornment based on model of the precedence effect, IEEE. Trans. Instrumentation and Measurement, 46, (1997), 842-846 [8] Omologo, M., and Svaizer, P.,: Acoustic source location in noisy and reverberant environment using CSP analysis, Proc. ICASSP’96, Atlanta, GA, (1996), 921-924 [9] Wang, H., and Chu, P., Voice source localization for automatic camera pointing system in videoconferencing, ICASSP-97., 1, (1997), 187 -190 [10] Suyama, K., Takahashi, K., and Hirabayashi, R.,: A robust technique for sound source localization in consideration of room capacity, presented at Proc. IEEE Workshop on Appl. Sig. Proc. to Audio and Acoust., New Paltz, NY, (2001) 63-66 [11] Weng, J., and Guentchev, K. Y.,: Three-dimensional sound localization form a compact non-coplanar array of microphones using three-based learning, J. Acoust. Soc. Am., 110, (2001) 311-323 [12] Griebel, S. M., and Brandstein, M. S.,: Microphone array source localization using realizable delay vectors, IEEE Workshop on Appl. Sig. Proc. to Audio and Acoust., New Paltz, NY, (2001). 71-74 [13] Yao, K., Hudson, R. E., Reed, C. W., Chen, D., and Lorenzelli, F.,: Blind beamforming on a randomly distributed sensor array system, IEEE J. Selected areas in communications, 16 (1998) 1555–1567 [14] Reed, C.W., Hudson, R., and Yao,K.,: Direct joint source localization and propagation speed estimation. In Proc. ICASSP’99, Phoenix, AZ, (1999) 1169-1172 [15] Haykin, S.,: Array Signal Processing, Englewood-Cliffs, NJ: Prentice-Hall, (1985) [16] Taff, L. G.,: Target localization from bearings-only observations, IEEE Trans. Aerosp. Electron., 3, issue 1, (1997) 2-10 [17] Oshman, Y., and Davidson, P., Optimization of observer trajectories for bearings-only target localization, IEEE Trans. Aerosp. Electron., 35, issue 3, (1999), 892-902 [18] Kaplan, K. M., Le, Q., and Molnar, P.,: Maximum likelihood methods for bearings-only target localization, Proc IEEE ICASSP, 5, (2001), 3001-3004 [19] Carter G. C.,: Coherence and Time Delay Estimation, IEEE Press, 1993. [20] Special issue on time-delay estimation, IEEE Trans. ASSP 29, (1981) [21] Smith, J.O., and Abel, J.S.,: Closed form least square source location estimation from range difference measurements. IEEE Trans. ASSP 35 (1987) 1661–1669
19
[22] D. Li, and Y. H. Hu,: Energy Based Collaborative Source Localization Using Acoustic Micro-Sensor Array, J. EUROSIP Applied Signal Processing, 2002 (to appear) [23] Kinsler, L.E., et al.,: Fundamentals of Acoustics. NY, NY: John, Wiley, and Sons, Inc., 1982 [24] K. Srodecki,: Evaluation of the reverberation decay quality in rooms using the autocorrelation function and the cepstrum analysis, Acustica, vol. 80, pp. 216-25, 1994. [25] Y. L. Li, M. J. White, and S. J. Franke,: New fast ?eld programs for anisotropic sound propagation through an atmosphere with a wind velocity pro?le, Journal of the Acoustical Society of America, vol. 95, pp. 718-26, 1994. [26] E. M. Salomons,: Downwind propagation of sound in an atmosphere with a realistic sound-speed pro?le: a semianalytical ray model, Journal of the Acoustical Society of America, vol. 95, pp. 2425-36, 1994. [27] T. Watanabe and S. Yamada,: Sound attenuation through absorption by vegetation, Journal of the Acoustical Society of Japan, vol. 17, pp. 175-82, 1996.
Xiaohong Sheng Xiaohong Sheng received B. Eng degree from Zhejiang University, China, and M. Eng degree from NanYang Technological University, Singapore in 1998. From 1998 to 1999, she worked for the Audio group at R&D, Creative Technology Ltd, Singapore. She received MSEE degree at the department of electrical and computer engineering, University of WisconsinMadison in 2001. She is currently a Ph.D. Candidate in the same department. Her research interests include Collaborative Signal Processing in wireless sensor network, Multi-Media Signal Processing and Statistical Signal Processing in Communication.
Yu Hen Hu (M’83-SM’88-F’99) Yu Hen Hu received BSEE degree from National Taiwan University, Taiwan, ROC in 1976. He received MS and PHD degree both in electrical engineering from University of Southern California, Los Angeles, CA in 1980 and 1982 respectively. Currently, he is a professor at the electrical and computer engineering department of the University of Wisconsin - Madison, WI, USA. Previously, he has been with the electrical engineering department of the Southern Methodist University, Dallas, TX. Dr. Hu’s research interests include multi-media signal processing, design methodology and implementation of signal processing algorithms and systems, and neural network signal processing. He published more than 200 journal and conference papers, edited two books Programmable Digital Signal Processors, and Handbook of Neural Network Signal Processing. Dr. Hu is a fellow of IEEE. He served as associate editors for IEEE Transactions on Signal Processing, IEEE Signal Processing letters, Journal of VLSI signal processing, European Journal of Applied Signal Processing. He served as secretary of IEEE signal processing society, board of governors of IEEE neural network council, chair of IEEE signal processing society, neural network signal processing technical committee.