奇异点保留滤波降噪的周跳探测方法

摘要

本发明公开了一种奇异点保留滤波降噪的周跳探测方法，包括如下步骤：Step1、建立MW组合观测序列；Step2、初始化；Step3、剔除MW组合观测序列中的野值；Step4：对序列进行一阶差分；Step5：对差分序列进行极值中值滤波；Step6：滑动平均降噪，得到历元i处的估计值和其方差的估计值；Step7：根据周跳判定准则判定是否周跳，并记录周跳历元和跳变量；Step8：滑动平均窗，直到滑动到最后一历元。本发明可以实现对MW组合降噪的同时保持周跳信号不被滤除，从而达到提高探测精度的效果，解决MW组合由于噪声过大导致的一到二周小周跳无法精准探测的问题。

权利要求

1.一种奇异点保留滤波降噪的周跳探测方法，其特征在于，包括如下步骤：Step1：建立MW组合观测序列  ，为第i历元MW组合观测值；

Step2：初始化，选取合适的滑动平均窗长L和极值中值滤波窗长l，并满足关系；

Step3、通过中值滤波器剔除MW组合观测序列中的野值；

Step4：对序列：

进行一阶差分得到：

其中 ， ，k是整数；

Step5：对Step4中得到的两个一阶差分序列分别进行窗长为l的极值中值滤波，滤波后得到序列：其中 为极值中值滤波后的 序列；

Step6：计算窗长内各历元 的值：，

其中，是整数；

得到历元i处MW组合量的估计值 和其方差的估计值 ：；

Step7：若 ，则计算 和 ：，

其中 为：

，

并且，若 ，其中const1是一个常数，则判定历元i处存在周跳，并记录周跳历元i和跳变量 ；

Step8：若i后存在 个观测值，则令 并返回Step4，直到滑动到最后一历元，否则结束循环。

2.根据权利要求1所述的奇异点保留滤波降噪的周跳探测方法，其特征在于：在步骤Step1中，由RINEX观测文件读取获得伪距与载波相位数据，做线性组合得到MW组合观测序列 。

3.根据权利要求2所述的奇异点保留滤波降噪的周跳探测方法，其特征在于：在步骤Step3中，通过一个窗长为3的中值滤波器剔除MW组合观测序列中的野值。

4.根据权利要求1或2或3所述的奇异点保留滤波降噪的周跳探测方法，其特征在于：在步骤Step5中，所述极值中值滤波过程如下：若 是长度为T的待滤波序列，为极值中值滤波窗长，令 ，则：Step5.1：搜索出集合 ：中的最大值 、最小值 和中值 ，并令，

其中 是一个给定的阈值；

Step5.2：按照公式：更新 ；

Step5.3：若 ，令 返回Step5.1。

5.根据权利要求4所述的奇异点保留滤波降噪的周跳探测方法，其特征在于：在Step5.2中，更新 时，若 则取 。

说明书

奇异点保留滤波降噪的周跳探测方法

技术领域

[0001] 本发明属于GNSS观测数据处理的技术领域，涉及一种GNSS精密单点定位预处理技术，具体涉及一种奇异点保留滤波降噪的周跳探测方法。

背景技术

[0002] 目前，随着GNSS（Global Navigation Satellite System，全球导航卫星系统）在各种领域中的应用和发展，不同的用户需要不同精度的定位服务。在测绘、大气、地震监测等领域要求的服务精度越来越高。相比于伪距观测量，载波相位观测量具有精度高、噪声低等优点，被广泛的应用于高精度定位定轨服务中。但由于载波相位本身存在周跳、整周模糊度不固定等问题，使得在高精度定位定轨服务中，必须首先实现周跳的精准探测和模糊度的解算。

[0003] 周跳探测是模糊度解算的前提，只有做到精准的周跳探测才能正确的参数化模糊度，构造正确的法方程求解模糊度，提高GNSS精密定位的精度和可靠性。当GNSS接收机锁相环失锁就会产生周跳，在低信噪比以及卫星仰角过低时都有可能发生周跳。特别的，对于我国北斗卫星导航系统，GEO（Geosynchronous orbit，地球同步轨道）卫星有着卫星轨道高、卫星几何差、需要频繁进行机动调整等原因导致周跳偏多。特别的对于一些低仰角（15°-20°）的GEO卫星，更是存在频繁的一周小周跳，经典的TurboEdit方法不能准确可靠地探测该周跳。

[0004] MW组合作为周跳探测经典组合量，有着不含几何量和电离层残差误差影响的优点，但噪声偏大不利于小周跳的探测。

[0005] TurboEdit方法采用MW组合作为主要探测量，对于某一历元  采用如下公式进行估计和方差更新：

[0006]

[0007] 其中 为历元 的MW组合估计值， 为该历元方差估计量。当 历元的实际量减去历元的估计量的绝对值大于 时，首先判断其是否为野值，若否则判定为周跳并记录 历元，随后从新初始化更新 和 。实际上该更新过程就是一个滑动平均过程，一旦某历元存在小周跳未被检测出，就会导致方差更新发散，最终导致该历元以后的小周跳根本无法被探测出；另外在初始化过程中若出现周跳，将严重影响后续周跳探测的精度。所以，该方法仅适用于周跳不频繁出现的情况。同时，受伪距噪声的影响，该方法几乎无法探测一到二周的小周跳。

发明内容

[0008] 本发明的目的在于避免现有技术中的不足而提供一种奇异点保留滤波降噪的周跳探测方法，其可以实现对MW组合降噪的同时保持周跳信号不被滤除，从而达到提高探测精度的效果，解决MW组合由于噪声过大导致的一到二周小周跳无法精准探测的问题。

[0009] 本发明的目的通过以下技术方案实现：

[0010] 提供一种奇异点保留滤波降噪的周跳探测方法，包括如下步骤：

[0011] Step1、建立MW组合观测序列 ， 为第i历元MW组合观测值；

[0012] Step2、初始化，选取合适的滑动平均窗长 和极值中值滤波窗长 ，并满足关系；

[0013] Step3、通过中值滤波器剔除MW组合观测序列中的野值；

[0014] Step4：对序列：

[0015]

[0016] 进行一阶差分得到：

[0017]

[0018] 其中 ， ，是整数；

[0019] Step5：对Step4中得到的两个一阶差分序列分别进行窗长为的极值中值滤波，滤波后得到序列：

[0020]

[0021] 其中 为极值中值滤波后的 序列；

[0022] Step6：计算窗长 内各历元 的值：

[0023] ，

[0024] 其中，是整数；

[0025] 得到历元 处MW组合量的估计值 和其方差的估计值 ：

[0026] ；

[0027] Step7：若 ，则计算 和 ：

[0028] ，

[0029] 其中 为：

[0030] ，

[0031] 并且，若 ，其中const1是一个常数，则判定历元 处存在周跳，并记录周跳历元 和跳变量 ；

[0032] Step8：若 后存在 个观测值，则令 并返回Step4，直到滑动到最后一历元，否则结束循环。

[0033] 作为进一步的改进，在步骤Step1中，由RINEX观测文件读取获得伪距与载波相位数据，做线性组合得到MW组合观测序列 。

[0034] 作为进一步的改进，在步骤Step3中，通过一个窗长为3的中值滤波器剔除MW组合观测序列中的野值。

[0035] 作为进一步的改进，在步骤Step5中，所述极值中值滤波过程如下：

[0036] 若 是长度为T的待滤波序列，为极值中值滤波窗长，令 ，则：

[0037] Step5.1：搜索出集合 ：

[0038]

[0039] 中的最大值 、最小值 和中值 ，并令

[0040] ，

[0041] 其中 是一个给定的阈值；

[0042] Step5.2：按照公式：

[0043]

[0044] 更新 ；

[0045] Step5.3：若 ，令 返回Step5.1。

[0046] 作为进一步的改进，由于MW组合的一阶差分序列几乎是零均值的，故在Step5.2中，更新 时，若 则取 。

[0047] 通过本发明的方法，可以得到保留周跳的降噪序列 和其方差估计序列 ，若利用该序列和Step7给出的准则进行周跳探测就可大大提高探测精度，实现利用MW组合探测频繁小周跳的目的。本发明也可以序贯的进行降噪并完成周跳探测，在延迟 个历元后给出周跳探测结果。本发明还可以解决TurboEdit方法无法探测频繁小周跳的问题，提高MW组合探测周跳的能力，可以序贯的进行周跳探测，并能应用在实时精密单点定位数据预处理工作中，提供正确的周跳历元，方便模糊的的估算，提高定位精度。

附图说明

[0048] 利用附图对本发明作进一步说明，但附图中的实施例不构成对本发明的任何限制，对于本领域的普通技术人员，在不付出创造性劳动的前提下，还可以根据以下附图获得其它的附图。

[0049] 图1是本发明奇异点保留滤波降噪的周跳探测方法的流程图。

[0050] 图2是本发明的效果图。

具体实施方式

[0051] 为了使本领域的技术人员更好地理解本发明的技术方案，下面结合附图和具体实施例对本发明作进一步详细的描述，需要说明的是，在不冲突的情况下，本申请的实施例及实施例中的特征可以相互组合。

[0052] 如图1所示，本发明实施例提供的一种奇异点保留滤波降噪的周跳探测方法，包括如下步骤：

[0053] Step1、由RINEX（Receiver Independent Exchange Format，即与接收机无关的交换格式）观测文件读取获得伪距与载波相位数据，做线性组合得到MW组合观测序列 ，为第i历元MW组合观测值：

[0054] Step2、选取合适的滑动平均窗长 和极值中值滤波窗长 ，满足关系。

[0055] Step3、初始化，对滑窗L内的数据，通过一个窗长为3的中值滤波器剔除MW组合观测序列中的野值。按照下述Step4到Step7计算第 历元（滑动平均窗L中心历元点）点的 和 ，并令 。

[0056] Step4、对序列

[0057]

[0058] 进行一阶差分得到：

[0059] ，

[0060] 其中 。

[0061] Step5、进行上述两个一阶差分序列分别进行窗长为 的极值中值滤波，滤波后得到序列：

[0062] ，

[0063] 其中 为极值中值滤波后的 序列。

[0064] Step6、计算窗 内各历元 的值：

[0065]

[0066] 得到历元 处MW组合量的估计值 和其方差的估计值 ：

[0067] 。

[0068] Step7、若 ，则计算 和 ：

[0069] ，

[0070] 其中 为：

[0071] ，

[0072] 并且若 （const1是一个常数），则判定历元 处存在周跳，并记录周跳历元 和跳变量 。

[0073] Step8、若后存在 个观测值，则令 并返回Step4；滑窗L向后移动一个历元，直到滑动到最后一历元。

[0074] 通过上述步骤，可以得到保留周跳的降噪序列 和其方差估计序列 ，若利用该序列和Step7给出的准则进行周跳探测就可大大提高探测精度，实现利用MW组合探测频繁小周跳的目的。本发明也可以序贯的进行降噪并完成周跳探测，在延迟 个历元后给出周跳探测结果。本发明还可以解决TurboEdit方法无法探测频繁小周跳的问题，提高MW组合探测周跳的能力，可以序贯的进行周跳探测，并能应用在实时精密单点定位数据预处理工作中，提供正确的周跳历元，方便模糊的的估算，提高定位精度。

[0075] 作为进一步优选的实施方式，在步骤Step5中，所述极值中值滤波过程如下：

[0076] 若 是长度为T的待滤波序列，为极值中值滤波窗长，令 ，则：

[0077] Step5.1：搜索出集合 ：

[0078]

[0079] 中的最大值 、最小值 和中值 ，并令

[0080] ，

[0081] 其中 是一个给定的阈值；

[0082] Step5.2：按照公式：

[0083]

[0084] 更新 ；

[0085] Step5.3：若 ，令 返回Step5.1。

[0086] 作为进一步优选的实施方式，由于MW组合的一阶差分序列几乎是零均值的，故在Step5.2中，更新 时，若 则取 。

[0087] 如图2所示，图2中的上图为Jfng站C01卫星MW组合周跳探测仿真效果图，图中黄色点图为MW组合原始数据点，蓝色实线为新算法滤波降噪后的MW组合平滑曲线，红色虚线为阈值曲线，黑色三角点为新算法探测得到的发生周跳的突变点。从图中可以看出全部的15个周跳点均被探测出来，特别的对于4到5小时内加入的几个比较频繁的一周小周跳也被探测出来。图2中的下图为Jfng站C05卫星MW组合周跳探测实测图，该卫星仰角低，经常发生频繁一周小周跳。可以看出，大部分小周跳都被探测出，误探数量4漏探数量8，效果远高于TurboEdit方法探测效果。

[0088] 上面的描述中阐述了很多具体细节以便于充分理解本发明，但是，本发明还可以采用其他不同于在此描述的其他方式来实施，因此，不能理解为对本发明保护范围的限制。

[0089] 总之，本发明虽然列举了上述优选实施方式，但是应该说明，虽然本领域的技术人员可以进行各种变化和改型，除非这样的变化和改型偏离了本发明的范围，否则都应该包括在本发明的保护范围内。