基于BDS/GPS 數(shù)據(jù)的中國(guó)及周邊地區(qū)電離層暴層析成像研究

尹萍 寧澤浩 閆曉鵬

（中國(guó)民航大學(xué)電子信息與自動(dòng)化學(xué)院,天津 300300）

引 言

電離層是地球表面60 km 以上的大氣由于受太陽(yáng)輻射的照射處于部分電離或完全電離狀態(tài)的等離子體區(qū)域，主要由電子、正離子和中性分子及原子構(gòu)成.而電離層的不規(guī)則運(yùn)動(dòng)會(huì)對(duì)人造衛(wèi)星、飛機(jī)的電子設(shè)備造成信號(hào)中斷和信號(hào)誤碼等影響.隨著GNSS 理論的快速發(fā)展，電離層層析成像(computerized ionospheric tomography,CIT)技術(shù)展現(xiàn)出探測(cè)電離層電子密度三維分布和變化方面的優(yōu)勢(shì)，其主要依靠區(qū)域地面接收機(jī)接收到的衛(wèi)星數(shù)據(jù)，計(jì)算得到傳播路徑上的總電子含量(total electron content,TEC)，并利用相關(guān)的層析算法，解算出區(qū)域內(nèi)的電子密度.

1998 年，Austen 等人在國(guó)際上提出了CIT 的思想，同時(shí)利用海軍導(dǎo)航衛(wèi)星系統(tǒng)(navy navigation satellite system,NNSS)數(shù)據(jù)重建了臺(tái)站平面與衛(wèi)星的二維電子密度剖面[1]；然而，利用NNSS 數(shù)據(jù)進(jìn)行的層析成像只能獲得電子密度隨緯度-高度變化的二維分布圖像，而在實(shí)際應(yīng)用與研究中掌握電離層三維時(shí)空分布至關(guān)重要[2].由于GPS 具有軌道高、服務(wù)范圍廣等特點(diǎn)，基于GPS 的CIT 技術(shù)得到了快速發(fā)展，從數(shù)據(jù)源方面有效改善了利用NNSS 數(shù)據(jù)僅能獲得電子密度二維分布的局限性.2005 年徐繼生等人利用GPS 地面臺(tái)網(wǎng)和掩星觀測(cè)數(shù)據(jù)，結(jié)合CIT技術(shù)重建了電離層大尺度結(jié)構(gòu)的時(shí)空演變[3].2007 年Jin 等人利用韓國(guó)GPS 網(wǎng)數(shù)據(jù)，重建了2003 年韓國(guó)上空電離層電子密度分布，并通過(guò)IRI-2001 模型給出的F2層峰值密度NmF2、TEC 和另一個(gè)獨(dú)立的電離層探測(cè)儀數(shù)據(jù)對(duì)比進(jìn)一步驗(yàn)證了該方法的有效性[4].2011年楊劍等人利用GPS 觀測(cè)數(shù)據(jù)對(duì)高度60～960 km的電離層電子密度分布進(jìn)行了三維層析重建，研究地震前電離層異常情況[5].2012 年陳必焰將神經(jīng)網(wǎng)絡(luò)的方法應(yīng)用到了基于GPS 的CIT 技術(shù)并探測(cè)日本地震前的電離層異常，分析了“3·11”地震前后區(qū)域上空電離層電子密度時(shí)空變化[6].

太陽(yáng)爆發(fā)所輻射的能量使地球磁場(chǎng)產(chǎn)生劇烈擾動(dòng)，稱為“磁暴”.磁暴的產(chǎn)生會(huì)使電離層的電子產(chǎn)生不規(guī)則的變化，導(dǎo)致電離層TEC 和電子密度的異常，包括電離層閃爍、電離層暴以及電離層行擾(travelling ionospheric disturbance,TID)等等，利用CIT 技術(shù)研究磁暴期間電子密度變化成為國(guó)內(nèi)外學(xué)者研究的重點(diǎn).2000 年徐繼生等人利用CIT 技術(shù)重建電離層電子密度二維剖面，研究了中強(qiáng)度磁暴對(duì)低緯度區(qū)域電離層的影響[7].2005 年Yizengaw 等人利用CIT 技術(shù)反演了2001-03-31 磁暴期間150°E 經(jīng)度面內(nèi)南半球電離層電子密度的變化[8].2014 年姚宜斌等人利用CIT 技術(shù)分別反演并分析了2003—2006年間3 次不同類型強(qiáng)磁暴發(fā)生期間歐洲區(qū)域上空電離層電子密度的三維時(shí)空分布以及電離層擾動(dòng)結(jié)構(gòu)的演化過(guò)程[9].2017 年趙海山等人利用CIT 技術(shù)對(duì)歐洲區(qū)域磁暴前和磁暴期間的峰值電子密度變化進(jìn)行了分析[10].

然而，由于GPS 主要采用單一的中軌軌道(medium earth orbit,MEO)衛(wèi)星，該數(shù)據(jù)容易受到衛(wèi)星運(yùn)動(dòng)的影響[11]，而B(niǎo)DS 中的地球同步軌道(geosynchronous orbit,GEO)衛(wèi)星由于自身的靜地特性，可提供更好的數(shù)據(jù)集研究經(jīng)度(時(shí)間)上的變化，使觀測(cè)更可靠[12].同時(shí)，將GNSS 多星座的數(shù)據(jù)應(yīng)用于CIT，相較于單星座的GPS，有更多的衛(wèi)星射線穿過(guò)目標(biāo)區(qū)域電離層，進(jìn)一步提高了CIT 反演的精度[13].

本文在英國(guó)巴斯大學(xué)多儀器數(shù)據(jù)分析系統(tǒng)(multi instrument data analysis system,MIDAS)CIT 算法[14]基礎(chǔ)上，利用IGS 官網(wǎng)提供的GPS 和BDS 觀測(cè)數(shù)據(jù)對(duì)2017-09-07—09 中國(guó)及其周邊地區(qū)出現(xiàn)的電離層暴進(jìn)行層析成像.借助中國(guó)地區(qū)兩個(gè)獨(dú)立的測(cè)高儀站點(diǎn)數(shù)據(jù)，評(píng)估單、雙星座層析算法獲取的F2層峰值密度NmF2、F2層峰值高度hmF2和TEC 等電離層重要參數(shù)，以討論BDS 數(shù)據(jù)在層析方法上的可用性和可靠性.在此基礎(chǔ)上，對(duì)2017-09-08 兩次強(qiáng)擾動(dòng)時(shí)段中國(guó)及周邊低緯地區(qū)的層析結(jié)果進(jìn)一步研究，觀測(cè)出赤道等離子體泡(equatorial plasma bubble，EPB)現(xiàn)象，并分析暴時(shí)電離層時(shí)空分布特性.

1 電離層暴特征分析

由于兩次日冕物質(zhì)拋射(coronal mass ejections,CMEs)于2017-09-07T23:30UT 和09-08T12:00UT 先后經(jīng)過(guò)地球，產(chǎn)生了一場(chǎng)強(qiáng)地磁暴[15].通過(guò)分析從NASA 數(shù)據(jù)庫(kù)獲取的相關(guān)物理參數(shù)和地磁數(shù)據(jù)，選取發(fā)生在2017-09-08 的電離層暴為主要研究對(duì)象.為了說(shuō)明該電離層暴的特點(diǎn)，圖1 給出了2017-09-07—09 的太陽(yáng)風(fēng)速變化vsw，行星際磁場(chǎng)Bz分量IMFBz，太陽(yáng)風(fēng)動(dòng)壓Psw，極光電子急流AE 指數(shù)，地磁指數(shù)Dst 和Kp.

圖1 2017-09-07—09 地球物理數(shù)據(jù)變化Fig.1 Variations in geophysical data on Sep.7-9,2017

伴隨著第一次CMEs，vsw在2017-09-07T23:30UT由473.4 km/s 增長(zhǎng)至599.7 km/s；IMFBz在20:00UT向南轉(zhuǎn)并到達(dá)-10 nT 且保持2 h，然后在23:30UT 再次南轉(zhuǎn)到達(dá)最小值-31.2 nT；Psw由1.31 nPa 增長(zhǎng)至4.35 nPa；Dst 指數(shù)處于-144 nT，Kp 值為8，說(shuō)明地磁活動(dòng)已到達(dá)首個(gè)擾動(dòng)階段.伴隨著第二次CMEs，vsw從2017-09-08T11:30UT 開(kāi)始持 續(xù)下降；IMFBz在11:55UT 瞬間南轉(zhuǎn)降至-17.4 nT；Psw從8.5 nPa 下降至2.9 nPa；同時(shí)Kp 出現(xiàn)最大值(超過(guò)8)，Dst 值仍處于-100 nT 以下，表明該時(shí)刻已到達(dá)第二個(gè)擾動(dòng)階段.另外，AE 在2017-09-08T00:00UT 和 13:00—14:00UT 均超過(guò)2 000 nT，表明由磁暴產(chǎn)生的兩次擾動(dòng)均伴隨著強(qiáng)烈的極光活動(dòng).

從圖1 2017-09-07—09 地球物理數(shù)據(jù)變化可以看出，Dst 和Kp 值的變化都說(shuō)明了地磁擾動(dòng)最強(qiáng)烈的時(shí)段主要發(fā)生在2017-09-08.為了體現(xiàn)暴時(shí)電離層的異常特性，重點(diǎn)研究處于較平靜期的09-07 和處于磁暴期的09-08 的電離層時(shí)空分布.

2 觀測(cè)數(shù)據(jù)與層析成像技術(shù)

2.1 觀測(cè)數(shù)據(jù)

作為本文電離層層析算法的輸入數(shù)據(jù)，圖2 給出了中國(guó)及周邊地區(qū)所用站點(diǎn)的分布以及該地區(qū)上空電離層穿刺點(diǎn)(ionospheric pierce point,IPP)的覆蓋情況.層析所用到的中國(guó)及周邊地區(qū)GPS 站點(diǎn)分布和BDS 站點(diǎn)分布如圖2(a)所示，黑色代表可采集GPS 觀測(cè)數(shù)據(jù)的地面站，紅色代表可同時(shí)采集GPS和BDS 觀測(cè)數(shù)據(jù)的地面站，全部合計(jì)26 個(gè)站點(diǎn)(部分站點(diǎn)位置重合)，其中8 個(gè)站可采集BDS 數(shù)據(jù)，已覆蓋中國(guó)、日本、韓國(guó)、泰國(guó)、越南等國(guó)家.但相比于GPS，該地區(qū)可采集BDS 數(shù)據(jù)的地面站數(shù)量較為有限，且在30～55°N 的站點(diǎn)分布稀疏，從數(shù)量、分布均勻程度上看，GPS 地面站更具有優(yōu)勢(shì).圖2(b)給出了2017-09-08T01:00—02:00UT 中國(guó)及周邊地區(qū)在假定電離層薄殼高度為350 km 時(shí)，沿地面站與GPS/BDS 衛(wèi)星路徑上IPP 的軌跡，即輸入數(shù)據(jù)的覆蓋范圍，共使用BDS 衛(wèi)星14 顆，C01～C05 為GEO衛(wèi)星，C06～C10、C13 為傾斜地球同步軌道(inclined geosynchronous orbit，IGSO)衛(wèi)星，C11、C12、C14 為MEO 衛(wèi)星[16]，其中“5GEO+6IGSO”用來(lái)保證服務(wù)區(qū)域內(nèi)的性能要求，3 顆MEO 衛(wèi)星用來(lái)實(shí)現(xiàn)擴(kuò)大地區(qū)覆蓋(01:00—02:00UT 期間MEO 衛(wèi)星未經(jīng)過(guò)該地區(qū)上空).從圖2(b)可以看出GPS IPP 覆蓋區(qū)域較大，其中東部更為集中.BDS IPP 覆蓋較集中在5～30°N和100～130°E.

同時(shí)，本文借助中國(guó)地區(qū)的兩個(gè)測(cè)高儀數(shù)據(jù)來(lái)驗(yàn)證層析成像結(jié)果的精度，一個(gè)是位于武漢地區(qū)的WU430(位置：30.5°N,114.4°E)；另一個(gè)是位于三亞地區(qū)的SA418(位置：18.34°N,109.42°E)，其位置如圖2圓圈所示.測(cè)高儀數(shù)據(jù)從Digital Ionogram Database官網(wǎng)下載，地址：http://ulcar.uml.edu/ DIDBase/[17].

圖2 中國(guó)及周邊地區(qū)GPS/BDS 站點(diǎn)分布及2017-09-08T01:00—02:00UT 在350 km 上空處的IPP 覆蓋Fig.2 Distribution of GPS/BDS ground stations and IPP tracks at 350 km over China and adjacent areas during 01:00—02:00UT on Sep.8,2017

2.2 數(shù)據(jù)預(yù)處理

首先將電離層看成高度為h的薄層模型[18]，式(1)為利用時(shí)間間隔為30 s 的GPS/BDS 觀測(cè)數(shù)據(jù)計(jì)算電波視線方向上的斜TEC[13](slant TEC,STEC)：

式中：f1、f2為信號(hào)頻率；λ1、λ2為波長(zhǎng)；L1、L2為載波相位.GPS選用f1= 1 575.42 MHz，f2=1 227.60 MHz，BDS 選用f1= 1 561.098 MHz，f2=1 207.14 MHz.為得到精確的電離層STEC，必須消除由GPS、BDS 衛(wèi)星和接收機(jī)帶來(lái)的硬件延遲[19].針對(duì)不同系統(tǒng)，在式(1)的基礎(chǔ)上分別采用以下兩種方法來(lái)計(jì)算硬件延遲.

針對(duì)GPS，記GPS 衛(wèi)星和接收機(jī)硬件延遲之和為QGPS，式(2)為QGPS與STEC 和垂直TEC(vertical TEC,VTEC)的關(guān)系[19]：

F為映射函數(shù)，

式中：R為地球半徑；h為電離層質(zhì)心高度；E為衛(wèi)星仰角.

假設(shè) VTEC(t)是關(guān)于時(shí)間和電離層IPP(λ ,φ)的函數(shù)：

式中：a(t)、b(t)、c(t)為 系數(shù)；λ0、φ0分別為觀測(cè)站經(jīng)度、緯度.將式(5)代入式(2)，采用最小二乘法估算硬件延遲QGPS.

針對(duì)BDS，記BDS 衛(wèi)星和接收機(jī)硬件延遲之和為QBDS，利用式(6)[20]和IGS 官網(wǎng)公布的硬件延遲偏差改正/差分碼偏差(differential code biases,DCB)文件[21]對(duì)QBDS進(jìn)行解算[22]：

式中：c為光速；Δt1為DCB文件中總?cè)簳r(shí)延(total group delay,TGD)和頻間偏置(inter-frequency bias,IFB)的和.

2.3 多星座數(shù)據(jù)層析成像技術(shù)

本文選用的層析方法是由英國(guó)巴斯大學(xué)Invert小組開(kāi)發(fā)的MIDAS 層析技術(shù)[14].傳統(tǒng)的單星座層析方法是將每條GPS 衛(wèi)星信號(hào)射線上的STEC 進(jìn)一步采樣(采樣時(shí)間間隔10 min)，再作為層析反演的輸入，數(shù)學(xué)公式為

式中：zGPS代表GPS 系統(tǒng)下的 STEC值；H代表衛(wèi)星射線穿過(guò)格網(wǎng)時(shí)截距組成的系數(shù)矩陣；x代表所求電子密度.

雙星座層析方法是先對(duì)每條GPS、BDS 衛(wèi)星信號(hào)射線上的STEC 進(jìn)行硬件延遲修正，再進(jìn)一步采樣(采樣時(shí)間間隔10 min)，然后將整體作為層析反演的輸入，數(shù)學(xué)公式為

又可寫為

式中：zGPS+BDS代表雙系統(tǒng)下的 STEC值 ；PGPS、PBDS分別為GPS、BDS 所求 STEC時(shí) 的權(quán)值，本文設(shè)為1；zGPS、zBDS分別為GPS、BDS 所求STEC 值.

由于待反演區(qū)域地面GNSS 接收機(jī)和衛(wèi)星數(shù)量有限且分布不均勻，很多像素內(nèi)沒(méi)有GPS、BDS 衛(wèi)星射線穿過(guò)導(dǎo)致H矩陣嚴(yán)重秩虧.為克服這種問(wèn)題，在MIDAS 中使用映射矩陣，式(8)可轉(zhuǎn)化為

式中，y代表變換基的解.MIDAS 采用了多種函數(shù)來(lái)構(gòu)造映射矩陣M[23]，本文采用Chapman 模型生成二階經(jīng)驗(yàn)正交函數(shù)(empirical orthogonal functions,EOF)來(lái)構(gòu)建垂向的映射矩陣，獲得電子密度的解為

選取電離層重構(gòu)區(qū)域?yàn)?5～55°N和70～140°E，經(jīng)緯度分辨率為4°×4°，高度分辨率為40 km，時(shí)間分辨率為15 min.考慮到多徑效應(yīng)和本地噪聲的影響以及可用的觀測(cè)數(shù)據(jù)量，在電離層TEC 的計(jì)算和分析中選取衛(wèi)星仰角大于5°的觀測(cè)數(shù)據(jù).

此外，將雙星座層析結(jié)果和單星座層析結(jié)果進(jìn)行比對(duì)來(lái)說(shuō)明融合BDS 數(shù)據(jù)的電離層層析技術(shù)在電離層監(jiān)測(cè)上的改善情況.為方便描述，在后面敘述中將單星座簡(jiǎn)稱為‘GPS’，雙星座簡(jiǎn)稱為‘GPS+BDS’ .

3 結(jié)果與討論

3.1 NmF2 層析結(jié)果分析

如圖3 所示，選取電離層相對(duì)平靜的2017-09-07、2017-09-09 和暴時(shí)的2017-09-08，用中國(guó)區(qū)域中低緯度SA418 站點(diǎn)、WU430 站點(diǎn)處兩個(gè)測(cè)高儀的NmF2分別與GPS 和GPS+BDS 層析結(jié)果進(jìn)行對(duì)比，討論利用GPS+BDS 層析技術(shù)在電離層平靜期和暴時(shí)反演電子密度的精度和可靠性.測(cè)高儀數(shù)據(jù)的時(shí)間間隔為15 min，用紅色‘○’表示.

圖3 2017-09-07—09 兩站點(diǎn)處GPS、GPS+BDS 層析結(jié)果和測(cè)高儀的NmF2 對(duì)比Fig.3 Comparison of NmF2 between GPS,GPS+BDS and ionosonde measurements at SA418 station and WU430 station on Sep.7-9,2017

較平靜期圖3(a)SA418 站點(diǎn)的層析結(jié)果，GPS 層析結(jié)果在4:00—7:00UT 期間波動(dòng)劇烈，且與SA418測(cè)高儀數(shù)據(jù)變化趨勢(shì)相差較大，而GPS+BDS 層析結(jié)果在此期間在變化趨勢(shì)上更接近測(cè)高儀數(shù)據(jù).圖3(b)中，兩種層析結(jié)果顯示都與WU430 測(cè)高儀觀測(cè)數(shù)據(jù)相差不大，且在2:00—11:00UT、22:00—0:00UT 期間幾乎重合(由于測(cè)高儀原因，目前無(wú)法對(duì)比其他時(shí)間段實(shí)測(cè)數(shù)據(jù)).觀察暴時(shí)的圖3(c)和(d)，圖3(c)中兩種層析結(jié)果的NmF2變化趨勢(shì)與SA418 測(cè)高儀觀測(cè)數(shù)據(jù)大體一致，均能反映當(dāng)日NmF2的變化，7:00—11:00UT 期間，GPS+BDS 層析結(jié)果的變化趨勢(shì)更接近SA418 測(cè)高儀數(shù)據(jù)；圖3(d)中，兩種層析結(jié)果在趨勢(shì)上仍保持一致.觀察恢復(fù)期的圖3(e)和(f)，圖3(e)與圖3(a)相似，在4:00—7:00UT 期間GPS層析結(jié)果波動(dòng)較大；而圖3(f)中部分時(shí)間點(diǎn)如09-09T8:00UT 的GPS+BDS 層析結(jié)果更接近WU430 測(cè)高儀的變化趨勢(shì).

通過(guò)對(duì)SA418 和WU430 測(cè)高儀三天數(shù)據(jù)進(jìn)行比較發(fā)現(xiàn)，WU430 測(cè)高儀NmF2整體偏小.由圖3(c)和(d)中的兩個(gè)測(cè)高儀得知，在09-08T01:00—02:00UT期間，由于受第一個(gè)擾動(dòng)影響，SA418 測(cè)高儀NmF2與同時(shí)間段的09-07 和09-09 相比，存在明顯的突變(升高和降低)；同時(shí)，在09-08T02:00—07:00UT 期間，WU430 測(cè)高儀NmF2較同時(shí)間段的09-07 和09-09 也表現(xiàn)出同樣的升高；在09-08T13:00—20:00UT 期間，受第二個(gè)擾動(dòng)影響，SA418 測(cè)高儀NmF2觀察到明顯的劇烈擾動(dòng).同時(shí)，從兩種層析結(jié)果上均可以看到，WU430 測(cè)高儀NmF2低于09-07 的同時(shí)間段，這可能與后文分析得到的EPB 有關(guān).

表1 為2017-09-07—09 SA418 和WU430兩個(gè)站點(diǎn)處GPS、GPS+BDS 層析NmF2和測(cè)高儀數(shù)據(jù)的統(tǒng)計(jì)量均值和均方根誤差(root mean square error,RMSE)對(duì)比.整體來(lái)看，BDS 數(shù)據(jù)對(duì)三天中個(gè)別時(shí)段的GPS+BDS 層析結(jié)果略有優(yōu)化，這可能是因?yàn)镚PS 系統(tǒng)中均為MEO 衛(wèi)星，而B(niǎo)DS 中的GEO 和IGSO 衛(wèi)星在地球上空的運(yùn)動(dòng)范圍相對(duì)固定，其中五顆GEO 衛(wèi)星在緯圈方向間隔1 000 km 左右且衛(wèi)星的IPP 幾乎保持靜止，為IPP 提供連續(xù)觀測(cè)數(shù)據(jù)的同時(shí)不受電離層時(shí)空變化等其他多種效應(yīng)的影響[24].

表1 2017-09-07—09 兩站點(diǎn)處GPS、GPS+BDS 層析結(jié)果和測(cè)高儀的NmF2 比較Tab.1 Comparison of NmF2 errors between GPS,GPS+BDS tomographic results and ionosonde data at SA418 and WU430 on Sep.7-9,2017 1011 el·m-3

此外，由圖2 還觀察到，SA418 所在地區(qū)缺少可接收BDS 數(shù)據(jù)的站點(diǎn)，同時(shí)BDS IPP 覆蓋效果略差于WU430，導(dǎo)致SA418 處誤差明顯高于WU430：均值最高相差1.5×1011el·m-3，最低相差0.7×1011el·m-3；RMSE值最高相差2.3×1011el·m-3，最低相差1.0×1011el·m-3.說(shuō)明站點(diǎn)分布的均勻性和IPP 的覆蓋程度影響了層析算法的精度.在具有相同站點(diǎn)數(shù)的情況下，GPS+BDS 層析結(jié)果更精確，進(jìn)一步驗(yàn)證了在電離層相對(duì)平靜期和暴時(shí)兩種狀態(tài)下，融入BDS 數(shù)據(jù)的層析反演結(jié)果更精確.

3.2 hmF2 層析結(jié)果分析

為討論GPS+BDS 層析技術(shù)對(duì)hmF2的影響，同時(shí)考慮到測(cè)高儀數(shù)據(jù)完好性問(wèn)題，圖4 給出了2017-09-08 SA418 和WU430 兩個(gè)站點(diǎn)處GPS、GPS+BDS 為數(shù)據(jù)源層析反演得到的hmF2與測(cè)高儀的hmF2對(duì)比結(jié)果.SA418 測(cè)高儀hmF2時(shí)間間隔為7.5 min，WU430測(cè)高儀hmF2時(shí)間間隔為15 min，站點(diǎn)觀測(cè)數(shù)據(jù)在某些時(shí)刻缺失，用零值代替.

圖4(a)中 由于IMFBz在2017-09-08T11:55UT期間進(jìn)行了大幅度南轉(zhuǎn)，可以觀察到SA418 測(cè)高儀hmF2顯著增加：在11:00—12:00UT 期間急劇上升了近100 km，并在12:00 達(dá)到該時(shí)段的最大值，而在之后的12:00—18:00UT 存在明顯高度上的起伏，幅度范圍在180 km左右.而GPS、GPS+BDS 層析結(jié)果在hmF2上的高度分布和變化趨勢(shì)并沒(méi)有太大差異.圖4(b)與(a)相反，GPS+BDS 層析結(jié)果中hmF2在垂直結(jié)構(gòu)上表現(xiàn)出明顯差異，尤其在6:00—12:00UT，在垂直高度上與GPS 層析結(jié)果相差近100 km；在9:00—12:00UT 有逐步升高的趨勢(shì)，與WU430 測(cè)高儀的觀測(cè)結(jié)果相符.

圖4 2017-09-08 兩站點(diǎn)處GPS、GPS+BDS 層析結(jié)果和測(cè)高儀的hmF2 對(duì)比Fig.4 Comparisons of hmF2 between GPS,GPS+BDS and ionosonde measurements at SA418 station and WU430 station on Sep.8,2017

表2 所示為2017-09-08 兩站點(diǎn)處GPS,GPS+BDS層析結(jié)果和測(cè)高儀的hmF2誤差比較.可以看出，GPS 與GPS+BDS 層析結(jié)果和SA418 測(cè)高儀之間hmF2的誤差幾乎持平，而WU430 站點(diǎn)GPS+BDS 層析結(jié)果的誤差更小一些.結(jié)合圖4(b)和表2 分析推斷：GPS+BDS 層析結(jié)果可以改進(jìn)暴時(shí)的hmF2精度，而B(niǎo)DS 數(shù)據(jù)更密集的地區(qū)(如中國(guó)武漢地區(qū))得到的層析結(jié)果的hmF2精度更高.

表2 2017-09-08 兩站點(diǎn)處GPS、GPS+BDS 層析結(jié)果和測(cè)高儀的hmF2 比較Tab.2 Comparison of hmF2 errors between GPS,GPS+BDS tomographic results and ionosonde data at SA418 and WU430 on Sep.8,2017 km

圖5 為2017-09-07 —08T02:00UT 和15:00UT 114°E經(jīng)度線上電子密度隨緯度與高度分布情況的GPS+BDS 層析結(jié)果.紅色虛線交點(diǎn)處為WU430 所在地(30°N，114°E)上空400 km 處位置.其中，2017-09-08T02:00UT 處于第一個(gè)電離層擾動(dòng)時(shí)間段，該時(shí)間段內(nèi)114°E 經(jīng)度線位置上的GPS IPP 和BDS IPP 覆蓋最為集中；2017-09-08T 15:00UT 處于第二個(gè)擾動(dòng)時(shí)間段.

圖5 2017-09-07—08 不同時(shí)間時(shí)刻114°E 經(jīng)度鏈上電子密度隨緯度與高度變化的GPS+BDS 層析結(jié)果Fig.5 Variations of electron density with latitude and height along 114°E longitude at 02:00UT and 15:00UT on Sep.7-8,2017

圖5(a)為2017-09-07T02:00UT 電子密度隨緯度和高度分布圖，結(jié)合圖1 可知該時(shí)刻Kp 指數(shù)為2.7，地磁處于較平靜階段.電子密度沿緯度變化較為平滑，電子密度峰值處于370 km 附近.隨著緯度的增加，電子密度值逐漸減小，由圖2(b)可知，在42～45°N 沒(méi)有BDS 數(shù)據(jù)，同時(shí)GPS 數(shù)據(jù)稀疏，所以無(wú)法得到其正確的電離層電子密度.由于該時(shí)刻IMFBz南轉(zhuǎn)，瞬間穿透電場(chǎng)(prompt penetration electric field,PPEF)導(dǎo)致東向滲透電場(chǎng)增加，使赤道噴泉效應(yīng)(equatorial fountain effect)增強(qiáng)[25]，高度上存在明顯下降，hmF2由15°N 的370 km 降至35°N 的320 km左右.

圖5(b)為2017-09-08T02:00UT 電子密度分布圖，該時(shí)刻Kp 指數(shù)為8，已經(jīng)達(dá)到了強(qiáng)磁暴的等級(jí).可以看出WU430 測(cè)高儀的hmF2由09-07 的345～385 km 提升至09-08 的385～425 km，其電子密度在09-08 的第一個(gè)電離層擾動(dòng)時(shí)間段強(qiáng)于同時(shí)段的09-07，由 09-07 的 4×1011el·m-3提升至 09-08 的9×1011el·m-3，由此可以明確本次電離層暴為正暴[26].此外，不同緯度地區(qū)電離層暴相有所不同，圖5(a)～(b)中，低緯地區(qū)的NmF2顯著降低，這可能與赤道異常結(jié)構(gòu)的動(dòng)力學(xué)變化有關(guān)[7].在赤道異常峰區(qū)以北的較高緯區(qū)(地理緯度22°N 以北)，主要表現(xiàn)為正暴相；在赤道異常峰區(qū)，主要表現(xiàn)為負(fù)暴相.

圖5(c)和5(d)為2017-09-07—08T15:00UT 電子密度分布圖，該時(shí)刻Kp 指數(shù)分別為2 和7.3，09-08 該時(shí)刻再次處于強(qiáng)擾動(dòng)階段(第二個(gè)電離層擾動(dòng)時(shí)間段).觀察發(fā)現(xiàn)，由于電離層電子密度受太陽(yáng)活動(dòng)的影響，NmF2普遍較??；此外，與09-08T02:00UT 變化情況相反，由于武漢附近受EPB 影響[15]，09-08T15:00UT的NmF2較09-07 有所下降，由09-07 的5×1011el·m-3降至09-08 的3×1011el·m-3，與圖3(d)的變化一致；高度上的變化同樣有所體現(xiàn)，由09-07 的370 km 提升至09-08 的400 km.

3.3 TEC 層析結(jié)果分析

為了討論BDS 數(shù)據(jù)對(duì)電離層層析算法反演TEC 的影響，本節(jié)對(duì)比中國(guó)及周邊低緯地區(qū)2017-09-08T14:00UT(22:00LT)GPS+BDS 和GPS重構(gòu)TEC 分布圖像，如圖6 所示.通過(guò)前面圖3 的分析可知2017-09-08T13:30—14:30UT 期間存在幅度較大的電離層TEC 擾動(dòng).

由文獻(xiàn)[15]可知，2017-09-08T12:45—17:00UT在赤道附近地區(qū)伴有明顯的EPB 現(xiàn)象，并以電子密度耗竭的形式出現(xiàn).由于電離層擾動(dòng)和周圍東向的極化電場(chǎng)同時(shí)促進(jìn)了EPB 緯度方向的發(fā)展[27]，圖6(a)和(b)均在當(dāng)?shù)匾归g15～45°N 和80～120°E 地區(qū)附近出現(xiàn)了EPB，TEC 值大約為5 TECU.圖6(a)中發(fā)現(xiàn)EPB 能以較完好的兩個(gè)單獨(dú)的平行結(jié)構(gòu)出現(xiàn)，TEC 沿地磁場(chǎng)線從東南到西北的變化與文獻(xiàn)[15]中TEC 變化一致；而6(b)中低緯地區(qū)(90～120°E,15～25°N)的TEC 無(wú)法呈現(xiàn)連續(xù)的耗竭結(jié)構(gòu).結(jié)合圖2(a)的BDS IPP 分析可得，BDS 數(shù)據(jù)的融合應(yīng)用彌補(bǔ)了GPS 層析結(jié)果的不足，進(jìn)一步提高了中國(guó)及周邊低緯地區(qū)電離層TEC 異常結(jié)構(gòu)的反演精度.由于篇幅所限，有關(guān)EPB 的三維分析及其周邊的閃爍現(xiàn)象將在另一篇文章中詳細(xì)分析.

圖6 中國(guó)及周邊地區(qū)2017-09-08T14:00UT 的TEC 分布圖像Fig.6 TEC maps over China and adjacent areas at 14:00UT on Sep.8,2017

4 結(jié) 論

本文借助一種基于GPS 數(shù)據(jù)的CIT 技術(shù)，并融合BDS 系統(tǒng)中的三種類型衛(wèi)星數(shù)據(jù)對(duì)2017-09-07—09 電離層暴期間15～55°N、70～140°E 范圍內(nèi)的中國(guó)及周邊地區(qū)上空電離層進(jìn)行重構(gòu).通過(guò)對(duì)地球物理數(shù)據(jù)和地磁指數(shù)數(shù)據(jù)研究發(fā)現(xiàn)，2017-09-08 電離層在01:00—02:00UT 和12:00—17:00UT 分別發(fā)生了兩次擾動(dòng).本文采用以GPS 和GPS+BDS 為數(shù)據(jù)源的CIT 結(jié)果，分別與位于武漢(WU430 站)和三亞(SA418 站)地區(qū)的兩個(gè)測(cè)高儀的NmF2數(shù)據(jù)進(jìn)行對(duì)比，來(lái)討論融合BDS 數(shù)據(jù)的CIT 技術(shù)在電離層平靜期和擾動(dòng)期的適用性和層析結(jié)果的精度.通過(guò)分析得到：GPS+BDS CIT 算法在電離層平靜期和暴時(shí)均適用，并觀察到兩次擾動(dòng)期間NmF2趨勢(shì)均有明顯變化，在地面數(shù)據(jù)覆蓋相對(duì)充分的情況下，同GPS 層析結(jié)果相比，GPS+BDS 層析算法得到的NmF2更接近測(cè)高儀的變化趨勢(shì)且誤差更小，初步得到了該方法在NmF2精度上的提高.

其次，利用GPS、GPS+BDS 層析算法獲取SA418和WU430 站上空F2層峰值高度hmF2數(shù)據(jù)，并與兩站觀測(cè)到的hmF2分析對(duì)比，可知BDS 數(shù)據(jù)對(duì)暴時(shí)電離層垂向分布精度可進(jìn)一步改進(jìn)，并通過(guò)SA418 站和WU430 站的Mean 和RMSE 對(duì)比發(fā)現(xiàn)，其改善程度依然受BDS IPP 分布的影響，在BDS IPP 分布更密集的地區(qū)，得到的層析結(jié)果中hmF2的精度更高.此外，利用GPS+BDS 層析結(jié)果進(jìn)一步分析了2017-09-07—08 各自02:00UT 和15:00UT 沿114°E 經(jīng)度鏈上電子密度隨緯度與高度的變化，進(jìn)一步驗(yàn)證了本次磁暴為正暴，而WU430 站15:00UT 電子密度的下降可能與等離子體泡有關(guān)；觀察到兩次擾動(dòng)均對(duì)該地的hmF2有一定的提高.

最后，通過(guò)對(duì)中國(guó)及鄰近地區(qū)上空GPS+BDS 反演與GPS 反演的TEC 分布圖發(fā)現(xiàn)，在2017-09-08T14:00UT 中低緯度地區(qū)存在沿地磁場(chǎng)線向西北方向分布的等離子體泡，而B(niǎo)DS 數(shù)據(jù)可優(yōu)化中低緯地區(qū)電離層分布結(jié)構(gòu)，進(jìn)一步驗(yàn)證了融合BDS 數(shù)據(jù)的CIT技術(shù)對(duì)反演電離層TEC 異常結(jié)構(gòu)的提高.BDS 目前已成為電離層探測(cè)的重要手段之一，但通過(guò)中國(guó)及周邊地區(qū)內(nèi)IGS 站和IPP 的對(duì)比來(lái)看，BDS 衛(wèi)星的數(shù)量以及IGS 站比較有限，對(duì)區(qū)域內(nèi)反演所改善的范圍造成一定的影響.在本文中，改善范圍更偏重中國(guó)東部和東南部.本文僅使用了2017 年BDS 衛(wèi)星和IGS 站，數(shù)據(jù)量仍有許多問(wèn)題需要解決.下一步考慮使用近期BDS 數(shù)據(jù)進(jìn)行反演，并考慮如何進(jìn)一步提高四維CIT 方法的空間分辨率.

致謝本文GPS/BDS 觀測(cè)數(shù)據(jù)從國(guó)際GNSS 服務(wù)(international GNSS service,IGS)網(wǎng)站獲取.地球物理數(shù)據(jù)分別從https://omniweb.gsfc.nasa.gov/form/omni_min.html 和 http://spidr.ionosonde.net/spidr 網(wǎng)站獲取，測(cè)高儀數(shù)據(jù)來(lái)自Digital Ionogram Database官網(wǎng)http://giro.uml.edu/didbase/scaled.php，CODE數(shù)據(jù)從http://ftp.aiub.unibe.ch/CODE/網(wǎng)站獲取.層析軟件 MIDAS 由英國(guó)巴斯大學(xué) Mitchell 教授提供，作者在此深表感謝.