CN102890275B - 一种基于交叉熵最小化的压缩感知雷达成像方法 - Google Patents

Abstract

本发明公开了一种基于交叉熵最小化的压缩感知雷达成像方法，包括：1)根据雷达平台与目标场景的几何关系模型建立压缩感知雷达的回波与目标场景之间观测方程，初始化雷达系统参数和迭代参数；2)根据观测方程更新迭代参数和目标场景，计算场景目标雷达散射截面积幅度的直方图分布；3)根据步骤2)的直方图分布计算目标雷达散射截面积幅度中满足交叉熵最小化的阈值；4)用步骤3)中计算得到的阈值采用硬阈值方法更新目标场景，判断重建场景是否满足特定的终止条件，如满足，则输出目标场景，否则重复步骤2～4重建场景满足条件。本发明有效解决目标场景中的杂波和系统观测噪声对稀疏场景重建的影响，实现压缩感知雷达成像。

Description

一种基于交叉熵最小化的压缩感知雷达成像方法

技术领域

本发明涉及雷达成像技术领域，是一种基于交叉熵最小化的压缩感知雷达成像方法，适用于稀疏场景和变换稀疏场景的目标重建。

背景技术

雷达成像技术是从接收到的回波重建场景目标分布的逆滤波问题。近一二十年来，雷达成像技术得到了长足的发展，为环境监测、资源勘察，目标识别、战场态势评估提供了新的途径。随着雷达成像技术的发展以及其应用需求的推动，传统的基于奈奎斯特采样定理的雷达成像数据采样方式最终导致系统采集的数据速率急剧增大。

近几年，出现了一种新的采样理论-压缩采样或称压缩感知(CS，compressive sampling or compressed sensing)，该方法是在采样过程的同时实现信号压缩，即以低于奈奎斯特速率的采样率进行采样，并能以极高的准确率恢复出原信号(参见《Compressed sensing》，刊于IEEETransactions on Information Theory，vol.52，Apr.2006，pp.1289-1306)。

基于压缩感知的信号处理可以降低数据的采样率。利用压缩感知技术对原始雷达回波数据进行成像处理的一个重要前提是目标场景具有稀疏特性。对于稀疏场景信号x，可以对回波信号进行欠采样(即以低于奈奎斯特速率进行采样)，然后利用凸优化的方法重建目标场景。在实现时，通过引入M×N(K＜M＜N)的采样矩阵Φ和雷达系统矩阵Ψ，得到雷达回波可以表示为：

y＝ΦΨx

由于M＜N，y即为采样并压缩后的信号。在接收端，通过凸优化的方法对目标场景重建。目前，已有很多文献提出了对信号进行重建的方法，如Basis Pursuit  算法(参见《Matching pursuit with time-frequencydictionaries》，刊于IEEE Transactions on Signal Processing，vol.41，Dec.1993，pp.3397-3415)、Orthogonal Matching Pursuit算法(参见《Signal Recovery from Random Measurements via Orthogonal MatchingPursuit》，刊于IEEE Transactions on Information Theory，vol.53，Dec.2007，pp.4655-4666)等等。

基于压缩感知的信号处理可以降低数据的采样率。基于这个理论，Baraniuk将CS应用到雷达的一维成像和二维成像中(参见《CompressiveRadar Imaging》，刊于IEEE Radar Conference，Boston，MA，USA，Apr.17-20，2007，128-133)，Herman引入Alltop序列作为二维CS成像雷达系统的发射脉冲信号(参见《High-Resolution Radar via Compressed Sensing》，刊于IEEE Transactions on Signal Processing，vol.57，June 2009，pp.2275-2284)。

目前的重建算法用于对压缩感知雷达目标场景进行重建时，一般是通过对目标场景的稀疏目标数，或者对目标场景的噪声水平进行估计，然后通过正则项控制重建精度，这些方法在实际应用中由于先验信息不足导致实用性较差，研究一种适用于实际压缩感知雷达成像的算法是从事压缩感知雷达研究的学者们十分关注的一个问题。

发明内容

本发明的目的是公开一种基于交叉熵最小化的压缩感知雷达成像方法，在雷达成像处理中，利用交叉熵最小化的方法区别图像中的目标和噪声，达到了利用很少回波信号实现对目标场景的准确成像，消除了压缩感知雷达成像中杂波和噪声对成像结果的影响，克服了目前常用的稀疏重建算法处理压缩感知雷达成像时对目标场景中的噪声和杂波处理存在的不足。

为达到上述目的，本发明的技术解决方案是：

一种基于交叉熵最小化的压缩感知雷达成像方法，其包括步骤：

(1)根据雷达平台与目标场景的几何关系模型，建立压缩感知雷达的回波与目标场景之间的观测方程，并初始化雷达系统参数和迭代参数；

(2)根据观测方程更新迭代参数和目标场景，并计算场景目标雷达散射截面积幅度的直方图分布；

(3)根据目标场景的目标雷达散射截面积幅度的直方图分布，计算目标雷达散射截面积幅度中满足交叉熵最小化的阈值；

(4)利用步骤(3)中计算得到的阈值，采用硬阈值的方法更新目标场景，并判断重建场景是否满足特定的终止条件，如果满足条件，则输出目标场景，否则，重复步骤(2)～(4)直到重建场景满足条件。

所述的基于交叉熵最小化的压缩感知雷达成像方法，其所述步骤(1)中，根据雷达平台与目标场景的几何关系模型建立的压缩感知雷达的回波与目标场景之间的观测方程为：

y＝ΦΨx+n＝Θx+n

其中，Ψ为系统矩阵，Φ为观测矩阵，y为获取的采样时间序列，x为具有稀疏特性的观测目标场景，Θ为测量矩阵，n为噪声向量；

初始化雷达系统参数和迭代参数：雷达系统参数的初始化是根据雷达系统设计时的参数进行初始化，目标雷达散射截面初始化为任意非零向量，迭代步数初始化为0，其他参数以初始化为特定值，然后根据观测方程建立重建目标函数f(x)为：

f(x)＝||y-Θx||²

并计算目标函数f(x)的梯度的最小莱布尼兹常数L(f)＝2λ_max(Θ^TΘ)，其中，Θ^T表示Θ的转置，λ_max(Θ^TΘ)表示取矩阵Θ^TΘ的最大特征值。

所述的基于交叉熵最小化的压缩感知雷达成像方法，其所述步骤(2)中，根据观测方程更新迭代参数和目标场景，并计算场景目标雷达散射截面积幅度的直方图分布，分为两步：

1)根据第k步迭代的场景目标雷达散射截面向量x^k，采用梯度法更新迭代场景目标雷达散射截面向量u^k为：

$u^k=x^k-\frac{1}{L}\▿f\left(x^k\right)$

其中，

表示目标函数f(x)的梯度在x＝x^k处的值；

2)根据更新迭代的场景目标雷达散射截面向量u^k，计算出场景目标雷达散射截面向量u^k幅度的直方图分布，具体步骤如下：

A、计算场景目标雷达散射截面向量u^k中模值最大值的元素与模值最小的元素

的差值：

$d^k=u_{\max }^k-u_{\min }^k$

B、将场景目标雷达散射截面向量u^k按模值从小到大分成L个集合，即：

$u_{\min }^k~u_{\min }^k+d^k/L,\·\·\·,u_{\min }^k+(L-1)d^k/L~u_{\max }^k$

然后计算出场景目标雷达散射截面向量u^k模值分布的统计直方图h。

所述的基于交叉熵最小化的压缩感知雷达成像方法，其所述步骤(3)中，根据目标场景的目标雷达散射截面积幅度的直方图分布，计算目标雷达散射截面积幅度中满足交叉熵最小化的阈值，是利用迭代的方法对步骤(2)中计算出的目标场景的直方图分布h进行迭代处理，计算出满足交叉熵最小化的最优阈值μ。

所述的基于交叉熵最小化的压缩感知雷达成像方法，其所述步骤(4)中，利用步骤(3)中计算得到的最优的阈值μ，采用硬阈值的方法更新目标场景，并判断重建场景是否满足特定的终止条件，如果目标场景满足条件，则输出目标场景，否则重复步骤(2)～(4)直到重建场景满足条件，输出目标场景。

所述的基于交叉熵最小化的压缩感知雷达成像方法，其所述终止条件，根据需要设定，最终是实现对场景目标雷达散射截面成像。

本发明方法与现有技术相比优点在于：

(1)本发明充分利用基于交叉熵最小化的理论可以很好的实现目标和背景的分离，将交叉熵最小化的理论引入到现在的压缩感知雷达成像处理中，可以有效的抑制回波中的杂波和背景噪声；

(2)本发明通过将交叉熵最小化的方法和迭代阈值算法结合在一起进行雷达成像处理，可以加快收敛速度，提高重建精度。

附图说明

图1为本发明的基于交叉熵最小化的压缩感知雷达成像方法的流程图；

图2为利用本发明的基于交叉熵最小化的压缩感知雷达成像方法重建的仿真结果图；其中：

图2a为真实的目标场景；

图2b为处理结果示意图；

图3为基于交叉熵最小化的压缩感知雷达成像方法与之前的迭代阈值算法重建精度比较图，其中，纵坐标的为目标函数f(x)＝||y-Θx||²。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚，下面结合附图和实施例仿真情况对本发明作进一步的详细描述。

本发明的一种基于交叉熵最小化的压缩感知雷达成像方法，利用基于交叉熵最小化的理论可以很好的将稀疏信号处理中的将稀疏信号与噪声分离的思想，将交叉熵最小化理论应用到压缩感知雷达成像处理中，提出了一种基于交叉熵最小化的压缩感知雷达成像处理算法。如图1所示，为本发明的一种基于交叉熵最小化的压缩感知雷达成像方法实施流程图，具体包含以下4个步骤：

(1)根据雷达平台与目标场景的几何关系模型建立的压缩感知雷达的回波与目标场景向量之间的观测方程：

y＝ΦΨx+n＝Θx+n                     (1)

其中y为获取的采样时间序列；Φ为观测矩阵，由回波信号的采样方式决定；Ψ为系统矩阵，主要由系统的观测几何决定；x为具有稀疏特性的观测目标场景向量，Θ为由系统矩阵和观测矩阵共同构成的测量矩阵，n为噪声向量.

根据前面构造的观测方程，初始化雷达系统参数和迭代参数：雷达系统参数的初始化可以根据雷达系统设计时的参数进行初始化；对迭代参数进行初始化时，目标场景可以初始化为任意非零向量x⁰，迭代步数k初始化为0，其它一些迭代参数可以初始化为任意值，然后根据观测方程建立重建目标函数f(x)为：

f(x)＝||y-Θx||²                         (2)

并计算目标函数f(x)的梯度的最小莱布尼兹常数L(f)＝2λ_max(Θ^TΘ)，其中Θ^T表示Θ的转置，λ_max(Θ^TΘ)表示取矩阵Θ^TΘ的最大特征值；

(2)根据观测方程更新迭代参数和目标场景，并计算场景目标雷达散射截面积幅度的直方图分布的实现过程分为两步：

1)根据第k步迭代的场景目标雷达散射截面向量x^k，采用梯度法更新迭代场景目标雷达散射截面向量u^k为：

$u^k=x^k-\frac{1}{L}\▿f\left(x^k\right)---\left(3\right)$

其中表示目标函数f(x)的梯度在x＝x^k处的值.

2)根据更新迭代的场景目标雷达散射截面向量u^k，计算出场景目标雷达散射截面积u^k的幅度直方图分布，具体步骤如下：

a.计算场景目标雷达散射截面向量向量u^k中模值最大值的元素

与模值最小的元素的差值：

$d^k=u_{\max }^k-u_{\min }^k---\left(4\right)$

b.将场景目标雷达散射截面向量向量u^k按模值从小到大分成L个集合，即：

$u_{\min }^k~u_{\min }^k+d^k/L,\·\·\·,u_{\min }^k+(L-1)d^k/L~u_{\max }^k---\left(5\right)$

然后计算出场景目标雷达散射截面向量u^k幅度分布的统计直方图h；

(3)根据目标场景的雷达散射截面积向量u^k的幅度的直方图分布计算目标雷达散射截面积幅度中满足交叉熵最小化的阈值μ^k的实现过程如下。

设有L个灰度等级的目标场景的雷达散射截面积向量u^k的幅度的直方图h的前景和背景在直方图中的零阶矩和一阶矩分别为：

$m_{0a}\left(t\right)=\underset{i=1}{\overset{t-1}{\mathrm{\Σ}}}h\left(i\right),$ $m_{0b}\left(t\right)=\underset{i=t}{\overset{L}{\mathrm{\Σ}}}h\left(i\right),$

(6)

$m_{1a}\left(t\right)=\underset{i=1}{\overset{t-1}{\mathrm{\Σ}}}i\·h\left(i\right),$ $m_{1b}\left(t\right)=\underset{i=t}{\overset{L}{\mathrm{\Σ}}}i\·h\left(i\right).$

它们的比例定义为：

${\mathrm{\μ}}_a\left(t\right)=\frac{m_{1a}\left(t\right)}{m_{0a}\left(t\right)},$ ${\mathrm{\μ}}_b\left(t\right)=\frac{m_{1b}\left(t\right)}{m_{0b}\left(t\right)}---\left(7\right)$

定义评判准则函数为：

η(t)＝-m_1a(t)·log(μ_a(t))-m_1b(t)·log(μ_b(t))                        (8)

最优的阈值t_op是通过最小化(8)式求得的，即：

$t_{\mathrm{op}}=\underset{t}{\arg \min }\mathrm{\η}\left(t\right)---\left(9\right)$

通过计算可以得到最优的阈值为：

$t=\frac{{\mathrm{\μ}}_b\left(t\right)-{\mathrm{\μ}}_a\left(t\right)}{\log \left({\mathrm{\μ}}_b\left(t\right)\right)-\log \left({\mathrm{\μ}}_a\left(t\right)\right)}---\left(10\right)$

Li等人给出了一个如下在离散情况下通过迭代计算最优阈值的方法：

$t_{n+1}=\mathrm{round}\left(\frac{{\mathrm{\μ}}_b\left(t_n\right)-{\mathrm{\μ}}_a\left(t_n\right)}{\log \left({\mathrm{\μ}}_b\left(t_n\right)\right)-\log \left({\mathrm{\μ}}_a\left(t_n\right)\right)}\right)---\left(11\right)$

给一个初始假设的t₀直到迭代收敛，收敛的条件是t_n+1＝t_n，选择t_opt＝t_n，这里round(x)函数表示与x最接近的整数。(参见《An Iterative Algorithmfor Minimum Cross Entropy Thresholding》，刊于Pattern RecognitionLetters，vol.19，March 1998，pp.771-776)

那么满足交叉熵最小化的目标雷达散射截面积幅度阈值μ^k通过下式选取：

${\mathrm{\μ}}^k=u_{\min }^k+(t_{\mathrm{opt}}-1)d^k/L---\left(12\right)$

(4)利用步骤(3)中计算得到的阈值μ^k采用硬阈值的方法更新目标场景向量并判断重建场景向量是否满足特定的终止条件步骤如下：

a.更新迭代步数，k→k+1，利用步骤(3)中算得到的阈值μ^k采用如下公式中的方法更新目标场景：

$x^{k+1}=\mathrm{hard}(u^k,{\mathrm{\μ}}^k)$

$=\left\{\begin{array}{cc}0& u_i^k<{\mathrm{\&mu;}}^k\\ u_i^k& u_i^k\&GreaterEqual;{\mathrm{\&mu;}}^k\end{array}\right.,i=1,\&CenterDot;\&CenterDot;\&CenterDot;,n---\left(13\right)$

b.判断步骤a中得到的更新重建场景向量x^k+1是否满足的终止条件，这里关于终止条件的选择主要有两种方法：

一种是满足目标函数的在x＝x^k+1的值小于特定的值tol(这里tol指预先设定好的很小的量)，即：

f(x^k+1)＝||y-Θx^k+1||²≤tol                            (14)

第二种是上一步的迭代的结果与本次迭代的结果的相对变化量ζ小于特定值tol(这里tol指预先设定好的很小的量)，即：

$\mathrm{\&zeta;}=\frac{{\left|\right|x^k-x^{k+1}\left|\right|}^2}{{\left|\right|x^k\left|\right|}^2}\&le;\mathrm{tol}---\left(15\right)$

选取上面的一种条件作为终止条件是，如果目标场景向量满足条件则输出目标场景向量，否则重复步骤(2)～(4)直到重建场景向量满足条件，最终实现对场景目标雷达散射截面成像。

下面通过仿真对本发明方法进行验证。实验中仿真的是一个机载平台的雷达回波数据，发射信号为线性调频脉冲信号，实验中用到的一些参数为：载频为10GHz，脉冲宽度为1μs，机载速度为150m/s，斜距为15000m，场景的分辨率为：方位向的分辨率2m，和距离向的分辨率为2m，按照奈奎斯特采样时，目标场景的大小选取的是31×31，场景中包括30个幅度和相位均随机的目标信号和一些幅度微弱的杂波信号，真实的目标场景如图2(a)所示，并在回波信号中加入信噪比为20dB的随机噪声，采用本算法进行成像处理是只随机采集10％的基于奈奎斯特采样采集的回波数据，利用该算法可以准确重建目标场景，处理结果如图2(b)所示。

选取与前面仿真式样相同的目标场景和雷达参数，分别对存在20dB噪声情况的回波信号在只随机采集10％的基于奈奎斯特采样采集的回波数据的10％的情况下采用我们这里提出的基于交叉熵最小化的迭代算法CEM-IHT与之前经典的IST算法以及改进的IST算法IST_opt比较，进行100次蒙特卡洛仿真后，发现在相同的迭代次数下，CEM-IHT具有较高的精度，具体比较情况如图3所示，所以基于交叉熵最小化的压缩感知雷达成像方法是一种可用于实际雷达成像处理的算法。

以上所述，仅为本发明中的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉该技术的人在本发明所揭露的技术范围内，可理解想到的变换或替换，都应涵盖在本发明包含的范围之内，因此，本发明的保护范围以权利要求书的保护范围为准。

Claims (6)

1.一种基于交叉熵最小化的压缩感知雷达成像方法，其特征在于：包括步骤：

(1)根据雷达平台与目标场景的几何关系模型，建立压缩感知雷达的回波与目标场景之间的观测方程，并初始化雷达系统参数和迭代参数；

(2)根据观测方程更新迭代参数和目标场景，并计算场景目标雷达散射截面积幅度的直方图分布；

(3)根据目标场景的目标雷达散射截面积幅度的直方图分布，利用迭代方法计算目标雷达散射截面积幅度中满足交叉熵最小化的阈值；

(4)利用步骤(3)中计算得到的阈值，采用硬阈值的方法更新目标场景，并判断重建场景是否满足特定的终止条件，如果满足条件，则输出目标场景，否则，重复步骤(2)～(4)直到重建场景满足条件。

2.根据权利要求1所述的基于交叉熵最小化的压缩感知雷达成像方法，其特征在于：所述步骤(1)中，根据雷达平台与目标场景的几何关系模型建立的压缩感知雷达的回波与目标场景之间的观测方程为：

y＝ФΨx+n＝Θx+n

其中，y为获取的采样时间序列，Φ为观测矩阵，Ψ为系统矩阵，x为具有稀疏特性的观测目标场景向量，Θ为测量矩阵，n为噪声向量；

初始化雷达系统参数和迭代参数：雷达系统参数的初始化是根据雷达系统设计时的参数进行初始化，目标雷达散射截面初始化为任意非零向量，迭代步数初始化为0，其他参数以初始化为特定值，然后根据观测方程建立重建目标函数f(x)为：

f(x)＝||y-Θx||²

并计算目标函数f(x)的梯度的最小莱布尼兹常数L(f)＝2λ_max(Θ^TΘ)，其中Θ^T表示Θ的转置，λ_max(Θ^TΘ)表示取矩阵Θ^TΘ的最大特征值。

3.根据权利要求1所述的基于交叉熵最小化的压缩感知雷达成像方法，其特征在于：所述步骤(2)中，根据观测方程更新迭代参数和目标场景，并计算场景目标雷达散射截面积幅度的直方图分布，分为两步：

1)根据第k步迭代的场景目标雷达散射截面向量x^k，采用梯度法更新迭代场景目标雷达散射截面向量u^k为：

$u^k=x^k-\frac{1}{L}\&dtri;f\left(x^k\right)$

其中，

表示目标函数f(x)的梯度在x＝x^k处的值；

2)根据更新迭代的结果，计算出场景目标雷达散射截面向量u^k的直方图分布，具体步骤如下：

A、计算场景目标雷达散射截面向量u^k中模值最大的元素

与模值最小的元素

的差值d^k，即：

$d^k=u_{\max }^k-u_{\min }^k$

B、将场景目标雷达散射截面向量u^k按模值从小到大分成L个集合，即：

$u_{\min }^k~u_{\min }^k+d^k/L,...,u_{\min }^k+(L-1)d^k/L~u_{\max }^k$

然后计算出场景目标雷达散射截面向量u^k模值分布的统计直方图h。

4.根据权利要求1所述的基于交叉熵最小化的压缩感知雷达成像方法，其特征在于：所述步骤(3)中，根据目标场景的目标雷达散射截面积幅度的直方图分布，利用迭代方法计算目标雷达散射截面积幅度中满足交叉熵最小化的阈值，是对步骤(2)中计算出的场景目标雷达散射截面积幅度的直方图分布进行迭代处理，计算出目标雷达散射截面积幅度中满足交叉熵最小化的阈值μ。

5.根据权利要求1所述的基于交叉熵最小化的压缩感知雷达成像方法，其特征在于：所述步骤(4)中，利用步骤(3)中计算得到的阈值，采用硬阈值的方法更新目标场景，并判断重建场景是否满足特定的终止条件，如果目标场景满足条件，则输出目标场景，否则重复步骤(2)～(4)直到重建场景满足条件，输出目标场景。

6.根据权利要求5所述的基于交叉熵最小化的压缩感知雷达成像方法，其特征在于：所述终止条件，根据需要设定，最终是实现对场景目标雷达散射截面成像。

Images