基于光度曲线对空间高轨小尺寸目标特征识别方法与流程

本发明属于光学观测目标特征识别领域，涉及空间高轨航天器特征识别的方法，特别涉及一种基于光度曲线实时最优地对空间高轨小尺寸目标特征识别与运动估计的收敛方法。

背景技术：

空间高轨小尺寸目标指的是高轨航天器，对其特征识别是当前空间态势感知的难点和重点之一，例如空间非合作目标的特征识别，对空间失效或异常合作目标的特征获取与诊断。其难点在于目标尺寸小，而且距离地面数万公里，常用雷达受到探测能力以及噪声的限制，无法实现对高轨目标的探测；而光学传感器受到高轨目标距离的限制无法获得目标的图像，仅可以获得目标亮度变化的光度曲线数据，既有设备难以直接观测到目标特征。由于现有观测手段的限制，近年来出现一种新的途径，基于目标对太阳光线反射的信号与目标的运动以及众多特征相关，因而基于目标的光度观测成为一种可行的空间高轨目标特征识别的方法，典型的光度观测数据曲线如图1所示，其单位为星等。

基于光度曲线获取目标的特征较早出现在对空间星体的分析，这种方式主要特征是对光度观测数据采用离线处理以获得目标特征，相比于传统方法离线反演的途径分析目标特性，这种离线处理的方法也被用于对高轨目标进行特征分析。考虑到空间目标相对于空间星体距离地面要近很多，所以这种离线处理的方法存在的主要缺点是目标的特征与其运动相耦合，目标特征受到其运动的影响，由此导致其光度观测曲线特征受到目标运动的影响，在观测数据中产生耦合效应，在进行目标特征识别的同时还需要进行相关运动信息的估计。同时高轨目标存在的多种特征之间存在耦合关系，对这些特征进行识别是基于离线分析的方法难以实现的。而且由于目标在运动之中离线处理的方法难以实时获取目标的特征，也难以实时获得目标的运动信息。

相对于离线处理方法用于高轨目标特征识别的缺点，近年来基于光度观测数据实时处理的方式受到重视，其特征是采用递推滤波技术实时获得目标运动信息与特征信息。在公开发表的文献中目前可实时获取目标特征信息仅有的工作由纽约州立大学johnl.crassidis教授与合作者开展，他们2014年在aiaajgcd发表的文献表明，在运动学模型的基础上设计模型集，可以通过多模型方法实现对目标位置、速度、姿态与形状、尺寸等联合估计。其基础是对目标运动信息的估计和对目标特性的认识，基于对光压的缺点是算法难以收敛，在此基础上实现对目标多种特征的有效识别则更难实现。对目标运动信息估计的工作代表性的有2009年美国空军实验室在aiaajgcd发表论文基于光度观测数据对垂直与轴向旋转的三级火箭发动机残骸的姿态角估计，申请者在2017年在《光学学报》发表的文章中在随机建模基础上对一般的卫星目标进行的姿态角与角速度联合估计。

原有方法对高轨目标特征识别方法建立在运动学模型基础上，光压项对算法性能至关重要，但是目标的动力学特性中光压项对目标轨道和姿态的影响很微弱，而且在较短时间段内难以体现其对目标轨迹和姿态产生影响，对于这种实时递推算法难以对这种光压项产生有效估计，从而导致难以有效估计目标的位置和速度，进而影响了对目标特征的识别。

技术实现要素：

本发明的目的是提供一种基于光度曲线对空间高轨小尺寸目标特征识别方法，以克服现有方法的不足，本发明的方法可以快速实现算法收敛并实现特征识别。

本发明通过以下技术方案来实现上述目的：

一种基于光度曲线对空间高轨小尺寸目标特征识别方法，包括以下步骤：

步骤1、计算在地面站观测高轨目标对太阳光反射的光度信息，考虑太阳、卫星、地面站之间的相对运动关系，并结合目标的各个面对太阳光线的反射以及姿态运动对观测数据的影响，计算得到对目标的光度观测数据；

步骤2、构造目标特征识别的动态模型，在目标姿态估计基础上开展对目标特征识别，动态模型由目标的姿态运动模型构成；

步骤3、建立由n个模型构成的模型集m＝{m1,m2,…,mn}，在动态模型的基础上基于模型集中目标可能特征构成的模型分别构造该多模型算法的成员滤波器，根据观测数据进行滤波，得到每个滤波器在当前时刻k的输出；

步骤4、计算在当前时刻每个成员滤波器的后验概率密度，第l个模型对应的模型概率

步骤5、对步骤2的结果递推计算第个模型当前时刻的模型概率；

步骤6、根据步骤3的递推式，多模型算法运行之后，根据具有最小似然函数或者最小量测协方差的模型会有最大的模型概率，判断与真实模型最接近或者为真实模型，完成了模型识别；

步骤7、通过多模型自适应算法识别出目标特征，以及实现目标运动状态的融合估计，最后的估计结果为每个滤波器估计结果的加权平均。

步骤1中，空间目标由有限个平面组成，观测站测到的光度大小为各个面反射到地面光测站光度大小的总和；光度值mapp的计算模型为

其中，csun,vis＝455w/m²为可见光照射到单位平方米目标上面的功率，fobs为太阳经卫星反射至地面的辐射量；

其中，asc为辐射面积，为卫星指向太阳的单位向量，为卫星指向地面观测站的单位向量，robs为卫星到地面观测站的距离，fsun为太阳光辐射至卫星表面的辐射量；

为地球惯性坐标系下卫星表面的单位法向量，基于星体坐标系中的法向量通过姿态矩阵a(q)的转换得到，此处四元数q是表示目标姿态的变量；

其中ρtotal为双向反射分布函数部分，与目标材料性能有关，针对brdf建模，采用phong模型，假定光反射量主要由漫反射部分ρdiff和镜面反射部分ρspec构成：

ρtotal＝ρdiff+ρspec。

步骤2中，只考虑目标的姿态运动和目标特征两方面因素，构造的动态模型为：

其中，q为姿态向量，ω为角速度向量，w1,w2为模型噪声。

步骤4中，由多模型方法第个模型对应的模型概率为：

式中，为第个滤波器条件下，k时刻观测值的条件概率；yk为以往k个时刻观测值的集合，即{y0,y1,…,yk}；为第个滤波器的信息，为第滤波器的协方差；

上式中计算为：

步骤5中，第个模型当前时刻的模型概率为：

步骤6中，在k时刻状态估计结果为每个滤波器估计结果的加权融合：

式中为第l个滤波器的状态估计的结果，此处状态向量x即代表目标运动信息，也代表目标的特征信息，此时状态向量的协方差为：

其中，为第l个滤波器的状态协方差。

相对于现有技术，本发明具有以下技术效果：

该发明的有益效果是在姿态估计的基础上实现特征识别，无需考虑轨道动力学模型，因而通过使用比已有方法更简单的动态模型，通过分析目标的动力学特征与光度曲线特征，重新确定了目标的特征与目标运动之间的耦合关系，通过对高轨目标动力学特征的分析，采用了在目标姿态估计基础上建立对目标特征识别的技术路线，本发明与现有技术主要区别在于不使用下面的引力、光压等作用下的轨道模型：

通过本发明的方案不仅实现了对目标形状、尺寸、材料参数等特征有效识别，也实现了对目标的姿态角估计。基于这种更为简单的动态模型构造多模型估计算法实现对高轨目标的特征识别，不仅确保了算法的收敛能力、降低了计算量，更重要是避免现有算法难以收敛的不足，较好确保了算法对目标特征有效识别能力，是实现对目标特征识别与姿态同步最优估计的有效解决方案。

附图说明

图1对空间高轨小尺寸目标观测的典型光度数据曲线

图2算法框架图；

图3角速度误差图；

图4姿态误差图；

图5形状识别结果图；

图6多模型下尺寸估计的姿态与角速度估计结果，(a)角速度误差；(b)姿态误差；

图7多模型下尺寸估计结果；

图8多模型下材料参数估计的姿态、角速度、四元数估计结果，(a)角速度误差；(b)姿态误差；

图9多模型下材料参数估计结果；

图10多模型下尺寸与材料参数识别的姿态与角速度估计结果；(a)角速度误差；(b)姿态误差；

图11多模型下尺寸与材料参数识别结果。

具体实施方式

下面结合具体实施例和附图对本发明作进一步说明：

在多模型框架下基于光度观测通实现对空间高轨小尺寸目标特征识别与运动估计的快速收敛，完整步骤如下：

步骤1以空间高轨航天器为目标，计算在地面观测时候目标对太阳光反射后地面可以观测到的光度信息，考虑太阳、卫星、地面观测站之间的相对运动关系，并结合目标形状与姿态运动对观测值的影响，根据公示可以计算得到光度观测曲线。

空间目标由有限个平面组成，观测站测到的光度大小为各个面反射到地面光测站光度大小的总和。光度值mapp的计算模型为

其中csun,vis＝455w/m²为可见光照射到单位平方米目标上面的功率，fobs为太阳经卫星反射至地面的辐射量。

其中，asc为辐射面积，为卫星指向太阳的单位向量，为卫星指向地面观测站的单位向量，robs为卫星到地面观测站的距离，fsun为太阳光辐射至卫星表面的辐射量。

为地球惯性坐标系下卫星表面的单位法向量，基于星体坐标系中的法向量通过姿态矩阵a(q)的转换得到，如式所示，此处四元数q是表示目标姿态的变量。

其中ρtotal为双向反射分布函数部分，与目标材料性能有关。针对brdf建模，本文采用phong模型，假定光反射量主要由漫反射部分ρdiff和镜面反射部分ρspec构成：

ρtotal＝ρdiff+ρspec

步骤2、构造特征识别的动态模型，此处是本发明与前人发明不同的地方，采取了在目标姿态估计基础上建立对目标特征识别的技术路线，动态模型为：

其中，q为姿态向量，ω为角速度向量，w1,w2为模型噪声。

本发明与现有技术主要区别在于不使用下面的引力、光压等作用下的轨道模型：

其中r为卫星到地心的距离向量，为引力向量，asrp为光压向量。通过本发明的方案不仅实现了对目标形状、尺寸、材料参数等特征有效识别，也实现了对目标的姿态角估计。

图2是通过多模型方法实现对目标特征识别的算法框架图，针对不确定性特征，基于目标可能的模型构造模型集，该图中基于模型集建立的这些成员滤波器为一个并行框架，该算法的输出通过右端的估计器状态表示，其中目姿态与特征二者共同构成了状态向量x，下标k表示在当前时刻实时输出的状态向量，pk|k为多模型滤波器状态估计的协方差，p(mi|yk)表示当前时刻每个模型在算法中所占据的概率，所有模型概率之和为1。

步骤3、建立由n个模型构成模型集m＝{m1,m2,…,mn}，此处的模型集可根据先验信息(形状、尺寸、材料)构造，其中包含了目标可能的模型；基于模型集中每个模型分别构造该多模型算法的成员滤波器，根据观测数据进行滤波，得到每个滤波器在当前时刻k的输出。

步骤4、计算在当前时刻每个成员滤波器的后验概率密度，第个模型对应的模型概率：

式中，为第个滤波器条件下，k时刻观测值的条件概率。yk为以往k个时刻观测值的集合，即{y0,y1,…,yk}。为第个滤波器的信息，为第滤波器的协方差。

上式中计算为：

步骤5、根据多模型方法，目标可能的特征由所对应模型的概率所代表，概率越高表示这种特征的可能性越高，对步骤2的结果递推计算第l个模型当前时刻的模型概率：

步骤6、根据步骤3的递推式，多模型算法运行之后，具有最小似然函数或者最小量测协方差的模型会有最大的模型概率，即是与真实模型最接近或者就是真实模型，完成了模型识别。

步骤7、多模型自适应算法不仅能识别出目标特征外，还可以实现目标运动状态的估计，算法在k时刻状态估计结果为每个滤波器估计结果的融合：

式中为第l个滤波器的状态估计的结果，此处x即包含目标运动信息，也包含目标特征信息，此时状态向量的协方差为：

其中，为第l个滤波器的状态协方差。

实施例1

结合对一个四棱柱卫星形状识别的具体例子对本发明做进一步描述，实现基于光度曲线对空间高轨小尺寸目标特征识别算法的快速收敛。基本的仿真环境设定为：目标选地球同步轨道卫星，轨道倾角为30°，卫星在经度120°处；地面观测站设在东经108°、北纬34°，高度0km处；除了可计算数据或者易于获得性数据，其他数据都由stk给出。除了涉及到模型集设计需要改变对应参数以外，一般情况下每个面的面积为60m²；rdiff＝0.4，rspec＝0.5。仿真开始时间为2015年5月22日05:00:00ut，仿真结束时间为2015年5月22日07:00:00ut。

根据本发明算法，特征识别的步骤如下：

第一步，根据给定的目标轨道与姿态信息计算得到观察时段内的光度观测数据；

第二步，根据下面公式建立目标的姿态运动学模型，其运动的不确定性由模型的噪声水平决定，该值的大小受到目标高轨动力学特征决定；

第三步，建立算法的模型集

目标可能的特征通常可以通过先验信息获知，模型集建立在先验信息的基础上，模型集代表了目标可能的特征，例如在步骤3中，对于形状不确定棱柱状目标，其模型集由四个可能形状对象的法向量所代表。

通过多模型自适应估计办法实现形状识别，而典型棱柱形目标可能的形状为四棱柱、六棱柱、八棱柱、十二棱柱，所以模型集的设计就是结合目标的尺寸、材料参数等参数的设计四个可能对象构成模型集。形状建模主要是定义模型下每个面的法向量；除了与z轴垂直的两个面容易确定外，主要是其他面的法向量的定义，以下分别为四棱柱、六棱柱、八棱柱、十二棱柱形状模型的各个面的法向量，其中四棱柱的各个面法向量为：

六棱柱的各个面法向量为

八棱柱的各个面法向量为

十二棱柱各个面法向量为

此处，

则本例中算法模型集由上述四种可能的成员模型构成m＝{m1,m2,m3,m4}。

第四步，计算在以往k个观测值条件下，每个单一滤波器的后验概率密度，第个模型对应的模型概率

第五步，基于2结果计算第个模型在k时刻的模型概率。

第六步，根据步骤3的递推式，多模型算法运行之后，具有最小似然函数或者最小量测协方差的模型会有最大的模型概率，即是与真实模型最接近或者就是真实模型，完成了模型识别。

第七步，通过多模型自适应算法识别出目标特征，以及实现目标运动状态的融合估计，系统最后的估计结果为每个滤波器估计结果的加权平均。

在上例仿真环境下，光度数据来自于四棱柱形状模型，即四棱柱为真实模型。图3，图4分别为角速度估计误差和姿态估计误差，误差都收敛到3σ范围内；初始情况下，各个模型的模型概率都设置为1/4，每个模型的权重在多模型算法初始时刻都相同，图5为多模型估计所得到的各个模型的概率图，在刚开始0～40时刻，模型二的概率最大，表现为是最可能的真实模型，但是之后很快真实模型即四棱柱的概率收敛到1，其他模型的概率收敛到0，实现了对形状的识别。

实施例2

尺寸估计：

通过多模型自适应估计办法实现同一形状模型下的不同尺寸的估计，所以模型集的设计就是关于尺寸参数的设计，其他参数都与单模型下保持一致。形状选为四棱柱，每个模型集由5个可能尺寸构成，尺寸参数用形状中六个面的面积来表征，真实模型为a＝[606060606060]，此将考虑真实模型包含在模型集中。

真实模型包含在模型集中的模型集设计为：

模型一：a＝[606060606060]

模型二：a＝[60601206012060]

模型三：a＝[60606012060120]

模型四：a＝[303030303030]

模型五：a＝[120120120120120120]

在仿真当中，光度数据来自于模型一，即模型一为真实模型；图6分别为角速度估计误差和姿态估计误差，其估计误差都收敛到3σ范围内；图7为多模型估计所得到的各个模型的概率图，在刚开始0～80时刻，模型三和模型一的概率最大，表现为是最可能的真实模型，但是之后很快模型一的概率收敛到1，其他模型包括模型三的概率收敛到0，实现了对尺寸的识别。

实施例3

材料参数估计：

通过多模型自适应估计办法实现同一形状模型下的材料参数的估计，所以模型集的设计就是关于材料参数的设计，其他参数都与单模型下保持一致。形状模型选为四棱柱形状，模型集设计5个尺寸模型，因为本文选用a-s模型进行量测建模，所以材料参数选择a-s模型中的参数，此处选择nv和nu进行模型集设计。

模型一：nu＝[101010101010]，nv＝[101010101010]

模型二：nu＝[101010101010]，nv＝[555555]

模型三：nu＝[555555]，nv＝[101010101010]

模型四：nu＝[555555]，nv＝[555555]

模型五：nu＝[202020202020]，nv＝[555555]

在仿真环境中，光度数据来自于模型一，即模型一为真实模型，模型一的材料参数为nu＝[101010101010]、nv＝[101010101010]；图8为角速度估计误差、姿态估计误差，误差都收敛到3σ范围内；初始情况下，各个模型的模型概率都设置为1/5，每个模型的权重在多模型算法初始时刻都相同，图9为多模型估计所得到的各个模型的概率图，在刚开始0～50时刻，模型五和模型一的概率最大，表现为是最可能的真实模型，之后在50～370时刻之间，模型五和模型一的概率交叉增大，但在380时刻之后模型一的概率逐渐收敛到1，其他模型包括模型五的概率收敛到0，实现了对材料参数的估计。

实施例4

尺寸与材料参数估计：

通过多模型自适应估计办法实现同一形状模型下，尺寸和材料参数的同步估计，所以模型集的设计就是关于材料参数和尺寸参数的组合设计，其他参数都与单模型下保持一致。形状模型选为四棱柱形状，模型集设计10个尺寸和材料参数的组合模型，模型集设计办法采用交叉排列组合，尺寸模型使用上文真实模型包含在模型集中的尺寸设计结果，五个尺寸模型分别定义为modelsize_1、modelsize_2、modelsize_3、modelsize_4、modelsize_5；材料参数模型设计采用上文设计结果，五个材料参数分别定义为modelmaterial_1、modelmaterial_2、modelmaterial_3、modelmaterial_4、modelmaterial_5；所以尺寸和材料参数组合的模型集设计为：

模型一：modelsize_1+modelmaterial_1

模型二：modelsize_1+modelmaterial_4

模型三：modelsize_2+modelmaterial_2

模型四：modelsize_2+modelmaterial_5

模型五：modelsize_3+modelmaterial_3

模型六：modelsize_3+modelmaterial_1

模型七：modelsize_4+modelmaterial_4

模型八：modelsize_4+modelmaterial_2

模型九：modelsize_5+modelmaterial_5

模型十：modelsize_5+modelmaterial_3

基于上述模型集实现对未知目标姿态估计与材料参数识别的结果分别见图10和图11。

本发明的保护范围并不限于上述的实施例，对于本领域的普通技术人员来说，倘若对本发明进行的各种改动和变形属于本发明权利要求及等同技术范围内，则本发明的意图也包含这些改动和变形在内。