0 引言

近年来，导航定位、资源遥感、海洋气象等卫星系统以得天独厚的空间优势获取了丰富的地球信息资源，产生了巨大的经济和社会效益。然而随着太空探索的持续升温，人类愈发频繁的航天活动让整个空间日益拥挤，数以千万计的空间碎片威胁着每一个航天器的运行安全。作为空间态势信息的重要组成，海量、高速运动的空间目标的监管问题成为各国所面临的共同挑战^[1]。这其中，空间目标索引的优劣将直接影响数据查询和检索的速度，一定程度上决定了数据处理的效率。

目前国内外提及空间目标时空索引的文献相对较少，考虑到空间目标同样属于运动目标，可以为其提供参考。从实现原理上，运动目标时空索引主要包含3类：1) R树及其变形树^[2-4]，基于最小外包矩形建立，可直接对占据一定范围的空间对象进行索引，然而随着空间维数的增加，死空间和重叠问题较为严重；2) 四叉树及其变形树^[5]，依据分布密度无重叠地将空间递归划分为大小相等的4个象限，可有效索引多维数据，但树深的差异可能导致空间利用率的急剧下降；3) 网格结构及其变形算法^[6-7]，将研究区域划分为固定的网格，其中记录所包含的空间实体，查询操作时首先定位网格位置，进而对其中包含的对象进行检索，但其对线或面对象进行索引时同样存在冗余问题。

上述索引方法的提出均具有一定的针对性，就前两类常用方法来说，空间目标数量大、运行速度高，索引的快速更新将使其结构变得高维和复杂，进而大大降低查询和检索的效率。考虑到在惯性空间中目标受轨道动力学约束，运行相对稳定并非毫无规律，因此本文建立了基于轨道约束的空间目标球面网格索引，并据此提出区域查询的具体应用方案。

1 空间目标球面网格索引构建原理

空间目标受轨道动力学约束，以地球为焦点作椭圆运动，其绝对位置虽然随时间高速变化，但其运行轨道一定时间内在惯性空间保持相对稳定。本文提出的网格索引方法是在惯性空间，将地球的外层空间进行球面网格剖分，计算各空间目标的轨道所经过的网格，并将该目标编号记录到对应的网格中，实现对高速运动目标的稳定索引，为空间目标的各类查询打下基础。需要说明的是，目前球面网格剖分方法很多，包括经纬网格剖分^[8]、正多面体网格剖分^[9-10]、Voronoi网格剖分^[11-12]和球面等分剖分^[13]等几大类，考虑到实用、高效的原则，本文采用球面等间隔经纬网格剖分方法进行索引构建。具体构建步骤如下。

1) 在地心惯性系——J2000球坐标系下对地球外层空间进行球面等经纬网格剖分，如图 1所示，球面网格依照赤经、赤纬方向划分，赤经方向从春分点开始，逆时针为正，顺时针为负，实际应用中为了处理的方便，由赤经-180°依次对网格进行编号；赤纬方向从南极开始向北极计算网格编号，每个网格可以用grid[col][row]来表示，对应的赤经赤纬范围为[α_col, α_col+1, β_row, β_row+1]，其中：

$ \left\{ \begin{array}{l} {\alpha _{col}} = \Delta \alpha \times col\\ {\alpha _{col + 1}} = \Delta \alpha \times {\rm{(}}col + 1{\rm{)}}\\ {\beta _{row}} = \Delta \beta \times row\\ {\beta _{row + 1}} = \Delta \beta \times {\rm{(}}row + 1{\rm{)}} \end{array} \right. $

(1)

其中：Δα为网格赤经方向的间隔，Δβ为网格赤纬方向的间隔。

2) 依次计算目标所经过的所有网格编号，具体的计算方法是：依照一定的步长预报一个周期内空间目标的位置，计算这些位置所落的空间网格。

依据轨道预报模型计算空间目标在给定时刻J2000坐标系下的空间直角坐标，依照式(2) 将坐标换算到对应的天球坐标系下，换算过程中需要注意象限的判断及分母为0时的处理。

$ \left\{ \begin{array}{l} R{\rm{ = }}\sqrt {{X^2} + {Y^2} + {Z^2}} \\ \alpha = \arctan \left( {Y/X} \right)\\ \beta = \arctan \left( {Z/\sqrt {{X^2} + {Y^2}} } \right) \end{array} \right. $

(2)

其中，α的取值范围是[－π, π]，β的取值范围是[－π/2, π/2]，空间目标所在网格编号的计算公式为：

$ \left\{ \begin{array}{l} col = \left\lfloor {\left( {\alpha + {\rm{\pi }}} \right)/\Delta \alpha } \right\rfloor \\ row = \left\lfloor {\left( {\beta + {\rm{\pi }}/2} \right)/\Delta \beta } \right\rfloor \end{array} \right. $

(3)

其中，$\left\lfloor \cdot \right\rfloor $为向下取整符号。

3) 每个空间网格中维护一个目标列表，如果某一个空间目标在其预报周期内有点落在该网格中，则将该目标编号记录在网格的目标列表中。目标列表的形式如下：

$ grid{\rm{[}}col{\rm{][}}row{\rm{]}}: < num, I{D_1}, I{D_2}, \ldots, I{D_n} > $

(4)

其中grid[col][row]表示编号为(col, row)的网格，〈num, ID₁, ID₂, …, ID_n〉为经过对应网格的目标列表，num为目标数，ID₁, ID₂, …, ID_n为详细的目标编号。

4) 索引结构更新。由于受到地球扁率、大气阻尼等摄动因素的影响，空间目标的轨道会发生变化，为此需要对空间目标索引结构进行更新，方法是建立一个空间目标编号Hash表，Hash表的结构如式(5) 所示，Hash表记录了对应的空间目标所经过的网格，索引结构更新时，首先获取需要更新目标所经过的所有网格编号，然后在所有经过网格中删除对应的目标编号，最后按照步骤2) 和步骤3) 重新将该目标插入索引结构，完成更新。

$ ID:gri{d_1}, gri{d_2}, \ldots, gri{d_n} $

(5)

其中：ID表示目标编号；grid₁, grid₂, …, grid_n为该目标所穿过的网格。

2 基于索引的区域查询应用方案

查询是目标索引研究的主要目的，查询的性能是评价索引的重要指标。针对空间目标的查询应用主要有区域查询、K最近邻查询、聚集查询、连续查询和密度查询等类型^[14]。其中，区域查询应用最为广泛，可解算空间目标经过某区域的时间或探寻指定的时间段内经过某区域的空间目标，完成空间目标的过境分析预报。

空间目标区域查询的一般过程，即传统方法是：在指定的时间范围内，按照一定步长逐个预报目标的空间位置，通过判断星下点与区域的关系确定其进入和离开指定区域的时间。SKYMAP Pro^[15]是由英国SKYMAP公司开发的一款业界公认的天文仿真分析软件，采用如上方法根据用户输入的空间目标轨道数据完成过境分析，其更新效率将在后文中加以分析。

本文提出了基于球面网格索引的区域查询应用方案，同时用于检验该索引方法的效率。下面给出指定时刻，该区域查询方案的具体步骤，如图 2所示。

1) 查询时段离散化，要查询一段时间内指定区域的空间目标过境情况，需要将查询时间区间按照一定的步长进行离散化，将连续的时间转变为一系列查询时间片。

2) 区域惯性空间边界计算，对于给定的区域，计算其在查询时刻所对应的惯性系下的区域边界。假设区域范围的边界点在大地坐标系下的坐标为(L₁, B₁), (L₂, B₂), …, (L_n, B_n)，则首先将其换算到J2000坐标系下的空间直角坐标，然后依照式(2) 计算出其对应的J2000地心球坐标系下的坐标(α₁, β₁), (α₂, β₂), …, (α_n, β_n)。

3) 交叉网格确定，计算查询区域与空间球面索引网格的相交情况，可以用如下的方法进行计算：

① 首先将查询区域球面网格映射到(α, β)参数坐标系中，变为二维平面网格，如图 3所示。

② 求取查询区域的外包围矩形，如图 3所示，外包围矩形的4个角点坐标分别是：(α_min, β_min), (α_max, β_min), (α_max, β_max), (α_min, β_max)，其中：

$ \left\{ \begin{array}{l} {\alpha _{{\rm{min}}}} = {\rm{min\{ }}{\alpha _1}, {\alpha _2}, \ldots, \alpha {}_n{\rm{\} }}\\ {\alpha _{{\rm{max}}}} = {\rm{max\{ }}{\alpha _1}, {\alpha _2}, \ldots, \alpha {}_n{\rm{\} }}\\ {\beta _{{\rm{min}}}} = {\rm{min\{ }}{\beta _1}, {\beta _2}, \ldots, {\beta _n}{\rm{\} }}\\ {\beta _{{\rm{max}}}} = {\rm{max\{ }}{\beta _1}, {\beta _2}, \ldots, {\beta _n}{\rm{\} }} \end{array} \right. $

(6)

③ 依次判断查询区域外包围矩形内的网格是否与查询区域相交，记录下相交的网格编号grid′₁, grid′₂, …, grid′_k。图 3中浅灰色填充部分表示查询区域外包围矩阵内的网格，深灰色填充部分为与查询区域相交的网格。

④ 查询区域跨边界问题的处理。如果查询区域跨越边界，则需要按照边界将查询区域切分成多个子查询区域，分别计算各个子查询区域与网格的相交情况，图 4是几种切分示意图。

4) 网格内目标位置计算与判断，取出相交的网格中的目标编号，在给定的时刻，根据轨道预报模型计算这些目标的星下点坐标，并判断这些目标的星下点是否在查询区域内。

5) 递归查询，根据时间步长确定下一查询时刻，按步骤1) 到步骤3) 循环执行，直至时间超界结束。

3 实验验证与分析

为了验证所提方法的优越性，本文开展了相关实验，电脑硬件环境为Intel i5 CPU处理器，4GB内存，NVIDIA NVS 4200M显卡(1GB内存)，操作系统为Windows 7平台，编程环境为Qt5.4.1。实验数据主体为STK网站公布的13864条空间目标TLE双行根数，利用SGP4/SDP4预报模型进行轨道预报。

如第1章所述，目前运动目标的时空索引主要包含R树、四叉树和网格索引3类方法，在惯性系下若采用前两种方法进行区域查询，则在查询时间段内的每一个离散时刻均需要重新构建索引，众所周知索引的意义在于通过对象的提前编排来提高检索的效率，其构建的时长远大于利用索引进行查询的时长，如此而言，针对空间目标的区域查询时长将随查询时长线性增长且十分耗时，因而理论分析即可排除前两种方法应用的可能性，仅仅需要比较传统方法(见第2章)，即未构建索引直接判断目标星下点和区域位置关系，与本文所提方法之间的优劣。

从原理上讲，影响传统方法和本文所提方法进行空间目标区域查询效率的因素包括查询时长、目标数和查询区域大小。其中，查询时长对于两种方法的影响是一致的，不再讨论；由于要逐目标判断与查询区域的关系，因此传统方法的效率主要受限于目标数，而与区域大小几无关系；本文方法由于提前将所有空间目标的轨道计算工作置于索引构建中，所以在实际的区域查询时，不需要对所有目标进行遍历，只需判断与查询区域相交的网格中的目标是否满足条件即可，查询效率得以提升，因此，本文首先对传统方法用时与空间目标数之间的关系进行分析实验，进而针对不同区域大小对所提方法的时间效率开展验证，最后在固定区域下测试了本文查询方法的准确性。

3.1 传统方法与目标数关系分析

传统方法采用逐目标验证其星下点位置与区域关系的思路进行区域查询分析，实验起止时间为2013-06-03T12:00:00至2013-06-04T12:00:00，步长为60s，范围为经度[50, 115]，纬度[25, 30]的矩形区域进行查询。表 1列出了不同空间目标数下的分析用时，其中第1列给出了实验分析中不同的空间目标数；第2列给出了对应数量下采用传统方法进行查询所消耗的时间；第3列给出了单个目标平均查询用时。图 5展示了传统方法分析用时和空间目标数之间的关系，并绘制了拟合曲线。

从实验结果可以看出，随着空间目标数的增加，利用传统方法进行区域查询用时所消耗的时间呈线性增长，平均每个目标耗时均稳定在0.09ms左右，从而印证了传统方法耗时与空间目标数之间的线性正相关性。

利用传统方法对目前在轨工作的1300多个空间目标进行24h内的过境预报时，需要消耗超过85ms的时间；如果算上目前能观测的空间碎片，空间目标超过17000，这时消耗的时间将达到1530ms。可以预见，随着人类太空活动的加剧和观测能力的提升，供分析空间目标的数必将不断增加，采用传统方法进行查询的速度显然无法满足需求。

3.2 时间效率验证与分析

由前面分析可知，查询区域的面积不同时，基于空间目标时空索引的区域查询所消耗的时间也不同，为此本文固定查询时间，变化查询区域来测试基于空间目标时空索引的区域查询所用时间，实验中起止查询时间段为2013-06-08T15:32:23至2013-06-08T15:43:23，步长为60s。表 2列出了查询区域范围及用时，图 6展示了两种方法区域查询用时和查询区域所占网格数目的关系。

从实验结果可知，当查询区域所占网格数增加时，基于空间目标时空索引的区域查询所用的时间呈线性增长，而传统方法则恒定在1250ms上下，实验结果与前文理论分析相符，本文所提方法受查询面积影响较大，但是当所占网格数小于2750时，所用的查询时间仍然小于传统方法。实际应用当中，查询区域面积一般远小于2750个网格数，此时基于球面网格索引的区域查询具有很大的优势。

3.3 准确性分析

为了测试所提方法的准确性，本文进一步选取时间段为2013-06-08T15:32:23至2013-06-08T15:43:23，步长为60s，范围为经度[50, 115]，纬度[25, 30]的矩形区域进行实验。表 3列出了传统方法和本文方法所查询的过境空间目标具体数，其中第1列时刻指查询时离散化的时间片，这里以儒略日的形式给出，第2列给出了方法1在对应时刻所查询出的过境空间目标数，第3列给出了方法2在对应时刻所查询出的过境空间目标数，第4列给出了两种方法所查询到的过境目标数差值。

从实验结果可以看出，本文方法查询得到的过境目标数和传统方法的相同，此外通过逐目标人工比对，两种方法的查询目标结果同样一致，一般可以将传统方法查询的结果作为基准值，因此该实验结果验证了基于空间目标时空索引的区域查询方案准确性。

4 结语

本文针对海量、高速空间目标的快速索引问题，通过惯性空间内目标位置的预先计算，构建基于轨道约束的空间目标球面网格索引。在此基础上，提出相应的空间目标区域查询应用方案，支持一定时长范围内的目标的区域过境分析。实验结果表明传统的查询方法耗时与空间目标数线性正相关，因此无法满足海量空间目标的查询需求；与其相比，对于一般的区域范围，本文方法在保证准确率的前提下可以大幅提高时间效率。后续将考虑开展基于空间目标球面网格索引的K近邻查询、空间碎片预警等方面的应用研究。