遥感卫星对区域目标可见性的快速计算方法
鄂智博, 李俊峰     
清华大学 航天航空学院, 北京 100084
摘要:遥感卫星对区域目标可见时间窗口的快速计算,对卫星观测规划任务具有重要意义。该文提出了一种快速算法,首先根据卫星视场类型以及与地面点的几何关系,给出判断卫星对地面点可见性的方法;然后通过离散时间步长结合二分搜索的方法,精确计算卫星对单个离散点的可见时间窗口;最后根据所有离散边界点的可见窗口最终得到对区域目标的可见时间窗口。仿真结果表明:该算法具有很高的效率和精度,与STK(satellite tool kit)计算得到的结果偏差在0.01 s量级。
关键词遥感卫星    区域目标    可见性    快速计算方法    
Fast simulation algorithm for area target visibility using remote sensing satellites
E Zhibo, LI Junfeng     
School of Aerospace Engineering, Tsinghua University, Beijing 100084
Abstract: Fast calculations of area target visibility are important for satellite mission scheduling. This paper presents a fast simulation algorithm to compute the visibility of area targets by remote sensing satellites. A method was developed to judge the visibility of a ground point from a satellite based on the satellite's view field and the geometric relation between the satellite and the ground point. Then, the in-view period of a single point was obtained using large time step discretization and a binary search. The in-view period of all the discrete boundary points was then used to determine the visible time window for an area target. Numerical results show that the algorithm is efficient and accurate with the error between the algorithm and satellite tool kit (STK) predictions of aboutt 0.01 s.
Key words: remote sensing satellites     area target     visibility     fast simulation algorithm    

遥感卫星是执行对地面观测任务的空间资源,观测范围大,可以获得高质量的地表图像。在气象预报、灾害监测及军事侦察中都有重要应用[1-3]。卫星经过地面目标时,在视场约束下可以观测目标的时间范围称为可见时间窗口。它是卫星对地面目标执行观测任务的基础,任何观测任务都必须在可见窗口内执行[4-6],时间窗口计算具有重要的意义。

目前,计算卫星观测区域目标的可见时间窗口还没有高效的通用算法。跟踪传播法是计算卫星可见性的传统方法[7],包括工程中普遍应用的卫星工具包STK(satellite tool kit)均采用该算法。跟踪传播法利用较小的离散时间步长,预报卫星轨道位置,在每一个离散时间步长判断卫星对目标点的可见性,进而得到卫星对目标的可见时间窗口。该方法精度高,可以计算多种要求下的可见性,但效率很低,计算较慢,很难适应工程实际中仿真效率要求。在保证计算精度下,研究快速计算方法,提高可见窗口计算效率具有很大价值。

关于卫星可见性快速计算,已经有一些初步的研究。文[8]提出了真近点角迭代补偿的方法来计算椭圆轨道卫星的可见性,从而避免求解Kepler方程,提高了计算效率。文[9-10]提出了一种考虑地球引力长期摄动的卫星可见性计算方法,先在二体情况下粗略搜索计算卫星的大致的过顶时间,再考虑复杂引力精确搜索提高精度。文[11-12]利用插值函数拟合原可见性函数,将卫星可见性求解转化为线性方程求根问题,提高了计算效率。文[13]分析了轨道摄动对可见时间窗口计算的影响。然而,这些卫星可见性计算方法只适用于点目标,无法计算区域目标的时间窗口。目前,针对区域目标可见性快速计算的研究较少。文[14]基于可视扩展区域与轨道过滤,提出一种针对区域目标的可见时间窗口快速计算方法,适用于二体轨道,也无法考虑卫星侧视情况。文[15]基于Poincaré映射法来求解可见时间窗口,可以得到解析表达式,但因需要数值积分导致效率较低。文[16]提出一种基于迭代修正的时间窗口计算方法,首先求解不考虑自转的可见时间窗口,然后根据自转特征对计算结果进行迭代修正。可见,目前区域目标快速计算方法采用的轨道模型简单,很难适用于工程实际的复杂情况,并且计算精度还需进一步提高。因此,研究高效、精确、适应性强的区域目标可见窗口的快速计算具有重要意义。

本文针对区域目标,提出了一种快速计算卫星可见性的算法。首先针对圆锥视场和矩形视场,分别给出了卫星对地面点可见性的判别条件;然后先利用大步长初步确定开始与结束时刻的范围,再利用二分搜索精确计算开始与结束时刻,进而得到区域目标的可见时间窗口。

1 卫星与地面点几何分析

卫星对地面目标的观测过程中,卫星与地面目标的相对位置随时间不断变化。只有当地面目标位于卫星视场的俯仰与滚动方向的可视范围内,才认为卫星对该地面点目标可见,如图 1所示。

图 1 卫星的可视范围

卫星的可视范围跟搭载的相机视场类型有关,本文就2种最常用的视场类型给出可见性判据。

1.1 卫星与地面点的位置关系

卫星在一个轨道周期内,大部分时间与地面目标不在地球的同侧,此时卫星对地面目标是不可见的。因此,可先判断地面目标是否与卫星在地球同侧,再判断其是否位于卫星可视范围内。

卫星与地球位置如图 2所示,若地面点T在地固系中的位置RRe为地球半径,H为卫星到地心的距离,卫星某时刻t相对地固系的位置r(t)和速度v(t),则r(t)与R的夹角大于θmax,则卫星对地面目标不可见。若r(t)与R的夹角小于θmax,则需进一步计算地面目标是否位于卫星可视范围内。

图 2 卫星与地球位置

计算卫星指向地面点的角度时,需要将所有位置量在本体系下表示出来。卫星与地面点的位置关系如图 3所示。

图 3 卫星与地面点相对位置关系

卫星和地面点在J2000惯性系的位置矢量分别为rSatrpoint,其相对位置

$ \Delta \mathit{\pmb{r}} = {\mathit{\pmb{r}}_{{\rm{point}}}} - {\mathit{\pmb{r}}_{{\rm{Sat}}}}. $ (1)

J2000惯性系到卫星轨道坐标系的转换矩阵可表示为

$ \begin{array}{l} \;\;\;\;\;\;\;\;\;{\mathit{\pmb{L}}_{oi}} = {\mathit{\pmb{L}}_y}\left( { - \left( {\frac{\pi }{2} + u} \right)} \right){\mathit{\pmb{L}}_x}\left( {\left( { - \left( {\frac{\pi }{2} - i} \right)} \right){\mathit{\pmb{L}}_z}\left( \Omega \right)} \right) = \\ \left[ \begin{array}{l} - {\rm{cos\Omega sin}}{u} - {\rm{sin\Omega cos}}{i}{\rm{cos}}{u}\;\;\;\;{\rm{cos\Omega cos}}{i}{\rm{cos}}{u} - {\rm{sin\Omega sin}}{u}\;\;\;\;\;\;\;\;{\rm{cos}}{u}{\rm{sin}}{i}\\ \;\;\;\;\;\;\;\;\;\;\; - {\rm{sin\Omega sin}}{i}\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;{\rm{cos\Omega sin}}{i}\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;{-\rm{cos}}{i}{\rm{}}\\ \;\;{\rm{sin\Omega cos}}{i}{\rm{sin}}{u}{-\rm{cos\Omega cos}}{u}\;\;\;\;\;{-\rm{sin\Omega cos}}{u}{-\rm{cos\Omega cos}}{i}{\rm{sin}}{u}\;\;\;\;\;\;\;{-\rm{sin}}{i}{\rm{sin}}{u} \end{array} \right]{\rm{.}} \end{array} $ (2)

其中LxLyLz分别为绕各自轴的坐标转换矩阵。则卫星与地面点的矢量差Δr在卫星轨道坐标系中可表示为

$ \Delta {\mathit{\pmb{r}}_{{\rm{orbit}}}} = {\mathit{\pmb{L}}_{oi}}\cdot\Delta {\mathit{\pmb{r}}}. $ (3)

便可根据Δrorbit与轨道系z轴的几何关系确定地面点是否位于卫星可视范围内。

1.2 圆锥视场

圆锥视场的可视范围为一个圆锥形,如图 4所示。

图 4 圆锥视场

其中z坐标轴为卫星光轴。由图可以看出,卫星与地面点的矢量差Δr与卫星光轴夹角θ关系为

$ \theta = {\rm{arccos}}\frac{{\Delta {\mathit{\pmb{r}}}\cdot{\mathit{\pmb{e}}_z}}}{{\left\| {\Delta {\mathit{\pmb{r}}}} \right\|}}. $ (4)

小于圆锥视场半角时可见,大于圆锥视场半角时不可见。

1.3 矩形视场

矩形视场的可视范围为一个矩形,如图 5所示。

图 5 矩形视场

由图可见,矩形视场在滚动和俯仰方向的可视范围不同。因此,只有当卫星与地面点的矢量差ΔrxOz平面内的投影Δrxz,与在yOz平面内的投影Δryz分别小于相应方向的矩形视场半角时才为可见,否则不可见,即

$ {\theta _{xz}} = {\rm{arccos}}\frac{{\Delta {\mathit{\pmb{r}}_{xz}}\cdot{\mathit{\pmb{e}}_z}}}{{\left\| {\Delta {\mathit{\pmb{r}}_{xz}}} \right\|}} < {\theta _v}, $ (5)
$ {\rm{ }}{\theta _{yz}} = {\rm{arccos}}\frac{{\Delta {\mathit{\pmb{r}}_{yz}}\cdot{\mathit{\pmb{e}}_z}}}{{\left\| {\Delta {\mathit{\pmb{r}}_{yz}}} \right\|}} < {\theta _h}. $ (6)
2 卫星对地面点的可见性计算

本文采用了粗略计算与二分搜索结合的方法计算卫星对地面点的可见性。

2.1 粗略计算

粗略计算过程中,离散时间步长较大。目的是快速确定可见开始时刻与可见结束时刻所在的大致时间区间,主要步骤如下。

步骤1  已知观测任务的起始与结束时刻,将仿真时段划分为若干子区间段。设M0为卫星每次经过轨道升交角距(argument of latitude)等于90°的点。从开始时刻到卫星第1次过M0为第1段,用(tS1, tE1)表示,从卫星第1次过M0到第2次过M0为第2段,用(tS2, tE2)表示,依次划分,最后一次过M0到仿真结束时刻为最后一段。

步骤2  将航天器轨道推进的整个过程离散为m等份,即m+1个离散点(取100 s为步长,记相应时间序列为ti, i=1, 2, …, m+1)。

步骤3  计算每个离散时间点ti,卫星对地面点的滚动角和俯仰角,若在可视范围内,则认为卫星在该离散时间点ti对点目标可见。

步骤4  ti可见,ti-1不可见时,则可见窗口开始时刻位于ti-1ti之间;当ti-1可见,ti不可见时时,则可见窗口结束时刻位于ti-1ti之间。

2.2 精确计算

通过粗略计算,可见时间窗口开始与结束时刻的大致区间已经计算出来。利用二分搜索,可以高效地计算出精确的可见时间窗口。具体步骤如下。

步骤1  若可见窗口开始时刻位于时间区间(t0, tf)内,计算中间时刻tm的卫星对地面点的可见性。若tm可见,则把其值赋给tf;若tm不可见,则把其值赋给t0

步骤2  重复步骤1,直至tft0 < εε为给定的误差限。

步骤3  可见开始时刻为$t = \frac{1}{2}\left( {{t_0} + {t_f}} \right)$

由此便求得可见开始时刻,可见结束时刻的计算思路与之相同。

3 卫星对区域目标的可见性计算

计算卫星对区域目标的可见性难度要远远大于对单个地面点的可见性,主要体现在只要区域目标任意部分在卫星的视场范围内,则记为卫星对该区域目标可见。因此,传统的跟踪传播法在计算卫星对区域目标的可见性时,大多会将区域目标离散成网格区域。计算卫星对网格节点可见性,从而得到对整个区域目标的可见性。但是,这种方法计算量大,计算效率低。本文对此进行了改进,在保证精度的前提下,提升计算效率。具体步骤如下。

步骤1  已知观测任务的起始与结束时刻,按卫星每次过轨道升交角距等于90°的点M0的时刻对仿真时段进行分割,将仿真时段划分为若干小段。开始时刻到卫星第1次过纬度最大点为第1小段,用(tS1, tE1)表示,从卫星第1次过M0到第2次过M0为第2段,用(tS2, tE2)表示,依次划分,最后一次过M0到仿真结束时刻为最后一段。

步骤2  A, B为区域目标的某一边上的端点,在A, B之间取若干离散点,记为N1, N2, …, NN,分别计算卫星对各个离散点的可见时间窗口。如果对每个离散点Ni卫星在该周期内都无时间窗口,则表明卫星在该周期内对该区域无时间窗口。

步骤3  在一个仿真时段内,对其中有时间窗口的点,记为[tSNi, tENi],则该边在该周期内的时间窗口为:

$ t{S_{{\rm{AB}}}} = {\rm{min}}\left\{ {t{S_{{N_i}}}} \right\},$ (7)
$ t{E_{{\rm{AB}}}} = {\rm{min}}\left\{ {t{E_{{N_i}}}} \right\}. $ (8)

步骤4  对区域目标,由N条边组成,分别记作E1, E2, …, EN。如果对每条边Ei,卫星在该周期内都无时间窗口,则表明卫星在该周期内对该区域无时间窗口。对其中有时间窗口的边,记为[tSEi, tEEi],则区域目标在该周期内的时间窗口为:

$ t{S_{{\rm{Area}}}} = {\rm{min}}\left\{ {t{S_{{E_i}}}} \right\},$ (9)
$ t{E_{{\rm{Area}}}} = {\rm{max}}\left\{ {t{E_{{E_i}}}} \right\}. $ (10)

步骤5  由于区域目标边界的可见时间窗口是由该边离散点的时间窗口统计得到的,因此计算结果的精度与离散点的数目有关。为了保证计算效率,同时提高可见时间窗口的计算精度。在得到区域目标的可见开始时刻与结束时刻后,需要对其进行修正。即在[tSArea-Δtε, tSAreatε]与[tSArea-Δtε, tSAreatε]范围内,提高边界离散点数目,重新计算可见开始时刻与结束时刻。

步骤6  最后对每个周期内的时间窗口进行汇总,则可以得到卫星对区域目标在整个仿真时段的时间窗口。

4 数值模拟

为说明本文所述方法的有效性和高效性,本文计算了携带圆锥视场相机的卫星与携带矩形视场相机的卫星对区域目标的可见时间窗口,并与STK的结果进行比较。

仿真起始时刻取为2020年12月18日00:00:00 (UTCG),仿真结束时刻为2020年12月19日00:00:00(UTCG),仿真时间为1 d。仿真环境为4 GHz CPU的计算机,编程语言为C语言。卫星的轨道参数如表 1所示,其中卫星搭载的2种载荷,矩形视场水平和垂直半角均为30°;圆锥视场,圆锥半角为30°。粗略搜索步长设为100 s,精确时间步长设为0.000 1 s,跟踪传播算法的时间步长设为0.1 s。区域目标为多边形,顶点经纬度如表 2所示。

表 1 卫星轨道参数
参数 参数值
半长轴/km 7 128.14
离心率 0.0
倾角/(°) 19.925
升交点赤经/(°) 0.0
近地点幅角/(°) 219.484
平近点角/(°) 326.698

表 2 区域目标顶点经纬度
顶点 经度/(°) 纬度/(°)
1 110.0 22.0
2 110.0 12.0
3 118.0 12.0
4 118.0 22.0

两颗卫星对区域目标的可见性计算所需CPU时间见表 3,对区域目标的可见时间窗口具体信息见表 4。本文算法可见性计算结果与STK计算结果偏差对比如图 6所示。

表 3 可见性计算所需CPU时间
视场类型 快速算法/s 传播算法/s
圆锥 1.093 214.983
矩形 1.406 289.893

表 4 可见时间窗口信息
视场类型 编号 可见窗口开始时刻 持续时间/s
圆锥 1 2020-12-18 03:16:24.410 260.006
2 2020-12-18 05:03:11.299 283.231
3 2020-12-18 06:50:09.743 275.930
4 2020-12-18 08:36:58.416 279.174
5 2020-12-18 10:23:58.834 290.094
矩形 1 2020-12-18 03:16:06.085 298.151
2 2020-12-18 05:02:57.545 310.715
3 2020-12-18 06:50:04.612 282.821
4 2020-12-18 08:36:53.304 295.551
5 2020-12-18 10:23:42.731 307.789

图 6 计算结果与STK结果偏差

表 3可见,本文提出的快速算法比跟踪传播算法有明显的优势。在保证精度的前提下,快速算法计算圆锥视场相机观测区域目标的时间窗口所需CPU时间为0.245 s,而精度为0.1 s的跟踪传播法所需时间为214.983 s,快速算法所需时间仅为跟踪传播算法的0.114%,计算矩形视场可见性所需时间也仅为跟踪传播算法的0.144 9%。因此,相比于跟踪传播算法,快速算法的精度更高,效率更高。

表 4可知,携带圆锥视场相机的卫星1 d内有5个轨道周期能观测到区域目标,携带矩形视场相机的卫星也有5个轨道周期能观测到区域目标。从可见窗口持续时间可知,携带圆锥视场相机的卫星对区域目标仅有1 500 s的可见时长,仅占整个规划时间段的1.73%。因此,绝大部分时间卫星对区域目标不可见。采用大时间步长粗略计算,可以有效减少计算量,提高计算效率。

图 6可以看出,圆锥视场相机观测区域目标时间窗口长度误差不超过0.248 s,平均绝对误差为0.104 5 s,相对误差小于0.089 8%。矩形视场相机观测区域目标窗口长度误差不大于0.606 s,平均绝对误差为0.372 s,相对误差小于1.973%。可以看出,本文算法精度很高,满足绝大多数工程需求。偏差主要来自卫星轨道计算误差,如果能提高轨道计算精度,则可见时间窗口计算可以进一步提高。

图 6中还可以看出,矩形视场相机观测区域目标时间窗口长度与STK计算结果偏差较大。主要原因在于区域目标边界离散点数目较少,当矩形视场的角点先进入或离开区域时,可见开始时刻与可见结束时刻计算不准确,如图 7所示。该项误差可以通过局部加密区域边界离散点数目来进行修正,修正后的计算结果与STK结果偏差如图 8所示。

图 7 矩形视场与区域目标位置关系

图 8 修正计算结果与STK结果偏差

图 8可以看出,修正后的圆锥视场相机观测区域目标时间窗口长度误差不超过0.032 s,平均绝对误差为0.012 2 s,相对误差小于0.011%。矩形视场相机观测区域目标窗口长度误差不大于0.069 s,平均绝对误差为0.032 2 s,相对误差小于0.024 4%。可以看出,修正后的结果精度相对于修正前提高了近10倍,与STK计算结果偏差在0.01 s量级,可以满足绝大多数工程需求。

为了证明结论的一般性,本文计算了10组不同形状区域目标的可见性,与STK计算结果进行对比,并将相对误差绘制于图 9。由图 9可以看出,圆锥视场相机观测区域目标时间窗口长度相对误差不超过0.216 8%,平均相对误差为0.067 9%。矩形视场相机观测区域目标时间窗口长度相对误差不超过0.317 8%,平均相对误差为0.015 9%。可见,本文算法具有较高精度。