CN103679780A - 一种空间目标实时模拟方法 - Google Patents

Abstract

本发明涉及一种空间目标实时模拟方法，属于航天仿真技术领域。本发明首先构建可视化的星空背景，然后结合空间目标的轨道根数，利用可编程处理器完成位置、大小和颜色的高速解算和绘制，最后结合视点位置和用户裁剪信息，完成目标方向空间态势的生成。本发明充分利用了GPU的并行处理和高精度浮点运算特性，加速效果明显，实时性好，为空间态势表达中目标的实时模拟提供可靠的手段。

Description

一种空间目标实时模拟方法

技术领域

本发明涉及一种空间目标实时模拟方法，属于航天仿真技术领域。

背景技术

随着人类探索太空的步伐越来越快，作为直观展示空间目标分布、空间环境状态的重要手段，基于航天仿真技术的空间态势表达受到了越来越多的关注。相较于地面态势的表达，空间态势表达最大的特点在于所要描述的时空范围非常之大，由此包含的空间目标的数量级别更是数以万计，这对于目标的实时模拟是一个巨大的挑战。

目前，国内外已有多家单位对空间态势表达进行了一定的探索，但更多的集中于星空背景、特定或少量空间目标的仿真模拟，对空间目标的实时模拟研究较少。OpenGL作为图形硬件的一种软件接口，是目前最流行的3D图形和模型库之一，可以辅助完成空间目标的一般模拟。然而，对于空间目标的实时模拟，采用OpenGL固定功能管线处理数据无法满足速度上的要求。图形处理器（GPU）能够从硬件上支持T&L（Transform and Lighting），具有并行性好和浮点运算能力强的特点，通过图形硬件编程（如图1所示）可以分担GPU大部分的T&L运算，从而大幅加快3D渲染的速度。随着GPU的应用范围进一步扩大，一些通用的计算任务通过着色语言转移至GPU运行，这为空间目标的实时位置解算和绘制渲染提供了新的思路，

发明内容

本发明的目的是提供一种空间目标实时模拟方法，以解决采用OpenGL固定功能管线处理方式来实时模拟空间目标导致时速度无法满足的问题。

本发明为解决上述技术问题而提供一种空间目标实时模拟方法，该实时模拟方法的步骤如下：

1）计算恒星位置、大小和亮度，并根据恒星像点的成像模型，绘制星空背景；

2）根据空间目标轨道根数利用GPU实时解算空间目标的位置，并分别利用GPU中的顶点着色器和片段着色器设置空间目标的大小和颜色；

3）结合视点位置和用户裁剪信息，将空间目标渲染成目标相应方向的空间态势。

所述步骤2）中对空间目标位置实时解算及属性设置的具体步骤如下；

A.建立轨道定向四元数与轨道角元素之间的关系，利用四元数对空间目标的轨道姿态进行描述；

B.将空间目标内置属性中的轨道根数与顶点位置、表面法线和纹理坐标进行关联；

C.利用OpenGL着色语言对顶点着色器进行编程，实现对空间目标位置的解算，并设置其大小；

D.对空间目标的图元进行装配和光栅化处理，将图元分解为更小的片元，并为片元分配不同的透明度，结合透明度变化对目标颜色进行设置。

所述步骤A中轨道定向四元数与轨道角元素之间的关系为：

$q_x=\cos \frac{i}{2}\cos \frac{\mathrm{\Ω}+f}{2},q_y=\sin \frac{i}{2}\cos \frac{\mathrm{\Ω}-f}{2}$

$q_z=\sin \frac{i}{2}\sin \frac{\mathrm{\Ω}-f}{2},q_w=\cos \frac{i}{2}\sin \frac{\mathrm{\Ω}+f}{2}$

其中q_x,q_y,q_z,q_w为轨道姿态的四元数，Ω为升交点赤经，i为轨道倾角，σ＝ω+f为升交点角距，ω为近地点角距，f为真近点角。

所述步骤C中利用轨道根数解算空间目标的天球坐标系坐标的过程如下：

a.由历元平近点角M₀、平均角速度n和当前历元时间time，计算平近点角，

M＝M₀+time*n；

b.由偏心率e，依据开普勒方程M＝E-e sin E，迭代求解偏近点角E；

c.由长半轴a，计算目标在轨道平面直角坐标系的坐标，

$\left(\begin{array}{c}x\\ y\end{array}\right)=a\left(\begin{array}{c}\cos E-e\\ \sqrt{(1-e^2)\sin E}\end{array}\right);$

d.利用轨道姿态四元数计算目标在天球坐标系中的坐标，

${\left(\begin{array}{c}X\\ Y\\ Z\end{array}\right)}_{\mathrm{CS}}=\left(\begin{array}{c}{\mathrm{xq}}_x^2+{\mathrm{yq}}_x{q}_y+z{q}_x{q}_z+a_0{q}_w-a_1{q}_z+a_2{q}_y\\ {\mathrm{xq}}_x{q}_y+{\mathrm{yq}}_y^2+{\mathrm{zq}}_y{q}_z+a_0{q}_z+a_1{q}_w-a_2{q}_x\\ {\mathrm{xq}}_x{q}_z+{\mathrm{yq}}_y{q}_z+{\mathrm{zq}}_z^2-a_0{q}_y+a_1{q}_x+a_2{q}_z\end{array}\right)$

其中 $\left(\begin{array}{c}{a}_0\\ a_1\\ a_2\end{array}\right)=\left(\begin{array}{c}{\mathrm{xq}}_w-{\mathrm{yq}}_z\\ {\mathrm{xq}}_z-{\mathrm{yq}}_w\\ {\mathrm{yq}}_x-{\mathrm{xq}}_y\end{array}\right)$ 为过渡矩阵，q_x,q_y,q_z,q_w为轨道姿态四元数。

所述步骤D中为片元分配的透明度模型如下：

$\mathrm{opacity}=1-\frac{{(x-x_c)}^2+{(y-y_c)}^2}{\max \{{(x-x_c)}^2+{(y-y_c)}^2\}}$

其中，(x,y)为片元坐标，(x_c,y_c)为图元中心坐标。

本发明的有益效果是:本发明首先构建可视化的星空背景，然后结合空间目标的轨道根数，利用GPU实时解算空间目标的位置，并利用GPU中的顶点着色器和片段着色器设置空间目标的大小和颜色，最后结合视点位置和用户裁剪信息，将空间目标渲染成目标相应方向的空间态势，从而完成目标空间态势的生成。本发明充分利用了GPU的并行处理和高精度浮点运算特性，加速效果明显，实时性好，为空间态势表达中目标的实时模拟提供了可靠的手段。

附图说明

图1是可编程图形硬件流水线；

图2是本发明空间目标实时模拟方法的流程图；

图3是空间目标位置解算及属性设定示意图；

图4是轨道六根数示意图；

图5是目标片元透明度修改前后对比图；

图6是着色器创建流程图；

图7是近地空间不同类型航天器的实时模拟结果图；

图8是太阳系小行星带的模拟结果图。

具体实施方式

下面结合附图对本发明的具体实施方式作进一步的说明。

OpenGL作为图形硬件的一种软件接口，是目前最流行的3D图形和模型库之一，可以辅助完成空间目标的一般模拟，但是，对于空间目标的实时模拟，采用OpenGL固定功能管线处理数据无法满足速度上的要求。图形处理器（GPU）能够从硬件上支持T&L（Transform and Lighting），具有并行性好和浮点运算能力强的特点，通过图形硬件编程（如图1所示）可以分担GPU大部分的T&L运算，从而大幅加快3D渲染的速度。随着GPU的应用范围进一步扩大，一些通用的计算任务通过着色语言转移至GPU运行，这为空间目标的实时位置解算和绘制渲染提供了新的思路，为此，本发明基于上述思路，提供了一种空间目标的实时模拟方法，该实时模拟方法的步骤如下：

1.计算恒星位置、大小和亮度，并根据恒星像点的成像模型，绘制星空背景；

2.根据空间目标轨道根数利用GPU实时解算空间目标的位置，并设置空间目标的大小和颜色；

3.结合视点位置和用户裁剪信息，将空间目标渲染成目标相应方向的空间态势。

下面针对上述三个步骤进行详细说明。

1.星空背景绘制

在空间态势表达中建立准确的三维星空背景，是建立逼真空间场景，进行目标实时模拟的基础。本发明通过相关的模型计算出恒星位置、大小、亮度，结合恒星像点的成像模型，利用OpenGL绘制出位置准确、亮度分明的星空，为空间目标的实时模拟提供真实、可靠的背景。该过程具体包括恒星位置模型的计算、恒星亮度模型的计算、恒星大小模型的计算和恒星像点成像模型的计算。

恒星位置模型的计算过程如下：

恒星相对太阳系的运动满足一定的规律，已知恒星的历元平位置和自行即可利用公式（1）推算任意时刻恒星的平位置，利用公式（2）进而得到恒星t时刻在地心赤道坐标系中的位置。

$\left\{\begin{array}{c}{\mathrm{\α}}_t={\mathrm{\α}}_0+{\mathrm{\μ}}_{\mathrm{\α}}(t-t_0)\\ {\mathrm{\δ}}_t={\mathrm{\δ}}_0+{\mathrm{\μ}}_{\mathrm{\δ}}(t-t_0)\end{array}\right.---\left(1\right)$

其中，μ_α和μ_δ分别是恒星在赤纬和赤经方向的自行，可以从依巴谷星表中读取得到。

由恒星的赤经、赤纬坐标(α_t,δ_t)转换到地心直角坐标(x,y,z)的计算公式为：

$\left\{\begin{array}{c}x=r\cos {\mathrm{\α}}_t\cos {\mathrm{\δ}}_t\\ y=r\cos {\mathrm{\δ}}_t\sin {\mathrm{\α}}_t\\ z=r\sin {\mathrm{\δ}}_t\end{array}\right.---\left(2\right)$

其中，r为地球到恒星的距离，单位为光年，利用星表中的秒差距可以换算得到。

恒星亮度模型的计算过程如下：

恒星的亮度无法按照真实情况进行定量表现，考虑星等和视觉两个方面两个因素通过线性插值的方法来进行定性模拟。

B＝B₀-k·G   （3）

其中，G为恒星星等，B₀和k的设定首先满足计算后的亮度值位于OpenGL颜色信息限定范围0.0～1.0之间，其次根据对视觉亮度的需求进一步调整。例如，星表的最大星等为13，设想将恒星亮度控制在0.5～0.1之间，可以设定B₀＝1.0，k＝0.038。

恒星大小模型的计算过程如下：

对于恒星大小的模拟同样无法按照真实情况定量表现，采用的方法等同于亮度模型，即通过线型插值模拟不同星等的恒星的大小。

R＝R₀-l·G   （4）

其中，G为恒星星等，恒星大小的绘制在OpenGL中没有刚性的限制，可以根据视觉感受自行设置R₀和l。例如，星表的最大星等为13，设想将恒星大小控制在1～3个像素之间，可以设定R₀＝3，l＝0.1534。

恒星像点成像模型的计算过程为：

无穷远处的恒星可以看作具有一定光谱特征的点光源，因此模拟星图的像点为光学***的点扩散函数，理想情况下其能量分布满足二维高斯分布，

I(x,y)＝(φ/2πσ²)exp[-(x-x₀)²/2σ²-(y-y₀)²/2σ²]   （5）

其中φ为恒星的辐射功率；(x₀,y₀)为光斑中心位置坐标。

2.空间目标位置实时解算及属性设置

该过程主要包括利用四元数描述目标轨道、将轨道根数与内置属性进行变量关联、目标位置解算和目标颜色设置。

利用四元数描述目标轨道的过程如下：

对于轨道的描述通常采用轨道椭圆长半轴、偏心率、轨道倾角等六个参数（图4所示），其中长半轴和偏心率可以确定轨道的大小和形状，真近点角可以确定卫星在轨道上的瞬时位置，轨道倾角、升交点赤经和近地点角距则用来确定轨道面的姿态。传统的对于姿态方位的描述多采用欧拉角，然而会导致方程中出现复杂繁琐的三角函数和附加的奇异点，若使用四元数作为描述轨道姿态的参数，就可以避开这些缺点。

因此，数据准备阶段将建立轨道定向四元数与轨道角元素的关联，假设描述轨道姿态的四元数为q_x,q_y,q_z,q_w，其与轨道角元素的关系模型如下所示：

$q_x=\cos \frac{i}{2}\cos \frac{\mathrm{\Ω}+f}{2},q_y=\sin \frac{i}{2}\cos \frac{\mathrm{\Ω}-f}{2}$

                                       （6）

$q_z=\sin \frac{i}{2}\sin \frac{\mathrm{\Ω}-f}{2},q_w=\cos \frac{i}{2}\sin \frac{\mathrm{\Ω}+f}{2}$

其中，Ω为升交点赤经，i为轨道倾角，σ＝ω+f为升交点角距，ω为近地点角距，f为真近点角。

将轨道根数与内置属性变量关联的过程为：

轨道根数指的是描述天体（空间目标）在其轨道运行状态的一组参数，目前GPU仅仅可以处理顶点位置、法向量、纹理坐标等图元信息，因此将其应用于通用科学计算的前提是，将参与运算的相关参数存入顶点缓冲区等待调用，因此，本发明将轨道根数值存入顶点缓冲区，各参数与内置属性变量的对应关系如表1所示。

表1轨道根数与内置属性关联关系

通过关联，可以通过标准的OpenGL顶点数据接口将轨道根数数据从应用程序发送到顶点着色器程序，在OpenGL中，所有的几何物体最终都描述成一组有序的顶点。

利用轨道根数解算目标的天球坐标系坐标的步骤如下：

①由历元平近点角M₀、平均角速度n和当前历元时间time，计算平近点角

M＝M₀+time*n   （7）

②依据开普勒方程M＝E-e sin E，迭代求解偏近点角E；

③计算目标在轨道平面直角坐标系的坐标

$\left(\begin{array}{c}x\\ y\end{array}\right)=a\left(\begin{array}{c}\cos E-e\\ \sqrt{(1-e^2)\sin E}\end{array}\right)---\left(8\right)$

④利用轨道姿态四元数计算目标在天球坐标系中的坐标。

${\left(\begin{array}{c}X\\ Y\\ Z\end{array}\right)}_{\mathrm{CS}}=\left(\begin{array}{c}{\mathrm{xq}}_x^2+{\mathrm{yq}}_x{q}_y+z{q}_x{q}_z+a_0{q}_w-a_1{q}_z+a_2{q}_y\\ {\mathrm{xq}}_x{q}_y+{\mathrm{yq}}_y^2+{\mathrm{zq}}_y{q}_z+a_0{q}_z+a_1{q}_w-a_2{q}_x\\ {\mathrm{xq}}_x{q}_z+{\mathrm{yq}}_y{q}_z+{\mathrm{zq}}_z^2-a_0{q}_y+a_1{q}_x+a_2{q}_z\end{array}\right)---\left(9\right)$

其中 $\left(\begin{array}{c}{a}_0\\ a_1\\ a_2\end{array}\right)=\left(\begin{array}{c}{\mathrm{xq}}_w-{\mathrm{yq}}_z\\ {\mathrm{xq}}_z-{\mathrm{yq}}_w\\ {\mathrm{yq}}_x-{\mathrm{xq}}_y\end{array}\right)$ 为过渡矩阵，q_x,q_y,q_z,q_w为轨道姿态四元数。

上述过程需要通过GLSL将其翻译为着色器可以识别的“语言”，然后交由GPU进行高速并行运算，GLSL有许多特性与C++和Java相似，对应以上步骤的关键代码如下所示：

float sma=gl_Vertex.x;//长半轴

float ecc=gl_Vertex.y;//偏心率

float M0=gl_Vertex.z;//历元平近点角

float nu=gl_Normal.x;//平均角速度

float M=M0+time*nu;//计算平近点角

float E=M;

for(int i=0;i<4;i+=1)

E=M+ecc*sin(E);//迭代循环4次求偏近点角

vec3position=vec3(sma*(cos(E)-ecc),sma*(sin(E)*sqrt(1.0-ecc*ecc)),0.0);\n"//卫星在轨道平面直角坐标系的坐标；

gl_Position=gl_ModelViewProjectionMatrix*vec4(position,1.0);//模型视图投影变换后输出坐标；此外，程序利用自定义的uniform变量将目标大小“传入”顶点着色器程序，通过对gl_PointSize赋值完成对目标大小的设置。

利用GLSL编程对目标颜色进行设置

处理完顶点后，与各个顶点相关联的全部属性就会完全得以确定，经过图元装配和光栅化后，图元将被分解为更小的片元，这些片元对应于目标帧缓冲区的像素。而片元处理器是一个处理片元值及相关联数据的可编程单元，可以用来执行传统的图形操作，本发明通过GLSL编程实现对目标颜色的设置。

如图5左图所示，栅格化后的图元形成的片元组合为方形分布，而真实的成像模型为圆形分布，且由内及外颜色渐浅，为此可以采用公式（10）为片元分配不同的透明度来实现。

$\mathrm{opacity}=1-\frac{{(x-x_c)}^2+{(y-y_c)}^2}{\max \{{(x-x_c)}^2+{(y-y_c)}^2\}}$

其中，(x,y)为片元坐标，(x_c,y_c)为图元中心坐标，如图5右图所示，片元透明度按照与中心片元距离的平方呈线性变化，图元整体由内及外透明度逐渐降低，颜色渐浅，符合目标成像模型。程序利用自定义的uniform变量将目标颜色“传入”片元着色器程序，结合透明度变化通过对gl_FragColor赋值完成对目标颜色的设置。

创建GLSL着色器对象并把它们链接起来创建可执行着色器程序所需要的步骤如图6所示。

①创建着色器对象；

②成功编译着色器代码；

③创建一个着色器程序，并将适当的着色器对象链接到程序中；

④使用着色器进行顶点或片段处理。

经过这样一个步骤后才能启动图形硬件进行加速处理，并最终高效地完成坐标解算和属性设置。

3.空间目标的实时模拟

本发明运用OpenGL三维图形库和Qt编程语言实现了目标的实时模拟绘制，首先利用现有的成熟技术绘制出具有准确时空关系和亮度等级的星空背景，然后根据空间目标轨道根数利用GPU实时解算空间目标的位置，设置大小及颜色，如图1所示，在管线绘制的过程中，图元装配、光栅化、深度测试等操作仍然由固定功能的管线来完成。因此，还需要根据视点位置完成场景的裁剪，最终将目标投影到高性能的显示器上，完成空间目标模拟。

本发明用于空间态势表达的空间目标实时模拟方法的实验验证：

为了验证本发明提出的目标实时模拟方法的可行性和适用性，以“空间态势表达”为背景，对银河系小行星带以及近地空间不同类型航天器进行实时模拟试验验证：在以本发明方法为核心的模拟器***中输入时间、设置目标大小及颜色、调整视点位置等参数后，由工作站生成相应的光学信号并显示在高分辨率显示器上，其结果如下：图7为近地空间近万颗不同类型航天器的实时模拟结果，试验对不同类别的航天器目标进行了颜色和大小的区分。图8是太阳系小行星带的模拟结果，增大仿真步长，小行星目标的位置随时间不断发生变化，展示流畅。

从试验结果可以看出，本发明用于空间态势表达的空间目标实时模拟方法速度快、可靠性强、效果逼真，充分利用了GPU并行、高效的计算特性，能够为未来空间态势的表达提供实时的目标模拟。

Claims (5)

1.一种空间目标实时模拟方法，其特征在于，该模拟方法的步骤如下：

1）计算恒星位置、大小和亮度，并根据恒星像点的成像模型，绘制星空背景；

2）根据空间目标轨道根数利用GPU实时解算空间目标的位置，并分别利用GPU中的顶点着色器和片段着色器设置空间目标的大小和颜色；

3）结合视点位置和用户裁剪信息，完成目标方向的空间态势生成。

2.根据权利要求1所述的空间目标实时模拟方法，其特征在于，所述步骤2）中对空间目标位置实时解算及属性设置的具体步骤如下；

A.建立轨道定向四元数与轨道角元素之间的关系，利用四元数对空间目标的轨道姿态进行描述；

B.将空间目标内置属性中的轨道根数与顶点位置、表面法线和纹理坐标进行关联；

C.利用OpenGL着色语言对顶点着色器进行编程，实现对空间目标位置的解算，并设置其大小；

D.对空间目标的图元进行装配和光栅化处理，将图元分解为更小的片元，并为片元分配不同的透明度，结合透明度变化对目标颜色进行设置。

3.根据权利要求2所述的空间目标实时模拟方法，其特征在于，所述步骤A中轨道定向四元数与轨道角元素之间的关系为：

$q_x=\cos \frac{i}{2}\cos \frac{\mathrm{\&Omega;}+f}{2},q_y=\sin \frac{i}{2}\cos \frac{\mathrm{\&Omega;}-f}{2}$

$q_z=\sin \frac{i}{2}\sin \frac{\mathrm{\&Omega;}-f}{2},q_w=\cos \frac{i}{2}\sin \frac{\mathrm{\&Omega;}+f}{2}$

其中q_x,q_y,q_z,q_w为轨道姿态的四元数，Ω为升交点赤经，i为轨道倾角，σ＝ω+f为升交点角距，ω为近地点角距，f为真近点角。

4.根据权利要求2所述的空间目标实时模拟方法，其特征在于，所述步骤C中利用轨道根数解算空间目标的天球坐标系坐标的过程如下：

a.由历元平近点角M₀、平均角速度n和当前历元时间time，计算平近点角，

M＝M₀+time*n；

b.由偏心率e，依据开普勒方程M＝E-e sin E，迭代求解偏近点角E；

c.由长半轴a，计算目标在轨道平面直角坐标系的坐标，

$\left(\begin{array}{c}x\\ y\end{array}\right)=a\left(\begin{array}{c}\cos E-e\\ \sqrt{(1-e^2)\sin E}\end{array}\right);$

d.利用轨道姿态四元数计算目标在天球坐标系中的坐标，

${\left(\begin{array}{c}X\\ Y\\ Z\end{array}\right)}_{\mathrm{CS}}=\left(\begin{array}{c}{\mathrm{xq}}_x^2+{\mathrm{yq}}_x{q}_y+z{q}_x{q}_z+a_0{q}_w-a_1{q}_z+a_2{q}_y\\ {\mathrm{xq}}_x{q}_y+{\mathrm{yq}}_y^2+{\mathrm{zq}}_y{q}_z+a_0{q}_z+a_1{q}_w-a_2{q}_x\\ {\mathrm{xq}}_x{q}_z+{\mathrm{yq}}_y{q}_z+{\mathrm{zq}}_z^2-a_0{q}_y+a_1{q}_x+a_2{q}_z\end{array}\right)$

其中 $\left(\begin{array}{c}{a}_0\\ a_1\\ a_2\end{array}\right)=\left(\begin{array}{c}{\mathrm{xq}}_w-{\mathrm{yq}}_z\\ {\mathrm{xq}}_z-{\mathrm{yq}}_w\\ {\mathrm{yq}}_x-{\mathrm{xq}}_y\end{array}\right)$ 为过渡矩阵，q_x,q_y,q_z,q_w为轨道姿态四元数。

5.根据权利要求2所述的空间目标实时模拟方法，其特征在于，所述步骤D中为片元分配的透明度模型如下：

$\mathrm{opacity}=1-\frac{{(x-x_c)}^2+{(y-y_c)}^2}{\max \{{(x-x_c)}^2+{(y-y_c)}^2\}}$

其中，(x,y)为片元坐标，(x_c,y_c)为图元中心坐标。

Images