尹萍 纵瑞龙
（中国民航大学电子信息与自动化学院, 天津 300300）
地基增强系统(Ground-based Augmentation Systems, GBAS)是一种基于L1频率GNSS信号的精密进近和着陆导航系统[1]。当电离层出现异常时，其产生的TEC空间梯度会导致GBAS产生几米或者十几米不同范围的误差。GBAS进近服务D型(GBAS approach service type D, GAST-D)可以支持Ⅲ类飞机的精密进近和着陆，GBAS在引导飞机进近的过程中需要极高的稳定性。电离层是影响GBAS安全运行的最大误差源，因此国际标准要求地面和机载完整性监测仪在实时操作中监测电离层异常[2]。目前，很多国家一直在研究他们所在区域的电离层异常，并建立了自己国家的电离层威胁模型。美国局域增强系统(Local Area Augmentation System, LAAS)的电离层威胁模型是基于确定美国大陆(CONUS)观测到的极端电离层延迟梯度建立的，其最大斜电离层延迟梯度为413 mm/km[3]。Minchan Kim等人估计了2000—2004年韩国地区的异常电离层梯度，并为该地区定义了一个电离层威胁模型，其TEC空间梯度上限为160 mm/km[4]。Emilien等人在欧洲地区进行了类似的研究，在该地区检测到的最大梯度为403 mm/km[5]，其梯度值要低于CONUS。上述三个地区主要位于中纬度地区。在巴西检测到最大梯度约为850 mm/km[6]，这可能与低纬地区出现的等离子体泡有关。
等离子体泡是指在赤道和低纬地区夜间电离层经常观察到的等离子体密度耗竭现象[7]。由离子扰动引起的瑞利泰勒不稳定机制，被普遍认为是等离子体泡出现的机理。它是一种局部电离层耗竭，其边缘有陡峭的梯度[8]。Satio等人利用三维等离子体泡模型研究了等离子体泡对GBAS的影响，仿真结果表明单个等离子体泡可造成垂直定位误差超过10 m[9]。因此，研究等离子体泡产生的TEC空间梯度是非常重要的。
文献[10]首先采用时间步长法，计算出电离层梯度参数估计值，然后基于膨胀电离层梯度的高斯分布函数获取优化后的电离层参数标准值，发明出一种电离层梯度参数的确定方法、装置及系统。王志鹏等人利用短基线双站法计算出了广东省的电离层异常产生的TEC空间梯度[11]。Saito等人为了估计电离层延迟梯度，使用由Fujita等人开发的单频载波和码辅助技术(single-frequency carrier-based and code-aided technique, SF-CBCA)，并进一步细化梯度监测，研究了日本石垣岛观测到的等离子体泡有关的电离层延迟的TEC空间梯度特性，并应用在低纬地区GBAS的定位误差修正[12]。
随着GNSS理论的快速发展，电离层层析技术(computerized ionospheric tomography, CIT)展现出探测电离层电子密度三维分布和变化方面的优势。2004年Yin等人利用多仪器数据分析系统(multiinstrument data analysis system, MIDAS)层析反演分析了美国地区2000-06风暴时的电离层电子密度分布图像，得到了较好的结果[13]；
2019年尹萍等人利用多分辨率层析技术研究了电离层暴时大尺度电离层行扰(large-scale travelling ionospheric disturbance,LSTID)的特征[14]。根据文献[15]记载，2017-09-08T12:45—17:00UT这个时间段内，在中国及周边地区有明显的等离子体泡现象出现。
本文在英国巴斯大学MIDAS CIT算法[16]基础上，利用IGS官网提供的GNSS数据对2017-09-08中国及周边中低纬地区出现的电离层暴进行层析成像，并对出现的等离子体泡利用短基线双站法和时间步长法对其产生的TEC空间梯度进行分析。结果表明，在CIT的辅助下，能够较好地分析等离子体泡产生的TEC空间梯度特性。
1.1 多分辨率层析成像技术
本文选用MIDAS多分辨率层析技术反演出电离层电子密度的三维分布图像，层析输入数据是采用双频载波相位进行计算，然后经过预处理（诸如相位平滑伪距等）和对TEC进一步采样（时间采样率5 min）得到的结果。MIDAS单分辨率层析技术通常只使用一次反演来重建大分辨率电离层结构，然而MIDAS多分辨率层析成像可以通过多次反演的方式实现较高精度的层析结果。
本文采用多次反演的层析算法[17]。首先设定反演级数level=[n,m]，第一次反演获取较低分辨率的层析图像，第二次反演是在第一次反演结果的基础上将反演级数在空间上以2的n次幂倍降低、时间上以2的m次幂倍降低，这样可以获得较高分辨率的层析图像。例如，反演级数在空间和时间上都选为1次，即level=[1, 1]，如果第一次空间分辨率为4°×4°、时间间隔为30 min，那么第二次反演将获得空间分辨率为2°×2°、时间间隔为15 min的电离层层析图像，以此类推。
本文中层析的格网划分为15~35°N、100~124°E，每一格网区域间隔经纬度相差4°。本文在多分辨率层析算法中，空间分辨率参数设定为8°×8°，时间间隔设置为40 min，反演级数level设置为[3,3]。即多分辨率层析反演结果的空间分辨率从8°×8°到1°×1°，时间间隔从40 min到5 min，TEC重构所选用的数据均是层析1°×1°的结果。
1.2 TEC空间梯度计算
1.2.1 短基线双站法计算TEC空间梯度
短基线双站法计算TEC空间梯度模型如图1所示，其算法原理为：将S1接收站和S2接收站观测到同一颗卫星K的电离层延迟值作差，然后除以两个接收机之间的基线距离∆S得到TEC空间梯度。在利用短基线双站法计算TEC空间梯度时，卫星和接收机会带来一定的硬件偏差，影响计算出来的TEC空间梯度的准确性，需要对接收机和卫星的硬件偏差进行修正。
图1 短基线双站法和时间步长法模型Fig.1 Short baseline station-pairs method and time-step method
首先，将电离层看成高度为350 km的薄壳模型[18]，式(1)为利用时间采样率为30 s的GPS观测数据计算的电波视线方向上的斜向TEC(slant TEC,STEC)[19]。
式中：f1、f2为信号频率，GPS选用f1=1 575.42 MHz，f2=1 227.60 MHz；
L1、L2为载波相位；
λ1、λ2为波长。
计算出来的STEC中包含一些硬件偏差，为了使其结果更加准确，参照文献[20]进行了一系列的硬件偏差修正。
为了计算电离层穿刺点(ionospheric pierce point,IPP)处的电离层延迟值，需要将STEC转换为垂向TEC(vertical TEC, VTEC)，计算VTEC的倾斜因子表达式[21]如下：
式中：Re为地球半径；
θ为卫星仰角；
h为电离层薄壳模型高度，此处选定350 km。由此可以求出VTEC和STEC之间的转换关系：
电离层延迟值计算公式为
式中，f=1575.42MHz。本文所说的TEC均指VTEC.
利用短基线双站法计算TEC空间梯度的公式为
式中，Iion,S1、Iion,S2分别为两个接收机在IPP处的电离层延迟值。
1.2.2 时间步长法计算TEC空间梯度
时间步长法计算TEC空间梯度模型如图1所示，其基本原理是利用一个IGS接收站观测同一颗卫星在历元t1时刻的电离层延迟值和在历元t2时刻的电离层延迟值进行差值计算，然后除以两个不同历元时刻IPP之间的距离。时间步长法计算TEC空间梯度的表达式为
式中：Iion(t1)、Iion(t2)分别为两个不同历元时刻IPP处的电离层延迟值；
dIPP为两个历元时刻IPP之间的距离。
与短基线双站法相比，时间步长法与GBAS地面站和飞机进近场景相似度很小，所以此方法与GBAS联系起来并不直观。尽管如此，因为时间步长法捕获了与GBAS地面站接收机和飞机接收机相同的电离层效应，所以利用时间步长法计算TEC空间梯度时，卫星和接收机的硬件偏差可以被消除。
时间步长法还可以通过改变两个历元之间的时间间隔，得到所需的IPP之间的距离，时间间隔一般取100 s左右。
本文中主要是用短基线双站法计算TEC空间梯度，而时间步长法主要用来对计算结果进行验证。这是由于短基线双站法在计算TEC空间梯度时，卫星和接收机产生的硬件偏差对计算TEC空间梯度存在干扰。虽然本文对卫星和接收机产生的硬件偏差进行了修正，但为了确保修正的准确性，使用了两种方法并进行比较。
1.3 数据选取
1.3.1 多分辨率层析输入数据
图2给出了2017-09-08T13:30—14:30UT中国及周边中低纬地区（15~35°N，100~124°E）在假定电离层薄壳模型高度为350 km时，沿地面站与GPS/BDS卫星路径上电离层IPP的轨迹，即输入数据的覆盖范围。可以看出：GPS IPP覆盖区域较大，主要集中在19~27°N和108~120°E；
BDS IPP覆盖范围主要集中在19~27°N和108~116°E。BDS IPP分布能够弥补GPS IPP分布空白部分，为电离层层析反演算法提供了更多的数据源。
图2 2017-09-08T13:30—14:30UT中国及周边中低纬地区350 km上空处IPP分布Fig.2 IPP tracks over 350 km above middle and low latitude areas in China and its surrounding areas at 13:30—14:30UT on September 8, 2017
1.3.2 短基线双站法/时间步长法选取的数据
为研究中国及周边中低纬地区出现的等离子体泡产生的TEC空间梯度，从香港地区筛选11个GNSS接收站来研究等离子体泡产生的TEC空间梯度。
香港地区位于22.1~22.4°N、113.8~114.2°E，处于中国低纬地区，能够用来分析中国及周边中低纬地区等离子体泡产生的TEC空间梯度。
本次研究采用2017-09-08从香港连续运行参考站(Continuously Operating Reference Station, CORS)网记录的数据进行分析。本研究中计算TEC空间梯度所用到的TEC均由香港接收站提供时间采样率为30 s的双频GPS观测数据计算出。图3展示了香港地区11个GNSS接收站以及香港国际机场的位置。
图3 香港11个GNSS接收站分布和国际机场位置Fig.3 Distribution of eleven GNSS receive stations and location of the International Airport in Hong Kong
由文献[15]可知，2017-09-08T12:45—17:00UT在赤道附近地区伴有明显的等离子体泡出现，下面针对出现的等离子体泡现象对其所产生的TEC空间梯度进行分析。
2.1 TEC层析图像分布
利用多分辨率层析成像技术对2017-09-08T13:10UT、13:15UT、13:20UT、13:55UT、14:00UT、14:05UT六个时刻进行层析结果分析，其结果如图4所示。
图4 2017-09-08六个时刻的TEC分布图像Fig.4 TEC distribution image of six times on September 8, 2017
利用香港地区提供的GNSS接收站观测数据进行大量的实验，结果发现最大梯度时刻主要出现在13:20UT、14:05UT左右。为了更加清晰地了解异常TEC空间梯度产生的原因，选取13:10—13:20UT和13:55—14:05UT两个时段的TEC层析图像进行分析，这两个时段的TEC分布图像覆盖了出现最大梯度的时刻，可以看出有明显的等离子体泡现象出现。图4中红色箭头所指为等离子体泡异常结构，其内部TEC与背景电离层TEC相差5~15 TECU[15,22]。
2.2 等离子体泡空间结构特征
选取14:05UT层析电子密度三维分布图像来分析等离子体泡空间结构特征。首先，层析算法设置的电离层格网在垂向上为40 km间隔，高度范围为100~820 km，共计19层。其次，为了简化本节只选用了其中最具代表性的三层来展现电子密度的分布，分别为底部电离层(140 km)、峰值电离层(380 km)、顶部高度层(780 km)，结果如图5所示。注意：为了清晰展示等离子体泡，三幅图上的电子密度色棒范围做了相应调整。
图5 2017-09-08T14:05UT不同高度上空处的电子密度分布结果Fig.5 Electron density distribution at different altitudes at 14:05UT on September 8, 2017
从图5(a)可以明显地看出在140 km处电子密度出现耗竭，如红色箭头指示的区域。同理，从图5(b)和图5(c)也看到了同区域的电子密度耗竭，因此可以直观地看出等离子体泡的电子密度耗竭是从电离层底部延伸到电离层顶部区域。由于等离子体泡内部的耗竭在等离子体泡边缘处会产生一定的TEC空间梯度，接下来研究等离子体泡边缘处产生的TEC空间梯度大小。
2.3 TEC空间梯度结果分析
从香港地区11个GNSS接收站中任选两个接收站组成一个站对，利用短基线双站法解算TEC空间梯度，由于不同接收机之间观测到的不同卫星IPP轨迹不尽相同，其观测到不同卫星产生的TEC空间梯度值也不一样。从中筛选出较大的TEC空间梯度值，紧接着进行分析是否由于接收机产生的故障导致大的TEC空间梯度值，然后剔除接收信号质量不好的接收机。在表1中统计了所有GNSS接收站与所有可见卫星之间产生的最大TEC空间梯度值及其出现时刻。
表1 2017-09-08香港地区11个站对的最大TEC空间梯度值及其出现时刻Tab.1 The maximum TEC spatial gradient value and its occurrence times of eleven stations in Hong Kong on September 8, 2017
2.3.1 HKPC和HKMW站对观测G32号卫星异常TEC空间梯度分析
为了分析表1中香港地区11个站对观测到最大TEC空间梯度产生的原因，选取表1中最大TEC空间梯度所对应的站对(HKPC和HKMW)进行分析。HKPC和HKMW两接收站之间的基线距离为4.82 km。
为了更加直观地分析两个接收机观测同一颗卫星之间产生的电离层延迟偏差，图6展现了HKPC和HKMW接收站观测G32号卫星随时间变化的双频电离层延迟值。
图6 2017-09-08两站观测G32号卫星的电离层延迟曲线Fig.6 The ionospheric delay curve of G32 satellite observed by the two stations on September 8, 2017
图7中蓝色曲线为HKPC和HKMW接收站观测G32号卫星的站间TEC空间梯度值随时间变化曲线，黑色曲线为利用时间步长法(选择HKMW接收站进行观测)接收站观测G32号卫星不同历元下IPP之间的TEC空间梯度值随时间变化曲线，为了使相邻IPP之间的距离分布更加接近HKPC和HKMW接收站之间的基线距离，设置时间抽样间隔为60 s，使得IPP之间的分布距离大约为4.93 km。
图7 2017-09-08两站观测G32号卫星的TEC空间梯度和仰角曲线Fig.7 The TEC spatial gradient and elevation curve of the G32 satellite observed by the two stations on September 8, 2017
图7中用红色圆圈表示HKPC接收站观测到G32号卫星随时间变化的仰角，绿色三角表示HKMW接收站观测到G32号卫星随时间变化的仰角。从图7可以看出两个站观测G32号卫星随时间变化的仰角几乎是重合的，其原因是由于两个接收站的地理位置分布很近，所以当观测到同一颗卫星时随时间变化的仰角几乎是一样的。
由图7还可以看出：使用短基线双站法计算出的TEC空间梯度结果在14:06UT左右出现了最大峰值，其TEC空间梯度高达133.16 mm/km；
利用时间步长法(HKMW)计算的结果与短基线双站法计算的结果在趋势上以及最大峰值上相吻合。为了更加清晰地了解异常TEC空间梯度在此时刻出现的原因，我们利用多分辨率层析成像进行辅助，重构出14:05UT时刻的TEC分布图像，如图8所示。
图8 2017-09-08T14:05UT多分辨率层析法重构的TEC分布图像Fig.8 TEC distribution image reconstructed by multiresolution tomography at 14:05UT on September 8, 2017
图8中红色圆点代表HKPC站点位置，粉色圆点代表HKMW站点位置，由于HKPC和HKMW接收站之间的基线距离为4.82 km，分布相对较近，而在图8中1°大约为108 km的距离，远大于这两个接收站之间的基线距离，所以图8中两个接收站(用不同颜色圆点代表)的位置在背景图上看起来几乎是重叠在一起的。
图8中红色三角标记为HKPC站在14:05UT观测到G32号卫星的IPP位置，粉色三角标记为HKMW站在14:05UT观测到G32号卫星的IPP位置，可以看出在14:05UT时刻，HKPC和HKMW接收站观测到G32号卫星的IPP处于等离子体泡边缘，结合图7可以推断14:06UT出现的峰值是由等离子体泡边缘引起的。
2.3.2 HKSS和HKST站对观测G10号卫星异常TEC空间梯度分析
为了验证观测其他卫星出现大的TEC空间梯度是否由等离子体泡边缘引起的，从表1中选取TEC空间梯度最小的两个接收站用于分析。选取了HKSS和HKST接收站去观测G10号卫星，HKSS和HKST两站之间的基线距离为9.62 km。
为了更加直观地分析两个接收机观测同一颗卫星之间产生的电离层延迟偏差，图9展示了HKSS和HKST接收站观测G10号卫星随时间变化的双频电离层延迟值。
图9 2017-09-08两站观测G10号卫星的电离层延迟曲线Fig.9 The ionospheric delay curve of G10 satellite observed by the two stations on September 8, 2017
图10中蓝色曲线为HKSS和HKST接收站观测G10号卫星的站间TEC空间梯度值随时间的变化，黑色曲线为利用时间步长法(选择HKST接收站进行观测)接收站观测G10号卫星不同历元下IPP之间TEC空间梯度值随时间的变化，为了使IPP之间的距离分布更加接近HKSS和HKST接收站之间的基线距离，设置时间抽样间隔为120 s，使得IPP之间的分布距离大约为9.48 km。
图10 2017-09-08两站观测G10号卫星的TEC空间梯度和仰角曲线Fig.10 The TEC spatial gradient and elevation curve of the G10 satellite observed by the two stations on September 8, 2017
图10中用红色圆圈表示HKSS接收站观测到G10号卫星随时间变化的仰角，绿色三角表示HKST接收站观测到G10号卫星随时间变化的仰角。从图10可以看出两个站观测G10号卫星随时间变化的仰角几乎是重合的，其原因同图7中仰角重合所作的解释。
由图10还可以看出：使用短基线双站法计算出的TEC空间梯度结果在13:20UT左右出现了最大峰值，其TEC空间梯度高达99.72 mm/km；
利用时间步长法(HKST)计算的结果与短基线双站法计算的结果在趋势上以及最大峰值上相吻合。下面利用多分辨率电离层层析成像重构13:20UT的TEC分布图像，分析梯度产生原因，结果如图11所示。
图11 2017-09-08T13:20UT多分辨率层析法重构的TEC分布图像Fig.11 TEC distribution image reconstructed by multiresolution tomography at 13:20UT on September 8, 2017
图11中红色圆点代表HKSS站点位置，粉色圆点代表HKST站点位置，由于HKSS和HKST接收站之间的基线距离为9.62 km，分布相对较近，在图中1°大约为108 km的距离，远大于这两个接收站之间的基线距离，所以在图11中两个接收站(用不同颜色圆点表示)的位置在背景图上看起来几乎是重叠在一起的。
图11中红色三角标记为HKSS站在13:20UT观测到G10号卫星的IPP位置，粉色三角标记为HKST站在13:20UT观测到G10号卫星的IPP位置。可以看出在13:20UT时刻，HKSS和HKST接收站观测到G10号卫星的IPP也处于等离子体泡边缘，由此推断产生的峰值梯度是由等离子体泡边缘造成的。
本文借助一种基于GPS+BDS的多分辨率层析成像技术，采用多次反演的层析算法对2017-09-08电离层暴期间15~35°N、100~124°E范围内的中国及周边地区上空电离层进行重构。结果发现在13:10—13:20UT和13:55—14:05UT两个时段内伴有等离子体泡现象出现。针对中国及周边中低纬地区出现的等离子体泡现象，本文从香港CORS网获取多个GNSS接收站数据对出现的等离子体泡产生的TEC空间梯度特性分析。
首先把香港地区11个GNSS接收站两两进行配置计算出观测到不同卫星产生的TEC空间梯度并进行统计筛选，然后选择比较大的TEC空间梯度值。通过统计，发现出现在香港空域附近的等离子体泡产生的TEC空间梯度最大为133.16 mm/km，且结合短基线双站法和时间步长法得到了很好的验证，在CIT辅助下，更能直观地分析出现较大TEC空间梯度出现的原因。
本文运用多分辨率层析成像技术进行辅助，采用短基线双站法和时间步长法能够较为准确地判定等离子体泡边缘处的TEC空间梯度，这为GBAS引导的飞机着陆的安全性提供了很好的电离层威胁模型参数，同时为消除卫星定位误差奠定了基础。此外，由于低纬地区IGS数量有限，而短基线双站法需要分布密集的IGS站，但本文只选择香港地区的GNSS接收站进行分析，所以数据量仍有许多问题需要解决。下一步考虑使用BDS的中圆轨道卫星结合GPS卫星对等离子体泡进行全方位研究。
致谢：本文香港地区GPS观测数据从香港政府网（ftp://ftp.geodetic.gov.hk/）网站上获取，精密星历(SP3)文件从(http://ftp.aiub.unibe.ch/CODE_MGEX/CODE/)网站上获取。硬件偏差DCB文件从(http://ftp.aiub.unibe.ch/CODE)网站上获取。本文所使用的MIDAS层析成像软件是由英国巴斯大学Mitchell教授提供的，作者在此深感表谢。
猜你喜欢 接收站层析电离层 LNG接收站的发展趋势煤气与热力(2021年9期)2021-11-06一种电离层TEC格点预测模型空间科学学报(2021年4期)2021-08-30LNG接收站扩建工程低压泵国产化应用管理煤气与热力(2021年7期)2021-08-23Kalman滤波估算电离层延迟的一种优化方法空间科学学报(2021年2期)2021-07-21犬细小病毒量子点免疫层析试纸条的研制中国动物传染病学报(2021年3期)2021-07-21LNG接收站工程项目设计进度的控制化工管理(2021年7期)2021-05-13电离层对中高轨SAR影响机理研究雷达学报(2017年6期)2017-03-26A族链球菌胶体金免疫层析试纸条的制备及应用分子影像学杂志(2015年3期)2015-12-04Linux Shell语言在电离层解算中的应用导航定位学报(2015年2期)2015-06-05LNG接收站贸易交接在线取样技术天然气与石油(2015年2期)2015-02-28
