一种地震同震地表破裂分段信息识别方法和系统

技术领域

本发明属于地震勘探技领域，尤其涉及一种地震同震地表破裂分段信息识别方法和系统。

背景技术

同震地表破裂是地震发生时发震断层滑动面错动地表，并与之形成的交线。复杂构造背景下，地震的孕育、发生机制复杂，因此引起地震的发震断层往往是分段的，且每一段的几何特征、滑动情况、空间分布均有不同。合成孔径雷达(Synthetic Aperture Radar，SAR)能够从卫星对地观测的角度为震源参数反演提供数据支撑。在利用InSAR数据作为约束条件进行震源参数反演时，首先，需要根据参考资料(地表破裂迹线、余震空间分布、形变场的梯度或不连续性等)，预设断层的组合模式，即单段断层模式或多段组合模式，建立初始的同震破裂几何模型，而后，再以SAR观测数据或全球导航卫星系统(Global NavigationSatellite System，GNSS)测量数据作为约束，进行同震破裂的滑动分布反演，进而获取地震的震源参数。

合理预设断层组合方案是关系震源参数模拟准确度的重要问题，该问题仍没有得到很好的解决。一个重要的限制因素就是缺乏宏观尺度上高效率、低成本、高可靠度的断层迹线(断层面与地表面交线)分段、走向信息提取技术。这严重制约了对震源发生、发展过程的认识与研究水平。

以2021年5月22日我国青海玛多7.4级地震为例，当前的研究现状是：玛多地震的发震断层具有复杂的几何形态，沿主破裂面的两端均表现出分支破裂特征【文献(Wang,W.,Fang,L.,Wu,J.et al.Aftershock sequence relocation of the 2021MS7.4 MaduoEarthquake,Qinghai,China[J],Science China.Earth Sciences,2021,64(8),1371-1380)；文献(He,K.,Wen,Y.,Xu,C.et al.Fault geometry and slip distribution ofthe 2021Mw 7.4Maduo,China,earthquake inferred from InSAR measurements andrelocated aftershocks[J],Seismological Research Letters,2022,93(1),8-20)文献(王未来,房立华,吴建平等.2021年青海玛多MS7.4地震序列精定位研究[J].中国科学:地球科学,2021,51(07):1193-1202)】，并且断层深部比浅部的几何形态与运动学特征更为复杂【文献(杨君妍,孙文科,洪顺英等.2021年青海玛多7.4级地震的同震变形分析[J],地球物理学报,2021,64(8),2671-2683)】。以往研究采用的单一断层模型难以精准地揭示九寨沟地震复杂的断层几何形状和分段破裂过程。国内外针对玛多地震发震断层地表破裂的研究众多，但对玛多地震发震构造的认识尚不统一，尤其是对主断层如何分段的问题上，国内外不同研究组的观点各异：文献(Hong,S.,Liu,M.,Liu,T.et al.Fault Source Model andStress Changes of the 2021Mw 7.4Maduo Earthquake,China,Constrained by InSARand GPS Measurements[J],Bulletin of the Seismological Society of America,2022.)认为玛多地震的发震断层不是简单的单断层破裂模式，而是由一条长约156km的西北走向断层和一条长约24km的东南走向断层组成；文献(Wang,S.,Song,C.,Li,S.etal.Resolving co-and early post-seismic slip variations of the 2021MW 7.4Maduoearthquake in east Bayan Har block with a block-wide distributed deformationmode from satellite synthetic aperture radar data[J],Earth and PlanetaryPhysics,2022,6,106-120)及文献(Zhang,Q.,Wu,Y.,Guo,N.et al.Research ondeformation characteristics of the 2021Qinghai Maduo MS7.4 earthquake throughcoseismic dislocation inversion[J],Advances in Space Research,2022.)认为玛多地震的发震断层至少可分为2段；文献(Tong,X.,Xu,X.,Chen,S.Coseismic Slip Model ofthe 2021Maduo Earthquake,China from Sentinel-1InSAR Observation[J],RemoteSensing,2022,14(3),436)和文献(Chen,K.,Avouac,J.P.,Geng,J.,et al.The 2021Mw7.4Madoi earthquake:an archetype bilateral slip-pulse rupture arrested at asplay fault[J],Geophysical Research Letters,2022,49,e2021GL095243)认为发震断层可分为4段，而文献(Zhao,D.,Qu,C.,Chen,H.,et al.Tectonic and geometric controlon fault kinematics of the2021Mw7.3 Maduo(China)earthquake inferredfrominterseismic,coseismic,and postseismic InSAR observations[J],GeophysicalResearch Letters,2021,48(18),e2021GL095417)以及文献(Chen,H.,Qu,C.,Zhao,D.,etal.Rupture Kinematics and Coseismic Slip Model of the 2021Mw 7.3Maduo(China)Earthquake:Implications for the Seismic Hazard of the Kunlun Fault[J],RemoteSensing,2021,13(16),3327)则采用5段式分段模拟方案进行了震源参数反演；进一步地，文献(He,L.,Feng,G.,Wu,X.,et al.Coseismic and Early Postseismic Slip Models ofthe 2021Mw 7.4Maduo Earthquake(Western China)Estimated by Space-BasedGeodetic Data[J],Geophysical Research Letters,2021,48(24),e2021GL095860)、文献(Liu,C.,Bai,L.,Hong,S.,et al.Coseismic deformation of the 2021Mw7.4 Maduoearthquake from joint inversion of InSAR,GPS,and teleseismic data[J],Earthquake Science,2021,34)以及文献(祝爱玉,王永哲,李永华等.基于InSAR地表形变约束的玛多Ms7.4地震孕育发生机理数值模拟研究[J],地球物理学报,2021,64(12),4548-4561)在反演震源参数时将主断层分为6段；文献(Jin,Z.,and Fialko,Y.Coseismic andEarly Postseismic Deformation Due to the 2021M7.4 Maduo(China)Earthquake[J],Geophysical Research Letters,2021,48(21),e2021GL095213)则采用了7段式分段方案，而文献(Xu,L.,Chen,Q.,Zhao,J.,et al.An Integrated Approach for Mapping Three-Dimensional CoSeismic Displacement Fields from Sentinel-1TOPSData Based onDInSAR,POT,MAI and BOI Techniques:Application to the 2021Mw 7.4MaduoEarthquake[J],Remote Sensing,2021,13(23),4847)则认为8段式的分段方案能得到更为精确的玛多地震震源参数。此外，实地勘察方面，文献(潘家伟,白明坤,李超等.2021年5月22日青海玛多Ms7.4地震地表破裂带及发震构造[J],地质学报,2021,95(6),1655-1670)指出，经过现场实地勘察，发现玛多地震的地表破裂至少存在8段分支；文献(盖海龙,姚生海,杨丽萍等.青海玛多“5·22”MS7.4级地震的同震地表破裂特征、成因及意义[J],地质力学学报,2021,27(6),899-912)通过地表调查、高分辨率光学图像解译、无人机低空摄影测量手段结合，得到了6段典型同震破裂资料。GNSS观测与震源模拟方面，文献(王迪晋,王东振,赵斌等.2021年青海玛多Mw 7.4地震GNSS同震形变场及其断层滑动分布[J],地球物理学报,2022,65(2),537-551)以及文献(Wang,M.,Wang,F.,Jiang,X.et al.GPS determinedcoseismic slip of the2021Mw7.4 Maduo,China,earthquake and its tectonicimplication[J],Geophysical Journal International,2022,228(3),2048-2055)利用GNSS数据集约束同震滑动模型，基于4段式断层破裂模型，得到了玛多地震发震断层的滑动特征。

虽然SAR的振幅、相位信息已被广泛应用于量测、解译和理解地震形变，但是，在利用InSAR观测资料作为约束的震源参数分段模拟研究中，存在一个制约震源参数反演精度的痛点问题，即难以确定复杂断层的分段方案(如断层应该分几段？每一段长度、走向如何？等)，这给震源参数模拟研究带来了较大的困难。虽然无人机航拍摄影技术能够提供较为精细的地表破裂信息，但其适用场景较小；人员现场勘察也能够提供较为精准、可靠的地表破裂资料，但该方法往往具有滞后性且经济成本和危险系数较高。星载SAR大地测量技术虽然可以全天时、全天候、大面积、低成本地对地震灾害做出快速响应，但目前其主要的应用领域是地表形变测量。

准确获知断层的几何特征和运动学参数对于深入理解构造运动和孕震过程都是一件意义重大的基础性工作，也是一直以来难以解决的问题。目前国内外研究者在复杂构造地震震源参数反演研究中，普遍存在以下共识：在震源参数模拟时，特别是在复杂构造背景下，地震发震断层通常不是单段、均一的。每一段的几何形态与滑动特征各异，而常规技术手段难以客观、高效、低成本地为制定断层分段方案提供决策依据，为震源参数反演研究带来了困难，直接影响了震源机制解的精度，成为严重制约卫星大地测量震源机制反演研究与应用水平的“卡脖子”技术瓶颈。

发明内容

为解决上述技术问题，本发明提出一种地震同震地表破裂分段信息识别方法的技术方案，以解决上述技术问题。

本发明第一方面公开了一种地震同震地表破裂分段信息识别方法，所述方法包括：

步骤S1、选取一幅震后SAR图像和两幅震前SAR图像，并将两幅震前SAR图像组成震前干涉对，将一幅震前SAR图像和一幅震后SAR图像组成同震干涉对；获取与震后SAR图像和震前SAR图像拍摄日期在预定义时间范围内的高分辨率图像；获取不小于震后SAR图像和震前SAR图像覆盖范围的数字高程模型数据；

步骤S2、应用所述震前干涉对和同震干涉对，得到震前干涉图和同震干涉图；对所述震前干涉图和同震干涉图分别进行距离向和方位向多视处理后，再去除平地相位及地形相位，得到震前差分干涉图和同震差分干涉图；应用滤波后的所述震前差分干涉图和同震差分干涉图，计算震前局部相位方差和同震局部相位方差；应用所述震前局部相位方差和同震局部相位方差，计算震前干涉相关系数和同震干涉相关系数，并进行地理编码；

步骤S3、通过地理编码后的震前干涉相关系数和同震干涉相关系数，得到升轨干涉相关系数差值图和降轨干涉相关系数差值图；基于所述升轨干涉相关系数差值图与降轨干涉相关系数差值图、所述高分辨率图像和所述数字高程模型数据，通过目视解译，提取断层的段数、各段长度、走向以及空间位置。

根据本发明第一方面的方法，在所述步骤S2中，所述应用所述震前干涉对和同震干涉对，得到震前干涉图和同震干涉图的方法包括：

将两幅震前SAR图像配准到相同坐标范围内；将震后SAR图像配准到震前SAR图像的坐标范围内；

对配准后的成震前干涉对和同震干涉对进行复数共轭相乘运算，得到震前干涉图和同震干涉图。

根据本发明第一方面的方法，在所述步骤S2中，所述去除平地相位及地形相位的方法包括：

采用精密轨道数据，模拟得到因地球曲率导致的平地相位；

采用所述数字高程模型数据，模拟得到因地形导致的地形相位；

从距离向和方位向多视处理后的震前干涉图和同震干涉图中去除所述平地相位及地形相位。

根据本发明第一方面的方法，在所述步骤S2中，所述应用滤波后的所述震前差分干涉图和同震差分干涉图，计算震前局部相位方差和同震局部相位方差的方法包括：

对滤波后的所述震前差分干涉图和同震差分干涉图，采用5×5像素的窗口，对相位概率分布函数进行积分运算，计算所述窗口内滤波后的所述震前差分干涉图和同震差分干涉图的局部相位方差，得到震前局部相位方差和同震局部相位方差。

根据本发明第一方面的方法，在所述步骤S2中，所述应用所述震前局部相位方差和同震局部相位方差，计算震前干涉相关系数和同震干涉相关系数的方法包括：

其中，γ

根据本发明第一方面的方法，在所述步骤S2中，所述进行地理编码的方法包括：

基于世界大地测量系统-84地理坐标系对前干涉相关系数和同震干涉相关系数进行地理编码，得到地理编码后的震前干涉相关系数和地理编码后的同震干涉相关系数。

根据本发明第一方面的方法，在所述步骤S3中，所述基于所述升轨干涉相关系数差值图与降轨干涉相关系数差值图、所述高分辨率图像和所述数字高程模型数据，通过目视解译，提取断层的段数、各段长度、走向以及空间位置的方法包括：

以完全覆盖研究区为标准，将所述升轨干涉相关系数差值图与降轨干涉相关系数差值图裁剪至相同地理范围，得到裁剪后的升轨干涉相关系数差值图和裁剪后的降轨干涉相关系数差值图；

将所述高分辨率图像和所述数字高程模型数据分别裁剪至与裁剪后的升轨干涉相关系数差值图和裁剪后的降轨干涉相关系数差值图相同的地理范围；

以裁剪后的升轨干涉相关系数差值图和裁剪后的降轨干涉相关系数差值图中的负值作为指标，分别提取震区范围内的直线型负值轨迹，作为初步的地表破裂迹线，形成初步的地表破裂迹线草图，记作升轨迹线和降轨迹线；

提取升轨迹线和降轨迹线中震区范围内直线型负值轨迹的并集，结合裁剪后的高分辨率图像和数字高程模型数据，通过目视解译，去除由于山脊岩土松动滑落等原因导致的直线型负值干涉相关系数迹线，最终得到地震发震断层的地表破裂迹线图，进而获知段数、各段长度、走向以及空间位置。

本发明第二方面公开了一种地震同震地表破裂分段信息识别系统，所述系统包括：

第一处理模块，被配置为，选取一幅震后SAR图像和两幅震前SAR图像，并将两幅震前SAR图像组成震前干涉对，将一幅震前SAR图像和一幅震后SAR图像组成同震干涉对；获取与震后SAR图像和震前SAR图像拍摄日期在预定义时间范围内的高分辨率图像；获取不小于震后SAR图像和震前SAR图像覆盖范围的数字高程模型数据；

第二处理模块，被配置为，应用所述震前干涉对和同震干涉对，得到震前干涉图和同震干涉图；对所述震前干涉图和同震干涉图分别进行距离向和方位向多视处理后，再去除平地相位及地形相位，得到震前差分干涉图和同震差分干涉图；应用滤波后的所述震前差分干涉图和同震差分干涉图，计算震前局部相位方差和同震局部相位方差；应用所述震前局部相位方差和同震局部相位方差，计算震前干涉相关系数和同震干涉相关系数，并进行地理编码；

第三处理模块，被配置为，通过地理编码后的震前干涉相关系数和同震干涉相关系数，得到升轨干涉相关系数差值图和降轨干涉相关系数差值图；基于所述升轨干涉相关系数差值图与降轨干涉相关系数差值图、所述高分辨率图像和所述数字高程模型数据，通过目视解译，提取断层的段数、各段长度、走向以及空间位置。

本发明第三方面公开了一种电子设备。电子设备包括存储器和处理器，存储器存储有计算机程序，处理器执行计算机程序时，实现本公开第一方面中任一项的一种地震同震地表破裂分段信息识别方法中的步骤。

本发明第四方面公开了一种计算机可读存储介质。计算机可读存储介质上存储有计算机程序，计算机程序被处理器执行时，实现本公开第一方面中任一项的一种地震同震地表破裂分段信息识别方法中的步骤。

综上，本发明提出的方案能够客观、快速、低成本地改进当前研究时效性较低、主观性较强、人力财力成本较大的问题，进而有益于提高震源参数反演的精度，为认识、理解复杂构造条件下地震的孕育发生机理提供可靠的大地测量结果与模型基础。

附图说明

为了更清楚地说明本发明具体实施方式或现有技术中的技术方案，下面将对具体实施方式或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图是本发明的一些实施方式，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1为根据本发明实施例的一种地震同震地表破裂分段信息识别方法的流程图；

图2为根据本发明实施例的Sentinel-1卫星SAR数据升轨、降轨数据处理得到的2021年青海玛多7.4级地震干涉相关系数图；

图3为根据本发明实施例的Sentinel-1数据处理得到的玛多地震升轨、降轨干涉相关系数差分结果及高分辨率图像和数字高程模型数据；

图4为根据本发明实施例的2021年5月22日玛多7.4级地震发震断层的地表破裂迹线图；

图5为根据本发明实施例的一种地震同震地表破裂分段信息识别系统的结构图；

图6为根据本发明实施例的一种电子设备的结构图。

具体实施方式

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例只是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

本发明第一方面公开了一种地震同震地表破裂分段信息识别方法。图1为根据本发明实施例的一种地震同震地表破裂分段信息识别方法的流程图，如图1所示，所述方法包括：

步骤S1、选取一幅震后SAR图像和两幅震前SAR图像，并将两幅震前SAR图像组成震前干涉对，将一幅震前SAR图像和一幅震后SAR图像组成同震干涉对；获取与震后SAR图像和震前SAR图像拍摄日期在预定义时间范围内的高分辨率图像；获取不小于震后SAR图像和震前SAR图像覆盖范围的数字高程模型数据；

步骤S2、应用所述震前干涉对和同震干涉对，得到震前干涉图和同震干涉图；对所述震前干涉图和同震干涉图分别进行距离向和方位向多视处理后，再去除平地相位及地形相位，得到震前差分干涉图和同震差分干涉图；应用滤波后的所述震前差分干涉图和同震差分干涉图，计算震前局部相位方差和同震局部相位方差；应用所述震前局部相位方差和同震局部相位方差，计算震前干涉相关系数和同震干涉相关系数，并进行地理编码；

步骤S3、通过地理编码后的震前干涉相关系数和同震干涉相关系数，得到升轨干涉相关系数差值图和降轨干涉相关系数差值图；基于所述升轨干涉相关系数差值图与降轨干涉相关系数差值图、所述高分辨率图像和所述数字高程模型数据，通过目视解译，提取断层的段数、各段长度、走向以及空间位置。

在步骤S1，选取一幅震后SAR图像和两幅震前SAR图像，并将两幅震前SAR图像组成震前干涉对，将一幅震前SAR图像和一幅震后SAR图像组成同震干涉对；获取与震后SAR图像和震前SAR图像拍摄日期在预定义时间范围内的高分辨率图像；获取不小于震后SAR图像和震前SAR图像覆盖范围的数字高程模型数据。

具体地，SAR图像选取是进一步处理得到干涉相关系数的数据基础，由欧洲空间局ESA公开数据库，检索下载与SAR数据同期的轨道文件；高分辨率图像用于目视解译排除由于降雪、人类活动等非断层活动原因导致的干涉相关系数的下降。

在步骤S2，选取一幅震后SAR图像和两幅震前SAR图像，并将两幅震前SAR图像组成震前干涉对，将一幅震前SAR图像和一幅震后SAR图像组成同震干涉对；获取与震后SAR图像和震前SAR图像拍摄日期在预定义时间范围内的高分辨率图像；获取不小于震后SAR图像和震前SAR图像覆盖范围的数字高程模型数据。

在一些实施例中，在所述步骤S2中，所述应用所述震前干涉对和同震干涉对，得到震前干涉图和同震干涉图的方法包括：

将两幅震前SAR图像配准到相同坐标范围内；将震后SAR图像配准到震前SAR图像的坐标范围内；

对配准后的成震前干涉对和同震干涉对进行复数共轭相乘运算，得到震前干涉图和同震干涉图。

所述去除平地相位及地形相位的方法包括：

采用欧洲空间局ESA公开发布的精密轨道数据，模拟得到因地球曲率导致的平地相位；

采用美国地质调查局USGS公开发布的所述数字高程模型数据，模拟得到因地形导致的地形相位；

从距离向和方位向多视处理后的震前干涉图和同震干涉图中去除所述平地相位及地形相位。

所述应用滤波后的所述震前差分干涉图和同震差分干涉图，计算震前局部相位方差和同震局部相位方差的方法包括：

对滤波后的所述震前差分干涉图和同震差分干涉图，采用5×5像素的窗口，对相位概率分布函数进行积分运算，计算所述窗口内滤波后的所述震前差分干涉图和同震差分干涉图的局部相位方差，得到震前局部相位方差和同震局部相位方差。

所述应用所述震前局部相位方差和同震局部相位方差，计算震前干涉相关系数和同震干涉相关系数的方法包括：

其中，γ

所述进行地理编码的方法包括：

基于世界大地测量系统-84地理坐标系对前干涉相关系数和同震干涉相关系数进行地理编码，得到地理编码后的震前干涉相关系数和地理编码后的同震干涉相关系数。

在步骤S3，通过地理编码后的震前干涉相关系数和同震干涉相关系数，得到升轨干涉相关系数差值图和降轨干涉相关系数差值图；基于所述升轨干涉相关系数差值图与降轨干涉相关系数差值图、所述高分辨率图像和所述数字高程模型数据，通过目视解译，提取断层的段数、各段长度、走向以及空间位置。

在一些实施例中，在所述步骤S3中，所述基于所述升轨干涉相关系数差值图与降轨干涉相关系数差值图、所述高分辨率图像和所述数字高程模型数据，通过目视解译，提取断层的段数、各段长度、走向以及空间位置的方法包括：

以完全覆盖研究区为标准，将所述升轨干涉相关系数差值图与降轨干涉相关系数差值图裁剪至相同地理范围，得到裁剪后的升轨干涉相关系数差值图和裁剪后的降轨干涉相关系数差值图；

将所述高分辨率图像和所述数字高程模型数据分别裁剪至与裁剪后的升轨干涉相关系数差值图和裁剪后的降轨干涉相关系数差值图相同的地理范围；

以裁剪后的升轨干涉相关系数差值图和裁剪后的降轨干涉相关系数差值图中的负值作为指标，分别提取震区范围内的直线型负值轨迹，作为初步的地表破裂迹线，形成初步的地表破裂迹线草图，记作升轨迹线和降轨迹线；

提取升轨迹线和降轨迹线中震区范围内直线型负值轨迹的并集，结合裁剪后的高分辨率图像和数字高程模型数据，通过目视解译，去除由于山脊岩土松动滑落等原因导致的直线型负值干涉相关系数迹线，最终得到地震发震断层的地表破裂迹线图，进而获知段数、各段长度、走向以及空间位置。

实施例

1.数据选取与干涉组合构建。

1.1数据选取

1.1.1SAR图像选取

SAR数据为Sentinel-1数据，数据格式为单视复数。升轨(Ascending)SAR数据选择2021年5月8日，2021年5月20日、2021年5月26日，分别记为A20210508、A20210520、A20210526。降轨(Descending)SAR数据选择2021年5月8日，2021年5月20日、2021年5月26日，分别记为D20210508、D20210520、D20210526。

1.1.2高分辨率图像选取

选取2020年12月31日Google Earth光学影像。

1.1.3数字高程模型数据DEM获取

美国USGS公开发布的SRTM 30m分辨率DEM下载，数据格式为Net-CDF 4网格。DEM的范围需完全覆盖SAR图像的范围。

1.1.4SAR数据轨道文件获取

欧空局ESA公开发布的Sentinel-1精密星历轨道文件下载，数据后缀名为.EOF。要求与SAR数据的成像时刻同期。

1.2干涉组合构建

分别构建升轨干涉对、降轨干涉对。具体地，升轨震前干涉对为：A20210508-A20210520；升轨同震干涉对为：A20210520-A20210526；降轨震前干涉对为；D20210508-D20210520；降轨同震干涉对为D20210520-D20210526。

2.干涉相关系数计算

2.1、配准。

震前干涉对配准时，将A20210520配准到A20210508坐标范围内；同震干涉对配准时，将A20210526配准到A20210520坐标范围内。

2.2、生成干涉图。

对步骤2.1所得配准后的震前干涉对A20210508-A20210520、配准后的同震干涉对A20210520-A20210526进行复数共轭相乘运算。

2.3、多视。

对步骤2.2所得震前干涉图、同震干涉图分别进行距离向、方位向多视处理多视比例为8(距离向)：2(方位向)。

2.4、去除平地相位及地形相位。

对步骤2.3所得多视后的震前干涉图、多视后的同震干涉图，采用欧洲空间局ESA公开发布的精密轨道数据、采用美国地质调查局USGS公开发布的DEM，分别去除平地相位、地形相位。

2.5、滤波。

采用Goldstein方法对步骤2.4所得震前差分干涉图和同震差分干涉图进行滤波，得到滤波后的震前差分干涉图和同震差分干涉图。

2.6、局部相位方差计算。

对步骤2.5所得的震前差分干涉图和同震差分干涉图采用5×5像素的窗口，对PDF进行积分运算，计算该窗口内干涉图的局部相位方差σ

2.7、干涉相关系数计算。

基于步骤2.6所得的局部相位方差σ

2.8、地理编码。

对步骤2.7所得的震前干涉相关系数和同震干涉相关系数，基于WGS-84地理坐标系进行地理编码，得到地理编码后的升轨震前干涉相关系数如图2中的a所示、降轨震前干涉相关系数，如图2中的c所示，和地理编码后的升轨同震干涉相关系数如图2中的b所示、降轨同震干涉相关系数如图2中的d所示。

3.地表破裂迹线特征参量特征提取

3.1、用步骤2.8所得地理编码后的升轨同震干涉相关系数如图2中的b所示，减去地理编码后的升轨震前干涉相关系数如图2中的a所示，得到升轨干涉相关系数差值图；用地理编码后的降轨同震干涉相关系数。如图图2中的d所示，减去地理编码后的降轨震前干涉相关系数，如图2中的c所示，得到降轨干涉相关系数差值图。

3.2、基于步骤3.1所得的干涉相关系数差值图和步骤1.1.2所得的高分辨率图像，如图3中的c所示，通过目视解译，提取断层的段数、各段长度、走向以及空间位置。

3.2.1以完全覆盖研究区为标准，将步骤3.1所得的升轨干涉相关系数差值图和降轨干涉相关系数差值图裁剪至相同地理范围，得到裁剪后的升轨干涉相关系数差值图，如图3中的b所示，和裁剪后的降轨干涉相关系数差值图，如图3中的d所示。

3.2.2将步骤1.1.2所得高分辨率图像、步骤2.4所用DEM分别裁剪至与图3中的b和图3中的d相同的范围。得到裁剪后的DEM，如图3中的a所示，和裁剪后的高分辨率图像，如图3中的c所示。

3.2.3以图3中的b和图3中的d中的负值作为指标，分别提取震区范围内的直线型负值轨迹，作为初步的地表破裂迹线，形成初步的地表破裂迹线草图。

3.2.4提取步骤3.2.3所得地表破裂迹线草图中震区范围内直线型负值轨迹的并集，结合步骤3.2.2所得裁剪后的高分辨率图像、DEM数据，通过目视解译，去除由于山脊岩土松动滑落等原因导致的直线型负值干涉相关系数迹线，最终得到地震发震断层的地表破裂迹线图，如图4所示，进而获知段数、各段长度、走向以及空间位置。

综上，本发明提出的方案能够将InSAR技术拓展应用到复杂构造背景下地表破裂精准识别与几何特征精细刻画上，是传统InSAR技术在地震震源参数反演领域的应用创新。由实际应用需求作牵引，瞄准目前基于InSAR资料的震源机制反演研究中存在的共性需求和痛点问题，即：难以确定发震断层分几段、在何处分段以及各段断层走向的问题，有针对性地予以解决，进而实现高精度、精细震源参数反演。从充分挖掘利用SAR相位干涉信息这一角度入手，突破了传统技术在断层地表破裂信息识别应用研究领域普遍存在的效率低、成本高、主观性强、缺乏参考资料等局限性，提出了一种经济、高效、客观、直观的解决方案。

本发明第二方面公开了一种地震同震地表破裂分段信息识别系统。图5为根据本发明实施例的一种地震同震地表破裂分段信息识别系统的结构图；如图5所示，所述系统100包括：

第一处理模块101，被配置为，选取一幅震后SAR图像和两幅震前SAR图像，并将两幅震前SAR图像组成震前干涉对，将一幅震前SAR图像和一幅震后SAR图像组成同震干涉对；获取与震后SAR图像和震前SAR图像拍摄日期在预定义时间范围内的高分辨率图像；获取不小于震后SAR图像和震前SAR图像覆盖范围的数字高程模型数据；

第二处理模块102，被配置为，应用所述震前干涉对和同震干涉对，得到震前干涉图和同震干涉图；对所述震前干涉图和同震干涉图分别进行距离向和方位向多视处理后，再去除平地相位及地形相位，得到震前差分干涉图和同震差分干涉图；应用滤波后的所述震前差分干涉图和同震差分干涉图，计算震前局部相位方差和同震局部相位方差；应用所述震前局部相位方差和同震局部相位方差，计算震前干涉相关系数和同震干涉相关系数，并进行地理编码；

第三处理模块103，被配置为，通过地理编码后的震前干涉相关系数和同震干涉相关系数，得到升轨干涉相关系数差值图和降轨干涉相关系数差值图；基于所述升轨干涉相关系数差值图与降轨干涉相关系数差值图、所述高分辨率图像和所述数字高程模型数据，通过目视解译，提取断层的段数、各段长度、走向以及空间位置。

根据本发明第二方面的系统，所述第一处理模块101具体被配置为，SAR图像选取是进一步处理得到干涉相关系数的数据基础，由欧洲空间局ESA公开数据库，检索下载与SAR数据同期的轨道文件；高分辨率图像用于目视解译排除由于降雪、人类活动等非断层活动原因导致的干涉相关系数的下降。

根据本发明第二方面的系统，所述第二处理模块102具体被配置为，所述应用所述震前干涉对和同震干涉对，得到震前干涉图和同震干涉图的方法包括：

将两幅震前SAR图像配准到相同坐标范围内；将震后SAR图像配准到震前SAR图像的坐标范围内；

对配准后的成震前干涉对和同震干涉对进行复数共轭相乘运算，得到震前干涉图和同震干涉图。

所述去除平地相位及地形相位的方法包括：

采用欧洲空间局ESA公开发布的精密轨道数据，模拟得到因地球曲率导致的平地相位；

采用美国地质调查局USGS公开发布的所述数字高程模型数据，模拟得到因地形导致的地形相位；

从距离向和方位向多视处理后的震前干涉图和同震干涉图中去除所述平地相位及地形相位。

所述应用滤波后的所述震前差分干涉图和同震差分干涉图，计算震前局部相位方差和同震局部相位方差的方法包括：

对滤波后的所述震前差分干涉图和同震差分干涉图，采用5×5像素的窗口，对相位概率分布函数进行积分运算，计算所述窗口内滤波后的所述震前差分干涉图和同震差分干涉图的局部相位方差，得到震前局部相位方差和同震局部相位方差。

所述应用所述震前局部相位方差和同震局部相位方差，计算震前干涉相关系数和同震干涉相关系数的方法包括：

其中，γ

所述进行地理编码的方法包括：

基于世界大地测量系统-84地理坐标系对前干涉相关系数和同震干涉相关系数进行地理编码，得到地理编码后的震前干涉相关系数和地理编码后的同震干涉相关系数。

根据本发明第二方面的系统，所述第三处理模块103具体被配置为，所述基于所述升轨干涉相关系数差值图与降轨干涉相关系数差值图、所述高分辨率图像和所述数字高程模型数据，通过目视解译，提取断层的段数、各段长度、走向以及空间位置的方法包括：

以完全覆盖研究区为标准，将所述升轨干涉相关系数差值图与降轨干涉相关系数差值图裁剪至相同地理范围，得到裁剪后的升轨干涉相关系数差值图和裁剪后的降轨干涉相关系数差值图；

将所述高分辨率图像和所述数字高程模型数据分别裁剪至与裁剪后的升轨干涉相关系数差值图和裁剪后的降轨干涉相关系数差值图相同的地理范围；

以裁剪后的升轨干涉相关系数差值图和裁剪后的降轨干涉相关系数差值图中的负值作为指标，分别提取震区范围内的直线型负值轨迹，作为初步的地表破裂迹线，形成初步的地表破裂迹线草图，记作升轨迹线和降轨迹线；

提取升轨迹线和降轨迹线中震区范围内直线型负值轨迹的并集，结合裁剪后的高分辨率图像和数字高程模型数据，通过目视解译，去除由于山脊岩土松动滑落等原因导致的直线型负值干涉相关系数迹线，最终得到地震发震断层的地表破裂迹线图，进而获知段数、各段长度、走向以及空间位置。

本发明第三方面公开了一种电子设备。电子设备包括存储器和处理器，存储器存储有计算机程序，处理器执行计算机程序时，实现本发明公开第一方面中任一项的一种地震同震地表破裂分段信息识别方法中的步骤。

图6为根据本发明实施例的一种电子设备的结构图，如图6所示，电子设备包括通过系统总线连接的处理器、存储器、通信接口、显示屏和输入装置。其中，该电子设备的处理器用于提供计算和控制能力。该电子设备的存储器包括非易失性存储介质、内存储器。该非易失性存储介质存储有操作系统和计算机程序。该内存储器为非易失性存储介质中的操作系统和计算机程序的运行提供环境。该电子设备的通信接口用于与外部的终端进行有线或无线方式的通信，无线方式可通过WIFI、运营商网络、近场通信(NFC)或其他技术实现。该电子设备的显示屏可以是液晶显示屏或者电子墨水显示屏，该电子设备的输入装置可以是显示屏上覆盖的触摸层，也可以是电子设备外壳上设置的按键、轨迹球或触控板，还可以是外接的键盘、触控板或鼠标等。

本领域技术人员可以理解，图6中示出的结构，仅仅是与本公开的技术方案相关的部分的结构图，并不构成对本申请方案所应用于其上的电子设备的限定，具体的电子设备可以包括比图中所示更多或更少的部件，或者组合某些部件，或者具有不同的部件布置。

本发明第四方面公开了一种计算机可读存储介质。计算机可读存储介质上存储有计算机程序，计算机程序被处理器执行时，实现本发明公开第一方面中任一项的一种地震同震地表破裂分段信息识别方法中的步骤。

请注意，以上实施例的各技术特征可以进行任意的组合，为使描述简洁，未对上述实施例中的各个技术特征所有可能的组合都进行描述，然而，只要这些技术特征的组合不存在矛盾，都应当认为是本说明书记载的范围。以上实施例仅表达了本申请的几种实施方式，其描述较为具体和详细，但并不能因此而理解为对发明专利范围的限制。应当指出的是，对于本领域的普通技术人员来说，在不脱离本申请构思的前提下，还可以做出若干变形和改进，这些都属于本申请的保护范围。因此，本申请专利的保护范围应以所附权利要求为准。