趙建磊, 李海陽,*

(1. 國(guó)防科技大學(xué)空天科學(xué)學(xué)院, 湖南 長(zhǎng)沙 410073;2. 空天任務(wù)智能規(guī)劃與仿真湖南省重點(diǎn)實(shí)驗(yàn)室, 湖南 長(zhǎng)沙 410073)

0 引 言

隨著人類探索地球外層空間活動(dòng)的不斷深入,空間飛行器尤其是非合作空間飛行器的數(shù)量日益增長(zhǎng)。出于對(duì)有限空間資源的爭(zhēng)奪,惡意抵近服務(wù)如干擾、偵照、抓捕等行為逐漸增加,飛行器在軌安全受到嚴(yán)重威脅。通常情況下,飛行器根據(jù)實(shí)際任務(wù)需要對(duì)自身進(jìn)行的軌道機(jī)動(dòng)、軌道維持、抵近觀測(cè)、姿態(tài)控制等操作的信息是不公開的。在此背景下,非合作目標(biāo)的軌道機(jī)動(dòng)信息的感知、檢測(cè)和識(shí)別對(duì)及時(shí)發(fā)現(xiàn)非合作飛行器的惡意抵近服務(wù)意圖,維護(hù)在軌飛行器的安全顯得尤為重要。

國(guó)內(nèi)外學(xué)者從不同角度對(duì)軌道機(jī)動(dòng)檢測(cè)和機(jī)動(dòng)信息識(shí)別等問題進(jìn)行研究。Kelecy等人對(duì)脈沖小推力作用下的定軌與機(jī)動(dòng)檢測(cè)進(jìn)行研究,對(duì)比分析了最小二乘估計(jì)方法與卡爾曼濾波方法的適用性。Lemmens等人分別提出了基于兩行軌道要素(two line elements, TLE)一致性檢驗(yàn)的低軌衛(wèi)星機(jī)動(dòng)檢測(cè)方法和支持近實(shí)時(shí)機(jī)動(dòng)檢測(cè)的TLE時(shí)間序列分析方法。Pittelkau和Acquesta等人分別提出了機(jī)動(dòng)檢測(cè)中避免歷史數(shù)據(jù)中噪聲干擾的新方法,其中前者采用中值濾波器,而后者則采用被稱為Box具有噪聲的基于密度的聚類方法(density-based spatial clustering of applications with noise, DBSCAN)。Bai等人利用衛(wèi)星的歷史軌道信息采用聚類的方法進(jìn)行了機(jī)動(dòng)檢測(cè),分析了均值聚類、層次聚類和模糊聚類在機(jī)動(dòng)檢測(cè)中的應(yīng)用。Yu等人通過分析機(jī)動(dòng)前后的攝動(dòng),利用相對(duì)動(dòng)力學(xué)方程推導(dǎo)出兩種機(jī)動(dòng)模型,提出了基于軌道動(dòng)力學(xué)模型的機(jī)動(dòng)檢測(cè)方法。Daniel等人提出了一種新的基于馬氏距離的二元假設(shè)檢測(cè)方法,用于表征和檢測(cè)航天器異常。Bao等人提出了一種基于貝葉斯理論和多模型的濾波器,使用最小均方誤差方法來適應(yīng)機(jī)動(dòng)目標(biāo)可能進(jìn)行的多種運(yùn)動(dòng),以此進(jìn)行聯(lián)合檢測(cè)和弱目標(biāo)跟蹤。Clark等人通過將TLE分析和狀態(tài)傳播這兩種眾所周知的檢測(cè)方法融合成混合檢測(cè)方法來提高檢測(cè)的計(jì)算速度和精度提高,提高了空間態(tài)勢(shì)感知(space situational awareness, SSA)算法的性能。Justin等人利用光聲特征檢測(cè)方案實(shí)現(xiàn)了航天器近實(shí)時(shí)機(jī)動(dòng)檢測(cè)和參數(shù)估計(jì),增強(qiáng)了對(duì)小機(jī)動(dòng)事件的識(shí)別能力。Jiang等人探討了跟蹤執(zhí)行未知連續(xù)機(jī)動(dòng)的衛(wèi)星問題,提出一種估計(jì)衛(wèi)星狀態(tài)和機(jī)動(dòng)加速度的魯棒擴(kuò)展卡爾曼濾波器,與無偏最小方差輸入和狀態(tài)估計(jì)方法相比,該方法在具有相同的位置估計(jì)精度的同時(shí)機(jī)動(dòng)期間速度估計(jì)誤差降低了5倍左右。Christopher從目標(biāo)歷史TLE記錄中獲得的數(shù)據(jù)訓(xùn)練神經(jīng)網(wǎng)絡(luò),用于發(fā)現(xiàn)新的TLE記錄中的異常軌道機(jī)動(dòng)、錯(cuò)誤標(biāo)記或其他原因引起的意外偏差,該方法能夠避免依賴于個(gè)別目標(biāo)軌道狀態(tài)或特定的傳播模型,提高了檢測(cè)的準(zhǔn)確性。Li等人利用目標(biāo)的歷史TLE數(shù)據(jù),分別用統(tǒng)計(jì)的方法和自組織映射聚類的方法對(duì)LEO軌道衛(wèi)星的機(jī)動(dòng)檢測(cè)問題進(jìn)行了研究。崔紅正等人提出針對(duì)不同推力非合作目標(biāo)的機(jī)動(dòng)檢測(cè)方法,結(jié)合地基與天基觀測(cè)數(shù)據(jù),可以滿足多數(shù)非合作目標(biāo)軌道機(jī)動(dòng)檢測(cè)需求。王慶瑞等人將距離變化率作為特征量,提出一種基于概率判決模型的自適應(yīng)軌道機(jī)動(dòng)檢測(cè)方法,該方法基于Neyman-Pearson準(zhǔn)則求解序列數(shù)據(jù)的變軌判決門限。

以上方法中基于歷史軌道信息的機(jī)動(dòng)檢測(cè)方法只能判斷是否發(fā)生機(jī)動(dòng),不能對(duì)機(jī)動(dòng)參數(shù)進(jìn)行有效的定量識(shí)別,而對(duì)于機(jī)動(dòng)信息的識(shí)別,大多是基于天基或地基平臺(tái)對(duì)目標(biāo)的連續(xù)測(cè)角或測(cè)距追蹤,結(jié)合動(dòng)力學(xué)模型和濾波方法進(jìn)行機(jī)動(dòng)軌道定軌及機(jī)動(dòng)識(shí)別。但對(duì)于多數(shù)非合作目標(biāo),難以對(duì)其進(jìn)行連續(xù)追蹤觀測(cè)并且動(dòng)力學(xué)模型是未知的,只能通過稀疏短弧段觀測(cè)數(shù)據(jù)在無機(jī)動(dòng)假設(shè)下定軌得到若干個(gè)離散的軌道狀態(tài)。而當(dāng)衛(wèi)星發(fā)生軌道機(jī)動(dòng)時(shí),通過軌道測(cè)量得到的衛(wèi)星軌道根數(shù)相比于不發(fā)生軌道機(jī)動(dòng)時(shí)的軌道根數(shù),通常含有較大的誤差,且在短時(shí)間內(nèi)有限條軌道根數(shù)無法通過濾波等手段處理誤差。

針對(duì)以上問題,本文基于稀疏短弧段觀測(cè)數(shù)據(jù)定軌得到的有限條軌道測(cè)量狀態(tài),對(duì)非合作飛行器的機(jī)動(dòng)參數(shù)進(jìn)行識(shí)別。首先建立基于稀疏軌道根數(shù)的飛行器軌道反演模型,然后基于軌道反演對(duì)飛行器機(jī)動(dòng)參數(shù)進(jìn)行識(shí)別,最后通過算例對(duì)本文的方法進(jìn)行校驗(yàn)。

1 非合作飛行器軌道反演模型

1.1 基本假設(shè)

在本文的研究中,做以下假設(shè)。

(1) 當(dāng)變軌推力作用時(shí)間遠(yuǎn)遠(yuǎn)小于變軌前后的軌道周期時(shí),可以假設(shè)推力隨時(shí)間變化的函數(shù)是脈沖函數(shù),使得脈沖矢量與原推力產(chǎn)生的沖量一樣,因此本文中涉及的變軌沖量均考慮為脈沖沖量。

(2) 軌道測(cè)量?jī)H存在隨機(jī)誤差,不存在系統(tǒng)誤差。

(3) 以美國(guó)戰(zhàn)略司令部發(fā)布的TLE中軌道機(jī)動(dòng)后恢復(fù)時(shí)間為例,空間飛行器在發(fā)生機(jī)動(dòng)后,通過測(cè)定軌得到穩(wěn)定的TLE軌道通常需要2～7天,本文假設(shè)飛行器在得到穩(wěn)定的軌道之前不會(huì)施加下一次軌道機(jī)動(dòng)。

(4) 在非合作飛行器的一個(gè)軌道周期內(nèi),至少能得到一段短弧段觀測(cè)數(shù)據(jù),結(jié)合假設(shè)(3),即對(duì)于軌道周期小于2天的非合作飛行器,相鄰短弧段觀測(cè)數(shù)據(jù)之間至多考慮存在一次軌道機(jī)動(dòng)。

1.2 沖量作用下飛行器軌道預(yù)報(bào)模型

在二體假設(shè)下,天體為質(zhì)點(diǎn),空間飛行器只受中心引力作用。而實(shí)際飛行任務(wù)中,中心天體一般為非理想球體,因此飛行器還會(huì)受中心天體非球形引力、大氣阻力、太陽光壓、第三體引力等各種攝動(dòng)因素的影響。

在J2000地心慣性坐標(biāo)系中的飛行器動(dòng)力學(xué)模型如下所示:

(1)

當(dāng)飛行器的初始運(yùn)動(dòng)狀態(tài)可知,在脈沖機(jī)動(dòng)作用下,飛行器在任意時(shí)刻的軌道狀態(tài)表示如下:

(2)

1.3 飛行器軌道反演

(3)

定義軌道反演的狀態(tài)偏差為在所有測(cè)量點(diǎn)時(shí)刻(=0,1,…,),反演軌道狀態(tài)與測(cè)量狀態(tài)的加權(quán)偏差和,即

(4)

(5)

=DU/VU

(6)

(7)

式中: ,(=,,)為對(duì)應(yīng)于沿軌道不同方向的測(cè)量誤差的比例系數(shù)。

狀態(tài)偏差反映了通過軌道測(cè)量數(shù)據(jù)反演的軌道與測(cè)量軌道的接近程度,狀態(tài)偏差越小,說明反演的軌道與測(cè)量軌道更接近,同時(shí)與測(cè)量軌道的離散狀態(tài)相比,通過測(cè)量狀態(tài)反演的軌道是一條連續(xù)的軌道,減小了測(cè)量數(shù)據(jù)中的隨機(jī)誤差的影響,反演得到的軌道更加接近真實(shí)軌道。因此,可以通過最優(yōu)化方法使?fàn)顟B(tài)偏差最小,以此反演非合作飛行器的軌道。

2 基于軌道反演的機(jī)動(dòng)參數(shù)識(shí)別方法

2.1 機(jī)動(dòng)參數(shù)識(shí)別優(yōu)化模型

(1) 優(yōu)化變量

由式(3)和式(4)可知,狀態(tài)偏差與反演軌道的初始軌道狀態(tài)和中間機(jī)動(dòng)參數(shù)有關(guān),即

(8)

因此,將優(yōu)化變量選為反演軌道初始狀態(tài)和中間機(jī)動(dòng),即

(9)

(2) 優(yōu)化變量邊界約束

(10)

式中:,,,,,分別為飛行器軌道的半長(zhǎng)軸、偏心率、軌道傾角、升交點(diǎn)赤經(jīng)、近地點(diǎn)角距和真近點(diǎn)角;和分別為根據(jù)飛行器軌道類型設(shè)定的半長(zhǎng)軸邊界。

(3) 目標(biāo)函數(shù)

目標(biāo)函數(shù)取為狀態(tài)偏差最小,即

(11)

值得注意的是,由于在設(shè)計(jì)變量中引入的中間機(jī)動(dòng)通常是多于實(shí)際機(jī)動(dòng)的,多余的機(jī)動(dòng)在優(yōu)化結(jié)果中為小值,因此在優(yōu)化流程的最后進(jìn)行一次小機(jī)動(dòng)值的修正,即若存在設(shè)計(jì)變量|Δ|(=1,2,…,)的數(shù)值小于一定門限,認(rèn)為在該時(shí)刻未發(fā)生機(jī)動(dòng)。

(12)

2.2 “全局尋優(yōu)+局部修正”的雙重優(yōu)化求解

上述問題可以歸結(jié)為非線性規(guī)劃問題,關(guān)于這種非線性規(guī)劃問題的求解,有很多算法可以使用,如差分進(jìn)化(differential evolution, DE)算法、序列二次規(guī)劃(sequential quadratic programming, SQP)算法等。此處選取DE算法和SQP算法,其中DE算法全局搜索能力較好,但效率不高,尤其當(dāng)軌道外推使用高精度模型進(jìn)行時(shí),會(huì)花費(fèi)大量時(shí)間;序列二次規(guī)劃算法局部尋優(yōu)效果好,但對(duì)初值比較敏感。因此,在使用高精度軌道外推模型求解時(shí),為了提高搜索效率以及避免初值陷入局部最優(yōu),提出一種“全局尋優(yōu)+局部修正”的雙重優(yōu)化算法進(jìn)行計(jì)算,計(jì)算流程如圖1所示。

圖1 雙重優(yōu)化求解流程

在二體條件下,通過DE算法尋找滿足約束的全局最優(yōu)解;

以步驟1中最優(yōu)解為初值,利用SQP算在高精度模型條件下,尋找滿足約束的最優(yōu)解;

剔除小機(jī)動(dòng)值,重新進(jìn)行一遍SQP優(yōu)化流程。

3 仿真校驗(yàn)

本節(jié)對(duì)上述的飛行器軌道反演模型及機(jī)動(dòng)參數(shù)識(shí)別方法進(jìn)行仿真校驗(yàn)。仿真在計(jì)算機(jī)處理器為Intel i7-7700 CPU 3.60 GHz,運(yùn)行內(nèi)存為16 GB的C++環(huán)境下進(jìn)行。

3.1 仿真參數(shù)設(shè)定

地球同步靜止軌道由于具有獨(dú)特的高軌和地球同步特性,成為了通信、氣象、導(dǎo)彈預(yù)警、廣播電視和數(shù)據(jù)中繼等重要衛(wèi)星的首選軌道,但同樣也將面臨嚴(yán)峻的惡意服務(wù)考驗(yàn)。因此,對(duì)地球同步靜止軌道衛(wèi)星機(jī)動(dòng)信息的感知、檢測(cè)和識(shí)別顯得尤為重要。本節(jié)以地球同步靜止軌道為例,對(duì)本文提出的模型和方法進(jìn)行校驗(yàn)。

表1 仿真初始軌道根數(shù)

表2 標(biāo)稱軌道中間機(jī)動(dòng)參數(shù)

對(duì)上述標(biāo)稱軌道每隔7 h選取該時(shí)刻的軌道狀態(tài)作為測(cè)量樣本,并施加不同的均值為0的高斯隨機(jī)定軌測(cè)量誤差,以此模擬實(shí)際情況下稀疏短弧段定軌得到的有限條離散軌道測(cè)量數(shù)據(jù)。不同工況下的定軌誤差如表3所示。

表3 定軌誤差標(biāo)準(zhǔn)差

3.2 算例分析

使用本文描述的機(jī)動(dòng)參數(shù)識(shí)別方法,不同測(cè)量點(diǎn)時(shí)刻對(duì)應(yīng)狀態(tài)偏差的權(quán)重因子均取為1,優(yōu)化目標(biāo)函數(shù)取為式(7),目標(biāo)函數(shù)中沿軌道不同方向的測(cè)量誤差的比例系數(shù), (=,,)取值如下所示：

(13)

其中，取為015 m/s,DE算法種群規(guī)模為240,代數(shù)為300,SQP算法最大迭代代數(shù)為50。

每種工況分別進(jìn)行50次蒙特卡羅仿真,機(jī)動(dòng)參數(shù)識(shí)別結(jié)果與標(biāo)稱機(jī)動(dòng)參數(shù)的偏差如表4所示。從結(jié)果來看,定軌誤差越大,則機(jī)動(dòng)參數(shù)的識(shí)別誤差越大。

表4 機(jī)動(dòng)參數(shù)識(shí)別結(jié)果

為進(jìn)一步說明雙重優(yōu)化算法的收斂情況及機(jī)動(dòng)參數(shù)識(shí)別效果,選取工況2中的一次仿真,初始軌道測(cè)量數(shù)據(jù)如表5所示。雙重優(yōu)化算法的迭代收斂情況如圖2所示,仿真顯示雙重優(yōu)化算法中DE算法在100代左右基本收斂,SQP算法在20代左右基本收斂。

圖2 雙重優(yōu)化算法收斂情況

表5 離散軌道測(cè)量數(shù)據(jù)

在表5中，定軌結(jié)果以經(jīng)典軌道根數(shù)的形式給出,括號(hào)中數(shù)據(jù)依次代表半長(zhǎng)軸、偏心率、軌道傾角、升交點(diǎn)赤經(jīng)、近地點(diǎn)角距、真近點(diǎn)角,其中半長(zhǎng)軸的單位為m,角度的單位為(°)。

同時(shí),為比較說明本文優(yōu)化算法的優(yōu)勢(shì),分別用DE算法和SQP算法對(duì)算例進(jìn)行求解。表6給出了不同算法求解的最大機(jī)動(dòng)識(shí)別誤差和計(jì)算時(shí)間,比較可知,SQP算法機(jī)動(dòng)識(shí)別誤差較大,原因是SQP算法對(duì)初值較為敏感,隨機(jī)給定初值容易使算法陷入局部收斂;DE算法計(jì)算時(shí)間較長(zhǎng),原因在于直接使用高精度軌道外推模型進(jìn)行全局搜索時(shí),由于計(jì)算量驟增,導(dǎo)致計(jì)算效率極低。相對(duì)而言,本文提出的“全局尋優(yōu)+局部修正”的雙重優(yōu)化算法對(duì)本文的模型適應(yīng)性很好。

表6 不同算法機(jī)動(dòng)識(shí)別誤差及計(jì)算時(shí)間

經(jīng)SQP高精度修正后的飛行器機(jī)動(dòng)參數(shù)反演結(jié)果如表7和表8所示,分別給出了剔除小機(jī)動(dòng)前后的機(jī)動(dòng)參數(shù)識(shí)別情況?？梢钥闯?剔除小機(jī)動(dòng)后的識(shí)別結(jié)果與不剔除小機(jī)動(dòng)相比結(jié)果更加接近標(biāo)稱軌道施加的機(jī)動(dòng),機(jī)動(dòng)識(shí)別的速度誤差在0.08 m/s以內(nèi),時(shí)間誤差在3 min以內(nèi)。

表7 飛行器機(jī)動(dòng)參數(shù)識(shí)別結(jié)果

表8 飛行器機(jī)動(dòng)參數(shù)識(shí)別結(jié)果(剔除小機(jī)動(dòng))

為進(jìn)一步說明反演軌道與實(shí)際軌道的偏差,圖3給出了在仿真時(shí)段反演得到的軌道與標(biāo)稱軌道的位置偏差及軌道跡向、法向和徑向偏差。可以看出,在整個(gè)仿真時(shí)段,反演軌道與標(biāo)稱軌道的最大位置偏差不超過1 200 m。

圖3 反演軌道與標(biāo)稱軌道的偏差

4 結(jié) 論

本文在稀疏軌道信息條件下,針對(duì)非合作飛行器軌道機(jī)動(dòng)參數(shù)的識(shí)別問題,探討了一種可行的思路。建立了飛行器軌道反演模型及基于軌道反演的機(jī)動(dòng)參數(shù)識(shí)別方法,并提出相應(yīng)的“全局尋優(yōu)+局部修正”的雙重優(yōu)化求解算法。仿真結(jié)果顯示,該方法能夠在稀疏軌道信息的條件下,對(duì)非合作飛行器的機(jī)動(dòng)參數(shù)進(jìn)行較為有效的識(shí)別;空間飛行器定軌誤差越大,則機(jī)動(dòng)參數(shù)的識(shí)別誤差越大;工況2中算例機(jī)動(dòng)參數(shù)識(shí)別結(jié)果顯示,剔除小機(jī)動(dòng)后的識(shí)別結(jié)果與不剔除小機(jī)動(dòng)相比結(jié)果偏差更小,反演軌道與標(biāo)稱軌道的位置偏差不超過1 200 m,機(jī)動(dòng)識(shí)別的速度誤差在0.08 m/s以內(nèi),時(shí)間誤差在3 min以內(nèi)。