张晴晴,王剑,刘瑞华

(中国民航大学 电子信息与自动化学院,天津300300)

0 引 言

大圆航线是指球面上两点之间的最短航线,飞机在远距离飞行过程中按照大圆航线航行[1-2].目前,大地主题的解算主要有高斯平均引数公式、史赖伯公式、巴乌曼投影法、贝塞尔公式以及文森特(Vincenty)嵌入系数法等[3-4]. 其中,高斯平均引数公式和史赖伯公式通常用于短距离大地主题解算、巴乌曼投影法广泛用于中距离大地主题解算、贝塞尔公式适用于中长距离解算,1975年T.Vincenty以贝塞尔公式为基础,推导出嵌套系数公式也可用于中长距离解算.

随着全球导航卫星系统的发展,飞机在飞行过程中可以脱离地面导航台的约束,实现点到点的飞行[5]. 但飞机在飞行过程中可能会出现可见卫星数少于5颗的情况,导致接收机自主完好性监测(RAIM)空洞,对飞行安全造成威胁. 因此,需要对航路上的各时空点进行RAIM预测,而精确推算某一时刻的飞机位置是进行RAIM预测的前提. 《使用全球定位系统(GPS)的机载辅助导航设备的最低运行性能标准》(RTCA DO-208)建议了机载卫星导航设备飞行性能的测试条件和测试方法. 在进行飞行性能测试时需要进行航道选择、路径计算等. 飞机飞行过程中大部分时间都处于巡航阶段[6],因此本文采用RTCA DO-208推荐的Vincenty嵌入系数法计算巡航阶段的航程和航向,然后依据航程和飞行速度推算过点时间,从而得到精确的四维航迹,为飞机进行精准的RAIM预测提供数据支持,进而更好地保障飞机飞行安全.

1 大地主题正反算算法原理

1.1 大地主题反算

已知WGS-84椭球上两个点的大地纬度和经度,找到两点之间最短路径的距离和方位称为大地主题反算.

WGS-84椭球航线的模型如图1所示,起点经纬度坐标为P1(L1,B1),终点坐标为P2(L2,B2).α1是起点经线按顺时针方向到大椭圆航线夹角,称为起点方位角,α2是终点经线按顺时针方向到大椭圆航线夹角,称为终点方位角.

由起点P1和终点P2计算起点方位角α1(单位:(°))、终点方位角α2(单位:(°))、起点地到终点的距离s(单位:m)即为大地主题反算的内容.

首先将起点和终点的纬度从度转换成弧度记为φ1、φ2,然后计算经度差:

ΔL=(π/180)(L2-L1).

(1)

接着对P1、P2点的纬度进行归化处理:

β1=arctan[(1-f)tanφ1].

(2)

β2=arctan[(1-f)tanφ2].

(3)

式中,f=1/298.257为WGS-84椭球的扁率.

设初始迭代条件λ0=ΔL,迭代公式为式(4),直到满足条件|λk+1-λk|<ε时终止迭代,其中ε=10-12度[7].

λk+1= ΔL+(1-C)fsinα{σ+Csinσ×

[cos2σm+Ccosσ(-1+2cos22σm)],

(4)

sinσ= (cosβ2sinλk)2+(cosβ1sinβ2-

cosσ= sinβ1sinβ2+cosβ1cosβ2cosλk,

式中:α为大椭圆航线在赤道的方位角;σ为起点与终点间的球面角距;σm为大椭圆航线与赤道的交点到大椭圆航线中点的球面角距.

最后分别计算距离s、起点和终点各自的方位角α1、α2[7]:

s=bA(σ-Δσ),

(5)

-768+u2(320-175u2)},

4sin2σ)(-3-4cos22σm)cos2σm]}，

sinβ1cosβ2cosλk+1).

(6)

cosβ1sinβ2cosλk+1).

(7)

式中:a=6378137 m、b=6356752.3142 m分别为WGS-84椭球的长半轴和短半轴.

1.2 大地主题正算

已知两点之间的距离和方位以及起点经纬度,计算终点经纬度称为大地主题正算.

设椭球上两点之间的距离为s,起点坐标为P1(L1,B1),起点到终点的方位为α1,计算终点经纬度(L2,B2)的公式为[8]:

B2=arctan

180/π.

(8)

180/π.

(9)

迭代公式为

(10)

4sin2σ)(-3-4cos22σm)cos 2σm]}.

2σm=2σ1+σ

=2arctan(tan(β1)/cos(α1))+σ.

式中,参数β1、b、A、σm与大地主题反算算法中含义相同. 对式(10)进行迭代计算,满足|Δσ|<ε时停止计算,ε=10-12,初始迭代时σ=s/bA.

2 大椭圆航线计算及验证

大椭圆航线的验证分为正、反验证. DO-208中给出了7组典型位置点,如表1所示,各位置点参考结果如表2所示[7].

表1 示例点经纬度

表2 示例点参考结果

根据表1中7组位置点经纬度,按照Vincenty大地主题反算,可计算大椭圆航线的距离和起点方位角以及终点方位角. 反算结果与表2中参考结果之差如图2所示,可以看出:起点方位角和终点方位角在第4个示例点的误差最大,分别为(-4.9784×10-5)°和(4.9744×10-5)°,其他示例点的误差接近于0;各个测试点的距离误差均在±5×10-3m的范围内.

根据7个示例点起始经纬度以及示例点间的距离和起点方位角,进行正算验证,得到终点经纬度如表3所示.

表3 大椭圆航线正算输出示例比较

从正算结果看出:第3个测试点经度差的绝对值较大,纬度差的绝对值只在第4个测试点偏大.总体而言,各计算结果误差的绝对值接近于零,各测试点计算结果的误差均在(±5.56×10-8)°的范围内.

经过上述验证,根据大椭圆航线反算得到的航线距离和航向以及根据大椭圆航线正算得到终点经纬度与RTCA DO-208中给出的参考结果误差极小,因此可以根据大椭圆航线正反算推算飞行过程中各航迹点的坐标.

3 航迹推算及验证

3.1 航迹推算

航迹推算主要包括航空器飞行轨迹推算和航空器过航路点时间推算[9-10]. 飞机从起点到终点的飞行路径并不是严格按大椭圆航线规划,因为飞机在一些地形复杂、气象条件恶劣、空管限制的区域必须绕行,如高峰、雷雨云、禁飞区等. 在整个飞行过程还会途经若干航路点,所以飞行过程中采用分段大圆航线飞行.

本文采用在大椭圆航线上进行内插,这里的内插指在两个航路点之间按照指定的距离间隔取点.长距离靠近大圆航线缩短航程从而提高航行效率,然后在各内插点之间按照恒向线航行,避免了大椭圆航线航向时刻变化的问题,方便导航.

根据起点、航路点以及终点的经纬度坐标,利用大地主题反算得到起始点到下一个航路点的距离和初始的起点方位. 然后通过设置飞行速度和时间间隔唯一确定飞行距离以及内插点的个数. 根据大地主题正算求得第一个内插点的经度和纬度,将该内插点作为起点,按照大地主题反算更新方位角,再通过大地主题正算得到下一内插点的经纬度,直至求出起点到下一个航路点间所有的内插点. 然后以此类推,得到起点到终点的所有内插点,各内插点的坐标即为所推算的航迹点的坐标. 另外,时间间隔和速度需根据不同的飞行阶段设定. 具体算法流程如图3所示.

3.2 仿真验证

以上海虹桥机场(ZSSS)到成都双流国际机场(ZUUU)实际飞行过程为例,已知巡航阶段途经的六个航路点以及过点时间. 飞机巡航高度为9174 m,巡航速度为700 km/h,在各途经航路点之间等时间间隔进行内插,推算航行过程中的航迹点坐标. 设时间间隔为5 min,推算出相邻两个途经航路点之间的航空器位置,推算结果如图4和图5所示.

仿真得到经过各途经航路点的时间如表4所示.在巡航阶段,飞机实际飞行速度并不是固定不变的,在仿真中使用的速度是整个巡航阶段的平均速度,所以会导致过点时间有所不同. 虽然仿真得到的过点时间与实际过点时间有差异,但两者时间相差最大为47 s,能够较为精确地对飞行航迹进行预测.

表4 过点时间

4 结束语

本文通过对基于WGS-84椭球的航迹点坐标推算算法的分析,以上海虹桥机场到成都双流国际机场巡航阶段实际飞行过程为例,仿真计算了飞机经过各航路点的时间. 结果显示:实际过点时间与不考虑实际飞行过程中风速、流量管制以及偏航等因素影响的情况下得到的仿真过点时间相差均小于1 min,最小仅为9 s,最大为47 s,从而验证了该航迹推算算法的可靠性和可用性.