[发明专利]基于硬阈值OMP的线阵SAR稀疏成像方法

权利要求书

1.一种基于硬阈值OMP的线阵SAR稀疏成像方法，其特征是它包括以下步骤：

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

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

根据步骤1中初始化的平台速度矢量线阵天线各阵元初始位置矢量和雷达系统的脉冲重复频率PRF，采用公式P‾n(l)=P‾n(0)+V‾·l/PRF,n=1,2,···,N,l=1,2,···,K,]]>计算得到第n个线阵天线阵元在第l个方位向慢时刻的位置矢量，记为其中N为步骤1中线阵天线阵元总数，K为步骤1的方位向慢时刻总数；采用公式n＝1,2,…,N，l＝1,2,…,K，m＝1,2,…,M，计算得到在第l个方位向慢时刻线阵SAR观测场景目标空间Ω中第m个单元格到第n个线阵天线阵元的距离，记为其中||·||₂表示定义2中的向量L2范数，为步骤2中初始化得到观测场景目标空间Ω中第m个单元格的坐标矢量，M为步骤2中初始化的场景目标空间Ω中单元格总数；采用公式τnm(l)=2·R(P‾n(l),P‾m)/C,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)为观测场景目标空间Ω中所有单元格回波的总数，采用定义12中传统的合成孔径雷达原始回波仿真方法产生提供，近似表示公式可表示为s(t,l,n)=Σm=1Mαm·exp[-j·2·π·fc·τnm(l)]exp{j·π·fdr·[t-τnm(l)]2},]]>其中Σ(·)表示求和运算符号，exp(·)表示e指数运算符号，f_c为步骤1初始化得到的雷达工作中心频率，f_dr为步骤1初始化得到的发射信号调频斜率，α_m为步骤2初始化得到的场景目标空间中第m个单元格的散射系数，t为距离向的第t个快时刻，j为虚数单位，即-1的开根值，π为圆周率；将所有线阵SAR原始回波信号s(t,l,n)按顺序排列组成向量，记为回波信号向量S，回波信号向量S由O行1列组成，其中O＝T·K·N,T为步骤1中初始化的距离向快时刻总数，K为步骤1初始化的方位向慢时刻总数，N为步骤1初始化的线阵天线的阵元总数；采用公式φ_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＝1,2,…,O，计算得到观测场景目标空间Ω中第m个单元格在回波信号向量S第i个元素信号对应的时延函数，记为φ_i(m)；令矩阵A为线阵SAR原始回波信号向量S与观测场景目标空间Ω中所有单元格散射系数向量α之间的测量矩阵，测量矩阵A由线阵SAR观测场景目标空间Ω中所有单元格对应的时延函数构成，具体表达式为

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

步骤4、设定硬阈值正交匹配追踪(简称OMP)算法的初始参数：

初始化硬阈值OMP算法的参数包括：硬阈值OMP算法重构迭代处理的最大迭代次数，记做MaxIter；目标散射系数最大最小对比度门限，记做η₀；目标散射射系数变化率门限，记做β₀；重构残余误差门限，记做ε₀；观测场景目标散射系数向量α的初始迭代值记为α⁽⁰⁾，一般α⁽⁰⁾的值选择为α⁽⁰⁾＝0或者α⁽⁰⁾＝A^HS，其中A为步骤3中得到的线阵SAR测量矩阵，S为步骤3中得到的线阵SAR原始回波信号向量，上标H表示共轭转置运算符号；重构残余误差的初始迭代值记为r⁽⁰⁾，一般r⁽⁰⁾的值选择为r⁽⁰⁾＝S；索引集合的初始迭代值，记为Ξ⁽⁰⁾，一般Ξ⁽⁰⁾的值选择为其中表示空集；k表示硬阈值OMP算法中的第k迭代次数，k的初始值设置为k＝0，且k的取值范围是从0到MaxIter；

步骤5、寻找测量矩阵与重构残余误差的最大相干项：

采用公式计算得到硬阈值OMP算法第k次迭代过程中测量矩阵与重构残余误差的最大相干项，记为I_m，其中表示求取满足括号中最大值时对应自变量m的最优值，||·||₂为向量的L2范数，Α_m为测量矩阵Α中的第m列，A为步骤3中得到的线阵SAR测量矩阵，右上角T为转置运算符号，r^(k-1)为算法第k-1次迭代过程中得到的重构残余误差，k表示硬阈值OMP算法中的第k迭代次数；若k＝1，r^(k-1)的值为步骤4中得到的初始迭代值r⁽⁰⁾，否则r^(k-1)通过硬阈值OMP算法第k-1次迭代过程中步骤10的重构残余误差计算提供；

步骤6、更新索引集合：

采用公式Ξ^(k)＝Ξ^(k-1)∪I_m计算得到硬阈值OMP算法第k次迭代过程中的索引集合，记为Ξ^(k)，其中Ξ^(k-1)为硬阈值OMP算法第k-1次迭代过程中得到的索引集合，I_m为步骤5得到的测量矩阵与重构残余误差的最大相干项，∪表示并集运算符号，k表示硬阈值OMP算法中的第k迭代次数；若k＝1，Ξ^(k-1)的值为步骤4中得到的初始值Ξ⁽⁰⁾；k＝2时，Ξ^(k-1)的值为硬阈值OMP算法第₁次迭代过程中得到的索引集合Ξ⁽¹⁾；k＝3时，Ξ^(k-1)的值为硬阈值OMP算法第2次迭代过程中得到的索引集合Ξ⁽²⁾；依次类推，在第k次迭代中且k＞1时，Ξ^(k-1)为硬阈值OMP算法第k-1次迭代过程中得到的索引集合；

步骤7、重构观测场景目标空间的散射系数：

采用公式和计算得到硬阈值OMP算法第k次迭代过程中观测场景目标空间的散射系数向量，记为α^(k)，其中为观测场景目标空间散射系数向量α^(k)中由索引集合Ξ^(k)对应的元素向量，Ξ^(k)为步骤6中硬阈值OMP算法第k次迭代过程中得到的索引集合，为观测场景目标空间散射系数向量α^(k)中索引集合对应的元素组成的向量，为索引集合Ξ^(k)在观测场景目标空间Ω的补集，为矩阵的伪逆矩阵，右上角为矩阵伪逆运算符号，为测量矩阵A中由索引集合对应的列组成的矩阵，A为步骤3中得到的线阵SAR测量矩阵，S为步骤3中得到的线阵SAR回波信号向量，∪表示并集运算符号，k表示硬阈值OMP算法中的第k迭代次数；

步骤8、计算目标散射系数最大最小对比度和目标散射系数变化率：

采用公式计算得到硬阈值OMP算法第k次迭代过程中的目标散射系数最大最小对比度，记为其中其中为步骤7中得到的观测场景目标空间散射系数向量α^(k)中由索引集合Ξ^(k)对应的元素向量，α^(k)为硬阈值OMP算法第k次迭代得到的观测场景目标空间的散射系数向量，Ξ^(k)为步骤6中第k次迭代过程中得到的索引集合，k表示硬阈值OMP算法中的第k迭代次数，max|·|表示求取向量中元素的最大绝对值，min|·|表示求取向量中元素的最小绝对值；采用公式计算得到第k次迭代过程中的目标散射系数变化率，记为其中α^(k)为硬阈值OMP算法第k次迭代得到观测场景目标空间的散射系数向量，α^(k-1)为硬阈值OMP算法第k-1次迭代得到观测场景目标空间的散射系数向量，α^(k-2)为硬阈值OMP算法第k-2次迭代得到观测场景目标空间的散射系数向量，||·||₂为向量L2范数；

步骤9、迭代终止判定：

如果且则执行步骤10，否则硬阈值OMP算法终止迭代，此刻硬阈值OMP算法第k次迭代得到的散射系数向量值α^(k)即为观测场景目标空间Ω最终的散射系数向量，其中为步骤8中得到的硬阈值OMP算法第k次迭代目标散射系数最大最小对比度，为步骤8中得到的硬阈值OMP算法第k次迭代目标散射系数变化率，η₀为步骤4中初始化得到的目标散射系数最大最小对比度门限；β₀为步骤4中初始化得到的目标散射射系数变化率门限，k为硬阈值OMP算法中的第k迭代次数；

步骤10、计算重构残余误差：

采用公式r^(k)＝S-Αα^(k)计算得到硬阈值OMP算法第k次迭代过程中的重构残余误差，记为r^(k)，其中S为步骤3中得到的线阵SAR回波信号向量，A为步骤3中得到的线阵SAR测量矩阵，α^(k)为步骤7中得到的硬阈值OMP算法第k次迭代的观测场景目标空间散射系数向量，k表示硬阈值OMP算法中的第k迭代次数；

步骤11、迭代终止判定：

如果若r^(k)≥ε₀且k≤MaxIter，则k←k+1，则执行步骤5至步骤9，否则终止迭代，此刻硬阈值OMP算法第k次迭代得到的散射系数向量值α^(k)即为线阵SAR观测场景目标空间Ω最终的散射系数向量，其中r^(k)为步骤10中得到的硬阈值OMP算法第k次迭代的重构残余误差，k表示硬阈值OMP算法中的第k迭代次数，MaxIter为步骤4中初始化得到的硬阈值OMP算法重构迭代处理的最大迭代次数；最后将观测场景目标空间散射系数向量α^(k)转换成三维矩阵形式，得到线阵SAR观测场景目标空间Ω的三维成像结果。