CN103439693B - 一种线阵sar稀疏重构成像与相位误差校正方法 - Google Patents

Abstract

本发明公开了一种线阵三维SAR稀疏重构成像以及模型相位误差校正方法，它是针对实际线阵三维SAR成像场景空间中散射目标稀疏的特征，建立了线阵SAR原始回波信号与场景目标空间中散射系数的线性测量矩阵，同时考虑了线阵三维SAR实际测量中的存在的相位误差因素，通过迭代优化处理估计相位误差矩阵与场景目标的散射相位向量，实现线阵SAR三维空间中稀疏目标成像处理，抑制了相位误差对成像的影响，提高了线阵SAR稀疏重构成像的稳定性和线阵SAR的成像精度。本发明可以应用于合成孔径雷达成像，地球遥感等领域。

Description

一种线阵SAR稀疏重构成像与相位误差校正方法

技术领域：

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

背景技术：

线阵三维合成孔径雷达(Linear array SAR)成像系统是将线阵天线安装在运动平台上的一种新型合成孔径雷达系统。线阵三维SAR通过发射带宽信号，利用线阵天线以及运动平台的运动合成二维虚拟阵列，从而可以实现场景目标的三维成像。与传统合成孔径雷达(SAR)系统相比，线阵三维SAR具有以下独特优点：1)三维成像能力，提供更丰富的目标信息，更有利于目标特征提取和识别；2)多模式工作能力，除了传统侧视模式，还可工作于下视模式和前视模式，克服了传统SAR成像的盲区；3)复杂地形成像，可对城市、高山、峡谷和森林等复杂地形进行三维成像。线阵三维SAR作为一项新的空间对地观测和侦察手段，无论在民用还是军事应用领域都有着更为广泛的发展空间。线阵三维SAR成像系统可用于对城市地区进行三维测绘；对复杂地形区域的地形测绘以便于进行对地观测、灾害预测与灾害评估、协助飞行器的导航及盲降；对隐藏在深峡谷的军事目标进行监视和侦察搜索，获取目标三维图像，提高目标识别能力等。线阵三维SAR已经成了目前合成孔径雷达成像技术领域的研究热点。

当今线阵三维SAR成像算法主要基于匹配滤波理论，如三维距离-多普勒(RD)算法和三维后向投影(BP)算法，见参考文献“G.Fornaro,F.Serafino,and F.Soldovieri.Three-dimensional focusing with multipass SAR data.IEEE Trans.Geosci.Remote Sens,vol.41,no.3,pp.507–517,Mar.2003.”和“Shi Jun,Zhang Xiaoling,Yang Jianyu,Wang yinbo.Surface-Tracing-Based LASAR 3-D Imaging Method viaMultiresolution Approximation.IEEE Trans.Geosci.Remote Sens,vol.46,no.11,pp.3719–3730,Nov.2008.”。根据匹配滤波理论，若平台高度、雷达入射角和雷达系统工作中心频率固定之后，线阵三维SAR系统的切航迹向分辨率取决于线阵天线长度。线阵天线长度越长，切航迹向分辨率就越高。但在实际中，由于载体(飞机或卫星)系统尺寸有限，线阵天线长度往往只有几米或者几十米，难以实现高精度成像，制约了线阵三维SAR在某些高分辨雷达成像领域的应用。如何提高成像精度成为当今线阵三维SAR技术的关键问题。

线阵三维SAR成像是从原始回波信号中重构出目标散射系数，该成像过程本质上是一个逆散射问题求解过程。如果能够建立线阵三维SAR回波信号和目标散射系数的线阵测量模型，线阵三维SAR成像问题就可等效为三维空间目标散射系数的线性方程逆求解问题。在线阵SAR成像的三维场景空间中，由于大多数区域不包含散射点(如，空气)或散射点被其他散射点遮挡而无法被入射波束照射(如，地下目标)，线阵SAR三维图像往往表现出典型的空间稀疏特征。因此线阵三维SAR成像问题可以进一步转化为稀疏目标的散射系数的估计与重构问题，在成像处理过程中只需要估计出包含稀疏目标空间单元的散射系数，并不需要估计场景目标空间所有单元的散射系数。压缩传感稀疏信号重构的主要思想是建立测量信号与被测量信号的线性测量模型，通过求解数值最优化问题准确重构出原始信号。压缩传感理论详见参考文献“D.Donoho.Compressed sensing.IEEE Trans.Inf.Theory,vol.52,no.4,pp:1289-1306,April 2006.”。压缩传感稀疏重构方法对信号测量模型的精确度要求很高，如果信号测量模型存在误差或者不精确，压缩传感稀疏重构方法的重构精确度就会严重退化，有时甚至会出现错误的结果。在线阵三维SAR回波数据实际采集过程中，由于平台运动误差、时间同步误差和大气传播效应等测量不确定性因素，不可避免的会在线阵三维SAR原始信号测量模型中引入相位误差。因此，若利用稀疏重构方法进行成像，线阵三维SAR信号测量模型不可能被精确已知。

因此，为了获得良好的线阵三维SAR稀疏重构成像效果，我们必须考虑回波信号测量模型相位误差对线阵三维SAR成像的影响，并在重构过程中对模型相位误差进行估计与补偿。根据本人了解，当前还没有出现考虑回波信号测量模型相位误差的线阵三维SAR稀疏重构成像方法。

发明内容：

为了提高线阵三维SAR成像精度，并抑制测量模型相位误差对线阵三维SAR成像的影响，本发明结合线阵三维SAR场景目标的稀疏特征，提供了一种基于线性测量模型的线阵三维SAR稀疏重构成像以及模型相位误差校正方法，该方法的思路是：利用线阵三维SAR系统参数、场景目标空间参数与原始回波信号的相互关系，建立线阵三维SAR回波信号与三维场景目标散射系数之间带有相位误差的线性测量模型，然后基于该信号线性测量模型，利用多次迭代估计方法完成测量模型相位误差估计校正和场景目标散射系数重构；在每一次迭代过程中，先利用最小均方误差准则方法对模型相位误差进行估计，然后对线性测量模型进行补偿相位误差，再利用压缩传感稀疏重构算法对三维场景空间中的稀疏目标进行成像。该方法的特点是：在成像中考虑了线阵三维SAR信号测量模型相位误差的影响，通过构建带有相位误差的线阵三维SAR回波信号线性测量模型，利用多次迭代估计和最优化方法实现对线阵三维SAR模型相位误差和稀疏目标散射系数的估计。该方法的优势是：适用于带有未知模型相位误差的线阵三维SAR稀疏成像，无需线阵三维SAR相位误差分布的先验信息，可以对任意相位误差进行估计和校正，并且在欠采样测量情况下也可以良好地估计出线阵三维SAR相位误差和稀疏目标散射系数。

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

定义1、线阵三维合成孔径雷达(Linear array SAR,LASAR)

线阵三维合成孔径雷达成像系统是将线性阵列天线固定于运动平台上并与平台运动方向垂直，结合运动平台的运动以合成二维平面阵列，并进行三维成像的一种合成孔径雷达系统。

定义2、稀疏信号

如果一个离散信号中非零值的个数远小于信号本身的长度，则该信号可认为是稀疏的。设X＝[x₁,x₂,…,x_N]^T为N个离散信号组成的列向量，其中x₁表示向量X中的第1个元素，x₂表示向量X中的第2个元素，x_N表示向量X中的第N个元素，，右上角T为转置运算符号。如果向量X中仅有K个非零元素或远大于零的元素，则X是K稀疏向量，的值定义为信号向量X的稀疏度。

定义3、稀疏重构

稀疏重构是建立测量信号的线性测量模型，并利用稀疏重构算法由测量值和测量矩阵重构出原始稀疏信号的过程。稀疏重构方法主要基于正则化理论，稀疏重构算法包括：最小范数算法、匹配追踪系列算法、最小全变分算法、迭代阈值算法等。

定义4、线阵三维SAR慢时刻与快时刻

线阵三维SAR运动平台飞过一个方位向合成孔径长度所需要的时间称为慢时间，雷达系统以一定时间长度的重复周期发射接收脉冲，因此慢时间可以表示为一个以脉冲重复周期为步长的离散化时间变量，其中每一个脉冲重复周期离散时间变量值为一个慢时刻。快时刻是指在一个脉冲重复周期内，距离向采样回波信号的时间间隔变量。

详见文献“合成孔径雷达成像原理”，皮亦鸣等编著，电子科技大学出版社出版。

定义5、线阵三维SAR场景目标空间

线阵SAR场景空间是指现实空间中所有待观测的场景目标点的集合。在不同空间坐标系下有不同的表示，但一旦坐标系确立以后其表示是唯一的。一般情况下为了方便成像，线阵SAR场景目标空间取为地面坐标系。本发明中用以下数学关系表示场景目标空间Ω：

其中和表示构成场景目标空间Ω的地表正交坐标基，分别表示水平横向、水平纵向和垂直地表的高度向。为场景目标空间的单元位置向量，x、y和z分别表示该点的水平横向、水平纵向和高度向坐标,表示实数域。

定义6、线阵三维SAR成像空间

线阵三维SAR成像空间是指将场景空间中的散射点投影到切航迹向－沿航迹向－距离向的三维空间坐标系，该空间由线阵SAR成像空间中的三个相互正交的坐标基确定。本发明中用以下数学关系表示成像空间M：

其中和表示构成成像空间M的正交坐标基，分别表示切航迹向、沿航迹向和距离向。为成像空间中的待观测点向量，u、v和z分别表示该点的切航迹向、沿航迹向和距离向坐标。

定义7、线阵三维SAR系统的传统理论分辨率

线阵三维SAR系统的传统理论分辨率是指利用基于匹配滤波理论和雷达系统模糊函数计算得到线阵三维SAR系统在距离向、方位向和切航迹向的分辨率。对于收发共用天线，线阵三维SAR距离向的分辨率记为ρ_r，近似表达式为方位向的分辨率记为ρ_a，近似表达式为切航迹向的分辨率记为ρ_c，近似表达式为详见参考文献“Shi,J.,et al.,APC trajectorydesign for One-Active linear-array three-dimensional imaging SAR,IEEETransactions on Geoscience and Remote Sensing,Vol.48,No.3,pp:1470-1486,2010”。

定义8、范数

设X是数域上线性空间，表示复数域，若它满足如下性质：||X||≥0，且||X||＝0仅有X＝0；||aX||＝|a|||X||，a为任意常数；||X₁+X₂||≤||X₁||+||X₂||，则称||X||为X空间上的范数(norm)。如果X＝[x₁,x₂,…,x_N]^T为N×1维离散信号，向量X的LP范数表达式为L1范数表达式为L2范数表达式为详见文献“矩阵理论”，黄廷祝等编著，高等教育出版社出版。

定义9、最小均方误差准则方法

若以Z表示原始信号向量，表示为Z估计信号向量，E为原始信号和估计信号向量之间的均方误差

$E={\left|\right|Z-\hat{Z}\left|\right|}_2$

均方误差最小准则方法就是选择一组估计向量使得均方误差E最小的过程。

$\hat{Z}=\arg \underset{\hat{Z}}{\min }\left(E\right)=\arg \underset{\hat{Z}}{\min }\left({\left|\right|Z-\hat{Z}\left|\right|}_2\right)$

其中，||·||₂表示向量L2范数，表示求取满足条件的自变量的最小值。

定义10、对角矩阵与单位矩阵

对角矩阵首先是行数和列数相等的方阵，同时矩阵的主对角元素不全为零，矩阵非主对角元素全为零。单位矩阵指的是主对角元素全为1，非主对角元素全为0的对角矩阵。

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

合成孔径雷达原始回波仿真方法是指基于合成孔径雷达成像原理仿真出一定系统参数条件下具有SAR回波信号特性的原始信号的方法，详细内容可参考文献：“InSAR回波信号与系统仿真研究”，张剑琦，哈尔滨工业大学硕士论文。

本发明提供的一种基于线性模型的线阵SAR稀疏重构成像与相位误差校正方法，它包括以下步骤：

步骤一、初始化线阵三维SAR系统参数：

初始化线阵三维SAR系统参数包括：平台速度矢量，记做线阵天线各阵元初始位置矢量，记做其中n为天线各阵元序号，为自然数，n＝1,2,...,N，N为线阵天线的阵元总数；线阵天线长度，记做L；相邻天线阵元之间的间距，记做d；雷达工作中心频率，记做f_c；雷达发射基带信号的信号带宽，记做B_r；雷达发射信号脉冲宽度，记做T_P；发射信号调频斜率，记做f_dr；雷达接收波门持续宽度，记做T_o；雷达接收系统的采样频率，记做f_s；雷达发射系统的脉冲重复频率，记做PRF；雷达接收系统接收波门相对于发射信号发散波门的延迟，记做T_D；天线在方位向的有效孔径长度，记做D_a；光在空气中的转播速度，记做C；距离向快时刻，记做t，t＝1,2,…,T，T为距离向快时刻总数；方位向慢时刻，记做l，l＝1,2,…,K，K为方位向慢时刻总数；上述参数均为线阵三维SAR系统标准参数，雷达系统的脉冲重复时间PRI，天线在方位向的有效孔径长度D_a在线阵三维SAR系统设计过程中已经确定；平台速度矢量及线阵天线各阵元初始位置矢量在线阵三维成像合成孔径雷达观测方案设计中已经确定；根据线阵三维SAR成像系统方案和观测方案，线阵三维SAR成像方法需要的初始化成像系统参数均为已知。

步骤二、初始化线阵三维SAR的场景目标空间参数：

初始化线阵三维SAR场景目标空间参数包括：以雷达波束照射场区域地平面和垂直于该地平面向上的单位向量所构成的空间直角坐标作为线阵三维SAR的场景目标空间Ω；将场景目标空间Ω均匀划分成大小相等的立体单元网格，单元网格在水平横向、水平纵向和高度向边长分别记为d_x、d_y和d_z，单元网格大小择为线阵三维SAR系统对应的传统理论分辨率；场景目标空间Ω中第m个单元格的坐标矢量，记做m表示场景目标空间Ω中第m个单元格，m＝1,2,…,M,M为场景目标空间Ω中单元格总数；场景目标空间Ω中单元格的散射系数向量，记做α，向量α由M行1列组成；散射系数向量α中第m个元素的散射系数，记做σ_m；场景目标空间Ω在线阵三维SAR成像方案设计中已经确定。

步骤三、建立回波信号与场景目标散射系数的线性测量矩阵：

根据步骤一中初始化的平台速度矢量线阵天线各阵元初始位置矢量和雷达系统的脉冲重复频率PRF，采用公式n＝1,2,…,N，l＝1,2,…,K，计算得到第n个线阵天线阵元在第l个方位向慢时刻的位置矢量，记为其中N为步骤一中线阵天线阵元总数，K为步骤一的方位向慢时刻总数；采用公式 $R({\stackrel{\‾}{P}}_n\left(l\right),{\stackrel{\‾}{P}}_m)={\left|\right|{\stackrel{\‾}{P}}_n\left(l\right)-{\stackrel{\‾}{P}}_m\left|\right|}_2,n=\mathrm{1,2},...,N,$ $l=\mathrm{1,2},...,K,m=\mathrm{1,2},...,M,$ 计算得到在第l个方位向慢时刻线阵三维SAR场景目标空间Ω中第m个单元格到第n个线阵天线阵元的距离，记为其中||·||₂表示向量L2范数，为步骤二中初始化得到场景目标空间Ω中第m个单元格的坐标矢量，M为步骤二中初始化的场景目标空间Ω中单元格总数；采用公式 ${\mathrm{\τ}}_{\mathrm{nm}}\left(l\right)=2\·R({\stackrel{\‾}{P}}_n\left(l\right),{\stackrel{\‾}{P}}_m)/C,n=\mathrm{1,2},...,N,l=\mathrm{1,2},...,K,m=\mathrm{1,2},...,M,$ 计算得到在第l个方位向慢时刻线阵三维SAR场景目标空间Ω中第m个单元格到第n个线阵天线阵元的时间延时，记为τ_nm(l)，其中C为步骤一中初始化得到的光在空气中的传播速度；在第l个方位向慢时刻和第t个距离向快时刻中线阵三维SAR第n个线阵天线阵元的原始回波数据记为s(t,l,n)，t＝1,2,…,T，l＝1,2,…,K，n＝1,2,…,N，其中T为步骤一中初始化的距离向快时刻总数；在线阵三维SAR实际观测中，s(t,l,n)可由数据接收机得到；若在仿真过程中，s(t,l,n)为场景目标空间Ω中所有单元格回波的总数，采用合成孔径雷达原始回波仿真方法产生得到，近似表示公式可以表示为 $s(t,l,n)=\underset{m=1}{\overset{M}{\mathrm{\Σ}}}{\mathrm{\σ}}_m\·\exp [-j\·2\·\mathrm{\π}\·f_c\·{\mathrm{\τ}}_{\mathrm{nm}}\left(l\right)]\exp \{j\·\mathrm{\π}\·f_{\mathrm{dr}}\·{[t-{\mathrm{\τ}}_{\mathrm{nm}}\left(l\right)]}^2\},$ 其中Σ(·)表示求和运算符号，exp(·)表示e指数运算符号，f_c为步骤一初始化得到的雷达中心频率，f_dr为步骤一初始化得到的发射信号调频斜率，σ_m为步骤二初始化得到的场景目标空间中第m个单元格的散射系数，j为虚数单位，π为圆周率；将所有线阵三维SAR原始回波信号s(t,l,n)按顺序排列组成向量，记为回波信号向量S，回波信号向量S由X行1列组成，其中X＝T·K·N；

采用公式 ${\mathrm{\φ}}_i\left(m\right)=\exp [-j\·2\·\mathrm{\π}\·f_c\·{\mathrm{\τ}}_{\mathrm{nm}}\left(l\right)]\exp \{j\·\mathrm{\π}\·f_{\mathrm{dr}}\·{[t-{\mathrm{\τ}}_{\mathrm{nm}}\left(l\right)]}^2\},$ $t=\mathrm{1,2},...,T,l=\mathrm{1,2},...,K,n=\mathrm{1,2},...,N,m=\mathrm{1,2},...,M,i=[(t-1)K+l-1]N+n,$ 计算得到场景目标空间Ω中第m个单元格在回波信号向量S第i个元素信号对应的时延函数，记为φ_i(m)，其中i的取值范围为i＝1,2,…,X；令矩阵A为线阵三维SAR原始回波信号与场景目标空间所有单元格的线性测量矩阵，测量矩阵A由场景目标空间Ω中所有单元格的时延函数构成，具体表达式为

其中，φ₁(1)为场景目标空间Ω中第1个单元格在回波信号向量S第1个元素信号对应的时延函数，φ₁(2)为场景目标空间Ω中第2个单元格在回波信号向量S第1个元素信号对应的时延函数，φ₁(M)为场景目标空间Ω中第M个单元格在回波信号向量S第1个元素信号对应的时延函数，φ₂(1)为场景目标空间Ω中第1个单元格在回波信号向量S第2个元素信号对应的时延函数，φ₂(2)为场景目标空间Ω中第2个单元格在回波信号向量S第2个元素信号对应的时延函数，φ₂(M)为场景目标空间Ω中第M个单元格在回波信号向量S第2个元素信号对应的时延函数，φ_X(1)为场景目标空间Ω中第1个单元格在回波信号向量S第X个元素信号对应的时延函数，φ_X(2)为场景目标空间Ω中第2个单元格在回波信号向量S第X个元素信号对应的时延函数，φ_X(M)为场景目标空间Ω中第M个单元格在回波信号向量S第X个元素信号对应的时延函数，φ₁(1),φ₁(2),…,φ₁(M)分别为场景目标空间Ω中第1,2,…,M个单元格在回波信号向量S第1个元素信号对应的时延函数向量，φ₂(1),φ₂(2),…,φ₂(M)分别为场景目标空间Ω中第1,2,…,M个单元格在回波信号向量S第2个元素信号对应的时延函数向量，φ_X(1),φ_X(2),…,φ_X(M)分别为场景目标空间Ω中第1,2,…,M个单元格在回波信号向量S第X个元素信号对应的时延函数向量；线阵三维SAR的线性测量矩阵A为X行M列的二维矩阵。

步骤四、设定稀疏重构算法的初始参数：

初始化稀疏重构算法的参数包括：稀疏重构处理的最大迭代次数，记做MaxIter；稀疏过程处理中的迭代终止阈值，记做ε；线阵三维SAR线阵测量模型中的相位误差矩阵，记为R，相位误差矩阵R为X行X列的对角二维矩阵，即矩阵R非对角元素全部为0；相位误差矩阵R的初始迭代值记为R⁽⁰⁾，R⁽⁰⁾选择为单位矩阵，即R⁽⁰⁾＝I_X，其中矩阵I_X表示为X行X列的二维单位矩阵；场景目标散射系数向量α的初始迭代值记为α⁽⁰⁾，α⁽⁰⁾的值选择为α⁽⁰⁾＝A^HS，其中A为步骤三中得到的线阵三维SAR线性测量矩阵，S为步骤三中得到的线阵三维SAR回波信号向量；k表示稀疏重构算法中的第k迭代次数，k的初始值设置为k＝0。

步骤五、模型相位误差矩阵估计：

在稀疏重构算法对线阵三维SAR原始回波信号处理过程中，场景目标散射系数向量α的第k次迭代值记为α^(k)，相位误差矩阵R的第k次迭代值记为R^(k)，其中k为步骤四中稀疏重构算法的第k迭代次数；令迭代次数k加1，若k＝1，α^(k-1)的值为步骤四中得到的初始迭代值α⁽⁰⁾，否则α^(k-1)通过步骤六目标散射系数估计得到；采用公式计算得到相位误差矩阵R的第k次迭代值R^(k)，其中表示求取满足条件的自变量矩阵R的最小值，||·||₂为向量L2范数，A为步骤三中得到的线阵三维SAR线性测量矩阵，S为步骤三中得到的线阵三维SAR回波信号向量。

步骤六、目标散射系数稀疏重构：

若k＝1，R^(k-1)的值为步骤四中得到的初始迭代值R⁽⁰⁾，否则R^(k-1)通过步骤五模型相位误差矩阵估计得到；采用公式计算得到场景目标散射系数向量α的第k次迭代值α^(k)，其中表示求取满足条件的自变量矩阵α的最小值，||·||₂为向量L1范数，||·||₂为向量L2范数，A为步骤三中得到的线阵三维SAR线性测量矩阵，S为步骤三中得到的线阵三维SAR回波信号向量，R^(k)为步骤五中得到的相位误差矩阵的第k次迭代值。

步骤七、迭代判定与迭代终止：

如果同时满足k≤MaxIter与两个条件，则重复步骤五和步骤六，其中k为稀疏重构算法的第k次迭代次数，MaxIter为步骤四中初始化得到的稀疏重构处理的最大迭代次数，ε为步骤四初始化得到的稀疏重构处理的迭代终止阈值，α^(k)为步骤五中得到的场景目标散射系数向量α的第k次迭代值，α^(k-1)为步骤五中得到的场景目标散射系数向量α的第k-1次迭代值；如果不满足k≤MaxIter与任意一个条件，则稀疏重构算法停止迭代，第k-1次迭代值α^(k-1)即为最终的场景目标散射系数向量，得到对应场景目标空间Ω中所有单元格的散射系数值，第k-1次迭代值R^(k-1)即为最终的模型相位误差矩阵，完成基于线阵测量模型的线阵三维SAR系数成像以及相位误差校正方法；最后将散射系数向量α^(k-1)转换成三维矩阵形式，得到线阵SAR场景目标空间Ω的三维成像结果。

本发明的创新点在于针对实际线阵三维SAR成像场景空间中散射目标稀疏的特征，建立了线阵SAR原始回波信号与场景目标空间中散射系数的线性测量矩阵，将线阵三维SAR成像问题转换为稀疏目标信号的稀疏重构问题，同时考虑了线阵三维SAR实际测量中的存在的相位误差因素，建立了线阵三维SAR测量模型的相位误差模型，最后利用线阵三维SAR线性测量模型，通过迭代优化处理估计相位误差矩阵与场景目标的散射相位向量，实现线阵SAR三维空间中稀疏目标成像处理，抑制了相位误差对成像的影响，提高了线阵SAR稀疏重构成像的稳定性，并提高了线阵SAR的成像精度。

本发明的优点在于适用于带有未知模型相位误差的线阵三维SAR稀疏成像，在成像处理过程中线阵三维SAR稀疏成像和误差估计同时完成，抑制模型相位误差对线阵三维SAR稀疏成像的影响，提高了线阵三维SAR的稀疏重构成像性能；只需要建立线阵三维SAR回波信号与场景目标空间的线阵测量矩阵，不需要相位误差分布的先验信息，在稀疏测量样本情况下也适用。本发明可以应用于合成孔径雷达成像，地球遥感等领域。

附图说明：

图1为线阵三维SAR成像几何关系图。

其中，线阵天线长度为L，PRI表示线阵SAR发射信号的脉冲重复时间，为线阵天线中第n个阵元在方位向第l个慢时刻的位置矢量，x_n、y_n和z_n分别表示线阵天线中第n个阵元的水平横向、水平纵向和高度向坐标；表示场景中的第m空间单元格的位置向量，为在慢时刻l时第m个单元格到线阵天线第n阵元的距离，x、y和z分别表示该点的水平横向、水平纵向和高度向坐标。

图2为本发明所提供方法的流程示意框图。

图3为本发明具体实施方式采用的线阵三维SAR系统仿真参数表。

图4为通过本发明提供方法得到的9点目标线阵三维SAR成像结果。

原始仿真场景目标为9个点，散射系数都为1，它们坐标位置分别为[0,0,0]m，[2,2,2]m，[2,-2,2]m，[-2,2,2]m，[-2,-2,2]m，[2,2,-2]m，[2,-2,-2]m，[-2,2,-2]m，[-2,-2,-2]m。图4中为空间中9个散射点的仿真成像结果，分别用字母1,2,3,4,5,6,7,8,9表示，其中横坐标为切航迹方向，纵坐标为沿航迹方向，垂直坐标为高度向。从图4中可以看出，本发明提供的方法可以实现线阵三维SAR仿真中9个点目标的稀疏重构成像处理，并且在成像结果没有出现高旁瓣。

具体实施方式

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

步骤一、选择仿真所需的线阵三维SAR系统参数，其值如图3所示，利用计算机仿真产生线阵三维SAR的原始回波数据。

初始化线阵三维SAR成像系统参数包括：运动平台速度矢量线阵天线的阵元总数N＝201，线阵天线各个阵元的初始位置矢量其中n为第n个天线阵元序号，n＝1,2,…,N，N＝201，线阵天线长度L＝3m，线阵天线相邻阵元之间的间距d＝0.015m，雷达中心频率f_c＝10GHz,雷达发射机基带信号的信号宽度B_r＝300MHz，雷达发射信号脉冲宽度T_P＝10^-6s，雷达发射信号调频斜率f_dr＝3×10¹⁴Hz/s，雷达接收系统的采样频率f_s＝500MHz，雷达系统发射信号的脉冲重复频率PRF＝500Hz，发射信号脉冲重复时PFI＝1/PRF＝2ms，线阵天线在方位向的有效孔径长度D_a＝1m，光在空气中的传播速度C＝3×10⁸m/s，距离向快时刻总数T＝500，距离向快时刻序列t＝1,2,…,500，方位向慢时刻总数K＝150，方位向慢时刻序列l＝1,2,…,150。

步骤二、初始化线阵三维SAR的场景目标空间参数：

以雷达波束照射场区域地平面和垂直于该地平面向上的单位向量所构成的空间直角坐标作为线阵三维SAR的场景目标空间Ω。初始化场景目标空间Ω的大小为128×128×128像素，场景目标空间Ω的中心坐标位置位于[0,0,0]，每一个单元网格在水平横向、水平纵向和高度向边长为d_x＝d_y＝d_z＝0.5m，计算可以得到场景目标空间的总单元格数M＝2097152，场景目标空间Ω中每一个单元格的位置为 ${\stackrel{\‾}{P}}_m=[(x^{\′}-64)\·0.5,(y^{\′}-64)\·0.5,(z^{\′}-64)\·0.5]m,$ 其中x′＝1,2,…,128，y′＝1,2,…,128，z′＝1,2,…,128，m＝[(x′-1)·128+y′-1]·128+z′。为场景目标空间Ω中第m个单元格的位置矢量，m表示场景目标空间Ω中第m个单元格，m＝1,2,…,M，M＝2097152。在场景目标空间Ω中加入仿真点目标，点目标数个数为9个，它们散射系数值均为1，坐标位置分别为[0,0,0]、[2,2,2]、[2,-2,2]、[-2,2,2]、[-2,-2,2]、[2,2,-2]、[2,-2,-2]、[-2,2,-2]、[-2,-2,-2]，单位均为m；场景目标空间Ω中没有包含点目标单元格的散射系数设置为0。将场景目标空间Ω中所有单元格的目标散射系数顺序排列组成散射向量α，散射系数向量α在线阵SAR三维成像观测方案仿真中已经确定。场景目标散射系数向量α在仿真过程中为已知参数，向量α由M行1列组成，σ_m为向量α中对应场景目标空间Ω中第m个单元格的散射系数值。在本仿真场景目标空间中，只有包含点目标的9个单元格散射系数值σ设置为1，其余单元格的散射系数都为0。

步骤三、建立回波信号与场景目标散射系数的线性观测矩阵与相位误差矩阵：

采用公式计算得到第n个线阵天线阵元在第l个方位向慢时刻的位置矢量其中n表示天线阵元序号n＝1,2,…,N，N＝201，l表示方位向慢时刻序号l＝1,2,…,K，K＝150，运动平台初始位置矢量 ${\stackrel{\‾}{P}}_n\left(0\right)=[0.015\·(n-101),\mathrm{0,1000}]m,$ 运动平台速度矢量 $\stackrel{\‾}{V}=\left[\mathrm{0,100,0}\right]m/s,$ 脉冲重复频率PRF＝500Hz。采用公式计算得到在第l个方位向慢时刻线阵三维SAR场景目标空间Ω中第m个单元格到第n个线阵天线阵元的距离其中n＝1,2,…,N，N＝201，l＝1,2,…,K，K＝150，m表示场景目标空间Ω中第m个单元格，m＝1,2,…,M，M＝2097152，||·||₂表示定义8中的向量L2范数，为步骤二中初始化得到场景目标空间Ω中第m个单元格的坐标矢量。采用公式计算得到在第l个方位向慢时刻线阵三维SAR场景目标空间Ω中第m个单元格到第n个线阵天线阵元的时间延时τ_nm(l)，其中n＝1,2,…,N，l＝1,2,…,K，m＝1,2,…,M，C＝3×10⁸m/s。采用公式和合成孔径雷达原始回波仿真方法得到场景目标空间Ω的近似线阵三维SAR原始信号回波s(t,l,n)，其中t表示距离向快时刻序号，t＝1,2,…,T，T＝500，l＝1,2,…,K，n＝1,2,…,N，Σ(·)表示求和运算符号，exp(·)表示e指数运算符号，雷达中心频率f_c＝10GHz，雷达发射信号调频斜率f_dr＝3×10¹⁴Hz/s，σ_m为步骤二初始化得到的场景目标空间中第m个单元格的散射系数，j为虚数单位(即-1的开根值)，π＝3.14159。在原始回波信号s(t,l,n)的每一个数据中加入均值为0，误差范围从-3到3的随机相位误差。将所有线阵三维SAR原始回波信号s(t,l,n)按顺序排列组成回波信号向量S，回波信号向量S由X行1列组成，其中X＝T·K·N＝15075000。采用公式φ_i(m)＝exp[-j·2·π·f_c·τ_nm(l)]exp{j·π·f_dr·[t-τ_nm(l)]²}，t＝1,2,…,T，T＝500，l＝1,2,…,K，K＝150，n＝1,2,…,N，N＝201，m＝1,2,…,M，M＝2097152，i＝[(t-1)K+l-1]N+n，计算得到场景目标空间Ω中第m个单元格在回波信号向量S第i个元素信号对应的时延函数φ_i(m)，其中i的取值范围为i＝1,2,…,X，X＝15075000。

采用矩阵表达公式

计算得到线阵三维SAR原始回波信号与场景目标空间所有单元格的线性测量矩阵A，线性测量矩阵A为X行M列的二维矩阵。

步骤四、设定重构算法的初始参数：

稀疏重构处理的最大迭代次数MaxIter＝100，稀疏过程处理中的迭代终止阈值ε＝0.001；线阵三维SAR线阵测量模型中的相位误差矩阵，记为R，相位误差矩阵R为X行X列的对角二维矩阵，X＝15075000，即矩阵R非对角元素全部为0。相位误差矩阵R的初始迭代值记为R⁽⁰⁾，一般R⁽⁰⁾选择为单位矩阵，即R⁽⁰⁾＝I_X，其中矩阵I_X表示为X行X列的二维单位矩阵。场景目标散射系数向量α的初始迭代值记为α⁽⁰⁾，一般α⁽⁰⁾的值选择为α⁽⁰⁾＝A^HS，其中A为步骤三中得到的线阵三维SAR线性测量矩阵，S为步骤三中得到的线阵三维SAR回波信号向量。k表示稀疏重构算法中的第k迭代次数，k的初始值设置为k＝0。

步骤五、模型相位误差矩阵估计：

在稀疏重构算法对线阵三维SAR原始回波信号处理过程中，场景目标散射系数向量α的第k次迭代值记为α^(k)，相位误差矩阵R的第k次迭代值记为R^(k)，其中k为步骤四中稀疏重构算法的第k迭代次数。令迭代次数k加1，若k＝1，α^(k-1)的值为步骤四中得到的初始迭代值α⁽⁰⁾，否则α^(k-1)通过步骤六目标散射系数估计得到。采用公式计算得到相位误差矩阵R的第k次迭代值R^(k)，其中表示求取满足条件的自变量矩阵R的最小值，||·||₂为向量L2范数，A为步骤三中得到的线阵三维SAR线性测量矩阵，S为步骤三中得到的线阵三维SAR回波信号向量。

步骤六、目标散射系数稀疏重构：

若k＝1，R^(k-1)的值为步骤四中得到的初始迭代值R⁽⁰⁾，否则R^(k-1)通过步骤五模型相位误差矩阵估计得到。采用公式计算得到场景目标散射系数向量α的第k次迭代值α^(k)，其中表示求取满足条件的自变量矩阵α的最小值，||·||₂为向量L1范数，||·||₂为向量L2范数，A为步骤三中得到的线阵三维SAR线性测量矩阵，S为步骤三中得到的线阵三维SAR回波信号向量，R^(k)为步骤五中得到的相位误差矩阵的第k次迭代值。

步骤七、迭代判定与迭代终止：

对迭代次数k和重构误差进行判定，如果同时满足k≤MaxIter与两个条件，其中，||·||₂为向量L2范数，MaxIter＝100，ε＝0.001，则重复步骤五和步骤六，α^(k)为步骤六中得到的场景目标散射系数向量α的第k次迭代值，如果迭代次数k大于1，α^(k-1)为步骤六中得到的场景目标散射系数向量α的第k-1次迭代值，如果k等于1，α^(k-1)为步骤四中初始化得到场景目标散射系数向量α⁽⁰⁾；如果不满足k≤MaxIter与任意一个条件，则稀疏重构算法停止迭代，此刻稀疏重构算法第k-1次迭代值α^(k-1)即为最终的场景目标散射系数向量，得到对应场景目标空间Ω中所有单元格的散射系数值，第k-1次迭代值R^(k-1)即为最终的模型相位误差矩阵，完成基于线阵测量模型的线阵三维SAR系数成像以及相位误差校正方法。最后将最终得到的散射系数向量α^(k-1)转换成128×128×128三维矩阵形式，得到线阵三维SAR场景目标空间Ω的三维成像结果。

通过本发明具体实施方式可以看出，本发明通过建立线阵三维SAR原始回波信号与场景目标空间散射系数的线性测量模型，并结合线阵三维SAR场景目标空间的稀疏特征，将线阵三维SAR成像处理过程转换成为线性方程组的稀疏求解过程。另外，在构造线性测量模型时，本发明还考虑了实际线阵三维SAR系统测量不精确引起的模型相位误差对稀疏重构成像的影响，并基于线阵三维SAR信号的线性测量模型实现了目标散射系数稀疏成像和相位误差校正，提高了线阵三维SAR切航迹向成像精度。本发明提供了一种基于线性测量模型的线阵三维SAR稀疏成像和相位误差校正方法，该方法结合线阵三维SAR系统参数和目标场景参数，建立回波信号与稀疏目标散射系数的线性测量模型以及相位误差模型，然后基于信号观测模型和相位误差模型，利用最小均方误差准则方法对模型相位误差进行估计与校正，利用稀疏重构算法对三维场景空间中的稀疏目标进行成像。

Claims (1)

1.一种线阵SAR稀疏重构成像与相位误差校正方法，其特征是该方法包括如下步骤：

步骤1、初始化线阵三维SAR系统参数：

初始化线阵三维SAR系统参数包括：平台速度矢量，记做线阵天线各阵元初始位置矢量，记做其中n为天线各阵元序号，为自然数，n＝1,2,...,N，N为线阵天线的阵元总数；线阵天线长度，记做L；相邻天线阵元之间的间距，记做d；雷达工作中心频率，记做f_c；雷达发射基带信号的信号带宽，记做B_r；雷达发射信号脉冲宽度，记做T_P；发射信号调频斜率，记做f_dr；雷达接收波门持续宽度，记做T_o；雷达接收系统的采样频率，记做f_s；雷达发射系统的脉冲重复频率，记做PRF；雷达接收系统接收波门相对于发射信号发散波门的延迟，记做T_D；天线在方位向的有效孔径长度，记做D_a；光在空气中的转播速度，记做C；距离向快时刻，记做t，t＝1,2,…,T，T为距离向快时刻总数；方位向慢时刻，记做l，l＝1,2,…,K，K为方位向慢时刻总数；上述参数均为线阵三维SAR系统标准参数，雷达系统的脉冲重复时间PRI，天线在方位向的有效孔径长度D_a在线阵三维SAR系统设计过程中已经确定；平台速度矢量V及线阵天线各阵元初始位置矢量在线阵三维成像合成孔径雷达观测方案设计中已经确定；根据线阵三维SAR成像系统方案和观测方案，线阵三维SAR成像方法需要的初始化成像系统参数均为已知；

步骤2、初始化线阵三维SAR的场景目标空间参数：

初始化线阵三维SAR场景目标空间参数包括：以雷达波束照射场区域地平面和垂直于该地平面向上的单位向量所构成的空间直角坐标作为线阵三维SAR的场景目标空间Ω；将场景目标空间Ω均匀划分成大小相等的立体单元网格，单元网格在水平横向、水平纵向和高度向边长分别记为d_x、d_y和d_z，单元网格大小选择为线阵三维SAR系统对应的传统理论分辨率；场景目标空间Ω中第m个单元格的坐标矢量，记做m表示场景目标空间Ω中第m个单元格，m＝1,2,…,M,M为场景目标空间Ω中单元格总数；场景目标空间Ω中单元格的散射系数向量，记做α，向量α由M行1列组成；散射系数向量α中第m个元素的散射系数，记做σ_m；场景目标空间Ω在线阵三维SAR成像方案设计中已经确定；

步骤3、建立回波信号与场景目标散射系数的线性测量矩阵：

根据步骤1中初始化的平台速度矢量线阵天线各阵元初始位置矢量和雷达系 统的脉冲重复频率PRF，采用公式n＝1,2,…,N，l＝1,2,…,K，计算得到第n个线阵天线阵元在第l个方位向慢时刻的位置矢量，记为其中N为步骤一中线阵天线阵元总数，K为步骤1的方位向慢时刻总数；采用公式n＝1,2,…,N，l＝1,2,…,K，m＝1,2,…,M，计算得到在第l个方位向慢时刻线阵三维SAR场景目标空间Ω中第m个单元格到第n个线阵天线阵元的距离，记为其中||·||₂表示向量L2范数，为步骤2中初始化得到场景目标空间Ω中第m个单元格的坐标矢量，M为步骤2中初始化的场景目标空间Ω中单元格总数；采用公式n＝1,2,…,N，l＝1,2,…,K，m＝1,2,…,M，计算得到在第l个方位向慢时刻线阵三维SAR场景目标空间Ω中第m个单元格到第n个线阵天线阵元的时间延时，记为τ_nm(l)，其中C为步骤1中初始化得到的光在空气中的传播速度；在第l个方位向慢时刻和第t个距离向快时刻中线阵三维SAR第n个线阵天线阵元的原始回波数据记为s(t,l,n)，t＝1,2,…,T，l＝1,2,…,K，n＝1,2,…,N，其中T为步骤1中初始化的距离向快时刻总数；在线阵三维SAR实际观测中，s(t,l,n)可由数据接收机得到；若在仿真过程中，s(t,l,n)为场景目标空间Ω中所有单元格回波的总数，采用合成孔径雷达原始回波仿真方法产生得到，近似表示公式可以表示为其中Σ(·)表示求和运算符号，exp(·)表示e指数运算符号，f_c为步骤1初始化得到的雷达中心频率，f_dr为步骤1初始化得到的发射信号调频斜率，σ_m为步骤2初始化得到的场景目标空间中第m个单元格的散射系数，j为虚数单位，π为圆周率；将所有线阵三维SAR原始回波信号s(t,l,n)按顺序排列组成向量，记为回波信号向量S，回波信号向量S由X行1列组成，其中X＝T·K·N；采用公式φ_i(m)＝exp[-j·2·π·f_c·τ_nm(l)]exp{j·π·f_dr·[t-τ_nm(l)]²}，t＝1,2,…,T，l＝1,2,…,K，n＝1,2,…,N，m＝1,2,…,M，i＝[(t-1)K+l-1]N+n，计算得到场景目标空间Ω中第m个单元格在回波信号向量S第i个元素信号对应的时延函数，记为φ_i(m)，其中i的取值范围为i＝1,2,…,X；令矩阵A为线阵三维SAR原始回波信号与场景目标空间所有单元 格的线性测量矩阵，测量矩阵A由场景目标空间Ω中所有单元格的时延函数构成，具体表达式为

其中，φ₁(1)为场景目标空间Ω中第1个单元格在回波信号向量S第1个元素信号对应的时延函数，φ₁(2)为场景目标空间Ω中第2个单元格在回波信号向量S第1个元素信号对应的时延函数，φ₁(M)为场景目标空间Ω中第M个单元格在回波信号向量S第1个元素信号对应的时延函数，φ₂(1)为场景目标空间Ω中第1个单元格在回波信号向量S第2个元素信号对应的时延函数，φ₂(2)为场景目标空间Ω中第2个单元格在回波信号向量S第2个元素信号对应的时延函数，φ₂(M)为场景目标空间Ω中第M个单元格在回波信号向量S第2个元素信号对应的时延函数，φ_X(1)为场景目标空间Ω中第1个单元格在回波信号向量S第X个元素信号对应的时延函数，φ_X(2)为场景目标空间Ω中第2个单元格在回波信号向量S第X个元素信号对应的时延函数，φ_X(M)为场景目标空间Ω中第M个单元格在回波信号向量S第X个元素信号对应的时延函数，φ₁(1),φ₁(2),…,φ₁(M)分别为场景目标空间Ω中第1,2,…,M个单元格在回波信号向量S第1个元素信号对应的时延函数向量，φ₂(1),φ₂(2),…,φ₂(M)分别为场景目标空间Ω中第1,2,…,M个单元格在回波信号向量S第2个元素信号对应的时延函数向量，φ_X(1),φ_X(2),…,φ_X(M)分别为场景目标空间Ω中第1,2,…,M个单元格在回波信号向量S第X个元素信号对应的时延函数向量；线阵三维SAR的线性测量矩阵A为X行M列的二维矩阵；

步骤4、设定稀疏重构算法的初始参数：

初始化稀疏重构算法的参数包括：稀疏重构处理的最大迭代次数，记做MaxIter；稀疏过程处理中的迭代终止阈值，记做ε；线阵三维SAR线阵测量模型中的相位误差矩阵，记为R，相位误差矩阵R为X行X列的对角二维矩阵，即矩阵R非对角元素全部为0；相位误差矩阵R的初始迭代值记为R⁽⁰⁾，R⁽⁰⁾选择为单位矩阵，即R⁽⁰⁾＝I_X，其中矩阵I_X表示为X行X列的二维单位矩阵；场景目标散射系数向量α的初始迭代值记为α⁽⁰⁾，α⁽⁰⁾的值选择为 α⁽⁰⁾＝A^HS，其中A为步骤3中得到的线阵三维SAR线性测量矩阵，S为步骤3中得到的线阵三维SAR回波信号向量；k表示稀疏重构算法中的第k迭代次数，k的初始值设置为k＝0；

步骤5、模型相位误差矩阵估计：

在稀疏重构算法对线阵三维SAR原始回波信号处理过程中，场景目标散射系数向量α的第k次迭代值记为α^(k)，相位误差矩阵R的第k次迭代值记为R^(k)，其中k为步骤4中稀疏重构算法的第k迭代次数；令迭代次数k加1，若k＝1，α^(k-1)的值为步骤4中得到的初始迭代值α⁽⁰⁾，否则α^(k-1)通过步骤6目标散射系数估计得到；采用公式计算得到相位误差矩阵R的第k次迭代值R^(k)，其中表示求取满足条件的自变量矩阵R的最小值，||·||₂为向量L2范数，A为步骤3中得到的线阵三维SAR线性测量矩阵，S为步骤3中得到的线阵三维SAR回波信号向量；

步骤6、目标散射系数稀疏重构：

若k＝1，R^(k-1)的值为步骤4中得到的初始迭代值R⁽⁰⁾，否则R^(k-1)通过步骤5模型相位误差矩阵估计得到；采用公式计算得到场景目标散射系数向量α的第k次迭代值α^(k)，其中表示求取满足条件的自变量矩阵α的最小值，||·||₂为向量L1范数，||·||₂为向量L2范数，A为步骤3中得到的线阵三维SAR线性测量矩阵，S为步骤3中得到的线阵三维SAR回波信号向量，R^(k)为步骤5中得到的相位误差矩阵的第k次迭代值；

步骤7、迭代判定与迭代终止：

如果同时满足k≤MaxIter与两个条件，则重复步骤5和步骤6，其中k为稀疏重构算法的第k次迭代次数，MaxIter为步骤4中初始化得到的稀疏重构处理的最大迭代次数，ε为步骤4初始化得到的稀疏重构处理的迭代终止阈值，α^(k)为步骤5中得到的场景目标散射系数向量α的第k次迭代值，α^(k-1)为步骤5中得到的场景目标散射系数向量α的第k-1次迭代值；如果不满足k≤MaxIter与任意一个条件，则稀疏重构算法停 止迭代，第k-1次迭代值α^(k-1)即为最终的场景目标散射系数向量，得到对应场景目标空间Ω中所有单元格的散射系数值，第k-1次迭代值R^(k-1)即为最终的模型相位误差矩阵，完成基于线阵测量模型的线阵三维SAR系数成像以及相位误差校正方法；最后将散射系数向量α^(k-1)转换成三维矩阵形式，得到线阵SAR场景目标空间Ω的三维成像结果。