朱重儒，朱立東

(電子科技大學(xué)通信抗干擾技術(shù)國(guó)家級(jí)重點(diǎn)實(shí)驗(yàn)室，成都 611731)

0 引言

無(wú)源衛(wèi)星定位技術(shù)具有定位范圍廣、偵查方式隱蔽以及易于實(shí)現(xiàn)等特點(diǎn)，在民用和軍事領(lǐng)域都發(fā)揮著極大的作用[1-2]。對(duì)于低軌道衛(wèi)星定位，多星座體制定位依賴于多顆衛(wèi)星協(xié)同定位，可采集更多輻射源信息，定位精度有所提升，但與此同時(shí)，多星座體制的時(shí)差定位對(duì)星間時(shí)鐘同步要求嚴(yán)格，時(shí)差信息對(duì)定位精度影響較大，給實(shí)際定位應(yīng)用增加了難度，單顆衛(wèi)星相對(duì)多星座體制，成本低廉，定位原理簡(jiǎn)單，避開了復(fù)雜多變的組網(wǎng)設(shè)計(jì)，且多普勒信息容易得到，因此單星定位技術(shù)仍有較好的發(fā)展前景[3]。

對(duì)于來(lái)自地面的未知載頻信號(hào)，需要進(jìn)行載頻估計(jì)，目前常用的載頻估計(jì)方法主要是硬件法和軟件法[4]。硬件法方面，不同的工程項(xiàng)目會(huì)設(shè)計(jì)出不同電路圖，硬件法適用性差且成本高昂，不適合用于定位系統(tǒng)中的載頻估計(jì)環(huán)節(jié)。對(duì)于軟件法，算法設(shè)計(jì)更加靈活方便，常見的算法有傅里葉算法、線性預(yù)測(cè)法、最小二乘法、Chirp-Z變換法等。傅里葉算法頻率估計(jì)精度低，且傅里葉變換計(jì)算量大，線性預(yù)測(cè)法更適合具有線性關(guān)系的參數(shù)估計(jì)模型[5]，最小二乘法初值的選擇具有任意性，載頻估計(jì)結(jié)果受初值選取的影響較大[6]，Chirp-Z變換法要求確定待估計(jì)頻率的范圍，以對(duì)頻譜進(jìn)行細(xì)化。受計(jì)算復(fù)雜度和測(cè)頻精度的相互制約，衛(wèi)星對(duì)地面輻射源進(jìn)行精確實(shí)時(shí)的定位得不到保證，如何平衡計(jì)算復(fù)雜度和測(cè)頻精度，是單星定位中亟需解決的問(wèn)題。

在單星定位實(shí)現(xiàn)過(guò)程中，需要采集衛(wèi)星不同時(shí)刻的星歷信息和接收機(jī)信號(hào)的參數(shù)，數(shù)據(jù)量大。針對(duì)最小二乘載頻估計(jì)法具有計(jì)算復(fù)雜度低、測(cè)頻精度不高的特點(diǎn)，本文提出一種基于加權(quán)最小二乘載頻估計(jì)的單星定位算法，在保證較小計(jì)算復(fù)雜度的同時(shí)，提高測(cè)頻精度，使得定位點(diǎn)更加靠近真實(shí)的位置。另外，本文分別研究了通信信號(hào)載頻、衛(wèi)星接收信號(hào)時(shí)長(zhǎng)兩方面對(duì)單星定位結(jié)果的影響，選擇合適的信號(hào)載頻和接收信號(hào)時(shí)長(zhǎng)，可提高定位精度。

1 單星定位原理

通常通信信號(hào)和雷達(dá)信號(hào)均被認(rèn)為是窄帶信號(hào)，本文假設(shè)地面輻射源的發(fā)射信號(hào)和衛(wèi)星接收機(jī)接收信號(hào)均為如式(1)所示的窄帶信號(hào)[7]：

s(t)=A(t)·exp(j·2πf0t)

(1)

其中，f0是發(fā)射機(jī)信號(hào)的載頻，A(t)為發(fā)射信號(hào)的包絡(luò)。對(duì)于窄帶信號(hào)的傳輸，可以知道，當(dāng)衛(wèi)星處于運(yùn)動(dòng)狀態(tài)時(shí)，接收機(jī)的接收信號(hào)會(huì)受到多普勒效應(yīng)的影響，接收信號(hào)的載頻等于發(fā)射信號(hào)的載頻疊加多普勒頻率，由信號(hào)模型可知，接收信號(hào)的瞬時(shí)頻率為：

f=f0+fd

(2)

其中，fd表示多普勒頻率，且多普勒頻率與衛(wèi)星和輻射源的徑向速度滿足式(3)的關(guān)系式，c表示電磁波傳輸速度，大小為2997 924 58 m/s。

(3)

r′(t)是衛(wèi)星和目標(biāo)輻射源之間徑向距離的導(dǎo)數(shù)，即二者的徑向速度的大小。徑向速度又可以用相對(duì)位置矢量和相對(duì)速度矢量表示，如式(4)所示[8]：

(4)

其中，T表示矢量的轉(zhuǎn)置，v是衛(wèi)星與目標(biāo)之間的相對(duì)速度矢量，r是衛(wèi)星與目標(biāo)之間的相對(duì)位置矢量，|r|是衛(wèi)星與目標(biāo)之間的距離大小。將式(3)和式(4)代入式(2)，則信號(hào)頻率可表示為：

(5)

根據(jù)單星無(wú)源定位的運(yùn)動(dòng)模型，在第i個(gè)時(shí)刻觀測(cè)到的多普勒頻率可表示為[9]

(6)

其中，ζi是第i個(gè)時(shí)刻的測(cè)頻誤差，是服從獨(dú)立高斯分布的白噪聲成分，F(xiàn)i(x,y,z)，可視為低軌衛(wèi)星與輻射源在相對(duì)運(yùn)動(dòng)過(guò)程中的多普勒效應(yīng)。當(dāng)衛(wèi)星過(guò)境輻射源上空時(shí)，可在N個(gè)時(shí)刻接收到輻射源的信號(hào)，N個(gè)時(shí)刻的載頻測(cè)量值用矩陣形式表示如式(7)所示：

φ=f0·F+ζ

(7)

其中，φ代表N個(gè)時(shí)刻衛(wèi)星接收機(jī)接收到的信號(hào)頻率的矢量，F(xiàn)是可以表征多普勒效應(yīng)的矢量，ζ表示信號(hào)接收時(shí)受到的干擾噪聲。在式(7)中，只有信號(hào)載頻f0是未知，由載頻估計(jì)方法計(jì)算得到fk,l，建立如式(8)的代價(jià)函數(shù)模型[10]：

J(Lk,Bi)=‖φ-fk,l·F‖2

(8)

選擇使得式(8)最小的(Lk,Bl)作為輻射源位置的解算結(jié)果，本文仿真將式(8)的代價(jià)函數(shù)取倒數(shù)，數(shù)據(jù)顯示更加直觀。

2 加權(quán)最小二乘估計(jì)原理

最小二乘法的核心思想是通過(guò)已有數(shù)據(jù)建立自變量和因變量之間的最佳對(duì)應(yīng)關(guān)系，在建立關(guān)系的過(guò)程中優(yōu)化模型參數(shù)。對(duì)于確定的目標(biāo)函數(shù)模型，自變量和因變量是已知的變量，系數(shù)或者常量是帶估計(jì)的參數(shù)，根據(jù)最小二乘計(jì)算公式建立目標(biāo)函數(shù)模型，進(jìn)行加權(quán)處理并計(jì)算估計(jì)值[11]。由式(6)和式(7)可構(gòu)建系統(tǒng)模型為：

Y=f0·X+ζ

(9)

其中，Y為模型的輸出，即衛(wèi)星接收信號(hào)載頻φ，X為模型的輸入，即X=F(x,y,z)，包含衛(wèi)星星歷等信息，ζ為高斯白噪聲，即接收信號(hào)載頻誤差，f0即為本文估計(jì)的輻射源發(fā)射信號(hào)載頻，根據(jù)最小二乘法理論，可得到f0粗略的估計(jì)值為：

f′=(XTX)-1XTY

(10)

估計(jì)得到的參數(shù)f′再代入模型中得到的輸出為Y′，估計(jì)系統(tǒng)輸出Y′與實(shí)際輸出Y存在誤差e，也可以被稱為殘差，殘差的表達(dá)式表示為：

e=Y-Y′=Y-f′X

(11)

隨機(jī)噪聲誤差項(xiàng)的方差ζ與系統(tǒng)輸入的解釋變量X存在相關(guān)性，上面的最小二乘模型存在異方差，為了使得f′為無(wú)偏估計(jì)值，通過(guò)殘差e構(gòu)造誤差權(quán)重矩陣W。

(12)

其中，權(quán)重函數(shù)為wi=tanh(e)/e， 則新的加權(quán)最小二乘參數(shù)估計(jì)可表示為[1]：

f′=(XTWX)-1XTWY

(13)

3 基于加權(quán)最小二乘載頻估計(jì)的單星定位算法

采用通常使用的WGS-84地球模型，對(duì)于地球表面的目標(biāo)輻射源而言，至少需要兩個(gè)定位曲面就可以確定其位置坐標(biāo)。當(dāng)衛(wèi)星觀測(cè)了足夠次數(shù)的信號(hào)頻率后，在星下足夠大的區(qū)域內(nèi)以經(jīng)緯度進(jìn)行網(wǎng)格劃分，定義二維網(wǎng)格點(diǎn)集∑。令k,l分別表示網(wǎng)格的縱橫計(jì)數(shù)，對(duì)每個(gè)網(wǎng)格點(diǎn)(Lk,Bl)∑，結(jié)合星歷信息中的衛(wèi)星位置坐標(biāo)、衛(wèi)星速度和加速度以及衛(wèi)星接收機(jī)接收到的信號(hào)頻率，由式(7)和式(14)計(jì)算f0的加權(quán)最小二乘解為：

fk,l=(FTWF)-1FTWφ

(14)

對(duì)每一個(gè)求得的載頻計(jì)算式(15)代價(jià)函數(shù)，選擇使得最小的(Lk,Bl)作為輻射源位置的最優(yōu)估計(jì)。

J(Lk,Bl)=‖φ-fk,l·F‖2

(15)

基于加權(quán)最小二乘載頻估計(jì)的單星定位算法的具體步驟如表1所示。

表1 單星定位算法步驟

4 實(shí)驗(yàn)仿真

設(shè)置地面站輻射源經(jīng)緯度坐標(biāo)為(118.5°,22.5°)，發(fā)射信號(hào)頻率為1.3 GHz。衛(wèi)星上的接收機(jī)只能在地面站輻射源的過(guò)境區(qū)域接收到信號(hào)，于是在表示接收機(jī)信號(hào)頻率的星歷中找到與估計(jì)載頻最接近的時(shí)刻，以這個(gè)時(shí)刻為中心，在前29 s和后30 s時(shí)間段內(nèi)，取對(duì)應(yīng)的衛(wèi)星位置、速度、加速度以及接收機(jī)接收頻率作為采樣數(shù)據(jù)，且對(duì)這60 s的接收信號(hào)頻率疊加服從高斯分布的測(cè)頻噪聲。根據(jù)表1中的步驟進(jìn)行仿真。

在STK軟件中，設(shè)置好衛(wèi)星運(yùn)動(dòng)場(chǎng)景，分別導(dǎo)出衛(wèi)星圍繞地球一周的坐標(biāo)、速度和加速度，以及衛(wèi)星接收機(jī)上的接收信號(hào)頻率。發(fā)送信號(hào)和接收信號(hào)頻率如圖1所示，當(dāng)衛(wèi)星與地面站發(fā)生相對(duì)運(yùn)動(dòng)時(shí)，接收頻率的大小發(fā)生變化。查文獻(xiàn)[12]可知，對(duì)于1.3 GHz信號(hào)頻率穩(wěn)定度是10-8，誤差在十幾赫茲量級(jí)，于是在接收端疊加均值為0，方差為15 Hz的正態(tài)隨機(jī)噪聲[12]。顯然，在取得的60 s時(shí)間內(nèi)，衛(wèi)星經(jīng)歷了靠近地面站再遠(yuǎn)離地面站的過(guò)程。

圖1 發(fā)送信號(hào)和接收信號(hào)頻率

通過(guò)衛(wèi)星的坐標(biāo)求出該時(shí)間段的星下點(diǎn)軌跡，找到星下點(diǎn)軌跡的中心點(diǎn)，經(jīng)緯度變化幅度為±5°，確定第一次網(wǎng)格搜索的經(jīng)緯度范圍，經(jīng)度范圍為[112.6169 122.6169]，緯度范圍為[17.7030 27.7030]，網(wǎng)格大小為0.1*0.1，第一次搜索得到兩個(gè)定位結(jié)果，再分別對(duì)兩個(gè)定位結(jié)果的領(lǐng)域進(jìn)行精搜索，經(jīng)緯度變化幅度為±0.1，網(wǎng)格大小為0.01*0.01。在每一個(gè)網(wǎng)格點(diǎn)對(duì)信號(hào)載頻進(jìn)行估計(jì)。前面已經(jīng)得知，載頻估計(jì)的精度影響著輻射源估計(jì)的精度，最小二乘法估計(jì)得到的信號(hào)載頻沒有加權(quán)最小二乘法估計(jì)得到的載頻精確，二者載頻估計(jì)的誤差如圖2所示，本文提出的方法載頻估計(jì)誤差保持在-60 Hz到60 Hz之間，而最小二乘法估計(jì)的載頻誤差在-160 Hz到220 Hz之間，與前面的結(jié)論一致。

圖2 載頻估計(jì)誤差

定位結(jié)果如圖3所示，黑色的虛線表示衛(wèi)星在60 s內(nèi)的星下點(diǎn)軌跡，藍(lán)色三角表示地面站輻射源實(shí)際坐標(biāo)，其余6個(gè)點(diǎn)表示定位結(jié)果坐標(biāo)，其中3個(gè)是第一次網(wǎng)格尋優(yōu)時(shí)的定位結(jié)果，另外3個(gè)是第二次網(wǎng)格尋優(yōu)定位結(jié)果，根據(jù)先驗(yàn)信息，舍棄模糊值。由仿真顯示，對(duì)于實(shí)際輻射源經(jīng)緯度坐標(biāo)(118.5°,22.5°)，定位結(jié)果坐標(biāo)為(118.4915°,22.4826°)，定位誤差為2.449 km。

圖3 定位仿真結(jié)果

第一次網(wǎng)格尋優(yōu)時(shí)，每個(gè)網(wǎng)格點(diǎn)對(duì)應(yīng)1-8式中代價(jià)函數(shù)的值如圖4所示。由衛(wèi)星運(yùn)動(dòng)場(chǎng)景知道，星下點(diǎn)軌跡確定的網(wǎng)格坐標(biāo)，在定位中總會(huì)存在著兩個(gè)對(duì)稱的點(diǎn)使得代價(jià)函數(shù)的值達(dá)到最大(最小)，最優(yōu)值對(duì)應(yīng)的經(jīng)緯度即為定位點(diǎn)。

圖4 粗搜索代價(jià)函數(shù)仿真結(jié)果

由第一次代價(jià)函數(shù)尋優(yōu)得到的三個(gè)定位點(diǎn)經(jīng)緯度坐標(biāo)分別是[118.4896°,22.5255°]、[117.8896°,22.4255°]、[118.5896°,22.5255°]，分別在其鄰域進(jìn)行第二次代價(jià)函數(shù)尋優(yōu)的仿真結(jié)果如圖5所示，(a),(b),(c)分別代表第二次網(wǎng)格搜索后代價(jià)函數(shù)的三維圖，第二次代價(jià)函數(shù)尋優(yōu)得到的三個(gè)定位點(diǎn)經(jīng)緯度坐標(biāo)分別是[118.4972°,22.5013°]、[118.4968°,22.4862°]、[118.5014°,22.5130°]，比第一次搜索的定位結(jié)果更加準(zhǔn)確。

圖5 精搜索代價(jià)函數(shù)仿真結(jié)果

5 定位精度影響因素分析

5.1 信號(hào)載頻對(duì)定位精度的影響

為研究輻射源信號(hào)載頻對(duì)定位精度的影響，現(xiàn)設(shè)置輻射源位置，衛(wèi)星軌道高度等參數(shù)不變，分別設(shè)定信號(hào)載頻為1.3、2.6、3.4(單位：GHz)，其多普勒頻率和多普勒頻率變化率如圖6所示。

圖6 多普勒頻率及其變化率

用相同的網(wǎng)格尋優(yōu)對(duì)輻射源進(jìn)行定位，且在同一信號(hào)載頻時(shí)，解算時(shí)設(shè)置載頻均方誤差分別為10 Hz、30 Hz、60 Hz、120 Hz、180 Hz、240 Hz、360 Hz、480 Hz，定位結(jié)果如圖7所示。

圖7 信號(hào)載頻對(duì)定位結(jié)果的影響

通過(guò)圖6可以發(fā)現(xiàn)當(dāng)信號(hào)載頻越小時(shí)，衛(wèi)星與地面站的相對(duì)運(yùn)動(dòng)速度減小，多普勒頻率效應(yīng)減弱，且多普勒頻率變化率減小。通過(guò)圖7可以發(fā)現(xiàn)信號(hào)載頻越小，定位誤差越大。

5.2 接收信號(hào)時(shí)長(zhǎng)對(duì)定位精度的影響

為研究信號(hào)時(shí)長(zhǎng)對(duì)定位精度的影響，現(xiàn)設(shè)置輻射源位置，衛(wèi)星軌道高度等參數(shù)不變，分別設(shè)定信號(hào)時(shí)長(zhǎng)為20、40、60、80(單位：s)，該時(shí)長(zhǎng)指以星下點(diǎn)軌跡中心為對(duì)稱的時(shí)間范圍。不同信號(hào)時(shí)長(zhǎng)對(duì)應(yīng)的定位誤差如圖8所示。

圖8 接收信號(hào)時(shí)長(zhǎng)對(duì)定位結(jié)果的影響

通過(guò)圖8可以發(fā)現(xiàn)當(dāng)信號(hào)時(shí)長(zhǎng)越長(zhǎng)，定位誤差越小；相反，若信號(hào)時(shí)長(zhǎng)越短，定位誤差越大。由于信號(hào)時(shí)長(zhǎng)越長(zhǎng)，接收機(jī)能夠獲得頻率信息越多，最小二乘法估計(jì)的載頻更加準(zhǔn)確，因此能夠提高定位精度。

6 結(jié)束語(yǔ)

本文將輻射源定位問(wèn)題轉(zhuǎn)化為發(fā)射信號(hào)載頻估計(jì)問(wèn)題，提出一種基于加權(quán)最小二乘載頻估計(jì)的單星定位算法，通過(guò)理論推導(dǎo)和仿真驗(yàn)證，可以發(fā)現(xiàn)加權(quán)最小二乘法估計(jì)得到的載頻精度更高，從而輻射源實(shí)際的坐標(biāo)漂移更小，定位誤差更小。此外也分析了信號(hào)載頻、信號(hào)時(shí)長(zhǎng)等因素對(duì)定位結(jié)果的影響，發(fā)現(xiàn)軌道傾角對(duì)定位結(jié)果影響較小，軌道高度越低，發(fā)射信號(hào)載頻越大，接收信號(hào)時(shí)長(zhǎng)越長(zhǎng)，定位精度都有一定程度的提高，為低軌衛(wèi)星的軌道設(shè)計(jì)和后續(xù)算法改進(jìn)提供了新的思路。