CN104833970B - 一种适用于存在遮挡效应的旋转目标稀疏多普勒成像方法 - Google Patents

Abstract

本发明公开了一种适用于存在遮挡效应的旋转目标稀疏多普勒成像方法，包括如下步骤：首先，建立存在遮挡效应的旋转目标的时变雷达回波模型；对存在遮挡效应的旋转目标的时变雷达回波模型进行稀疏信号表示；然后，根据时变雷达回波模型的稀疏信号表示形式，将对存在遮挡效应的旋转目标的成像过程转化为2维Fused Lasso模型的稀疏优化过程，并利用凸优化工具包CVX求解2维Fused Lasso模型的稀疏优化过程，得到稀疏后向散射强度矩阵B的估计矩阵；最后，取估计矩阵中各行的最大值构成等效后向散射强度向量，即实现成像。

Description

一种适用于存在遮挡效应的旋转目标稀疏多普勒成像方法

技术领域

本发明属于雷达成像技术领域，涉及一种适用于存在遮挡效应的旋转目标稀疏多普勒成像方法，可用于对存在遮挡效应的旋转目标进行成像。

背景技术

多普勒成像方法能够提供目标在距离维和横向距离维的精确信息，是实现对3维旋转目标进行2维成像的一种常用方法。针对不同的旋转目标设计有效的多普勒成像方法在实际应用中具有重要意义。例如，对太空碎片进行多普勒成像是空间碎片监测系统中的关键技术之一。目前，随着人类空间活动的日益频繁，空间碎片的数量也在与日俱增，已成为了外太空中极具危害的一种污染源。它们不仅拥有庞大的数量，还具有很高的运动速度，通常在近地轨道上能够达到10km/s，从而极易与在轨飞行器发生碰撞，导致大灾难的发生。有效的太空碎片多普勒成像在防撞指示、碰撞程度预测以及碎片的回收清理等工作中均具有重要的指导意义。

在多普勒成像过程中，目标的雷达回波在整个观测时间段内都被看作位于单个距离单元内，并且被一系列的雷达发射脉冲序列采样。根据散射中心模型理论，多普勒成像的关键在于估计目标上主要散射中心的位置和后向散射强度，并且估计精度需要尽可能地高。仅用多普勒信息的传统单快拍成像方法，其假设在整个成像过程中观测目标上所有散射中心始终对于雷达是可见的。但在实际中，当目标的自旋轴与雷达视线不一致时，会导致目标上的部分散射中心在部分观测间隔内对于雷达是不可见的，即遮挡效应，此时的遮挡效应相当于对目标在方位向进行了非均匀采样，这个因素在已有成像算法中没有被涉及到。

一般而言，雷达成像的过程可看作是由测量场景的电磁散射信号来反推场景信息的过程。根据雷达成像的数学模型，可以将雷达成像的过程看成是一个逆问题。由于遮挡效应的存在，测量数据的信息量通常是有限的，无法满足人们对场景信息详细程度的要求。所以，上述逆问题是病态的，无法通过经典的最小二乘估计方法得以很好地解决。传统的多普勒成像方法通过对最小二乘估计中不可逆或不稳定的部分作近似处理来实现成像，但代价是具有一定的主瓣宽度和旁瓣效应，即造成场景图像中的细节信息有所损失。

发明内容

本发明的目的在于克服上述已有技术的不足，提出一种适用于存在遮挡效应的旋转目标稀疏多普勒成像方法，可实现对存在遮挡效应的旋转目标进行成像，并且能够自适应地估计遮挡区域。

本发明的实现思路是：将对存在遮挡效应的旋转目标进行多普勒成像的过程看作是求解一个逆问题的过程，利用目标上的散射中心是稀疏分布的(观测场景具有一定的稀疏性)以及目标上各个散射中心的后向散射现象随时间缓慢变化(慢变)这两个假设，在多普勒成像过程中增加相应的稀疏约束，将病态的逆问题转化成良态的2维Fused Lasso模型来求解，从而实现成像。

为达到上述技术目的，本发明采用以下技术方案予以实现。

一种适用于存在遮挡效应的旋转目标稀疏多普勒成像方法，其特征在于，包括以下步骤：

步骤1，建立存在遮挡效应的旋转目标的时变雷达回波模型；

步骤2，对存在遮挡效应的旋转目标的时变雷达回波模型进行稀疏信号表示；

步骤3，根据时变雷达回波模型的稀疏信号表示形式，将对存在遮挡效应的旋转目标的成像过程转化为2维Fused Lasso模型的稀疏优化过程；

步骤4，利用凸优化工具包CVX来求解2维Fused Lasso模型的稀疏优化过程，得到稀疏后向散射强度矩阵B的估计矩阵

步骤5，取估计矩阵中各行的最大值构成等效后向散射强度向量，即实现成像。

本发明的有益效果是：将对存在遮挡效应的旋转目标进行多普勒成像的过程看作是求解逆问题的过程，利用目标上的散射中心是稀疏分布的以及目标上各个散射中心的后向散射现象具有慢变特性这两个假设，在多普勒成像过程中引入稀疏约束，将病态的逆问题转化成良态的2维Fused Lasso模型来求解，从而实现对存在遮挡效应的旋转目标的多普勒成像，同时自适应地估计目标的遮挡区域。

附图说明

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

图1是本发明的流程图。

图2是用本发明对不同遮挡比例下的旋转目标的成像结果图，其中：

图2a是本发明对遮挡比例为30％的旋转目标的成像结果图；

图2b是本发明对遮挡比例为60％的旋转目标的成像结果图；

图2c是本发明对遮挡比例为90％的旋转目标的成像结果图；

图2d是本发明对遮挡比例为95％的旋转目标的成像结果图。

图3是本发明对不同信噪比下遮挡比例为50％的旋转目标的成像结果图，其中：

图3a是本发明对20dB信噪比下遮挡比例为50％的旋转目标的成像结果图；

图3b是本发明对15dB信噪比下遮挡比例为50％的旋转目标的成像结果图；

图3c是本发明对10dB信噪比下遮挡比例为50％的旋转目标的成像结果图。

图4是电磁仿真的目标模型示意图。

图5是电磁仿真目标的一维高分辨距离像随时间的变化图。

图6是本发明对电磁仿真数据的成像结果图。

具体实施方式

参照图1，本发明的适用于存在遮挡效应的旋转目标稀疏多普勒成像方法，包括以下步骤：

步骤1，建立存在遮挡效应的旋转目标的时变雷达回波模型。

假设旋转目标由K个散射中心构成，当宽带雷达带宽的距离分辨率大于目标尺寸时，旋转目标上各个散射中心的雷达子回波都被看作是来自于同一个距离单元。未考虑遮挡效应时，旋转目标的各次雷达回波均为目标上K个散射中心的雷达子回波的线性叠加。因此，未考虑遮挡效应的旋转目标的第n次雷达回波s(n)可表示为：

其中，A_k是第k个散射中心的后向散射系数，λ是雷达发射波波长，ω为旋转目标的旋转角速度，F_r是雷达脉冲重复频率，r_k和θ_k分别为第k个散射中心的半径和起始角度，N为雷达回波次数。

考虑到噪声的影响，根据式(1)可以将未考虑遮挡效应的旋转目标的N次雷达回波表示为

s＝Φα+w (2)

其中，s＝[s(0),s(1),...,s(N-1)]^T表示未考虑遮挡效应的旋转目标的雷达回波向量，上标T为转置，Φ为N×K维矩阵，该矩阵中第n+1行第k列元素为(Φ)_(n+1)k＝exp[-j4πr_ksin(θ_k+n(ω/F_r))/λ]，α为K维后向散射强度向量，α＝[A₁,A₂,...,A_K]^T，w为复高斯噪声。

对于未考虑遮挡效应的情况，各个时刻的雷达回波均为目标上所有散射中心雷达子回波的线性叠加。实际上，当目标作旋转运动时，目标体上部分区域被遮挡而且被遮挡的区域是在不停地变换着的。具体来说，在不同的观测时刻，目标上的散射中心轮流进入遮挡区。进入遮挡区的散射中心对目标的雷达回波是没有贡献的，即这些散射中心的雷达子回波为零。因此，每个时刻的雷达回波是动态变化的，是未被遮挡的散射中心的雷达子回波的线性叠加。对应地，式(2)中的后向散射强度向量α随时间是动态变化的，各个时刻的后向散射强度向量α中仅有部分元素为非零值。此时，遮挡效应相当于对目标在方位向进行了非均匀采样。根据存在遮挡效应时后向散射强度向量α的动态变化特性及式(2)，可得到存在遮挡效应的旋转目标在第t时刻的雷达回波s_o(t)为：

s_o(t)＝φ_tα_t+w_t,t＝0,...,N-1 (3)

其中，φ_t为矩阵Φ的第t+1行，α_t为第t时刻的后向散射强度向量，w_t为第t时刻的复高斯噪声。

令s_o＝[s_o(0),s_o(1),…,s_o(N-1)]^T，s_o表示存在遮挡效应的旋转目标的雷达回波向量，后向散射强度矩阵A由各个时刻的后向散射强度向量构成，A＝[α₀,α₁,…,α_N-1]，其维度为K×N。由式(3)，可以得到存在遮挡效应的旋转目标的时变雷达回波模型为：

s_o＝Φ⊙A+w (4)

Φ为N×K维矩阵，该矩阵中第n+1行第k列元素为(Φ)_(n+1)k＝exp[-j4πr_ksin(θ_k+n(ω/F_r))/λ]，“⊙”定义为向量-向量乘积，即Φ⊙A＝[φ₀α₀,φ₁α₁,…,φ_N-1α_N-1]^T，w为复高斯噪声。

步骤2，对存在遮挡效应的旋转目标的时变雷达回波模型进行稀疏信号表示。

其中，稀疏信号表示所需的空间字典通过离散网格点模型来构建，稀疏性通过假设存在遮挡效应的旋转目标上的散射中心是稀疏分布的来引入。

假设旋转目标上的散射中心分布在观测区域C内，其半径的观测窗为[r_L,r_H]，初始角度的观测窗为[0,2π)。利用离散网格点模型将观测区域C划分成M个空间分割区，每个空间分割区中至多存在一个散射中心。

观测区域C记为：

其中，为第m个空间分割区的半径，i₁＝1,…,N_r，N_r为等间隔划分半径观测窗的数目，Δr为半径划分间隔；为第m个空间分割区的初始角度，i₂＝1,…,N_θ，N_θ为等间隔划分初始角度观测窗的数目，Δθ为初始角度划分间隔；m为第m个空间分割区的标号，m＝(i₁-1)×N_θ+i₂；M＝N_r×N_θ是空间分割区的总个数。

首先，利用来自可能存在散射中心的M个空间分割区的雷达子回波向量来构造空间字典Ψ。其中，第m个可能存在的散射中心(位于第m个空间分割区)的N维雷达子回波向量可构成空间字典Ψ的第m个原子s_m为

然后，将空间字典Ψ的M个原子依次排列，空间字典Ψ为

Ψ＝[s₀,…,s_m,…,s_M]_N×M (7)

假设N次雷达回波来自M个空间分割区的K个散射中心(K＜＜M)，即目标上的散射中心是稀疏分布的，可根据式(4)，对存在遮挡效应的旋转目标的时变雷达回波模型进行稀疏信号表示，具有如下形式：

s_o＝Ψ⊙B+w (8)

其中，稀疏后向散射强度矩阵B由各时刻的稀疏后向散射强度向量β_t,t＝0,…,N-1构成，B＝[β₀,β₁,...,β_N-1]，维度为N×M；稀疏后向散射强度向量β_t为M维列向量，稀疏度为K(即存在K个非零元素)，稀疏后向散射强度向量β_t中的非零元素表示第t时刻在观测区域C相应的空间分割区存在散射中心，零元素则表示第t时刻在观测区域C相应的空间分割区不存在散射中心；将空间字典Ψ表示为行向量构成形式，即其中为Ψ的第t+1行行向量，则w为复高斯噪声。

步骤3，根据时变雷达回波模型的稀疏信号表示形式，将对存在遮挡效应的旋转目标的成像过程转化为2维Fused Lasso模型的稀疏优化过程。

根据时变雷达回波模型的稀疏信号表示形式，可将对存在遮挡效应的旋转目标的成像过程看作是一个逆问题，在最小二乘估计的基础上引入两个稀疏约束，第一个稀疏约束表示目标上的散射中心是稀疏分布的，即任一观测时刻的稀疏后向散射强度向量β在空间维上具有稀疏性；第二个稀疏约束表示目标上各散射中心的稀疏后向散射强度向量β具有慢变特性，即任一观测时刻的稀疏后向散射强度向量β在时间维上具有分段恒定特性。整个过程可通过2维Fused Lasso模型来描述：

其中，稀疏后向散射强度矩阵B由各时刻的稀疏后向散射强度向量β_t,t＝0,…,N-1构成，B_m,l表示第(l-1)时刻第m个可能散射中心(位于半径为初始角度为的空间分割区)的稀疏后向散射强度，τ₁和τ₂分别为任一时刻的稀疏后向散射强度向量β的稀疏度控制参数及光滑性控制参数。

步骤4，利用凸优化工具包CVX来求解2维Fused Lasso模型的稀疏优化过程，得到稀疏后向散射强度矩阵B的估计矩阵

具体地，利用凸优化工具包CVX来求解式(9)所表示的2维Fused Lasso模型的稀疏优化问题，得到稀疏后向散射强度矩阵B的估计矩阵

步骤5，取估计矩阵中各行的最大值构成等效后向散射强度向量，即实现成像。其中，等效后向散射强度向量中非零元素的位置表征了目标上各散射中心的在观测区域C的空间位置，等效后向散射强度向量中非零元素的值则表征了目标上各散射中心的后向散射强度。

另外，通过对估计矩阵中各非零行的所有相邻两个跳跃点间(由照射区进入遮挡区，再由遮挡区进入照射区)的零值或近零值进行计数，然后关于估计矩阵的所有非零行的该计数值做平均(记为)，可估计得到遮挡区域的范围为

其中，ω为目标的自旋角速度，F_r为雷达脉冲重复频率，即采样频率。

本发明的效果通过以下数值仿真实验和电磁仿真实验作进一步说明：

(1)数值仿真实验参数设置

雷达参数设置如下：脉冲重复频率为300Hz,波长为0.03m。观测目标包含四个各向同性的独立散射中心，其空间位置分别为(0.03m,0.25π),(0.05m,0.75π),(0.015m,1.2π)，(0.08m,1.75π)，坐标为极坐标表达方式，其中，空间位置的第一个元素表示散射中心所处位置的半径值，空间位置的第二个元素表示散射中心所处位置的初始角度；各独立散射中心对应的后向散射强度分别为0.3，0.5，0.8，1.0；目标的自旋角频率为10πrad/s，信噪比定义为单次回波的平均能量与噪声功率之比。

(2)数值仿真实验内容

1)检验本发明在不同遮挡区域尺寸下的成像性能。如图2所示，其中，图2a是本发明对遮挡比例为30％的旋转目标的成像结果图，图2b是本发明对遮挡比例为60％的旋转目标的成像结果图，图2c是本发明对遮挡比例为90％的旋转目标的成像结果图，图2d是本发明对遮挡比例为95％的旋转目标的成像结果图。上述图中各采用极坐标表示形式，方位角代表初始角度，逆时针方向为正方向，同心圆代表半径。

2)检验本发明在不同信噪比下的成像性能，如图3所示，其中，图3a,图3b,图3c分别是本发明对20dB，15dB,及10dB信噪比下遮挡比例为50％的旋转目标的成像结果图。

(3)数值仿真实验结果分析

从图2可看出，本发明针对不同遮挡比例的旋转目标均能较好地成像。当遮挡比例达到90％，本发明仍能准确地估计目标上四个散射中心的位置，但是稀疏后向散射强度的估计精度随遮挡范围增大而有所降低。当遮挡区域达到95％时，成像平面出现了一些虚假散射中心，这是因为此时的雷达回波数据缺失严重，无法分辨真实的散射中心。从稀疏重构的角度来看，当观测样本十分有限时，精确重构将无法得到保证。

从图3可看出，随着信噪比的提高，本发明的成像性能逐渐改善。当信噪比为20dB时，四个散射中心能够完全分辨，且稀疏后向散射强度估计值基本接近真值；当信噪比降为15dB时，仍能分辨四个散射中心，但是稀疏后向散射强度估计值的精度有所下降；当信噪比为10dB时，出现了虚假散射中心，其强度与最弱的散射中心基本持平。这是因为当噪声足够强时，噪声将淹没弱散射中心，从而导致本发明无法将弱散射中心从噪声中分离出来。这与稀疏理论中含噪观测会恶化重构性能的事实相符合。

(4)电磁仿真实验参数设置

本实验采用专业3维电磁仿真软件CST Studio Suite 2011产生雷达观测样本。目标采用简单的锥台模型，如图4所示，材质为理想良导体。其中，锥台的上底面半径为0.1m，下底面半径为0.3m，高度为0.3m,表面刻有5个半球形凹槽，凹槽的空间坐标分别为(0.14m,0),(0.18m,0.5π),(0.20m,π),(0.12m,1.25π)及(0.24m,1.5π)，坐标为极坐标表达方式，空间坐标的第一个元素表示散射中心所处位置的半径值，空间坐标的第二个元素表示散射中心所处位置的初始方位角，凹槽深度为0.01m。仿真的雷达参数设置如下：载频为10.01GHz，脉冲重复频率为360Hz，驻留时间为0.5s。目标的自旋频率为2r/s。

(5)电磁仿真实验内容

参照图5，为了说明自旋目标的遮挡效应是真实存在的，通过产生一组宽带数据来获得目标的一维高分辨距离像。图5中，横坐标为方向角，单位为(度)，纵坐标为径向距离，坐标为(m)。图5展示了根据电磁仿真数据获得的一维距离像随时间的变化情况，从图5中可清晰地看到5个不完整的正弦曲线，表明遮挡效应确实是显著存在的。实验中信噪比设定为15dB。

(6)电磁仿真实验结果分析：

参照图6，是本发明对电磁仿真数据的成像结果图。从图6可看出，本发明能够从含噪观测样本中准确分辨出5个凹槽。另外，根据式(10)可估计得到照射区比例大概在26.4％(通过计算得到)，估计值较真实值要稍偏低些，主要原因在于来自照射区与遮挡区的边界地带的雷达回波相对来说比较弱，从而在边界地带很难将五个散射中心分辨开来。

显然，本领域的技术人员可以对本发明进行各种改动和变型而不脱离本发明的精神和范围。这样，倘若本发明的这些修改和变型属于本发明权利要求及其等同技术的范围之内，则本发明也意图包含这些改动和变型在内。

Claims (3)

1.一种适用于存在遮挡效应的旋转目标稀疏多普勒成像方法，其特征在于，包括如下步骤：

步骤1，建立存在遮挡效应的旋转目标的时变雷达回波模型；

步骤2，对存在遮挡效应的旋转目标的时变雷达回波模型进行稀疏信号表示；

步骤3，根据时变雷达回波模型的稀疏信号表示形式，将对存在遮挡效应的旋转目标的成像过程转化为2维Fused Lasso模型的稀疏优化过程；

其中，所述2维Fused Lasso模型描述为：

$\begin{array}{cc}s.t.& \underset{m=1}{\overset{M}{\Σ}}\underset{l=1}{\overset{N}{\Σ}}\left|B_{m,l}\right|\≤{\mathrm{\τ}}_1\end{array}$

$\underset{m=1}{\overset{M}{\Σ}}\underset{l=1}{\overset{N-1}{\Σ}}|B_{m,l+1}-B_{m,l}|\≤{\mathrm{\τ}}_2$

其中，稀疏后向散射强度矩阵B由各时刻的稀疏后向散射强度向量β_t，t＝0，…，N-1构成；B_m，l表示第(l-1)时刻第m个可能散射中心的稀疏后向散射强度；τ₁和τ₂分别为任一时刻的稀疏后向散射强度向量β的稀疏度控制参数及光滑性控制参数；

步骤4，利用凸优化工具包CVX来求解2维Fused Lasso模型的稀疏优化过程，得到稀疏后向散射强度矩阵B的估计矩阵

步骤5，取估计矩阵中各行的最大值构成等效后向散射强度向量，即实现成像。

2.如权利要求1所述的适用于存在遮挡效应的旋转目标稀疏多普勒成像方法，其特征在于，所述步骤1的具体子步骤为：

假设旋转目标由K个散射中心构成，当宽带雷达带宽的距离分辨率大于目标尺寸时，旋转目标上各个散射中心的雷达子回波都被看作是来自于同一个距离单元；因此，未考虑遮挡效应的旋转目标的第n次雷达回波s(n)可表示为

$\begin{array}{c}s\left(n\right)=\underset{k=1}{\overset{K}{\Σ}}{A}_k\exp \[-j4{\mathrm{\πr}}_{k}sin({\mathrm{\θ}}_k+n(\mathrm{\ω}/F_r\left)\right)/\mathrm{\λ}\]\\ n=0,1,\mathrm{...},N-1\end{array}---\left(1\right)$

其中，A_k是第k个散射中心的后向散射系数，λ是雷达发射波波长，ω为旋转目标的旋转角速度，F_r是雷达脉冲重复频率，r_k和θ_k分别为第k个散射中心的半径和起始角度，N为雷达回波次数；

考虑到噪声的影响，根据式(1)可以将未考虑遮挡效应的旋转目标的N次雷达回波表示为

s＝Φα+w (2)

其中，s＝[s(0)，s(1)，...，s(N-1)]^T表示未考虑遮挡效应的旋转目标的雷达回波向量，上标T为转置，Φ为N×K维矩阵，该矩阵中第n+1行第k列元素为(Φ)_(n+1)k＝exp[-j4πr_ksin(θ_k+n(ω/F_r))/λ]，α为K维后向散射强度向量，α＝[A₁，A₂，...，A_K]^T，w为复高斯噪声；

存在遮挡效应的旋转目标在第t时刻的雷达回波s_o(t)为：

s_o(t)＝φ_tα_t+w_t，t＝0，...，N-1 (3)

其中，φ_t为矩阵Φ的第t+1行，α_t为第t时刻的后向散射强度向量，w_t为第t时刻的复高斯噪声；

令s_o＝[s_o(0)，s_o(1)，…，s_o(N-1)]^T，s_o表示存在遮挡效应的旋转目标的雷达回波向量，后向散射强度矩阵A由各个时刻的后向散射强度向量构成，A＝[α₀，α₁，…，α_N-1]，其维度为K×N；由式(3)，可以得到存在遮挡效应的旋转目标的时变雷达回波模型为：

s_o＝Φ⊙A+w (4)

Φ为N×K维矩阵，该矩阵中第n+1行第k列元素为(Φ)_(n+1)k＝exp[-j4πr_ksin(θ_k+n(ω/F_r))/λ]，“⊙”定义为向量-向量乘积，即Φ⊙A＝[φ₀α₀，φ₁α₁，…，φ_N-1α_N-1]^T，w为复高斯噪声。

3.如权利要求1所述的适用于存在遮挡效应的旋转目标稀疏多普勒成像方法，其特征在于，所述步骤2的具体子步骤为：

假设旋转目标上的散射中心分布在观测区域C内，其半径的观测窗为[r_L，r_H]，初始角度的观测窗为[0，2π)；利用离散网格点模型将观测区域C划分成M个空间分割区，每个空间分割区中至多存在一个散射中心；

再利用来自可能存在散射中心的M个空间分割区的雷达子回波向量来构造空间字典Ψ；其中，第m个可能存在的散射中心的N维雷达回波向量可构成空间字典Ψ的第m个原子s_m；然后，将空间字典Ψ的M个原子依次排列，空间字典Ψ为

Ψ＝[s₀，…，s_m，…，s_M]_N×M (7)

最后，假设N次雷达回波来自M个空间分割区的K个散射中心，即目标上的散射中心是稀疏分布的，对存在遮挡效应的旋转目标的时变雷达回波模型进行稀疏信号表示，具有如下形式：

s_o＝Ψ⊙B+w (8)

其中，稀疏后向散射强度矩阵B由各时刻的稀疏后向散射强度向量β_t，t＝0，…，N-1构成，B＝[β₀，β₁，...，β_N-1]，维度为N×M；稀疏后向散射强度向量β_t为M维列向量，稀疏度为K，稀疏后向散射强度向量β_t中的非零元素表示第t时刻在观测区域C相应的空间分割区存在散射中心，零元素则表示第t时刻在观测区域C相应的空间分割区不存在散射中心；将空间字典Ψ表示为行向量构成形式，即其中为Ψ的第t+1行行向量，则w为复高斯噪声。