基于轨道根数预测的星载导航接收机定位授时方法和装置

技术领域

本发明涉及星载导航接收机的技术领域，尤其是涉及一种基于轨道根数预测的星载导航接收机定位授时方法和装置。

背景技术

全球导航卫星系统(Global Navigation Satellite System，GNSS)已经成为人们获取位置和授时信息的基本手段。不单单是地面用户需要利用全球导航卫星系统提供授时导航定位服务，太空中其它功能的卫星同样如此，比如遥感卫星、通信卫星等。星载导航接收机可以为宿主卫星提供精确的位置、速度和时间信息，是导航领域重要的技术课题和发展方向。

中高轨道卫星也可以使用星载导航接收机进行定位授时。中高轨道卫星的星载导航接收机相较于低轨导航接收机和传统地面导航接收机面临着更大的挑战。中高轨道卫星常见的轨道有20000km轨道、36400km轨道和大椭圆轨道等。在这些高度下，星载导航接收机天线需要安装在宿主卫星对地面接收全球导航卫星发射的主瓣未被地球遮挡部分和旁瓣信号部分，以实现导航定位的目的。事实上，中高轨道卫星以接收导航卫星的旁瓣信号为主，导航卫星旁瓣信号非常弱，而且中高轨道卫星处于高速运动状态。因此，中高轨卫星的星载导航接收机比利用导航卫星主瓣信号的低轨和地面导航接收机需要更加先进的定位授时方法，以确保定位授时结果的精度和稳定性，以满足波束定向卫星通信等载荷的应用要求。

发明内容

有鉴于此，本发明的目的在于提供一种基于轨道根数预测的星载导航接收机定位授时方法和装置，充分的利用卫星在轨运动规律对导航状态进行预测，从而改善导航定位和授时结果的精度，以更好的满足中高轨波束定向通信载荷等应用对在轨定位和授时精度的要求。

第一方面，本发明实施例提供了一种基于轨道根数预测的星载导航接收机定位授时方法，应用于星载导航接收机，包括：确定历史历元时刻所述星载导航接收机所搭载卫星的历史状态量，其中，所述历史状态量包括以下至少之一：历史位置矢量、历史速度矢量、历史钟差和历史频差；基于所述搭载卫星的轨道根数和所述历史状态量，计算所述搭载卫星的目标预测量，其中，所述目标预测量包括：位置预测量和速度预测量；基于所述目标预测量确定所述搭载卫星的第一导航观测量，其中，所述第一导航观测量包括以下至少之一：第一伪距和第一载波相位；对所述搭载卫星的所述第一导航观测量进行计算，得到所述搭载卫星在当前历元时刻的当前状态量，其中，所述当前状态量包括以下至少之一：当前位置矢量、当前速度矢量、钟差和频差；基于所述当前状态量确定所述搭载卫星的卫星授时时刻；并将所述卫星授时时刻和所述当前状态量作为所述搭载卫星的导航结果。

进一步地，确定历史历元时刻所述星载导航接收机所搭载卫星的历史状态量包括：利用迭代最小二乘法获取所述历史历元时刻所对应的历史状态量。

进一步地，基于所述搭载卫星的轨道根数和所述历史状态量，计算所述搭载卫星的目标预测量包括：根据所述历史状态量建立线性观测方程，根据所述线性观测方程建立误差方程；通过迭代算法和最小二乘平差法求解所述误差方程，得到所述轨道根数的估计量；根据所述轨道根数的估计量和所述线性观测方程，得到所述搭载卫星的目标预测量，其中，所述目标预测量包括以下至少之一：位置预测量和速度预测量。

进一步地，基于所述目标预测量确定所述搭载卫星的第一导航观测量包括：在当前历元时刻采集导航卫星的观测信息，并基于所述观测信息获取所述搭载卫星的第二导航观测量，其中，所述第二导航观测量包括以下至少之一：第二伪距和第二载波相位，所述导航卫星为用于对所述搭载卫星提供导航信号的卫星；对所述目标预测量进行计算，得到预测信息；所述预测信息包括以下至少之一：所述搭载卫星的第三伪距和所述搭载卫星的第三载波相位；基于所述第二导航观测量和所述预测信息确定所述搭载卫星的第一导航观测量。

进一步地，基于所述第二导航观测量和所述预测信息确定所述搭载卫星的第一导航观测量包括：获取预设伪距门限范围和预设载波相位门限范围；基于所述第二伪距和所述第三伪距之间的差值以及所述预设伪距门限范围，在所述第二伪距和所述第三伪距中确定所述第一伪距；基于所述第二载波相位和第三载波相位之间的差值以及所述预设载波相位门限范围，在所述第二载波相位和所述第三载波相位中确定所述第一载波相位；将确定出的第一伪距和第一载波相位确定为第一导航观测量。

进一步地，对所述搭载卫星的所述第一导航观测量进行计算，得到所述搭载卫星在当前历元时刻的当前状态量包括：对所述第一导航观测量进行差分变换，得到差分变换结果；基于所述差分变换结果构建非线性滤波解算方程，并通过卡尔曼滤波法求解所述解算方程，并将结算结果作为所述搭载卫星在当前历元时刻的当前状态量。

进一步地，所述当前历元时刻为星载导航接收机的秒脉冲的上升沿；基于所述当前状态量确定所述搭载卫星的卫星授时时刻包括：对所述当前状态量中的钟差和频差进行质量改善处理，得到改善之后的钟差和改善之后的频差；根据所述改善之后的钟差和所述改善之后的频差调整所述秒脉冲的相位和频率，得到调整之后的秒脉冲，以使所述调整之后的秒脉冲的上升沿与标准时间整秒起始时刻同步；将所述调整之后的秒脉冲的上升沿确定为所述搭载卫星的卫星授时时刻。

进一步地，所述方法还包括：将所述卫星授时时刻广播给所述搭载卫星上的目标电子系统；在下一个卫星授时时刻到来之前，将所述当前状态量广播给所述目标电子系统。

进一步地，所述方法还包括：采用最小二乘残差平方和法对所述当前状态量进行自主完好性校验；对所述当前状态量进行合理性校验。

第二方面，本发明实施例还提供了一种基于轨道根数预测的星载导航接收机定位授时装置，包括：初始状态量获取模块、第一位置速度解算模块、观测量变换模块、第二位置速度解算模块、钟差和频差估计模块，其中，所述初始状态量获取模块，用于获取历史历元时刻所述星载导航接收机所搭载卫星的历史状态量，其中，所述历史状态量包括以下至少之一：历史位置矢量、历史速度矢量、历史钟差和历史频差；所述第一位置速度解算模块，用于基于所述搭载卫星的轨道根数和所述历史状态量，计算所述搭载卫星的目标预测量，其中，所述目标预测量包括：位置预测量和速度预测量；所述观测量变换模块，用于基于所述目标预测量确定所述搭载卫星的第一导航观测量，其中，所述第一导航观测量包括以下至少之一：第一伪距和第一载波相位；所述第二位置速度解算模块，用于对所述搭载卫星的所述第一导航观测量进行计算，得到所述搭载卫星的当前状态量，其中，所述当前状态量包括以下至少之一：当前位置矢量、当前速度矢量、钟差和频差；所述钟差和频差估计模块，用于基于所述当前状态量确定所述搭载卫星的卫星授时时刻；并将所述卫星授时时刻和所述当前状态量作为导航结果。

本发明实施例提供了一种基于轨道根数预测的星载导航接收机定位授时方法和装置，利用星载导航接收机所搭载卫星的历史状态量和轨道根数模型求解卫星当前的位置和速度矢量，以确保位置和速度矢量的精度，再进一步根据已经求解的位置和速度矢量求解钟差和频差。对比已有方法，本发明的方法可更加充分的利用卫星在轨运动规律对导航状态进行预测，从而改善导航定位和授时结果的精度，以更好的满足中高轨波束定向通信载荷等应用对在轨定位和授时精度的要求。

本发明的其他特征和优点将在随后的说明书中阐述，并且，部分地从说明书中变得显而易见，或者通过实施本发明而了解。本发明的目的和其他优点在说明书、权利要求书以及附图中所特别指出的结构来实现和获得。

为使本发明的上述目的、特征和优点能更明显易懂，下文特举较佳实施例，并配合所附附图，作详细说明如下。

附图说明

为了更清楚地说明本发明具体实施方式或现有技术中的技术方案，下面将对具体实施方式或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图是本发明的一些实施方式，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1为本发明实施例提供的中高轨导航接收机接收导航卫星的主瓣未被地球遮挡部分和旁瓣信号示意图；

图2为本发明实施例提供的一种基于轨道根数预测的星载导航接收机定位授时方法流程图；

图3为本发明实施例提供的一种基于目标预测量确定搭载卫星的第一导航观测量的方法流程图；

图4为本发明实施例提供的另一种基于轨道根数预测的星载导航接收机定位授时方法流程图；

图5为本发明实施例提供的位置和速度预测方法的模型示意图；

图6为本发明实施例提供的星载导航接收机定位授时方法中频差质量改善处理方法的流程图；

图7为本发明实施例提供的星载导航接收机定位授时方法中钟差质量改善处理方法的流程图；

图8为本发明实施例提供的一种基于轨道根数预测的星载导航接收机定位授时装置的示意图；

图9为本发明实施例提供的一种授时型星载导航接收机以及星载导航应用系统的功能连接框图；

图10为本发明实施例方法与现有技术的对比图。

图标：

100-授时型星载导航接收机；102-天线；104-射频前端模块；106-捕获跟踪模块；108-观测量获取模块；110-观测量变换模块；112-位置和速度解算模块；114-钟差和频差估计模块；116-钟差和频差质量改善模块；118-遥测模块；120-秒脉冲模块；10-星载导航应用系统；200-星上计算机；300-测控通信模块；400-被授时载荷；500-地面通信系统；600-星上数据总线；700-秒脉冲数据总线。

具体实施方式

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合附图对本发明的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

在本发明的描述中，需要说明的是，术语“中心”、“上”、“下”、“左”、“右”、“竖直”、“水平”、“内”、“外”等指示的方位或位置关系为基于附图所示的方位或位置关系，仅是为了便于描述本发明和简化描述，而不是指示或暗示所指的装置或元件必须具有特定的方位、以特定的方位构造和操作，因此不能理解为对本发明的限制。此外，术语“第一”、“第二”、“第三”仅用于描述目的，而不能理解为指示或暗示相对重要性。

另外，在本发明实施例的描述中，除非另有明确的规定和限定，术语“安装”、“相连”、“连接”应做广义理解，例如，可以是固定连接，也可以是可拆卸连接，或一体地连接；可以是机械连接，也可以是电连接；可以是直接相连，也可以通过中间媒介间接相连，可以是两个元件内部的连通。对于本领域的普通技术人员而言，可以具体情况理解上述术语在本发明中的具体含义。

图1是中高轨导航接收机接收导航卫星的主瓣未被地球遮挡部分和旁瓣信号示意图。可以看到，图中的宿主卫星因为轨道高度高于导航卫星，只能通过接收导航卫星的主瓣信号未被地球遮挡部分(图中给出的张角为7.45°)和旁瓣信号进行导航定位。

全球导航卫星系统主要包括美国的全球定位系统(Global Positioning System，GPS)、中国的北斗导航卫星系统(BeiDou Navigation Satellite System，BDS)、俄罗斯的全球导航卫星系统(Global Navigation Satellite System，GLONASS)和欧洲的伽利略卫星导航系统(Galileo Satellite navigation system，Galileo)。为改善定位结果，中高轨道星载导航接收机通常联合接收多种卫星导航信号。在我国，星载导航接收机最常见的体制为联合接收中国的北斗导航卫星系统信号和美国的全球定位系统信号。

需要说明的是，卫星的定位授时包括，星载导航接收机获取当前历元宿主卫星的位置矢量、速度矢量，称作定位；星载导航接收机输出秒脉冲及其上升沿所对应的世界协调时给宿主卫星，称作授时。

下面以GPS和加北斗双模接收机为例讲述本发明，本领域的一般技术人员都能将本实施例扩展到更多模的导航接收机。

实施例一：

根据本发明实施例，提供了一种基于轨道根数预测的星载导航接收机定位授时方法实施例，需要说明的是，在附图的流程图示出的步骤可以在诸如一组计算机可执行指令的计算机系统中执行，并且，虽然在流程图中示出了逻辑顺序，但是在某些情况下，可以以不同于此处的顺序执行所示出或描述的步骤。

图2是根据本发明实施例的一种基于轨道根数预测的星载导航接收机定位授时方法，应用于星载导航接收机，如图2所示，该方法包括如下步骤：

步骤S102，确定历史历元时刻星载导航接收机所搭载卫星的历史状态量，其中，历史状态量包括以下至少之一：历史位置矢量、历史速度矢量、历史钟差和历史频差。

在本发明中，不考虑宿主卫星(即，星载导航接收机所搭载卫星)自转和尺寸的影响，位置矢量、速度矢量既是宿主卫星在历元时刻的轨道位置和运动速度矢量，也是星载导航接收机的轨道位置和速度。星载导航接收机还会解算出当前历元自身的时钟与导航卫星时钟的频率误差(称为频差)和与世界协调时的差(称为钟差)。星载导航接收机会根据钟差和频差调整下一个历元时刻即秒脉冲的上升沿，使其最靠近世界协调时整秒。

其中，星载导航接收机一般会在最靠近世界协调时整秒的时刻采集导航观测量，进行导航解算，而这个时刻称为历元。该历元时刻也是星载导航接收机的秒脉冲的上升沿。

需要说明的是，在本方法实施例中，确定历史历元时刻星载导航接收机所搭载卫星的历史状态量所用方法为迭代最小二乘法。

步骤S104，基于搭载卫星的轨道根数和历史状态量，计算搭载卫星的目标预测量，其中，目标预测量包括：位置预测量和速度预测量。

其中，轨道根数是描述卫星在轨运动的一系列参数，本领域的一般技术人员都能理解。常见的轨道根数有北美航空司令部的两行根数、16根数和26根数等。两行根数通常用作卫星轨道预报、16根数通常在导航卫星中使用，而更复杂的根数比如26根数通常用于更加准确的描述低轨卫星的在轨运动。

面向中高轨应用场景，本实施例使用如下16根数作为轨道根数：

其中，这些参数含义如下表所示或者请参考导航卫星的星历定义：

步骤S106，基于目标预测量确定搭载卫星的第一导航观测量，其中，第一导航观测量包括以下至少之一：第一伪距和第一载波相位。

步骤S108，对搭载卫星的第一导航观测量进行计算，得到搭载卫星在当前历元时刻的当前状态量，其中，当前状态量包括以下至少之一：当前位置矢量、当前速度矢量、钟差和频差。

步骤S110，基于当前状态量确定搭载卫星的卫星授时时刻；并将卫星授时时刻和当前状态量作为搭载卫星的导航结果。

本发明实施例提供了一种基于轨道根数预测的星载导航接收机定位授时方法，利用星载导航接收机所搭载卫星的历史状态量和轨道根数模型求解卫星当前的位置和速度矢量，以确保位置和速度矢量的精度，再进一步根据已经求解的位置和速度矢量求解钟差和频差。对比已有方法，本发明的方法可更加充分的利用卫星在轨运动规律对导航状态进行预测，从而改善导航定位和授时结果的精度，以更好的满足中高轨波束定向通信载荷等应用对在轨定位和授时精度的要求。

具体地，在上述步骤S104中，基于搭载卫星的轨道根数和历史状态量，计算搭载卫星的目标预测量，具体包括如下步骤：

步骤S1041，根据历史状态量建立线性观测方程，根据线性观测方程建立误差方程；

步骤S1042，通过迭代算法和最小二乘平差法求解所述误差方程，得到轨道根数的估计量；

步骤S1043，根据轨道根数的估计量和线性观测方程，得到搭载卫星的目标预测量，其中，目标预测量包括以下至少之一：位置预测量和速度预测量。

可选地，图3为本发明实施例提供的一种基于目标预测量确定搭载卫星的第一导航观测量的方法流程图，如图3所示，上述步骤S106具体包括下步骤：

步骤S301，在当前历元时刻采集导航卫星的观测信息，并基于观测信息获取搭载卫星的第二导航观测量，其中，第二导航观测量包括以下至少之一：第二伪距和第二载波相位，导航卫星为用于对搭载卫星提供导航信号的卫星；

步骤S302，对目标预测量进行计算，得到预测信息；预测信息包括以下至少之一：搭载卫星的第三伪距和搭载卫星的第三载波相位；

步骤S303，基于第二导航观测量和预测信息确定搭载卫星的第一导航观测量。

通过以上描述可知，在本发明实施例中，通过搭载卫星的轨道根数和历史状态量，能够确定搭载卫星的位置预测量和速度预测量，依此对导航卫星的观测量进行量化估计，并且根据预测信息对导航观测量进行筛选，从而能够改善导航定位和授时结果的精度。

具体地，步骤S303包括如下步骤：

步骤S3031，获取预设伪距门限范围和预设载波相位门限范围；

其中，预设伪距门限范围和预设载波相位门限范围可根据需要来限定，优选地，预设伪距门限范围可以为3米至30米，更为优选地，预设伪距门限范围为10米。预设载波相位门限范围可以为0.1米至10米。优选地，预设载波相位门限范围为1米。

步骤S3032，基于第二伪距和第三伪距之间的差值以及预设伪距门限范围，在第二伪距和第三伪距中确定所述第一伪距；

例如，获取预设伪距门限范围为10m，如果第二伪距和第三伪距之间的差值超过了预设伪距门限范围，即，第二伪距和第三伪距之间的差值大于10m，则将第三伪距确定为第一伪距；如果第二伪距和第三伪距之间的差值没有超过预设伪距门限范围，即，第二伪距和第三伪距之间的差值小于等于10m，则将第二伪距确定为第一伪距。

步骤S3033，基于第二载波相位和第三载波相位之间的差值以及预设载波相位门限范围，在第二载波相位第三载波相位中确定第一载波相位；

例如，获取预设载波相位门限范围为1m，如果第二载波相位和第三载波相位之间的差值超过了预设载波相位门限范围，即，第二载波相位和第三载波相位之间的差值大于1m，则将第三载波相位确定为第一载波相位；如果第二载波相位和第三载波相位之间的差值没有超过预设载波相位门限范围，即，第二载波相位和第三载波相位之间的差值小于等于1m，则将第二载波相位确定为第一载波相位。

步骤S3034，将确定出的第一伪距和第一载波相位确定为第一导航观测量。

具体地，步骤S108包括如下步骤：

步骤S1081，对第一导航观测量进行差分变换，得到差分变换结果；

步骤S1082，基于差分变换结果构建非线性滤波解算方程，并通过卡尔曼滤波法求解解算方程，并将结算结果作为搭载卫星在当前历元时刻的当前状态量。

具体地，步骤S110包括如下步骤：

步骤S1101，对所述当前状态量中的钟差和频差进行质量改善处理，得到改善之后的钟差和改善之后的频差；

步骤S1102，根据所述改善之后的钟差和所述改善之后的频差调整所述秒脉冲的相位和频率，得到调整之后的秒脉冲，以使所述调整之后的秒脉冲的上升沿与标准时间整秒起始时刻同步；其中，上述当前历元时刻即为星载导航接收机的秒脉冲的上升沿；

步骤S1103，将所述调整之后的秒脉冲的上升沿确定为所述搭载卫星的卫星授时时刻。

通过以上描述可知，本发明实施例提供的基于轨道根数预测的星载导航接收机定位授时方法，通过差分变换将导航观测量中与时钟有关的部分消除，利用变换后的结果和轨道根数模型求解卫星当前的位置和速度矢量，以确保位置和速度矢量的精度，再进一步根据已经求解的位置和速度矢量求解钟差和频差。对比已有方法，本发明的方法可更加充分的利用卫星在轨运动规律对导航状态进行预测，从而改善导航定位和授时结果的精度，以更好的满足中高轨波束定向通信载荷等应用对在轨定位和授时精度的要求。

实施例二：

根据本发明实施例，提供了另一种基于轨道根数预测的星载导航接收机定位授时方法实施例，需要说明的是，在附图的流程图示出的步骤可以在诸如一组计算机可执行指令的计算机系统中执行，并且，虽然在流程图中示出了逻辑顺序，但是在某些情况下，可以以不同于此处的顺序执行所示出或描述的步骤。

图4是根据本发明实施例的另一种基于轨道根数预测的星载导航接收机定位授时方法，应用于星载导航接收机，如图4所示，该方法包括如下步骤：

S1，星载导航接收机进行状态量初始化；

S2，星载导航接收机在当前历元时刻采集导航观测量(即上述实施例一种的第二导航观测量，本实施例中均称为导航观测量)；

S3，基于轨道根数和多个历史状态量进行计算当前历元位置(即上述实施例一中的位置预测量，本实施例中均称为当前历元位置)和速度预测量；

S4，进行异常观测量剔除以获得保留伪距(即上述实施例一中的第一伪距，本实施例中均称为保留伪距)和保留载波相位(即上述实施例一中的第一载波相位，本实施例中均称为保留载波相位)；

S5，通过差分变换消除导航观测量中与时钟有关的部分；

S6，利用所述变换后的导航观测量和非线性滤波解算当前历元的位置和速度；

S7，根据所述当前历元的状态量求解钟差和频差；

S8，对所述当前历元状态量进行自主完好性校验和合理性校验；

S9，对所述频差进行质量改善处理；

Sa，对所述钟差进行质量改善处理；

Sb，调整下一历元时刻，使其与世界协调时整秒时刻的同步误差最小；

Sc，输出导航定位结果，回到步骤S2。

上述步骤S1进一步包括：

S21，设定初始化历元数m0；

S22，比较当前历元数k和所述初始化历元数m0；

S23，在当前历元时刻采集导航观测量；

S24，当k<m0时，进行初始化，具体包括：使用迭代最小二乘法获得当前历元状的位置、速度、钟差以及频差；

S25，根据所述钟差以及频差调整历元时刻，使其与世界协调时整秒时刻的同步误差最小；

S26，当k≥m0时，初始化结束；

在本实施例中，为了16个轨道根数以及初始化首个与世界协调时整秒同步的历元时刻，取所述初始化历元数m0至少为16。在初始化完成后，初始参考时刻变成已知数，这时只需要至少15个历元状态量就可以求解出新的轨道根数。

图5为本发明实施例提供的位置和速度预测方法的模型示意图。在上述步骤S3中，星载导航接收机基于轨道根数和多个历史状态量进行计算当前历元位置和速度预测量。如图5所示，选择多个历史状态量计算出轨道根数，轨道根数从本质上描述了宿主卫星的运动轨迹。基于轨道根数，可以计算出稳定的位置和速度预测量。选择的多个历史历元状态量一方面要能够代表一个适当长度的宿主卫星轨道弧段，比如宿主卫星π/6地心张角的轨道弧段，另一方面，要选择数量适当的历元状态量数目，确保降低噪声和异常的影响准确稳定的获得轨道根数。

在本发明实施例的一个可选实施方式中，步骤S3进一步包括以下步骤：

S31，设定初始化单位权误差相对变化量ε；

S32，比较当前单位权误差相对变化量ε_i和所述初始化单位权误差相对变化量ε；

S33，当ε_i≥ε时，进行迭代，具体包括：根据得到的线性观测方程建立误差方程，使用最小二乘平差法获得当前历元轨道根数修正量和当前历元轨道根数估计量；

S34，当ε_i<ε时，根据所述当前轨道根数估计量，得到当前历元位置和速度预测量。

上述步骤S33进一步包括：

S331，选择预设数量历史历元，将这些历元的位置矢量作为观测列向量；

S332，根据轨道根数建立线性观测方程建立误差方程；

S333，根据所述误差方程，利用最小二乘原理获得当前历元轨道根数修正量；

S334，根据所述修正量，获得当前历元轨道根数估计量。

在上述步骤S331中，选择足够数量历史历元的方法，包括以下步骤：

S3311，确定备选历史历元所在的轨道弧段，一般为包括当前历元的π/6地心张角轨道弧段；

S3312，剔除备选历史历元所在的轨道弧段中几何因子最大的1/3历元；

S3313，在剩余历元中按照时间间隔均匀的选择m个；

在上述步骤S331中，所述观测列向量满足如下关系：

其中p_k＝[p_k,x,p_k.y,p_k,z],k＝1,…,m为所述历元时刻的轨道位置。

在上述步骤S332中，所述的线性观测方程满足如下关系：

上述线性观测方程是由s(t)＝s(R,t)在R_i处泰勒展开并舍去二阶和二阶以上小量后的结果。其中，R为所述待求轨道根数中除初始参考时刻外的15个参数，具体地，R_i为估值在第i次迭代后的结果，δR_i为所述待求轨道根数中除初始参考时刻外的15个参数的修正量，为轨道位置对所述轨道根数的偏导数。

在WGS-84地心地固直角坐标系中，所述轨道根数和轨道位置[p_k,x,p_k,y,p_k,z]满足如下关系：

其中，上述关系式中：

为地球自转角速度常数值；

r＝a(1-e cos E)+C_rssin(2φ)+C_rccos(2φ)；

u＝φ+C_ussin(2φ)+C_uccos(2φ)；

E可通过迭代方法求解，E_j＝M+esin(E_j-1)，GM为万有引力相干常数，j为迭代次数，通常地，迭代2至3次即可得到十分精确的结果。

在上述步骤S332中，所述的误差方程为:

其中，L＝s(t)-s(R_i,t)，分别为x,y,z三个方向分量的雅克比矩阵。

在上述步骤S333中，所述当前历元轨道根数修正量满足如下关系：

上述方程是由最小二乘原理得到的。其中，P为所述星载接收机位置观测向量的权矩阵。

在上述步骤S334中，所述当前历元轨道根数估计量满足如下关系：

在上述步骤S34中，所述轨道根数和所述当前历元状态量X′_k中接收机的位置和速度满足如下关系：

上述位置表达式与上述步骤S332中的关系相同，上述速度表达式中，

其中：

上述步骤S4进一步包括：

S41，根据所述当前历元位置和速度预测量计算预测伪距(即上述实施例一中的第三伪距)和预测载波相位(即上述实施例一种的第三载波相位)；

S42，比较所述观测伪距(即上述实施例一中的第三伪距)和预测伪距以及比较所述观测载波相位(即上述实施例一中的第三载波相位)和预测载波相位来确定所述保留伪距和保留载波相位。

在上述步骤S41中，对于当前历元k，有伪距和载波相位满足如下关系：

其中，i＝1,…,n_s为接收到的导航卫星序号，为第i颗导航卫星的位置，为当前历元接收机的位置,为第i颗导航卫星的速度，为当前历元接收机的速度，为未知接收机钟差，为未知接收机频差，为伪距高斯噪声，为相位高斯噪声。

在上述步骤S42中，通过比较预测值(预测伪距和预测载波相位)与观测值(观测伪距和观测载波相位)，将与预测值相差较大的观测值剔除。具体地，剔除过程包括以下步骤：

S421，设定伪距门限范围(即上述实施例一中的预设伪距门限范围)以及载波相位门限范围(即上述实施例一中的预设载波相位门限范围)；

S422，如果观测伪距和预测伪距的差值超出所述伪距门限范围，则选取所述预测伪距作为所述保留伪距，否则选取所述观测伪距作为所述保留伪距；同时，如果所述观测的载波相位和预测载波相位的差值超出所述载波相位门限值，则选取所述预测载波相位作为所述保留载波相位，否则，选取所述观测的载波相位作为所述保留载波相位。

所述伪距门限范围和所述载波相位门限范围可根据需要来限定，优选地，所述伪距门限范围可以为3米至30米，更为优选地，所述伪距门限范围为10米。所述载波相位门限范围可以为0.1米至10米。优选地，所述载波相位门限范围为1米。

上述步骤S5进一步包括：

S51，根据伪距定位原理构建所述GPS和北斗联合解算方程；

S52，通过差分变换得到变换后的导航观测量矩阵，在变换后的导航观测量矩阵中，与时钟有关的部分被消除；

在上述步骤S51中，所述解算方程为：

ρ＝r+δt_u-δt+I+T+ε_ρ

其中，ρ为保留伪距，r为接收机与卫星的距离，δt_u为接收机钟差，δt为卫星钟差，I为电离层延时，T为对流层延时，ε_ρ为伪距测量误差。其中，时间参量δt、I、T、ε_ρ均以该时间参量乘以光速后的距离量表示。进一步地，可将所述解算方程简化为校正后的伪距观测方程：

ρ_c＝r+δt_u

其中，ρ_c为校正后的伪距。

在步骤S52中，可将所述伪距信息和载波相位信息写成矩阵形式如下：

其中，为接收到的所有GPS卫星伪距，n_gps为接收到GPS卫星总数，为接收到的所有GPS卫星载波相位，为接收到的所有BDS卫星伪距，n_bd为接收到BDS卫星总数，为接收到的所有BDS卫星载波相位。

为了消除伪距和载波相位表达式中的钟差和频差，我们分别构造(n_gps-1)×n_gps阶差分映射矩阵D_gps＝[I_gps,-1_gps]和(n_bd-1)×n_bd阶差分映射矩阵D_bd＝[I_bd,-1_bd]，其中I_gps为(n_gps-1)×(n_gps-1)阶单位矩阵，I_bd为(n_bd-1)×(n_bd-1)阶单位矩阵，-1_gps和-1_bd分别为全1列向量。我们构造总差分映射矩阵D，有

那么通过差分映射矩阵变换后的观测量就只与接收机和导航卫星的位置和速度有关，与时钟相关的部分被消除，得到变换后的结果Z_k，有

第二个差分变换是利用历史历元载波周期φ_i,k对载波相位进行平滑，有

其中，λ为导航信号载波的波长，这个变换相当于一个z变换如下的FIR滤波器：

在上述步骤S6中，我们定义新的状态为所述星载导航接收机第k个历元的位置和速度矢量。

所述非线性滤波法可为Cubature卡尔曼滤波法(CKF)。

具体地，该离散非线性系统可表示为

其中，Y_k为当前历元观测量(接收机位置和速度),当前历元随机系统噪声W_k～N(0,Q)，当前历元随机观测噪声V_k～N(0,R)，f(J_k)为系统状态方程，h(J_k)为系统观测方程。则CKF状态更新方法如下：

对于状态更新，假设k-1历元的后验概率密度已知，对状态误差协方差矩阵进行Cholesky分解，有

那么Cubature采样点可计算如下，对i＝1,2,…,2n，n为状态向量维数：

其中，传播到第k历元为

第k历元的状态预测估计为

状态误差协方差预测估计值为

对于观测更新，对误差协方差P_k|k-1进行Cholesky分解，有

那么Cubature采样点可计算如下:

传播Cubature采样点，有

Y_i,k|k-1＝h(J_i,k|k-1)

第k历元的观测预测估计为

观测误差自协方差预测估计值为

观测误差互协方差预测估计值为

估计卡尔曼增益：

最终状态估计值为

第k历元状态协方差为

在上述步骤S7中，所述钟差和频差满足如下关系：

其中，LSE为最小二乘估计。

上述步骤S8进一步包括：

S81，对所述当前历元的状态量X_k进行自主完好性校验，若校验通过，则采纳X_k作为当前历元的状态量，否则，采纳X′_k作为当前历元的状态量；

S82，合理性校验：若发生连续自主完好性校验不通过或定位结果明显不合理，则回到步骤S1，重新开始状态初始化。

在上述步骤S81中，可采用最小二乘残差自主完好性校验方法，即最小残差平方和法。其中残差向量b的每一分量定义为星载导航接收机对相应卫星的校正后伪距值ρ_c减去几何距离预测值r与星载导航接收机钟差δt_u之和：

b＝ρ_c-r-δt_u

而残差平方和εSSE为标量，定义为残差向量与其转置的乘积：

ε_SSE＝b^Tb

ε_SSE的大小体现着各个测量值之间的一致性。由于定位解算中有四个独立有效的控制方程，因此ε_SSE在理论上服从自由度为N-4的χ²分布，其中N为有效卫星的个数。根据χ²分布概率密度函数，可以由设定的误警率P_fa确定ε_SSE的门限值T_SSE，当ε_SSE超过门限值T_SSE时则认为各测量值之间不一致，也就是其中某些测量值出现了错误，即表示校验不通过。

可见卫星数越多且卫星几何分布越好时，自主完好性校验算法的可靠性就越高。优选地，当同一时刻星载导航接收机可见卫星数大于等于5时，所述自主完好性校验才能够有效进行。

在上述步骤S82中，可设定自主完好性校验连续不通过最大历元数N_SSE，当连续N_SSE个历元发生自主完好性校验不通过时，或定位结果明显与常识或预期不符时，说明合理性校验不通过，即定位解算部分功能出现错误，此时重新返回步骤S31，开始初始化过程，以保证定位功能可以正常有效运行。

此外，若未出现合理性校验不通过的现象，可令k＝k+1，循环进入步骤S3，进行下一历元的解算方程构建和接收机状态解算工作。

可选地，如图6所示，在上述步骤S9中，对各个历元的频差进行质量改善处理，具体包括以下步骤：

S91，设定初始化历元数N_f；

S92，比较当前历元数k和所述初始化历元数N_f；

S93，当k<N_f时，设定经过质量改善处理的当前历元的频差等于当前历元的频差df_k，即：以及

S94，当k≥N_f时，设定经过质量改善处理的当前历元的频差其中，α为滤波系数，且0<α≤1。

在上述步骤S91中，所述初始化历元数N_f可根据需要自行设定，如可为10个历元、15个历元。

在上述步骤S93中，当前历元数k小于所述初始化历元数N_f时，即当已经历的历元数未达到给定的初始化历元数N_f时，不对当前历元的频差进行质量改善处理。

在上述步骤S94中，当当前历元数k大于等于所述初始化历元数N_f时，才对当前历元的频差进行质量改善处理。通过加权平均的方式可使所述星载导航接收机频差的变化得到平滑，从而可减小随机噪声误差的影响以及提高所述星载导航接收机频差的精度。

可选地，如图7所示，在上述步骤Sa中，对所述钟差进行质量改善处理，具体包括以下步骤：

Sa1，设定初始化历元数N_t以及钟差限定量Δ_max；

Sa2，比较当前历元数k和初始化历元数N_t；

Sa3，当k<N_t时，设定经过质量改善处理的当前历元的钟差等于当前历元的钟差dt_k，即：以及

Sa4，当k≥N_t时，判断当前历元的dt_k与差值是否超出最大范围[-Δ_max,Δ_max]；

Sa5，若超过了所述最大范围，且成立，则令

Sa6，若超过了所述最大范围，且成立，则令以及

Sa7，若未超出所述最大范围，则令

在上述步骤SA1，所述钟差限定量Δ_max用于在k≥N_t时，将当前历元的钟差dt_k与经过改善处理的钟差的差值限定在所述最大范围内。所述钟差限定量Δ_max的取值范围为100纳秒至1微秒。

在上述步骤Sa3中，当已经历的历元数k未达到所述设定的初始化历元数N_t时，不对当前历元的频差进行质量改善处理。

在上述步骤Sa5-Sa7中，通过将当前历元的钟差与当前历元质量改善处理后的频差的差值限定在给定的最大范围[-Δ_max,Δ_max]内，对超出范围的值，以范围边界取代之。该方法保证钟差与频差的差值在给定范围内，从而提高钟差精度。

在Sc中，输出导航定位结果，其具体步骤为：

步骤Sc1，所述星载导航接收机可将所述调整后的秒脉冲通过秒脉冲广播总线广播给星上其它电子系统。所述秒脉冲广播总线优选为多点低电压差分总线，该多点低电压差分总线可以使所述秒脉冲以低电压差分信号传输，从而可实现几百Mbps的传输速率，并且低电压差分信号的低压幅和低电流驱动输出可进一步降低噪声和减少功耗。此外，采用多点低电压差分总线传输还可减小电磁干扰。

步骤Sc2，所述星载导航接收机将所述位置、速度、时间转换为满足卫星总线设计要求的信息格式，并通过卫星数据总线广播给卫星上其它电子系统。其中，所述信息格式包括精确到秒的时间、速度、位置、精度几何因子以及校验和。所述校验和是对精确到秒的时间、速度、位置、精度几何因子的校验结果。

步骤Sc3，所述星载导航接收机生成和更新与时间相关的遥测量，并在星上计算机请求对应遥测量时给出应答的方法，其中星载导航接收机生成和更新的与时间相关的遥测量包括：升交点地方时，降交点地方时，以及所述{时间，位置，速度}组成的三元组。

升交点为星载导航接收机沿轨道从南向北运动时与赤道面的焦点，该点地方时即为升交点地方时。降交点为星载导航接收机沿轨道从北向南运动时与赤道面的焦点，该点地方时即为降交点地方时。所述时间、位置、速度三元组表示星载导航接收机经过解算和自主完好性、合理性校验后得到的接收机状态解算结果。

在上述步骤Sc3中，获得所述升交点地方时和降交点地方时的方法为：用地心地固坐标系(Earth-Fixed Coordinate System，简称ECEF)表示解算结果中的位置，当当前历元最终解算结果的位置的Z轴解算结果与上一历元相比符号发生变化时，进行如下判断：如果Z轴速度为正值，则判断此时的时间为升交点地方时；否则，认为此时的时间为降交点地方时。

当当前历元最终解算结果的位置的Z轴解算结果与上一历元相比符号发生变化时，说明星载导航接收机在该历元通过了轨道与赤道面交点，此时若Z轴速度为正值，根据地心地固坐标系定义，可知星载导航接收机在该历元的运动方向为由南向北运动，从而此时的时间为升交点地方时；若Z轴速度为负值，根据地心地固坐标系定义，可知星载导航接收机在该历元的运动方向为由北向南运动，从而此时的时间为降交点地方时。

在上述步骤Sc3中，更新所述遥测量中所述{时间、位置、速度}三元组的方法为：

对之前缓存的{时间，位置，速度}三元组的遥测状态进行检测，只有当之前缓存的{时间，位置，速度}三元组已经被星上计算机遥测完成后，才用当前历元的{时间，位置，速度}三元组替代之前缓存的{时间，位置，速度}三元组，否则不进行更新。

通过以上描述可知，本发明实施例提供的基于轨道根数预测的星载导航接收机定位授时方法，通过差分变换将导航观测量中与时钟有关的部分消除，利用变换后的结果和轨道根数模型求解卫星当前的位置和速度矢量，以确保位置和速度矢量的精度，再进一步根据已经求解的位置和速度矢量求解钟差和频差。对比已有方法，本发明的方法可更加充分的利用卫星在轨运动规律对导航状态进行预测，从而改善导航定位和授时结果的精度，以更好的满足中高轨波束定向通信载荷等应用对在轨定位和授时精度的要求。

实施例三：

本发明实施例还提供了一种基于轨道根数预测的星载导航接收机定位授时装置，该装置主要用于执行本发明实施例上述内容所提供的一种基于轨道根数预测的星载导航接收机定位授时方法，以下对本发明实施例提供的一种基于轨道根数预测的星载导航接收机定位授时装置做具体介绍。

图8是根据本发明实施例的一种基于轨道根数预测的星载导航接收机定位授时装置的示意图，如图8所示，该装置主要包括：初始状态量获取模块81、第一位置速度解算模块82、观测量变换模块83、第二位置速度解算模块84、钟差和频差估计模块85，其中，

初始状态量获取模块81，用于获取历史历元时刻星载导航接收机所搭载卫星的历史状态量，其中，历史状态量包括以下至少之一：历史位置矢量、历史速度矢量、历史钟差和历史频差；

第一位置速度解算模块82，用于基于搭载卫星的轨道根数和历史状态量，计算搭载卫星的目标预测量，其中，目标预测量包括：位置预测量和速度预测量；

观测量变换模块83，用于基于目标预测量确定所述搭载卫星的第一导航观测量，其中，第一导航观测量包括以下至少之一：第一伪距和第一载波相位；

第二位置速度解算模块84，用于对搭载卫星的第一导航观测量进行计算，得到搭载卫星的当前状态量，其中，当前状态量包括以下至少之一：当前位置矢量、当前速度矢量、钟差和频差；

钟差和频差估计模块85，用于基于当前状态量确定搭载卫星的卫星授时时刻；并将卫星授时时刻和当前状态量作为导航结果。

所属领域的技术人员可以清楚地了解到，为描述的方便和简洁，上述描述的系统和装置的具体工作过程，可以参考前述方法实施例中的对应过程，在此不再赘述。

可选地，如图9所示，本发明实施例还提供了一种授时型星载导航接收机100以实现星载导航接收机定位授时方法。

如图9所示，该授时型星载导航接收机100包括：天线102、射频前端模块104、捕获跟踪模块106、观测量获取模块108、观测量变换模块110、位置和速度解算模块112、钟差和频差估计模块114、钟差和频差质量改善模块116、遥测模块118以及秒脉冲模块120，其中，

天线102，用于接收卫星导航信号；

射频前端模块104，用于接收并处理卫星导航信号，其中，处理操作包括以下至少之一：滤波、变频以及模数采样；变频可以为下变频；

捕获跟踪模块106，用于接收经射频前端模块104处理后的卫星导航信号，并捕获和跟踪目标卫星的导航信号；

观测量获取模块108，用于接收跟踪后的目标卫星的导航信号以及秒脉冲，并在秒脉冲的上升沿采集该跟踪后的目标卫星的导航信号，获得观测量；

观测量变换模块110，用于接收观测量，进行异常观测量剔除以获得保留观测量，并对保留观测量进行差分变换，获得变换后的观测量；其中，保留观测量包括以下至少之一：保留伪距和保留载波相位；

位置速度解算模块112，用于基于轨道根数和多个历史状态量进行计算当前历元位置和速度预测量，输出给观测量变换模块110，同时接收来自观测量变换模块110的变换后的观测量，构建解算方程，并通过解算获得星载导航接收机的位置、速度；

钟差和频差估计模块114，用于根据星载导航接收机的位置和速度进行最小二乘估计，获得时间、钟差和频差；

钟差和频差质量改善模块116，用于对钟差和频差估计模块114解算获得的频差和钟差进行质量改善处理；

遥测模块118，用于生成和更新与时间相关的遥测量，并在星上计算机请求对应遥测量时给出应答；

秒脉冲模块120，用于产生所述秒脉冲，并根据质量改善处理后的频差和钟差调节秒脉冲的频率和相位，使秒脉冲的上升沿与标准时间整秒起始时刻同步。

本发明实施例所提供的装置，其实现原理及产生的技术效果和前述方法实施例相同，为简要描述，装置实施例部分未提及之处，可参考前述方法实施例中相应内容。

如图9所示，本发明还进一步提供一种星载导航应用系统10，包括：授时型星载导航接收机100、星上计算机200、测控通信模块300、被授时载荷400以及地面通信系统500，其中，星上计算机200通过星上数据总线600与遥测模块118相连接，被授时载荷400通过秒脉冲数据总线700与秒脉冲模块120相连接，

测控通信模块300，用于在星上计算机200与地面通信系统500之间实现指令与数据传输；

星上计算机200，用于执行卫星上日常管理和记录功能，包括对授时型星载导航接收机100的日常管理和记录；

被授时载荷400，用于为卫星上需要授时的信息提供专门应用载荷，如被授时载荷400可以为移动通信载荷；

地面通信系统500，用于为地面人员对卫星进行远程控制和跟踪提供输入输出界面。

可选地，地面通信系统500可对授时型星载导航接收机100进行远程控制和跟踪。

通过以上描述可知，本发明实施例提供的基于轨道根数预测的星载导航接收机定位授时方法和装置，通过差分变换将导航观测量中与时钟有关的部分消除，利用变换后的结果和轨道根数模型求解卫星当前的位置和速度矢量以确保位置和速度矢量的精度，再进一步根据已经求解的位置和速度矢量求解钟差和频差。对比已有方法，本发明的方法可更加充分的利用卫星在轨运动规律对导航状态进行预测，从而改善导航定位和授时结果的精度，以更好的满足中高轨波束定向通信载荷等应用对在轨定位和授时精度的要求。

图10为本发明实施例方法与现有技术的对比图。星载导航接收机定位授时方法在文献中已有较多方法，它们之间的本质区别在于对历史状态的利用程度和处理细节的不同。本发明利用轨道根数实现了π/6左右弧度的利用，兼顾了实时授时和较高精度实时定位解算的需求，特别适合因为灵敏度受限可用导航卫星数不足或者轨道高度太高导致导航卫星信号太弱的场景。

对比已有卡尔曼滤波+J2摄动累方法，本发明的方法可更加充分的利用卫星在轨运动规律对导航状态进行预测，从而改善导航定位和授时结果的精度，以更好的满足中高轨波束定向通信载荷等应用对在轨定位和授时精度的要求。

另外，本领域技术人员还可以在本发明精神内做其它变化，当然，这些依据本发明精神所做的变化，都应包含在本发明所要求保护的范围之内。

除非另外具体说明，否则在这些实施例中阐述的部件和步骤的相对步骤、数字表达式和数值并不限制本发明的范围。

在这里示出和描述的所有示例中，任何具体值应被解释为仅仅是示例性的，而不是作为限制，因此，示例性实施例的其他示例可以具有不同的值。

应注意到：相似的标号和字母在下面的附图中表示类似项，因此，一旦某一项在一个附图中被定义，则在随后的附图中不需要对其进行进一步定义和解释。

附图中的流程图和框图显示了根据本发明的多个实施例的系统、方法和计算机程序产品的可能实现的体系架构、功能和操作。在这点上，流程图或框图中的每个方框可以代表一个模块、程序段或代码的一部分，所述模块、程序段或代码的一部分包含一个或多个用于实现规定的逻辑功能的可执行指令。也应当注意，在有些作为替换的实现中，方框中所标注的功能也可以以不同于附图中所标注的顺序发生。例如，两个连续的方框实际上可以基本并行地执行，它们有时也可以按相反的顺序执行，这依所涉及的功能而定。也要注意的是，框图和/或流程图中的每个方框、以及框图和/或流程图中的方框的组合，可以用执行规定的功能或动作的专用的基于硬件的系统来实现，或者可以用专用硬件与计算机指令的组合来实现。

在本申请所提供的几个实施例中，应该理解到，所揭露的系统、装置和方法，可以通过其它的方式实现。以上所描述的装置实施例仅仅是示意性的，例如，所述单元的划分，仅仅为一种逻辑功能划分，实际实现时可以有另外的划分方式，又例如，多个单元或组件可以结合或者可以集成到另一个系统，或一些特征可以忽略，或不执行。另一点，所显示或讨论的相互之间的耦合或直接耦合或通信连接可以是通过一些通信接口，装置或单元的间接耦合或通信连接，可以是电性，机械或其它的形式。

所述作为分离部件说明的单元可以是或者也可以不是物理上分开的，作为单元显示的部件可以是或者也可以不是物理单元，即可以位于一个地方，或者也可以分布到多个网络单元上。可以根据实际的需要选择其中的部分或者全部单元来实现本实施例方案的目的。

另外，在本发明各个实施例中的各功能单元可以集成在一个处理单元中，也可以是各个单元单独物理存在，也可以两个或两个以上单元集成在一个单元中。

最后应说明的是：以上所述实施例，仅为本发明的具体实施方式，用以说明本发明的技术方案，而非对其限制，本发明的保护范围并不局限于此，尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，其依然可以对前述实施例所记载的技术方案进行修改或可轻易想到变化，或者对其中部分技术特征进行等同替换；而这些修改、变化或者替换，并不使相应技术方案的本质脱离本发明实施例技术方案的精神和范围，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应所述以权利要求的保护范围为准。