慣性/天文組合導(dǎo)航中的時空坐標(biāo)轉(zhuǎn)換及其誤差分析

嚴(yán)恭敏，戴晨杰，陳若彤

（西北工業(yè)大學(xué) 自動化學(xué)院，西安 710072）

慣性/天文組合導(dǎo)航是一種重要的高精度自主導(dǎo)航方法，可應(yīng)用于艦船、飛機(jī)、導(dǎo)彈和衛(wèi)星等多種運(yùn)載平臺[1-4]。天文導(dǎo)航基于星體的光電信號進(jìn)行姿態(tài)測量，其絕對測量精度高，可達(dá)3～5″。天文導(dǎo)航輔助能有效抑制慣性導(dǎo)航的陀螺漂移，從而減緩慣導(dǎo)姿態(tài)、方位和定位誤差的積累。比如，應(yīng)用于美國B-2“幽靈”隱形轟炸機(jī)上的NAS-26 慣性/天文導(dǎo)航系統(tǒng)，其純慣性導(dǎo)航定位誤差為0.5 n mile/h，而若采取組合導(dǎo)航方式在10 小時內(nèi)的定位誤差可小于0.2 n mile[5,6]。

在慣性/天文組合導(dǎo)航中，慣性導(dǎo)航的姿態(tài)、速度和定位參數(shù)一般是在地球坐標(biāo)系（比如WGS84 和ITRF2000，簡記為e 系，前者主要應(yīng)用于導(dǎo)航領(lǐng)域而后者應(yīng)用于大地測量，在導(dǎo)航中可忽略兩者之間的微小差異）上描述和求解的，而天文導(dǎo)航的星體信號觀測是在天球坐標(biāo)系（比如J2000，簡記為i 系）上表示的，要對這兩種分導(dǎo)航系統(tǒng)進(jìn)行信息比對和融合，就必然涉及到導(dǎo)航信息在兩類坐標(biāo)系之間轉(zhuǎn)換的問題。J2000 協(xié)議天球坐標(biāo)系（協(xié)議慣性坐標(biāo)系）的X軸和Z軸分別指向J2000.0 時刻的平春分點(diǎn)與平天極，Y軸與X、Z構(gòu)成右手直角坐標(biāo)系，其中J2000.0 時刻為地球時（Terrestrial Time,TT）的2000年1月1日12:00 時[7]。在導(dǎo)航應(yīng)用中可取J2000 的原點(diǎn)為地球質(zhì)心。在導(dǎo)航專業(yè)中大家都比較熟悉WGS84 的定義，這里就不再具體給出。進(jìn)行組合導(dǎo)航信息轉(zhuǎn)換時，常常是將天文導(dǎo)航的姿態(tài)測量值從i 系轉(zhuǎn)換至e 系，若借助基礎(chǔ)天文標(biāo)準(zhǔn)庫（Standards Of Fundamental Astronomy,SOFA）可很容易實現(xiàn)該轉(zhuǎn)換[8]。

SOFA 是國際天文聯(lián)合會（International Astronomical Union,IAU）贊助的開源項目，旨在為天文計算提供權(quán)威的算法程序，避免了許多不必要的重復(fù)工作，使科技人員將主要精力集中在解決問題和思路創(chuàng)新上。SOFA 程序集實現(xiàn)了各類時間（世界時、歷書時、國際原子時、地球時、太陽系質(zhì)心力學(xué)時、協(xié)調(diào)世界時和儒略日）轉(zhuǎn)換，以及各類坐標(biāo)系（歲差、章動、視恒星時、極移）旋轉(zhuǎn)轉(zhuǎn)換，其中很多坐標(biāo)系轉(zhuǎn)換又有多種版本，比如歲差就有IAU1976、IAU2000、IAU2006等不同模型。SOFA的坐標(biāo)系轉(zhuǎn)換精度極高，可達(dá)到毫角秒（milli-arcsecond,mas）甚至微角秒（μas）級精度。但是，SOFA 主要面向天文學(xué)計算和空間大地測量學(xué)應(yīng)用[7-10]，它涉及較多的天文學(xué)知識，特別是章動參數(shù)計算比較煩瑣，對普通導(dǎo)航技術(shù)人員而言難以快速理解、掌握和直接借用。本文主要從導(dǎo)航專業(yè)的實用角度出發(fā)，簡要介紹了天文學(xué)中常用的時間系統(tǒng)概念和空間參考系轉(zhuǎn)換知識，著重關(guān)注它們之間的相互聯(lián)系，對IAU 提供的歲差、章動和視恒星時計算模型進(jìn)行簡化，給出了一套簡單有效的天文導(dǎo)航信息轉(zhuǎn)換算法，能夠滿足慣性/天文導(dǎo)航中角秒級的天文測姿精度需求。

1 時間系統(tǒng)

1.1 世界時

自十七世紀(jì)以來，天文學(xué)家們普遍以地球自轉(zhuǎn)作為時間尺度。以地球自轉(zhuǎn)為基礎(chǔ)的時間計量系統(tǒng)稱為世界時（Universal Time,UT）。1820年法國科學(xué)院正式提出秒長的定義：全年中所有真太陽日的平均長度的1/86400 為1 秒，即后來所謂的平太陽秒或世界時秒。人們曾認(rèn)為地球自轉(zhuǎn)速率是均勻的，因而UT 曾被誤認(rèn)為是一種均勻時，1960年前作為基本時間計量系統(tǒng)被廣泛應(yīng)用。后來，隨著科技的發(fā)展，人們發(fā)現(xiàn)地球自轉(zhuǎn)的速率并不均勻，存在長期減緩和短周期變化（穩(wěn)定度約為10-8量級），且自轉(zhuǎn)軸在地球內(nèi)部的位置是不斷變化的（即極移現(xiàn)象），情況比較復(fù)雜。

UT 中又有三種不同的概念：UT0、UT1 和UT2。其中UT0 是由全球分布的多個觀測站觀測恒星的視運(yùn)動而確定的，UT1 是UT0 加上極移修正得到的，UT2 是UT1 加上地球自轉(zhuǎn)的季節(jié)性變化改正得到的。文中用于描述地球自轉(zhuǎn)角度的世界時均指UT1。

1.2 歷書時

與地球自轉(zhuǎn)相比，行星繞日公轉(zhuǎn)的周期要穩(wěn)定得多，1958年國際天文學(xué)聯(lián)合會（IAU）決議從1960年起采用歷書時（Ephemeris Time,ET）取代UT 作為描述天體運(yùn)動、編制天體歷書中所采用的時間系統(tǒng)。ET 的基本時間單位——?dú)v書秒等于1899年12月31日12:00 時起算的回歸年長度的1/31556925.9747，1900年1月1日12:00 時為ET 的起點(diǎn)。ET 的觀測精度隨著觀測設(shè)備和技術(shù)的提高而逐漸提高，從早期的ET0、ET1 發(fā)展到后來的ET2，利用連續(xù)幾年的天文觀測可達(dá)到10-10量級的時間穩(wěn)定度，高于以地球自轉(zhuǎn)為基準(zhǔn)的時間系統(tǒng)UT，但由于ET 的精度依然有限，且不能實時給出精確的ET 時間，所以它在使用上仍有一定的局限性。

1967年起用原子時代替歷書時作為基本的時間計量系統(tǒng)，但在天文歷表上仍用歷書時。1984年起天文計算和歷表上所用的時間單位，也都以原子時秒為基礎(chǔ)了。

1.3 國際原子時

1967年，第13 屆國際計量大會（法文Conference Generale des Poids et Mesures,CGPM）通過國際單位秒（SI 秒）決議[11,12]：位于海平面上的銫原子133（Cs133）基態(tài)的兩個超精細(xì)能級間在零磁場中躍遷振蕩9192631770 個周期所持續(xù)的時間為一個原子時秒。CGPM 在1971年正式確定，國際原子時的秒長與1900.0 歷元時刻的ET 秒長相等，取1958年1月1日0 時0 分0 秒（UT1）的瞬間作為同年同月同日0時0 分0 秒國際原子時（International Atomic Time,TAI）。但事后發(fā)現(xiàn)，在該瞬間兩者相差0.0039 s，這一差值就作為歷史事實而保留了下來。ET 的秒長是通過對長周期的天文觀測再等分而得到的，為了提高精度需要進(jìn)行長時間的觀測；而TAI 的秒長是靠原子鐘來復(fù)現(xiàn)的，原子躍遷頻率的周期很短，通過對短周期的觀測再倍乘后便可得到高精度的秒長，容易測定和應(yīng)用。早期原子鐘的穩(wěn)定度為10-10，而目前原子鐘的穩(wěn)定度已達(dá)10-16甚至更高。

1.4 地球時與太陽系質(zhì)心力學(xué)時

ET 是太陽質(zhì)心坐標(biāo)系中的一種均勻時間尺度，它是牛頓運(yùn)動方程中的獨(dú)立變量。在牛頓力學(xué)中，時間是與空間的位置及能量無關(guān)的一個獨(dú)立變量，它既可以作為行星繞日運(yùn)行運(yùn)動學(xué)方程中的自變量，也可以作為衛(wèi)星繞地球運(yùn)動的運(yùn)動學(xué)方程中的自變量。隨著現(xiàn)代觀測技術(shù)和計時精度的不斷改善，牛頓力學(xué)經(jīng)典理論與觀測結(jié)果之間的矛盾開始顯現(xiàn)。1976年，IAU引入相對論時間尺度，給出了地球質(zhì)心力學(xué)時（Terrestrial Dynamic Time,TDT,1991年改稱地球時TT)和太陽系質(zhì)心力學(xué)時。TT 主要是用于解算繞地球質(zhì)心旋轉(zhuǎn)的天體（如人造衛(wèi)星）的運(yùn)動方程，并編制其星歷時所用的一種時間系統(tǒng)；而太陽系質(zhì)心力學(xué)時（Barycentric Dynamic Time,TDB）是用以解算繞太陽系質(zhì)心旋轉(zhuǎn)的天體（如行星）的運(yùn)動方程，并編制其星歷時所用的時間系統(tǒng)。TT 是連續(xù)且均勻的時間系統(tǒng)，本質(zhì)上是基于原子物理所定義的原子時系統(tǒng)，可以視為是一種在大地水準(zhǔn)面上實現(xiàn)的與SI 秒相一致的理想化的原子時。TT 的起始?xì)v元為1971年1月1日0時0 分32.184 秒（UT1），其中32.184 秒為起始?xì)v元時刻與TAI 的差值，這樣就可使TT 與過去使用的ET相銜接，只需把過去歷表中的ET 改為TT 便可繼續(xù)使用。

TT 和TDB 之間的差別主要是由地月系繞日運(yùn)動的引力勢變化相對論效應(yīng)引起的，兩者間的關(guān)系為：

式中，M為地月系質(zhì)心繞日軌道的平近點(diǎn)角，變化周期為一年；O2為可忽略的微小量。由式(1)可見，TT和TDB 之間不存在長期趨勢項，最大差異僅約1.7 ms。

1.5 協(xié)調(diào)世界時

在確定國際原子時TAI 起點(diǎn)之后，由于地球自轉(zhuǎn)速率不均勻，即世界時UT1 不均勻，其秒長將變得越來越長，導(dǎo)致UT1 與TAI 之間的時差逐年累積。協(xié)調(diào)世界時（Coordinated Universal Time,UTC）兼顧了對UT1 時刻和TAI 秒長兩者的需要，UTC 秒長與TAI秒長一致，在時刻上則要求盡量接近于UT1，是介于TAI 和UT1 之間一種折中的時間系統(tǒng)。UTC 采用閏秒方法來解決TAI 與UT1 之間的累積偏差，使得UTC和UT1 之間的差值總是保持在±0.9 s 之內(nèi)。1972年CGPM 決定啟用UTC，從1958年到1971年末UT1已落后于TAI 將近10 s，從1972年至今（2021年12月）UTC 又進(jìn)行了27 次閏秒（最近一次發(fā)生在2017年1月1日），因此目前有差值dTAI=TAI-UTC=37s。

1979年12月，UTC 取代UT1 作為無線電通信中的標(biāo)準(zhǔn)時間，UTC 也是世界上大多數(shù)國家的法定時間[11,12]，許多傳感器輸出的時標(biāo)信息都是基于UTC 時間的，標(biāo)準(zhǔn)慣導(dǎo)系統(tǒng)和天文導(dǎo)航系統(tǒng)也輸出UTC 時間。UTC 與UT1、TT 之間的轉(zhuǎn)換關(guān)系為：

其中，dUT1=UT1-UTC和dTAI=TAI-UTC 的具體數(shù)值可從國際地球自轉(zhuǎn)服務(wù)（www.iers.org）網(wǎng)站上查閱[13]，n即為UTC 的閏秒數(shù)，目前取值為37。

1.6 儒略日

儒略日（Julian Day,JD）是一種不涉及年、月等概念的長期連續(xù)的記日法，在天文學(xué)中有著廣泛的應(yīng)用，計算跨越多年的兩個時刻之間的間隔時采用這種方法特別方便。由年-月-日-時轉(zhuǎn)換為儒略日的計算式為[10]：

式中，常數(shù)1721013.5 為公元1年1月1日0 時的儒略日；Y、M、D分別為公歷中的年、月、日數(shù)，h為小時數(shù)（分鐘和秒轉(zhuǎn)化為小時的小數(shù)部分），可為UT1、TT 或TDB 等時間系統(tǒng)之一，需根據(jù)實際情況選用，Y年-M月-D日-h時時刻可簡單記為t（Y,M,D,h）時刻；int()為取整符號。特別的，儒略日的起點(diǎn)時間（JD=0日）為公元前4713年1月1日12:00 時；J2000.0 起始?xì)v元（2000年1月1日12:00 時）對應(yīng)的儒略日為JD=2451545.0。

另外，儒略世紀(jì)（Julian Century）是天文學(xué)中常用的一種長間隔時間單位，兩個儒略日J(rèn)D1和JD0之間的儒略世紀(jì)數(shù)T為：

2 空間坐標(biāo)變換矩陣

首先給出基本坐標(biāo)變換矩陣的定義：如有一個坐標(biāo)系（記為a 系）繞其OX軸（或OY軸、或OZ軸）旋轉(zhuǎn)α角度得另一坐標(biāo)系（記為b 系），則從a 系到b系的坐標(biāo)變換矩陣可分別表示為：

式中，sα，cα分別表示sinα，cosα。后面將(α)分別簡寫為Ri(α),(i=X,Y,Z)。

2.1 歲差變換矩陣

由于太陽、月球及行星對地球上赤道隆起部分的作用力矩而導(dǎo)致赤道平面的進(jìn)動，稱為赤道歲差，其角速率每年約西移50.39″；行星對地月質(zhì)心繞日公轉(zhuǎn)黃道面影響，造成春分點(diǎn)移動，稱為黃道歲差，其角速率每年約東移0.1″，而且還會使黃赤交角產(chǎn)生變化。綜合赤道歲差和黃道歲差，IAU1976 歲差模型給出黃經(jīng)歲差ψ和黃赤交角ε計算公式，分別為：

若取TTDB=0.5世紀(jì)數(shù)（即從2000年到2050年），式(7)右端的項最大取值僅為 1.11161′≈0.27′，忽略式(7)右端0.1″量級及更微小的項，近似有：

后文將遵循同樣的規(guī)則：如果t時刻歷元直到2050年某項誤差的影響為0.1″量級或更小，則直接對其進(jìn)行忽略，以簡化計算。

根據(jù)式(7)中的兩個參數(shù)推導(dǎo)，IAU1976 歲差模型采用三次旋轉(zhuǎn)的方式，計算從J2000 協(xié)議天球坐標(biāo)系（i 系）至t時刻的平天球坐標(biāo)系（記為P 系）的坐標(biāo)變換矩陣（即歲差矩陣），為：

其中，三個歲差角分別為：

IAU1976 歲差模型在2003年之前曾被廣泛使用，后來IAU2000 歲差模型在IAU1976 的基礎(chǔ)上簡單地對黃經(jīng)歲差的速率和交角歲差的速率進(jìn)行了修正，2006年后又啟用了最新的IAU2006 歲差模型，其精度可達(dá)1 μas 量級[10]。針對本文研究的角秒級的慣性/天文導(dǎo)航應(yīng)用，IAU1976 模型的精度就足夠了。

2.2 章動變換矩陣

IAU2000A 章動序列模型包含678 個日月章動項和687 個行星章動項，其計算精度優(yōu)于0.2 mas，計算比較復(fù)雜，而IAU2000B 章動序列只包含77 個日月章動項和一些行星章動偏差項，它與IAU2000A 之差不超過1 mas[7]。IAU1980 章動模型更為簡單，下面以它為例作介紹。

IAU1980 章動模型亦采用三次旋轉(zhuǎn)的方式，計算從t時刻的平天球坐標(biāo)系（P 系）至t時刻的瞬時天球坐標(biāo)系（記為N 系）的坐標(biāo)變換矩陣（即章動矩陣），為：

其中，Δψ、Δε分別為黃經(jīng)章動和黃赤交角章動，計算公式為：

式中，l是月球的平近點(diǎn)角，l′是太陽的平近點(diǎn)角，D是太陽到月球的平角距，Ω是月球升交點(diǎn)的平黃經(jīng)，L是月球的平黃經(jīng)，F(xiàn)=L-Ω，以上各參數(shù)的具體表達(dá)式如式(13)所示。

式(12)中，Ai、Ai′為黃經(jīng)章動振幅及其變化率，Bi、Bi′為交角章動振幅及其變化率；N1i、N2i、...，N5i為計算幅角βi時用到的一些系數(shù)，其取值為0、±1、±2、±4。前6 項振幅較大的系數(shù)摘錄如表1。

當(dāng)TTDB＜0.5世紀(jì)數(shù)時，表1中所有變化率系數(shù)Ai′、Bi′的影響均很小，且Ai、Bi中絕對值小于0.1者的影響也很小，因此僅保留絕對值大于0.1 者（表1中粗斜體數(shù)據(jù)），這時式(12)可以直接展開，近似為式(14)所示的少數(shù)幾項之和，極大降低了計算量。

表1 IAU1980 章動序列系數(shù)表（按Ai 絕對值由大至小排列）[10]Tab.1 IAU1980 nutation series listed according to abs(Ai)

2.3 視恒星時變換矩陣

相對于慣性空間的地球自轉(zhuǎn)角度（Earth Rotation Angle,ERA），可用當(dāng)?shù)刈游缇€相對于天球上的非旋轉(zhuǎn)原點(diǎn)（Non-Rotational Origin,NRO）的角運(yùn)動來度量。格林尼治恒星時（Greenwich Sidereal Time,GST）以春分點(diǎn)通過格林尼治子午圈為參考起點(diǎn)，在數(shù)值上等于春分點(diǎn)相對于格林尼治子午圈的時角，格林尼治恒星時字面上雖帶“時”字，其實給出的是0～2π 之間的角度值。由于章動的影響，春分點(diǎn)有平春分點(diǎn)和真春分點(diǎn)之分，GST 也有格林尼治平恒星時角（Greenwich Mean Sidereal Time,GMST）和格林尼治視恒星時角（Greenwich Apparent Sidereal Time,GAST）之分。

在IAU2000 模型中，ERAEθ為

其中，dJDUT1=JDUT1-JD00，為J2000.0參考?xì)v元(JD00=2451545.0)至t時刻歷元JDUT1的儒略日間隔，注意計算JDUT1的時間應(yīng)該使用UT1 時間。式(15)的計算與歲差及章動無關(guān)。由式(15)還可知，地球轉(zhuǎn)速的不均勻即意味著世界時UT1 的不均勻。

在地球自轉(zhuǎn)角度的基礎(chǔ)上考慮歲差的影響，GMSTMθ的計算公式為：

注意到，式(16)右端的多項式θψε近似表示黃經(jīng)總歲差在天赤道上的投影，即近似有θψε=ψ· cosε。

進(jìn)一步考慮章動的影響，得GASTθA的計算公式為：

其中，Δψ· cosε為黃經(jīng)章動在天赤道上的投影，即赤經(jīng)章動。這表明，GAST 是瞬時春分點(diǎn)與格林尼治零度經(jīng)線間的夾角，其變化主要取決于地球自轉(zhuǎn)角，同時與歲差和章動導(dǎo)致瞬時春分點(diǎn)的移動有關(guān)。

最后，將瞬時地球坐標(biāo)系（記為E 系）繞其Z 軸旋轉(zhuǎn) Aθ-角度后便得瞬時天球坐標(biāo)系（N 系），該變換矩陣即為視恒星時變換矩陣

2.4 極移變換矩陣

由于地球表面上的物質(zhì)運(yùn)動（如洋流等）以及地球內(nèi)部的物質(zhì)運(yùn)動（如地幔對流等），地球自轉(zhuǎn)軸在其內(nèi)部會緩慢變化，從而引起地球瞬時極點(diǎn)在地面上位置的不斷移動，導(dǎo)致瞬時地球坐標(biāo)系與協(xié)議地球坐標(biāo)系之間的差異。極移變換矩陣為：

其中，e 為協(xié)議地球坐標(biāo)系（比如WGS84 坐標(biāo)系），xp和yp是瞬時極相對于國際協(xié)議原點(diǎn)（Conventional International Origin,CIO）的角度值，其方向分別是從CIO 指向0°經(jīng)線和270°經(jīng)線方向，兩者大小一般小于0.5″，即時數(shù)值變化可在www.iers.org網(wǎng)站查得[13]。

3 簡化計算與誤差分析

前面已對IAU 模型各計算公式進(jìn)行了近似和簡化，忽略了0.1″及更小量級的項。下面將對時間系統(tǒng)類型（UTC、UT1、TDB 和TT）的選用問題做進(jìn)一步的分析和簡化。

如果時間輸入?yún)?shù)UT1 相差1 s，根據(jù)ERA 計算式(15)求誤差，可得：

該值近似是地球在1 s 內(nèi)轉(zhuǎn)過的角度（15.04″）。因此，為了使ERA 的計算誤差達(dá)到0.1″量級，需要UT1 的輸入精度達(dá)到0.01 s；若要ERA 的計算誤差為1″量級，則UT1 的精度只需達(dá)到0.1 s 即可。

在GMST 的計算式(16)中，若以UTC 代替TT（目前TT-UTC=37+32.184=69.184s），則由θψε引起的誤差僅為：

由此可見，以UTC 代替TT 引起的GMST 計算誤差極小，完全可忽略不計。

同理，在歲差參數(shù)式(10)中，以UTC 代替TT，其對歲差矩陣的影響也可忽略，即有：

最后，在章動計算式(14) 中，若僅對Δψ≈A1sin（N51Ω）作誤差估計，以UTC 代替TT，可得其引起的誤差為：類似的，可計算得式(14)中其它項的誤差，也都非常小。

上述分析表明，最主要的時間誤差來源于ERA的計算，其UT1 輸入精度至少應(yīng)該達(dá)到0.1 s，而其他天文參數(shù)計算中的時間TT 或TDB 均可以直接使用UTC 時間代替，對導(dǎo)航應(yīng)用而言沒有任何區(qū)別。

更進(jìn)一步，在式(21)中，若將69.184 s 擴(kuò)大為1天（24×3600 s）跨度，則有：

同理，在1 天跨度下，式(22)(23)的誤差也都不會大于0.1″量級。因此，除ERA 須實時計算外，其他天文參數(shù)的計算可以每天僅更新計算一次，便可滿足角秒級天文導(dǎo)航的精度需求。

4 仿真驗證

以J2000.0 歷元為時間中點(diǎn)，給出公元1900年～2100年兩百年間的歲差（Precession）、章動（Nutation）及格林尼治視恒星時（GAST）計算結(jié)果，采樣間隔為36.525 天，并在固定間隔中加入均值為0、均方差為12 天的正態(tài)分布隨機(jī)抖動。

對比仿真結(jié)果如圖1所示。圖1(a)(c)(e)為SOFA算法程序（調(diào)用函數(shù)iauPn06a 和iauGst06a）給出的結(jié)果，將作為參考真值使用。其中圖1(a)為矩陣RPi按“3-1-2”方式給出的俯仰、橫滾和方位三個歲差歐拉角，分別記為θP，γP，ψP，圖1(c)為矩陣RNP按“3-1-2”方式給出的俯仰、橫滾和方位三個章動歐拉角，分別記為θN，γN，ψN，圖1(e)為GAST角度θA。圖1(b)(d)(f)為本文簡化算法（除計算ERA 外，其他公式中的時間輸入?yún)?shù)均加入均值為0、均方差為0.33 天的正態(tài)分布隨機(jī)抖動）與SOFA 算法參考基準(zhǔn)之間的歲差歐拉角誤差、章動歐拉角誤差和GAST 誤差。由誤差圖1(b)(d)(f)可見，三個章動歐拉角誤差基本上小于2″；歲差俯仰誤差角δθP很小、歲差橫滾誤差角δγP基本小于0.5″；只有歲差方位誤差角δψP和GAST 角度誤差δψA在兩端稍大些但兩者取值大致互為相反數(shù)。

圖1 歲差、章動、GAST 歐拉角及其誤差曲線Fig.1 Precession,nutation,GAST Euler angles and the corresponding errors

根據(jù)2 節(jié)的各項坐標(biāo)系變換，綜合可得從協(xié)議天球坐標(biāo)系（J2000,i 系）到協(xié)議地球坐標(biāo)系（WGS84,e系）的坐標(biāo)系轉(zhuǎn)換矩陣[14,15]，為：

式(25)轉(zhuǎn)換過程示意圖參見圖2。

圖2 從i 系到e 系的坐標(biāo)系轉(zhuǎn)換過程Fig.2 Flowchart of coordinate system transformation from i-frame to e-frame

若不考慮極移，即令xp=yp=0（REe=I3×3），使用圖1中的仿真數(shù)據(jù)并按式(25)計算從i 系到e 系轉(zhuǎn)換矩陣Rei，同樣按“3-1-2”方式給出i 系至e 系的歐拉角，圖3是本文算法與SOFA 算法參考基準(zhǔn)之間的歐拉角誤差。

圖3 從J2000 到WGS84 轉(zhuǎn)換（Re i）的歐拉角誤差Fig.3 Euler angle errors of coordinate system transformation from J2000 to WGS84 (Re i)

由圖3可見，前述圖1中誤差δψP和δψA之間大部分相互抵消了，從而使得這里Rei的三個歐拉角誤差基本都小于0.5″。這說明，對于實際導(dǎo)航應(yīng)用，本文所提簡化算法的誤差在2050年之前基本可小于0.2″。

5 結(jié) 論

在慣性/天文組合導(dǎo)航中，需要將在J2000 天球坐標(biāo)系（i 系）中測量的天文傳感器信息轉(zhuǎn)換到WGS84地球坐標(biāo)系（e 系），在專業(yè)的天文學(xué)中其精確的轉(zhuǎn)換過程比較煩瑣，特別是章動矩陣的計算量較大。本文研究表明，除地球自轉(zhuǎn)角的計算需要比較精確的UT1 時間（0.01 s 量級）外，其他天文參數(shù)計算中的時間輸入均可直接使用UTC 時間代替，甚至在一天之內(nèi)僅需更新計算一次從i 系到e 系的坐標(biāo)轉(zhuǎn)換，都不會引起明顯的坐標(biāo)轉(zhuǎn)換誤差，能夠滿足角秒級的天文傳感器測姿精度需求。仿真結(jié)果顯示，與SOFA 參考基準(zhǔn)相比，所提簡化算法的轉(zhuǎn)換精度達(dá)到了0.2″。