基于空变运动误差分析的微波光子超高分辨SAR成像方法

陈潇翔 邢孟道

(西安电子科技大学雷达信号处理国家重点实验室 西安 710071)

(西安电子科技大学信息感知技术协同创新中心 西安 710071)

1 引言

SAR具有全天时，全天候的能力，目前已经广泛应用于军事与民用的各个领域。超高分辨一直是SAR成像研究的一个重要方向，超高分辨图像能够表征更加精细的目标，是实现目标识别重要基础。传统机载SAR经过几十年的技术发展与革新，其成像分辨率已从原先的米级发展到亚米级[1—3]。为实现更高的成像分辨率，需要雷达具备发射超宽带线性调频信号的能力，传统微波雷达系统中的电子器件在高频频段性能受限，难以实现超宽带线性调频信号的产生与模数转换，已成为限制超高分辨SAR成像的主要因素[4,5]。近年来提出的微波光子雷达结合光子技术和微波技术，利用光倍频技术实现超宽带雷达波形产生，解决了传统雷达在超大带宽线性度不佳的问题，突破了超高分辨SAR的技术瓶颈。2017年，我国某单位成功研制了10 GHz的微波光子超宽带雷达系统[6]，并进行了车载实验，其分辨率相较于近几年的SAR分辨率提升了近一个数量级。随着分辨率的大幅提升，给SAR成像也带来了新的难题，主要体现在运动误差的估计与空变补偿。

在超高分辨领域，针对运动误差的估计与空变补偿，已经提出了一些有效的方法。现有的运动误差估计方法可以分参数化估计与非参数化估计两类[7]，参数化估计方法的代表为采用MD(Map Drift)方法估计飞机的法向加速度[8]，非参数化估计的代表为采用PGA(Phase Gradient Autofocus)估计点目标的剩余相位[9]。由于非参数化估计方法具有更高的估计精度，因此在超高分辨领域更多的是采用PGA等方法估计点目标的剩余相位，进而反演得到精确的运动误差。为了提升PGA的估计精度与鲁棒性，学者们已经对PGA进行了多次改进，例如文献[10]提出的SPGA(Squint Phase Gradient Autofocus)方法，通过补偿点目标的高阶固有相位，进而提升了估计精度。文献[11]提出了采用MD−PGA的方法，通过MD减少了PGA估计相位在多孔径相位拼接中的误差，提升全孔径的估计精度与鲁棒性。在运动误差的空变补偿方面，文献[12]提出采用CZT变换进行包络的距离空变校正，能够在RCMC(Range Curve Migration Correction)之前完成包络的空变补偿，其方法在0.05 m的机载SAR数据上获得良好的成像效果。文献[13]在求解距离空变误差补偿量时摒弃了1阶泰勒展开近似，在3维误差精确已知时对距离空变误差的补偿量进行了精确求解。文献[14]提出了一种SFT−SATA(Scaled Fourier Transform SubAperture Topography−and Aperture−dependent)算法，该算法在两步MOCO(MOtion COmpensation)后进行运动误差的方位空变补偿，解决了运动误差方位空变对方位聚焦的影响。针对非平地场景，文献[15]通过在场景中设置4个合作观测点，结合精确的DEM(Digit−al Elevation Model)与惯导技术，采用高度分层处理，最后通过CLEAN[16]技术实现多图像融合，实现对货船的高精度成像，分辨率达到0.11 m。

上述运动误差估计与空变补偿算法，在传统SAR情形下均能获得较好的成像效果，不过这些方法在微波光子所能涵盖的超高分辨领域，例如厘米级甚至毫米级，其有效性仍需重新考量。原因在于，在估计方面，MD与PGA等误差估计方法均基于相位误差非空变特性，然而超高分辨引起的空变特性十分严重，如果不先去除误差空变的影响将严重影响估计结果。在补偿方面，基于中心波束平面近似补偿的传统补偿算法没有考虑方位空变对包络的影响，导致徙动校正无法拉直包络，进一步影响方位空变的相位补偿。运动误差的补偿总地来说可以分为两个步骤进行，第1步为惯导数据的补偿，第2步为基于数据估计的误差补偿。由于惯导精度受限，甚至有些情况下没有惯导系统，使得第1步惯导数据补偿无法实现运动误差的全部补偿，即第1步补偿后必定会存在剩余的3维运动误差。同时传统的惯导数据补偿是基于中心波束平面近似实现的，其补偿精度进一步降低，因此第2步基于数据估计的误差补偿是必要的。基于数据估计的误差补偿首先需要通过误差估计获得运动误差，然而传统误差估计算法基于相位非空变假设，其估计结果为1维向量，即在每个方位时刻仅存在一个误差相位值，通过对多距离块分块估计可以获得近似的波束平面误差结果，因此第2步补偿在现有技术下只能进行波束中心平面近似补偿。综上所述，现有的补偿算法可以归结为中心波束平面近似下的补偿，因此，开展关于中心波束平面近似假设下的补偿精度研究，从而确定其可行范围，能够为实现微波光子超高分辨成像提供依据。因此，本文首先分析了中心波束平面补偿算法的剩余空变误差，通过评估包络的剩余空变量与相位的剩余空变量，提出了运动误差空变影响判定准则。接着针对微波光子SAR系统条件满足准则时提出采用频域或时频域方法进行成像，对不满足的情况，采用BP算法进行成像[17]，并对BP算法的成像原理进行介绍说明，分析了采用BP算法进行成像的优势与劣势。然后采用机载SAR飞行参数，通过微波光子SAR成像参数与传统SAR成像参数的对比分析与成像，验证所提判定准则与成像方法的有效性，同时也说明了运动误差空变对微波光子超高分辨SAR成像的影响。最后，对录取的车载10 GHz微波光子超高分辨SAR实测数据进行分析与成像处理，实验结果表明所提方法的有效性。

本文安排如下，第2节建立几何与信号模型；第3节分析中心波束平面补偿算法的剩余空变误差，提出运动误差空变影响判定准则；第4节提出基于空变运动误差分析的微波光子超高分辨SAR成像方法；第5节对所提方法进行了仿真验证；第6节给出车载10 GHz微波光子超高分辨SAR实测数据的分析与成像结果；第7节对本文进行了总结与展望。

2 信号建模

首先，描述雷达录取数据过程的几何关系。为简化下面的分析，假设雷达进行直线平飞运动，并以0时刻雷达所在位置建立3维坐标系，其几何模型如图1所示。在图1中，实线为理想的雷达运动轨迹，虚线表示实际的雷达运动轨迹，两条线在O点相交，v为理想的雷达运动速度，H为飞行高度，用(vta，0，0)表示理想情况下ta时刻的雷达坐标，其中，ta表示方位时刻。向量(ta)=[ex(ta)，ey(ta)，ez(ta)]为雷达在ta时刻的运动误差向量，P为场景中的任意点目标，其坐标为(x0，y0，z0)。

图1 数据录取几何模型Fig.1 Data acquisition geometric model

下面介绍回波数据的录取与预处理过程，在走停模型假设下方位时刻ta可以看成以PRF 进行间隔过距离傅里叶变换、去RVP项与tr=fr/γ等价替换，其中，fr表示Dechirp接收信号的频域表示，tr表示直采等效距离快时间，转换为直采信号的形式，通过转换之后的P点的2维时域信号可以表示为采样的离散点，回波数据可看成是一个2维矩阵。由于微波光子超高分辨SAR系统距离向发射超宽带信号，如果直接对回波数据进行接收，按照奈奎斯特采样定理，采样率极高，现有AD与存储技术无法实现，因此需采用Dechirp方式进行回波接收。为了更好地与传统SAR成像方法结合，本文采用文献[18]的方法，将Dechirp接收的回波信号，依次通

其中，win 表示方位窗函数，¸表示波长，B表示等效距离带宽，tn表示雷达波束中心照射到目标(x0，y0，z0)的方位时间，R(ta;x0，y0，z0)为P点的斜距历程。同时，虽然微波光子可以实现超宽带信号的发射，但由于系统还不够完善，录取的回波数据直接脉压时，无法达到理想效果，因此有必要先对脉压后的数据进行一次距离向的PGA，估计剩余相位误差，距离向PGA的实现过程与方位向类似，只要把数据转置后进行PGA相位估计即可。

3 剩余空变误差分析

使在3维运动误差精确已知时，也无法实现运动误差的精确补偿。需要解释的是，这里的精确补偿指的并不是完全补偿，假设A,B,C为场景中的任意3点，基于中心波束平面补偿的传统补偿算法对B,C用A点的运动误差进行补偿，从而导致补偿不精

运动误差空变补偿是超高分辨成像的关键问题，直接决定了图像的聚焦质量。传统的运动误差补偿算法基于中心波束平面近似，这些补偿方法即确。运动误差对成像的影响随成像分辨率的提升更加突出，传统SAR系统由于成像分辨率在米级或者亚米级，误差的空变对成像的影响往往可以忽略，然而微波光子SAR系统能够实现厘米甚至毫米量级的地物观测，分辨率提升一个量级以上，中心波束平面补偿算法是否能适用于微波光子SAR成像需要进行重新评估，从而判断传统运动补偿算法能否实现微波光子SAR超高分辨成像。本文通过计算中心波束平面补偿算法的剩余空变误差，从而界定传统基于中心波束平面的补偿算法在微波光子SAR成像中的适用范围，目前最优的中心波束平面补偿算法为文献[12]所提方法，其采用调频Z变换(Chirp−Z Transform, CZT)去除了距离包络的线性空变，并对距离空变的相位进行了补偿，其方法已在0.05 m分辨的机载SAR实测数据中取得良好的成像效果。

在理想情况下，按照图1的运动方式，用Echoi，ta(tr)表示ta时刻录取的通过距离脉压之后的数据，假设点P为ta时刻照射到的场景中的任一点目标，其与雷达的瞬时距离为R，则点P的回波将以sinc 函数的形式落在Echoi，ta(tr)中的位置上，并附带一个的相位。当存在运动误差时，飞机会偏离理想的运动轨迹，此时雷达录取的回波表示为Echor，ta(tr)，由于Echoi，ta(tr)/=Echor，ta(tr)，即相应时刻录取的回波数据存在差异，所以需要进行运动补偿才能获得理想的成像结果。传统SAR运动补偿算法对回波包络与相位进行分开补偿，首先通过包络平移与变标操作实现包络空变补偿，接着通过相位补偿，使Echor，ta(tr)与Echoi，ta(tr)尽可能相等，包络空变补偿与相位补偿的顺序可以互换，因为空变的包络补偿并不会引起数据的相位变化，因此不会对剩余空变量分析造成影响，本文采用先包络后相位的空变补偿方式进行剩余空变误差分析。在3维运动误差已知时，由于中心波束平面补偿算法仅在各自的时刻ta处补偿各自的误差，与相邻时刻无关，因此中心波束平面补偿算法的剩余空变误差与时间ta无关，可以仅分析雷达一个时刻的剩余空变误差，为简化分析，设定雷达处于(0，0，0)的位置，场景中的点目标表示为(x，y，z=-H)，运动误差表示为(ex，ey，ez)。

首先，计算中心波束平面补偿算法的包络补偿量。中心波束平面补偿算法的包络补偿量通过中心波束平面近似求得，即忽略方位空变性，因此可以设定场景中的点为(x=0，y，z=-H)。理想情况下，点目标的回波时间可以表示为，实际由于运动误差的存在其位置表示为，包络校正通过对trr进行变换实现，即构建一个函数映射tri=f(trr)。然而，如果直接通过对f(trr)在trr=2Rs/c处进行泰勒展开近似，将无法满足校正后tri在2Rs/c不存在运动误差。因此，为了实现tri在2Rs/c不存在运动误差，构建函数trr=g(tri)，通过tri与trr相等，得到g(tri)的表达式为

由于CZT只能实现线性变标，因此对式(2)在tri=2Rs/c处进行1阶泰勒展开，可以得到g1(tri)的表达式，如下所示

进而可以得到f(trr)的表达式为

从式(4)可以得出，非空变的包络补偿量为

下面计算场景中任意点的剩余包络误差，设目标P的坐标为(x0，y0，z=-H)，则目标P的回波时刻trr0表示为

通过函数tri=f(trr)的包络校正，其最终的聚焦位置表示为f(trr0)，则包络的剩余误差为

剩余相位误差的计算方式与包络剩余误差计算方式类似，都是通过波束中心平面近似进行计算，即由于相位补偿是直接相乘，不需要CZT变标实现，因此无需对相位进行线性拟合。同时，在进行实际相位补偿时，由于包络校正时改变了包络的位置，但未改变相应的相位值，用式(7)来表示包络校正之后的回波信号

其中，(x，y)为满足的场景中所有点。所以进行实际相位补偿时需要考虑包络变动的影响，但这并不影响剩余空变误差进行分析。

每个trr处需要补偿的相位量，可以用式(8)表示

其中，h(trr)表示trr对应的点目标在无运动误差时的位置，tri=trr可得h(trr)的表达式为

trr的每个值可以用g1(tri)进行求解，因此补偿的相位为求解剩余相位时只需采用式(8)，因为式(8)给出了实际录取位置处于trr时刻的所有点目标补偿的相位量，式(10)只是在实现补偿时才需要，与分析剩余误差无关，因此剩余的相位误差为根据SAR成像理论的基本理论，当包络误差小于1/4个距离单元时，可以忽略包络误差的影响，当剩余相位误差小于时，可以忽略相位空变对成像的影响。其中距离的分辨单元用表示，因此包络空变的判定准则如式(13)的Ruleprofile所示。

相位空变的判定准则求解存在很大的不同，因为影响方位聚焦的剩余相位为全孔径的累积误差，不由单一时刻的误差决定。在求解上述剩余相位量时，没有考虑方位时间的影响，随着方位时间的变化，原本在波束边缘的点，会有一个从波束边缘变化到波束中心再变化到波束边缘的过程。因此，通过使得y0等于一个恒定的值y00,x0选取为全孔径范围内的所有点，同时在每个不同的x0处随机调整运动误差(ex，ey，ez)的大小与方向，通过式(11)可以估计方位向运动误差的影响，求得的剩余相位误差历程为

从而可以得到相位空变的判定准则，如式(13)的Rulephase所示。综上，运动误差空变影响的判定准则如式(13)所示。当满足判定条件时，可以认为采用中心波束平面补偿后的剩余运动空变误差不会对成像造成影响，这为微波光子SAR实现超高分辨成像提供了选择成像方法的依据。

在进行实际计算时，因为误差(ex，ey，ez)是随机产生的，为增加鲁棒性需要进行多次φ(x0)的生成，并取平均作为φ(x0)的估计结果。同时由于2次项是影响方位聚焦的主要因素，且线性相位不影响方位聚焦，因此需要对φ(x0)去掉线性项，并用2次函数拟合作为最终的计算结果。

4 成像方法

式(13)中的判定准则给出了中心波束平面补偿算法的可行范围，当设计的微波光子SAR系统满足判定条件时，可以先采用中心波束平面补偿算法进行运动误差空变补偿，然后采用频域校正方法实现徙动校正，例如RMA算法，也可以采用CS, RD等时频域校正方法，最后通过匹配滤波实现方位聚焦实现成像，为表述方便，称其为第1类成像算法，其流程图如图2中左侧框所示。

图2 算法流程图Fig.2 Algorithm flowchart

当设计的微波光子SAR系统不满足式(13)的关系时，剩余的空变将对成像造成影响，无法采用上述方法获得良好的聚焦效果。BP算法是一种时域匹配算法，通过对特定轨迹回波数据进行相干积累获得无近似的高分辨图像，理论上BP算法可以实现对任意轨迹任意分辨率的SAR成像。因此当设计的微波光子SAR系统不满足式(13)的关系时，可以采用BP算法获得超高分辨的精聚焦图像。选取地面作为成像平面，BP算法的成像过程可以表示为

然而由于其计算复杂度相较于上述成像方法高出许多，运算效率较低。为了实现更快速的BP成像，文献[19]提出了FFBP算法，通过递归拼接，实现BP成像的效果并大大降低了运算效率。实际处理时，由于惯导精度不一定能满足成像需要，但其精度往往可以满足包络的分辨需求，因此采用BP成像后只会存在一定的方位向散焦，可以对BP成像后的结果采用非参数化估计方法，例如PGA或者最小熵，对局部图像进行剩余误差估计与补偿，实现图像精聚焦，称其为第2类成像算法，其算法流程图如图2中右侧框所示。通过结合空变影响判定准则与在不同条件下的成像算法，本文提出了基于空变运动误差分析的微波光子超高分辨SAR成像方法，能够针对不同的微波光子SAR系统，在满足成像精度的同时，选取更加快速的成像方式，其算法流程图如图2所示。

5 仿真分析与成像

为了分析所提方法的有效性。在同一组载机飞行参数下，分别设计一组传统SAR成像参数与微波光子SAR成像参数，用空变影响判定准则进行分析，说明运动误差空变对微波光子超高分辨SAR成像的影响。采用第1类与第2类成像算法，分别对传统SAR成像参数与微波光子SAR成像参数进行点仿真成像，验证所提成像方法的有效性。载机飞行参数如表1所示，传统SAR与微波光子SAR的成像参数如表2所示，仿真点选取为场景中心点，这是因为雷达直线飞行进行地域观测时存在方位平移不变性，场景中同一距离单元上的任意点目标波束照射的历程保持一致，即其回波会覆盖斜距历程的一片范围，因此足以表示误差空变对成像的影响。

表1 载机飞行参数Tab.1 Flight parameters

表2 SAR成像参数Tab.2 SAR imaging parameters

图3 传统SAR参数分析Fig.3 Traditional SAR parameter analysis

图4 微波光子SAR参数分析Fig.4 Microwave photonic−based SAR parameter analysis

参数分析的仿真结果如图3与图4所示，图3为传统SAR参数分析结果，图3(a)为添加的3维运动误差，图3(b)为剩余包络空变单元数，图3(c)为剩余空变相位。从图3(b)中可以看出，包络的最大偏移量为0.05个距离单元，远小于1个距离单元，因此剩余的包络空变可以忽略。在图3(c)中，通过2次拟合后，最大的剩余相位差为0.05 rad，其偏差远小于，因此剩余空变相位不会对成像造成影响。图4为微波光子SAR参数分析结果，图4(a)为添加的3维运动误差，误差大小与传统SAR一致，图4(b)为剩余包络空变单元数，图4(c)为剩余空变相位。从图4(b)中可以看出，最大的包络偏移量高达5.5个距离单元，剩余的包络误差空变严重。在图4(c)中，相位的剩余量在4 rad左右，因此剩余的相位空变也将对成像造成影响。

对比图3与图4的结果可以发现，在相同载机飞行情况下，传统雷达系统的剩余空变量不会对成像造成严重影响，而对微波光子超高分辨雷达的影响较大。为了更加直观地体现剩余空变误差对成像的影响，对数据进行成像处理，图5为传统SAR参数成像结果，图6为微波光子SAR参数成像结果。其中，图5(a)与图6(a)表示成像时添加的3维运动误差，仿真选取的运动误差为某次机载实际飞行的惯导参数，由于传统雷达系统参数的成像分辨率较低，相应的合成孔径点数远小于微波光子SAR，所以对于同一个实测机载运动误差，通过spline拟合获得各自的运动误差。图5(b)与图6(b)为采用第1类成像算法获得的成像结果的等高线图，图5(b)存在一定的散焦，但散焦不严重，而图6(b)已经完全散焦，可以看出运动误差的剩余空变对传统SAR参数成像影响不大，而对微波光子SAR成像影响十分严重。图5(c)与图6(c)为采用第2类成像算法获得的成像结果的等高线图，从图中可以看出第2类成像算法均能获得良好的聚焦效果，说明第2类成像算法可以用于超高分辨的微波光子SAR成像。

图5 传统SAR参数成像结果Fig.5 Traditional SAR imaging results

图6 微波光子SAR参数成像结果Fig.6 Microwave photonic−based SAR imaging results

虽然第2类成像算法的聚焦性能较好，然而其计算复杂度较高，下面分析两类成像算法计算量上的差异。决定算法计算效率的主要步骤为FFT与插值，因此分析的时候仅需要对FFT与插值引起的计算量进行分析。首先指出N点FFT的复数乘法次数为log2(N2/2)，复数加法次数为log2N2，同样CZT的实现采用FFT实现，这里认为计算量与FFT一致，插值均采用16点插值实现，其计算量为16N次复数加法与16N次复数乘法。假设原始数据的大小为N×N的矩阵，当实现第1类补偿成像时，需要先通过一次距离向的CZT变换，同时实现2维匹配滤波成像，至少需要两次距离向的FFT与两次方位向的FFT，当采用计算复杂度最高的RMA算法进行徙动校正时，会多加一次全数据的插值操作，因此第1类算法的最高计算复杂度为5Nlog2(N2/2)+16N2次复数乘法与5Nlog2N2+16N2次复数加法。第2类算法采用FFBP实现，当成像结果为N×N的图像时，其每一步需要2N2次插值操作，供需log2N步操作，其计算复杂度为32N2log2N次的复数加法与乘法。可以直观地看出第1类算法的复杂度在量级，而第2类算法的复杂度在量级，其复杂度大于第1类算法，因此在可以采用第1类算法时尽可能地避免采用第2类成像算法。

6 实测数据分析与成像

2017年我国某单位成功研制了带宽10 GHz的微波光子超宽带雷达，由于微波光子雷达仍处于研制阶段，没有进行机载实测挂飞，仅进行了车载SAR成像实验，其参数如表3所示。采用本文所提方法进行空变影响判断。由于车载实验时汽车行驶平稳，运动误差相对较小，因此设定最大3维运动误差偏移量为5 cm。通过计算可以得到全孔径波束为40 m，因此可以设定成像场景的大小为半径20 m的圆形区域。

图7 10 GHz车载微波光子雷达参数分析Fig.7 10 GHz microwave photonic−based SAR parameter analysis

图8 雷峰塔微波光子雷达成像结果Fig.8 10 GHz microwave photonic−based SAR imaging results

表3 车载微波光子雷达系统参数Tab.3 SAR system parameters

图7(a)为所设定的3维运动误差。图7(b)为计算得到的剩余包络空变单元数，图中标记了最大的偏移单元，从图中可以看出，最大偏移量在0.18个距离单元，因此补偿后其剩余包络空变不会对成像造成影响。图7(c)为计算得到的剩余相位误差历程，从图中可以看出，最大的剩余空变相位有0.7 rad，处于边界条件，即剩余的相位方位空变对成像影响较小，可根据实际成像的效果进行判断。采用本文所提方法得到的实测数据处理结果如图8所示，其中图8(a)为第1类成像方法的结果图，图8(b)为第2类成像方法的结果图，图8(c)为成像结果的光学对比图，从结果可以看出其在满足成像条件下两类成像方法的成像结果相似，且成像效果较好。图9为成像结果的局部放大图，放大区域为图8中红色框所示部分，可以看出，第2类成像算法略微好于第1类成像算法的结果，与本文仿真分析结果一致，验证了所提方法的有效性。

图9 雷峰塔微波光子雷达成像结果局部方法图Fig.9 Typical area of 10 GHz microwave photonic−based SAR imaging results

7 结束语

针对运动误差空变对实现微波光子雷达超高分辨SAR成像的影响，本文提出了一种基于空变运动误差分析的超高分辨成像方法，为微波光子SAR实现超高分辨成像提供了选择成像方法的依据。本文首先分析了中心波束平面补偿算法的剩余空变误差，通过评估包络的剩余空变量与相位的剩余空变量，提出了运动误差空变影响判定准则。接着针对微波光子SAR系统条件满足准则时提出采用频域或时频域方法进行成像，对不满足的情况，采用BP算法进行成像。然后采用机载SAR飞行参数，通过微波光子SAR成像参数与传统SAR成像参数的对比分析与成像，验证所提判定准则与成像方法的有效性，同时也说明了运动误差空变对微波光子超高分辨SAR成像的影响。最后，对录取的车载10 GHz微波光子超高分辨SAR实测数据进行分析与成像处理，实验结果表明所提方法的有效性。