CN105242285A - 一种基于卫星通信的无人机导航数据被欺骗识别方法 - Google Patents

Abstract

本发明公开了一种基于卫星通信的无人机导航数据被欺骗识别方法，包括步骤一：建立与无人机相关的地心坐标系、机体坐标系统、大地坐标系；步骤二：获得坐标系之间的变换关系；步骤三：根据无人机天线跟踪误差电压跟踪卫星得到卫星在机体坐标下的方位角和俯仰角；步骤四：根据惯性导航测量的机体的方位角、俯仰角、横滚角，计算卫星在飞机所在点大地坐标系统下的大地方位角和俯仰角；步骤五：根据卫星在飞机所在点大地坐标系统下的大地方位角和俯仰角，反推出飞机所在的经度和纬度θ_Lat；本方法利用已有的卫星通信系统就可以进行导航数据被欺骗进行识别，无需要增加其它辅助设备。

Description

一种基于卫星通信的无人机导航数据被欺骗识别方法

技术领域

本发明涉及一种基于卫星通信的无人机导航数据被欺骗识别方法，属于无人机导航抗欺骗技术领域。

背景技术

无人机属于智能设备，由于没有人在无人机上，地面完全靠接收机载的导航数据来判定无人机的具体位置，常用的无人机采用惯性导航+卫星导航相结合的方式进行导航。可以利用惯性导航的不易受干扰和卫星导航的高精度的特点，提供高精度、稳定的位置导航数据。如果机载的导航数据受到欺骗式干扰，如干扰机载的卫星导航数据，当干扰采用渐进式干扰时，即使是采用了组合导航的惯导也无法区分是由于无人机的飞行导致的位置变化，还是由于受到了欺骗干扰导致的位置偏差。导航系统受到干扰后，地面将无法知道飞机的准确位置，可能地面人员看到飞机是往本场着陆飞行，实际上飞机已经远离本场降落到其它机场，在地面上如何识别导航是否被欺骗是必须要解决的问题。

发明内容

本发明的目的是为了解决上述问题，判断在地面识别导航数据是否受到欺骗，从而采取相应的措施，是一种基于卫星通信的无人机导航数据被欺骗识别方法。

一种基于卫星通信的无人机导航数据被欺骗识别方法，包括以下几个步骤：

步骤一：建立与无人机相关的地心坐标系、机体坐标系统、大地坐标系；

步骤二：获得坐标系之间的变换关系；

步骤三：根据无人机天线跟踪误差电压跟踪卫星得到卫星在机体坐标下的方位角和俯仰角；

步骤四：根据惯性导航测量的机体的方位角、俯仰角、横滚角，计算卫星在飞机所在点大地坐标系统下的大地方位角和俯仰角；

步骤五：根据卫星在飞机所在点大地坐标系统下的大地方位角和俯仰角，反推出飞机所在的经度和纬度θ_Lat；

本发明的优点在于：

(1)本方法可以利用已有的卫星通信系统就可以进行导航数据被欺骗进行识别，无需要增加其它辅助设备；

(2)本方法对位置误差进行了分析，并通过实际飞行数据验证：通过卫星通信推算的导航数据平均误差为27.76km，精度满足导航数据欺骗识别要求。

附图说明

图1是本发明的方法流程图；

图2是本发明的地心坐标系示意图；

图3是本发明的大地坐标系示意图；

图4是本发明的机体坐标示意图；

图5是本发明的地心坐标系到大地坐标系之间的变换关系示意图；

图6是本发明的某无人机飞行试验中实际经度和推算经度示意图；

图7是本发明的某无人机飞行试验中实际纬度和推算纬度示意图；

图8是本发明的某无人机飞行试验中总误差比较示意图；

具体实施方式

下面将结合附图和实施例对本发明作进一步的详细说明。

本发明是一种基于卫星通信的无人机导航数据被欺骗识别方法，如图1所示，包括以下几个步骤：

步骤一：建立与无人机相关的地心坐标系、机体坐标系统、大地坐标系；

设飞机的经度和纬度分别为θ_LatJ，卫星的经度地球半径R_e，飞机到地面的高度为h，卫星到地面的高度为H，卫星到飞机的距离为D。

设地心坐标系为O(x_e,y_e,z_e)，如图2所示，其中，O_e原点为地球的中心，x_e轴过地球中心指向真北方向，y_e轴过地球中心指向正东方向，在赤道平面内，z_e轴与x_ey_e轴成右手关系。

设大地坐标系为O(x_g,y_g,z_g)，如图3所示，其中，O原点为地球质心，x_g轴过地球质心水平指向真北方向，y_g轴过地球质心，指向正东方向，z_g轴过地球质心，指向地心方向，与x_gy_g轴成右手关系，其中x_goy_g为水平面，x_goz_g为垂直面。

设机体坐标系为O(x_j,y_j,z_j)，如图4所示，其中，O原点为机体重心，x_j轴过质心的机体滚动轴线，飞行方向为正(横滚轴)；y_j轴在机体机翼水平面内，指向右边，z_j轴在飞机水平对称平面内，与x_jy_j轴成右手关系。

步骤二：获得坐标系之间的变换关系；

机体相对大地坐标系下的角度定义为：

偏航角为机体ox_j轴在大地坐标系x_goy_g水平面的投影与ox_g的夹角，自真北起顺时针方向为正，与ox_j构成右手旋，范围[0,360°]；

俯仰角θ为机体ox_j轴与大地坐标系x_goy_g水平面内的投影之间的夹角，以机头抬高为正，与oy_j构成右手旋，范围[-90°,90°]；

横滚角φ为机体oz_j轴与大地坐标系x_goz_g垂直面之间的夹角，右侧下降为正，与ox_j构成右手旋，范围[-90°,90°]。

图5为本发明涉及的地心坐标系到大地坐标系之间的变换关系示意图。

大地坐标与机体坐标之间的转换关系为：

$x_j=T_g^j{x}_g---\left(1\right)$

其中， $x_j=\left[\begin{array}{c}{x}_j\\ y_j\\ z_j\end{array}\right]$ 表示机体坐标下的一个向量， $x_g=\left[\begin{array}{c}{x}_g\\ y_g\\ z_g\end{array}\right]$ 表示大地坐标坐标下的一个向量，为从大地坐标系到机体坐标系的旋转矩阵，具体为：

机体坐标与大地坐标之间的转换关系为：

$x_g=T_j^g{x}_j---\left(3\right)$

其中：为从机体坐标系到大地坐标系的旋转矩阵，具体为：

飞机所在点大地坐标与地心坐标之间的转换关系为：

$x_g=T_e^g{x}_e---\left(5\right)$

其中： $x_e=\left[\begin{array}{c}{x}_e\\ y_e\\ z_e\end{array}\right]$ 表示地心坐标下的一个向量，飞机所在的点在大地坐标系位置的经度和纬度分别为θ_Lat，地球半径为R，飞机的高度为h，为从地心坐标系到飞机所在点的大地坐标系的旋转矩阵。

步骤三：根据无人机天线跟踪误差电压跟踪卫星得到卫星在机体坐标下的方位角和俯仰角；

设卫星在大地坐标系中的方位、俯仰角为A_gS和E_gS，卫星在机体坐标系中方位、俯仰角为A_jS和E_jS，飞机姿态角为偏航角俯仰角θ、横滚角φ，D为卫星到飞机之间的距离，卫星在机体坐标下的方位角、俯仰角为A_jS和E_jS：

$x_{gS}=\left[\begin{array}{c}D\cos E_{gS}\cos A_{gS}\\ D\cos E_{gS}\sin A_{gS}\\ Dsin{E}_{gS}\end{array}\right]---\left(7\right)$

$x_{jS}=\left[\begin{array}{c}D\cos E_{jS}\cos A_{jS}\\ D\cos E_{jS}\sin A_{jS}\\ -D\sin E_{jS}\end{array}\right]---\left(8\right)$

其中，X_gS为卫星在飞机所在点大地坐标系下的坐标，X_jS为卫星在机体坐标系下的坐标，将(7)(8)式代入(1)式，得：

通过上式可以计算出卫星在机体坐标系下方位角A_jS和俯仰角E_jS。

$A_{jS}={\mathrm{tg}}^{-1}\frac{y_{jS}}{x_{jS}},E_{jS}={\tan }^{-1}\frac{z_{jS}}{\sqrt{{x^2}_{jS}+{y^2}_{jS}}}---\left(10\right)$

步骤四：根据惯性导航测量的机体的方位角、俯仰角、横滚角，计算卫星在飞机所在点大地坐标系统下的大地方位角和俯仰角；

已知卫星在机体坐标系中方位、俯仰角为A_jS和E_jS及飞机的偏航角俯仰角θ、横滚角φ，将(7)(8)式代入(3)式，得：

通过上式可以计算出卫星在飞机所在点大地坐标系统下的大地方位角A_gS和俯仰角E_gS。

$A_{gS}={\mathrm{tg}}^{-1}\frac{y_{gS}}{x_{gS}},E_{gS}=-{\tan }^{-1}\frac{z_{gS}}{\sqrt{{x^2}_{gS}+{y^2}_{gS}}}---\left(12\right)$

步骤五：根据卫星在飞机所在点大地坐标系统下的大地方位角和俯仰角，反推出飞机所在的经度和纬度θ_Lat；

卫星在飞机所在点大地坐标系统下的大地方位角A_gS和俯仰角E_gS与飞机的坐标位置之间的关系：

假设卫星的经度为卫星到地面的高度为H，则卫星在地心坐标系下的坐标X_eS为：

卫星在飞机所在点的大地坐标系下的坐标X_gS为：

$x_{gS}=\left[\begin{array}{c}{x}_{gS}\\ y_{gS}\\ z_{gS}\end{array}\right]---\left(14\right)$

飞机位置点在地心坐标系下的坐标A为：

由飞机所在点大地坐标与地心坐标之间的转换关系得：

$x_{gS}=T_e^g(x_e-A)=T_e^g{x}_e-T_e^{g}A---\left(16\right)$

由上式得：因为： $A_{gS}={\mathrm{tg}}^{-1}\frac{y_{gS}}{x_{gS}},E_{gS}=-{\tan }^{-1}\frac{z_{gS}}{\sqrt{{x^2}_{gS}+{y^2}_{gS}}}$

得到：

由已得知的A_gS、E_gS、反推出飞机经度、纬度和θ_Lat：

设由式(20)化简得：

$A_{gS}={\mathrm{tg}}^{-1}\frac{tan\mathrm{\λ}}{{\mathrm{sin\θ}}_{Lat}}={\sin }^{-1}\frac{sin\mathrm{\λ}}{\sqrt{1-{{\mathrm{cos\θ}}^2}_{Lat}{\cos }^2\mathrm{\λ}}}---\left(21\right)$

$E_{gS}={\tan }^{-1}\frac{{\mathrm{cos\θ}}_{Lat}cos\mathrm{\λ}-(R+h)/(R+H)}{\sqrt{1-{{\mathrm{cos\θ}}^2}_{Lat}{\cos }^2\mathrm{\λ}}}---\left(22\right)$

由公式(21)简化得到： $\sqrt{1-{{\mathrm{cos\θ}}^2}_{Lat}{\cos }^2\mathrm{\λ}}=\frac{sin\mathrm{\λ}}{{\mathrm{sinA}}_{gS}}---\left(23\right)$

设 $\frac{sin\mathrm{\λ}}{{\mathrm{sinA}}_{gS}}=T,(R+h)/(R+H)=H^{\′}$

则公式(23)简化为：cosθ² _Latcos²λ＝1-T²

同理由公式(22)简化得到：(tan²E_gS+1)T²+2H'tanE_gST+(H'²-1)＝0(24)

则公式(24)简化为： $T=\frac{-H^{\′}.{\mathrm{tgE}}_{gs}\±\sqrt{{\mathrm{tg}}^2{E}_{gs}-H^{\′2}+1}}{({\mathrm{tg}}^2{E}_{gs}+1)}$

由公式(23)(24)结合计算得到飞机经度和纬度θ_Lat：

图6为本发明系统在某无人机飞行试验中实际经度和推算经度示意图；

其中，横轴代表时间，单位为秒；纵轴代表经度，单位为度；平滑轨迹为飞机实际经度轨迹；点画线曲线为本方法推算出经度轨迹。

图7为本发明系统在某无人机飞行试验中实际纬度和推算纬度示意图；

其中，横轴代表时间，单位为秒；纵轴代表纬度，单位为度；平滑轨迹为飞机实际纬度轨迹；点画线曲线为本方法推算出纬度轨迹。

步骤六：通过实际飞行数据，根据计算出的飞机的经纬度和导航下传的经纬度进行比较，验证本方法导航数据欺骗识别技术的有效性。

由于机载天线跟踪误差以及天线波束有一定宽度，即使天线主波束没有正对准卫星时，也可以进行信息传输，这种情况将会影响A_jS和E_jS的精度；另外一方面由于无人机平台的惯性测量部件测量的飞机机体方位角俯仰角θ、横滚角φ存在误差，也可能导致推算出的无人机位置存在误差。以某无人机实际飞行数据为例，将无人机导航系统给出的位置数据和根据卫星的大地角数据推算出的飞机的位置数据进行了比对实验。推算出的位置误差平均值27.76km，标准差:10.74，最大值：80km。实验表明：当导航未受到干扰时，推算出的经纬度在实际经纬度上下变化，如果受到欺骗后，得到的误差将会是单向变化。

图8为本发明系统在某无人机飞行试验中总误差比较示意图，其中，横轴代表时间，单位为秒；纵轴代表总误差，单位为千米；如图8所示曲线为本方法计算出的飞机的经纬度和导航下传的经纬度进行比较得到的总误差曲线。

Claims (1)

1.一种基于卫星通信的无人机导航数据被欺骗识别方法，包括以下几个步骤：

步骤一：建立与无人机相关的地心坐标系、机体坐标系统、大地坐标系；

设飞机经度为飞机纬度为θ_LatJ，卫星经度为地球半径为R_e，飞机到地面的高度为h，卫星到地面的高度为H，卫星到飞机的距离为D；

设地心坐标系为O(x_e,y_e,z_e)，原点O_e为地球的中心，x_e轴过地球中心指向真北方向，y_e轴过地球中心指向正东方向，在赤道平面内，z_e轴与x_ey_e轴成右手关系；

设大地坐标系为O(x_g,y_g,z_g)，原点O为地球质心，x_g轴过地球质心水平指向真北方向，y_g轴过地球质心，指向正东方向，z_g轴过地球质心，指向地心方向，与x_gy_g轴成右手关系，其中x_goy_g为水平面，x_goz_g为垂直面；

设机体坐标系为O(x_j,y_j,z_j)，原点O为机体重心，x_j轴过质心的机体滚动轴线，飞行方向为正；y_j轴在机体机翼水平面内，指向右边，z_j轴在飞机水平对称平面内，与x_jy_j轴成右手关系；

步骤二：获得坐标系之间的变换关系；

机体相对大地坐标系下的角度定义为：

偏航角为机体ox_j轴在大地坐标系x_goy_g水平面的投影与ox_g的夹角，自真北起顺时针方向为正，与ox_j构成右手旋，范围[0,360°]；

俯仰角θ为机体ox_j轴与大地坐标系x_goy_g水平面内的投影之间的夹角，以机头抬高为正，与oy_j构成右手旋，范围[-90°,90°]；

横滚角φ为机体oz_j轴与大地坐标系x_goz_g垂直面之间的夹角，右侧下降为正，与ox_j构成右手旋，范围[-90°,90°]；

大地坐标与机体坐标之间的转换关系为：

$x_j=T_g^j{x}_g---\left(1\right)$

其中， $x_j=\left[\begin{array}{c}{x}_j\\ y_j\\ z_j\end{array}\right]$ 表示机体坐标下的一个向量， $x_g=\left[\begin{array}{c}{x}_g\\ y_g\\ z_g\end{array}\right]$ 表示大地坐标坐标下的一个向量，为从大地坐标系到机体坐标系的旋转矩阵，具体为：

机体坐标与大地坐标之间的转换关系为：

$x_g=T_j^g{x}_j---\left(3\right)$

其中：为从机体坐标系到大地坐标系的旋转矩阵，具体为：

飞机所在点大地坐标与地心坐标之间的转换关系为：

$x_g=T_e^g{x}_e---\left(5\right)$

其中： $x_e=\left[\begin{array}{c}{x}_e\\ y_e\\ z_e\end{array}\right]$ 表示地心坐标下的一个向量，飞机所在的点在大地坐标系位置的经度和纬度分别为θ_Lat，地球半径为R，飞机的高度为h，为从地心坐标系到飞机所在点的大地坐标系的旋转矩阵；

步骤三：根据无人机天线跟踪误差电压跟踪卫星得到卫星在机体坐标下的方位角和俯仰角；

设卫星在大地坐标系中的方位、俯仰角为A_gS和E_gS，卫星在机体坐标系中方位、俯仰角为A_jS和E_jS，飞机姿态角为偏航角俯仰角θ、横滚角φ，D为卫星到飞机之间的距离，卫星在机体坐标下的方位角、俯仰角为A_jS和E_jS：

$x_{gS}=\left[\begin{array}{c}D\cos E_{gS}\cos A_{gS}\\ {\mathrm{DcosE}}_{gS}{\mathrm{sinA}}_{gS}\\ {\mathrm{DsinE}}_{gS}\end{array}\right]---\left(7\right)$

$x_{jS}=\left[\begin{array}{c}D\cos E_{jS}\cos A_{jS}\\ {\mathrm{DcosE}}_{jS}{\mathrm{sinA}}_{jS}\\ -{\mathrm{DsinE}}_{jS}\end{array}\right]---\left(8\right)$

其中，X_gS为卫星在飞机所在点大地坐标系下的坐标，X_jS为卫星在机体坐标系下的坐标，将(7)(8)式代入(1)式，得：

通过上式计算出卫星在机体坐标系下方位角A_jS和俯仰角E_jS；

$A_{jS}={\mathrm{tg}}^{-1}\frac{y_{jS}}{x_{jS}},E_{jS}={\tan }^{-1}\frac{z_{jS}}{\sqrt{{x^2}_{jS}+{y^2}_{jS}}}---\left(10\right)$

步骤四：根据惯性导航测量的机体的方位角、俯仰角、横滚角，计算卫星在飞机所在点大地坐标系统下的大地方位角和俯仰角；

已知卫星在机体坐标系中方位、俯仰角为A_jS和E_jS及飞机的偏航角俯仰角θ、横滚角φ，将(7)(8)式代入(3)式，得：

通过上式计算出卫星在飞机所在点大地坐标系统下的大地方位角A_gS和俯仰角E_gS；

$A_{gS}={\mathrm{tg}}^{-1}\frac{y_{gS}}{x_{gS}},E_{gS}=-{\tan }^{-1}\frac{z_{gS}}{\sqrt{{x^2}_{gS}+{y^2}_{gS}}}---\left(12\right)$

步骤五：根据卫星在飞机所在点大地坐标系统下的大地方位角和俯仰角，反推出飞机所在的经度和纬度θ_Lat；

卫星在飞机所在点大地坐标系统下的大地方位角A_gS和俯仰角E_gS与飞机的坐标位置之间的关系：

假设卫星的经度为卫星到地面的高度为H，则卫星在地心坐标系下的坐标X_eS为：

卫星在飞机所在点的大地坐标系下的坐标X_gS为：

$x_{gS}=\left[\begin{array}{c}{x}_{gS}\\ y_{gS}\\ z_{gS}\end{array}\right]---\left(14\right)$

飞机位置点在地心坐标系下的坐标A为：

由飞机所在点大地坐标与地心坐标之间的转换关系得：

$x_{gS}=T_e^g(x_e-A)=T_e^g{x}_e-T_e^{g}A---\left(16\right)$

由上式得：

因为： $A_{gS}={\mathrm{tg}}^{-1}\frac{y_{gS}}{x_{gS}},E_{gS}=-{\tan }^{-1}\frac{z_{gS}}{\sqrt{{x^2}_{gS}+{y^2}_{gS}}}$

得到：

由已得知的A_gS、E_gS、反推出飞机经度、纬度和θ_Lat：

设由式(20)化简得：

$A_{gS}={\mathrm{tg}}^{-1}\frac{tan\mathrm{\λ}}{{\mathrm{sin\θ}}_{Lat}}={\sin }^{-1}\frac{sin\mathrm{\λ}}{\sqrt{1-{{\mathrm{cos\θ}}^2}_{Lat}{\cos }^2\mathrm{\λ}}}---\left(21\right)$

$E_{gS}={\tan }^{-1}\frac{{\mathrm{cos\θ}}_{Lat}cos\mathrm{\λ}-(R+h)/(R+H)}{\sqrt{1-{{\mathrm{cos\θ}}^2}_{Lat}{\cos }^2\mathrm{\λ}}}---\left(22\right)$

由公式(21)简化得到： $\sqrt{1-{{\mathrm{cos\θ}}^2}_{Lat}{\cos }^2\mathrm{\λ}}=\frac{sin\mathrm{\λ}}{{\mathrm{sinA}}_{gS}}---\left(23\right)$

设 $\frac{sin\mathrm{\λ}}{{\mathrm{sinA}}_{gS}}=T,(R+h)/(R+H)=H^{\′}$

则公式(23)简化为：cosθ² _Latcos²λ＝1-T²

同理由公式(22)简化得到：(tan²E_gS+1)T²+2H'tanE_gST+(H'²-1)＝0(24)

则公式(24)简化为： $T=\frac{-H^{\′}.{\mathrm{tgE}}_{gs}\±\sqrt{{\mathrm{tg}}^2{E}_{gs}-H^{\′2}+1}}{({\mathrm{tg}}^2{E}_{gs}+1)}$

由公式(23)(24)结合计算得到飞机经度和纬度θ_Lat：

。