CN105116422A - 一种星载高动态gnss接收机及其导航方法 - Google Patents

Abstract

本发明公开了一种星载高动态GNSS接收机及其导航方法，该接收机基于ARM+FPGA，包括顺次相连的射频前端处理模块、基带信号数字处理模块和定位解算模块，导航方法如下：射频前端处理模块对卫星信号进行放大、变频、滤波和模数转换，最终得到数字中频信号；基带信号数字处理模块对数字中频信号进行捕获、跟踪、位同步和帧同步处理，得到导航测量值和导航电文；定位解算模块利用导航测量值和导航电文进行定位解算，最终得到用户的卫星信息，获取定位观测值。本发明采用ARM+FPGA的架构，具有很强的灵活性，能够加载不同的算法，且FPGA的并行处理结构能够保证GNSS接收机的实时性。

Description

一种星载高动态GNSS接收机及其导航方法

技术领域

本发明涉及卫星导航技术领域，特别是一种基于ARM+FPGA的星载高动态GNSS接收机及其导航方法。

背景技术

GNSS全球卫星导航系统具有全球、全天候、连续三维定位和导航能力的系统。GNSS系统虽然主要为地面应用而设计，但也给航天应用带来了显著的便利。一台低成本星载GNSS接收机可以提供许多关键航天传感器的功能，包括定轨、姿态和姿态变化率测量以及精密授时。采用GNSS技术的航天器比现有同等性能的系统在价格、功率、质量和可靠性上有显著改善。GNSS的航天应用对星载GNSS接收机的动态性要求很高，同时由于宇宙中的恶劣环境，使得GNSS接收机在太空中的应用于地面和航空等的应用不一样，因此在设计星载GNSS接收机要考虑很多特点。

卫星在发射和升空阶段要经历很高的动态，加速度会超过10g，卫星在轨运行时速度将超过7km/s，现有GNSS接收机未充分考虑高动态对GNSS信号的捕获和跟踪的影响，且存在结构复杂、兼容性差、成本高的问题。

发明内容

本发明的目的在于提供一种能够为低轨卫星提供高精度导航定位以及精确授时的星载高动态GNSS接收机及其导航方法。

实现本发明目的的技术解决方案为：一种星载高动态GNSS接收机，该接收机基于ARM+FPGA，包括顺次相连的射频前端处理模块、基带信号数字处理模块和定位解算模块：

射频前端处理模块，对卫星信号进行放大、变频、滤波和模数转换，最终得到数字中频信号；

基带信号数字处理模块，对数字中频信号进行捕获、跟踪、位同步和帧同步处理，得到导航测量值和导航电文；

定位解算模块，利用导航测量值和导航电文进行定位解算，最终得到用户的卫星信息，获取定位观测值。

一种基于所述星载高动态GNSS接收机的导航方法，包括以下步骤：

步骤1，射频前端处理模块对卫星信号进行放大、变频、滤波和模数转换，最终得到数字中频信号；

步骤2，基带信号数字处理模块对数字中频信号进行捕获、跟踪、位同步和帧同步处理，得到导航测量值和导航电文；

步骤3，定位解算模块利用导航测量值和导航电文进行定位解算，最终得到用户的卫星信息，获取定位观测值。

本发明与现有技术相比，其显著优点是：(1)采用ARM+FPGA结构设计，成本较低，并且能很方便的对算法进行改进和更新，适合工程上使用；(2)为通用型GNSS接收机，可以兼容GPS信号与北斗信号，兼容性好；(3)采用合适的器件，在体积、功耗、性能之间取得最佳平衡，搭建符合实际应用要求的样机，适合微纳卫星上使用。

附图说明

图1是本发明星载高动态GNSS接收机及其导航方法的信号处理流程图。

图2是本发明星载高动态GNSS接收机的系统硬件结构图。

图3是本发明星载高动态GNSS接收机的射频前端处理模块的结构图。

图4是本发明星载高动态GNSS接收机及其导航方法的捕获电路框图。

图5是本发明星载高动态GNSS接收机及其导航方法的跟踪环路框图。

图6是本发明星载高动态GNSS接收机及其导航方法的二阶环路数字滤波器方框图。

图7是本发明星载高动态GNSS接收机及其导航方法的由二阶锁频环辅助的三阶锁相环的滤波器方框图。

具体实施方式

下面结合附图及具体实施例对本发明作进一步详细说明。

结合图1～3，本发明星载高动态GNSS接收机，该接收机基于ARM+FPGA，包括顺次相连的射频前端处理模块、基带信号数字处理模块和定位解算模块：射频前端处理模块，对卫星信号进行放大、变频、滤波和模数转换，最终得到数字中频信号；基带信号数字处理模块，对数字中频信号进行捕获、跟踪、位同步和帧同步处理，得到导航测量值和导航电文；定位解算模块，利用导航测量值和导航电文进行定位解算，最终得到用户的卫星信息，获取定位观测值。

优选地，所述射频前端处理模块采用型号为MAX2769的GNSS接收机芯片。所述基带信号数字处理模块采用Altera公司的FPGA芯片EP3C55F484CBN。所述定位解算模块采用ARM芯片STM32F407VET6，主频168MHz，内存192KB，Cortex-M4F核带硬件除法和浮点计算。如图2所示，具体硬件说明如下：

1、射频前端处理模块选用MAX2769，该芯片提供2种配置方式，一种是通过3线SPI接口，为8至10号引脚，通过CS，SCLK，SDATA按照一定的时序对片内寄存器进行配置，这种方式需要外部主机编程来实现，但是灵活性很好。另外一种方式是将8至10引脚作为配置选择引脚，通过对各引脚进行拉高和拉低，来选择8种既定的典型配置，此方式实现方便，但是灵活性较差。8至10引脚是作为SPI串口输入还是作为配置选择引脚是由26脚PGM接不同电平来决定的。本发明采用第一种方案，用FPGA来设计一个SPI接口，并通过此接口向MAX2769写入控制字来实现MAX2769的配置。

2、基带信号数字处理模块采用Altera公司的FPGA芯片EP3C55F484CBN，通过SPI总线接收MAX2769产生的数字中频信号，并提供给通道相关器处理，累加器锁存I/Q信号后触发累加中断；TIC锁存器锁存相关量并触发TIC中断，同时输出PPS秒脉冲。

FPGA在导航系统中主要用于GNSS信号处理与控制。其中用于信号处理的相关器是参照GP2021设计的。FPGA中目前实现了32个通道，通道中的GPS相关器的资源占用为750LE/通道，在不倍频纯使用LE来实现相关器的同类代码中，已经非常优化。

3、定位解算模块采用ARM芯片STM32F407VET6，主频168MHz，内存192KB，Cortex-M4F核带硬件除法和浮点计算。

4、UART接口电路。本发明中ARM对外有两个串口，由于工作电压的不同，需要配置电压转换芯片才能与标准串行接口设备连接。MAX232芯片是美信(MAXIM)公司专为RS-232标准串口设计的单电源电平转换芯片，使用+5v单电源供电；其低功耗关断模式可以将功耗减小到5uW以内，器件特别适合电池供电系统。

5、电源电路。本发明的电压设计采用5V电压输入。经线性稳压源ams1117-3.3芯片稳压后，提供系统所需的3.3V电源。

结合图1～2，本发明为星载高动态GNSS接收机的导航方法，包括以下步骤：

步骤1，射频前端处理模块对卫星信号进行放大、变频、滤波和模数转换，最终得到数字中频信号；结构如图3所示，具体如下：

射频前端处理模块接收可见卫星的信号，将该信号经过放大后，再与射频模块中的振荡器产生的正弦波本振信号进行混频而下变频成中频信号，对该中频信号进行滤波放大，最后由AD芯片采样，将中频信号离散化成数字中频信号。

步骤2，基带信号数字处理模块对数字中频信号进行捕获、跟踪、位同步和帧同步处理，得到导航测量值和导航电文；基带信号数字处理包括捕获、跟踪、位同步和帧同步这4个过程，具体如下：

(2.1)基带信号数字处理模块对数字中频信号进行捕获

GNSS接收机的基带信号数字处理模块复制出一定频率的载波和某颗卫星一定码相位的伪码，然后将其与GNSS中频信号相乘，并对相乘结果进行时间为T_coh的积分。积分以后同相和正交支路上的结果的表达式分别为：

Q(n)＝aD(n)R(τ)sinc(f_eT_coh)sinφ_e

I(n)＝aD(n)R(τ)sinc(f_eT_coh)cosφ_e

式中，a为信号的幅值，D(n)为数据码，R(τ)为中频信号中的C/A码与本地接收机复制出的伪码之间相关性，sinc为辛格函数，f_e为中频信号的中心频率与本地接收机复制出的载波的频差，T_coh为相干积分时间，φ_e为中频信号与本地接收机复制出的载波之间的初始相位差。

积分以后再对两路结果取平方和，其结果为：

V²＝I²(n)+Q²(n)＝(aR(τ)sinc(f_eT_coh))²

式中，V为非相干积分幅值，I、Q分别为I、Q路信号。

由上式可知，在GNSS接收机的基带信号数字处理模块复制的载波与GNSS信号频率一致且本地复制出的伪码与GNSS信号中的信号一致时，检测量V会出现一个峰值。因此，信号的捕获是一个3维搜索的过程，捕获控制模块不断的改变3个参量：卫星号、载波频率、当前卫星标号下的伪码的码相位，以此来产生各种频率的本地载波和各颗卫星的不同码相位的本地伪码。将其与GNSS中频输入信号进行相关并取平方和后，对其大小进行判定。信号捕获电路如图4所示。

本发明采用的三维搜索的方案如下：先锁定某颗卫星的编号和载波频点，对该卫星的伪码的各个码相位进行遍历，然后再更新载波频点，在新的频点下对码相位进行遍历，如此反复，直到遍历完所有的载波频点后再更新卫星的编号。本发明伪码的搜索步进为半码相位，载波频率的搜索步进为400HZ。

本发明采用的捕获判定方法如下：在每一次频点下进行一次捕获判定，读取不同半码相位的2046*2个相关值，并计算得到2046个非相干积分幅值，找出这一组非相干积分幅值的最大值并计算这一组非相干积分幅值的平均值，算出峰均比，若大于某门限则认为捕获成功，否则更新频点重新进行搜索。门限值计算公式为：

$V_t={\mathrm{\σ}}_n\sqrt{-2{\mathrm{lnP}}_{fa}}$

式中，V_t为门限值，σ_n为给定噪声信号功率，P_fa为虚警率。本发明门限值设置为1.44×10⁶。

(2.2)跟踪处理

本发明采用的跟踪环由跟踪伪码的码跟踪环和跟踪载波的载波跟踪环组成。整个跟踪环路的结构框图如图5所示。

(2.2.1)码跟踪环

GNSS接收机基带信号处理模块中的码跟踪环是一种延时锁定环。伪码发生器通过延时电路，复制出三个伪码序列，其中超前码(E码)和滞后码(L码)分别是即时码(P码)超前和滞后半个码片得到的扩频码，将其与中频信号相关以后，码环通过对E和L两路的相关结果进行鉴相和滤波，并将滤波结果反馈到伪码发生器的NCO(NumericalControlledOscillator，数控振荡器)中。以此来实现本地伪码与接收信号中的伪码的完全对齐。

玛跟踪环路采用延迟锁定环鉴别器进行检相，超前、及时、滞后支路上的相关值幅值分别为：

$E=\sqrt{I_E^2+Q_E^2}$

$P=\sqrt{I_P^2+Q_P^2}$

$L=\sqrt{I_L^2+Q_L^2}$

式中，E、P、L分别为超前、即时、滞后支路上的相关值幅值，I_E、I_P、I_L分别为超前、即时、滞后支路上的I路信号，Q_E、Q_P、Q_L分别为超前、即时、滞后支路上的Q路信号。

伪码自相关函数主峰是一个对称的三角形，如果即时码与接收信号伪码相位上保持一致，E和L相等；如果即时码与接收信号伪码相位不一致，则E和L不等，根据两者之间的差异可以反映出即时码和接收信号扩频码的相位差值。本发明采用的鉴相器为非相干超前减滞后幅值法：

${\mathrm{\δ}}_{cp}=\frac{1}{2}\frac{E-L}{E+L}$

式中，δ_cp为码相位差异，E、L分别为超前、滞后支路上的相关值幅值。

然后将鉴相结果进行滤波，码跟踪环使用2阶环路数字滤波器，如图6所示，其传递函数为：

$F\left(s\right)=\frac{1}{K}(a_2{\mathrm{\ω}}_n+\frac{{\mathrm{\ω}}_n^2}{s})$

式中，K为环路增益，a₂为滤波参数，ω_n为特征频率。

则码跟踪环的系统函数为

$H\left(s\right)=\frac{KF\left(s\right)}{s+KF\left(s\right)}=\frac{a_2{\mathrm{\ω}}_{n}s+{\mathrm{\ω}}_n^2}{s^2+a_2{\mathrm{\ω}}_{n}s+{\mathrm{\ω}}_n^2}$

其中各个参数的值为a₂＝1.414，滤波器中的ω_n是由其对应滤波器的噪声带宽决定的，其中B_L＝0.53ω_n。本发明中码环滤波器噪声带宽为1Hz。

(2.2.2)载波跟踪环

载波跟踪环由锁频环和锁相环两种环路构成，载波跟踪环路是对相关结果进行鉴相、鉴频和滤波，并将滤波结果反馈给载波NCO，以此来最终锁定载波的频率和相位。

(a)锁相环：跟踪环路P路就是与即时的伪码进行相关，所以相关以后的结果与捕获处理的结果相同。

$\frac{Q}{I}=\frac{aD\left(n\right)R\left(\mathrm{\τ}\right)\sin c\left(f_e{T}_{coh}\right){\mathrm{sin\φ}}_e}{aD\left(n\right)R\left(\mathrm{\τ}\right)\sin c\left(f_e{T}_{coh}\right){\mathrm{cos\φ}}_e}=\tan \left({\mathrm{\φ}}_e\right)$

载波跟踪环路的鉴相器的鉴相公式为：

${\mathrm{\φ}}_e={\tan }^{-1}\left(\frac{Q}{I}\right)$

式中，φ_e为鉴相器输出误差，Q、I分别为Q、I路信号。

最后将鉴相结果进行滤波后，反馈给载波NCO以此调整载波环路。

(b)锁频环：由于捕获时频率搜索步进为400Hz，得到的载波频率估计值与实际值相差较大，而锁相环虽然精度高但是牵引速度慢，因此需要锁频环路来快速的跟踪信号。锁频环主要是对信号的频率进行锁定，并且其动态范围比所锁相环更宽，能够快速锁定输入信号。

锁频环路的鉴频除了用到了P路当前时刻的相关结果，还用到了上一次P路的相关结果，若设aD(n)R(τ)sinc(f_eT_coh)为A(n)，则：

$I\left(n\right)=\frac{1}{2}A\left(n\right)cos\left({\mathrm{\φ}}_e\right(n\left)\right)$

$Q\left(n\right)=\frac{1}{2}A\left(n\right)sin\left({\mathrm{\φ}}_e\right(n\left)\right)$

其中φ_e(n)是第n个历元的本地载波与输入信号的相位差，定义点乘P_dot和叉乘P_cross分别为：

$\begin{array}{c}{P}_{dot}=I(n-1)I\left(n\right)+Q(n-1)Q\left(n\right)\\ =\frac{1}{4}A(n-1)A\left(n\right)\left(\cos \right({\mathrm{\φ}}_e\left(n\right)\left)\cos \right({\mathrm{\φ}}_e\left(n-1\right))+\sin ({\mathrm{\φ}}_e\left(n\right))\sin \left({\mathrm{\φ}}_e\right(n-1\left)\right))\\ \frac{1}{4}A(n-1)A\left(n\right)\left({\mathrm{\φ}}_e\right(n)-{\mathrm{\φ}}_e(n-1\left)\right)\end{array}$

$\begin{array}{c}{P}_{cross}=I(n-1)Q\left(n\right)+Q(n-1)I\left(n\right)\\ =\frac{1}{4}A(n-1)A\left(n\right)\left(\sin \right({\mathrm{\φ}}_e\left(n\right)\left)\cos \right({\mathrm{\φ}}_e\left(n-1\right))-\cos ({\mathrm{\φ}}_e\left(n\right))\sin \left({\mathrm{\φ}}_e\right(n-1\left)\right))\\ \frac{1}{4}A(n-1)A\left(n\right)\sin \left({\mathrm{\φ}}_e\right(n)-{\mathrm{\φ}}_e(n-1\left)\right)\end{array}$

因此，本发明采用的载波跟踪环路锁频环的鉴频器的鉴频公式为：

${\mathrm{\φ}}_e\left(n\right)=\frac{arctan2(P_{cross},P_{dot})}{t\left(n\right)-t(n-1)}$

式中，φ_e(n)为鉴频器输出误差，t为采样时间。

本发明使用锁频环辅助锁相环的载波跟踪环路结构。在捕获刚成功并切换到跟踪状态的时候，由于捕获得到的粗略值精度不够，锁频环在载波跟踪环中起主导作用，它将快速的牵入信号，当本地载波的频率与GNSS信号的频率较为接近后，锁相环在载波跟踪环中起主导作用，使本地信号更加精确的与输入信号同步。载波跟踪环使用二阶锁频环辅助三阶锁相环的滤波器，如图7所示。其中二阶锁频环的传递函数为：

$F\left(s\right)=\frac{1}{K}(a_2{\mathrm{\ω}}_n+\frac{{\mathrm{\ω}}_n^2}{s})$

式中，K为环路增益，a₂为滤波参数，ω_n为特征频率。

三阶锁相环的传递函数为：

$F\left(s\right)=\frac{1}{K}(b_3{\mathrm{\ω}}_n+\frac{a_3{\mathrm{\ω}}_n^2}{s}+\frac{{\mathrm{\ω}}_n^3}{s^2})$

式中，K为环路增益，a₃、b₃为滤波参数，ω_n为特征频率。

其中各个参数的值为a₃＝1.1，b₃＝2.4，滤波器中的ω_n是由其对应滤波器的噪声带宽决定的，其中B_L＝0.7845ω_n。本发明中载波环路滤波器噪声带宽为18Hz。

(2.3)位同步处理

若跟踪成功，即本地的载波与伪码分别与GNSS中频信号中的载波和伪码完全对齐，此时P路中的同相支路的相关结果包含了导航比特的全部信息。但是跟踪环的积分值是每一个伪码周期即1ms输出一个，而一个导航比特持续时长20ms。跟踪以后还需要找出导航比特的跳变边界，以此边界来将每20个跟踪环输出结果组成一个导航比特。

本发明采用的位同步算法如下：将跟踪环的输出结果进行循环编号，以T为周期，T为正整数20，并建立20个计数器存储两个相邻跟踪结果的跳变次数，如编号2的跟踪结果与编号3的跟踪结果发生信号跳变，则计数器2的计数值加1。一共设立2个门限值N₁和N₂，其中N₁＞N₂，当20个计数器中有一个的值到达N₁时，检查其他19个计数器，若其值全部都小于N₂，则认定为位同步成功，值为N₁的计数器所代表的跳变边界就是导航比特的跳变边界，否则认定同步失败，清零所有计数器重新进行位同步。

(2.4)帧同步处理

帧同步的目的有两点，第一是找到每个子帧的起始位置，以正确划分导航电文中的30比特长度的字；第二是确定是否存在由于180度相位模糊导致的导航比特反相。

本发明中使用的帧同步算法如下：

(2.4.1)将导航比特与子帧的同步码做相关运算，当相关结果为8或-8时，认为找到帧同步码的起始位置，进入步骤(2.4.2)；

(2.4.2)按步骤(2.4.1)中得到的起始位置将导航比特流进行划分，每30比特组成一个字，并对每个字进行奇偶校验，若正确，则进入步骤(2.4.3)，否则返回步骤(2.4.1)；

(2.4.3)提取导航数据中的周内时、子帧号与上一次操作在该步骤中获得的对应值比较，正确则进入步骤(2.4.4)，否则返回步骤(2.4.1)；

(2.4.4)对步骤(2.4.1)中找到的起始位置的后一个子帧的前8个比特进行判断，如果依然为子帧同步头则重复步骤(2.4.2)～(2.4.4)，否则重返回步骤(2.4.1)。

步骤3，定位解算模块利用导航测量值和导航电文进行定位解算，最终得到用户的卫星信息，获取定位观测值，具体为：

(3.1)伪距的测量

伪距采用到达时间测距原理，即用卫星信号的发射时间减去接收机接收时间最后乘以光速得到，卫星发射时间通过下式求出：

$t^s=TOW+(30w+b)\×0.020+(c+\frac{CP}{1023})\×0.001\left(s\right)$

式中，t为卫星发射时间，TOW为交接字中的周内时，w表示在当前子帧中接收机已经接收到完整字的个数，b表示在当前字中已经接收到完整的数据比特数量，c表示在当前比特中已经接收的完整伪码周期个数，CP表示当前时刻的伪码相位；

(3.2)卫星位置解算

卫星是在既定轨道上绕着地球运动的，它的位置是时间t的函数，由星历中包含的卫星开普勒轨道参数，计算出卫星轨道的极坐标方程：

$\left\{\begin{array}{c}r=a_s(1-e_s\cos E)\\ v=arctan\left(\frac{\sqrt{1-e_s^2}\sin E}{\cos E-e_s}\right)\end{array}\right.$

式中，(r，v)为卫星所在位置的极坐标，a_s为卫星轨道的长半径，e_s为卫星轨道偏心率，E为卫星轨道偏近点角；

(3.3)接收机位置解算

若设卫星i的坐标为(x_i,y_i,z_i)，接收机到该卫星的伪距为ρ_i，接收机的坐标(x_u,y_u,z_u)，卫星时钟与接收机本地时钟钟差为δt_u，则有伪距ρ_i的公式：

${\mathrm{\ρ}}_i=\sqrt{{(x_i-x_u)}^2+{(y_i-y_u)}^2+{(z_i-z_u)}^2}+{\mathrm{\δt}}_u$

其中，卫星的位置(x_i,y_i,z_i)和卫星与接收机的伪距为ρ_i都是已知量，可以通过导航电文中的信息求得。接收机的坐标(x_u,y_u,z_u)和钟差δt_u为未知量，若接收机能获取4颗以上的卫星的导航电文，就可以列出四个上述方程，从而解算出接收机的位置。由于方程组是非线性的，本发明采用牛顿迭代及其线性化方法对方程组进行求解，其具体步骤如下：

(3.3.1)设置方程初始解，迭代前给方程组的4个未知数设定一个初始值，初始值的设置分为两种情况：若是首次定位，则全部设为0；若已经成功定位，则将上一次的结果设置为本次迭代的初始值；

(3.3.2)线性化方程组，对伪距ρ_i的公式进行泰勒展开，得：

${\mathrm{\ρ}}_i-{\mathrm{\ρ}}_i(x,y,z,{\mathrm{\δt}}_u)=\frac{\∂{\mathrm{\ρ}}_i}{\∂x}\·\mathrm{\Δ}x+\frac{\∂{\mathrm{\ρ}}_i}{\∂y}\·\mathrm{\Δ}y+\frac{\∂{\mathrm{\ρ}}_i}{\∂z}\·\mathrm{\Δ}z+1\·{\mathrm{\δt}}_u$

式中，Δx、Δy、Δz、Δδt_u为最小二乘法的解；

其中：

$\frac{\∂{\mathrm{\ρ}}_i}{\∂x}=\frac{-(x_i-x)}{\sqrt{{(x_i-x_u)}^2+{(y_i-y_u)}^2+{(z_i-z_u)}^2}}=\frac{-(x_i-x)}{r_i}$

$\frac{\∂{\mathrm{\ρ}}_i}{\∂y}=\frac{-(y_i-y)}{\sqrt{{(x_i-x_u)}^2+{(y_i-y_u)}^2+{(z_i-z_u)}^2}}=\frac{-(y_i-y)}{r_i}$

$\frac{\∂{\mathrm{\ρ}}_i}{\∂z}=\frac{-(z_i-z)}{\sqrt{{(x_i-x_u)}^2+{(y_i-y_u)}^2+{(z_i-z_u)}^2}}=\frac{-(z_i-z)}{r_i}$

将上式写成矩阵形式可得：

$G\left[\begin{array}{c}\mathrm{\Δ}x\\ \mathrm{\Δ}y\\ \mathrm{\Δ}z\\ \mathrm{\Δ}\mathrm{\δ}{t}_u\end{array}\right]=b$

其中

$G=\left[\begin{array}{cccc}\frac{\∂{\mathrm{\ρ}}_1}{\∂x}{|}_{x_{k-1}}& \frac{\∂{\mathrm{\ρ}}_1}{\∂y}{|}_{y_{k-1}}& \frac{\∂{\mathrm{\ρ}}_1}{\∂z}{|}_{z_{k-1}}& 1\\ \frac{\∂{\mathrm{\ρ}}_2}{\∂x}{|}_{x_{k-1}}& \frac{\∂{\mathrm{\ρ}}_2}{\∂y}{|}_{y_{k-1}}& \frac{\∂{\mathrm{\ρ}}_2}{\∂z}{|}_{z_{k-1}}& 1\\ \frac{\∂{\mathrm{\ρ}}_3}{\∂x}{|}_{x_{k-1}}& \frac{\∂{\mathrm{\ρ}}_3}{\∂y}{|}_{y_{k-1}}& \frac{\∂{\mathrm{\ρ}}_3}{\∂z}{|}_{z_{k-1}}& 1\\ \frac{\∂{\mathrm{\ρ}}_4}{\∂x}{|}_{x_{k-1}}& \frac{\∂{\mathrm{\ρ}}_4}{\∂y}{|}_{y_{k-1}}& \frac{\∂{\mathrm{\ρ}}_4}{\∂z}{|}_{z_{k-1}}& 1\end{array}\right]$

$b=\left[\begin{array}{c}{\mathrm{\ρ}}_1-r_1(k-1)-\mathrm{\δ}{t}_{u,k-1}\\ {\mathrm{\ρ}}_2-r_2(k-1)-{\mathrm{\δt}}_{u,k-1}\\ {\mathrm{\ρ}}_3-r_3(k-1)-{\mathrm{\δt}}_{u,k-1}\\ {\mathrm{\ρ}}_4-r_4(k-1)-{\mathrm{\δt}}_{u,k-1}\end{array}\right]$

其中，δt_u,k-1表示第k-1次迭代求出的钟差，r_i(k-1)表示第k-1次迭代求出的接收机与对应卫星的距离，k＝1表示步骤(3.3.1)中设置的初始值；

(3.3.3)利用最小二乘法公式求解方程组：

$\left[\begin{array}{c}\mathrm{\Δ}x\\ \mathrm{\Δ}y\\ \mathrm{\Δ}z\\ \mathrm{\Δ}\mathrm{\δ}{t}_u\end{array}\right]={\left(G^{T}G\right)}^{-1}{G}^{T}b$

(3.3.4)更新非线性方程组的根：

$\left[\begin{array}{c}{x}_k\\ y_k\\ z_k\\ {\mathrm{\δt}}_{u,k}\end{array}\right]=\left[\begin{array}{c}{x}_{k-1}\\ y_{k-1}\\ z_{k-1}\\ {\mathrm{\δt}}_{u,k-1}\end{array}\right]+\left[\begin{array}{c}\mathrm{\Δ}x\\ \mathrm{\Δ}y\\ \mathrm{\Δ}z\\ \mathrm{\Δ}\mathrm{\δ}{t}_u\end{array}\right]$

(3.3.5)判断牛顿迭代收敛性：每次迭代，步骤(3.3.3)的结果会逐渐减小，当矢量长度值小于门限的时候，说明方程组的解已经收敛，则停止迭代，否则重返(3.3.2)；本发明判断方式为检查此次计算得到的位移向量Δx的长度||Δx||是否小于预先设定的门限值0.001；最后一次迭代步骤(3.3.4)的值即为接收机的位置坐标和时钟钟差。一般情况下3至5次迭代即可收敛。

需要说明的是，在定位计算前，普通接收机需要设置卫星仰角滤角，仰角滤角是一个门限值，任何低于这个滤角值的卫星都将被“过滤”掉而不用于定位计算中。一般来说低仰角卫星信号的大气延时校正误差可能很大，并且它的多路径效应又可能很严重，因而通常认为低仰角卫星对改善定位精度的益处抵不上它所带来的较大测量误差和定位误差的坏处。然而对于星载GNSS接收机来说，低仰角的卫星信号仍能参与计算，且接收机捕获到的卫星仰角可以低于水平线。综合考虑对流层和电离层的高度，本发明星载GNSS接收机捕获处理中的仰角滤角设置为-20°。

综上，本发明采用ARM+FPGA的架构，具有很强的灵活性，能够加载不同的算法，且FPGA的并行处理结构能够保证GNSS接收机的实时性。本发明将基带处理中的计算复杂度低但是处理频率较高的相关操作划分给FPGA，将后续处理频率较低但是计算复杂度高的运算划分给ARM。

Claims (9)

1.一种星载高动态GNSS接收机，其特征在于，该接收机基于ARM+FPGA，包括顺次相连的射频前端处理模块、基带信号数字处理模块和定位解算模块：

射频前端处理模块，对卫星信号进行放大、变频、滤波和模数转换，最终得到数字中频信号；

基带信号数字处理模块，对数字中频信号进行捕获、跟踪、位同步和帧同步处理，得到导航测量值和导航电文；

定位解算模块，利用导航测量值和导航电文进行定位解算，最终得到用户的卫星信息，获取定位观测值。

2.根据权利要求1所述的星载高动态GNSS接收机，其特征在于，所述射频前端处理模块采用型号为MAX2769的GNSS接收机芯片。

3.根据权利要求1所述的星载高动态GNSS接收机，其特征在于，所述基带信号数字处理模块采用Altera公司的FPGA芯片EP3C55F484CBN。

4.根据权利要求1所述的星载高动态GNSS接收机，其特征在于，所述定位解算模块采用ARM芯片STM32F407VET6，主频168MHz，内存192KB，Cortex-M4F核带硬件除法和浮点计算。

5.一种基于权利要求1所述星载高动态GNSS接收机的导航方法，其特征在于，包括以下步骤：

步骤1，射频前端处理模块对卫星信号进行放大、变频、滤波和模数转换，最终得到数字中频信号；

步骤2，基带信号数字处理模块对数字中频信号进行捕获、跟踪、位同步和帧同步处理，得到导航测量值和导航电文；

步骤3，定位解算模块利用导航测量值和导航电文进行定位解算，最终得到用户的卫星信息，获取定位观测值。

6.根据权利要求5所述的星载高动态GNSS接收机的导航方法，其特征在于，步骤1所述射频前端处理模块对卫星信号进行放大、变频、滤波和模数转换，最终得到数字中频信号，具体为：射频前端处理模块接收可见卫星的信号，将该信号经过放大后，再与射频模块中的振荡器产生的正弦波本振信号进行混频而下变频成中频信号，对该中频信号进行滤波放大，最后由AD芯片采样，将中频信号离散化成数字中频信号。

7.根据权利要求5所述的星载高动态GNSS接收机的导航方法，其特征在于，步骤2所述基带信号数字处理模块对数字中频信号进行捕获、跟踪、位同步和帧同步处理，得到导航测量值和导航电文，具体为：

(2.1)基带信号数字处理模块对数字中频信号进行捕获

采用三维搜索，具体如下：先锁定一颗卫星的编号和载波频点，对该卫星的伪码的各个码相位进行遍历，然后再更新载波频点，在新的频点下对码相位进行遍历，如此反复，直到遍历完所有的载波频点后再更新卫星的编号；

采用的捕获判定方法如下：在每一次频点下进行一次捕获判定，读取不同半码相位的2046*2个相关值，并计算得到2046个非相干积分幅值，找出该组非相干积分幅值的最大值并计算该组非相干积分幅值的平均值，算出峰均比，若大于门限值则认为捕获成功，否则更新频点重新进行搜索，其中门限值V_t计算公式为：

$V_t={\mathrm{\σ}}_n\sqrt{-2{\mathrm{lnP}}_{fa}}$

式中，V_t为门限值，σ_n为给定噪声信号功率，P_fa为虚警率；

(2.2)跟踪处理

(2.2.1)码跟踪环

采用的鉴相器为非相干超前减滞后幅值法：

${\mathrm{\δ}}_{cp}=\frac{1}{2}\frac{E-L}{E+L}$

式中，δ_cp为码相位差异，E、L分别为超前、滞后支路上的相关值幅值；

然后将鉴相结果进行滤波，码跟踪环使用2阶环路数字滤波器，传递函数F(s)为：

$F\left(s\right)=\frac{1}{K}(a_2{\mathrm{\ω}}_n+\frac{{\mathrm{\ω}}_n^2}{s})$

式中，K为环路增益，a₂为滤波参数，ω_n为特征频率；

则码跟踪环的系统函数为

$H\left(s\right)=\frac{KF\left(s\right)}{s+KF\left(s\right)}=\frac{a_2{\mathrm{\ω}}_{n}s+{\mathrm{\ω}}_n^2}{s^2+a_2{\mathrm{\ω}}_{n}s+{\mathrm{\ω}}_n^2}$

(2.2.2)载波跟踪环

(a)锁相环：本发明采用的载波跟踪环路锁相环的鉴相器的鉴相公式为：

${\mathrm{\φ}}_e={\tan }^{-1}\left(\frac{Q}{I}\right)$

式中，φ_e为鉴相器输出误差，Q、I分别为Q、I路信号；

(b)锁频环：本发明采用的载波跟踪环路锁频环的鉴频器的鉴频公式为：

${\mathrm{\φ}}_e\left(n\right)=\frac{arctan2(P_{cross},P_{dot})}{t\left(n\right)-t(n-1)}$

式中，P_dot＝I(n-1)I(n)+Q(n-1)Q(n)，P_cross＝I(n-1)Q(n)+Q(n-1)I(n)，φ_e(n)为鉴频器输出误差，t为采样时间；

本发明中载波跟踪环使用二阶锁频环辅助三阶锁相环的结构，其中二阶锁频环的传递函数为：

$F\left(s\right)=\frac{1}{K}(a_2{\mathrm{\ω}}_n+\frac{{\mathrm{\ω}}_n^2}{s})$

式中，K为环路增益，a₂为滤波参数，ω_n为特征频率；

三阶锁相环的传递函数为：

$F\left(s\right)=\frac{1}{K}(b_3{\mathrm{\ω}}_n+\frac{a_3{\mathrm{\ω}}_n^2}{s}+\frac{{\mathrm{\ω}}_n^3}{s^2})$

式中，K为环路增益，a₃、b₃为滤波参数，ω_n为特征频率；

(2.3)位同步处理

本发明采用的位同步算法如下：将跟踪环的输出结果进行循环编号，以T为周期，T为正整数，并建立T个计数器存储两个相邻跟踪结果的跳变次数，设立2个门限值N₁和N₂，其中N₁＞N₂，当T个计数器中有一个的值到达N₁时，检查其他T-1个计数器：若值全部都小于N₂，则认定为位同步成功，值为N₁的计数器所代表的跳变边界就是导航比特的跳变边界；否则认定同步失败，清零所有计数器重新进行位同步；

(2.4)帧同步处理

(2.4.1)将导航比特与子帧的同步码做相关运算，当相关结果为8或-8时，认为找到帧同步码的起始位置，进入步骤(2.4.2)；

(2.4.2)按步骤(2.4.1)中得到的起始位置将导航比特流进行划分，每30比特组成一个字，并对每个字进行奇偶校验，若正确，则进入步骤(2.4.3)，否则返回步骤(2.4.1)；

(2.4.3)提取导航数据中的周内时、子帧号与上一次操作在该步骤中获得的对应值比较，正确则进入步骤(2.4.4)，否则返回步骤(2.4.1)；

(2.4.4)对步骤(2.4.1)中找到的起始位置的后一个子帧的前8个比特进行判断，如果依然为子帧同步头则重复步骤(2.4.2)～(2.4.4)，否则重返回步骤(2.4.1)。

8.根据权利要求5所述的星载高动态GNSS接收机的导航方法，其特征在于，步骤3所述定位解算模块利用导航测量值和导航电文进行定位解算，最终得到用户的卫星信息，获取定位观测值，具体为：

(3.1)伪距的测量

伪距采用到达时间测距原理，即用卫星信号的发射时间减去接收机接收时间最后乘以光速得到，卫星发射时间通过下式求出：

$t^s=TOW+(30w+b)\×0.020+(c+\frac{CP}{1023})\×0.001\left(s\right)$

式中，t为卫星发射时间，TOW为交接字中的周内时，w表示在当前子帧中接收机已经接收到完整字的个数，b表示在当前字中已经接收到完整的数据比特数量，c表示在当前比特中已经接收的完整伪码周期个数，CP表示当前时刻的伪码相位；

(3.2)卫星位置解算

卫星是在既定轨道上绕着地球运动的，它的位置是时间t的函数，由星历中包含的卫星开普勒轨道参数，计算出卫星轨道的极坐标方程：

$\left\{\begin{array}{c}r=a_s(1-e_s\cos E)\\ v=arctan\left(\frac{\sqrt{1-e_s^2}\sin E}{\cos E-e_s}\right)\end{array}\right.$

式中，(r，v)为卫星所在位置的极坐标，a_s为卫星轨道的长半径，e_s为卫星轨道偏心率，E为卫星轨道偏近点角；

(3.3)接收机位置解算

若设卫星i的坐标为(x_i,y_i,z_i)，接收机到该卫星的伪距为ρ_i，接收机的坐标(x_u,y_u,z_u)，卫星时钟与接收机本地时钟钟差为δt_u，则有伪距ρ_i的公式：

${\mathrm{\ρ}}_i=\sqrt{{(x_i-x_u)}^2+{(y_i-y_u)}^2+{(z_i-z_u)}^2}+{\mathrm{\δt}}_u$

若接收机能获取4颗以上的卫星的导航电文，则列出四个上述方程得到方程组，从而解算出接收机的位置。

9.根据权利要求8所述的星载高动态GNSS接收机的导航方法，其特征在于，步骤(3.3)所述方程组是非线性的，采用牛顿迭代及线性化方法对方程组进行求解，具体步骤如下：

(3.3.1)设置方程初始解，迭代前给方程组的4个未知数设定一个初始值，初始值的设置分为两种情况：若是首次定位，则全部设为0；若已经成功定位，则将上一次的结果设置为本次迭代的初始值；

(3.3.2)线性化方程组，对伪距ρ_i的公式进行泰勒展开，得：

${\mathrm{\ρ}}_i-{\mathrm{\ρ}}_i(x,y,z,{\mathrm{\δt}}_u)=\frac{\∂{\mathrm{\ρ}}_i}{\∂x}\·\mathrm{\Δ}x+\frac{\∂{\mathrm{\ρ}}_i}{\∂y}\·\mathrm{\Δ}y+\frac{\∂{\mathrm{\ρ}}_i}{\∂z}\·\mathrm{\Δ}z+1\·{\mathrm{\Δ\δt}}_u$

式中，Δx、Δy、Δz、Δδt_u为最小二乘法的解；

其中：

$\frac{\∂{\mathrm{\ρ}}_i}{\∂x}=\frac{-(x_i-x)}{\sqrt{{(x_i-x_u)}^2+{(y_i-y_u)}^2+{(z_i-z_u)}^2}}=\frac{-(x_i-x)}{r_i}$

$\frac{\∂{\mathrm{\ρ}}_i}{\∂y}=\frac{-(y_i-y)}{\sqrt{{(x_i-x_u)}^2+{(y_i-y_u)}^2+{(z_i-z_u)}^2}}=\frac{-(y_i-y)}{r_i}$

$\frac{\∂{\mathrm{\ρ}}_i}{\∂z}=\frac{-(z_i-z)}{\sqrt{{(x_i-x_u)}^2+{(y_i-y_u)}^2+{(z_i-z_u)}^2}}=\frac{-(z_i-z)}{r_i}$

将上式写成矩阵形式得：

$G\left[\begin{array}{c}\mathrm{\Δ}x\\ \mathrm{\Δ}y\\ \mathrm{\Δ}z\\ \mathrm{\Δ}\mathrm{\δ}{t}_u\end{array}\right]=b$

其中：

$G=\left[\begin{array}{cccc}\frac{\∂{\mathrm{\ρ}}_1}{\∂x}{|}_{x_{k-1}}& \frac{\∂{\mathrm{\ρ}}_1}{\∂y}{|}_{y_{k-1}}& \frac{\∂{\mathrm{\ρ}}_1}{\∂z}{|}_{z_{k-1}}& 1\\ \frac{\∂{\mathrm{\ρ}}_2}{\∂x}{|}_{x_{k-1}}& \frac{\∂{\mathrm{\ρ}}_2}{\∂y}{|}_{y_{k-1}}& \frac{\∂{\mathrm{\ρ}}_2}{\∂z}{|}_{z_{k-1}}& 1\\ \frac{\∂{\mathrm{\ρ}}_3}{\∂x}{|}_{x_{k-1}}& \frac{\∂{\mathrm{\ρ}}_3}{\∂y}{|}_{y_{k-1}}& \frac{\∂{\mathrm{\ρ}}_3}{\∂z}{|}_{z_{k-1}}& 1\\ \frac{\∂{\mathrm{\ρ}}_4}{\∂x}{|}_{x_{k-1}}& \frac{\∂{\mathrm{\ρ}}_4}{\∂y}{|}_{y_{k-1}}& \frac{\∂{\mathrm{\ρ}}_4}{\∂z}{|}_{z_{k-1}}& 1\end{array}\right]$

$b=\left[\begin{array}{c}{\mathrm{\ρ}}_1-r_1(k-1)-\mathrm{\δ}{t}_{u,k-1}\\ {\mathrm{\ρ}}_2-r_2(k-1)-{\mathrm{\δt}}_{u,k-1}\\ {\mathrm{\ρ}}_3-r_3(k-1)-{\mathrm{\δt}}_{u,k-1}\\ {\mathrm{\ρ}}_4-r_4(k-1)-{\mathrm{\δt}}_{u,k-1}\end{array}\right]$

其中，δt_u,k-1表示第k-1次迭代求出的钟差，r_i(k-1)表示第k-1次迭代求出的接收机与对应卫星的距离，k＝1表示步骤(3.3.1)中设置的初始值；

(3.3.3)利用最小二乘法公式求解方程组：

$\left[\begin{array}{c}\mathrm{\Δ}x\\ \mathrm{\Δ}y\\ \mathrm{\Δ}z\\ \mathrm{\Δ}\mathrm{\δ}{t}_u\end{array}\right]={\left(G^{T}G\right)}^{-1}{G}^{T}b$

(3.3.4)更新非线性方程组的根：

$\left[\begin{array}{c}{x}_k\\ y_k\\ z_k\\ {\mathrm{\δt}}_{u,k}\end{array}\right]=\left[\begin{array}{c}{x}_{k-1}\\ y_{k-1}\\ z_{k-1}\\ {\mathrm{\δt}}_{u,k-1}\end{array}\right]+\left[\begin{array}{c}\mathrm{\Δ}x\\ \mathrm{\Δ}y\\ \mathrm{\Δ}z\\ \mathrm{\Δ}\mathrm{\δ}{t}_u\end{array}\right]$

(3.3.5)判断牛顿迭代收敛性：当矢量长度值小于门限的时候，说明方程组的解已经收敛，则停止迭代，否则重返(3.3.2)；最后一次迭代步骤(3.3.4)的值即为接收机的位置坐标和时钟钟差。