黃曉斌， 張 燕， 魯 力， 肖 銳

(空軍預(yù)警學(xué)院， 湖北武漢 430019)

0 引言

雷達(dá)的優(yōu)點(diǎn)是白天黑夜均能探測(cè)遠(yuǎn)距離的目標(biāo)，且不受霧、云和雨的阻擋，具有全天候、全天時(shí)的特點(diǎn)，并有一定的穿透能力。雷達(dá)已廣泛應(yīng)用于社會(huì)經(jīng)濟(jì)發(fā)展(如氣象預(yù)報(bào)、資源探測(cè)、環(huán)境監(jiān)測(cè)等)和科學(xué)研究(天體研究、大氣物理、電離層結(jié)構(gòu)研究等)中，但其應(yīng)用最廣的還是軍事領(lǐng)域，如防空警戒、反導(dǎo)預(yù)警、空間監(jiān)視、指揮引導(dǎo)等。

空間目標(biāo)監(jiān)視雷達(dá)是空間監(jiān)視領(lǐng)域的骨干裝備，與常規(guī)防空雷達(dá)相比，主要增加了空間目標(biāo)軌道確定功能，軌道改進(jìn)實(shí)現(xiàn)該功能的核心模塊，軌道改進(jìn)模塊的研制涉及非常復(fù)雜的航天動(dòng)力學(xué)知識(shí)，這對(duì)從事雷達(dá)系統(tǒng)設(shè)計(jì)的非軌道力學(xué)專業(yè)的研究人員來說帶來極大困難。雖有一些開源定軌軟件可供借鑒，但這些軟件主要用于處理光學(xué)觀測(cè)數(shù)據(jù)，并不適合處理雷達(dá)觀測(cè)數(shù)據(jù)。為此，作者著手開發(fā)RadarOrbDet函數(shù)庫來輔助非軌道力學(xué)專業(yè)的人員進(jìn)行空間目標(biāo)監(jiān)視雷達(dá)的系統(tǒng)設(shè)計(jì)，目前已完成坐標(biāo)轉(zhuǎn)換模塊、初軌確定模塊和攝動(dòng)分析模塊的開發(fā)，本文研究軌道改進(jìn)模塊的開發(fā)。

軌道改進(jìn)是空間目標(biāo)軌道確定的重點(diǎn)和難點(diǎn)，它涉及的內(nèi)容較多，包括攝動(dòng)分析、偏微分方程的變分法分析、最小二乘優(yōu)化等理論，本文聚焦空間目標(biāo)雷達(dá)定軌中的批處理軌道改進(jìn)問題，著重介紹了最小二乘軌道改進(jìn)原理，講解了模型中的偏導(dǎo)數(shù)求解方法，描述了程序?qū)崿F(xiàn)流程，最后用ODTK軟件驗(yàn)證了本文軟件設(shè)計(jì)的有效性。

1 最小二乘軌道改進(jìn)原理

最小二乘軌道確定的基本原理是找到一條軌道，使得理論觀測(cè)值和實(shí)際觀測(cè)值之間的殘差平方和(即損耗函數(shù))最小，即求解，使得式(1)的值最小。

()=[-()][-()]

(1)

式中:=[();()]是衛(wèi)星在時(shí)刻的位置和速度構(gòu)成的狀態(tài)參量，由通過軌道預(yù)報(bào)可以得到一條參考軌道；=[;;…;]是次雷達(dá)實(shí)際觀測(cè)矢量，每次觀測(cè)矢量包含距離、方位和俯仰三個(gè)分量(,,)；()是由參考軌道計(jì)算得到的理論觀測(cè)值。

由于是一個(gè)非線性函數(shù)，式(1)描述的是一個(gè)非線性最小二乘問題，求解異常困難，通常采用泰勒展開將其轉(zhuǎn)化為一個(gè)線性最小二乘問題，然后通過迭代方式求解它的線性最小二乘近似解，表達(dá)式為

(2)

(3)

2 模型中的偏導(dǎo)數(shù)求解

式(2)的求解關(guān)鍵在于偏導(dǎo)數(shù)雅克比矩陣的計(jì)算，不失一般性，以某一時(shí)刻的雷達(dá)觀測(cè)量為例。為書寫方便，除時(shí)刻外，省略其他時(shí)刻索引。

利用求導(dǎo)鏈?zhǔn)椒▌t，轉(zhuǎn)換為

(4)

式中，矩陣是當(dāng)前時(shí)刻觀測(cè)量對(duì)當(dāng)前狀態(tài)參量的偏導(dǎo)數(shù)矩陣，是當(dāng)前時(shí)刻狀態(tài)參量對(duì)初始時(shí)刻狀態(tài)參量的偏導(dǎo)數(shù)矩陣，也稱為狀態(tài)誤差傳遞矩陣。

2.1 觀測(cè)偏導(dǎo)數(shù)矩陣求解

在忽略光行差等因素后，雷達(dá)觀測(cè)量只與衛(wèi)星當(dāng)前的位置有關(guān)，而與其速度無關(guān)。因此，矩陣可拆分成

(5)

式中，是衛(wèi)星在慣性坐標(biāo)系(ECI坐標(biāo)系)下的位置矢量。

雷達(dá)觀測(cè)是基于站心地平坐標(biāo)進(jìn)行的(如南-東-天坐標(biāo)系，SEZ)，直接計(jì)算觀測(cè)量對(duì)ECI坐標(biāo)系下的位置矢量的偏導(dǎo)數(shù)比較困難，因此，矩陣的計(jì)算需進(jìn)一步運(yùn)用鏈?zhǔn)角髮?dǎo)法則進(jìn)行分解：

(6)

式中，是衛(wèi)星在SEZ坐標(biāo)系下的直角坐標(biāo)，用(,,)表示；是雷達(dá)對(duì)衛(wèi)星的觀測(cè)矢量在地固坐標(biāo)系(ECEF坐標(biāo)系)下的表示；是衛(wèi)星在ECEF坐標(biāo)系的位置矢量。

注意到

=-

(7)

式中是雷達(dá)在ECEF坐標(biāo)系下的位置矢量，是個(gè)常矢量。

因此

(8)

根據(jù)坐標(biāo)變換有

(9)

根據(jù)式(9)有

(10)

同理有

(11)

將式(8)、(10)和(11)代入式(6)有

(12)

在坐標(biāo)系下，由直角坐標(biāo)和極坐標(biāo)的關(guān)系，易得

(13)

將式(13)代入式(12)，有

(14)

計(jì)算式(14)中的兩個(gè)矩陣時(shí)，都需要在當(dāng)前觀測(cè)時(shí)刻進(jìn)行。最后，將式(14)代入式(5)后，可得矩陣。

2.2 狀態(tài)誤差傳遞矩陣求解

設(shè)衛(wèi)星的運(yùn)動(dòng)方程為

(15)

則

(16)

其中，矩陣為

(17)

由于采用迭代方式進(jìn)行最小二乘問題求解，因此對(duì)矩陣的求解精度要求不高，從而對(duì)矩陣的求解精度也要求不高。因此，矩陣中的加速度可以只考慮地球中心引力和低階的非球形引力攝動(dòng)，以加快計(jì)算機(jī)求解速度。

在只考慮地球引力(包含非球形攝動(dòng))情況下，加速度只與位置有關(guān)，而與速度無關(guān)，并且它的計(jì)算通常是在ECEF坐標(biāo)系下進(jìn)行，具體表達(dá)式可以參考文獻(xiàn)[18]，篇幅所限，這里不再贅述。

最后還需將ECEF坐標(biāo)系下的偏導(dǎo)數(shù)矩陣轉(zhuǎn)換到ECI坐標(biāo)系下。在忽略科里奧利力和離心力情況下，有如下關(guān)系式：

(18)

計(jì)算得到當(dāng)前時(shí)刻的矩陣后，采用數(shù)值微分方程方法求解。

3 程序設(shè)計(jì)

圖1是最小二乘軌道改進(jìn)算法的程序設(shè)計(jì)流程圖。

圖1 最小二乘軌道改進(jìn)算法的程序設(shè)計(jì)流程圖

程序分為三層循環(huán)，第一層循環(huán)是數(shù)值微分方程求解的內(nèi)部循環(huán)，輸出結(jié)果是第時(shí)刻的參考軌道和狀態(tài)誤差傳遞矩陣；第二層循環(huán)是遍歷每個(gè)觀測(cè)值，根據(jù)最小二乘原理求解單次軌道改進(jìn)量；第三層循環(huán)是執(zhí)行多次最小二乘軌道改進(jìn)，使最終結(jié)果趨近于實(shí)際軌道狀態(tài)。

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

4.1 接口調(diào)用示例

采用C語言按照第3部分的程序設(shè)計(jì)原理實(shí)現(xiàn)了最小二乘軌道改進(jìn)的接口函數(shù)rodLeastSquare，并集成到RadarOrbDet庫中。rodAccelHarmonic接口函數(shù)的軟件調(diào)用方法如圖2所示。

圖2 最小二乘軌道改進(jìn)接口函數(shù)說明

4.2 實(shí)驗(yàn)比較

在STK11.2中仿真一顆衛(wèi)星，其歷元軌道參數(shù)設(shè)置如表1所示，運(yùn)動(dòng)模式設(shè)置為HPOP(高精度軌道預(yù)報(bào))，考慮21階非球形攝動(dòng)。雷達(dá)站址的經(jīng)緯高為(120°,30°,0)。選擇2021-09-07 19:15:01至2021-09-07 19:15:10間隔1秒的10點(diǎn)數(shù)據(jù)，對(duì)距離、方位和俯仰分別加上100 m、01°和01°的觀測(cè)誤差。

表1 衛(wèi)星軌道根數(shù)

采用本文軟件和AGI公司開發(fā)的一款先進(jìn)的軌道確定與分析軟件——Orbit Determination Tool Kit(ODTK)進(jìn)行定軌比較，ODTK與本文軟件計(jì)算結(jié)果如表2所示，采用式(19)所描述的位置和速度相對(duì)誤差進(jìn)行定量比較。

(19)

從表2可以看出，定軌的位置誤差為0.059 6%，速度誤差為0.35%，表明本文軟件算法是有效的，其誤差主要來源為計(jì)算程序的坐標(biāo)與時(shí)間轉(zhuǎn)換、計(jì)算截?cái)嗾`差等。

表2 定軌結(jié)果對(duì)比表(ODTK與本文軟件)

5 結(jié)束語

本文針對(duì)空間目標(biāo)雷達(dá)定軌中的批處理軌道改進(jìn)問題，介紹了最小二乘軌道改進(jìn)的原理，描述了模型中的偏導(dǎo)數(shù)求解，給出了程序設(shè)計(jì)思路，仿真表明了軟件設(shè)計(jì)的有效性。此外，本文開發(fā)的RadarOrbDet函數(shù)庫可以很好地輔助空間目標(biāo)監(jiān)視雷達(dá)系統(tǒng)設(shè)計(jì)人員進(jìn)行定軌指標(biāo)的快速設(shè)計(jì)和驗(yàn)證。