CN113484862B - 一种自适应的高分宽幅sar清晰重构成像方法 - Google Patents
一种自适应的高分宽幅sar清晰重构成像方法 Download PDFInfo
- Publication number
- CN113484862B CN113484862B CN202110889943.4A CN202110889943A CN113484862B CN 113484862 B CN113484862 B CN 113484862B CN 202110889943 A CN202110889943 A CN 202110889943A CN 113484862 B CN113484862 B CN 113484862B
- Authority
- CN
- China
- Prior art keywords
- channel
- time
- sampling
- echo
- azimuth
- 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.)
- Active
Links
- 238000003384 imaging method Methods 0.000 title claims abstract description 94
- 230000003044 adaptive effect Effects 0.000 title claims abstract description 22
- 238000000034 method Methods 0.000 claims abstract description 66
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 30
- 230000033001 locomotion Effects 0.000 claims abstract description 16
- 238000002939 conjugate gradient method Methods 0.000 claims abstract description 13
- 238000005070 sampling Methods 0.000 claims description 115
- 239000013598 vector Substances 0.000 claims description 19
- 230000006835 compression Effects 0.000 claims description 18
- 238000007906 compression Methods 0.000 claims description 18
- 239000011159 matrix material Substances 0.000 claims description 14
- 238000012545 processing Methods 0.000 claims description 13
- 238000002592 echocardiography Methods 0.000 claims description 7
- 238000010845 search algorithm Methods 0.000 claims description 3
- 238000009825 accumulation Methods 0.000 claims description 2
- 238000005457 optimization Methods 0.000 abstract description 12
- 238000013441 quality evaluation Methods 0.000 abstract 1
- 238000005516 engineering process Methods 0.000 description 6
- 238000005259 measurement Methods 0.000 description 6
- 230000000737 periodic effect Effects 0.000 description 4
- 208000004350 Strabismus Diseases 0.000 description 2
- 238000004458 analytical method Methods 0.000 description 2
- 230000005540 biological transmission Effects 0.000 description 2
- 230000021615 conjugation Effects 0.000 description 2
- 238000001914 filtration Methods 0.000 description 2
- 238000003672 processing method Methods 0.000 description 2
- 238000004088 simulation Methods 0.000 description 2
- 238000013459 approach Methods 0.000 description 1
- 238000004364 calculation method Methods 0.000 description 1
- 230000008094 contradictory effect Effects 0.000 description 1
- 238000005314 correlation function Methods 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 230000005672 electromagnetic field Effects 0.000 description 1
- 238000004519 manufacturing process Methods 0.000 description 1
- 238000013507 mapping Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000000926 separation method Methods 0.000 description 1
- 238000001228 spectrum Methods 0.000 description 1
- 238000002945 steepest descent method Methods 0.000 description 1
- 230000001629 suppression Effects 0.000 description 1
- 238000003325 tomography Methods 0.000 description 1
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)作为一种具有全天时、全天候、信息量丰富的遥感成像技术，已成为当今对地观测的重要手段。高分辨率和宽幅成像(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天线平台之间的距离历史，然后将对应的时域回波信号进行相干累加，从而重构出每个像素的散射系数其中，σ(m)表示第m个散射点的后向散射系数，m＝1,2,…,M，M表示图像空间划分的总网格点数；s(r(n,m),n)表示第m个散射点在第n个慢时刻经插值后得到的回波值r(n,m)表示第m个散射点在第n个慢时刻到雷达的双程斜距；ω为发射信号的载波频率；N表示一个合成孔径内的慢时刻总数。
通常，使用均匀采样回波的成像场景中第m个像素的BP成像分量可以表示为详见文献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图像，其图像熵可以表示为：
其中，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梯度表示为或gradf，其中/>表示向量微分算子。
对于向量，相对于n×1向量x的梯度算子记作定义为/>
对于矩阵，标量函数f(A)相对于m×n实矩阵A的梯度为一m×n矩阵，简称梯度矩阵，定义为
详细内容可参考文献“微积分”，傅英定，谢云荪主编，高等教育出版社出版。
定义9、标准共轭梯度方法
共轭梯度法是一种重要的求解无约束最优化问题的方法。其基本思想是把共轭性与最速下降方法相结合，利用已知点处的梯度构造一组共轭方向，并沿这组方向进行搜索，求出目标函数的极小点。Fletcher-Reeves共轭梯度法计算步骤如下：
(1)给定初始点x(1)，允许误差ε>0，置
(2)若则停止计算；否则，做一维搜索，求λj，满足
令y(j+1)＝y(j)+λjd(j).
(3)如果j<n，则进行步骤(4)；否则，进行步骤(5).
(4)令其中
置j:＝j+1，转步骤(2).
(5)令x(k+1)＝y(n+1)，y(1)＝x(k+1)，置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准则的数学表达式为：
详细内容可参考文献“最优化理论与算法”，陈宝林等编著，清华大学出版社出版。
定义11、复数取模运算
对于复数、复向量A，其模值
详细内容可参考文献“矩阵理论”，黄廷祝，钟守铭，李正良编著，高等教育出版社出版。
定义12、虚数单位j
虚数单位j定义为二次方程j2+1＝0的两个根的一个。
详细内容可参考文献“电磁场与电磁波”，谢处方主编，高等教育出版社出版。
定义13、偏导数运算
多元函数偏导数的定义是：一个多变量的函数(或称多元函数)，对其中一个变量(导数)微分，而保持其他变量恒定(相对于全导数，在其中所有变量都允许变化)。函数f关于变量x的偏导数写为fx'或写为
详细内容可参考文献“微积分”，傅英定，谢云荪主编，高等教育出版社出版。
定义14、矩阵范数
范数是具有“长度”概念的函数。在线性代数、泛函分析及相关的数学领域，是一个函数，其为向量空间内的所有向量赋予非零的正长度或大小。矩阵范数(matrix norm)亦译矩阵模是数学中矩阵论、线性代数、泛函分析等领域中常见的基本概念，是将一定的矩阵空间建立为赋范向量空间时为矩阵装备的范数。在n维复空间中，最常见范数为：
详细内容可参考文献“矩阵理论”，黄廷祝，钟守铭，李正良编著，高等教育出版社出版。
本发明提供的一种自适应的高分宽幅SAR清晰重构成像方法，它包括以下步骤：
步骤1、初始化自适应的高分宽幅SAR清晰重构成像所需要的多通道SAR系统参数：
初始化自适应的高分宽幅SAR清晰重构成像所需要的多通道SAR系统参数，包括：多通道SAR的观测空间为地面三维坐标系，记为X-Y-Z，其中X表示水平面横轴，Y表示水平面纵轴，Z表示水平垂直轴；多通道SAR的载波频率fcw，信号带宽Bs，距离向采样率fsr，方位向带宽Br；多通道SAR成像系统沿着Y轴平行方向进行匀速直线运动，平台速度矢量为平台飞行高度h；多通道SAR的通道数目为N，记第一个发射信号的通道为Tx，它同时也是接收回波的第一个通道Rx1；其余的通道与第一个通道同时接收回波，依次记为Rx2,Rx3,...,RxN；雷达中心频率fc，雷达系统发射信号的脉冲重复频率PRF，发射信号脉冲重复时间/>
设每个通道之间的距离为d，则第n个通道与第一个通道之间距离的基线长度可以表示为xn＝(n-1)·d。由定义3中等效相位中心(EPC)原理可知，通道的等效相位中心之间的间隔为d/2。当PRF，v和d满足时，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为两个相邻通道之间的采样时间间隔；
由于存在通道间隔引起的基线误差，设考虑误差后的通道间隔测量值为记满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间为/>其中，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、计算由非均匀采样回波重构的均匀采样回波：
采用公式计算重构系数Ψ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为两个相邻通道之间的采样时间间隔；/>为满足均匀时间采样的第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中脉冲重复频率，∏为累乘符号；
再采用公式计算得到由非均匀采样的信号重构获得的均匀采样的回波信号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为两个相邻通道之间的采样时间间隔；/>为满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间，k′代表脉冲重复周期，取值k'＝1,2,...,K，n′代表通道，取值1,2,...,N；L表示插值周期，/>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。
采用公式计算通道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个像素的回波的双程延迟的距离向快时间，为满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间，k′代表脉冲重复周期，取值k'＝1,2,...,K，n′代表通道，取值1,2,...,N；j为定义12中的虚数单位；
步骤6、计算成像场景中第m个散射点在重构出的均匀采样通道n′中的BP成像分量：
将步骤4中计算出的由非均匀采样的信号重构获得的均匀采样的信号sr(τ(tk',n',m),tk',n')，采用公式计算成像场景中第m个散射点在重构出的均匀采样通道n′中的BP成像分量，记为In′(m)，其中Ψln为步骤5中的l·n维度向量的重构系数，Iln(m)为步骤5中计算得到的通道n中的子图像，l＝k-k'，m取值为1,2,...,M，M为步骤2中定义的总网格数；
步骤7、计算重构出的均匀采样回波的BP成像：
再采用公式计算第m个散射点的成像结果，其中，In'(m)是步骤6中计算得到的成像场景中第m个散射点的在重构出的均匀采样通道n′中的BP成像分量。
步骤8、计算图像熵目标函数：
步骤7中计算的由步骤6中的公式/>则与重构系数Ψ相关的第m个散射点的成像结果表达式为/>
采用j＝1,2,...,J代替n′ln，其中J＝N2L，考虑了不为负数，且要保证重构时取满原始回波ln个方位数；则In′ln(m)可写为Ij(m)，Ψn′ln记为一维向量Ψ＝[ψ1,ψ2,...,ψj,...,ψJ]，
成像场景中第m个散射点的BP成像分量记为记/>为重构系数的变量估计，则图像熵目标函数可记为/>其中，/>m为成像场景中第m个散射点，/>为与重构系数的变量估计/>相关的第m个散射点的成像结果，M为步骤2中所声明的总网格数，|·|为定义11中的复数取模运算；
步骤9、计算图像熵目标函数的梯度下降方向：
采用公式计算目标函数图像熵的梯度下降方向/>其中，/>是定义8中的梯度运算，/>为定义13中的偏导数运算，/>为步骤8中定义的重构系数的变量估计，j＝1,2,...,J，J＝N2L为步骤8中定义；
步骤10、设定基于共轭梯度法的最优求解算法的初始参数：
初始化基于共轭梯度法的最优求解算法参数包括：最大迭代次数，记为Nitermax；重构出的子图像，记第m个散射点的BP成像后子图像为Ij(m)，j＝1,2,...,J，将步骤8中定义的重构系数的初始迭代值记为/>搜索方向的迭代初始值记为d(0)，d(0)的值选择为目标函数图像熵的负梯度方向，/>其中/>为定义8中的梯度运算；q表示Armijo算法中的第q次迭代次数，设q的初始值q＝0；第q次迭代的搜索方向记为d(q)；在满足q<Nitermax条件下，重复执行步骤11-步骤13，第q次迭代过程中将会更新/>值为/>且q＝q+1；直到q＝Nitermax时，输出/>作为最优重构系数估计结果，记为/>
步骤11、采用Armijo线性搜索算法计算搜索步长：
求解公式计算出当前第q次迭代的搜索步长λ(q)，其中，/>为所求目标函数，min为求解其最小值，d(q)为步骤10中定义的第q次迭代的搜索方向，/>为步骤10中定义的第n次迭代过程中更新的/>值，ε(·)为定义7中的图像熵；
步骤12、计算下次迭代的重构系数
当q<Nitermax-1时，采用公式计算下次迭代用的重构系数其中，λ(q)为步骤11中计算出的第q次迭代的搜索步长，/>为步骤10中定义的第q次迭代过程中更新的/>值，q为步骤10中定义的第q次迭代，Nitermax为步骤10中定义的最大迭代次数；当q＝Nitermax-1时，记最优重构系数估计结果/>返回步骤10，输出/>作为最优重构系数估计结果，其中，q为步骤10中定义的第q次迭代，Nitermax为步骤10中定义的最大迭代次数，/>为步骤10中定义的第q次迭代过程更新后的/>值；
步骤13、更新搜索方向：
采用公式计算更新后的搜索方向d(q+1)，并更新q＝q+1，其中/> 为步骤10中定义的第q次迭代过程更新后的/>值，/>为定义8中的梯度方向，||·||为定义14中的矩阵范数；
步骤14、计算重构后的清晰图像：
采用步骤8中推导的公式其中，j为步骤8中定义的代替n′ln的变量，j＝1,2,...,J，J＝N2L，将步骤10中计算得到的最优重构系数估计结果中的/>分别作为/>中的ψj带入公式/>计算第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轴平行方向进行匀速直线运动，平台速度矢量为平台飞行高度h＝20Km；多通道SAR的通道数目为N＝4，记第一个发射信号的通道为Tx，它同时也是接收回波的第一个通道Rx1；其余的通道与第一个通道同时接收回波，依次记为Rx2,Rx3,...,RxN；雷达中心频率fc＝10GHz，雷达系统发射信号的脉冲重复频率PRF＝700Hz，发射信号脉冲重复时间/>
设每个通道之间的距离为d＝1m，由于存在通道间隔引起的基线误差，设考虑误差后的通道间隔测量值为则第n个通道与第一个通道之间距离的基线长度可以表示为/>由等效相位中心(EPC)原理可知，通道的等效相位中心之间的间隔为当PRF，v和d满足/>时，EPC是均匀分布的；当EPC非均匀分布时，方位向上的第n个通道在第k次脉冲重复周期的采样时间点记为tk,n＝kTr+tn,k＝1,2,...,K，其中，/>为脉冲重复时间，K＝1024为方位向回波的采样周期总数，tn＝(n-1)·Tb是第n个通道的采样时间偏差，即方位向空间上同时接收信号的N个通道可以分别等效为第一个发射/接收通道(Tx/Rx1)在(n-1)·Tb时间后运动的位置，Tb＝0.263ms为两个相邻通道之间的采样时间间隔。
由于存在通道间隔引起的基线误差，设考虑误差后的通道间隔测量值为记满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间为/>其中，k′代表脉冲重复周期，取值k'＝1,2,...,K，K＝1024为方位向回波的采样周期总数，为脉冲重复时间，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、计算由非均匀采样回波重构的均匀采样回波：
采用公式计算重构系数Ψ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为两个相邻通道之间的采样时间间隔；/>为满足均匀时间采样的第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中脉冲重复频率，∏为累乘符号；
再采用公式计算得到由非均匀采样的信号重构获得的均匀采样的回波信号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为两个相邻通道之间的采样时间间隔；/>为满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间，k′代表脉冲重复周期，取值k'＝1,2,...,K，n′代表通道，取值1,2,...,N；L＝6表示插值周期，/>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。采用公式计算通道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个像素的回波的双程延迟的距离向快时间，/>为满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间，k′代表脉冲重复周期，取值k'＝1,2,...,K，n′代表通道，取值1,2,...,N；j为定义12中的虚数单位；
步骤6、计算成像场景中第m个散射点在重构出的均匀采样通道n′中的BP成像分量：
将步骤4中计算出的由非均匀采样的信号重构获得的均匀采样的信号sr(τ(tk',n',m),tk',n')，采用公式计算成像场景中第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成像：
再采用公式计算第m个散射点的成像结果，其中，In'(m)是步骤6中计算得到的成像场景中第m个散射点的在重构出的均匀采样通道n′中的BP成像分量。
步骤8、计算图像熵目标函数：
步骤7中计算的由步骤6中的公式/>则与重构系数Ψ相关的第m个散射点的成像结果表达式可以写为/>采用j＝1,2,...,J代替n′ln，其中J＝N2L＝112，考虑了不为负数，且要保证重构时取满原始回波ln个方位数；则In′ln(m)可写为Ij(m)，Ψn′ln可写为一维向量Ψ＝[ψ1,ψ2,...,ψj,...,ψJ]，成像场景中第m个散射点的BP成像分量可写为/>记/>为重构系数的变量估计，则图像熵目标函数可记为/>其中，/> m为成像场景中第m个散射点，/>为与重构系数的变量估计/>相关的第m个散射点的成像结果，M＝5000000为步骤2中所声明的总网格数，|·|为定义11中的复数取模运算；
步骤9、计算图像熵目标函数的梯度下降方向：
采用公式计算目标函数图像熵的梯度下降方向/>其中，/>是定义8中的梯度运算，/>为定义13中的偏导数运算，/>为步骤8中定义的重构系数的变量估计，j＝1,2,...,J，J＝N2L＝112为步骤8中定义；
步骤10、设定基于共轭梯度法的最优求解算法的初始参数：
初始化基于共轭梯度法的最优求解算法参数包括：最大迭代次数，记为Nitermax＝10；重构出的子图像，记第m个散射点的BP成像后子图像为Ij(m)，j＝1,2,...,J，J＝112为步骤8中定义，将步骤8中定义的重构系数的初始迭代值记为/>根据/>计算初始迭代值/>搜索方向的迭代初始值记为d(0)，d(0)的值选择为目标函数图像熵的负梯度方向，/>其中/>为定义8中的梯度运算；q表示Armijo算法中的第q次迭代次数，设q的初始值q＝0；第q次迭代的搜索方向记为d(q)；在满足q<Nitermax条件下，重复执行步骤11-步骤13，第q次迭代过程中将会更新/>值为/>且q＝q+1；直到q＝Nitermax时，输出/>作为最优重构系数估计结果，记为/>
步骤11、采用Armijo线性搜索算法计算搜索步长：
求解公式计算出当前第q次迭代的搜索步长λ(q)，其中，/>为所求目标函数，min为求解其最小值，d(q)为步骤10中定义的第q次迭代的搜索方向，/>为步骤10中定义的第n次迭代过程中更新的/>值，ε(·)为定义7中的图像熵；
步骤12、计算下次迭代的重构系数
当q<Nitermax-1时，采用公式计算下次迭代用的重构系数其中，λ(q)为步骤11中计算出的第q次迭代的搜索步长，/>为步骤10中定义的第q次迭代过程中更新的/>值，q为步骤10中定义的第q次迭代，Nitermax＝10为步骤10中定义的最大迭代次数；当q＝Nitermax-1时，记最优重构系数估计结果/>返回步骤10，输出/>作为最优重构系数估计结果，其中，q为步骤10中定义的第q次迭代，Nitermax为步骤10中定义的最大迭代次数，/>为步骤10中定义的第q次迭代过程更新后的/>值；
步骤13、更新搜索方向：
采用公式计算更新后的搜索方向d(q+1)，并更新q＝q+1，其中/> 为步骤10中定义的第q次迭代过程更新后的/>值，/>为定义8中的梯度方向，||·||为定义14中的矩阵范数；
步骤14、计算重构后的清晰图像：
采用步骤8中推导的公式其中，j为步骤8中定义的代替n′ln的变量，j＝1,2,...,J，J＝N2L＝112为步骤8中定义，将步骤10中计算得到的最优重构系数估计结果/>中的/>分别作为/>中的ψj带入公式/>计算第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轴平行方向进行匀速直线运动，平台速度矢量为平台飞行高度h；多通道SAR的通道数目为N，记第一个发射信号的通道为Tx，它同时也是接收回波的第一个通道Rx1；其余的通道与第一个通道同时接收回波，依次记为Rx2,Rx3,...,RxN；雷达中心频率fc，雷达系统发射信号的脉冲重复频率PRF，发射信号脉冲重复时间/>
设每个通道之间的距离为d，则第n个通道与第一个通道之间距离的基线长度可以表示为xn＝(n-1)·d；通道的等效相位中心之间的间隔为d/2；当PRF，v和d满足时，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为两个相邻通道之间的采样时间间隔；
由于存在通道间隔引起的基线误差，设考虑误差后的通道间隔测量值为记满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间为/>其中，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、计算由非均匀采样回波重构的均匀采样回波：
采用公式计算重构系数Ψ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为两个相邻通道之间的采样时间间隔；/>为满足均匀时间采样的第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为脉冲重复频率，Π为累乘符号；
再采用公式计算得到由非均匀采样的信号重构获得的均匀采样的回波信号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为两个相邻通道之间的采样时间间隔；/>为满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间，k′代表脉冲重复周期，取值k'＝1,2,...,K，n′代表通道，取值1,2,...,N；L表示插值周期，/>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；
采用公式计算通道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个像素的回波的双程延迟的距离向快时间，为满足均匀时间采样的第n′通道在第k′次脉冲重复周期的采样时间，k′代表脉冲重复周期，取值k'＝1,2,...,K，n′代表通道，取值1,2,...,N；j为虚数单位；
步骤6、计算成像场景中第m个散射点在重构出的均匀采样通道n′中的BP成像分量：
将步骤4中计算出的由非均匀采样的信号重构获得的均匀采样的信号sr(τ(tk',n',m),tk',n')，采用公式计算成像场景中第m个散射点在重构出的均匀采样通道n′中的BP成像分量，记为In′(m)，其中Ψln为步骤5中的l·n维度向量的重构系数，Iln(m)为步骤5中计算得到的通道n中的子图像，l＝k-k'，m取值为1,2,...,M，M为步骤2中定义的总网格数；
步骤7、计算重构出的均匀采样回波的BP成像：
再采用公式计算第m个散射点的成像结果，其中，In'(m)是步骤6中计算得到的成像场景中第m个散射点的在重构出的均匀采样通道n′中的BP成像分量；
步骤8、计算图像熵目标函数：
步骤7中计算的由步骤6中的公式/>则与重构系数Ψ相关的第m个散射点的成像结果表达式为/>
采用j＝1,2,...,J代替n′ln，其中J＝N2L，考虑了不为负数，且要保证重构时取满原始回波ln个方位数；则In′ln(m)可写为Ij(m)，Ψn′ln记为一维向量Ψ＝[ψ1,ψ2,...,ψj,...,ψJ]，
成像场景中第m个散射点的BP成像分量记为记/>为重构系数的变量估计，则图像熵目标函数可记为/>其中，/> m为成像场景中第m个散射点，/>为与重构系数的变量估计/>相关的第m个散射点的成像结果，M为步骤2中所声明的总网格数，|·|为复数取模运算；
步骤9、计算图像熵目标函数的梯度下降方向：
采用公式计算目标函数图像熵的梯度下降方向/>其中，是梯度运算，/>为偏导数运算，/>为步骤8中重构系数的变量估计，j＝1,2,...,J，J＝N2L为步骤8中定义；
步骤10、设定基于共轭梯度法的最优求解算法的初始参数：
初始化基于共轭梯度法的最优求解算法参数包括：最大迭代次数，记为Nitermax；重构出的子图像，记第m个散射点的BP成像后子图像为Ij(m)，j＝1,2,...,J，将步骤8中定义的重构系数的初始迭代值记为/>搜索方向的迭代初始值记为d(0)，d(0)的值选择为目标函数图像熵的负梯度方向，/>其中/>为梯度运算；q表示Armijo算法中的第q次迭代次数，设q的初始值q＝0；第q次迭代的搜索方向记为d(q)；
在满足q＜Nitermax条件下，重复执行步骤11-步骤13，第q次迭代过程中将会更新值为且q＝q+1；直到q＝Nitermax时，输出/>作为最优重构系数估计结果，记为/>
步骤11、采用Armijo线性搜索算法计算搜索步长：
求解公式计算出当前第q次迭代的搜索步长λ(q)，其中，/>为所求目标函数，min为求解其最小值，d(q)为步骤10中定义的第q次迭代的搜索方向，/>为步骤10中定义的第n次迭代过程中更新的/>值，ε(·)为图像熵；
步骤12、计算下次迭代的重构系数
当q＜Nitermax-1时，采用公式计算下次迭代用的重构系数/>其中，λ(q)为步骤11中计算出的第q次迭代的搜索步长，/>为步骤10中定义的第q次迭代过程中更新的/>值，q为步骤10中定义的第q次迭代，Nitermax为步骤10中定义的最大迭代次数；
当q＝Nitermax-1时，记最优重构系数估计结果返回步骤10，输出/>作为最优重构系数估计结果，其中，q为步骤10中定义的第q次迭代，Nitermax为步骤10中定义的最大迭代次数，/>为步骤10中定义的第q次迭代过程更新后的/>值；
步骤13、更新搜索方向：
采用公式计算更新后的搜索方向d(q+1)，并更新q＝q+1，其中 为步骤10中定义的第q次迭代过程更新后的/>值，/>为梯度方向，||·||为矩阵范数；
步骤14、计算重构后的清晰图像：
采用步骤8中推导的公式其中，j为步骤8中定义的代替n′ln的变量，j＝1,2,...,J，J＝N2L，将步骤10中计算得到的最优重构系数估计结果/>中的ψ1,ψ2,...,ψj,...,ψJ，分别作为/>中的ψj带入公式/>计算第m个散射点重构后获得的清晰图像IΨ(m)，其中，Ij(m)为步骤8中定义的用j代替n′ln后的第m个散射点在满足均匀时间采样的第n′通道的图像In′ln(m)，n′代表通道，取值1,2,...,N；
经过以上步骤后就得到第m个散射点重构后获得的清晰图像IΨ(m)。
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 CN113484862A (zh) | 2021-10-08 |
| CN113484862B true 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) |
Families Citing this family (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115436939B (zh) * | 2022-07-13 | 2024-11-12 | 电子科技大学 | 一种用于多通道高分宽幅sar的多级分辨快速成像方法 |
| CN117849800B (zh) * | 2024-03-07 | 2024-05-14 | 中国科学院空天信息创新研究院 | 多角度sar图像序列生成方法、装置、设备及存储介质 |
Citations (15)
| 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图像目标提取方法 |
-
2021
- 2021-08-04 CN CN202110889943.4A patent/CN113484862B/zh active Active
Patent Citations (15)
| 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)
| Title |
|---|
| "Unambiguous reconstruction for multichannel nonuniform sampling SAR signal based on image fusion";Zhou L等;《IEEE Access》(第第8期期);第71558-71571页 * |
| "一种边缘增强的高分辨率遥感影像目标检测方法";韩玲 等;《第七届高分辨率对地观测学术年会论文集》;第220-231页 * |
| "基于自适应滤波的DPC-MABSAR方位向信号重建";陈倩 等;《电子与信息学报》;第第34卷卷(第第6期期);第1331-1336页 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN113484862A (zh) | 2021-10-08 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| US7212149B2 (en) | System, method and computer program product for detecting and tracking a moving ground target having a single phase center antenna | |
| CN107229048B (zh) | 一种高分宽幅sar动目标速度估计与成像方法 | |
| EP3364212A1 (en) | A method and an apparatus for computer-assisted processing of sar raw data | |
| CN111505639B (zh) | 基于变重频采样模式的合成孔径雷达宽幅稀疏成像方法 | |
| CN105699969B (zh) | 基于广义高斯约束的最大后验估计角超分辨成像方法 | |
| CN113484862B (zh) | 一种自适应的高分宽幅sar清晰重构成像方法 | |
| CN110568434B (zh) | 一种多通道匀加速sar动目标二维速度估计方法 | |
| CN105137425B (zh) | 基于卷积反演原理的扫描雷达前视角超分辨方法 | |
| CN113534151B (zh) | 基于离网稀疏贝叶斯学习的双频段isar成像方法 | |
| CN106291543A (zh) | 一种运动平台扫描雷达超分辨成像方法 | |
| CN101710173B (zh) | 一种星机联合双基地合成孔径雷达时域成像方法 | |
| Ahmad et al. | A beamforming approach to stepped-frequency synthetic aperture through-the-wall radar imaging | |
| Chen et al. | Forward looking imaging of airborne multichannel radar based on modified iaa | |
| CN107024684B (zh) | 一种空间高速运动目标干涉式三维成像方法 | |
| CN114325700B (zh) | 一种星载多通道sar动目标成像方法 | |
| CN114089333B (zh) | 一种基于直升机平台的sar振动误差估计及补偿方法 | |
| CN109143235B (zh) | 一种双基前视合成孔径雷达地面运动目标检测方法 | |
| CN112684446B (zh) | 基于最小熵准则的Bi-ISAR横向定标与畸变校正方法 | |
| CN116400353A (zh) | 一种用于多通道雷达前视成像的实合成孔径成像方法 | |
| US12181573B2 (en) | Radar processing method and system | |
| CN110780273B (zh) | 一种混合正则化方位超分辨成像方法 | |
| CN113759371A (zh) | 一种多通道sar复图像域相位和基线误差联合估计方法 | |
| CN115436939B (zh) | 一种用于多通道高分宽幅sar的多级分辨快速成像方法 | |
| CN111239733B (zh) | 一种基于极化方式的短时移动正交信号宽幅成像方法 | |
| CN111025292A (zh) | 一种四发mimo雷达发射波形设计方法 |
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 |