歐 明，甄衛(wèi)民，徐繼生，劉 鈍，於 曉

(1.武漢大學電子信息學院,湖北 武漢 430079;2.中國電波傳播研究所,山東 青島 266107)

0 引 言

電離層層析成像(CIT)是隨著衛(wèi)星無線電探測技術(shù)的發(fā)展而興起的一種電離層測量技術(shù),它利用星載無線電信標設(shè)備發(fā)射的無線電波從外部掃描電離層,根據(jù)電離層對波作用的積分效應(低維投影效應)確定目標介質(zhì)的高維分布。1986年Austen等[1]在國際上首次提出CIT的設(shè)想,利用該設(shè)想并結(jié)合極軌衛(wèi)星的觀測數(shù)據(jù),可重構(gòu)出電離層電子密度的二維分布。由于電離層層析技術(shù)適合監(jiān)測大尺度的電離層電子密度的時空分布,且建造和運行費用低廉,這一新設(shè)想立即引發(fā)了國際上相繼開展電離層CT成像的實驗和理論研究[2-3]。

隨著上世紀90年代以來全球衛(wèi)星導航系統(tǒng)的興起,基于全球定位系統(tǒng)GPS的電離層層析技術(shù)逐漸發(fā)展起來。自從Kunitsyn等[4]在國際上首次證實了基于高軌衛(wèi)星系統(tǒng)的電離層層析技術(shù)的可行性以來,許多研究者開展了基于GPS觀測的電離層層析成像技術(shù)研究。Hansen等[5]最早給出了函數(shù)基電離層層析公式,并利用隨機反演方法和WAAS系統(tǒng)的GPS數(shù)據(jù)反演獲得電離層電子密度分布。Ruffini等[6]利用相關(guān)函數(shù)法對全球電離層進行了層析成像,Howe等[7]借助于Kalman濾波方法仿真重構(gòu)了四維電離層結(jié)構(gòu)。Yizengaw等[8]利用澳大利亞140°E子午線附近五個IGS站的實測數(shù)據(jù)進行了電離層反演研究。在國內(nèi),徐繼生[9]、鄒玉華[10]、聞德保[11],樂新安[12]等也開展了相關(guān)研究并取得許多有意義的成果。

由于地基GPS電離層層析成像觀測本身的局限,如信號傳播路徑主要是垂直方向、觀測視角有限等,基于地基GPS的CIT只有融合一定的先驗知識才能獲得唯一解。其中常用的辦法是直接引入背景電離層模型參數(shù)作為電子密度先驗信息對成像算法進行約束,再利用行作用技術(shù)(RAM)對方程進行迭代求解,以獲取一個穩(wěn)定的電子密度解[9],但此時電離層層析成像的垂直分辨率依然受限于背景電離層模型的精度,層析成像反演得到的電子密度峰值高度與背景模型相比并無明顯改善。參考文獻[13]開始探討利用其它電離層探測數(shù)據(jù),如垂測、斜測的數(shù)據(jù)輔助地基GPS進行電離層層析成像,取得了較好的層析效果。

電離層掩星探測技術(shù)是一種伴隨GPS技術(shù)發(fā)展起來的可用于長期穩(wěn)定測量從地面至800 km高空電離層電子密度`的新技術(shù)。通過低軌衛(wèi)星接收GPS信號,隨著GPS和低軌衛(wèi)星的相對運動,信號將完成從電離層頂部直至地面的掃描過程,完成一次掩星的探測過程。掩星觀測具有高精度、高垂直分辨率、全天候觀測的特點,與地基GPS恰好能形成優(yōu)勢互補,可有效提升電離層層析成像的反演精度和垂直分辨率。地基GPS與掩星聯(lián)合開展電離層層析成像作為一個發(fā)展前景非常廣闊的新領(lǐng)域,國內(nèi)外已有很多學者開展了探索研究[14-17]。氣象/電離層及氣候衛(wèi)星探測系統(tǒng)(COSMIC)由6顆低軌小衛(wèi)星組成,每天能在全球范圍內(nèi)提供大約2 000 次掩星電離層資料,非常適于聯(lián)合地基GPS進行層析成像研究[18]。

本文試圖在定量的角度上對地基GPS與掩星聯(lián)合的電離層層析成像方法性能進行分析。選擇中國23個地基GPS站與COSMIC星座掩星事件作為仿真輸入,開展基于地基GPS與掩星聯(lián)合電離層層析成像方法研究。仿真結(jié)果表明:地基GPS與掩星聯(lián)合電離層層析成像方法獲取的電離層F2層峰值電子密度NmF2、F2層峰值高度hmF2及總電子含量(TEC)在平均反演誤差和均方根誤差方面均小于地基GPS方法,特別是在hmF2上表現(xiàn)最為明顯。

1 CIT原理與方法

CIT作為一個反演問題,它通過一系列GPS衛(wèi)星和接收機(包括地面GPS接收機和低軌衛(wèi)星掩星接收機)間無線電信號傳播路徑上的積分總電子含量TEC測量來重構(gòu)區(qū)域內(nèi)未知電離層電子密度分布[1]。接收機所獲得的TEC可以表示為沿信號傳播路徑上電子密度的積分,有

(1)

式中: TECi為總電子含量; Ne(r)為電子密度,隨時間和空間而變化; S為地面接收機至衛(wèi)星的視線路徑。

根據(jù)離散反演理論的方法,式(1)可以變換為下列線性方程組的問題

(2)

式中: G為地基GPS觀測的路徑數(shù)目；R為掩星觀測的路徑數(shù)目；N為網(wǎng)格數(shù)；向量d由地基GPS或掩星接收機的TEC觀測數(shù)據(jù)TECi組成; A為信號傳播路徑在離散化網(wǎng)格中的截距; X代表電子密度值;而e為測量與離散化后引入的誤差。將矩陣進行組合,式(2)可以簡化為以下形式

d=AX+e.

(3)

(4)

2 數(shù)值仿真及結(jié)果分析

2.1 仿真場景

選擇中國區(qū)域的北京、長春、德令哈、廣州、哈爾濱、海拉爾、昆明、拉薩、瀘州、瓊中、上海、綏陽、泰安、塔什、烏魯木齊、武漢、烏什、西安、下關(guān)、廈門、西寧、鹽池、鄭州等23個地基GPS觀測站,聯(lián)合COSMIC星座(包括6顆衛(wèi)星)的掩星測量進行CIT數(shù)值仿真。

設(shè)定緯度 13～55°N,經(jīng)度 70～140°E,高度100～800 km范圍內(nèi)的區(qū)域為CIT的區(qū)域,考慮到電離層的變化特征和計算機的性能的限制,將層析區(qū)域內(nèi)劃分為緯度間隔 1°,經(jīng)度間隔1°,高度間隔20 km。每次成像選取23個地基GPS臺站和COSMIC衛(wèi)星觀測的2 h內(nèi)的GPS射線進行,GPS衛(wèi)星和COSMIC衛(wèi)星位置由兩行軌道參數(shù)(TLE)計算獲得。

對于地基GPS,一般接收機的采樣間隔為30 s,由于衛(wèi)星軌道較高,GPS衛(wèi)星相對某一個固定接收機的運動角速度很小(周期約12 h),傳播路徑對應層析成像區(qū)域的網(wǎng)格而言基本沒有變化[12],為節(jié)省存儲空間和減小計算量,取0.5 h作為觀測間隔[12],觀測仰角的下限設(shè)置為15°;對于掩星而言,由于每次可觀測的掩星事件的時間僅約1～5 min,因此設(shè)置掩星的采樣間隔為15 s.

圖1(a)示出了2011年3月10日參與層析成像反演的地基GPS接收機位置分布及接收機-GPS衛(wèi)星間射線的穿刺點(450 km高度)在2 h內(nèi)的運動軌跡;圖1(b)示出了COSMIC掩星接收機與GPS衛(wèi)星間射線碰撞點的位置分布。從圖1可以看到地基觀測射線較為密集,基本覆蓋了層析成像區(qū)域,掩星在中國西南部觀測數(shù)目較多,東南部區(qū)域較少。

圖1 地基GPS穿刺點和掩星碰撞點分布(a)觀測站地理頒布及穿刺點軌跡(圓點為觀測點)；(b)COSMIC掩星觀測碰撞點分布

2.2 仿真流程

利用地基GPS和掩星進行電離層層析成像,首先應判斷GPS衛(wèi)星對于地基GPS或掩星接收機而言是否可見;接著對接收機與GPS衛(wèi)星間的信號傳播路徑是否“穿越”設(shè)定的反演區(qū)域內(nèi)進行判定,若滿足條件,在完成對所有地基和掩星接收機的遍歷后,利用接收機與GPS的位置的幾何關(guān)系構(gòu)建反演矩陣。

模擬層析成像區(qū)域內(nèi)的電離層電子密度的“真實”分布采用意大利國際中心理論物理(ICTP)與奧地利Graz大學提出的NeQuick經(jīng)驗電離層模型[19]進行計算,模型輸入的太陽輻射通量F10.7指數(shù)設(shè)定為125,通過式(1)的積分方程,可模擬出各地基GPS接收機和COSMIC掩星接收機觀測到的電離層TEC數(shù)據(jù)。為模擬測量噪聲的影響,所有TEC數(shù)據(jù)加入了約1～3 TECU的隨機誤差。

選取美國空軍研究實驗室研發(fā)的參數(shù)化電離層模型(PIM)作為電離層層析成像的背景電離層模型[20]用于MART算法的迭代初值,由于PIM模型與NeQuick模型在建模方法上各有不同,這樣做能減小背景模型對成像方法精度的影響。

最后，根據(jù)仿真得到的TEC數(shù)據(jù)進行層析成像,獲取電離層電子密度信息，具體的仿真流程如圖2所示。

圖2 電離層層析成像仿真流程圖

2.3 仿真結(jié)果與分析

為驗證算法穩(wěn)定性和精度,對算法反演的所有網(wǎng)格點的電子密度進行反演精度評估，分別對層析成像獲得的電子密度反演誤差進行評估,定義:

ΔNe=Ne_tomo-Ne_real，

(5)

(6)

(7)

式中:Ne_tomo表示層析反演(或背景模型)計算得到的電子密度;Ne_real表示真實的電子密度;N表示所有的反演算例個數(shù)。并對各算法的反演誤差進行統(tǒng)計分析。誤差分布柱狀圖如圖3所示。圖中PIM代表背景電離層模型,GPS代表地基GPS方法,GPS+RO代表聯(lián)合方法,由分析結(jié)果可以看出,背景模型、地基GPS、聯(lián)合方法反演得到的電離層電子密度誤差基本符合正態(tài)分布特征。其中背景模型誤差最大,地基GPS方法次之,聯(lián)合方法最小。

圖3 電子密度層析反演誤差比較

計算各成像方法的ΔNe絕對平均誤差AE(ΔNe):其中背景模型為1.2×1011el.m-3,地基GPS方法為0.8×1011el.m-3,聯(lián)合方法為0.5×1011el.m-3; 計算ΔNe均方根誤差RMS(ΔNe):背景模型為1.9×1011el.m-3,地基GPS方法為1.4×1011el.m-3,聯(lián)合方法為1.0×1011el.m-3,相對精度改善了近40%。聯(lián)合方法比僅利用地基GPS進行電離層層析成像在反演精度和算法穩(wěn)定性上均表現(xiàn)更優(yōu)。

F2層峰值電子密度NmF2以及F2層峰值高度hmF2是反映電離層電子密度剖面形狀的重要參量,特別是hmF2反演精度常選擇作為驗證CIT算法垂直分辨率的重要參考量[9-10]。為實現(xiàn)對地基GPS方法以及聯(lián)合方法的精確性和有效性的驗證,利用重構(gòu)的電子密度值計算出所有反演區(qū)域上空背景模型(PIM)、聯(lián)合方法計算出每個反演時刻對應的峰值電子密度NmF2和hmF2,定義:

ΔNmF2=NmF2_tomo-NmF2_real,

(8)

ΔhmF2=hmF2_tomo-hmF2_real,

(9)

其中: ΔNmF2表示F2層峰值電子密度絕對反演誤差;NmF2_tomo為層析反演的F2層峰值電子密度;NmF2_real為真實的F2層峰值電子密度; ΔhmF2表示F2層峰值高度絕對反演誤差;hmF2_tomo為層析反演的F2層峰值高度;hmF2_real為真實的F2層峰值高度,同樣參照式(6)和式(7)計算以上兩個量的平均誤差和均方根誤差,分析結(jié)果如表1和表2所示,表中“PIM”代表背景電離層模型,“GPS”代表地基GPS方法,“GPS+RO”代表聯(lián)合方法。

表1電離層層析成像反演ΔNmF2性能比較單位:×1011/(el.m-3)

層析時刻(UTC)平均誤差PIM GPSGPS 均方根誤差PIM GPSGPS+RO+RO 00∶00-02∶001.61.30.91.91.31.1 02∶00-04∶003.52.11.31.21.40.9 04∶00-06∶005.52.41.71.92.21.5 06∶00-08∶006.02.01.41.92.31.8 08∶00-10∶004.61.31.01.51.71.2 10∶00-12∶002.90.91.01.51.31.4 12∶00-14∶002.50.70.71.40.90.9 14∶00-16∶001.60.60.51.40.70.8 16∶00-18∶001.20.40.31.40.60.6 18∶00-20∶001.20.30.31.40.40.4 20∶00-22∶003.30.50.31.80.70.4 22∶00-24∶003.40.80.41.71.00.6 日均值3.11.10.81.61.21.0

表2 電離層層析成像得到的ΔhmF2性能比較 (單位∶km)

從表1可以看出,ΔNmF2誤差分布基本遵循白天較大,夜間較小的分布特征。從全天來看,日平均誤差背景模型為3.1×1011el.m-3;地基GPS方法為1.1×1011el.m-3,相比背景模型而言反演精度有了明顯提高;而聯(lián)合方法為0.8×1011el.m-3,精度比僅利用地基GPS有了近25%的改善。在反演穩(wěn)定性方面,日均方根誤差背景模型為1.6×1011el.m-3;地基GPS方法為1.2×1011el.m-3,而聯(lián)合方法為1.0×1011el.m-3,同樣是聯(lián)合方法更為出色。

對于F2層峰值電子密度反演精度方面,從表2可以看出,背景模型ΔhmF2誤差分布與時間變化并無呈現(xiàn)明顯規(guī)律,這與PIM模型與NeQuick模型的時間錯位有一定關(guān)系。從全天來看,日平均誤差背景模型為27.7 km;地基GPS方法為24.7 km,相比較背景模型而言反演精度改善有限;而聯(lián)合方法為15.2 km,精度比僅利用地基GPS有了近50%的改善。同樣的,聯(lián)合方法在反演的穩(wěn)定性方面也更好,誤差均方根下降到不到20 km.必須指出的是,不管是背景模型PIM或NeQuick模型,它們輸出的hmF2一般在250～400 km范圍內(nèi)變化,并不會像NmF2一樣,不同模型在不同時間上會出現(xiàn)一個甚至數(shù)個數(shù)量級的差別,因此地基GPS與掩星聯(lián)合方法在hmF2參量的改進程度上并沒有像NmF2那么大,這是由于采用經(jīng)驗電離層模型仿真而引起的,與方法本身無關(guān)。

同樣,利用電離層層析成像方法獲取了區(qū)域電離層電子密度分布后,根據(jù)電子密度隨路徑積分得到電離層TEC的原理,可以重構(gòu)出中國區(qū)域的垂直TEC分布,以此可評估層析成像方法的輸出與輸入的TEC數(shù)據(jù)之間的“符合”程度。定義垂直TEC絕對重構(gòu)誤差為ΔTEC=|TEC_tomo-TEC_real|,其中TEC_tomo為層析成像重構(gòu)出的電離層TEC值,TEC_real為真實的電子密度重構(gòu)出的電離層TEC值。以UT00∶00-02∶00為例,分析結(jié)果如圖4所示。

從圖4可以看出,背景模型的電離層TEC誤差分布呈南高北低的特征,基本符合我國南方處于電離層磁低緯而北部屬于磁中緯的特點。利用地基GPS進行層析成像后,電離層TEC誤差有了非常明顯的下降,特別是在我國北緯25°以南的區(qū)域,誤差從7～10 TECU下降為不足～2 TECU,東北部區(qū)域的電離層TEC誤差也有明顯下降。而聯(lián)合方法表現(xiàn)則更為出色,反演的中國區(qū)域上空的電離層TEC誤差有了更加明顯的下降;從圖2(b)可以看出,由于較多掩星射線路徑穿越了西南區(qū)域,聯(lián)合方法在西南區(qū)域的TEC重構(gòu)誤差相比地基GPS改善更為明顯。

對所有12個時段層析成像的電離層TEC誤差進行統(tǒng)計分析的結(jié)果如表3所示,從表中可以看出,ΔTEC變化特征與NmF2非常類似,背景模型TEC誤差在白天偏大,夜間偏小。在日平均誤差方面,背景模型為6.7 TECU,地基GPS方法為1.3 TECU,聯(lián)合方法為1.0 TECU;均方根誤差,背景模型日平均為3.5 TECU,地基GPS方法為1.8 TECU,而聯(lián)合方法為1.6 TECU.從分析結(jié)果來看,聯(lián)合方法在TEC重構(gòu)精度方面同樣比地基GPS方法更優(yōu)。

圖4 電離層TEC反演誤差比較(a)背景模型ΔTEC分布(b)地基GPS方法ΔTEC分布(c)地基GPS與掩星聯(lián)合方法ΔTEC分布

表3電離層ΔTEC誤差比較(單位:TECU)

層析時刻(UTC)平均誤差PIM GPSGPS 均方根誤差PIM GPSGPS+RO+RO 00∶00-02∶004.81.51.04.92.31.8 02∶00-04∶005.40.90.84.41.51.4 04∶00-06∶009.81.81.63.02.12.1 06∶00-08∶0011.42.11.82.02.52.4 08∶00-10∶009.21.71.53.22.32.2 10∶00-12∶006.81.31.33.31.91.9 12∶00-14∶006.31.01.03.41.81.7 14∶00-16∶004.00.90.93.41.71.7 16∶00-18∶002.80.70.73.11.31.3 18∶00-20∶002.70.60.43.30.90.7 20∶00-22∶007.81.00.53.81.60.8 22∶00-24∶009.61.40.93.61.91.6 日均值6.71.31.03.51.81.6

3 結(jié)束語

掩星技術(shù)是未來天基電離層探測領(lǐng)域最具發(fā)展前景的領(lǐng)域之一,由于其能夠?qū)﹄婋x層進行垂直方向的掃描,獲得非常好的電離層探測垂直分辨率,與地基GPS電離層探測聯(lián)合,能有效的提升電離層層析成像的電子密度反演精度。本文采用實際的GPS和COSMIC衛(wèi)星星歷,選擇NeQuick模型模擬電離層電子密度的“真實”分布,對地基GPS與掩星聯(lián)合的電離層層析成像方法進行了研究。利用中國區(qū)域23個地基GPS站及COSMIC星座進行的仿真表明,地基GPS與聯(lián)合掩星電離層層析成像,比僅利用地基GPS方法,在反演的電子密度和電子總含量方面,均表現(xiàn)出明顯的優(yōu)勢。特別是在hmF2的反演精度方面,相比地基GPS,精度提高了50%左右,有效的提升了電離層層析成像的垂直分辨率水平。

隨著全球其他國家的衛(wèi)星導航系統(tǒng)的發(fā)展,利用地基GPS和掩星進行電離層觀測將會有更好的空間覆蓋性,這有利用電離層層析成像的精度和分辨率的進一步的提升。

[1]AUSTEN J R, FRANKE S J, LIU C H,etal. Application of computerized tomography technique to ionospheric research[C]//URSI and COSPAR International Beacon Satellite Symposium on Radio Beacon Contribution to the Study of Ionization and Dynamics of the Ionosphere and to Corrections to Geodesy and Technical Workshop, 1986.

[2]RAYMUND T D, AUSTEN J R, FRANKE S J,etal. Application of computerized tomography to the investigation of ionospheric structures [J]. Radio Science, 1990, 25(5): 771-789.

[3]PRYSE S,KERSLEY E.A preliminary experimental test of ionospheric tomography [J]. Journal of Atomospheric and Solar-Terrestrial Physics,1992(54): 1007-1012.

[4]KUNITSYN V E, ANDREEVA E S, RAZINKOV O G. Possibilities of the near-space environment radiotomography [J].Radio Science, 1997, 32 (5)：1953-1963.

[5]HANSEN A J, WALTER T, ENGE P. Ionospheric correction using tomography[C]//In Proceeding of Institute of Navigation ION GPS-97, Kasas City, Missouri, USA, September, 16-19, 1997: 249-260.

[6]RUFFINI G, FLORES A, RIUS A. GPS tomography of the ionospheric electron content with a correlation functional[J].IEEE Transactions on Geoscience and Remote Sensing,1998,36(1):143-153.

[7]HOWE B M, RUNCIMAN K, SECAN J A. Tomography of Ionospheric: Four-dimensional simulations[J]. Radio Science,1998,33(1):109-128.

[8]YIZEGAE E, MOLDWIN M B, DYSON P L,etal, Using tomography of GPS TEC to routinely determine ionospheric average electron density profiles[J]. Journal of Atomospheric and Solar-Terrestrial Physics, 2007, 69(3):314-321.

[9]徐繼生, 馬淑英, 吳雄斌，等. 一次中強磁暴期間低緯電離層響應的CT成像[J].地球物理學報, 2000, 43(2): 145-151.

[10]鄒玉華. GPS地面臺網(wǎng)和掩星觀測結(jié)合的時變?nèi)S電離層層析[D].武漢:武漢大學, 2004.

[11]聞德保.基于GPS的電離層層析成像算法及其應用研究[D]. 中科院測量與地球物理研究所, 2007.

[12]樂新安,萬衛(wèi)星,劉立波,等.基于Gauss-Markov卡爾曼濾波的電離層數(shù)值同化現(xiàn)報預報系統(tǒng)的構(gòu)建——以中國及周邊地區(qū)為例的觀測系統(tǒng)模擬試驗[J]. 地球物理學報,2010,53(4):787-795,

[13]ZHAO H S, XU Z W, WU J,etal. Ionospheric tomography by combining vertical and oblique ionograms with TEC retrieved from tri-band beacon [J]. Journal of Geophysical Research, 2010: 115.

[14]LI Hui, YUAN Yunbin, LI Zhishen,etal. Ionospheric electron concentration imaging using combination of LEO satellite data with ground-gased GPS observations over China[J]. IEEE Transactions on Geo-science and Remote Sensing, 2012,50(5):1728-1734.

[15]MITCHELL C N, CANNON P S, SPENCER P S. Multi-instrument data analysis system imaging of the ionosphere [R]. The United States Air Force European Office of Aerospace Research and Development, 2002: 1-30.

[16]BUST G S, GARNER T W, GAUSSIRAN II T L. Ionospheric data assimilation three-dimensional (IDA3D): A global, multisensor, electron density specification algorithm[J].Journal of Geophysical Research,2004,109,A11312.doi:10.1029/2003 JAo10234.

[17]YUE Xinan, WILLIAM S S, YU Chenglin,etal. Data assimilation retrieval of electron density profiles from radio occultation measurements [J]. IEEE Transactions on Geoscience and Remote Sensing,2011,116,A03317.doi:10.1029/2010J A015980,

[18]徐賢勝, 洪振杰, 郭 鵬,等, 掩星電離層資料反演以及結(jié)果驗證[J]. 物理學報, 2010: 59(3), 2163-2168.

[19]NAVA B, RADICELLA S M, AZPILICUETA F,etal. Data ingestion into NeQuick 2.[J]. Radio Science, 2011,46(6):1-8.

[20]DANIELL R E, BROWN L D,ANDERSON D N,etal. Parameterized ionospheric model: A global ionospheric parameterization based on first principles models [J]. Radio Science, 1995, 30(5): 1499-1510.