CN103064071B - 基于稀疏分解的雷达目标属性散射中心特征提取方法 - Google Patents

Abstract

本发明公开一种基于稀疏分解的雷达目标属性散射中心特征提取方法，主要解决现有基于对雷达图像进行图像分割提取属性散射中心方法的模型失配、特征易丢失以及参数估计精度较低的问题。其实现过程是：首先利用噪声样本建立散射中心强度门限，在雷达图像中进行强散射中心检测，并确定散射中心参数的取值集合；然后根据属性散射中心模型，利用坐标轮回下降技术构建超分辨字典，通过求解0范数优化问题得到目标属性散射中心参数超分辨估计集合，根据散射中心参数集合提取目标及其重要部件的几何尺寸特征。本发明能有效提取目标属性散射中心、超分辨估计散射中心参数、精确估计目标及其重要部件几何尺寸，可用于雷达目标分类识别。

Description

基于稀疏分解的雷达目标属性散射中心特征提取方法

技术领域

本发明属于雷达技术领域，涉及一种雷达目标属性散射中心特征提取方法，可用于估计目标及其重要部件的几何尺寸，为目标分类识别提供重要的特征信息。

背景技术

雷达成像技术是20世纪50年代发展起来的，雷达图像是目标的二维散射图。传统雷达成像是以点散射模型为基础的，该模型只包含目标散射点位置信息，但仅利用目标散射点的位置信息构建的识别特征并不能完备表征雷达图像中目标的本质属性。在光学区，扩展目标的高频电磁散射响应可以用一组独立分布的散射体，或称散射中心的电磁散射响应之和近似表示。目标的散射中心主要产生于目标的边缘、拐点、棱角及尖端等不连续点部位，代表了目标的精细物理结构，所以散射中心模型能够更贴切地描述目标属性，也在雷达目标识别领域有着重要的应用。

基于几何绕射理论和物理光学理论，1999年Michael J.Gerry和Lee C.Potter提出了一个适用于合成孔径雷达(SAR)的参数化模型—属性散射中心模型，见[M.J.Gerry，L.C.Potter，I.J.Gupta，and A.van der Merwe，A parametric model for syntheticaperture radar measurements[J].IEEE Transactions on Antennas and Propagation，1999，Vol.47，NO.7，pp.1179-1188]。属性散射中心模型用一组参数描述每个散射中心的位置、形状、方向以及幅度等，这些属性都是关系目标的重要信息；同点散射模型相比，属性散射中心模型包含了更丰富的可用于目标分类识别的特征。

目标散射中心的特征提取本质上是一个从目标回波数据中估计各个散射中心参数的过程。由于属性散射中心模型结构复杂以及参数维数较高，增加了模型参数估计的复杂性。现有方法通过对以点散射模型为基础得到的雷达图像进行图像分割，得到阶数较低的目标散射区或者是孤立的散射中心，利用近似最大似然方法估计目标的属性散射中心参数。首先由于现有方法利用根据点散射模型得到的雷达图像提取属性散射中心，所以存在模型失配问题；其次由于该类方法是以图像分割为基础的，所以这类方法要求图像质量较高，此外当目标某些部件散射强度较弱时，通过图像分割方法是难以正确检测的，这就导致目标的一些重要特征容易丢失。除此以外由于该类方法的优化问题为非凸问题，并且存在很多局部最小解，所以存在参数初始化、模型阶数选择和散射中心结构类别判别的问题，使得最终参数估计精度较低。

发明内容

本发明的目的在于提出一种基于稀疏分解的雷达目标属性散射中心特征提取的方法，以解决现有方法中存在的模型失配、特征易丢失和参数估计精度较低问题。

本发明是这样实现的：

一．技术思路

雷达回波中，目标散射场绝大部分能量仅由少量强散射中心贡献，说明雷达回波在属性散射中心的参数空间具有很强的稀疏性。考虑属性散射中心参数空间维数较高，导致构造的冗余字典维数远远大于我们可处理的维数，本发明的方法基于稀疏分解理论，并利用坐标轮回下降技术降低问题维数并构建超分辨字典，通过求解0范数优化问题提取属性散射中心特征；建立散射中心强度门限并在雷达图像中进行强散射中心检测，结合先验信息确定散射中心参数的取值集合；根据属性散射中心模型，利用坐标轮回下降法技术构建超分辨字典，得到属性散射中心参数估计集合与目标及其重要部件几何尺寸特征。

二．技术方案

本发明的实现步骤包括如下：

1)根据雷达图像中噪声设定散射中心强度门限ξ，将雷达图像中强度大于ξ的散射中心确定为强散射中心，根据检测到的强散射中心坐标(x，y)确定坐标参数x，y的取值范围，由先验信息确定散射中心长度L取值范围，设定分布式散射中心方位角最终确定散射中心参数的集合

2)根据所述的参数集合构建字典

2a)将字典的参数集合离散化，即将相邻坐标参数x的间隔设为一个距离分辨单元长度ρ_r，将相邻坐标参数y的间隔设为一个方位分辨单元长度ρ_a，将相邻散射中心长度参数L的间隔设为ρ_a；

2b)根据属性散射中心模型，产生对应不同参数的原子d_i(f，φ)

$d_i(f,\mathrm{\φ})=\exp \left(\frac{-j4\mathrm{\πf}}{c}(x_i\cos \mathrm{\φ}+y_i\sin \mathrm{\φ})\right)\·$

$\sin c\left(\frac{2\mathrm{\&pi;f}}{c}{L}_i\sin (\mathrm{\&phi;}-\stackrel{\&OverBar;}{{\mathrm{\&phi;}}_i})\right)---<1>$

i＝1，…N₀

其中i表示原子序号，N₀表示初步估计(x，y，L)时字典的原子个数，exp(·)为自然指数函数，sinc(·)为辛克函数，为集合离散化后的第i组参数，f为雷达发射信号频率，φ为雷达波束方位角，c为光速；

将不同原子d_i(f，φ)列向量化，并对原子向量进行能量归一化，构建字典为：

$D_0(x,y,L|{\widehat{\mathrm{\&theta;}}}_0)=[{\hat{d}}_1,...,{\hat{d}}_i,...,{\hat{d}}_{N_0}]$

${\hat{d}}_i=\frac{\mathrm{vec}\left(d_i(f,\mathrm{\&phi;})\right)}{{\left|\right|\mathrm{vec}\left(d_i(f,\mathrm{\&phi;})\right)\left|\right|}_2}---<2>$

其中vec(·)表示列向量化操作，||·||₂为2范数算子；

3)利用正交匹配追踪OMP方法求解如式<3>的0范数优化问题，由σ₁更新散射中心参数集合 为散射中心参数(x，y,L)的初步估计集合；

${\mathrm{\&sigma;}}_1=\arg \underset{\mathrm{\&sigma;}}{\min }{\left|\right|\mathrm{\&sigma;}\left|\right|}_0,$ $s.t.{\left|\right|s-D_0(x,y,L|{\widehat{\mathrm{\&theta;}}}_0)\&CenterDot;\mathrm{\&sigma;}\left|\right|}_2\&le;\mathrm{\&epsiv;}---<3>$

其中，σ为待优化散射系数向量，σ₁为优化的散射系数向量，s为频域观测信号矩阵的列向量化，||·||₀为0范数算子，ε为能量误差约束因子，根据目标支撑区占整幅雷达图像的能量比确定；

4)根据所述的参数集合构建字典

4a)根据先验信息确定散射中心长度L取值范围，将相邻散射中心长度参数L的间隔设为ρ_a，得到字典的散射中心长度参数L的离散取值；根据雷达回波数据录取的方位角域确定分布式散射中心方位角取值范围，将相邻分布式散射中心方位角参数的间隔设为得到字典的分布式散射中心方位角参数的离散取值；散射中心坐标参数

4b)根据离散化的散射中心长度参数L与分布式散射中心方位角由式<1>产生原子向量，并对原子进行能量归一化，构建字典为：

$D_1(L,\stackrel{\&OverBar;}{\mathrm{\&phi;}}|{\widehat{\mathrm{\&theta;}}}_1)=[{\hat{d}}_1,...,{\hat{d}}_i,...,{\hat{d}}_{N_1}]---<4>$

其中N₁表示估计参数时字典的原子个数；

5)利用正交匹配追踪OMP方法求解如式<5>的0范数优化问题，由σ₂更新散射中心参数集合 为参数的估值集合；

${\mathrm{\&sigma;}}_2=\arg \underset{\mathrm{\&sigma;}}{\min }{\left|\right|\mathrm{\&sigma;}\left|\right|}_0,$ $s.t.{\left|\right|s-D_1(L,\stackrel{\&OverBar;}{\mathrm{\&phi;}}|{\widehat{\mathrm{\&theta;}}}_1)\&CenterDot;\mathrm{\&sigma;}\left|\right|}_2\&le;\mathrm{\&epsiv;}---<5>$

其中σ₂为优化的散射系数向量；

6)根据所述的参数集合构建超分辨字典

6a)在以中的坐标参数x为中心，长度为ρ_r的邻域内取值，将相邻参数x的间隔设为ρ_r/N_s，得到字典的坐标参数x的离散取值，N_s为超分辨倍数，一般N_s=2，4，8，...；在以中的y为中心，长度为ρ_a的邻域内取值，将相邻参数y的间隔设为ρ_a/N_s，得到字典的坐标参数y的离散取值；在以中的L为中心，长度为ρ_a的邻域内取值，将相邻参数L的间隔设为ρ_a/N_s，得到字典的散射中心长度参数L的离散取值；分布式散射中心方位角参数

6b)根据离散化的字典参数x，y，L由式<1>产生原子向量，并对原子进行能量归一化，构建超分辨字典

$D_2(x,y,L|{\widehat{\mathrm{\&theta;}}}_2)=[{\hat{d}}_1,...,{\hat{d}}_i,...,{\hat{d}}_{N_2}]---<6>$

其中，N₂表示超分辨字典的原子个数；

7)利用正交匹配追踪OMP方法求解如式<7>的0范数优化问题，由σ₃更新散射中心参数集合得到(x，y,L)的超分辨估计集合与目标及其重要部件几何尺寸特征；

${\mathrm{\&sigma;}}_3=\arg \underset{\mathrm{\&sigma;}}{\min }{\left|\right|\mathrm{\&sigma;}\left|\right|}_0,$ $s.t.{\left|\right|s-D_2(x,y,L|{\widehat{\mathrm{\&theta;}}}_2)\&CenterDot;\mathrm{\&sigma;}\left|\right|}_2\&le;\mathrm{\&epsiv;}---<7>$

其中σ₃为优化的散射系数向量。

本发明与现有技术相比具有以下优点：

(1)现有提取属性散射中心特征的方法存在模型失配问题，并且由于图像分割会导致目标的一些重要散射中心易于丢失。本发明提出的基于稀疏分解的雷达目标属性散射中心特征提取方法从频域观测数据出发，通过求解0范数优化问题得到频域观测数据的稀疏表示，可以有效提取目标散射中心。

(2)现有方法基于对雷达图像分割的结果，进行参数初始化和模型阶数选择，通过近似最大似然方法估计雷达图像子区域中的属性散射中心的参数，这类方法估计的参数精度较低。本发明方法通过利用坐标轮回下降技术构建超分辨字典，对频域观测数据进行稀疏分解，得到雷达目标属性散射中心参数的超分辨估计。

(3)由于本发明方法能够有效提取目标散射中心，并得到散射中心参数的超分辨估计，所以本发明根据提取的散射中心参数集合，可以得到目标及其重要部件的精确几何尺寸特征。

附图说明

图1是本发明的流程图；

图2是属性散射中心模型中频率依赖因子α=1与频率依赖因子α=0时的距离频率域数据变化波形；

图3是属性散射中心模型中方位依赖因子γ=10^-11与方位依赖因子γ=0时的方位频率域数据变化波形；

图4是雷达目标T72坦克的原始雷达图像；

图5是通过本发明方法对图4重构的雷达目标T72坦克的雷达图像；

图6是原始雷达目标D7铲车雷达图像；

图7是通过本发明方法对图6重构的雷达目标D7铲车雷达图像。

具体实施方式

一、技术原理

传统雷达图像是以点散射模型为基础得到的，该模型只包含目标散射点位置信息，但仅利用目标散射点的位置信息构建的识别特征并不能完备表征雷达图像中目标的本质属性。在光学区，扩展目标的高频电磁散射响应可以用一组独立分布的散射体，或称散射中心的电磁散射响应之和近似表示。属性散射中心模型用一组参数描述每个散射中心的位置、形状、方向以及幅度等，这些属性都是关系目标的重要信息；同点散射模型相比，属性散射中心模型包含了更丰富的可用于目标分类识别的特征。

根据属性散射中心模型可知，目标中第i个散射中心频率-方位二维回波信号为：

$E_i(f,\mathrm{\&phi;};{\mathrm{\&theta;}}_i)=A_i\&CenterDot;{\left(j\frac{f}{f_c}\right)}^{{\mathrm{\&alpha;}}_i}\exp \left(\frac{-j4\mathrm{\&pi;f}}{c}(x_i\cos \mathrm{\&phi;}+y_i\sin \mathrm{\&phi;})\right)\&CenterDot;---<8>$

$\sin c\left(\frac{2\mathrm{\&pi;f}}{c}{L}_i\sin (\mathrm{\&phi;}-\stackrel{\&OverBar;}{{\mathrm{\&theta;}}_i})\right)\exp (-2\mathrm{\&pi;f}{\mathrm{\&gamma;}}_i\sin \mathrm{\&phi;})$

其中，i表示散射中心序号，f为雷达发射信号频率，φ为雷达方位角，exp(·)为自然指数函数，sin c(·)为辛克函数，c为光速，θ_i表示第i个散射中心的参数向量Ai为散射中心的散射强度，x_i为距离维坐标，y_i为方位维坐标，L_i为分布式散射中心方位维的长度，为分布式散射中心的方位角，α_i为频率依赖因子，一般α_i∈{-1，-0.5,0,0.5,1}，γ_i为局部式散射中心的方位依赖因子。

由各个散射中心的回波信号之和，即可构成目标频率-方位二维回波信号：

$E(f,\mathrm{\&phi;};\mathrm{\&theta;})=\underset{i=1}{\overset{M}{\mathrm{\&Sigma;}}}{E}_i(f,\mathrm{\&phi;};{\mathrm{\&theta;}}_i),$ ${\mathrm{\&theta;}}^T=[{\mathrm{\&theta;}}_1^T,...,{\mathrm{\&theta;}}_i,...,{\mathrm{\&theta;}}_M^T]---<9>$

其中i表示散射中心序号，E_i(f，φ；θ_i)为第i个散射中心的回波信号，M为散射中心个数，θ表示M个散射中心参数矩阵，(·)^T表示转置操作；

将式<9>的目标回波信号用矩阵形式表达，其表达式为：

s=D(θ)σ+n    <10>

其中s为目标回波信号E(f，φ；θ)的列向量化，D(θ)为散射中心参数矩阵θ对应的字典，σ为散射系数向量，n表示高斯白噪声。

在雷达回波中，由于目标散射场的绝大部分能量仅由少量强散射中心贡献，故说明雷达回波在属性散射中心的参数空间具有很强的稀疏性。根据雷达回波在属性散射中心的参数空间的稀疏性，通过求解0范数优化问题即可得到观测数据s的稀疏表示，以及散射中心参数估计集合，通过散射中心参数集合可以估计目标及其重要部件的几何尺寸特征。

由于属性散射中心参数空间维数较高，导致构造的冗余字典维数远远大于需要处理的维数，因此可通过分析属性散射中心模型中的频率依赖因子α、方位依赖因子γ对回波信号的影响，对属性散射中心模型进行简化，降低参数空间维数。

当雷达的中心频率f_c=9.6GHz，带宽B=591MHz，其它散射中心参数为0时，仅仅改变频率依赖因子α或者方位依赖因子γ，而频率依赖因子α与方位依赖因子γ对频域回波数据有如下影响：

a)由于雷达相对带宽较小，频率依赖因子α使散射中心的散射强度在距离频率维可发生微小变化，如图2所示，其中图2(a)中频率依赖因子α=1，图2(b)中频率依赖因子α=0。

b)由于雷达方位角φ变化范围较小，局部式散射中心方位依赖因子γ的变化能使局部式散射中心的散射强度在方位频率上稍有变化，如图3所示，其中图3(a)中方位依赖因子γ=10^-11，图3(b)中方位依赖因子γ=0。

通过以上分析可知，频率依赖因子α与方位依赖因子γ对目标频域回波数据影响很小，因此估计参数时可以忽略频率依赖因子α与方位依赖因子γ，即认为频率依赖因子α=0，方位依赖因子γ=0，则原始属性散射中心模型式<8>可简化为：

${\hat{E}}_i(f,\mathrm{\&phi;};{\mathrm{\&theta;}}_i)=A_i\&CenterDot;\exp \left(\frac{-j4\mathrm{\&pi;f}}{c}(x_i\cos \mathrm{\&phi;}+y_i\sin \mathrm{\&phi;})\right)\&CenterDot;$ $<11>$

$\sin c\left(\frac{2\mathrm{\&pi;f}}{c}{L}_i\sin (\mathrm{\&phi;}-\stackrel{\&OverBar;}{{\mathrm{\&theta;}}_i})\right)$

根据简化后的属性散射中心模型，基于稀疏分解，利用坐标轮回下降技术，可以得到目标散射中心参数的超分辨估计，通过散射中心参数集合可以估计目标及其重要部件的几何尺寸特征。

二、实现步骤

参照图1，本发明的实现步骤如下：

步骤1，建立散射中心强度门限ξ进行强散射中心检测，确定散射中心参数取值范围。

1a)雷达成像：

由式<11>可知，雷达回波信号在频率-方位极坐标空间下均匀采样，雷达成像将回波信号重采样或者插值，得到在笛卡尔坐标系(f_x，f_y)下均匀采样的回波信号E(f_x，f_y)，其中f_x为距离频率，f_y为方位频率；

f_x＝fcos(φ)，f_y＝fsin(φ)        <12>

其中，φ为雷达方位角；

对回波信号E(f_x，f_y)做二维逆傅里叶变换，得到目标的雷达图像；

1b)强散射中心检测，确定散射中心参数取值范围：

根据雷达图像中噪声设定散射中心强度门限ξ，将雷达图像中强度大于ξ的散射中心确定为强散射中心，根据检测到的强散射中心坐标(x，y)确定坐标参数x，y的取值范围，由先验信息确定散射中心长度L取值范围，设定分布式散射中心方位角最终确定散射中心参数的集合

例如，当雷达图像数据为公开MSTAR实测数据中方位角为80.774185°俯仰角为15°的T72坦克数据时，根据雷达图像中噪声设定散射中心强度门限ξ＝0.28，确定的散射中心距离坐标参数x的取值范围为[11.1675,14.467]米，方位坐标参数y的取值范围为[13.40625,4.875]米，散射中心长度L取值范围为[0，4.0625]米，分布式散射中心方位角

步骤2，根据所述的参数集合构建字典

2a)将字典的参数集合离散化，即将相邻坐标参数x的间隔设为一个距离分辨单元长度ρ_r，将相邻坐标参数y的间隔设为一个方位分辨单元长度ρ_a，将相邻散射中心长度参数L的间隔设为ρ_a；

例如，当雷达图像数据为公开MSTAR实测数据中方位角为80.774185°俯仰角为15°的T72坦克数据时，根据雷达图像中噪声设定散射中心强度门限ξ＝0.28，确定的散射中心距离坐标参数x的取值范围为[11.1675，14.467]米，将相邻坐标参数x间隔为0.2538米，方位坐标参数y的取值范围为[13.40625,4.875]米，将相邻坐标参数y的间隔设为0.203125米，散射中心长度L取值范围为[0，4.0625]米，将相邻散射中心长度参数L的间隔设为0.203125米，分布式散射中心方位角

2b)根据属性散射中心模型，产生对应不同参数的原子d_i(f，φ)：

$d_i(f,\mathrm{\&phi;})=\exp \left(\frac{-j4\mathrm{\&pi;f}}{c}(x_i\cos \mathrm{\&phi;}+y_i\sin \mathrm{\&phi;})\right)\&CenterDot;$

$\sin c\left(\frac{2\mathrm{\&pi;f}}{c}{L}_i\sin (\mathrm{\&phi;}-\stackrel{\&OverBar;}{{\mathrm{\&phi;}}_i})\right)---<13>$

i＝1，…N₀

其中，f为雷达发射信号频率，φ为雷达方位角，exp(·)为自然指数函数，sin c(·)为辛克函数，c为光速，i表示原子序号，N₀表示初步估计(x，y,L)时字典的原子个数，为集合离散化后的第i组参数。

2c)将不同原子d_i(f，φ)列向量化，并对原子向量进行能量归一化，构建字典

$D_0(x,y,L|{\widehat{\mathrm{\&theta;}}}_0)=[{\hat{d}}_1,...,{\hat{d}}_i,...,{\hat{d}}_{N_0}]$

${\hat{d}}_i=\frac{\mathrm{vec}\left(d_i(f,\mathrm{\&phi;})\right)}{{\left|\right|\mathrm{vec}\left(d_i(f,\mathrm{\&phi;})\right)\left|\right|}_2}$ <14>

i=1，…N₀

其中，i表示原子序号，N₀表示初步估计(x，y,L)时字典的原子个数，vec(·)表示列向量化操作，||·||₂为2范数算子。

步骤3，利用正交匹配追踪OMP求解如式<15>的0范数优化问题：

${\mathrm{\&sigma;}}_1=\arg \underset{\mathrm{\&sigma;}}{\min }{\left|\right|\mathrm{\&sigma;}\left|\right|}_0,$ $s.t.{\left|\right|s-D_0(x,y,L|{\widehat{\mathrm{\&theta;}}}_0)\&CenterDot;\mathrm{\&sigma;}\left|\right|}_2\&le;\mathrm{\&epsiv;}---<15>$

其中σ为待优化散射系数向量，σ₁为优化的散射系数向量，s为频域观测信号矩阵的列向量化，||·||₀为0范数算子，ε为能量误差约束因子，其根据目标支撑区占整幅雷达图像的能量比确定。

本步骤的求解步骤如下：

3a)输入字典

3b)初始化，即将信号余量r初始化为频域观测信号矩阵的列向量s，将散射系数向量σ初始化为0，将原子索引向量a初始化为空，设初始迭代次数k=1，并开始迭代；

3c)计算字典与信号余量r的相关系数向量C

$C=D_0{(x,y,L|{\widehat{\mathrm{\&theta;}}}_0)}^H\&CenterDot;r---<16>$

其中(·)^H表示共轭转置，相关系数向量C中最大元素的序号为a_k，可得原子索引向量a为：

a=[a₁，…，a_k]          <17>

3d)利用最小二乘方法计算散射系数向量σ

其中表示伪逆，D₀(：，a)为字典D₀中列号为a_i的k个列向量构成的矩阵i＝1，…，k，k为当前迭代次数；

更新余量r：

r＝s-D₀(：，a)·σ(a)      <19>

其中σ(a)为散射系数向量σ中序号为a_i的k个元素构成的列向量i＝1，…，k；

3e)计算重构能量比η_k；

${\mathrm{\&eta;}}_k=\frac{\mathrm{\&Sigma;}{\left|\right|s_k\left|\right|}_2^2}{\mathrm{\&Sigma;}{\left|\right|s\left|\right|}_2^2}---<20>$

s_k=D₀(：，a)·σ(a)

其中s_k为重构信号，表示2范数的平方；

3f)如果η_k-η_k-1≤δ停止迭代，σ₁=σ，其中δ为相邻迭代能量比之差的门限，取δ=0.001；否则k=k+1并转至步骤3c)。

步骤4，由优化的散射系数向量σ₁更新散射中心参数集合 为散射中心参数(x，y,L)的初步估计集合。

步骤5，根据所述的参数集合构建字典

5a)根据先验信息确定散射中心长度L取值范围，将相邻散射中心长度参数L的间隔设为ρ_a，得到字典的散射中心长度参数L的离散取值；根据雷达回波数据录取的方位角域确定分布式散射中心方位角取值范围，将相邻分布式散射中心方位角参数的间隔设为得到字典的分布式散射中心方位角参数的离散取值；散射中心坐标参数

例如，当雷达图像数据为公开MSTAR实测数据中方位角为80.774185°俯仰角为15°的T72坦克数据时，散射中心长度L取值范围为[0，4.0625]米，将相邻散射中心长度参数L的间隔设为0.203125米，分布式散射中心方位角的取值范围为[-3°，3°]，将相邻分布式散射中心方位角参数的间隔设为0.25°。

5b)根据离散化的散射中心长度参数L与分布式散射中心方位角由式<13>产生原子向量，并对原子向量进行能量归一化，构建字典为：

$D_1(L,\stackrel{\&OverBar;}{\mathrm{\&phi;}}|{\widehat{\mathrm{\&theta;}}}_1)=[{\hat{d}}_1,...,{\hat{d}}_i,...,{\hat{d}}_{N_1}]---<21>$

其中，i表示原子序号，i＝1，…N₁，N₁表示估计参数时字典的原子个数。

步骤6，利用正交匹配追踪OMP求解如式<22>的0范数优化问题：

${\mathrm{\&sigma;}}_2=\arg \underset{\mathrm{\&sigma;}}{\min }{\left|\right|\mathrm{\&sigma;}\left|\right|}_0,$ $s.t.{\left|\right|s-D_1(L,\stackrel{\&OverBar;}{\mathrm{\&phi;}}|{\widehat{\mathrm{\&theta;}}}_1)\&CenterDot;\mathrm{\&sigma;}\left|\right|}_2\&le;\mathrm{\&epsiv;}---<22>$

其中σ₂为优化的散射系数向量，其求解步骤如下：

6a)输入字典

6b)初始化，即将信号余量r初始化为频域观测信号矩阵的列向量s，将散射系数向量σ初始化为0，将原子索引向量a初始化为空，设初始迭代次数k=1，并开始迭代；

6c)计算字典与信号余量r的相关系数向量C

$C=D_1{(L,\stackrel{\&OverBar;}{\mathrm{\&phi;}}|{\widehat{\mathrm{\&theta;}}}_1)}^H\&CenterDot;r---<23>$

其中(·)^H表示共轭转置，相关系数向量C中最大元素的序号为a_k，可得原子索引向量a为：

a=[a₁，…，a_k]        <24>

6d)利用最小二乘方法计算散射系数向量σ

其中表示伪逆，D₁(：，a)为字典D₁中列号为a_i的k个列向量构成的矩阵i＝1，…，k，k为当前迭代次数；

更新余量r：

r=s-D₁(：，a)·σ(a)       <26>

其中σ(a)为散射系数向量σ中序号为a_i的k个元素构成的列向量i＝1，…，k；

6e)计算重构能量比η_k；

${\mathrm{\&eta;}}_k=\frac{\mathrm{\&Sigma;}{\left|\right|s_k\left|\right|}_2^2}{\mathrm{\&Sigma;}{\left|\right|s\left|\right|}_2^2}---<27>$

s_k=D₁(：，a)·σ(a)

其中s_k为重构信号，表示2范数的平方；

6f)如果η_k-η_k-1≤δ停止迭代，σ₂=σ，其中δ为相邻迭代能量比之差的门限，取δ=0.001；否则k=k+1并转至步骤6c)。

步骤7，由优化的散射系数向量σ₂更新散射中心参数集合 为参数的估值集合。

步骤8，根据所述的参数集合构建超分辨字典

8a)确定字典的距离维坐标参数x离散取值：在以参数集合中的距离维坐标参数x为中心，长度为ρ_r的邻域内取值，将相邻参数x的间隔设为ρ_r/N_s，得到字典的距离维坐标参数x的离散取值，其中N_s为超分辨倍数，一般N_s=2，4，8，...；

8b)确定字典的方位维坐标参数y离散取值：

在以参数集合中的方位维坐标参数y为中心，长度为ρ_a的邻域内取值，将相邻参数y的间隔设为ρ_a/N_s，得到字典的方位维坐标参数y的离散取值；

8c)确定字典的散射中心长度参数L离散取值：

在以参数集合中的散射中心长度参数L为中心，长度为ρ_a的邻域内取值，将相邻参数L的间隔设为ρ_a/N_s，得到字典的散射中心长度参数L的离散取值；

8d)确定字典的分布式散射中心方位角参数离散取值：

以参数集合中的分布式散射中心方位角参数为字典的分布式散射中心方位角参数

8e)根据上述步骤确定的距离维坐标参数x，方位维坐标参数y，散射中心长度参数L，分布式散射中心方位角参数离散取值，由式<13>产生与离散化的字典参数组合对应的原子向量，并对原子进行能量归一化，构建超分辨字典为：

$D_2(x,y,L|{\widehat{\mathrm{\&theta;}}}_2)=[{\hat{d}}_1,...,{\hat{d}}_i,...,{\hat{d}}_{N_2}]---<28>$

其中，i表示原子序号，i＝1...N₂，N₂表示超分辨字典的原子个数。

步骤9，利用正交匹配追踪OMP方法求解如式<29>的0范数优化问题：

${\mathrm{\&sigma;}}_3=\arg \underset{\mathrm{\&sigma;}}{\min }{\left|\right|\mathrm{\&sigma;}\left|\right|}_0,$ $s.t.{\left|\right|s-D_2(x,y,L|{\widehat{\mathrm{\&theta;}}}_2)\&CenterDot;\mathrm{\&sigma;}\left|\right|}_2\&le;\mathrm{\&epsiv;}---<29>$

其中σ₃为优化的散射系数向量，其求解步骤如下：

9a)输入字典

9b)初始化，即将信号余量r初始化为频域观测信号矩阵的列向量s，将散射系数向量σ初始化为0，将原子索引向量a初始化为空，设初始迭代次数k=1，并开始迭代；

9c)计算字典与信号余量r的相关系数向量C

$C=D_2{(x,y,L|{\widehat{\mathrm{\&theta;}}}_2)}^H\&CenterDot;r---<30>$

其中(·)^H表示共轭转置，相关系数向量C中最大元素的序号为a_k，得到原子索引向量a为：

a=[a₁，…，a_k]      <31>

9d)利用最小二乘方法计算散射系数向量σ

其中表示伪逆，D₂(：，a)为字典D₂中列号为a_i的k个列向量构成的矩阵i＝1，…，k，k为当前迭代次数；

更新余量r：

r=s-D₂(：，a)·σ(a)         <33>

其中σ(a)为散射系数向量σ中序号为a_i的k个元素构成的列向量i＝1，…，k；

9e)计算重构能量比η_k：

${\mathrm{\&eta;}}_k=\frac{\mathrm{\&Sigma;}{\left|\right|s_k\left|\right|}_2^2}{\mathrm{\&Sigma;}{\left|\right|s\left|\right|}_2^2}---<34>$

s_k=D₂(：，a)·σ(a)

其中s_k为重构信号，表示2范数的平方；

9f)如果η_k-η_k-1≤δ停止迭代，σ3=σ其中δ为相邻迭代能量比之差的门限，取δ=0.001；否则k=k+1并转至步骤9c)。

步骤10，由优化的散射系数向量σ₃更新散射中心参数集合得到(x，y,L)的超分辨估计集合根据散射中心参数集合估计目标及其重要部件几何尺寸特征。

例如，当雷达图像数据为公开MSTAR实测数据中方位角为80.774185°俯仰角为15°的T72坦克数据时，根据最终提取的散射中心及其参数集合估计T72坦克的几何尺寸为：T72坦克炮筒长度为6.4米，坦克体长度为6.8米，T72坦克的总长度为9.4米，T72坦克的宽度估值为3.8米。

本发明的效果通过以下对实测数据的实验进一步说明：

1)实验场景：

实验所用的实测数据为公开的MSTAR数据集中方位角为80.774185°俯仰角为15°的T72坦克的合成孔径雷达SAR数据以及方位角为0.307442°俯仰角为17°的D7铲车的SAR数据，雷达的中心频率f_c=9.6GHz，带宽B=591MHz。T72坦克真实几何尺寸为：车体长6.41m，车宽3.52m，炮筒长6.155m，炮向前时炮筒伸出长度3.035m。D7铲车真实几何尺寸为：车宽2.56m，车长4.2m，车高3.294m，履带板宽度0.56m，履带中心距1.981m，接地履带长度2.72m。2)实验内容：

2a)对于方位角为80.774185°俯仰角为15°的T72坦克MSTAR数据，利用坐标轮回下降法构建超分辨字典，通过求解0范数优化问题得到散射中心参数估计集合。图4为原始T72坦克SAR图像，根据提取的散射中心及其参数重构的T72坦克图像如图5所示，由图4与图5可知本发明能够准确重构目标图像。

2b)对于方位角为0.307442°俯仰角为17°的D7铲车MSTAR数据，利用坐标轮回下降法构建超分辨字典，通过求解0范数优化问题得到散射中心参数估计集合。图6为原始D7铲车SAR图像，根据提取的散射中心及其参数重构的T72坦克图像如图7所示，由图6与图7可知本发明能够准确重构目标图像。

3)实验结果分析：

表1本发明方法处理T72坦克数据收敛过程

表2本发明方法处理D7铲车数据收敛过程

表1以及表2分别给出了本发明的方法处理T72坦克和D7铲车数据过程中重构能量比、字典原子个数以及散射中心个数的收敛情况。由于步骤3假定而仅仅对参数(x，y,L)进行初步估计，所以步骤3的重构能量较低；当步骤6在步骤3的基础上考虑分布式散射中心方位角时，重构能量比得到提高，提取的散射中心个数降低；步骤9中字典的参数分辨率增加，重构能量比得到提高，散射中心个数也进一步减少，字典原子个数得到降低，这是由于步骤3至步骤9之间利用了坐标轮回下降技术。表1、表2中数据收敛过程表明本发明方法能够有效提取散射中心。

根据最终提取的散射中心及其参数可以估计T72坦克与D7铲车的几何尺寸。

表3用于估计T72坦克炮筒长度的散射中心参数

表3表明，根据散射中心的坐标和散射中心长度估计坦克炮筒长度为6.4米。

表4用于估计T72坦克体长度的散射中心参数

表4表明，根据散射中心的坐标参数和散射中心长度估计坦克体长度为6.8米。

由目标支撑区可以得到T72坦克的总长度估值为9.4米；根据T72坦克上边缘的散射中心至炮筒散射中心的距离为T72坦克宽度的一半，得到T72坦克的宽度估值为3.8米。

表5用于估计D7铲车铲刀长度的散射中心参数

表5表明，根据散射中心的坐标参数和散射中心长度估计D7铲车铲刀长为3米。

表6用于估计D7铲车宽度的散射中心参数

表6表明，根据散射中心的坐标参数和散射中心长度估计D7铲车宽度为2.6米。

表7用于估计D7铲车履带宽度的散射中心参数

表7表明，根据散射中心的坐标参数和散射中心长度估计D7铲车履带宽度为0.6米。

由提取的履带散射中心可以估计D7铲车的接地履带长度为2.6米；由目标支撑区得到D7铲车车长估值为4.31米。

结合T72坦克与D7铲车的真实几何尺寸，可知本发明方法估计的目标及其重要部件的几何尺寸误差在6％以内，实验结果说明本发明方法能够有效地准确提取目标属性散射中心特征和精确估计目标及其重要部件的几何尺寸。

Claims (4)

1.一种基于稀疏分解的雷达目标属性散射中心特征提取方法，包括如下步骤：

1）根据雷达图像中噪声设定散射中心强度门限ξ，将雷达图像中强度大于ξ的散射中心确定为强散射中心，根据检测到的强散射中心坐标(x,y)确定坐标参数x,y的取值范围，由先验信息确定散射中心长度L取值范围，设定分布式散射中心方位角最终确定散射中心参数的集合

2）根据所述的参数集合构建字典

2a）将字典的参数集合离散化，即将相邻坐标参数x的间隔设为一个距离分辨单元长度ρ_r，将相邻坐标参数y的间隔设为一个方位分辨单元长度ρ_a，将相邻散射中心长度参数L的间隔设为ρ_a；

2b）根据属性散射中心模型，产生对应不同参数的原子d_i(f,φ)

$\begin{array}{c}{d}_i(f,\mathrm{\&phi;})=\exp \left(\frac{-j4\mathrm{\&pi;f}}{\text{c}}(x_i\cos \mathrm{\&phi;}+y_i\sin \mathrm{\&phi;})\right)\&CenterDot;\\ \sin c\left(\frac{2\mathrm{\&pi;f}}{c}{L}_i\sin (\mathrm{\&phi;}-{\stackrel{\&OverBar;}{\mathrm{\&phi;}}}_i)\right)\\ i=1,\&CenterDot;\&CenterDot;\&CenterDot;N_0\end{array}---<1>$

其中i表示原子序号，N₀表示初步估计(x,y,L)时字典的原子个数，exp(·)为自然指数函数，sinc(·)为辛克函数，为集合离散化后的第i组参数，f为雷达发射信号频率，φ为雷达波束方位角，c为光速；

将不同原子d_i(f,φ)列向量化，并对原子向量进行能量归一化，构建字典 $D_0(x,y,L|{\widehat{\mathrm{\&theta;}}}_0)$ 为：

$\begin{array}{c}{D}_0(x,y,L|{\widehat{\mathrm{\&theta;}}}_0)=[{\hat{d}}_1,\&CenterDot;\&CenterDot;\&CenterDot;,{\hat{d}}_i,\&CenterDot;\&CenterDot;\&CenterDot;,{\hat{d}}_{N_0}]\\ {\hat{d}}_i=\frac{\mathrm{vec}\left(d_i\left(\mathrm{fm\&phi;}\right)\right)}{\left|\right|\mathrm{vec}\left(d_i(f,\mathrm{\&phi;})\right)|{|}_2}\end{array}---<2>$

其中vec(·)表示列向量化操作，||·||₂为2范数算子；

3）利用正交匹配追踪OMP方法求解如式<3>的0范数优化问题，由σ₁更新散射中心参数集合为散射中心参数(x,y,L)的初步估计集合；

${\mathrm{\&sigma;}}_1=\arg \underset{\mathrm{\&sigma;}}{\min }{\left|\right|\mathrm{\&sigma;}\left|\right|}_0,s.t.{\left|\right|s-D_0(x,y,L|{\widehat{\mathrm{\&theta;}}}_0)\&CenterDot;\mathrm{\&sigma;}||}_2\&le;\mathrm{\&epsiv;}---<3>$

其中，σ为待优化散射系数向量，σ₁为优化的散射系数向量，s为频域观测信号矩阵的列向量化，||·||₀为0范数算子，ε为能量误差约束因子，根据目标支撑区占整幅雷达图像的能量比确定；

4）根据所述的参数集合构建字典

4a）根据先验信息确定散射中心长度L取值范围，将相邻散射中心长度参数L的间隔设为ρ_a，得到字典的散射中心长度参数L的离散取值；根据雷达回波数据录取的方位角域确定分布式散射中心方位角取值范围，将相邻分布式散射中心方位角参数的间隔设为得到字典的分布式散射中心方位角参数的离散取值；散射中心坐标参数

4b）根据离散化的散射中心长度参数L与分布式散射中心方位角由式<1>产生原子向量，并对原子进行能量归一化，构建字典

$D_1(L,\stackrel{\&OverBar;}{\mathrm{\&theta;}}|{\widehat{\mathrm{\&theta;}}}_1)=[{\hat{d}}_1,\&CenterDot;\&CenterDot;\&CenterDot;,{\hat{d}}_i,\&CenterDot;\&CenterDot;\&CenterDot;,{\hat{d}}_{N_0}]---<4>$

其中N₁表示估计参数时字典的原子个数；

5）利用正交匹配追踪OMP方法求解如式<5>的0范数优化问题，由σ₂更新散射中心参数集合为参数的估值集合；

${\mathrm{\&sigma;}}_2=\arg \underset{\mathrm{\&sigma;}}{\min }{\left|\right|\mathrm{\&sigma;}\left|\right|}_0,s.t.{\left|\right|s-D_1(L,\stackrel{\&OverBar;}{\mathrm{\&phi;}}|{\widehat{\mathrm{\&theta;}}}_1)\&CenterDot;\mathrm{\&sigma;}\left|\right|}_2\&le;\mathrm{\&epsiv;}---<5>$

其中σ₂为优化的散射系数向量；

6）根据所述的参数集合构建超分辨字典

6a）在以中的坐标参数x为中心，长度为ρ_r的邻域内取值，将相邻参数x的间隔设为ρ_r/N_s，得到字典的坐标参数x的离散取值，N_s为超分辨倍数，N_s＝2,4,8,...；在以中的y为中心，长度为ρ_a的邻域内取值，将相邻参数y的间隔设为ρ_a/N_s，得到字典的坐标参数y的离散取值；在以中的L为中心，长度为ρ_a的邻域内取值，将相邻参数L的间隔设为ρ_a/N_s，得到字典的散射中心长度参数L的离散取值；分布式散射中心方位角参数

6b）根据离散化的字典参数x,y,L由式<1>产生原子向量，并对原子进行能量归一化，构建超分辨字典

$D_2(x,y,L|{\widehat{\mathrm{\&theta;}}}_2)=[{\hat{d}}_1,\&CenterDot;\&CenterDot;\&CenterDot;,{\hat{d}}_i,\&CenterDot;\&CenterDot;\&CenterDot;,{\hat{d}}_{N_2}]---<6>$

其中，N₂表示超分辨字典的原子个数；

7）利用正交匹配追踪OMP方法求解如式<7>的0范数优化问题，由σ₃更新散射中心参数集合得到(x,y,L)的超分辨估计集合与目标及其重要部件几何尺寸特征；

${\mathrm{\&sigma;}}_3=\arg \underset{\mathrm{\&sigma;}}{\min }{\left|\right|\mathrm{\&sigma;}\left|\right|}_0,s.t.{\left|\right|s-D_2(x,y,L|{\widehat{\mathrm{\&theta;}}}_2)\&CenterDot;\mathrm{\&sigma;}\left|\right|}_2\&le;\mathrm{\&epsiv;}---<7>$

其中σ₃为优化的散射系数向量。

2.根据权利要求1所述的方法，其中步骤3所述的利用正交匹配追踪OMP方法求解如式<3>的0范数优化问题，按如下步骤进行：

3a）输入字典 $D_0(x,y,L|{\widehat{\mathrm{\&theta;}}}_0);$

3b）初始化，即将信号余量r初始化为频域观测信号矩阵的列向量s，将散射系数向量σ初始化为0，将原子索引向量a初始化为空，设初始迭代次数k＝1，并开始迭代；

3c）计算字典与信号余量r的相关系数向量C

$C=D_0{(x,y,L|{\widehat{\mathrm{\&theta;}}}_0)}^H\&CenterDot;r---<8>$

其中(·)^H表示共轭转置，相关系数向量C中最大元素的序号为a_k，可得原子索引向量a为：

a＝[a₁,…,a_k]             <9>

3d）利用最小二乘方法计算散射系数向量σ

                     <10>

其中表示伪逆，D₀(:,a)为字典D₀中列号为a_i的k个列向量构成的矩阵i=1,…,k，k为当前迭代次数；

更新余量r：

r=s-D₀(:,a)·σ(a)            <11>

其中σ(a)为散射系数向量σ中序号为a_i的k个元素构成的列向量i=1,…,k；

3e）计算重构能量比η_k；

$\begin{array}{c}{\mathrm{\&eta;}}_k=\frac{\mathrm{\&Sigma;}{\left|\right|s_k\left|\right|}_2^2}{\mathrm{\&Sigma;}{\left|\right|s\left|\right|}_2^2}\\ s_k=D_0(:,a)\&CenterDot;\mathrm{\&sigma;}\left(a\right)\end{array}---<12>$

其中s_k为重构信号，表示2范数的平方；

3f）如果η_k-η_k-1≤δ停止迭代，其中δ为相邻迭代能量比之差的门限，取δ＝0.001；否则k＝k+1并转至步骤3c）。

3.根据权利要求1所述的方法，其中步骤5所述的利用正交匹配追踪OMP方法求解如式<5>的0范数优化问题，按如下步骤进行：

5a）输入字典 $D_1(L,\stackrel{\&OverBar;}{\mathrm{\&theta;}}|{\widehat{\mathrm{\&theta;}}}_1);$

5b）初始化，即将信号余量r初始化为频域观测信号矩阵的列向量s，将散射系数向量σ初始化为0，将原子索引向量a初始化为空，设初始迭代次数k＝1，并开始迭代；

5c）计算字典与信号余量r的相关系数向量C

$C=D_1{(L,\stackrel{\&OverBar;}{\mathrm{\&theta;}}|{\widehat{\mathrm{\&theta;}}}_1)}^H\&CenterDot;r---<13>$

其中(·)^H表示共轭转置，相关系数向量C中最大元素的序号为a_k，可得原子索引向量a为：

a＝[a₁,…,a_k]             <14>

5d）利用最小二乘方法计算散射系数向量σ

                    <15>

其中表示伪逆，D₁(:,a)为字典D₁中列号为a_i的k个列向量构成的矩阵i=1,…,k，k为当前迭代次数；

更新余量r：

r=s-D₁(:,a)·σ(a)            <16>

其中σ(a)为散射系数向量σ中序号为a_i的k个元素构成的列向量i=1,…,k；

5e）计算重构能量比η_k；

$\begin{array}{c}{\mathrm{\&eta;}}_k=\frac{\mathrm{\&Sigma;}{\left|\right|s_k\left|\right|}_2^2}{\mathrm{\&Sigma;}{\left|\right|s\left|\right|}_2^2}\\ s_k=D_1(:,a)\&CenterDot;\mathrm{\&sigma;}\left(a\right)\end{array}---<17>$

其中s_k为重构信号，表示2范数的平方；

5f）如果η_k-η_k-1≤δ停止迭代，其中δ为相邻迭代能量比之差的门限，取δ＝0.001；否则k＝k+1并转至步骤5c)。

4.根据权利要求1所述的方法，其中步骤7所述的利用正交匹配追踪OMP方法求解如式<7>的0范数优化问题，按如下步骤进行：

7a）输入字典 $D_2(x,y,L|{\widehat{\mathrm{\&theta;}}}_2);$

7b）初始化，即将信号余量r初始化为频域观测信号矩阵的列向量s，将散射系数向量σ初始化为0，将原子索引向量a初始化为空，设初始迭代次数k＝1，并开始迭代；

7c）计算字典与信号余量r的相关系数向量C

$C=D_2{(x,y,L|{\widehat{\mathrm{\&theta;}}}_2)}^H\&CenterDot;r---<18>$

其中(·)^H表示共轭转置，相关系数向量C中最大元素的序号为a_k，得到原子索引向量a为：

a＝[a₁,…,a_k]           <19>

7d）利用最小二乘方法计算散射系数向量σ

             <20>

其中表示伪逆，D₂(:,a)为字典D₂中列号为a_i的k个列向量构成的矩阵i=1,…,k，k为当前迭代次数；

更新余量r：

r=s-D₂(:,a)·σ(a)             <21>

其中σ(a)为散射系数向量σ中序号为a_i的k个元素构成的列向量i=1,…,k；

7e）计算重构能量比η_k：

$\begin{array}{c}{\mathrm{\&eta;}}_k=\frac{\mathrm{\&Sigma;}{\left|\right|s_k\left|\right|}_2^2}{\mathrm{\&Sigma;}{\left|\right|s\left|\right|}_2^2}\\ s_k=D_2(:,a)\&CenterDot;\mathrm{\&sigma;}\left(a\right)\end{array}---<22>$

其中s_k为重构信号，表示2范数的平方；

7f）如果η_k-η_k-1≤δ停止迭代，其中δ为相邻迭代能量比之差的门限，取δ＝0.001；否则k＝k+1并转至步骤7c)。