CN104111458A - 基于双重稀疏约束的压缩感知合成孔径雷达成像方法 - Google Patents

Abstract

本发明公开了一种基于双重稀疏约束的压缩感知合成孔径雷达成像的方法；主要解决降采样时传统合成孔径雷达成像方法不能得到完整SAR图像及时间损耗大的问题。其实现步骤是：1.雷达照射地面稀疏场景D，获得该场景的回波信号r；2.构造场景的方位向基矩阵A，场景的距离向基矩阵B；3.根据方位、距离向基矩阵获得方位向测量矩阵Θ，距离向测量矩阵Ω和回波测量矩阵S；4.根据方位向、距离向测量矩阵，回波测量矩阵构造拉格朗日函数f(Y,U)；5.用交替迭代乘子法求解步骤4中的拉格朗日函数f(Y,U)，得到场景D的图像。本发明能实现对降采样时合成孔径雷达的高分辨成像，且成像速度快，可用于大面积地形测绘和制图学。

Description

基于双重稀疏约束的压缩感知合成孔径雷达成像方法

技术领域

本发明属于雷达技术领域，具体地说是一种能够有效地处理大场景的压缩感知合成孔径雷达CS-SAR成像方法，可用于大面积地形测绘，制图学，全天候全球侦察。

背景技术

合成孔径雷达SAR是一种具有高分辨力的成像雷达，具有全天时，全天候成像能力，在军事和民用方面得到了广泛的应用。SAR的高分辨，在径向距离上依靠宽带信号，几百兆赫的频带可将距离分辨单元缩小到亚米级；在方位上则依靠雷达平台运动，等效地在空间形成很长的线性阵列，并将各次回波存储作合成的阵列处理。

Shannon-Nyquist采样定理指出，如果以不低于两倍信号带宽的频率对信号进行采样，再通过一定的处理，如利用低通滤波器从取样信号频谱中得到原信号的频谱在进行反变换，就可以精确重构该信号。传统合成孔径雷达成像算法，例如Range-Doppler算法，ChirpScaling算法，根据这一定理利用大于两倍带宽的采样频率对回波信号进行采样。较高的采样频率对接收端的A/D(Analog-to-digital)转换器来说无疑是一个巨大的挑战，同时大量的数据也增加了存储和传输的负担。

近几年兴起的压缩感知CS理论指出，可以利用非自适应的线性观测方法获得少量的观测向量，通过求解一个优化问题来高概率地重构出稀疏或者可压缩的信号，其核心思想是基于混叠采样特征解耦信号重建。根据这一理论，Jungang Yang等人在文献“Jungang Yang,John Thompson,Xiaotao Huang,Tian Jin,and Zhimin Zhou.Segmented Reconstruction forCompressed Sensing SAR Imaging.IEEE Transactions on Geoscience and Remote Sensing,Vol.51,NO.7,July2013.”提出了一种CS-SAR成像算法，将场景向量化并构造观测矩阵，再利用接收回波进行重构的方法。尽管CS-SAR成像算法可以在低于Shannon-Nyquist采样频率下进行成像，但是这种方法往往需要花费大量的时间，不能达到实时性的成像处理，且成像有效性较差，使得其只能处理较小的场景。

发明内容

本发明的目的在于针对上述已有技术的不足，提出一种基于双重稀疏约束的压缩感知合成孔径雷达成像方法，以实现直接对二维稀疏场景进行重构，在保持高分辨率的同时减少成像时间。

为实现上述目的，本发明的技术方案包括如下步骤：

(1)雷达载机沿着航向前行，不断向地面稀疏场景D发射线性调频脉冲信号其中为快时间，t_a为慢时间；

(2)雷达向地面稀疏场景D发射脉冲的同时接收该场景D的回波脉冲，得到回波信号为

$r(t_a,\hat{t})=\underset{i=1}{\overset{P}{\mathrm{\Σ}}}\underset{j=1}{\overset{Q}{\mathrm{\Σ}}}{G}_{\mathrm{ij}}{w}_a(t_a-\frac{x_i}{v})\exp (-j\frac{2\mathrm{\π}{v}^2}{\mathrm{\λ}{R}_0}{(t_a-\frac{x_i}{v})}^2)w_r(\hat{t}-\frac{2R_j}{C})\exp \left(\mathrm{j\π\γ}{(\hat{t}-\frac{2R_j}{C})}^2\right)\exp (-j\frac{4\mathrm{\π}{f}_c}{C}{R}_j)$

式中，f_c为载频，γ为调频斜率，v为载机速度，C为电磁波传播速度，w_a(·)为方位窗函数，w_r(·)为线性调频脉冲信号时间窗函数，P为地面稀疏场景D在方位向的离散网格个数，Q为地面稀疏场景D在距离向的离散网格个数，i为地面稀疏场景D在方位向的第i个离散网格，j为地面稀疏场景D在距离向的第j个离散网格，G_ij为(i,j)处散射点T的散射系数，R₀为雷达到地面稀疏场景D中心的垂直斜距，R_j为雷达到点目标T的垂直斜距，x_i为点目标T在三维空间中的x轴坐标；

(3)对回波信号进行二维离散采样，得到如下矩阵形式：

其中，M为方位向发射脉冲的数目，N为每个脉冲内，距离向采样点数；为回波信号在第(m,n)个采样时刻的样本值，其表示式如下：

$r(t_{a,m},{\hat{t}}_n)=\underset{i=1}{\overset{P}{\mathrm{\Σ}}}\underset{j=1}{\overset{Q}{\mathrm{\Σ}}}{G}_{\mathrm{ij}}{w}_a(t_{a,m}-\frac{x_i}{v})\exp (-j\frac{2\mathrm{\π}{v}^2}{\mathrm{\λ}{R}_0}{(t_{a,m}-\frac{x_i}{v})}^2)w_r({\hat{t}}_n-\frac{2R_j}{C})\exp \left(\mathrm{j\π\γ}{({\hat{t}}_n-\frac{2R_j}{C})}^2\right)\exp (-j\frac{4\mathrm{\π}{R}_j}{\mathrm{\λ}})$

式中，λ为载波波长；

(4)构造方位向基矩阵A：

其中， $a(t_{a,m},i)=w_a(t_{a,m}-\frac{x_i}{v})\exp (-j\frac{2\mathrm{\π}{v}^2}{\mathrm{\λ}{R}_0}{(t_{a,m}-\frac{x_i}{v})}^2);$

(5)构造距离向基矩阵B：

其中， $b(j,{\hat{t}}_n)=w_r({\hat{t}}_n-\frac{2R_j}{C})\exp \left(\mathrm{j\π\γ}{({\hat{t}}_n-\frac{2R_j}{C})}^2\right)\exp (-j\frac{4\mathrm{\π}{R}_j}{\mathrm{\λ}});$

(6)根据随机观测矩阵Φ_a,Φ_b，方位向基矩阵A，距离向基矩阵B，回波矩阵r，获得方位向测量矩阵Θ，距离向测量矩阵Ω和回波测量矩阵S

Θ＝Φ_aA

Ω＝BΦ_b，

S＝Φ_arΦ_b

(7)根据距离向测量矩阵Ω，地面稀疏场景D的散射系数矩阵G，得到分离变量矩阵Y和中间变量矩阵U：

Y＝GΩ

U＝G^H，

式中，H表示共轭转置；

(8)利用步骤(6)和步骤(7)得到的参数，构造拉格朗日函数f(Y,U)：

$\underset{Y,U}{\min }f(Y,U)=\{{\left|\right|Y\left|\right|}_{\mathrm{2,1}}+\mathrm{\α}{\left|\right|U\left|\right|}_{{2,l}_p}^{l_p}+\frac{{\mathrm{\β}}_1}{2}{\left|\right|\mathrm{\ΘY}-S\left|\right|}_F^2+\frac{{\mathrm{\β}}_2}{2}{\left|\right|{\mathrm{\Ω}}^{H}U-Y^H\left|\right|}_F^2\},$

其中表示求函数最小值的运算符，稀疏p范数参数0＜l_p≤1，α为正则参数，β₁为距离向惩罚因子，β₂为方位向惩罚因子，||·||_F表示求矩阵的Frobenius范数；||Y||_2,1表示先对Y的各行向量求2范数，再将所得的2范数结果构成的列向量求1范数；表示先对U的各行向量求2范数，再将所得的2范数结果构成列向量求l_p范数；

(9)用交替方向乘子法对步骤(8)中的拉格朗日函数f(Y,U)进行交替迭代求解，得到中间变量矩阵U的最终迭代结果U_*；

(10)对步骤(9)得到中间变量矩阵U的最终迭代结果U_*取共轭转置，得到地面稀疏场景的散射系数矩阵G＝(U_*)^H，再对该散射系数矩阵G取模值，得到地面稀疏场景D的图像。

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

第一，本发明通过对场景内散射点的散射系数矩阵及中间变量矩阵加以混合范数稀疏约束，使用二维矩阵稀疏重构的计算成像方法，突破了传统采样数据的限制，实现了合成孔径雷达的高分辨成像。

第二，本发明直接对二维稀疏场景图像进行重构成像，避免了场景图像向量化的重构方法产生的观测矩阵过大，计算量急剧增加的问题，具有大场景快速成像，小场景实时成像的优点，对成像速度要求高的应用场合有明显的优势。

附图说明

图1是本发明的实现流程图；

图2是本发明仿真使用的圆环稀疏场景的图像；

图3是本发明仿真使用的点目标稀疏场景的图像；

图4是用本发明方法重构出的圆环稀疏场景的图像。

图5是用本发明方法重构出的点目标稀疏场景的误差与时间关系图。

图6是用传统CS-SAR算法重构出的点目标稀疏场景的误差与时间关系图。

具体实施方法

下面结合附图对本发明做进一步的详细描述。

参照图1，本发明的具体实施步骤如下：

步骤1：获得离散二维回波矩阵r。

(1a)SAR雷达载机沿着航向前行，不断向地面稀疏场景D发射线性调频脉冲信号：

$s(t_a,\hat{t})=p\left(\hat{t}\right)\exp \left(j2\mathrm{\π}{f}_{c}t\right)$

其中f_c为载频，t为全时间，为快时间，t_a为慢时间，这三个时间的关系为 $p\left(\hat{t}\right)=w_r\left(\hat{t}\right)\exp \left(\mathrm{j\π\γ}{\hat{t}}^2\right)$ 为发射信号复包络，γ为调频斜率， $w_r\left(\hat{t}\right)=\mathrm{rect}(\hat{t}/T_r),$ 表示LFM信号的时间窗函数，T_r是脉冲持续的时间。

(1b)雷达向地面稀疏场景D发射脉冲的同时接收该场景D的回波脉冲，得到回波信号为

$r(t_a,\hat{t})=\underset{i=1}{\overset{P}{\mathrm{\Σ}}}\underset{j=1}{\overset{Q}{\mathrm{\Σ}}}{G}_{\mathrm{ij}}{w}_a(t_a-\frac{x_i}{v})\exp (-j\frac{2\mathrm{\π}{v}^2}{\mathrm{\λ}{R}_0}{(t_a-\frac{x_i}{v})}^2)w_r(\hat{t}-\frac{2R_j}{C})\exp \left(\mathrm{j\π\γ}{(\hat{t}-\frac{2R_j}{C})}^2\right)\exp (-j\frac{4\mathrm{\π}{f}_c}{C}{R}_j)$

其中，f_c为载频，γ为调频斜率，v为载机速度，C为电磁波传播速度，w_a(·)为方位窗函数，w_r(·)为线性调频脉冲信号时间窗函数，P为地面稀疏场景D在方位向的离散网格个数，Q为地面稀疏场景D在距离向的离散网格个数，i为地面稀疏场景D在方位向的第i个离散网格，j为地面稀疏场景D在距离向的第j个离散网格，G_ij为(i,j)处散射点T的散射系数，R₀为雷达到地面稀疏场景D中心的垂直斜距，R_j为雷达到点目标T的垂直斜距，x_i为点目标T在三维空间中的x轴坐标；

(1c)将步骤(1b)接收到稀疏场景D的回波进行二维离散采样，写成如下矩阵形式：

其中M为方位向发射脉冲的数目，N为每个脉冲内距离向采样点数；表示回波的第(m,n)个采样时刻的样本值，其表示式如下：

$r(t_{a,m},{\hat{t}}_n)=\underset{i=1}{\overset{P}{\mathrm{\Σ}}}\underset{j=1}{\overset{Q}{\mathrm{\Σ}}}{G}_{\mathrm{ij}}{w}_a(t_{a,m}-\frac{x_i}{v})\exp (-j\frac{2\mathrm{\π}{v}^2}{\mathrm{\λ}{R}_0}{(t_{a,m}-\frac{x_i}{v})}^2)w_r({\hat{t}}_n-\frac{2R_j}{C})\exp \left(\mathrm{j\π\γ}{({\hat{t}}_n-\frac{2R_j}{C})}^2\right)\exp (-j\frac{4\mathrm{\π}{R}_j}{\mathrm{\λ}})$

式中λ为载波波长。

步骤2：构造场景的方位向基矩阵A，场景的距离向基矩阵B。

(2a)构造场景的方位向基矩阵A：

其中 $a(t_{a,m},i)=w_a(t_{a,m}-\frac{x_i}{v})\exp (-j\frac{2\mathrm{\π}{v}^2}{\mathrm{\λ}{R}_0}{(t_{a,m}-\frac{x_i}{v})}^2),$ w_a(·)为方位窗函数，v为载机速度，λ为载波波长，M为方位向发射脉冲的数目，P为地面稀疏场景D在方位向的离散网格个数，R₀为雷达到地面稀疏场景D中心的垂直斜距，R_j为雷达到点目标T的垂直斜距，x_i为点目标T在三维空间中的x轴坐标；

(2b)构造场景的距离向基矩阵B：

其中 $b(j,{\hat{t}}_n)=w_r({\hat{t}}_n-\frac{2R_j}{C})\exp \left(\mathrm{j\π\γ}{({\hat{t}}_n-\frac{2R_j}{C})}^2\right)\exp (-j\frac{4\mathrm{\π}{R}_j}{\mathrm{\λ}}),$ w_r(·)为线性调频脉冲信号时间窗函数，C为电磁波传播速度，λ为载波波长，N为每个脉冲内距离向采样点数，Q为地面稀疏场景D在距离向的离散网格个数。

步骤3：根据步骤1和步骤2中的矩阵变量，获得方位向测量矩阵Θ，距离向测量矩阵Ω和回波测量矩阵S。

(3a)根据方位向随机观测矩阵Φ_a，方位向基矩阵A，获得方位向测量矩阵Θ：

Θ＝Φ_aA；

(3b)根据距离向随机观测矩阵Φ_b，距离向基矩阵B，获得距离向测量矩阵Ω：

Ω＝BΦ_b；

(3c)根据方位向随机观测矩阵Φ_a，距离向随机观测矩阵Φ_b，回波矩阵r，获得回波测量矩阵S：

S＝Φ_arΦ_b。

步骤4：根据步骤3中的矩阵变量构造拉格朗日函数f(Y,U)。

(4a)根据距离向测量矩阵Ω，地面稀疏场景D的散射系数矩阵G，得到分离变量矩阵Y和中间变量矩阵U：

Y＝GΩ

U＝G^H，

式中，H表示共轭转置；

(4b)对分离变量矩阵Y和中间变量矩阵U加以稀疏约束，得到如下双重稀疏约束函数：

$\underset{Y,U}{\min }\{{\left|\right|Y\left|\right|}_{\mathrm{2,1}}+\mathrm{\α}{\left|\right|U\left|\right|}_{2,l_p}^{l_p}\}$

$s.t.\left\{\begin{array}{c}S=\mathrm{\ΘY}+E\\ Y^H={\mathrm{\Ω}}^{H}U\end{array}\right.,$

其中表示求函数最小值的运算符，稀疏p范数参数0＜l_p≤1，α为正则参数，Ω距离向测量矩阵，Θ为方位向测量矩阵，S为回波测量矩阵，E是与S同维度的噪声矩阵，H表示共轭转置；||Y||_2,1表示先对Y的各行向量求2范数，再将所得的2范数结果构成的列向量求1范数；表示先对U的各行向量求2范数，再将所得的2范数结果构成列向量求l_p范数；

(4c)根据(4b)中的约束函数构造如下拉格朗日函数：

$\underset{Y,U}{\min }f(Y,U)=\{{\left|\right|Y\left|\right|}_{\mathrm{2,1}}+\mathrm{\α}{\left|\right|U\left|\right|}_{2,l_p}^{l_p}+\frac{{\mathrm{\β}}_1}{2}{\left|\right|\mathrm{\ΘY}-S\left|\right|}_F^2+\frac{{\mathrm{\β}}_2}{2}{\left|\right|{\mathrm{\Ω}}^{H}U-Y^H\left|\right|}_F^2\}$

其中β₁为距离向惩罚因子，β₂为方位向惩罚因子，||·||_F表示求矩阵的Frobenius范数。

步骤5：用交替迭代乘子法求解步骤4中的拉格朗日函数f(Y,U)，得到成像后的目标场景的图像。

(5a)初始化中间变量矩阵U为全1矩阵，分离变量矩阵Y为全1矩阵，设置正则化参数α＞0，距离向惩罚因子β₁＞0，方位向惩罚因子β₂＞0，稀疏p范数参数0＜l_p≤1，迭代步数k＝0，迭代终止条件ε＝10^-6；

(5b)构造中间变量对角矩阵Π(U_k)和分离变量对角矩阵Π(Y_k)：

$\mathrm{\Π}\left(U_k\right)=\mathrm{diag}\left({({\left|\right|U_k(j,:)\left|\right|}_2^2+\mathrm{\δ})}^{\frac{l_p}{2}-1}\right),$

$\mathrm{\Π}\left(Y_k\right)=\mathrm{diag}\left({({\left|\right|Y_k(i,:)\left|\right|}_2^2+\mathrm{\δ})}^{-\frac{1}{2}}\right),$

式中diag(·)表示将向量排列成对角矩阵的形式，||·||₂表示求向量的2范数，修正因子δ＞0，U_k为第k次迭代后中间变量矩阵U的值，Y_k为第k次迭代后分离变量矩阵Y的值；U_k(j,:)为第k次迭代后中间变量矩阵U_k的第j行构成的向量，j＝1,2,…,Q，Q为地面稀疏场景D在距离向的离散网格个数；Y_k(i,:)为第k次迭代后分离变量矩阵Y_k的第i行构成的向量，i＝1,2,…，P，P为地面稀疏场景D在方位向的离散网格个数；

(5c)按照如下公式更新中间变量矩阵U，得到第k+1次迭代后中间变量矩阵U_k+1：

$U_{k+1}={\mathrm{\Π}}^{-1}\left(U_k\right)\mathrm{\Ω}{({\mathrm{\Ω}}^H{\mathrm{\Π}}^{-1}\left(U_k\right)\mathrm{\Ω}+\frac{l_p}{{\mathrm{\β}}_2}I)}^{-1}{\left(Y_k\right)}^H,$

式中，Ω为距离向测量矩阵，I为单位矩阵，Π^-1(·)表示对角矩阵Π(·)的逆矩阵；

(5d)按照如下公式更新分离变量矩阵Y，得到第k+1次迭代后分离变量矩阵Y_k+1：

式中，Θ为方位向测量矩阵，S为回波测量矩阵，为矩阵Moore–Penrose的逆；

(5e)根据第k次迭代后中间变量矩阵U_k，第k+1次迭代后中间变量矩阵U_k+1，得到均方误差M：

$M=\frac{{\left|\right|U_{k+1}-U_k\left|\right|}_F^2}{{\left|\right|U_k\left|\right|}_F^2},$

式中，||·||_F表示求矩阵的Frobenius范数；

(5f)判断均方误差M≤ε是否成立，若成立执行步骤(9g)；否则令k＝k+1，返回步骤(9b)继续迭代运行，其中ε为迭代终止条件；

(5g)令U_*＝U_k+1，输出中间变量矩阵U的最终迭代结果U_*；

(5h)对步骤(5g)得到中间变量矩阵U的最终迭代结果U_*取共轭转置，得到地面稀疏场景的散射系数矩阵G＝(U_*)^H，再对该散射系数矩阵G取模值，得到地面稀疏场景D的图像。

本发明的效果可以通过下述仿真实验加以说明：

1.仿真条件

(1a)运行平台配置：

CPU：Intel(R)Core(TM)2Duo CPU E45002.20GHz；

内存：2GB(亿能DDR2667MHZ)

操作系统：Windows7旗舰版32位SP1操作系统；

仿真软件：MATLAB R(2011b)。

(1b)仿真参数设置

发射信号采用线性调频信号，发射信号参数以及实验仿真参数设置如表1所示。圆环稀疏场景的方位向网格个数P＝201，距离向网格个数Q＝201。

表1 发射信号参数以及实验仿真参数

2.仿真内容与结果

仿真1，根据表1的仿真参数，用合成孔径雷达对图2所示的圆环稀疏场景，进行探测并获取回波，用本发明方法重构出圆环稀疏场景的图像，如图4所示。

仿真2，根据表1的仿真参数，用合成孔径雷达对如图3所示的点目标稀疏场景，进行探测并获得回波，用现有的将二维场景向量化的CS-SAR算法和本发明方法对点目标稀疏场景进行稀疏重构，并记录这两种不同成像方法对场景进行重构的误差与时间关系，结果如图5和图6所示。其中图5是现有的将二维场景向量化的CS-SAR算法的仿真结果，图6是本发明方法的仿真结果。

3.仿真结果分析

由重构出的圆环稀疏场景的图像4可以看出，本发明方法实现了对稀疏场景的高分辨重构；

由图5和图6可以看出，相对于现有的将二维场景向量化的CS-SAR算法，本发明的方法大大加速了对场景的成像时间，实现了小场景实时成像，大场景快速成像的目的，具有显著的实用性。

Claims (4)

1.一种基于双重稀疏约束的压缩感知合成孔径雷达成像方法，包括如下步骤：

(1)雷达载机沿着航向前行，不断向地面稀疏场景D发射线性调频脉冲信号其中为快时间，t_a为慢时间；

(2)雷达向地面稀疏场景D发射脉冲的同时接收该场景D的回波脉冲，得到回波信号为

$r(t_a,\hat{t})=\underset{i=1}{\overset{P}{\mathrm{\Σ}}}\underset{j=1}{\overset{Q}{\mathrm{\Σ}}}{G}_{\mathrm{ij}}{w}_a(t_a-\frac{x_i}{v})\exp (-j\frac{2\mathrm{\π}{v}^2}{\mathrm{\λ}{R}_0}{(t_a-\frac{x_i}{v})}^2)w_r(\hat{t}-\frac{2R_j}{C})\exp \left(\mathrm{j\π\γ}{(\hat{t}-\frac{2R_j}{C})}^2\right)\exp (-j\frac{4\mathrm{\π}{f}_c}{C}{R}_j)$

式中，f_c为载频，γ为调频斜率，v为载机速度，C为电磁波传播速度，w_a(·)为方位窗函数，w_r(·)为线性调频脉冲信号时间窗函数，P为地面稀疏场景D在方位向的离散网格个数，Q为地面稀疏场景D在距离向的离散网格个数，i为地面稀疏场景D在方位向的第i个离散网格，j为地面稀疏场景D在距离向的第j个离散网格，G_ij为(i,j)处散射点T的散射系数，R₀为雷达到地面稀疏场景D中心的垂直斜距，R_j为雷达到点目标T的垂直斜距，x_i为点目标T在三维空间中的x轴坐标；

(3)对回波信号进行二维离散采样，得到如下矩阵形式：

其中，M为方位向发射脉冲的数目，N为每个脉冲内，距离向采样点数；为回波信号在第(m,n)个采样时刻的样本值，其表示式如下：

$r(t_{a,m},{\hat{t}}_n)=\underset{i=1}{\overset{P}{\mathrm{\Σ}}}\underset{j=1}{\overset{Q}{\mathrm{\Σ}}}{G}_{\mathrm{ij}}{w}_a(t_{a,m}-\frac{x_i}{v})\exp (-j\frac{2\mathrm{\π}{v}^2}{\mathrm{\λ}{R}_0}{(t_{a,m}-\frac{x_i}{v})}^2)w_r({\hat{t}}_n-\frac{2R_j}{C})\exp \left(\mathrm{j\π\γ}{({\hat{t}}_n-\frac{2R_j}{C})}^2\right)\exp (-j\frac{4\mathrm{\π}{R}_j}{\mathrm{\λ}})$

式中，λ为载波波长；

(4)构造场景的方位向基矩阵A：

其中， $a(t_{a,m},i)=w_a(t_{a,m}-\frac{x_i}{v})\exp (-j\frac{2\mathrm{\π}{v}^2}{\mathrm{\λ}{R}_0}{(t_{a,m}-\frac{x_i}{v})}^2);$

(5)构造场景的距离向基矩阵B：

其中， $b(j,{\hat{t}}_n)=w_r({\hat{t}}_n-\frac{2R_j}{C})\exp \left(\mathrm{j\π\γ}{({\hat{t}}_n-\frac{2R_j}{C})}^2\right)\exp (-j\frac{4\mathrm{\π}{R}_j}{\mathrm{\λ}});$

(6)根据方位向随机观测矩阵Φ_a，距离向随机观测矩阵Φ_b，方位向基矩阵A，距离向基矩阵B，回波矩阵r，获得方位向测量矩阵Θ，距离向测量矩阵Ω和回波测量矩阵S：

Θ＝Φ_aA

Ω＝BΦ_b，

S＝Φ_arΦ_b

(7)根据距离向测量矩阵Ω，地面稀疏场景D的散射系数矩阵G，得到分离变量矩阵Y和中间变量矩阵U：

Y＝GΩ

U＝G^H，

式中，H表示共轭转置；

(8)利用步骤(6)和步骤(7)得到的参数，构造拉格朗日函数f(Y,U)：

$\underset{Y,U}{\min }f(Y,U)=\{{\left|\right|Y\left|\right|}_{\mathrm{2,1}}+\mathrm{\α}{\left|\right|U\left|\right|}_{{2,l}_p}^{l_p}+\frac{{\mathrm{\β}}_1}{2}{\left|\right|\mathrm{\ΘY}-S\left|\right|}_F^2+\frac{{\mathrm{\β}}_2}{2}{\left|\right|{\mathrm{\Ω}}^{H}U-Y^H\left|\right|}_F^2\},$

其中表示求函数最小值的运算符，稀疏p范数参数0＜l_p≤1，α为正则参数，β₁为距离向惩罚因子，β₂为方位向惩罚因子，||·||_F表示求矩阵的Frobenius范数；||Y||_2,1表示先对Y的各行向量求2范数，再将所得的2范数结果构成的列向量求1范数；表示先对U的各行向量求2范数，再将所得的2范数结果构成列向量求l_p范数；

(9)用交替方向乘子法对步骤(8)中的拉格朗日函数f(Y,U)进行交替迭代求解，得到中间变量矩阵U的最终迭代结果U_*；

(10)对步骤(9)得到中间变量矩阵U的最终迭代结果U_*取共轭转置，得到地面稀疏场景的散射系数矩阵G＝(U_*)^H，再对该散射系数矩阵G取模值，得到地面稀疏场景D的图像。

2.根据权利要求1所述的基于双重稀疏约束的压缩感知合成孔径雷达成像方法，其中所述步骤(1)线性调频脉冲信号表示如下：

$s(t_a,\hat{t})=p\left(\hat{t}\right)\exp \left(j2\mathrm{\π}{f}_{c}t\right),$

式中f_c为载频，t为全时间，为快时间，t_a为慢时间，这三个时间的关系为 $p\left(\hat{t}\right)=w_r\left(\hat{t}\right)\exp \left(\mathrm{j\π\γ}{\hat{t}}^2\right)$ 为发射信号复包络， $w_r\left(\hat{t}\right)=\mathrm{rect}(\hat{t}/T_r)$ 表示线性调频脉冲信号的时间窗函数，T_r是脉冲持续的时间，γ为调频斜率。

3.根据权利要求1所述的基于双重稀疏约束的压缩感知合成孔径雷达成像方法，其中所述步骤(8)构造拉格朗日函数f(Y,U)，按如下步骤进行：

(8a)对分离变量矩阵Y和中间变量矩阵U加以稀疏约束，得到如下双重稀疏约束函数：

$\underset{Y,U}{\min }\{{\left|\right|Y\left|\right|}_{\mathrm{2,1}}+\mathrm{\α}{\left|\right|U\left|\right|}_{2,l_p}^{l_p}\}$

$s.t.\left\{\begin{array}{c}S=\mathrm{\ΘY}+E\\ Y^H={\mathrm{\Ω}}^{H}U\end{array}\right.,$

其中表示求函数最小值的运算符，稀疏p范数参数0＜l_p≤1，α为正则参数，Ω距离向测量矩阵，Θ为方位向测量矩阵，S为回波测量矩阵，E是与S同维度的噪声矩阵，H表示共轭转置；||Y||_2,1表示先对Y的各行向量求2范数，再将所得的2范数结果构成的列向量求1范数；表示先对U的各行向量求2范数，再将所得的2范数结果构成列向量求l_p范数；

(8b)根据(8a)中的约束函数构造如下拉格朗日函数

$\underset{Y,U}{\min }f(Y,U)=\{{\left|\right|Y\left|\right|}_{\mathrm{2,1}}+\mathrm{\α}{\left|\right|U\left|\right|}_{2,l_p}^{l_p}+\frac{{\mathrm{\β}}_1}{2}{\left|\right|\mathrm{\ΘY}-S\left|\right|}_F^2+\frac{{\mathrm{\β}}_2}{2}{\left|\right|{\mathrm{\Ω}}^{H}U-Y^H\left|\right|}_F^2\}$

其中β₁为距离向惩罚因子，β₂为方位向惩罚因子，||·||_F表示求矩阵的Frobenius范数。

4.根据权利要求1所述的基于双重稀疏约束的压缩感知合成孔径雷达成像方法，其中所述步骤(9)用交替方向乘子法对步骤(8)中的拉格朗日函数f(Y,U)进行交替迭代求解，按如下步骤进行：

(9a)初始化中间变量矩阵U为全1矩阵，分离变量矩阵Y为全1矩阵，设置正则化参数α＞0，距离向惩罚因子β₁＞0，方位向惩罚因子β₂＞0，稀疏p范数参数0＜l_p≤1，迭代步数k＝0，迭代终止条件ε＝10^-6；

(9b)构造中间变量对角矩阵Π(U_k)和分离变量对角矩阵Π(Y_k)：

$\mathrm{\Π}\left(U_k\right)=\mathrm{diag}\left({({\left|\right|U_k(j,:)\left|\right|}_2^2+\mathrm{\δ})}^{\frac{l_p}{2}-1}\right),$

$\mathrm{\Π}\left(Y_k\right)=\mathrm{diag}\left({({\left|\right|Y_k(i,:)\left|\right|}_2^2+\mathrm{\δ})}^{-\frac{1}{2}}\right),$

式中diag(·)表示将向量排列成对角矩阵的形式，||·||₂表示求向量的2范数，修正因子δ＞0，U_k为第k次迭代后中间变量矩阵U的值，Y_k为第k次迭代后分离变量矩阵Y的值；U_k(j,:)为第k次迭代后中间变量矩阵U_k的第j行构成的向量，j＝1,2,…,Q，Q为地面稀疏场景D在距离向的离散网格个数；Y_k(i,:)为第k次迭代后分离变量矩阵Y_k的第i行构成的向量，i＝1,2,…,P，P为地面稀疏场景D在方位向的离散网格个数；

(9c)按照如下公式更新中间变量矩阵U，得到第k+1次迭代后中间变量矩阵U_k+1：

$U_{k+1}={\mathrm{\Π}}^{-1}\left(U_k\right)\mathrm{\Ω}{({\mathrm{\Ω}}^H{\mathrm{\Π}}^{-1}\left(U_k\right)\mathrm{\Ω}+\frac{l_p}{{\mathrm{\β}}_2}I)}^{-1}{\left(Y_k\right)}^H,$

式中，Ω为距离向测量矩阵，I为单位矩阵，Π^-1(·)表示对角矩阵Π(·)的逆；

(9d)按照如下公式更新分离变量矩阵Y，得到第k+1次迭代后分离变量矩阵Y_k+1：

式中，Θ为方位向测量矩阵，S为回波测量矩阵，为矩阵Moore–Penrose的逆；

(9e)根据第k次迭代后中间变量矩阵U_k，第k+1次迭代后中间变量矩阵U_k+1，得到均方误差M：

$M=\frac{{\left|\right|U_{k+1}-U_k\left|\right|}_F^2}{{\left|\right|U_k\left|\right|}_F^2},$

式中，||·||_F表示求矩阵的Frobenius范数；

(9f)判断均方误差M≤ε是否成立，若成立执行步骤(9g)；否则令k＝k+1，返回步骤(9b)继续迭代运行，其中ε为迭代终止条件；

(9g)令U_*＝U_k+1，输出中间变量矩阵U的最终迭代结果U_*。