CN113484862A - 一种自适应的高分宽幅sar清晰重构成像方法 - Google Patents

一种自适应的高分宽幅sar清晰重构成像方法 Download PDF

Info

Publication number
CN113484862A
CN113484862A CN202110889943.4A CN202110889943A CN113484862A CN 113484862 A CN113484862 A CN 113484862A CN 202110889943 A CN202110889943 A CN 202110889943A CN 113484862 A CN113484862 A CN 113484862A
Authority
CN
China
Prior art keywords
time
channel
sampling
echo
reconstruction
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.)
Granted
Application number
CN202110889943.4A
Other languages
English (en)
Other versions
CN113484862B (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 CN202110889943.4A priority Critical patent/CN113484862B/zh
Publication of CN113484862A publication Critical patent/CN113484862A/zh
Application granted granted Critical
Publication of CN113484862B publication Critical patent/CN113484862B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO 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/00—Systems 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/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G01S13/9021—SAR image post-processing techniques
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO 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/00—Systems 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/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G01S13/9094—Theoretical aspects

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Physics & Mathematics (AREA)
  • Electromagnetism (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • General Physics & Mathematics (AREA)
  • Radar Systems Or Details Thereof (AREA)

Abstract

本发明公开了一种自适应的高分宽幅SAR清晰重构成像方法，它通过使用将SAR平台非线性运动完全补偿的IDR方法生成用于重构的子图像，然后基于重构图像质量的原则建立优化模型，估计清晰成像的重构系数。用图像熵作为图像质量评判标准，清晰的图像有低的图像熵；用共轭梯度法解决优化问题，即通过迭代使重构图像的图像熵最小化。最后通过估计的重构系数对子图像进行加权，以获得清晰的重构图像。本发明既保留了IDR算法可以精确补偿SAR平台的非线性运动的优势，又优化了IDR算法中的子图像加权求和中的权重，本发明可以应用于合成孔径雷达成像，尤其是HRWS SAR成像和地球遥感等领域。

Description

一种自适应的高分宽幅SAR清晰重构成像方法
技术领域
本发明属于雷达技术领域，它特别涉及了高分宽幅合成孔径雷达(SAR)成像技术领域。
背景技术
合成孔径雷达(SAR)作为一种具有全天时、全天候、信息量丰富的遥感成像技术，已成为当今对地观测的重要手段。高分辨率和宽幅成像(HRWS)已成为合成孔径雷达(SAR)技术的重要趋势。但对于传统的单通道SAR，需要较低的脉冲重复频率(PRF)才能获得距离向的较为宽幅的图像，这与方位向上为避免模糊而需要的较高RPF的相矛盾。多通道SAR技术是解决这一矛盾的有效方法。在多通道SAR中，第一个通道以低PRF发射信号，以确保较宽的距离向扫描范围，同时，所有通道接收回波信号以增加方位向分辨率。但是，多通道HRWSSAR系统通常会导致方位向周期性非均匀采样，这仍会导致方位向模糊。目前用于抑制方位向模糊的方法主要有：滤波器组重建算法，时空自适应处理方法，频谱重建算法，图像域重构方法(IDR)，详见文献N.Gebert,G.Krieger,and A.Moreira,“Digital beamforming onreceive:Techniques and optimization strategies for high-resolution wide-swathsar imaging,”Aerospace and Electronic Systems IEEE Transactions on,vol.45,no.2,pp.564–592,2009.。但在多通道SAR系统中存在机械制造误差和测量误差导致的基线误差，严重降低了上述重建方法的性能。为解决这个问题，一种方法是估计基线误差。在基于干涉相位的多通道GMTI-SAR系统中，基线误差估计方法得到了广泛研究。然而由于多普勒模糊，用于GMTI的方法可能是无效的HRWS SAR基线误差估计。S.Zhang提出了一种基于相关函数的HRWS SAR基线估计方法，但这可能会增加噪声分量并产生许多交叉项。Huang提出了一种采用迭代自适应方法的沿线基线估算方法，尽管此方法对幅度和相位误差具有鲁棒性，但由于空间相关运算，在杂波噪声比相对较低的情况下，该方法可能无效。
发明内容：
针对具有基线误差的HRWS SAR成像，本发明提出了一种自适应清晰重构(AUR)方法，这是一种无需基线估计即可根据图像质量自适应重构图像的新方法。为了进行这种自适应重构，本发明提出的方法基于有图像加权求和特点、适用于弯曲轨迹SAR重构的IDR重构方法。首先，使用可以将SAR平台非线性运动完全补偿的IDR方法生成用于重构的子图像，然后基于重构图像质量的原则建立优化模型，估计清晰成像的重构系数。用图像熵作为图像质量评判标准，清晰的图像有低的图像熵；用共轭梯度法解决优化问题，即通过迭代使重构图像的图像熵最小化。最后通过估计的重构系数对子图像进行加权，以获得清晰的重构图像。该方法的特点是：1)AUR是一种基于图像质量自适应重构图像的新方法，而不是在图像重构之前估算基线；2)通过使用IDR，AUR可以利用基线误差执行曲线轨迹SAR重构，从而补偿非线性运动以生成子图像。
为了方便描述本发明的内容，首先作以下术语定义：
定义1、脉冲重复频率(PRF)
PRF(pulse repetition frequency)即脉冲重复频率，简称重频，每秒钟发射的脉冲数目，是脉冲重复间隔(pulse repetition interval,PRI)的倒数。脉冲重复间隔就是一个脉冲和下一个脉冲之间的时间间隔。周期性的脉冲重复间隔时间为脉冲重复周期。详细内容可参考文献“雷达成像技术”，保铮等编著，电子工业出版社出版。
定义2、多通道高分宽幅(HRWS)SAR
为了获得方位的高分辨，传统的单通道合成孔径雷达系统在满足方位天线面积最小的同时还必须满足高脉冲重复频率(PRF)要求。但设置有限高的PRF是为了获得距离宽测绘带而需要克服距离模糊的必要条件，即单通道合成孔径雷达系统中，高分辨和宽测绘带是一对矛盾的性能指标。多通道合成孔径雷达系统的应用可以很好地解决这个问题，它结合数字波束形成技术(Digital Beam-Forming，DBF)可高效地实现高分辨宽测绘带对地观测，在地面动目标检测、高分宽幅成像及抑制干扰方面有显著优势，详见文献：郭晓江.高分宽幅SAR成像处理关键技术研究[D].上海交通大学,2018.
定义3、等效相位中心(EPC)原理
对于多通道SAR系统，可以看作是在位于发射机和接收机中间的等效相位中心(EPC)上发送和接收信号的SAR，详见文献G.Krieger,N.Gebert,and A.Moreira,“Unambiguous sar signal reconstruction from nonuniform displaced phase centersampling,”IEEE Geoscience and Remote Sensing Letters,vol.1,no.4,pp.260–264,2004.
定义4、标准的合成孔径雷达(SAR)距离向压缩方法
SAR距离向压缩方法是指利用合成孔径雷达发射参数进行匹配滤波的过程，主要包括：采用距离向参考信号对合成孔径雷达距离向回波信号进行匹配滤波压缩和距离走动补偿得到合成孔径雷达的成像。详细内容可参考文献“雷达成像技术”，保铮等编著，电子工业出版社出版。
定义5、标准图像域重构(IDR)方法
IDR方法是指基于SAR时域成像方法，利用重构系数从非均匀采样的回波信号重构均匀采样信号，已被广泛研究以抑制方位角模糊。详见文献N.Gebert,G.Krieger,andA.Moreira,“Digital beamforming on receive:Techniques and optimizationstrategies for high-resolution wide-swath SAR imaging,”Aerospace andElectronic Systems IEEE Transactions on,vol.45,no.2,pp.564–592,2009.
定义6、标准后向投影(BP)算法
BP算法是从层析成像技术中借鉴的处理方法，该方法通过将回波数据逐点投影到图像空间各像素，实现各散射点能量的积累。与其他方法相比，BP方法成像原理简单，具有运动补偿精度高，无几何失真，适用于大场景、大斜视、任意模式和任意轨迹的SAR成像等特点。基本思想是通过计算成像区域内每一像素点到孔径长度内SAR天线平台之间的距离历史，然后将对应的时域回波信号进行相干累加，从而重构出每个像素的散射系数
Figure BDA0003195541020000031
其中，σ(m)表示第m个散射点的后向散射系数，m＝1,2,...,M，M表示图像空间划分的总网格点数；s(r(n,m),n)表示第m个散射点在第n个慢时刻经插值后得到的回波值r(n,m)表示第m个散射点在第n个慢时刻到雷达的双程斜距；ω为发射信号的载波频率；N表示一个合成孔径内的慢时刻总数。
通常，使用均匀采样回波的成像场景中第m个像素的BP成像分量可以表示为
Figure BDA0003195541020000032
详见文献S.Jun,Z.Xiaoling,Y.Jianyu,and W.Chen,“Apc trajectory design for one-active linear-array three-dimensional imaging sar,”IEEE Transactions on Geoscience and Remote Sensing,vol.48,pp.1470–1486,March 2010.
定义7、图像熵
图像熵表示为图像灰度级集合的比特平均数，描述了图像信源的平均信息量。
对于BP成像后的SAR图像，其图像熵可以表示为：
Figure BDA0003195541020000033
其中，
Figure BDA0003195541020000034
I(m)表示的值，M表示图像空间划分的总网格点数。
当SAR图像质量越好，图像熵值越小。准确的重构系数估计会产生清晰的图像，从而最小化了图像的熵，详见文献：T.Zeng,R.Wang,and F.Li,“Sar image autofocusutilizing minimum-entropy criterion,”IEEE Geoscience and Remote SensingLetters,vol.10,no.6,pp.1552–1556,2013.。
定义8、梯度
对于多元函数，在点P上的梯度不是零向量，则它的方向是这个函数在P上最大增长的方向，而它的量是在这个方向上的增长率。
对于标量函数f梯度表示为
Figure BDA0003195541020000041
或gradf，其中
Figure BDA0003195541020000042
表示向量微分算子。
对于向量，相对于n×1向量x的梯度算子记作
Figure BDA0003195541020000043
定义为
Figure BDA0003195541020000044
对于矩阵，标量函数f(A)相对于m×n实矩阵A的梯度为一m×n矩阵，简称梯度矩阵，定义为
Figure BDA0003195541020000045
详细内容可参考文献“微积分”，傅英定，谢云荪主编，高等教育出版社出版。
定义9、标准共轭梯度方法
共轭梯度法是一种重要的求解无约束最优化问题的方法。其基本思想是把共轭性与最速下降方法相结合，利用已知点处的梯度构造一组共轭方向，并沿这组方向进行搜索，求出目标函数的极小点。Fletcher-Reeves共轭梯度法计算步骤如下：
(1)给定初始点x(1)，允许误差ε＞0，置
Figure BDA0003195541020000046
(2)若
Figure BDA0003195541020000047
则停止计算；否则，做一维搜索，求λj，满足
Figure BDA0003195541020000048
令y(j+1)＝y(j)+λjd(j).
(3)如果j＜n，则进行步骤(4)；否则，进行步骤(5).
(4)令
Figure BDA0003195541020000049
其中
Figure BDA00031955410200000410
置j:＝j+1，转步骤(2).
(5)令x(k+1)＝y(n+1)，y(1)＝x(k+1)，
Figure BDA0003195541020000051
置j＝1，k:＝k+1，转步骤(2).
详细内容可参考文献“最优化理论与算法”，陈宝林等编著，清华大学出版社出版。
定义10、Armijo线性搜索方法
一维线性搜索是指在得到点x(k)后，按某种规则确定一个方向d(k)，再从x(k)出发，沿方向d(k)在直线(或射线)上求目标函数的极小点，从而得到x(k)的后继点x(k+1)，重复以上做法，直到求得问题的解。
Armijo-Goldstein准则(Armijo准则)是不精确的一维搜索的一大准则，为了能使算法收敛(求最优解)。即要使我们的不精确的一维搜索的步长满足一定的规则，使之后的求最优解的过程不至于因为步长过大或者过小而不收敛。
Armijo准则的核心思想有两个：①目标函数值应该有足够的下降；②一维搜索的步长α不应该太小。
Armijo准则的数学表达式为：
Figure BDA0003195541020000052
Figure BDA0003195541020000053
详细内容可参考文献“最优化理论与算法”，陈宝林等编著，清华大学出版社出版。
定义11、复数取模运算
对于复数、复向量A，其模值
Figure BDA0003195541020000054
详细内容可参考文献“矩阵理论”，黄廷祝，钟守铭，李正良编著，高等教育出版社出版。
定义12、虚数单位j
虚数单位j定义为二次方程j2+1＝0的两个根的一个。
详细内容可参考文献“电磁场与电磁波”，谢处方主编，高等教育出版社出版。
定义13、偏导数运算
多元函数偏导数的定义是：一个多变量的函数(或称多元函数)，对其中一个变量(导数)微分，而保持其他变量恒定(相对于全导数，在其中所有变量都允许变化)。函数f关于变量x的偏导数写为fx'或写为
Figure BDA0003195541020000055
详细内容可参考文献“微积分”，傅英定，谢云荪主编，高等教育出版社出版。
定义14、矩阵范数
范数是具有“长度”概念的函数。在线性代数、泛函分析及相关的数学领域，是一个函数，其为向量空间内的所有向量赋予非零的正长度或大小。矩阵范数(matrix norm)亦译矩阵模是数学中矩阵论、线性代数、泛函分析等领域中常见的基本概念，是将一定的矩阵空间建立为赋范向量空间时为矩阵装备的范数。在n维复空间中，最常见范数为：
Figure BDA0003195541020000061
详细内容可参考文献“矩阵理论”，黄廷祝，钟守铭，李正良编著，高等教育出版社出版。
本发明提供的一种自适应的高分宽幅SAR清晰重构成像方法，它包括以下步骤：
步骤1、初始化自适应的高分宽幅SAR清晰重构成像所需要的多通道SAR系统参数：
初始化自适应的高分宽幅SAR清晰重构成像所需要的多通道SAR系统参数，包括：多通道SAR的观测空间为地面三维坐标系，记为X-Y-Z，其中X表示水平面横轴，Y表示水平面纵轴，Z表示水平垂直轴；多通道SAR的载波频率fcw，信号带宽Bs，距离向采样率fsr，方位向带宽Br；多通道SAR成像系统沿着Y轴平行方向进行匀速直线运动，平台速度矢量为V；平台飞行高度h；多通道SAR的通道数目为N，记第一个发射信号的通道为Tx，它同时也是接收回波的第一个通道Rx1；其余的通道与第一个通道同时接收回波，依次记为Rx2,Rx3,...,RxN；雷达中心频率fc，雷达系统发射信号的脉冲重复频率PRF，发射信号脉冲重复时间
Figure BDA0003195541020000062
设每个通道之间的距离为d，则第n个通道与第一个通道(参考通道)之间距离的基线长度可以表示为xn＝(n-1)·d。由定义3中等效相位中心(EPC)原理可知，通道的等效相位中心之间的间隔为d/2。当PRF，v和d满足
Figure BDA0003195541020000063
时，EPC是均匀分布的。当EPC非均匀分布时，方位向上的第n个通道在第k次脉冲重复周期的采样时间点记为tk,n＝kTr+tn,k＝1,2,...,K，其中，Tr为脉冲重复时间，K为方位向回波的采样周期总数，tn＝(n-1)·Tb是第n个通道的采样时间偏差，即方位向空间上同时接收信号的N个通道可以分别等效为第一个发射/接收通道(Tx/Rx1)在(n-1)·Tb时间后运动的位置，Tb为两个相邻通道之间的采样时间间隔；
由于存在通道间隔引起的基线误差，设考虑误差后的通道间隔测量值为
Figure BDA0003195541020000064
记满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间为
Figure BDA0003195541020000071
其中，k′代表脉冲重复周期，取值k'＝1,2,...,K，n′代表通道，取值1,2,...,N；
记距离向压缩后回波为sr(τ(t,m),t)，其中t为方位向慢时间，τ为在方位向慢时间t时的从EPC到第m个像素的回波的双程延迟的距离向快时间，m是图像域的第m个像素。则回波的非均匀时间采样可以表示为sr(τ(tk,n,m),tk,n)，其中，tk,n为方位向快时间，τ(tk,n,m)为在方位向慢时间tk,n时的从EPC到第m个像素的回波的双程延迟的距离向快时间；回波的均匀时间采样可以表示为sr(τ(tk′,n′,m),tk′,n′)，其中，tk′,n′为方位向快时间，τ(tk′,n′,m)为在方位向慢时间tk′,n′时的从EPC到第m个像素的回波的双程延迟的距离向快时间；
步骤2、初始化自适应的高分宽幅SAR清晰重构成像的观测场景目标空间参数：
初始化自适应的高分宽幅SAR清晰重构成像的场景空间参数，包括：以地平面作为IDR成像方法BP算法的投影空间，将图像空间划分为大小相等的单元网格点，根据处理区域选择对应数据段，数据段起始位置为天线波束进入处理区域时刻到天线波束滑出处理区域时刻，其总长度为处理区域长度加上合成孔径长度。初始化距离向像素数Nr，方位向像素数量Na，距离向网格分辨率dr，方位向网格分辨率da，总网格数M个；
步骤3、生成原始回波数据，并进行距离向脉冲压缩，得到脉冲压缩后的回波数据：
生成多通道SAR的原始回波数据，采用定义4中标准合成孔径雷达回波数据距离向脉冲压缩方法对原始回波数据进行距离向脉冲压缩，得到距离向压缩后的多通道SAR回波数据，记作sr(τ(t,m),t)，即步骤1中的距离向压缩后的回波，其中t为方位向慢时间，τ为在方位向慢时间t时的从EPC到第m个像素的回波的双程延迟的距离向快时间，m是图像域的第m个像素；
步骤4、计算由非均匀采样回波重构的均匀采样回波：
采用公式
Figure BDA0003195541020000072
计算重构系数Ψk,n(tk′,n′)，其中，N为通道数目，tk,n＝kTr+tn,k＝1,2,...,K是方位向上的第n个通道在第k次脉冲重复周期的采样时间点，Tr为脉冲重复时间，K为方位向回波的采样周期总数，tn＝(n-1)·Tb是第n个通道的采样时间偏差，即方位向空间上同时接收信号的N个通道可以分别等效为第一个发射/接收通道(Tx/Rx1)在(n-1)·Tb时间后运动的位置，Tb为两个相邻通道之间的采样时间间隔；
Figure BDA0003195541020000081
为满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间，k′代表脉冲重复周期，取值k'＝1,2,...,K，n′代表通道，取值1,2,...,N；tn＝(n-1)·Tb，n＝1,2,...,N是第n个通道的采样时间偏差；tq为第q个通道中的采样时间偏差，tq＝(q-1)·Tb，q＝1,2,...,N；π为圆周率，PRF为定义1中脉冲重复频率，Π为累乘符号；
再采用公式
Figure BDA0003195541020000082
计算得到由非均匀采样的信号重构获得的均匀采样的回波信号sr(τ(tk',n',m),tk',n')，其中N为通道总数目，sr(τ(tk,n,m),tk,n)是步骤1中的回波的非均匀时间采样，sr(τ(tk′,n′,m),tk′,n′)为步骤1中的回波的均匀时间采样；tk,n＝kTr+tn,k＝1,2,...,K是方位向上的第n个通道在第k次脉冲重复周期的采样时间点，Tr为脉冲重复时间，K为方位向回波的采样周期总数，tn＝(n-1)·Tb是第n个通道的采样时间偏差，即方位向空间上同时接收信号的N个通道可以分别等效为第一个发射/接收通道(Tx/Rx1)在(n-1)·Tb时间后运动的位置，Tb为两个相邻通道之间的采样时间间隔；
Figure BDA0003195541020000083
为满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间，k′代表脉冲重复周期，取值k'＝1,2,...,K，n′代表通道，取值1,2,...,N；L表示插值周期，
Figure BDA0003195541020000085
sr(τ(tk,n,m),tk,n)为步骤1中定义的回波的非均匀时间采样，其中，tk,n为方位向快时间，τ(tk,n,m)为在方位向慢时间tk,n时的从EPC到第m个像素的回波的双程延迟的距离向快时间；sr(τ(tk′,n′,m),tk′,n′)为步骤1中定义的回波的均匀时间采样，其中，tk′,n′为方位向快时间，τ(tk′,n′,m)为在方位向慢时间tk′,n′时的从EPC到第m个像素的回波的双程延迟的距离向快时间；Ψk,n(tk′,n′)为步骤4声明的重构系数；
步骤5、计算第m个散射点在通道n中的子图像：
记k'与k的差值为l，k'＝k-l，则步骤4中计算的重构系数Ψk,n(tk′,n′)可以记作与tk',n'无关的l·n维度向量的重构系数Ψln。
采用公式
Figure BDA0003195541020000084
计算通道n中的子图像，其中，ω是发射信号中心角频率，m是图像域的第m个像素，sr(τ(tk,n,m),tk,n)是步骤1中的回波的非均匀时间采样，τ(tk,n,m)为在方位向慢时间tk,n时的从EPC到第m个像素的回波的双程延迟的距离向快时间，tk,n＝kTr+tn,k＝1,2,...,K是方位向上的第n个通道在第k次脉冲重复周期的采样时间点，Tr为脉冲重复时间，K为方位向回波的采样周期总数；τ(tk′,n′,m)为在方位向慢时间tk′,n′时，从EPC到第m个像素的回波的双程延迟的距离向快时间，
Figure BDA0003195541020000091
为满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间，k′代表脉冲重复周期，取值k'＝1,2,...,K，n′代表通道，取值1,2,...,N；j为定义12中的虚数单位；
步骤6、计算成像场景中第m个散射点在重构出的均匀采样通道n′中的BP成像分量：
将步骤4中计算出的由非均匀采样的信号重构获得的均匀采样的信号sr(τ(tk',n',m),tk',n')，采用公式
Figure BDA0003195541020000092
计算成像场景中第m个散射点在重构出的均匀采样通道n′中的BP成像分量，记为In′(m)，其中Ψln为步骤5中的l·n维度向量的重构系数，Iln(m)为步骤5中计算得到的通道n中的子图像，l＝k-k'，m取值为1,2,...,M，M为步骤2中定义的总网格数；
步骤7、计算重构出的均匀采样回波的BP成像：
再采用公式
Figure BDA0003195541020000093
计算第m个散射点的成像结果，其中，In'(m)是步骤6中计算得到的成像场景中第m个散射点的在重构出的均匀采样通道n′中的BP成像分量。
步骤8、计算图像熵目标函数：
步骤7中计算的
Figure BDA0003195541020000094
由步骤6中的公式
Figure BDA0003195541020000095
则与重构系数Ψ相关的第m个散射点的成像结果表达式为
Figure BDA0003195541020000096
采用j＝1,2,...,J代替n′ln，其中J＝N2L，k表示相干积累时回波的那个序号，考虑了不为负数，且要保证重构时取满原始回波ln个方位数；则In′ln(m)可写为Ij(m)，Ψn′ln记为一维向量Ψ＝[ψ1,ψ2,...,ψj,...,ψJ]，
成像场景中第m个散射点的BP成像分量记为
Figure BDA0003195541020000097
记
Figure BDA0003195541020000098
为重构系数的变量估计，则图像熵目标函数可记为
Figure BDA0003195541020000099
其中，
Figure BDA00031955410200000910
m为成像场景中第m个散射点，
Figure BDA0003195541020000101
为与重构系数的变量估计
Figure BDA0003195541020000102
相关的第m个散射点的成像结果，M为步骤2中所声明的总网格数，|·|为定义11中的复数取模运算；
步骤9、计算图像熵目标函数的梯度下降方向：
采用公式
Figure BDA0003195541020000103
计算目标函数图像熵的梯度下降方向
Figure BDA0003195541020000104
其中，
Figure BDA0003195541020000105
是定义8中的梯度运算，
Figure BDA0003195541020000106
为定义13中的偏导数运算，
Figure BDA0003195541020000107
为步骤8中定义的重构系数的变量估计，j＝1,2,...,J，J＝N2L为步骤8中定义；
步骤10、设定基于共轭梯度法的最优求解算法的初始参数：
初始化基于共轭梯度法的最优求解算法参数包括：最大迭代次数，记为Nitermax；重构出的子图像，记第m个散射点的BP成像后子图像为Ij(m)，j＝1,2,...,J，将步骤8中定义的重构系数
Figure BDA0003195541020000108
的初始迭代值记为
Figure BDA0003195541020000109
搜索方向的迭代初始值记为d(0)，d(0)的值选择为目标函数图像熵的负梯度方向，
Figure BDA00031955410200001010
其中
Figure BDA00031955410200001011
为定义8中的梯度运算；n表示Armijo算法中的第n次迭代次数，设n的初始值n＝0；第n次迭代的搜索方向记为d(n)；
在满足n＜Nitermax条件下，重复执行步骤11-步骤13，第n次迭代过程中将会更新
Figure BDA00031955410200001012
值为
Figure BDA00031955410200001013
且n＝n+1；直到n＝Nitermax时，输出
Figure BDA00031955410200001014
作为最优重构系数估计结果，记为
Figure BDA00031955410200001015
步骤11、采用Armijo线性搜索算法计算搜索步长：
求解公式
Figure BDA00031955410200001016
计算出当前第n次迭代的搜索步长λ(n)，其中，
Figure BDA00031955410200001017
为所求目标函数，min为求解其最小值，d(n)为步骤10中定义的第n次迭代的搜索方向，
Figure BDA00031955410200001018
为步骤10中定义的第n次迭代过程中更新的
Figure BDA00031955410200001019
值，ε(·)为定义7中的图像熵；
步骤12、计算下次迭代的重构系数
Figure BDA00031955410200001020
当n＜Nitermax-1时，采用公式
Figure BDA00031955410200001021
计算下次迭代用的重构系数
Figure BDA00031955410200001022
其中，λ(n)为步骤11中计算出的第n次迭代的搜索步长，
Figure BDA00031955410200001023
为步骤10中定义的第n次迭代过程中更新的
Figure BDA00031955410200001024
值，n为步骤10中定义的第n次迭代，Nitermax为步骤10中定义的最大迭代次数；
当n＝Nitermax-1时，记最优重构系数估计结果
Figure BDA0003195541020000111
返回步骤10，输出
Figure BDA0003195541020000112
作为最优重构系数估计结果，其中，n为步骤10中定义的第n次迭代，Nitermax为步骤10中定义的最大迭代次数，
Figure BDA0003195541020000113
为步骤10中定义的第n次迭代过程更新后的
Figure BDA0003195541020000114
值；
步骤13、更新搜索方向：
采用公式
Figure BDA0003195541020000115
计算更新后的搜索方向d(n+1)，并更新n＝n+1，其中
Figure BDA0003195541020000116
为步骤10中定义的第n次迭代过程更新后的
Figure BDA0003195541020000117
值，
Figure BDA0003195541020000118
为定义8中的梯度方向，||·||为定义14中的矩阵范数；
步骤14、计算重构后的清晰图像：
采用步骤8中推导的公式
Figure BDA0003195541020000119
其中，j为步骤8中定义的代替n′ln的变量，j＝1,2,...,J，J＝N2L，将步骤10中计算得到的最优重构系数估计结果
Figure BDA00031955410200001110
中的ψ1,ψ2,...,ψj,...,ψJ，分别作为
Figure BDA00031955410200001111
中的ψj带入公式
Figure BDA00031955410200001112
计算第m个散射点重构后获得的清晰图像IΨ(m)，其中，Ij(m)为步骤8中定义的用j代替n′ln后的第m个散射点在满足均匀时间采样的第n′通道的图像In′ln(m)，n′代表通道，取值1,2,...,N；
经过以上步骤后就得到本发明的第m个散射点重构后获得的清晰图像IΨ(m)。
本发明的实质是在IDR方法的基础上提出了一种自适应的高分宽幅SAR清晰重构成像(AUR)方法。该方法既保留了IDR算法可以精确补偿SAR平台的非线性运动的优势，又优化了IDR算法中的子图像加权求和中的权重。AUR方法用图像熵作为衡量标准，用共轭梯度法估计清晰成像的重构系数的最优化，代替传统IDR算法中的直接子图像加权求和，从而获得清晰的重构图像。根据本人了解，当前还没有出现自适应的高分宽幅SAR清晰重构成像方法像方法；
本发明的创新点在于针对HRWS SAR成像中的基线误差，基于IDR算法中子图像加权求和的特征，建立了图像熵模型衡量图像的清晰度，并提供了基于Armijo算法的子图像加权系数优化方法，通过迭代使重构图像的图像熵最小化的，能自适应地调整IDR方法重构成像结果；
本发明的优点在于利用图像熵模型衡量图像的清晰度进行重构系数的优化，是一种基于图像质量自适应重构图像的方法，不需要在IDR算法重构图像之前计算基线；保留了IDR方法可以补偿非直线运动的优点，同样适用于曲线轨迹SAR重构；基于IDR方法，可以抑制多通道SAR系统在方位上形成周期性的非均匀采样导致的图像模糊。本发明可以应用于合成孔径雷达成像，尤其是HRWS SAR成像和地球遥感等领域。
附图说明
图1为本发明所提供的一种自适应的高分宽幅SAR清晰重构成像方法的处理流程示意框图。
图2为本发明具体实施方式采用的多通道SAR系统仿真参数表。
具体实施方式
本发明主要采用仿真实验地方法进行验证，所有步骤、结论都在MATLAB上验证正确。具体实施步骤如下：
步骤1、初始化自适应的高分宽幅SAR清晰重构成像所需要的多通道SAR系统参数：
初始化自适应的高分宽幅SAR清晰重构成像所需要的多通道SAR系统参数，包括：多通道SAR的观测空间为地面三维坐标系，记为X-Y-Z，其中X表示水平面横轴，Y表示水平面纵轴，Z表示水平垂直轴；多通道SAR的载波频率fcw＝9.6GHz，信号带宽Bs＝150MHz，距离向采样率fsr＝210MHz，方位向带宽Br＝1900Hz；多通道SAR成像系统沿着Y轴平行方向进行匀速直线运动，平台速度矢量为
Figure BDA0003195541020000121
平台飞行高度h＝20Km；多通道SAR的通道数目为N＝4，记第一个发射信号的通道为Tx，它同时也是接收回波的第一个通道Rx1；其余的通道与第一个通道同时接收回波，依次记为Rx2,Rx3,...,RxN；雷达中心频率fc＝10GHz，雷达系统发射信号的脉冲重复频率PRF＝700Hz，发射信号脉冲重复时间
Figure BDA0003195541020000122
设每个通道之间的距离为d＝1m，由于存在通道间隔引起的基线误差，设考虑误差后的通道间隔测量值为
Figure BDA0003195541020000123
则第n个通道与第一个通道(参考通道)之间距离的基线长度可以表示为
Figure BDA0003195541020000124
由等效相位中心(EPC)原理可知，通道的等效相位中心之间的间隔为
Figure BDA0003195541020000131
当PRF，v和d满足
Figure BDA0003195541020000132
时，EPC是均匀分布的；当EPC非均匀分布时，方位向上的第n个通道在第k次脉冲重复周期的采样时间点记为tk,n＝kTr+tn,k＝1,2,...,K，其中，
Figure BDA0003195541020000133
为脉冲重复时间，K＝1024为方位向回波的采样周期总数，tn＝(n-1)·Tb是第n个通道的采样时间偏差，即方位向空间上同时接收信号的N个通道可以分别等效为第一个发射/接收通道(Tx/Rx1)在(n-1)·Tb时间后运动的位置，Tb＝0.263ms为两个相邻通道之间的采样时间间隔。
由于存在通道间隔引起的基线误差，设考虑误差后的通道间隔测量值为
Figure BDA0003195541020000134
记满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间为
Figure BDA0003195541020000135
其中，k′代表脉冲重复周期，取值k'＝1,2,...,K，K＝1024为方位向回波的采样周期总数，
Figure BDA0003195541020000136
为脉冲重复时间，n′代表通道，取值n′＝1,2,...,N，其中N＝4为上文声明的多通道SAR的通道数目；
记距离向压缩后回波为sr(τ(t,m),t)，其中t为方位向慢时间，τ为在方位向慢时间t时的从EPC到第m个像素的回波的双程延迟的距离向快时间，m是图像域的第m个像素。则回波的非均匀时间采样可以表示为sr(τ(tk,n,m),tk,n)，其中，tk,n为方位向快时间，τ(tk,n,m)为在方位向慢时间tk,n时的从EPC到第m个像素的回波的双程延迟的距离向快时间；回波的均匀时间采样可以表示为sr(τ(tk′,n′,m),tk′,n′)，其中，tk′,n′为方位向快时间，τ(tk′,n′,m)为在方位向慢时间tk′,n′时的从EPC到第m个像素的回波的双程延迟的距离向快时间；
步骤2、初始化自适应的高分宽幅SAR清晰重构成像的观测场景目标空间参数：
初始化自适应的高分宽幅SAR清晰重构成像的场景空间参数，包括：以地平面作为IDR成像方法BP算法的投影空间，将图像空间划分为大小相等的单元网格点，根据处理区域选择对应数据段，数据段起始位置为天线波束进入处理区域时刻到天线波束滑出处理区域时刻，其总长度为处理区域长度加上合成孔径长度。初始化距离向像素数Nr＝1000，方位向像素数量Na＝5000，距离向网格分辨率dr＝1m，方位向网格分辨率da＝1m，总网格数M＝5000000个；
步骤3、生成原始回波数据，并进行距离向脉冲压缩，得到脉冲压缩后的回波数据：
生成多通道SAR的原始回波数据，采用定义4中标准合成孔径雷达回波数据距离向脉冲压缩方法对原始回波数据进行距离向脉冲压缩，得到距离向压缩后的多通道SAR回波数据，记作sr(τ(t,m),t)，即步骤1中声明的距离向压缩后的回波，其中t为方位向慢时间，τ为在方位向慢时间t时的从EPC到第m个像素的回波的双程延迟的距离向快时间，m是图像域的第m个像素；
步骤4、计算由非均匀采样回波重构的均匀采样回波：
采用公式
Figure BDA0003195541020000141
计算重构系数Ψk,n(tk′,n′)，其中，N＝4为步骤1中声明的通道数目，tk,n＝kTr+tn,k＝1,2,...,K是方位向上的第n个通道在第k次脉冲重复周期的采样时间点，Tr为脉冲重复时间，K＝1024为方位向回波的采样周期总数，tn＝(n-1)·Tb是第n个通道的采样时间偏差，即方位向空间上同时接收信号的N个通道可以分别等效为第一个发射/接收通道(Tx/Rx1)在(n-1)·Tb时间后运动的位置，Tb＝0.263ms为两个相邻通道之间的采样时间间隔；
Figure BDA0003195541020000142
为满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间，k′代表脉冲重复周期，取值k'＝1,2,...,K，n′代表通道，取值1,2,...,N；tn＝(n-1)·Tb，n＝1,2,...,N是第n个通道的采样时间偏差；tq为第q个通道中的采样时间偏差，tq＝(q-1)·Tb，q＝1,2,...,N；π为圆周率，PRF＝700Hz为定义1中脉冲重复频率，Π为累乘符号；
再采用公式
Figure BDA0003195541020000143
计算得到由非均匀采样的信号重构获得的均匀采样的回波信号sr(τ(tk',n',m),tk',n')，其中N＝4为通道总数目，sr(τ(tk,n,m),tk,n)是步骤1中声明的回波的非均匀时间采样，sr(τ(tk′,n′,m),tk′,n′)为步骤1中声明的回波的均匀时间采样；tk,n＝kTr+tn,k＝1,2,...,K是方位向上的第n个通道在第k次脉冲重复周期的采样时间点，Tr为脉冲重复时间，K＝1024为方位向回波的采样周期总数，tn＝(n-1)·Tb是第n个通道的采样时间偏差，即方位向空间上同时接收信号的N个通道可以分别等效为第一个发射/接收通道(Tx/Rx1)在(n-1)·Tb时间后运动的位置，Tb＝0.263ms为两个相邻通道之间的采样时间间隔；
Figure BDA0003195541020000151
为满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间，k′代表脉冲重复周期，取值k'＝1,2,...,K，n′代表通道，取值1,2,...,N；L＝6表示插值周期，
Figure BDA0003195541020000152
sr(τ(tk,n,m),tk,n)为步骤1中定义的回波的非均匀时间采样，其中，tk,n为方位向快时间，τ(tk,n,m)为在方位向慢时间tk,n时的从EPC到第m个像素的回波的双程延迟的距离向快时间；sr(τ(tk′,n′,m),tk′,n′)为步骤1中定义的回波的均匀时间采样，其中，tk′,n′为方位向快时间，τ(tk′,n′,m)为在方位向慢时间tk′,n′时的从EPC到第m个像素的回波的双程延迟的距离向快时间；Ψk,n(tk′,n′)为步骤4声明的重构系数；
步骤5、计算第m个散射点在通道n中的子图像：
记k'与k的差值为l，k'＝k-l，则步骤4中计算的重构系数Ψk,n(tk′,n′)可以记作与tk',n'无关的l·n维度向量的重构系数Ψln。采用公式
Figure BDA0003195541020000153
计算通道n中的子图像，其中，ω是发射信号中心角频率，ω＝2πf＝19.2πG(rad/s)，m是图像域的第m个像素，sr(τ(tk,n,m),tk,n)是步骤1中声明的回波的非均匀时间采样，τ(tk,n,m)为在方位向慢时间tk,n时的从EPC到第m个像素的回波的双程延迟的距离向快时间，tk,n＝kTr+tn,k＝1,2,...,K是方位向上的第n个通道在第k次脉冲重复周期的采样时间点，Tr为脉冲重复时间，K为方位向回波的采样周期总数；τ(tk′,n′,m)为在方位向慢时间tk′,n′时，从EPC到第m个像素的回波的双程延迟的距离向快时间，
Figure BDA0003195541020000154
为满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间，k′代表脉冲重复周期，取值k'＝1,2,...,K，n′代表通道，取值1,2,...,N；j为定义12中的虚数单位；
步骤6、计算成像场景中第m个散射点在重构出的均匀采样通道n′中的BP成像分量：
将步骤4中计算出的由非均匀采样的信号重构获得的均匀采样的信号sr(τ(tk',n',m),tk',n')，采用公式
Figure BDA0003195541020000155
计算成像场景中第m个散射点在重构出的均匀采样通道n′中的BP成像分量，记为In′(m)，其中Ψln为步骤5中声明的l·n维度向量的重构系数，Iln(m)为步骤5中计算得到的通道n中的子图像，l＝k-k'，m取值为1,2,...,M，M＝5000000为步骤2中所声明的总网格数；
步骤7、计算重构出的均匀采样回波的BP成像：
再采用公式
Figure BDA0003195541020000161
计算第m个散射点的成像结果，其中，In'(m)是步骤6中计算得到的成像场景中第m个散射点的在重构出的均匀采样通道n′中的BP成像分量。
步骤8、计算图像熵目标函数：
步骤7中计算的
Figure BDA0003195541020000162
由步骤6中的公式
Figure BDA0003195541020000163
则与重构系数Ψ相关的第m个散射点的成像结果表达式可以写为
Figure BDA0003195541020000164
采用j＝1,2,...,J代替n′ln，其中J＝N2L＝112，k表示相干积累时回波的那个序号，考虑了不为负数，且要保证重构时取满原始回波ln个方位数；则In′ln(m)可写为Ij(m)，Ψn′ln可写为一维向量Ψ＝[ψ1,ψ2,...,ψj,...,ψJ]，成像场景中第m个散射点的BP成像分量可写为
Figure BDA0003195541020000165
记
Figure BDA0003195541020000166
为重构系数的变量估计，则图像熵目标函数可记为
Figure BDA0003195541020000167
其中，
Figure BDA0003195541020000168
m为成像场景中第m个散射点，
Figure BDA0003195541020000169
为与重构系数的变量估计
Figure BDA00031955410200001610
相关的第m个散射点的成像结果，M＝5000000为步骤2中所声明的总网格数，|·|为定义11中的复数取模运算；
步骤9、计算图像熵目标函数的梯度下降方向：
采用公式
Figure BDA00031955410200001611
计算目标函数图像熵的梯度下降方向
Figure BDA00031955410200001612
其中，
Figure BDA00031955410200001613
是定义8中的梯度运算，
Figure BDA00031955410200001614
为定义13中的偏导数运算，
Figure BDA00031955410200001615
为步骤8中定义的重构系数的变量估计，j＝1,2,...,J，J＝N2L＝112为步骤8中定义；
步骤10、设定基于共轭梯度法的最优求解算法的初始参数：
初始化基于共轭梯度法的最优求解算法参数包括：最大迭代次数，记为Nitermax＝10；重构出的子图像，记第m个散射点的BP成像后子图像为Ij(m)，j＝1,2,...,J，J＝112为步骤8中定义，将步骤8中定义的重构系数
Figure BDA00031955410200001616
的初始迭代值记为
Figure BDA00031955410200001617
根据
Figure BDA00031955410200001618
计算初始迭代值
Figure BDA00031955410200001619
搜索方向的迭代初始值记为d(0)，d(0)的值选择为目标函数图像熵的负梯度方向，
Figure BDA0003195541020000171
其中
Figure BDA0003195541020000172
为定义8中的梯度运算；n表示Armijo算法中的第n次迭代次数，设n的初始值n＝0；第n次迭代的搜索方向记为d(n)；在满足n＜Nitermax条件下，重复执行步骤11-步骤13，第n次迭代过程中将会更新
Figure BDA0003195541020000173
值为
Figure BDA0003195541020000174
且n＝n+1；直到n＝Nitermax时，输出
Figure BDA0003195541020000175
作为最优重构系数估计结果，记为
Figure BDA0003195541020000176
步骤11、采用Armijo线性搜索算法计算搜索步长：
求解公式
Figure BDA0003195541020000177
计算出当前第n次迭代的搜索步长λ(n)，其中，
Figure BDA0003195541020000178
为所求目标函数，min为求解其最小值，d(n)为步骤10中定义的第n次迭代的搜索方向，
Figure BDA0003195541020000179
为步骤10中定义的第n次迭代过程中更新的
Figure BDA00031955410200001710
值，ε(·)为定义7中的图像熵；
步骤12、计算下次迭代的重构系数
Figure BDA00031955410200001711
当n＜Nitermax-1时，采用公式
Figure BDA00031955410200001712
计算下次迭代用的重构系数
Figure BDA00031955410200001713
其中，λ(n)为步骤11中计算出的第n次迭代的搜索步长，
Figure BDA00031955410200001714
为步骤10中定义的第n次迭代过程中更新的
Figure BDA00031955410200001715
值，n为步骤10中定义的第n次迭代，Nitermax＝10为步骤10中定义的最大迭代次数；当n＝Nitermax-1时，记最优重构系数估计结果
Figure BDA00031955410200001716
返回步骤10，输出
Figure BDA00031955410200001717
作为最优重构系数估计结果，其中，n为步骤10中定义的第n次迭代，Nitermax为步骤10中定义的最大迭代次数，
Figure BDA00031955410200001718
为步骤10中定义的第n次迭代过程更新后的
Figure BDA00031955410200001719
值；
步骤13、更新搜索方向：
采用公式
Figure BDA00031955410200001720
计算更新后的搜索方向d(n+1)，并更新n＝n+1，其中
Figure BDA00031955410200001721
为步骤10中定义的第n次迭代过程更新后的
Figure BDA00031955410200001722
值，
Figure BDA00031955410200001723
为定义8中的梯度方向，||·||为定义14中的矩阵范数；
步骤14、计算重构后的清晰图像：
采用步骤8中推导的公式
Figure BDA00031955410200001724
其中，j为步骤8中定义的代替n′ln的变量，j＝1,2,...,J，J＝N2L＝112为步骤8中定义，将步骤10中计算得到的最优重构系数估计结果
Figure BDA0003195541020000181
中的ψ1,ψ2,...,ψj,...,ψJ，分别作为
Figure BDA0003195541020000182
中的ψj带入公式
Figure BDA0003195541020000183
计算第m个散射点重构后获得的清晰图像IΨ(m)，其中，Ij(m)为步骤8中定义的用j代替n′ln后的第m个散射点在满足均匀时间采样的第n′通道的图像In′ln(m)，n′代表通道，取值1,2,...,N，N＝4为步骤1中定义。
至此，我们得到自适应的高分宽幅SAR清晰重构成像，整个方法结束。
通过本发明具体实施方式可以看出，本发明通过的IDR方法成像原理简单、运动补偿精度高、无几何失真、适用于大场景、大斜视、任意模式、任意轨迹的优势，利用AUR方法用图像熵作为衡量标准，采用共轭梯度法估计清晰成像的重构系数的最优化，代替传统IDR算法中的直接子图像加权求和，从而获得清晰的重构图像，相对于传统更便于系统实现，不需要在IDR算法重构图像之前计算基线；保留了IDR方法可以补偿非直线运动的优点，同样适用于曲线轨迹SAR重构；并且相对于传统IDR方法，本方法更能有效，可以抑制多通道SAR系统在方位上形成周期性的非均匀采样导致的图像模糊，提高成像质量。

Claims (1)

1.一种自适应的高分宽幅SAR清晰重构成像方法，其特征是它包括以下步骤：
步骤1、初始化自适应的高分宽幅SAR清晰重构成像所需要的多通道SAR系统参数：
初始化自适应的高分宽幅SAR清晰重构成像所需要的多通道SAR系统参数，包括：多通道SAR的观测空间为地面三维坐标系，记为X-Y-Z，其中X表示水平面横轴，Y表示水平面纵轴，Z表示水平垂直轴；多通道SAR的载波频率fcw，信号带宽Bs，距离向采样率fsr，方位向带宽Br；多通道SAR成像系统沿着Y轴平行方向进行匀速直线运动，平台速度矢量为
Figure FDA0003195541010000011
平台飞行高度h；多通道SAR的通道数目为N，记第一个发射信号的通道为Tx，它同时也是接收回波的第一个通道Rx1；其余的通道与第一个通道同时接收回波，依次记为Rx2,Rx3,...,RxN；雷达中心频率fc，雷达系统发射信号的脉冲重复频率PRF，发射信号脉冲重复时间
Figure FDA0003195541010000012
设每个通道之间的距离为d，则第n个通道与第一个通道(参考通道)之间距离的基线长度可以表示为xn＝(n-1)·d；通道的等效相位中心之间的间隔为d/2；当PRF，v和d满足
Figure FDA0003195541010000013
EPC是均匀分布的；当EPC非均匀分布时，方位向上的第n个通道在第k次脉冲重复周期的采样时间点记为tk,n＝kTr+tn,k＝1,2,...,K，其中，Tr为脉冲重复时间，K为方位向回波的采样周期总数，tn＝(n-1)·Tb是第n个通道的采样时间偏差，即方位向空间上同时接收信号的N个通道可以分别等效为第一个发射/接收通道(Tx/Rx1)在(n-1)·Tb时间后运动的位置，Tb为两个相邻通道之间的采样时间间隔；
由于存在通道间隔引起的基线误差，设考虑误差后的通道间隔测量值为
Figure FDA0003195541010000015
记满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间为
Figure FDA0003195541010000014
其中，k′代表脉冲重复周期，取值k'＝1,2,...,K，n′代表通道，取值1,2,...,N；
记距离向压缩后回波为sr(τ(t,m),t)，其中t为方位向慢时间，τ为在方位向慢时间t时的从EPC到第m个像素的回波的双程延迟的距离向快时间，m是图像域的第m个像素；则回波的非均匀时间采样可以表示为sr(τ(tk,n,m),tk,n)，其中，tk,n为方位向快时间，τ(tk,n,m)为在方位向慢时间tk,n时的从EPC到第m个像素的回波的双程延迟的距离向快时间；回波的均匀时间采样可以表示为sr(τ(tk′,n′,m),tk′,n′)，其中，tk′,n′为方位向快时间，τ(tk′,n′,m)为在方位向慢时间tk′,n′时的从EPC到第m个像素的回波的双程延迟的距离向快时间；
步骤2、初始化自适应的高分宽幅SAR清晰重构成像的观测场景目标空间参数：
初始化自适应的高分宽幅SAR清晰重构成像的场景空间参数，包括：以地平面作为IDR成像方法BP算法的投影空间，将图像空间划分为大小相等的单元网格点，根据处理区域选择对应数据段，数据段起始位置为天线波束进入处理区域时刻到天线波束滑出处理区域时刻，其总长度为处理区域长度加上合成孔径长度；初始化距离向像素数Nr，方位向像素数量Na，距离向网格分辨率dr，方位向网格分辨率da，总网格数M个；
步骤3、生成原始回波数据，并进行距离向脉冲压缩，得到脉冲压缩后的回波数据：
生成多通道SAR的原始回波数据，采用标准合成孔径雷达回波数据距离向脉冲压缩方法对原始回波数据进行距离向脉冲压缩，得到距离向压缩后的多通道SAR回波数据，记作sr(τ(t,m),t)，即步骤1中的距离向压缩后的回波，其中t为方位向慢时间，τ为在方位向慢时间t时的从EPC到第m个像素的回波的双程延迟的距离向快时间，m是图像域的第m个像素；
步骤4、计算由非均匀采样回波重构的均匀采样回波：
采用公式
Figure FDA0003195541010000021
计算重构系数Ψk,n(tk′,n′)，其中，N为通道数目，tk,n＝kTr+tn,k＝1,2,...,K是方位向上的第n个通道在第k次脉冲重复周期的采样时间点，Tr为脉冲重复时间，K为方位向回波的采样周期总数，tn＝(n-1)·Tb是第n个通道的采样时间偏差，即方位向空间上同时接收信号的N个通道可以分别等效为第一个发射/接收通道(Tx/Rx1)在(n-1)·Tb时间后运动的位置，Tb为两个相邻通道之间的采样时间间隔；
Figure FDA0003195541010000022
为满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间，k′代表脉冲重复周期，取值k'＝1,2,...,K，n′代表通道，取值1,2,...,N；tn＝(n-1)·Tb，n＝1,2,...,N是第n个通道的采样时间偏差；tq为第q个通道中的采样时间偏差，tq＝(q-1)·Tb，q＝1,2,...,N；π为圆周率，PRF为脉冲重复频率，Π为累乘符号；
再采用公式
Figure FDA0003195541010000023
计算得到由非均匀采样的信号重构获得的均匀采样的回波信号sr(τ(tk',n',m),tk',n')，其中N为通道总数目，sr(τ(tk,n,m),tk,n)是步骤1中的回波的非均匀时间采样，sr(τ(tk′,n′,m),tk′,n′)为步骤1中的回波的均匀时间采样；tk,n＝kTr+tn,k＝1,2,...,K是方位向上的第n个通道在第k次脉冲重复周期的采样时间点，Tr为脉冲重复时间，K为方位向回波的采样周期总数，tn＝(n-1)·Tb是第n个通道的采样时间偏差，即方位向空间上同时接收信号的N个通道可以分别等效为第一个发射/接收通道(Tx/Rx1)在(n-1)·Tb时间后运动的位置，Tb为两个相邻通道之间的采样时间间隔；
Figure FDA0003195541010000031
为满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间，k′代表脉冲重复周期，取值k'＝1,2,...,K，n′代表通道，取值1,2,...,N；L表示插值周期，
Figure FDA0003195541010000032
Figure FDA0003195541010000033
为步骤1中回波的非均匀时间采样，其中，tk,n为方位向快时间，τ(tk,n,m)为在方位向慢时间tk,n时的从EPC到第m个像素的回波的双程延迟的距离向快时间；sr(τ(tk′,n′,m),tk′,n′)为步骤1中回波的均匀时间采样，其中，tk′,n′为方位向快时间，τ(tk′,n′,m)为在方位向慢时间tk′,n′时的从EPC到第m个像素的回波的双程延迟的距离向快时间；Ψk,n(tk′,n′)为步骤4声明的重构系数；
步骤5、计算第m个散射点在通道n中的子图像：
记k'与k的差值为l，k'＝k-l，则步骤4中计算的重构系数Ψk,n(tk′,n′)可以记作与tk',n'无关的l·n维度向量的重构系数Ψln；
Figure FDA0003195541010000034
Figure FDA0003195541010000035
计算通道n中的子图像，其中，ω是发射信号中心角频率，m是图像域的第m个像素，sr(τ(tk,n,m),tk,n)是步骤1中的回波的非均匀时间采样，τ(tk,n,m)为在方位向慢时间tk,n时的从EPC到第m个像素的回波的双程延迟的距离向快时间，tk,n＝kTr+tn,k＝1,2,...,K是方位向上的第n个通道在第k次脉冲重复周期的采样时间点，Tr为脉冲重复时间，K为方位向回波的采样周期总数；τ(tk′,n′,m)为在方位向慢时间tk′,n′时，从EPC到第m个像素的回波的双程延迟的距离向快时间，
Figure FDA0003195541010000036
为满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间，k′代表脉冲重复周期，取值k'＝1,2,...,K，n′代表通道，取值1,2,...,N；j为虚数单位；
步骤6、计算成像场景中第m个散射点在重构出的均匀采样通道n′中的BP成像分量：
将步骤4中计算出的由非均匀采样的信号重构获得的均匀采样的信号sr(τ(tk',n',m),tk',n')，采用公式
Figure FDA0003195541010000041
计算成像场景中第m个散射点在重构出的均匀采样通道n′中的BP成像分量，记为In′(m)，其中Ψln为步骤5中的l·n维度向量的重构系数，Iln(m)为步骤5中计算得到的通道n中的子图像，l＝k-k'，m取值为1,2,...,M，M为步骤2中定义的总网格数；
步骤7、计算重构出的均匀采样回波的BP成像：
再采用公式
Figure FDA0003195541010000042
计算第m个散射点的成像结果，其中，In'(m)是步骤6中计算得到的成像场景中第m个散射点的在重构出的均匀采样通道n′中的BP成像分量；
步骤8、计算图像熵目标函数：
步骤7中计算的
Figure FDA0003195541010000043
由步骤6中的公式
Figure FDA0003195541010000044
则与重构系数Ψ相关的第m个散射点的成像结果表达式为
Figure FDA0003195541010000045
采用j＝1,2,...,J代替n′ln，其中J＝N2L，k表示相干积累时回波的那个序号，考虑了不为负数，且要保证重构时取满原始回波ln个方位数；则In′ln(m)可写为Ij(m)，Ψn′ln记为一维向量Ψ＝[ψ1,ψ2,...,ψj,...,ψJ]，
成像场景中第m个散射点的BP成像分量记为
Figure FDA0003195541010000046
记
Figure FDA0003195541010000047
为重构系数的变量估计，则图像熵目标函数可记为
Figure FDA0003195541010000048
其中，
Figure FDA0003195541010000049
Figure FDA00031955410100000417
m为成像场景中第m个散射点，
Figure FDA00031955410100000410
为与重构系数的变量估计
Figure FDA00031955410100000411
相关的第m个散射点的成像结果，M为步骤2中所声明的总网格数，|·|为复数取模运算；
步骤9、计算图像熵目标函数的梯度下降方向：
采用公式
Figure FDA00031955410100000412
计算目标函数图像熵的梯度下降方向
Figure FDA00031955410100000413
其中，
Figure FDA00031955410100000414
是梯度运算，
Figure FDA00031955410100000415
为偏导数运算，
Figure FDA00031955410100000416
为步骤8中重构系数的变量估计，j＝1,2,...,J，J＝N2L为步骤8中定义；
步骤10、设定基于共轭梯度法的最优求解算法的初始参数：
初始化基于共轭梯度法的最优求解算法参数包括：最大迭代次数，记为Nitermax；重构出的子图像，记第m个散射点的BP成像后子图像为Ij(m)，j＝1,2,...,J，将步骤8中定义的重构系数
Figure FDA0003195541010000051
的初始迭代值记为
Figure FDA0003195541010000052
搜索方向的迭代初始值记为d(0)，d(0)的值选择为目标函数图像熵的负梯度方向，
Figure FDA0003195541010000053
其中
Figure FDA0003195541010000054
为梯度运算；n表示Armijo算法中的第n次迭代次数，设n的初始值n＝0；第n次迭代的搜索方向记为d(n)；
在满足n＜Nitermax条件下，重复执行步骤11-步骤13，第n次迭代过程中将会更新
Figure FDA0003195541010000055
值为
Figure FDA0003195541010000056
且n＝n+1；直到n＝Nitermax时，输出
Figure FDA0003195541010000057
作为最优重构系数估计结果，记为
Figure FDA0003195541010000058
步骤11、采用Armijo线性搜索算法计算搜索步长：
求解公式
Figure FDA0003195541010000059
计算出当前第n次迭代的搜索步长λ(n)，其中，
Figure FDA00031955410100000510
为所求目标函数，min为求解其最小值，d(n)为步骤10中定义的第n次迭代的搜索方向，
Figure FDA00031955410100000511
为步骤10中定义的第n次迭代过程中更新的
Figure FDA00031955410100000512
值，ε(·)为图像熵；
步骤12、计算下次迭代的重构系数
Figure FDA00031955410100000513
当n＜Nitermax-1时，采用公式
Figure FDA00031955410100000514
计算下次迭代用的重构系数
Figure FDA00031955410100000515
其中，λ(n)为步骤11中计算出的第n次迭代的搜索步长，
Figure FDA00031955410100000516
为步骤10中定义的第n次迭代过程中更新的
Figure FDA00031955410100000517
值，n为步骤10中定义的第n次迭代，Nitermax为步骤10中定义的最大迭代次数；
当n＝Nitermax-1时，记最优重构系数估计结果
Figure FDA00031955410100000518
返回步骤10，输出
Figure FDA00031955410100000519
作为最优重构系数估计结果，其中，n为步骤10中定义的第n次迭代，Nitermax为步骤10中定义的最大迭代次数，
Figure FDA00031955410100000520
为步骤10中定义的第n次迭代过程更新后的
Figure FDA00031955410100000521
值；
步骤13、更新搜索方向：
采用公式
Figure FDA00031955410100000522
计算更新后的搜索方向d(n+1)，并更新n＝n+1，其中
Figure FDA00031955410100000523
Figure FDA00031955410100000524
为步骤10中定义的第n次迭代过程更新后的
Figure FDA00031955410100000525
值，
Figure FDA00031955410100000526
为梯度方向，||·||为矩阵范数；
步骤14、计算重构后的清晰图像：
采用步骤8中推导的公式
Figure FDA0003195541010000061
其中，j为步骤8中定义的代替n′ln的变量，j＝1,2,...,J，J＝N2L，将步骤10中计算得到的最优重构系数估计结果
Figure FDA0003195541010000062
中的ψ1,ψ2,...,ψj,...,ψJ，分别作为
Figure FDA0003195541010000063
中的ψj带入公式
Figure FDA0003195541010000064
计算第m个散射点重构后获得的清晰图像IΨ(m)，其中，Ij(m)为步骤8中定义的用j代替n′ln后的第m个散射点在满足均匀时间采样的第n′通道的图像In′ln(m)，n′代表通道，取值1,2,...,N；
经过以上步骤后就得到本发明的第m个散射点重构后获得的清晰图像IΨ(m)。
CN202110889943.4A 2021-08-04 2021-08-04 一种自适应的高分宽幅sar清晰重构成像方法 Active CN113484862B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110889943.4A CN113484862B (zh) 2021-08-04 2021-08-04 一种自适应的高分宽幅sar清晰重构成像方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110889943.4A CN113484862B (zh) 2021-08-04 2021-08-04 一种自适应的高分宽幅sar清晰重构成像方法

Publications (2)

Publication Number Publication Date
CN113484862A true CN113484862A (zh) 2021-10-08
CN113484862B CN113484862B (zh) 2023-10-17

Family

ID=77944458

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110889943.4A Active CN113484862B (zh) 2021-08-04 2021-08-04 一种自适应的高分宽幅sar清晰重构成像方法

Country Status (1)

Country Link
CN (1) CN113484862B (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117849800A (zh) * 2024-03-07 2024-04-09 中国科学院空天信息创新研究院 多角度sar图像序列生成方法、装置、设备及存储介质

Citations (15)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102854505A (zh) * 2012-09-10 2013-01-02 电子科技大学 一种加权稀疏驱动自聚焦sar成像方法
CN103713288A (zh) * 2013-12-31 2014-04-09 电子科技大学 基于迭代最小化稀疏贝叶斯重构线阵sar成像方法
CN104111458A (zh) * 2014-07-29 2014-10-22 西安电子科技大学 基于双重稀疏约束的压缩感知合成孔径雷达成像方法
CN104391295A (zh) * 2014-09-02 2015-03-04 电子科技大学 一种图像熵最优的压缩传感sar稀疏自聚焦成像方法
CN107015225A (zh) * 2017-03-22 2017-08-04 电子科技大学 一种基于自聚焦的sar平台初始高度误差估计方法
CN107037429A (zh) * 2017-04-17 2017-08-11 电子科技大学 基于门限梯度追踪算法的线阵sar三维成像方法
CN107229048A (zh) * 2017-06-06 2017-10-03 电子科技大学 一种高分宽幅sar动目标速度估计与成像方法
CN107748362A (zh) * 2017-10-10 2018-03-02 电子科技大学 一种基于最大锐度的线阵sar快速自聚焦成像方法
CN109061642A (zh) * 2018-07-13 2018-12-21 电子科技大学 一种贝叶斯迭代重加权稀疏自聚焦阵列sar成像方法
CN110109101A (zh) * 2019-04-04 2019-08-09 电子科技大学 一种基于自适应阈值的压缩感知三维sar成像方法
CN110109107A (zh) * 2019-04-24 2019-08-09 电子科技大学 一种合成孔径雷达频域bp算法的运动误差补偿方法
CN110135267A (zh) * 2019-04-17 2019-08-16 电子科技大学 一种大场景sar图像细微目标检测方法
CN111145337A (zh) * 2019-12-13 2020-05-12 电子科技大学 基于分辨率逼近的快速稀疏重构的线阵sar三维成像方法
CN111797717A (zh) * 2020-06-17 2020-10-20 电子科技大学 一种高速高精度的sar图像船只检测方法
CN112837331A (zh) * 2021-03-08 2021-05-25 电子科技大学 一种基于自适应形态重建模糊三维sar图像目标提取方法

Patent Citations (15)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102854505A (zh) * 2012-09-10 2013-01-02 电子科技大学 一种加权稀疏驱动自聚焦sar成像方法
CN103713288A (zh) * 2013-12-31 2014-04-09 电子科技大学 基于迭代最小化稀疏贝叶斯重构线阵sar成像方法
CN104111458A (zh) * 2014-07-29 2014-10-22 西安电子科技大学 基于双重稀疏约束的压缩感知合成孔径雷达成像方法
CN104391295A (zh) * 2014-09-02 2015-03-04 电子科技大学 一种图像熵最优的压缩传感sar稀疏自聚焦成像方法
CN107015225A (zh) * 2017-03-22 2017-08-04 电子科技大学 一种基于自聚焦的sar平台初始高度误差估计方法
CN107037429A (zh) * 2017-04-17 2017-08-11 电子科技大学 基于门限梯度追踪算法的线阵sar三维成像方法
CN107229048A (zh) * 2017-06-06 2017-10-03 电子科技大学 一种高分宽幅sar动目标速度估计与成像方法
CN107748362A (zh) * 2017-10-10 2018-03-02 电子科技大学 一种基于最大锐度的线阵sar快速自聚焦成像方法
CN109061642A (zh) * 2018-07-13 2018-12-21 电子科技大学 一种贝叶斯迭代重加权稀疏自聚焦阵列sar成像方法
CN110109101A (zh) * 2019-04-04 2019-08-09 电子科技大学 一种基于自适应阈值的压缩感知三维sar成像方法
CN110135267A (zh) * 2019-04-17 2019-08-16 电子科技大学 一种大场景sar图像细微目标检测方法
CN110109107A (zh) * 2019-04-24 2019-08-09 电子科技大学 一种合成孔径雷达频域bp算法的运动误差补偿方法
CN111145337A (zh) * 2019-12-13 2020-05-12 电子科技大学 基于分辨率逼近的快速稀疏重构的线阵sar三维成像方法
CN111797717A (zh) * 2020-06-17 2020-10-20 电子科技大学 一种高速高精度的sar图像船只检测方法
CN112837331A (zh) * 2021-03-08 2021-05-25 电子科技大学 一种基于自适应形态重建模糊三维sar图像目标提取方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
ZHOU L等: ""Unambiguous reconstruction for multichannel nonuniform sampling SAR signal based on image fusion"", 《IEEE ACCESS》, no. 8, pages 71558 - 71571, XP011785638, DOI: 10.1109/ACCESS.2020.2987196 *
陈倩 等: ""基于自适应滤波的DPC-MABSAR方位向信号重建"", 《电子与信息学报》, vol. 34, no. 6, pages 1331 - 1336 *
韩玲 等: ""一种边缘增强的高分辨率遥感影像目标检测方法"", 《第七届高分辨率对地观测学术年会论文集》, pages 220 - 231 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117849800A (zh) * 2024-03-07 2024-04-09 中国科学院空天信息创新研究院 多角度sar图像序列生成方法、装置、设备及存储介质
CN117849800B (zh) * 2024-03-07 2024-05-14 中国科学院空天信息创新研究院 多角度sar图像序列生成方法、装置、设备及存储介质

Also Published As

Publication number Publication date
CN113484862B (zh) 2023-10-17

Similar Documents

Publication Publication Date Title
CN107132535B (zh) 基于变分贝叶斯学习算法的isar稀疏频带成像方法
CN103091674B9 (zh) 基于hrrp序列的空间目标高分辨成像方法
CN107229048B (zh) 一种高分宽幅sar动目标速度估计与成像方法
US6037892A (en) Method for automatic focusing of radar or sonar imaging systems using high-order measurements
CN111505639B (zh) 基于变重频采样模式的合成孔径雷达宽幅稀疏成像方法
CN112099007B (zh) 适用于非理想天线方向图的方位向多通道sar模糊抑制方法
CN109507666B (zh) 基于离网变分贝叶斯算法的isar稀疏频带成像方法
CN104020471A (zh) 一种基于分块处理的sar实时成像方法及系统
Callow et al. Stripmap phase gradient autofocus
CN111781595B (zh) 基于匹配搜索和多普勒解模糊的复杂机动群目标成像方法
CN110244303A (zh) 基于sbl-admm的稀疏孔径isar成像方法
Peng et al. Inverse synthetic aperture radar rotation velocity estimation based on phase slope difference of two prominent scatterers
CN111722227B (zh) 基于近似观测矩阵的聚束sar压缩感知成像方法
CN112147608A (zh) 一种快速高斯网格化非均匀fft穿墙成像雷达bp方法
CN112114313A (zh) 结合正交匹配追踪算法的isar稀疏采样成像方法
CN108845318B (zh) 基于Relax算法的星载高分宽幅成像方法
CN113608218B (zh) 一种基于后向投影原理的频域干涉相位稀疏重构方法
EP2817655A2 (en) Systems and methods for image sharpening
CN113484862A (zh) 一种自适应的高分宽幅sar清晰重构成像方法
CN110133646A (zh) 基于nlcs成像的双基前视sar的多通道两脉冲杂波对消方法
Xin et al. ISAR imaging of target with complex motion associated with the fractional Fourier transform
CN108931770B (zh) 基于多维贝塔过程线性回归的isar成像方法
CN114780911B (zh) 一种基于深度学习的海洋宽测绘带距离解模糊方法
CN113985407B (zh) 一种基于解耦原子范数最小化的高精度多带融合方法
Dong et al. High-resolution and wide-swath imaging of spaceborne SAR via random PRF variation constrained by the coverage diagram

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant