基于DInSAR技术的2021年比如M_%28S%296.1地震同震形变场提取与分析.pdf

SEISMOLOGICAL AND GEOMAGNETICOBSERVATION AND RESEARCH第 44 卷 第 2 期2023 年 4 月Vol.44 No.2Apr.2023地 震 地 磁 观 测 与 研 究doi:10.3969/j.issn.1003-3246.2023.02.017第一作者简介：桑杰卓玛（1998），女，毕业于长安大学地质工程与测绘学院测绘工程专业，助理工程师，从事地震监测预报相关工作。E-mail：基金项目：中国地震局地震监测、预报、科研三结合课题（项目编号：3JH-202201032）本文收到日期：2022-12-15基于DInSAR技术的2021年比如MS 6.1 地震同震形变场提取与分析桑杰卓玛 索朗南杰 次仁多吉 格桑卓玛（中国拉萨 850000 西藏自治区地震局）摘要 2021 年 3 月 19 日西藏自治区那曲市比如县发生 MS 6.1 地震，因西藏地区地球物理台站少且分布稀疏，无法较好地获取地震发生前后地球物理场变化信息。基于Sentinel-1 卫星升、降轨 SAR 数据及数字高程模型 SRTM3，利用合成孔径雷达差分干涉测量技术（DInSAR）技术，获取震区地震同震形变场信息，判断地表形变场的变化。结果表明：本次地震由正断兼走滑的断层运动造成，升、降轨沿雷达视线向的最大沉降量分别为 2.4 cm（升轨）和 2.8 cm（降轨），最大抬升量分别为 0.6 cm（升轨）和 3.5 cm（降轨）。关键词 地震；DInSAR；地震同震形变场；Sentinel-1 卫星0 引言据中国地震台网测定，2021 年 3 月 19 日 14 时 11 分西藏自治区那曲市比如县发生 MS 6.1 地震（31.94 N，92.74 E），震源深度 10 km。此次地震发生在羌塘地块中部安多 错那地堑东侧附近，距安多南断裂约 16 km（袁伏全等，2021）。由于西藏地区地球物理台站数量少、分布稀疏，且观测手段单一，无法获取地震发生前后较高质量的地球物理场变化信息。研究区域地处高原，植被稀疏、建筑少、人口密度低，适宜采用卫星遥感空对地监测手段进行地表形变监测。合成孔径雷达（SAR）采用主动式微波成像模式，可以全天时、全天候工作且不会受到云、雨的影响，其中合成孔径雷达差分干涉测量（Differential Interferometric Synthetic Aperture Radar，DInSAR）技术通过卫星系统进行空间对地观测，不受时间、空间以及天气状况的限制，具有全天候观测、地表穿透性强等特点，不仅测量精度达毫米级，还可识别潜在的地表沉降信息，在地震形变监测中具有重要作用，目前被广泛用于同震形变场的研究。张艳梅等（2007）认为，地表形变是地震发生的直接表观现象，利用雷达差分干涉技术测量地表形变，对于地震学研究具有重要意义；刘学武等（2011）认为 DInSAR 技术可有效监测快速、激烈的地表形变，并利用该技术获取了玉树 MW 6.9 地震的同震形变场；杨红磊等（2014）利用升降轨 SAR 干涉对，联合 DInSAR 技术和 MAI 技术，采用最小二乘方法求解三维地表形变场；张明等（2017）利用两景 Sentinel-1A 影像，经 DInSAR 技术处理得到135第 2 期2016 年青海门源 MW 5.9 地震的同震形变位移，通过两步法位错反演模型，分析确定门源地震的主要影响区域及断层的滑动分布；宫熙雯等（2021）为获取真实的地表形变，采用联合DInSAR技术和MAI技术的地震三维形变场反演方法，以2017年九寨沟MS 7.0地震为例，融合升降轨 Sentinel-1A 数据，联合解算获取研究区地表在竖直、南北、东西向上的形变，构建地震三维同震形变场。本研究基于 DInSAR 技术，利用研究区域两景升降轨 Sentinel-1 数据和研究区域精细DEM，获取 2021 年 3 月 19 日西藏比如 MS 6.1 地震的地震同震形变场，判断地表形变场的变化。本研究使用的精细 DEM 为 SRTM3，于 2000 年 2 月 11 日由美国发射的“奋进”号航天飞机上搭载的 SRTM 系统完成。SRTM（Shuttle Radar Topography Mission，航天飞机雷达地形测绘使命)由美国太空总署（NASA）和国防部国家测绘局（NIMA）以及德国与意大利航天机构联合测量，数据覆盖全球陆地面积 80%以上，包含中国全境。本文使用SRTM 两种精度中的 3 arc-seconds。1 DInSAR 基本原理DInSAR 是利用同一地区不同时相的合成孔径雷达（SAR）影像，结合获取研究区域精细 DEM 消除干涉图中地形因素的影响，以检测地表微小形变的技术。张威等（2014）研究发现，在地面沉降监测中应用 DInSAR，可以提高监测资料的准确性和可靠性。DInSAR 充分利用 SAR 技术的回波相位信息，将测量精度提高到毫米级，可在不同微波频段、不同极化状态下进行形变识别，是其他测量技术无法比拟的，具有全天候观测、地表穿透性强等特点，可以识别潜在地表沉降信息，目前被广泛运用在地震形变、火山运动、冰川运动、地面沉降以及山体滑坡等监测中。基于 DInSAR 数据，采用二轨法、三轨法和四轨法获取地表形变信息。本研究基于 DInSAR 的二轨法进行形变数据处理，该方法由 Massonnet 等（1993）提出，测量原理见图 1 所示，使用 2 幅 SAR 影像组成一个干涉对，即将研究区域形变之前的主影像和形变之后的辅影像组成生成包含地表信息的干涉图。利用 SRTM3 反演地形相位，并从干涉相位中去除地形相位信息，获取具备形变信息的地形干涉图。2 研究区域及数据源2.1 研究区域2021 年比如 MS 6.1 地震震中位于羌塘块体内部的安多 错那盆地。该盆地是青藏高原中部的一个断陷盆地，总体呈NNE走向，具有东北部高、西南部低的特点。实地调查表明，该盆地由边界正断层控制，其主要边界断裂之一的安多南缘断裂距此次地震震中最近（约图 1 基于 DInSAR 二轨法的形变测量原理Fig.1 Deformation measurement principle based on DInSAR two-track methodS1S2BRP2P1R+R0B2hB1B2B1桑杰卓玛等：基于 DInSAR 技术的 2021 年比如 MS 6.1 地震同震形变场提取与分析地震地磁观测与研究44 卷 13616 km），该断裂为正断层，长约 40 km，走向 N50 60 E，倾向 NW，在全新世存在明显垂直活动，约 2.8 百万年以来平均垂直活动速率最小约 0.19 mm/a（吴中海等，2005）。田雷等（2021）对此次比如地震分析发现，不同研究机构给出的震源机制解差异较大，但矩震级及震源深度结果较为一致，其中矩震级约 MW 5.8，震源深度约 11 km，综合分析认为，该地震属于正断兼走滑型地震。李玉丽等（2022）认为发震区域地质构造复杂，既有NW 走向右旋性质的走滑断裂，也存在 NE 走向的左旋走滑断裂，还发育了近 NS 向、以正断性质为主的张性断裂（图 2）。由于震中位于高海拔地区，地形复杂，地质调查研究程度相对较低，相关研究人员对于发震构造与发震机制存在不同认识。图 2 震中位置示意Fig.2 Schematic diagram of epicenter location3233 N3195 E94939291比如6.1级地震断层f7f6f5f4f3f1f2 断 裂f1 赤布张错扎萨断裂f2 安多色哇断裂f3 错那西南断裂f4 崩错断裂f5 安多北缘断裂f6 安多南缘断裂f7 怒江断裂带2.2 数据源哨兵 1 号（Sentinel-1）卫星于 2014 年 4 月发射成功，卫星轨道高度 693 km，重访周期 12 d。Sentinel-1 是由欧洲航天局和欧洲委员会基于欧洲航天局哥白尼计划（Global Monitoring for Environment and Security，GMES）发射的地球观测卫星，其应用范围大，具备全天时天候的雷达成像系统，具有多种成像模式，如条带模式（SW）、干涉宽幅模式（IW）、超宽辐模式（EW）、波模式（WV），并实现了单极化和双极化等多种极化方式，不受气候及昼夜影响。由于干涉宽幅模式（IW）一般运用于地表形变监测及地形信息提取等（表 1），故本研究使用 Sentinel-1A 卫星 IW 成像模式的斜距单视复数产品（SLC），空间分辨率为 5 m20 m，极化方式为 VV，并选取地震发生前后 4 景升降轨 SAR 数据用于同震形变场反演。SAR 影像基本信息见表 2。137第 2 期表 1 Sentinel-1 卫星的工作模式及应用领域Table 1 Basic working-pattern conditions of Sentinel-1 satellite模式入射角/()幅宽/km地面分辨率极化方式应用领域条带模式2045805 m5 mHH+HV、HH、VH+VV、VV提供应急管理支援等干涉宽幅模式29462505 m20 mHH+HV、HH、VH+VV、VV地表形变监测、地形信息提取等超宽幅模式194740020 m40 mHH+HV、HH、VH+VV、VV需大范围覆盖和短重访周期区域等波模式2235，353820205 m5 mHH、VV海洋参数的获取表 2 SAR 影像基本参数Table 2 SAR imagery basic parameters序号获取时间模式极化方式地面分辨率空间基线/m雷达入射角/()12022-03-07降轨VV5 m20 m-43.70839.58122022-03-1932022-03-12升轨VV5 m20 m-25.90939.57742022-03-243 同震形变场获取及分析3.1 技术原理及处理过程差分干涉图中存在一些因各种相干因素影响而产生的噪声相位，干涉相位条纹越清晰表示测量结果越精确，反之则精确度越低。为使测量结果更为精确，需采用滤波方法提高干涉相位信噪比和数据处理精度。本研究选用Goldstein滤波方法，通过计算局部条纹频谱，设计基于条纹频谱的自适应滤波，对差分干涉图进行滤波处理，得到图 3。Goldstein 滤波方法存在滤波器可变，干涉条纹清晰度高，可减少由时间、空间基线引起的失相干噪声。图 3 滤波后的差分干涉图（a）降轨；（b）升轨Fig.3 Filtered differential interferogram 0(a)(b)200 000400 000kmN0200 000400 000kmN经滤波处理后，需采用相位解缠方法，将缠绕的干涉相位恢复到真实值。相位解缠方法有区域增长法（Region Growing）、最小费用流法（Minimum Cost Flow）、德洛内 MCF（Delaunay MCF）等。其中，最小费用流法在大面积低相干限制增长导致相位解缠困难时，解缠效果较好，本研究即选择该方法进行相位解缠，设置相干性阈值为0.3，所得结果见图4。桑杰卓玛等：基于 DInSAR 技术的 2021 年比如 MS 6.1 地震同震形变场提取与分析地震地磁观测与研究44 卷 138图 5 比如 MS 6.1 地震 LOS 向同震形变场(a)基于降轨数据的同震形变场；(b)基于升轨数据的同震形变场Fig.5 LOS direction coseismic deformation field of the Biru MS 6.1 earthquake 93 100E32 100N32 0031 50092 50092 300N县界形变量km0302010405-0.900.30.81.11.6-2.83.5N县界形变量km0302010405-1.200.20.61.01.5-2.4-0.793 100E32 100N32 0031 50092 50092 300相位解缠后为排除轨道参数不正确导致的轨道误差，需进行轨道精炼后将相位转化为形变，并进行地理编码得到雷达视线向（LOS）的形变结果。图 4 相位解缠图（a）降轨；（b）升轨Fig.4 Phase unwrapping diagram0200 000400 000kmN高：5.4低：-7.70(a)(b)200 000400 000kmN高：10.2低：-0.43.2 同震形变场分析基于差分干涉方法，使用二轨法获取比如 MS 6.1 地震的同震形变场。利用比如原震区SAR 影像数据和外部 DEM SRTM3 数据得到干涉图，将升轨和降轨干涉图进行去除平地效应和地形影响的处理，经滤波及相干系数计算，得到升轨和降轨的相位解缠图，据成像几何及轨道参数把相位转换为高程值完成相位转形变的操作，地理编码后得到LOS向的形变，即传感器方向上的形变。选用哨兵 1 号（Sentinel-1）卫星 2021 年 3 月 7 日、19 日降轨道数据和 3 月 12 日、24 日升轨道数据作为二轨法 DInSAR 的主从影像，选择 VV 极化方式，将制图分辨率设定为 25 m，进行 DInSAR 处理，结果见图 5。形变结果是一个单波段的灰度图像（为使结果更加清晰直观，文中提供彩色图件，见图 5），像元值即为形变量，正值代表朝传感器方向（LOS）运动，即为抬升（表示地面隆升），139第 2 期负值代表远离传感器方向（LOS）运动，即为沉降（表示地面沉降）。由图 5 可见，基于升轨和降轨数据获取的比如 MS 6.1 地震同震形变结果有所差异，应为不同轨道侧视角造成。从降轨结果可以看出，沿 LOS 方向，地震引发的降轨最大抬升量为 3.5 cm，最大沉降量为 2.8 cm，震区存在 2 个方向相反的形变区域，其中北西侧地面沉降，东南侧地面抬升；从升轨结果看，沿 LOS 方向，地震区域主要表现为沉降趋势，其中地面最大沉降量为 2.4 cm，最大抬升量为 0.6 cm，且地震区域北西侧地表存在沉降。降轨结果表现出较为明显的正断层形变特征。地震发生在安多盆地南缘断裂，断裂呈EW 走向。地震震源机制解反演结果显示，此次地震破裂面为正断兼具走滑性质的断层面。结合断层运动与升、降轨 LOS 向沉降、抬升之间的关系，判定本次 MS 6.1 地震同震形变场由断层兼具走滑运动造成。4 结论 震源机制解显示，2021 年 3 月 19 日 14 时 11 分西藏自治区那曲市比如县 MS 6.1 地震破裂面呈正断兼走滑型性质。由于西藏地区地球物理台站分布稀疏且观测手段单一，利用 DInSAR 提取地震同震形变场，可将测量精度提高到毫米级，并可从不同微波频段、极化状态下识别形变。基于升、降轨 Sentinel-1 数据及数字高程模型 SRTM3，采用 DInSAR技术获取 2021 年比如 MS 6.1 地震同震形变场，通过分析雷达视线向（LOS）的方向，判断震区地面的沉降或抬升，从而得出地面形变量，结果发现，主体形变区沿 NW 向展布，升、降轨沿 LOS 向最大沉降量为 2.4 cm（升轨）和 2.8 cm（降轨），最大抬升量为 0.6 cm（升轨）和 3.5 cm（降轨），且整体以沉降为主。将该结果与断层运动相结合，判定本次形变场由断层兼具走滑运动造成。哨兵 1 号（Sentinel-1）数据由欧州航空局（European Space Agency，ESA）免费提供，在此表示感谢。参考文献宫熙雯，张正加，陈启浩，等.融合 DInSAR 和 MAI 技术反演九寨沟地震三维形变场 J.测绘工程，2021，30（5）：49-57.田雷，张小涛，解孟雨，等.2021 年 3 月 19 日西藏比如 MS 6.1 地震总结 J.地震地磁观测与研究，2021，42（2）：42-57.李玉丽，李启雷.2021 年西藏比如 MS 6.1 地震序列发震构造分析 J.地震研究，2022，45（1）：54-65.刘学武，胡波，莫玉娟，等.DInSAR 技术监测玉树地震同震形变场的研究 J.测绘工程，2011，20（6）：38-39.吴中海，赵希涛，吴珍汉，等.西藏安多-错那湖地堑的第四纪地质、断裂活动及其运动学特征分析J.第四纪研究，2005，25（4）：490-502.杨红磊，彭军还.基于 DInSAR 和 MAI 技术揭示地震三维形变场 J.地球物理学进展，2014，29（6）：2580-2586.袁伏全，黄浩，孙世瑞，等.2021 年西藏比如 MS 6.1 地震震源机制解与发震构造分析 J.地震地磁观测与研究，2021，42（Z1）：186-188.张明，高涵，牛玉芬，等.2016 年门源地震震源机制 DInSAR 同震形变反演 J.地球物理学进展，2017，32（3）：1 089-1 094.张威，吴仍武，刘勇.D-InSAR 技术在地面变形监测中的应用 J.科技信息，2014，（1）：79.张艳梅，江在森，李斐，等.基于DInSAR的2001年昆仑山口西8.1级大地震形变探测J.大地测量与地球动力学，2007，27（3）：12-17.Massonnet D,Rossi M,Carmona C,et al.The displacement field of the Landers earthquake mapped by radar interferometryJ.Nature,1993,364(6 433):138-142.桑杰卓玛等：基于 DInSAR 技术的 2021 年比如 MS 6.1 地震同震形变场提取与分析地震地磁观测与研究44 卷 140Extraction and analysis of co-seismic deformation field of 2021 Biru MS 6.1 earthquake based on DInSAR technologySangjiezhuoma，Suolangnanjie，Cirenduoji and Gesangzhuoma(Earthquake Agency of Tibet Autonomous Region,Lhasa 850000,China)AbstractOn March 19,2021,an earthquake with MS 6.1 occurred in Nagqu City,Tibet Autonomous Region.Due to the few geophysical stations in Tibet and their sparse distribution,it is impossible to obtain the change information in the geophysical field before and after the earthquake.Based on Sentinel-1 satellite ascending and descending SAR data and digital elevation model SRTM3,the coseismic deformation field information of earthquake area is obtained by using differential synthetic aperture radar interferometry(DInSAR)technology,and the change of surface deformation field is judged.The results show that the earthquake was caused by normal fault and strike-slip fault movement.The maximum subsidence of the ascending and descending tracks along the radar line of sight are 2.4 cm(ascending track)and 2.8 cm(descending track)respectively,and the maximum uplift are 0.6 cm(ascending track)and 3.5 cm(descending track)respectively.Keywords:earthquake，DInSAR，earthquake coseismic deformation field，Sentinel-1 satellite