CN113671498A - 基于低秩与双重稀疏矩阵分解的sar射频干扰抑制方法 - Google Patents

Abstract

本发明属于SAR信号处理技术领域，公开了一种基于低秩与双重稀疏矩阵分解的SAR射频干扰抑制方法。在对SAR回波信号时频特征分析的基础上，提出RFI时频矩阵具有低秩与稀疏特性，结合目标回波信号的稀疏性，建立联合低秩与双重稀疏特性约束的SAR回波信号分离模型。采用交替迭代投影策略将该多重约束的信号分离优化问题转化为两个优化求解子问题：干扰低秩重构与信号稀疏恢复。针对干扰的低秩重构问题本采用双边随机投影策略进行低秩估计，并结合硬阈值方法实现RFI与目标回波信号的稀疏求解。本发明对RFI进行了精细化约束，使得模型更加精确，提高了模型对于数据的拟合精度。

Description

基于低秩与双重稀疏矩阵分解的SAR射频干扰抑制方法

技术领域

本发明属于雷达信号处理技术领域，更进一步涉及合成孔径雷达(SAR)信号处理领域中的基于低秩与双重稀疏矩阵分解的射频干扰抑制方法。本发明可用于对SAR回波信号中的射频干扰(RFI)进行抑制，恢复遮盖的目标信息，获得高质量的目标成像效果，显著增强SAR图像解译能力。

背景技术

SAR具有全天时，全天候，远作用距离、大测绘带宽、高分辨成像等的优势，在目标检测与识别、资源勘探与地质测绘等领域都具有广阔的应用前景。然而，随着电磁辐射源的急剧增加与频谱资源的重叠利用，SAR面临的电磁环境日益复杂，并且SAR系统为了实现高分辨成像需要发射大带宽信号，使其不可避免受到同频段RFI的影响。RFI可分为窄带干扰(NBI)和宽带干扰(WBI)。首先，RFI的存在会显著降低SAR回波信号的信噪比，其次，RFI会影响SAR系统成像参数的估计精度，使成像结果散焦。因此，研究有效的RFI抑制算法具有重要的应用价值。

R.T.Lord等人在文献“Radio frequency interference suppression appliedtosynthetic aperture radar data”(General Assembly and Scientific Symposium ofthe International Union ofRadio Science，2005年)提出了基于频域陷波的干扰抑制方法，该方法简单且易于实践。但是该方法会存在潜在的目标信号损失问题。

Tao等人在文献“Wideband Interference Mitigation in High-ResolutionAirborne Synthetic Aperture Radar Data”(IEEE Transactions on Geoscience&Remote Sensing，第54卷第1期，2016年1月)中通过将WBI抑制问题转化为瞬时频谱的NBI抑制问题，提出了基于瞬时频谱特征子空间滤波(ISEF)的干扰抑制方法。然而，该方法一般适用于统计特性缓变或者不变的信号，当面对统计特性随时间变化较大的回波信号时干扰抑制结果将不够理想。

Su等人在文献“Narrow-band interference suppression via RPCA-basedsignal separation in time-frequency domain，”(IEEE Journal of Selected Topicsin Applied Earth Observations and Remote Sensing，第10卷第11期，2017年八月)中基于对SAR回波信号时频谱图低秩与稀疏成分的分析，利用GoDec(Go Decomposition)算法进行矩阵分解实现NBI与目标回波信号的分离。然而，该方法仅利用了干扰与信号的部分特征，模型不够精确，会产生一定的重构误差。

发明内容

针对上述存在的问题，本发明的目的在于提供一种基于低秩与双重稀疏矩阵分解的SAR射频干扰抑制方法，在对SAR回波信号时频分析的基础上，提出RFI时频矩阵低秩且稀疏，以及目标回波信号时频稀疏的假设，利用信号分量的时频域特征差异构建基于矩阵分解的RFI与目标回波信号分离模型，最后采用交替迭代策略实现RFI的精确重构与目标回波信号的低损恢复。本发明对RFI进行了精细化约束，使得模型更加精确，提高了模型对于数据的拟合精度。

本发明的基本思路是：首先对单次含RFI回波信号利用短时傅里叶变换转换至时频域，并将其表征为RFI，目标回波信号以及加性噪声的叠加形式。在对SAR回波信号时频特征分析的基础上，提出RFI时频矩阵具有低秩与稀疏特性，结合目标回波信号的稀疏性，建立联合低秩与双重稀疏特性约束的SAR回波信号分离模型。采用交替迭代投影策略将该多重约束的信号分离优化问题转化为两个优化求解子问题：干扰低秩重构与信号稀疏恢复。针对干扰的低秩重构问题本发明采用双边随机投影策略进行低秩估计，并结合硬阈值方法实现RFI与目标回波信号的稀疏求解，交替迭代直至原信号分离优化问题误差收敛。然后，使用重构RFI进行相消实现目标回波信号的低损恢复。

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

一种基于低秩与双重稀疏矩阵分解的SAR射频干扰抑制方法，包括以下步骤：

步骤1，建立SAR沿方位向的单次含RFI回波模型，利用短时傅里叶变换将单次含RFI回波信号转换至时频域，得到时频域回波信号；

步骤2，基于时频域回波信号的时频特征分析，提出RFI时频矩阵具有低秩与稀疏特性，并结合目标回波信号时频稀疏特征，构建联合低秩与稀疏约束的RFI与目标回波信号分离模型；

步骤3，采用交替迭代投影策略将所述联合低秩与稀疏约束的RFI与目标回波信号分离模型分解为RFI重构与目标回波信号稀疏恢复两个子问题的求解；通过双边随机投影策略结合硬阈值投影实现对RFI的重构；

步骤4，从时频域回波信号中消去重构的RFI，得到时频域目标回波信号的稀疏恢复，对时频域目标回波信号的稀疏恢复进行傅里叶逆变换，得到射频干扰抑制后的目标回波信号。

与现有技术相比，本发明的有益效果为：

本发明通过短时傅里叶变换技术，充分利用SAR回波的时序与频谱特征，并在此基础上，引入RFI的时频低秩与稀疏特征，结合目标回波信号的稀疏分布假设，构建了精确的RFI与目标回波信号分离模型。使用交替迭代投影策略简化RFI与目标回波信号的联合优化问题，并结合双边随机投影方法与硬阈值投影方法实现RFI的精确重构与目标回波信号的低损恢复。本发明进一步提高了RFI抑制模型精度，有效实现RFI抑制，并显著降低SAR目标回波信号重构误差。

附图说明

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

图1为本发明提供的基于低秩与双重稀疏矩阵分解的SAR射频干扰抑制方法的实现流程图；

图2为实测含NBI和WBI的SAR回波数据时频表征结果。其中，(a)和(b)分别为NBI和WBI的时频谱图，(c)和(d)分别为(a)和(b)的特征值分布，(e)和(f)分别为(a)和(b)的幅值分布；

图3为无干扰单次SAR回波数据的时频分析结果。其中，(a)为单次SAR回波时频图，(b)为该次SAR回波的时频矩阵幅值分布结果；

图4为X波段中国某机载SAR实测数据干扰抑制后成像结果对比。其中，(a)为原始SAR图像，(b)为GoDec算法干扰抑制处理后的成像结果，(c)为本发明方法干扰抑制处理后的成像结果；

图5为Sentinel-1B星载SAR的舰船场景实测数据干扰抑制后成像结果对比。其中，(a)为原始SAR图像，(b)为GoDec算法干扰抑制处理后的成像结果，(c)为本发明方法干扰抑制处理后的结果。

具体实施方式

下面将结合实施例对本发明的实施方案进行详细描述，但是本领域的技术人员将会理解，下列实施例仅用于说明本发明，而不应视为限制本发明的范围。

参考图1，本发明提供的一种基于低秩与双重稀疏矩阵分解的SAR射频干扰抑制方法，包括以下步骤：

步骤1，建立SAR沿方位向的单次含RFI回波模型，利用短时傅里叶变换将单次含RFI回波信号转换至时频域，得到时频域回波信号；

SAR沿方位向单次回波可以表示为信号、干扰与噪声的线性叠加形式：

s(k)＝x(k)+i(k)+n(k)

其中，k表示距离快拍序号，s(·)、x(·)、i(·)、n(·)分别表示原始回波信号、目标回波信号、射频干扰以及噪声信号。

由于射频干扰与目标回波信号在时域和频域特征耦合严重，且本发明为了充分利用SAR回波的多维域信息，利用短时傅里叶变换将单次方位回波映射至距离时频域，得时频域回波信号：

STFT_X＝STFT_I+STFT_S+STFT_N

其中，STFT_X、STFT_I、STFT_S与STFT_N分别表示原始回波信号的时频矩阵、射频干扰的时频矩阵、目标回波信号的时频矩阵以及加性噪声的时频矩阵。

步骤2，基于时频域回波信号的时频特征分析，提出RFI时频矩阵具有低秩与稀疏特性，并结合目标回波信号时频稀疏特征，构建联合低秩与稀疏约束的RFI与目标回波信号分离模型；

子步骤2a，分析RFI与目标回波信号时频低秩与稀疏特性；

图2(a)和图2(b)分别是NBI和WBI的时频谱图，横坐标为时间单元，纵坐标为频率单元。由于RFI的能量远高于目标回波信号能量，相应的时频谱图中较亮区域即对应着RFI，可以看到RFI仅占有限的时频单元，且在时频域具有聚集特性，可以认为RFI时频矩阵具有稀疏和低秩特性。图2(c)和图2(d)分别为NBI和WBI的特征值从大到小分布图，横坐标为特征值索引，纵坐标为特征值大小。根据特征值大小与信号分量功率的关联可以发现RFI对应的大特征值占较少部分，而其余信号分量对应的小特征值占大部分，因此说明了RFI时频矩阵相对于回波时频矩阵具有低秩性。图2(e)和图2(f)分别为图2(a)和图2(b)所示时频谱图的幅值分布结果，横坐标为元素下标索引值，纵坐标为幅值，可以看到该时频谱图的幅度分布具有明显的转折，其中RFI对应的大幅值仅占所有元素的较小比例，说明了RFI具有稀疏性。与定性分析结果一致，RFI时频矩阵为低秩且稀疏的矩阵。

图3(a)是无干扰SAR回波的时频谱图，横坐标为时间单元，纵坐标为频率单元。目标信号相对比于噪声具有更强的回波功率，因此无干扰SAR回波时频谱图中较大的幅值对应着目标回波信号，可以看到目标回波信号的时频分布异于噪声，其占据的时频单元数量相对较少。图3(b)是图3(a)所示时频谱图的幅值分布图，可以看到其幅值分布曲线具有明显的转折，且目标回波信号对应的大幅值仅占较少比例，说明目标回波信号具有稀疏性，其时频矩阵为稀疏矩阵。

子步骤2b，构建联合低秩与双重稀疏约束的RFI与目标回波信号分离模型；

不失一般性，假设噪声服从复高斯分布。高斯噪声作为RFI与目标回波信号的重构误差，本发明给出基于L2范数的最小化误差优化函数：

其中，||·||_F表示求L2范数。

由于STFT_I为低秩、稀疏矩阵，STFT_S为稀疏矩阵的性质，对上述最小误差优化函数加入约束条件，得到联合低秩与双重稀疏约束的RFI与目标回波信号分离模型：

其中，rank(STFT_I)与card(STFT_I)分别为RFI时频矩阵的秩与基数，r与k₁分别为RFI时频矩阵的秩与基数的阈值；card(STFT_S)为目标回波信号时频矩阵的基数，k₂为card(STFT_S)的阈值。r、k₁、k₂为预设常数，在具体操作过程中根据结果进行相应调整。

步骤3，采用交替迭代投影策略将所述联合低秩与稀疏约束的RFI与目标回波信号分离模型分解为RFI重构与目标回波信号稀疏恢复两个子问题的求解；通过双边随机投影策略结合硬阈值投影实现对RFI的重构；

交替迭代求解联合低秩与稀疏约束的RFI与目标回波信号分离模型，即求解RFI与目标回波信号的联合优化问题：

将步骤2b中的优化问题分解为RFI重构与目标回波稀疏恢复的两个子问题，并进行交替迭代更新。其中第t次迭代更新公式为：

其中，

表示第t+1次迭代更新的STFT_I结果，

表示第t+1次迭代更新的STFT_S结果。

通常，对于每次迭代过程中RFI时频矩阵的重构可以通过奇异值阈值方法结合逐元素的硬阈值投影求解：

其中，

表示子集，可以看出，

先根据

的前r个奇异向量更新，再根据其中数值从大到小排列前k₁个元素组成的非零子集更新，而

根据

中数值从大到小排列前k₂个元素组成的非零子集更新。svd表示奇异值分解操作；Λ是奇异值矩阵，为对角阵，其对角线上的元素为λ_i′，λ_i′是第i′个奇异值，U、V为酉矩阵；

按照奇异值阈值方法，每一步迭代都需要进行一次全矩阵的奇异值分解处理，而奇异值分解运算量与矩阵规模相关，当矩阵较大时采用奇异值阈值方法会消耗大量的时间。

针对奇异值分解的计算耗时问题，本发明采用双边随机投影方法进行每次迭代过程中的RFI重构，双边随机投影算法可以在保证估计精度的前提下显著降低每次迭代的运算量。其定义如下，对于大小为m×n的稠密矩阵X，其r阶近似矩阵为：

其中Y₁＝XA₁，Y₂＝X^TA₂，A₁和A₂分别是n×r和m×r阶的随机矩阵。

因此，射频干扰时频矩阵的重构问题求解过程可以更新为：

其中，L′是X₁的近似矩阵，

表示近似矩阵L′中满足Ω₁的元素组成的子集，

表示Ω₁的补集；

表示近似矩阵L′中满足Ω₁的元素，|Ω₁|≤k₁表示Ω₁的长度不大于k₁，即表示L′中前k₁个最大的元素；

同理，目标回波信号时频矩阵的稀疏求解为：

其中，L″是X₂的近似矩阵，

表示近似矩阵L″中满足Ω₂的元素组成的子集，

表示Ω₂的补集；

表示近似矩阵L″中满足Ω₂的元素，|Ω₂|≤k₂表示Ω₂的长度不大于k₂；

按照上述射频干扰重构与目标回波信号稀疏恢复的求解过程，进行交替迭代估计，直至RFI重构结果收敛即可获得RFI的重构结果。

步骤4，从时频域回波信号中消去重构的RFI，得到时频域目标回波信号的稀疏恢复，对时频域目标回波信号的稀疏恢复进行傅里叶逆变换，得到射频干扰抑制后的目标回波信号。

为了减少重构误差带来的信号损失，本发明采用相消策略最大化保留目标信号：

其中，

表示射频干扰抑制后的目标回波信号，

表示干扰时频矩阵的最终估计结果，ISTFT表示逆STFT变换。

仿真实验

下面通过仿真实验对本发明的效果做进一步说明。

使用本发明方法对SAR实测回波数据进行处理，获得SAR成像结果，并进行定性与定量的分析，验证本发明对RFI抑制的有效性。

(1)对干扰抑制后的SAR回波数据进行成像处理；

按照本发明方法对SAR回波数据进行逐脉冲的干扰抑制处理，并经过SAR成像处理即可获得无干扰且聚焦效果良好的SAR图像。如附图4、附图5所示。

下面结合附图4和附图5对干扰抑制效果进行定性分析。

图4(a)为受到NBI影响且未进行干扰抑制处理的原始SAR数据直接成像结果，可以看出，成像结果中存在由于NBI形成的沿距离向亮线，其遮挡目标建筑物和地面场景，严重影响SAR图像解译；图4(b)为利用GoDec算法进行干扰抑制处理后的SAR成像结果，可以看出，NBI得到有效抑制，干扰遮盖下的目标场景得到恢复，但是SAR图像在细节处由于信号损失产生了一定程度的图像散焦，如矩形框中图像所示，其给图像解译带来了挑战；图4(c)为利用本发明所提算法进行干扰抑制处理后的SAR成像结果，可以看处，其在干扰抑制的同时信号损失更少，图像聚焦度更优。

图5(a)为受到WBI影响且未进行干扰抑制处理的原始SAR数据直接成像结果，可以看出，成像结果中存在大片由于WBI形成的亮线，该图像中舰船目标几乎全部被WBI覆盖，严重影响SAR目标识别；图5(b)为利用GoDec算法进行干扰抑制处理后的SAR成像结果，可以看出，WBI得到有效抑制，舰船目标场景得到恢复，但由于信号损失，舰船目标模糊；图5(c)为利用本发明所提算法进行干扰抑制处理后的SAR成像结果，可以看出，其在干扰抑制的同时舰船目标聚焦效果良好，图像相对更加清晰。

(2)对干扰效果进行定量评估。

为了对本发明所提算法的干扰抑制效果进行定量评估，采用具有代表性的平均梯度与梯度差作为评估指标。

平均梯度(Average Gradient，AG)代表着图像的锐利度，表征了图像细节与纹理的丰富度，数学定义为

其中

知

分别表示图像(m，n)处的灰度值在垂直方向和水平方向上的梯度。AG值越大，代表着SAR图像场景的边缘特征越清晰。I表示SAR成像结果图，其大小为M×N。

灰度差(Gray level difference，GLD)用来表示边缘的灰度变化，具体数学表征为

GLD值越大，表示SAR图像目标细节信息越丰富。

表1和表2分别是对两组实测数据SAR图像的干扰抑制效果的评估结果。

表1窄带干扰抑制效果评估

表2宽带干扰抑制效果评估

上述结果表明：GoDec算法和本发明所提算法皆可实现射频干扰的有效抑制，但是本发明所提方法在NBI与WBI实测数据上相较于GoDec算法具有更大的AG与GLD指标，表明采用本发明所提方法进行干扰抑制处理后的SAR成像结果对比度更高，目标场景轮廓与舰船目标更为清晰，与定性分析结果一致，进一步说明本发明的射频干扰抑制效果较好。

虽然，本说明书中已经用一般性说明及具体实施方案对本发明作了详尽的描述，但在本发明基础上，可以对之作一些修改或改进，这对本领域技术人员而言是显而易见的。因此，在不偏离本发明精神的基础上所做的这些修改或改进，均属于本发明要求保护的范围。

Claims (6)

1.一种基于低秩与双重稀疏矩阵分解的SAR射频干扰抑制方法，其特征在于，包括以下步骤：

步骤1，建立SAR沿方位向的单次含RFI回波模型，利用短时傅里叶变换将单次含RFI回波信号转换至时频域，得到时频域回波信号；

步骤2，基于时频域回波信号的时频特征分析，提出RFI时频矩阵具有低秩与稀疏特性，并结合目标回波信号时频稀疏特征，构建联合低秩与稀疏约束的RFI与目标回波信号分离模型；

步骤3，采用交替迭代投影策略将所述联合低秩与稀疏约束的RFI与目标回波信号分离模型分解为RFI重构与目标回波信号稀疏恢复两个子问题的求解；通过双边随机投影策略结合硬阈值投影实现对RFI的重构；

步骤4，从时频域回波信号中消去重构的RFI，得到时频域目标回波信号的稀疏恢复，对时频域目标回波信号的稀疏恢复进行傅里叶逆变换，得到射频干扰抑制后的目标回波信号。

2.根据权利要求1所述的基于低秩与双重稀疏矩阵分解的SAR射频干扰抑制方法，其特征在于，步骤1中，所述单次含RFI回波模型为：

s(k)＝x(k)+i(k)+n(k)

其中，k表示距离快拍序号，s(·)、x(·)、i(.)、n(.)分别表示原始回波信号、目标回波信号、射频干扰以及噪声信号；

所述时频域回波信号的表达式为：

STFT_X＝STFT_I+STFT_S+STFT_N

其中，STFT_X、STFT_I、STFT_S与STFT_N分别表示原始回波信号的时频矩阵、射频干扰的时频矩阵、目标回波信号的时频矩阵以及加性噪声的时频矩阵。

3.根据权利要求2所述的基于低秩与双重稀疏矩阵分解的SAR射频干扰抑制方法，其特征在于，所述构建联合低秩与稀疏约束的RFI与目标回波信号分离模型，具体为：

首先，高斯噪声作为RFI与目标回波信号的重构误差，设噪声服从复高斯分布，则给出基于L2范数的最小化误差优化函数：

其中，||·||_F表示求L2范数；

由于STFT_I为低秩、稀疏矩阵，STFT_S为稀疏矩阵，对上述最小误差优化函数加入约束条件，得到联合低秩与双重稀疏约束的RFI与目标回波信号分离模型：

其中，rank(STFT_I)与card(STFT_I)分别为RFI时频矩阵的秩与基数，r与k₁分别为RFI时频矩阵的秩与基数的阈值；card(STFT_S)为目标回波信号时频矩阵的基数，k₂为card(STFT_S)的阈值。

4.根据权利要求3所述的基于低秩与双重稀疏矩阵分解的SAR射频干扰抑制方法，其特征在于，所述采用交替迭代投影策略将所述联合低秩与稀疏约束的RFI与目标回波信号分离模型分解为RFI重构与目标回波信号稀疏恢复两个子问题的求解，具体为：

将分离模型分解为RFI重构与目标回波稀疏恢复两个子问题，并进行交替迭代更新，其中第t次迭代更新公式为：

其中，

表示第t+1次迭代更新的STFT_I结果，

表示第t+1次迭代更新的STFT_S结果；

每次迭代过程中的RFI重构问题通过双边随机投影算法进行求解；

按照射频干扰重构与目标回波信号稀疏恢复的求解过程，进行交替迭代估计，直至RFI重构结果收敛即可获得RFI的重构结果。

5.根据权利要求4所述的基于低秩与双重稀疏矩阵分解的SAR射频干扰抑制方法，其特征在于，所述双边随机投影算法的求解过程为：

首先，定义对于大小为m×n的稠密矩阵X，其r阶近似矩阵为：

其中，Y₁＝XA₁，Y₂＝X^TA₂，A₁和A₂分别是n×r和m×r阶的随机矩阵；

据此，将射频干扰时频矩阵的重构问题求解过程更新为：

其中，L′是X₁的近似矩阵，

表示近似矩阵L′中满足Ω₁的元素组成的子集，

表示Ω₁的补集；

表示近似矩阵L′中满足Ω₁的元素，|Ω₁|≤k₁表示Ω₁的长度不大于k₁，即表示L′中前k₁个最大的元素；|·|表示取模；上标T表示矩阵转置，

表示

的奇异值分解的前r个奇异向量，λ_i′表示第i′个奇异值；

同理，目标回波信号时频矩阵的稀疏求解为：

其中，L″是X₂的近似矩阵，

表示近似矩阵L″中满足Ω₂的元素组成的子集，

表示Ω₂的补集；

表示近似矩阵L″中满足Ω₂的元素，|Ω₂|≤k₂表示Ω₂的长度不大于k₂。

6.根据权利要求1所述的基于低秩与双重稀疏矩阵分解的SAR射频干扰抑制方法，其特征在于，所述射频干扰抑制后的目标回波信号为：

其中，

表示干扰时频矩阵的重构结果，ISTFT表示逆STFT变换。

Images