PSI和CRI联合算法用于苏通大桥基础沉降监测

张宁宁，陈志坚，陈元俊，边 磊

（1.华电电力科学研究院，杭州 310030；2.河海大学 地球科学与工程学院，南京 210098；3.福建省建江水利水电咨询有限公司，福州 350001；4.山东电力工程咨询院有限公司，济南 250013）

1 引 言

苏通大桥桥位区河床松软覆盖层厚约300 m，主塔群桩基础由131根直径2.8 m和2.5 m、长为117 m （北塔）和114 m（南塔）的钻孔灌注桩组成，桩径比仅为2.29，群桩基础沉降问题复杂、涉及面广，加之缺少就近而稳固的参照基准点，致使变形观测工作量大，技术难度高。由于传统的精密大地测量方法外业工作量大、作业时间长，其观测结果和精度受环境因素尤其是潮位影响非常大，从而降低了观测结果的可用性。故综合采用精密微压传感器、剖面沉降和高精度静力水准技术对施工过程中的群桩基础沉降和差异沉降进行了跟踪观测。然而，限于传感器测点（或剖面）数量和观测范围的局限性，有必要辅助采用大尺度的PSI观测技术。

由合成孔径雷达干涉测量（interferometry synthetic aperture radar，即InSAR）扩展而来的差分干涉测量（differential interferometry synthetic aperture radar，D-InSAR）在观测地表的微小形变时，其精度已达到毫米级[1－3]，是目前惟一的一种同时具有全天候、低成本、大覆盖和高精度优点的极具应用潜力的地表形变监测手段。然而，在实际应用中，尤其在植被覆盖地区，常规D-InSAR技术严重受时间、空间失相干的影响，同时由于大气环境对雷达波传播过程的多重干扰，致使D-InSAR技术的应用范围受到限制，而且极大地影响了观测精度[4]。为了弥补常规D-InSAR的局限性，国外学者逐步发展了两种利用高相干性散射体来监测地面微小形变的方法：一种是由意大利Tele-Rilevamento Europa公司 Ferretti等[5－6]提出的永久散射体干涉测量方法（permanent scatterers InSAR，简称PSI）；一种是由德国波茨坦地质中心Xia等[7]提出的人工角反射器（corner reflectors InSAR，简称CRI）方法。已有研究表明，PSI和CRI联合解算方法可以用来监测城市建筑物和大型工程（如大坝、桥梁、管线）的安全性[8－9]。

为了提高研究区的相干性，进而更好地得到桥位区形变结果，本研究借助了大桥自身的永久散射体（PS）特性，同时在重点监测部位设置了人工角反射器（CR），利用PSI和CRI联合解算可提取出桥位区的 PS点，得到不同施工阶段索塔群桩基础的沉降情况。实践证明，随着合成孔径雷达（SAR）数据的积累，借助PSI和CRI联合算法技术可以对大桥的安全性进行实时的观测。

2 PSI和CRI联合算法

2.1 永久散射体（PS）

所谓永久散射体，是指在相当长的时间内对雷达波始终保持强反射特性，而且不因时间和气候的变化而改变的地面物体。通常是尺寸和形状并不统一且不规则的人工地物，如桥梁、金属塔、堤岸等，或者空旷地区正对着雷达波束入射方向的天然地物，如裸露山峰、巨石等。这些地物具有强的回波信号，在图像上表现为一系列亮点或一定形状的亮线[10]。图1为苏通大桥研究区域的SAR影像，已经建成的苏通大桥以亮线形式表现出来。

2.2 人工角反射器（CR）

人工角反射器是利用导电性能和导磁性能良好、电容率大的金属材料制成的二面角或三面角形状，且表面为实体或网状的一种人工几何体。三面CR在雷达图像上是相应于3个面交点的1个亮点，它可以大范围地反射雷达波信号，利用角反射效应将接收到的雷达波束经过几次反射后，形成很强的回波信号，在所获得的SAR图像中呈现了亮度很强的十字星状亮斑，能够很明显地从背景地物中识别出来[11－12]。图 2为苏通大桥南主墩承台上架设的CR（它的表面由铝板做成）。

与PS不同的是，CR的几何形状、尺寸、结构可根据需要人为地控制，因此，其雷达散射截面（radar cross section，简称RCS）可以精确计算出来。

图1 苏通大桥在SAR图像中的PS特征Fig.1 PS characteristic of Sutong Bridge in SAR image

图2 研究区域架设的CRFig.2 CR installed in the study area

RCS是用来度量目标在雷达波照射下所产生的回波强度大小的。RCS常以 m2或 dBm2为单位（σ(dBm2)=10lg[σ(m2)]）。要使三面CR的RCS最大，CR顶点到开口面的中心之间的连线必需与入射波方向平行，其最大RCS的计算如下[13]：式中：L为CR开口的边长；λ为波长。由上式可知，在波长一定的情况下，CR尺寸越大，则RCS越大。例如对C波段而言，边长为1 m和1.2 m的CR的散射截面分别为31.3 dBm2和34.4 dBm2。因此，需要根据工程的具体情况来选择CR的边长。

散射截面的公式表明，散射截面与材料特性没有显式关系。但应选择反射能力强的材料，反射效率越高，实现最大散射截面的可能越大。在各种金属中，铝的反射性仅次于银和金。故在制作时采用3 mm厚铝板外加1 mm厚镀锌铁皮的双层结构。

2.3 PSI和CRI联合解算算法

PSI与CRI联合解算的基本思想是：将覆盖同一地区的多幅（一般>20幅）SAR影像按成像时间先后排序，以选取的公共主影像为基准，将其余所有影像都配准并取样到主影像空间，这样同地区的K+ 1幅SAR影像总共可形成K个干涉对。将所有干涉对逐一进行相位差分处理，可得到K幅干涉相位图。然后借助PSI与CRI联合解算算法提取SAR影像中的 PS目标，将后续处理和分析焦点集中于这些高信噪比的PS点上，PS上相应的不确定性如地形数据误差和大气延迟误差等可以通过一定的算法进行分离，从而提高 PS上地表形变参数的估计精度；其他低信噪比像素集的形变信息可通过内插方法计算得到[14－16]。

在差分干涉图中，PSI与CRI联合解算的差分干涉相位也可表达为式（2），完成差分处理后，便生成时间序列的多时相差分干涉纹图。

式中：φdefo为地表形变相位；φatmos为大气相位；φoffset为随机误差相位；φdem_error为高程误差和外部DEM误差引起的残差相位；φnoise为噪声残差相位。上述各项误差利用传统的 D-InSAR无法进行有效地分离，而通过分析生成的多时相差分干涉纹图中每个像元的幅度稳定性系数，先提取出经过长时间的时间间隔仍具有较好相干性的像元作为PS，然后研究 PS点的相位变化，从这些点的干涉相位中去除 φatmos、φdem_error和视线方向目标物体的偏移值，分析得到可信度高的形变测量值，最终生成地表形变的平均偏移率图，以此来监测微小的基础沉降。

在这一过程中，大气效应的去除、PS点的识别以及其形变速率的求取是非常关键和复杂的步骤。

（1）大气效应的估计与去除。分析表明，确定了 φatmos和φdem_error后，便能够计算式（2）中的各相位分量，从而完成相位解缠。虽然大气效应对每一景SAR影像表现了一种很强的去相关，但通过对长时间序列的多景SAR影像的综合分析，大气效应的影响可以被估计并去除。

φatmos必须用相位纠正来估计和调整，这就需要在监测区域附近选择一个稳定的区域作为大气误差估计的初始值。如果给出一个假定的数字高程模型（DEM）误差，就能够计算出每一景相位中关联的误差 φdem_error；因为每一幅图像获取时间是已知的，所以如果给出一个假定的形变模型，则φdefo的相关值就能计算出来。因此，假定φdefo和 φdem_error的估值存在，那么受噪声相位影响的φatmos+φoffset可以被计算出来。这是因为假设φatmos+φoffset对于所选择的稳定研究区域是常数，而φnoise也可以通过对整个地区进行均值处理来消去。同时，去除大气扰动相位和DEM误差相位之后，PS的时间相干性值能达到一个很高的值，并且这个值是可信的。

（2）PS点视线方向的移动速度的求取。在求取φdefo时，假设地表线性形变速度为常数v，式（2）中的φdefo即为

式中： Ti为时间基线。

对PS点差分干涉相位构建相位系统为

式中：a为常数相位值；d为回波信号相位；pξ和pη为线性相位分量沿方位向ξ和距离向η的坡度值；B为垂直基线距；Δq为地面高程残差相位值；T为时间基线距；E为大气、噪声及PS点非线性动态等残差相位项。

研究表明，在得出Δφi后，利用多幅图像，对式（4）迭代求解，即可求得v，进而得到φdefo。目前，PSI技术只适用于小范围的研究，在PS点上，只要大气效应贡献值估计出来并被去除，就可以计算出毫米级的地表形变。

3 PSI技术在桥位区的应用

3.1 CR的布设

桥梁本身虽具有很好的 PS特性，但对于水域微波会产生镜面反射，雷达天线接收不到回波信号，图像会呈暗黑色调。为了增强桥位区的相干性，在主桥的2个主墩、近塔北辅助墩和辅桥北主墩的桩基础承台上分别安装了边长为1.2、1.2、1、1 m的CR，同时在北岸施工码头以及北岸和南岸的稳定区域分别安装了直角边长为1.2、1、1.2 m的CR。

3.2 桥位区的SAR数据情况

根据苏通大桥基础工程的施工时间表，对欧洲太空局数据进行了查询，最终决定选用的数据类型为Envisat的ASAR SLC影像的降轨数据，影像的成像模式为Image Mode，入射角模式为IS2（入射角 22°），极化方式为 VV。影像覆盖范围的经纬度为 120°25′～121°43′E，31°30′～32°33′N，覆盖面积为 100 km×100 km，桥位区的中心经纬度为120°59′E 和31°47′N，研究范围为 5 km×30 km。

已接收到ESA提供的数据20景，获取日期分别为：2003-11-09、2004-03-28、2004-05-02、2004-08-15、2004-10-24、2006-01-22、2006-02-26、2006-05-07、2006-08-20、2006-09-24、2006-10-29、2006-12-03、2007-03-18、2007-05-27、2007-08-05、2007-09-09、2007-12-23、2008-03-02、2008-07-20、2008-08-24和2008-12-07（后面11景是架设了CR之后的）。从架设CR之后的影像可以清晰地识别出3个CR点，图3为2007-03-18的SAR影像中显示的CR点。

图3 CR在SAR影像中的特征Fig.3 Feature of CR in SAR image

利用EV-InSAR软件对20景架设CR前后的数据两两进行InSAR处理，由这20景数据的基线参数表可以看出，2004-03-28获取的数据与其他数据的基线距最优，因此，在处理过程中，将其作为主图像，其余的 19景数据作为幅图像，采用 SRTM DEM作为外部DEM，借助EV-InSAR中的CTM模块进行D-InSAR处理。最后通过PSI和CRI联合解算技术提取PS点并对其进行形变量计算。

3.3 桥位区的形变分析

由相干系数较高的像对可以得到桥位区较好的差分干涉。其中，图4为2004-03-28至2008-08-24期间的相干性图，其相干性值最大可达0.96，图5为相应的差分干涉图。

为了提取备受关注的桥墩处形变量，本文依据PSI和CRI联合解算方法，选用“双重阈值”PS自动探测算法中的相干系数和振幅离差2个约束条件相结合的方法[17]，搜索该研究区域内的稳定散射目标进行 PS点的选取。成功提取出 2004-03-28至2008-08-24期间的PS点共计1124个，图6为计算得到的PS点形变速率图（颜色深浅表示PS点形变速度大小，负值表示下沉，背景为平均后向散射图像）。

图4 2004-03-28至2008-08-24相干性图Fig.4 Coherent map from 2004-03-28 to 2008-08-24

图5 2004-03-28至2008-08-24差分干涉图Fig.5 Differential interferogram from 2004-03-28 to 2008-08-24

图6 2004-03-28至2008-08-24期间PS点形变速度Fig.6 Deformation speed of PS point from 2004-03-28 to 2008-08-24

显然，在施工过程中，大桥主塔墩的沉降问题最突出，是监测的重点。本文以北索塔为研究对象，承台施工完成后就可以提取到大桥主墩处的PS点，根据获取的PS-1形变情况，对照施工情况，对承台、索塔和桥面铺装等施工过程中产生的群桩基础沉降量进行分析，分析结果如图7所示。

从图可以看出，在2004-03-28至2008-8-24期间，处于苏通大桥主墩处的PS-1形变速度为14.8 mm/a，形变量在0～-78.4 mm之间，形变如此之大的原因在于，该期间苏通大桥正处于建设过程中，随着施加荷载的增大，主墩的沉降也在不断增加。2008年7月大桥正式通车，截止到大桥正常运营，群桩基础的沉降量达到-77.6 mm。

图7 PS-1点形变量分析Fig.7 Deformation analysis of point PS-1

4 结果验证

为了验证该结果的可靠性，结合工程实例，利用ABAQUS计算软件对苏通大桥主塔超大型群桩基础进行了三维非线性有限元分析，建立群桩基础的有限元模型。

4.1 计算范围选取

在群桩基础沉降计算的有限元分析中采用完整的计算模型。承台模型尺寸与实际一致，其平面尺寸为102.7 m×48.1 m。根据有限元计算经验，计算范围取基础平面尺寸的2.5～4.0倍，故土体横桥向方向取400 m，纵桥向方向取200 m，竖直方向的深度通过试算确定，以不影响桩端沉降为原则，竖直方向取174 m。

4.2 计算参数取值

多年的研究表明，特别是在变形相对小的情况下，Duncan-Chang非线性弹性模型模拟桩基工程的变形性状效果较好[18－19]。它通过调整弹性系数近似地考虑了土体的弹塑性变形的特性，并用于增量计算，能反映应力路径对变形的影响。

苏通大桥桥位区位于长江下游潮汐河段，河床覆盖层深厚。桥位区勘察资料显示：桥位区共分为22个工程地质层。桩持力层范围内土层以粉细砂、中粒砂和粉质黏土为主。在现有的计算条件下，22个工程地质层是不能一一考虑的，必须进行概化，土层的概化按相似的原则进行，概化后的材料参数（见表1）可按土层厚度进行加权平均。

表1 土层材料参数Table1 Material parameters of soil layers

桩截面是圆形，利用周长等效将其简化为八角形。桩在计算时不考虑钢护筒和桩身变径的影响，桩和承台的本构模型采用线弹性模型。材料参数采用混凝土与钢筋的等效参数。桩的弹性模量E=35.6 GPa，泊松比λ=0.167，重度γ=24.5 kN/m3，k、n、Rf、G、F、D、Ks为邓-肯张模型中的参数，其中，k、n、Rf均是反映土体弹性模量随应力变化的非线性参数，G、F、D是反映泊松比随应力变化的非线性参数，Ks是反映土体卸载的参数；承台各部分参数按照体积等效原则确定，材料参数见表2。所建概化模型如图8所示。

表2 承台各部分材料参数Table2 Material parameters of pile cap

图8 群桩基础模型剖分网格Fig.8 Mesh generation of pile group foundation model

4.3 计算结果

为了模拟整个施工过程中的沉降情况，数值计算按照以下几个工况进行：桩基施工、承台浇注、索塔施工、箱梁吊装和桥面铺装。承台分5层浇注，作用于基桩的荷载约为1615900 kN，索塔分69层浇筑作用到承台顶面的荷载约为785180 kN，箱梁荷载约为641100 kN，桥面铺装荷载为318627 kN。

根据以上工况分别计算了群桩基础的沉降情况，其中，图9为从桩基施工到通车运营期间的承台沉降云图。

图9 大桥通车运营后承台沉降云图Fig.9 Settlement nephogram of pile cap in the operation of bridge

对比图7、9的计算结果可知，从图9沉降云图中提取到PS-1点的沉降量为-81.2 mm（负值表示下沉），比通过PSI和CRI联合解算得到的沉降量要大，两者绝对误差为3.6 mm，相对误差为4.64%。分析计算值偏大的原因主要有以下几点：①未考虑钢护筒的影响；②未考虑桩底后注浆对桩端土体的影响；③未考虑液态混凝土对桩周土体的挤密作用。如果在计算过程中充分考虑以上几点，两者的结果会有较好的一致性。

5 结 论

超大型深水群桩基础的沉降性状是一个极为复杂的问题，针对这一工程难题，本文进行了一定的研究，得出以下结论：

（1）借助大桥自身的永久散射特性，加上桥位区架设的CR大大增加了整个监测区的相干性，运用PSI和CRI联合算法提取出了整个桥位区大桥从群桩基础施工到通车运营后的沉降量。

（2）与有限元计算结果的对比表明，利用PSI和CRI联合算法进行大型桥梁地基基础沉降监测，可以获得有效可靠的观测结果。

（3）尽管由于本研究中获取的SAR数据有限，可能会影响观测结果的精度，但随着监测区SAR数据的不断累积以及PSI和CRI联合解算算法的不断改进，苏通大桥群桩基础沉降量也将更加准确地观测，为大桥在运营中的安全性评价提供了重要的依据。

该研究成果已在苏通大桥得到很好的验证，解决了处于复杂环境中的类似建筑物基础沉降观测难的问题，具有很好的工程应用前景。

[1]王超,张红,刘智. 星载合成孔径雷达干涉测量[M].北京: 科学出版社,2002: 88－95.

[2]STREZZI T,WEGMULLER U,TOSI L,et al. Land subsidence monitoring with differential SAR interferometry[J]. Photogrammetric Engineering & Remote Sensing,2001,67(11): 1261－1270.

[3]王桂杰,谢谟文,邱骋,等. D-INSAR技术在大范围滑坡监测中的应用[J]. 岩土力学,2010,31(4): 1337－1344.WANG Gui-jie,XIE Mo-wen,QIU Cheng,et al.Application of D-INSAR technique to landslide monitoring[J]. Rock and Soil Mechanics,2010,31(4):1337－1344.

[4]MASSONNET D,RABAUTE T. Radar interferometry:Limits and potential[J]. IEEE Transactions on Geoscience and Remote Sensing,1993,31(2): 455－464.

[5]FERRETTI A,PRATI C,ROCCA F. Nonlinear subsidence rate estimation using permanent scatterers in differential SAR interferometry[J]. IEEE Transactions on Geoscience and Remote Sensing,2000,38(5): 2202－2212.

[6]FERRETTI A,PRATI C,ROCCA F. Permanent scatterers in SAR interferometry[J]. IEEE Transactions on Geoscience and Remote Sensing,2001,39(1): 8－20.

[7]XIA Ye,KAUFMANN H,GUO Xiao-fang. Differential SAR interferometry using corner reflectors[C]//IEEE 2002 International Geoscience and Remote Sensing Symposium (IGARSS02). Toronto: [s. n.],2002: 1243－1246.

[8]李德仁,廖明生,王艳. 永久散射体雷达干涉测量技术[J]. 武汉大学学报(信息科学版),2004,29(8): 664－667.LI De-ren,LIAO Mingsheng,WANG Yan. Progress of permanent scatterer interferometry[J]. Geomatics and Information Science of Wuhan University,2004,29(8):664－667.

[9]ALLIEVI J,AMBROSI C,CERIANI M,et al.Monitoring slow mass movements with the permanent scatterers technique[C]//IEEE 2003 International Geoscience and Remote Sensing Symposium(IGARSS03). Toulouse: [s. n.],2003: 215－217.

[10]COLESANT C,FERRETTI A,PRATI C,et al.Monitoring landslides and tectonic motions with the permanent scatterers technique[J]. Engineering Geology,2003,68(1－2): 3－14.

[11]XIA Y. INSAR activities in central Asia using mobile SAR receiving station[C]//IEEE 2001 International Geoscience and Remote Sensing Symposium(IGARSS01). Sydney: [s. n.],2001: 407－409.

[12]XIA Y,KAUFMANN H,GUO X F. Landslide monitoring in the Three Gorges area using D-InSAR and corner reflectors[J]. Photogrammetric Engineering and Remote Sensing,2004,70(10): 1167－1172.

[13]赵英时. 遥感应用分析原理与方法[M]. 北京: 科学出版社,2003: 134－139.

[14]傅文学,田庆久,郭小方,等. PS技术及其在地表形变监测中的应用现状与发展[J]. 地球科学进展,2006,21(11): 1193－1198.FU Wen-xue,TIAN Qing-jiu,GUO Xiao-fang,et al. The permanent scatterers technique and its application tomonitoring ground deformation[J]. Advances in Earth Science,2006,21(11): 1193－1198.

[15]CHUL J H,MIN K D. Observing coalmining subsidence from JERS-1 permanent scatterer analysis[C]//IEEE 2005 International Geoscience and Remote Sensing Symposium (IGARSS05). Seoul: [s. n.],2005: 4578－4581.

[16]COLESANTI C,MOUELIC S L,BENNANI M,et al.Detection of mining related ground instabilities using the permanent scatterers technique—A case study in the east of France[J]. International Journal of Remote Sensing,2005,26(1): 201－207.

[17]罗小军,黄丁发,刘国祥. 时序差分雷达干涉中永久散射体的自动探测[J]. 西南交通大学学报,2007,42(4):414－418.LUO Xiao-jun,HUANG Ding-fa,LIU Guo-xiang.Automated detection of permanent scatterers in time serial differential radar interferometry[J]. Journal of Southwest Jiaotong University,2007,42(4): 414－418.

[18]魏鸿,唐伯明,董元帅. 基于Biot理论的沙井地基固结沉降分析[J]. 重庆交通大学学报(自然科学版),2009,28(5): 887－891.WEI Hong,TANG Bo-ming,DONG Yuan-shuai.Consolidation settlement of sand drain subgrade based on Biot theory[J]. Journal of Chongqing Jiaotong University(Natural Science),2009,28(5): 887－891.

[19]HE Xu-hui,YU Zhi-wu,CHEN Zheng-qing. Finite element model updating of existing steel bridge based on structural health monitoring[J]. Journal of Central South University (Science and Technology),2008,15:399－403.