一種導航衛星實時精密定軌的方法

專利名稱：一種導航衛星實時精密定軌的方法
技術領域：
本發明涉及一種衛星導航技術，具體的說是一種導航衛星實時精密定軌的方法。
背景技術：
導航衛星精密定軌作為衛星導航系統的關鍵技術，一直受到國內外眾多研究 機構的重視，而成為相關研究領域的前沿課題。自從1994年國際GNSS地球動力 學服務組織成立以來，精密GPS軌道與鐘差產品在精度、時間延遲上不斷的改進。 目前，IGS提供的事后精密星歷軌道精度優于5cm，鐘差精度優于0.1ns，可滿 足高精度用戶的事后定位需求。隨著精密單點定位技術應用的廣泛深入，已有越 來越多的GPS用戶對導航星座軌道的實時性和精度提出了更高要求。為促進GPS 高精度實時應用領域的發展，IGS組織、美國的JPL、歐洲的BKG、 GFZ、 ESA等 知名研究機構都在開展相關研究。
導航衛星實時軌道改進最早在廣域差分系統中提出，早期通過幾何法對廣播 星歷計算的衛星位置進行實時定軌，雖然此方法簡單、易于實現，也能夠在一定 程度上改善軌道的精度，但由于沒有考慮衛星軌道的運動特性，其噪聲大而且易 受觀測的幾何形狀和數據丟失的影響，以致軌道預報的精度難以保證，因此，我 們對于高精度的導航衛星實時定軌應充分考慮衛星的動力學信息。導航衛星定軌 是一個大型復雜的數據處理系統，涉及到復雜的動力學模型和觀測模型。由于導 航衛星實時精密軌道生成過程中數據傳輸、數據處理以及實時軌道產品發播都存 在時間延遲，實時用戶需求的衛星位置一般通過實時更新的精密衛星軌道初值和 力模型參數向外積分衛星運動方程得到。因而，實時軌道的精度與衛星軌道初值 和力模型參數精度、軌道預報模型精度、外推時間長短三者相關。又由于GNSS 衛星軌道預報模型及方法的改進空間很小，因此導航衛星實時定軌技術的關鍵在 于快速精確地更新軌道初值與力模型參數。
目前，實時精密定軌參數一般采用線性濾波的方法，并且逐個歷元更新衛星 軌道。對于高度在20000多公里的導航衛星，在非機動及無故障的正常情況下， 衛星軌道變化是平滑的。因此，短時間的軌道預報精度很高。實踐表明，基于l 天對IGS事后精密軌道進行擬合結果并外推1小時軌道的精度達到5cm左右，與 估計軌道精度基本相當。在采用逐歷元濾波過程中，為提高解算速度以滿足實時 需求，系統中估計參數的個數與觀測值的個數要求盡量少，這就會影響到系統模型以及參數解的精度。因而，對于實時定軌而言，是否有必要采用逐歷元估計的 方法是值得進一步研究的。
一般而言，濾波方法不可避免地存在初始化、數值穩定性、數據飽和以及處 理速度等棘手問題。初始化是指觀測數據累積到一定的時間跨度后參數估計才能 提供穩定可靠的結果。在實時定軌中這個問題非常突出，因為初始化觀測時間要 大于12小時，初始化所需要的時間至少等于處理12小時數據的計算時間。同時 濾波中穩健的實時數據質量控制是非常困難的問題，由此產生的濾波發散將導致 估計系統崩潰，系統需要重新初始化。
美國JPL的RTG軟件采用平方根信息濾波的方法實現實時精密定軌，也同樣 存在以上這些方面的問題，其參與實時解算的測站數量在70個左右，定軌精度 在20cm左右，甚至差于IGS的快速軌道IGU預報部分精度；并且實時定軌系統 需要多個備份系統支持，確保在處理系統崩潰的狀態下迅速恢復提供用戶實時軌 道。

發明內容
本發明的目的就是為了克服上述背景技術的不足之處，提供一種導航衛星實 時精密定軌的方法，它能有效解決實時定軌中初始化、數值穩定性、數據飽和以 及處理速度等方面的問題。
為了實現上述目的，本發明提供的一種導航衛星精密定軌的方法，該方法包
括以下步驟
(1) 將衛星軌道按觀測時間分割成連續等時長的短弧段，設定每個短弧段
的時長為T小時，設定處理軌道數據的滑動數據窗口長度為N個短弧段時長為
(2) 處理當前短弧段，即滑動數據窗口內最后一個短弧段上獲得的實時數 據，生成該弧段上的短弧法方程，其短弧法方程的公式為a^々&,=『i4，其中^a 為該弧段的軌道初始位置與力模型初值相應的改正參數，iV^為法方程矩陣，k 為短弧段標記；
(3) 將滑動數據窗口內前N-1個短弧段生成的短弧法方程進行合并；其短
弧合并利用的公式如下
兩個相連弧段的軌道合并公式
<formula>formula see original document page 5</formula>W個相連弧段的軌道合并公式
<formula>formula see original document page 6</formula> 式2
其中，O(UJ是狀態轉移矩陣；(D描述了函數關系模型上的連接關系，Pw 描述了隨機模型統計特性上的連接，體現了前后弧段對應參數映射關系的約束程
度；
根據短弧軌道間參數類型的不同相連弧段間連接時①與&的設置分為以下 三類，
(a) 對于前后弧段相同或具有確定性的狀態轉移關系的參數，在軌道連接 處是相同的，設尺為無窮大實現，
(b) 對于前后兩個弧段不相關的分段參數或白噪聲參數，設0) = 0， Pw=0
實現，
(c) 對于隨機過程噪聲參數，通過其過程噪聲的標準差給定相應的^值實
現；
(4)將當前短弧段的法方程與由前N-1個短弧段法方程合并成的一個法方 程進行合并；
(5) 求解由整個滑動數據窗口內短弧合并而成的法方程，得到衛星軌道初 值和力模型參數，通過軌道積分得到下一時段的預報軌道，作為實時軌道發送給 用戶；
(6) 將滑動窗口向下移動一個短弧段，返回步驟(2)。在上述技術方案中，步驟(1)中所述的滑動數據窗口的長度即N*T的值 不超過24小時。
在上述技術方案中，步驟(2)中所述的短弧段上獲得的實時數據包括衛星 初始狀態、力模型參數、地球自轉參數、大氣參數、測站坐標。
在上述技術方案中，步驟(3)中所述的短弧合并公式推導過程如下
(A)衛星運動狀態連接，在導航衛星運動的實際過程中，在正常情況下， 衛星運動軌跡和速度是連續的。如圖2所示，在相鄰的兩個短弧的交界處即同一 時刻，衛星的狀態與受到的攝動力是相同的，可用以下的關聯條件表示
式3
其中，^(U、 4(4+,)、 ^是第A弧段^,時刻的衛星位置、速度與力模型參數； 4+1ft+1)、 4+1(U、《w是第A: + l弧段L,時刻的衛星位置、速度與力模型參數。
將上述關聯條件通過兩個相連弧段各自的軌道初值與力模型參數可表示為 如下式4， 5。
<formula>formula see original document page 7</formula>
其中，
^U。、 4+1(U。、 &+1(4+1)。為"1弧段4+1時刻的軌道與力模型參數初值； &U、 &U、 、+1("+1)為"1弧段~+1時刻軌道初始位置與力模型初
值相應的改正參數；
"U。、 4(U。、 A(U。為A弧段U寸刻軌道與力模型參數初值；
5 ("，54ft) ， ^Aft)為6弧段4時刻軌道初始位置與力模型初值相應的改正參數；
O(UJ是狀態轉移矩陣。
(B) 觀測模型參數連接，與衛星連續運動相同，導航衛星的觀測在正常情 況下也是連續的，由此相連兩個弧段的觀測模型參數也是關聯的。觀測模型參數 包括對流層參數、衛星鐘差參數、接收機參數、地球自轉參數、地面站坐標、模 糊度參數等，這些待估參數可歸類為隨機過程噪聲參數與狀態參數。對于狀態參 數同樣可以建立與上述衛星運動狀態連接類似的關聯方程，只是狀態轉移矩陣具 有不同的定義。對于隨機過程噪聲參數， 一般可以通過以下公式表示
HA+^ 式6 由于v^可以認為是一獨立隨機過程噪聲，并服從正態分布，上式可以轉換為
<formula>formula see original document page 8</formula> 式7
其中Q為相應過程噪聲的標準差，考慮到^為零均值，^A可設為零，A、 A+1 為隨機過程參數，M^^為映射矩陣，這樣隨機過程噪聲參數與狀態參數一樣建 立類似的關聯方程。
(C) 由此將上述兩種數學模型統一起來表述，將不同的參數都用義表示，這 樣關聯條件則表示為
<formula>formula see original document page 8</formula> 式8
約束方程表述為
<formula>formula see original document page 8</formula>
式9
對于短弧軌道間不同參數的連接或合并，可以通過設置參數前后弧段的映射 因子o與動態噪聲&來實現。實際上，o描述了函數關系模型上的連接關系，A
描述了隨機模型統計特性上的連接，體現了前后弧段對應參數映射關系的約束程 度，況為觀測值殘差，以下討論不同類型參數相連弧段間連接時O與A的設置。
對于前后弧段相同或具有確定性的狀態轉移關系的參數，如衛星狀態參數、測站坐標參數、連續的模糊度參數等，在軌道連接處是相同的，可以設^為無 窮大實現；
對于前后兩個弧段不相關的分段參數或白噪聲參數，如衛星鐘差、接收機鐘 差等，可以設0 = 0， &=0實現；
對于隨機過程噪聲參數，可以通過其過程噪聲的標準差給定相應的&值實現。
(D)基于相連兩個弧段各自的短弧法方程信息，與邊界約束條件方程，推 導短弧合并的統一數學模型。
假設相連弧段+1對應的法方程為
<formula>formula see original document page 9</formula>式10
邊界約束條件方程為
<formula>formula see original document page 9</formula> 式11
將上式io假設為先驗約束，式ii作為虛擬觀測方程，則軌道合并的問題就 可以轉換為附有先驗信息的平差問題，采用最小方差估計，要求滿足性能函數
J(x)為最小
<formula>formula see original document page 9</formula> 式12
則^M-0可推導得到 <formula>formula see original document page 9</formula>簡化表達，設<formula>formula see original document page 9</formula>上式13可簡化為:<formula>formula see original document page 10</formula>式1
以上推導了兩個相連弧段的軌道合并公式，下面直接給出一般性的w個弧段
的軌道合并公式
<formula>formula see original document page 10</formula>
式2
相比較現有技術而言，運用本方法實現導航衛星實時精密定軌具有如下優 點 一、處理速度快，由于法方程合并中大量的參數都可以預先消去，如除了最 后一個弧段外的其它弧段上的軌道參數、模糊度參數等，這樣，合并的計算量就 可以大大減小，從而提高數據處理速度。二、快速初始化， 一旦實時系統需要重 新啟動，此方法只需處理當前短弧段上的觀測數據，然后與己經存在的法方程合 并就可以實時更新軌道，從而實現快速初始化。三、數值穩定，由于使用滑動數 據窗口的方法，每次軌道更新只用當前窗口中的數據，陳舊的觀測信息被去除， 從而可以完全克服"數據飽和"導致的動力學模型參數估計精度降低甚至錯誤的 問題。四、能有效實現數據質量的控制，大大提高實時定軌的精度。


圖1為本發明一種導航衛星實時精密定軌的方法中滑動數據窗口定軌示意
圖2為本發明方法中短弧連接示意圖。
圖3為本發明方法中實施例試算每小時GPS實時軌道與IGS事后軌道比較 結果圖。
具體實施方式
下面結合附圖對本發明作進一步的描述。
如圖1所示， 一種導航衛星實時精密定軌的方法，主要是將實時定軌通過 建立短弧法方程和短弧法方程合并兩個并行進程實現，該方法包括以下步驟
(1) 將衛星軌道按觀測時間分割成連續等時長的短弧段，設定每個短弧段 的時長為T小時，設定處理軌道數據的滑動數據窗口長度為N個短弧段時長為 N+T小時，其N+T的值不超過24小時；
(2) 處理當前短弧段，即滑動數據窗口內最后一個短弧段上獲得的實時數
據，建立該弧段上包括衛星初始狀態、力模型參數、地球自轉參數、大氣參數、
測站坐標的短弧法方程，其短弧法方程的公式為7V^&^ =『^，其中5A為該弧 段的軌道初始位置與力模型初值相應的改正參數，7V^為法方程矩陣，k為短弧 段標記；
(3)短弧法方程合并進程將前N-1個弧段生成的法方程進行合并；其短弧 合并利用的公式如下 兩個相連弧段的軌道合并公式
—①^
一尸O
「化-
式i
"個弧段的軌道合并公式
-尸①
(2,1)<formula>formula see original document page 11</formula>
式2
其中，O(UJ是狀態轉移矩陣，①描述了函數關系模型上的連接關系，Pw 描述了隨機模型統計特性上的連接，體現了前后弧段對應參數映射關系的約束程度；
根據短弧軌道間參數類型的不同相連弧段間違接時0與^的設置分為以下 三類，
(a) 對于前后弧段相同或具有確定性的狀態轉移關系的參數，如衛星狀態 參數、測站坐標參數、連續的模糊度參數等，在軌道連接處是相同的，設尺為 無窮大實現；
(b) 對于前后兩個弧段不相關的分段參數或白噪聲參數，如衛星鐘差、接 收機鐘差等，設0 = 0， &=0實現；
(c) 對于隨機過程噪聲參數，通過其過程噪聲的標準差給定相應的&值實
現；
(4) 將當前弧段的法方程與由前N-1個弧段法方程合并成的一個法方程合
并；
(5) 求解由整個滑動窗口中短弧合并而成的法方程，得到衛星軌道初值和 力模型參數，通過軌道積分得到下一時段的預報軌道，作為實時軌道發送給用戶；
(6) 將滑動窗口向下移動一個弧段，返回步驟(2)。 本實施例采用全球70個IGS站2006年197天至203天一周的觀測數據，
所選參考站均提供1小時的準實時觀測數據文件，進行每小時更新一次的實時精 密軌道試算。
具體實施方式
如下通過最小二乘方法對每個小時的觀測文件形成 法方程，實時定軌的滑動窗口長度設為24小時，前23小時為先驗信息，第24 小時實時估計軌道參數，基于估計的實時軌道參數外推1小時的軌道作為提供用 戶的實時軌道，以下分析中的實時軌道均為此外推軌道，從197天開始往前遞推 獲得198至203天的實時軌道，197天軌道作為實時系統啟動的初始信息。
圖3給出了以每小時為單位6天實時軌道與IGS事后精密軌道比較統計的 平均RMS值的時間序列，從圖中可以看出90%以上的RMS值小于6cm， 98% 以上的RMS值小于8cm，最大RMS值在10cm以內，絕大部分衛星平均RMS 值在4-8cm之間，由此試算結果表明使用本方法每小時更新的實時定軌與IGS 事后精密軌道比較，精度達到5cm左右。
權利要求
1. 一種導航衛星實時精密定軌的方法，其特征在于該方法包括以下步驟(1)將衛星軌道按觀測時間分割成連續等時長的短弧段，設定每個短弧段的時長為T小時，設定處理軌道數據的滑動數據窗口長度為N個短弧段時長為N*T小時；(2)處理當前短弧段，即滑動數據窗口內最后一個短弧段上獲得的實時數據，生成該弧段上的短弧法方程，其短弧法方程的公式為Nk，kδxk＝Wkk，其中δxk為該弧段的軌道初始位置與力模型初值相應的改正參數，Nk，k為法方程矩陣，k為短弧段標記；(3)將滑動數據窗口內前N-1個短弧段生成的短弧法方程進行合并；其短弧合并利用的公式如下兩個相連弧段的軌道合并公式式1n個相連弧段的軌道合并公式 式2其中，Φ(tk+1，tk)是狀態轉移矩陣；Φ描述了函數關系模型上的連接關系，Pw描述了隨機模型統計特性上的連接，體現了前后弧段對應參數映射關系的約束程度；根據短弧軌道間參數類型的不同相連弧段間連接時Φ與Pw的設置分為以下三類，(a)對于前后弧段相同或具有確定性的狀態轉移關系的參數，在軌道連接處是相同的，設Pw為無窮大實現，(b)對于前后兩個弧段不相關的分段參數或白噪聲參數，設Φ＝0，Pw＝0實現，(c)對于隨機過程噪聲參數，通過其過程噪聲的標準差給定相應的Pw值實現；(4)將當前短弧段的法方程與由前N-1個短弧段法方程合并成的一個法方程進行合并；(5)求解由整個滑動數據窗口內短弧合并而成的法方程，得到衛星軌道初值和力模型參數，通過軌道積分得到下一時段的預報軌道，作為實時軌道發送給用戶；(6)將滑動窗口向下移動一個短弧段，返回步驟(2)。
2. 根據權利要求1所述的一種導航衛星實時精密定軌的方法，其特征是 步驟(1)中所述的滑動數據窗口的長度即N^^T的值不超過24小時。
3. 根據權利要求1所述的一種導航衛星實時精密定軌的方法，其特征是 步驟(2)中所述的短弧段上獲得的實時數據包括衛星初始狀態、力模型參數、 地球自轉參數、大氣參數、測站坐標。
全文摘要
本發明公開了一種基于滑動數據窗口的導航衛星實時精密定軌的方法，該方法的主要思想是在對狀態X_k，包括衛星初始狀態、力模型參數、地球自轉參數、大氣參數、測站坐標等進行估計時，只利用離k時刻最近的前N個量測，而把之前其余的量測完全甩掉，這里N是預先規定的滑動窗口長度，再通過建立，合并，求解法方程等步驟獲得精確的第N弧段的衛星軌道初值和力模型參數，通過軌道積分得到下一時段的預報軌道，作為實時軌道發送給用戶。采用本方法進行衛星定軌具有數據處理速度快，軌道精度高，數值穩定等優點。
文檔編號G01S1/04GK101435863SQ20081023742
公開日2009年5月20日 申請日期2008年12月25日 優先權日2008年12月25日
發明者劉經南, 唐衛明, 闖 施, 樓益棟, 葛茂榮, 趙齊樂 申請人:武漢大學