- Email: [email protected]
Refining inaccurate sensor positions using target at unknown location
Signal Processing 92 (2012) 2097–2104
Contents lists available at SciVerse ScienceDirect
Signal Processing journal homepage: www.elsevier.com/locate/sigpro
Refining inaccurate sensor positions using target at unknown location M. Sun, K.C. Ho n Department of Electrical and Computer Engineering, University of Missouri, Columbia, MO 65211, USA
a r t i c l e i n f o
abstract
Article history: Received 17 August 2011 Received in revised form 19 January 2012 Accepted 26 January 2012 Available online 5 February 2012
This paper addresses the problem of improving the receiver positions in a sensor array using the positioning measurements of targets at unknown locations. Using an efficient estimate of the target locations, a computationally efficient estimator is proposed to refine the sensor positions using the same set of measurements from the targets. The considered measurements include TOA and TDOA. The proposed estimator has closedform and is able to reach the CRLB accuracy under small Gaussian noise which is supported by theoretical analysis and simulation studies. & 2012 Elsevier B.V. All rights reserved.
Keywords: Localization Sensor position errors Time difference of arrival (TDOA) Time of arrival (TOA)
1. Introduction Many engineering applications in practice requires the localization of a signal source, such as in radar and sonar [1,2], wireless communications [3] and more recently sensor networks [4,5]. The most common form of localization is to use a number of sensors to acquire the emitted source signal, compute the positioning measurements such as time of arrival (TOA) or time difference of arrival (TDOA) of the received signals at different sensors and estimate the source location from the positioning measurements. Traditionally, the sensor positions are assumed to be exactly known. Modern localization system often deploys the sensors randomly or uses dynamic sensors. In some applications such as sensor networks, quite often the sensors are without GPS to save costs and the sensor positions are not known precisely. Inaccurate sensor positions lead to degradation in localization performance if the sensor position errors are ignored [6,7]. A number of algorithms became
n
Corresponding author. Tel.: þ573 882 8023; fax: þ 573 882 0397. E-mail addresses: [email protected] (M. Sun), [email protected] (K.C. Ho). 0165-1684/$ - see front matter & 2012 Elsevier B.V. All rights reserved. doi:10.1016/j.sigpro.2012.01.022
available for localization problems with inexact sensor positions in literatures [7–11]. Methods [7,11] have been proposed to improve the source location estimate by taking into account the statistical behavior of the sensor position errors. However, relatively little emphasis has been allocated to improve the sensor positions so that they will be able to identify better the location of another source later. In this work, we are interested in improving the sensor positions given a set of positioning measurements from a single or multiple sources at unknown locations. A direct approach would be to jointly estimate the source (target) and sensor positions using the maximum likelihood (ML) method [12–14]. The non-linear nature of the problem makes it difficult to solve and it requires iterations with good initial guesses to obtain the ML estimate. Instead of solving the ML function directly, Yang and Ho [8] recently proposed an algebraic closed-form solution that jointly estimates the target and sensor positions. The target positions are able to achieve the Cramer Rao low bound (CRLB) accuracy, however, the sensor positions cannot. Another approach is to estimate the target locations first and then use the target location estimate to refine the sensor positions. This approach enables the target positions to be identified quickly for immediate actions with
2098
M. Sun, K.C. Ho / Signal Processing 92 (2012) 2097–2104
much less computations. Indeed, the solution in [7] estimates the target position only, which, through the design of a proper weighting matrix to account for the sensor position error, can reach the CRLB performance. After locating the target, the refinement of sensor positions can be followed later off-line when the sensing system is less congested to ease the computation burden. While the sensor position improvement stage in [15] by considering the unknown target as a calibration source can be applied to the current problem, it is not able to achieve the CRLB accuracy. In addition, the method in [15] assumes the position error of the target and the measurement noise are independent, which is clearly not the case for the current problem at hand because the target position is derived from the same set of measurements. This paper proposes a solution to achieve the CRLB performance in refining sensor positions. It takes the estimated target locations and the same set of measurements to improve the sensor positions. The proposed estimator is closed-form and does not require iterations or initialization. Perhaps the most important aspect of the estimator is that it reaches the CRLB accuracy for small Gaussian noise before thresholding effect occurs. This is proven analytically under small noise assumption and the performance is supported by simulations. We derive the proposed estimator for commonly used observations in localization that include TOA and TDOA. The measurement errors considered are Gaussian for ease of CRLB derivation and performance comparison, although the proposed estimator can be applied without Gaussian noise assumption. Nevertheless, the Gaussian measurement model has been validated under realistic environments in [16] and [17] and has been extensively used in literatures [7–9,11–17]. The development focuses on the scenario of a single target. The extension of the estimator to multiple target situation is also provided in which the performance can be shown to reach the CRLB accuracy as well. The proposed algorithm will be useful in applications such as sensor networks, where the sensors do not have GPS to identify their accurate positions, the computational resources are limited in locating targets and the refinement in sensor positions results in performance improvement of subsequent targets. In the following, Section 2 formulates the sensor position refinement problem mathematically. Section 3 develops the new estimator in closed-form. Section 4 conducts performance analysis with respect to the CRLB accuracy. The proposed solution is extended to multiple target scenario in Section 5. Section 6 presents the simulation results and Section 7 concludes the paper. 2. Problem formulation We shall consider 2D localization in this study for ease of illustration. There are M stationary sensors whose known but erroneous positions are si ¼ soi þ Dsi ,i ¼ 1; 2, . . . ,M, where soi ¼ ½xoi ,yoi T is the true position and Dsi is the error component. The collection of all the sensor positions forms s ¼ ½sT1 ,sT2 , . . . ,sTM T ¼ so þ Ds, where Ds ¼ ½DsT1 , DsT2 , . . . , DsTM T . The sensor position error Ds is assumed to be a zero-mean Gaussian random vector with covariance matrix Qs. Besides the sensors, the target appears at the position
uo ¼ ½xo ,yo T that is not known. The observations to locate the target are TOAs or TDOAs and they are represented by the measurement vector m ¼ mo þ nm
ð1Þ
o
where m is the noise free value. nm is the measurement noise and it is modeled as a zero-mean Gaussian random vector with covariance matrix Qm. The cross-covariance between nm and Ds is Qms. Qms in most cases has all zero elements. We include it here because the proposed method works even if correlation exists between the measurement noise and sensor position errors. For the case of TOA, mo ¼ ½r o1 ,r o2 , . . . ,r oM T and r oi is related to the unknown target position uo by r oi ¼ Juo soi J,
i ¼ 1; 2, . . . ,M
ð2Þ o
and JnJ represents the 2-norm. For TDOA, m . . . ,r oM1 T and r oi1 ¼ r oi r o1 ,
i ¼ 2; 3, . . . ,M
¼ ½r o21 ,r o31 , ð3Þ
where sensor 1 denotes the reference sensor in obtaining the TDOAs. Note that we have multiplied the TOA and TDOA with the known signal propagation speed when forming m. We shall suppose an estimate of the target position uo has been obtained using m that reaches CRLB accuracy. Several methods such as those in [7,11] are able to achieve this purpose. They do not need joint estimation of the sensor and target positions and can produce the target position estimate very fast. The requirement to achieve the CRLB is that the noise in the positioning measurements is small. Based on the target position estimate, the objective is to improve the sensor position s as much as possible using the same measurement vector m. 3. Proposed solution Instead of finding the sensor position so directly, the proposed solution estimates the error component Ds to correct s. The sensor position error Ds is a random vector. Hence we are estimating a particular realization of Ds. The set of pseudo-linear equations is composed of two parts: one is related to the positioning measurements and the other is from the probability distribution of the noisy sensor positions. By exploring the efficient estimate of the target position, a closed-form solution can be directly formulated. 3.1. Pseudo-linear equations 3.1.1. TOA Squaring both sides of TOA measurement r i ¼ r oi þ Dr i yields 2r oi Dr i C r 2i Juo soi J2
ð4Þ
where (2) has been used and ðDr i Þ2 has been ignored. The true sensor position soi is not available and only their noisy values are known. After substituting soi ¼ si Dsi , we can rewrite (4) as r oi Dr i C 12ðr 2i Juo si J2 Þðuo si ÞT Dsi
ð5Þ
M. Sun, K.C. Ho / Signal Processing 92 (2012) 2097–2104
when the second order error term DsTi Dsi is ignored. The matrix form of (5) is BT nm CqoT PoT Ds
ð6Þ
where BT ¼ diagðr o1 ,r o2 , . . . ,r oM Þ
ð7Þ
and diagða, b, . . . , gÞ denotes a block diagonal matrix with diagonal blocks a, b, . . . , g. On the right hand side, the M 1 vector qoT is qoT ¼ 12½r 21 Juo s1 J2 ,r 22 Juo s2 J2 , . . . ,r 2M Juo sM J2 T
ð8Þ
and PoT is an M 2M block-diagonal matrix defined as PoT ¼ diagððuo s1 ÞT ,ðuo s2 ÞT , . . . ,ðuo sM ÞT Þ:
3.1.2. TDOA Expressing the TDOA measurement ri1 as r i1 ¼ r oi1 þ Dr i1 and applying (3) yield r i1 þ r o1 ¼ r oi þ Dr i1 . Squaring both sides and using (2), we have o oT o 2r oi Dr i1 Cr 2i1 þ 2r i1 r o1 þ soT 1 s1 si si
i ¼ 2; 3, . . . ,M
ð10Þ
where ðDr i1 Þ2 has been neglected. We shall express soi as soi ¼ si Dsi . Different from the TOA case, we have the extra unknown variable r o1 that is dependent on uo and so1 . Using the Taylor series expansion around s1 and retaining up to linear error term, we have r o1 ¼ Juo so1 J C r~ o1 þ qTuo ,s1 Ds1 r~ o1
where the subscript m, which can be T or D, represents the measurement type. We shall use the common form (18) in subsequent development. Eq. (18) is constructed solely from the measurements. The statistical knowledge of the sensor position error Ds is exploited next. We follow the approach in [18] where the prior is utilized and by stacking (18) with ns ¼ Ds. Ds on the right-hand side is treated as the unknown parameter to be estimated and on the left-hand side ns is interpreted as a zero-mean random vector of position error with covariance matrix Qs. The resulting composite pseudo-linear equation is Bn C qo Po Ds
ð9Þ
Eq. (6) is interpreted as a pseudo-linear equation with respect to the unknown Ds since the true target position uo that appears in qoT and PoT is also not known. The lefthand side corresponds to the equation error due to the measurement noise nm.
þ 2ðsoi so1 ÞT uo ,
2099
ð11Þ
o
ð19Þ T T n ¼ ½nTm ,nTs T , qo ¼ ½qoT m , 02M1 o
where B ¼ diagðBm ,I2M2M Þ, T T and Po ¼ ½PoT and Ds m ,I2M2M . The unknowns are u (2Mþ2 variables) while Bn is the equation error. 3.2. Computing the sensor position estimate Defining the cost function to be minimized as J ¼ ðqo Po DsÞT Rðqo Po DsÞ
ð20Þ
where R is the weighting matrix obtained by taking the inverse of the equation error covariance matrix [19], i.e. R ¼ E½ðBnÞðBnÞT 1 ¼ ðBQBT Þ1 , " # Q m Q ms T Q ¼ E½nn ¼ Q Tms Q s
ð21Þ
R is always invertible as long as the target is not colocated with a sensor and the measurement noise covariance matrix Q is full rank, which is expected to be the case in practice. We shall denote the solution that minimizes (20) as u and Ds^ . Obviously, they satisfy @J ¼0 ð22Þ @Ds
¼ Ju s1 J and qa,b ¼ ðabÞ=JabJ represents where the unit vector pointing from b to a. As a result, (10) reduces to
Evaluating the above gradient and simplifying, we have
r oi Dr i1 Cqoi þ ðr i1 qTuo ,s1 þ ðuo s1 ÞT ÞDs1 ðuo si ÞT Dsi
Ds^ ¼ ðPT RPÞ1 PT Rq
ð12Þ
ð23Þ o
where qoi ¼
u, Ds^
T o T T 1 2 ~o 2ðr i1 þ 2r i1 r 1 þs1 s1 si si þ 2ðsi s1 Þ u Þ
ð13Þ
Stacking (12) for i¼2,3,y,M together gives the matrix equation BD nm C qoD PoD Ds
ð14Þ
where BD ¼ diagðr o2 ,r o3 , . . . ,r oM Þ
ð15Þ
qoD ¼ ½qo2 ,qo3 , . . . ,qoM T
ð16Þ
and the (i 1)th row of the (M 1) 2M matrix o
½ðr i1 quo ,s1 þ ðu s1 ÞÞ i ¼ 2; 3, . . . ,M
T
0T2ði2Þ1
o
ðu si Þ
T
PoD
is
o
where P and q are P and q with u replaced by u. Eq. (23) indicates that a closed-form solution of Ds can be directly obtained when u is available, which can come from another estimator as long as it gives the efficient target position estimate. u can be represented as u ¼ uo þ Du and its covariance is covðuÞ ¼ CRLBðuo Þ, where Du is the estimation error. Finding Ds using (23) requires the weighting matrix R. The matrix B in the weighting matrix R depends on the true sensor–source distances and they are replaced by the computed distances using the noisy sensor position s and the source estimate u. After Ds^ is obtained from (23), the estimate of so is s^ ¼ sDs^
ð24Þ
0T2ðMiÞ1 , ð17Þ
Regardless of the measurement type as TOA or TDOA, the solution equation for Ds has the generic form Bm nm Cqom Pom Ds
o
ð18Þ
4. Performance analysis In this Section, the performance of the proposed algorithm is evaluated and compared with the CRLB accuracy. The CRLB of the sensor positions is derived in Appendix A
2100
M. Sun, K.C. Ho / Signal Processing 92 (2012) 2097–2104
and it is given in (A.9). We first express qo in (19), by substituting uo ¼ uDu, as
of the proposed estimator is approximately efficient under C1 and C2
qo CqLDu
covðs^ Þ C CRLBðso Þ
ð25Þ
where the second order error terms have been ignored and L ¼ ½LTm ,0T2M2 T where 0 is used to denote a matrix of all zeros. For the case of TOA, using the definition of qoT in (8) gives Lm ¼ LT ¼ ½ðus1 Þ ðus2 Þ ðusM ÞT
ð26Þ
For the case of TDOA, applying the approximation that r~ o1 ¼ Juo s1 J C r~ 1 qTu,s1 Du where r~ 1 ¼ Jus1 J, we have from the definition of the elements of qoD in (16) Lm ¼ LD ¼ ½r 21 qu,s1 þðs2 s1 Þ r 31 qu,s1 þðs3 s1 Þ r M1 qu,s1 þ ðsM s1 ÞT
ð27Þ
Eq. (19) can now be rewritten, by substituting (25) and Po ¼ PDP, as q CPDsþ Bn þ LDu
ð28Þ
where the second order error term DPDs has been ignored. After substituting (28) in (23), the estimation error of Ds^ can be expressed as e9Ds^ Ds ¼ Z^
1
PT RBn þ Z^
1
T Y^ Du
ð29Þ
Y^ ¼ LT RP
ð30Þ
We now evaluate the covariance matrix of the proposed estimator. Eq. (24) can be expressed as s^ ¼ so þ DsDs^ . Hence the covariance matrix is covðs^ Þ ¼ covðDs^ Þ ¼ E½eeT
ð31Þ
From (29) and noting that R ¼ ðBQBT Þ1 , we have 1 1 T 1 1 T ^ Z^ 1 covðs^ Þ ¼ Z^ þ Z^ Y^ covðuÞY^ Z^ þ Z^ ðY^ T þ TT YÞ
ð32Þ o
where covðuÞ ¼ CRLBðu Þ and T9E½Du nT BRP
ð33Þ
We would like to prove that (32) and (A.9) are nearly identical under the following two small noise conditions C1: 9nm,i 9=r oi C 0; C2: JDsi J=r oi C 0, for i ¼ 1; 2, . . . ,M
The sensor position estimator and its performance analysis presented so far only use measurements from a single target. We shall extend the developed algorithm to the multiple target scenario and validate that the obtained sensor position estimate remains to be nearly efficient under C1 and C2. Suppose there are K targets whose true positions oT oT T uo ¼ ½uoT 1 ,u2 , . . . ,uK are not known. The measurement vector becomes m ¼ ½mT1 ,mT2 , . . . ,mTK T ¼ mo þ nm , where nm ¼ ½nTm,1 ,nTm,2 , . . . ,nTm,K T . Given the efficient estimate of the target position vector u ¼ uo þ Du obtained from the measurement vector m, we are interested to refine the sensor positions. Recalling (18) is the equation for the measurements from a single target, we rewrite it here for target i as Bm,i nm,i C qom,i Pom,i Ds
ð36Þ
B1 P C P~
Bm nm Cqom Pom Ds
ð34Þ
As a result of (34) and using (21), we have Z^ C Z and Y^ C Y from (30) and (A.4). Also, we have shown in Appendix C that if (34) is valid, T C 0. Thus, using an estimate u with CRLB accuracy for the proposed algorithm, its covariance matrix (32) is nearly identical to (A.9) and the performance
ð37Þ
where Bm ¼ diagðBm,1 ,Bm,2 , . . . ,Bm,K Þ. On the right-hand oT oT T and Pom ¼ ½PoT side of (37), qom ¼ ½qoT m,1 , m,1 ,qm,2 , . . . ,qm,K oT T PoT m,2 , . . . ,Pm,K . Combining (37) and the statistical knowledge of the sensor position error yields the pseudo-linear equation having the same form as (19). Once the pseudo-linear equation (19) is formed, the procedure to compute the sensor position estimate is straightforward and is exactly the same as the one developed in Section 3.B. The performance analysis follows the same as in Section 4. The only difference is that L in (25) is given by L ¼ ½LTm ,0T2M2K T and Lm ¼ diagðLm,1 ,Lm,2 , . . . ,Lm,K Þ. Lm,i is defined as in (26) for TOA or (27) for TDOA with u replaced by ui . From Appendix C, under C1 and C2, we have for i¼1, 2,y,K, o o B1 m,i L m,i C ð@mi =@ui Þ,
where nm,i is the ith element of nm in (1). C1 states the measurement noise is small compared to the distance between the target and the sensor. C2 means that the 2norm of the ith sensor position error is much smaller than the distance between the target and the sensor. Both conditions are satisfied if the noise is sufficiently small or the target is far away from the sensor array. Appendix B shows that under C1 and C2, ~ B1 L CL,
5. Multiple target scenario
Stacking (36) for i¼1,2,y,K gives the equation of all measurements,
where Z^ ¼ PT RP,
ð35Þ
o o B1 m,i Pm,i C ð@mi =@s Þ
ð38Þ
and (34) is valid. As a result, (35) holds true for the case of multiple targets as well. 6. Simulations The localization system in the simulation study is composed of M ¼4 sensors and their true positions are so1 ¼ ½8,16T ,so2 ¼ ½28,40T ,so3 ¼ ½20,8T and so4 ¼ ½0; 8T as shown in Fig. 1. The measurement noise covariance matrix is Q m ¼ 0:001IMM for TOA. For TDOA, the covariance matrix is Qm ¼0.001J where J is an (M 1) (M 1) matrix with diagonal elements equal to 1 and all other elements 12. This form of covariance assumes ML estimation of TDOAs whose accuracy is independent of the choice of reference sensor [20]. The measurement noise and sensor
M. Sun, K.C. Ho / Signal Processing 92 (2012) 2097–2104
20
2101
25
Targets Sensors
Proposed method Before refinement CRLB
20 10log (Sensor position MSE)
10
y
0
−10
−20
−30
−40 −20
15
10
5
0
−5 −10
0
10
20
30
40
50
−5
x
position errors are independent so that Qms ¼0. The sensor positions are not known at the beginning and they are obtained through three emitting sources at known positions ao1 ¼ ½19:4,2:8T , ao2 ¼ ½1:8,4:8T and ao3 ¼ ½5; 17T from TOA measurements. We also assume the sensors are active and can communicate with each other so that TOAs between so1 and so4 and between so2 and so3 are also available to identify their positions. For illustration purpose, the TOA noise values to obtain the sensor positions are independent and identically distributed zero-mean Gaussian random variables with noise power s2r . The ML estimator is applied to obtain the sensor position vector s using these TOAs. The covariance matrix Qs is set to the associated CRLB, where the true sensor positions in the bound are replaced by the estimated values. The initial sensor position vector therefore has error and the amount of sensor position errors is proportional to s2r . At a given signal-to-noise ratio (SNR), the corresponding noise power is obtained using (2) of [4] with the effective signal bandwidth set to 20 MHz. The localization process proceeds as follows: (i) obtain the sensor position vector s using the ML method based on the TOAs from ao1 ,ao2 ,ao3 and among the sensors; (ii) localize the targets using the TOA/TDOA measurements between the targets and the sensors, where the method in [11] is used as an example to produce the efficient target position estimate u; (iii) apply the proposed algorithm using the estimated target locations and the same set of measurements to improve the sensor position vector s. P l^ o 2 The performance measure is MSEðs^ Þ ¼ N l ¼ 1 J s s J =N, o l^ where s is the estimate of s at ensemble l and N ¼2000 is the number of ensemble runs. We first examine the performance of the proposed method using a single target located at unknown position uo1 ¼ ½4,24T . The results are shown in Fig. 2 for TOA and Fig. 3 for TDOA. Each figure plots the MSEðs^ Þ and the trace of CRLBðso Þ as a function of SNR, before and after applying the proposed sensor refinement method. The proposed estimator can improve the sensor locations by about 3 dB for TOA and 1 dB TDOA. Furthermore, it is able to achieve the CRLB accuracy before the thresholding effect occurs. In Fig. 3, the sensor position estimate from [8] is also plotted for comparison. The MSE of the sensor position
5
10
15
20
SNR (dB)
Fig. 2. Performance of the proposed method in improving the sensor positions using TOA measurements from an unknown target whose true location is uo1 . 20 Proposed method Solution in [8] Before refinement CRLB
15 10log (Sensor position MSE)
Fig. 1. Sensor network configuration.
0
10
5
0
−5
0
5
10
15
20
SNR (dB)
Fig. 3. Performance of the proposed method in improving the sensor positions using TDOA measurements from an unknown target whose true location is uo1 .
estimate almost overlaps with the variance of the sensor positions before estimation. This means that the solution in [8] is not able to improve the sensor positions much as expected. It is observed though that the proposed method has lower noise threshold than [8] before the performance deviates from the bound. The computational complexity of the two methods for the case of TDOA is also examined. Both the solution in [8] and the proposed estimator involve the inverse of a matrix that accounts for most of the computation burden. The complexity of the method in [8] is dominated by its first stage WLS computation that jointly estimates the source and sensor positions. The size of the matrix that needs to be inverted for 2D localization is 3Kþ2M, where 2K is for the K source coordinates, another K is for K auxiliary variables one for each source and 2M comes from the estimation of M sensor positions. However, the proposed algorithm only requires to estimate the 2M sensor positions and the inverse matrix has a size of 2M only. The positions of the K sources can be
M. Sun, K.C. Ho / Signal Processing 92 (2012) 2097–2104
estimated with a matrix inverse of size 3K [11]. Therefore the computation complexity is Oðð3K þ2MÞ3 Þ for the method in [8] and Oðð3KÞ3 þ ð2MÞ3 Þ for the proposed estimator. The proposed estimator has lower complexity than [8]. The reduction in sensor position errors can be represented by the concentration ellipses. Fig. 4 gives the concentration ellipses of the sensor positions before (dashed line) and after (solid line) improvement using TOAs measurements where SNR is set to be 4.55 dB. Figs. 5 and 6 present the results using another target at uo2 ¼ ½49; 19T . Note that different from the target at uo1 , the target at uo2 is outside of the area enclosed by the sensors. Compared to the case that target is at uo1 , the improvement in sensor positions is less. Nevertheless, the proposed method remains to reach the CRLB accuracy. The performance of the proposed method using multiple targets is also evaluated. Figs. 7 and 8 present the results with two targets at uo1 ¼ ½4,24T and uo2 ¼ ½49; 19T . The measurements are TOA and they are independent from the two targets. Compared to the single target case, the
20
10
5
0
−5
0
5
20
25
y
−10
−20
−30
Proposed method Before refinement CRLB
20 10log (Sensor position MSE)
0
−40
15
10
5
0
−5
−10
0
10
20
30
−5
40
0
5
10
15
20
SNR (dB)
x
Fig. 4. Sensor position improvement illustration using uncertainty ellipse.
Fig. 7. Performance of the proposed method in improving the sensor positions using TOA measurements from two unknown targets whose true locations are uo1 and uo2 .
25
20
Proposed method Before refinement CRLB
20
Proposed method Solution in [8] Before refinement CRLB
15 10log (Sensor position MSE)
10log (Sensor position MSE)
15
Fig. 6. Performance of the proposed method in improving the sensor positions using TDOA measurements from an unknown target whose true location is uo2 .
Target Sensors Before refinement After refinement
10
15
10
5
10
5
0
0
−5
10 SNR (dB)
20
−50 −20
Proposed method Solution in [8] Before refinement CRLB
15 10log (Sensor position MSE)
2102
−5
0
5
10
15
20
SNR (dB)
Fig. 5. Performance of the proposed method in improving the sensor positions using TOA measurements from an unknown target whose true location is uo2 .
−5
0
5
10
15
20
SNR (dB)
Fig. 8. Performance of the proposed method in improving the sensor positions using TDOA measurements from two unknown targets whose true locations are uo1 and uo2 .
M. Sun, K.C. Ho / Signal Processing 92 (2012) 2097–2104
25 Target position MSE, without refining sensor positions Target position CRLB, without refining sensor positions Target position MSE, with refining sensor positions Target position CRLB, with refining sensor positions
10log (Target position MSE)
20
15
2103
where c is a constant not depending on ho . We can obtain from (A.1) @ ln pðx; ho Þ o ~ T Rðxx ~ ¼ ½L~ P Þ @ho
ðA:2Þ
where 2 3 @mo o ~L9 @x ¼ 4 @uo 5, @uo 0
10
2 3 @mo @xo 4 ~ o P9 o ¼ @s 5, @s I
2M2
5
R~ ¼ Q 1
2M2M
ðA:3Þ 0
−5 −5
0
5
10
15
20
Hence, the Fisher information matrix (FIM) is, after postmultiplying (A.2) by its transpose and taking expectation, # 2 " T T X Y @ ln pðx; ho Þ L~ R~ L~ L~ R~ P~ ¼ FIMðho Þ9E ¼ T T YT Z @ho @hoT P~ R~ L~ P~ R~ P~
SNR (dB)
ðA:4Þ
Fig. 9. Localization accuracy of a new subsequent target with and without improving the sensor positions from a previously observed target at uo1 , TOA measurements.
improvement in sensor positions is larger as expected. The proposed method remains to reach the CRLB accuracy. The main purpose of improving sensor positions is to increase the localization accuracy in the next round of target positioning. To illustrate this, Fig. 9 shows the localization result of a new target at ½36; 39T with (crosses) and without (squares) improving the sensor positions from the single target at uo1 ¼ ½4,24T that is present earlier. The measurement is TOA and the method in [11] is used to locate the new target. It is clear that there is about 5 dB performance improvement when using the refined sensor positions. We would like to clarify that the proposed algorithm is intended for use when the SNR is high. When the SNR is low, simulation shows its performance deviates from the CRLB.
The partial derivatives are given below. For TOA, the ith row of ð@mo =@uo Þ is ð@r oi =@uo ÞT ¼ qTuo ,so , i
i ¼ 1; 2, . . . ,M
ðA:5Þ
and ð@mo =@so Þ is a block diagonal matrix given by ð@mo =@so Þ ¼ diagðqTuo ,so , qTuo ,so , . . . , qTuo ,so Þ 1
ðA:6Þ
M
2
For TDOA, the (i 1)th rows of ð@mo =@uo Þ and ð@mo =@so Þ are ð@r oi1 =@uo ÞT ¼ ½quo ,so quo ,so T , 1
i
i ¼ 2; 3, . . . ,M
ðA:7Þ
ð@r oi1 =@so ÞT ¼ ½qTuo ,so 0T2ði2Þ1 qTuo ,so 0T2ðMiÞ1 , 1
i
i ¼ 2; 3, . . . ,M
ðA:8Þ
Invoking the block matrix inversion formula in (A.4), the CRLBðso Þ is CRLBðso Þ ¼ ðZYT X1 YÞ1 ¼ Z1 þZ1 YT CRLBðuo ÞYZ1 ðA:9Þ
7. Conclusions
where CRLBðu Þ ¼ ðXYZ
In this paper, an explicit solution was proposed to refine the sensor positions using TOA or TDOA measurements from a single or multiple targets at unknown locations. Suppose that an efficient estimate of target locations is obtained by some techniques, the proposed method refines the sensor positions using the same set of measurements from the unknown targets. This approach enables the target positions to be identified quickly first with much less computations while leaving the improvement of sensor positions at a later time when computation resources are available. Performance analysis and simulations showed that the proposed method is able to achieve the CRLB accuracy under small Gaussian noise.
T 1
Y Þ
.
Appendix B. Proof of (34) In (34), the definitions of B, P and L are under (19) and (25). We begin with the proof of (34) for TOA. Bm and Pm are BT, LT and PT defined below (6) and LT is given in (26). After substituting their expressions and simplifying yield us us2 usM T 1 B1 L ¼ ðB:1Þ o o T T r r ro 1
2
M
T T T o o o B1 T PT ¼ diagððus1 Þ =r 1 ,ðus2 Þ =r 2 , . . . ,ðusM Þ =r M Þ
ðB:2Þ From C1 and C2 and based on the form of Du in (C.3), JDuJ=r oi C 0. Hence under both conditions,
Appendix A. Derivation of the CRLB From the noise model, the logarithm of the density function of the data vector x ¼ ½mT ,sT T parameterized on the unknown vector ho ¼ ½uoT ,soT T is ln pðx; ho Þ ¼ c12ðxxo ÞT Q 1 ðxxo Þ
1
o
ðA:1Þ
ðusi Þ=r oi ¼ ðuo soi þ DuDsi Þ=r oi C quo ,so , i
i ¼ 1; 2, . . . ,M ðB:3Þ
Putting (B.3) back into (B.1) and (B.2) and comparing with (A.5) and (A.6) validate (34).
2104
M. Sun, K.C. Ho / Signal Processing 92 (2012) 2097–2104
For TDOA, after using BD and PD given below (14) and LD in (27), the (i 1)th row of B1 D L D is
References
ðr i1 qTu,s1 þ ðsi s1 ÞT Þ=r oi ,
[1] F. Ahmad, M. Amin, Noncoherent approach to through-the-wall radar localization, IEEE Transactions on Aerospace and Electronic Systems 42 (2006) 1405–1419. [2] S. Coraluppi, Multistatic sonar localization, IEEE Journal of Oceanic Engineering 31 (2006) 964–974. [3] G. Sun, J. Chen, W. Guo, K.J.R. Liu, Signal processing techniques in network-aided positioning: a survey of state-of-the-art positioning designs, IEEE Signal Processing Magazine 22 (2005) 12–23. [4] S. Gezici, T. Zhi, G. Giannakis, H. Kobayashi, A. Molisch, H. Poor, Z. Sahinoglu, Localization via ultra-wideband radios: a look at positioning aspects for future sensor networks, IEEE Signal Processing Magazine 22 (2005) 70–84. [5] T. Li, A. Ekpenyong, Y.-F. Huang, Source localization and tracking using distributed asynchronous sensors, IEEE Transactions on Signal Processing 54 (2006) 3991–4003. [6] A. Ferreol, P. Larzabal, M. Viberg, On the asymptotic performance analysis of subspace DOA estimation in the presence of modeling errors: case of MUSIC, IEEE Transactions on Signal Processing 54 (2006) 907–920. [7] K.C. Ho, X. Lu, L. Kovavisaruch, Source localization using TDOA and FDOA measurements in the presence of receiver location errors: analysis and solution, IEEE Transactions on Signal Processing 55 (2007) 684–696. [8] L. Yang, K.C. Ho, An approximately efficient TDOA localization algorithm in closed-form for locating multiple disjoint sources with erroneous sensor positions, IEEE Transactions on Signal Processing 57 (2009) 4598–4615. [9] K.W.K. Lui, W.-K. Ma, H.C. So, F.K.W. Chan, Semidefinite programming approach to sensor network node localization with anchor position uncertainty, in: Proceedings of the IEEE International Conference on Acoustics, Speech, and Signal Processing (ICASSP’09), Taipei, Taiwan, 2009, pp. 2245–2248. [10] W. Chiu, B. Chen, C. Yang, Robust relative location estimation in wireless sensor networks with inexact position problems, IEEE Transactions on Mobile Computing 99 (2011) 1. [11] M. Sun, K.C. Ho, An asymptotically efficient estimator for TDOA and FDOA positioning of multiple disjoint sources in the presence of sensor location uncertainties, IEEE Transactions on Signal Processing 59 (2011) 3434–3440. [12] Y. Rockah, P.M. Schultheiss, Array shape calibration using sources in unknown location. Part I: far-field source, IEEE Transactions on Acoustics, Speech, and Signal Processing ASSP-35 (1987) 286–299. [13] Y. Rockah, P.M. Schultheiss, Array shape calibration using sources in unknown location. Part II: near-field source and estimator implementation, IEEE Transactions on Acoustics, Speech, and Signal Processing ASSP-35 (1987) 724–735. [14] A. Weiss, B. Friedlander, Array shape calibration using sources in unknown locations: a maximum likelihood approach, IEEE Transactions on Acoustics, Speech, and Signal Processing ASSP-37 (1989) 1958–1966. [15] L. Yang, K.C. Ho, Alleviating sensor position error in source localization using calibration emitters at inaccurate locations, IEEE Transactions on Signal Processing 58 (2010) 67–83. [16] N. Patwari, J.N. Ash, S. Kyperountas, A.O. Hero III, R.L. Moses, N.S. Correal, Locating the nodes: cooperative localization in wireless sensor networks, IEEE Signal Processing Magazine 22 (2005) 54–69. [17] B. Alavi, K. Pahlavan, Modeling of the TOA-based distance measurement error using UWB indoor radio measurements, IEEE Communications Letters 10 (2006) 275–277. [18] H.W. Sorenson, Parameter Estimation: Principles and Problems, Marcel Dekker, New York, 1980. [19] S.M. Kay, Fundamentals of Statistical Signal Processing, Estimation Theory, Prentice Hall, Englewood Cliffs, NJ, 1993. [20] K.C. Ho, W. Xu, An accurate algebraic solution for moving source location using TDOA and FDOA measurements, IEEE Transactions on Signal Processing 52 (2004) 2453–2463.
i ¼ 2; 3, . . . ,M
ðB:4Þ
and the (i 1)th row of B1 D PD is ½ðr i1 qTu,s1 þ ðus1 ÞT Þ=r oi 0T2ði2Þ1 ðusi ÞT =r oi 0T2ðMiÞ1 , i ¼ 2; 3, . . . ,M
ðB:5Þ
It can be easily verified that under C1 and C2, r~ 1 C r o1 and r i1 C r oi1 . Using these approximations and (B.3) yield ðr i1 qTu,s1 þ ðsi s1 ÞT Þ=r oi C ðr oi1 ðus1 ÞT þ r o1 ðsi s1 ÞT Þ=ðr oi r o1 Þ C qTuo ,so qTuo ,so 1
ðB:6aÞ
i
ðr i1 qTu,s1 þ ðus1 ÞT Þ=r oi C ðr oi1 ðus1 ÞT þ r o1 ðus1 ÞT Þ=ðr oi r o1 Þ CqTuo ,so
1
ðB:6bÞ
Applying (B.6) to (B.4) and (B.5) and (B.5) and comparing with (A.7) and (A.8), we have the results in (34). Appendix C. Value of T under conditions C1 and C2 We shall express Du in terms of n ¼ ½nTm ,nTs T . According to the CRLB theorem [19], an unbiased estimator may be found that attains the bound in that covðhÞ ¼ FIMðho Þ1 if and only if ð@ ln pðx; ho Þ=@ho Þ ¼ FIMðho ÞðgðxÞho Þ and h ¼ gðxÞ is the efficient estimate of ho . Thus, for an estimator that yields the CRLB accuracy, the estimation error must be of the form
Dh ¼ gðxÞho ¼ FIMðho Þ1 ð@ ln pðx; ho Þ=@ho Þ
ðC:1Þ
Note that in our case, u is assumed to be an efficient estimate of uo obtained from existing literatures, such as [8,11]. Therefore, it corresponds to the first two elements of gðxÞ with covðuÞ ¼ CRLBðuo Þ. FIMðho Þ1 in (C.1) can be expressed as, when applying the block matrix inversion formula to (A.4), " # CRLBðuo ÞYZ1 CRLBðuo Þ FIMðho Þ1 ¼ ðC:2Þ 1 T Z Y CRLBðuo Þ CRLBðso Þ where CRLBðuo Þ is given below (A.9). Using (C.2) and (A.2) in (C.1), we have from the first two elements ~ T Rn ~ Du ¼ CRLBðuo Þ½I YZ1 ½L~ P
ðC:3Þ
Hence ~ T E½DunT ¼ CRLBðuo Þ½I YZ1 ½L~ P
ðC:4Þ
1 has been used. Under C1 and C2 where E½nnT ¼ Q ¼ R~ ~ 1 P ¼ R~ P. ~ and substituting (21) and (34), BRP ¼ BB1 RB As a result, (33) becomes " # T L~ R~ P~ T CCRLBðuo Þ½I YZ1 ¼0 ðC:5Þ T P~ R~ P~
where the definitions of Y and Z in (A.4) have been used.
Contents lists available at SciVerse ScienceDirect
Signal Processing journal homepage: www.elsevier.com/locate/sigpro
Refining inaccurate sensor positions using target at unknown location M. Sun, K.C. Ho n Department of Electrical and Computer Engineering, University of Missouri, Columbia, MO 65211, USA
a r t i c l e i n f o
abstract
Article history: Received 17 August 2011 Received in revised form 19 January 2012 Accepted 26 January 2012 Available online 5 February 2012
This paper addresses the problem of improving the receiver positions in a sensor array using the positioning measurements of targets at unknown locations. Using an efficient estimate of the target locations, a computationally efficient estimator is proposed to refine the sensor positions using the same set of measurements from the targets. The considered measurements include TOA and TDOA. The proposed estimator has closedform and is able to reach the CRLB accuracy under small Gaussian noise which is supported by theoretical analysis and simulation studies. & 2012 Elsevier B.V. All rights reserved.
Keywords: Localization Sensor position errors Time difference of arrival (TDOA) Time of arrival (TOA)
1. Introduction Many engineering applications in practice requires the localization of a signal source, such as in radar and sonar [1,2], wireless communications [3] and more recently sensor networks [4,5]. The most common form of localization is to use a number of sensors to acquire the emitted source signal, compute the positioning measurements such as time of arrival (TOA) or time difference of arrival (TDOA) of the received signals at different sensors and estimate the source location from the positioning measurements. Traditionally, the sensor positions are assumed to be exactly known. Modern localization system often deploys the sensors randomly or uses dynamic sensors. In some applications such as sensor networks, quite often the sensors are without GPS to save costs and the sensor positions are not known precisely. Inaccurate sensor positions lead to degradation in localization performance if the sensor position errors are ignored [6,7]. A number of algorithms became
n
Corresponding author. Tel.: þ573 882 8023; fax: þ 573 882 0397. E-mail addresses: [email protected] (M. Sun), [email protected] (K.C. Ho). 0165-1684/$ - see front matter & 2012 Elsevier B.V. All rights reserved. doi:10.1016/j.sigpro.2012.01.022
available for localization problems with inexact sensor positions in literatures [7–11]. Methods [7,11] have been proposed to improve the source location estimate by taking into account the statistical behavior of the sensor position errors. However, relatively little emphasis has been allocated to improve the sensor positions so that they will be able to identify better the location of another source later. In this work, we are interested in improving the sensor positions given a set of positioning measurements from a single or multiple sources at unknown locations. A direct approach would be to jointly estimate the source (target) and sensor positions using the maximum likelihood (ML) method [12–14]. The non-linear nature of the problem makes it difficult to solve and it requires iterations with good initial guesses to obtain the ML estimate. Instead of solving the ML function directly, Yang and Ho [8] recently proposed an algebraic closed-form solution that jointly estimates the target and sensor positions. The target positions are able to achieve the Cramer Rao low bound (CRLB) accuracy, however, the sensor positions cannot. Another approach is to estimate the target locations first and then use the target location estimate to refine the sensor positions. This approach enables the target positions to be identified quickly for immediate actions with
2098
M. Sun, K.C. Ho / Signal Processing 92 (2012) 2097–2104
much less computations. Indeed, the solution in [7] estimates the target position only, which, through the design of a proper weighting matrix to account for the sensor position error, can reach the CRLB performance. After locating the target, the refinement of sensor positions can be followed later off-line when the sensing system is less congested to ease the computation burden. While the sensor position improvement stage in [15] by considering the unknown target as a calibration source can be applied to the current problem, it is not able to achieve the CRLB accuracy. In addition, the method in [15] assumes the position error of the target and the measurement noise are independent, which is clearly not the case for the current problem at hand because the target position is derived from the same set of measurements. This paper proposes a solution to achieve the CRLB performance in refining sensor positions. It takes the estimated target locations and the same set of measurements to improve the sensor positions. The proposed estimator is closed-form and does not require iterations or initialization. Perhaps the most important aspect of the estimator is that it reaches the CRLB accuracy for small Gaussian noise before thresholding effect occurs. This is proven analytically under small noise assumption and the performance is supported by simulations. We derive the proposed estimator for commonly used observations in localization that include TOA and TDOA. The measurement errors considered are Gaussian for ease of CRLB derivation and performance comparison, although the proposed estimator can be applied without Gaussian noise assumption. Nevertheless, the Gaussian measurement model has been validated under realistic environments in [16] and [17] and has been extensively used in literatures [7–9,11–17]. The development focuses on the scenario of a single target. The extension of the estimator to multiple target situation is also provided in which the performance can be shown to reach the CRLB accuracy as well. The proposed algorithm will be useful in applications such as sensor networks, where the sensors do not have GPS to identify their accurate positions, the computational resources are limited in locating targets and the refinement in sensor positions results in performance improvement of subsequent targets. In the following, Section 2 formulates the sensor position refinement problem mathematically. Section 3 develops the new estimator in closed-form. Section 4 conducts performance analysis with respect to the CRLB accuracy. The proposed solution is extended to multiple target scenario in Section 5. Section 6 presents the simulation results and Section 7 concludes the paper. 2. Problem formulation We shall consider 2D localization in this study for ease of illustration. There are M stationary sensors whose known but erroneous positions are si ¼ soi þ Dsi ,i ¼ 1; 2, . . . ,M, where soi ¼ ½xoi ,yoi T is the true position and Dsi is the error component. The collection of all the sensor positions forms s ¼ ½sT1 ,sT2 , . . . ,sTM T ¼ so þ Ds, where Ds ¼ ½DsT1 , DsT2 , . . . , DsTM T . The sensor position error Ds is assumed to be a zero-mean Gaussian random vector with covariance matrix Qs. Besides the sensors, the target appears at the position
uo ¼ ½xo ,yo T that is not known. The observations to locate the target are TOAs or TDOAs and they are represented by the measurement vector m ¼ mo þ nm
ð1Þ
o
where m is the noise free value. nm is the measurement noise and it is modeled as a zero-mean Gaussian random vector with covariance matrix Qm. The cross-covariance between nm and Ds is Qms. Qms in most cases has all zero elements. We include it here because the proposed method works even if correlation exists between the measurement noise and sensor position errors. For the case of TOA, mo ¼ ½r o1 ,r o2 , . . . ,r oM T and r oi is related to the unknown target position uo by r oi ¼ Juo soi J,
i ¼ 1; 2, . . . ,M
ð2Þ o
and JnJ represents the 2-norm. For TDOA, m . . . ,r oM1 T and r oi1 ¼ r oi r o1 ,
i ¼ 2; 3, . . . ,M
¼ ½r o21 ,r o31 , ð3Þ
where sensor 1 denotes the reference sensor in obtaining the TDOAs. Note that we have multiplied the TOA and TDOA with the known signal propagation speed when forming m. We shall suppose an estimate of the target position uo has been obtained using m that reaches CRLB accuracy. Several methods such as those in [7,11] are able to achieve this purpose. They do not need joint estimation of the sensor and target positions and can produce the target position estimate very fast. The requirement to achieve the CRLB is that the noise in the positioning measurements is small. Based on the target position estimate, the objective is to improve the sensor position s as much as possible using the same measurement vector m. 3. Proposed solution Instead of finding the sensor position so directly, the proposed solution estimates the error component Ds to correct s. The sensor position error Ds is a random vector. Hence we are estimating a particular realization of Ds. The set of pseudo-linear equations is composed of two parts: one is related to the positioning measurements and the other is from the probability distribution of the noisy sensor positions. By exploring the efficient estimate of the target position, a closed-form solution can be directly formulated. 3.1. Pseudo-linear equations 3.1.1. TOA Squaring both sides of TOA measurement r i ¼ r oi þ Dr i yields 2r oi Dr i C r 2i Juo soi J2
ð4Þ
where (2) has been used and ðDr i Þ2 has been ignored. The true sensor position soi is not available and only their noisy values are known. After substituting soi ¼ si Dsi , we can rewrite (4) as r oi Dr i C 12ðr 2i Juo si J2 Þðuo si ÞT Dsi
ð5Þ
M. Sun, K.C. Ho / Signal Processing 92 (2012) 2097–2104
when the second order error term DsTi Dsi is ignored. The matrix form of (5) is BT nm CqoT PoT Ds
ð6Þ
where BT ¼ diagðr o1 ,r o2 , . . . ,r oM Þ
ð7Þ
and diagða, b, . . . , gÞ denotes a block diagonal matrix with diagonal blocks a, b, . . . , g. On the right hand side, the M 1 vector qoT is qoT ¼ 12½r 21 Juo s1 J2 ,r 22 Juo s2 J2 , . . . ,r 2M Juo sM J2 T
ð8Þ
and PoT is an M 2M block-diagonal matrix defined as PoT ¼ diagððuo s1 ÞT ,ðuo s2 ÞT , . . . ,ðuo sM ÞT Þ:
3.1.2. TDOA Expressing the TDOA measurement ri1 as r i1 ¼ r oi1 þ Dr i1 and applying (3) yield r i1 þ r o1 ¼ r oi þ Dr i1 . Squaring both sides and using (2), we have o oT o 2r oi Dr i1 Cr 2i1 þ 2r i1 r o1 þ soT 1 s1 si si
i ¼ 2; 3, . . . ,M
ð10Þ
where ðDr i1 Þ2 has been neglected. We shall express soi as soi ¼ si Dsi . Different from the TOA case, we have the extra unknown variable r o1 that is dependent on uo and so1 . Using the Taylor series expansion around s1 and retaining up to linear error term, we have r o1 ¼ Juo so1 J C r~ o1 þ qTuo ,s1 Ds1 r~ o1
where the subscript m, which can be T or D, represents the measurement type. We shall use the common form (18) in subsequent development. Eq. (18) is constructed solely from the measurements. The statistical knowledge of the sensor position error Ds is exploited next. We follow the approach in [18] where the prior is utilized and by stacking (18) with ns ¼ Ds. Ds on the right-hand side is treated as the unknown parameter to be estimated and on the left-hand side ns is interpreted as a zero-mean random vector of position error with covariance matrix Qs. The resulting composite pseudo-linear equation is Bn C qo Po Ds
ð9Þ
Eq. (6) is interpreted as a pseudo-linear equation with respect to the unknown Ds since the true target position uo that appears in qoT and PoT is also not known. The lefthand side corresponds to the equation error due to the measurement noise nm.
þ 2ðsoi so1 ÞT uo ,
2099
ð11Þ
o
ð19Þ T T n ¼ ½nTm ,nTs T , qo ¼ ½qoT m , 02M1 o
where B ¼ diagðBm ,I2M2M Þ, T T and Po ¼ ½PoT and Ds m ,I2M2M . The unknowns are u (2Mþ2 variables) while Bn is the equation error. 3.2. Computing the sensor position estimate Defining the cost function to be minimized as J ¼ ðqo Po DsÞT Rðqo Po DsÞ
ð20Þ
where R is the weighting matrix obtained by taking the inverse of the equation error covariance matrix [19], i.e. R ¼ E½ðBnÞðBnÞT 1 ¼ ðBQBT Þ1 , " # Q m Q ms T Q ¼ E½nn ¼ Q Tms Q s
ð21Þ
R is always invertible as long as the target is not colocated with a sensor and the measurement noise covariance matrix Q is full rank, which is expected to be the case in practice. We shall denote the solution that minimizes (20) as u and Ds^ . Obviously, they satisfy @J ¼0 ð22Þ @Ds
¼ Ju s1 J and qa,b ¼ ðabÞ=JabJ represents where the unit vector pointing from b to a. As a result, (10) reduces to
Evaluating the above gradient and simplifying, we have
r oi Dr i1 Cqoi þ ðr i1 qTuo ,s1 þ ðuo s1 ÞT ÞDs1 ðuo si ÞT Dsi
Ds^ ¼ ðPT RPÞ1 PT Rq
ð12Þ
ð23Þ o
where qoi ¼
u, Ds^
T o T T 1 2 ~o 2ðr i1 þ 2r i1 r 1 þs1 s1 si si þ 2ðsi s1 Þ u Þ
ð13Þ
Stacking (12) for i¼2,3,y,M together gives the matrix equation BD nm C qoD PoD Ds
ð14Þ
where BD ¼ diagðr o2 ,r o3 , . . . ,r oM Þ
ð15Þ
qoD ¼ ½qo2 ,qo3 , . . . ,qoM T
ð16Þ
and the (i 1)th row of the (M 1) 2M matrix o
½ðr i1 quo ,s1 þ ðu s1 ÞÞ i ¼ 2; 3, . . . ,M
T
0T2ði2Þ1
o
ðu si Þ
T
PoD
is
o
where P and q are P and q with u replaced by u. Eq. (23) indicates that a closed-form solution of Ds can be directly obtained when u is available, which can come from another estimator as long as it gives the efficient target position estimate. u can be represented as u ¼ uo þ Du and its covariance is covðuÞ ¼ CRLBðuo Þ, where Du is the estimation error. Finding Ds using (23) requires the weighting matrix R. The matrix B in the weighting matrix R depends on the true sensor–source distances and they are replaced by the computed distances using the noisy sensor position s and the source estimate u. After Ds^ is obtained from (23), the estimate of so is s^ ¼ sDs^
ð24Þ
0T2ðMiÞ1 , ð17Þ
Regardless of the measurement type as TOA or TDOA, the solution equation for Ds has the generic form Bm nm Cqom Pom Ds
o
ð18Þ
4. Performance analysis In this Section, the performance of the proposed algorithm is evaluated and compared with the CRLB accuracy. The CRLB of the sensor positions is derived in Appendix A
2100
M. Sun, K.C. Ho / Signal Processing 92 (2012) 2097–2104
and it is given in (A.9). We first express qo in (19), by substituting uo ¼ uDu, as
of the proposed estimator is approximately efficient under C1 and C2
qo CqLDu
covðs^ Þ C CRLBðso Þ
ð25Þ
where the second order error terms have been ignored and L ¼ ½LTm ,0T2M2 T where 0 is used to denote a matrix of all zeros. For the case of TOA, using the definition of qoT in (8) gives Lm ¼ LT ¼ ½ðus1 Þ ðus2 Þ ðusM ÞT
ð26Þ
For the case of TDOA, applying the approximation that r~ o1 ¼ Juo s1 J C r~ 1 qTu,s1 Du where r~ 1 ¼ Jus1 J, we have from the definition of the elements of qoD in (16) Lm ¼ LD ¼ ½r 21 qu,s1 þðs2 s1 Þ r 31 qu,s1 þðs3 s1 Þ r M1 qu,s1 þ ðsM s1 ÞT
ð27Þ
Eq. (19) can now be rewritten, by substituting (25) and Po ¼ PDP, as q CPDsþ Bn þ LDu
ð28Þ
where the second order error term DPDs has been ignored. After substituting (28) in (23), the estimation error of Ds^ can be expressed as e9Ds^ Ds ¼ Z^
1
PT RBn þ Z^
1
T Y^ Du
ð29Þ
Y^ ¼ LT RP
ð30Þ
We now evaluate the covariance matrix of the proposed estimator. Eq. (24) can be expressed as s^ ¼ so þ DsDs^ . Hence the covariance matrix is covðs^ Þ ¼ covðDs^ Þ ¼ E½eeT
ð31Þ
From (29) and noting that R ¼ ðBQBT Þ1 , we have 1 1 T 1 1 T ^ Z^ 1 covðs^ Þ ¼ Z^ þ Z^ Y^ covðuÞY^ Z^ þ Z^ ðY^ T þ TT YÞ
ð32Þ o
where covðuÞ ¼ CRLBðu Þ and T9E½Du nT BRP
ð33Þ
We would like to prove that (32) and (A.9) are nearly identical under the following two small noise conditions C1: 9nm,i 9=r oi C 0; C2: JDsi J=r oi C 0, for i ¼ 1; 2, . . . ,M
The sensor position estimator and its performance analysis presented so far only use measurements from a single target. We shall extend the developed algorithm to the multiple target scenario and validate that the obtained sensor position estimate remains to be nearly efficient under C1 and C2. Suppose there are K targets whose true positions oT oT T uo ¼ ½uoT 1 ,u2 , . . . ,uK are not known. The measurement vector becomes m ¼ ½mT1 ,mT2 , . . . ,mTK T ¼ mo þ nm , where nm ¼ ½nTm,1 ,nTm,2 , . . . ,nTm,K T . Given the efficient estimate of the target position vector u ¼ uo þ Du obtained from the measurement vector m, we are interested to refine the sensor positions. Recalling (18) is the equation for the measurements from a single target, we rewrite it here for target i as Bm,i nm,i C qom,i Pom,i Ds
ð36Þ
B1 P C P~
Bm nm Cqom Pom Ds
ð34Þ
As a result of (34) and using (21), we have Z^ C Z and Y^ C Y from (30) and (A.4). Also, we have shown in Appendix C that if (34) is valid, T C 0. Thus, using an estimate u with CRLB accuracy for the proposed algorithm, its covariance matrix (32) is nearly identical to (A.9) and the performance
ð37Þ
where Bm ¼ diagðBm,1 ,Bm,2 , . . . ,Bm,K Þ. On the right-hand oT oT T and Pom ¼ ½PoT side of (37), qom ¼ ½qoT m,1 , m,1 ,qm,2 , . . . ,qm,K oT T PoT m,2 , . . . ,Pm,K . Combining (37) and the statistical knowledge of the sensor position error yields the pseudo-linear equation having the same form as (19). Once the pseudo-linear equation (19) is formed, the procedure to compute the sensor position estimate is straightforward and is exactly the same as the one developed in Section 3.B. The performance analysis follows the same as in Section 4. The only difference is that L in (25) is given by L ¼ ½LTm ,0T2M2K T and Lm ¼ diagðLm,1 ,Lm,2 , . . . ,Lm,K Þ. Lm,i is defined as in (26) for TOA or (27) for TDOA with u replaced by ui . From Appendix C, under C1 and C2, we have for i¼1, 2,y,K, o o B1 m,i L m,i C ð@mi =@ui Þ,
where nm,i is the ith element of nm in (1). C1 states the measurement noise is small compared to the distance between the target and the sensor. C2 means that the 2norm of the ith sensor position error is much smaller than the distance between the target and the sensor. Both conditions are satisfied if the noise is sufficiently small or the target is far away from the sensor array. Appendix B shows that under C1 and C2, ~ B1 L CL,
5. Multiple target scenario
Stacking (36) for i¼1,2,y,K gives the equation of all measurements,
where Z^ ¼ PT RP,
ð35Þ
o o B1 m,i Pm,i C ð@mi =@s Þ
ð38Þ
and (34) is valid. As a result, (35) holds true for the case of multiple targets as well. 6. Simulations The localization system in the simulation study is composed of M ¼4 sensors and their true positions are so1 ¼ ½8,16T ,so2 ¼ ½28,40T ,so3 ¼ ½20,8T and so4 ¼ ½0; 8T as shown in Fig. 1. The measurement noise covariance matrix is Q m ¼ 0:001IMM for TOA. For TDOA, the covariance matrix is Qm ¼0.001J where J is an (M 1) (M 1) matrix with diagonal elements equal to 1 and all other elements 12. This form of covariance assumes ML estimation of TDOAs whose accuracy is independent of the choice of reference sensor [20]. The measurement noise and sensor
M. Sun, K.C. Ho / Signal Processing 92 (2012) 2097–2104
20
2101
25
Targets Sensors
Proposed method Before refinement CRLB
20 10log (Sensor position MSE)
10
y
0
−10
−20
−30
−40 −20
15
10
5
0
−5 −10
0
10
20
30
40
50
−5
x
position errors are independent so that Qms ¼0. The sensor positions are not known at the beginning and they are obtained through three emitting sources at known positions ao1 ¼ ½19:4,2:8T , ao2 ¼ ½1:8,4:8T and ao3 ¼ ½5; 17T from TOA measurements. We also assume the sensors are active and can communicate with each other so that TOAs between so1 and so4 and between so2 and so3 are also available to identify their positions. For illustration purpose, the TOA noise values to obtain the sensor positions are independent and identically distributed zero-mean Gaussian random variables with noise power s2r . The ML estimator is applied to obtain the sensor position vector s using these TOAs. The covariance matrix Qs is set to the associated CRLB, where the true sensor positions in the bound are replaced by the estimated values. The initial sensor position vector therefore has error and the amount of sensor position errors is proportional to s2r . At a given signal-to-noise ratio (SNR), the corresponding noise power is obtained using (2) of [4] with the effective signal bandwidth set to 20 MHz. The localization process proceeds as follows: (i) obtain the sensor position vector s using the ML method based on the TOAs from ao1 ,ao2 ,ao3 and among the sensors; (ii) localize the targets using the TOA/TDOA measurements between the targets and the sensors, where the method in [11] is used as an example to produce the efficient target position estimate u; (iii) apply the proposed algorithm using the estimated target locations and the same set of measurements to improve the sensor position vector s. P l^ o 2 The performance measure is MSEðs^ Þ ¼ N l ¼ 1 J s s J =N, o l^ where s is the estimate of s at ensemble l and N ¼2000 is the number of ensemble runs. We first examine the performance of the proposed method using a single target located at unknown position uo1 ¼ ½4,24T . The results are shown in Fig. 2 for TOA and Fig. 3 for TDOA. Each figure plots the MSEðs^ Þ and the trace of CRLBðso Þ as a function of SNR, before and after applying the proposed sensor refinement method. The proposed estimator can improve the sensor locations by about 3 dB for TOA and 1 dB TDOA. Furthermore, it is able to achieve the CRLB accuracy before the thresholding effect occurs. In Fig. 3, the sensor position estimate from [8] is also plotted for comparison. The MSE of the sensor position
5
10
15
20
SNR (dB)
Fig. 2. Performance of the proposed method in improving the sensor positions using TOA measurements from an unknown target whose true location is uo1 . 20 Proposed method Solution in [8] Before refinement CRLB
15 10log (Sensor position MSE)
Fig. 1. Sensor network configuration.
0
10
5
0
−5
0
5
10
15
20
SNR (dB)
Fig. 3. Performance of the proposed method in improving the sensor positions using TDOA measurements from an unknown target whose true location is uo1 .
estimate almost overlaps with the variance of the sensor positions before estimation. This means that the solution in [8] is not able to improve the sensor positions much as expected. It is observed though that the proposed method has lower noise threshold than [8] before the performance deviates from the bound. The computational complexity of the two methods for the case of TDOA is also examined. Both the solution in [8] and the proposed estimator involve the inverse of a matrix that accounts for most of the computation burden. The complexity of the method in [8] is dominated by its first stage WLS computation that jointly estimates the source and sensor positions. The size of the matrix that needs to be inverted for 2D localization is 3Kþ2M, where 2K is for the K source coordinates, another K is for K auxiliary variables one for each source and 2M comes from the estimation of M sensor positions. However, the proposed algorithm only requires to estimate the 2M sensor positions and the inverse matrix has a size of 2M only. The positions of the K sources can be
M. Sun, K.C. Ho / Signal Processing 92 (2012) 2097–2104
estimated with a matrix inverse of size 3K [11]. Therefore the computation complexity is Oðð3K þ2MÞ3 Þ for the method in [8] and Oðð3KÞ3 þ ð2MÞ3 Þ for the proposed estimator. The proposed estimator has lower complexity than [8]. The reduction in sensor position errors can be represented by the concentration ellipses. Fig. 4 gives the concentration ellipses of the sensor positions before (dashed line) and after (solid line) improvement using TOAs measurements where SNR is set to be 4.55 dB. Figs. 5 and 6 present the results using another target at uo2 ¼ ½49; 19T . Note that different from the target at uo1 , the target at uo2 is outside of the area enclosed by the sensors. Compared to the case that target is at uo1 , the improvement in sensor positions is less. Nevertheless, the proposed method remains to reach the CRLB accuracy. The performance of the proposed method using multiple targets is also evaluated. Figs. 7 and 8 present the results with two targets at uo1 ¼ ½4,24T and uo2 ¼ ½49; 19T . The measurements are TOA and they are independent from the two targets. Compared to the single target case, the
20
10
5
0
−5
0
5
20
25
y
−10
−20
−30
Proposed method Before refinement CRLB
20 10log (Sensor position MSE)
0
−40
15
10
5
0
−5
−10
0
10
20
30
−5
40
0
5
10
15
20
SNR (dB)
x
Fig. 4. Sensor position improvement illustration using uncertainty ellipse.
Fig. 7. Performance of the proposed method in improving the sensor positions using TOA measurements from two unknown targets whose true locations are uo1 and uo2 .
25
20
Proposed method Before refinement CRLB
20
Proposed method Solution in [8] Before refinement CRLB
15 10log (Sensor position MSE)
10log (Sensor position MSE)
15
Fig. 6. Performance of the proposed method in improving the sensor positions using TDOA measurements from an unknown target whose true location is uo2 .
Target Sensors Before refinement After refinement
10
15
10
5
10
5
0
0
−5
10 SNR (dB)
20
−50 −20
Proposed method Solution in [8] Before refinement CRLB
15 10log (Sensor position MSE)
2102
−5
0
5
10
15
20
SNR (dB)
Fig. 5. Performance of the proposed method in improving the sensor positions using TOA measurements from an unknown target whose true location is uo2 .
−5
0
5
10
15
20
SNR (dB)
Fig. 8. Performance of the proposed method in improving the sensor positions using TDOA measurements from two unknown targets whose true locations are uo1 and uo2 .
M. Sun, K.C. Ho / Signal Processing 92 (2012) 2097–2104
25 Target position MSE, without refining sensor positions Target position CRLB, without refining sensor positions Target position MSE, with refining sensor positions Target position CRLB, with refining sensor positions
10log (Target position MSE)
20
15
2103
where c is a constant not depending on ho . We can obtain from (A.1) @ ln pðx; ho Þ o ~ T Rðxx ~ ¼ ½L~ P Þ @ho
ðA:2Þ
where 2 3 @mo o ~L9 @x ¼ 4 @uo 5, @uo 0
10
2 3 @mo @xo 4 ~ o P9 o ¼ @s 5, @s I
2M2
5
R~ ¼ Q 1
2M2M
ðA:3Þ 0
−5 −5
0
5
10
15
20
Hence, the Fisher information matrix (FIM) is, after postmultiplying (A.2) by its transpose and taking expectation, # 2 " T T X Y @ ln pðx; ho Þ L~ R~ L~ L~ R~ P~ ¼ FIMðho Þ9E ¼ T T YT Z @ho @hoT P~ R~ L~ P~ R~ P~
SNR (dB)
ðA:4Þ
Fig. 9. Localization accuracy of a new subsequent target with and without improving the sensor positions from a previously observed target at uo1 , TOA measurements.
improvement in sensor positions is larger as expected. The proposed method remains to reach the CRLB accuracy. The main purpose of improving sensor positions is to increase the localization accuracy in the next round of target positioning. To illustrate this, Fig. 9 shows the localization result of a new target at ½36; 39T with (crosses) and without (squares) improving the sensor positions from the single target at uo1 ¼ ½4,24T that is present earlier. The measurement is TOA and the method in [11] is used to locate the new target. It is clear that there is about 5 dB performance improvement when using the refined sensor positions. We would like to clarify that the proposed algorithm is intended for use when the SNR is high. When the SNR is low, simulation shows its performance deviates from the CRLB.
The partial derivatives are given below. For TOA, the ith row of ð@mo =@uo Þ is ð@r oi =@uo ÞT ¼ qTuo ,so , i
i ¼ 1; 2, . . . ,M
ðA:5Þ
and ð@mo =@so Þ is a block diagonal matrix given by ð@mo =@so Þ ¼ diagðqTuo ,so , qTuo ,so , . . . , qTuo ,so Þ 1
ðA:6Þ
M
2
For TDOA, the (i 1)th rows of ð@mo =@uo Þ and ð@mo =@so Þ are ð@r oi1 =@uo ÞT ¼ ½quo ,so quo ,so T , 1
i
i ¼ 2; 3, . . . ,M
ðA:7Þ
ð@r oi1 =@so ÞT ¼ ½qTuo ,so 0T2ði2Þ1 qTuo ,so 0T2ðMiÞ1 , 1
i
i ¼ 2; 3, . . . ,M
ðA:8Þ
Invoking the block matrix inversion formula in (A.4), the CRLBðso Þ is CRLBðso Þ ¼ ðZYT X1 YÞ1 ¼ Z1 þZ1 YT CRLBðuo ÞYZ1 ðA:9Þ
7. Conclusions
where CRLBðu Þ ¼ ðXYZ
In this paper, an explicit solution was proposed to refine the sensor positions using TOA or TDOA measurements from a single or multiple targets at unknown locations. Suppose that an efficient estimate of target locations is obtained by some techniques, the proposed method refines the sensor positions using the same set of measurements from the unknown targets. This approach enables the target positions to be identified quickly first with much less computations while leaving the improvement of sensor positions at a later time when computation resources are available. Performance analysis and simulations showed that the proposed method is able to achieve the CRLB accuracy under small Gaussian noise.
T 1
Y Þ
.
Appendix B. Proof of (34) In (34), the definitions of B, P and L are under (19) and (25). We begin with the proof of (34) for TOA. Bm and Pm are BT, LT and PT defined below (6) and LT is given in (26). After substituting their expressions and simplifying yield us us2 usM T 1 B1 L ¼ ðB:1Þ o o T T r r ro 1
2
M
T T T o o o B1 T PT ¼ diagððus1 Þ =r 1 ,ðus2 Þ =r 2 , . . . ,ðusM Þ =r M Þ
ðB:2Þ From C1 and C2 and based on the form of Du in (C.3), JDuJ=r oi C 0. Hence under both conditions,
Appendix A. Derivation of the CRLB From the noise model, the logarithm of the density function of the data vector x ¼ ½mT ,sT T parameterized on the unknown vector ho ¼ ½uoT ,soT T is ln pðx; ho Þ ¼ c12ðxxo ÞT Q 1 ðxxo Þ
1
o
ðA:1Þ
ðusi Þ=r oi ¼ ðuo soi þ DuDsi Þ=r oi C quo ,so , i
i ¼ 1; 2, . . . ,M ðB:3Þ
Putting (B.3) back into (B.1) and (B.2) and comparing with (A.5) and (A.6) validate (34).
2104
M. Sun, K.C. Ho / Signal Processing 92 (2012) 2097–2104
For TDOA, after using BD and PD given below (14) and LD in (27), the (i 1)th row of B1 D L D is
References
ðr i1 qTu,s1 þ ðsi s1 ÞT Þ=r oi ,
[1] F. Ahmad, M. Amin, Noncoherent approach to through-the-wall radar localization, IEEE Transactions on Aerospace and Electronic Systems 42 (2006) 1405–1419. [2] S. Coraluppi, Multistatic sonar localization, IEEE Journal of Oceanic Engineering 31 (2006) 964–974. [3] G. Sun, J. Chen, W. Guo, K.J.R. Liu, Signal processing techniques in network-aided positioning: a survey of state-of-the-art positioning designs, IEEE Signal Processing Magazine 22 (2005) 12–23. [4] S. Gezici, T. Zhi, G. Giannakis, H. Kobayashi, A. Molisch, H. Poor, Z. Sahinoglu, Localization via ultra-wideband radios: a look at positioning aspects for future sensor networks, IEEE Signal Processing Magazine 22 (2005) 70–84. [5] T. Li, A. Ekpenyong, Y.-F. Huang, Source localization and tracking using distributed asynchronous sensors, IEEE Transactions on Signal Processing 54 (2006) 3991–4003. [6] A. Ferreol, P. Larzabal, M. Viberg, On the asymptotic performance analysis of subspace DOA estimation in the presence of modeling errors: case of MUSIC, IEEE Transactions on Signal Processing 54 (2006) 907–920. [7] K.C. Ho, X. Lu, L. Kovavisaruch, Source localization using TDOA and FDOA measurements in the presence of receiver location errors: analysis and solution, IEEE Transactions on Signal Processing 55 (2007) 684–696. [8] L. Yang, K.C. Ho, An approximately efficient TDOA localization algorithm in closed-form for locating multiple disjoint sources with erroneous sensor positions, IEEE Transactions on Signal Processing 57 (2009) 4598–4615. [9] K.W.K. Lui, W.-K. Ma, H.C. So, F.K.W. Chan, Semidefinite programming approach to sensor network node localization with anchor position uncertainty, in: Proceedings of the IEEE International Conference on Acoustics, Speech, and Signal Processing (ICASSP’09), Taipei, Taiwan, 2009, pp. 2245–2248. [10] W. Chiu, B. Chen, C. Yang, Robust relative location estimation in wireless sensor networks with inexact position problems, IEEE Transactions on Mobile Computing 99 (2011) 1. [11] M. Sun, K.C. Ho, An asymptotically efficient estimator for TDOA and FDOA positioning of multiple disjoint sources in the presence of sensor location uncertainties, IEEE Transactions on Signal Processing 59 (2011) 3434–3440. [12] Y. Rockah, P.M. Schultheiss, Array shape calibration using sources in unknown location. Part I: far-field source, IEEE Transactions on Acoustics, Speech, and Signal Processing ASSP-35 (1987) 286–299. [13] Y. Rockah, P.M. Schultheiss, Array shape calibration using sources in unknown location. Part II: near-field source and estimator implementation, IEEE Transactions on Acoustics, Speech, and Signal Processing ASSP-35 (1987) 724–735. [14] A. Weiss, B. Friedlander, Array shape calibration using sources in unknown locations: a maximum likelihood approach, IEEE Transactions on Acoustics, Speech, and Signal Processing ASSP-37 (1989) 1958–1966. [15] L. Yang, K.C. Ho, Alleviating sensor position error in source localization using calibration emitters at inaccurate locations, IEEE Transactions on Signal Processing 58 (2010) 67–83. [16] N. Patwari, J.N. Ash, S. Kyperountas, A.O. Hero III, R.L. Moses, N.S. Correal, Locating the nodes: cooperative localization in wireless sensor networks, IEEE Signal Processing Magazine 22 (2005) 54–69. [17] B. Alavi, K. Pahlavan, Modeling of the TOA-based distance measurement error using UWB indoor radio measurements, IEEE Communications Letters 10 (2006) 275–277. [18] H.W. Sorenson, Parameter Estimation: Principles and Problems, Marcel Dekker, New York, 1980. [19] S.M. Kay, Fundamentals of Statistical Signal Processing, Estimation Theory, Prentice Hall, Englewood Cliffs, NJ, 1993. [20] K.C. Ho, W. Xu, An accurate algebraic solution for moving source location using TDOA and FDOA measurements, IEEE Transactions on Signal Processing 52 (2004) 2453–2463.
i ¼ 2; 3, . . . ,M
ðB:4Þ
and the (i 1)th row of B1 D PD is ½ðr i1 qTu,s1 þ ðus1 ÞT Þ=r oi 0T2ði2Þ1 ðusi ÞT =r oi 0T2ðMiÞ1 , i ¼ 2; 3, . . . ,M
ðB:5Þ
It can be easily verified that under C1 and C2, r~ 1 C r o1 and r i1 C r oi1 . Using these approximations and (B.3) yield ðr i1 qTu,s1 þ ðsi s1 ÞT Þ=r oi C ðr oi1 ðus1 ÞT þ r o1 ðsi s1 ÞT Þ=ðr oi r o1 Þ C qTuo ,so qTuo ,so 1
ðB:6aÞ
i
ðr i1 qTu,s1 þ ðus1 ÞT Þ=r oi C ðr oi1 ðus1 ÞT þ r o1 ðus1 ÞT Þ=ðr oi r o1 Þ CqTuo ,so
1
ðB:6bÞ
Applying (B.6) to (B.4) and (B.5) and (B.5) and comparing with (A.7) and (A.8), we have the results in (34). Appendix C. Value of T under conditions C1 and C2 We shall express Du in terms of n ¼ ½nTm ,nTs T . According to the CRLB theorem [19], an unbiased estimator may be found that attains the bound in that covðhÞ ¼ FIMðho Þ1 if and only if ð@ ln pðx; ho Þ=@ho Þ ¼ FIMðho ÞðgðxÞho Þ and h ¼ gðxÞ is the efficient estimate of ho . Thus, for an estimator that yields the CRLB accuracy, the estimation error must be of the form
Dh ¼ gðxÞho ¼ FIMðho Þ1 ð@ ln pðx; ho Þ=@ho Þ
ðC:1Þ
Note that in our case, u is assumed to be an efficient estimate of uo obtained from existing literatures, such as [8,11]. Therefore, it corresponds to the first two elements of gðxÞ with covðuÞ ¼ CRLBðuo Þ. FIMðho Þ1 in (C.1) can be expressed as, when applying the block matrix inversion formula to (A.4), " # CRLBðuo ÞYZ1 CRLBðuo Þ FIMðho Þ1 ¼ ðC:2Þ 1 T Z Y CRLBðuo Þ CRLBðso Þ where CRLBðuo Þ is given below (A.9). Using (C.2) and (A.2) in (C.1), we have from the first two elements ~ T Rn ~ Du ¼ CRLBðuo Þ½I YZ1 ½L~ P
ðC:3Þ
Hence ~ T E½DunT ¼ CRLBðuo Þ½I YZ1 ½L~ P
ðC:4Þ
1 has been used. Under C1 and C2 where E½nnT ¼ Q ¼ R~ ~ 1 P ¼ R~ P. ~ and substituting (21) and (34), BRP ¼ BB1 RB As a result, (33) becomes " # T L~ R~ P~ T CCRLBðuo Þ½I YZ1 ¼0 ðC:5Þ T P~ R~ P~
where the definitions of Y and Z in (A.4) have been used.