劉原華,孫毅,牛新亮
(1.西安郵電大學(xué) 通信與信息工程學(xué)院,陜西 西安 710121;2.中國(guó)空間技術(shù)研究院西安分院,陜西 西安710100)
作為一種新型的微波遙感技術(shù)、全球衛(wèi)星導(dǎo)航系統(tǒng)反射(GNSS-R)技術(shù)具有豐富的信號(hào)源,低成本、低功耗、全天時(shí)、全天候等優(yōu)點(diǎn),它可以獲得各種全球海洋動(dòng)力學(xué)和環(huán)境參數(shù)如海面高度(SSH),海面風(fēng)場(chǎng)、海冰和海水鹽度. 此外,GNSS-R衛(wèi)星平臺(tái)的成本和有效載荷相對(duì)較低,結(jié)合專用指向天底天線、軌道和星座的參數(shù)優(yōu)化設(shè)計(jì),可形成6~8顆衛(wèi)星星座,完成多項(xiàng)任務(wù),具備優(yōu)良性能.
在已知接收機(jī)坐標(biāo)情況下,時(shí)間延遲和反射信號(hào)路徑長(zhǎng)度的誤差主要來(lái)自海洋表面的鏡面反射點(diǎn)定位誤差和大氣信號(hào)傳輸誤差. 此外,雙基雷達(dá)系統(tǒng)的輸出信噪比與信號(hào)傳播路徑有關(guān),而信號(hào)傳播路徑又與鏡面反射點(diǎn)的位置直接相關(guān). 在海面遙感應(yīng)用中,鏡面反射點(diǎn)的精確位置信息也可以用來(lái)確定反射信號(hào)是來(lái)自海面、陸地還是冰面. 鏡面反射點(diǎn)作為參考點(diǎn)GNSS-R信號(hào)的反射幾何及相關(guān)參數(shù),其定位誤差影響GNSS-R遙感產(chǎn)品的準(zhǔn)確性[1].
鏡面反射點(diǎn)定位方法可分為物理方法和幾何方法. 物理方法以接收信號(hào)的處理為基礎(chǔ),幾何方法以鏡面反射點(diǎn)滿足的幾何條件為基礎(chǔ)進(jìn)行研究. 現(xiàn)有的鏡面反射點(diǎn)搜索方法包括黃金分割法[2]、二分法[3]、Gleason算法[4]、S.C.WU[5]、C.Wagner[6]、橢球面算法[7]、AdaGrad自適應(yīng)算法[8]等. Gleason算法精度較高,但是迭代次數(shù)較多. 黃金分割法迭代較少,算法經(jīng)典,但仍然存在一定的提升空間. 二分法模型簡(jiǎn)單,易于計(jì)算,不過(guò)由于采用的是角度作為迭代依據(jù),因此存在一定的誤差. 本文在分析基本幾何反射關(guān)系的基礎(chǔ)上,提出了2種鏡面反射點(diǎn)搜索算法,進(jìn)一步改進(jìn)鏡點(diǎn)搜索的算法.
首先,提出了基于斐波那契法的鏡點(diǎn)搜索算法,并和黃金分割法進(jìn)行比較,從數(shù)論的角度證明了其算法的優(yōu)越性. 另外提出了一種基于求導(dǎo)極值的搜索算法,該方法省去了不斷迭代的過(guò)程,從而優(yōu)化了計(jì)算. 仿真結(jié)果表明,這種算法能夠計(jì)算出鏡點(diǎn)的位置,對(duì)鏡點(diǎn)的搜索具有一定的意義.
在討論直射信號(hào)與反射信號(hào)的關(guān)系中,引入了鏡面反射點(diǎn)的概念[9],即從反射區(qū)域反射的反射信號(hào)中路徑最短的理論反射點(diǎn). 其數(shù)學(xué)表達(dá)式為
f=|R-S|+|T-S|=
(1)
式中:R為接收機(jī)坐標(biāo);T為衛(wèi)星坐標(biāo);S為球面坐標(biāo);其大致幾何關(guān)系如圖1所示,圖中R為接收機(jī);T為衛(wèi)星;O為地心;Re為地球半徑,在RT線段中一點(diǎn)為P,P所對(duì)應(yīng)的星下點(diǎn)為S,R的星下點(diǎn)為a,T的星下點(diǎn)為b,ar和at分別表示s到R、T之間的仰角.
鏡面反射點(diǎn)有以下三個(gè)特征:
1)位于橢球面弧段ab上;
2)入射與出射角相等,即ar=at;
3)到R、T點(diǎn)距離和為最短.
圖1 基本幾何關(guān)系
斐波那契法與黃金分割法類似,也是用于單峰函數(shù). 在計(jì)算過(guò)程中,也是第一次迭代需要2個(gè)試探點(diǎn),與黃金分割法的主要區(qū)別在于區(qū)間長(zhǎng)度縮短比例不是常數(shù),而是斐波那契數(shù)[10]. 由于斐波那契數(shù)列自身規(guī)律,迭代時(shí)相較黃金分割法有一定提升,因此帶來(lái)了計(jì)算次數(shù)的減少.
設(shè)有數(shù)列{Fk},滿足條件:
(2)
可將斐波那契數(shù)列列表如表1所示.
表1 斐波那契數(shù)列
k0123456…Fk11235813…
其通項(xiàng):
(3)
斐波那契在迭代中計(jì)算試探點(diǎn)的公式:
(4)
其具體計(jì)算步驟如下:
1)給定初始區(qū)間[R,T]和最終區(qū)間長(zhǎng)度L=1 m,求計(jì)算函數(shù)值的次數(shù)n. 使
Fn≥|R-T|/L.
(5)
置k=1,計(jì)算試探點(diǎn)λ1和μ1:
(6)
(7)
計(jì)算λ1和μ1所對(duì)應(yīng)的星下點(diǎn)到RT距離和f(λ1)和f(μ1).
2)若f(λ1)>f(uk),則轉(zhuǎn)步驟3);若f(λk)≤f(λk),則轉(zhuǎn)步驟4).
3)令ak+1=λk,bk+1=bk,λk+1=μk
計(jì)算試探點(diǎn)μk+1:
(8)
若k=n-2,則轉(zhuǎn)步驟6);否則,計(jì)算函數(shù)值f(μk+1),則轉(zhuǎn)步驟5).
4) 令ak+1=ak,bk+1=μk,μk+1=λk,計(jì)算試探點(diǎn)λk+1:
(9)
若k=n-2,則轉(zhuǎn)步驟6);否則,計(jì)算f(λk+1),轉(zhuǎn)步驟5).
5)置k=k+1,轉(zhuǎn)步驟2).
在相同n次迭代計(jì)算時(shí),用黃金分割法,最終區(qū)間長(zhǎng)度為
(10)
用斐波那契法,最終區(qū)間長(zhǎng)度為
(11)
由式(10)~(11)可知:
(12)
圖2 計(jì)算流程圖
在WGS84坐標(biāo)系下,當(dāng)R(xr,yr,zr)、T(xt,yt,zt)點(diǎn)位置給定之后. 假設(shè)變量K∈[0,1],令P=R+k·(T-R)=(xr+k·(xt-xr),yr+k·(yt-yr),zr+k·(zt-zr)),當(dāng)K從0連續(xù)變化到1時(shí),P就取遍了RT線段中的所有點(diǎn),即RT線段中的點(diǎn)均可用變量K和常數(shù)來(lái)表示,如圖3所示.
圖3 變量K表示RT線段
2)由映射算法知此緯度所對(duì)應(yīng)的球徑ρ
凳子上擺了本書(shū)。書(shū)打開(kāi),封皮朝上,我看到了兩個(gè)字:詩(shī)集。我很詫異,柴米河,沒(méi)幾人識(shí)字,想不到別呦呦竟愛(ài)看書(shū),這一點(diǎn)倒和我一樣。我愛(ài)看談狐說(shuō)鬼的書(shū),詩(shī),幾乎不看。
那么P所對(duì)應(yīng)的星下點(diǎn)S=p×ρ=
3)R、T為已知數(shù),即可建立星下點(diǎn)S到R、T的距離和函數(shù):
f(k)=
其圖像如圖4所示.
圖4 函數(shù)大致圖像
4)此時(shí),鏡面反射點(diǎn)的搜索問(wèn)題,也就轉(zhuǎn)化函數(shù)f(k)的極值問(wèn)題. 利用MATLAB求導(dǎo)函數(shù)對(duì)函數(shù)f(k)進(jìn)行一階求導(dǎo),得到f′(k). 再利用solve解方程函數(shù)求解方程.
f′(k)=0.
求解出K值之后,回代入P的坐標(biāo),即可求得所對(duì)應(yīng)的星下點(diǎn)S,此點(diǎn)即為鏡面反射點(diǎn).
傳統(tǒng)做法大部分基于嘗試法,在RT區(qū)間內(nèi)不斷進(jìn)行迭代,根據(jù)距離最小或者角度相等的原則,不斷縮小搜索區(qū)間,最終區(qū)間小于一定長(zhǎng)度即可得到目標(biāo)點(diǎn). 而這樣的做法在函數(shù)形式已知的情況下較為冗余,增加了算法編程的復(fù)雜度.
本文做法基于距離和函數(shù)的極值進(jìn)行求導(dǎo),相較二分法等基于角度進(jìn)行不斷迭代的做法,避免了由于徑向和法向不統(tǒng)一而帶來(lái)的問(wèn)題. 同時(shí),本文算法只需要一步求導(dǎo)、一步解方程,相較黃金分割法,又省去了不斷迭代的步驟,因此優(yōu)化了計(jì)算. 對(duì)鏡點(diǎn)的搜索具有一定參考意義.
本文采取導(dǎo)航模擬器生成全球衛(wèi)星導(dǎo)航系統(tǒng)(GNSS)接收機(jī)和衛(wèi)星在WGS84坐標(biāo)系下的位置信息,接收機(jī)R的坐標(biāo)為[-1747836.86, 4316984.79,-5091148.03],衛(wèi)星T的坐標(biāo)為[18524153.94,1465929.77,-19218107.19]. 鏡面反射點(diǎn)的參考坐標(biāo)為[-648295.67,3698445.399,-5138472]. 基于MATLAB對(duì)本文所提出的方法進(jìn)行四次仿真,閾值選取1 m. 通過(guò)分析四種方法的計(jì)算結(jié)果、迭代次數(shù)、計(jì)算時(shí)間,比較說(shuō)明本文提出的兩種方法的正確性,優(yōu)越性. 由于篇幅原因,在此僅列出一次仿真的結(jié)果,如表2所示.剩余仿真結(jié)果以柱狀圖的形式展示如圖5和圖6所示.
表2 鏡面反射點(diǎn)位置估計(jì)仿真
圖5 和參考點(diǎn)之間的距離
圖6 四次仿真的計(jì)算時(shí)間
如表2所示,極值求導(dǎo)做法距離參考點(diǎn)的距離是最短的,而黃金分割法和斐波那契法則次之,二分法的距離最遠(yuǎn).
進(jìn)行多次仿真后,圖5顯示前三種方法的結(jié)果基本都在10 m以內(nèi),極值求導(dǎo)的距離是最小的,二分法則偏離參考點(diǎn)50 m左右. 斐波那契法偏離的距離要小于黃金分割法但是多于極值求導(dǎo). 其主要原因在于二分法收斂原則與其他方法不同.
從四次計(jì)算的時(shí)間結(jié)果來(lái)看,極值求導(dǎo)的時(shí)間少于黃金分割法,幾近于斐波那契.而二分法的計(jì)算時(shí)間是最少的,其根本原因在于二分法模型簡(jiǎn)單,每一次計(jì)算只需要求解其中間的一個(gè)點(diǎn)的RT距離和,而黃金分割法和斐波那契法相對(duì)復(fù)雜一些. 同時(shí),極值求導(dǎo)的時(shí)間略微少于黃金分割法,近似等于斐波那契法. 綜上分析,斐波那契法在偏移距離和計(jì)算時(shí)間上對(duì)傳統(tǒng)做法略有提升. 而極值求導(dǎo)在偏移距離上有一定幅度的提升,在計(jì)算時(shí)間上和傳統(tǒng)做法相近.
隨著GNSS的迅速發(fā)展,針對(duì)鏡點(diǎn)搜索的需要,本文提出2種鏡面反射點(diǎn)搜索算法. 首先提出了斐波那契算法,從數(shù)論角度證明了該方法相較黃金分割法有一定提升,進(jìn)一步減少了迭代次數(shù). 另外,本文提出了基于極值求導(dǎo)的搜索方法,相較傳統(tǒng)做法,它省去了不斷迭代的過(guò)程,采用極值求導(dǎo)的思想計(jì)算出鏡點(diǎn),具有一定的優(yōu)越性. 最后采用導(dǎo)航模擬器生成的數(shù)據(jù)進(jìn)行仿真,結(jié)果表明,兩種方法均可以有效計(jì)算出鏡面反射點(diǎn)的位置,具有一定的借鑒意義.