CN106548475A - 一种适用于空间非合作自旋目标运动轨迹的预测方法 - Google Patents

Abstract

本发明涉及一种适用于空间非合作自旋目标运动轨迹的预测方法，首先获得相机坐标系下目标位置坐标，由已知的相机姿态，求解惯性坐标系下目标位置坐标。构造NAR神经网络，通过BPTT算法对神经网络进行训练，偏差收敛后对目标位置值给出预测输出，以保证系统鲁棒性与采样连续性。位置信息同时用于计算姿态变换四元数，同时根据估计参数通过运动学和动力学方程计算姿态变换四元数，有最小二乘回归估计参数值，通过方程计算预测结果。该方法用于空间在轨服务任务中，快速准确的跟踪目标轨迹并在参数收敛后，获得长时间的目标轨迹预测信。

Description

一种适用于空间非合作自旋目标运动轨迹的预测方法

技术领域

本发明属于信息处理与信息融合领域，涉及一种适用于空间非合作自旋目标运动轨迹的预测方法。

背景技术

空间智能机器人具有视觉辨识与运动控制系统，可完成对合作或非合作目标的在轨捕获与释放，执行如轨道垃圾清理、航天器在轨维护等空间任务。为完成对目标的在轨捕获任务，在接近与停靠目标的过程中，目标相对空间智能机器人往往存在一定的线速度与角速度，由于对信息处理以及位姿解算需要耗费一定时间，因此实时追踪目标运动轨迹并对运动目标的轨迹进行预测是十分必要的，也是开展在轨服务的仿真与试验验证基础。

前期的研究工作以及地面验证试验表明，在目标位置坐标获得过程中，可能出现的相机抖动、光照变化、目标漂浮及遮挡等影响，造成采样时间较短，采样过程不连续及采样噪声较大等问题。因此，设计基于信息测量，对运动目标的轨迹进行追踪与预测具有较强鲁棒性的方法便十分必要。

已有的方法包括双卡尔曼滤波、牛顿拉格朗日插值等，双卡尔曼滤波因为每一个迭代周期都需要对惯性参量重新估算，增加了系统的运算负荷，降低了效率。牛顿拉格朗日插值等采用外插法进行预测，这使得预测精度很难有效提高。

发明内容

要解决的技术问题

为了避免现有技术的不足之处，本发明提出一种适用于空间非合作自旋目标运动轨迹的预测方法，针对具备视觉观测系统的空间智能机器人对观测目标轨迹的跟踪与预测的问题。

技术方案

一种适用于空间非合作自旋目标运动轨迹的预测方法，其特征在于步骤如下：

步骤1：由视觉相机获得的非合作目标特征点相对测量平台的位置坐标，对此位置坐标以及已知的测量平台的位姿信息进行坐标变换解算目标的惯性系坐标；

步骤2：构造非线性递归神经网络NAR，以步骤1解算出的目标惯性系坐标，对NAR模型进行训练，得到目标位置坐标短期预测值；

设计最小二乘参数估计器LSE，对非合作目标的运动状态以及动力学模型参数进行估计，在获得动力学模型参数估计值的基础上，通过动力学模型运算实现目标运动轨迹预测；

当出现目标位置因振动及遮挡时，采用NAR网络预测值作为最小二乘参数估计器LSE的输入，在参数估计收敛后，实现目标运动轨迹预测。

所述步骤2构造非线性递归神经网络NAR以及的方法为学习算法的方法：

步骤a、设计网络结构：在多层感知器引入反馈连接，反馈连接有一组联系单元context，记忆网络过去的状态，并在下一时刻连同网络输入一起作为隐含层输入；

步骤b、设计递归网络学习算法：设计截断历时反向传播学习算法BPTT(h)，计算最小化的目标函数定义为瞬时误差的平方和：

其中，E(k)为k时刻误差平方和，e_j(k)为神经元j在k时刻输出值与真实值偏差。

计算BPTT(h)对神经元j的局部梯度：

对于j∈A，其中δ_j(l)为第l层神经元j梯度，net_j(l)为第l层神经元j的输出值。

计算过程如下：

第一步，对区间(k-h,k)内各时刻数据，按信号流逐层进行前向计算，保存输入数据和网络权值；

第二步，应用上述前向计算结果，执行逐层反向传播过程，计算公式为

其中，w_kj(l)为k时刻l层神经元j的权值。

上述计算反向进行到时刻k-h+1时，所计算步数等于h；

第三步，当反向传播计算到k-h+1时，对神经元j的权值进行如下调整

其中，η为学习率，Δw_ij(k)为权值调整量，x_i(l-1)为输入数据；

第四步，通过训练好的网络，对目标位置进行循环短时预测。

所述步骤2设计最小二乘参数估计器LSE，对非合作目标的运动状态以及动力学模型参数进行估计，实现目标运动轨迹预测的方法为；

1、根据惯性系下非合作目标的位置坐标P(x_w,y_w,z_w,t-1)，P(x_w,y_w,z_w,t)，估算姿态变化：

其中nn_out为神经网络预测输出，S为转换闸值。

2、由P(x_w,y_w,z_w,t-1)，P(x_w,y_w,z_w,t)计算变换矩阵D

其中：S_jk＝P_j(x_w,y_w,z_w,t-1)P_k(x_w,y_w,z_w,t)

3、构造二次型函数F(q)＝q^TDq

其中q为四元数描述的目标姿态，采用最大特征值法求姿态估计值：

4、采用最小二乘辨识系统参量，求得辨识参量k

其中q为四元数，ω为角速度；

采用忽略外力作用下系统动力学方程：

因动力学公式中惯性参量无法单独辨识，所以定义：

其中，惯性参量满足以下条件I_xx+I_yy＞I_zz,I_yy+I_zz＞I_xx,I_zz+I_xx＞I_yy，则改写为

其中

取待辨识参量：k＝[p_x p_y p_z ω₁(0) ω₂(0) ω₃(0)]^T

在给出k初值后，采用系统运动学方程和忽略外力作用下系统动力学方程构造姿态偏差：

其中，q_v＝[q₁ q₂ q₃]，为四元数向量部分

非线性最优解为：

5、由辨识出的系统特征参量k，进行系统轨迹预测

在已知k条件下，可积分求解下式：

计算旋转矩阵估计值：

由旋转矩阵估计值计算

有益效果

本发明提出的一种适用于空间非合作自旋目标运动轨迹的预测方法，利用数据的局部性与全局性，针对每一时刻检测获得的数据，采用最小二乘估计对参数进行非线性回归，获得参数的估计值，进而通过运动学与动力学方程运算得到目标轨迹预测值。针对数据的全局性，构造非线性递归神经网络，通过大量数据与在线训练，获得数据内在的周期性规律，解决参数估计过程中因噪声所产生的数据丢失问题。该方法用于空间在轨服务任务中，快速准确的跟踪目标轨迹并在参数收敛后，获得长时间的目标轨迹预测信。

本发明具体拥有以下优点：

1.预测效果好，结果准确性高

本发明基于信息融合技术，结合了最小二乘估计法对轨迹预测的长期性与非线性递归神经网络对轨迹预测的快速性，解决了拉格朗日插值法等因对信息内在联系未建模造成的预测偏差较大问题。

2.预测时间长

本发明通过对系统参数的辨识获得了系统的内部特性，根据系统运动学与动力学方程，可对目标运动轨迹做出长时间的准确预测。

3.具有良好的鲁棒性，能够应用于多种复杂场景，抗噪声性能比较强

本发明针对位置数据获取过程中可能存在的振动、遮挡等问题，创新的提出了融合非线性递归神经网络方法，该方法的短期准确预测能力，使得数据获取过程不会因噪声而产生数据中断、采样间隔随机等问题，增强了系统的鲁棒性与适应性。

4.计算量小，计算速度快，能够满足实时性要求

本发明中对数据的处理采用了对微分方程的一阶近似，数据存储与计算同步，短时与长时预测并行，在满足精度要求的前提下，大大减小了计算量，提高了运行速度。

附图说明

图1为本发明的算法流程意图

图2为NAR神经网络结构图

图3为参数辨识结果

图4为位置预测结果

具体实施方式

现结合实施例、附图对本发明作进一步描述：

如图1所示，本发明实施例提供的应用于空间绳系机器人的抓捕点附近障碍物的检测方法包括：

步骤一：相机坐标系下目标位置解算

设相机坐标系中目标位置为P(x_c,y_c,z_c),目标在相机坐标系中成像位置(u,v)，则有以下公式成立：

式中(u₀,v₀)为图像坐标系的中心，(a_x,a_y)为有效焦距长度。

步骤二：惯性坐标系下目标位置解算

待求惯性系下坐标P(x_w,y_w,z_w)便可以表示为：

其中旋转矩阵R、平移矩阵T分别为：

T＝[t_x,t_y,t_z]^T (3)

将公式(3)，(1)代入(2)，可得到

其中B＝[x_w,y_w,z_w,1]^T，待求惯性系下的坐标P(xw,yw,zw)可通过对相机坐标系中成像位置(u,v)的实时检测及公式(4)解算获得，坐标点的集合便构成目标运动轨迹。

步骤三：NAR网络结构设计

为使神经网络具有预测感知能力，在神经网络设计过程中对多层感知器引入了反馈连接，反馈连接有一组联系单元(context)，它记忆网络过去的状态，并在下一时刻连同网络输入一起作为隐含层输入，这种网络的反馈连接规定了网络的内部状态和记忆形式，使得其输入不仅依赖当前输入，也与过去的输入有关，这一性质使得部分反馈网络具有动态记忆与预测感知能力。根据预测问题的性质，网络结构选择为内部时延的反馈递归网络。具体网络连接与结构见附图2

步骤四：NAR学习算法设计

选择递归网络学习算法，由于只在一个固定的时间序列内存储相关的输入数据和网络状态历史记录，选择截断历时反向传播学习算法(BPTT(h)),该时间序列的长度h称截断深度，任何h时刻之前的信息均不需要存储。最小化的目标函数定义为瞬时误差的平方和。

BPTT(h)对神经元j的局部梯度定义为

对于j∈A有k-h＜l≤k

该算法按下述步骤进行：

1)对区间(k-h,k)内各时刻数据，按信号流逐层进行前向计算，保存输入数据和网络权值。

2)应用上述前向计算结果，执行逐层反向传播过程，计算公式为

上述计算反向进行到时刻k-h+1时，所计算步数等于h。

3)当反向传播计算到k-h+1时，对神经元j的权值进行如下调整

其中，η为学习率。

步骤五：预测结果输出

Y_t+1＝f(X_c(t),X_c(t-1),…X_c(t-m),U_t) (9)

式中，m为隐层最大时延值，Xc为隐层反馈；Ut为系统输入。

步骤六：计算(t-1,t)姿态变化

根据P(x_w,y_w,z_m,t-1)，P(x_w,y_w,z_m,t)，考虑到振动、遮挡等目标位置不可测状态，估算姿态变化：

其中nn_out为神经网络预测输出。

由P(x_w，y_w，z_w，t-1)，P(x_w，y_w，z_w，t)计算变换矩阵D

其中

S_jk＝P_j(t-1)P_k(t) (12)

构造二次型函数

F(q)＝q^TDq (13)

其中q为四元数描述的目标姿态，采用最大特征值法求姿态估计值：

步骤七：由参数估计姿态变化

系统运动学方程：

忽略外力作用下系统动力学方程：

其中q为四元数，ω为角速度。

因动力学公式中惯性参量无法单独辨识，所以作如下定义：

其中，惯性参量满足以下条件

I_xx+I_yy＞I_zz,I_yy+I_zz＞I_xx,I_zz+I_xx＞I_yy

则公式可改写为

其中

步骤八:辨识特征参数

取待辨识参量

k＝[p_x p_y p_z ω₁(0) ω₂(0) ω₃(0)]^T (19)

在给出K初值后，由公式(15)(16)构造姿态偏差：

其中，q_v＝[q₁ q₂ q₃]，

非线性最优解定义为：

步骤九：运动学动力学解算预测值

由公式(15)，(16)，在已知k条件下，可积分求解

计算旋转矩阵估计值：

可由旋转矩阵估计值计算P(x_w,y_w,z_m,t+1)

至此，便可获得轨迹预测值。

Claims (3)

1.一种适用于空间非合作自旋目标运动轨迹的预测方法，其特征在于步骤如下：

步骤1：由视觉相机获得的非合作目标特征点相对测量平台的位置坐标，对此位置坐标以及已知的测量平台的位姿信息进行坐标变换解算目标的惯性系坐标；

步骤2：构造非线性递归神经网络NAR，以步骤1解算出的目标惯性系坐标，对NAR模型进行训练，得到目标位置坐标短期预测值；

设计最小二乘参数估计器LSE，对非合作目标的运动状态以及动力学模型参数进行估计，在获得动力学模型参数估计值的基础上，通过动力学模型运算实现目标运动轨迹预测；

当出现目标位置因振动及遮挡时，采用NAR网络预测值作为最小二乘参数估计器LSE的输入，在参数估计收敛后，实现目标运动轨迹预测。

2.根据权利要求1所述适用于空间非合作自旋目标运动轨迹的预测方法，其特征在于：所述步骤2构造非线性递归神经网络NAR以及的方法为学习算法的方法：

步骤a、设计网络结构：在多层感知器引入反馈连接，反馈连接有一组联系单元context，记忆网络过去的状态，并在下一时刻连同网络输入一起作为隐含层输入；

步骤b、设计递归网络学习算法：设计截断历时反向传播学习算法BPTT(h)，计算最小化的目标函数定义为瞬时误差的平方和：

$E\left(k\right)=\frac{1}{2}\underset{j\∈A}{\Σ}{e}_j^2\left(k\right)$

其中，E(k)为k时刻误差平方和，e_j(k)为神经元j在k时刻输出值与真实值偏差。

计算BPTT(h)对神经元j的局部梯度：

${\mathrm{\δ}}_j\left(l\right)=-\frac{\∂E\left(l\right)}{\∂{\mathrm{net}}_j\left(l\right)}$

对于j∈A，其中δ_j(l)为第l层神经元j梯度，net_j(l)为第l层神经元j的输出值。

计算过程如下：

第一步，对区间(k-h,k)内各时刻数据，按信号流逐层进行前向计算，保存输入数据和网络权值；

第二步，应用上述前向计算结果，执行逐层反向传播过程，计算公式为

${\mathrm{\&delta;}}_j\left(l\right)=\left\{\begin{array}{cc}{f}^{\&prime;}\left({\mathrm{net}}_j\right(l\left)\right)e_j\left(l\right)& l=k\\ f^{\&prime;}\left({\mathrm{net}}_j\right(l\left)\right)\underset{k\&Element;A}{\&Sigma;}{w}_{kj}\left(l\right){\mathrm{\&delta;}}_k(l+1)& k-h<l<k\end{array}\right.$

其中，w_kj(l)为k时刻l层神经元j的权值。

上述计算反向进行到时刻k-h+1时，所计算步数等于h；

第三步，当反向传播计算到k-h+1时，对神经元j的权值进行如下调整

${\mathrm{\&Delta;w}}_{ij}\left(k\right)=\mathrm{\&eta;}\underset{l=k-h+1}{\overset{k}{\&Sigma;}}{\mathrm{\&delta;}}_j\left(l\right)x_i(l-1)$

其中，η为学习率，Δw_ij(k)为权值调整量，x_i(l-1)为输入数据；

第四步，通过训练好的网络，对目标位置进行循环短时预测。

3.根据权利要求1所述适用于空间非合作自旋目标运动轨迹的预测方法，其特征在于：所述步骤2设计最小二乘参数估计器LSE，对非合作目标的运动状态以及动力学模型参数进行估计，实现目标运动轨迹预测的方法为；

1、根据惯性系下非合作目标的位置坐标P(x_w,y_w,z_w,t-1)，P(x_w,y_w,z_w,t)，估算姿态变化：

$P(x_w,y_w,z_w,t)=\left\{\begin{array}{cc}P(x_w,y_w,z_w)& \left|\right|P(x_w,y_w,z_w,t)-P(x_w,y_w,z_w,t-1)\left|\right|<S\\ {\mathrm{nn}}_{out}& else\end{array}\right\}$

其中nn_out为神经网络预测输出，S为转换闸值。

2、由P(x_w,y_w,z_w,t-1)，P(x_w,y_w,z_w,t)计算变换矩阵D

$D=\left[\begin{array}{cccc}{S}_{11}+S_{22}+S_{33}& S_{23}-S_{32}& S_{31}-S_{13}& S_{12}-S_{21}\\ S_{23}-S_{32}& S_{11}-S_{22}-S_{33}& S_{12}+S_{21}& S_{31}+S_{13}\\ S_{31}-S_{13}& S_{12}+S_{21}& -S_{11}+S_{22}-S_{33}& S_{23}+S_{32}\\ S_{12}-S_{21}& S_{31}+S_{13}& S_{23}+S_{32}& -S_{11}-S_{22}+S_{33}\end{array}\right]$

其中：S_jk＝P_j(x_w,y_w,z_w,t-1)P_k(x_w,y_w,z_w,t)

3、构造二次型函数F(q)＝q^TDq

其中q为四元数描述的目标姿态，采用最大特征值法求姿态估计值：

$\hat{q}=\arg \underset{\left|\right|q\left|\right|=1}{max}F\left(q\right)$

4、采用最小二乘辨识系统参量，求得辨识参量k

$\stackrel{\&CenterDot;}{q}=\left[\begin{array}{c}{\stackrel{\&CenterDot;}{q}}_0\\ {\stackrel{\&CenterDot;}{q}}_1\\ {\stackrel{\&CenterDot;}{q}}_2\\ {\stackrel{\&CenterDot;}{q}}_3\end{array}\right]=\frac{1}{2}\left[\begin{array}{c}-{\mathrm{\&omega;}}_1{q}_1-{\mathrm{\&omega;}}_2{q}_2-{\mathrm{\&omega;}}_3{q}_3\\ {\mathrm{\&omega;}}_1{q}_0+{\mathrm{\&omega;}}_3{q}_2-{\mathrm{\&omega;}}_2{q}_3\\ {\mathrm{\&omega;}}_2{q}_0-{\mathrm{\&omega;}}_3{q}_1+{\mathrm{\&omega;}}_1{q}_3\\ {\mathrm{\&omega;}}_3{q}_0+{\mathrm{\&omega;}}_2{q}_1-{\mathrm{\&omega;}}_1{q}_2\end{array}\right]$

其中q为四元数，ω为角速度；

采用忽略外力作用下系统动力学方程：

因动力学公式中惯性参量无法单独辨识，所以定义：

$p_x=\frac{I_{yy}-I_{zz}}{I_{xx}},p_y=\frac{I_{zz}-I_{xx}}{I_{yy}},p_z=\frac{I_{xx}-I_{yy}}{I_{zz}}$

其中，惯性参量满足以下条件I_xx+I_yy＞I_zz,I_yy+I_zz＞I_xx,I_zz+I_xx＞I_yy，则改写为

其中

取待辨识参量：k＝[p_x p_y p_z ω₁(0) ω₂(0) ω₃(0)]^T

在给出k初值后，采用系统运动学方程和忽略外力作用下系统动力学方程构造姿态偏差：

其中，q_v＝[q₁ q₂ q₃]，为四元数向量部分

非线性最优解为：

5、由辨识出的系统特征参量k，进行系统轨迹预测

在已知k条件下，可积分求解下式：

$\mathrm{\&omega;}(t,t+1)=-{\&Integral;}_t^{t+1}{I}^{-1}(\mathrm{\&omega;}\&times;I\mathrm{\&omega;})dt$

$q(t,t+1)={\&Integral;}_t^{t+1}\stackrel{\&CenterDot;}{q}dt$

计算旋转矩阵估计值：

$\hat{R}(t,t+1)=\left[\begin{array}{ccc}{{\hat{q}}_0}^2+{{\hat{q}}_1}^2-{{\hat{q}}_2}^2-{{\hat{q}}_3}^2& 2({\hat{q}}_1{\hat{q}}_2-{\hat{q}}_0{\hat{q}}_3)& 2({\hat{q}}_1{\hat{q}}_3+{\hat{q}}_0{\hat{q}}_2)\\ 2({\hat{q}}_2{\hat{q}}_1+{\hat{q}}_0{\hat{q}}_3)& ({{\hat{q}}_0}^2-{{\hat{q}}_1}^2+{{\hat{q}}_2}^2-{{\hat{q}}_3}^2)& 2({\hat{q}}_2{\hat{q}}_3-{\hat{q}}_0{\hat{q}}_1)\\ 2({\hat{q}}_3{\hat{q}}_1-{\hat{q}}_0{\hat{q}}_2)& 2({\hat{q}}_3{\hat{q}}_2+{\hat{q}}_0{\hat{q}}_1)& ({{\hat{q}}_0}^2-{{\hat{q}}_1}^2-{{\hat{q}}_2}^2+{{\hat{q}}_3}^2)\end{array}\right]$

由旋转矩阵估计值计算

$\hat{P}(x_w,y_w,z_w,t+1)=\hat{R}(t,t+1)*P(x_w,y_w,z_w,t).$