基於修改線性規劃法及偏差修正法之無線感測網路定位

(1)

國立

交

通

大

學

電信工程研究所

碩

士

論

文

基於修改線性規劃法及偏差修正法之無線

感測網路定位

Modified Linear Programming and Bias

Correction in Wireless Sensor Networks

研究生：邱祺祥

指導教授：謝世福教授

(2)

2

基於修改線性規劃法及偏差修正法之無線感測網路

定位

Modified Linear Programming and bias correction in

Wireless Sensor Networks

研究生：邱祺祥 Student：C.H. Chiu

指導教授：謝世福 Advisor：S. F. Hsieh

國立交通大學

電信工程研究所

碩士論文

A Thesis

Submitted to Department of Electrical Engineering Institute of Communications Engineering

National Chiao Tung University In Partial Fulfillment of the Requirements

For the Degree of Master of Science

In

Communication Engineering

Hsinchu 2012

Hsinchu, Taiwan, Republic of China

(3)

i

基於修改線性規劃法及偏差修正法之無線感測

網路定位

學生:邱祺祥指導教授:謝世福

國立交通大學電信工程研究所

中文摘要

在無線感測網路發展快速的現今，其應用的範圍十分廣泛，其中又以定位為一個熱門的議題。定位可以應用在室內以及戶外區域，藉由基地台接收手機或是待測物的訊號去推測出之間的距離，來求得最終待測物的位置。而室內或是大都會定位系統中，最容易造成較大定位誤差的干擾為 NLOS(非直視性)情形，由於手機待測物及基地台之間有障礙物，因此使得量測的距離有過長的現象。NLOS 效應之減緩可以分成兩種方向：知道統計特性及不知道統計特性。不知道 NLOS 統計特性情況下則是藉由非線性最小方差將 NLOS 偏差當作未知數估計，或是藉由一組粗估位置尋求殘餘量中帶有 NLOS 偏差的訊息，並且加以估計 NLOS 偏差量來修正開始粗估的位置，過去也有提出利用線性規劃法利用 LOS 量測建立目標方程式，NLOS 量測建立可行域來輔助定位。本論文將著重於不知道 NLOS 統計特性的情形，修改過去所提出的基於泰勒線性化並且帶限制式的 NLOS 偏差修正方法，藉由不考慮 NLOS 限制式，能有較短的求解時間以換取定位準確度，以及用泰勒線性以及新增變數法來修改線性規劃法中基於相減法的目標方程式，並且進一步提出基於泰勒線性化法之可行域之建立，改善定位準確度，最後則是將此可行域之建立，套用在現有線性規畫法之延伸，最後將藉由電腦模擬驗證我們所提出之方法。

(4)

ii

Modified Linear Programming and bias correction in

Wireless Sensor Networks

Student: C. H. Chiu Advisor：S. F. Hsieh

Department of Communication Engineering

National Chiao Tung University

Abstract

Localization in wireless sensor network in a Non-Line-of-Sight environment, , has become an important issue recent years. NLOS (Non0Line-of-Sight) is one of the most crucial

biases when localizing the mobile at indoor or urban area. It can cause large errors in our solution of localization because the measurements are much longer than the actual distance.

There are two common ways to mitigate the NLOS situation: known and unknown statistics. In the past, there are some localization algorithms with unknown statistics. Use nonlinear

weighted LS treating NLOS biases into unknowns with constraints. Or roughly localize the mobile pretending all the BS are LOS, and use residual concept to estimate the NLOS bias.

Then correct the “rough” position. Linear programming use LOS measurements constructing cost function, and use NLOS measurements build feasible region to help positioning. This

thesis will focus on the unknown statistic properties. By using different linearization methods and without constraints to modify the bias correction method. Shortening the computation

time in exchange for positioning accuracy. And a significant improvement can be achieved by modifying the cost function of the linear programming method and the establishment of the

(5)

iii

致謝

轉眼間，碩士生涯即將告一段落，回想當初剛進來交大時懵懵懂懂的我，更有迷失在做研究的路途，一直到今天能完成這本論文，有太多太多的人需要感謝。首先謝謝我的家人，在我漫漫求學生涯中經濟上的無虞並且不時的給予我支持與鼓勵。感謝我的指導教授，用一個嚴謹的態度教導我去面對研究上的困難。感謝我最好的朋友們，俊佑、梅爾許傑巽、張爸、jamie 大大、言異常、冠杰、阿邦、台妹及美月能讓我在新竹的研究生活中多了很多不同，也感謝假掰銘、瘋婷、賭懶民、資順、布萊德秉璁、韋廷、玉莘這些不在新竹的好朋友們一起瘋過的時間，更要感謝依蓁在這些時間，一起完成碩士學業，並互相鼓勵以及走過很多的人生旅途，也謝謝吳爸吳媽的照顧，讓我在新竹有著另一個家的感覺。回首這些日子，點點滴滴浮現在心中，研究生涯只占了人生的一小部分，這幾年中讓我對於自己審視了很多不足，不久後將踏入人生的另一個階段，脫離學生生活，相信更是一個不一樣的旅程，有著更多的困難與挑戰，想起一句話，人生的路途上我們一定會遇到困難和挫折，沒有人會注意你怎麼摔倒，但一定會有人注意你怎麼爬起來，希望未來也能莫忘初衷，這些年所學習到的一切，所遇過的挫折，都能幫助我在未來人生中有著更堅定的心智去面對所有。

(6)

iv

中文摘要

………i

Abstract

………..…ii

致謝

...……….iii

圖片目錄

……….………vi

表格目錄

……….…..viii

1. 簡介

………...1

2. 定位演算法

………...……….4 2.1 量測……….4 2.2 直視性情況下之 TOA 定位演算法………..…….5 2.2.1 最大概似機率估計………...……….…….6 2.2.2 線性最小方差估計…………...………7 2.2.2a 相減法………7 2.2.2b 泰勒線性化法………8 2.2.2c 新增變數法...10 2.2.3 利用一般化最小方差法作權重 C(N,2)相減法...11

3. NLOS 情況下利用偏差修正法之定位

…... .16 3.1 權重非線性最小方差……...16 3.2 帶限制式之 NLOS 偏差修正法...…...17

(7)

v 3.3 修改 NLOS 偏差修正法…...21

4.NLOS 情況下利用修改線性規劃法之定位

...25 4.1 現有線性規劃法介紹...25 4.2 利用其他線性化方法建立線性規劃目標方程式...30 4.2a 利用泰勒線性化法建立目標方程式...30 4.2b 利用新增變數法建立目標方程式...33 4.3 利用泰勒線性化建立可行域……... 34 4.4 現有線性規劃法延伸...42

5.電腦模擬

...48 5.1 利用 GLS 作權重 C(N,2)相減法...48 5.2 修改 NLOS 修正偏差法...50 5.3 線性規劃法之比較...53 5.3.1 修改目標方程式之線性規劃法...53 5.3.2 利用泰勒線性化形成可行域之比較...55 5.4 利用人工 NLOS 擴張線性規劃法比較...58

6.結論

...61

參考文獻

...62

(8)

vi

圖片目錄

圖 1.1 無線感測網路定位系統示意圖...1 圖 2.1 LOS 情況之 TOA 量測示意圖...5 圖 3.1 手機及基地台間 NLOS 發生情況示意圖...16 圖 3.2 帶限制式之 NLOS 偏差估計修正流程圖...22 圖 3.3 修改 NLOS 偏差修正法之流程圖...23 圖 4.1 矩形放鬆之可行域示意圖...28 圖 4.2 泰勒線性化之對於圓之幾何關係...35 圖 4.3 泰勒線性化之可行域示意圖...35 圖 4.4 原始，矩形放寬，泰勒放寬之可行域示意圖比較之一...36 圖 4.5 兩個 NLOS 基地台情況下封閉區間可能比較...38 圖 4.6 原始，矩形放寬，泰勒放寬之可行域示意圖比較之二...38 圖 4.7 泰勒線性化之非封閉可行域示意圖...39 圖 4.8 原始，矩形放寬，泰勒放寬之可行域示意圖比較之三...39 圖 4.9 矩形放寬結合泰勒放寬之可行域示意圖...40 圖 4.10 線性規劃目標方程式及可行域之流程圖...42 圖 4.11 加上人工 NLOS 後放大的 LOS 示意圖...43 圖 4.12 加上一組人工 NLOS 後可行域形成示意圖...43 圖 4.13 矩形放寬加上一組人工 NLOS 可行域之變化...44 圖 4.14 泰勒放寬加上一組人工 NLOS 可行域之變化...44 圖 5.1 C(N,2)相減法與利用一般化最小方差作權重之 C(N,2)相減法比較...49 圖 5.2 帶限制式 NLOS 偏差修正跟修改過後偏差修正方法比較...51 圖 5.3 改變 NLOS 個數針對定位誤差<5m 比例之比較...52 圖 5.4 模擬環境之基地台分布...54

(9)

vii 圖 5.5 不同線性化建立之線性規劃法比較...55 圖 5.6 不同可行域之線性規劃法比較...56 圖 5.7 固定 NLOS 統計特性改變雜訊之 MSE 比較...57 圖 5.8 固定雜訊統計特性改變 NLOS 之 MSE 比較...58 圖 5.9 矩形及泰勒之擴增比較...59 圖 5.10 矩形及泰勒之擴增比較 (泰勒差) ...60

(10)

viii

表格目錄

(11)

1

第一章

簡介

在最近幾年，隨著無線網路的蓬勃發展，許多無線網路之應用[1]也成為熱門的話題，包含了醫療，搜索及商用等，目前被廣泛使用的無線定位系統為

GPS(global position system)，GPS 利用了人造衛星以及待測者之間的資訊來找出待測者的位置，而在一般使用者中，利用 GPS 需要較大的代價，因此目前現有的研究皆是如何從少數的 GPS 系統或式已知基地台位置，藉由跟未知手機或是待測者的位置距離關係量測，來找出手機或待測者的最終估計。，以圖 1.1 為例 BS 1、2、…N 為移動式基地台，MS 1、2…M 為手機。圖 1.1 無線感測網路定位系統示意圖較常見的是利用 TOA(Time of Arrival)[2]方法來量測手機跟基地台之間的距離， . . . .

(12)

2 利用手機訊號至基地台間的飛行時間乘上速度將可以得到我們估測的距離。而定位系統中，手機以及基地台之間分成兩種狀態分別為:LOS(Line-of-Sight)及 NLOS(Non-Line-of-Sight)，當手機待測物及基地台之間沒有障礙物的情形下，一般而言皆可得到一個夠準的定位結果，不過在都會區或是室內定位系統中， NLOS 為一個常見的現象，使得我們原本的距離量測造成過長的效應，因此直接影響了我們最後定位的結果變差。過去所提出的 NLOS 減緩方法分成兩種，一種是利用 NLOS 的統計特性，結合雜訊以及 NLOS 偏差之機率分布函數[4][19]接著利用最大概似估計法求解，也有利用多次接收手機及基地台間之訊號，並將已知 NLOS 偏差統計特性加以修正量測來估計[20]。另一種則是不需知道 NLOS 統計特性，過去有利用加上幾何限制式的最佳化方法[4]，也有假設各量測為 LOS 先做粗估，接著用估計後的殘餘量作權重[21]，或是利用此殘餘量估計 NLOS 偏差並且加以修正粗估[5]，也有利用線性規劃法並且加上幾何的限制式[6]，由 LOS 量測數據來做定位目標方程式，NLOS 量測數據來建立可行域。本文將著重於 NLOS 發生情況下不須知道統計特性之下的減緩，利用幾種常見的線性化方法[1]，相減法，新增變數法，泰勒線性化法來修改過去所提出的方法，帶有限制式之修正偏差法用不帶限制式取代，可以較短時間的位置估計換取一些不準度的代價，以及改寫現有的相減線性規劃法之目標方程式，利用泰勒以及新增變數法形成，並且提出不同於現有的矩形放鬆，利用泰勒線性化建立可行域，藉此得到更貼近原本可行域的區間，已達到更精準的定位結果。接下來介紹論文的整體架構，第二章為在 LOS 情況下之定位演算法，而著重於線性化方法之介紹，並且提出一個利用 GLS 針對 C(N,2)相減法的權重問題，避開了原本無法求得 covariance 反矩陣的問題。第三章中則為 NLOS 情況下之定位演算法，針對過去提出帶限制式偏差修正法提出修改。第四章則針對現有的線性規劃法，由目標方程式及可行域兩著之建立的角度利用不同的線性化方法來做

(13)

3

修改，並且跟現有的方法做比較以及結合，在第五章電腦模擬中將驗證我們提出

(14)

4

第二章

定位演算法

2.1 量測

無線感知網路定位系統中，為了求得手機的位置，首先我們需要基地台量測手機以及基地台之間的量測資訊，本節將簡介其中兩種量測模型， 1. TOA (Time-of-Arrival) [3]：藉由量測手機發出的訊號到基地台之間的飛行時間，其中飛行時間表示第 i 個基地台的時間量測。因此，手機跟第 i 個基地台之間的量測距離就可以從時間以及訊號傳輸速度的相乘關係式求得。如下式 (2.1) 一般而言，訊號傳輸過程中，TOA 量測會遭受到雜訊的干擾，甚至是障礙物的干擾造成過長的距離量測，此量測模型如下 (2.2)

其中，是手機跟基地台之間的真實距離，是 AWGN (additive white

Gaussian noise) 雜訊表示為



2



~ 0,

j j

n N  ，則是障礙物造成大於零的

NLOS 偏差，而當時，則為 LOS 情況。

2. TDOA (Time-difference-of Arriaval) [3]：

藉由兩個不同基地台所量測到的時間量測差值，可以得到兩個基地台之間的距離量測差，此量測 model 如下

(2.3)

其中，表示第 i 個及第 j 個基地台之間的距離差，TDOA 量測定位演算法[3]

(15)

5

2.2 直視性情況下之 TOA 定位演算法

在無線感測網路定位系統中，當手機與基地台間沒有障礙物的干擾時，我們稱 此情況為 LOS (Line-of-Sight)，每一個已知基地台位置跟未知手機位置之間的真實距離為，如下圖示圖 2.1 LOS 情況之 TOA 量測示意圖假設共有 N 個基地台，訊號經由空氣中傳播後基地台分別可得到 N 組 TOA 量測，回顧 2.1 節中所提 (2.2)，當時 i i i i r  - + nx xi d + n , i = 1,2, ,N(2.4) 其中表示從第 i 個基地台量測的 TOA 得到的量測距離，則是真實距離，



x y



 x 為手機座標位置，x_i 



x_i y_i



為第 i 個基地台座標位置，n 則是高斯_i 雜訊



2



0, _i N  。

(16)

6

2.2.1 最大概似估計法

藉由我們距離量測的雜訊模型中，我們可以將第 i 組基地台所收到的距離量測的機率密度 pdf (probability density function) 表示如下式









2 2 2 1 | e x p 2 2 i i i i i r p r    _{ }          x x x (2.5) 假設每個基地台接收量測時所接受的雜訊干擾皆為獨立，我們可以將概似函數寫成如下式













2 1 2 2 1 1 1 1 , , , | | exp 2 2 N N i i N i _N N i i i i i r p r r r p r        _{ }         







x x x x (2.6) 因此，最大概似估計法想找一組 x 能最大化上式的 pdf，





2 1 1 1 ˆ a r g m a x e x p 2 2 N i i N N i i i i r           _{ }   _ _  _  _    _ _    





x x x x (2.7) 在(2.7)中，最大化等同於我們最小化其中的相加項次，因此我們也可以將之改寫成





2 2 1 ˆ a r g m i n 2 N i i i i r    _{ }     _ _   



 x x x x (2.8) (2.8)中表示了最大概似估計法[1]中的最佳解，也是一種 WLS (Weighted Least Square)解，可以表式成





2 1 ˆ a r g m i nN i i i i r      _   _ 



 x x x x (2.9) 其中，_i 為每一組基地台量測根據統計特性上的訊息我們給予不同的重要性乘積，藉此將較可靠的量測給予較重的權重，而較不可靠的量測則反之，以得到更

(17)

7

精準的解。如果我們不去考量基地台量測的統計特性，則(2.9)式則會變成一個單

純的 NLS (Nonlinear Least Square)，_i皆等於一。

上述的求解方法，皆為非線性的方程式，求解過程需要付出較多的運算複雜度，可以藉由疊代最小方差 [7]方法來求解，但是還是需要較高的運算代價，因此接下來將介紹一些常見的線性化方法，藉此降低我們運算複雜度。

2.2.2 線性最小方差估計法

[1][15] 常見的線性化方法有三種，包含了相減法，新增變數法，還有泰勒線性化法，將在本節內做介紹。

2.2.2.a 相減法

相減法[1][15]是常見的線性化方法之一，首先我們需要將 (2.4) 做左右平方的動作，可以得到 2 2 2 2 i i i i i r  d  d n n (2.10) 上列 (2.10) 中手機及基地台之真實距離d_i 



xx_i

 

2 yy_i



2 ，可以改寫成 2 2 2 2 2 2

2

i i i i i i i i

r



x



y



x



y



xx



yy



d n



n

(2.11) 其中， 2 2 x y 是未知數非線性項次，而 2 2d n_i _in_i 是距離影響誤差及誤差平方項次，為了消除我們的非線性項次，我們找一組基地台 k 的量測資訊當作參考基地台來跟第 i 組基地台量測做相減，如下 2 2 2 2 2 2

2

k k k k k k k k

r



x



y



x



y



xx



yy



d n



n

(2.12) 將(2.11)和(2.12)相減後可得到









2 2 2 2 2x x_k x_i 2y y_k y_i r_i r_k R_k  R_i 2d n_{i i}2d n_k _k n_i n_k , ik (2.13)

(18)

8 其中 2 2 i i i R x y ， 2 2 k k k R x y _。我們先假設 k=1，取第一組基地台當作參考基地台，(2.13)為一組線性化關係式，我們將(2.13)改寫成矩陣型式 D D



D



H x

y

n

_(2.14) 其中



 





 





 



2 2 1 2 1 2 1 2 1 2 2 1 3 1 3 1 3 1 2 2 1 1 1 1

2 ,

,

N N D D N N N N

x

y

r

R

x

y

r

R

x

y

r

R









 







_





_{ }

_































_

 



_





_

_





H

y

D H 矩陣內表示基地台座標差矩陣，y 包含已知量測及基地台座標平方相加減向_D 量，n 表示誤差向量。 _D 一般而言，我們不知道(2.14)誤差向量y ，上述的矩陣方程式則可以寫成 _D



Hx

b

可以利用 WLLS [15](Weighted Linear Least-Square)求手機位置的解如下





1 T T T D D D

ˆ



 _D

x

H W H

H W y

(2.15) 其中 W 是用來做權重的矩陣，是從 n 向量的共變數矩陣的反矩陣求得





1 D D T

E









_

_

W

n n

(2.16)

2.2.2.b 泰勒級數線性化法

泰勒線性化法我們先將(2.4)中的d 寫成 _i



 



2 2

( , )

i i i i

d



F x y



x



x



y



y

(2.17) 利用泰勒級數展開將(2.17)改寫成

(19)

9



,





0

,

0





0

,

0



T i i i

F x y



F x y

 



_

F x y



_

 

高階項次

假設我們的參考點



x y₀, ₀



和真實手機位置夠接近，則高階項次可以被忽略，而 (2.17)可被線性化成下式

 

,



0

,

0





0





0



i i i i

F

F x y

x

y

x

y















(2.18) 其中







 



2 2 0

,

0 0 0 ,0 i i i i

F x y



x



x



y



y



d

0 0 ,0 ,0

,

i i i i i i

F

x

F

y

x

d

y

d



_





_





將上列式子整理可得 0 0 0 0 ,0 0 0 ,0 ,0 ,0 ,0 i i i i i i i i i i

x

y

x

y

x

y

r

d

x

y

d



_





_





_





_





 



















































(2.19) 將(2.19)整理成矩陣型式

H

T S

x



y

T S

+

n

S

T



H

T S



x

y

T S (2.20) 其中 0 1 0 1 1,0 1,0 1 1,0 0 2 0 2 2 2,0 2,0 2,0 ,0 0 0 ,0 ,0 , TS TS N N N N N N x x y y d d r d x x y y r d d d r d x x y y d d                       _ _ _ _       _ _ _         H y 0 0 ,0 ,0 0 0 ,0 ,0 i i i i i i x x y y d d x y d d  _   _   __ __ __ __     可得解為

(20)

10





1 TS TS

ˆ

T T TS TS TS TS TS 



x

H W H

H W y

(2.21) 其中誤差之權重為誤差共變數矩陣之反矩陣 _TS



_TS _TST



1 E    _ _ W n n 。藉由泰勒級數線性化方法，假設我們有夠準的初始值參考點去做泰勒級數展開，可以藉此得到一個夠準的解，並且避開了高運算複雜度的非線性疊代之求解過程。

2.2.2.c 新增變數法

最後一個也是另一種常見的線性化方法為新增變數法[1]，我們將在以下做介 紹整理，首先由(2.17)平方



 



2 2 2 2 i i i i

r



d



x



x



y



y

(2.22) 展開後 2 2 2 2 2

2

i i i i i

r



x



y



x



y



xx



yy

(2.23) 令非線性項次 2 2 Rx y 可改寫整理成 2 2 2

2

i i i i i

r



x



y

 

x x



y y



R

(2.24) 接著將之改寫成矩陣型式

H p

A



y

A

+

n

A



H p

A



y

A (2.25) 其中 2 2 2 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2

2

1

2

1 ,

2

1

A A N N N N N

x

y

r

x

y

x

y

r

x

y

x

y

r

x

y



















_



_











_

_











_



_













H

y

[ ]

x y R

T



p

(21)

11 可得解為





1 A A

ˆ

T T A A A A A 



p

H W H

H W y

(2.26) 其中誤差之權重為誤差共變數矩陣之反矩陣 _A



_A _AT



1 E    _ _ W n n 。從上述式子中，我們可以看到此方法的關鍵在於將非線性項次 2 2 x y 令其成為一組新的變數 R，增加方程式一個自由度，也避開了我們非線性求解的過程。

2.2.3 利用 GLS 作權重 C(N,2)相減法

由 2.2.2a 中，我們介紹了相減法線性化，其中採取了一組基地台當作我們的參考基地台，藉此用其他基地台量測去跟此參考基地台量測相減去除非線性項次，過去有人也利用了不只一組基地台去當作參考基地台來做相減法[8]，在此稱做 C(N,2)相減法，仿造(2.13)中我們將 N 組基地台之間，兩兩相減，其中 i > k，因此能將之整理成以下的矩陣型式

A x



p

(2.27) 其中 1 2 1 2 1 1 2 3 2 3 1 1

2

N N N N N N

x

y

x

y

x

y

x

_

x

y

_

y





















 

_

























A

2 2 2 1 2 1 2 2 2 2 1 1 2 2 3 1 3 1 2 2 1 1

,

N N i i i N N N N

r

R

r

R

x

y

r

R

r

_

R

_



 















_{ }

_













 



























p

(22)

12 藉由線性最小方差估計可以求得 1

ˆ (



T

)

 T

x

A A

A p

_(2.28) 值得一提的是，C(N,2)相減法跟 2.2.2a 的相減法不同在於，我們多加了 C(N,2) – (N-1)組相減關係式，雖然新增出來的這些關係式皆可由原本 N-1 組內線性組合而成，但我們利用多出來的幾組能平均掉參考基地台選取不同造成的效應，可以得到較原本好一點的結果。然而，我們在此也想仿照原本相減法中，利用誤差 Cov 矩陣特性，去做權重，將較可靠的量測組給予較重的權重，反之則給予較小權重，藉此能得到更可靠的解答，因此我們將討論其中誤差擾動向量部分如下

p



p

c



p

n(2.29) 其中 C(N,2) 1 c   p 常數項及 C(N,2) 1 n   p 雜訊擾動向量分別為 2 2 2 1 2 1 2 2 1 1 2 2 3 1 3 1 2 2 1 1 N N c N N N N d d R R d d R R d d R R d d _ R R _           _ _ _                     p 2 2 2 2 1 1 2 1 2 2 1 1 1 2 2 3 3 2 2 3 2 2 2 1 1 1 2 2 2 2 2 2 2 2 N N N n N N N N N N d n d n n n d n d n n n d n d n n n d n d _n _ n n _           _ _ _                     p 由p 中可以看出，每一列的雜訊擾動項次會隨著我們手機以及各基地台之間的_n 距離d 增加時而同時受到影響。 _i 我們試著想要利用上述 (2.29) 中的訊息，藉由 WLS 來求解，寫成下式 1 1 1

ˆ (



T 

)

 T 

x

A C A

A C p

_(2.30)

(23)

13 其中 C(N,2) C(N,2) ( _n) [ _n _nT] Cov E     C p p p 然而，在 (2.30) 中我們發現 C 並非是一個滿秩矩陣，因為其中雜訊擾動向量之間並非每一個元素皆為線性獨立， (2.29)中以四個基地台量測為例 2 2 2 2 1 1 2 1 2 2 3 3 1 1 3 1 2 2 4 4 1 1 4 1 2 2 3 3 2 2 3 2 2 2 4 4 2 2 4 2 2 2 4 4 3 3 4 3 2 2 2 2 2 2 2 2 2 2 2 2 n d n d n n n d n d n n n d n d n n n d n d n n n d n d n n n d n d n n n       _ _ _     _ _ _           _ _ _           p 可以看到p 中的四到六列皆可從一到三列中線性組合得到，因此誤差擾動向量 _n 將其外積後，會得到一組 det=0 的 C 矩陣，因此我們無法求得其反矩陣，也無法由此方法去達到我們想要做的權重動作，我們轉而尋求另一種方法，一般化最小方差法 (Generalize Least-Square)[9]，想要藉由 GLS 方法去避開我們需要利用 C 的反矩陣才能得到權重效果，此方法的求解過程介紹如下由(2.27)中可改寫成

A x + p

n

= p

(2.31) 其中，p 的共變數矩陣為一個對稱並且帶有正定特性，我們可以將之分解成_n T  C BB ，其中B C(N,2) C(N,2) 然而為了避開 C 反矩陣才能做權重的動作，利用[9] 所提出的一般化最小方差法形式求解將可以改寫如下

min

T Ax+Bv=p

v v

(2.32) C(N,2) 1  v 以下將簡略介紹 GLS 法的求解步驟首先我們先將  C(N,2) 2 A 矩陣做 QR 分解

(24)

14 1 2 2 C(N,2)-2

,

[

]

T



 



 

 

R

Q A

Q

0

(2.33) 其中  C(N,2) C(N,2) Q 為正交矩陣，  2 2 R 。 因此，(2.32)中的 Ax + Bv = p 可以寫成

(

)

T



T

Q Ax + Bv

Q p

_(2.33) 也等同於 1 1 1 2 2 2 T T T T T T









 



_

_



_

_

 

 









R

Q B Q B

Q p

x

v

0 Q B Q B

Q p

(2.34) 接下來我們想找一組  C(N,2) C(N,2) Z 使得 2







1



2 T



Q BZ

S

0 Z

Z

(2.35) 其中S C(N,2) C(N,2)-2 為一組下三角矩陣，Z 為正交矩陣， C(N,2) C(N,2)-2  1   Z ， C(N,2) 2 2   Z 。再將(2.34)改寫成 1 1 1 1 1 2 2 2 T T T T T T



 







 

_

_



 







 

 



 



R

Q BZ

Z v

Q p

x

0 S

0 Z v

Q p

(2.36) 我們可以由上述(2.36)注意到，求解方程式中的”下半部”只跟 v 有關係，因此 我們可以藉此求 v，如下 1 2 T



T

SZ v

Q p

(2.37) 然而，因為上列(2.37)中未知數大於我們的觀察向量，因此我們利用僞逆反矩陣

方法求出一組最小範數解(minimum norm solution)的 v

† 1 2

ˆ



(

T

)

T

v

SZ

Q b

(2.38)

(25)

15 有了(2.38)中的 ˆv 後，我們將之代入(2.36)中的”上半部”，藉此去求得手機位置， 整理過後如下 1 1

ˆ

T T





Rx

Q p Q Bv

(2.39) 而手機最終估計位置則可從上式中求得





† 1 1

ˆ



T



T

ˆ

x

R

Q p Q Bv

(2.40) 此方法值得一提的部分在於，利用改寫成一般化最小方差的過程，我們不僅僅 避開了原本由於無法求得 C 的反矩陣而無法藉此去將較可靠的量測數據以及較 不可靠的分別做較大權重以及較小權重的動作，更因此可以求得更精準的解。本章節主要專注在於在 TOA 量測模型下，直視性(LOS)情況下的定位演算法介紹，簡介了三種線性化最小方差，更進一步再相減法線性化之 C(N,2)的權重問題，利用一般化最小方差法避開了原本共變數矩陣 C 不是滿秩而無法得到反矩陣來作權重。

(26)

16

第三章

NLOS 情況下利用偏差修正法之定位

不同於第二章所介紹的 LOS 情況下定位演算法，本章節將針對手機以及基地台之間遭受 NLOS (Non-Line-of-Sight)情況下之定位演算法，在都會區或是室內空間中，手機及基地台之間的距離量測，容易遭受障礙物干擾，造成過長的距離估計，這些過長的距離量測使得我們手機估計位置將造成更大的誤差，NLOS 情況示意圖如下圖 3.1 手機及基地台間 NLOS 發生情況示意圖

3.1 權重非線性最小方差

本節將討論 NLOS 偏差將之當作未知變數[4]，寫成權重非線性最小方差如下 2 , 1

ˆ = a r g m i n

(

( )

)

N i i i i i

r

d

b









x b

x

(3.1) (3.1)中相加項次表示著我們對應每一組基地台量測，位置估計的殘餘量，_i為權

(27)

17 重係數，想找出一組手機位置以及分別對應上的 NLOS 偏差使得我們的殘餘量達到最小，因為b 有 N 組，因此求解過程中，未知數大過我們觀察向量，過去_i 的文章[4]提出了增加了 NLOS 限制式，給予每組 NLOS 偏差一個範圍區間限制，藉此縮小我們可行域，如下

b

L

 

b

U (3.2) 其中，NLOS 偏差上界 b 由手機以及基地台之間接收到的距離量測之幾何關係_U [10]可得到





, min , , j=1,2, ,N, j i i=1,2, ,N U i i j i j b  r  r d  (3.3) 其中，d_{i j}_, 為第 i 組及第 j 組基地台之間的真實距離，而 NLOS 偏差的下界由於 NLOS 為大於零，因此我們令b_L0 。 由上述(3.1)當作目標方程式對應(3.2)的限制式所形成的可行域，我們會得到一組非線性帶有限制式的方程式 2 , 1

ˆ = arg min

(

( )

)

L U N i i i i i subject to

r

d

b



  





x b b b b

x

(3.4)

可藉由 Sequential Quadratic Programming [11]來解這組帶限制式之最佳化問題。

3.2 帶限制式之偏差修正法[5]

由 3.1 中我們可以看到，(3.1)WNLS 式子中將手機位置估計以及 NLOS 偏差估計當成未知數，並且藉由 NLOS 偏差的限制式，去將各基地台量測的殘餘量相加達到最小，本節將介紹一種方法分成三步驟，首先假設各基地台皆為 LOS 並且利用泰勒線性化之 LS 粗估求得的解，用殘餘量 l 概念，利用內點法求 NLOS 偏差效應，最後將估出的 NLOS 偏差修正粗估的手機位置，有別於上一節所介紹之 NLOS 減緩方法，以下將做簡單介紹。

(28)

18 首先，遭受 NLOS 干擾的量測模型可由 (2.2) 寫成以下矩陣形式

r



d x

( )

 

n b

(3.5) 其中 r 表示距離量測向量，d(x)表示手機及基地台真實距離之位置函數，b 則為 我們 NLOS 偏差向量。由 (3.5) 中形成的非線性目標方程式將之改寫成二次式矩陣形式如下

J

  

x



r d x



( )



b



T

R

1



r d x



( )



b



(3.6) 其中R為雜訊之共變數矩陣。由 (3.6) 中我們將此非線性之方程式利用 2.2.2.b 所介紹之泰勒線性化法，並加以改寫可以得到以下之矩陣形式

y

T S



H

T S

x

 

n

(3.7)

b

其中H 及_TS y 從 2.2.2.b 中之 (2.20) 可得。 _TS 而由上述 (3.7) 中之泰勒線性化後之量測模型我們將之改寫成

J

  

T S T S



T

R

1



T S T S









x

y

H x b

y

H x b

(3.8) 為了最小化我們 (3.8) 的目標方程式，將之對 x 取一階微分並令之為零，可以得到我們 LS 解如下



₁



1 ₁



₁



1 ₁

ˆ

T T T T TS

R

TS TS

R

TS TS

R

TS TS

R

     





x

H

y

H

b

_(3.9) 上式 (3.9) 中如果我們先忽略了 NLOS 偏差項次 b，則會是一組較不準的手機位 置估測





1 1 1 T T TS

R

TS TS

R

TS   



x

H

y

_(3.10) 其中x為假設各基地台皆為 LOS，是一組假設 NLOS 干擾不存在的位置估計。 接著加以考慮 NLOS 偏差 b，(3.9)可以改寫成

x

ˆ

 

x Vb

(3.11) 其中 V 為一組對應 b 之修正矩陣，由 (3.9) 中可得到

(29)

19





1 1 1 T T TS

R

TS TS

R

  

 

V

H

(3.12) 實際情況下，b 是一組未知向量無法從距離量測中直接得知，但可以被估計，首 先，利用第一步驟粗估的x其殘餘量如下

z



y

_{T S}



H x

_{T S} (3.13) 將 (3.7) 及 (3.10) 代入上式 (3.13) 中可以得到

z





H

T S

x

 

n



b



T S

H



ˆ

x





(3.14)

V b

因為我們想利用此殘餘量關係去找出 b 的未知估計向量，因此將(3.14)整理後可 以得到

z





I



H

T S



V b



H

T S





ˆ

x



x



(3.15)

n

其中 I 為單位矩陣。 令S I H V ，誤差向量_TS wHTS



x x ˆ



n ，寫成

z



Sb



w

(3.16) (3.16)中，我們可以看到 z 以及 S 皆可求得，w 之共變數矩陣也可以求如下列敘 述

Q

w



E





ww

T







H P H

TS x TST



R

(3.17) 其中









₁



1

ˆ

T _T TS TS

E





R

 



_



_



x

P

x x x x

H

有了上述的 (3.16) 式，加上 (3.2) 式中 b 的限制式，增加 b 估計的可行域資訊， 可以將之寫成









1

ˆ

_{arg min}

L U T subject to   





_w



b b b b

b

z Sb

Q

z Sb

_(3.18) 由上式 (3.18) 求出 ˆb 後代入 (3.11) 中，修正原本的粗估計x

x

ˆ

 

x Vb

ˆ

(3.19)

(30)

20 最後我們需要將上述方法遞回，直到最新一次 ˆk x 跟上一次的 _ˆk1 x 誤差小於一個門檻，當作我們最終估計，反之則將 ˆk x 再次代入 (3.5) 當作泰勒線性化參考點， 並再次做估計。流程圖如下分成三個區塊，第一個區塊為假設無 NLOS 之粗估，接著是 NLOS 偏差估計，最後則是藉由第二區塊的偏差估計，來修正我們第一區塊的粗估解。

Yes

ˆ

k



k



k k

x

V b

 

T k T k k k TS TS

E





R



_

_





w x

Q

ww

H P

H

 





1 1 T k k k TS

R

TS  



x

P

H

 





1

 

1 1 T T k k k TS

R

TS TS

R

  

 

V

H

k k k TS

 

S

I

H V

k k TS



TS



y

H x n

k



k



k

z

S b w

( )



 

r

d x

n b

1 1

ˆ

(

)

k k k k TS

d

TS  

 



y

r

x

H x

1 1

ˆ

(

)

k k k k TS

d

TS  

 



y



 

r

x



H x

 

1 1 1 T T k k k k k TS

R

TS TS

R

TS   



x

H

y

k k k k TS TS





z

y

H x





1





ˆ

_{arg min}

L U T k k k k k k subject to   





_w



b b b b

b

z

S b

Q

z

S b

1

ˆ

k

x

1

ˆ

k



ˆ

k





x

No

(31)

21

3.3 修改 NLOS 偏差修正法

然而，上一節所介紹的方法需要有著較精準的泰勒初始值參考點，才能使得我 們遞回估計較快達到夠準的最終解，並且每次迴圈在估計 b 時，都需要求一組帶 限制式之最佳化問題，需要較高的運算複雜度，因此我們為了減少運算上的複雜度，我們將上述方法稍作修改。首先，在第一區塊中，我們利用 2.2.2.c 所介紹過的新增變數法替代泰勒線性化法以求假設無 NLOS 之粗估，如

y

A



H p n

A



(3.20) 其中 2 2 2 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 1 2 2 1 , , [ ] 2 2 1 T A A N N N N N x y r x y x y r x y x y R x y r x y           _ _  _ _              _ _          H y p 接著在第二區塊中，則擷取了原本利用泰勒線性化法求得此粗估後，利用殘餘量概念得到的 NLOS 偏差估計關係式，把我們第一區塊所求得的粗估解代入，並且將 NLOS 偏差估計的限制式拿掉，避開了需要帶限制式之最佳化問題，我們 NLOS 偏差估計的目標方程式如下

ˆ a r g m i n





1





T _





_w



b

z

S b

Q

z

S b

_(3.21) 其中和 (3.18) 中不同就在於，我們少了限制式。 其中值得一提的是，藉由 S 和H 相乘如下 _A







₁



1 ₁

0

TS TS TS T T TS TS TS

R

TS TS

R

TS   

 







S H

I

H V H

H

(3.22)

上式(3.22)中，由 Sylvester’s inequality theorem[12]可知道 det(S)=0，因此 S 並非 一組滿秩矩陣，我們尋求最小範數解(minimum norm solution)來求 NLOS 偏差，

(32)

22 利用權重偽逆反矩陣可以得到下式

b

ˆ



S z

† (3.23) 其中 1 1 2

_,

2  



_w



_w

S

Q S

z

Q z

接著則是利用我們所求出的 NLOS 偏差估計，進一步去修正我們的距離量測 r， 如下

r

 

r b

ˆ

(3.24) 接下來再回到第一步驟將此修正過後的距離量測利用新增變數法求出下一次的粗估位置估計進入下一次修正迴圈，直到跟上一次的估計誤差小於一組設定門檻。兩種方法之間的方塊圖比較如下 1. 帶限制式之 NLOS 偏差估計修正流程圖 圖 3.2 帶限制式之 NLOS 偏差估計修正流程圖

ˆ

k

x

1

ˆ

k



ˆ

k





x

k

x

b

ˆ

k

TS linearize

Bias estimation

With constrain

correction

1

ˆ

k

x

Yes No

(33)

23 2. 修改 NLOS 偏差修正法之流程圖 圖 3.3 修改 NLOS 偏差修正法之流程圖由上述兩個方塊圖我們最後整理出兩種方法差別在於，第一種利用泰勒線性化須要先有一組展開參考點，而是新增變數法不需要的，而第二塊區塊裡面則是將原本泰勒線性化法所求粗估解形成的殘餘量，利用帶限制式最佳化法來求解，而我們的方法不同於，少了限制式並且直接利用最小範數解求出 NLOS 偏差估計，最後，不同於第一個方法利用 NLOS 偏差估計去修正我們粗估解，而是用來修正我們的距離量測藉以減緩 NLOS 效應，使其更接近 LOS 的距離量測，然後在將此修正過後的量測再次利用新增變數法求解，以下為此方法完整步驟流程

ˆ

k

b

ˆ

k

x

r

N o 1

ˆ

k



ˆ

k





x

Yes

新增變數

Bias estimation weighted Pseudo-inv 1

ˆ

k



k



k

r

b

新增變數

0

ˆx

1

ˆ

k

x

(34)

24

Yes

No

1

ˆ

k



ˆ

k





x





1





ˆ

_{arg min}

T k



k



k k  k



k w b

b

z

S b

Q

z

S b

 

T k T k k k TS TS

E





R



_

_





w x

Q

ww

H P

H

 





1 1 T k k k TS

R

TS  



x

P

H

 





1

 

1 1 T T k k k TS

R

TS TS

R

  

 

V

H

k k k TS

 

S

I

H V

k k k TS TS





z

y

H x

0 A



A



y

H p n

k



k



k

z

S b w

( )



 

r

d x

n b

1 1

ˆ

(

)

k k k k TS

d

TS  

 



y

r

x

H x

1 1

ˆ

(

)

k k k k TS

d

TS  

 



y

r

x

H x





1 0 1 1 0

ˆ

T T A

R

A A

R

A   



x

H

y

k A



A



y

H p n





1 1 1

ˆ

k T T k A

R

A A

R

A   



x

H

y

1

ˆ

k k k k A 







r

b

y

(35)

25

第四章

NLOS 情況下利用修改線性規劃法之定

位

本章節將討論 NLOS 效應減緩的其他方法，首先，假設手機對應於每個基地

台之間 LOS 及 NLOS 狀態能先判斷[17][18](identification)，前一章所討論的方法不需要事先知道 LOS 及 NLOS 狀態，而是利用假設各基地台及手機間皆 LOS 的

概念粗估一個手機位置，接著估計 NLOS 偏差藉以修正粗估的手機位置，本章節則是利用已知 LOS 狀態的距離量測建立我們的線性目標方程式，而 NLOS 狀態的距離量測則是用來建立線性限制式來形成一個可行域限制我們搜索區間，在 4.1 中將先介紹其中一種線性規畫法以減緩 NLOS 效應之方法，4.2 及 4.3 將針對 4.1 的方法基礎，由目標方程式及可行域的角度分別作修改。

4.1 現有線性規劃法介紹

本節將介紹利用目前被提出過利用線性規劃法來減緩 NLOS 效應之方法[6]，首先，假設手機以及各基地台之間 LOS 以及 NLOS 狀態能被正確辨識，則由(2.2) 中 TOA 的模型可以寫成

LOS

:

r

Lj



d ( )

Lj

x



n

j

,

j



1,

,

m

L (4.1a)

NLOS

:

r

Nj



d ( )

Nj

x

 

n

j

b

j

,

j



1,

,

m

N (4.1b)

其中，下標 L 及 N 分別代表 LOS 及 NLOS 之量測，m 及_L m 分別表示 LOS 及_N

NLOS 量測個數，雜訊 n 為 AWGN。

首先考慮 LOS 基地台量測，由 (4.1a) 中我們假設一共有m 個 LOS 量測，利_L

(36)

26 2 2 ,1 ,2 ,1 ,2 ,2 ,1 ,2 ,1 2 2 ,1 , ,1 , , ,1 , ,1 2 2 ,2 ,3 ,2 ,3 ,3 ,2 ,3 ,2 2 2 ,( 1) , ,( 1) , ₂ , ,( 1)

1 ,

2

L L L L L L L L L L N L L N L N L L N L L L L L L L L L L N L N L N L N _M L N L N

x

y

r

R

x

y

r

R

x

y

r

R

x

_

x

y

_

y

r

_ 































_































A

p

, ,( 1) ₁

,

L N L N _M

R

_ 











































2 2 , , , L i L i L i

R



x



y

_(4.2) 將上述矩陣形式中令 M = C(N,2)，則此 M 組方程式中，寫成通式如下

a x

i j



b y

i j



, i ,

c

i j



j

1 ,

m

L

,

(4.3)

i < j

其中 , ,

,

, , ij L i L j ij L i L j

a



x



x

b



y



y

2 2 , , , ,

2

L j L i L j L i ij

r

R

c







接著由上述 (4.3) 中我們定義殘餘量為

e

i j



a x

i j



b y

i j



, i ,

c

i j



j

1 ,

m

,

L (4.4)

i < j

而手機位置最終估計的最小化目標方程式就是將所有 (4.4) 中殘餘量平方和，如下 2 ,

ˆ a r g m i n

_ij i j j i

e





 

x

_(4.5) 在這裡我們不將上述 (4.5) 中之位置估計利用 (2.28) 中的 LLS 來求解，由(4.5) 中，我們將變數

e

ij視為一組不帶限制式的鬆弛變數[13]，試著將上式中由

e

ij的非線性目標方程式轉換成線性規劃式，尋求線性規劃之解法，首先我們定義 , ij i j j i

Z

e



 

(4.6)

(37)

27 接著我們將令

e

i j

e

i j

e

,

i j

e

i j

e

,

i j

0

   







_(4.7) 並將 (4.6) 中改寫成





, i j i j i j j i

Z

e

 

e





 

(4.8) 在此特別說明，由[13]中可以知道，上式 (4.8) 中的最佳解使得我們目標方程式 Z 能達到最小的時候，



e e

ij

,

ij



  將有一個是等於

e

ij而另一個則為零。接著，由 (4.4) 中可以將限制式寫成

,

1, 2,

,

;

ij ij ij ij ij L

a x b y e



 



e





c

i j



m i



j

_(4.9) 由上式 (4.9) 中可觀察到我們一共有 2M 個非負的鬆弛變數以及[x y]一共(2M+2) 個未知數，定義 T 2 M + 2

基於修改線性規劃法及偏差修正法之無線感測網路定位

國立

交

通

大

學

電信工程研究所

碩

士

論

文