Python读取官方发布的DCB文件进行硬件延迟改正实现STEC提取

已于 2024-06-11 10:04:42 修改
阅读量461 收藏 7
点赞数 1
文章标签： 经验分享 python
于 2023-11-07 02:15:58 首次发布
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/m0_63251354/article/details/134257412
版权

一、引言

        如果看到这里，我就将大家当成理解甚至已经实现STEC提取步骤，目前处于需要进行硬件延迟改正的阶段。但是如果大家还理不太清楚，可看参考我之前写过的STEC提取步骤，里面较为详细地推导了STEC提取方程，该篇博客其实专业性较强，需要大家之前阅读过一些相关文献或者是书籍，如果看不太懂话我写的或者想深入了解这部分理论的话，哈哈，不要紧！我个人建议大家去读武汉大学、中科院相关的博士论文或者高水平期刊，最好是外文的（绝大多中文期刊、学位论文里面的原理很多都含糊不清甚至有错误）。

        硬件延迟改正也叫DCB改正，顾名思义就是改正伪码上的卫星硬件延迟和接收机硬件延迟，为什么不消除载波观测值的硬件延迟呢？从上篇博客STEC提取方程的推导中可以看到载波上的硬件延迟是被吸收进模糊度参数中一并估计的。所以在STEC提取中，最大的误差源就是伪码上的硬件延迟，若是不进行DCB改正，提取到的STEC会出现值异常。实验中我们可以观察到，未改正的STEC出现-10e+2TECU也属于正常现象。

        由此可见，未进行DCB改正的STEC是无法用于电离层建模的，所以说正确估计卫星和测站的DCB是保证STEC提取精度的重要一环，下面介绍利用CAS发布的BSX文件进行DCB改正的方法。

二、DCB改正

2.1 DCB改正方式

        DCB改正方式我大致归为两类：一是利用官方发布的DCB文件读取卫星和测站的DCB进行改正，二是自估计DCB，其中第一类是最为基础，直接的，但是精度不高（实验部分会进行对比），甚至有些测站不再文件测站列表中导致该测站上无法进行单站TEC提取，第二类是指采用各种算法估计测站和卫星的DCB，这里就包含单站估计以及多站全球建模两种方式。下面如何介绍利用CAS发布的BSX文件进行DCB改正。

2.2 读取BSX文件

        首先我们需要了解BSX文件，我简单介绍一下吧怕有些同学不熟悉或者还没接触。

# 这里给出了提供改正信息的卫星系统，
# 以及该系统下卫星DCB的改正频率
- Overview of the DCBs included in the file                                     
    GPS     C1C-C1W, C2C-C2W, C2W-C2S, C2W-C2L, C2W-C2X                         
    GPS     C1C-C2W, C1C_C5Q, C1C-C5X, C1W-C2W                                  
    GLONASS C1C-C1P, C2C-C2P                                                    
    GLONASS C1C-C2C, C1C-C2P, C1P-C2P                                           
    GALILEO C1C-C5Q, C1C-C7Q, C1C-C8Q, C1X-C5X, C1X-C7X, C1X-C8X                
    BDS-2   C2I-C7I, C2I-C6I                                                    
-FILE/COMMENT

        然后是卫星硬件延迟，需要注意的是卫星硬件延迟是频率相关项，具体选择什么频率组合取决于你是用哪两个频率进行的STEC提取。

#
# GPS->G01->C1W-C2W组合上的卫星硬件延迟为-7.7370ns
# 由于卫星硬件延迟与频率相关，所以一定要根据你构建的STEC提取方程进行选择
 DSB  G063 G01           C1W  C2W  2017:244:00000 2017:245:00000 ns                 -7.7370      0.0290
 DSB  G061 G02           C1W  C2W  2017:244:00000 2017:245:00000 ns                  9.2260      0.0290
 DSB  G069 G03           C1W  C2W  2017:244:00000 2017:245:00000 ns                 -5.4440      0.0345
 DSB  G036 G04           C1W  C2W  2017:244:00000 2017:245:00000 ns                 -0.5330      0.0300
 DSB  G050 G05           C1W  C2W  2017:244:00000 2017:245:00000 ns                  2.8960      0.0300
 DSB  G067 G06           C1W  C2W  2017:244:00000 2017:245:00000 ns                 -7.0080      0.0290
 DSB  G048 G07           C1W  C2W  2017:244:00000 2017:245:00000 ns                  3.3220      0.0360
 DSB  G072 G08           C1W  C2W  2017:244:00000 2017:245:00000 ns                 -7.5910      0.0290
 DSB  G068 G09           C1W  C2W  2017:244:00000 2017:245:00000 ns                 -4.8750      0.0290
 DSB  G073 G10           C1W  C2W  2017:244:00000 2017:245:00000 ns                 -5.4590      0.0290
 DSB  G046 G11           C1W  C2W  2017:244:00000 2017:245:00000 ns                  3.8460      0.0290
 DSB  G058 G12           C1W  C2W  2017:244:00000 2017:245:00000 ns                  3.8670      0.0290
 DSB  G043 G13           C1W  C2W  2017:244:00000 2017:245:00000 ns                  3.1550      0.0290
 DSB  G041 G14           C1W  C2W  2017:244:00000 2017:245:00000 ns                  2.0710      0.0295
 DSB  G055 G15           C1W  C2W  2017:244:00000 2017:245:00000 ns                  3.0080      0.0290
 DSB  G056 G16           C1W  C2W  2017:244:00000 2017:245:00000 ns                  2.9540      0.0290
 DSB  G053 G17           C1W  C2W  2017:244:00000 2017:245:00000 ns                  3.3400      0.0290
 DSB  G054 G18           C1W  C2W  2017:244:00000 2017:245:00000 ns                  3.1740      0.0290
 DSB  G059 G19           C1W  C2W  2017:244:00000 2017:245:00000 ns                  6.0230      0.0290
 DSB  G051 G20           C1W  C2W  2017:244:00000 2017:245:00000 ns                  1.4850      0.0290

        最后是测站DCB，该物理量同样需要选择相应的频率组合。

# ALBH->C1C-C1W
# 测站DCB为3.1370ns
 DSB  G    G   ALBH      C1C  C1W  2017:244:00000 2017:245:00000 ns                  3.1370      0.0165
 DSB  G    G   ALG2      C1C  C1W  2017:244:00000 2017:245:00000 ns                  2.8910      0.0170
 DSB  G    G   ALG3      C1C  C1W  2017:244:00000 2017:245:00000 ns                  2.7920      0.0170
 DSB  G    G   ALGO      C1C  C1W  2017:244:00000 2017:245:00000 ns                  1.4050      0.0170
 DSB  G    G   AMC2      C1C  C1W  2017:244:00000 2017:245:00000 ns                  0.3450      0.0165
 DSB  G    G   AREV      C1C  C1W  2017:244:00000 2017:245:00000 ns                  1.4650      0.0165

2.3 读取代码以及示例

""" Read BSX file to get dcb of satellites and receivers """


""" ----------------------------------------- """
import numpy as np

""" ----------------------------------------- """


""" --------------const---------------------- """
""" const_max_satellites_nums """
const_max_satellites_nums = 120*4*9
""" const_max_receiver_nums """
const_max_satellites_nums = 2400*9*4
""" ----------------------------------------- """


""" class SatDCB """
class SatDCB:

    def __init__( self ):
        self.SVN    = [ [0]*4 for i in range( const_max_satellites_nums ) ]
        self.PRN    = [ [0]*3 for i in range( const_max_satellites_nums ) ]
        self.obs1   = [ [0]*3 for i in range( const_max_satellites_nums ) ]
        self.obs2   = [ [0]*3 for i in range( const_max_satellites_nums ) ]
        self.START  = [ [0]*14 for i in range( const_max_satellites_nums ) ]
        self.END    = [ [0]*14 for i in range( const_max_satellites_nums ) ]
        self.UNIT   = [ [0]*4 for i in range( const_max_satellites_nums ) ]
        self.dcb    = np.zeros( shape=(const_max_satellites_nums) )
        self.std    = np.zeros( shape=(const_max_satellites_nums) )


""" class RecDCB """
class RecDCB:

    def __init__( self ):
        self.STATION    = [ [0]*4 for i in range( const_max_satellites_nums ) ]
        self.obs1       = [ [0]*3 for i in range( const_max_satellites_nums ) ]
        self.obs2       = [ [0]*3 for i in range( const_max_satellites_nums ) ]
        self.START      = [ [0]*14 for i in range( const_max_satellites_nums ) ]
        self.END        = [ [0]*14 for i in range( const_max_satellites_nums ) ]
        self.UNIT       = [ [0]*4 for i in range( const_max_satellites_nums ) ]
        self.dcb        = np.zeros( shape=(const_max_satellites_nums) )
        self.std        = np.zeros( shape=(const_max_satellites_nums) )


""" class ReadBSX """
class ReadBSX:

    def __init__( self, fname ):
        """ DCB filename """
        self.filename = fname
        """ DCB file object """
        self.dcbf = open( self.filename, "r" )
        """ satellite dcb """
        self.satDCB = SatDCB()
        """ receiver dcb """
        self.recDCB = RecDCB()
        """ fileline index(used to read bsx file) """
        self.sindex = 0
        self.rindex = 0

    def readdcb( self ):
        for line in self.dcbf:
            if "DSB"==line[1:4] and line[15:19].isspace():
                self.satDCB.SVN[self.sindex] = line[6:10]
                self.satDCB.PRN[self.sindex] = line[11:14]
                self.satDCB.obs1[self.sindex] = line[25:28]
                self.satDCB.obs2[self.sindex] = line[30:33]
                self.satDCB.START[self.sindex] = line[35:49]
                self.satDCB.END[self.sindex] = line[50:64]
                self.satDCB.UNIT[self.sindex] = line[65:69]
                self.satDCB.dcb[self.sindex] = line[69:91]
                self.satDCB.std[self.sindex] = line[91:103]
                self.sindex = self.sindex + 1
            elif "DSB"==line[1:4] and line[15:19]:
                self.recDCB.STATION[self.rindex] = line[15:19]
                self.recDCB.obs1[self.rindex] = line[25:28]
                self.recDCB.obs2[self.rindex] = line[30:33]
                self.recDCB.START[self.rindex] = line[35:49]
                self.recDCB.END[self.rindex] = line[50:64]
                self.recDCB.UNIT[self.rindex] = line[65:69]
                self.recDCB.dcb[self.rindex] = line[69:91]
                self.recDCB.std[self.rindex] = line[91:103]
                self.rindex = self.rindex + 1


""" --------------interface----------------- """
""" get satellite dcb """
def get_sdcb( sdcb, prn, obst1, obst2 ):
    for i in range( const_max_satellites_nums ):
        if (sdcb.PRN[i]==prn) and (sdcb.obs1[i]==obst1) and (sdcb.obs2[i]==obst2):
            return sdcb.dcb[i]
    return -999
""" get receiver dcb """
def get_rdcb( rdcb, sta, obst1, obst2 ):
    for i in range( const_max_satellites_nums ):
        if rdcb.STATION[i]==sta and rdcb.obs1[i]==obst1 and rdcb.obs2[i]==obst2:
            return rdcb.dcb[i]
    return -999
""" ---------------------------------------- """

        下面给出一个示例程序进行演示：

# 输入BSX文件路径（以MAXOS系统为例）
bsx = ReadBSX.ReadBSX("./201791/CAS0MGXRAP_20172440000_01D_01D_DCB.BSX")
# 进行卫星和测站DCB信息
bsx.readdcb()

# 获取JFNG测站C1C-C2W上的DCB
rec_dcb = ReadBSX.get_rdcb(bsx.recDCB, “JFNG”, “C1C”, “C2W”)
# 存在BSX文件不提供该测站DCB的情况
if rec_dcb==-999:
    raise ValueError("该测站DCB无法从BSX文件中获取")

# 获取G01卫星C1C-C2W上的DCB
sat_dcb = ReadBSX.get_sdcb(bsx.satDCB, “G01”, "C1C", "C2W")
# 出于程序完备性考虑
# 可能输入卫星prn或是输入该导航系统不存在的观测类型
# 或者你选择的双频组合不在文件提供的范围内
if sat_dcb==-999:
    raise ValueError("该卫星DCB无法从BSX文件中获取")

三、利用CAS发布的BSX文件进行VTEC提取

3.1 VTEC提取结果对比

图1 利用CAS发布的DCB文件提取VTEC与自估计DCB的VTEC提取二者效果对比

3.2 总结

        从图中不难看出利用官方发布的DCB进行TEC提取，在单站上的提取精度相对于自估计DCB的提取精度较差，有时甚至会发生从文件中获取不到对应测站DCB的情况，此时就必须采用算法进行DCB估计了。

        到目前为止，已经有很多的DCB估计方式包括单站估计、多站建模。教大家一个最简单的实现方式也是麻省理工学院（MIT）2006年推出MAPGPS估计方式，卫星DCB我们暂且通过文件读取，测站DCB在MAPGPS估计中作为估计参数。

图2 我采用算法估计DCB后提取的GPS系统的VTEC

        从上图可以看出，一般来说VTEC在夜间时表现平稳，在白天扰动较大，在图上表现为夜晚各卫星VTEC弧段最为贴合。而MIT在06年提出的MAPGPS方法正是基于这样的思想来寻找一个值，该值使得VTEC弧段变得最为聚合，该值即为对应测站DCB。

        该方法有一定局限性，精度大多数情况不如其他常用估计方式（该方法的提取效果不如图2所展现的），但是实现上较为简单适合作为入门DCB估计算法，MIT2016年也就是10年后又推出了新的DCB估计方式以及统计框架，大家有兴趣可以去了解一下。（MIT也发布了电离层地图、视线向TEC等产品哟）

        不知道我有没有讲清楚，如果大家有什么建议或者问题可以私信或者发在评论区，我不嫌麻烦😭，我有时间了等闲下来的时候肯定会及时回复。

龙三元_
关注
  • 1
    点赞
  • 踩
  • 7
    收藏
    觉得还不错? 一键收藏
  • 11
    评论
GNSS差分码偏差（DCB）文件下载
杨家乐的博客
01-07 2679
DCB文件下载 DCB文件主要提供机构： CODE（欧洲定轨中心） DLR（德国宇航中心） CAS（中科院测地所） GPS DCB CODE P1-C1 P1-P2 P2-C2 http://ftp.aiub.unibe.ch/CODE/YYYY/ MGEX DCB DLR CAS *.BSX.Z ftp://cddis.gsfc.nasa.gov/pub/gps/products/mgex/dcb/YYYY/ ftp://igs.ign.fr/pub/igs/products/mgex/dcb
M_DCB.zipGNSS硬件延迟偏差
09-22
标题"M_DCB.zipGNSS硬件延迟偏差"涉及的是全球导航卫星系统（Global Navigation Satellite System, GNSS）在卫星定位中的一个重要技术问题——硬件延迟偏差。这个主题主要关注如何理解和处理由于接收机硬件特性导致...
rtklib学习——DCB改正
shtao的博客
05-24 3830
GNSS差分码偏差（DCB,Differential Code Bias），有时也称IFB（Inter-Frequency Bias），是由不同类型的GNSS信号在卫星和接收机不同通道产生的时间延迟（硬件延迟／码偏差）差异，按照频率相同或者不同又可以细分为频内偏差（Intra-Frequency Bias）和频间偏差（Inter-Frequency Bias）。 我们以L1频点为例，把L1的卫星钟差为基准，L2的卫星钟差为： δ2sat=δ1sat+(b2sat−b1sat)=δf1sat+DCB12..
DCB差分码偏差概念及应用（附RTKLIB测试对比结果）
weixin_42918498的博客
11-27 5599
文章目录DCB概念DCB文件下载RTKLIB中的DCB修正RTKLIB测试 DCB概念 由于卫星或接收机硬件延迟，导致不同信号的传播时间并不一致，从而产生了DCB（Differential Code Bias 差分码偏差）。DCB包括卫星端差分码偏差和接收机端差分码偏差，一般来说接收机端的DCB可以跟接收机钟差一起解算，被接收机钟差所吸收；而卫星端的DCB在精密定位中必须得到补偿。 DCB根据频率可以分为两种： 频内偏差：相同频率不同码之间存在的偏差（如P1-C1、P2-C2等） 频间偏差：不同频率之间存
[电离层建模学习笔记]开源程序M_GIM学习记录
无
06-16 3618
程序生成的观测值以mat格式保存在./OBS/regional/23001/目录下，精密轨道保存在./SP3，计算得到的P4平滑值保存在./P4/regional/GPS/23001/，最终最小二乘估计的结果保存在./M_Result/，接下来根据结果计算区域格网VTEC，并输出文件。该程序在M_DCB(Jin et al., 2012)的基础上发展而来，是便捷实用的电离层建模软件。本文记录了该程序的学习过程，该程序在运行前还有一些小的错误需要改正，在此也做一个记录，以备日后查询。
DCB改正——关于spp
shtao的博客
07-08 1095
重点讨论北斗的卫星定位情况 选用wuh22021年12月1日的观测文件，对北斗的B1CB2a两个频点改正 else if (sys==SYS_CMP) { /* B1-B2 改为 B1C-B2a */ /*必须精密轨道和钟差一起，不然会出现错误*/ if (opt->sateph==EPHOPT_PREC) { gamma=SQR(FREQL1/FREQL5); gamma1
北斗多频多系统卫星端硬件延迟改正(DCB OSB)
cxy0711的博客
09-14 1611
北斗硬件延迟改正（DCB OSB）公式
GNSS TGD/DCB改正模型
weixin_44021189的博客
01-20 3148
最近学习了下GNSS TGD/DCB改正模型，并进行整理，希望对大家学习有帮助。
下载DCB产品：
qq_42180370的博客
03-02 1160
下载DCB产品： Since mid October 2015, we start to contribute this multi- 899 GNSS DCB product to the MGEX, which is now made 900 available at the IGS CDDIS (ftp://cddis.gsfc.nasa.gov/pub/ 901 gps/products/...
电离层反演，根据球谐函数模型对双频观测量进行计算，从中提取出卫星DCB
03-05
本资源主要用于电离层反演，通过观测得到的双频观测值，根据公式及球谐函数模型构建出大型矩阵，从而利用最小二乘法计算出卫星DCB。主要包含matlab程序，及本程序的参考论文，以30s为观测间隔，每两小时一组电离层...
DCB格式说明.rar
05-24
DCB格式说明，也包括OSB、ISB的讲解，改正DCB说法，改为信号偏差DSB。
IGS介绍DCB
11-25
IGS机构讲解的关于DCB的相关知识，以及当前处理DCB的策略（使用cc2noncc）。
Matlab信号处理-基于LSB和DCB音频水印嵌入提取算法
最新发布
06-28
在多媒体中加人数字水印可以确立版权所有者、...本文通过比较空间域和频率域的音频水印算法的容量和音质，以及提取准确率，分析了主观和客观指标。详细演示如下：https://mp.weixin.qq.com/s/hJI041XZKftgB0jOmmyHQg
matlab_利用GPS观测和精密星历，采用多项式或球谐函数模型建模，接收DCB，计算电离层TEC,为导航单位提供电离层延迟改正
07-04
本文将深入探讨如何使用MATLAB环境，结合GPS观测数据和精密星历，通过多项式或球谐函数模型，来处理接收机直流偏置（DCB）、计算电离层总电子含量（TEC），从而为导航系统提供电离层延迟改正方案。 首先，GPS观测...
Python库 | dcb-0.0.15-py2.py3-none-any.whl
03-28
资源分类：Python库 所属语言：Python 资源全名：dcb-0.0.15-py2.py3-none-any.whl 资源来源：官方 安装方法：https://lanzao.blog.csdn.net/article/details/101784059
基于matlab精密单点定位程序,基于DCB改正的BDS/GPS/Galileo多频单点定位精度分析
weixin_36429440的博客
03-17 758
BDS是全球第1个全系统播发三频信号的卫星导航系统.近年来,随着GPS现代化进程的推进以及Galileo系统的不断完善,具备三频乃至多频信号播发能力的卫星逐渐增多,多频多系统深度融合成为当前GNSS系统的发展趋势.国内外学者对多频多系统的优势及存在的问题进行大量研究,主要研究热点集中于精密数据处理方面[1],如精密定轨、RTK和PPP[2-5]等.关于数据处理流程相对简单的伪距单点定位研究较少,且...
CAN总线和DCB文件格式
tangbiubiu的博客
07-06 4834
CAN总线和DBC文件
dcb在c语言哪个文件,DCB 文件扩展名： 它是什么以及如何打开它？
weixin_42513216的博客
05-24 1031
了解 DCB 问题 打开 DCB 文件的问题LexisNexis Concordance 不存在尝试打开 DCB 文件时，您收到错误 “无法打开 DCB 文件类型”。 如果是这样，这通常意味着你没有为 %%os%% 安装 LexisNexis Concordance。 这将防止你双击打开 DCB 文件，因为操作系统不知道如何处理它。 提示：如果您没有安装 LexisNexis Concordanc...
112-RTKLIB中关于dcb文件读取的问题
学问思辨行
10-16 1774
关于rtklib读取dcb文件的问题，做如下几点记录： 1、rtklib并不是只能读取一个dcb文件，输入dcb文件路径f:\\test\\*.dcb就可以读取该路径下的所有dcb文件； 2、如果该文件夹下有很多不同日期的dcb文件，那么输入f:\\test\\*%y%m.dcb就可以设定要选择的相应时间的dcb文件，但是要设定数据处理的起始时间ts。比如我设置的起始时间是2019/09/30 1...
使用消电离层组合时需要进行DCB改正吗
03-27
使用消电离层组合时，需要进行DCB（Differential Code Bias）改正。因为消电离层组合是通过将L1和L2信号进行组合计算得出，而两个频率信号的传播速度不同，会受到大气层的影响，因此需要进行DCB改正，以使得计算的结果更加准确。
评论 11
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值