一种多分辨率的电离层层析方法与流程

本发明属于电离层探测技术领域，尤其是涉及一种多分辨率的电离层层析方法。

背景技术：

电离层层析作为一种全天候、大范围的电离层探测技术，具有费用低、操作简单、探测范围广等诸多优势，对于电离层不同尺度结构变化及全球电离层环境监测具有重要意义。但是由于陆基gnss数据覆盖的局限性以及地面站分布不均匀，导致层析反演数据源缺乏且分布不均匀，从而限制了电离层结构的重构精度，尤其是对中小尺度结构的重构产生严重制约。

近年来，随着一些地区陆基gnss地面站的增加和多星座gnss的建设，为电离层成像算法的发展提供了良好的契机。赵海山等人针对传统乘法代数重构算法(mart)迭代精度不高的问题,提出一种自适应电离层层析新算法，一方面，该算法根据射线穿越像素点的截距和电子密度值的综合影响,合理地分配迭代差值；另一方面,提出一种与电子密度值相关的自适应松弛因子,有效克服传播噪声对电子密度反演的影响。汤俊等人针对电离层电子密度重构问题,提出了一种综合利用总变差最小化与乘法代数重构算法的电离层层析成像算法。该算法对反演模型的参数施加总变差约束,以提高反演过程的稳定性和结果的精确性。明峰等人提出一种附加尺度因子的电离层层析方法，将实测stec与迭代中当前电子密度计算的stec的比值作为基本尺度,依据电子密度的空间相关性和网格中心到观测射线的距离,对没有信号穿过的网格也计算出对应尺度进行调整,使其向真值逐渐逼近,提高反演精度。以上这些方法在一定程度上提高了电离层层析反演的精度，但这些层析方法是基于一次反演来重构电离层结构，导致电离层层析反演精度不能达到一个较理想的状态。

技术实现要素：

有鉴于此，本发明旨在提出一种多分辨率的电离层层析方法，以实现多分辨率的电离层结构层析成像。

为达到上述目的，本发明的技术方案是这样实现的：

一种多分辨率的电离层层析方法，包括以下步骤：

(1)将反演区域的电离层格网化，将每一个格网分为三维的体积像素，穿过格网的信号路径上电子密度积分即为tec，并设置反演的时间函数；

(2)对gnss观测数据预处理，计算步骤(1)所需的tec；

(3)对步骤(1)中穿过格网的信号路径创建经验正交函数；

(4)对步骤(3)得到的经验正交函数构建补偿映射矩阵m，并对矩阵做正则化处理；

(5)对步骤(4)得到正则化处理后的矩阵和步骤(2)计算得到的tec构建正规方程，并估计方程的解，得到方程最优解，即一次反演结果；

(6)设定level的值，并再次进行步骤(3)、(4)的操作；

(7)根据步骤(5)得到的一次反演结果作为输入数据和步骤(6)得到处理后的矩阵再次构建正规方程，估计方程的解，得到二次反演的最优解。

进一步的，所述步骤(1)具体包括：第一次反演用setgrid函数将格网设置为较粗略的4°×4°，并将时间间隔设置为30min。

进一步的，所述步骤(2)具体包括：对观测数据进行周跳和整周模糊度检测，然后利用载波相位观测法计算穿过电离层格网的信号路径上的tec。

进一步的，所述步骤(3)具体包括：利用chapman函数或国际参考模型(iri)从一组电离层电子密度的径向分布来创建经验正交函数(eof)。

进一步的，所述步骤(4)具体包括：在反演过程中，电子密度被假定为给定的一组三维eof的线性组合，即映射矩阵使用一组eof来映射电子密度的径向变化，经度和纬度映射电子密度的水平方向，反演的时间窗映射电子密度的时间轴，然后分别在水平方向、径向、时间上对利用映射矩阵求出的电子密度进行二次平滑处理。

进一步的，所述步骤(5)具体包括：使用最小残差法将最小范数残差解引入到正规方程求取最优解或者二次规划法求方程最优解。

进一步的，所述步骤(6)具体包括：设置level＝[2,1]，level代表空间和时间分辨率，二次反演是在第一次反演基础上将空间和时间分辨率以2的次幂倍降低，即在第二次反演中电离层格网是1°×1°，时间为15min。

进一步的，所述步骤(7)具体包括：将一次反演结果内插到由level设定的反演格网及时间函数里，然后构建正规方程，使用最小残差法或二次规划法求解此次正规方程的最优解。

本发明的另一目的在于提出一种多分辨率的电离层层析装置，具体方案是这样实现的：

一种多分辨率的电离层层析装置，包括：

用于将反演区域的电离层格网化，将每一个格网分为三维的体积像素，穿过格网的信号路径上电子密度积分即为tec，并设置反演的时间函数的电离层网格化模块；

用于对gnss观测数据预处理，计算电离层网格化模块所需的tec的tec计算模块；

用于对电离层网格化模块中穿过格网的信号路径创建经验正交函数的经验正交函数创建模块；

用于对经验正交函数创建模块得到的经验正交函数构建补偿映射矩阵m，并对矩阵做正则化处理的经验正交函数修正模块；

用于对经验正交函数修正模块得到正则化处理后的矩阵和tec计算模块计算得到的tec构建正规方程，并估计方程的解，得到方程最优解，即一次反演结果的一次反演结果模块；

用于设定level的值，并再次利用进行经验正交函数创建模块、经验正交函数修正模块进行经验正交函数创建、经验正交函数修正操作的数据处理模块；

用于根据一次反演结果模块得到的一次反演结果作为输入数据和数据处理模块得到处理后的矩阵再次构建正规方程，估计方程的解，得到二次反演的最优解的二次反演结果模块。

相对于现有技术，本发明所述的一种多分辨率的电离层层析方法及系统具有以下优势：

(1)本发明对电离层结构进行二次重构，较大程度地提高了电离层层析成像的精度；

(2)本发明通过对level的设定，可以实现不同较小分辨率的电离层层析结果。

附图说明

构成本发明的一部分的附图用来提供对本发明的进一步理解，本发明的示意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定。在附图中：

图1为本发明提供的多分辨率电离层层析技术方法的流程图；

图2为电离层格网化结构示意图；

图3(a)为mhj45站点2015年3月16日4°×4°的no-multi和multi与ionos对比图；

图3(b)为mhj45站点2015年3月16日1°×1°的mhj45no-multi和multi与ionos对比图；

图4(a)为if843站点2015年3月16日4°×4°的no-multi与multi和ionos对比图；

图4(b)为if843站点2015年3月16日1°×1°的no-multi与multi和ionos对比图。

具体实施方式

需要说明的是，在不冲突的情况下，本发明中的实施例及实施例中的特征可以相互组合。

下面将参考附图并结合实施例来详细说明本发明。

如图1所示，本发明提供的多分辨率电离层层析方法按顺序进行的下列步骤：

(1)将反演区域的电离层格网化，将每一个格网分为三维的体积像素，穿过格网的信号路径上电子密度积分即为tec，并设置反演的时间函数；

使用setgrid函数将电离层格网化，然后设卫星与接收机传输路径为h，该射线穿过每一体素时在垂直方向的截距记为h1,h2,...,hi...假设每个网格中的电子密度为定值x1,x2,...,xi...则单位格网内tec值为hixi，则一条传输路径上的tec值可认为是射线穿过网格内的截距与该网格内的电子密度的乘积之和，即tec＝h1x1+h2x2+…+hixi+…；则所有射线上电子含量为hx＝z，x即为所求的电离层电子密度。并将一天的反演时间段设置成30min的时间间隔，从而进行4°×4°的第一次反演。

(2)对gnss观测数据预处理，计算步骤(1)所需的tec；

具体步骤如下：

gps接收机的rinex文件包含以下可观测量。

其中p1和p2是从精密码中获得的伪距，l1和l2代表信号的载波相位。p0是电离层自由伪距，n和λ分别代表整周模糊度和载波波长。ε代表接收机和卫星硬件误差分量。

比较相位和码距的可观测量之间的不同，可以得到关于tec的两个方程

综上所述，tec值可通过式2.5中的噪声项计算，也可通过式2.6中与其相关联的整数模糊度的偏置项来计算。卫星在可见范围内(除了在周跳期间，如较大的跳变或突变时)，其整周模糊度恒定，发生大周跳时的偏置量等于式2.5和式2.6求差的加权平均，并且该权重与信噪比相关。相比利用伪距计算tec值，载波相位计算的tec精度更高，，因此通常采用载波相位法来计算穿过电离层格网的信号路径上的tec。

(3)解决步骤(1)中穿过格网的路径不完备问题，创建经验正交函数；

(a)chapman函数创建经验正交函数(eof)：

电子密度分布chapman函数一般方程为：

其中χ是太阳天顶角；n是电子密度；k是一个有值的常数，表示chapmanβ层时，k＝1，表示chapmanα层时，z由下式给出：

其中h和hm分别代表高度和峰值密度高；h是标高。我们通常是通过chapman的β函数构建正交函数，一般需要以下参数：白天260km～320km范围内的峰值高度，晚上260km～360km范围内的峰值高度，以及60km～90km的标高。此外，还需输入一定经纬度范围内的电离层格网即setgrid函数和level设定的格网，以获取电子密度反演重构算法所需的三维空间经验正交函数集。

(b)iri模型创建经验正交函数(eof)：

国际参考电离层(iri)是地球电离层等离子体参数的一个标准规范。iri模型综合了全球范围的电离层探测声纳、散射雷达等众多电离层探测仪器的观测数据,描述了50km～1500km范围内电离层密度和温度的月平均变化。用iri模型系数创建经验正交函数，一般需要以下参数：所需经验正交函数个数，日期以及指定的电离层空间格网。具体步骤如下：

1、通过iri模型获取测区中心位置且与观测数据同时期的电离层剖面时间序列；

2、求取这个时间序列的协方差矩阵；

3、对此协方差矩阵进行特征值分解从而求取相应的eof。

(4)对步骤(2)得到经验正交函数构建补偿映射矩阵m，并对矩阵做正则化处理。

在反演过程中，电子密度被假定为经验正交函数(eof)的线性组合，即映射矩阵使用步骤(3)计算的经验正交函数(eof)来映射电子密度的径向变化，经度和纬度映射电子密度的水平方向，反演的时间窗映射电子密度的时间轴，其中经纬度间隔及时间窗分别按照步骤(1)设置，形成映射矩阵m，即方程hx＝z变为(hm)x＝z，然后引入正则化矩阵r，在水平、径向、时间上对方程的解进行二次平滑。

(5)对步骤(4)得到处理后的矩阵和步骤(2)计算得到的tec构建正规方程，并估计方程的解，得到方程最优解，即一次反演结果。

具体步骤：

(a)构建正规方程h^th+cr^tr＝h^tz，c是用户定义的正则化常数

(b)利用最小残差法计算方程的解

最小残差方法(minres)试图将最小范数残差解x引入到线性方程组，例如：ax＝b，在最小残差法里要求n阶系数矩阵a必须对称但不必正定，等式右侧向量b的行或列的长度必须为n。

(c)利用二次规划法(quadprog)计算方程的解

二次规划法(quadprog)求解问题的形式：

其中，hx≤b

这样可以求得使上述问题最小化的向量x，为使该问题有一个有限的最小值解，需要保证h为正定矩阵。quadprog算法需要两个或两个以上的输入，例如上限和下限，等式和不等式约束等，以对反演问题的解施加约束。为了使变换基中的逆问题求解为非负，quadprog算法通过不等式mx≥0来对反演问题施加约束，这保证了eof系数的向量为正数集。然后，最后的解是基于原始基，最后求得的电子密度剖面分布为：x＝mx

选择这两种方法中任意一种求出方程最优解，即为一次反演结果。

(6)设定level的值，并再次进行步骤(3)、(4)的操作；

设置level＝[2,1]，即反演区域的电离层格网变为4°×4°/2^2，反演时间间隔为30min/2^1，创建经验正交函数及构建映射矩阵，并对矩阵正则化处理同步骤(3)、(4)。

(7)根据步骤(5)得到的一次反演结果作为输入数据和步骤(6)得到处理后的正则矩阵再次构建正规方程，估计方程的解，得到二次反演的最优解。

再次使用步骤(5)的方法构造正规方程，然后根据最小残差法或二次规划法求解方程最优解，得到二次反演结果。

本发明提供的多分辨率层析方法的效果可以通过以下实验结果进一步说明，并对本发明方法计算出的结果和测高仪数据进行比较分析。

实验数据描述：gps观测数据由unavco机构提供，精密星历数据由igs官网下载，测高仪数据来自digitalionogramdatabase官网，其中一个测高仪位于美国西部的idahonationallab(代码:if843，位置：北纬43.81°，西经112.68°)，其当地时lt＝世界时ut-7；另一个位于美国东部的millstonehill(代码:mhj45，位置：北纬42.6°，西经71.5°)，其当地时lt＝ut-5。

表一给出了本发明方法和单分辨率方法与测高仪数据对比结果，从统计结果来看，本发明方法计算结果相对于测高仪数据来说精度较高；另外本发明方法的层析成像分辨率较高。

表一multi，no-multi层析结果和测高仪的nmf2误差比较

以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。