基于無線傳感器網絡時鐘同步的定位算法研究

摘 要： 討論了在無線傳感器定位節點（文中稱為錨節點）位置已知條件下，基于信號到達時間差（TDOA）的時間同步算法。根據所有傳感器定位節點彼此之間周期性地廣播信號，測量信號到達時間（TOA）來估計時間頻率的偏移。當信源發送信號，由于傳感器鐘差造成的初始TDOA誤差可以通過一系列補償計算來減少誤差。在一個加性高斯白噪聲模型里，通過卡爾曼濾波提高鐘差和頻差的估計精度。在此基礎上改善標簽節點的定位精度。

關鍵詞： 無線傳感器網絡； 卡爾曼濾波； 定位； 到達時間差； 到達時間； 鐘差估計

中圖分類號： TN911?34 文獻標識碼： A 文章編號： 1004?373X（2015）07?0023?05

0 引 言

隨著無線電技術的發展，人們在利用信息上已經不再滿足于相互之間的通信和聯絡，對目標進行有效定位和跟蹤，獲取有用信息加以利用，成為人們更進一步的技術追求和現實需要。伴隨著物聯網的發展，無線傳感器網絡逐漸成為其概念應用的重要技術支撐。傳感器節點必須明確自身位置才能詳細說明“在什么位置或區域發生了特定事件”，實現對外部目標的定位和追蹤[1]，從而為人們獲知事態，做出相應決策提供依據。諸多新近的定位技術，較為準確的定位方法是基于傳感器定位節點的測距方法。這些方法通常使用TOA、TOF或者TDOA進行定位。在基于TOA/TDOA的無線傳感器網絡定位中，都會根據各自系統需要獲得相對準確的時間同步從而進行定位。

在諸多文獻中無線傳感器的時間同步方法有很多。其中引用較多的方法是參考廣播同步算法（Reference Broadcast Synchronization）[2]、傳感器時間同步協議（TPSN）[3]、泛洪時間同步協議（FTSP）[4]。另外，Noh等人提出基于雙向消息交換數學評估時間同步方法的性能[5]。他們提到的方法被用來解決不同系統中存在的問題。但是，兩個連續不斷的TOA測量估計時間同步仍沒有達到滿意的效果。

本文使用一個補償算法來補償鐘差以達到提高定位精度的目的。無線傳感器網絡中的錨節點彼此之間周期性地廣播信號或者數據包。每個錨節點利用自身接收器的本地時鐘將信號的TOA信息記錄下來。骨干網絡的服務器收集所有的時間戳來估計計算每對錨節點之間的頻差。采用卡爾曼濾波進一步提高鐘差和頻差的精確性。

1 模型建立

TDOA定位技術是通過檢測信號到達兩個錨節點的時間差來確定標簽節點（即需定位節點）的位置[6]。當標簽節點向多個錨節點發送信號時，由于錨節點距離標簽節點不一樣，不同的錨節點接收信號有一個時間差。而信號波在同一個介質中的傳輸速率相同，因此可以計算出標簽節點到達兩個錨節點的距離差[7]，由雙曲線計算原理得出定位位置。

假設一個有[N]個錨節點的無線傳感器網絡，其錨節點位置為：[xi，]其中[i=1，2，…，N，]一個標簽節點的位置設為[x0。][N]個錨節點的內部時鐘的標稱頻率為[f0，]未知時鐘偏移為[εi，]每個錨節點有一個未知的起始時間[Ωi。]同時將這種錨節點周期性的廣播信號稱為“眨眼”信號。

第[i]個錨節點的時鐘實際頻率可以表示為：

[fi=f0（1+εi）] （1）

假設當前發送廣播包為錨節點[p，]在錨節點[j]處測量的TOA時間可以表示為：

[tj，p（k）=kTs+djp1+εj1+εp+Ωj+nj，p（k）] （2）

式中：[djp=xj-xpc]為兩個錨節點之間的傳播時延；[Ωj，j=1，2，…，N]為各個錨節點的起始時間；[Ts]為前后兩個廣播包發送的時間間隔；[nj，p（k）]為TOA的測量噪聲，可以看作加性高斯白噪聲（AWGN）。經過泰勒級數展開有：

[1+εj1+εp≈1+εj-εp] （3）

TOA可以近似表示為：

[tj，p（k）=kTs+djp1+εj-εp+Ωj+nj，p（k）] （4）

令[ejp=εj-εp]為兩個錨節點之間的頻差，由于晶振制作工藝和水平等原因其數值通常介于點典型值[10-6]和[10-4]之間。對于無線傳感器網絡，錨節點之間的距離通常在300 m以內，[djpejp]則為一個遠小于[10-10]的時間數值，對于定位結果的影響不大，這里可以略去。本地測量的TOA可以進一步簡化為：

[tj，p（k）≈kTs+kTsejp+djp+Ωj+nj，p（k）] （5）

則標簽節點到錨節點[i]和錨節點[j]之間的TDOA為：

[ti，p（k）-tj，p（k）=kTs+kTseip+dip+Ωi+ni，p（k）-kTs+kTsejp+djp+Ωj+nj，p（k）=kTseij+dip-djp+Ωi-Ωj+ni（k）-nj（k） ] （6）

2 估計頻差和鐘差

通過對連續兩次錨節點廣播獲得的TOA進行分析計算，可獲得錨節點鐘差和頻差的估計值。每個錨節點按照時間間隔[T]發送信號，其他錨節點根據他們各自的時鐘接收信號并記錄TOA時間。連續兩次“眨眼”信號的TDOA測量過程如圖1所示。

圖1 連續兩次“眨眼”測量的TOA信息

[ti，p（k+1）-tj，p（k+1）-（dip-djp）=ti，q（k）-tj，q（k）-（diq-djq）+T?eij（k）+ni（k）-nj（k）] （7）

令[xij（k）]為錨節點[q]發送信號測量的鐘差：

[xij（k）=ti，q（k）-tj，q（k）-diq-djq] （8）

[xij（k+1）=xij（k）+T?eij（k）+ni（k）-nj（k）] （9）

式中：[eij（k）]為錨節點之間存在的頻差；[ni（k）]為測量噪聲，可以看作高斯白噪聲。錨節點之間的鐘差實際上是由頻差隨時間累積造成的，在實際系統中錨節點之間的頻差可能會隨著時間而發生偏移，因此本文采用卡爾曼濾波對鐘差和頻差進行濾波估計。

由于[eij（k）]是晶振的相對頻率偏移，這里可以把隨機頻率偏移看成是一個方差為[σ2f]的加性高斯白噪聲。

錨節點[i]和[j]之間的頻偏可表示為：

[eij（k+1）=eij（k）+nf，i（k）-nf，j（k）] （10）

式中：[nf，i（k）]為錨節點晶振的隨機頻率擾動，可以看作AWGN。讓[?ij（k）]為錨節點[q]發送信號的本地真實鐘差，真實鐘差可表示為：

[?ij（k+1）=?ij（k）+T?eij（k）] （11）

假設[sil（k）]為狀態向量，狀態轉移矩陣可以用[Φ]來表示：

[sil（k）=?il（k）eil（k），Φil=1T01，]

[nil（k）=T?nf，i（k）-T?nf，l（k）nf，i（k）-nf，l（k）]

得到時鐘偏移狀態轉移方程為：

[si1（k+1）=Φi1si1（k）+ni1（k）， i=2，3，…，N] （12）

將[N-1]個方程合并：

[s（k+1）=Φs（k）+Γv（k）] （13）

其中：

[s（k）=s21（k）s31（k）?sN1（k）2（N-1）×1]

[Φ=Φ21Φ31?ΦN1]

[Γ=-TT0…0-110…?-T0T…0-101…0?????-T00…T-100…12（N-1）×N]

[v（k）=nf，1（k）nf，2（k）…nf，N（k）T]

其中[v（k）～Ν（0，Qv），][Qv]是頻率波動向量的協方差矩陣，結合式（9）和式（11），鐘差觀測方程為：

[xi1（k）=?i1（k）+ni（k）-n1（k），i=2，3，…，N] （14）

設：

[H=1000010?????????10（N-1）×2（N-1）]

[u（k）=n1（k）n2（k）…nN（k）T]

[C=-110…0-101??????0-10…01（N-1）×N]

鐘差觀測方程可以轉變為矩陣形式：

[y（k）=Hs（k）+Cu（k）] （15）

其中[u（k）～Ν（0，Qu）。]

卡爾曼濾波過程如下：

（1） 預測

[s（k+1|k）=Φs（k|k）] （16）

（2） 預測誤差方差陣

[Ps（k+1|k）=Φ?Ps（k|k）?ΦT+ΓQv（k）ΓT] （17）

（3） 卡爾曼增益

[K（k+1）=Ps（k+1|k）HT[H?Ps（k+1|k）?HT+CQu（k+1）CT]-1] （18）

（4） 濾波

[s（k+1|k+1）=s（k+1|k）+K（k+1）[y（k+1）-H?s（k+1|k）]] （19）

（5） 濾波誤差方差陣

[Ps（k+1|k+1）=[I-K（k+1）?H]?Ps（k+1|k）] （20）

3 標簽定位

3.1 TDOA測量補償

由于標簽根據各自周期發射信號，一些類似于CSMA/CA這樣的介質訪問方法，可以用于避免錨節點“眨眼”和標簽發送信號的沖突。他們發送信號的機制類似圖2所示。

圖2 標簽和錨節點測量的TOA信息

標簽節點到達兩個錨節點的TDOA為：

[ψij（kT+τ0）=ti，0（kT+τ0）-tj，0（kT+τ0）=?ij（kT+τ0）+τ0eij（k）+di，0-dj，0+ni，0（kT+τ0）-nj，0（kT+τ0）] （21）

其中[τ0]是最后錨節點“眨眼”的耗費時間。減去最初TDOA，則得到獲得補償的TDOA：

[zij（kT+τ0）=ψij（kT+τ0）-?ij（k）-τ0eij（k）] （22）

[zij（kT+τ0）=?ij（kT+τ0）-?ij（k）-τ0eij（k）+di，0-d1，0+ni，0（kT+τ0）-n1，0（kT+τ0）] （23）

對于所有[N]個錨節點，有[N-1]個獨立的TDOA補償變量[zi1（kT+τ0），i=2，3，…，N：]

[zi1（kT+τ0）=1τ0??i1（k）-?i1（k）ei1（k）-ei1（k）+-τ0τ0?nf，1（kT+τ0）nf，i（kT+τ0）+di，0-d1，0+-11?n1，0（kT+τ0）ni，0（kT+τ0）] （24）

則標簽節點的TDOA為：

[z=f（x0）+w] （25）

其中：

[z=z21（kT+τ0），…，zN1（kT+τ0）T] （26）

[f（x0）=1cx0-x2-x0-x1??x0-xi-x0-x1?x0-xN-x0-x1] （27）

[w=Dse（k）+Cτv0（kT+τ0）+Cu0（kT+τ0）] （28）

[se（k）=s（k）-s（k）] （29）

[D=1τ00……0001τ00…?……………?00…01τ0（N-1）×2（N-1）] （30）

[Cτ=-τ0τ00…0-τ00τ0??????0-τ00…0τ0（N-1）×N] （31）

由于[se（k）]和[u0（kT+τ0）]是各自獨立的，所以[w]的協方差矩陣為：

[Qw=E（wwT）=CQuCT+CτQvCTτ+DM（kk）DT] （32）

3.2 最大似然估計

似然方程為：

[p（z，x0）=1（2π）MN/2det（Qw）12?exp-12z-f（x0）TQ-1wz-f（x0）] （33）

[x0]的最大似然估計為：

[x0=argminz-f（x0）TQ-1wz-f（x0）] （34）

假定在高斯白噪聲下，最大似然估計有一個最小二乘估算，但是因為函數的非線性，總體而言公式（34）閉合形式的解并不存在，因此需要數值極小化。在這里，使用一個連續的線性算法[8]，概述如下：

（1） 設[x0（m）]為估計量經過[m]次迭代的結果，并有[x0=x0（m）+Δ（m），]對[f（x0）]在[x0（m）]附近進行線性化，得到：

[f（x0）≈f（x0（m））+G（x0（m））Δ（m）] （35）

其中[G（x0（m））]為Jacobian矩陣：

[G（x0）=?f（x0）?x0] （36）

將式（35）代入式（34），可以解決最小線性化問題：

[Δ（m）=GT（x0（m））Q-1wG（x0（m））-1?GT（x0（m））Q-1wz-f（x0（m））] （37）

其中：

[G（x0）=1crT2（x0）-rT1（x0），…，rTN（x0）-rT1（x0）T] （38）

[ri（x0）]定義了一個單位長度的矢量：

[ri（x0）=x0-xix0-xi，i=1，2，…，N] （39）

（2） [（m+1）]次迭代估計為：

[x0（m+1）=x0（m）+Δ（m）] （40）

注意到當[GT（x0（m））Q-1wG（x0（m））-1]很大的時候，迭代的方法會停在某個局部最小值并且不會收斂。兩步過程，即從一個粗柵格搜索并持續一個迭代的過程可以被用來尋找全局最小值。通過公式（40）不斷進行迭代的過程，如果前后兩次迭代的結果差距大于0.1 mm，則繼續進行迭代過程，直到二者的差距小于0.1 mm時，迭代結束。當前所得到的[x0]即為標簽定位的結果。

標簽位置估計的CRLB是Fisher信息矩陣（FIM）的轉置：

[J（x0）=E??x0lnp（s；x0）??x0lnp（s；x0）T] （41）

因為高斯噪聲，公式（41）變為：

[J（x0）=?f（x0）?x0TQ-1w?f（x0）?x0=GT（x0）Q-1wG（x0）] （42）

假定[J-1（x0）]存在，[x0]任意元素的方差為[[J-1（x0）]nn，]標簽位置估計的CRLB為：

[σ2x，t=n=12[J-1（x0）]nn] （43）

4 仿真結果

在本節，通過仿真對時鐘頻率偏移的補償算法進行測試。定義一個二維空間，四個錨節點的位置如圖3所示。錨節點記下其他錨節點或者標簽發送來的TOA數據包的時戳，并且將這些信息發送給定位服務器，由它進行時鐘偏移估計計算。

圖3 錨節點位置設置

4.1 時鐘差異估計結果

通過蒙特卡洛仿真實現文中提及的頻率偏移估計算法的評估，如圖4所示。四個錨節點的頻率各自顯示出10，15，20，30的周期性正弦變化。[{σe（1，j），j=2，3，4}]是三個各自獨立的[{e1j，j=2，3，4}]頻率偏移平方根方差。均方根誤差結果是平均超過10 000個獨立噪聲實現。

通過仿真，隨著數據量的增加，隨著之前數據的代入，誤差估計效果將進一步提高。

4.2 標簽定位結果

在圖5中，隨著不同的[T]周期均方根誤差和CRLB顯示了標簽定位結果，首先能得出如果[T]值取得足夠小，RMSE和CRLB就會減小，這是由于時鐘漂移的變化遠小于“眨眼”周期。雖然在“眨眼”周期[T]很小時，能得到更加精確的標簽定位估計，但是由于隨之而來的錨節點“眨眼”產生的大量的信息和數據，從考慮系統穩定角度需要限制定位標簽的數量。所以一個適當的[T]值是必要的。另一方面，當TOA測量標準偏差在[σf]之下，RMSE會遠大于CRLB。這是因為TOA測量標準差[σn]小于[σf]時，由于隨機的時鐘波動而引起的鐘差占主導因素。

圖4 對TOA測量時鐘頻率偏移RMSE估計，

T=0.1，0.2，0.4，0.8

圖5 在不同錨節點眨眼間隔下，標簽位置估計

的CRLB與RMSE比較

5 結 語

本文提出了一個新的無線傳感器TDOA時間同步補償算法。所有傳感器或者錨節點都裝備有相同標稱頻率的石英時鐘。實際上錨節點彼此異步工作，錨節點間的頻率偏移通過錨節點“眨眼”估計。卡爾曼濾波被用來提高估計性能并減少計算量。這種方法簡單并且極大方便異步錨節點的無線傳感器網絡。

參考文獻

[1] VARESE C， RABAEY J M， BEUTEL J. Locationing in distri?buted ad?hoc wireless sensor network [C]// Proceedings of the 2001 IEEE International Conference on Acoustics， Speech， and Signal. Salt Lake： IEEE Signal Processing Society， 2001， 4： 2037?2040.

[2] ELSON J， GIROD L， ESTRIN D. Fine?grained network time synchronization using reference broadcasts [C]// Proceedings of 5th Symposium on Operating System Design and Implementation， Boston， MA： [s.n.]， 2002： 147?163.

[3] GANERIWAL S， KUMAR R， SRIVASTAVA M B. Timing?syn protocol for sensor networks [C]// Proceedings of SenSys. Los Angeles， CA： [s.n.]， 2003： 138?149.

[4] MAROTI M， KUSY B， SIMON G， et al. The ?ooding time synchronization protocol [C]// Proceedings of SenSys. Baltimore， MD： [s.n.]， 2004： 39?49.

[5] NOH K， CHAUDHARI Q M， SERPEDIN E， et al. Novel clock phase offset and skew estimation using two?way timing message exchanges for wireless sensor networks [J]. IEEE Transactions on Communications， 2007， 55（4）： 766?777.

[6] 許秀蘭.無線傳感器網絡距離無關定位算法研究[J].計算機與現代化，2009（10）：6?9.

[7] 張健，李鷗.改進的無線傳感器網絡TDOA定位算法[J].計算機工程與應用，2010，46（16）：117?120.

[8] KAY S M. Fundamentals of statistical signal processing， Volume I： estimation theory [M]. [S.l.]： Prentice Hall， 1993.