單站區域實時電離層異常改正擬合模型研究*

孫 帆，寧一鵬，邢建平*，王 森，代培培

(1.山東大學微電子學院，山東 濟南 250100；2.山東建筑大學測繪地理信息學院，山東 濟南 250101)

電離層位于地面上空，影響頻率信號傳播，是導航定位結果出現偏差的重要原因，制約著導航定位精度，對天基、地基系統的穩定運行及可靠性，航天器在電離層的運行安全，無線電波的傳輸，以及全球導航衛星系統具有很大影響，例如我國在電離層高發期多次觀測到電離層閃爍導致全球衛星導航系統(Global NavigaTIon Satellite System，GNSS)信號失鎖[1]，鑒于此，研究精確度高、實時性優良的電離層模型在民用、軍事等領域的應用非常必要且具有重要的作用。

目前電離層模型研究進展可分為三類。第一類為建立電離層延遲經驗模型，如國際GPS 服務(International GNSS Service，IGS)提供的全球電離層模型Klobuchar[2]、Abdus Salam 發布的NeQuick 模型[3]等經驗模型，第二類為計算得到的數學擬合模型，如文獻[4-5]分別基于球諧函數和球冠諧函數擬合方法建立區域電離層(Total Electron Content，TEC)模型；周晨等[6]采用雙頻載波偽距平滑方法建立動態電離層模型；第三類為利用神經網絡、時間序列等方法得到預測模型，如陳春，吳振森等[7]提出通過神經網絡ANN 方法提前1 h 預測電離層TEC 模型，謝劭峰等[8]利用IGS 中心提供的數據采用Holt-Winters’加法模型和乘法模型建立TEC 預報模型。

以上為傳統電離層建模方法，如文獻[9]依據最新發布數據建立格網點利用時間序列法得到所求位置點具有的電離層數據發布延遲誤差，誤差較大且不能反映電離層異常時的情況；文獻[10]利用偽距解算得到實時電離層延遲，受限于所用衛星偽距數據質量及數目，且高精度計算方法需考慮整周模糊度的消除，計算較為繁瑣。

本文將發布的格網數據進行預測建模，結合多頻載波偽距法經電離層異常條件判定改正，建立改正擬合模型。通過預測模型降低因數據質量不良造成的偽距模型誤差，通過實時偽距模型提高預測模型與時間的關聯度，經異常判定建立模型，有效提高模型實時性及模型擬合精度。本文對比分析建立模型與多頻載波偽距法、Holt-Winters’對比，驗證模型的正確性及可行性。

1 電離層模型

1.1 單站區域電離層球諧級數模型

1.1.1 載波相位平滑偽距原理

基于GNSS 觀測數據的偽距觀測值建立電離層模型，針對偽距法建立TEC 模型精確度低的問題，采用載波相位平滑偽距法可有效平滑偽距誤差，削弱觀測噪聲，提高偽距觀測值精度。GNSS 系統[11]雙頻偽距觀測方程和載波相位觀測方程如下:

式中:j表示衛星編號，k表示頻率編號，ρ表示衛星到觀測站接收機的偽距觀測值，L表示載波相位觀測值，c表示光速，dt1表示接收機鐘差，dt2表示衛星鐘差，T與I分別表示信號傳播路徑對流層與電離層延遲誤差，br 和bs 分別表示接收機與衛星硬件延遲誤差；ε表示偽距碼觀測噪聲，λ表示波長，N表示載波相位整周模糊度[12]。

式中:f1和f2分別表示和Lj2 載波頻率，由于頻間差分碼偏差(Differential Code Bias，DCB)參數通常一天之內幾乎無變化，代入IGS 發布的DCB 文件中觀測站接收機和衛星的頻間DCB1、DCB2 參數值，得到信號路徑上的TEC 值(Slant Total Electron Content，STEC)[13]。

1.1.2 穿刺點經緯度計算

式中:φPP、λPP表示衛星到觀測站的穿刺點坐標，φ、λ表示觀測點坐標，H表示電離層高度值(一般情況為300 km～450 km，本文選用450 km)，RE表示地球平均半徑，為6 371 km，El、Az 分別表示衛星相對于觀測站的高度角和方位角，z和z＇分別表示觀測站和穿刺點的天頂角[14]。

圖1 所示為穿刺點原理圖，可見衛星穿刺點坐標與地面坐標正上空不是同一個點。

圖1 穿刺點原理圖

1.1.3 VTEC 計算

擬合兩個組合得到的測點上空電離層(Vertical Total Electron Content，VTEC)值，若雙頻時P1＝1，P2＝0，或者P1＝0，P2＝1，通過計算得到的穿刺點經緯度坐標，得到投影函數，從而計算得到VTEC[12]。

將STEC 通過單層模型投影函數轉化為觀測站上空垂直電離層VTEC，轉化公式為:

式中:F(z)表示穿刺點處的投影函數，z表示接收器天頂角，RE表示地球半徑，Hion表示薄層高度。

1.1.4 球諧函數模型

經過式(5)計算后，將不同VTEC 值及其穿刺點代入到電離層球諧函數模型，表示如下:

式中:Cnm，Snm為球諧函數的參數，Pnm(sinβ)為締合Legendre 函數，nmax為模型階數，ms為過穿刺點的經線與過地心和太陽連線的經線構成的夾角，β為地磁緯度。

1.2 Holt-Winters’加法預測模型

選取時間序列法中的周期季節性時間序列法更符合電離層變化趨勢，提高電離層預測模型準確度。使用歐洲定軌中心(CODE)發布的前10 天數據，預測得到當天數據，構建電離層預測模型。

加法模型初始值計算公式如下:

式中:Xt表示t時刻觀測值，It表示t時刻季節成分，bi表示t時刻趨勢成分，Ft+m表示m時刻趨勢成分，m表示預測的時刻數，L表示季節長度，α、β、γ表示平滑參數[9]。

1.3 擬合改正模型

1.3.1 初步擬合模型

基于雙頻載波相位平滑偽距法得到的球諧函數模型受數據質量影響，精確度不穩定，圍繞真實值跳動較大，造成誤差較大問題；依據CODE 發布數據進行的時間序列法得到的預測值模型，因其實時性關聯不夠緊密，存在延遲問題，造成預測模型雖然整體符合真實值變化，但當電離層出現異常擾動變化時，該模型無法及時預測，從而產生誤差?；诖耍疚膶⒍哌M行改正擬合，使其兩種方法優勢互補，彌補傳統模型的不足。

通過加權最小二乘法對10 天內3 萬歷元數據處理，初步得到確定系數a，b。

1.3.2 電離層異常判定標準

初步擬合模型在精準度上較兩種傳統模型略有提升，但在電離層出現異常擾動情況時，受時間序列法原理制約，誤差仍然較大，拉低初步擬合模型整體精準度，針對這一問題，本文使用電離層異常判斷標準判定是否發生異常，為后面異常情況下改正擬合提供依據。

電離層發生異常時，增加使用載波相位平滑偽距法計算得到的球諧函數模型在擬合模型中占比，降低預測得到的時間序列模型占比。異常符合判斷時刻前后5 天的VTEC 值服從均值為μ標準差為σ的正態分布[15]，異常判定上下界表示如下:

C1iup表示判定上界，C2idown表示判定下界，下標i表示時刻，VTECmedian表示預測值中值，VTEC(Projection)表示預測模型的VTEC。

2 算例分析

2.1 建模方法

考慮太陽磁暴和地理位置，本文選取MIZU 參考站，時間為2013 年1 月10 日(BRDC0100)到20日(BRDC0200)，中心經緯度范圍為10°×10°以內，半徑＜600 km[9]的GPS 數據作為分析數據構建模型。使用載波相位平滑偽距法提高偽距觀測精度，衛星截止仰角設為15°，投影函數采用標準投影函數，時間采樣率為30 s，電離層薄層高度設置為450 km[16]。對比分析載波相位平滑偽距法得到的球諧函數模型[11]，利用CODE 發布數據構建的周期季節性時間序列法Holt-Winters 加法預測模型，以及電離層異常改正擬合后的建立模型進行區域VTEC 建模，以0.15°間隔劃分網格點，求解球諧級函數得到模型系數，反代入網格點，計算得出各網格點處的VTEC 值，預測值法使用CODE 發布的數據利用Holt-Winters 預測得到各網格點處的VTEC 值，基于以上兩種模型進行最小二乘法擬合得到初步擬合模型，在初步擬合模型基礎上通過異常判斷改正異常點，構建改正擬合模型。

2.2 結果與分析

2.2.1 初步加權擬合

通過最小二乘法如表1 所示得到未發生異常時刻時載波相位平滑偽距法與預測模型的加權系數，當a＝0.642，b＝0.358 時，擬合模型最貼近真實值。

表1 時間序列法模型系數及其對應誤差

2.2.2 異常判定

由于觀測站所處區域在1 月10 日至20 日均無太陽磁暴現象，但17/18/20 均有小規模地震發生，考慮其對電離層擾動影響，分別分析地震前后VTEC 值，取當前時刻前后5 天的預測模型VTEC(或foF2)中值和標準差作為背景值，以2 倍標準差作為判斷標準，如果載波相位平滑偽距法的VTEC計算值超出背景值上下邊界，則認為該時刻電離層出現異常，經電離層異常閾值判定，增加載波相位平滑偽距法所占權重，降低預測模型所占權重，通過最小二乘法修正異常時刻系數a＝0.67，b＝0.33，異常修正后均方差(Mean Square Error，RMSE)提高0.5。

如圖2 所示，c1、c2為異常判定的上界和下界，y3計算值為載波相位平滑偽距法建立的球諧級函數模型得到的測站點上空VTEC 值，當計算值超出上下閾值邊界時，判定電離層發生異常擾動。對此時的加權擬合值進行重新加權，增大計算值的所占比重。

圖2 計算值上下邊界異常判定

2.2.3 異常判定改正前后的模型與真實值對比

相比初步擬合模型，改正加權系數比重后，如圖3所示，在28-34、53-56、148-153、155-167 時刻精度均有提高，尤其在53-56 時刻提升最為明顯。

圖3 異常點處改正前后與真實值的對比

2.2.4 模型比較

如圖4 所示，y1代表VTEC 真實值，y2代表VTEC初步擬合值。預測模型在192 個時刻內59 個時刻精度較高，相比于真實值，誤差低于0.2VTEC；32 個時刻預測值與真實值相差較大，誤差高于2VTEC，其他時刻計算值與真實值誤差均保持在(0.2，2)以內，整體日RMSE 為22.441。由此可知，預測模型總體滿足2VTEC 精度，且精度較高的時刻較多，但某些時刻存在較大波動，且當遇到由于電離層擾動造成的真實值突變時，預測模型往往無法有效預測，使得預測模型產生較大誤差，若能有效探測電離層是否發生異常，并進行改正，可大大提高預測模型精度。

圖4 預測值與真實值對比

如圖5 所示，載波相位平滑偽距法得到的球諧函數模型在192 個時刻內共39 個時刻精度較高，相比于真實值，誤差低于0.2VTEC；10 個時刻預測值與真實值相差較大，誤差高于2VTEC，其他時刻預測值與真實值誤差均保持在(0.2，2)以內，日RMSE在14.618。由此可知，球諧函數模型，誤差較大時刻遠遠少于預測值模型，整體可以描述真實值變化情況。但受數據影響造，與真實值間的差值波動較大，高精度時刻較少。

圖5 球諧模型值與真實值

如圖6 所示，初步擬合值模型在192 個時刻內41 個時刻精度較高，相比于真實值，誤差低于0.2VTEC，3 個時刻精度較低，與真實值誤差高于2VTEC，其他時刻預測值與真實值誤差均保持在(0.2，2)以內，整體日RMSE 為9.672。由此可知，初步擬合模型，中和了以上兩種模型的優點，彌補了兩種模型的缺點，與預測值模型相比，低精度時刻變少，與球諧函數模型相比高精度時刻增多，并且總體擬合精度提高33.83%。

圖6 初步擬合模型與真實值對比

如圖7 所示，異常判斷改正后的加權擬合模型在192 個時刻內61 個時刻精度較高，相比于真實值，誤差低于0.2VTEC，5 個時刻精度較低，相比于真實值，誤差高于2VTEC，其他時刻預測值與真實值誤差均保持在(0.2，2)以內，日RMSE 在9.111。由此可見，改正后的擬合模型更好地將以上兩種模型優點繼承下來，在保留高精度時刻不變甚至略有增加的情況下，低精度時刻略有降低，且總體擬合模型精度提高了37.67%。

圖7 異常判定后擬合模型與真實值對比

表2 顯示加權擬合后的模型在精度上優于僅使用載波相位平滑偽距以及季節性時間序列法建立的球諧模型以及預測模型；隨著利用數據的增多，載波相位平滑偽距法和季節性時間序列法的日RMSE 值均在降低，而異常改正后的擬合模型在精度上更高，且在發生地震等影響電離層，使其發生異常擾動時，該種方法精確度相比較為判斷改正的精度更高。

表2 10 日內4 種方法與真實值的均方根誤差

測站及其周邊接收機對20 日24 h 內上空經度135°到145°，北緯35°到45°的電離層誤差累積渲染圖像如圖8～圖11 所示，初步加權擬合后的模型在精度上已經在大部分區域優于僅使用載波相位平滑偽距以及季節性時間序列法；而經過異常改正后的擬合模型更是對初步擬合模型誤差較大區域進行了明顯的改善。

圖8 預測值模型

圖9 球諧函數模型

圖10 初步擬合模型

圖11 異常判定后擬合模型

3 分析與總結

本文提出了一種電離層異常改正擬合模型，對比分析區域電離層建模的2 種傳統模型，選取10 天內3 萬個歷元數據作為分析數據，以0.15°為單位的網格點地圖進行日誤差渲染，研究結果如下:

①載波相位平滑偽距法構建球諧函數模型與預測模型經異常判定條件后，建立的改正擬合模型，可有效提高預測模型實時性，減輕載波相位平滑偽距法構建球諧函數模型的不穩定性，提高電離層模型精度，支撐電離層建模研究。

②文中的試驗是基于日本測站MIZU 的周邊區域10 天內的數據來分析驗證的，沒有對其他不同地區測站和更長時間段的區域進行分析。并且該種模型加權改正系數受地理分布位置、磁暴、地震和海嘯等電離層擾動發生頻率規模差異和觀測數據質量和密集程度等多種因素的影響，未來為了達到一個更為詳細的試驗結論和更精細的模型，需要在時間長度上、電離層異常擾動的標準、不同測站地區上進行更加完善的試驗。