CN102645651B - 一种sar层析超分辨成像方法 - Google Patents

Abstract

本发明公开了一种SAR层析超分辨成像方法，本发的方法明首先对SAR层析各基线得到的二维成像结果进行配准处理，逐像素按照基线顺序构建高度向信号，并进行调频校正，构造冗余矩阵；根据高度向散射系数的稀疏性特点，将高度向成像问题建模为稀疏信号的重构问题，利用加权最小范数的约束条件，通过迭代求其稀疏解，实现目标高度维成像。将本发明的方法应用于SAR层析高度维成像，解决了同一区域的二维SAR图像数量少和(/或)轨迹分布不均匀且只有很少部分的轨迹间隔(基线)满足奈奎斯特采样定理的问题，同时本发明的方法能使目标能量更集中，即提高了分辨率，解决了基线数目少的问题。

Description

一种SAR层析超分辨成像方法

技术领域

本发明属于电子信号处理技术领域，涉及空间遥感和空地观测信息处理技术领域，特别涉及机载合成孔径雷达(Synthetic Aperture Radar，SAR)层析三维成像技术。

背景技术

合成孔径雷达是二十世纪雷达技术发展的重要里程碑，其利用雷达回波信号的相关性，累积雷达运动过程中回波信号的多普勒频移，从而在雷达的运动方向上合成等效的雷达孔径，实现方位向的高分辨成像。因SAR采用主动式工作模式，对自身发射电磁波的反射回波进行成像处理，不受光照、温度等外界环境的限制，可实现全天时、全天候的区域监测成像，且对植被、沙漠覆盖等介质具有穿透能力而在灾害评估、环境监测、海洋观测、资源勘查、植被监测、测绘和军事等领域具有广泛的应用前景。

然而，传统SAR是对三维场景的二维成像，且由于距离向成像位置由目标与雷达的距离决定，SAR成像结果实际上是具有相同距离的所有散射体散射系数的迭加，因此在二维SAR成像中存在着圆柱对称模糊(即成像结果无法分辨具有相同距离的不同散射体)、叠岩模糊(雷达接收到目标上部反射回波先于下部反射回波)等问题，难以满足越来越高的成像精度和复杂环境侦察的要求。而问题的根源在于二维SAR成像的本质是三维场景空间到二维成像平面沿距离为半径的圆柱面的投影。因此，解决问题的最直接方法就是获得目标的第三维信息，避免三维空间到二维空间平面投影的信息损失。

SAR层析(Tomography SAR，TomoSAR)三维成像是传统二维SAR成像的三维扩展，通过将从一定数量满足一定分布规律的雷达运动轨迹上获得的二维SAR图像联合，在二维SAR成像平面的高度向再次利用孔径合成理论建立等效孔径获得第三维的高分辨能力，并与方位向、距离向的高分辨成像结合实现散射体分布的三维成像。SAR层析三维成像具有对目标三维空间直接定位、分辨的三维成像能力，是解决二维SAR成像圆柱对称模糊的最直接方法之一，大大地拓展了SAR的应用潜力，特别是在体散射结构(例如：森林、城市等区域)的分析领域。

SAR层析三维成像技术因其真正的三维成像能力，受到越来越多的关注。其中高度维的成像主要采用基于谱分析的方法(例如直接傅里叶变换方法、MUSIC，ESPRIT方法等)，首先对均匀采样、满足奈奎斯特采样定理的高度向信号进行相位校正，再通过谱分析的方法获得目标的高程和幅度，然而在实际条件下，SAR层析高度维成像将面对两个主要问题：一是基线分布不均匀，有时候采样间隔过大与奈奎斯特定理相违背而使成像结果恶化；另一个是受限于实际条件，基线数目非常少，合成孔径长度太短而导致高度维分辨率不高。

发明内容

本发明的目的是为了解决现有的SAR层析超分辨成像方法存在的上述问题，提出了一种SAR层析超分辨成像方法。

本发明详细技术方案如下：一种SAR层析超分辨成像方法，具体包括如下步骤：

步骤1：对SAR层析各基线得到的二维成像结果进行配准处理，逐像素按照基线顺序构建高度向信号，并进行调频校正，构造冗余矩阵D：

其中，

r₀为雷达轨迹与成像目标的最短距离，v_n＝n·Δv为目标成像高度坐标，Δv为高度向采样间隔，z_m为雷达所在高度；λ表示雷达发射信号的波长，M表示基线数目，N为高度向采样点数；

步骤2：初始化源向量α₀＝D⁺s，迭代次数k＝1，迭代终止误差ε，其中，α＝[b₀，b₁，…，b_i，…，b_N-1]^T为高度向成像散射系数矢量，α₀为α的初始化向量，s为观测信号；D⁺表示冗余矩阵D的广义逆矩阵，b_i表示第i+1个成像目标的高度向成像散射系数；

步骤3：更新权对角矩阵W_k满足q_k＝W_k ^-1α_k，W_k为第k次迭代结果，p的取值范围(0，1)，a_k-1[i]表示源向量α_k-1的第i个分量；

步骤4：估计使目标函数

最小的中间向量q_k，其中η是与噪声水平有关的正则化参数；

步骤5：更新源向量α_k＝W_kq_k，如果||α_k-α_k-1||≤ε，迭代终止，α_k即为最终所求高度像；否则k＝k+1，跳转至步骤3。

本发明的有益效果：本发的方法明首先对SAR层析各基线得到的二维成像结果进行配准处理，逐像素按照基线顺序构建高度向信号，并进行距离、调频校正，构造冗余矩阵；根据高度向散射系数的稀疏性特点，将高度向成像问题建模为稀疏信号的重构问题，利用加权最小范数的约束条件，通过迭代求其稀疏解，实现目标高度维成像。将该发明的方法应用于SAR层析高度维成像，解决了同一区域的二维SAR图像数量少和(/或)轨迹分布不均匀且只有很少部分的轨迹间隔(基线)满足奈奎斯特采样定理的问题，同时本发明的方法通过将观测系统及信号的先验信息引入观测信号的基函数空间，获得了非均匀基线分布的处理能力，同时由于冗余基函数具有根据观测系统和信号特性自适应调节的能力，从而使观测信号在基函数域具有超分辨特性，能使目标能量更集中，解决了基线数目少的问题。本发明的方法能够降低了成像的成本和飞行风险，拓展了SAR层析三维成像在机载平台的应用范围。

附图说明

图1为本发明的SAR层析超分辨成像方法的流程示意图。

图2为本发明的SAR层析成像几何模型示意图。

具体实施方式

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

本发明的SAR层析超分辨成像方法，如图1所示，包括以下步骤：

步骤1：对SAR层析各基线得到的二维成像结果进行配准处理，逐像素按照基线顺序构建高度向信号，并进行调频校正，构造冗余矩阵D：

其中，

r₀为雷达轨迹与成像目标的最短距离，v_n＝n·Δv为目标成像高度坐标，Δv为高度向采样间隔，z_m为雷达所在高度；λ表示雷达发射信号的波长，M表示基线数目，N为高度向采样点数；

机载SAR层析三维成像几何模型如图2所示。x，r，z分别表示方位向、斜距向以及高度向，v为垂直于x-r平面的法线方向，各条航迹间隔为d(等间距)，基线长度总和L＝(M-1)d，相当于高度向合成孔径长度。雷达轨迹与散射体的最短距离为r₀，由于r₀＞＞L/2，可以忽略z方向与v方向的微小差异，认为二者方向相同。为方便讨论，假定SAR系统的点扩展函数为标准二维狄拉克函数，则对于单点目标(r₀，v₀)，第m航迹雷达所获得的高度向回波信号可以表示为：

$S\left(z\right)=b(r_0,v_0)\exp [-j\frac{2\mathrm{\π}}{\mathrm{\λ}}2R(r_0,v_0)]---\left(1\right)$

其中，z为雷达当前高度，b(r₀，v₀)表示该点的复反射系数，λ为波长，R(r₀，v₀)为目标到雷达距离：

$R(r_0,v_0)=\sqrt{r_0^2+{(z-v_0)}^2}\≈r_0+\frac{{(z-v_0)}^2}{2r_0}---\left(2\right)$

将(2)代入(1)中，然后对其进行距离和调频校正，得：

$S^{\′}\left(z\right)=b(r_0,v_0)\exp [-j\frac{2\mathrm{\π}}{\mathrm{\λ}{r}_0}({r_0}^2-2z{v}_0)]---\left(3\right)$

对(3)式作傅里叶变换，将回波信号从空域z变换至空频域v，便可得到点(r₀，v₀)的高度像：

$I\left(v\right)=b(r_0,v_0)\exp [-j\frac{2\mathrm{\π}}{\mathrm{\λ}{r}_0}{v_0}^2)L\sin c\left[\frac{2\mathrm{\πL}}{\mathrm{\λ}{r}_0}(v-v_0)\right]---\left(4\right)$

式(1)表示的是单点目标高度向回波信号，真实情况下，各航迹二维成像结果位于(x₀，r₀)处的像素可表示为该分辨单元内高度向上所有目标的积分形式：

$S\left(z\right)=\∫b(r_0,v)\exp [-j\frac{2\mathrm{\π}}{\mathrm{\λ}}2R(r_0,v_0)]\mathrm{dv}---\left(5\right)$

对式(5)进行距离、调频校正，并且将其近似表示为过采样离散形式：

$S_m^{\′}=\underset{n=0}{\overset{N-1}{\mathrm{\Σ}}}{b}_m\·\exp [-j\frac{2\mathrm{\π}}{\mathrm{\λ}{r}_0}({v_n}^2-2z_m{v}_n)]---\left(6\right)$

其中，b_n为高度向第n个采样点的复散射系数，目标成像高度坐标v_n＝n·Δv，z_m为航迹所在高度。考虑到实际成像中系统观测噪声的存在，定义ω＝[ω₁，ω₂，…，ω_M]^T是独立同分布的零均值高斯白噪声，将式(6)表示为加噪的矩阵形式：

s＝D×α+ω    (7)

其中s＝[S′₁+ω₁，S′₂+ω₂，…，S′_M+ω_M]^T表示观测信号。

为系统观测矩阵，即冗余矩阵，其中，元素

α＝[b₀，b₁，…，b_i，…，b_N-1]^T即为高度向成像散射系数矢量。

可以看出，SAR层析高度维成像就是对式(7)的求解过程，即从观测数据s中反演散射系数矢量α。然而对于稀疏基线系统，满足M＜＜N，方程(7)的求解变成一个欠定问题，解的个数不唯一。对于高频雷达，SAR层析高度维信号可等效于少数几个散射中心反射的回波叠加形式，相当于高度维信号满足稀疏条件。因此，将稀疏性先验引入高度维聚焦，便可获得式(7)的唯一解。

步骤2：初始化源向量α₀＝D⁺s，迭代次数k＝1，迭代终止误差ε，其中，α＝[b₀，b₁，…，b_i，…，b_N-1]^T为高度向成像散射系数矢量，α₀为α的初始化向量，s为观测信号；D⁺表示冗余矩阵D的广义逆矩阵，b_i表示第i+1个成像目标的高度向成像散射系数，迭代终止误差ε的取值范围为(0，0.003)。

步骤3：更新权对角矩阵

W_k满足q_k＝W_k ^-1α_k，W_k为第k次迭代结果，p的取值范围(0，1)，α_k-1[i]表示源向量α_k-1的第i个分量；这里的p可以理解为l_p范数中的p。

步骤4：估计使目标函数

最小的中间向量q_k，其中η是与噪声水平有关的正则化参数；

步骤5：更新源向量α_k＝W_kq_k，如果||α_k-α_k-1||≤ε，迭代终止，α_k即为最终所求高度像；否则k＝k+1，跳转至步骤3。这里，ε为迭代终止误差，迭代终止误差ε的取值范围为(0，0.003)。

需要说明的是步骤2至步骤5是在加权最小范数的约束条件下，通过迭代求其稀疏解，实现目标成像。

首先采用尺度变换W^-1α＝q，将无噪声情况的式(7)转化为如下优化模型

$\left\{\begin{array}{c}\min {\left|\right|q\left|\right|}_2^2\\ s.t:\mathrm{DWq}=S\end{array}\right.---\left(8\right)$

该模型的迭代格式如下：

α_k＝W_kq_k＝W_k(DW_k)⁺S    (9)

其第k次迭代中的线性变换因子W_k为前一次迭代结果α_k-1元素组成的对角矩阵，记为W_k＝(diag(α_k-1))。

噪声环境下，在以上模型的基础上，引入l_p模(p≤1)作为问题的稀疏度量函数，利用拉格朗日乘数法将(8)转化为无约束模型：

$\min J\left(\mathrm{\α}\right)={\left|\right|\mathrm{D\α}-s\left|\right|}_2^2+\mathrm{\η}{\left|\right|\mathrm{\α}\left|\right|}_p^p---\left(10\right)$

式中，η是与噪声水平有关的正则化参数。

根据最优化理论，利用牛顿迭代法求解式(10)，整理得用于噪声环境的模型(10)的迭代格式：

α_k＝W_k ²D^H(DW_k ²D^H+ηI)^-1s    (11)

其中， $W_k=\mathrm{diag}\left\{{\left|a_{k-1}\right[i\left]\right|}^{1-\frac{p}{2}}\right\}.$

观察式(11)可知，每一次迭代同样可以等效于先针对中间变量q求其最小范数解，再由α＝Wq重构α的过程。又因为在第k步迭代中的权更新矩阵W_k是与上一次迭代结果α_k-1有关的对角矩阵，所以源向量α_k的第i个分量可以简化表示成：

a_k[i]＝W_k[i，i]·q_k[i]    i＝1，2，…，N    (12)

这里，a_k[i]对应着α＝[b₀，b₁，…，b_i，…，b_N-1]^T中的b_i-1，

可见，如果中间向量q_k的大量元素为零，即q_k满足稀疏性，那么源向量α_k也是稀疏的。因此可以利用直接求解q的稀疏估计来迫使α更快更有效的向真实稀疏源位置收敛，将式(10)的无约束模型便转化为：

$\min J\left(\mathrm{\α}\right)={\left|\right|\mathrm{DWq}-s\left|\right|}_2^2+\mathrm{\η}{\left|\right|q\left|\right|}_p^p---\left(13\right)$

式(13)便是本发明方法的无约束模型。对(13)的求解可以通过较传统牛顿法更为有效的拟牛顿算法(BFGS)来直接估计，正则化参数η的选择应该是在模型拟合度和稀疏性度量之间进行折衷，当选取

时可以有效地去除噪声的干扰，其中σ表示系统噪声的标准差。

采用机载雷达E-SAR的参数对本发明所述的成像方法进行成像仿真。设雷达单次飞行的高度为7000m，仰视角θ＝45°，工作波长0.0313m，PRF为1.2kHz，脉冲宽度8μs，线性调频信号调频率为12.5MHz/μs，方位向波束宽度17°，俯仰向波束宽度30°，基线间距为5m，共8条基线。对各基线接收的回波信号进行二维SAR成像处理后，将所有成像结果按基线顺序进行配准处理。逐像素按照基线顺序构建高度向信号，并进行调频校正。按照1m采样间隔对成像高度向空间进行采样，结合SAR层析系统的频率、高度、基线分布构建观测矩阵，最后将高度向信号表示为信号稀疏表示的形式，按照本发明中的SAR层析超分辨成像方法求解步骤即可得到该方位-距离单元内散射体的高度向成像，通过合理的选择正则化参数，此方法能够在获得良好成像结果的同时，对信号中的噪声起到很好的抑制作用。

本领域的普通技术人员将会意识到，这里所述的实施例是为了帮助读者理解本发明的原理，应被理解为本发明的保护范围并不局限于这样的特别陈述和实施例。本领域的普通技术人员可以根据本发明公开的这些技术启示做出各种不脱离本发明实质的其它各种具体变形和组合，这些变形和组合仍然在本发明的保护范围内。

Claims (2)

1.一种SAR层析超分辨成像方法，其特征在于，具体包括如下步骤：

步骤1：对SAR层析各基线得到的二维成像结果进行配准处理，逐像素按照基线顺序构建高度向信号，并进行调频校正，构造冗余矩阵D：

其中，

r₀为雷达轨迹与成像目标的最短距离，v_n=n·Δv为目标成像高度坐标，Δv为高度向采样间隔，z_m为雷达所在高度；λ表示雷达发射信号的波长，M表示基线数目，N为高度向采样点数；

步骤2：初始化源向量α₀=D⁺s，迭代次数k=1，迭代终止误差ε，其中，

α=[b₀,b₁,…,b_i,…,b_N-1]^T为高度向成像散射系数矢量，α₀为α的初始化向量，s为观测信号；

D⁺表示冗余矩阵D的广义逆矩阵，b_i表示第i+1个成像目标的高度向成像散射系数；

步骤3：更新权对角矩阵

W_k满足q_k=W_k ^-1α_k，W_k为第k次迭代结果，p的取值范围（0，1），a_k-1[i]表示源向量α_k-1的第i个分量；

步骤4：估计使目标函数

最小的中间向量q_k，其中η是与噪声水平有关的正则化参数；

步骤5：更新源向量α_k=W_kq_k，如果‖α_k-α_k-1‖≤ε，迭代终止，α_k即为最终所求高度像；否则k=k+1，跳转至步骤3。

2.根据权利要求1所述的SAR层析超分辨成像方法，其特征在于，步骤2中所述的迭代终止误差ε的取值范围为(0,0.003)。

Images