卫星导航中的时间系统主要由各大卫星系统建立分为中国北斗导航系统的BDT、美国GPS系统的GPST、欧洲GALILEO系统的GST、俄罗斯GLONASS系统的GLONASST以及将各大系统统一化的协调世界时UTC时间系统所以在整个GNSS导航模组内部需要实现各大系统向UTC时间系统转换函数个别模组以GPS作为主系统主要依靠GPST来维持整个代码的时间算法逻辑。首先了解下UTC时间的定义和来源。世界时UT是以地球自转形成的时间系统其并不是一个严格均匀的时间系统由于地球自转速度变慢经过极移修正的世界时UT1在以每年1s的时间变慢。基于此引入高精度时间侧标标准国际原子时TAI其基于全球数百个高精度原子钟主要是铯原子钟的读数加权平均。原子钟的振荡频率极其稳定因此TAI是一个连续、均匀、高度稳定的时间尺度。它的秒长是固定不变的。将精密计数和真实观测结合在一起即将UT1和TAI结合在一起使用形成协调世界时UTC其起始历元为1972 年1 月1 日00 时00 分00 秒TAI当UTC和UT1差距超过±0.9s后UTC采用闰秒的方式加插1s使UTC尽量接近UT1这也是我们在GPST/BDT/GST中导航电文中会遇到的闰秒操作基础。闰秒区分正闰秒和负闰秒只不过当前发生的调整均为正闰秒闰秒的调整总是被安排在UTC时间6月30日或12月31日的的最后一分钟的最后一秒之后即在23:59:59后的时间是23:59:60该时间则为调整的闰秒也就是该1分钟包含61s。从1972年建立UTC时间起到2025年11月一共发生过27次闰秒具体时间如下表1 闰秒发生时间闰秒序号闰秒年份闰秒日期 (UTC)闰秒时分秒 (UTC)11972年6月30日23:59:6021972年12月31日23:59:6031973年12月31日23:59:6041974年12月31日23:59:6051975年12月31日23:59:6061976年12月31日23:59:6071977年12月31日23:59:6081978年12月31日23:59:6091979年12月31日23:59:60101981年6月30日23:59:60111982年6月30日23:59:60121983年6月30日23:59:60131985年6月30日23:59:60141987年12月31日23:59:60151989年12月31日23:59:60161990年12月31日23:59:60171992年6月30日23:59:60181993年6月30日23:59:60191994年6月30日23:59:60201995年12月31日23:59:60211997年6月30日23:59:60221998年12月31日23:59:60232005年12月31日23:59:60242008年12月31日23:59:60252012年6月30日23:59:60262015年6月30日23:59:60272016年12月31日23:59:60上一次闰秒发生在2016年年底到现在2025年已经8年多没有宣布调整闰秒的时间节点2022年曾有国际消息传出最晚于2035年停止引入闰秒的操作但未来具体闰秒如何处理还很未知。1. 各大时间系统介绍各大系统虽然初始时间点UTC不一致但是原理基本一致所以在此选择最熟悉的BDT做介绍和延伸。BDT 采用国际单位制SI秒为基本单位连续累计不闰秒起始历元为UTC时2006 年1 月1 日00 时00 分00 秒记为BDT0北斗系统的所有时间都是以此做为基准BDS的周计数bds_week以该时间所在周开始计数。BDT与UTC之间的闰秒信息在导航电文中播报自2006年到2015年发生过4次闰秒可以从表1中看到分别为闰秒序号2006年1月1日00:00:0012008年12月31日23:59:6022012年6月30日23:59:6032015年6月30日23:59:6042016年12月31日23:59:60即2017年1月1日后BDTUTC4s当前代码中在没有获取到BDS导航电文中的UTCmodel信息时直接采用的方式就是使用4s作为BDT和UTC的时间转换函数获取到有效的UTCmodel后则使用UTCmodel中的闰秒信息进行时间转换。1.1 BDS UTCmodel 解码BDS导航电文的UTCmodel参数主要包括1A0UTC32bitBDT相对于UTC的钟差2A1UTC24bitBDT相对于UTC的钟速3ΔtLS8bit新的闰秒生效前BDT相对于UTC的累积闰秒改正数4WNLSF8bit新的闰秒生效的周计数占8比特为DN对应的整周计数模256。WNLSF在模256之前和WN之差的绝对值不超过127。6DN8bit新的闰秒生效的周内日计数7ΔtLSF8bit新的闰秒生效后BDT相对于UTC的累积闰秒改正数。考虑到导航电文中存在误码情况并且UTCmodel直接影响时间维持所以在存储计算的时候需要判断闰秒信息是否有效1判断DN是否在0~6的范围内2判断ΔtLS和ΔtLS是否都为0并且差值是否满足≤1当前闰秒一次只会调整1s仿照RTKLIB中的UBLOX.c300bit frame存储在uint32_t bds_frame[10]里面基于BDS D1电文解码#define P2_30 9.313225746154785E-10 /* 2^-30 */ #define P2_50 8.881784197001252E-16 /* 2^-50 */ typedef struct { /* utc struct */ double a0; double a1; int8_t ls,lsf,dn,wnlsf; unsigned char check; } utc_model_t, *putc_model_t; void decode_bdsd1_utc(uint8_t *bds_frame, putc_model_t pbds_utc){ unsigned char fraid 0, pageid 0; uint16_t i 0; utc_model_t utc_temp; memset(utc_temp, 0, sizeof(utc_temp)); fraid getbitu (bds_frame,i 15, 3); pageid getbitu (bds_frame,i 43, 7); if(5 fraid 10 pageid ){ utc_temp.ls getbitu2(buff, i50, 2, i60,6); utc_temp.lsf getbitu (bds_frame, i66, 8); utc_temp.wnlsf getbitu (bds_frame, i74, 8); utc_temp.a0 getbits2(bds_frame, i90, 22, i120, 10)*P2_30 ; utc_temp.a1 getbits2(bds_frame, i130, 12, i150, 12)*P2_50; utc_temp.dn getbitu (bds_frame, i162, 8); if((0 utc_temp.dn) (utc_temp.dn 6)){ if((0 ! utc_temp.ls) (0 ! utc_temp.lsf)){ if(abs(utc_temp.ls - utc_temp.lsf) 1){ utc_temp.check 1; memcpy(pbds_utc, utc_temp, sieeof(utc_temp)); } } } } }1.2 BDS UTCmodel 使用UTCmodel的使用可以具体参见BDS ICD文件中的解释说明。假设当前时间的北斗周计数和周内秒分别为week_bds和tow_bds首先计算tow_bds在一周内的天计数dowt_bds[模86400](1) 当指示闰秒生效的周计数WNLSF和周内天计数DN还没到来时而且用户当前时刻tE处在DN2/3之前则UTC与BDT之间的变换关系为2若用户当前的系统时刻tE处在指示闰秒生效的周计数WNLSF和周内天计数DN2/3到DN5/4之间则UTC与BDT之间的变换关系为3当指示闰秒生效的周计数WNLSF和周内天计数DN已经过去且用户当前的系统时刻tE处在DN5/4之后则UTC与BDT之间的变换关系为代码实现方面需要考虑以下几个方面1和解码部分一样也需要判断重点参数的有效性2考虑bds_week%256的mod_week和wnlsf之间关系的判断3考虑到dn和跨周同时发生的情况所以在衡量DN2/3和DN5/4数值的时候需要将当前时间week和周内天综合形成bds_day来判断4考虑只有wnls的第dn天的23:59:59之后的那一秒才能闰秒形成235960所以需要给出闰秒标志。假设开发人员现在拿不到UTCmodel的相关参数这个时候需要原生GNSS模组给出闰秒信息比如上篇文章提到的UBX、和芯星通或者市面上可以遇到的其他GNSS模组均需要特定语句拿到闰秒但是有时候原生结果会晚一点时间输出所以本来闰秒为UTC235959~235960~000000就变成了UTC235959~000000~000000假设高频率输出的时候有可能也会出现闰秒时刻加载在000000中间该措施需要注意和衡量。最好的解决方法就是需要GNSS模组直接输出子帧信息算法内部解码获得UTCmodel再做时间处理。1.3 其他系统与BDT规定类似其他系统也有着时间初始零点的设置BDS导航电文中week计数有13bit最多可以计数到8192周当前时间范围内暂时可以不用考虑周计数翻转的问题。GPST的初始时间点是UTC时间1980年1月6日00 时00 分00 秒记为GPST0在GPS 导航电文中播发的week计数只有10bit其最多只能表示1024周所以在代码代码要做好GPS周翻转计数次数的存储防止在week翻转的时候出现时间维持异常的问题。GST的初始时间点是UTC时间1998年8月22日00 时00 分00 秒记为GST0Galileo导航电文中播发的week计数只有12bit可表示范围到4096周。GLONASS的时间跟随UTC时间相差3hour可以写成并且GLONSS的时间维持和其他系统都不一样可以直接参见TRTKLIB中的解码函数towtime2gpst(gpst2utc(geph-tof),week); todfmod(tow,86400.0); tow-tod; toftk_h*3600.0tk_m*60.0tk_s-10800.0; /* lt-utc */ if (toftod-43200.0) tof86400.0; else if (toftod43200.0) tof-86400.0; geph-tofutc2gpst(gpst2time(week,towtof)); toetb*900.0-10800.0; /* lt-utc */ if (toetod-43200.0) toe86400.0; else if (toetod43200.0) toe-86400.0; geph-toeutc2gpst(gpst2time(week,towtoe)); /* utc-gpst */1.4 RTKLIB 中相关时间转换1.4.1 统一时间格式在RTKLIB或者实际代码中会比较bdt1和bdt0的时间大小当week和tow都已知后会直接将其转换到以UTC0点对应的数值轴上即UTC0为0s后续分别基于该数值进行直接累加结构体形式如下//来源RTKLIB typedef struct { /* time struct */ time_t time; /* time (s) expressed by standard time_t */ double sec; /* fraction of second under 1 s */ } gtime_t; /* beidou time (bdt) to time --------------------------------------------------- * convert week and tow in beidou time (bdt) to gtime_t struct * args : int week I week number in bdt * double sec I time of week in bdt (s) * return : gtime_t struct *-----------------------------------------------------------------------------*/ extern gtime_t bdt2time(int week, double sec) { gtime_t tepoch2time(bdt0); if (sec-1E9||1E9sec) sec0.0; t.time86400*7*week(int)sec; t.secsec-(int)sec; return t; }1.4.2 各大系统时间转换由于从2016年起闰秒信息再未更新 所以大部分GPS和UTC之间的闰秒转换直接使用18sRTKLIB中有根据闰秒表和年月日判断选择闰秒的操作const static double gpst0[]{1980,1, 6,0,0,0}; /* gps time reference */ const static double gst0 []{1999,8,22,0,0,0}; /* galileo system time reference */ const static double bdt0 []{2006,1, 1,0,0,0}; /* beidou time reference */ static double leaps[MAXLEAPS1][7]{ /* leap seconds (y,m,d,h,m,s,utc-gpst) */ {2017,1,1,0,0,0,-18}, {2015,7,1,0,0,0,-17}, {2012,7,1,0,0,0,-16}, {2009,1,1,0,0,0,-15}, {2006,1,1,0,0,0,-14}, {1999,1,1,0,0,0,-13}, {1997,7,1,0,0,0,-12}, {1996,1,1,0,0,0,-11}, {1994,7,1,0,0,0,-10}, {1993,7,1,0,0,0, -9}, {1992,7,1,0,0,0, -8}, {1991,1,1,0,0,0, -7}, {1990,1,1,0,0,0, -6}, {1988,1,1,0,0,0, -5}, {1985,7,1,0,0,0, -4}, {1983,7,1,0,0,0, -3}, {1982,7,1,0,0,0, -2}, {1981,7,1,0,0,0, -1}, {0} }; extern gtime_t gpst2utc(gtime_t t) { gtime_t tu; int i; for (i0;leaps[i][0]0;i) { tutimeadd(t,leaps[i][6]); if (timediff(tu,epoch2time(leaps[i]))0.0) return tu; } return t; } /* utc to gpstime -------------------------------------------------------------- * convert utc to gpstime considering leap seconds * args : gtime_t t I time expressed in utc * return : time expressed in gpstime * notes : ignore slight time offset under 100 ns *-----------------------------------------------------------------------------*/ extern gtime_t utc2gpst(gtime_t t) { int i; for (i0;leaps[i][0]0;i) { if (timediff(t,epoch2time(leaps[i]))0.0) return timeadd(t,-leaps[i][6]); } return t; } /* gpstime to bdt -------------------------------------------------------------- * convert gpstime to bdt (beidou navigation satellite system time) * args : gtime_t t I time expressed in gpstime * return : time expressed in bdt * notes : ref [8] 3.3, 2006/1/1 00:00 BDT 2006/1/1 00:00 UTC * no leap seconds in BDT * ignore slight time offset under 100 ns *-----------------------------------------------------------------------------*/ extern gtime_t gpst2bdt(gtime_t t) { return timeadd(t,-14.0); } /* bdt to gpstime -------------------------------------------------------------- * convert bdt (beidou navigation satellite system time) to gpstime * args : gtime_t t I time expressed in bdt * return : time expressed in gpstime * notes : see gpst2bdt() *-----------------------------------------------------------------------------*/ extern gtime_t bdt2gpst(gtime_t t) { return timeadd(t,14.0); } /* galileo system time to time ------------------------------------------------- * convert week and tow in galileo system time (gst) to gtime_t struct * args : int week I week number in gst * double sec I time of week in gst (s) * return : gtime_t struct *-----------------------------------------------------------------------------*/ extern gtime_t gst2time(int week, double sec) { gtime_t tepoch2time(gst0); if (sec-1E9||1E9sec) sec0.0; t.time86400*7*week(int)sec; t.secsec-(int)sec; return t; } /* time to galileo system time ------------------------------------------------- * convert gtime_t struct to week and tow in galileo system time (gst) * args : gtime_t t I gtime_t struct * int *week IO week number in gst (NULL: no output) * return : time of week in gst (s) *-----------------------------------------------------------------------------*/ extern double time2gst(gtime_t t, int *week) { gtime_t t0epoch2time(gst0); time_t sect.time-t0.time; int w(int)(sec/(86400*7)); if (week) *weekw; return (double)(sec-w*86400*7)t.sec; }RTKLIB中将GPS作为主系统所有的时间系统都是使用GPST做转换即其他系统转换为gtime格式需要将BDT/GST先转gtime格式GPST再转为UTC时间。1.4.3 输出时间格式从上述来看基本还是各大系统的内部时间对外输出的时候一般需要转换为年月日时分秒的格式来看比如NMEA数据中的GGA/RMC语句、RINEX数据中的观测量时间可以直接使用RTKLIB中的转换函数time2epoch(。/* convert calendar day/time to time ------------------------------------------- * convert calendar day/time to gtime_t struct * args : double *ep I day/time {year,month,day,hour,min,sec} * return : gtime_t struct * notes : proper in 1970-2037 or 1970-2099 (64bit time_t) *-----------------------------------------------------------------------------*/ extern gtime_t epoch2time(const double *ep) { const int doy[]{1,32,60,91,121,152,182,213,244,274,305,335}; gtime_t time{0}; int days,sec,year(int)ep[0],mon(int)ep[1],day(int)ep[2]; if (year1970||2099year||mon1||12mon) return time; /* leap year if year%40 in 1901-2099 */ days(year-1970)*365(year-1969)/4doy[mon-1]day-2(year%40mon3?1:0); sec(int)floor(ep[5]); time.time(time_t)days*86400(int)ep[3]*3600(int)ep[4]*60sec; time.secep[5]-sec; return time; } /* time to calendar day/time --------------------------------------------------- * convert gtime_t struct to calendar day/time * args : gtime_t t I gtime_t struct * double *ep O day/time {year,month,day,hour,min,sec} * return : none * notes : proper in 1970-2037 or 1970-2099 (64bit time_t) *-----------------------------------------------------------------------------*/ extern void time2epoch(gtime_t t, double *ep) { const int mday[]{ /* # of days in a month */ 31,28,31,30,31,30,31,31,30,31,30,31,31,28,31,30,31,30,31,31,30,31,30,31, 31,29,31,30,31,30,31,31,30,31,30,31,31,28,31,30,31,30,31,31,30,31,30,31 }; int days,sec,mon,day; /* leap year if year%40 in 1901-2099 */ days(int)(t.time/86400); sec(int)(t.time-(time_t)days*86400); for (daydays%1461,mon0;mon48;mon) { if (daymday[mon]) day-mday[mon]; else break; } ep[0]1970days/1461*4mon/12; ep[1]mon%121; ep[2]day1; ep[3]sec/3600; ep[4]sec%3600/60; ep[5]sec%60t.sec; }2. 如何维持时间系统本篇文章区分自研GNSS模组以及GNSS算法/大融合算法对于时间系统的维持两者是不一样的维持思路。2.1 GNSS模组对于GNSS模组而言即就是在模组内部需要完成卫星信号的捕获跟踪、帧同步、计算观测量、维持时间系统的功能在维持时间系统模块的具体流程为1帧同步获得frame_tow基于frame_towestimate_pathtow来更新本地时间tow如果是多个卫星在同时在一个历元均帧同步成功则选择一个最高CNR的卫星做为初始化本地时间tow的参考卫星2解码星历获得周计数week不同系统一定要做区分week周计数出现的频率也要做部分。举例假设BDS系统只有2个卫星获得week判断week1和week2是否一致一致则更新rcv_week否则则不更新假设当前有3个及以上的卫星获得week后需要判断n个week中出现频率最多的week并且出现频率1/2的week才用来更新rcv_week;3根据各大系统之间的week计数关系做相互校验比如GPS_week和BDS_week相差1356GPS_week和Galileo_week相差1024周因为BDS_week的数据有效位比较长所以常用来转换计算GPS的周翻转次数。注意点1跨周的维持个别模组直接采用的方式是tow采用累加方式一直到604800.0s将tow翻转为0week1;2接收机调钟因为计算接收机钟差一直在累加常用的方式是在钟差累积到一定数值的时候直接调整此时baseband需要对于上报观测量间隔进行调整。3观测量调整在RINEX那一节也提到在输出伪距、载波、多普勒的时候将钟差、钟漂调整量修正可以参考RTCM协议中的计算公式2.2 GNSS 算法/融合算法对于GNSS算法即直接使用外部模组输出的RTCM/UBLOX数据做PVT/RTK计算/与INS融合等算法不同GNSS模组输出的时间信息有差异在此区分下RTCM和UBLOX数据流程。对于RTCM数据其观测量数据中只有tow时间没有week时间所以无法对齐时间所以在使用RTKLIB中的rtkcov解码RTCM数据的时候需要输入时间信息在代码中常用的方法是1基于RTCM星历获得week首先基于tow和eph_toe的差值做星历有限期的判断才能使用星历的week后续判断和GNSS模组一致不同系统一定要做区分week周计数出现的频率也要做部分。举例假设BDS系统只有2个卫星获得week判断week1和week2是否一致一致则更新rcv_week否则则不更新假设当前有3个及以上的卫星获得week后需要判断n个week中出现频率最多的week并且出现频率1/2n的week才用来更新rcv_week;2基于从GNSS模组定位信息中的时间信息week、tow和从星历中拿到的week做校验判断防止跨周出现异常。对于UBLOX数据观测量中有week和tow信息直接使用所增加的步骤就是从GNSS模组定位信息中的时间信息week、tow和从星历中拿到的week做校验判断防止跨周出现异常。注意点1整点对齐功能举例个别模组在输出观测量数据的时候并不是标准的5hz在0.0、0.2、0.4、0.6、0.8的时候输出而是在0.12、0.32、0.52、0.72、0.92的时刻输出所以要做要整点观测量数据的归算2跨周功能非常重要在代码内部维持的时候需要做重点关注