1.   引言

星载重轨干涉合成孔径雷达(Synthetic Aperture Radar, SAR)可以快速、大范围地进行地表形变测量，其技术已趋于成熟，正在向实用化的方向发展。而机载重轨干涉SAR可以灵活、机动地对点、线等小面积地物以及滑坡、地震等快速形变进行高精度探测。而且，对于低频SAR系统，机载重轨干涉SAR还能够对林下地形进行测绘。但是，与星载重轨干涉SAR相比，机载重轨干涉SAR仍处于试验阶段。

在机载重轨干涉SAR试验上，德国宇航局(DLR)做了很多开拓性的工作。2006年，Reigber等人^[1]利用TopoSAR系统获取的P波段重轨干涉SAR数据，在去除了残余运动误差的基础上，得到了精度为6 m的数字高程模型(Digital Elevation Model, DEM)。2018年又提出利用双频双基线的方式提取平原地区的高程，其精度约为14 cm^[2]。2003年和2006年，他们利用E-SAR系统在瑞士阿尔卑斯山地区采集了Aletsch冰河的L和P波段的差分干涉SAR(Differential Interferometric Synthetic Aperture Radar, DInSAR)数据，反演出了该地区冰河的2维和3维速度场，得到了L波段的精度比P波段好的结论^[3]。2013年，DLR还发表了一幅德国矿区的L波段机载DInSAR干涉结果^[4]，两次数据采集的时间间隔为6个月，但即使在植被覆盖地区，相干性也非常好。同时，DLR在长时序机载SAR差分干涉方面也进行了一些初步尝试^[5-8]。美国喷气推进实验室(Jet Propulsion Laboratory, JPL)在机载重轨干涉SAR方面也做了不少工作。他们详细分析了进行机载重轨干涉SAR的条件^[9,10]，并使用无人机SAR(Unpiloted Aerial Vehicle SAR, UAVSAR)系统开展了多次重轨飞行试验。2008年5月和7月，他们利用两次飞行的重轨数据得到加州Lost Hills地区由于石油开采而造成的地表形变图^[11]。同年，利用3次重复飞行的L波段机载SAR数据，探测加州San Andreas断层附近的滑坡状况^[12]。2009年5月和6月，他们又在格陵兰和冰岛采集了一系列L波段重轨干涉SAR数据，用于监测这些地区冰河和冰原的变化^[13]。意大利和巴西的研究者们利用OrbiSAR-1机载SAR系统也进行了多次X波段重轨干涉SAR试验^[14,15]，并利用2009年和2010年采集的多次重轨飞行的X和P波段数据反演地表的3维形变，用于评估管道沿线的地质灾害对石油/燃气管道的威胁^[16]。关于国外机载重轨干涉SAR详细的研究进展可以参见文献[17]。最近两年，为了回避频域成像算法固有的平地假设、波束中心造成的误差，有学者提出使用时域后向投影(Back Projection, BP)算法进行成像并进行机载重轨干涉SAR信息提取^[18,19]。

国内由于载机平台的航线控制精度难于满足干涉SAR的要求，在重轨干涉SAR飞行试验开展的数量和成果上与国外相比都存在一定的差距。2010年下半年，中科院电子所利用P波段全极化SAR系统进行了多次重复轨道飞行试验，考虑高程、孔径以及残余运动误差的影响，成功获取了试验区的DEM，与同时获取的X波段双天线干涉SAR系统获取的DEM相比，精度在6 m左右^[20]。近几年，国防科大也利用小型无人机载P波段SAR系统采集了贵州某地区的重轨干涉SAR数据，并提取了该地区的DEM，但由于缺少地面外定标，DEM精度未进行定量化评估。目前，国内机载重轨干涉SAR试验还停留在DEM提取阶段，由于尚未实现载机无人驾驶控制系统的改造，靠飞行员手动控制难以达到差分干涉SAR的基线要求等因素，因而鲜见好的机载差分干涉SAR试验结果。

2.   L波段机载重轨干涉SAR系统

本文介绍的L波段机载SAR系统由双极化有源天线阵面、低功率单元、综合处理单元、电源单元、显控单元、数据记录仪、高精度惯导等组成，该雷达系统的功能如下：低功率单元中的频率源模块产生基准频率，并由波形产生模块生成发射激励信号，进行预功率放大后，送入发射/接收(Transmitter/Receiver, T/R)组件进一步放大后，经天线阵面辐射至指定区域，回波信号经天线阵面接收后，通过T/R组件、接收机后，进行低噪声放大、下变频、正交解调和模数(Analog to Digital, AD)量化，最后将数据送入大容量记录仪存储记录。通过控制系统极化选择开关，交替发射H, V极化信号，双通道同时接收H,V极化数据，在两个脉冲周期内形成HH, HV, VH和VV 4种极化回波数据。对于此次飞行试验，L波段SAR系统搭载在Y-7飞机上，安装在机舱的左腹部，如图1所示。

图  1  L波段机载重轨干涉SAR飞行系统

Figure  1.  L-band SAR system for airborne repeat-pass interferometry

下载: 全尺寸图片 幻灯片

机载SAR重轨干涉实验中，由于环境温度的变化，会出现系统参数的漂移。为此，在飞行中对发射机功率、收发通道增益和相位等重要参数进行内定标。本系统的内定标分为3个方面：发射通道定标、接收通道定标以及参考定标。发射通道定标用于监测和补偿发射信号功率以及包括T/R组件在内的SAR系统发射通道增益和相位的变化；接收定标用于监测和补偿包括T/R组件在内的SAR系统接收通道增益和相位的变化；内定标测量得到的SAR系统收发增益和相位包含了内定标通道及其周边路径的影响，为了消除这种影响，还进行了参考定标。系统内定标环路如图2所示。

图  2  系统内定标环路

Figure  2.  Internal calibration loop of the SAR system

下载: 全尺寸图片 幻灯片

由于重轨干涉SAR对载机的位置、姿态的测量精度要求很高，因此该机载SAR系统配备了高精度的差分全球定位系统/惯性测量单元(Differential Globe Positioning System/Inertial Measurement Unit, DGPS/IMU)组合的惯导系统，用于精确测量载机的位置和姿态，其测量精度如表1所示。

表  1  组合惯导的精度

Table  1.  Precisions of the DGPS/IMU navigation system

	实时测量精度(1σ)	事后处理精度(1σ)
航向测量精度(°)	0.02	0.005
姿态测量精度(°)	0.005	0.002
定位测量精度(m)	3～5	0.02 (水平向)， 0.05 (高度向)
速度测量精度(m/s)	0.03	0.05

下载: 导出CSV 

| 显示表格

3.   机载飞行试验

此次飞行试验中雷达系统参数如表2所示。飞行试验区域位于陕西省渭南市，图3是重复轨道飞行的航线设计图，其中，在A至B段和C至D段雷达采集数据。选择黄河边的一块平坦地区作为定标场，在定标场内沿着距离向布置了10个边长为2.5 m的三面角反射器，其形状如图4所示，其中8个角反的开口朝AB航线方向，如图3中的标注，用于辐射定标和极化定标，两个角反的开口朝向CD航线方向，未在图3中标出。

表  2  L波段SAR系统参数

Table  2.  Parameters for L-band SAR system

系统参数	数值
波长	0.23 m
信号带宽	200 MHz
载机速度	～108 m/s
载机高度	～4600 m
天线长度(方位向，距离向)	1.35 m, 0.3 m
分辨率(方位向，距离向)	1 m, 1 m
采样间隔(方位向，距离向)	～0.6 m, 0.6 m
波束中心视角	50°
极化方式	全极化
雷达发射功率	1000 W
图像信噪比	30 dB
脉冲重复频率(PRF)	1400 Hz

下载: 导出CSV 

| 显示表格

图  3  飞行试验的航线设计图

Figure  3.  Route planning for the SAR campaign

下载: 全尺寸图片 幻灯片

图  4  布设于定标场内的角反射器

Figure  4.  Corner reflectors deployed in the calibrating area

下载: 全尺寸图片 幻灯片

载机绕定标场进行多次重复飞行，如图5所示。SAR系统只在图中直线区域获取数据，单个航迹1次获取的数据长度约为20 km。由于飞行过程中，航迹的控制完全依靠飞行员手动实时调整，所以其平行度难以精确保证，在单次数据获取期间飞行航迹间存在着交叉现象，如图6所示。为了获得较好的DEM精度，L波段的重轨干涉基线应选择在二三十米左右。根据此原则，同时考虑地面情况，在9圈轨迹中选择了较为平行的两段进行重轨干涉SAR处理，图6显示了所选择的两段轨迹的局部，其中红色表示主图像航迹，绿色代表辅图像航迹，蓝色代表未使用的航迹。与航迹对应的试验区1正好位于定标场附近，地势较为平坦，在CD段成像，图像中心位于图5中蓝色$*$号所标记区域，其L波段HH极化幅度图如图7所示。试验区2的数据也是在CD段获取，靠近航线的C端，如图3所示，图像中心位于图5中红色$*$号所标记的区域，该区域地表分布着大大小小的沟壑，如图8所示。表3列出了这两个试验区重轨干涉的基线、获取时间间隔以及成像场景大小。

表  3  试验区重轨数据的参数

Table  3.  Parameters for the repeat-pass SAR data

参数	试验区1	试验区2
基线长(m)	28.1811	10.5500
基线角(rad)	1.99282	3.05920
场景大小(方位向km×斜距向km)	5.0×7.8	3.5×7.8
主/辅图像数据获取时间	12:22/11:55	12:33/12:20

下载: 导出CSV 

| 显示表格

图  5  重轨干涉飞行试验的飞行航迹图

Figure  5.  Flight tracks for the repeat-pass interferometric campaign

下载: 全尺寸图片 幻灯片

图  6  试验区重轨干涉航迹(局部)

Figure  6.  The repeat-pass flight tracks for the test areas (partial)

下载: 全尺寸图片 幻灯片

图  7  试验区1的HH极化幅度图

Figure  7.  The amplitude image of HH polarization for the first test area

下载: 全尺寸图片 幻灯片

图  8  试验区2的HH极化幅度图

Figure  8.  The amplitude image of HH polarization for the second test area

下载: 全尺寸图片 幻灯片

4.   处理过程和方法

4.1   机载SAR信号处理方法

4.1.1   距离向幅相误差校正

假设发射通道定标信号、接收通道定标信号以及参考定标信号分别为${s_{\rm{t}}}(t)$, ${s_{\rm{r}}}(t)$和${s_{\rm{c}}}(t)$，通过式(1)可以计算出包含信号源、发射及接收通道误差在内的线性调频脉冲信号$R(t)$，提取该线性调频脉冲信号的幅度和相位误差来补偿距离向脉冲压缩参考函数，可以提高距离向压缩性能。

对于该L波段SAR系统，根据式(1)得到线性调频脉冲信号的频谱幅度如图9(a)所示，理论上该信号的频谱应为矩形，但误差的存在导致频谱的幅度出现波动。可以通过求频谱幅度的倒数均衡不同频率上的幅度。图9(b)是该线性调频信号中有效频谱范围内存在的3阶及以上相位误差。距离向频谱的幅度和相位误差可以在频域距离压缩时一起去除。对于该L波段机载SAR系统，未进行距离向幅相误差校正时的点目标压缩的结果如图10所示，可以看到距离向一侧副瓣升高，主瓣展宽，并且存在明显的成对回波。而经过距离向幅相误差校正后，点目标的距离向压缩接近理想情况，如图11所示。

图  9  L波段SAR系统线性调频脉冲信号的幅相误差

Figure  9.  The spectrum amplitude and phase errors of the chirp signal for the L-band SAR system

下载: 全尺寸图片 幻灯片

图  10  距离向幅相补偿前的点目标压缩结果

Figure  10.  The focused point image before amplitude and phase compensation

下载: 全尺寸图片 幻灯片

图  11  距离向幅相误差补偿后的点目标压缩结果

Figure  11.  The focused point image after amplitude and phase compensation

下载: 全尺寸图片 幻灯片

4.1.2   机载重轨干涉SAR成像处理方法

有关L波段实际数据成像，大多数采用CS方法^[21]。但本文在处理本系统实际数据时发现，在斜视角达到6°～7°时，传统ECS方法^[22]聚焦的效果很差。而$\omega$-k成像方法在2维频域通过Stolt插值来校正距离、方位耦合与距离频率、方位频率间的依赖关系，是一种更为精确的SAR成像方法，它可以实现宽孔径、大斜视角下的聚焦处理。为了与未来P波段SAR系统处理方法保持一致，该系统采用$\omega$-k成像方法。

实际的机载SAR成像处理需要进行运动补偿处理。$\omega$-k方法也能使用两级运动补偿方法，第2级运动补偿可以通过方位向逆压缩或利用文献[23]中的方法进行处理，两者是等价的。与ECS方法不同，对于$\omega$-k成像方法，2级运动补偿后的运动误差不仅随方位向频率变化，还随着距离向频率变化，即

式中，$\Delta x(r,{k_x},\omega )$是补偿量在方位向偏移的距离，${k_x}$是方位向波数，${\omega _0}$和$\omega$分别是距离向载波的中心频率和实际频率，$r$是目标的最近斜距，c是光速。为了得到高精度的SAR图像，在常规的2级运动补偿后还应进行3级运动补偿。该补偿需要在方位向压缩前的时域、对距离和方位向进行分块处理，分块的大小一般为64$\times$64像元。对于每个数据块，先作2维快速傅里叶变换到方位、距离向频域，再乘以式(3)所示的补偿相位

式中，${x_{\rm{c}}}$,${r_{\rm{c}}}$分别是数据块中心位置处的方位向位置和最近斜距，${r_{\rm{m}}}$是用于初级运动补偿时的参考斜距，${k_{{\rm{xi}}}}$,${\omega _{{\rm{rj}}}}$分别是方位向波数和距离向频率。式(3)中，所有的$\Delta r$表示的都是地面上的目标与载机的实际位置、参考位置间的距离之差。假设多普勒中心为0，那么，第1个$\Delta r$项是地面目标位置$({x_{\rm{c}}} + {x_{{\rm{shift}}}},{r_{\rm{c}}},h)$处与${x_{\rm{c}}} + \Delta x$处载机实际的位置和参考位置间的距离差，其中，

式(3)中，第2个$\Delta r$项是目标位置$({x_{\rm{c}}},{r_{\rm{c}}},{h_{\rm{m}}})$与${x_{\rm c}}$处的载机位置间的距离差，第3项是目标位置$({x_{\rm c}} + {x_{{\rm{shift}}}},{r_{\rm m}},h')$与${x_{\rm c}} + \Delta x$处的载机位置间的距离差，第4项是目标位置$({x_{\rm{c}}},{r_{\rm{m}}},{h_{\rm{m}}})$与${x_{\rm{c}}}$处的载机位置间的距离差。上述$\Delta r$是根据两点间的3维坐标计算的，$h$和$h'$分别是相应点目标位置处的高程，可以从外部DEM中获得。${h_{\rm{m}}}$是成像时的参考高程。因此，式(3)不仅补偿了与孔径有关的误差，同时也补偿了地形引起的相位误差。

机载重轨干涉SAR要求成像时两幅机载SAR图像的参考航迹平行，否则将在生成的干涉纹图中产生沿方位向分布的密集条纹，影响参数反演的正确性。因此，需要在进行主、辅图像成像之前，计算两者的参考轨迹。将所需处理的原始回波脉冲的起始、结束时间信息与后处理的组合惯导中的全球定位系统(Global Positioning System, GPS)时间信息精确对齐，利用组合惯导数据分别拟合该段时间内主辅图像的参考轨迹，同时计算它们的速度。最终的参考轨迹以两幅图像的平均航迹角作为航迹角、以两航迹的平均速度作为速度。然后分别进行$\omega$-k成像。实际轨道计算时，可以采用速度积分得到位置信息，也可以利用经纬度信息经过坐标转换变为大地坐标，从而得到位置信息，这里使用后者。机载重轨干涉SAR成像的处理流程如图12所示。

图  12  机载重轨干涉SAR成像处理流程

Figure  12.  The imaging diagram for airborne repeat-pass interferometric SAR

下载: 全尺寸图片 幻灯片

利用上述成像算法对该L波段SAR数据进行成像，对于重轨干涉，本文只使用了方位向总带宽的0.6倍进行成像。图13是成像区域中一个角反射器的成像结果，可以看到点目标压缩较为理想。

图  13  角反射器的成像结果

Figure  13.  The imaging results for the corner reflectors

下载: 全尺寸图片 幻灯片

4.1.3   残余运动误差的估计和补偿

经过地形和孔径补偿后的SAR图像中仍然存在着由于导航设备精度的限制而引起的载机位置的误差，即所谓的残余运动误差(Residual Motion Error, RME)。未补偿的误差会在生成的SAR图像中产生几何位移和相位误差。与机载双天线干涉SAR系统不同，机载重轨干涉由于两幅SAR图像中的残余运动误差相互独立，干涉计算时不能抵消，因此，由残余运动引起的相位误差会严重影响干涉产品的精度。文献中也提出了多种去除残余运动误差的方法^[24-29]，其中，对于重轨干涉来说，多斜视(MultiSquint technique, MS)方法^[29]最为有效和精确，但在实际处理过程中，主、辅图像间方位向的配准误差很大程度上决定了该方法的精度。

下面利用点目标多斜视方法(Multi-squint Technique with Point Targets, MTPT)^[27]的原理对MS的配准误差进行分析讨论。对于某一点目标，假设它在线性残余运动的影响下，在主、辅图像的两幅相邻子孔径图像m₁, m₂, s₁和s₂中的位置如图14所示。其中，${O_m}$, ${O_s}$是在理想轨迹下该点目标在主、辅图像中的方位位置，即目标的理想位置。它们之间存在$\Delta x$的配准误差。${O_{{m_1}}}$, ${O_{{m_2}}}$是该目标在主图像的第1, 2子孔径图像中的峰值位置，${O_{{s_1}}}$, ${O_{{s_2}}}$是它在辅图像的两子孔径图像中的相应峰值位置。不失一般性，可以认为在4个子孔径图像中，目标的峰值相位都是$\phi$。假设考虑的这6个方位位置都位于所有子孔径图像的共同主瓣内。那么，根据文献[27]可知，主图像两子孔径内残余运动的平均1阶导数可以通过测量目标理想位置${O_m}$处的相位差$\Delta {\phi _m}$获得。这里简单地将$\Delta {\phi _m}$等效为该1阶导数，如图14所示，图中的虚线分别表示第1,2子孔径图像中的方位向相位斜率，它与子孔径图像的多普勒中心直接相关^[30]。同样，辅图像两子孔径内残余运动的平均1阶导数$\Delta {\phi _s}$也可在理想目标位置${O_s}$处获得，如图14所示。另一方面，主、辅图像第1个子孔径的干涉相位$\Delta {\phi _{{\rm{interf1}}}}$也很容易在图14中表示出。需要注意的是，为了减少主、辅图像方位频谱去相干的影响，在成像时主、辅图像使用它们的平均多普勒中心和共同的方位向带宽成像，因此，主、辅图像对应的子孔径图像具有相同的多普勒中心，因而它们沿方位向的相位斜率是相同的。对于第1个子孔径，在该目标的共同主瓣内都可以正确测量主、辅图像的干涉相位，如图14中的阴影部分所示。同样，对于第2个子孔径图像也可在图中标出它们的干涉相位$\Delta {\phi _{{\rm{interf2}}}}$。

图  14  主、辅图像存在配准误差时，MS算法的原理示意图

Figure  14.  The principle diagram of MS when coregistration error exists between the master and slave images

下载: 全尺寸图片 幻灯片

对于MS算法来说，得到第1,2子孔径的干涉相位后，需要将它们相减，以获得主、辅图像残余运动误差差异的1阶导数，即$\Delta {\phi _{{\rm{interf1}}}} - \Delta {\phi _{{\rm{interf2}}}}$。而对于主、辅图像来说，它们的残余运动的1阶导数的差异可以用$\Delta {\phi _m} - \Delta {\phi _s}$表示。根据图14中所测相位之间的关系，容易得出

即由MS方法计算出来的结果与理想数值之间存在$\Delta \phi$的误差，而该误差是由方位向配准误差$\Delta x$引起的。由于方位向的相位斜率可用$2\pi {f_{{\rm{dc}}}}/v$表示^[30]，$v$是载机的速度，${f_{{\rm{dc}}}}$是图像的多普勒中心，那么根据图14中的关系，$\Delta \phi$可表示为

其中，${f_{{\rm{dc1}}}}$, ${f_{{\rm{dc2}}}}$分别是第1, 2子孔径的多普勒中心，${B_{{\rm{sep}}}}$是子孔径的间隔带宽。对于方位向所有的点目标，如果它们的配准误差都相同，那么只是在最终计算的1阶导数中引入了一个常数偏移，而这并不影响残余运动的估计。因为在MS方法中，为了消除该配准误差造成的相位偏置，减去了1阶导数的均值，使得最终的估计中常数项和1阶项待定。利用外部DEM信息可以很容易解算这两个未知数。但是，如果方位向点目标的配准误差并不完全相同，这会在残余运动估计中引入非线性误差。注意，这里的配准误差相对于理想方位位置而言。因此，对于成像在零多普勒频率的情况，只要主、辅图像的理想参考轨迹相互平行，那么对于所有方位位置上的点目标，配准误差都是一致的，与目标的本地视角无关。

而对于较大斜视角的情况下，目标通常成像在多普勒中心位置，即获取坐标系中(Acquisition Doppler, AD)。此时，由于目标在方位向会从零多普勒频率的位置上平移到多普勒中心频率位置。对于主图像，方位向平移量为

其中，${R_m}$是目标的最近斜距，$\beta$是斜视角。对于辅图像，如果成像的斜视角与主图像相同，那么，对于相同的目标在方位向平移时由于干涉基线的影响，目标的最近距离会有变化，其方位向平移量为

其中，$B$, $\alpha$分别是基线和基线角，${B_\parallel }$是基线在视线向的基线分量，$\theta$为视角。假设沿轨基线为0，那么本应该成像在相同方位位置的主、辅图像上的点目标，此时产生的位移为

如果不考虑这种情况，就会认为该位移是由残余运动误差引起的，从而得到错误的估计结果。通过式(9)可知，在大斜视角和大基线成像的情况下，由目标方位向位移引起的配准差异非常明显。在斜视成像的条件下，如果考虑目标的本地视角，主、辅图像方位向的配准偏移量在全图范围内都是不同的。

根据上述讨论，实际数据处理时，如果目标成像在零多普勒位置，由于主、辅图像的参考轨迹平行，所有目标的方位向配准误差都一致，因此，在进行主、辅图像配准时，所有方位向的配准偏移量都使用一个值，这里使用所有方位向配准偏移量的均值；而当成像在多普勒中心位置时，方位向配准偏移量要考虑式(9)的关系。实际数据处理时，为了避免引入不必要的处理误差，应尽量成像在零多普勒中心位置。由于MS方法只能估计出主、辅图像残余运动的差值，因此在补偿时并不能确定各自补偿多少，在实际处理时主图像不变，而将估计出的残余运动都补偿到辅图像中。

4.2   机载干涉处理方法

4.2.1   机载重轨干涉SAR观测模型

一般情况下，重轨干涉SAR观测几何(侧视/斜视)建模如图15所示。主、辅天线经过一定时间间隔对地面同一点${{P}}$进行重复观测，设主图像天线高度为$H$，空间基线及其倾角分别为$B$, $\alpha$, $B$位于YOZ平面，雷达视线在垂直于航迹的平面内的下视角为$\theta$，斜视角为$\gamma$，雷达波长为$\lambda$。在该干涉几何关系下见(图15)，两天线之间的干涉相位 $\varphi$为

图  15  机载重轨干涉SAR观测几何模型

Figure  15.  The observation geometry for airborne repeat-pass interferometric SAR

下载: 全尺寸图片 幻灯片

由于R₁已知，根据式(10)可求得R₂。由二面角定理可得基线B与${R_1}$之间的夹角$\beta$为

进而，由余弦定理得到主辅天线斜距差为

从而得到$\theta$。最后，由主天线观测构型，得到目标高度为

4.2.2   干涉处理方法

通常，干涉数据处理包括单视复(Single-Look Complex, SLC)图像频谱滤波、精配准、干涉图生成、基线估计、平地相位生成与去除、相位滤波、相位解缠、相位转高程及地理编码等^[31]，如图16所示。主辅图像间频域滤波能够降低距离向的频谱偏移去相干影响；图像精配准至少需要达到亚像元以上精度，以降低配准去相干影响；干涉基线主要根据精确的组合惯导参数及主辅图像系统参数进行计算，对于机载重轨运动平台，需要着重考虑时变基线和线性基线的估计与补偿，提高干涉反演的精度；去除平地相位后的干涉相位纹图经滤波和解缠后，得到高程相位；最后，根据干涉测量几何关系，逐点计算得到地面高程信息，经过地理定位及坐标转换等地理编码步骤，最终生成地理坐标系下的DEM产品。

图  16  干涉SAR数据处理步骤

Figure  16.  Data processing steps for interferometric SAR

下载: 全尺寸图片 幻灯片

需要指出的是，主、辅图像的精确配准影响干涉纹图的质量及相干性，在机载重轨干涉SAR测量中，由于平台的不稳定及惯性导航系统精度的局限性，图像对间存在严重扭曲和变形，数据配准及重采样处理需着重考虑^[32]。一方面，当系统处于斜视观测几何时多普勒中心频率不为0，那么，在图像精配准阶段则需要考虑非基带复数据的重采样处理，以正确恢复图像信息。另一方面，由于斜视SAR成像的斜距特性及重轨飞行状态和轨迹的复杂性，更容易引起复杂的图像几何畸变问题，给主、辅图像间的配准带来更多的困难。因此，在本文试验数据处理中，为了对主、辅图像严重的几何畸变关系进行精确拟合，这里采取分块处理的数据配准策略，每一分块配准区域内具有严格的局部畸变对应关系，从而能更精确地建立主、辅图像像素间的映射方程。

5.   试验结果和分析

该L波段SAR系统能够实时调整方位向波束，因此，对于该重轨干涉SAR数据，多普勒中心较小，成像时所有图像都成像在零多普勒位置上。由于缺少该地区精确的DEM，并且由于该区域非常平坦，成像过程中没有使用外部DEM，仅使用了该地区的平均高程。平均高程由角反射器位置上实地测量获得。重轨干涉仅对HH极化数据进行处理。由于该试验区地表较为稳定，地面不存在明显的形变，可以使用MS方法进行残余运动的估计和补偿。MS方法采用4个子孔径计算，当视线向误差的最大值小于0.1 rad或两次循环的计算结果的最大差异小于0.1 rad时停止循环。

5.1   试验区1

图17是试验区1主、辅图像实际轨迹与参考轨迹之间水平向(y方向)和垂直向(z方向)的偏移量，可以看到，在这段轨迹内飞行比较平稳，两者相差不超过5 m。试验区1的相干系数和干涉纹图如图18所示。图19是该试验区相干系数的统计。可以看到，由于该区域两图像获取的时间间隔只有半小时，基线也不长，因此，该区域的相干系数普遍大于0.8。图20是从试验区1中估计出的残余运动误差的水平分量和垂直分量，经过3次循环，MS算法收敛，如图20所示。将3次循环得到的残余运动相加，得到两个方向上总的残余运动${E_y}$, ${E_z}$，根据式(14)计算出残余运动误差在不同距离处的视线向误差$E{}_{{\rm{los}}}$，

图  17  试验区1主、辅图像载机实际航迹与参考轨迹之间的偏移量

Figure  17.  The displacement between the actual and referenced tracks for the first test area

下载: 全尺寸图片 幻灯片

图  18  试验区1的相干系数及干涉纹图

Figure  18.  The coherence and interferometry map for the first test area

下载: 全尺寸图片 幻灯片

图  19  试验区1的相干系数统计

Figure  19.  The coherence plot for the test area

下载: 全尺寸图片 幻灯片

图  20  利用MS方法估计的主辅图像干涉条纹中的高阶残余运动误差

Figure  20.  The high-order residual motion errors derived from the interferometric fringes with MS

下载: 全尺寸图片 幻灯片

其中，$\theta$是本地视角，$\pm$分别表示右、左侧视的情况。对于试验区1，不同距离处视线向残余运动误差如图20(c)所示。图21是残余运动补偿前后，去平地后的干涉纹图的比较。可以看出，补偿前沿方位向存在着明显的相位波动，而且该相位波动与图20(c)的波动趋势相符；而经过残余运动补偿后，该相位波动消除了。图22是该区域残余运动补偿后生成的高程图。图23是残余运动补偿前后反演的高程值的差异。去除图中由河流区域带来的去相干影响，该区域残余运动补偿前后高程值的差异约在–5.8～0.9 m之间。通过6个地面点测量值(其分布如图22所示)与反演的高程值比较，得到该区域的DEM精度约为3.40 m，如表4所示。其中，P5, P6两点为角反点，可以精确测量其地面经纬度信息用于平面定位精度评测，由实测信息P5($34^\circ31'25.8316''$, $109^\circ 30'32.6706''$), P6($34^\circ 31'24.2423''$, $109^\circ30'$$14.9873''$)及平面定位后${\rm P}5'$($34^\circ 31'25.83''$, $109^\circ 30'$$32.42''$), ${\rm P}6'$($34^\circ 31'24.24''$, $109^\circ 30'14.79''$)，得到定位精度分别为5.0 m和6.4 m。

表  4  试验区1干涉反演的DEM精度评价

Table  4.  Evaluation of the derived DEM for the first test area

序号	位置(列,行)	测量值(m)	结果值(m)	误差(m)
1	(2808,2524)	311.0	316.27	–5.27
2	(6324,2380)	312.6	311.92	0.68
3	(7980,2900)	313.5	310.54	2.96
4	(4284,1420)	318.3	317.04	1.26
5	(5200,1092)	315.1	319.39	–4.29
6	(5696, 620)	316.2	314.64	1.56
标准方差	3.40 m

下载: 导出CSV 

| 显示表格

图  21  残余运动误差补偿前后的干涉纹图(去平地效应之后)

Figure  21.  The interferometry map before and after RME compensation for the first test area (after flat-earth removing)

下载: 全尺寸图片 幻灯片

图  22  残余运动补偿后生成的高程图

Figure  22.  The derived DEM after RME compensation

下载: 全尺寸图片 幻灯片

图  23  残余运动补偿前后高程的差异

Figure  23.  The DEM difference map before and after RME compensation

下载: 全尺寸图片 幻灯片

5.2   试验区2

所选择的试验区2存在冲积沟地形，起伏较为明显，如图3和图8所示。由于事先未在该区域布置角反射器，因此干涉过程中的外定标点和高程验证点都是在飞行后去实地补测的。该试验区的相干系数及其统计如图24所示，大部分地区的相干系数都能达到0.8左右。图25是试验区2估计出的残余运动的水平和垂直分量，不同距离处总的残余运动误差如图25(c)所示。残余运动补偿前后该区域的干涉纹图如图26所示，观察图26(a)中沿方位向波动的相位，它与图25(c)中的计算结果颇为一致。经过残余运动补偿后，该相位波动消除了。反演的该区域的DEM如图27所示。通过实测的8个地面点(其分布如图26(b)所示)与反演的高程值比较，得到该区域的DEM精度约为2.85 m，如表5所示。其中，部分测量点的高程误差较大，可能与该成像区域沟壑较多，某些目标相干系数不高有关，同时也与测量点的地面坡度、定位精度相关。

表  5  试验区2干涉反演的DEM精度评价

Table  5.  Evaluation of the derived DEM for the second test area

序号	位置(列,行)	测量值(m)	结果值(m)	误差(m)
1	(10334, 3012)	543.0	544.44	1.44
2	(10392, 2930)	542.0	536.88	–5.12
3	(10398, 1276)	540.4	543.95	3.55
4	(7306, 2620)	509.8	508.65	–1.14
5	(5746, 3970)	515.2	512.98	–2.22
6	(2526, 2276)	463.3	465.24	1.94
7	(2558, 2274)	471.9	468.91	–2.99
8	(1768, 2368)	528.7	528.29	–0.41
标准方差	2.85 m

下载: 导出CSV 

| 显示表格

图  24  试验区2的相干系数

Figure  24.  The coherence for the second test area

下载: 全尺寸图片 幻灯片

图  25  利用MS方法估计的主辅图像干涉条纹中的高阶残余运动误差

Figure  25.  The high-order RME derived from the interferometric fringes with MS

下载: 全尺寸图片 幻灯片

图  26  试验区2残余运动补偿前后的相位图(相位解缠后)

Figure  26.  The interferometry map before and after RME compensation for the second test area (after phase unwrapping)

下载: 全尺寸图片 幻灯片

图  27  试验区2反演的DEM的3维显示图

Figure  27.  3-dimentinal display of the derived DEM for the second test area

下载: 全尺寸图片 幻灯片

6.   结束语

本文在详细研究了L波段重轨干涉SAR系统和试验的基础上，重点分析讨论了本系统重轨干涉SAR处理的过程和方法，尤其是残余运动的估计和补偿，最后对实际采集的L波段重轨干涉SAR数据进行处理和分析。试验结果表明，此次重轨干涉飞行试验获取的数据具有较高的相干性，利用地面实测点目标进行比对，两试验区反演的DEM精度能够分别达到3.40 m和2.85 m。需要指出的是，由于缺少试验区域高精度的DEM作为参考，精度分析只能在有限的地面实测点上进行，而由于实测点的分析结果并不能完全反映整个试验区域的情形，因此，这里的精度分析并不全面，重轨干涉全面的精度分析还需要进一步试验。

本次飞行试验没有采集到有效的差分干涉SAR像对，而且试验区域内也没有明显的发生形变的区域，因此，此次飞行没能进行差分干涉SAR数据分析处理。未来，在条件允许时会进行差分干涉SAR试验，甚至在飞行航迹控制较好的情况下，本系统可用于层析SAR研究。另外，由于本系统是全极化SAR系统，极化干涉SAR研究也是本系统未来应用的一个重要方向。