朱燕丽，周 鹏✉，张振华，张 晰，张 杰，王 影

（1 中国石油大学（华东）海洋与空间信息学院 青岛 266580 2 北京遥测技术研究所 北京 100076 3 自然资源部第一海洋研究所 青岛 266061）

引言

逆合成孔径雷达ISAR 利用非合作运动目标相对于雷达的转角进行高分辨成像。它能够全天候、全天时地工作，在军事和民用领域具有重要的应用价值[1]。在ISAR 系统中，距离向的高分辨是通过大时间带宽积的波形实现的，方位向的高分辨一般通过长相干积累时间内雷达观测目标的转角而获得。然而，在现实世界中，观察环境可能很复杂。可能会有强噪声源污染雷达回波，导致信噪比SNR 低。也可能很难获得连续较长的积累时间，导致方位分辨率较低，同时造成脉冲积累后的信噪比相对较低。当ISAR系统的方位分辨率较低时，将严重降低ISAR 图像质量，对目标的分类与识别精度造成非常不利的影响，因此研究ISAR 系统的超分辨率成像算法是非常必要的。

目前，研究ISAR 超分辨率成像的思路主要基于压缩感知CS（Compressive Sensing）理论[2]。基于CS的成像算法以其结构简单、图像质量高、计算复杂度低的优势，被认为是一种优秀的数据处理算法。如果信号是“稀疏的”，它就可以从很少的观测数据中被精确地重构出来。CS的框架主要分为三个部分：测量矩阵、稀疏表示、重构算法。其中重构算法的选取对重构精度和重构速度都有很大的影响。常用的重构算法包括凸优化算法、贪婪算法和统计学习算法等几大类。目前常用的线性LBI 算法[3]、OMP算法[4]、SBL 算法[5]分别是上述三类算法的典型代表。然而，上述算法在低信噪比的情况下都存在着重构性能不佳、计算复杂度高等不足。

近年来，快速分裂Bregman 迭代算法FSBI（Fast Split Bregman Iteration）因其抗噪性能强的优点，被广泛用于解决大规模优化问题，在压缩感知领域取得了良好的效果[6]。受此启发，为实现低信噪比下的ISAR 快速超分辨成像，本文提出了一种基于快速分裂Bregman 迭代的ISAR 超分辨成像算法。首先，在正则化框架下引入L1范数表示目标的稀疏先验信息，将方位分辨率提高问题转化为正则化问题；其次，为实现快速重构，将托普利兹矩阵的低位移秩特征[7]和Gohberg-Semencul 表示[8]相结合，实现了矩阵的快速求逆。通过上述两方面的措施，既能利用分裂Bregman 迭代算法在低信噪比条件下的重构能力[9]，又能保证快速成像。最后，通过仿真和实测数据处理验证了本文提出算法的性能。

1 压缩感知ISAR 成像模型的建立

经过脉冲压缩和平动补偿后的某一距离单元的回波可以表示为：

式中，K是散射点的数量，Ak为第k个散射点的回波幅度，tm为慢时间，MT为方位向时间长度，M为脉冲个数，fk为第k个散射点相对于参考点的多普勒频率，n(tm)表示加性噪声。

为实现方位向的超分辨率成像，分别定义多普勒采样间隔和时间采样间隔为Δfd和Δt，脉冲重复频率为fr，M=TM/Δt。离散化慢时间序列表示为tm=[1:M]T·Δt，Δt=1/fr，多普勒域为[-fr/2,fr/2]，各多普勒单元的散射点的强度用向量σ表示为σ=，σq（q=0,1,···,Q-1）表示位于第q个多普勒单元内的散射系数，=qΔfd-fr/2。只有当多普勒单元内有散射点时，σq才不为0。由于包含散射点的多普勒单元的数量很少，σ可以视为稀疏的。

式中，s为方位向信号向量，n为噪声向量。

在压缩感知理论中，对稀疏信号s的观测不是直接测量信号s本身，而是通过非相关测量将信号s投影到一组低维的观测向量，即有：

式中，y为观测向量，Φ为测量矩阵，M′为观测向量中元素的个数。由于M′＜M，直接求解式（3）是一个欠定问题，所以无法直接从y中解出信号s。将式（2）代入式（3）中，有：

式中，=A ΦΨ为感知矩阵，E=Φn为噪声矩阵。可以通过求解式（4）得到稀疏系数向量σ，从而最终重构出稀疏信号s。

可以将式（4）转化为如下优化问题：

式（5）可以利用传统的线性规划方法进行求解，但是噪声的扰动可能导致估算值远离真实值。因此将式（5）转化为以下正则化形式：

式中，μ是正则化参数，f(σ) 表示正则化项。

通常，对正则项的选择取决于目标先验信息。可以引入合理的正则项来提高分辨率。在ISAR 成像场景中，通常只关注强散射点目标。在稀疏假设下，将式（6）转化为以下形式：

通过正则化项‖σ1‖ 控制的稀疏性，通过-‖Aσ y控制的误差，同时保证了解的稀疏性和准确性，而且能较好地克服噪声的影响。

由于L1范数和L2范数都是凸函数，可以通过凸优化方法对式（7）求解。然而L1范数是不可微的，所以求解它具有挑战性。

2 快速分裂Bregman 迭代ISAR 超分辨率成像算法

在本节中，利用托普利兹矩阵的低位移秩特征和Gohberg-Semencul 定理对分裂Bregman 迭代SBI（Split Bregman Iteration）算法进行加速，以实现L1正则化问题的快速求解。

2.1 分裂Bregman 迭代算法

为了求解式（7），首先引入辅助变量g来代替σ，式（7）变为：

将上面的约束问题进一步转化为下面的非约束优化问题：

式中，γ是大于0的常数。表面上，式（8）和式（9）看起来是多余的，然而，它们能极大地降低计算复杂度并获得更有效的迭代策略。

另一个重要的概念是Bregman 距离。基于凸函数J(·)上的点u和w的Bregman 距离定义为[10]：

式中，向量∂ J(u) 是J(u)的梯度。

式（11）中的最小化问题可以转化为对2 个子问题σk+1和gk+1分别进行优化，即有：

对子问题σk+1的求解可以采用多种方法。当问题规模不大时，对式（12）第1 个式子中的σ进行求导并令导数等于0，可以得到：

式（13）中涉及到求逆运算，计算复杂度较高。为简化计算，可以使用Gauss Seidel 迭代法来求解[11]。子问题gk+1可直接利用软阈值法得到最优解，其表达式如下：

其中：

式中，sign()表示符号函数。在求解式（14）的过程中，式（12）中的第2 个式子可直接对b进行更新。

2.2 快速迭代计算算法

虽然用SBI 求解L1正则化问题需要三个步骤，但主要的计算复杂度通常来自于σ的求解，即式（13）。g和b的求解，只涉及到一些简单的加减运算，与式（13）相比，它们的计算复杂性可以忽略不计。因此，这里主要分析式（13）的计算复杂度。经分析后可知，其计算复杂度是O(N3)。

实际上，可以通过降低矩阵求逆引起的高计算复杂度来降低式（13）的计算复杂度。首先，将式（13）改写为：

其中，F=μATA+λI，dk=μATy+γ（g k-bk）。从矩阵A和I的结构来看，可以看出矩阵F具有托普利兹结构。因此，F-1计算复杂度的降低可以通过使用合适的Gohberg-Semencul 表示来解决。并且，式（16）的计算可以通过使用快速托普利兹向量乘法更有效地实现。

基于Gohberg-Semencul 表示的加速策略可以实现系数矩阵的快速求逆，大大降低了矩阵求逆的计算复杂度。

首先，自回归系数a=（1,a2,···,aM）T和预测误差e可以通过Levinson-Durbin 算法从Yule-Walker AR方程获得[12]，即有：

基于Gohberg-Semencul 表示，托普利兹矩阵的求逆可以表示为[13]：

其中：

因此，式(16)的计算可以转化为：

由于U和V的特殊结构，UUHdk和VVHdk的计算可以通过快速傅里叶变换FFT（Fast Fourier Transform）和矩阵截取来实现。

为了实现快速计算，构造如下矩阵：

U可以通过截取U1的前M行获得，UH可以通过截取U2的第M行到第2M-1 行获得。因此，U与向量的乘积可视为矩阵U1与向量的乘积结果的第1 到第M行，即U和向量FFT的结果的第1 到第M个元素，UH与向量相乘也是同理。同样地，VVHdk也可以通过两次FFT 和矩阵截取来计算得到。

2.3 计算复杂度分析

本文FSBI 算法的计算复杂度分析如下。假设迭代次数为K，首先，利用Levinson-Durbin 算法计算自回归系数a和预测误差e，计算复杂度为O((N-1)2)。式（25）的解可以通过四个托普利兹向量运算来实现，计算复杂度为7φ(2N)+φ(N)+O(4N)，dk的计算复杂度为φ(N)+O(4N)。因此，经过加速后的算法的计算复杂度为O(K((N-1)2+4N))+7Kφ(2N)+φ(N)。而SBI算法的计算复杂度为O((K+1)N3+5KN2+3KN+NlogN)+φ(N)。计算复杂度从O(N3)降低到O(N2)，计算复杂度大大降低。

2.4 正则化参数的选取

通常采用L-curve曲线法来确定正则化参数[14]。图1给出利用L-curve曲线确定正则化参数的示意图，其中横坐标表示正则化项残差的范数，纵坐标表示解的范数，曲线上的拐点代表了最佳正则化参数的位置。

在计算式（12）的过程中，需要确定μ和γ两个正则化参数。在第一次迭代中，参数γ的初值是根据经验确定的，然后通过L-curve 曲线法获得正则化参数μ的初值。具体而言，是把作为的函数并作曲线，曲线上的点与正则化参数μ一一对应。通过曲线的拐点处的值获得对应μ的值。同理，参数γ可以通过式（12）的第二个式子确定，L-curve 曲线法绘制作为的函数，通过曲线的拐点处的值获得对应μ的值

3 实验

首先，针对一维信号，采用蒙特卡洛模拟的方法，评估并比较了OMP 算法、SBI（Split Bregman Iteration）算法、LBI 算法、本文的FSBI 算法在不同信噪比条件下的稀疏信号重构性能；接下来，针对短积累时间和数据随机缺失两种情况，利用仿真的方法生成了低信噪比条件下匀速转动船只的ISAR 模拟回波。评估并比较了RD 算法、OMP 算法、SBI 算法、本文FSBI 算法的成像结果与图像性能指标；最后，针对实测的船只ISAR 回波数据进行了超分辨率成像处理，评估并比较了OMP 算法、SBI 算法、LBI 算法、FSBI 算法的成像结果与图像性能指标。对于实测船只ISAR 回波数据的处理，同样分为了短积累时间和数据随机缺失两种情况。

3.1 一维稀疏信号重构的性能分析

假设原始信号x是一个1 024×1的随机稀疏向量，其中有50 个非零值。经乘以测量矩阵后降维为300×1的观测向量，通过添加高斯白噪声使信噪比从-5 dB 到20 dB 变化。对于每种信噪比水平，通过200次的蒙特卡洛模拟实验，计算重构值与真实值之间的均方根相对误差RMSRE（Root Mean Square Relative Error）。这里的RMSRE 定义为：

式中，i为元素序号，N为向量中的元素个数，xri为重构值，xti为真实值。

图2 给出了OMP 算法、SBI 算法、LBI 算法、本文FSBI 算法的RMSRE 随信噪比水平的变化。需说明的是，由于FSBI 算法和SBI 算法的原理相同，二者的RMSRE 曲线相同。从图中可以看出，本文FSBI 算法和SBI 算法在各种信噪比水平下的RMSRE值均最小。

表1 给出了四种算法的平均运行时间，其中LBI 算法的运行时间最长，接下来依次是SBI 算法、FSBI 算法和OMP 算法。因此，与现有常用算法相比，本文FSBI 算法可以在可接受的运行时间内获得最佳的稀疏重建性能。

表1 不同算法平均运行时间的比较Table 1 Comparison of the average running time for the four algorithms

3.2 仿真数据的实验结果

针对短积累时间和数据随机缺失两种情况，利用仿真的方法生成了低信噪比条件下匀速转动船只的ISAR 模拟回波。对OMP 算法、SBI 算法、LBI 算法和FSBI 算法的成像结果和成像性能指标进行了评估比较。仿真中使用的主要参数如表2 所示。仿真过程中使用的船只散射点模型如图3 所示。针对短时间积累的情况，通过加入高斯白噪声，使得距离压缩后的回波信噪比为-5 dB 和0 dB。针对数据随机缺失的情况，数据缺失比例为25%。

表2 仿真过程中使用的主要参数Table 2 Main parameters used in simulation

首先，测试了短积累时间情况下利用几种不同算法获得的ISAR 图像。图4 和图5 显示了该情况下的ISAR 回波和利用几种不同算法获得的ISAR 图像。其中，图4（a）为模拟得到的ISAR 距离压缩回波，图4（b）、（c）、（d）、（e）分别为利用RD 算法、OMP 算法、LBI 算法、本文FSBI 算法获得的ISAR 成像结果。因为FSBI 算法和SBI 算法的原理相同，获得的ISAR 图像相同。它们之间的唯一区别是计算复杂性。因此，图中结果仅列出了FSBI 算法的成像结果。从图中可以观察到，OMP、LBI、SBI、FSBI 算法等几种基于压缩感知的成像算法获得了比RD 算法更好的聚焦图像。对比几种压缩感知算法获得的图像发现，OMP 算法产生的图像虚假点相对较多，LBI 算法产生的图像较暗淡，而FSBI 算法生成的ISAR 图像获得了最清晰的船只结构。为了进一步定量比较四种算法成像结果的质量和计算复杂度，计算了图像熵并记录了运行时间。表3 列出了四种算法的图像熵和运行时间。从表中可以看出，在SNR=-5 dB的条件下，FSBI 算法和SBI 算法的图像熵最小，在SNR=0 dB的条件下，LBI 算法的熵值最小，FSBI 算法与LBI 算法的熵值接近。与此同时，FSBI 算法的运行时间也相对较短。因而取得了最佳的成像质量。

表3 短积累时间情况下不同算法成像结果的图像熵和运行时间比较Table 3 Comparison of image entropy and running time of four algorithms in the case of short accumulation time

接下来，测试了数据随机缺失情况下利用几种不同算法获得的ISAR 图像，设置数据缺失比例为25%。图6 显示了该情况下的ISAR 距离压缩回波和利用几种不同算法获得的ISAR 图像。其中，图6（a）为模拟得到的ISAR 距离压缩回波，在这种情况下，由于采样网格不连续，RD 算法失效。图6（b）、（c）、（d）分别为利用OMP 算法、LBI 算法、FSBI 算法获得的ISAR 成像结果。同样计算了几幅图像的图像熵和运行时间。表4 列出了对应的计算结果。对比计算结果可发现，FSBI 算法的图像熵值最低、运行时间相对较短。综合而言，FSBI 算法的成像性能最优。

表4 数据随机缺失情况下不同算法成像结果的图像熵和运行时间比较Table 4 Comparison of image entropy and running time of four algorithms in the case of random loss of data

3.3 实测数据的实验结果

本节利用舰船ISAR 实测数据来进一步验证本文FSBI 算法的有效性。雷达信号的中心频率为9.25 GHz，带宽为500 MHz，脉冲重复频率为100 Hz。完整的雷达回波包含256 个脉冲，每个脉冲由256 个采样值组成。针对短时间积累的情况，通过加入高斯白噪声，使得距离压缩后的回波信噪比为0 dB和5 dB。针对数据随机缺失的情况，数据缺失比例为20%。

图7 和图8 显示了短时间积累情况下的实测数据的ISAR 回波和利用几种不同算法获得的ISAR 图像。其中，两图中的子图（a）为实测数据的ISAR 回波，两图中的子图（c）、（d）、（e）分别为利用OMP 算法、LBI 算法、FSBI 算法获得的ISAR 成像结果。计算了几幅图像的图像熵并记录了运行时间。表5 列出了相应的图像熵和运行时间。对比计算结果可发现，FSBI 算法的图像熵值相对较低，因而成像质量最优且运行时间相对较短。

表5 短时间积累情况下不同算法对实测数据成像结果的指标值Table 5 Index values of imaging results of the real data corresponding to four algorithms in the case of short accumulation time

图9 显示了数据随机缺失情况下利用几种不同算法获得的ISAR 图像。数据缺失比例为20%。其中，图9（a）为ISAR 距离压缩回波，图9（b）、（c）、（d）分别为利用OMP 方法算法、LBI 算法、FSBI算法获得的ISAR 成像结果。同样计算了几幅图像的图像熵和运行时间。表6 列出了对应的计算结果。对比计算结果可发现，所提本文方法算法（FSBI）的图像熵值最低、运行时间相对较短。因而，综合而言，FSBI 算法的成像质量性能最优。

表6 数据随机缺失情况下不同算法对实测数据成像结果的指标值Table 6 Index values of imaging results of the real data corresponding to four algorithms in the case of random loss of data

4 结束语

针对低信噪比下的ISAR 超分辨率成像问题，提出了一种基于快速分裂Bregman 迭代的ISAR 超分辨率成像算法。所提出的FSBI 算法的主要优势为：充分利用托普利兹矩阵的低位移秩特征和Gohberg-Semencul 表示实现了矩阵的快速求逆，从而既能利用SBI 算法在低信噪比条件下的重构能力，又实现了快速成像。通过对仿真和实测数据进行处理，验证了所提出的FSBI 算法相比于现有常用算法能够在低信噪比的情况下具有更佳的聚焦性能，且运行时间相对较短。