低轨卫星（LEO）定位之多普勒定位原理及推导

shannanshouyin

已于 2025-03-22 16:38:18 修改

阅读量1.1k

点赞数 15

分类专栏：卫星导航定位文章标签：信号处理算法矩阵线性代数

于 2025-03-19 14:59:07 首次发布

本文链接：https://blog.csdn.net/shannanshouyin/article/details/146368697

版权

卫星导航定位专栏收录该内容

3 篇文章

订阅专栏

在低轨卫星（LEO）定位时，采用多普勒定位方法具有显著优势，主要体现在以下几个方面：

1、定位精度高：低轨卫星离地面较近，信号传播路径短，信号强度相对较强，多普勒频移现象更为明显。通过精确测量卫星信号的多普勒频移，结合卫星的精确轨道信息以及相关的定位算法，可以实现较高精度的定位。

2、定位相应速度快：多普勒定位方法基于对卫星信号频率变化的实时监测。由于低轨卫星的运动速度相对较快，其与接收机之间的相对运动所产生的多普勒频移变化迅速，接收机能够快速捕捉到这些变化并进行处理。

3、抗干扰能力强：多普勒定位方法在复杂电磁环境和恶劣天气条件下具有较强的抗干扰能力。其定位原理主要依赖于对卫星信号频率特征的分析，而不是像一些其他定位方法那样对信号的幅度、相位等特征有高度依赖。因此，在面对电磁干扰、多径效应以及恶劣天气导致的信号衰减等情况时，多普勒定位方法能够通过对频率信息的稳健提取和处理，较为准确地计算出目标位置，保证定位系统的可靠性和稳定性，确保在各种复杂环境下都能为用户提供连续、准确的定位服务。

本文将详细介绍多普勒定位的方法并进行公式推导。

一、低轨卫星多普勒频移介绍

在LEO卫星运行中，LEO运动速度非常大，可达几十千米每秒，因此地面接收机与运动的LEO卫星会存在径向的相对运动，接收到的信号频率与真实信号频率不再相同，这种现象叫做多普勒效应，引起的频率改变称为多普勒频移。卫星和接收机在地心地固坐标系（ECEF）下，对于多普勒频移有如下公式：

$f_d = f_r - f_0 = \frac{\mathbf{v_r - v^s}}{c}\cdot\frac{\mathbf{x^s - x_u}}{\|\mathbf{x^s - x_u}\|}\cdot f^s\, \, \, \, \, \, \, \, (1)$

式中， $c$ 表示光速， $f_d$ 为多普勒频移量， $f_0$ 为地面接收机接收到的信号频率， $f^s$ 为卫星发射信号的原始频率， $\mathbf{v}_\mathbf{u}=\left[v_{x\ },v_y,\left.v_z\right]\right.^T$ 、 $\mathbf{v}^\mathbf{s}=\left[v_x^s\ ,v_y^s,\left.v_z^s\right]\right.^T$ 分别为地面接收机运动速度和卫星运动速度； $\mathbf{x}_\mathbf{u}=\ \left[x_u,y_u,\left.z_u\right]\right.^T$ 、 $\mathbf{x}^\mathbf{s}=\left[x^s,y^s,\left.z^s\right]\right.^T$ 分别表示地面接收机和卫星的位置，后续的公式计算推导均在ECEF坐标系下，不再重复说明。

令 $r^s=\|\mathbf{x^s - x_u}\|, \ \mathbf{x^s_u= x^s - x_u} , \ \mathbf{v^s_u = v^s - v_u}$ 可将式（1）简化为：

$f_d = -\frac{f_0}{c}\cdot\frac{\boldsymbol{v}^s_u\cdot\boldsymbol{x}^s_u}{r^s}\, \, \, \, \, \, \, \, (2)$

其中 $r^s$ 表示地面接收机和卫星之间的几何距离， $\mathbf{x}_\mathbf{r}^\mathbf{s}$ 表示二者的相对位置矢量， $\mathbf{v}_\mathbf{r}^\mathbf{s}$ 表示二者的相对速度矢量。

二、多普勒定位模型建立

在三维空间内，当地面接收机可以接收到卫星信号并测量出多普勒频移时，以该卫星位置矢量为顶点，卫星和地面接收机连线为母线构成的圆锥面上任意一点的多普勒频移值相等，该圆锥曲即为等多普勒频率面，其与地球面相交可以得到等多普勒线。如图1所示。

在 $t_1$ 时刻,设卫星位置为 $\mathbf{x}^\mathbf{s}=\left[x^s,y^s,\left.z^s\right]\right.^T$ ，速度为 $\mathbf{v}^\mathbf{s}=\left[v_x^s\ ,v_y^s,\left.v_z^s\right]\right.^T$ ,地面接收机位置为 $\mathbf{x}_\mathbf{u}=\ \left[x,y,\left.z\right]\right.^T$ ，速度为 $\mathbf{v}_\mathbf{u}=\left[v_{x\ },v_y,\left.v_z\right]\right.^T$ ，二者相对位置矢量为 $\mathbf{x}_\mathbf{u}^\mathbf{s}=\mathbf{x}^\mathbf{s}-\mathbf{x}_\mathbf{u}=\left[x^s\right.-x,y^s-y,z^s-\left.z\right]$ ，相对速度矢量为 $\mathbf{v}_\mathbf{u}^\mathbf{s}={\mathbf{v}^\mathbf{s}-\mathbf{v}}_\mathbf{u}=\left[v_x^s-v_{x\ },v_y^s-v_y,v_z^s-v_z\right]=[v_{xu}^s,v_{yu}^s,v_{zu}^s]$ ， ${\dot{\delta}}_{t_u}$ 为接收机的钟漂，在不考虑其他误差时，卫星的多普勒频移观测方程可表示为:

$f_d^s\left(t_1\right)=\frac{f_0}{c}\ \bullet\ \frac{v_{ux}\left(x-x^s\right)+v_{uy}\left(y-y^s\right)+v_{uz}\left(z-z^s\right)}{\sqrt{\left(x^s-x\right)^2+\left(y^s-y\right)^2+\left(z^s-z\right)^2}}\ +\ {\dot{\delta}}_{t_u}\, \, \, \, \, \, \, \, (3)$

对于地面的接收机，卫星的位置、速度可由星历提供，信号频率以及多普勒频移值可由接收机自主测得，不考虑频率测量误的情况下，方程（3）仅存在接收机坐标 $\left[x,y,\left.z\right]\right.^T$ 和接收机钟漂 ${\dot{\delta}}_{t_u}$ 四项未知数,组成了解向量为 $\left[x,y,z,\left.{\dot{\delta}}_{t_u}\right]\right.^T$ 。当累计观测到三颗卫星时，三颗卫星的地面等多普勒线交汇处可确定接收机的唯一位置，如图2所示。

由于存在四个未知量，当累计观测四颗及以上卫星时，即可求解方程组全部未知数。

在累计观测卫星数量n≥4情况下，可列出如下多普勒定位方程组：

$\begin{cases} f_d^1(\boldsymbol{u}) = \frac{f_0}{c} \cdot \frac{v_{ux}^{1}(x - x^1)+v_{uy}^{1}(y - y^1)+v_{uz}^{1}(z - z^1)}{\sqrt{(x^1 - x)^2+(y^1 - y)^2+(z^1 - z)^2}} + \dot{\delta_{t_u}}\\ f_d^2(\boldsymbol{u}) = \frac{f_0}{c} \cdot \frac{v_{ux}^{2}(x - x^2)+v_{uy}^{2}(y - y^2)+v_{uz}^{2}(z - z^2)}{\sqrt{(x^2 - x)^2+(y^2 - y)^2+(z^2 - z)^2}} + \dot{\delta_{t_u}}\\ \qquad\qquad\vdots\\ f_d^n(\boldsymbol{u}) = \frac{f_0}{c} \cdot \frac{v_{ux}^{n}(x - x^n)+v_{uy}^{n}(y - y^n)+v_{uz}^{n}(z - z^n)}{\sqrt{(x^n - x)^2+(y^n - y)^2+(z^n - z)^2}} + \dot{\delta_{t_u}} \end{cases}\, \, \, \, \, \, \, \, (4)$

三、多普勒定位求解

在解算多普勒定位方程组时，需要先给用户一个初始估计解，再利用牛顿迭代法进行解算，最后将收敛值作为最终定位结果。可以选取地球质心 $[0,0,0]$ 或者附近的参考站坐标作为解的坐标部分，钟漂部分则取0即可。

对于多普勒定位方程，假设用户初始坐标为 $\mathbf{x}_\mathbf{0}=[x_0,y_0,z_0]^T$ ,接收机初始钟漂为 ${\dot{\delta}}_{t_0}$ ,组成初始解为 $\mathbf{u}_\mathbf{0}=[x_0,y_0,z_0,\dot{\delta _{t_0}}]^T$ 可将式（3）进行泰勒一阶展开线性化为下式：

$f_d^s\left(\mathbf{u}\right)=f_d^s\left(\mathbf{u}_\mathbf{0}\right)+\left.\frac{\partial f}{\partial x}\right|_{x=x_0}\left(x-x_0\right)+\left.\frac{\partial f}{\partial y}\right|_{y=y_0}\left(y-y_0\right)\\+\left.\frac{\partial f}{\partial z}\right|_{z=z_0}\left(z-z_0\right)+{\dot{\delta}}_{t_u}-{\dot{\delta}}_{t_0}\, \, \, \, \, \, \, \, (5)$

其中，

$\left.\frac{\partial f}{\partial x}\right|_{x = x_0} = \frac{f_0}{c} \cdot \left[\frac{v_x}{\|\mathbf{x}_0 - \mathbf{x}^s\|} - (x_0 - x^s) \cdot \frac{\mathbf{v}^s \cdot (\mathbf{x}_0 - \mathbf{x}^s)}{\|\mathbf{x}_0 - \mathbf{x}^s\|^3}\right],\\ \left.\frac{\partial f}{\partial y}\right|_{y = y_0} = \frac{f_0}{c} \cdot \left[\frac{v_y}{\|\mathbf{x}_0 - \mathbf{x}^s\|} - (y_0 - y^s) \cdot \frac{\mathbf{v}^s \cdot (\mathbf{x}_0 - \mathbf{x}^s)}{\|\mathbf{x}_0 - \mathbf{x}^s\|^3}\right],\\ \left.\frac{\partial f}{\partial z}\right|_{z = z_0} = \frac{f_0}{c} \cdot \left[\frac{v_z}{\|\mathbf{x}_0 - \mathbf{x}^s\|} - (x_0 - x^s) \cdot \frac{\mathbf{v}^s \cdot (\mathbf{x}_0 - \mathbf{x}^s)}{\|\mathbf{x}_0 - \mathbf{x}^s\|^3}\right]\, \, \, \, \, \, \, \, (6)$

将式（4）每一个方程进行如式（5）的线性变化，可转化为矩阵表达式：

$\mathbf{G}\cdot\Delta\boldsymbol{x}_{k - 1}=\Delta f\, \, \, \, \, \, \, \, (7)$

其中，

$\mathbf{G}= \begin{bmatrix} I_x^{(1)}(\boldsymbol{x}_{k - 1}) & I_y^{(1)}(\boldsymbol{x}_{k - 1}) & I_z^{(1)}(\boldsymbol{x}_{k - 1}) & 1\\ I_x^{(2)}(\boldsymbol{x}_{k - 1}) & I_y^{(2)}(\boldsymbol{x}_{k - 1}) & I_z^{(2)}(\boldsymbol{x}_{k - 1}) & 1\\ \cdots & \cdots & \cdots & \cdots\\ I_x^{(n)}(\boldsymbol{x}_{k - 1}) & I_y^{(n)}(\boldsymbol{x}_{k - 1}) & I_z^{(n)}(\boldsymbol{x}_{k - 1}) & 1 \end{bmatrix}$

$\Delta\boldsymbol{x}_{k - 1}= \begin{bmatrix} x_k - x_{k - 1}\\ y_k - y_{k - 1}\\ z_k - z_{k - 1}\\ \delta_{t_k} - \delta_{t_{k - 1}} \end{bmatrix}$

式（6）中，k为迭代次数， $\mathbf{u}_{\mathbf{k}-\mathbf{1}}$ 为第k-1次迭代解算 $\mathbf{u}$ 的估计值， $\Delta f$ 为接收机接收卫星信号的多普勒频移值与预测值的插值向量， $\mathbf{G}$ 为雅可比矩阵，由式（6）可知其内部参数如下式：

$I_x^{(s)}(\boldsymbol{x}_{k - 1}) = \frac{f_0}{c} \cdot \left[\frac{v_x^s - v_x}{\|\boldsymbol{x}_{k - 1} - \boldsymbol{x}^s\|} - (x_{k - 1} - x^s) \cdot \frac{\boldsymbol{v}_u^s \cdot \boldsymbol{x}_{k - 1}^s}{\|\boldsymbol{x}_{k - 1} - \boldsymbol{x}^s\|^3}\right]\, \, \, \, \, \, \, \, (8)$

由式（8）可知，矩阵 $\mathbf{G}$ 只与各颗卫星相对于接收机的几何位置及速度相关。

对于式（7），可根据最小二乘法原理求解 $\Delta x$ ，求解得：

$\Delta\boldsymbol{x}_{k - 1}=(\mathbf{G}^T\mathbf{G})^{-1}\mathbf{G}^T\cdot\Delta f\, \, \, \, \, \, \, \, \, (9)$

关于线性方程组的求解问题可以看博主的另一篇文章：

最小二乘法求解线性方程组推导（行空间、列空间、投影矩阵）-CSDN博客

将解算出的 $\Delta \mathbf{x_{k-1}}$ 对上一次迭代的估计解 $\mathbf{x}_{\mathbf{k}-\mathbf{1}}$ 进行修正，则更新后的估计解为 $\mathbf{x}_\mathbf{k}=\mathbf{x}_{\mathbf{k}-\mathbf{1}}+\Delta \mathbf{x_{k-1}}$ ，并带入下一次迭代运算，直到 $\Delta \mathbf{x_{k-1}}$ 小于设定的阈值，则退出迭代，将此时的 $\mathbf{x}_\mathbf{k}$ 作为最终的用户三维坐标解及接收机的钟漂解，算法流程图如图3所示。