黃曉斌， 張 燕， 肖 銳

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

0 引言

空間目標(biāo)軌道確定是空間目標(biāo)監(jiān)視雷達(dá)的一項(xiàng)重要任務(wù)，主要的數(shù)據(jù)處理步驟（如圖 1 所示）是：①坐標(biāo)轉(zhuǎn)換；②初軌確定；③軌道改進(jìn)；④編目。在數(shù)據(jù)處理過(guò)程中，涉及到大量的天文學(xué)和數(shù)學(xué)知識(shí)。雖然國(guó)內(nèi)外已有一些開(kāi)源定軌軟件可供借鑒，如 OrbFit[1]、ORSA[2]、NEOPROP[3]、Novas[4]、ODTBX[5]等，但這些軟件都是針對(duì)天體軌道問(wèn)題而研發(fā)的，且主要是處理光學(xué)觀測(cè)數(shù)據(jù)，并不適合處理雷達(dá)觀測(cè)數(shù)據(jù)。為此，本文開(kāi)發(fā)了 1 套雷達(dá)定軌函數(shù)庫(kù)（RadarOrbDet），能夠有助于雷達(dá)技戰(zhàn)指標(biāo)的快速實(shí)現(xiàn)。RadarOrbDet 庫(kù)按照上述4 個(gè)處理步驟開(kāi)發(fā)，本文主要介紹其坐標(biāo)轉(zhuǎn)換模塊。

本文通過(guò)討論雷達(dá)定軌中涉及的基本坐標(biāo)系、坐標(biāo)轉(zhuǎn)換的數(shù)學(xué)原理，給出程序設(shè)計(jì)思路，最后利用衛(wèi)星工具包（satellite tool kit,STK）軟件驗(yàn)證開(kāi)發(fā)模塊的有效性。

圖1 雷達(dá)定軌數(shù)據(jù)處理流程

1 坐標(biāo)系簡(jiǎn)介

雷達(dá)觀測(cè)基于地球坐標(biāo)系，空間目標(biāo)軌道是基于天球坐標(biāo)系，這就涉及到地球坐標(biāo)系與天球坐標(biāo)系之間的轉(zhuǎn)換。

協(xié)議天球坐標(biāo)系由國(guó)際天文學(xué)聯(lián)合會(huì)（The International Astronomical Union, IAU）和國(guó)際地球自轉(zhuǎn)和參考系服務(wù)組織（The International Earth Rotation and Reference Systems Service, IERS）發(fā)布，目前采用的是國(guó)際天球參考系（international celestial reference frame, ICRS）。依據(jù)坐標(biāo)原點(diǎn)的不同，ICRS 可分為太陽(yáng)系質(zhì)心天球參考系（barycentric celestial reference system, BCRS）和地球質(zhì)心天球參考系（geocentric celestial reference system, GCRS）。BCRS 用于計(jì)算行星的運(yùn)動(dòng)軌道及編制星表；GCRS用于計(jì)算衛(wèi)星軌道及編制衛(wèi)星星歷。ICRS 由國(guó)際天球參考框架（international celestial reference frame,ICRF）來(lái)實(shí)現(xiàn)。在1997 年IAU 第23 屆大會(huì)上，通過(guò)并決定自1998-01-01 起，在天文研究、空間探測(cè)、大地測(cè)量以及地球動(dòng)力學(xué)等領(lǐng)域中采用ICRS[6]。

協(xié)議地球坐標(biāo)系由國(guó)際地球參考系（international terrestrial reference system, ITRS）實(shí)現(xiàn)。GCRS 是1 個(gè)相當(dāng)好的準(zhǔn)慣性系，衛(wèi)星的軌道計(jì)算一般都是在 GCRS 中進(jìn)行。這就必須涉及到GCRS 與ITRS 間的坐標(biāo)轉(zhuǎn)換問(wèn)題[7]。

2 坐標(biāo)轉(zhuǎn)換原理

在雷達(dá)定軌中，需要將在站心地平坐標(biāo)系下的雷達(dá)觀測(cè)數(shù)據(jù)轉(zhuǎn)換到GCRS 坐標(biāo)系中。

2.1 站心地平坐標(biāo)系Xh與ITRS 坐標(biāo)系XGO的轉(zhuǎn)換

站心地平坐標(biāo)系Xh與 ITRS 坐標(biāo)系XGO的定義如表1 所示。

表1 站心地平坐標(biāo)系與ITRS 坐標(biāo)系的定義

站心地平坐標(biāo)系與ITRS 坐標(biāo)系之間的轉(zhuǎn)換公式為

2.2 ITRS 與GCRS 坐標(biāo)系的轉(zhuǎn)換

ITRS 與GCRS 的轉(zhuǎn)換早期是基于春分點(diǎn)的。目前，IERS 2010 年建議使用基于無(wú)旋轉(zhuǎn)原點(diǎn)（nonrotating origin, NRO）的轉(zhuǎn)換方法[9]?；?IAU 2006/2000A-CIO 模型的轉(zhuǎn)換流程[9]如圖2 所示。

圖2 “IAU 2006/2000A-CIO based”坐標(biāo)轉(zhuǎn)換流程

轉(zhuǎn)換過(guò)程中涉及到2 個(gè)中間坐標(biāo)系：地球中間坐標(biāo)系（terrestrial intermediate reference system, TIRS）和天球中間坐標(biāo)系（celestial intermediate reference system, TIRS）它們的定義見(jiàn)表2。

表2 TIRS 坐標(biāo)系與CIRS 坐標(biāo)系的定義

表 2 中：TIP（terrestrial intermediate pole）為地球瞬時(shí)自轉(zhuǎn)軸；TIO（terrestrial intermediate origin）為地球中間零點(diǎn)；CIP（celestial intermediate pole）為天球瞬時(shí)自轉(zhuǎn)軸；CIO（celestial intermediate origin）為天球中間零點(diǎn)。

在t時(shí)刻，ITRS 和 GCRS 的轉(zhuǎn)換是 2 個(gè) 3 維直角坐標(biāo)系間的轉(zhuǎn)換，可以寫(xiě)成

式中：M(t)、RCIO(t)和W(t)分別是由于 CIP 在 GCRS中的運(yùn)動(dòng)（歲差章動(dòng)）、地球的自轉(zhuǎn)以及 CIP 在ITRS 中的運(yùn)動(dòng)（極移）引起的旋轉(zhuǎn)矩陣。它們的具體表達(dá)式見(jiàn)文獻(xiàn)[9]。

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

圖3 給出了坐標(biāo)轉(zhuǎn)換的程序設(shè)計(jì)流程圖，其轉(zhuǎn)換過(guò)程分為 3 個(gè)部分：①將目標(biāo)在測(cè)站地平坐標(biāo)系的極坐標(biāo)轉(zhuǎn)換為ITRS 坐標(biāo)系下的直角坐標(biāo)；②根據(jù)轉(zhuǎn)換時(shí)刻以及地球定向參數(shù)（Earth orientation parameters，EOP）數(shù)據(jù)形成 ITRS 到 GCRS 的轉(zhuǎn)換矩陣；③將轉(zhuǎn)換矩陣與ITRS 下的直角坐標(biāo)相乘得到GCRS 下的直角坐標(biāo)。轉(zhuǎn)換流程中涉及一些時(shí)間系統(tǒng)的轉(zhuǎn)換可參見(jiàn)文獻(xiàn)[7]。圖3 中的“TT 時(shí)間”表示地球時(shí)（Terrestrial Time）。

圖3 坐標(biāo)轉(zhuǎn)換程序設(shè)計(jì)流程

在轉(zhuǎn)換過(guò)程中有如下幾個(gè)問(wèn)題應(yīng)該注意：

1）按照式（1）的要求，應(yīng)該輸入站點(diǎn)的天文經(jīng)緯度和地理經(jīng)緯度。前者用來(lái)進(jìn)行坐標(biāo)旋轉(zhuǎn)，后者用來(lái)計(jì)算站點(diǎn)的大地直角坐標(biāo)。但通常只給出站點(diǎn)的地理經(jīng)緯度，在進(jìn)行坐標(biāo)旋轉(zhuǎn)時(shí)，用地理經(jīng)緯度替代天文經(jīng)緯度，2 者雖有細(xì)微的差別，但在絕大多數(shù)場(chǎng)合可忽略不計(jì)，如果需要極高坐標(biāo)轉(zhuǎn)換精度時(shí)應(yīng)加以區(qū)分。

2）轉(zhuǎn)換過(guò)程中需要地球定向參數(shù)（Earth orientation parameters, EOP）數(shù)據(jù)，它可從網(wǎng)站下載獲得[10]，要注意表中數(shù)據(jù)的有效時(shí)段和數(shù)據(jù)條目的時(shí)間間隔（通常為1 d）。因此，在給定轉(zhuǎn)換時(shí)刻時(shí)，通常需要依據(jù)該表插值來(lái)獲取更為精確的數(shù)據(jù)，而當(dāng)轉(zhuǎn)換時(shí)刻超出EOP 數(shù)據(jù)表有效時(shí)間范圍時(shí)，一般可取最接近該轉(zhuǎn)換時(shí)刻的數(shù)據(jù)來(lái)替代，但當(dāng)超出時(shí)間范圍間隔過(guò)大時(shí)，其誤差應(yīng)加以考慮。

3）當(dāng)需要極高坐標(biāo)轉(zhuǎn)換精度時(shí)，在用EOP 數(shù)據(jù)表插值計(jì)算極移和世界協(xié)調(diào)時(shí)（coordinated universal time，UTC）與 1 類世界時(shí)（universal time，UT1）參數(shù)時(shí)，需要額外考慮海潮和固體潮的修正[7]，但同時(shí)會(huì)導(dǎo)致轉(zhuǎn)換速度的下降。

4 軟件接口調(diào)用及有效性驗(yàn)證

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

RadarOrbDet 庫(kù)坐標(biāo)轉(zhuǎn)換功能的使用方法如下：

步驟①：打開(kāi)VS2015，新建控制臺(tái)應(yīng)用程序，并命名為 radarorddet_Test，輸出模式設(shè)為 64 位Debug，具體方法可參考文獻(xiàn)[11]。

步驟②：從百度網(wǎng)盤(pán)下載RadarOrbDet 開(kāi)發(fā)包[12]，并拷貝到工程目錄中，如圖4 所示。

圖4 測(cè)試工程目錄結(jié)構(gòu)

步驟③：在radarorbdet_Test.cpp 文件中輸入如圖5 所示的代碼。

圖5 函數(shù)接口調(diào)用示例

4.2 接口函數(shù)有效性驗(yàn)證

在 STK 11.2 版中仿真 1 顆衛(wèi)星，其軌道參數(shù)設(shè)置如圖 6 所示，測(cè)站位置參數(shù)設(shè)置如圖 7所示。

圖6 衛(wèi)星軌道參數(shù)設(shè)置

圖7 測(cè)站位置參數(shù)設(shè)置

仿真時(shí)間段設(shè)在2016-09-15 UTCG（格林尼治協(xié)調(diào)世界時(shí)） 04:00:00 至 2016-09-16 UTCG 04:00:00（STK 沒(méi)有最新的EOP 數(shù)據(jù)，故設(shè)此時(shí)間段保證STK 自帶的EOP 數(shù)據(jù)有效）。利用STK 的報(bào)表功能輸出測(cè)站對(duì)衛(wèi)星的觀測(cè)時(shí)刻、觀測(cè)距離、方位和俯仰數(shù)據(jù)，同時(shí)輸出STK 內(nèi)部計(jì)算所得的衛(wèi)星GCRS 坐標(biāo)，坐標(biāo)數(shù)據(jù)的采樣時(shí)間段為2016-09-15 UTCG04:00:00 至 2016-09-15UTCG 04:02:00，間隔1 s，總共121 點(diǎn)數(shù)據(jù)，圖8(a)顯示了距離測(cè)量數(shù)據(jù)，圖8(b)顯示了方位測(cè)量數(shù)據(jù)，圖8(c)顯示了俯仰測(cè)量數(shù)據(jù)。

利用 RadarOrbDet 軟件包計(jì)算 GCRS 坐標(biāo)（XROD,YROD,ZROD），采用式（3）描述的 GCRS 坐標(biāo)相對(duì)誤差與STK 的輸出進(jìn)行比較，圖9 顯示了相對(duì)誤差曲線（實(shí)線）。從圖9 中可以看出相對(duì)誤差在 1×10-7水平，經(jīng)分析表明，誤差的主要原因是STK 報(bào)表以文本方式輸出數(shù)據(jù)，導(dǎo)致數(shù)據(jù)產(chǎn)生精度截?cái)嗾`差，由此可知軟件包算法是有效的。此外，在圖9 中還繪制了俯仰角曲線（虛線），可以看出當(dāng)俯仰角越小時(shí)，轉(zhuǎn)換的誤差越大，其原因是俯仰角越小時(shí)，截?cái)嗾`差對(duì)坐標(biāo)轉(zhuǎn)換誤差的影響越大。

圖8 測(cè)站觀測(cè)的距離、方位和俯仰數(shù)據(jù)

圖9 GCRS 坐標(biāo)轉(zhuǎn)換相對(duì)誤差

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

本文介紹了 RadarOrbDet 庫(kù)中坐標(biāo)轉(zhuǎn)換模塊的基本數(shù)學(xué)原理，給出了程序設(shè)計(jì)思路，描述了調(diào)用方法，最后用STK 驗(yàn)證了算法的有效性。下一步將繼續(xù)介紹初軌確定模塊的功能。