CN107576961B - 一种互质降采样间歇合成孔径雷达稀疏成像方法 - Google Patents

Abstract

本发明公开了一种互质降采样间歇合成孔径雷达稀疏成像方法，它是结合互质降采样方式及压缩感知稀疏重构方法，先采用互质的两个整数对间歇合成孔径雷达中方位向脉冲数进行降采样间歇数据获取，实现合成孔径雷达的脉冲互质降采样间歇，然后分别对每一个质数得到合成孔径雷达数据进行压缩感知稀疏成像，得到两幅合成孔径雷达图像后，再进行栅瓣融合处理，得到最终的间歇合成孔径雷达图像，从而实现互质降采样间歇合成孔径雷达高精度成像。本发明具有相对于传统随机间歇方式更便于系统实现、并且相对于传统插值重采样成像算法提高间歇SAR稀疏成像质量、抑制图像中旁瓣水平等特点。

Description

一种互质降采样间歇合成孔径雷达稀疏成像方法

技术领域

本发明属于雷达技术领域，它特别涉及合成孔径雷达(SAR)成像技术领域。

背景技术

合成孔径雷达(SAR)作为一种全天时、全天候、信息量丰富的遥感成像技术，已成为对地观测的重要手段，在地形图像生成、目标探测与侦察、国土资源勘查和自然灾害监测等国民经济与军事领域发挥重要作用。近几年来，随着雷达硬件技术的快速发展及合成孔径雷达系统的任务多样化，实现多模式成像、具备多功能能力已成为现代合成孔径雷达系统设计的迫切需求。间歇合成孔径雷达(Interrupted SAR)是目前实现多模式、多功能合成孔径雷达系统的热点技术之一，其主要采用发射脉冲间隙工作方式，在雷达工作过程中部分时段的脉冲用于实现雷达测角、探测及跟踪等雷达功能，然后利用间歇技术使得部分时段的脉冲用于实现合成孔径雷达成像观测能力，详见文献“Stojanovic I,Novak L,Karl WC.Interrupted SAR persistent surveillance via group sparse reconstruction ofmultipass data[J].IEEE Transactions on Aerospace and Electronic Systems,2014,50(2):987-100”。因此，相对传统合成孔径雷达系统，间隙合成孔径雷达可具备目前多种雷达体制的功能，具有更重要、广泛的应用价值。目前，如何设计收发脉冲间歇方式及降采样条件下成像方法成为间歇合成孔径雷达应用面临的两个关键问题。

目前，间隙合成孔径雷达的脉冲间歇方式主要有凹縫间歇、周期性间隙、随机间歇等三种，凹縫间歇是指在雷达观测孔径中有一段连续的脉冲用于实现合成孔径雷达成像之外的其它雷达功能，周期性间隙指在雷达观测孔径中脉冲数被周期性用于实现合成孔径雷达成像之外的其它雷达功能，随机性间隙指在雷达观测孔径中脉冲数被随机性用于实现合成孔径雷达成像之外的其它雷达功能，详见文献“Tough R J A,Ward K D,Shepherd PW.Modelling and enhancement of interrupted synthetic aperture radar imagery[J].IET Radar,Sonar&Navigation,2007,1(6):407-413”。相对随机间歇，凹縫间歇和周期性间隙更便于雷达硬件设计实现。为了保证降采样条件下间隙合成孔径雷达的成像性能，目前主要采用插值重采样方法对丢失的间歇数据进行补偿，再通过传统合成孔径雷达成像算法进行成像处理，详见文献“Tough R J A,Ward K D,Shepherd P W,et al.TheEnhancement of Interrupted SAR Imagery[C],3rd EMRS DTC Technical Conference,2006”。但是，当间歇脉冲数过多时，插值及重采样方法难以精确进行数据拟合重建，导致间歇合成孔径雷达难以保持原来的成像质量。

近几年，压缩感知稀疏重构作为一种新提出的信号处理理论，突破了传统Nyquist采样定理约束，可利用远低于Nyquist采样率精确重构原始稀疏信号，在降采样条件下提高雷达性能等方面有着巨大的应用潜力，详见参考文献“D.L.Donoho.Compressedsensing.IEEE Transactions on Information Theory,2006,52(4):1289-1306”。目前，多种基于压缩感知稀疏重构理论的稀疏成像算法已经被应用于合成孔径雷达稀疏成像，详见参考文献“Candes E,Romberg J.Sparsity and incoherence in compressive sampling[J].Inverse problems,2007,23(3):969”。对于脉冲数据丢失的间歇合成孔径雷达系统，压缩感知稀疏重构可通过建立线性测量模型及L1范数最小化重构，从而可抑制间歇合成孔径雷达脉冲数据丢失的旁瓣栅瓣影响，详见文献“Vu D,Xu L,Xue M,et al.Nonparametricmissing sample spectral analysis and its applications to interrupted SAR[J].IEEE Journal of Selected Topics in Signal Processing,2012,6(1):1-14”。另外，针对间歇合成孔径雷达的不规则脉冲数据重构，也有学者提出一种基于稀疏正则化的稀疏重构方法，详见文献“Storm K,Murthy V,Selesnick I,et al.Sparsity-based methods forinterrupted radar data reconstruction[C].IEEE Radar Conference,2012,0107-0111”，该方法在离散傅立叶变换域采用最小二乘正则化方法可较好提高间歇合成孔径雷达成像质量。相对传统间歇合成孔径雷达插值重采样方法，压缩感知稀疏重构算法虽然可抑制旁瓣和提高间歇合成孔径雷达图像质量，但是压缩感知稀疏重构算法主要适用于凹縫间歇和随机间歇方式，而对于周期性间隙降采样的数据成像还会存在严重的栅瓣影响。

互质分布设计是阵列天线阵元布设中降低栅瓣水平的一种有效方法，其主要利用两个互质的倍数对天线阵元间距进行布设，然后通过适当加权处理即可实现稀疏阵列天线高性能的波束形成和波达角估计，详见文献“Qin S,Zhang Y D,Amin M G.Generalizedcoprime array configurations for direction-of-arrival estimation[J].IEEETransactions on Signal Processing,2015,63(6):1377-1390”。因具有抑制成像栅瓣水平的能力并且互质作为一种便于实现的周期降采样方法，互质脉冲降采样也被视为一种新的合成孔径雷达数据获取模式，应用合成孔径雷达海面舰船成像，详见文献“Di MartinoG,Iodice A.Coprime Synthetic Aperture Radar(CopSAR):A New Acquisition Modefor Maritime Surveillance[J].IEEE Transactions on Geoscience and RemoteSensing,2015,53(6):3110-3123”。但是，该方法中互质降采样合成孔径雷达还是局限于传统成像方法，成像图像中强目标高旁瓣仍难以抑制，影响合成孔径雷达图像的质量。

根据压缩感知稀疏重构和互质降采样的特性及优势，可以将两者结合应用于间歇合成孔径雷达系统，利用互质降采样实现间歇合成孔径雷达的脉冲间歇方式及降低合成孔径雷达成像的脉冲数，并且结合压缩感知稀疏重构方法提高互质降采样条件下合成孔径雷达的成像质量，进一步提高间歇合成孔径雷达系统的成像性能。

发明内容

本发明提出了一种互质降采样间歇合成孔径雷达稀疏成像方法，该方法结合互质降采样方式及压缩感知稀疏重构方法，先采用互质的两个整数对间歇合成孔径雷达中方位向脉冲数进行降采样间歇数据获取，实现合成孔径雷达的脉冲互质降采样间歇，然后分别对每一个质数得到合成孔径雷达数据进行压缩感知稀疏成像，得到两幅合成孔径雷达图像后，再进行栅瓣融合处理，得到最终的间歇合成孔径雷达图像，从而实现互质降采样间歇合成孔径雷达高精度成像。本发明具有相对于传统随机间歇方式更便于系统实现、并且相对于传统插值重采样成像算法提高间歇SAR稀疏成像质量、抑制图像中旁瓣水平等特点。

为了方便描述本发明的内容，首先作以下术语定义：

定义1、合成孔径雷达(SAR)

合成孔径雷达(SAR)是利用雷达天线与目标的相对运动合成一个大的虚拟天线孔径，以实现雷达方位向的高分辨成像能力，再通过发射宽带信号及脉冲压缩技术得到雷达距离向的高分辨成像能力，从而最终实现观测目标的二维高精度成像的一种雷达成像技术，详见文献“合成孔径雷达成像原理”，皮亦鸣等编著，电子科技大学出版社出版。

定义2、合成孔径雷达慢时刻和快时刻

合成孔径雷达慢时间是指雷达平台飞过一个合成孔径所需要的时间。雷达系统以一定的重复周期发射接收脉冲，因此慢时间可以表示为一个以重复周期为步长的离散化时间变量，其中每一个离散时间变量值为一个慢时刻；合成孔径雷达快时间是指雷达发射接收脉冲的一个周期的时间。由于雷达接收回波是以采样率进行采样，则快时刻可以表示为一个离散化的时间变量，每一个离散变量值为一个快时刻；详见文献“保铮,刑孟道,王彤.雷达成像技术[M].电子工业出版社,2005”。

定义3、合成孔径雷达成像空间

合成孔径雷达投影成像空间是指在合成孔径雷达数据成像时选取的成像空间，合成孔径雷达成像需要将回波数据投影到该成像空间进行聚焦处理。一般的，合成孔径雷达成像成像空间选择为斜距平面坐标系或者水平地面坐标系。

定义4、标准合成孔径雷达原始回波仿真方法

标准合成孔径雷达原始回波仿真方法是指给定雷达系统参数、平台轨迹参数以及观测场景参数等所需的参数条件下，基于合成孔径雷达成像原理仿真得到具有合成孔径雷达回波信号特性的原始回波信号的方法，详细内容可参考文献：“郭晓宇.高分辨SAR射频目标仿真技术研究[D].南京理工大学,2008”。

定义5、标准压缩感知合成孔径雷达测量矩阵构造方法

压缩感知合成孔径雷达测量矩阵构造方法是指给定雷达系统参数、平台轨迹参数以及成像空间参数等所需的参数条件下，利用合成孔径雷达系统与成像空间中网格点的时延相位构造压缩感知测量矩阵，该测量矩阵与成像空间中目标真实散射系数向量相乘即可得到合成孔径雷达回波数据，详见文献“Wei S J,Zhang X L,Shi J,et al.Sparsereconstruction for SAR imaging based on compressed sensing[J].Progress InElectromagnetics Research,2010,109:63-81”。

定义5、标准压缩感知稀疏重构方法

标准压缩感知稀疏重构方法是指压缩感知理论中用于稀疏数据重构的经典算法，主要包含了基追踪、正交基追踪及贪婪追踪等稀疏重构算法，详见文献“Donoho DL.Compressed sensing[J].IEEE Transactions on information theory,2006,52(4):1289-1306”。

本发明提供的一种互质降采样间歇合成孔径雷达稀疏成像方法，它包括以下几个步骤：

步骤1、初始化SAR雷达系统参数和成像空间：

初始化SAR雷达系统参数，包括：雷达载波波长，记为λ，雷达中心频率，记为f₀，雷达发射信号带宽，记为B，雷达发射脉冲时宽，记为T_r，雷达采样频率，记为F_s，雷达波束中心入射角，记为θ，雷达脉冲重复频率，记为PRF，平台运动速度矢量，记为V_r，雷达系统距离向采样点数，记为N_r，雷达系统方位向采样点数，记为N_a，雷达系统天线初始位置，记为P_S(0)；上述参数中，光在空气中的传播速度c为常数，雷达载波波长λ、雷达中心频率f₀、雷达发射信号带宽B、雷达发射脉冲时宽T_r、雷达采样频率F_s、雷达波束中心入射角θ、雷达脉冲重复频率PRF、平台运动速度矢量V_r、雷达系统距离向采样点数N_r、雷达系统方位向采样点数N_a、雷达系统天线初始位置P_S(0)在SAR雷达成像观测方案设计中已经确定；

初始化SAR成像空间参数，包括：以SAR雷达波束照射场区域的水平地面坐标系作为SAR成像空间，成像空间中心坐标位置位于[0,0]，成像空间的横轴和纵轴分别记为X轴和Y轴，成像空间X轴和Y轴范围大小分别记为W_X和W_Y，将成像空间均匀划网格得到X轴和Y轴的网格点总数分别记为N_X和N_Y，成像空间X轴第i个和Y轴第j个网格点的位置，记为P_T(i,j)，i＝1,2,…,N_X，j＝1,2,…,N_Y，P_T(i,j)的值表示为其中i和j为自然数，i表示成像空间X轴第i个网格点，j表示成像空间Y轴第j个网格点；

步骤2、获取SAR雷达系统的原始回波数据和天线位置：

SAR雷达系统在距离向第t个快时刻方位向第k个慢时刻的原始回波数据，记为S(t,k)，t＝1,2,…,N_r，k＝1,2,…,N_a，其中t和k为自然数，t表示距离向第t个快时刻，k表示方位向第k个慢时刻，N_r为步骤1初始化得到的SAR雷达系统距离向采样点数，N_a为步骤1初始化得到的SAR雷达系统方位向采样点数；在SAR系统实际成像中，原始回波数据S(t,k)，t＝1,2,…,N_r，k＝1,2,…,N_a，由SAR雷达系统数据接收机提供；在仿真成像中，原始回波数据S(t,k)，t＝1,2,…,N_r，k＝1,2,…,N_a，由标准合成孔径雷达原始回波仿真方法得到；

采用公式计算得到方位向第k个慢时刻的雷达系统天线位置，记为P_S(k)，k＝1,2,…,N_a，其中P_S(0)为步骤1初始化得到的雷达系统天线初始位置，V_r为步骤1初始化得到的平台运动速度矢量，PRF为步骤1初始化得到的雷达脉冲重复频率；

步骤3、对SAR方位向脉冲进行互质降采样的数据抽取：

初始化互质的两个整数，分别记为M₁和M₂；

对于步骤2中得到的原始回波数据S(t,k)，t＝1,2,…,N_r，k＝1,2,…,N_a，按照质数M₁的倍数进行方向数据抽取，即方位向慢时刻k的值被M₁整除时提取出对应的方位向慢时刻k、距离向第t个快时刻方位向第k个慢时刻原始回波数据S(t,k)和雷达系统天线位置P_S(k)，得到质数M₁降采样后的SAR雷达系统方位向采样点总数，记为N₁，得到质数M₁降采样后的距离向第t个快时刻方位向第m₁个慢时刻的SAR原始回波数据，记为S₁(t,m₁)，t＝1,2,…,N_r，m₁＝1,2,…,N₁，得到质数M₁降采样后方位向第m₁个慢时刻的雷达系统天线位置，记为P₁(m₁)，m₁＝1,2,…,N₁，其中m₁为自然数，m₁表示质数M₁降采样后方位向第m₁个慢时刻，N₁的值为 表示取整运算符号；

对于步骤2中得到的原始回波数据S(t,k)，t＝1,2,…,N_r，k＝1,2,…,N_a，按照质数M₂的倍数进行方向数据抽取，即方位向慢时刻k的值被M₂整除时提取出对应的方位向慢时刻k、距离向第t个快时刻方位向第k个慢时刻原始回波数据S(t,k)和雷达系统天线位置P_S(k)，得到质数M₂降采样后的SAR雷达系统方位向采样点总数，记为N₂，得到质数M₂降采样后的距离向第t个快时刻方位向第m₂个慢时刻的SAR原始回波数据，记为S₂(t,m₂)，t＝1,2,…,N_r，m₂＝1,2,…,N₂，得到质数M₂降采样后方位向第m₂个慢时刻的雷达系统天线位置，记为P₂(m₂)，m₂＝1,2,…,N₂，其中m₂为自然数，m₂表示质数M₂降采样后方位向第m₂个慢时刻，N₂的值为

步骤4、对互质降采样后SAR数据分别建立线性测量模型：

将步骤3中得到的质数M₁降采样后的SAR原始回波数据S₁(t,m₁)，t＝1,2,…,N_r，m₁＝1,2,…,N₁，按先行后列的顺序组成一个原始回波数据列向量，记为y₁，即向量y₁的值表示为y₁＝[S₁(1,1),S₁(2,1),…,S₁(N_r,N₁)]^T，其中S₁(1,1)表示质数M₁降采样后距离向第1个快时刻方位向第1个慢时刻的SAR原始回波数据，S₁(2,1)表示质数M₁降采样后距离向第2个快时刻方位向第1个慢时刻的SAR原始回波数据，S₁(N_r,N₁)表示质数M₁降采样后距离向第N_r个快时刻方位向第N₁个慢时刻的SAR原始回波数据，上标T表示转置运算符号，N₁为步骤3中得到的质数M₁降采样后SAR雷达系统方位向采样点总数；

将步骤3中得到的质数M₂降采样后的SAR原始回波数据S₂(t,m₂)，t＝1,2,…,N_r，m₂＝1,2,…,N₂，按先行后列的顺序组成一个原始回波数据列向量，记为y₂，即向量y₂的值表示为y₂＝[S₂(1,1),S₂(2,1),…,S₂(N_r,N₂)]^T，其中S₂(1,1)表示质数M₂降采样后距离向第1个快时刻方位向第1个慢时刻的SAR原始回波数据，S₂(2,1)表示质数M₂降采样后距离向第2个快时刻方位向第1个慢时刻的SAR原始回波数据，S₂(N_r,N₂)表示质数M₂降采样后距离向第N_r个快时刻方位向第N₂个慢时刻的SAR原始回波数据，N₂为步骤3中得到的质数M₂降采样后SAR雷达系统方位向采样点总数；

利用步骤1初始化得到的SAR雷达系统参数、成像空间参数和步骤3得到的质数M₁降采样后方位向第m₁个慢时刻的雷达系统天线位置P₁(m₁)，m₁＝1,2,…,N₁，采用标准压缩感知SAR测量矩阵构造方法得到质数M₁降采样后SAR原始回波的测量矩阵，记为A₁，矩阵A₁的维数大小为N_rM₁×N_XN_Y，其中N_r为步骤1初始化得到的SAR雷达系统距离向采样点数，N_X和N_Y分别为步骤1初始化得到的成像空间X轴和Y轴的网格点总数；

利用步骤1初始化得到的SAR雷达系统参数、成像空间参数和步骤3得到的质数M₂降采样后方位向第m₂个慢时刻的雷达系统天线位置P₂(m₂)，m₂＝1,2,…,N₂，采用标准压缩感知SAR测量矩阵构造方法得到质数M₂降采样后SAR原始回波的测量矩阵，记为A₂，矩阵A₂的维数大小为N_rM₂×N_XN_Y；

步骤5、采用压缩感知稀疏重构方法进行稀疏成像：

采用标准压缩感知稀疏重构方法对最优化问题方程进行求解，得到质数M₁降采样后SAR回波数据稀疏重构成像的列向量，记为向量的维数为N_XN_Y×1，向量中第w个元素记为w＝1,2,…,N_XN_Y，其中w为自然数，w表示为一个向量中第w个元素，y₁为步骤4中得到的质数M₁降采样后回波数据组成的列向量数据，A₁为步骤4中得到的质数M₁降采样后SAR原始回波的测量矩阵，表示求取满足条件的自变量x的最小值，||·||₁为向量L1范数，||·||₂为向量L2范数；

采用标准压缩感知稀疏重构方法对最优化问题方程进行求解，得到质数M₂降采样后SAR回波数据稀疏重构成像的列向量，记为向量的维数为N_XN_Y×1，向量中第w个元素记为w＝1,2,…,N_XN_Y，其中y₂为步骤4中得到的质数M₂降采样后回波数据组成的列向量数据，A₂为步骤4中得到的质数M₂降采样后SAR原始回波的测量矩阵；

步骤6、对互质降采样得到稀疏成像结果进行栅瓣抑制融合：

初始化一个维数为N_XN_Y×1的向量，记为其中向量第w个元素记为w＝1,2,…,N_XN_Y；将步骤5得到的向量和向量逐个元素进行比较，若w＝1,2,…,N_XN_Y，则向量中第w个元素赋值为否则向量中第w个元素赋值为最后得到的向量即为互质降采样栅瓣抑制融合后的稀疏重构向量；

步骤7、重排得到最终的稀疏成像结果：

将步骤6得到的稀疏重构向量按照先行后列的顺序组成一个维数为N_X×N_Y的二维矩阵，记为I，I即为互质降采样间歇SAR稀疏成像的最终图像，其中N_X和N_Y分别为步骤1初始化得到的成像空间X轴和Y轴的网格点总数。

本发明的创新点结合互质降采样特性及压感感知稀疏重构方法，将互质降采样特性及压感感知稀疏重构应用于间歇合成孔径雷达成像，提出了一种互质降采样间歇合成孔径雷达稀疏成像方法，该方法采用互质的两个整数对间歇SAR方位向实现脉冲间歇及数据降采样提取，再利用压缩感知稀疏重构方法分别进行稀疏成像，然后通过栅瓣抑制融合处理，从而实现了互质降采样间歇SAR的高精度稀疏成像。

本发明的优点是结合互质降采样方法及压感感知稀疏重构方法优势，利用互质降采样不仅可以降低间歇SAR方位向的采样数据，相对于传统随机间歇方式更便于系统实现，并且相对于传统插值重采样成像算法提高间歇SAR稀疏成像质量，抑制图像中旁瓣水平。

附图说明

图1为发明所提供方法的流程示意框图；

图2为合成孔径雷达系统及成像的各仿真参数值。

具体实施方式

本发明主要采用仿真实验的方法进行验证，所有步骤、结论都在MATLABR2014b软件上验证正确。具体实施步骤如下：

步骤1、初始化SAR雷达系统参数和成像空间：

初始化SAR雷达系统参数，包括：雷达载波波长λ＝0.03m，雷达中心频率f₀＝10¹⁰Hz，雷达发射信号带宽B＝5×10⁸Hz，雷达发射脉冲时宽T_r＝10-⁵s，雷达采样频率F_s＝6×10⁸Hz，雷达波束中心入射角θ＝45°，雷达脉冲重复频率PRF＝1000Hz，平台运动速度矢量V_r＝[0,100,0]m/s，雷达系统距离向采样点数N_r＝2048，雷达系统方位向采样点数N_a＝2048，雷达系统天线初始位置P_S(0)＝[-100,0,6000]m；

初始化SAR成像空间参数，包括：以SAR雷达波束照射场区域的水平地面坐标系作为SAR成像空间，成像空间中心坐标位置位于[0,0]m，成像空间的横轴和纵轴分别记为X轴和Y轴，成像空间X轴和Y轴范围大小分别W_X＝100m和W_Y＝100m，将成像空间均匀划网格得到X轴和Y轴的网格点总数分别为N_X＝500和N_Y＝500，成像空间X轴第i个和Y轴第j个网格点的位置，记为P_T(i,j)，i＝1,2,…,N_X，j＝1,2,…,N_Y，P_T(i,j)的值表示为其中i和j为自然数，i表示成像空间X轴第i个网格点，j表示成像空间Y轴第j个网格点；

步骤2、获取SAR雷达系统的原始回波数据和天线位置：

SAR雷达系统在距离向第t个快时刻方位向第k个慢时刻的原始回波数据，记为S(t,k)，t＝1,2,…,N_r，k＝1,2,…,N_a，其中t和k为自然数，t表示距离向第t个快时刻，k表示方位向第k个慢时刻，N_r为步骤1初始化得到的SAR雷达系统距离向采样点数N_r＝2048，N_a为步骤1初始化得到的SAR雷达系统方位向采样点数N_a＝2048；在仿真成像中，原始回波数据S(t,k)，t＝1,2,…,N_r，k＝1,2,…,N_a，由标准合成孔径雷达原始回波仿真方法得到；

采用公式k＝1,2,…,N_a，计算得到方位向第k个慢时刻的雷达系统天线位置，记为P_S(k)，k＝1,2,…,N_a，其中P_S(0)为步骤1初始化得到的雷达系统天线初始位置P_S(0)＝[-100,0,6000]m，V_r为步骤1初始化得到的平台运动速度矢量V_r＝[0,100,0]m/s，PRF为步骤1初始化得到的雷达脉冲重复频率PRF＝1000Hz；

步骤3、对SAR方位向脉冲进行互质降采样的数据抽取：

初始化互质的两个整数M₁＝4和M₂＝5；

对于步骤2中得到的原始回波数据S(t,k)，t＝1,2,…,N_r，k＝1,2,…,N_a，按照质数M₁＝4的倍数进行方向数据抽取，即方位向慢时刻k的值被M₁＝4整除时提取出对应的方位向慢时刻k、距离向第t个快时刻方位向第k个慢时刻原始回波数据S(t,k)和雷达系统天线位置P_S(k)，得到质数M₁降采样后的SAR雷达系统方位向采样点总数，记为N₁，得到质数M₁降采样后的距离向第t个快时刻方位向第m₁个慢时刻的SAR原始回波数据，记为S₁(t,m₁)，t＝1,2,…,N_r，m₁＝1,2,…,N₁，得到质数M₁降采样后方位向第m₁个慢时刻的雷达系统天线位置，记为P₁(m₁)，m₁＝1,2,…,N₁，其中m₁为自然数，m₁表示质数M₁降采样后方位向第m₁个慢时刻，N₁的值为 表示取整运算符号；

对于步骤2中得到的原始回波数据S(t,k)，t＝1,2,…,N_r，k＝1,2,…,N_a，按照质数M₂＝5的倍数进行方向数据抽取，即方位向慢时刻k的值被M₂＝5整除时提取出对应的方位向慢时刻k、距离向第t个快时刻方位向第k个慢时刻原始回波数据S(t,k)和雷达系统天线位置P_S(k)，得到质数M₂降采样后的SAR雷达系统方位向采样点总数，记为N₂，得到质数M₂降采样后的距离向第t个快时刻方位向第m₂个慢时刻的SAR原始回波数据，记为S₂(t,m₂)，t＝1,2,…,N_r，m₂＝1,2,…,N₂，得到质数M₂降采样后方位向第m₂个慢时刻的雷达系统天线位置，记为P₂(m₂)，m₂＝1,2,…,N₂，其中m₂为自然数，m₂表示质数M₂降采样后方位向第m₂个慢时刻，N₂的值为

步骤4、对互质降采样后SAR数据分别建立线性测量模型：

将步骤3中得到的质数M₁＝4降采样后的SAR原始回波数据S₁(t,m₁)，t＝1,2,…,N_r，m₁＝1,2,…,N₁，按先行后列的顺序组成一个原始回波数据列向量，记为y₁，即向量y₁的值表示为y₁＝[S₁(1,1),S₁(2,1),…,S₁(N_r,N₁)]^T，其中S₁(1,1)表示质数M₁降采样后距离向第1个快时刻方位向第1个慢时刻的SAR原始回波数据，S₁(2,1)表示质数M₁降采样后距离向第2个快时刻方位向第1个慢时刻的SAR原始回波数据，S₁(N_r,N₁)表示质数M₁降采样后距离向第N_r个快时刻方位向第N₁个慢时刻的SAR原始回波数据，上标T表示转置运算符号，N₁为步骤3中得到的质数M₁降采样后SAR雷达系统方位向采样点总数N₁＝512；

将步骤3中得到的质数M₂＝5降采样后的SAR原始回波数据S₂(t,m₂)，t＝1,2,…,N_r，m₂＝1,2,…,N₂，按先行后列的顺序组成一个原始回波数据列向量，记为y₂，即向量y₂的值表示为y₂＝[S₂(1,1),S₂(2,1),…,S₂(N_r,N₂)]^T，其中S₂(1,1)表示质数M₂降采样后距离向第1个快时刻方位向第1个慢时刻的SAR原始回波数据，S₂(2,1)表示质数M₂降采样后距离向第2个快时刻方位向第1个慢时刻的SAR原始回波数据，S₂(N_r,N₂)表示质数M₂降采样后距离向第N_r个快时刻方位向第N₂个慢时刻的SAR原始回波数据，N₂为步骤3中得到的质数M₂降采样后SAR雷达系统方位向采样点总数N₂＝409；

利用步骤1初始化得到的SAR雷达系统参数、成像空间参数和步骤3得到的质数M₁＝4降采样后方位向第m₁个慢时刻的雷达系统天线位置P₁(m₁)，m₁＝1,2,…,N₁，采用标准压缩感知SAR测量矩阵构造方法得到质数M₁降采样后SAR原始回波的测量矩阵，记为A₁，矩阵A₁的维数大小为N_rM₁×N_XN_Y，其中N_r为步骤1初始化得到的SAR雷达系统距离向采样点数N_r＝2048，N_X和N_Y分别为步骤1初始化得到的成像空间X轴和Y轴的网格点总数N_X＝1000和N_Y＝1000；

利用步骤1初始化得到的SAR雷达系统参数、成像空间参数和步骤3得到的质数M₂＝5降采样后方位向第m₂个慢时刻的雷达系统天线位置P₂(m₂)，m₂＝1,2,…,N₂，采用标准压缩感知SAR测量矩阵构造方法得到质数M₂降采样后SAR原始回波的测量矩阵，记为A₂，矩阵A₂的维数大小为N_rM₂×N_XN_Y；

步骤5、采用压缩感知稀疏重构方法进行稀疏成像：

采用标准压缩感知稀疏重构方法对最优化问题进行求解，得到质数M₁降采样后SAR回波数据稀疏重构成像的列向量，记为向量的维数为N_XN_Y×1，向量中第w个元素记为w＝1,2,…,N_XN_Y，其中w为自然数，w表示为一个向量中第w个元素，y₁为步骤4中得到的质数M₁降采样后回波数据组成的列向量数据，A₁为步骤4中得到的质数M₁降采样后SAR原始回波的测量矩阵，表示求取满足条件的自变量x的最小值，||·||₁为向量L1范数，||·||₂为向量L2范数；

采用标准压缩感知稀疏重构方法对最优化问题进行求解，得到质数M₂降采样后SAR回波数据稀疏重构成像的列向量，记为向量的维数为N_XN_Y×1，向量中第w个元素记为w＝1,2,…,N_XN_Y，其中y₂为步骤4中得到的质数M₂降采样后回波数据组成的列向量数据，A₂为步骤4中得到的质数M₂降采样后SAR原始回波的测量矩阵；

步骤6、对互质降采样得到稀疏成像结果进行栅瓣抑制融合：

初始化一个维数为N_XN_Y×1的向量，记为其中向量第w个元素记为w＝1,2,…,N_XN_Y，N_X和N_Y分别为步骤1初始化得到的成像空间X轴和Y轴的网格点总数N_X＝1000和N_Y＝1000；将步骤5得到的向量和向量逐个元素进行比较，若向量元素w＝1,2,…,N_XN_Y，则向量中第w个元素赋值为否则向量中第w个元素赋值为最后得到的向量即为互质降采样栅瓣抑制融合后的稀疏重构向量；

步骤7、重排得到最终的稀疏成像结果：

将步骤6得到的稀疏重构向量按照先行后列的顺序组成一个维数为N_X×N_Y的二维矩阵，记为I，I即为互质降采样间歇SAR稀疏成像的最终图像，其中N_X和N_Y分别为步骤1初始化得到的成像空间X轴和Y轴的网格点总数N_X＝1000和N_Y＝1000。

Claims (1)

1.一种互质降采样间歇合成孔径雷达稀疏成像方法，其特征是它包括以下几个步骤：

步骤1、初始化SAR雷达系统参数和成像空间：

初始化SAR雷达系统参数，包括：雷达载波波长，记为λ，雷达中心频率，记为f₀，雷达发射信号带宽，记为B，雷达发射脉冲时宽，记为T_r，雷达采样频率，记为F_s，雷达波束中心入射角，记为θ，雷达脉冲重复频率，记为PRF，平台运动速度矢量，记为V_r，雷达系统距离向采样点数，记为N_r，雷达系统方位向采样点数，记为N_a，雷达系统天线初始位置，记为P_S(0)；上述参数中，光在空气中的传播速度c为常数，雷达载波波长λ、雷达中心频率f₀、雷达发射信号带宽B、雷达发射脉冲时宽T_r、雷达采样频率F_s、雷达波束中心入射角θ、雷达脉冲重复频率PRF、平台运动速度矢量V_r、雷达系统距离向采样点数N_r、雷达系统方位向采样点数N_a、雷达系统天线初始位置P_S(0)在SAR雷达成像观测方案设计中已经确定；

初始化SAR成像空间参数，包括：以SAR雷达波束照射场区域的水平地面坐标系作为SAR成像空间，成像空间中心坐标位置位于[0,0]，成像空间的横轴和纵轴分别记为X轴和Y轴，成像空间X轴和Y轴范围大小分别记为W_X和W_Y，将成像空间均匀划网格得到X轴和Y轴的网格点总数分别记为N_X和N_Y，成像空间X轴第i个和Y轴第j个网格点的位置，记为P_T(i,j)，i＝1,2,…,N_X，j＝1,2,…,N_Y，P_T(i,j)的值表示为其中i和j为自然数，i表示成像空间X轴第i个网格点，j表示成像空间Y轴第j个网格点；

步骤2、获取SAR雷达系统的原始回波数据和天线位置：

SAR雷达系统在距离向第t个快时刻方位向第k个慢时刻的原始回波数据，记为S(t,k)，t＝1,2,…,N_r，k＝1,2,…,N_a，其中t和k为自然数，t表示距离向第t个快时刻，k表示方位向第k个慢时刻，N_r为步骤1初始化得到的SAR雷达系统距离向采样点数，N_a为步骤1初始化得到的SAR雷达系统方位向采样点数；在SAR系统实际成像中，原始回波数据S(t,k)，t＝1,2,…,N_r，k＝1,2,…,N_a，由SAR雷达系统数据接收机提供；在仿真成像中，原始回波数据S(t,k)，t＝1,2,…,N_r，k＝1,2,…,N_a，由标准合成孔径雷达原始回波仿真方法得到；

采用公式计算得到方位向第k个慢时刻的雷达系统天线位置，记为P_S(k)，k＝1,2,…,N_a，其中P_S(0)为步骤1初始化得到的雷达系统天线初始位置，V_r为步骤1初始化得到的平台运动速度矢量，PRF为步骤1初始化得到的雷达脉冲重复频率；

步骤3、对SAR方位向脉冲进行互质降采样的数据抽取：

初始化互质的两个整数，分别记为M₁和M₂；

对于步骤2中得到的原始回波数据S(t,k)，t＝1,2,…,N_r，k＝1,2,…,N_a，按照质数M₁的倍数进行方向数据抽取，即方位向慢时刻k的值被M₁整除时提取出对应的方位向慢时刻k、距离向第t个快时刻方位向第k个慢时刻原始回波数据S(t,k)和雷达系统天线位置P_S(k)，得到质数M₁降采样后的SAR雷达系统方位向采样点总数，记为N₁，得到质数M₁降采样后的距离向第t个快时刻方位向第m₁个慢时刻的SAR原始回波数据，记为S₁(t,m₁)，t＝1,2,…,N_r，m₁＝1,2,…,N₁，得到质数M₁降采样后方位向第m₁个慢时刻的雷达系统天线位置，记为P₁(m₁)，m₁＝1,2,…,N₁，其中m₁为自然数，m₁表示质数M₁降采样后方位向第m₁个慢时刻，N₁的值为 表示取整运算符号；

对于步骤2中得到的原始回波数据S(t,k)，t＝1,2,…,N_r，k＝1,2,…,N_a，按照质数M₂的倍数进行方向数据抽取，即方位向慢时刻k的值被M₂整除时提取出对应的方位向慢时刻k、距离向第t个快时刻方位向第k个慢时刻原始回波数据S(t,k)和雷达系统天线位置P_S(k)，得到质数M₂降采样后的SAR雷达系统方位向采样点总数，记为N₂，得到质数M₂降采样后的距离向第t个快时刻方位向第m₂个慢时刻的SAR原始回波数据，记为S₂(t,m₂)，t＝1,2,…,N_r，m₂＝1,2,…,N₂，得到质数M₂降采样后方位向第m₂个慢时刻的雷达系统天线位置，记为P₂(m₂)，m₂＝1,2,…,N₂，其中m₂为自然数，m₂表示质数M₂降采样后方位向第m₂个慢时刻，N₂的值为

步骤4、对互质降采样后SAR数据分别建立线性测量模型：

将步骤3中得到的质数M₁降采样后的SAR原始回波数据S₁(t,m₁)，t＝1,2,…,N_r，m₁＝1,2,…,N₁，按先行后列的顺序组成一个原始回波数据列向量，记为y₁，即向量y₁的值表示为y₁＝[S₁(1,1),S₁(2,1),…,S₁(N_r,N₁)]^T，其中S₁(1,1)表示质数M₁降采样后距离向第1个快时刻方位向第1个慢时刻的SAR原始回波数据，S₁(2,1)表示质数M₁降采样后距离向第2个快时刻方位向第1个慢时刻的SAR原始回波数据，S₁(N_r,N₁)表示质数M₁降采样后距离向第N_r个快时刻方位向第N₁个慢时刻的SAR原始回波数据，上标T表示转置运算符号，N₁为步骤3中得到的质数M₁降采样后SAR雷达系统方位向采样点总数；

将步骤3中得到的质数M₂降采样后的SAR原始回波数据S₂(t,m₂)，t＝1,2,…,N_r，m₂＝1,2,…,N₂，按先行后列的顺序组成一个原始回波数据列向量，记为y₂，即向量y₂的值表示为y₂＝[S₂(1,1),S₂(2,1),…,S₂(N_r,N₂)]^T，其中S₂(1,1)表示质数M₂降采样后距离向第1个快时刻方位向第1个慢时刻的SAR原始回波数据，S₂(2,1)表示质数M₂降采样后距离向第2个快时刻方位向第1个慢时刻的SAR原始回波数据，S₂(N_r,N₂)表示质数M₂降采样后距离向第N_r个快时刻方位向第N₂个慢时刻的SAR原始回波数据，N₂为步骤3中得到的质数M₂降采样后SAR雷达系统方位向采样点总数；

利用步骤1初始化得到的SAR雷达系统参数、成像空间参数和步骤3得到的质数M₁降采样后方位向第m₁个慢时刻的雷达系统天线位置P₁(m₁)，m₁＝1,2,…,N₁，采用标准压缩感知SAR测量矩阵构造方法得到质数M₁降采样后SAR原始回波的测量矩阵，记为A₁，矩阵A₁的维数大小为N_rM₁×N_XN_Y，其中N_r为步骤1初始化得到的SAR雷达系统距离向采样点数，N_X和N_Y分别为步骤1初始化得到的成像空间X轴和Y轴的网格点总数；

利用步骤1初始化得到的SAR雷达系统参数、成像空间参数和步骤3得到的质数M₂降采样后方位向第m₂个慢时刻的雷达系统天线位置P₂(m₂)，m₂＝1,2,…,N₂，采用标准压缩感知SAR测量矩阵构造方法得到质数M₂降采样后SAR原始回波的测量矩阵，记为A₂，矩阵A₂的维数大小为N_rM₂×N_XN_Y；

步骤5、采用压缩感知稀疏重构方法进行稀疏成像：

采用标准压缩感知稀疏重构方法对最优化问题方程进行求解，得到质数M₁降采样后SAR回波数据稀疏重构成像的列向量，记为向量的维数为N_XN_Y×1，向量中第w个元素记为w＝1,2,…,N_XN_Y，其中w为自然数，w表示为一个向量中第w个元素，y₁为步骤4中得到的质数M₁降采样后回波数据组成的列向量数据，A₁为步骤4中得到的质数M₁降采样后SAR原始回波的测量矩阵，表示求取满足条件的自变量x的最小值，||·||₁为向量L1范数，||·||₂为向量L2范数；

采用标准压缩感知稀疏重构方法对最优化问题方程进行求解，得到质数M₂降采样后SAR回波数据稀疏重构成像的列向量，记为向量的维数为N_XN_Y×1，向量中第w个元素记为w＝1,2,…,N_XN_Y，其中y₂为步骤4中得到的质数M₂降采样后回波数据组成的列向量数据，A₂为步骤4中得到的质数M₂降采样后SAR原始回波的测量矩阵；

步骤6、对互质降采样得到稀疏成像结果进行栅瓣抑制融合：

初始化一个维数为N_XN_Y×1的向量，记为其中向量第w个元素记为w＝1,2,…,N_XN_Y；将步骤5得到的向量和向量逐个元素进行比较，若w＝1,2,…,N_XN_Y，则向量中第w个元素赋值为否则向量中第w个元素赋值为最后得到的向量即为互质降采样栅瓣抑制融合后的稀疏重构向量；

步骤7、重排得到最终的稀疏成像结果：

将步骤6得到的稀疏重构向量按照先行后列的顺序组成一个维数为N_X×N_Y的二维矩阵，记为I，I即为互质降采样间歇SAR稀疏成像的最终图像，其中N_X和N_Y分别为步骤1初始化得到的成像空间X轴和Y轴的网格点总数。