本發(fā)明屬于月球著陸器精密定位與深空探測領(lǐng)域，特別是涉及一種基于四程中繼跟蹤模式的遠(yuǎn)月面著陸器精密定位方法。

背景技術(shù)：

目前，深空探測是世界各國競相爭逐的熱點(diǎn)，進(jìn)行深空探測，可以進(jìn)一步解答地球如何起源與演變、行星和太陽系形成和演化的過程、人類是不是宇宙中唯一的智慧生命、地球的未來將如何等一系列問題，同時(shí)有利于人類積極開發(fā)和利用空間資源。在深空探測中，探測器的定軌定位是任務(wù)成敗的關(guān)鍵，也是各種科學(xué)任務(wù)順利進(jìn)行的前提。精密的探測器軌道是進(jìn)行地形地貌測繪的基礎(chǔ)數(shù)據(jù)，也可以用于行星重力場的解算，進(jìn)而反演行星的內(nèi)部構(gòu)造。

月球背面的著陸探測一直是國際上的熱點(diǎn)和難點(diǎn)。月球背面具有不同于月球正面的獨(dú)特地質(zhì)構(gòu)造，是研究地月起源的重要突破口；月球背面沒有地球電磁波干擾，是進(jìn)行低頻射電天文觀測的天然理想場所。但是由于月球自轉(zhuǎn)與公轉(zhuǎn)同步，月球背面著陸器無法與地球測站進(jìn)行直接通信，任務(wù)難度和風(fēng)險(xiǎn)比較大，目前對月球背面的就位探測仍舊是一項(xiàng)國際空白。因此，開展月球背面的就位探測具有重要的科學(xué)意義和工程意義。

對于近月面著陸器的定位問題，一般采用運(yùn)動(dòng)學(xué)統(tǒng)計(jì)定位方法，采用測量模式包括：雙程/三程測距、測速，VLBI時(shí)延/時(shí)延率模型等等。在“嫦娥三號(hào)”著陸器的精密定位中使用了運(yùn)動(dòng)學(xué)統(tǒng)計(jì)定位方法，綜合了測距、測速數(shù)據(jù)和VLBI時(shí)延、時(shí)延率數(shù)據(jù)，著陸器定位的絕對精度在10m左右(曹建峰,張宇,胡松杰,等.2016.嫦娥三號(hào)著陸器精確定位與精度分析.武漢大學(xué)學(xué)報(bào)(信息科學(xué)版),41(2):274-278.doi:10.13203/j.whugis20140123.)，“玉兔號(hào)”漫游器與著陸器的相對定位采用了同波束VLBI技術(shù)，相對位置精度可達(dá)到米級(jí)(黃勇,昌勝騏,李培佳,等.2014.“嫦娥三號(hào)”月球探測器的軌道確定和月面定位.科學(xué)通報(bào),59(23):2268-2277.)。

但是對于月球背面著陸器的精密定位問題，如圖1所示，圖1a為近月面示意圖，圖1b為遠(yuǎn)月面示意圖，由于月球自身的遮擋，遠(yuǎn)月面著陸器與地球深空測控站無法通視，傳統(tǒng)的直接測量模式，如“嫦娥三號(hào)”中的雙程、三程測距/測速，VLBI時(shí)延/時(shí)延率測量模式將不再適用，因此有必要采取新的跟蹤測量模式對遠(yuǎn)月面著陸器進(jìn)行定位。

在對月球著陸器定位的過程中，涉及到坐標(biāo)系的轉(zhuǎn)換。地面測站位于地球，一般用地球參考框架描述；月球著陸器需要用月球主軸坐標(biāo)系來表示；而信號(hào)在太空中的傳輸要考慮大天體引力時(shí)延以及相對論效應(yīng)的影響，最終均需要在太陽系質(zhì)心坐標(biāo)系下(BCRS)表示。為了得到高精度的定位結(jié)果，需要不同坐標(biāo)系之間進(jìn)行高精度的轉(zhuǎn)換。各坐標(biāo)系之間的轉(zhuǎn)換如圖2所示。

地面測站坐標(biāo)由協(xié)議地固系(ITRS)轉(zhuǎn)為J2000地球慣性坐標(biāo)系下的坐標(biāo)，此轉(zhuǎn)換需要?dú)q差旋轉(zhuǎn)、章動(dòng)旋轉(zhuǎn)、極移旋轉(zhuǎn)以及地球自轉(zhuǎn)旋轉(zhuǎn)。關(guān)于歲差、章動(dòng)模型，隨著觀測值精度的不斷提高，也在不斷地修正。目前IERS2010規(guī)范推薦采用IAU2006決議的歲差章動(dòng)模型。J2000地球慣性系坐標(biāo)經(jīng)過平移、洛倫茲變換可以轉(zhuǎn)到太陽系質(zhì)心坐標(biāo)系下。

同地固坐標(biāo)系統(tǒng)的轉(zhuǎn)換中用到歲差章動(dòng)和極移一樣，從月球主軸坐標(biāo)系(PA)到月心天球坐標(biāo)系，要考慮月球天平動(dòng)(Libration)的影響。目前用于月球物理天平動(dòng)研究的主要手段仍然是對月激光測距(LLR)，由月球自轉(zhuǎn)引起測距的變化，通過數(shù)值積分來求得月球自轉(zhuǎn)相應(yīng)的三個(gè)歐拉角，通過旋轉(zhuǎn)三個(gè)歐拉角即可轉(zhuǎn)換到月心天球坐標(biāo)系下，進(jìn)一步經(jīng)過平移、洛倫茲變換可以轉(zhuǎn)到太陽系質(zhì)心坐標(biāo)系(BCRS)下。

目前對月球背面的著陸一直是國際空白，月球背面的著陸最大的難題是克服信號(hào)被月球遮擋的影響，同時(shí)面臨無法實(shí)時(shí)測控跟蹤著陸器等難題。

技術(shù)實(shí)現(xiàn)要素：

本發(fā)明的目的是克服傳統(tǒng)直接跟蹤測量模式無法對遠(yuǎn)月面著陸器進(jìn)行定位的缺陷，提供一種基于四程中繼跟蹤模式的遠(yuǎn)月面著陸器精密定位方法，該方法借助于中繼星的中繼跟蹤測量，可以消除月球自身的遮擋，實(shí)現(xiàn)對遠(yuǎn)月面著陸器的精密定位。

為達(dá)到上述發(fā)明目的，本發(fā)明的技術(shù)方案提供一種基于四程中繼跟蹤模式的遠(yuǎn)月面著陸器精密定位方法，建立月球著陸器-軌道器的四程中繼跟蹤測量模式的觀測模型，基于觀測模型實(shí)現(xiàn)遠(yuǎn)月面著陸器精密定位，

所述四程中繼跟蹤測量模式，將月球軌道器作為中繼星，使用月球著陸器和月球軌道器之間的鏈路進(jìn)行測量，設(shè)地球跟蹤站T_i在時(shí)刻i發(fā)送一個(gè)上行的信號(hào)給月球軌道器，月球軌道器位置為S_j，時(shí)刻為j；經(jīng)過轉(zhuǎn)發(fā)，信號(hào)傳送給月面著陸器，相應(yīng)位置為L_k，時(shí)刻為k，L_k位于月球背面；月面著陸器再次轉(zhuǎn)發(fā)給在軌的月球軌道器，此時(shí)軌道器已經(jīng)從S_j的位置運(yùn)動(dòng)至S_m，時(shí)刻為m；之后該信號(hào)下行至地球跟蹤站T_n，時(shí)刻為n；設(shè)時(shí)刻n、m、k、j、i分別相應(yīng)的參與者記為T_n,S_m,L_k,S_j,T_i，則有以下4個(gè)幾何距離，

R₁＝|X(S_m)-X(T_n)|

R₂＝|X(S_m)-X(L_k)|

R₃＝|X(S_j)-X(L_k)|

R₄＝|X(S_j)-X(T_i)|

其中，X(j)為參與者j在太陽系質(zhì)心坐標(biāo)系下下的位置矢量；

四程中繼測距的觀測模型建立如下，

R＝(R₁+c·RLT_nm)+R₂+R₃+(R₄+c·RLT_ij)+

c·[TDB(i)-UTC(i)]-c·[TDB(n)-UTC(n)]

＝c·[UTC(n)-UTC(i)]

其中，RLT_nm為由參與者T_n發(fā)射信號(hào)到S_m的廣義相對論時(shí)延，RLT_ij為由參與者T_i發(fā)射信號(hào)到S_j的廣義相對論時(shí)延，UTC(n)和UTC(i)分別為時(shí)刻n和時(shí)刻i對應(yīng)的協(xié)調(diào)世界時(shí)時(shí)標(biāo)，R為四程中繼測距值，TDB(n)和TDB(i)分別為時(shí)刻n和時(shí)刻i對應(yīng)的太陽系質(zhì)心力學(xué)時(shí)，c為光速；

設(shè)在一個(gè)多普勒積分周期內(nèi)，起始時(shí)刻T_s和終止時(shí)刻T_e各自對應(yīng)的四程中繼測距值為R_s、R_e，建立四程中繼測速的觀測模型如下，

其中，T_c為一個(gè)多普勒積分周期，RR為四程中繼測速值。

而且，根據(jù)四程中繼跟蹤測量模式，將月球軌道器的軌道和遠(yuǎn)月面著陸器的位置同時(shí)視為待估參數(shù)，采用精密定軌的方式進(jìn)行解算，最后得到兩者的最佳估值。

而且，采用精密定軌的方式進(jìn)行解算時(shí)，

計(jì)算四程中繼測距值R對遠(yuǎn)月面著陸器坐標(biāo)X(L_k)的偏導(dǎo)數(shù)H_R公式如下，

按照上式計(jì)算積分起始時(shí)刻T_s和積分終止時(shí)刻T_e分別對應(yīng)的四程中繼測距值R對遠(yuǎn)月面著陸器坐標(biāo)X(L_k)的偏導(dǎo)數(shù)H_s和H_e，則四程中繼測速值RR對遠(yuǎn)月面著陸器坐標(biāo)X(L_k)的偏導(dǎo)數(shù)H_RR為

將結(jié)果作為觀測時(shí)刻的觀測值對狀態(tài)向量的偏導(dǎo)數(shù)。

而且，將探測器軌道跟蹤數(shù)據(jù)文件分為多個(gè)觀測弧段，將初始參考軌道x₀^ref＝[r₀v₀ p₀ q₀]^T中的中繼星的初始位置r₀、初始速度v₀以及力模型參數(shù)p₀作為局部參數(shù)處理，記作x₁，維數(shù)為m₁；參數(shù)q₀是遠(yuǎn)月面著陸器位置坐標(biāo)，將其作為全局參數(shù)處理，記作x₂，維數(shù)為3；采用精密定軌的方式進(jìn)行解算時(shí)，融合多個(gè)觀測弧段相應(yīng)法方程，將非線性的方程線性化，使用迭代法求解。

本發(fā)明的有益效果：本發(fā)明基于四程中繼跟蹤模式，利用一顆中繼星，實(shí)現(xiàn)對遠(yuǎn)月面著陸器的事后精密定位；由于基于四程中繼跟蹤模式建立了中繼星和遠(yuǎn)月面著陸器間的聯(lián)系，解決了地球深空站無法直接跟蹤測量遠(yuǎn)月面著陸器的問題。由于該方式對著陸器以及中繼星有很強(qiáng)的幾何約束，將使精密定軌、定位的精度得到極大地提高。此外，隨著此模式觀測數(shù)據(jù)的增加，可將用于支持對月球自轉(zhuǎn)歐拉角，月球K2勒夫數(shù)以及月球重力場等反演，得到更為精確的月球物理模型，進(jìn)而為了解月球的起源與演化，深入研究月球內(nèi)部構(gòu)造提供重要的科學(xué)數(shù)據(jù)。本發(fā)明還提出將遠(yuǎn)月球著陸器坐標(biāo)作為全局參數(shù)進(jìn)行解算，理論上更為嚴(yán)密，可得到高精度的遠(yuǎn)月面著陸器坐標(biāo)。因此，本發(fā)明具有重要的科研價(jià)值和市場價(jià)值，是我國領(lǐng)先世界的重要研究項(xiàng)目。

附圖說明

圖1為傳統(tǒng)的雙程/三程測距、測速及VLBI模型示意圖，其中圖1a為近月面示意圖，圖1b為遠(yuǎn)月面示意圖；

圖2為本發(fā)明涉及到的各類坐標(biāo)系之間的轉(zhuǎn)換示意圖；

圖3為本發(fā)明實(shí)施例的基于月球著陸器-軌道器的四程中繼跟蹤測量模型示意圖。

圖4為本發(fā)明實(shí)施例的解算流程示意圖。

具體實(shí)施方式

以下結(jié)合附圖和實(shí)施例對本發(fā)明技術(shù)方案進(jìn)行具體描述。

為便于實(shí)施參考起見，首先分別介紹本發(fā)明涉及的精密定軌原理和基于月球著陸器-軌道器的四程中繼跟蹤測量模式(4W L-O)的基本原理：

(1)精密定軌原理

本發(fā)明中涉及到的精密定軌理論不考慮探測器搭乘火箭被送入指定軌道的過程，而是進(jìn)入指定軌道后，通過探測器自身發(fā)射無線電信號(hào)，建立與地面(或空間站和衛(wèi)星)跟蹤測控系統(tǒng)的聯(lián)系，從而對原始軌道進(jìn)行調(diào)控和微分改進(jìn)。具體闡述如下：

第一，根據(jù)牛頓運(yùn)動(dòng)定律，探測器(本文指月球軌道器)進(jìn)入指定軌道(初始軌道)后，受到宇宙中各種力(天體引力、太陽光壓力等)的影響將會(huì)產(chǎn)生加速度，從而改變原有的初始運(yùn)動(dòng)狀態(tài)(初始軌道)。該過程可用一個(gè)常微分方程組描述，稱之為力模型。但由于太空環(huán)境復(fù)雜且人類目前的認(rèn)知程度有限，不可能用該方程準(zhǔn)確無誤地描述出探測器的實(shí)際受力情況，所以該方程組存在誤差即力模型存在誤差，具體體現(xiàn)在兩方面：第一，力模型中的已知參數(shù)不夠準(zhǔn)確；第二，實(shí)際受力中還存在未建模的部分。對于第二點(diǎn)，需要整個(gè)科學(xué)的推動(dòng)才能對受力有進(jìn)一步的認(rèn)識(shí)，所以力模型的誤差主要指第一點(diǎn)。

第二，如果不能對探測器保持聯(lián)系，任何探測器都會(huì)很快失效。所以，探測器將會(huì)主動(dòng)發(fā)射無線電信號(hào)，由地面(或空間站、衛(wèi)星)跟蹤站的跟蹤測量設(shè)備進(jìn)行接收。根據(jù)收發(fā)信號(hào)鏈路的幾何關(guān)系，可以建立相應(yīng)的觀測模型。信號(hào)在實(shí)際傳播過程中，會(huì)受到復(fù)雜空間環(huán)境、相對論效應(yīng)等的影響，該影響通常作為改正項(xiàng)融入觀測模型中，可以視為觀測模型的模型參數(shù)。同樣，由于對這些影響的認(rèn)知有限等原因，這些模型參數(shù)也存在誤差。同時(shí)，跟蹤站的觀測設(shè)備本身存在一定的誤差，所以獲取的觀測資料也是存在誤差的。

結(jié)合上述兩點(diǎn)，精密定軌實(shí)際上可以這樣描述：利用帶有大量帶有誤差的觀測數(shù)據(jù)和并非精確的運(yùn)動(dòng)學(xué)微分方程，使用統(tǒng)計(jì)學(xué)原理解算在某種意義下的探測器初始軌道和相關(guān)模型(力模型、觀測模型)參數(shù)的最優(yōu)估值的過程。

由于一般將觀測誤差視為正態(tài)分布的隨機(jī)誤差，所以精密定軌實(shí)質(zhì)上為基于最小二乘的參數(shù)估計(jì)問題，即：尋找初始軌道和相關(guān)模型(力模型、觀測模型)參數(shù)的估計(jì)值，使得利用估計(jì)值計(jì)算的觀測殘差平方和最小。

其中，第一點(diǎn)中所述力模型的形式如下：

式(1)是一個(gè)描述力模型的常微分方程組，其中：x(t)為隨時(shí)間t變化的m維狀態(tài)向量，包含探測器的位置r(t)、速度v(t)、力模型相關(guān)參數(shù)p以及與觀測模型相關(guān)的參數(shù)q，m表示以上這些參數(shù)的個(gè)數(shù)的總和；為隨時(shí)間變化的m維狀態(tài)向量的一階導(dǎo)數(shù)，f(t,x(t))為的具體函數(shù)形式，x₀為歷元初始時(shí)刻t₀不精確的初值x(t₀)，包含探測器的初始位置r₀(t₀)、初始速度v₀(t₀)、(以下簡寫為r₀和v₀)力模型相關(guān)初始參數(shù)p₀以及與觀測模型相關(guān)的參數(shù)q₀，如式(2)所示。

其中，

令z表示n維觀測向量，即：

z_i(t_i)表示第i個(gè)觀測量，即在觀測時(shí)間t_i獲取的觀測量，i＝1,2,…n，n表示觀測弧段內(nèi)的觀測次數(shù)；z_i(t_i)用于表示第i個(gè)觀測量與探測器狀態(tài)的關(guān)系如下：

z_i(t_i)＝g_i(t_i,x(t_i))+ε_i＝h_i(t_i,x₀)+ε_i (4)

式(4)中：g_i表示觀測量與觀測時(shí)間t_i及該時(shí)刻的探測器狀態(tài)x(t_i)的函數(shù)關(guān)系；h_i()表示觀測量與觀測時(shí)間t_i及歷元初始時(shí)刻t₀的航天器狀態(tài)的函數(shù)關(guān)系；ε_i為在觀測時(shí)間t_i由于誤差造成的理論觀測量與實(shí)際觀測量的偏差。式(5)寫為向量表達(dá)式如下：

z＝h(x₀)+ε (5)

其中，z,h(x₀),ε分別為式(4)中所述的z_i(t_i),h_i(t_i,x₀),ε_i的向量形式；一方面，h(x₀)是一個(gè)關(guān)于未知矢量x₀的非線性函數(shù)，利用最小二乘方法求解軌道確定問題非常復(fù)雜；另一方面，真實(shí)的矢量x₀雖然未知，但是實(shí)際給出的初始探測器狀態(tài)可以視為探測器真實(shí)初始狀態(tài)x₀的近似值，稱之為初始參考軌道，記為x₀^ref。在該點(diǎn)進(jìn)行泰勒展開，并忽略高階項(xiàng)，稱之為線性化，可得：

寫為誤差方程的形式，有：

ε＝Δz-HΔx₀ (7)

Δx₀＝x₀-x₀^ref (8)

Δz＝z-h(x₀^ref) (9)

式(6)至式(10)中：△x₀為真實(shí)初始軌道x₀與不精確的初始參考軌道x₀^ref的偏差；△z為實(shí)際觀測量z與利用該不精確的參考軌道x₀^ref計(jì)算的觀測量的殘差；H為利用歷元t₀時(shí)刻參考軌道初值(即初始參考軌道x₀^ref)計(jì)算的觀測量相對于歷元t₀時(shí)刻狀態(tài)向量的偏導(dǎo)數(shù)，x(t_i)為歷元t_i時(shí)刻的參考軌道，H可以利用鏈?zhǔn)椒▌t進(jìn)行展開成(10)式所示的形式。是t_i時(shí)刻的觀測量對狀態(tài)向量x(t_i)的偏導(dǎo)數(shù)。

通過線性化，精密定軌問題可以簡化為線性最小二乘問題。獲取狀態(tài)改正量△x₀的最優(yōu)估值△x₀^lsq，使得經(jīng)狀態(tài)改正△x₀后的目標(biāo)函數(shù)J(Δx₀)取最小值。

J(Δx₀)＝ε^Tε＝(Δz-HΔx₀)^T(Δz-HΔx₀) (11)

根據(jù)線性方程理論，如果H滿秩，則目標(biāo)函數(shù)取極小值的條件為

計(jì)算上述偏導(dǎo)數(shù)，可得法方程形式如下：

(H^TH)Δx₀＝H^TΔz (13)

解得線性最小二乘問題解Δx₀^lsq的一般形式：

Δx₀^lsq＝(H^TH)^-1(H^TΔz) (14)

矩陣H^TH是一個(gè)m維對稱矩陣，通常稱為法方程矩陣。

如果在軌道估計(jì)之前，在已有待估參數(shù)△x₀初始值的前提下，同時(shí)具備待估參數(shù)△x₀的驗(yàn)前信息

式(15)中，E(△x₀)表示對△x₀取期望，cov(△x₀-x₀)表示對某個(gè)隨機(jī)變量序列取方差-協(xié)方差矩陣；給出的先驗(yàn)信息即為待估參數(shù)△x₀的期望x₀和方差-協(xié)方差矩陣P₀。引入變量如下，

這樣式(11)所示的目標(biāo)函數(shù)就可以寫為

相應(yīng)的法方程可以寫為

其最小二乘解Δx₀^lsq為

則修正后的航天器狀態(tài)x₀^lsq為

x₀^lsq＝x₀^ref+Δx₀^lsq (20)

對于不同觀測類型，測量數(shù)據(jù)精度各有差異。而上述式中將觀測數(shù)據(jù)按相同觀測精度處理，這在實(shí)際應(yīng)用不盡合理，故需要對數(shù)據(jù)進(jìn)行加權(quán)，相應(yīng)的權(quán)矩陣可以寫為如下形式：

W＝diag(σ₁^-2,…,σ_n^-2) (21)

σ_i為互不相關(guān)的測量誤差，i＝1,2,…n。當(dāng)誤差相互關(guān)聯(lián)時(shí)，權(quán)矩陣為非對角矩陣。相應(yīng)的法方程的形式為：

式(18)的最小二乘解為如下形式：

當(dāng)不考慮先驗(yàn)信息時(shí)，P₀^-1＝0；相應(yīng)的法方程變?yōu)椋?/p>

(H^TWH)Δx₀＝H^TWΔz (24)

相應(yīng)的最小二乘解為：

Δx₀^lsq＝(H^TWH)^-1(H^TWΔz) (25)

其中，式中的H按照式(10)進(jìn)行計(jì)算；式中的△z按照式(9)計(jì)算。

由于對非線性方程進(jìn)行了線性化，為了保證最終解的準(zhǔn)確性，根據(jù)精密定軌理論，整個(gè)過程需要迭代進(jìn)行。即將解算得到的x₀^lsq重新作為x₀^ref，重復(fù)式(1)至式(25)，直至與上一次得到的x₀^lsq之差的絕對值小于設(shè)定的限差δ即可。

故本發(fā)明可概括為：結(jié)合一種新的觀測模型——四程中繼跟蹤測量模型，將月球軌道器(中繼星)的軌道和遠(yuǎn)月面著陸器的位置同時(shí)視為待估參數(shù)，按照精密定軌的流程進(jìn)行解算，最后得到兩者的最佳估值的過程。

(2)月球著陸器-軌道器的四程中繼跟蹤測量模式(4W L-O)原理

和傳統(tǒng)雙程模式相比，基于月球著陸器-軌道器的四程中繼跟蹤模式使用了著陸器和軌道器之間的鏈路。如圖3所示，地球跟蹤站T_i(時(shí)刻i)發(fā)送一個(gè)上行信號(hào)給月球軌道器(中繼星，下文軌道器均指中繼星)(位置為S_j，時(shí)刻j)，經(jīng)過轉(zhuǎn)發(fā)，該信號(hào)傳送給月面著陸器(位置為L_k，時(shí)刻k，L_k位于月球背面)，著陸器再次轉(zhuǎn)發(fā)給在軌的月球軌道器(此時(shí)軌道器已經(jīng)從S_j的位置運(yùn)動(dòng)至S_m，時(shí)刻為m)，之后該信號(hào)下行至地球跟蹤站T_n(時(shí)刻n)。上述過程就是基于月球著陸器-軌道器的四程中繼跟蹤模式信號(hào)收發(fā)的一個(gè)完整過程。在實(shí)際測量中，地球跟蹤站T_n可以通過積分多普勒的方式得到4W L-O的距離變率，和雙程測速的模式類似。對于4W L-O測距和測速的理論值C的計(jì)算，和傳統(tǒng)的雙程測距測速類似，采用“回溯法”，即按照T_n-＞S_m-＞L_k-＞S_j-＞T_i的路徑，為便于計(jì)算，可視為有五個(gè)參與者T_n,S_m,L_k,S_j,T_i。地球跟蹤站T_n在n時(shí)刻接收到下行信號(hào)并記錄，但m，k，j，i時(shí)刻的相關(guān)參與者的狀態(tài)向量是未知的，需要進(jìn)行光行時(shí)解算依次獲取。其中，光行時(shí)解算公式如下：

式(26)和式(27)中，c為光速；下標(biāo)a為該支路的信號(hào)發(fā)射方，下標(biāo)b為該支路的信號(hào)接收方；t_b表示該支路的信號(hào)接收方的接收時(shí)刻，t_a表示該支路的信號(hào)發(fā)射方的發(fā)射時(shí)刻；對于4W L-O模式的四條鏈路，a＝i時(shí)b＝j(luò)，a＝j(luò)時(shí)b＝k，a＝k時(shí)b＝m，a＝m時(shí)b＝n。式(26)中，r_ab是BCRS(Barycentric Celestial Reference System)坐標(biāo)系下該支路信號(hào)發(fā)射方和接收方的幾何距離，是該支路信號(hào)傳播過程中的牛頓光行時(shí)，RLT_ab是該支路信號(hào)傳播過程中的廣義相對論時(shí)延。式(27)中，μ_s是太陽引力常數(shù)，μ_B是行星、外行星系統(tǒng)或月球的引力常數(shù)。γ為后牛頓參數(shù)，和代表該支路信號(hào)發(fā)射方和接收方分別在對應(yīng)時(shí)刻的離太陽質(zhì)心的幾何距離，和代表該支路信號(hào)發(fā)射方和接收方在對應(yīng)時(shí)刻的離行星、外行星系統(tǒng)或月球質(zhì)心的幾何距離。其中和代表該支路信號(hào)接收方和發(fā)射方在對應(yīng)時(shí)刻的日心坐標(biāo)向量，和代表該支路信號(hào)接收方和發(fā)射方在對應(yīng)時(shí)刻的行星質(zhì)心坐標(biāo)向量。代表在日心坐標(biāo)系下該支路信號(hào)接收方和發(fā)射方的幾何距離，代表在行星質(zhì)心坐標(biāo)系下該支路信號(hào)接收方和發(fā)射方的幾何距離。

本發(fā)明利用提出的月球著陸器-軌道器的四程中繼跟蹤觀測模型，同時(shí)結(jié)合傳統(tǒng)的雙程、三程測距測速觀測模式，實(shí)現(xiàn)中繼星的精密定軌和遠(yuǎn)月面著陸器的精密定位。具體步驟如下：

根據(jù)式(2)，設(shè)與軌道器相關(guān)的初始狀態(tài)的近似值為：x₀^ref＝[r₀ v₀ p₀ q₀]^T。其中，x₀^ref為隨時(shí)間變化的m維狀態(tài)向量，包含軌道器的初始位置r₀、速度v₀、力模型相關(guān)參數(shù)p₀以及與觀測模型相關(guān)的參數(shù)q₀；其中q₀為遠(yuǎn)月面著陸器的初始近似位置坐標(biāo)同時(shí)，設(shè)置迭代收斂的限差δ，在設(shè)置迭代收斂的限差時(shí)，分別設(shè)置局部參數(shù)和全局參數(shù)迭代收斂的限差δ₁和δ₂。

①根據(jù)軌道器的實(shí)際受力情況，同時(shí)在給定的觀測站的觀測弧段內(nèi)，選定x₀^ref，并將式(1)所示的微分方程組進(jìn)行數(shù)值積分，解算出一條以x₀^ref為初始條件且并不精確的參考軌道x^ref(t_i)，其中x^ref(t_i)代表觀測弧段內(nèi)某時(shí)刻t_i的參考軌道x^ref，i＝1,2,…n。同時(shí)，式(10)所示的也可以由式(1)所示的微分方程組數(shù)值積分得來。具體實(shí)施時(shí)，初始值、微分方程組的具體細(xì)節(jié)以及限差可由本領(lǐng)域技術(shù)人員自行預(yù)先設(shè)定。

②建立月球著陸器-軌道器的四程中繼跟蹤測量模式(4W L-O)的觀測模型

本步驟相當(dāng)于計(jì)算式(6)中的h(x₀^ref)。

四程中繼測距(4W L-O R)：

將下行信號(hào)到達(dá)T_n的時(shí)刻n記錄為協(xié)調(diào)世界時(shí)時(shí)標(biāo)UTC(n)，協(xié)調(diào)世界時(shí)通常簡稱UTC，首先需要將該UTC時(shí)標(biāo)轉(zhuǎn)換為太陽系質(zhì)心力學(xué)時(shí)標(biāo)TDB，可得到下行信號(hào)到達(dá)T_n的時(shí)刻的太陽系質(zhì)心力學(xué)時(shí)標(biāo)TDB(n)。之后通過光行時(shí)解依次得到S_m,L_k,S_j,T_i處對應(yīng)的太陽系質(zhì)心力學(xué)時(shí)TDB(m)，TDB(k)，TDB(j)，TDB(i)，最后再次通過時(shí)間轉(zhuǎn)換得到T_i處的協(xié)調(diào)世界時(shí)時(shí)標(biāo)UTC(i)。在進(jìn)行各個(gè)鏈路的光行時(shí)解算時(shí)，如式(26)和(27)所示，軌道器與地面之間的鏈路考慮廣義相對論時(shí)延項(xiàng)，軌道器與著陸器之間的鏈路由于距離很近，只考慮牛頓光行時(shí)，不考慮廣義相對論時(shí)延項(xiàng)。通過此“回溯”，可以得到相應(yīng)時(shí)刻各個(gè)參與者的狀態(tài)向量。定義X(J)為參與者J在BCRS下的位置矢量，J＝T_n,S_m,L_k,S_j,T_i，那么4個(gè)幾何距離依次為：

R₁＝|X(S_m)-X(T_n)| (28)

R₂＝|X(S_m)-X(L_k)| (29)

R₃＝|X(S_j)-X(L_k)| (30)

R₄＝|X(S_j)-X(T_i)| (31)

以上可作為四程測距測量方程。

至此，四程中繼測距(4W L-O R)的觀測模型建立如下：

其中，RLT_nm和RLT_ij為式(26)和式(27)中的廣義相對論時(shí)延，即RLT_nm為由參與者T_n發(fā)射信號(hào)到S_m的廣義相對論時(shí)延，RLT_ij為由參與者T_i發(fā)射信號(hào)到S_j的廣義相對論時(shí)延；UTC(n)和UTC(i)分別為上述T_n相應(yīng)時(shí)刻n和T_i相應(yīng)時(shí)刻i對應(yīng)的協(xié)調(diào)世界時(shí)時(shí)標(biāo)，TDB(n)和TDB(i)分別為時(shí)刻n和時(shí)刻i對應(yīng)的太陽系質(zhì)心力學(xué)時(shí)，R為四程中繼測距值，即基于月球著陸器-軌道器的四程中繼測距理論計(jì)算值。

四程中繼測速(4W L-O RR)：

基于月球著陸器-軌道器的四程中繼測速理論值按照積分多普勒的形式可表達(dá)為：在一個(gè)多普勒積分周期內(nèi)，起始時(shí)刻T_s和終止時(shí)刻T_e各自對應(yīng)的四程中繼測距值R_s、R_e之差相對于時(shí)間的變化率。首先將積分起始時(shí)刻T_s視為下行信號(hào)的到達(dá)時(shí)刻T_n，按照式(32)所示的觀測模型計(jì)算四程中繼測距值R_s，然后將積分終止時(shí)刻T_e同樣視為下行信號(hào)的到達(dá)時(shí)刻T_n，按照式(32)所示的觀測模型計(jì)算四程中繼測距值R_e。建立四程中繼測速(4W L-O RR)的觀測模型為：

上式中T_c為一個(gè)多普勒積分周期，等于T_e-T_s，RR為四程中繼測速值，即基于月球著陸器-軌道器的四程中繼測速理論值。

③計(jì)算觀測模型所得四程中繼測距觀測量(包括四程中繼測距值R、四程中繼測速值RR)對遠(yuǎn)月面著陸器坐標(biāo)的偏導(dǎo)數(shù)

本步驟相當(dāng)于計(jì)算式(10)中的

四程中繼測距值R對遠(yuǎn)月面著陸器坐標(biāo)X(L_k)的偏導(dǎo)數(shù)H_R：

按照式(34)計(jì)算積分起始時(shí)刻T_s和積分終止時(shí)刻T_e分別對應(yīng)的四程中繼測距值R對遠(yuǎn)月面著陸器坐標(biāo)X(L_k)的偏導(dǎo)數(shù)H_s和H_e，則四程中繼測速值RR對遠(yuǎn)月面著陸器坐標(biāo)X(L_k)的偏導(dǎo)數(shù)H_RR為：

計(jì)算所得偏導(dǎo)數(shù)H_RR即偏導(dǎo)數(shù)

④輸入實(shí)際觀測值序列(即n維觀測向量z)，按照式(9)計(jì)算Δz，并按照式(10)計(jì)算H矩陣。

⑤結(jié)合傳統(tǒng)的雙程/三程測距、測速等測量模型，結(jié)合實(shí)際測量精度，按照式(21)確定權(quán)矩陣W；其中，W也可以由本領(lǐng)域技術(shù)人員根據(jù)不同的算法(如：赫爾模特方差分量估計(jì)法)預(yù)先確定。

⑥生成如式(24)所示的法方程，根據(jù)式(25)解算出改正量Δx₀^lsq，根據(jù)式(20)計(jì)算出最佳估值，即修正后的航天器狀態(tài)x₀^lsq。其中，x₀^lsq向量中包含有遠(yuǎn)月面著陸器的位置坐標(biāo)。

⑦將解算出的x₀^lsq作為①中的x₀^ref，重復(fù)①～⑥，直至兩次得到的Δx₀^lsq之差的絕對值小于預(yù)設(shè)的限差δ時(shí)停止迭代，得到最終的最佳估值x₀^lsq。

本發(fā)明進(jìn)一步提出，在實(shí)際情況中，由于地球測站可視范圍有限，探測器變軌等原因，不能保證永遠(yuǎn)有連續(xù)且有效的觀測跟蹤數(shù)據(jù)，所以探測器軌道跟蹤數(shù)據(jù)文件可分為多個(gè)觀測弧段；同時(shí)，根據(jù)待估參數(shù)本身的性質(zhì)，一般將x₀^ref＝[r₀ v₀ p₀ q₀]^T中的r₀v₀p₀作為局部參數(shù)處理，記作x₁，維數(shù)為m₁，其中局部參數(shù)的特點(diǎn)是已知其每個(gè)弧段的初始近似值x₁^ref(i)，并解算每個(gè)弧段的參數(shù)的改正值Δx₁^lsq(i)，最終得到每個(gè)弧段的參數(shù)的最佳估值x₁^lsq(i)，其中，i＝n,n-1,…,1，代表第i個(gè)弧段；對于參數(shù)q₀，本發(fā)明中指遠(yuǎn)月面著陸器位置坐標(biāo)X(L_k)，將其作為全局參數(shù)處理，記作x₂，維數(shù)為3，其中全局參數(shù)的特點(diǎn)是與弧段無關(guān)，解算時(shí)融合多個(gè)弧段獲得最終的全局參數(shù)的改正值Δx₂^lsq以及最佳估值x₂^lsq。具體闡述如下：

對于每一個(gè)弧段，一般取1-2天的弧長。輸入某個(gè)弧段的初始近似值x₀^ref，綜合大量的四程中繼跟蹤測量值和傳統(tǒng)的雙程/三程測距、測速等測量值，對中繼星(軌道器)進(jìn)行精密定軌，在不考慮參數(shù)先驗(yàn)信息的前提下，按照步驟①～⑥，根據(jù)四程中繼跟蹤測量模式，針對單弧段中繼星精密定軌生成對著陸器坐標(biāo)的偏導(dǎo)數(shù)，生成如式(24)所示的單弧段法方程，為了區(qū)分局部參數(shù)和全局參數(shù)，將法方程矩陣按照局部參數(shù)和全局參數(shù)的維數(shù)分塊，可得如下形式：

Δx₁和Δx₂分別代表待估的局部參數(shù)(包括中繼星初軌、力模型參數(shù)，中繼星初軌此處為位置和速度)的改正值和全局參數(shù)(即遠(yuǎn)月面著陸器的位置坐標(biāo))的改正值，A₁₁對應(yīng)為局部參數(shù)的法方程矩陣，具體形式對應(yīng)式(24)，為m₁×m₁維的矩陣A₂₂對應(yīng)為全局參數(shù)的法方程矩陣，具體形式對應(yīng)式(24)，為3×3的矩陣A₂₁和A₁₂為矩陣分塊后的左下部分和右上部分；b₁和b₂按照式(24)所示對應(yīng)為m₁×1維的矩陣和3×1維的矩陣

為了更有效的解算待估參數(shù)，聯(lián)合多個(gè)弧段按照步驟①～⑥形成法方程，設(shè)有N個(gè)弧段可得到類似式(36)的一個(gè)序列的弧段法方程：

其中上標(biāo)(1)～(N)表示弧段編號(hào)1～N。式(37)中各符號(hào)的具體含義與式(36)相同。

融合多個(gè)弧段的法方程，形成如下融合法方程：

對上式中的聯(lián)合矩陣每個(gè)弧段部分采取矩陣行變換，可得到如下形式：

式(38)和式(39)中的各項(xiàng)含義與式(36)相同，其中上標(biāo)(1)～(N)表示弧段編號(hào)1～N。I代表單位矩陣，從下至上解算式(39)，首先解得全局參數(shù)向量的改正值Δx₂^lsq，即遠(yuǎn)月面著陸器坐標(biāo)的改正值，按照式(20)加到遠(yuǎn)月面著陸器坐標(biāo)的初始近似值x₂^ref得到最優(yōu)估值x₂^lsq。解算出Δx₂^lsq后，對式(39)從下至上依次回代，即可按順序解算出局部參數(shù)向量按照式(20)加到局部參數(shù)初始近似值x₁^ref(n),x₁^ref(n-1),…,x₁^ref(2),x₁^ref(1)，即可得各個(gè)弧段局部參數(shù)對應(yīng)的最優(yōu)估值x₁^lsq(n),x₁^lsq(n-1),…,x₁^lsq(2),x₁^lsq(1)。

由于上述過程基于將非線性的方程線性化，為更準(zhǔn)確的解算參數(shù)，數(shù)學(xué)上通常將線性化后的非線性方程使用迭代法求解。參見圖4，實(shí)施例提供的具體流程闡述如下：

第一，對于第一個(gè)弧段，按照①～⑥形成如式(36)所示的法方程。

第二，對于第2～N個(gè)弧段，按照①～⑥與第一個(gè)弧段按照上一步驟形成的法方程一起形成如式(37)所示的法方程。

第三，按照式(37)和式(38)融合多弧段法方程，形成如式(39)所示的形式；即生成弧段(1)～(N)對中繼星初軌、力模型參數(shù)(局部參數(shù))的偏導(dǎo)數(shù)，生成對著陸器坐標(biāo)(全局參數(shù))的偏導(dǎo)數(shù)，可生成弧段(1)～(N)的法方程，開始融合弧段(1)～(N)的法方程進(jìn)行著陸器坐標(biāo)(全局參數(shù))解算。

第四，解算出全局參數(shù)向量的第一次的改正值Δx₂^lsq，按照式(20)加到遠(yuǎn)月面著陸器坐標(biāo)的初始近似值x₂^ref得到第一次的最優(yōu)估值x₂^lsq。

第五，將當(dāng)前的最優(yōu)估值x₂^lsq重新作為初始近似值x₂^ref，并按照①～⑥的步驟更新第三步中的法方程中的矩陣形成新的多弧段融合的法方程。其中，i＝N,N-1,…1，表示第i個(gè)弧段。第一次執(zhí)行第五步時(shí)，當(dāng)前的最優(yōu)估值即采用第四步中第一次得到的全局參數(shù)的最優(yōu)估值x₂^lsq。

第六，根據(jù)新的法方程解算出新的全局參數(shù)向量的改正值Δx₂^lsq，按照式(20)加到遠(yuǎn)月面著陸器坐標(biāo)的初始近似值x₂^ref得到新的最優(yōu)估值x₂^lsq。

第七，重復(fù)第五步和第六步，直至當(dāng)次迭代和上一次迭代得到的Δx₂^lsq之差的絕對值小于預(yù)設(shè)的相應(yīng)限差δ₂停止迭代，得到最終的全局參數(shù)的最優(yōu)估值x₂^lsq。具體實(shí)施時(shí)，本領(lǐng)域技術(shù)人員可自行根據(jù)精度需要設(shè)置限差δ₂，可將流程設(shè)計(jì)為判斷是否滿足著陸器坐標(biāo)收斂條件，是則完成著陸器坐標(biāo)(全局參數(shù))的最后一次迭代，獲取著陸器坐標(biāo)(全局參數(shù))精密坐標(biāo)，否則返回進(jìn)行全局參數(shù)迭代。

第八，解算出最新的Δx₂^lsq后，按照從下到上的次序回代到由以上步驟形成的最新的多弧段融合的法方程中，形式如式(39)所示，依次解算出第一次的局部參數(shù)的改正值按照式(20)加到局部參數(shù)初始近似值x₁^ref(n),x₁^ref(n-1),…,x₁^ref(2),x₁^ref(1)，即可得各個(gè)弧段局部參數(shù)對應(yīng)的第一次的最優(yōu)估值x₁^lsq(n),x₁^lsq(n-1),…,x₁^lsq(2),x₁^lsq(1)。

第九，將第一次得到的局部參數(shù)的最優(yōu)估值x₁^lsq(i)重新作為初始近似值x₁^ref(i)，按照①～⑥的步驟更新法方程中的矩陣形成新的多弧段融合的法方程。其中，i＝N,N-1,…1，表示第i個(gè)弧段。

第十，根據(jù)新的法方程解算出新的局部參數(shù)向量的改正值按照式(20)加到各個(gè)弧段局部參數(shù)的初始近似值x₁^ref(i)得到新的最優(yōu)估值x₁^lsq(i)，其中，i＝N,N-1,…1，表示第i個(gè)弧段。

第十一，重復(fù)第九和第十步，直至當(dāng)次迭代和上次迭代得到的之差的絕對值小于預(yù)設(shè)的相應(yīng)限差δ₁停止迭代，得到最終局部參數(shù)的最優(yōu)估值x₁^lsq(i)，其中，i＝N,N-1,…1，表示第i個(gè)弧段。具體實(shí)施時(shí)，本領(lǐng)域技術(shù)人員可自行根據(jù)精度需要設(shè)置限差δ₁，可將流程設(shè)計(jì)為判斷是否滿足局部參數(shù)收斂條件，是則完成中繼星初軌、力模型參數(shù)(局部參數(shù))的最后一次迭代，獲取中繼星初軌、力模型參數(shù)估值，否則返回進(jìn)行局部參數(shù)迭代。

至此，根據(jù)以上流程，得到最終穩(wěn)定的中繼星精密軌道和遠(yuǎn)月面著陸器的精密坐標(biāo)。具體實(shí)施時(shí)，可采用計(jì)算機(jī)軟件技術(shù)實(shí)現(xiàn)自動(dòng)流程運(yùn)行。

以上所述，僅為本發(fā)明較佳的具體實(shí)施方式，本發(fā)明的保護(hù)范圍不限于此，任何熟悉本技術(shù)領(lǐng)域的技術(shù)人員在本發(fā)明披露的技術(shù)范圍內(nèi)，可顯而易見地得到的技術(shù)方案的簡單變化或等效替換均落入本發(fā)明的保護(hù)范圍內(nèi)。