基于最短路徑的GNSS-R鏡面反射點算法

劉原華, 賀立慶, 牛新亮

(1.西安郵電大學 通信與信息工程學院， 陜西 西安 710121; 2.中國空間技術研究院 西安分院， 陜西 西安 710000)

全球導航衛(wèi)星系統(tǒng)反射信號(global navigation satellite systems-reflectometry, GNSS-R)技術利用導航衛(wèi)星偽隨機測距碼或者載波的反射信號，提取目標反射面特征[1]。與基于后向散射單基工作模式的衛(wèi)星散射計和雷達高度計相比，GNSS-R技術采用擴頻通信技術，擁有大量信號源，被廣泛應用于海面測高[2]、海面測風[3-6]、海洋鹽度[7]和土壤濕度[8]探測等領域。

在處理反射信號時，GNSS-R遙感器搜索和捕獲信號的參考中心為鏡面反射點，精確定位鏡面反射點具有重要的意義[1]。定位GNSS-R鏡面反射點的算法有很多，最先被提出的是Gleason算法[9]，該算法采用單步定長線性迭代法，對鏡面反射點進行粗估算后，利用校正增益系數不斷收斂校正，使得新舊反射點之差達到指定的公差，從而得到精確的鏡面反射點。但Gleason算法收斂次數多，算法效率低。針對這些缺點，文[10]提出了角平分線變步長算法[10]，該算法在Gleason算法的基礎上通過改進校正增益系數，降低了收斂次數。但是，收斂次數仍較多，算法效率低。為了進一步提高算法效率，Wu S C算法[11]利用幾何關系得到初值鏡面反射點，不斷判斷入射角和反射角的大小并進行收斂校正，最終得到精確的鏡面反射點。但Wu S C算法未考慮橢球面上徑向和法向的差異，導致定位誤差。針對這一缺陷，線段二分法[12]通過不斷地二分發(fā)射機和接收機構成的線段來獲得滿足條件的鏡面反射點。但是線段二分法在判斷收斂條件時，需要不斷地調整入射角和反射角的大小并進行迭代計算。而用角度當作收斂判斷條件時，GNSS-R幾何關系中的入射角和反射角與斯涅爾定律中定義的入射角和反射角不同，導致定位存在一定誤差。

為了提高GNSS-R鏡面反射點的測量準確度與算法效率，本文采用全球定位系統(tǒng)(global position system, GPS)反射信號為研究對象，擬從GNSS-R鏡面反射點的定義[1]出發(fā)，利用反射路徑長度最短作為收斂判斷依據，迭代計算出鏡面反射點的位置信息。

1 鏡面反射點幾何關系及其定義

1.1 鏡面反射點的幾何關系

GNSS-R鏡面反射點的幾何關系模型如圖1所示。其中S表示鏡面反射點的位置，T為衛(wèi)星發(fā)射機的位置，R為衛(wèi)星接收機的位置，Rt和RR分別表示發(fā)射機和接收機到鏡面反射點的距離，θT表示發(fā)射機相對于鏡面反射點與地球面法線之間的夾角(入射角)，θR為衛(wèi)星接收機相對于鏡面反射點與地球面法線之間的夾角(反射角)，Re表示地球半徑。

圖1 GNSS-R鏡面反射點的幾何關系

1.2 鏡面反射點定義及其特征

GNSS-R鏡面反射點是GPS衛(wèi)星發(fā)射機經反射點到接收機的路徑最短的點[1]。令GNSS衛(wèi)星發(fā)射機在WGS-84參考系中的位置坐標為(xT,yT,zT)；GNSS衛(wèi)星接收機在WGS-84參考系中的位置坐標為(xR,yR,zR)；鏡面反射點在WGS-84坐標系中的位置坐標為(xS,yS,zS)。

反射路徑的長度fS的數學表達式為

當且僅當fS取最小值時，(xS,yS,zS)為鏡面反射點的位置坐標。

GNSS-R鏡面反射點有3個特征[13]，即鏡面反射點必定位于地球表面；鏡面反射點是衛(wèi)星發(fā)射機經反射點到接收機路徑最短的點；鏡面反射點滿足斯涅爾定律，即在鏡面反射點處，其入射角等于反射角。

2 基于最短路徑的鏡面反射點定位

2.1 最短路徑的原理

在迭代計算時，本文縮小收斂區(qū)間的判斷依據是比較反射路徑的長度，最終在給定的精度范圍內找到最短反射路徑對應的鏡面反射點的位置信息。

算法迭代計算的幾何關系如圖2所示，衛(wèi)星發(fā)射機T位置坐標表示為向量a=[xT,yT,zT],衛(wèi)星接收機R位置坐標表示為向量b=[xR,yR,zR],Sa,Sb分別為T,R 2點的星下點，∠TOR的角平分線與線段TR交點為P，其星下點為SP，C點為線段TR的二分點，其星下點SC。

圖2 最短路徑法迭代計算的幾何關系

為了避免采用單個點P做判斷出現的錯判舍棄區(qū)間錯誤，在應用最短路徑法確定鏡面反射點的迭代計算時，需要在線段TR中插入點P和C。通過判斷點P和C的位置，判斷反射路徑的長度，舍棄較長的反射路徑所處的區(qū)間，從而保證較短的反射路徑所處的區(qū)間被留下，提高算法的準確性。

2.2 最短路徑法的步驟

最短路徑法的步驟如下。

步驟1確定初始搜索區(qū)間為線段TR，設定精度ε為1米。

步驟2計算衛(wèi)星發(fā)射機T與接收機R的星下點Sa和Sb的坐標向量，其表達式分別為

步驟3求∠TOR的角平分線與線段TR的交點P及線段TR的中點C的位置坐標。并且，分別求得P,C 2點星下點Sp,SC的位置坐標。

步驟4若Sa與Sb之間的距離小于精度ε，則鏡面反射點設定為SP，結束；否則轉步驟5。

步驟5分別計算發(fā)射機位置T和接收機位置R到星下點位置SP與SC的反射路徑長度

fSP=|TSP|+|RSP|,fSC=|TSC|+|RSC|。

步驟6判斷SP與SC的位置。若Sa與SP之間的距離小于Sa與SC之間的距離，則轉步驟7；否則轉步驟8。

步驟7若fSP

步驟8若fSP≤fSC，則將C點設定為T點，新的搜索區(qū)間為線段CR，轉步驟2；否則將P點設定為R點，新的搜索區(qū)間設為線段TP，轉步驟2。

最短路徑法的流程如圖3所示。

圖3 最短路徑法的流程

3 仿真結果及分析

使用來自英國UK-DMC衛(wèi)星于2004年11月16日測到的海洋數據中的3組數據來驗證，數據組1發(fā)射機的位置坐標為(1 366 322.127 94，-23 127 631.879 20，-11 982 698.916 85)，接收機的位置坐標為(335 492.491 51，-6 460 107.391 44，-2 813 071.858 19)；數據組2中發(fā)射機的位置坐標為(-11 178 791.991 29，-13 160 191.204 99，20 341 528.127 54)，接收機的位置坐標為(-4 069 896.703 38，-3 583 236.963 73，4 527 639.271 75)；數據組3中發(fā)射機的位置坐標為(-11 130 828.076 797，-13 167 090.702 929，20 363 562.269 411)，接收機的位置坐標為(-4 163 790.851 951 254，-3 618 216.540 713 221，4 413 520.900 955 853)。采用WGS-84地球橢圓模型中的地球半徑值作為地球半徑。利用MATLAB軟件仿真。

3.1 算法的正確性

由GNSS-R鏡面反射點的3個特征可知，鏡面反射點的位置必定在地球表面。仿真數據組1計算出的大地坐標系下鏡面反射點位置信息為南緯23°64′37″，西經87°02′71″，高度為0；數據2的位置信息為北緯40°05′03″，西經138°63′84″，高度為0；數據組3的位置信息為北緯38°85′12″，西經139°01′03″，高度為0。滿足了GNSS-R鏡面反射點的位置必須在地球表面這一特征。

根據定義，GNSS-R鏡面反射點是衛(wèi)星發(fā)射機經反射點到接收機路徑最短的點。數據組1、2、3的反射路徑長度曲線分別如圖4所示。

圖4 反射路徑長度曲線

可以看出，數據組1、2、3分別在迭代21、23、23次時，反射路徑的長度達到最小并趨于穩(wěn)定，從而滿足了反射路徑最短要求。

根據斯涅爾定律，在GNSS-R鏡面反射點處，其入射角等于反射角。數據組1、2、3的入射角與反射角絕對值之差曲線如圖5所示。

圖5 入射角與反射角的絕對值之差曲線

可以看出，最短路徑法的收斂速度較快，當迭代10次后，入射角與反射角的絕對值之差趨于零，即2角相等。滿足了GNSS-R鏡面反射點入射角等于反射角的要求。

以上表明所提出的算法滿足了GNSS-R鏡面反射點定義的3個特征，證明了算法的正確性。

3.2 算法的效率

為了驗證最短路徑法的效率，分別采用Gleason算法、線段二分法與最短路徑法計算鏡面反射點，不同算法確定的鏡面反射點位置坐標和迭代次數如表1所示。可以看出，Gleason算法的迭代次數最多，效率最低，而線段二分法和最短路徑法的迭代次數相當，效率相當。

表1 不同算法確定的鏡面反射點位置結果對比

最短路徑法的收斂迭代次數相較Gleason算法的4 600多次，效率提升了98%以上。其原因在于Gleason算法的收斂判定方法為單步定長線性迭代，收斂緩慢，甚至可能不收斂。而最短路徑法確定的收斂區(qū)間內有且僅有1個最短的反射路徑長度，必定收斂。所以，最短路徑法的迭代次數較Gleason算法少，算法效率較高。

3.3 算法的準確性

接下來驗證算法的準確性。由于Gleason算法的迭代次數在4 600次以上，算法效率較低，所以，僅比較算法效率相當的最短路徑法與線段二分法在迭代停止時，鏡面反射點處的反射路徑長度。反射路徑長度越短，說明算法的準確性越高。2種算法迭代停止時的反射路徑長度曲線如圖6所示。

(a) 數據組1

(b) 數據組2

(c) 數據組3圖6 兩種算法迭代停止時的反射路徑長度曲線

在迭代停止時，對數據組1采用最短路徑法和線段二分法確定的鏡面反射點處的反射路徑長度分別為1.358 574 761 106 m和1.358 846 186 106 m；對數據組2采用最短路徑法和線段二分法確定的鏡面反射點處的反射路徑長度分別為1.389 822 487 106 m和1.389 927 239 106 m；對數據組3采用最短路徑法和線段二分法確定的鏡面反射點處的反射路徑長度分別為1.389 594 236 106 m和1.389 852 134 106 m?？梢钥闯觯捎米疃搪窂椒ㄓ嬎愠龅溺R面反射點處的反射路徑長度均少于線段二分法，說明最短路徑法的準確性高于線段二分法。原因在于，線段二分法采用角度判定作為收斂判定方法，該判定方法的幾何關系中入射角與反射角是由地球球心與鏡面反射點的連線來確定，而不是由地球表面的平面法線確定的。在斯涅爾定律中，入射角、反射角是平面法線與入射線和反射線所形成的夾角，因此用角度作為收斂判斷條件存在一定的誤差。最短路徑法的判定方法是尋找最短反射路徑長度，而GNSS-R鏡面反射點的定義是衛(wèi)星發(fā)射機經反射點到接收機路徑最短的點。因此，最短路徑法滿足了GNSS-R鏡面反射點定義的要求，算法結果更加準確。

綜上所述，最短路徑法計算出的結果滿足了GNSS-R鏡面反射點的3個特征，證明了該算法預測GNSS-R鏡面反射點的正確性；且與Gleason算法和線段二分法相比，最短路徑法在算法效率與準確性上有一定的優(yōu)越性。

4 結語

基于GNSS-R的基本幾何關系與鏡面反射點定義，提出了一種用最短反射路徑作為收斂判斷條件來計算鏡面反射點位置信息的方法。該算法在衛(wèi)星發(fā)射機與接收機形成的線段上，插入該線段二分點以及由發(fā)射機、地球球心與接收機3點所形成角的角平分線和發(fā)射機與接收機形成線段的交點。計算并比較該兩點的星下點的反射路徑長度，在給定的精度范圍內找到最短反射路徑對應的鏡面反射點的位置坐標。利用MATLAB軟件仿真，結果顯示，最短路徑法滿足鏡面反射點的3個特征，可以正確定位鏡面反射點。算法的收斂迭代次數與線段二分法基本相當，遠低于Gleason算法。與線段二分法相比，算法的準確性較高。