一种图像熵最优的压缩传感sar稀疏自聚焦成像方法
【专利摘要】本发明公开了一种图像熵最优的压缩传感SAR稀疏自聚焦成像方法，它是针对SAR回波信号测量模型中方位向相位误差对压缩传感SAR成像的影响，以及压缩传感SAR成像模型中的未知相位误差估计与补偿问题，利用成像模型中方位向相位误差特性和稀疏目标特征，采用压缩传感SAR图像熵作为评价准则，在每一次迭代处理过程中利用SAR方位向回波与观测目标的关系估计方位向相位误差，然后对压缩传感成像模型进行相位误差补偿，接着再进行压缩传感SAR成像，并采用逐次迭代方法使得压缩传感SAR成像的图像熵最优，从而提高了压缩传感SAR成像质量。
【专利说明】-种图像熵最优的压缩传感SAR稀疏自聚焦成像方法
【技术领域】：
[0001] 本技术发明属于雷达【技术领域】，它特别涉及了合成孔径雷达（SAR)成像技术领 域。
【背景技术】：
[0002] 由于具有全天时、全天候和大场景观测等优势，合成孔径雷达（SAR)已成为当前 大面积地形测绘的一项重要遥感技术，在地形测绘、自然灾害监测和自然资源调查等领域 发挥越来越大的作用。压缩传感稀疏重构作为一种近几年新提出的信号处理理论，突破了 传统Nyquist采样定理约束，可利用远低于Nyquist采样率精确重构原始稀疏信号（详 见参考文献 "D. L. Donoho. Compressed sensing. IEEE Transactions on Information Theory, 2006, 52 (4) : 1289-1306"），在降低SAR系统采样率和提高成像质量等方面有着巨 大的应用潜力。因此，压缩传感SAR成像已经成为了 SAR领域的新兴热点课题。压缩传感 重构方法对信号测量模型的精确度要求很高，如果信号测量模型存在误差或者不精确，压 缩传感重构方法的重构精确度就会严重退化，有时甚至会出现错误的结果。对于压缩传感 SAR实际数据成像，影响其成像效果的一个重要因素是运动误差估计与补偿（详见参考文 献''Jihua Tian, Sun Jinping1Han Xiao, and Zhang Bingchen. Motion compensation for Compressive Sensing SAR imaging with autofocus. IEEE Conference on Industrial Electronics and Applications(ICIEA), 1564-1567,2011"）。在压缩传感 SAR成像时需要 利用每个方位向平台位置建立SAR信号观测模型。但是，由于平台运动等因素存在测量误 差，使得每个方位向回波信号中混入相位误差，导致压缩传感SAR成像质量下降甚至无法 成像。因此在压缩传感SAR实际成像过程中，必须先对SAR数据每个方位向的未知相位误 差进行有效估计与补偿，才可对SAR原始回波数据进行高精度压缩传感成像。
[0003] 熵作为描述信息不确定性的尺度，在信号处理领域有着广泛的应用。对于SAR 成像，散焦图像对比度较低，表明图像各像素取值概率接近，不确定性较大，而聚焦图像对 比度较高，说明各像素取值的确定性增加。因此，图像熵可作为衡量SAR成像图像质量优 劣的重要尺度。对基于传统成像方法的SAR数据处理中，图像熵优化已被用于实现自聚 焦成像的一个重要手段（详见参考文献1.1^，6.3.1^11，311(1]\1^附.4111:〇;1；'〇。118；[1^01；' ISAR images based on entropy minimization. IEEE Transactions on Aerospace and Electronic Systems, 1999, 35(4) : 1240-1252")。然而，基于压缩传感的SAR成像方法与传 统SAR成像方法在信号处理上有本质的差别，基于压缩传感的SAR成像方法可从欠采样回 波数据中精确重构观测目标，从而实现高精度成像，但是传统SAR成像方法很难实现欠采 样回波数据高精度成像。目前，基于图像熵最优的传统自聚焦方法是针对全采样条件下的 SAR回波数据进行处理，在欠采样情况下这类方法自聚焦性能会严重下降甚至失效。因此， 为了在压缩传感SAR成像中利用图像熵进行自聚焦成像，不能直接利用传统图像熵自聚焦 成像方法，必须要结合压缩传感SAR信号重构模型和重构结果进行处理。
【发明内容】
：
[0004] 为了解决压缩传感SAR稀疏成像过程中未知方位向相位误差的估计与补偿问题， 本发明根据压缩传感SAR成像模型中方位向相位误差特性和稀疏目标特征，结合压缩传感 稀疏重构方法以及SAR图像熵，提出了一种图像熵最优的压缩传感SAR稀疏自聚焦成像方 法，本发明利用成像模型中方位向相位误差特性和稀疏目标特征，利用压缩传感SAR图像 熵作为评价准则，提出了一种图像熵最优的压缩传感SAR稀疏自聚焦成像方法，该方法在 每一次迭代处理过程中利用SAR方位向回波与观测目标的关系估计方位向相位误差，然后 对压缩传感成像模型进行相位误差补偿，接着再进行压缩传感SAR成像，并采用逐次迭代 方法使得压缩传感SAR成像的图像熵最优，从而提高了压缩传感SAR成像质量。
[0005] 为了方便描述本发明的内容，首先作以下术语定义：
[0006] 定义1、稀疏信号
[0007] 如果一个离散信号中非零值的个数远小于信号本身的长度，则该信号可认为是稀 疏的。设X= [Xl，x2，…，xN]T为N个离散信号组成的列向量，其中&表示向量X中的第1 个元素，X2表示向量X中的第2个元素，xN表示向量X中的第N个元素，右上角正体符号T 为转置运算符号。如果向量X中仅有Ktl个元素非零或远大于零，则向量X定义为Ktl稀疏向 量。详见文献"S. Mallat. A Wavelet Tour of Signal Processing:The Sparse Way. Access Online via Elsevier，2008"。
[0008] 定义2、范数
[0009] 设X是数域C上线性空间，C表示复数域，若它满足如下性质：| |X| I彡0,且 |x| I = 0仅有 X = 0, I IaXl I = |a| I |X| I，a为任意常数，I |Xi+X2| I < I IX1I | + | Ix2I I， 则称I |X| I为X空间上的范数，U ? I I表示范数符号，其中X1和X2为X空间上的任意两 个值。对于定义1中的NXl维离散信号向量X = [X1, X2，…，xN]T，向量X的LP范数表达
【权利要求】
1. 一种图像熵最优的SAR稀疏自聚焦成像方法，其特征是它包括如下步骤： 步骤1、初始化SAR系统参数： 初始化SAR系统参数包括：平台速度矢量，记做V ;天线初始位置矢量，即方位向O慢 时刻天线位置，记做P(O);雷达工作中心频率，记做f。；雷达载频波长，记做λ ;雷达发射基 带信号的信号带宽，记做4 ;雷达发射信号脉冲宽度，记做Tp ;雷达发射信号的调频斜率，记 做;雷达脉冲重复频率，记做PRF ;方位向等效天线长度，记做Da;雷达接收系统的采样频 率，记做fs ;光在空气中的传播速度，记做C ;距离向快时刻序列，记做t，t = 1，2,…，NK，Nk 为距离向快时刻总数；方位向慢时刻序列，记做1，1 = 1，2,…，Na，Na为方位向慢时刻总数； 上述参数均为SAR系统标准参数，在SAR系统设计和观测过程中已经确定；根据SAR成像系 统方案和观测方案，SAR成像方法需要的初始化系统参数均为已知； 步骤2、初始化SAR成像空间参数以及获取原始回波信号： 初始化SAR的成像空间参数，包括：以雷达波束照射场区域地平面所构成的空间坐标 作为SAR的成像空间，该成像空间记为Ω ;平台到SAR成像空间中心的参考斜距，记做Rref ; 将成像空间Ω均匀划分成大小相等的平面单元格，也称为分辨单元，单元网格在水平横 向、水平纵向边长分别记为4和d y，单元网格大小选择为SAR系统传统理论成像分辨率或 SAR系统传统理论成像分辨率的二分之一；观测场景目标空间Ω中第m个单元格的坐标矢 量，记做Pm，m表示SAR成像空间Ω中第m个单元格，m = 1，2,…，M，M为成像空间Ω中的 单元格总数；SAR成像空间Ω中所有单元格的散射系数按位置顺序排列组成向量，记做α， 向量α由M行1列组成；散射系数向量α中第m个元素，记做a m，m = 1，2,…，Μ ; 根据步骤1中初始化的平台速度矢量V，天线初始位置矢量P(O)和雷达系统的脉冲重 复频率PRF，采用公式P(I) =P(O)+V ^VPRF, 1 = 1，2，···，Na，计算得到天线在第1个方位 向慢时刻的位置矢量，记为P(I)，Na为步骤1的方位向慢时刻总数；采用公式R(P(1)，P m) =I |p⑴-PmI |2，1 = 1，2,…，乂，!!! = 1，2, "·，Μ，计算得到在第1个方位向慢时刻SAR成像 空间Ω中第m个单元格到天线的距离，记为R(P(l)，Pm)，其中I I · I |2表示定义2中的向 量L2范数，Pm为初始化得到成像空间Ω中第m个单元格的坐标矢量，M为初始化的成像空 间 Ω 中单元格总数；采用公式 Tm(I) = 2*R(P(l)，Pm)/C，l = 1，2,…，NA，m= 1，2,…，M， 计算得到在第1个方位向慢时刻SAR成像空间Ω中第m个单元格到天线的时间延时，记为 τ m(l)，其中C为步骤1中初始化得到的光在空气中的传播速度； 第1个方位向慢时刻和第t个距离向快时刻中SAR天线的原始回波数据记为s (t，1)， t = 1，2,…，NK，1 = 1，2,…，Na，其中Nk为步骤1中初始化的距离向快时刻总数；在SAR实 际成像中，回波数据s (t, 1)，t = 1,2,…，乂，1 = 1,2,…，Na，可由SAR系统中雷达数据接 收机提供； 步骤3、建立SAR回波信号的线性测量模型： 将步骤2中获取所有SAR原始回波信号s(t，1)按顺序排列组成向量，记为回波信号向 量S，回波信号向量S由D行1列组成，其中D = Na ·ΝΚ，Nk为步骤1中初始化的距离向快时 刻总数，Na为步骤1中初始化的方位向慢时刻总数； 采用公式
1 = 1，2, "'Na,!!! = 1，2, "'IVUd = 1，2, ···，0,计算得到成像空间Ω中第m个单元格在回 波信号向量S第d个元素信号对应的时延函数，记为Φ d(m)，其中exp ( ·)表示e指数运算 符号，f。为步骤1初始化得到的雷达工作中心频率，为步骤1初始化得到的发射信号调 频斜率，τπ(1)为步骤2得到的在第1个方位向慢时刻SAR成像空间Ω中第 m个单元格到 天线的时间延时，t为距离向的第t个快时刻，j为虚数单位（即-1的开根值），π为圆周 率； SAR原始回波信号向量S与成像空间Ω中所有单元格散射系数向量α之间的测量矩 阵，记为A ;测量矩阵A由SAR成像空间Ω中所有单元格对应的时延函数构成，A为D行M 列的二维矩阵，具体表达式为
其中，￠:(1)为成像空间Ω中第1个单元格在回波信号向量S第1个元素信号对应 的时延函数，为成像空间Ω中第2个单元格在回波信号向量S第1个元素信号对应 的时延函数，Φ1 (M)为成像空间Ω中第M个单元格在回波信号向量S第1个元素信号对应 的时延函数，Φ2(1)为成像空间Ω中第1个单元格在回波信号向量S第2个元素信号对应 的时延函数，Φ 2 (2)为成像空间Ω中第2个单元格在回波信号向量S第2个元素信号对应 的时延函数，Φ2 (M)为成像空间Ω中第M个单元格在回波信号向量S第2个元素信号对应 的时延函数，￠:/1)为成像空间Ω中第1个单元格在回波信号向量S第D个元素信号对应 的时延函数，Φ Β(2)为成像空间Ω中第2个单元格在回波信号向量S第D个元素信号对应 的时延函数，Φγ/Μ)为成像空间Ω中第M个单元格在回波信号向量S第D个元素信号对 应的时延函数，Φ^Ι)，…，C^ 1(M)分别为成像空间Ω中第1，2,…，M个单元格在回 波信号向量S第1个元素信号对应的时延函数向量，Φ 2(1)，Φ2(2)，···，Φ2(Μ)分别为成像 空间Ω中第1，2，···，M个单元格在回波信号向量S第2个元素信号对应的时延函数向量， iK⑴，、⑵，…，ΦΒ(Μ)分别为成像空间Ω中第1，2,…,M个单兀格在回波信号向量S第 D个元素信号对应的时延函数向量； 步骤4、对SAR成像空间的散射系数向量进行初步估计： 采用公式Aw=Ah.S计算得到SAR成像空间Ω中散射系数向量的初始估计值，记做 ，其中Ah为矩阵A的共轭转置，A为步骤3中得到的SAR测量矩阵，上标H为共轭转置 运算符号，S为步骤3中得到的SAR原始回波信号向量； 步骤5、初始化稀疏自聚焦成像算法所需的参数： 初始化稀疏自聚焦成像算法所需参数，包括：成像处理迭代的最大迭代次数， 记做Maxiter ;成像处理第i次迭代过程中SAR成像模型中的相位误差向量，记为
，相位误差向量Φ ω为D行1列的向量，i为自然数，表示为成像处 理的第i次迭代，i = 1，2, "'Maxitei,元素为向量Φω中的第1个元素,元素为 向量Φ ω中的第2个元素，元素 ff为向量Φ ω中的第D个元素，右上角正体符号T为转 置运算符号；成像处理第i次迭代过程中第1个方位向慢时刻的相位误差值，记为C，I =1，2,…，Na ;成像处理迭代的阈值,记做τ ;采用公式w(1) = diag(exp(j · Φ (1))), i = 1，2,…，Maxiter，计算得到成像处理第i次迭代过程中SAR测量模型的初始相位误差矩阵 值，记做W(i)，i = 1，2, "^Maxiter,其中相位误差矩阵W(i)为0行D列的对角矩阵，diag( ·) 为将向量元素作为对角矩阵中对角元素的运算符号；成像处理第i次迭代过程中SAR成像 图像的图像熵，记为L (i)，i = 1，2,…，Maxiter ;在成像迭代处理之前，成像处理第O次迭代 中SAR测量模型相位误差矩阵的初值，记为Ww，成像处理第O次迭代中SAR成像图像熵的 初值，记为L w ; 步骤6、估计SAR方位向的相位误差向量： 包括以下步骤： 步骤6. 1、获取目标回波的响应向量： 采用合成孔径雷达标准距离压缩方法对SAR原始回波数据s (t，I)，t = 1，2,…，Νκ， 1 = 1，2，···，ΝΑ，进行距离压缩处理，得到距离压缩后的SAR回波数据，记做Si(t，l)，t = 1，2,…，乂，1 = 1，2,…，Na，其中s (t，1)为步骤2得到的第1个方位向慢时刻和第t个距 离向快时刻中SAR天线的原始回波数据； 采用公式
l，2，...，NA，m= 1，2，...，M，和公式 S1= [β(1，1)，β (1，2)，...，M1，M)]，1 = 1，2，...，Na， 计算得到第I个方位向目标回波的响应向量，记为S1，其中β (1，I)为第I个方位向慢时刻 SAR成像空间Ω中第1个单元格对应的β (l，m)，β (1，2)为第1个方位向慢时刻SAR成 像空间Ω中第2个单元格对应的β (l，m)，β (1，Μ)为第1个方位向慢时刻SAR成像空间 Ω中第M个单元格对应的β (l，m)，Ceil(·)表示上取整运算符号，Rref为步骤2初始化得 到的平台到SAR成像空间中心参考斜距，C为步骤1初始化得到的光在空气中的传播速度， fs为步骤1初始化得到的雷达接收系统的采样频率，R(P (I)，Pm)为步骤2中得到的第1个 方位向慢时刻SAR成像空间Ω中第m个单元格到天线的距离，λ为步骤1中初始化得到 的波长； 令q为自然数，q取值范围为q=l，2,…，ΝΑ;当。=1时，采用公式
计算得到第q个方位向慢时刻的数据 向量，记为X ;当q = Na时，采用公式
计算得到第q个方位向慢时刻的数据向量，记为X ;当I < q < Na时，采用公式
，计算得到第q个方位向慢时刻的数据向量，记 为X，向量X的第m个元素记为Xm,其中f为成像处理第i次迭代过程中第1个方位向慢 时刻的相位误差值，1 = 1，2, "·，ΝΑ，S1为第1个方位向目标回波的响应向量； 采用公式Y = Sq计算得到第q个方位向慢时刻的数据向量，记为Y，向量Y的第m个元 素记为y^，其中Sq为1等于q时第1个方位向目标回波的响应向量，即Sq = S1 ; 采用公式 1
I M，计算得到元素 fm;再采用公式F =[f\，f2，…，fM]计算得到的数据向量，记为F，其中卜I2为绝对值平方运算符号，4为!11 =1时对应的元素 fm，f2为m = 2时对应的元素 fm，fM为m = M时对应的元素 fm，xm为向量 X的第m个元素，ym为向量Y的第m个元素； 步骤6. 2、利用施密特正交化方法获取单位正交向量： 采用公式
，M，计算得到元素 am和 bm，其中Xm为步骤6. 1得到的向量X中第m个元素，ym为步骤6. 1得到的向量Y中第m个元 素；再采用a = Iia1, a2,…，aM]、b = Dd1, b2,…，bM]，m = 1，2,…，M，计算得到向量a和向量 b，S1为m = 1时对应的元素 am，a2为m = 2时对应的元素 am，aM为m = M时对应的元素 am, Id1为m = 1时对应的元素 bm, b2为m = 2时对应的元素 bm, bM为m = M时对应的元素 bm ; 利用施密特正交化方法对向量a和向量b进行标准正交化，得到向量a和向量b构成 的平面对应的单位正交向量
，其中flI为向量i的第1个元素，&为 向量S的第2个元素，I1为向量I的第1个元素，为向量β的第2个元素； 步骤6. 3、计算二次型矩阵及参数估计： 采用公式
- 算得到二次型矩阵R，其中i,为向量的第1个元素，52为向量s的第2个元素，$为向量 的第1个元素，为向量:6 I的第2个元素，s和S分别为步骤6. 2得到的单位正交向量； 采用公式Γ · Λ · Γτ = R对二次型矩阵R进行特征值分解，得到二次型矩阵R的特 征向量矩阵记为Γ，Γτ为矩阵Γ的转置，得到二次型矩阵R的特征值矩阵记为Λ ; 采用公式η = (Γ ·Χ(ΓΙ) ^IT1计算得到参数的估计值，记为η，其中Χ(ι为向量F在 向量a和向量b所张成的二维平面的垂足点，I为单位矩阵，IT1为二次型矩阵R的逆矩阵； 步骤6. 4、估计方位向相位误差： 采用公式σ = [a，br· U ?R+ir'xjf算得到的向量，记为〇,其中向量〇维数 为2X 1，向量〇中第1个元素记为σ i，向量σ中第2个元素记为〇 2, a和b为步骤6. 2 得到的向量，Π 为步骤6. 3得到的估计值向量，R为步骤6. 3得到的二次型矩阵，Xtl为向量 f在向量a和b所张成的二维平面Ω内的垂足点，I为单位矩阵，上标表示矩阵求逆运算 符号； 采用公式
丨计算更新成像处理第i次迭代过程中第q个方位向慢时刻 的相位误差向量，其中atan( ·)为求解正切函数反函数运算符号；将相位误差向量Φω中 第（q-Ι) · Νκ+1至q · Nk个元素…，<1,.的值全部赋值为 <，其中Nk为距离向快 时刻总数，hU为相位误差向量Φω中的第（Cr 1) ·ΝΚ+1个元素，为相位误差向量 Φω中的第q·乂个元素； 步骤6. 5、方位向相位误差值逐个估计： 对于所有方位向序号q，q = 1，2,…，Na，采用步骤6. 1到步骤6. 4逐个方位向直到估 计所有方位向相位误差，最终得到方位向的相位误差向量Φω ; 步骤7、压缩传感稀疏成像： 采用公式W(i) = diag (exp (j · Φ ω))计算得到成像处理第i次迭代过程中的相位误差 矩阵，记做W(i)，其中Φ ω为步骤6. 5得到的成像处理第i次迭代过程中的相位误差向量， diag( ·)为将向量元素作为对角矩阵对角元素的运算符号； 采用公式
和标准压缩传感稀疏重构方法计算得到 成像处理第i次迭代过程中SAR目标成像空间的散射系数向量，记做，》，其中
表 示求取满足括号中最小值时对应自变量α的最优值，向量S为步骤3得到的SAR数据回波 信号向量，||表示向量L2范数的平方运算符号，Μ · N1表示向量Ll范数运算符号； 步骤8、计算SAR图像熵： 采用公式
计算得到成像处理 第i次迭代过程中SAR成像图像熵，记做L(1)，i = 1，2, "^Maxiter,其中(?1为步骤8中得 至Ll的散射系数向量#(〇的第m个元素，i为成像处理的第i次迭代，i = 1，2,…，Maxiter，M 为步骤2初始化得到的SAR成像空间Ω中划分的单元格总数，Hg为向量L2范数的平方运 算符号，Σ ( ·)表示向量元素求和运算符号，l〇g2( ·)表示底数为2的对数运算符号； 步骤9、成像算法迭代条件判定： 若满足条件I I Lw-LqA I 12彡τ且i彡Maxiter，则令成像处理迭代次数i的值加1，然 后重复执行步骤6到步骤9,其中L(i)为成像处理第i次迭代过程中SAR成像图像熵， 为成像处理第i-Ι次迭代过程中SAR成像图像熵，τ为步骤5中初始化得到的算法迭代阈 值，Maxiter为步骤5中初始化得到的成像处理最大迭代次数，I I · I |2为向量L2范数运算 符号； 若满足条件I I L⑴-LqI I 12 < τ或者i > Maxiter，则成像处理第i次迭代过程得到 的散射系数向量和相位误差向量Φ (i)为最终的SAR成像散射系数向量和相位误差估计 结果。
【文档编号】G01S13/90GK104391295SQ201410442888
【公开日】2015年3月4日 申请日期:2014年9月2日 优先权日:2014年9月2日
【发明者】韦顺军, 张晓玲, 熊海进 申请人:电子科技大学