CN103713288B - 基于迭代最小化稀疏贝叶斯重构线阵sar成像方法 - Google Patents

基于迭代最小化稀疏贝叶斯重构线阵sar成像方法 Download PDF

Info

Publication number
CN103713288B
CN103713288B CN201310752685.0A CN201310752685A CN103713288B CN 103713288 B CN103713288 B CN 103713288B CN 201310752685 A CN201310752685 A CN 201310752685A CN 103713288 B CN103713288 B CN 103713288B
Authority
CN
China
Prior art keywords
linear array
scattering coefficient
object space
vector
scene object
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Expired - Fee Related
Application number
CN201310752685.0A
Other languages
English (en)
Other versions
CN103713288A (zh
Inventor
张晓玲
韦顺军
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
University of Electronic Science and Technology of China
Original Assignee
University of Electronic Science and Technology of China
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by University of Electronic Science and Technology of China filed Critical University of Electronic Science and Technology of China
Priority to CN201310752685.0A priority Critical patent/CN103713288B/zh
Publication of CN103713288A publication Critical patent/CN103713288A/zh
Application granted granted Critical
Publication of CN103713288B publication Critical patent/CN103713288B/zh
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88Radar or analogous systems specially adapted for specific applications
    • G01S13/89Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G01S13/9004SAR image acquisition techniques
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88Radar or analogous systems specially adapted for specific applications
    • G01S13/89Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G01S13/904SAR modes
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S7/00Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/02Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
    • G01S7/41Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section

Abstract

本发明提供了一种基于迭代最小化稀疏贝叶斯重构线阵SAR成像方法,它是基于线阵SAR原始回波信号测量模型的先验分布假设,假设线阵SAR观测场景目标空间中散射系数的先验概率密度函数服从复指数先验分布,线阵SAR原始回波信号的后验概率密度函数服从高斯随机分布,再利用贝叶斯准则构造线阵SAR观测场景目标空间中散射系数的重构代价函数,通过复指数分布参数最优化估计以及迭代最小化重构代价函数实现线阵SAR观测场景目标空间中的散射系数稀疏重构,从而提高了线阵SAR稀疏成像的性能。本发明可以应用于合成孔径雷达成像和地球遥感等领域。

Description

基于迭代最小化稀疏贝叶斯重构线阵SAR成像方法
技术领域:
本技术发明属于雷达技术领域,它特别涉及了合成孔径雷达(SAR)成像技术领域。
背景技术:
对比传统二维合成孔径雷达(SAR),三维合成孔径雷达(三维SAR)可以获得观测场景目标的三维成像,能够更加精细地描述观测场景中目标的几何和散射特征,更利于目标特征提取和识别,是SAR成像技术未来发展的必然趋势以及当前的研究热点。三维SAR成像的基本原理是通过单个或者多个天线的运动合成一个大的虚拟二维面阵天线,从而获得面阵平面内二维高分辨,再结合脉冲压缩技术获得雷达视线方向高分辨率,实现对观测场景目标的三维成像。线阵合成孔径雷达(Linear Array SAR,LASAR,简称线阵SAR)是近几年新提出的一种新型三维SAR成像技术,该技术在运动平台上放置一列线阵天线,结合载荷平台的运动合成二维虚拟面阵的三维SAR成像技术。与其它三维SAR成像技术(如圆周SAR和层析SAR等)相比,线阵SAR具有多模式工作能力,除了可工作于传统SAR侧视模式,还可工作于下视模式以及前视模式,在实际成像具有应用灵活及多功能等优势。目前线阵SAR成像方法主要基于匹配滤波(Matched Filter,MF)理论,如三维距离-多普勒(RD)算法和三维后向投影(BP)算法,见参考文献“G.Fornaro,F.Serafino,and F.Soldovieri.Three-dimensional Focusing withMultipass 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-DImaging Method via Multiresolution Approximation.IEEE Trans.Geosci.Remote Sens,Vol.46,No.11,pp.3719–3730,Nov.2008.”,该类算法在频域或时域上通过线阵SAR回波数据的相参积累获得观测场景目标的三维成像。虽然传统匹配滤波算法运算效率较高,但匹配滤波成像算法由于受线阵长度和分辨率瑞利准则限制,成像分辨率较低而且在稀疏线阵时存在较严重主瓣展宽和旁瓣干扰,成像质量不能满足高分辨成像应用的要求。研究获取高分辨或者超分辨能力成像算法成为了当前线阵SAR成像技术的一个迫切需要。
线阵SAR成像是从原始回波信号中重构出目标散射系数的过程,建立线阵SAR回波信号和观测场景目标空间散射系数的线阵测量模型后,线阵SAR成像问题就可等效为三维观测场景目标空间散射系数的线性方程逆求解问题。在线阵SAR成像的三维观测场景目标空间中,由于大多数区域不包含散射点(如,空气)或散射点被其他散射点遮挡而无法被入射波束照射(如,地下目标),线阵SAR三维观测场景目标往往表现出典型的空间稀疏特征。因此线阵SAR成像可以转化为三维观测场景中稀疏目标散射系数的估计与重构过程,在线阵SAR稀疏成像过程中只需要估计三维观测场景中稀疏目标分辨单元的散射系数,并不需要估计整个三维观测场景所有分辨单元的散射系数。
近年来,压缩感知稀疏重构理论成为信号处理领域的研究热点,将压缩感知稀疏重构理论应用于线阵SAR成像技术可以突破传统匹配滤波理论成像的分辨率限制,为提高线阵SAR成像精度提供了一种新的技术途径。在稀疏信号重构理论中,相对于传统贪婪稀疏重构算法,如匹配追踪(Matching Pursuit,MP)算法和正交匹配追踪(Orthogonal Matching Pursuit,OMP)算法,参考文献“J.A.Tropp.Greed is Good:Algorithmic Results for Sparse Approximation.IEEETransactions on Information Theory,Vol.50,No.10,pp.2231-2242,2004”和“J.A.Tropp,A.C.Gilbert.Signal Recovery from Random Measurements via Orthogonal Matching Pursuit.IEEETransactions on Information Theory,Vol.53,No.12,pp.4655-4666,2007”,稀疏贝叶斯重构方法通过选择信号测量模型不同的先验概率分布,可以更加灵活和精确地构造稀疏信号的重构模型或重构函数,从而提高稀疏信号的重构精度,另外还可获得估计信号的协方差矩阵,能够评估重构算法获得的稀疏信号估计结果的误差范围。贝叶斯压缩采样(Bayesian CompressedSensing,BCS)算法是稀疏贝叶斯稀疏重构方法中的经典算法,详见参考文献“Ji.S,Xue.Y,Carin.L.Bayesian Compressive Sensing.IEEE Transactions on Signal Processing,Vol.56,No.6,pp.2346-2356,2008”。但是,BCS算法基于稀疏测量信号服从高斯随机分布的假设前提,算法中待确定的参数较多,在线阵SAR稀疏成像时BCS多个参数合理选择困难,参数选择不当时重构精度下降。
发明内容:
为了提高线阵SAR稀疏成像的精度,本发明结合线阵SAR三维观测场景目标的稀疏特征以及回波测量模型的先验分布,利用贝叶斯准则和似然函数构造重构代价函数,提供了一种基于迭代最小化稀疏贝叶斯重构的线阵SAR成像方法。该方法的主要思路是:利用线阵SAR雷达系统参数、运动平台参数和观测场景目标的空间参数与原始回波信号的相互关系,建立线阵SAR原始回波信号与三维观测场景目标散射系数之间的线性测量模型,然后基于线阵SAR原始回波信号线性测量模型的先验分布假设,假设线阵SAR观测场景目标空间中散射系数的先验概率密度函数服从复指数先验分布,线阵SAR原始回波信号的后验概率密度函数服从高斯随机分布,再利用贝叶斯准则和似然函数构造线阵SAR观测场景目标空间中散射系数的重构代价函数,通过迭代最小化重构代价函数实现线阵SAR观测场景目标空间中的散射系数稀疏重构。该方法的特点是:1)基于线阵SAR原始回波信号测量模型的先验概率分布假设,合理选择线阵SAR观测场景目标空间的散射系数和原始回波信号的概率密度函数分布;2)结合贝叶斯准则和似然函数构造线阵SAR的重构代价函数;3)利用重构代价函数的迭代最小化方法获得线阵SAR观测场景目标空间散射系数的稀疏成像结果。
为了方便描述本发明的内容,首先作以下术语定义:
定义1、稀疏信号
如果一个离散信号中非零值的个数远小于信号本身的长度,则该信号可认为是稀疏的。设X=[x1,x2,…,xN]T为N个离散信号组成的列向量,其中x1表示向量X中的第1个元素,x2表示向量X中的第2个元素,xN表示向量X中的第N个元素,右上角T为转置运算符号。如果向量X中仅有K0个元素非零或远大于零,则向量X定义为K0稀疏向量。详见文献“S.Mallat.A Wavelet Tour of Signal Processing:The Sparse Way.Access Online via Elsevier,2008”。
定义2、范数
设X是数域上线性空间,表示复数域,若它满足如下性质:||X||≥0,且||X||=0仅有X=0,||aX||=|a|||X||,a为任意常数,||X1+X2||≤||X1||+||X2||,则称||X||为X空间上的范数,||·||表示范数符号,其中X1和X2为X空间上的任意两个值。对于定义1中的N×1维离散信号向量X=[x1,x2,…,xN]T,向量X的LP范数表达式为其中xi为向量X的第i个元素,|·|表示绝对值符号,Σ|·|表示绝对值求和符号,向量X的L1范数表达式为向量X的L2范数表达式为向量X的L0范数表达式为且xi≠0。详见文献“矩阵理论”,黄廷祝等编著,高等教育出版社出版。
定义3、信号线性测量模型
对于一个数字信号测量系统,假设定义1中的N×1维离散信号向量X=[x1,x2,…,xN]T为该数字信号测量系统需要测量的信号,向量Y=[y1,y2,…,yM]T为该数字信号测量系统输出的M×1维离散信号向量,其中y1表示向量Y中的第1个元素,y2表示向量Y中的第2个元素,yM表示向量Y中的第M个元素,右上角T为转置运算符号。该测量系统的信号线性测量模型是指测量信号Y与被测量信号X的关系可以表示为Y=AX,其中A为M×N矩阵,矩阵A称为测量系统中信号X的测量矩阵。
定义4、贝叶斯准则
贝叶斯准则用来描述两个随机事件条件概率之间的关系。对于随机事件A和B,贝叶斯准则指出随机事件A和B的条件概率满足以下等式:
Pr ( A | B ) = Pr ( A | B ) Pr ( B ) Pr ( A )
其中,Pr(A|B)为事件B发生后事情A发生的条件概率,亦称为在事件B条件下事件A的后验概率,Pr(B|A)为事件A发生后事件B发生的条件概率,亦称为在事件A条件下事件B的后验概率,Pr(A)为事件A的先验概率或者边缘概率,Pr(B)为事件B的先验概率或者边缘概率。详见参考文献“概率论与数理统计(第4版),盛骤、谢式千和潘承毅著,高等教育出版社”。
定义5、线阵合成孔径雷达(Linear array SAR,LASAR,线阵SAR)
线阵合成孔径雷达成像是将线性阵列天线固定于载荷运动平台上并与平台运动方向与垂直,结合运动平台的运动以合成二维平面阵列实现阵列平面维二维成像,再利用雷达波束向回波延时实现距离一维成像,从而实现观测目标三维成像的一种合成孔径雷达技术。
定义6、线阵SAR慢时刻与快时刻
线阵SAR运动平台飞过一个方位向合成孔径长度所需要的时间称为慢时间,雷达系统以一定时间长度的重复周期发射接收脉冲,因此慢时间可以表示为一个以脉冲重复周期为步长的离散化时间变量,其中每一个脉冲重复周期离散时间变量值为一个慢时刻。快时刻是指在一个脉冲重复周期内,距离向采样回波信号的时间间隔变量。详见文献“合成孔径雷达成像原理”,皮亦鸣等编著,电子科技大学出版社出版。
定义7、线阵SAR观测场景目标空间
线阵SAR观测场景目标空间是指现实空间中所有待观测场景目标散射点的集合。观测场景目标空间在不同空间坐标系下有不同的表示,但一旦坐标系确立以后其表示是唯一的。一般情况下为了方便成像,线阵SAR观测场景目标空间取为地面坐标系。本发明中线阵SAR观测场景目标空间记为Ω,用以下数学关系表示:
其中表示构成线阵SAR场景目标空间Ω的地表正交坐标基,分别表示水平横向、水平纵向和垂直地表的高度向,为场景目标空间中一个分辨单元位置向量,x、y和z分别表示该分布单元的水平横向、水平纵向和高度向坐标,表示实数域。
定义8、线阵SAR成像空间
线阵SAR成像空间是指将场景目标空间中的散射点投影到切航迹向-沿航迹向-距离向的三维空间坐标系,该空间由线阵SAR成像空间中的三个相互正交的坐标基确定。本发明中线阵SAR成像空间记为M,用以下数学关系表示:
其中表示构成线阵SAR成像空间M的正交坐标基,分别表示切航迹向、沿航迹向和距离向,为成像空间中的待观测点位置向量,u、v和z分别表示该点的切航迹向、沿航迹向和距离向坐标,表示实数域。
定义9、线阵SAR传统理论成像分辨率
线阵SAR传统理论成像分辨率是指利用经典匹配滤波理论成像算法得到线阵SAR系统在距离向、方位向和切航迹向的成像分辨率。对于收发共用线阵天线阵元,线阵SAR距离向的分辨率记为ρr,近似表达式为其中C为光在空气中的传播速度,Br为线阵SAR发射信号的带宽;方位向的分辨率记为ρa,近似表达式为其中Da为天线在方位向的真实孔径;切航迹向的分辨率记为ρc,近似表达式为其中λ为线阵SAR雷达载频波长,R0为线阵SAR平台到观测场景中心的参考斜距,L为线阵天线长度。详见参考文献“Shi.J,Zhang.X.L,et al.,APC Trajectory Design for One-Active Linear-arrayThree-dimensional Imaging SAR,IEEE Transactions on Geoscience and Remote Sensing,Vol.48,No.3,pp:1470-1486,2010”。
定义10、共轭梯度算法
共轭梯度算法是一种求解特定大规模线性方程组数值解的快速方法,它要求线性方程组的系数矩阵是正定共轭矩阵。共轭梯度算法主要思路是采用迭代逼近估计的方法,在每次迭代中利用已知点处的梯度构造一组共轭方向,并沿这组方向进行搜素,求出目标函数的极小值点。在大规模线性方程组求解时,共轭梯度算法相对于逆矩阵求解方法具有所需存储量少和计算效率快等优点。详见参考文献“M.R.Hestenes,E.Stiefel.Methods of Conjugate Gradientsfor Solving Linear Systems.Journal of Research of the National Nureau of Standard,Vol.49,No.6,pp.409-463,1952”。
定义11、Stein无偏风险估计量方法
Stein无偏风险估计量方法是一种模型参数自适应估计的经典方法,该方法的主要思想是通过估计量的均方误差最小化进行最优参数估计。详见参考文献“B.Efron,C.Morris.DataAnalysis Using Stein's Estimator and Its Generalizations.Journal of the American StatisticalAssociation,vol.70,No.350,pp.311-319,1975”。
定义12、合成孔径雷达原始回波仿真方法
合成孔径雷达原始回波仿真方法是指基于合成孔径雷达成像原理仿真出一定系统参数条件下具有合成孔径雷达回波信号特性的原始信号的方法,详细内容可参考文献:“InSAR回波信号与系统仿真研究”,张剑琦,哈尔滨工业大学硕士论文。
本发明提供的一种基于迭代最小化稀疏贝叶斯重构的线阵SAR成像方法,它包括以下步骤:
步骤1、初始化线阵SAR系统参数:
初始化线阵SAR系统参数包括:平台速度矢量,记做线阵天线各阵元初始位置矢量,记做其中n为线阵天线中第n个阵元的序号,n为自然数,n=1,2,...,N,N为线阵天线的阵元总数;线阵天线长度,记做L;雷达工作中心频率,记做fc;雷达载频波长,记做λ;雷达发射基带信号的信号带宽,记做Br;雷达发射信号脉冲宽度,记做TP;雷达发射信号的调频斜率,记做fdr;雷达接收波门持续宽度,记做To;雷达接收系统的采样频率,记做fs;雷达发射系统的脉冲重复频率,记做PRF;雷达系统的脉冲重复时间,记为PRI;雷达接收系统接收波门相对于发射信号发散波门的延迟,记做TD;天线在方位向的有效孔径长度,记做Da;光在空气中的传播速度,记做C;距离向快时刻,记做t,t为自然数,t=1,2,…,T,T为距离向快时刻总数;方位向慢时刻,记做l,l为自然数,l=1,2,…,K,K为方位向慢时刻总数;上述参数均为线阵SAR系统标准参数,其中线阵天线的阵元总数N,线阵天线长度L,相邻天线阵元之间的间距d,雷达中心频率fc,雷达载频波长λ,雷达发射基带信号的信号带宽Br,雷达发射信号脉冲宽度TP,雷达发射信号调频斜率fdr,雷达接收波门持续宽度To,雷达接收系统的采样频率fs,雷达系统的脉冲重复频率PRF,雷达系统的脉冲重复时间PRI和雷达接收系统接收波门相对于发射信号发散波门的延迟TD,天线在方位向的有效孔径长度Da在线阵SAR系统设计过程中已经确定;平台速度矢量及线阵天线各阵元初始位置矢量在线阵SAR观测方案设计中已经确定。根据线阵SAR成像系统方案和观测方案,线阵SAR成像方法需要的初始化系统参数均为已知。
步骤2、初始化线阵SAR观测场景目标空间参数:
初始化线阵SAR的观测场景目标空间参数,包括:以雷达波束照射场区域地平面和垂直于该地平面向上的单位向量所构成的空间直角坐标作为线阵SAR的观测场景目标空间Ω;将观测场景目标空间Ω均匀划分成大小相等的立体单元格(亦称为分辨单元),单元网格在水平横向、水平纵向和高度向边长分别记为dx、dy和dz,单元格大小一般选择为线阵SAR系统传统理论成像分辨率或该分辨率的二分之一;观测场景目标空间Ω中第m个单元格的坐标矢量,记做m表示观测场景目标空间Ω中第m个单元格,m为自然数,m=1,2,…,M,M为观测场景目标空间Ω中的单元格总数;观测场景目标空间Ω中所有单元格的散射系数按位置顺序排列组成向量,记做α,向量α由M行1列组成;散射系数向量α中第m个元素的散射系数,记做αm。观测场景目标空间Ω在线阵SAR成像方案设计中已经确定。
步骤3、建立线阵SAR原始回波信号与观测场景目标散射系数的线性测量矩阵:
采用公式n=1,2,…,N,l=1,2,…,K,计算得到第n个线阵天线阵元在第l个方位向慢时刻的位置矢量,记为其中为步骤1初始化得到线阵天线各阵元的初始位置矢量,为步骤1中初始化得到的平台速度矢量,PRF为步骤1初始化得到的雷达系统脉冲重复频率,n为线阵天线中第n个阵元的序号,n=1,2,…,N,N为步骤1中初始化得到的线阵天线阵元总数,l为方位向慢时刻中第l个慢时刻序号,l=1,2,…,K,K为步骤1初始化得到的线阵SAR方位向慢时刻总数。
采用公式n=1,2,…,N,l=1,2,…,K,m=1,2,…,M,计算得到在方位向第l个慢时刻线阵SAR观测场景目标空间Ω中第m个单元格到第n个线阵天线阵元的距离,记为其中||·||2表示向量的L2范数,为步骤2中初始化得到观测场景目标空间Ω中第m个单元格的坐标矢量,m表示观测场景目标空间Ω中第m个单元格,m=1,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为距离向的第t个快时刻,T为步骤1中初始化得到的距离向快时刻总数。在线阵SAR实际成像中,s(t,l,n)可由数据接收机提供;在仿真过程中,s(t,l,n)为观测场景目标空间Ω中所有单元格回波的总数,采用传统的合成孔径雷达原始回波仿真方法产生得到,线阵SAR回波信号s(t,l,n)的近似表示公式表示为 s ( t , l , n ) = Σ m = 1 M α m · exp [ - j · 2 · π · f c · τ n m ( l ) ] exp { j · π · f d r · [ t - τ n m ( l ) ] 2 } , 其中Σ(·)表示求和运算符号,exp(·)表示e指数运算符号,fc为步骤1初始化得到的雷达工作中心频率,fdr为步骤1初始化得到的发射信号调频斜率,αm为步骤2初始化得到的观测场景目标散射系数α中第m个元素的散射系数,t为距离向的第t个快时刻,j为虚数单位,π为圆周率。
将所有线阵SAR原始回波信号s(t,l,n)按顺序排列组成向量,记为回波信号向量S,回波信号向量S由O行1列组成,其中O=T·K·N,T为步骤1中初始化得到的距离向快时刻总数,K为步骤1初始化得到的方位向慢时刻总数,N为步骤1初始化得到的线阵天线阵元总数。
采用公式φi(m)=exp[-j·2·π·fc·τnm(l)]exp{j·π·fdr·[t-τnm(l)]2},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)为观测场景目标空间Ω中第1个单元格在回波信号向量S第1个元素信号对应的时延函数,φ1(2)为观测场景目标空间Ω中第2个单元格在回波信号向量S第1个元素信号对应的时延函数,φ1(M)为观测场景目标空间Ω中第M个单元格在回波信号向量S第1个元素信号对应的时延函数,φ2(1)为观测场景目标空间Ω中第1个单元格在回波信号向量S第2个元素信号对应的时延函数,φ2(2)为观测场景目标空间Ω中第2个单元格在回波信号向量S第2个元素信号对应的时延函数,φ2(M)为观测场景目标空间Ω中第M个单元格在回波信号向量S第2个元素信号对应的时延函数,φO(1)为观测场景目标空间Ω中第1个单元格在回波信号向量S第O个元素信号对应的时延函数,φO(2)为观测场景目标空间Ω中第2个单元格在回波信号向量S第O个元素信号对应的时延函数,φO(M)为观测场景目标空间Ω中第M个单元格在回波信号向量S第O个元素信号对应的时延函数,φ1(1),φ1(2),…,φ1(M)分别为观测场景目标空间Ω中第1,2,…,M个单元格在回波信号向量S第1个元素信号对应的时延函数向量,φ2(1),φ2(2),…,φ2(M)分别为观测场景目标空间Ω中第1,2,…,M个单元格在回波信号向量S第2个元素信号对应的时延函数向量,φO(1),φO(2),…,φO(M)分别为观测场景目标空间Ω中第1,2,…,M个单元格在回波信号向量S第O个元素信号对应的时延函数向量。线阵SAR的线性测量矩阵A为O行M列的二维矩阵。
步骤4、设定迭代最小化稀疏贝叶斯重构算法的初始参数:
初始化传统的迭代最小化稀疏贝叶斯重构算法的参数包括:算法重构迭代处理的最大迭代次数,记做MaxIter;重构残余误差门限,记做ε0;相邻迭代的目标散射系数变化率门限,记为ω0;观测场景目标空间中散射系数的复指数分布参数,记为η;观测场景目标空间中散射系数复指数分布参数η的初始迭代值,记为η(0);观测场景目标空间中散射系数复指数分布参数η的取值区间,记为[ηminmax],其中ηmin为参数区间的最小值,ηmax为参数区间的最大值;观测场景目标空间中散射系数向量范数的参数,记为p,p的值在0至1之间;观测场景目标空间中散射系数向量范数的平滑因子,记为c0;线阵SAR原始回波中噪声方差的初始迭代值记为β(0),一般β(0)的值选择为β(0)=S,观测场景目标空间中散射系数向量α的初始迭代值记为α(0),一般α(0)的值选择为α(0)=0或者α(0)=AHS,其中A为步骤3中得到的线阵SAR测量矩阵,S为步骤3中得到的线阵SAR原始回波信号向量,上标H表示共轭转置运算符号;重构残余误差的初始迭代值记为r(0),一般r(0)的值选择为r(0)=S;观测场景目标散射系数的对角矩阵,记为Λ,矩阵Λ的初始迭代值记为Λ(0),Λ(0)的值为其中diag(·)表示利用向量构造对角矩阵的运算符号,p为初始化得到的散射系数向量范数参数,α(0)为初始化得到的观测场景目标空间中的散射系数向量,c0为初始化得到的观测场景目标空间中散射系数向量范数平滑因子,|·|为取绝对值运算符号,k表示迭代最小化稀疏贝叶斯重构算法中的第k迭代次数,k的初始值设置为k=0,且k的取值范围是从0到MaxIter。
步骤5、设定原始回波信号和目标散射系数的概率分布以及构建算法的重构代价函数:
假设线阵SAR原始回波信号S的概率密度函数服从复高斯随机分布其中f(S|α,β)表示在观测场景目标散射系数和噪声方差条件下线阵SAR原始回波信号的概率密度函数,表示均值向量为Αα、协方差矩阵为βI的复高斯随机分布函数,S为步骤3中得到的线阵SAR原始回波信号向量,A为步骤3中得到的线阵SAR测量矩阵,α为待估计的线阵SAR观测场景目标空间中的散射系数向量,β为待估计的线阵SAR原始回波中的噪声方差,I为单位矩阵,符号~表示服从概率分布符号;
假设线阵SAR观测场景目标空间中每一个分辨单元的散射系数独立同分布且服从复指数分布,散射系数的先验概率密度函数表示为其中f(α)表示观测场景目标空间中散射系数α的先验概率密度函数,αm表示观测场景目标空间中散射系数α第m个元素的散射系数值,η(k-1)为迭代最小化稀疏贝叶斯重构算法第k-1次迭代过程中得到的观测场景目标空间中散射系数的指数分布参数,若k=1时,η(k-1)的值为步骤4中初始化得到的η(0),否则η(k-1)的值利用步骤9中指数分布最优估计得到,p为步骤4初始化得到的散射系数向量范数参数,exp(·)表示指数运算符号,表示元素1至元素M相乘运算符号。根据贝叶斯准则以及f(S|α,β)和f(α),采用似然函数构造重构代价函数其中J(α,β)为线阵SAR观测场景空间中散射系数向量α和原始回波信号中噪声方差β的重构代价函数,O为步骤3得到的线阵SAR原始回波信号向量S的维数,α为待估计的线阵SAR观测场景目标空间中的散射系数向量,β为待估计的线阵SAR原始回波中的噪声方差,η(k-1)为迭代最小化稀疏贝叶斯重构算法第k-1次迭代过程中得到的观测场景目标空间中散射系数的指数分布参数,若k=1时,η(k-1)的值为步骤4中初始化得到的η(0),否则η(k-1)的值利用步骤9中指数分布最优估计得到,A为步骤3中得到的线阵SAR测量矩阵,p为步骤4初始化得到的散射系数向量范数参数,||·||p表示向量的LP范数,表示向量L2范数的平方,ln表示以e为底数的对数运算符号,e的值约等于2.71828。
步骤6、观测场景目标空间的散射系数向量估计:
采用公式 α ( k ) = arg min α J ( α , β ( k - 1 ) ) ≈ ( A H A + η ( k - 1 ) β ( k - 1 ) Λ ( k - 1 ) ) - 1 A H S 和定义10中传统共轭梯度算法计算,得到迭代最小化稀疏贝叶斯重构算法第k次迭代过程中观测场景目标空间的散射系数向量,记为α(k),其中表示求取满足括号中最小值时对应自变量α的最优值,J(α,β(k-1))为重构代价函数J(α,β)中噪声方差β的值为β(k-1)时得到的代价函数,J(α,β)步骤5中得到的重构代价函数,A为步骤3中得到的线阵SAR测量矩阵,S为步骤3中得到的线阵SAR原始回波信号向量,η(k-1)为迭代最小化稀疏贝叶斯重构算法第k-1次迭代过程中得到的观测场景目标空间中散射系数的指数分布参数,若k=1时,η(k-1)的值为步骤4中初始化得到的η(0),否则η(k-1)的值利用步骤9中指数分布最优估计得到,Λ(k-1)为迭代最小化稀疏贝叶斯重构算法第k-1次迭代过程中得到的观测场景目标散射系数的对角矩阵,若k=1时,Λ(k-1)的值为步骤4中初始化得到的Λ(0),否则Λ(k-1)的值利用步骤8中观测场景目标散射系数的对角矩阵计算得到,β(k-1)为迭代最小化稀疏贝叶斯重构算法第k-1次迭代过程中得到的噪声方差,若k=1时,β(k-1)的值为步骤4中初始化得到的β(0),否则β(k-1)的值利用步骤7中观测场景目标散射系数的对角矩阵计算得到,上标H表示共轭转置运算符号,上标-1表示矩阵求逆运算符号。
步骤7、原始回波信号中噪声方差估计:
采用公式计算得到迭代最小化稀疏贝叶斯重构算法第k次迭代过程中的噪声方差,记为β(k),其中表示求取满足括号中最小值时对应自变量β的最优值,J(α(k),β)为步骤5中代价函数J(α,β)中散射系数向量α的值为α(k)时得到的代价函数,A为步骤3中得到的线阵SAR测量矩阵,S为步骤3中得到的线阵SAR原始回波信号向量,O为步骤3中初始化得到的原始回波信号向量S的维数,α(k)为步骤6中得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中观测场景目标空间的散射系数向量,表示向量L2范数的平方,k表示迭代最小化稀疏贝叶斯重构算法中的第k迭代次数。
步骤8、计算观测场景目标散射系数的对角矩阵:
采用公式计算得到迭代最小化稀疏贝叶斯重构算法第k次迭代过程中观测场景目标散射系数的对角矩阵,记为Λ(k),其中diag(·)表示利用向量构造对角矩阵的运算符号,p为步骤4中初始化得到的范数参数,α(k)为步骤6得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中观测场景目标的散射系数向量,c0为步骤4中初始化得到的范数平滑因子;k表示迭代最小化稀疏贝叶斯重构算法中的第k迭代次数。
步骤9、观测场景目标空间中散射系数复指数分布参数估计:
采用公式 η ( k ) = arg min η min ≤ η ≤ η max ( Oβ ( k ) + | | S - Aα ( k ) | | 2 2 + 2 β ( k ) t r a c e ( ( A H A + ηβ ( k ) Λ ( k ) ) - 1 A H ) ) 和定义11中传统Stein无偏风险估计量方法估计,得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中的复指数分布参数,记为η(k),其中表示在[ηminmax]区间上求取满足括号中最小值时对应自变量η的最优值,ηmin和ηmax分别为步骤4中初始化得到的目标散射系数复指数分布参数区间的最小值和最大值,A为步骤3中得到的线阵SAR测量矩阵,S为步骤3中得到的线阵SAR原始回波信号向量,O为步骤3中初始化得到的原始回波信号向量S的维数,α(k)为步骤6中得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中观测场景目标空间的散射系数向量,β(k)为步骤7中得到的第k次迭代过程中的噪声方差,Λ(k)为步骤8得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中的观测场景目标散射系数的对角矩阵,表示向量L2范数的平方,trace(·)为矩阵对角求和运算符号,上标H表示共轭转置运算符号;上标-1表示矩阵求逆运算符号,k表示迭代最小化稀疏贝叶斯重构算法中的第k迭代次数。
步骤10、计算重构残余误差和相邻迭代目标散射系数变化率:
采用公式r(k)=S-Αα(k)计算得到迭代最小化稀疏贝叶斯重构算法第k次迭代过程中的重构残余误差,记为r(k),其中S为步骤3中得到的线阵SAR回波信号向量,A为步骤3中得到的线阵SAR测量矩阵,α(k)为步骤6中得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中观测场景目标的散射系数向量。
采用公式计算得到第k次迭代过程中的相邻目标散射系数变化率,记为ω(k),其中α(k)为迭代最小化稀疏贝叶斯重构算法第k次迭代得到观测场景目标空间散射系数向量,α(k-1)为迭代最小化稀疏贝叶斯重构算法第k-1次迭代得到观测场景目标空间散射系数向量,||·||2为向量L2范数,k表示迭代最小化稀疏贝叶斯重构算法中的第k迭代次数。
步骤11、迭代终止判定:
如果若r(k)≥ε0且ω(k)≥ω0且k≤MaxIter,则k←k+1,执行步骤5至11,否则终止迭代,此刻第k次迭代得到的散射系数向量值α(k)即为观测场景目标空间Ω最终的散射系数向量,其中r(k)为步骤10中得到第k次迭代重构残余误差,ε0为步骤4中初始化得到的重构残余误差的门限,ω(k)为步骤10中得到第k次迭代过程中的相邻目标散射系数变化率,ω0为步骤4中初始化得到的相邻迭代目标散射系数变化率的门限,k表示迭代最小化稀疏贝叶斯重构算法中的第k迭代次数,MaxIter为步骤4中初始化得到的算法重构处理的最大迭代次数。
最后将观测场景目标空间散射系数向量α(k)转换成三维矩阵形式,得到线阵SAR观测场景目标空间Ω的三维成像结果。
本发明的创新点在于针对线阵SAR观测场景目标空间中散射目标稀疏的特征,提供了基于迭代最小化稀疏贝叶斯重构的线阵SAR稀疏成像方法,该算法基于线阵SAR原始回波信号测量模型的先验分布假设,假设线阵SAR观测场景目标空间中散射系数的先验概率密度函数服从复指数先验分布,线阵SAR原始回波信号的后验概率密度函数服从高斯随机分布,再利用贝叶斯准则构造线阵SAR观测场景目标空间中散射系数服的重构代价函数,通过复指数分布参数最优化估计以及迭代最小化重构代价函数实现线阵SAR观测场景目标空间中的散射系数稀疏重构。根据本人了解,当前还没有出现基于迭代最小化稀疏贝叶斯重构的线阵SAR成像方法。
本发明的优点在于利用线阵SAR原始回波信号观测模型的先验概率分布假设,利用贝叶斯准则和似然函数可以更加合理地构造适用于线阵SAR观测场景目标空间中散射系数服的重构代价函数,从而提高了线阵SAR稀疏成像的性能。本发明可以应用于合成孔径雷达成像和地球遥感等领域。
附图说明:
图1为线阵SAR成像几何关系图
其中,长平行四边形框表示线阵,黑点表示线阵阵元,线阵天线长度为L,K为方位向慢时刻总数,l为方位向第l个慢时刻,PRI表示线阵SAR发射信号的脉冲重复时间,为线阵天线中第n个阵元在方位向第l个慢时刻的位置矢量,xn(l)、yn(l)和zn(l)分别表示线阵天线中第n个阵元在方位向第l个慢时刻的水平横向、水平纵向和高度向坐标;表示观测场景目标空间中第m单元格的位置向量,为在方位向第l个慢时刻时观测场景目标空间中第m个单元格到线阵天线第n阵元的距离,x、y和z分别表示观测场景目标空间中水平横向、水平纵向和高度向坐标,0表示观测场景目标空间中坐标原点。
图2为本发明所提供的基于迭代最小化稀疏贝叶斯重构的线阵SAR稀疏成像方法的处理流程示意框图
图3为本发明具体实施方式采用的线阵SAR系统仿真参数表
具体实施方式
本发明主要采用仿真实验的方法进行验证,所有步骤和结论都在MATLABR2008b软件上验证正确。具体实施步骤如下:
步骤1、初始化仿真所需的线阵SAR系统参数:
初始化线阵SAR系统参数值如图3所示,包括:运动平台速度矢量线阵天线的阵元总数N=201,线阵天线各个阵元的初始位置矢量其中n为线阵天线第n个阵元的序号,n=1,2,…,N,N=201,线阵天线长度L=3m,线阵天线相邻阵元之间的间距d=0.015m,雷达中心频率fc=10GHz,雷达发射机基带信号的信号宽度Br=300MHz,雷达发射信号脉冲宽度TP=10-6s,雷达发射信号调频斜率fdr=3×1014Hz/s,雷达接收系统的采样频率fs=500MHz,雷达系统发射信号的脉冲重复频率PRF=500Hz,发射信号脉冲重复时间线阵天线在方位向的有效孔径长度Da=1m,光在空气中的传播速度C=3×108m/s,距离向快时刻总数T=256,距离向快时刻序列t=1,2,…,T,其中t为距离向第t个快时刻,方位向慢时刻总数K=256,方位向慢时刻序列l=1,2,…,256,其中l为方位向第l个慢时刻。
步骤2、初始化线阵SAR的观测场景目标空间参数:
以雷达波束照射场区域地平面和垂直于该地平面向上的单位向量所构成的空间直角坐标作为线阵SAR的观测场景目标空间Ω。初始化观测场景目标空间Ω的大小为128×128×1像素,观测场景目标空间Ω的中心坐标位置位于[0,0,0],每一个单元网格在水平横向、水平纵向和高度向边长为dx=dy=dz=0.5m,计算得到观测场景目标空间的单元格总数M=16384,观测场景目标空间Ω中每一个单元格的位置为其中x′=1,2,…,128,y′=1,2,…,128,m=(x′-1)·128+y′。为观测场景目标空间Ω中第m个单元格的位置矢量,m表示观测场景目标空间Ω中第m个单元格,m=1,2,…,M,M=16384。在观测场景目标空间Ω中加入仿真点目标散射体,点目标散射体数个数为9个,它们散射系数值均为1,坐标位置分别为[0,0,0]、[10,0,0]、[-10,0,0]、[0,10,0]、[0,-10,0]、[10,10,0]、[10,-10,0]、[-10,10,0]、[-10,-10,0],单位均为m;观测场景目标空间Ω中没有包含点目标单元格的散射系数设置为0。将观测场景目标空间Ω中所有单元格的目标散射系数按单元格位置顺序排列组成散射系数向量α。确定观测场景目标空间Ω所有单元散射系数后,散射系数向量α在线阵SAR三维成像观测仿真过程中就已经确定。场景目标散射系数向量α由M行1列组成,αm为向量α中对应场景目标空间Ω中第m个单元格的散射系数值。在本仿真观测场景目标空间Ω中,只有包含点散射目标的9个单元格散射系数值α设置为1,其余单元格的散射系数都为0。利用计算机合成孔径雷达原始回波仿真方法产生线阵SAR的原始回波数据。
步骤3、建立线阵SAR原始回波信号与场景目标散射系数的线性测量矩阵:
采用公式计算得到线阵天线第n个阵元在第l个方位向慢时刻的位置矢量其中n表示天线第n个阵元序号,n=1,2,…,N,N为步骤1初始化得到的线阵阵元总数N=201,l表示方位向第l个慢时刻序号,l=1,2,…,K,K为步骤1初始化得到的线阵SAR方位向慢时刻总数K=256,为步骤1初始化得到的运动平台初始位置矢量 为步骤1初始化得到的运动平台速度矢量PRF为步骤1初始化得到的脉冲重复频率PRF=500Hz。采用公式计算得到在方位向第l个慢时刻线阵SAR场景目标空间Ω中第m个单元格到线阵天线第n个阵元的距离其中n=1,2,…,N,N=201,l=1,2,…,K,K=256,m表示场景目标空间Ω中第m个单元格,m=1,2,…,M,M为步骤2初始化得到的观测场景目标空间总单元格数M=16384,||·||2表示向量L2范数,为步骤2中初始化得到观测场景目标空间Ω中第m个单元格的坐标矢量。采用公式计算得到在方位向第l个慢时刻线阵SAR场景目标空间Ω中第m个单元格到第n个线阵天线阵元的时间延时τnm(l),其中n=1,2,…,N,l=1,2,…,K,m=1,2,…,M,C为步骤1初始化得到的光在空气中的传播速度C=3×108m/s。采用公式 s ( t , l , n ) = Σ m = 1 M α m · exp [ - j · 2 · π · f c · τ n m ( l ) ] exp { j · π · f d r · [ t - τ n m ( l ) ] 2 } 和合成孔径雷达原始回波仿真方法得到观测场景目标空间Ω的近似线阵SAR原始信号回波s(t,l,n),其中t表示距离向快时刻序号,t=1,2,…,T,T=256,l=1,2,…,K,K=256,n=1,2,…,N,N=201,Σ(·)表示求和运算符号,exp(·)表示e指数运算符号,雷达中心频率fc=10GHz,雷达发射信号调频斜率fdr=3×1014Hz/s,αm为步骤2初始化得到的场景目标空间中第m个单元格的散射系数,j为虚数单位(即-1的开根值),π=3.14159。将所有线阵SAR原始回波信号s(t,l,n)按位置序号顺序排列组成回波信号向量S,回波信号向量S由O行1列组成,其中O=T·K·N=13172736。采用公式φi(m)=exp[-j·2·π·fc·τnm(l)]exp{j·π·fdr·[t-τnm(l)]2},t=1,2,…,T,T=256,l=1,2,…,K,K=256,n=1,2,…,N,N=201,m=1,2,…,M,M=16384,i=[(t-1)K+l-1]N+n,计算得到观测场景目标空间Ω中第m个单元格在回波信号向量S第i个元素信号对应的时延函数φi(m),其中i表示回波信号向量S中第i个元素序号,i的取值范围为i=1,2,…,O,O=13172736。
采用矩阵表达公式
计算得到线阵SAR原始回波信号与场景目标空间所有单元格的线性测量矩阵A,其中,φ1(1)为观测场景目标空间Ω中第1个单元格在回波信号向量S第1个元素信号对应的时延函数,φ1(2)为观测场景目标空间Ω中第2个单元格在回波信号向量S第1个元素信号对应的时延函数,φ1(M)为观测场景目标空间Ω中第M个单元格在回波信号向量S第1个元素信号对应的时延函数,φ2(1)为观测场景目标空间Ω中第1个单元格在回波信号向量S第2个元素信号对应的时延函数,φ2(2)为观测场景目标空间Ω中第2个单元格在回波信号向量S第2个元素信号对应的时延函数,φ2(M)为观测场景目标空间Ω中第M个单元格在回波信号向量S第2个元素信号对应的时延函数,φO(1)为观测场景目标空间Ω中第1个单元格在回波信号向量S第O个元素信号对应的时延函数,φO(2)为观测场景目标空间Ω中第2个单元格在回波信号向量S第O个元素信号对应的时延函数,φO(M)为观测场景目标空间Ω中第M个单元格在回波信号向量S第O个元素信号对应的时延函数,φ1(1),φ1(2),…,φ1(M)分别为观测场景目标空间Ω中第1,2,…,M个单元格在回波信号向量S第1个元素信号对应的时延函数向量,φ2(1),φ2(2),…,φ2(M)分别为观测场景目标空间Ω中第1,2,…,M个单元格在回波信号向量S第2个元素信号对应的时延函数向量,φO(1),φO(2),…,φO(M)分别为观测场景目标空间Ω中第1,2,…,M个单元格在回波信号向量S第O个元素信号对应的时延函数向量。线阵SAR的线性测量矩阵A为O行M列的二维矩阵,其中O=13172736,M=16384。
步骤4、设定迭代最小化稀疏贝叶斯重构算法的初始参数:
初始化迭代最小化稀疏贝叶斯重构算法的参数包括:算法重构迭代处理的最大迭代次数MaxIter=500,重构残余误差门限ε0=0.0001,相邻迭代的观测场景目标散射系数变化率门限ω0=0.1,观测场景目标散射系数复指数分布的参数η初始迭代值η(0)=1,复指数分布参数的取值区间[ηminmax],其中参数区间的最小值ηmin=0.001,参数区间的最大值ηmax=100;观测场景目标空间中散射系数向量范数的参数p=0.9,观测场景目标空间中散射系数向量范数的平滑因子c0=0.000001,线阵SAR原始回波中的噪声方差的初始迭代值β(0)=S,观测场景目标散射系数向量α的初始迭代值α(0)=AHS,其中A为步骤3中得到的线阵SAR测量矩阵,S为步骤3中得到的线阵SAR原始回波信号向量,上标H表示共轭转置运算符号;重构残余误差的初始迭代值r(0)=S;观测场景目标散射系数的对角矩阵Λ初始迭代值其中diag(·)表示利用向量构造对角矩阵的运算符号,p为初始化的范数参数p=0.9,α(0)为初始化的观测场景目标散射系数向量,c0为初始化的范数平滑因子c0=0.000001;k表示迭代最小化稀疏贝叶斯重构算法中的第k迭代次数,k的初始值设置为k=0,且k的取值范围是从0到MaxIter,MaxIter=500。
步骤5、设定原始回波信号和目标散射系数的概率分布以及构建算法的重构代价函数:
假设线阵SAR原始回波信号S的概率密度函数服从复高斯随机分布其中f(S|α,β)表示在观测场景目标散射系数和噪声方差条件下线阵SAR原始回波信号的概率密度函数,表示均值向量为Αα、协方差矩阵为βI的复高斯随机分布函数,S为步骤3中得到的线阵SAR原始回波信号向量,A为步骤3中得到的线阵SAR测量矩阵,α为待估计的线阵SAR观测场景目标空间中的散射系数向量,β为待估计的线阵SAR原始回波中的噪声方差,I为单位矩阵,符号~表示服从概率分布符号;假设线阵SAR观测场景目标空间中每一个分辨单元的散射系数独立同分布且服从复指数分布,散射系数的先验概率密度函数表示为其中f(α)表示观测场景目标空间中散射系数α的先验概率密度函数,αm表示观测场景目标空间中散射系数α第m个元素的散射系数值,η(k-1)为迭代最小化稀疏贝叶斯重构算法第k-1次迭代过程中得到的观测场景目标空间中散射系数的指数分布参数,若k=1时,η(k-1)的值为步骤4中初始化得到的η(0),否则η(k-1)的值利用步骤9中指数分布最优估计得到,p为步骤4初始化得到的散射系数向量范数参数p=0.9,exp(·)表示指数运算符号,表示元素1至元素M相乘运算符号。根据贝叶斯准则以及f(S|α,β)和f(α),利用似然函数构造重构代价函数其中J(α,β)为线阵SAR观测场景空间中散射系数向量α和原始回波信号中噪声方差β的重构代价函数,O为步骤3得到的线阵SAR原始回波信号向量S的维数,α为待估计的线阵SAR观测场景目标空间中的散射系数向量,β为待估计的线阵SAR原始回波中的噪声方差,η(k-1)为迭代最小化稀疏贝叶斯重构算法第k-1次迭代过程中得到的观测场景目标空间中散射系数的指数分布参数,若k=1时,η(k-1)的值为步骤4中初始化得到的η(0)=1,否则η(k-1)的值利用步骤9中指数分布最优估计得到,A为步骤3中得到的线阵SAR测量矩阵,||·||p表示向量的LP范数,表示向量L2范数的平方,ln表示以e为底数的对数运算符号,e的值约等于2.71828。
步骤6、观测场景目标空间的散射系数向量估计:
采用公式 α ( k ) = arg min α J ( α , β ( k - 1 ) ) ≈ ( A H A + η ( k - 1 ) β ( k - 1 ) Λ ( k - 1 ) ) - 1 A H S 和定义10中传统共轭梯度算法计算得到迭代最小化稀疏贝叶斯重构算法第k次迭代过程中观测场景目标空间的散射系数向量,记为α(k);其中表示求取满足括号中最小值时对应自变量α的最优值,J(α,β(k-1))为代价函数J(α,β)中噪声方差β的值为β(k-1)时得到的重构代价函数,J(α,β)为步骤5中得到的重构代价函数,A为步骤3中得到的线阵SAR测量矩阵,S为步骤3中得到的线阵SAR原始回波信号向量,η(k-1)为迭代最小化稀疏贝叶斯重构算法第k-1次迭代过程中得到的观测场景目标空间中散射系数的指数分布参数,若k=1时,η(k-1)的值为步骤4中初始化得到的η(0),否则η(k-1)的值利用步骤9中迭代最小化稀疏贝叶斯重构算法第k-1迭代过程的指数分布最优估计得到,Λ(k-1)为算法第k-1次迭代过程中得到的观测场景目标散射系数的对角矩阵,若k=1时,Λ(k-1)的值为步骤4中初始化得到的Λ(0),否则Λ(k-1)的值利用步骤8中迭代最小化稀疏贝叶斯重构算法第k-1迭代过程的观测场景目标散射系数的对角矩阵计算得到,β(k-1)为算法第k-1次迭代过程中得到的噪声方差,若k=1时,β(k-1)的值为步骤4中初始化得到的β(0),否则β(k-1)的值利用步骤7中迭代最小化稀疏贝叶斯重构算法第k-1迭代过程的噪声方差估计得到,上标H表示共轭转置运算符号,上标-1表示矩阵求逆运算符号。
步骤7、原始回波信号中噪声方差估计:
采用公式计算得到迭代最小化稀疏贝叶斯重构算法第k次迭代过程中的噪声方差,记为β(k),其中表示求取满足括号中最小值时对应自变量β的最优值,J(α(k),β)为代价函数J(α,β)中散射系数向量α的值为α(k)时得到的重构代价函数,J(α,β)为步骤5中得到的重构代价函数,A为步骤3中得到的线阵SAR测量矩阵,S为步骤3中得到的线阵SAR原始回波信号向量,O为步骤3中初始化得到的原始回波信号向量S的维数O=13172736,α(k)为步骤6中得到的观测场景目标空间的散射系数向量,表示向量L2范数的平方,k表示迭代最小化稀疏贝叶斯重构算法中的第k迭代次数。
步骤8、计算观测场景目标散射系数的对角矩阵:
采用公式计算得到迭代最小化稀疏贝叶斯重构算法第k次迭代过程中观测场景目标散射系数的对角矩阵,记为Λ(k),其中diag(·)表示利用向量构造对角矩阵的运算符号,p为步骤4中初始化得到的范数参数p=0.9,α(k)为步骤6得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中观测场景目标的散射系数向量,c0为步骤4中初始化得到的范数平滑因子c0=0.000001,|·|为取绝对值运算符号,k表示迭代最小化稀疏贝叶斯重构算法中的第k迭代次数。
步骤9、观测场景目标空间中散射系数复指数分布参数估计:
采用公式 η ( k ) = arg min η min ≤ η ≤ η max ( Oβ ( k ) + | | S - Aα ( k ) | | 2 2 + 2 β ( k ) t r a c e ( ( A H A + ηβ ( k ) Λ ( k ) ) - 1 A H ) ) 和定义11中传统Stein无偏风险估计量方法估计得到迭代最小化稀疏贝叶斯重构算法第k次迭代过程中的指数分布参数,记为η(k),其中表示在[ηminmax]区间上求取满足括号中最小值时对应自变量η的最优值,其中ηmin和ηmax分别为步骤4中初始化得到的目标散射系数复指数分布参数区间的最小值和最大值,ηmin=0.001,ηmax=100;A为步骤3中得到的线阵SAR测量矩阵,S为步骤3中得到的线阵SAR原始回波信号向量,O为步骤3中初始化得到的线阵SAR原始回波信号向量S的维数O=13172736,α(k)为步骤6中得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中观测场景目标空间的散射系数向量,β(k)为步骤7中得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中的噪声方差,Λ(k)为步骤8得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中的观测场景目标散射系数的对角矩阵,表示向量L2范数的平方,trace(·)为矩阵对角求和运算符号,上标H表示共轭转置运算符号;上标-1表示矩阵求逆运算符号,k表示迭代最小化稀疏贝叶斯重构算法中的第k迭代次数。
步骤10、计算重构残余误差和相邻迭代目标散射系数变化率:
采用公式r(k)=S-Αα(k)计算得到迭代最小化稀疏贝叶斯重构算法第k次迭代过程中的重构残余误差,记为r(k),其中S为步骤3中得到的线阵SAR回波信号向量,A为步骤3中得到的线阵SAR测量矩阵,α(k)为步骤6中得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中观测场景目标的散射系数向量。采用公式计算得到迭代最小化稀疏贝叶斯重构算法第k次迭代过程中的相邻目标散射系数变化率,记为ω(k),其中α(k)为迭代最小化稀疏贝叶斯重构算法第k次迭代得到观测场景目标空间散射系数向量,α(k-1)为迭代最小化稀疏贝叶斯重构算法第k-1次迭代得到观测场景目标空间散射系数向量,||·||2为向量L2范数。
步骤11、迭代终止判定:
如果若r(k)≥ε0且ω(k)≥ω0且k≤MaxIter,则k←k+1,执行步骤5至11,否则迭代最小化稀疏贝叶斯重构算法终止迭代过程,此刻第k次迭代得到的散射系数向量值α(k)即为线阵SAR观测场景目标空间Ω最终的散射系数向量,其中r(k)为步骤10中得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中的重构残余误差,k为迭代最小化稀疏贝叶斯重构算法第k迭代次数,ε0为步骤4中初始化得到的重构残余误差的门限ε0=0.0001,ω(k)为步骤10中得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中的相邻目标散射系数变化率,ω0为步骤4中初始化得到的相邻迭代目标散射系数变化率的门限ω0=0.1,k表示迭代最小化稀疏贝叶斯重构算法中的第k迭代次数,MaxIter为步骤4中初始化得到的算法重构处理的最大迭代次数MaxIter=500。最后将观测场景目标空间散射系数向量α(k)转换成三维矩阵形式,得到线阵SAR观测场景目标空间Ω的三维成像结果。
通过本发明具体实施方式可以看出,本发明通过建立线阵SAR原始回波信号与场景目标空间散射系数的线性测量模型,利用线阵SAR原始回波信号测量模型的先验分布,将线阵SAR稀疏成像过程转换成为稀疏贝叶斯重构的稀疏求解过程。
本发明提供了基于迭代最小化稀疏贝叶斯重构的线阵SAR稀疏成像方法,该算法基于线阵SAR原始回波信号测量模型的先验分布假设,假设线阵SAR观测场景目标空间中散射系数的先验概率密度函数服从复指数先验分布,线阵SAR原始回波信号的后验概率密度函数服从高斯随机分布,再利用贝叶斯准则构造线阵SAR的重构代价函数,通过复指数分布参数最优化估计以及迭代最小化重构代价函数实现线阵SAR观测场景目标空间中的散射系数稀疏重构,提高了SAR稀疏成像的质量。

Claims (1)

1.一种基于迭代最小化稀疏贝叶斯重构线阵SAR成像方法,其特征是它包括以下步骤:
步骤1、初始化线阵SAR系统参数:
初始化线阵SAR系统参数包括:平台速度矢量,记做线阵天线各阵元初始位置矢量,记做其中n为线阵天线中第n个阵元的序号,n为自然数,n=1,2,...,N,N为线阵天线的阵元总数;线阵天线长度,记做L;雷达工作中心频率,记做fc;雷达载频波长,记做λ;雷达发射基带信号的信号带宽,记做Br;雷达发射信号脉冲宽度,记做TP;雷达发射信号的调频斜率,记做fdr;雷达接收波门持续宽度,记做To;雷达接收系统的采样频率,记做fs;雷达发射系统的脉冲重复频率,记做PRF;雷达系统的脉冲重复时间,记为PRI;雷达接收系统接收波门相对于发射信号发散波门的延迟,记做TD;天线在方位向的有效孔径长度,记做Da;光在空气中的传播速度,记做C;距离向快时刻,记做t,t为自然数,t=1,2,…,T,T为距离向快时刻总数;方位向慢时刻,记做l,l为自然数,l=1,2,…,K,K为方位向慢时刻总数;上述参数均为线阵SAR系统标准参数,其中线阵天线的阵元总数N,线阵天线长度L,相邻天线阵元之间的间距d,雷达中心频率fc,雷达载频波长λ,雷达发射基带信号的信号带宽Br,雷达发射信号脉冲宽度TP,雷达发射信号调频斜率fdr,雷达接收波门持续宽度To,雷达接收系统的采样频率fs,雷达系统的脉冲重复频率PRF,雷达系统的脉冲重复时间PRI和雷达接收系统接收波门相对于发射信号发散波门的延迟TD,天线在方位向的有效孔径长度Da在线阵SAR系统设计过程中已经确定;平台速度矢量及线阵天线各阵元初始位置矢量在线阵SAR观测方案设计中已经确定;根据线阵SAR成像系统方案和观测方案,线阵SAR成像方法需要的初始化系统参数均为已知;
步骤2、初始化线阵SAR观测场景目标空间参数:
初始化线阵SAR的观测场景目标空间参数,包括:以雷达波束照射场区域地平面和垂直于该地平面向上的单位向量所构成的空间直角坐标作为线阵SAR的观测场景目标空间Ω;将观测场景目标空间Ω均匀划分成大小相等的立体单元格,亦称为分辨单元,单元网格在水平横向、水平纵向和高度向边长分别记为dx、dy和dz,单元格大小选择为线阵SAR系统传统理论成像分辨率或该分辨率的二分之一;观测场景目标空间Ω中第m个单元格的坐标矢量,记做m表示观测场景目标空间Ω中第m个单元格,m为自然数,m=1,2,…,M,M为观测场景目标空间Ω中的单元格总数;观测场景目标空间Ω中所有单元格的散射系数按位置顺序排列组成向量,记做α,向量α由M行1列组成;散射系数向量α中第m个元素的散射系数,记做αm;观测场景目标空间Ω在线阵SAR成像方案设计中已经确定;
步骤3、建立线阵SAR原始回波信号与观测场景目标散射系数的线性测量矩阵:
采用公式n=1,2,…,N,l=1,2,…,K,计算得到第n个线阵天线阵元在第l个方位向慢时刻的位置矢量,记为其中为步骤1初始化得到线阵天线各阵元的初始位置矢量,为步骤1中初始化得到的平台速度矢量,PRF为步骤1初始化得到的雷达系统脉冲重复频率,n为线阵天线中第n个阵元的序号,n=1,2,…,N,N为步骤1中初始化得到的线阵天线阵元总数,l为方位向慢时刻中第l个慢时刻序号,l=1,2,…,K,K为步骤1初始化得到的线阵SAR方位向慢时刻总数;
采用公式n=1,2,…,N,l=1,2,…,K,m=1,2,…,M,计算得到在方位向第l个慢时刻线阵SAR观测场景目标空间Ω中第m个单元格到第n个线阵天线阵元的距离,记为其中||·||2表示向量的L2范数,为步骤2中初始化得到观测场景目标空间Ω中第m个单元格的坐标矢量,m表示观测场景目标空间Ω中第m个单元格,m=1,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为距离向的第t个快时刻,T为步骤1中初始化得到的距离向快时刻总数;在线阵SAR实际成像中,s(t,l,n)由数据接收机提供;在仿真过程中,s(t,l,n)为观测场景目标空间Ω中所有单元格回波的总数,采用传统的合成孔径雷达原始回波仿真方法产生得到,线阵SAR回波信号s(t,l,n)的近似表示公式表示为 s ( t , l , n ) = Σ m = 1 M α m · exp [ - j · 2 · π · f c · τ n m ( l ) ] exp { j · π · f d r . [ t - τ n m ( l ) ] 2 } , 其中Σ(·)表示求和运算符号,exp(·)表示e指数运算符号,fc为步骤1初始化得到的雷达工作中心频率,fdr为步骤1初始化得到的发射信号调频斜率,αm为步骤2初始化得到的观测场景目标散射系数α中第m个元素的散射系数,t为距离向的第t个快时刻,j为虚数单位,π为圆周率;
将所有线阵SAR原始回波信号s(t,l,n)按顺序排列组成向量,记为回波信号向量S,回波信号向量S由O行1列组成,其中O=T·K·N,T为步骤1中初始化得到的距离向快时刻总数,K为步骤1初始化得到的方位向慢时刻总数,N为步骤1初始化得到的线阵天线阵元总数;
采用公式φi(m)=exp[-j·2·π·fc·τnm(l)]exp{j·π·fdr·[t-τnm(l)]2},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)为观测场景目标空间Ω中第1个单元格在回波信号向量S第1个元素信号对应的时延函数,φ1(2)为观测场景目标空间Ω中第2个单元格在回波信号向量S第1个元素信号对应的时延函数,φ1(M)为观测场景目标空间Ω中第M个单元格在回波信号向量S第1个元素信号对应的时延函数,φ2(1)为观测场景目标空间Ω中第1个单元格在回波信号向量S第2个元素信号对应的时延函数,φ2(2)为观测场景目标空间Ω中第2个单元格在回波信号向量S第2个元素信号对应的时延函数,φ2(M)为观测场景目标空间Ω中第M个单元格在回波信号向量S第2个元素信号对应的时延函数,φO(1)为观测场景目标空间Ω中第1个单元格在回波信号向量S第O个元素信号对应的时延函数,φO(2)为观测场景目标空间Ω中第2个单元格在回波信号向量S第O个元素信号对应的时延函数,φO(M)为观测场景目标空间Ω中第M个单元格在回波信号向量S第O个元素信号对应的时延函数,φ1(1),φ1(2),…,φ1(M)分别为观测场景目标空间Ω中第1,2,…,M个单元格在回波信号向量S第1个元素信号对应的时延函数向量,φ2(1),φ2(2),…,φ2(M)分别为观测场景目标空间Ω中第1,2,…,M个单元格在回波信号向量S第2个元素信号对应的时延函数向量,φO(1),φO(2),…,φO(M)分别为观测场景目标空间Ω中第1,2,…,M个单元格在回波信号向量S第O个元素信号对应的时延函数向量;线阵SAR的线性测量矩阵A为O行M列的二维矩阵;
步骤4、设定迭代最小化稀疏贝叶斯重构算法的初始参数:
初始化传统的迭代最小化稀疏贝叶斯重构算法的参数包括:算法重构迭代处理的最大迭代次数,记做MaxIter;重构残余误差门限,记做ε0;相邻迭代的目标散射系数变化率门限,记为ω0;观测场景目标空间中散射系数的复指数分布参数,记为η;观测场景目标空间中散射系数复指数分布参数η的初始迭代值,记为η(0);观测场景目标空间中散射系数复指数分布参数η的取值区间,记为[ηminmax],其中ηmin为参数区间的最小值,ηmax为参数区间的最大值;观测场景目标空间中散射系数向量范数的参数,记为p,p的值在0至1之间;观测场景目标空间中散射系数向量范数的平滑因子,记为c0;线阵SAR原始回波中噪声方差的初始迭代值记为β(0),β(0)的值选择为β(0)=S,观测场景目标空间中散射系数向量α的初始迭代值记为α(0),α(0)的值选择为α(0)=0或者α(0)=AHS,其中A为步骤3中得到的线阵SAR测量矩阵,S为步骤3中得到的线阵SAR原始回波信号向量,上标H表示共轭转置运算符号;重构残余误差的初始迭代值记为r(0),r(0)的值选择为r(0)=S;观测场景目标散射系数的对角矩阵,记为Λ,矩阵Λ的初始迭代值记为Λ(0),Λ(0)的值为其中diag(·)表示利用向量构造对角矩阵的运算符号,p为初始化得到的散射系数向量范数参数,α(0)为初始化得到的观测场景目标空间中的散射系数向量,c0为初始化得到的观测场景目标空间中散射系数向量范数平滑因子,|·|为取绝对值运算符号,k表示迭代最小化稀疏贝叶斯重构算法中的第k迭代次数,k的初始值设置为k=0,且k的取值范围是从0到MaxIter;
步骤5、设定原始回波信号和目标散射系数的概率分布以及构建算法的重构代价函数:
假设线阵SAR原始回波信号S的概率密度函数服从复高斯随机分布其中f(S|α,β)表示在观测场景目标散射系数和噪声方差条件下线阵SAR原始回波信号的概率密度函数,表示均值向量为Αα、协方差矩阵为βI的复高斯随机分布函数,S为步骤3中得到的线阵SAR原始回波信号向量,A为步骤3中得到的线阵SAR测量矩阵,α为待估计的线阵SAR观测场景目标空间中的散射系数向量,β为待估计的线阵SAR原始回波中的噪声方差,I为单位矩阵,符号~表示服从概率分布符号;
假设线阵SAR观测场景目标空间中每一个分辨单元的散射系数独立同分布且服从复指数分布,散射系数的先验概率密度函数表示为其中f(α)表示观测场景目标空间中散射系数α的先验概率密度函数,αm表示观测场景目标空间中散射系数α第m个元素的散射系数值,η(k-1)为迭代最小化稀疏贝叶斯重构算法第k-1次迭代过程中得到的观测场景目标空间中散射系数的指数分布参数,若k=1时,η(k-1)的值为步骤4中初始化得到的η(0),否则η(k-1)的值利用步骤9中指数分布最优估计得到,p为步骤4初始化得到的散射系数向量范数参数,exp(·)表示指数运算符号,表示元素1至元素M相乘运算符号;根据贝叶斯准则以及f(S|α,β)和f(α),采用似然函数构造重构代价函数其中J(α,β)为线阵SAR观测场景空间中散射系数向量α和原始回波信号中噪声方差β的重构代价函数,O为步骤3得到的线阵SAR原始回波信号向量S的维数,α为待估计的线阵SAR观测场景目标空间中的散射系数向量,β为待估计的线阵SAR原始回波中的噪声方差,η(k-1)为迭代最小化稀疏贝叶斯重构算法第k-1次迭代过程中得到的观测场景目标空间中散射系数的指数分布参数,若k=1时,η(k-1)的值为步骤4中初始化得到的η(0),否则η(k-1)的值利用步骤9中指数分布最优估计得到,A为步骤3中得到的线阵SAR测量矩阵,p为步骤4初始化得到的散射系数向量范数参数,||·||p表示向量的LP范数,表示向量L2范数的平方,ln表示以e为底数的对数运算符号,e的值约等于2.71828;
步骤6、观测场景目标空间的散射系数向量估计:
采用公式 α ( k ) = arg m i n α J ( α , β ( k - 1 ) ) ≈ ( A H A + η ( k - 1 ) β ( k - 1 ) Λ ( k - 1 ) ) - 1 A H S 和传统共轭梯度算法计算,得到迭代最小化稀疏贝叶斯重构算法第k次迭代过程中观测场景目标空间的散射系数向量,记为α(k),其中表示求取满足括号中最小值时对应自变量α的最优值,J(α,β(k-1))为重构代价函数J(α,β)中噪声方差β的值为β(k-1)时得到的代价函数,J(α,β)步骤5中得到的重构代价函数,A为步骤3中得到的线阵SAR测量矩阵,S为步骤3中得到的线阵SAR原始回波信号向量,η(k-1)为迭代最小化稀疏贝叶斯重构算法第k-1次迭代过程中得到的观测场景目标空间中散射系数的指数分布参数,若k=1时,η(k-1)的值为步骤4中初始化得到的η(0),否则η(k-1)的值利用步骤9中指数分布最优估计得到,Λ(k-1)为迭代最小化稀疏贝叶斯重构算法第k-1次迭代过程中得到的观测场景目标散射系数的对角矩阵,若k=1时,Λ(k-1)的值为步骤4中初始化得到的Λ(0),否则Λ(k-1)的值利用步骤8中观测场景目标散射系数的对角矩阵计算得到,β(k-1)为迭代最小化稀疏贝叶斯重构算法第k-1次迭代过程中得到的噪声方差,若k=1时,β(k-1)的值为步骤4中初始化得到的β(0),否则β(k-1)的值利用步骤7中观测场景目标散射系数的对角矩阵计算得到,上标H表示共轭转置运算符号,上标-1表示矩阵求逆运算符号;
步骤7、原始回波信号中噪声方差估计:
采用公式计算得到迭代最小化稀疏贝叶斯重构算法第k次迭代过程中的噪声方差,记为β(k),其中表示求取满足括号中最小值时对应自变量β的最优值,J(α(k),β)为步骤5中代价函数J(α,β)中散射系数向量α的值为α(k)时得到的代价函数,A为步骤3中得到的线阵SAR测量矩阵,S为步骤3中得到的线阵SAR原始回波信号向量,O为步骤3中初始化得到的原始回波信号向量S的维数,α(k)为步骤6中得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中观测场景目标空间的散射系数向量,表示向量L2范数的平方,k表示迭代最小化稀疏贝叶斯重构算法中的第k迭代次数;
步骤8、计算观测场景目标散射系数的对角矩阵:
采用公式计算得到迭代最小化稀疏贝叶斯重构算法第k次迭代过程中观测场景目标散射系数的对角矩阵,记为Λ(k),其中diag(·)表示利用向量构造对角矩阵的运算符号,p为步骤4中初始化得到的范数参数,α(k)为步骤6得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中观测场景目标的散射系数向量,c0为步骤4中初始化得到的范数平滑因子;k表示迭代最小化稀疏贝叶斯重构算法中的第k迭代次数;
步骤9、观测场景目标空间中散射系数复指数分布参数估计:
采用公式 η ( k ) = arg min η min ≤ η ≤ η max ( Oβ ( k ) + | | S - Aα ( k ) | | 2 2 + 2 β ( k ) t r a c e ( ( A H A + ηβ ( k ) Λ ( k ) ) - 1 A H ) ) 和传统Stein无偏风险估计量方法估计,得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中的复指数分布参数,记为η(k),其中表示在[ηminmax]区间上求取满足括号中最小值时对应自变量η的最优值,ηmin和ηmax分别为步骤4中初始化得到的目标散射系数复指数分布参数区间的最小值和最大值,A为步骤3中得到的线阵SAR测量矩阵,S为步骤3中得到的线阵SAR原始回波信号向量,O为步骤3中初始化得到的原始回波信号向量S的维数,α(k)为步骤6中得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中观测场景目标空间的散射系数向量,β(k)为步骤7中得到的第k次迭代过程中的噪声方差,Λ(k)为步骤8得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中的观测场景目标散射系数的对角矩阵,表示向量L2范数的平方,trace(·)为矩阵对角求和运算符号,上标H表示共轭转置运算符号;上标-1表示矩阵求逆运算符号,k表示迭代最小化稀疏贝叶斯重构算法中的第k迭代次数;
步骤10、计算重构残余误差和相邻迭代目标散射系数变化率:
采用公式r(k)=S-Αα(k)计算得到迭代最小化稀疏贝叶斯重构算法第k次迭代过程中的重构残余误差,记为r(k),其中S为步骤3中得到的线阵SAR回波信号向量,A为步骤3中得到的线阵SAR测量矩阵,α(k)为步骤6中得到的迭代最小化稀疏贝叶斯重构算法第k次迭代过程中观测场景目标的散射系数向量;
采用公式计算得到第k次迭代过程中的相邻目标散射系数变化率,记为ω(k),其中α(k)为迭代最小化稀疏贝叶斯重构算法第k次迭代得到观测场景目标空间散射系数向量,α(k-1)为迭代最小化稀疏贝叶斯重构算法第k-1次迭代得到观测场景目标空间散射系数向量,||·||2为向量L2范数,k表示迭代最小化稀疏贝叶斯重构算法中的第k迭代次数;
步骤11、迭代终止判定:
如果若r(k)≥ε0且ω(k)≥ω0且k≤MaxIter,则k←k+1,执行步骤5至11,否则终止迭代,此刻第k次迭代得到的散射系数向量值α(k)即为观测场景目标空间Ω最终的散射系数向量,其中r(k)为步骤10中得到第k次迭代重构残余误差,ε0为步骤4中初始化得到的重构残余误差的门限,ω(k)为步骤10中得到第k次迭代过程中的相邻目标散射系数变化率,ω0为步骤4中初始化得到的相邻迭代目标散射系数变化率的门限,k表示迭代最小化稀疏贝叶斯重构算法中的第k迭代次数,MaxIter为步骤4中初始化得到的算法重构处理的最大迭代次数;
最后将观测场景目标空间散射系数向量α(k)转换成三维矩阵形式,得到线阵SAR观测场景目标空间Ω的三维成像结果。
CN201310752685.0A 2013-12-31 2013-12-31 基于迭代最小化稀疏贝叶斯重构线阵sar成像方法 Expired - Fee Related CN103713288B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201310752685.0A CN103713288B (zh) 2013-12-31 2013-12-31 基于迭代最小化稀疏贝叶斯重构线阵sar成像方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201310752685.0A CN103713288B (zh) 2013-12-31 2013-12-31 基于迭代最小化稀疏贝叶斯重构线阵sar成像方法

Publications (2)

Publication Number Publication Date
CN103713288A CN103713288A (zh) 2014-04-09
CN103713288B true CN103713288B (zh) 2015-10-28

Family

ID=50406413

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201310752685.0A Expired - Fee Related CN103713288B (zh) 2013-12-31 2013-12-31 基于迭代最小化稀疏贝叶斯重构线阵sar成像方法

Country Status (1)

Country Link
CN (1) CN103713288B (zh)

Families Citing this family (29)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105487052B (zh) * 2015-12-08 2017-10-17 电子科技大学 基于低相干性的压缩感知lasar稀布线阵优化方法
CN105388476B (zh) * 2015-12-28 2017-12-29 河南工业大学 一种基于联合稀疏模型的层析sar成像方法
CN105842693B (zh) * 2016-03-23 2018-03-30 哈尔滨工业大学 一种基于压缩感知的双通道sar动目标检测的方法
CN106405548A (zh) * 2016-08-23 2017-02-15 西安电子科技大学 基于多任务贝叶斯压缩感知的逆合成孔径雷达成像方法
CN106443675B (zh) * 2016-09-28 2018-12-25 北京航空航天大学 一种基于压缩感知的层析sar盲信源估计方法
CN106682621B (zh) * 2016-12-28 2019-01-29 西安电子科技大学 基于贝叶斯网络的sar图像识别方法
CN107015066B (zh) * 2017-03-27 2019-06-21 电子科技大学 一种基于稀疏贝叶斯学习的天线阵列故障诊断方法
CN107064930B (zh) * 2017-03-29 2020-02-18 西安电子科技大学 基于gpu的雷达前视成像方法
CN107132535B (zh) * 2017-04-07 2019-12-10 西安电子科技大学 基于变分贝叶斯学习算法的isar稀疏频带成像方法
CN107037429B (zh) * 2017-04-17 2020-06-16 电子科技大学 基于门限梯度追踪算法的线阵sar三维成像方法
CN108388828A (zh) * 2017-07-13 2018-08-10 中国科学院遥感与数字地球研究所 一种综合多源遥感数据的滨海湿地土地覆盖信息提取方法
CN108226927B (zh) * 2017-12-14 2020-03-27 电子科技大学 基于加权迭代最小稀疏贝叶斯重构算法的sar成像方法
CN108776339B (zh) * 2018-03-29 2021-08-17 清华大学 基于块稀疏迭代阈值处理的单比特合成孔径雷达成像方法
CN109100718B (zh) * 2018-07-10 2019-05-28 中国人民解放军国防科技大学 基于贝叶斯学习的稀疏孔径isar自聚焦与横向定标方法
CN109547029B (zh) * 2018-11-05 2023-04-14 东南大学 一种基于组稀疏结构的自适应匹配追踪信号重构方法
CN109447921A (zh) * 2018-12-05 2019-03-08 重庆邮电大学 一种基于重构误差的图像测量矩阵优化方法
CN110109101A (zh) * 2019-04-04 2019-08-09 电子科技大学 一种基于自适应阈值的压缩感知三维sar成像方法
CN109959932B (zh) * 2019-04-08 2023-06-30 西安电子科技大学 基于下降段曲线轨迹的雷达前视三维成像方法
CN110161499B (zh) * 2019-05-09 2022-07-01 东南大学 改进的稀疏贝叶斯学习isar成像散射系数估计方法
CN110082761A (zh) * 2019-05-31 2019-08-02 电子科技大学 分布式外辐射源雷达成像方法
CN110133656B (zh) * 2019-06-06 2022-05-03 电子科技大学 一种基于互质阵列的分解与融合的三维sar稀疏成像方法
CN111679277B (zh) * 2020-05-28 2022-05-03 电子科技大学 一种基于sbrim算法的多基线层析sar三维成像方法
CN111948652B (zh) * 2020-07-17 2023-05-05 北京理工大学 一种基于深度学习的sar智能参数化超分辨成像方法
CN113556131B (zh) * 2021-07-21 2022-05-03 中国人民解放军国防科技大学 一种复数域多任务贝叶斯压缩感知方法
CN113671492A (zh) * 2021-07-29 2021-11-19 西安电子科技大学 一种面向机动平台前视成像的samp重构方法
CN113484862B (zh) * 2021-08-04 2023-10-17 电子科技大学 一种自适应的高分宽幅sar清晰重构成像方法
CN113866766B (zh) * 2021-09-29 2024-03-22 电子科技大学 一种基于近场三维成像的雷达散射截面积精确外推方法
CN116540203B (zh) * 2023-07-04 2023-09-22 西安电子科技大学 基于快速稀疏贝叶斯的宽带雷达高速目标的相参积累方法
CN116990771B (zh) * 2023-08-04 2024-03-22 小儒技术(深圳)有限公司 一种自动利用雷达测量淤泥深度方法及系统

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102313888A (zh) * 2010-06-29 2012-01-11 电子科技大学 一种基于压缩传感的线阵sar三维成像方法
CN103439693A (zh) * 2013-08-16 2013-12-11 电子科技大学 一种线阵sar稀疏重构成像与相位误差校正方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8681037B2 (en) * 2011-04-28 2014-03-25 Raytheon Company Performance model for synthetic aperture radar automatic target recognition and method thereof

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102313888A (zh) * 2010-06-29 2012-01-11 电子科技大学 一种基于压缩传感的线阵sar三维成像方法
CN103439693A (zh) * 2013-08-16 2013-12-11 电子科技大学 一种线阵sar稀疏重构成像与相位误差校正方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
基于压缩传感的线阵SAR三维成像方法研究;韦顺军等;《宇航学报》;20111130;第32卷(第11期);2403-2409 *

Also Published As

Publication number Publication date
CN103713288A (zh) 2014-04-09

Similar Documents

Publication Publication Date Title
CN103713288B (zh) 基于迭代最小化稀疏贝叶斯重构线阵sar成像方法
CN103698763B (zh) 基于硬阈值正交匹配追踪的线阵sar稀疏成像方法
CN103439693B (zh) 一种线阵sar稀疏重构成像与相位误差校正方法
CN104833973B (zh) 基于半正定规划的线阵sar后向投影自聚焦成像方法
CN107037429B (zh) 基于门限梯度追踪算法的线阵sar三维成像方法
CN104950305B (zh) 一种基于稀疏约束的实波束扫描雷达角超分辨成像方法
CN105699969B (zh) 基于广义高斯约束的最大后验估计角超分辨成像方法
CN104391295A (zh) 一种图像熵最优的压缩传感sar稀疏自聚焦成像方法
CN106405548A (zh) 基于多任务贝叶斯压缩感知的逆合成孔径雷达成像方法
CN110244303B (zh) 基于sbl-admm的稀疏孔径isar成像方法
CN109212527B (zh) 用于高频地波雷达的大孔径分布式多站目标定位方法
CN107271993A (zh) 一种基于最大后验的扫描雷达角超分辨成像方法
CN106680817A (zh) 一种实现前视雷达高分辨成像的方法
CN108226927A (zh) 基于加权迭代最小稀疏贝叶斯重构算法的sar成像方法
CN105137425A (zh) 基于卷积反演原理的扫描雷达前视角超分辨方法
CN105137424A (zh) 一种杂波背景下实波束扫描雷达角超分辨方法
CN105676190B (zh) 一种校正合成孔径雷达回波数据的方法和装置
CN104950297A (zh) 基于矩阵1范数拟合的阵元误差估计方法
CN107607945A (zh) 一种基于空间嵌入映射的扫描雷达前视成像方法
CN102788979A (zh) 一种基于后向投影InSAR成像配准的GPU实现方法
CN105487052A (zh) 基于低相干性的压缩感知lasar稀布线阵优化方法
CN104251991A (zh) 一种基于稀疏度估计的分维度阈值迭代稀疏微波成像方法
CN113608218A (zh) 一种基于后向投影原理的频域干涉相位稀疏重构方法
Vachon et al. Validation of along-track interferometric SAR measurements of ocean surface waves
CN109696672A (zh) 一种基于空间结构关联性的高分辨率穿墙雷达成像方法

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C14 Grant of patent or utility model
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20151028

Termination date: 20191231