一种改正InSAR大气延迟的新方法

傅拓 ，张书毕

(1.中国矿业大学环境与测绘学院，江苏徐州 221116;2.中国矿业大学江苏省测绘与国土信息工程重点实验室，江苏徐州 221116)

1 引言

InSAR技术是近年来迅速发展起来的极具应用价值的空间对地观测新技术，具有监测精度高、范围大、成本低、空间连续覆盖等优点，为高边坡地质灾害区监测提供了一种新型的监测方法。但由于地质灾害多发生在暴雨频发、地质地貌复杂的区域，特殊的地理位置与气候使得InSAR技术应用中受大气延迟的影响非常严重，大气延迟误差是InSAR技术的主要误差源，严重时还会导致InSAR结果的错误解译［1］，因此，必须要对InSAR观测进行大气延迟改正。

2 InSAR大气延迟改正的研究

自1994年 Massonnet等［2］在利用 InSAR研究1992年发生在美国南加州地区Laders地震时首次发现了InSAR干涉图中的大气效应后，众多学者对消除或减弱InSAR中的大气延迟的方法展开了研究。近年来，InSAR大气延迟改正的理论及方法取得了快速的发展，国内外提出了多种改进InSAR大气延迟的方法。以下介绍几种主流的方法:

2.1 利用MODIS/MERIS水汽产品

MODIS(Moderate Resolution Imaging Spectroradiometer)［3］的近红外水汽产品比世界上最密集的GPS网络——美国南加州综合GPS网(SCIGN)的密度还要高10倍以上。2005年国内许多学者提出融合MODIS和 GPS 数据的 InSAR 干涉图改正方法［4－5］。但MODIS容易受到云的影响，且水汽产品存在系统误差，需要进行改正。

MERIS(Medium Resolution Imaging Spectrometer Instrument)水汽产品与ASAR(Advanced Synthetic Aperture Radar)影像同时获取，且得到的PWV空间分辨率高(300 m×300 m)，精度也比MODIS 还要高［6～7］，但其缺点是受云的影响更大。

2.2 利用PSInSAR技术

常规DInSAR技术进行滑坡监测会受到相位失相关和大气延迟影响，Ferretti［8～9］提出了仅仅跟踪成像区域内雷达散射特性较为稳定的目标而放弃那些失相关严重的分辨单元的方法。这些目标(如建筑物的墙角、屋顶和裸露的岩石等)几乎不受失相关噪声影响，即使在多年时间间隔的干涉对中仍然保持较高的干涉相关性，把这些稳定的目标称之为永久散射体(permanentscatterers，PS)。PSInSAR 技术［10～12］能克服传统InSAR技术中的去相干性和大气效应等困难，但利用PSInSAR技术进行InSAR大气改正的缺点在于，PS点通常分布于市区等人工建筑较多的地方或无植被覆盖的山峰、山脊等有裸露巨石的地方，对于研究山体滑坡等地质灾害受到了很大的限制，要得到准确可靠的结果，对同一区域内的SAR影像数量也有较高要求，一般需要至少30幅影像。

2.3 利用地面GPS观测

GPS具有高精度、全天候、可连续估计大气延迟的优点，GPS气象学的发展为精确反演对流层延迟提供了保障，亦为估计与改正InSAR大气效应提供了一种新的方法，因而利用GPS进行InSAR大气延迟改正的方法得到了越来越多的应用。但SAR影像与GPS测站的空间分辨率差异很大，很大程度上限制了大气延迟分布图的估计精度，为了对SAR影像进行逐像素的大气层延迟改正，必须对GPS获取的大气延迟进行加密插值，Janssen等［13］对反距离加权平均(IDW)、克里格(Kriging)、三次样条3种不同的插值算法用InSAR大气改正的效果进行了比较，结果表明反距离加权平均插值和克里格插值法比三次样条插值算法效果更佳;李志伟［14］提出了一种新的大气校正的方法，通过计算出平均大气延迟，然后利用GPS探测到的大气延迟对平均延迟进行修正，弥补了SAR干涉图水平方向上各向异性的问题。但GPS改正InSAR大气延迟仍然存在限制，当GPS测站附近的气象观测缺失时，大气延迟估计精度会降低，GPS采集数据连续运行站点之间的间隔一般为几十千米不等，测站点过于分散会降低观测精度，增加观测站的密度会受到地理环境和运作成本等因素的限制，即使设立临时性、低成本的观测点，效果仍然非常有限。针对以上几种方法存在的问题，本文提出了一种新的改正InSAR大气延迟的方法。

3 GPS和NCEP FNL数据改正InSAR大气延迟

3.1 NCEP全球分析最终产品FNL(Final)

NCEP/NCAR再分析资料［15］是1991年开始实施的美国国家气象中心(NMC)气候资料同化系统(CDAS)的一个重要部分，其包含了丰富的气象观测资料，可用来精确估计气象参数并提高大气延迟分布图的精度。该系统由NASA资助，以NCAR以前为航空与交通部门开发的用于美国大陆的产品为基础，将卫星资料、模式结果与先进的人工智能技术结合起来，快速确定和预报风暴与潜在扰动气流。

NCEP FNL资料［16］比再分析资料具有更高的时间、空间分辨率，其同化了全球搜集到的几乎所有观测信息，是全球业务模型生成的最后产品，可以认为是NCEP所有产品中用于模式存档的最佳选择。其数据水平格距为1°，水平方向共360×181个格点，垂直方向分为地表层和从1000～10hPa共26个等压层，有海平面气压、表面温度、海冰等9个单层变量和位势高度、温度、东西风分量、南北风分量、相对湿度共5个全层次变量。从1999年7月30日至今的NCEP FNL产品可以免费获取，其空间分辨率为1°×1°，时间分辨率为 6h(http://dss.ucar.edu/datasets/ds083.2/)。NCEP FNL数据集以GRIB格式存储，同样可用GrADS软件进行读取。图1为逐6h再分析资料集的分布图。

图1 逐6h再分析资料集的分布图

3.2 基于“凝固流”假设的InSAR大气改正方法

“凝固流”假设由泰勒于1938年提出［17］，该假设认为在平均风速的驱使下，空气流在平移过程中处于“凝固”状态，即在某一固定点上观测到的时间上的延迟波动是由某一折射率场在平均速度为V的流的驱使下，经过该点时产生的空间波动引起的，“凝固流”假设的主要思想是空气流在平均风速的驱使下进行平流输送。

“凝固流”假设运用到空气中的三维湿折射率Nwet(x，y，z，t)的数学描述可以表示为:(x，y，z)表示空间的三维坐标，t表示时间，u，v分别表示水平方向的纬向风和经向风。

“凝固流”假设研域平均风速U→—空间上为常数，但随着时间变化，则指定空间位置(x，y，z)上的湿折射率可以认为是在之前时刻t，由另一空间位置(x－ut，yvt，z)上的湿折射率通过平均风速U→—=(u，v)的平流传输得到的。

3.3 平均风速的估计

为了利用GPS天顶湿延迟(ZWD)时间序列估计平均风速，假设研究区域内的GPS天顶湿延迟(ZWD)在空间上的统计特性为各向同性，且在时间上服从广义稳态随机过程［18］，则两个GPS测站上时间连续的ZWD间的互相关函数可表示为:

R=［RxRy］T表示 GPS 测站间的空间距离，ρ、ρ+R分别表示GPS测站的空间位置，△t表示ZWD获取的时间间隔，t、t+△t分别表示ZWD的获取时刻，V为平均风速。ZWD间的自相关函数可表示为:

连续时间上的二维空间ZWD场的互相关函数可以通过平均风速与其自相关函数进行转换，将式(2)线性化有:

∇表示Hamilton算子。

利用最小二乘原理，建立GPS ZWD序列的自相关函数模型与互相关之间的关系后，即根据式(4)解算出平均风速V。ZWD在连续观测时刻的互相关函数可表示为:

式中，S(Rk)={(ρa，ρb):ρa－ρb≈Rk};k 表示 GPS测站个数，Rk表示 ZWD 间的空间距离;ρa、ρb分别表示时刻的 ZWD 位置。

利用式(4)，平均风速可以用自相关函数进行估计:

式中，tm、ts分别表示主副SAR影像的获取时刻。

由于地面GPS测站的分辨率低，利用其插值得到的湿延迟分布不能真实反映湿延迟空间分布状况。基于泰勒提出的“凝固流”假设，Onn［19］将这种思想用于InSAR大气改正，并指出这种方法比仅用GPS观测时更为有效。其中，准确地估计研究区域的平均风速是至关重要的。平均风速的大小及方向决定了“扩展”的ZWD控制点的位置，决定了基于这些控制点插值得到的湿延迟分布图的精度，从而决定了该方法对In-SAR大气延迟效应的改正效果。

但是，将研究区域内ZWD的变化全部归因于平均风速的变化是不严密的，GPS ZWD的变化是多方面因素的综合结果，除了与风速有关外，还与研究区域的温度等其他气象因素有关。

3.4 改进的InSAR大气改正方法

针对上述问题，本文提出了融合GPS和NCEP全球分析最终产品FNL(Final)改正InSAR大气延迟的方法。本方法的思想就是:“以时间换空间”。

利用NCEP FNL对研究区域的平均风速矢量进行估计时，首先获取UTC时间00h，06h，12h和18h上的平均风速的u分量和v分量序列，然后将u分量和v分量序列插值到SAR影像过境时刻。在利用NCEP FNL资料获取了研究区域内的平均风速后，即可利用“凝固流”假设，通过式(1)来获取空间分辨率更高的湿延迟控制点。假定GPS天线的位置为(xg，yg)，则在时刻tni该GPS测站获取的ZWD可以平流输送至SAR影像获取时刻ti的新位置(x，y)。式(1)可表示为:

4 结语

针对GPS改正InSAR大气延迟时，地面GPS测站的低空间分辨率影响了其生成的大气湿延迟分布的精度，从而影响了GPS对InSAR大气的改正效果这一问题，本文利用NCEP全球操作分析最终产品FNL，通过时间和空间插值，获取研究区域平均风速大小和方向的方法，提高湿延迟控制点的空间分辨率和差分湿延迟分布图的精度。改进方法与在理论上证明了减少InSAR干涉图中的残余大气相位的可行性，这为以后InSAR在地质灾害区域进行山体滑坡等的监测起到了很重要的推进作用。

［1］Hanssen，R.F.，Radar interferometry:data interpretation and error analysis［M］.xviii，308pp.，Kluwer Academic，Dordrecht;Boston，2001

［2］Massonnet，D.，K.Feigl，M.Rossi，and F.Adragna，Radar interferometric mapping of deformation in the year after the Landers earthquake［J］.Nature，369，227 ～230，1994

［3］Gao，B.C.;Kaufman，Y.J.Water vapor retrievals using Mod-erate Resolution Imaging Spectro－radiometer(MODIS)near－infrared channels［J］.J.Geophys.Res.2003，108，4389～4398

［4］Li，Z.H.;Muller，J.－P.;Cross，P.;Fielding，E.J.Interferometric synthetic aperture radar(InSAR)atmosphericcorrection:GPS，ModerateResolutionImagingSpectroradiometer(MODIS)，and InSAR integration［J］.J.Geophys.Res.2005，110，B03410

［5］Li，Z.W.，2005.Modeling atmospheric effects on repeatpass InSAR measurements［D］.PhD Dissertation.The Hong Kong Polytechnic University，Hong Kong

［6］Zhenhong LI.Correction of Atmospheric Water Vapour Effects on Repeat－Pass SAR Interferometry Using GPS，MODIS and MERIS Data［D］.University College London，2005.5

［7］Li，Z.H.;Muller，J.P.;Cross，P.Tropospheric correction techniques in repeat－pass SAR interferometry［C］.Proceedings of the FRINGE 2003 workshop，ESA ESRIN，Frascati，Italy，1 ～5，December 2003

［8］FerrettiA，PratiC and RoccaF.Perm anentscatterersin SAR interferometry［J］.IEEE Transactions on Geoscience and Rem ote Sensing，2000，38(5):2202 ～ 2212

［9］FerrettiA，PratiC and Rocca F.Perm anentscatterers inSAR interferom etry［J］.IEEE Transactions on Geoscience and Rem ote Sensing，2001，39(1):8 ～20

［10］Ferretti，A.;Prati，C.;Rocca，F.Nonlinear Subsidence Rate Estimation Using Permanent Scatters in Differential SAR Interferometry［J］.IEEE T.Geosci.Remote Sens.2000，38，2202 ～2212

［11］Ferretti，A.;Prati，C.;Rocca，F.Permanent Scatters in SAR interferometry［J］.IEEE T.Geosci.Remote Sens.2001，39，8 ～20

［12］Hooper，A.;Zebker，H.;Segall，P.;Kampes，B.A new method for measuring deformation on volcanoes and other natural terrains using InSAR persistent scatterers［J］.Geophys.Res.Lett.2004，31，L23611

［13］Janssen，V.;Ge，L.L.;Rizos，C.Tropospheric correction to SAR interferometry from GPS observations［J］.GPS Solut.2004，8，140 ～ 151.China，November.1997［C］.［s.l.］:［s.n.］，1997

［14］Li，Z.W.;Ding，X.L，Huang，C.;Wadge，G.;Zheng，D.W.Modeling of atmospheric effects on InSAR measurements by incorporating terrain elevation information［J］.J.Atmos.Terr.Phys.2006，66，1189 ～1194

［15］Kalnay E，M Kanamit su，R Kistler，et al.The NECP/NCAR 40 year reanalysis project［J］.Bull.Amer.Meteor.Soc.，1996，77(3)，437 ～470

［16］邓伟，陈海波，马振升等.NCEP FNL全球分析资料的解码及其图形显示［J］.气象与环境科学，2009，32(3):79～82

［17］Taylor，G.I.，1938.The Spectrum of Turbulence［J］.P.Roy.Soc.Lond.A.Mat.，164，476 ～490

［18］Emardson，T.R.andWebb，F.H.Estimating the motion of atmospheric water vapor using the global positioning system.GPS Solutions［J］.2002，6(1－2):58 ～64

［19］Onn，F.Modeling water vapor using GPS with application to mitigating InSAR atmospheric distortions［D］.Ph.D dissertation，Stanford University，2006;P.176