CN113640791A - 一种基于距离和瞬时速度的空间目标三维姿态重构方法 - Google Patents

Abstract

本发明公开了一种基于距离和瞬时速度的空间目标三维姿态重构方法，包括：获取所述雷达接收到的回波信号，所述回波信号由空间目标散射至所述雷达；根据所述回波信号，确定所述空间目标的所有散射中心对应的一维距离序列；通过卡尔曼滤波对所述一维距离序列进行关联处理；根据所述回波信号，估计所述散射中心对应的瞬时速度；基于卡尔曼滤波对所述散射中心对应的多普勒频率进行关联处理；基于关联处理后的一维距离序列和瞬时速度，对散射中心进行三维重构处理。本发明能够实现对空间目标三维姿态的精确重构。

Description

一种基于距离和瞬时速度的空间目标三维姿态重构方法

技术领域

本发明属于雷达技术领域，具体涉及一种基于距离和瞬时速度的空间 目标三维姿态重构方法。

背景技术

空间目标三维姿态能够反映目标的形状以及运动情况，在目标识别中 起着重要作用，因此为了能够准确地得到空间目标的三维姿态信息，需要研 究有效的空间目标三维姿态重构方法。

现有的空间目标三维重构方案通常需要预先确定镜头焦距与空间目标 距离之间的几何关系，且受到图像分辨率、空间目标位置及尺寸参数的影响 较大，操作复杂，尤其是无法针对非平稳运动的空间目标进行聚焦，导致三 维重构准确率低，甚至无法进行三维重构。

洪灵在其发表的论文“基于稀疏重构的空间目标感知方法研究”(西安 电子科技大学，博士学位论文，2015)中提出了一种基于一维高分辨距离向 的空间目标三维姿态重构方法。该方法的实现步骤如下：(1)采用稀疏优化 的方法进行来径向距离估计；(2)采用动态规划的方法实现多散射中心径向 距离历程的关联；(3)基于因式分解进行三维姿态重构。该方法仍然存在的 不足之处是：(1)一维距离像的分辨率受制于雷达带宽；(2)目标散射中心 坐标与真实坐标相差一个任意旋转矩阵。

发明内容

为了解决现有技术中存在的上述问题，本发明提供了一种基于距离和 瞬时速度的空间目标三维姿态重构方法。本发明要解决的技术问题通过以 下技术方案实现：

一种基于距离和瞬时速度的空间目标三维姿态重构方法，应用于电子 设备，所述电子设备对应有多个雷达，所述方法包括：步骤1：获取所述雷 达接收到的回波信号，所述回波信号由空间目标散射至所述雷达；步骤2： 根据所述回波信号，确定所述空间目标的所有散射中心对应的一维距离序 列；步骤3：通过卡尔曼滤波对所述一维距离序列进行关联处理；步骤4： 根据所述回波信号，估计所述散射中心对应的瞬时速度；步骤5：基于卡尔 曼滤波对所述散射中心对应的多普勒频率进行关联处理；步骤6：基于关联 处理后的一维距离序列和瞬时速度，对散射中心进行三维重构处理。

本发明的有益效果：

本发明能够克服距离分辨率受制于信号带宽的问题，以及能够避免重 构出来的目标散射中心坐标与真实坐标相差一个任意旋转矩阵，能够实现 对空间目标三维姿态的精确重构。

以下将结合附图及实施例对本发明做进一步详细说明。

附图说明

图1是本发明实施例提供的一种基于距离和瞬时速度的空间目标三维 姿态重构方法流程示意图；

图2是本发明实施例提供的一种雷达与空间目标之间的方位示意图；

图3是本发明实施例提供的一种空间目标几何形状示意图；

图4是本发明实施例提供的一种基于三部雷达得到的一维距离历程以 及关联后的一维距离历程示意图；

图5是本发明实施例提供的一种基于三部雷达得到的瞬时频率历程以 及关联后的瞬时频率历程示意图；

图6是本发明实施例提供的一种基于t＝0.08s时4个散射中心在参考坐 标系中的重构结果示意图；

图7是本发明实施例提供的一种在不同信噪比条件下A、B两个散射 中心的重构性能比较示意图；

图8是本发明实施例提供的一种电子设备示意图。

具体实施方式

下面结合具体实施例对本发明做进一步详细的描述，但本发明的实施 方式不限于此。

实施例一

请参见图1，图1是本发明实施例提供的一种基于距离和瞬时速度的空 间目标三维姿态重构方法流程示意图，所述方法包括：

步骤1：获取所述雷达接收到的回波信号，所述回波信号由空间目标散 射至所述雷达。

所述电子设备节点指的是指经由通信设施向计算机输入程序和数据或 接收计算机输出结果的设备，广义讲所述可移动电子设备包括但不局限于 服务器、个人计算机、笔记电脑，无线电定位设备等，所述电子设备中部署 有数据处理软件，如matlab仿真软件，cst电磁仿真软件等，所述电子设备 能够获取雷达接收到的数据，并通过所述数据处理软件对数据进行分析处 理，以进行空间目标三维姿态的重构。本发明以个人计算机为例进行说明， 可以理解，任何可以实现本发明的电子设备，均属于本发明的保护范围。

所述雷达在接收回波信号之前会向空间目标发射信号，所述回波信号 即为所述空间目标基于发射的所述信号向雷达散射的信号。

参见图2是本发明实施例提供的一种雷达与空间目标之间的方位示意 图，参见图3是本发明实施例提供的一种空间目标几何形状示意图，本发 明对所述雷达的数量以及所述空间目标的几何形状不作具体限制，本发明 以三个雷达以及空间目标为三棱锥为例进行说明。所述空间目标中包括四 个稳定的散射中心，分别位于ABCD四个位置。

步骤2：根据所述回波信号，确定所述空间目标的所有散射中心对应的 一维距离序列。

可选的，所述步骤2包括：

步骤2-1：根据所述回波信号，确定雷达与空间目标对应的质心之间的 距离r₀。

本发明可以基于发射信号与接收到回波信号对应的时延，确定雷达与 空间目标对应的质心之间的距离。

步骤2-2：基于所述距离r₀，对所述回波信号进行解线性调频处理，以 对所述回波信号进行平动补偿，以得到一维距离像，表示为：

其中，

为快时间，t_k是慢时间，A_i为第i个散射中心的散射系数，T_p为 脉冲宽度，μ为调频率，τ₀＝2r₀/c，Δτ_i(t_k)＝τ_i(t_k)-τ₀，j是虚数单位。

步骤2-3：对所述信号的每个脉冲在快时间维进行迭代自适应处理，以 得到高分辨一维距离像，并沿着慢时间维将一维距离像拼接起来得到高分 辨一维距离像序列。

可选的，所述步骤2-3包括：

步骤2-31：选择每个脉冲对应的回波信号y_M×1作为观测信号。

步骤2-32：基于预设超分辨倍数K，从离散傅里叶基中构造观测矩阵 Φ_M×N，并基于预设迭代自适应次数，获取信号频谱估计值，表示为：

其中，φ_i为观测矩阵中的第i个基向量，y为观测信号。

步骤2-33：计算每一个超分辨点处的初始功率，表示为：

步骤2-34：基于所述初始功率，计算协方差矩阵，表示为：

S＝Φdiag(p)Φ^H。

步骤2-35：基于所述协方差矩阵，更新所述信号频谱估计值，表示为：

步骤2-36：判断当前迭代次数与所述预设迭代自适应次数的大小关系。

步骤2-37：当所述迭代次数小于所述预设迭代自适应次数时，重复执 行步骤2-34至步骤2-37；或者，当所述迭代次数等于所述预设迭代自适应 次数时，继续执行步骤2-4。

步骤2-4：提取所述高分辨一维距离像序列中峰值对应的距离维坐标， 以得到散射中心对应的一维距离序列。

步骤3：通过卡尔曼滤波对所述一维距离序列进行关联处理。

可选的，所述步骤3包括：

步骤3-1：基于上一时刻空间目标上散射中心在雷达视线上的投影距离， 获取当前时刻散射中心在雷达视线上的投影距离的预测值，表示为：

x_k|k-1＝A_k-1x_k-1|k-1，

其中，x_k|k-1是基于上一时刻散射中心在雷达视线上投影距离得到的当 前时刻的预测值，x_k-1|k-1是对上一时刻散射中心在雷达视线上的投影距离关 联的最优结果。

步骤3-2：基于上一时刻距离历程对应的方差，获取当前时刻的距离历 程预测值对应的方差，表示为：

其中，P_k|k-1是x_k|k-1对应的方差，P_k-1|k-1是x_k-1|k-1对应的方差，A_k-1为k时 刻的状态转移矩阵，为k时刻的状态噪声的方差。

步骤3-3：获取当前时刻的散射中心在雷达视线上投影距离的测量值；

步骤3-4：根据所述距离历程测量值与所述距离历程预测值，更新当前 时刻散射中心在雷达视线上投影距离的最优估计值，表示为：

x_k|k＝x_k|k-1+K_k(y_k-Hx_k|k-1)，

其中，K表示卡尔曼增益，

H_k为系统 观测矩阵。

步骤3-5：获取所述最优估计值对应的方差，表示为：

P_k|k＝(I-K_kH_k)P_k|k-1。

步骤3-6：重复步骤3-1至步骤3-5，直至获取到所有时刻散射中心对 应的距离历程。

步骤4：根据所述回波信号，估计所述散射中心对应的瞬时速度。

可选的，所述步骤4包括：

步骤4-1：建立空间目标回波信号对应的时变自回归模型，表示为：

其中，x(n)表示非平稳时间序列，a_k(n)(k＝1,2,3,...,p；n＝0,1,2,...,N-1) 表示时变自回归系数，p表示时变自回归模型的阶数，ε_n相互独立，且服从 N(0,σ²)。

所述时变自回归系数a_k(n)通过一组时变基函数的线性组合表示：

其中，q表示扩展基数，f_m(n)表示基函数， a_k,m(k＝1,2p,m3＝...,q为线性组合的系数矩阵，称为时变自回归模型 的时不变系数，该系数矩阵元素为常量。

所述x(n)的p阶时变自回归模型可表示为：

步骤4-2：对所述时变自回归模型中的线性时不变系数a_k,m进行求解。

可选的，所述步骤4-2包括：

步骤4-21：获取前后向的预测误差总和，表示为：

其中，前向预测误差和后向预测误差分别为：

步骤4-22：通过最小化ξ来求解a_k,m，令ξ关于

的导数为零，并对方程化简得：

其中，

步骤4-23：使用最小二乘法对化简后的所述方程进行求解，得时不变 系数为：

a＝(C^TC)^-1Cd，

其中，C是化简后的所述方程中等式左边第一个矩阵，d是等式右边的 矩阵。

步骤4-3：对功率谱密度函数P_xx(Ω,n)的极点进行求解，以得到瞬时速度 的估计值。

可选的，所述步骤4-3包括：

步骤4-31：根据时变自回归模型确定信号的功率谱密度函数，表示为：

步骤4-32：令z＝e^-jω，并使P_xx(Ω,n)的分母为0：

z^p+a₁(n)z^p-1+a₂z^p-2+...+a_p＝0。

步骤4-33：假设n时刻步骤4-32中公式的根为z_k(n)k＝1,2,3...,p，确 定瞬时频率分量，表示为：

f_k(n)＝angle(z_k(n))·Fs/2π，

其中，Fs为采样频率。

步骤4-34：根据瞬时频率分量和多普勒频率公式

确定瞬时速 度的估计值，表示为：

v_k(n)＝f_k(n)λ/2。

步骤5：基于卡尔曼滤波对所述散射中心对应的多普勒频率进行关联处 理。

可选的，所述步骤5包括：

步骤5-1：假设基于时变自回归模型得到了有n个频率成分信号的频率 估计结果，经过卡尔曼滤波得到系统状态后，定义X矩阵，X矩阵中的各元 素表示使用不同频率分量更新第i个频率分量后相邻两个时刻的状态差，X 矩阵表示为：

其中，

为用k时刻第i个频率分量的观测值

对第i个频率分量预 测状态

的更新结果。

步骤5-2：搜索矩阵X的最小值，并标记最小值的行和列为i和j，并将 第i个频率分量与下一时刻第j个频率分量的观测值进行匹配，以及将X 矩阵中的第i行和第j列置为无穷大。

步骤5-3：重复执行步骤5-2，直到矩阵X中的所有元素均为无穷大， 以得到频率分量的估计值。

步骤5-4：根据瞬时速度计算公式，计算得到多个散射中心对应的瞬时 速度序列。

步骤6：基于关联处理后的一维距离序列和瞬时速度，对散射中心进行 三维重构处理。

可选的，所述步骤6包括：

步骤6-1：根据空间目标上N个散射中心在雷达视线上的投影距离D_m， 获得第i个散射中心在连体坐标系中的初始坐标c_i，其中，设置收敛条件为：

k表示迭代次数，目标距离矩阵表示为：

其中，T表示总的观测时刻。

步骤6-2：根据雷达网中的M个雷达到第i个散射中心的距离 r_1i,…,r_mi,…,r_Mi和瞬时多普勒

计算目标上第i个散射中心在参 考坐标系中的坐标估计值

步骤6-3：基于所述坐标估计值，确定参考坐标系与连体坐标系对应的 旋转矩阵Q。

可选的，所述步骤6-3包括：

步骤6-31：确定在不同坐标系下的两组数据之间的旋转关系，表示为：

其中，加权矩阵W_i＝w_iI，SO(K)＝{Q∈R^K×K Q^TQ＝I det(Q)＝1}，K表示 坐标维度，对于三维模型K＝3，t表示连体坐标系到参考坐标系的平移向 量。

步骤6-32：对

和c_i进行加权平均得

并将

和c_i与加权平均值

分别作差，得到参考坐标系和连体坐标系中 的坐标估计误差

及

步骤6-33：对所述连体坐标系到参考坐标系的平移向量t求导，并令导 数为零，以确定t的解

步骤6-34：将所述所述连体坐标系到参考坐标系的平移向量的解和所 述坐标估计误差代入所述旋转关系对应的公式中，以得到代价函数，表示为：

其中，

最小化代价函数J等价于最大化

步骤6-35：令

的奇异值分解为UΣV^T，得到Q的最优解为：

Q＝Vdiag([1^T,det(VU^T)]^T)U^T，

其中，1为长度为2的单位向量，det(VU^T)保证了Q满足det(Q)＝1，将 Q的最优解代入

中得到平移向量t。此时，得到Q和t对应的估计 值，记为

步骤6-36：将第m个雷达与空间目标上第i个散射中心的距离r_mi和瞬时 速度

相乘，得到关于

及

的线性方程：

将Q和t对应的估计值

代入，得到

及

对应的初步估计值。

步骤6-37：设

的修正量为Q_δ，Δt，Δω，

将Q对应 的最终估计值表示为

其中，求解出的Q满足SO(K)条件，其他量 表示为：

步骤6-38：将Q_δ，Δt，Δω，

代入到

中，得到公式如下：

步骤6-39：对步骤6-38中的公式进行矢量化处理，以得到关于校正量 Q_δ，Δt，Δω，

对应的线性方程，并基于加权最小二乘法求解校正量β， Δt，Δω和

以及将所述校正量β，Δt，Δω和

与其对应的初步估计 量进行结合，以得到最终估计值。

步骤6-4：基于所述旋转矩阵Q，第m个雷达在参考坐标系中的坐标a_m， 以及空间目标质心在参考坐标系中的坐标p变换到连体坐标系，得到T个 时刻内雷达视线单位矢量在连体坐标系中的变化矩阵C_m：

C_m＝[l_m1,l_m2,…l_mT]，

其中，

为第m个雷达j时刻在连体坐标系中 的坐标，

为空间目标质心在连体坐标系中的坐标。

步骤6-5：对一维距离矩阵D_m以及雷达视线矩阵C_m之间的关系 D_m＝S_mC_m进行矢量化，并基于最小二乘法求解得到空间目标上散射中心在 连体坐标系中的三维坐标S_m。

步骤6-6：将M个雷达得到的S_m进行平均，以得到最终散射中心在连 体坐标系中的三维坐标矩阵

步骤6-7：重复步骤6-2至步骤6-6，直至

达到收敛条件，结束迭代 处理，并基于达到收敛条件的

和每个时刻的旋转矩阵Q_i，计算得到空间目 标在参考坐标系中每个时刻对应的N个散射中心的坐标。

综上，本发明利用迭代自适应IAA方法得到超分辨的一维距离序列， 既提高了距离分辨率，进而使得三维姿态重构的准确度提高，又避免了距离 分辨率受制于雷达带宽的问题。另外，本发明利用目标距离和瞬时速度信息 来对三维姿态进行重构，避免了仅利用一维距离序列基于SVD分解的两步 重构算法存在任意旋转矩阵的问题，提高了空间目标三维姿态重构的精确 度。

下面结合仿真实验对本发明的效果做进一步描述。

1.仿真实验条件。

本发明的仿真实验的硬件平台为：处理器为Inter(R)Core(TM)i5-4590 CPU主频为3.30GHZ，内存64GB。

本发明的仿真实验的软件平台为：CST STUDIO SUITE 2019和 MATLAB R2018b。

2.仿真内容与结果。

本发明的仿真实验采用本发明提出的算法与现有算法(基于SVD分解 的两步重构方法)对空间目标的三维姿态重构，得到两组三维坐标，并与目 标个散射中心真实坐标进行对比。

实验通过3部雷达对目标进行观测，其雷达视线在进动坐标系中的方 位角和俯仰角分别为[0°,20°]，[90°,20°]，[-45,15°]，参见图2和图3，目标 的运动方式为进动，目标自旋频率为0.5Hz，进动频率为1Hz，进动角为10°。 雷达参数设置：雷达发射脉冲重复频率为500Hz，观测时间为1s，扫频范 围为8～12GHz，极化方式为水平极化。

参见图4是本发明实施例提供的一种基于三部雷达得到的一维距离历 程以及关联后的一维距离历程示意图，图4(a)是第1部雷达的数据处理 得到目标的一维距离像及关联后的一维距离历程；图4(b)是第2部雷达 的数据处理得到目标的一维距离像及关联后的一维距离历程；图4(c)是 第3部雷达的数据处理得到目标的一维距离像及关联后的一维距离历程。

参见图5是本发明实施例提供的一种基于三部雷达得到的瞬时频率历 程以及关联后的瞬时频率历程示意图，其中，每种图形对应的曲线表示一个 散射中心的瞬时频率。图5(a)是第1部雷达数据处理得到的频率分量结 果及关联后的瞬时频率结果；图5(b)是第2部雷达数据处理得到的频率 分量结果及关联后的瞬时频率结果；图5(c)是第3部雷达数据处理得到 的频率分量结果及关联后的瞬时频率结果。

参见图6是本发明实施例提供的一种基于t＝0.08s时4个散射中心在参 考坐标系中的重构结果示意图，图6是t＝0.08s时目标散射中心在参考坐标 系中的重构结果，其中，重构方法1表示基于SVD分解的两步重构方法， 重构方法2表示本发明提出的算法，4个散射中心的三维重构结果如表1所 示。

参见图7是本发明实施例提供的一种在不同信噪比条件下A、B两个 散射中心的重构性能比较示意图，算法的性能用归一化均方误差来评估。在 每个信噪比下做100次试验，以散射中心A，B为例，得到散射中心A，B 重构结果随信噪比变化情况。

表1 t＝0.08s目标散射中心重构结果

	理论坐标/m	重构方法1/m	重构方法2/m
				散射点A	(0.0689,0.0368,0.4432)	(0.1204,0.1843,0.5844)	(0.0803,0.0506,0.4702)
散射点B	(-0.3797,0.3472,-0.2744)	(-0.5433,0.2108,-0.4624)	(-0.3592,0.3218,-0.2643)
				散射点C	(0.6099,0.5729,-0.4471)	(0.8041,0.4786,-0.6219)	(0.5912,0.5618,-0.4215)
散射点D	(0.0822,-0.6100,-0.2668)	(0.0677,-0.4953,-0.3411)	(0.0719,-0.5847,-0.2443)

3.仿真结果分析

从图4可以看出，不同方位的雷达观测到各散射中心的一维距离历程不 同，这是因为雷达视线不同，散射中心在雷达视线上的投影也不同；经过关 联后可以得到正确的一维距离历程。

从图5可以看出，不同方位的雷达观测到各散射中心的瞬时频率历程不 同，经过关联后可以得到正确的频率历程。

从图6和表1可以看出，相比于仅对一维距离序列基于SVD分解的两 步重构方法，本发明提出的算法三维重构结果与理论坐标更接近，即重构方 法的精度更高。

从图7可以看出，本发明所提出的算法及利用距离重构的误差都是随着 信噪比的增加而减小，其次利用距离进行重构的误差要比同时利用距离速 度进行重构的误差要大。

实施例二

本发明实施例还提供了一种电子设备，如图8所示，包括处理器801、 通信接口802、存储器803和通信总线804，其中，处理器801，通信接口 802，存储器803通过通信总线804完成相互间的通信，

存储器803，用于存放计算机程序；

处理器801，用于执行存储器803上所存放的程序时，实现如下步骤：

步骤1：获取所述雷达接收到的回波信号，所述回波信号由空间目标散 射至所述雷达；

步骤2：根据所述回波信号，确定所述空间目标的所有散射中心对应的 一维距离序列；

步骤3：通过卡尔曼滤波对所述一维距离序列进行关联处理；

步骤4：根据所述回波信号，估计所述散射中心对应的瞬时速度；

步骤5：基于卡尔曼滤波对所述散射中心对应的多普勒频率进行关联处 理；

步骤6：基于关联处理后的一维距离序列和瞬时速度，对散射中心进行 三维重构处理。

上述电子设备提到的通信总线可以是外设部件互连标准(Peripheral ComponentInterconnect，PCI)总线或扩展工业标准结构(Extended Industry StandardArchitecture，EISA)总线等。该通信总线可以分为地址总线、数据 总线、控制总线等。为便于表示，图中仅用一条粗线表示，但并不表示仅有 一根总线或一种类型的总线。

通信接口用于上述电子设备与其他设备之间的通信。

存储器可以包括随机存取存储器(Random Access Memory，RAM)，也 可以包括非易失性存储器(Non-Volatile Memory，NVM)，例如至少一个磁盘 存储器。可选的，存储器还可以是至少一个位于远离前述处理器的存储装置。

上述的处理器可以是通用处理器，包括中央处理器(Central Processing Unit，CPU)、网络处理器(Network Processor，NP)等；还可以是数字信号处 理器(DigitalSignal Processing，DSP)、专用集成电路(Application Specific Integrated Circuit，ASIC)、现场可编程门阵列(Field-Programmable Gate Array， FPGA)或者其他可编程逻辑器件、分立门或者晶体管逻辑器件、分立硬件组 件。

本发明实施例提供的方法可以应用于电子设备。具体的，该电子设备可 以为：台式计算机、便携式计算机、智能移动终端、服务器等。在此不作限 定，任何可以实现本发明的电子设备，均属于本发明的保护范围。

对于装置/电子设备/存储介质实施例而言，由于其基本相似于方法实施 例，所以描述的比较简单，相关之处参见方法实施例的部分说明即可。

需要说明的是，本发明实施例的装置、电子设备及存储介质分别是应用 上述一种基于距离和瞬时速度的空间目标三维姿态重构方法的装置、电子 设备及存储介质，则上述一种基于距离和瞬时速度的空间目标三维姿态重 构方法的所有实施例均适用于该装置、电子设备及存储介质，且均能达到相 同或相似的有益效果。

应用本发明实施例所提供的终端设备，可以展示专有名词和/或固定词 组供用户选择，进而减少用户输入时间，提高用户体验。

该终端设备以多种形式存在，包括但不限于：

(1)移动通信设备：这类设备的特点是具备移动通信功能，并且以提供 话音、数据通信为主要目标。这类终端包括：智能手机(例如iPhone)、多媒 体手机、功能性手机，以及低端手机等。

(2)超移动个人计算机设备：这类设备属于个人计算机的范畴，有计算 和处理功能，一般也具备移动上网特性。这类终端包括：PDA、MID和UMPC 设备等，例如iPad。

(3)便携式娱乐设备：这类设备可以显示和播放多媒体内容。该类设备 包括：音频、视频播放器(例如iPod)，掌上游戏机，电子书，以及智能玩具 和便携式车载导航设备。

(4)其他具有数据交互功能的电子装置。

在本发明的描述中，需要理解的是，术语“中心”、“纵向”、“横向”、 “长度”、“宽度”、“厚度”、“上”、“下”、“前”、“后”、“左”、 “右”、“竖直”、“水平”、“顶”、“底”、“内”、“外”、“顺时 针”、“逆时针”等指示的方位或位置关系为基于附图所示的方位或位置 关系，仅是为了便于描述本发明和简化描述，而不是指示或暗示所指的装置 或元件必须具有特定的方位、以特定的方位构造和操作，因此不能理解为对 本发明的限制。

在本说明书的描述中，参考术语“一个实施例”、“一些实施例”、 “示例”、“具体示例”、或“一些示例”等的描述意指结合该实施例或示 例描述的具体特征、结构、材料或者特点包含于本发明的至少一个实施例或 示例中。在本说明书中，对上述术语的示意性表述不必须针对的是相同的实 施例或示例。而且，描述的具体特征、结构、材料或者特点可以在任何的一 个或多个实施例或示例中以合适的方式结合。此外，本领域的技术人员可以 将本说明书中描述的不同实施例或示例进行接合和组合。

尽管在此结合各实施例对本申请进行了描述，然而，在实施所要求保护 的本申请过程中，本领域技术人员通过查看所述附图、公开内容、以及所附 权利要求书，可理解并实现所述公开实施例的其他变化。在权利要求中， “包括”(comprising)一词不排除其他组成部分或步骤，“一”或“一个” 不排除多个的情况。单个处理器或其他单元可以实现权利要求中列举的若 干项功能。相互不同的从属权利要求中记载了某些措施，但这并不表示这些措施不能组合起来产生良好的效果。

本领域技术人员应明白，本申请的实施例可提供为方法、装置(设备)、 或计算机程序产品。因此，本申请可采用完全硬件实施例、完全软件实施例、 或结合软件和硬件方面的实施例的形式，这里将它们都统称为“模块”或“系 统”。而且，本申请可采用在一个或多个其中包含有计算机可用程序代码的 计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等) 上实施的计算机程序产品的形式。计算机程序存储/分布在合适的介质中， 与其它硬件一起提供或作为硬件的一部分，也可以采用其他分布形式，如通 过Internet或其它有线或无线电信系统。

本申请是参照本申请实施例的方法、装置(设备)和计算机程序产品的流 程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方 框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的 结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理 机或其他可编程数据处理设备的处理器以产生一个机器，使得通过计算机 或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一 个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。

这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理 设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存 储器中的指令产生包括指令装置的制造品，该指令装置实现在流程图一个 流程或多个流程和/或方框图一个方框或多个方框中指定的功能。

这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上， 使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现 的处理，从而在计算机或其他可编程设备上执行的指令提供用于实现在流 程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的 步骤。

以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明， 不能认定本发明的具体实施只局限于这些说明。对于本发明所属技术领域 的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干简单 推演或替换，都应当视为属于本发明的保护范围。

Claims (10)

1.一种基于距离和瞬时速度的空间目标三维姿态重构方法，应用于电子设备，所述电子设备对应有多个雷达，其特征在于，所述方法包括：

步骤1：获取所述雷达接收到的回波信号，所述回波信号由空间目标散射至所述雷达；

步骤2：根据所述回波信号，确定所述空间目标的所有散射中心对应的一维距离序列；

步骤3：通过卡尔曼滤波对所述一维距离序列进行关联处理；

步骤4：根据所述回波信号，估计所述散射中心对应的瞬时速度；

步骤5：基于卡尔曼滤波对所述散射中心对应的多普勒频率进行关联处理；

步骤6：基于关联处理后的一维距离序列和瞬时速度，对散射中心进行三维重构处理。

2.根据权利要求1所述的方法，其特征在于，所述步骤2包括：

步骤2-1：根据所述回波信号，确定雷达与空间目标对应的质心之间的距离r₀；

步骤2-2：基于所述距离r₀，对所述回波信号进行解线性调频处理，以对所述回波信号进行平动补偿，以得到一维距离像，表示为：

其中，

为快时间，t_k是慢时间，A_i为第i个散射中心的散射系数，T_p为脉冲宽度，μ为调频率，τ₀＝2r₀/c，Δτ_i(t_k)＝τ_i(t_k)-τ₀，j是虚数单位；

步骤2-3：对所述信号的每个脉冲在快时间维进行迭代自适应处理，以得到高分辨一维距离像，并沿着慢时间维将一维距离像拼接起来得到高分辨一维距离像序列；

步骤2-4：提取所述高分辨一维距离像序列中峰值对应的距离维坐标，以得到散射中心对应的一维距离序列。

3.根据权利要求2所述的方法，其特征在于，所述步骤2-3包括：

步骤2-31：选择每个脉冲对应的回波信号y_M×1作为观测信号；

步骤2-32：基于预设超分辨倍数K，从离散傅里叶基中构造观测矩阵Φ_M×N，并基于预设迭代自适应次数，获取信号频谱估计值，表示为：

其中，φ_i为观测矩阵中的第i个基向量，y为观测信号；

步骤2-33：计算每一个超分辨点处的初始功率，表示为：

步骤2-34：基于所述初始功率，计算协方差矩阵，表示为：

S＝Φdiag(p)Φ^H；

步骤2-35：基于所述协方差矩阵，更新所述信号频谱估计值，表示为：

步骤2-36：判断当前迭代次数与所述预设迭代自适应次数的大小关系；

步骤2-37：当所述迭代次数小于所述预设迭代自适应次数时，重复执行步骤2-34至步骤2-37；或者，当所述迭代次数等于所述预设迭代自适应次数时，继续执行步骤2-4。

4.根据权利要求1所述的方法，其特征在于，所述步骤3包括：

步骤3-1：基于上一时刻空间目标上散射中心在雷达视线上的投影距离，获取当前时刻散射中心在雷达视线上的投影距离的预测值，表示为：

x_k|k-1＝A_k-1x_k-1|k-1，

其中，x_k|k-1是基于上一时刻散射中心在雷达视线上投影距离得到的当前时刻的预测值，x_k-1|k-1是对上一时刻散射中心在雷达视线上的投影距离关联的最优结果；

步骤3-2：基于上一时刻距离历程对应的方差，获取当前时刻的距离历程预测值对应的方差，表示为：

其中，P_k|k-1是x_k|k-1对应的方差，P_k-1|k-1是x_k-1|k-1对应的方差，A_k-1为k时刻的状态转移矩阵，为k时刻的状态噪声的方差；

步骤3-3：获取当前时刻的散射中心在雷达视线上投影距离的测量值；

步骤3-4：根据所述距离历程测量值与所述距离历程预测值，更新当前时刻散射中心在雷达视线上投影距离的最优估计值，表示为：

x_k|k＝x_k|k-1+K_k(y_k-Hx_k|k-1)，

其中，K表示卡尔曼增益，

H_k为系统观测矩阵；

步骤3-5：获取所述最优估计值对应的方差，表示为：

P_k|k＝(I-K_kH_k)P_k|k-1；

步骤3-6：重复步骤3-1至步骤3-5，直至获取到所有时刻散射中心对应的距离历程。

5.根据权利要求1所述的方法，其特征在于，所述步骤4包括：

步骤4-1：建立空间目标回波信号对应的时变自回归模型，表示为：

其中，x(n)表示非平稳时间序列，a_k(n)(k＝1,2,3,...,p；n＝0,1,2,...,N-1)表示时变自回归系数，p表示时变自回归模型的阶数，ε_n相互独立，且服从N(0,σ²)；

所述时变自回归系数a_k(n)通过一组时变基函数的线性组合表示：

其中，q表示扩展基数，f_m(n)表示基函数，a_k，m(k＝1，2，3...，p；m＝0，1，2...q)为线性组合的系数矩阵，称为时变自回归模型的时不变系数，该系数矩阵元素为常量；

所述x(n)的p阶时变自回归模型可表示为：

步骤4-2：对所述时变自回归模型中的线性时不变系数a_k,m进行求解；

步骤4-3：对功率谱密度函数P_xx(Ω,n)的极点进行求解，以得到瞬时速度的估计值。

6.根据权利要求5所述的方法，其特征在于，所述步骤4-2包括：

步骤4-21：获取前后向的预测误差总和，表示为：

其中，前向预测误差和后向预测误差分别为：

步骤4-22：通过最小化ξ来求解a_k,m，令ξ关于

的导数为零，并对方程化简得：

其中，

步骤4-23：使用最小二乘法对化简后的所述方程进行求解，得时不变系数为：

a＝(C^TC)^-1Cd，

其中，C是化简后的所述方程中等式左边第一个矩阵，d是等式右边的矩阵。

7.根据权利要求5所述的方法，其特征在于，所述步骤4-3包括：

步骤4-31：根据时变自回归模型确定信号的功率谱密度函数，表示为：

步骤4-32：令z＝e^-jω，并使P_xx(Ω,n)的分母为0：

z^p+a₁(n)z^p-1+a₂z^p-2+...+a_p＝0；

步骤4-33：假设n时刻步骤4-32中公式的根为z_k(n)k＝1,2,3...,p，确定瞬时频率分量，表示为：

f_k(n)＝angle(z_k(n))·Fs/2π，

其中，Fs为采样频率；

步骤4-34：根据瞬时频率分量和多普勒频率公式

确定瞬时速度的估计值，表示为：

v_k(n)＝f_k(n)λ/2。

8.根据权利要求1所述的方法，其特征在于，所述步骤5包括：

步骤5-1：假设基于时变自回归模型得到了有n个频率成分信号的频率估计结果，经过卡尔曼滤波得到系统状态后，定义X矩阵，X矩阵中的各元素表示使用不同频率分量更新第i个频率分量后相邻两个时刻的状态差，X矩阵表示为：

其中，

为用k时刻第i个频率分量的观测值

对第i个频率分量预测状态

的更新结果；

步骤5-2：搜索矩阵X的最小值，并标记最小值的行和列为i和j，并将第i个频率分量与下一时刻第j个频率分量的观测值进行匹配，以及将X矩阵中的第i行和第j列置为无穷大；

步骤5-3：重复执行步骤5-2，直到矩阵X中的所有元素均为无穷大，以得到频率分量的估计值；

步骤5-4：根据瞬时速度计算公式，计算得到多个散射中心对应的瞬时速度序列。

9.根据权利要求1所述的方法，其特征在于，所述步骤6包括：

步骤6-1：根据空间目标上N个散射中心在雷达视线上的投影距离D_m，获得第i个散射中心在连体坐标系中的初始坐标c_i，其中，设置收敛条件为：

k表示迭代次数，目标距离矩阵表示为：

其中，T表示总的观测时刻；

步骤6-2：根据雷达网中的M个雷达到第i个散射中心的距离r_1i,…,r_mi,…,r_Mi和瞬时多普勒

计算目标上第i个散射中心在参考坐标系中的坐标估计值

步骤6-3：基于所述坐标估计值，确定参考坐标系与连体坐标系对应的旋转矩阵Q；

步骤6-4：基于所述旋转矩阵Q，第m个雷达在参考坐标系中的坐标a_m，以及空间目标质心在参考坐标系中的坐标p变换到连体坐标系，得到T个时刻内雷达视线单位矢量在连体坐标系中的变化矩阵C_m：

C_m＝[l_m1,l_m2,…l_mT]，

其中，

为第m个雷达j时刻在连体坐标系中的坐标，

为空间目标质心在连体坐标系中的坐标；

步骤6-5：对一维距离矩阵D_m以及雷达视线矩阵C_m之间的关系D_m＝S_mC_m进行矢量化，并基于最小二乘法求解得到空间目标上散射中心在连体坐标系中的三维坐标S_m；

步骤6-6：将M个雷达得到的S_m进行平均，以得到最终散射中心在连体坐标系中的三维坐标矩阵

步骤6-7：重复步骤6-2至步骤6-6，直至

达到收敛条件，结束迭代处理，并基于达到收敛条件的

和每个时刻的旋转矩阵Q_i，计算得到空间目标在参考坐标系中每个时刻对应的N个散射中心的坐标。

10.根据权利要求9所述的方法，其特征在于，所述步骤6-3包括：

步骤6-31：确定在不同坐标系下的两组数据之间的旋转关系，表示为：

其中，加权矩阵W_i＝w_iI，SO(K)＝{Q∈R^K×K Q^TQ＝I det(Q)＝1}，K表示坐标维度，对于三维模型K＝3，t表示连体坐标系到参考坐标系的平移向量；

步骤6-32：对

和c_i进行加权平均得

并将

和c_i与加权平均值

分别作差，得到参考坐标系和连体坐标系中的坐标估计误差

及

步骤6-33：对所述连体坐标系到参考坐标系的平移向量t求导，并令导数为零，以确定t的解

步骤6-34：将所述所述连体坐标系到参考坐标系的平移向量的解和所述坐标估计误差代入所述旋转关系对应的公式中，以得到代价函数，表示为：

其中，

最小化代价函数J等价于最大化

步骤6-35：令

的奇异值分解为UΣV^T，得到Q的最优解为：

Q＝Vdiag([1^T,det(VU^T)]^T)U^T，

其中，1为长度为2的单位向量，det(VU^T)保证了Q满足det(Q)＝1，将Q的最优解代入

中得到平移向量t。此时，得到Q和t对应的估计值，记为

步骤6-36：将第m个雷达与空间目标上第i个散射中心的距离r_mi和瞬时速度

相乘，得到关于

及

的线性方程：

并将Q和t对应的估计值

代入，得到

及

对应的初步估计值；

步骤6-37：设

的修正量为Q_δ，Δt，Δω，

将Q对应的最终估计值表示为

其中，求解出的Q满足SO(K)条件，其他量表示为：

步骤6-38：将Q_δ，Δt，Δω，

代入到

中，得到公式如下：

步骤6-39：对步骤6-38中的公式进行矢量化处理，以得到关于校正量，Q_δ，Δt，Δω，

对应的线性方程，并基于加权最小二乘法求解校正量β，Δt，Δω和

以及将所述校正量β，Δt，Δω和

与其对应的初步估计量进行结合，以得到最终估计值。

Images