基于GNSS 30s采样频率数据的电离层相位闪烁因子构建方法

技术领域

本发明涉及一种电离层闪烁因子构建方法，尤其适用于一种监测电离层闪烁天文灾害使用的基于GNSS 30s采样频率数据的电离层相位闪烁因子构建方法。

技术背景

电离层闪烁是赤道和极区频发的一种天文灾害，是近年来空间科学面临的一大难题，给GNSS系统位置、导航与时间服务的稳定性带来巨大挑战，例如2003年的万圣节电离层闪烁事件造成GNSS导航卫星中断服务，严重干扰了电力、银行系统以及北极航线的航行安全，造成经济损失。由此可见，电离层闪烁灾害对人类社会可持续发展构成了重大威胁，科学认识电离层闪烁发生规律，监测及预报电离层闪烁天文灾害对于防灾减灾具有重大现实意义。然而监测电离层闪烁通常需要电离层闪烁监测接收机，简称ISMR，其通常运行在50Hz的采样频率，导致其需要较大的存储空间，且价格较贵，进而限制了ISMR的布站数量。当前，世界各国经过近30年的时间仅建设了约150个ISMR测站，无法形成对全球的全面覆盖，难以利用ISMR对赤道和极区地区的电离层闪烁形成完备监测。相比于ISMR，采样频率通常在30s的测地型接收机分布十分广泛，在全球布设超过1万个观测台站，基本可以实现对全球的覆盖。建立基于30s采样频率的电离层闪烁因子有利于全面的监测全球电离层闪烁的变化，同时测地型接收机测站的设立时间更久，可以提供更长时间序列的观测数据，有助于对电离层闪烁变化的长周期项监测，提高电离层闪烁模型的精度。但随着采样频率的降低，信号更容易受到其他观测误差的干扰，难以提取出其中的电离层闪烁误差，同时其对较小电离层闪烁监测的响应也会降低。因此如何基于30s采样频率的GNSS数据构建合理的电离层闪烁因子是实现全球电离层闪烁完备监测面临的瓶颈难题之一。

目前，不少学者围绕基于30s采样频率的GNSS观测数据构建电离层闪烁因子展开研究。如Pi等提出基于无几何组合观测值提出了总电子含量变化率指数，简称ROTI，可以较好减弱接收机钟差的影响，已经成为电离层闪烁监测领域应用最为广泛的标准闪烁因子之一；ROTI可以给出每个观测历元的闪烁因子，在此基础上Sanz等通过对一段时间内所有可视卫星获得的TEC值以卫星高度角定权的方式进行加权平均从而提出了另一种闪烁因子——TEC起伏，其可以用来预报发生闪烁的时间段和区域，已经被欧洲地球静止导航重叠服务选为表征电离层活动的参数之一；Juan等将无电离层组合残差的标准差作为监测电离层闪烁的因子并提出了对应的修正无电离层组合中的接收机钟差方法；Ahmed等利用1Hz的斜路径总电子含量和小波变换技术提出了一种类相位闪烁因子，该因子被认为可以有效的预警欧洲部分北极区域的电离层闪烁事件。但以上因子给出的均是闪烁在两个频点共同作用的结果，难以利用这些闪烁因子研究电离层闪烁对GNSS每个频点信号的影响，此外以上研究缺乏对30s采样频率观测数据的阈值研究，难以利用这些因子准确监测电离层闪烁的发生。由此可见，虽然利用上述研究成果可以使得测地型接收机应用到电离层监测当中，但是目前针对基于GNSS 30s采样频率观测数据的电离层相位闪烁因子构建及其经验阈值确定尚缺乏深入研究。

发明内容

针对上述技术问题，本发明提供了一种监测电离层闪烁天文灾害，将分布广泛的低采样频率测地型接收机引入电离层闪烁监测领域，解决了难以在全球范围内完备监测电离层闪烁天文灾害的技术难题。

1.为了实现上述技术目的，本发明的基于GNSS 30s采样频率数据的电离层相位闪烁因子构建方法，其特征在于步骤为：

S1、对30s采样频率的GNSS原始载波相位观测数据进行大地测量趋势分离，修正载波相位中的卫星和测站之间的几何距离以及固体潮、天线相位中心、卫星钟差、相位缠绕和对流层静力学延迟误差；

S2、利用宽项组合、无几何关系组合和无电离层组合对GNSS 30s采样频率载波相位数据进行周跳探测与修复；

S3、利用修复周跳后的载波相位数据进行精密单点定位PPP，并结合随机游走模型估计并改正对流层湿延迟；

S4、利用精密单点定位PPP和无电离层组合方法改正接收机钟差，获得仅受电离层影响的观测值残差，

S5、利用小波变换的方法从观测值残差中提取电离层闪烁信号，以5分钟作为滑动窗口对提取出的电离层闪烁信号做标准差，构建出闪烁因子；

S6、利用电离层闪烁监测接收机获取电离层相位闪烁因子的相关性，结合构建的闪烁因子与电离层相位闪烁因子的相关性反向优化电离层闪烁特征频带区间、Morse小波的对称参数与时间带宽积；

S7、将电离层闪烁监测接收机提供的相位闪烁因子作为参考值，利用偏斜分析方法确定所构建闪烁因子与参考值之间90％最接近0值的差值对应的数值上边界范围，将该边界值作为所构建闪烁因子的经验阈值；

S8、以太阳活动、地磁活动的物理参数和测站坐标作为自变量，以经验阈值作为因变量，利用多元非线性回归的方法构建出经验阈值模型。

具体的：对30s采样频率的GNSS原始载波相位数据进行大地测量趋势分离：首先利用国际GNSS服务提供的精密轨道文件改正卫星和测站之间的几何距离，然后利用二阶潮汐模型改正固体潮误差，利用IGS提供的天线相位模型改正卫星和接收机天线相位中心误差，利用IGS提供的精密钟差文件改正卫星钟误差，改正相位缠绕误差，利用UNB3m模型改正对流层静力学延迟误差。

对GNSS 30s采样频率载波相位数据进行周跳探测与修复方法包括：

利用宽相组合和无几何关系组合对载波相位数据进行周跳探测与初次修复；

若在强电离层闪烁的条件下修复1或2周的小周跳，则通过构建无电离层组合观测值的方法对经过初次修复的载波相位数据进一步修复周跳，从而消除99％的电离层误差，获得较为精准的周跳探测检测值。

对流层湿延迟改正包括以下两步：

首先利用随机游走模型对载波相位数据中的对流层湿延迟进行初步估计，然后利用静态精密单点定位PPP方法对对流层湿延迟做进一步估计，其中将经过大地测量趋势分离和周跳修复后的载波相位数据残差作为PPP的观测量，将测站坐标作为已知量，仅将接收机钟差和对流层湿延迟作为待求量，利用扩展卡尔曼Kalman滤波方法对PPP进行解算，获得对流层湿延迟的改正量。

对接收机钟差改正包括以下两步：

首先利用静态精密单点定位给出的接收机钟差值做初步的估计，然后利用经过大地测量趋势分离、周跳和对流层湿延迟修正之后的载波相位数据残差构建无电离层组合，由于该组合值在历元间做差后可以消除整周模糊度和电离层闪烁的影响，其残差仅包括历元简差分后的接收机钟差，因此只要对该残差在时间上积分后即可获得接收机钟差的细节部分。

利用小波变换方法构建电离层闪烁因子及经验参数确定：利用Morse小波基将残差变换到频率域，并对特征频带区间的信号进行小波逆变换，提取出电离层闪烁信号；设定5分钟的滑动窗口，并对窗口内的电离层闪烁信号取标准差，完成基于测地型接收机的电离层相位闪烁因子的构建。

30s采样频率中的电离层闪烁特征频带区间、Morse小波的对称参数和时间带宽积的确定方法如下：

根据以下要求确定30s采样频率中的电离层闪烁特征频带区间：首先将电离层闪烁接收机提供的电离相位闪烁因子作为参考值，由于对30s采样频率的GNSS载波相位数据进行小波变换后，其在小波时频谱图的频率上边界通常小于0.014Hz，因此可假设电离层闪烁对应的上边界频率位于0.010Hz与0.013Hz之间，下边界频率要小于上边界频率，以0.001Hz为最小的频率分辨率，可组成42组可能的频带区间上下边界组合；分析由每组频率组合所获得的闪烁因子与参考值的相关性，从而获得最高相关性的组上下边界频率组合即为30s采样频率的GNSS载波相位数据中的电离层闪烁特征的频带区间；

根据以下要求确定Morse小波的对称参数和时间带宽积：以1为步长选取3至25之间的对称参数，以20为步长选取20至1000之间的时间带宽积，同时需要满足数学表示条件：

利用电离层闪烁监测接收机提供的相位闪烁因子确定所构建闪烁因子的经验阈值的方法如下：

由于采样频率的降低会使得所估计闪烁因子的幅值普遍大于由1Hz采样频率或电离层闪烁接收机给出的闪烁因子，从而使通常的经验阈值无法用于由30s采样频率数据构建的电离层闪烁因子，为此需要确定传统的经验阈值0.2rad增量的专门的经验阈值，以保证所构建的闪烁因子对电离层闪烁的监测的准确性。

为使得所构建闪烁因子可准确探测到90％电离层闪烁，因此认为闪烁因子与电离层闪烁接收机提供的相位闪烁因子做差之后的差值在0值附近90％的因子是准确的，并将累计概率密度为90％所对应的上边界作为传统阈值的增量，通过无电离层闪烁发生时段所构建闪烁因子与参考值之间差值的分布规律确认；

利用偏斜分析，确定差值的可能符合的分布函数：若不存在偏斜，则差值符合高斯分布函数，若存在偏斜，则符合指数高斯分布函数，利用每三小时内的差值，确定分布函数中的参数，并据此确定90％概率密度函数的上边界所对应幅值，即为阈值增量δ；在确定阈值增量δ的基础上，所构建闪烁因子的阈值thres表示如下：

thres＝0.2+δ

其中0.2为利用1Hz或电离层闪烁接收机直接给出的闪烁因子的经验阈值，单位为弧度。

构建闪烁因子阈值模型具体方法如下：

所述阈值不是一个定值，会随测站、太阳活动、地磁的变化而不同，为了在无法提供参考值时，仍可利用所构建闪烁因子进行电离层闪烁监测，需要对阈值进行建模预报；以表征太阳活动的物理参数F10.7、太阳黑子数量和表征地磁活动的Kp、Dst参数以及测站的坐标、地方时作为自变量，以阈值作为因变量，利用多元非线性回归的方法，构建阈值的模型，据此可实现对阈值的预报。

有益效果：

相比于现有利用50Hz高采样频率的电离层闪烁监测接收机开展电离层闪烁监测，本方法引入30s采样频率的普通测地型接收机进行电离层闪烁监测，大幅度降低电离层闪烁监测的成本，结合现有普通测地型接收机分布广泛的特性，有利于对区域电离层闪烁的完备监测，具有广泛的实用性。

附图说明

图1为本发明基于GNSS 30s采样频率数据的电离层相位闪烁因子构建方法的流程示意图。

具体实施方式

下面将结合图和具体实施过程对本发明做进一步详细说明：

如图1所示，本发明的基于GNSS 30s采样频率数据的电离层相位闪烁因子构建方法，步骤为：

首先对30s采样频率的GNSS原始载波相位观测数据进行大地测量趋势分离，修正载波相位中的卫星和测站之间的几何距离以及固体潮、天线相位中心、卫星钟差、相位缠绕和对流层静力学延迟误差；利用宽项HMW组合、无几何关系组合和无电离层组合进行周跳探测与修复；利用无周跳影响的载波相位数据进行精密单点定位PPP，并结合随机游走模型估计并改正对流层湿延迟；利用PPP和无电离层组合方法改正接收机钟差，获得仅受电离层影响的观测值残差，进而利用小波变换的方法提取电离层闪烁信号及构建闪烁因子，并结合所构建闪烁因子与参考值的相关性反向优化电离层闪烁的特征频带区间及Morse小波的对称参数与时间带宽积。为确定与所构建闪烁因子匹配的阈值，本发明利用偏斜分析方法确定了在无闪烁发生时段，所构建闪烁因子与参考值之间差值的分布函数及其有关参数，进而确定90％最接近0值的差值对应的上边界范围，将其作为所构建电离层闪烁因子的经验阈值；为提高经验阈值的适用性，构建了经验阈值与多种物理参数的模型，实现了经验参数的预报，完成电离层闪烁因子的构建及其经验阈值的确定。

具体步骤为：

步骤1：对30s采样频率的GNSS原始载波相位观测数据进行大地测量趋势分离，包括利用国际GNSS服务IGS提供的精密轨道文件改正卫星和测站之间的几何距离，利用二阶潮汐模型改正固体潮误差，利用IGS提供的天线相位模型改正卫星和接收机天线相位中心误差，利用IGS提供的精密钟差文件改正卫星钟误差，利用Wu等人1993提出的方法改正相位缠绕误差，利用UNB3m模型改正对流层静力学延迟误差。

步骤2：周跳探测与修复主要包括如下两步：首先利用宽相HMW组合和无几何关系组合对载波相位数据进行周跳探测与修复，该方法在通常情况下可以修复大多数周跳，但难以修复1或2周的小周跳，尤其在较强电离层闪烁的条件下。为此需要对经过初次修复的载波相位数据进一步通过构建无电离层组合观测值的方法进一步修复周跳，该种方法可消除99％的电离层误差，可以提供较为精准的检测值，但因其会受到不敏感周跳组的影响，导致有少部分周跳无法探测到，通过采用了小波变换的方法来构建电离层闪烁因子从而克服少量周跳的不利影响。

步骤3：对流层湿延迟改正包括以下两步：首先利用随机游走模型对载波相位数据钟的对流层湿延迟进行初步估计，然后利用静态精密单点定位PPP方法对对流层湿延迟做进一步估计，其中将经过大地测量趋势分离和周跳修复后的载波相位数据残差作为PPP的观测量，将测站坐标作为已知量，仅将接收机钟差和对流层湿延迟作为待求量，利用扩展卡尔曼Kalman滤波方法对PPP进行解算，获得对流层湿延迟的改正量。

步骤4：接收机钟差改正包括以下两步：首先利用步骤3中静态精密单点定位给出的接收机钟差值做初步的估计，但在接收机存在钟跳Clock Jump的情况下，PPP的方法无法对接收机钟的细节部分进行准确估计。然后利用经过步骤1、2、3修正之后的载波相位数据残差构建无电离层组合，由于该组合值在历元间做差后可以消除整周模糊度和电离层闪烁的影响，其残差仅包括历元简差分后的接收机钟差，对该残差在时间上积分后即可获得接收机钟差的细节部分。

步骤5：利用小波变换方法构建电离层闪烁因子及经验参数确定；利用Morse小波基将步骤4提供的残差变换到频率域，并对特征频带区间的信号进行小波逆变换，提取出电离层闪烁信号；设定一定长度5min的平移窗口，并对窗口内的电离层闪烁信号取标准差，完成基于测地型接收机的电离层相位闪烁因子的构建。

Morse小波确定如下三个参数：30s采样频率数据中的电离层闪烁特征频带区间、Morse小波的对称参数和时间带宽积。将电离层闪烁接收机提供的相位闪烁因子作为参考值，确定以上参数的方法如下。

由于对30s采样频率的GNSS载波相位数据进行小波变换后，其在小波时频谱图的频率上边界通常小于0.014Hz，因此可假设电离层闪烁对应的上边界频率位于0.010Hz与0.013Hz之间，下边界频率要小于上边界频率，以0.001Hz为最小的频率分辨率，可组成42组可能的频带区间上下边界组合。分析由每组频率组合所获得的闪烁因子与参考值的相关性。可获得最高相关性的组上下边界频率组合即为30s采样频率的GNSS载波相位数据中的电离层闪烁对应的频率区间。

Morse小波的对称参数γ和时间带宽积β需要满足如下条件，

根据以上条件，以1为步长选取3至25之间的对称参数，以20为步长选取20至1000之间的时间带宽积，同时兼顾上式提供的条件，组成待选取最优参数组合的搜索区间。分析每组参数组合获得的闪烁因子与参考值的相关性，选取出相关性较高的几组参数组合以缩小收索区间；进一步通过分析在无电离层闪烁时段，由各组参数组合获得的闪烁因子与参考值之间的插值的分布，找出可提供较为集中在0值附近的插值的参数组合作为较优的参数组合。

步骤6：由30s采样频率数据构建的电离层闪烁因子的阈值的确定；采样频率的降低会使得所估计闪烁因子的幅值普遍大于由1Hz采样频率或电离层闪烁接收机给出的闪烁因子，这使得通常的经验阈值无法用于由30s采样频率数据构建的电离层闪烁因子，为此需要确定其专门的经验阈值，即相比传统的经验阈值0.2rad增量的部分，以保证所构建的闪烁因子对电离层闪烁的监测的准确性。具体确定方法如下。

为使得所构建闪烁因子可准确探测到90％电离层闪烁，我们认为与参考值做差之后的插值在0值附近90％的因子是准确的，并将累计概率密度为90％所对应的上边界作为传统阈值的增量。因此需要研究在无电离层闪烁发生时段，所构建闪烁因子与参考值之间差值的分布规律；利用偏斜分析，确定差值的可能符合的分布函数：若不存在偏斜，则符合高斯分布函数，若存在偏斜，则符合指数高斯分布函数。利用每三小时内的差值，确定分布函数中的参数，并据此确定90％概率密度函数的上边界所对应幅值，即为阈值增量，在确定阈值增量δ的基础上，所构建闪烁因子的阈值thres的数学表达式如下：

thres＝0.2+δ

其中0.2为利用1Hz或电离层闪烁接收机直接给出的闪烁因子的经验阈值，单位为弧度。

由于确定的阈值不是一个定值，会随测站、太阳活动、地磁的变化而不同，为了在无法提供参考值时，仍可利用所构建闪烁因子进行电离层闪烁监测，需要对阈值进行建模预报；以表征太阳活动的物理参数F10.7、太阳黑子数量和表征地磁活动的Kp、Dst参数以及测站的坐标、地方时作为自变量，以阈值作为因变量，利用多元非线性回归的方法，构建阈值的模型，据此可实现对阈值的预报。