光子探測激光測距的目標跟蹤建模與仿真

蔣 燕，孫志峰，周智睿，孫 通，饒 強，殷天峰

(國網湖北省電力有限公司信息通信公司, 武漢 430000)

0 引 言

由于空中各種噪聲的影響，傳統的激光雷達不能直接根據極其微弱的回波信號來測量目標的相對位置，因此在接收端需要有一個靈敏的探測器來接收微弱的回波信號。近年來，單光子計數方法利用弱光照射下光子探測器輸出電信號自然離散的特點，采用脈沖甄別技術和數字計數技術識別并提取出極其微弱的信號[1-2]。單光子探測技術信噪比高、抗漂移性好、時間穩定性好，便于計算機模擬仿真。但是到目前為止，關于單光子計數技術的文獻大多集中于單光子計數硬件系統的設計[3-4]，對于數據處理部分也僅局限于測量靜態目標的相對距離[5-8]，目前還未查到關于光子探測激光測距跟蹤運動目標的相關文獻。

本文以單光子技術測量運動目標的相對距離為背景，利用實驗得出的光子探測器的回波強度和噪聲強度經驗公式，對光子探測器接收的回波數據進行仿真并采用圖像處理的方法，從大量的回波數據中提取出目標的運動軌跡。

1 回波模型

某一飛行器正以一定的速度接近一個以恒定速度運動的目標。飛行器上裝有激光器來發射脈沖，同時有光子探測器來接收光子回波。在運動過程中，飛行器的運動狀態在一維空間中可以隨時改變(即可以加速、減速和勻速)。通過發射機接收到的光子，根據脈沖雷達測距公式可以生成時空二維回波事件的散點圖。

1.1 飛行器與目標之間的相對距離

本文通過分析目標相對飛行器的運動狀態建模和仿真，涉及到運動狀態(距離、速度和加速度)均為目標相對飛行器的運動狀態。假設Ri(i=1,2,3,…)為一個周期開始的初始距離；t1為發出的脈沖與目標相遇所用的時間；t2為脈沖遇到目標產生回波信號返回過程中被飛行器接收所用的時間；v0為初始速度；R為飛行器收到信號回波時與目標的相對距離；c為光速。

對于勻速運動，接收回波光子的時間計算表達式如下：

對于勻加速運動，a為加速度，接收回波光子的時間計算表達式如下：

計算完延遲時間后，根據測距公式即可計算出理想運動狀態下飛行器根據接收的信號回波測量出與目標之間的距離。

1.2 接收信號的回波光子強度與回波噪聲強度

激光測距時，影響系統性能的主要因素有：激光單脈沖的能量Et、目標反射率ρ、接收系統光學效率ηr、激光發射系統光學效率ηt、大氣透過率τ、光學系統口徑D、激光發散角α和探測器靈敏度等。激光測距對遠距離的點目標測距時，探測器上接收到的能量[9]表達式如下：

式中：Er為探測器接收到的功率；Ar為有效接收面積；At為目標雷達的截面積；θi為太陽入射方向與目標法向夾角。

當目標處于日照區時，目標表面的反射光會降低探測器的性能，此外，探測器本身的暗計數也會在一定程度上引入噪聲。上述噪聲信號會產生較高的虛警數據，給后續的距離數據提取帶來很大的壓力。因此，在系統設計中需要最大限度地降低背景噪聲對探測器的影響。假定視場角為θfov，太陽光譜照度為E，ηq為光子探測器的量子效率，濾光片帶寬為Δλ，H為普朗克常量，測距系統接收到的背景光噪聲功率P的表達式如下：

探測器上每秒接收到的平均光子事件次數為

式中：hν為532 nm單光子能量；ndark為探測器本身的暗計數；tgate為門控信號寬度；在激光采樣次數fpulse一定的情況下，噪聲事件的次數主要與濾波片帶寬、探測器象元視場角、望遠鏡口徑、探測器暗計數以及門控信號寬度有關。

1.3 回波概率

根據統計光學理論，激光測距中目標表面散射回來的激光回波信號在強度分布上是一個泊松分布[2]。因此，此處可以采用泊松分布來近似系統地探測模型中接收回波的個數。根據光子計數激光測距的原理，將探測的距離分成若干個距離柵格(每個距離柵格作為系統的觀測精度，在該精度內最多只能接收一個光子回波信號)，每個距離柵格對應激光探測過程中的一個時間段，該時間段內探測器能探測到光子回波信號的概率為

式中：Nsn(i)=Nn(i)+Ns(i)，Nsn(i)為第i個距離單元內回波光子信號與回波噪聲的總個數，Ns(i)為第i個距離單元內回波光子信號的平均個數，Nn(i)為第i個距離單元內回波噪聲的平均個數；Psn(i)為第i個距離單元內光子回波出現的概率。

1.4 生產回波事件

高噪聲、低信號強度是光子探測回波數據的特點。在回波事件中包含回波光子信號和回波噪聲。根據統計光學可知，在基于光子計數的激光測距中，探測器接收到的光子數目滿足負二項式分布[2]，由于負二項式分布的曲線是由一個個小階梯構成，因此在有限的時間段內，可以將探測器接收到的回波信號看作是均勻分布。

在一段時間間隔內，回波光子落在未知斜率的一段直線上(假定飛行器進行勻速運動)，背景噪聲在距離門內是隨機分布的，可以假定背景噪聲遵循均勻分布特征。根據式(6)分別計算出每個周期內回波光子和噪聲出現的總概率，然后按照每個周期的概率生成均勻分布的隨機數，滿足隨機數小于等于理論概率的事件記作1，不滿足的記作0。然后在每個周期內提取Psn=1處的回波作為接收器接收到的事件。每個周期內的Psn=1的事件都對應有出現回波，根據每個周期內回波出現的時間點，可以計算每個周期內接收到回波事件所對應的目標與飛行器間的距離，這些回波事件是一些雜亂的散點(即為Psn=1的回波事件)，如圖1所示，信號被淹沒在噪聲中。

圖1 200 ms內生成的回波事件Figure 1 Echo events generated within 200 ms

2 軌跡提取

在發射和接收脈沖之前，信號的接收者并不知道飛行器與目標的相對距離、速度以及加速度。因此，為了之后能夠精確地跟蹤目標相對于飛行器的運動軌跡，必須先積累一段時間得到飛行器初始運動狀態的先驗信息。

2.1 目標運動軌跡范圍的確定

在接收到的回波數據中，信號一般分布在散點密集的地方，噪聲則均勻分布在整個探測區域內。本文通過將N個周期內接收到的回波事件個數進行疊加找到回波數目最大的時刻(該時刻為接收機最后一次接收到回波的時刻)，由此初步估計接收機接收回波光子信號的位置范圍。

假設目標與飛行器的初始距離為15 km，飛行器以2 000 m/s勻速運動，積累的觀測周期數目為500，每個周期為1 ms，理論上出現回波數目最多的個數在90 000 ns左右。

確定好了目標在積累時間段末的大致位置后，以數目最多的回波所在位置為基準，向前后各取一定距離的門寬(前門限到基準位置要遠小于后門限到基準位置。為了不丟失信號，一般把后門限盡量設置大一些。)，得到篩選出來的回波事件。

2.2 飛行器與目標相對距離的提取

由于接收到的回波事件是一群密集的散點，如圖1所示，信號回波淹沒在噪聲[10]中，很難發現并提取出來。即使將每個周期的回波事件進行累加得到一條特征線(可以認為是隨時間變化的飛行器和目標的相對距離，因為信號集中在某一個區域出現，而噪聲只是均勻分布在整個區域)，要想在散亂的點中直接提取這條曲線還是非常困難的。本文先用回波數據進行預處理得到目標運動的模型，在已知模型的基礎上用隨機抽樣一致性(Random Sample Consensus，RANSAC)算法精確提取回波光子信號。

由于數據量比較大，本文首先采用常用的大數據處理方法萊特準則剔除部分無關緊要的數據。對于剩下的散點，將其分成若干個小區間后，通過統計每個區間的光子數目粗略提取目標軌跡的范圍。每個周期內信號出現的概率總是大于噪聲出現的概率，可以根據回波事件出現的個數提取回波光子。以橫坐標為觀測時間(時間間隔為一個周期)，縱坐標為觀測距離，建立平面坐標系。然后將這個平面分割成棋盤格類型的若干個小區間，對每個區間統計回波事件的個數Nij(i=1,2,…;j=1,2,…)。最后搜索并提取每列中Nij(j=1,2,3,…)最大的那個區間并用多項式對這些位置上的點初步擬合后得到目標與飛行器相對距離變化的軌跡方程，為下一步RANSAC算法提取回波光子提供初始的數學模型。

RANSAC[11-12]為一種穩健的模型參數估計算法。在一組包含異常數據的樣本中，通過與數據相關的數學模型，給出有效的樣本數據。RANSAC算法的假設是正確數據內點(inliers)可以被模型描述，而異常數據(outliers)則無法適應數學模型，可以認為是數據集中的噪聲。數據預處理后，可以將軌跡方程當作先驗數學模型，然后用RANSAC算法對模型進行優化并提取目標的相對運動軌跡。

在積累的時間段用RANSAC算法精確提取目標的相對運動軌跡后，就知道了飛行器相對目標的距離隨著時間的變化規律，為后面的跟蹤提供先驗信息。具體的跟蹤步驟如下：

(1) 根據積累時間段內目標的相對運動軌跡對后N(N的大小可調整)個周期預測目標的相對位置；

(2) 以每個相對位置為中心，取一定距離寬度的門構成波門；

(3) 取波門內距離波門中心最近的回波事件作為回波光子。如果波門內沒有光子就將波門中心作為回波光子的位置。

(4) 將N個周期內得到的回波光子與前M(M的個數可調)個周期的回波進行最小二乘擬合，更新得到新的相對運動軌跡。重復以上過程直至跟蹤完畢。

3 實驗結果與分析

3.1 預測目標的先驗信息

該仿真依據的實驗樣機的具體參數為：激光發散角α=1 mrad，空間目標雷達截面積為0.5 m2，激光波段為532 nm，Et=20 μJ,ρ=0.1，ηr=0.6，ηt=0.9，τ=1，D=60 mm。假設已知目標相對飛行器的初始距離為15 000 m，目標的相對速度為2 000 m/s，相對加速度為-10 m/s2，激光發射脈沖的周期T為1 ms，積累時間為500個周期。根據式(1)～(2)和式(6)計算出500個周期內飛行器接收回波光子的概率Psn(i)(i=1,2,3,…,500)的理論值，然后根據這500個概率按照1.4節所述的方式生成并篩選回波事件，如圖1所示。

按照2.1節所述方法初步提取目標的相對運動軌跡范圍后，按2.2節所述方法對回波散點數據進行預處理。將縱軸平均分成1 000份，橫軸平均分成500份，生成50萬個小區間Iij(i=1,2,…,1 000,j=1,2,…,500)。對每個區間內，統計該區域的回波個數為Nij(i=1,2,…,1 000;j=1,2,…,500)。提取每列Nij數目最大的區間作為回波光子所在的區域，如圖2所示。

圖2 均值濾波后的回波信號Figure 2 Echo signal after mean filtering

將這些點擬合得到目標運動軌跡方程，為下一步的RANSAC算法提供初始的數學模型。最后根據該初始模型，用RANSAC算法優化該模型并提取積累時間段內目標的相對運動軌跡，信號回波與目標理想運動軌跡的關系如圖3所示。

圖3 RANSAC算法提取目標運動軌跡Figure 3 Target trajectory extracted by RANSAC algorithm

3.2 根據先驗信息跟蹤運動的目標

跟蹤從500 ms開始，跟蹤周期500～510 ms的目標相對運動軌跡。圖4所示為RANSAC算法優化回波點云后得到的目標運動軌跡。

圖4 RANSAC算法預測的10個周期內目標運動軌跡Figure 4 Target trajectory predicted byRANSAC algorithm in 10 periods

圖4中紅圈為觀測到的脈沖，由圖可知，觀測到的距離與預測的距離有一定差異。一般觀測到的脈沖比較少，無法直接提取有用信息，需要使用部分舊脈沖作為補充重新擬合。此處取前面200周期(周期大小可以調整)即可，使用300～510周期的數據重新擬合，擬合后的結果如圖5所示。

圖5 300～510 ms回波點云重新更新的運動軌跡Figure 5 The trajectory within 300～510 ms after echo point cloud update

圖5中藍點是前一段時間積累的回波信息，紅圈是新觀測到的新回波信息。 圖中倒數第2個點就是一個噪聲脈沖，這時只需要將誤差比較大的個別點剔除，再擬合即可。下一個大周期的跟蹤按該方法執行，直到運動結束，最后得到的目標跟蹤軌跡與理想目標運動軌跡的誤差如圖6所示。

圖6 目標理想相對運動軌跡與跟蹤軌跡的誤差圖Figure 6 Error diagram of ideal relative Motion trajectory and tracking trajectory of target

由圖可知，越往后跟蹤誤差越大。尤其是當飛行器的運動狀態發生了突變，跟蹤精度會大大下降。對于目標突變狀況的跟蹤，還需要建立特定的模型來實現。

4 結束語

本文模擬了光子探測激光測距對慢速運動(相對光速而言)目標的跟蹤。根據實驗所得的經驗公式，用統計的方法生成了回波事件，然后將圖像處理的相關算法運用到大量回波數據中提取目標的軌跡。仿真后得出目標的軌跡信息與目標理論的運動軌跡誤差在5 m以內。光子探測激光測距能在較強的背景噪聲條件下獲得較高的定位精度，將其作為北斗導航系統的一個輔助定位模塊會大大提升北斗導航系統定位的精確度。