一种基于初至时差和波形叠加的微地震震源定位方法，步骤如下：S1：输入速度模型；S2：拾取并输入实际初至时间，以及读取地震数据；S3：依据S1中速度模型，计算可行解区域内所有网格点到每一个检波器的理论初至时间表；S4：构造走时残差函数Tr；S5：构造波形叠加函数Ews；S6：输入权重系数β，依据S4中走时残差函数Tr和S5中波形叠加函数Ews，构造改进目标函数；S7：通过网格搜索法寻找改进目标函数的最小值，对应的最优解即震源位置。本发明针对基于走时目标函数的定位方法对初至误差敏感的问题，结合初至时差和波形叠加信息构造改进的目标函数，可增强定位方法的抗噪能力，可提高反演方法的收敛性，从而可改善微地震震源定位的精度。
1.一种基于初至时差和波形叠加的微地震震源定位方法，其特征在于，包括以下步骤：S1：输入速度模型；
S3：依据步骤S1中速度模型，计算可行解区域内所有网格点到每一个检波器的理论初至时间表；
所述可行解区域内所有网格点是指以射孔点为中心，给定一个求解区域，然后对该区域进行网格化，每个网格点均为可能的震源点；
S4：基于步骤S2中输入的实际初至时间和步骤S3中计算的理论初至时间表，依据公式(1)和(2)构造走时残差函数Tr；
上式中：M、N分别是观测到的横波和纵波初至的个数； 和 分别是第j个检波器记录的横波和纵波初至时间， 和 分别是与之对应的第j个检波器的横波及纵波理论计算时间； 和 分别是第i个检波器记录的横波和纵波初至时间， 和 分别是与之对应的第i个检波器的横波及纵波理论计算时间；Tshift是观测及理论初至之间的一个恒定的漂移量，以解决实际资料监测中微地震震源发震时间未知的难题；
S5：将步骤S2中读入的地震数据、步骤S4中的走时残差函数对应微地震事件中第m个检波器记录的可识别纵波初至时间 和步骤S3中理论初至时间表结合，通过极性修正，然后依据公式(3)、(4)和(5)构造波形叠加函数Ews；
式中ΔP、ΔS分别是纵、横波滑动时窗样点数；dt是采样间隔；n是检波器个数；
分别是假设震源点到第i个检波器的纵、横波理论时间 与第m个检波器的纵波理论时间 之差， 分别是第i个检波器纵、横波振幅的符号系数，ui是地震数据第i道振幅值；
S6：输入权重系数β，依据步骤S4中走时残差函数Tr和步骤S5中波形叠加函数Ews，按照公式(6)构造改进目标函数；
S7：通过网格搜索法寻找改进目标函数的最小值，输出此时对应的最优解，即震源位置。
2.根据权利要求1所述的一种基于初至时差和波形叠加的微地震震源定位方法，其特征在于：在所述步骤S1中，速度模型可通过声波测井资料和地质分层资料建立初始速度模型，然后利用射孔资料进行速度模型校正。
3.根据权利要求1所述的一种基于初至时差和波形叠加的微地震震源定位方法，其特征在于：在所述步骤S5中，极性修正可通过符号系数 实现，即判断检波器记录中时窗内振幅符号，需要正负号保持一致，此处保证波形的同向叠加。
4.根据权利要求1所述的一种基于初至时差和波形叠加的微地震震源定位方法，其特征在于：在所述步骤S6中，权重系数β>0，大小与微地震资料的信噪比情况成反比。
5.根据权利要求1所述的一种基于初至时差和波形叠加的微地震震源定位方法，其特征在于：在所述步骤S7中，在网格搜索定位中，假设网格点是震源位置时，Tr为最小，Ews为最大，即-βEws为最小，故改进目标函数Q为最小值。
一种基于初至时差和波形叠加的微地震震源定位方法
技术领域
[0001] 本发明属于石油地球物理勘探开发技术领域，具体涉及一种基于初至时差和波形叠加的微地震震源定位方法。
背景技术
[0002] 水力压裂微地震监测技术是监测致密含油气和页岩气等非常规储层改造及效果的油藏地球物理技术，其应用结果可评价压裂效果、调整压裂设计和井网布置，并对下一步开发提供有效指导，从而提高非常规油气藏的产能。
[0003] 在水力压裂监测资料处理中，其核心部分是微地震震源定位。前人基本上是利用单一的走时信息或者波形信息来建立目标函数，可分为基于走时的射线追踪定位和基于波形的偏移定位方法。而在实际水力压裂实时监测作业中，采集的微地震资料可能存在信噪比不高、能量较弱和背景噪音较大等特点，极大地影响了初至拾取的精度。基于走时信息目标函数的定位方法对初至误差更为敏感，而基于波形信息目标函数的定位方法虽然不需要精确的初至时间，但是其计算量巨大，而导致定位效率不高。为此，亟需结合走时信息和波形信息建立目标函数，研发一种基于初至时差和波形叠加的微地震震源定位方法。
发明内容
[0004] 本发明的目的在于克服上述背景技术的不足，而提供一种基于初至时差和波形叠加的微地震震源定位方法，该方法利用目标函数中微地震资料的走时和波形信息同时约束定位，并且具有一定抗噪能力，可提高反演方法的收敛性，以改善微地震震源定位的精度。
[0005] 为实现上述目的，本发明提供的一种基于初至时差和波形叠加的微地震震源定位方法，针对任意一个微地震震源(事件)定位而言，包括以下步骤：
[0006] S1：输入速度模型；
[0007] S2：拾取并输入实际初至时间，以及读取地震数据；
[0008] S3：依据步骤S1中速度模型，计算可行解区域内所有网格点到每一个检波器的理论初至时间表；
[0009] S4：基于步骤S2中输入的实际初至时间和步骤S3中计算的理论初至时间表，依据公式(1)和(2)构造走时残差函数Tr；
[0010]
[0011]
[0012] 上式中：M、N分别是观测到的横波和纵波初至的个数； 和 分别是第j个检波器记录的横波和纵波初至时间， 和 分别是与之对应的第j个检波器的横波及纵波理论计算时间； 和 分别是第i个检波器记录的横波和纵波初至时间， 和 分别是与之对应的第i个检波器的横波及纵波理论计算时间；Tshift是观测及理论初至之间的一个恒定的漂移量，以解决实际资料监测中微地震震源发震时间未知的难题；
[0013] S5：将步骤S2中读入的地震数据、走时目标函数对应微地震事件中第m个检波器记录的可识别纵波初至时间 和步骤S3中理论初至时间表结合，通过极性修正，然后依据公式(3)、(4)和(5)构造波形叠加函数Ews；
[0014]
[0015]
[0016]
[0017] 式中ΔP、ΔS分别是纵横波滑动时窗样点数；dt是采样间隔；n是检波器个数；分别是假设震源点到第i个检波器的纵横波理论时间 与第m个检波器的
纵波理论时间 之差， 分别是第i个检波器纵横波振幅的符号系数，ui是地震数据第i道振幅值；
[0018] S6：输入权重系数β，依据步骤S4中走时残差函数Tr和步骤S5中波形叠加函数Ews，按照公式(6)构造改进目标函数；
[0019] Q＝Tr-βEws (6)
[0020] S7：通过网格搜索法寻找改进目标函数的最小值，输出此时对应的最优解，即震源位置。
[0021] 上述技术方案中，在所述步骤S1中，速度模型可通过声波测井资料和地质分层资料建立初始速度模型，然后利用射孔资料进行速度模型校正。
[0022] 上述技术方案中，在所述步骤S3中，可行解区域内所有网格点是指以射孔点为中心，给定一个求解区域，然后对该区域进行网格化，每个网格点均为可能的震源点。
[0023] 上述技术方案中，在所述步骤S5中，极性修正可通过符号系数 实现，即判断检波器记录中时窗内振幅符号，需要正负号保持一致，此处保证波形的同向叠加。
[0024] 上述技术方案中，在所述步骤S6中，权重系数β>0，大小与微地震资料的信噪比情况成反比。
[0025] 上述技术方案中，在所述步骤S7中，在网格搜索定位中，假设网格点是震源位置时，Tr为最小，Ews为最大(即-βEws为最小)，故改进目标函数Q为最小值。
[0026] 与现有技术相比，本发明具有如下优点：
[0027] 本发明利用改进目标函数中微地震资料的走时和波形信息同时约束定位，增强了定位方法的抗噪能力，提高了反演方法的收敛性，可改善微地震震源定位的精度。
附图说明
[0028] 图1是本发明的实施流程示意图；
[0029] 图2是基于均匀介质模型的二维微地震井中观测系统，横、纵坐标分别表示水平距离x和深度z；
[0030] 图3a是理论模型的z分量地震记录；
[0031] 图3b是加入随机噪音后理论模型的z分量地震记录；
[0032] 图4a是基于走时目标函数方法的定位结果图，横、纵坐标分别表示水平距离x和深度z；
[0033] 图4b是本发明方法的定位结果图，横、纵坐标分别表示水平距离x和深度z。
具体实施方式
[0034] 下面结合附图和实施例对本发明的实施方式做进一步的详细说明，但它们并不构成对本发明的限定，仅作举例而已。同时通过说明本发明的优点将变得更加清楚和容易理解。
[0035] 本实施例的一种基于初至时差和波形叠加的微地震震源定位方法，如图1所示，包括以下步骤：
[0036] S1：输入速度模型；
[0037] S2：拾取并输入实际初至时间，以及读取地震数据；
[0038] S3：依据步骤S1中速度模型，计算可行解区域内所有网格点到每一个检波器的理论初至时间表；
[0039] S4：基于步骤S2中输入的实际初至时间和步骤S3中计算的理论初至时间表，依据公式(1)和(2)构造走时残差函数Tr；
[0040]
[0041]
[0042] 上式中：M、N分别是观测到的横波和纵波初至的个数； 和 分别是第j个检波器记录的横波和纵波初至时间， 和 分别是与之对应的第j个检波器的横波及纵波理论计算时间； 和 分别是第i个检波器记录的横波和纵波初至时间， 和 分别是与之对应的第i个检波器的横波及纵波理论计算时间；Tshift是观测及理论初至之间的一个恒定的漂移量，以解决实际资料监测中微地震震源发震时间未知的难题；
[0043] S5：将步骤S2中读入的地震数据、走时目标函数对应微地震事件中第m个检波器记录的可识别纵波初至时间 和步骤S3中理论初至时间表结合，通过极性修正，然后依据公式(3)、(4)和(5)构造波形叠加函数Ews；
[0044]
[0045]
[0046]
[0047] 式中ΔP、ΔS分别是纵横波滑动时窗样点数；dt是采样间隔；n是检波器个数；分别是假设震源点到第i个检波器的纵横波理论时间 与第m个检波器的
纵波理论时间 之差， 分别是第i个检波器纵横波振幅的符号系数，ui是地震数据第i道振幅值；
[0048] S6：输入权重系数β，依据步骤S4中走时残差函数Tr和步骤S5中波形叠加函数Ews，按照公式(6)构造改进目标函数；
[0049] Q＝Tr-βEws (6)
[0050] S7：通过网格搜索法寻找改进目标函数的最小值，输出此时对应的最优解，即震源位置。
[0051] 上述技术方案中，在所述步骤S1中，速度模型可通过声波测井资料和地质分层资料建立初始速度模型，然后利用射孔资料进行速度模型校正。
[0052] 上述技术方案中，在所述步骤S3中，可行解区域内所有网格点是指以射孔点为中心，给定一个求解区域，然后对该区域进行网格化，每个网格点均为可能的震源点。
[0053] 上述技术方案中，在所述步骤S5中，用大致准确的纵波初至时间 代替了常规波形叠加方法中未知的震源发震时间。
[0054] 上述技术方案中，在所述步骤S5中，极性修正可通过符号系数 实现，即判断检波器记录中时窗内振幅符号，需要正负号保持一致，此处保证波形的同向叠加。
[0055] 上述技术方案中，在所述步骤S6中，权重系数β大小(β>0)取决于微地震资料的信噪比情况，如果信噪比较低，则波形叠加部分权重要大些，故设置β大些，反之设置β小些。
[0056] 上述技术方案中，在所述步骤S7中，在网格搜索定位中，假设网格点是震源位置时，Tr为最小，Ews为最大(即-βEws为最小)，故改进目标函数Q为最小值。
[0057] 实施例：首先建立基于均匀介质模型的二维微地震井中观测系统，纵波速度Vp＝4500m/s，横波速度Vs＝2500m/s，密度为2.425g/cm3。如图2所示，水平距离x＝2000m处监测井(直井)中10级检波器深度z分别位于1500-1590m，其道间距为10m。设置微地震震源发生的位置为(x，z)＝(2200，1570)m。然后采用100Hz的雷克子波进行二维弹性波动方程正演模拟，采样间隔为0.5ms，得到该理论模型的z分量地震记录(图3a)，图3b是对图3a加入随机噪音后理论模型的z分量地震记录。
[0058] 为了测试本发明方法对初至误差的稳定性，对10个检波器微地震事件初至分别加入200次[-2，2]ms随机误差。然后依据公式(1)和(2)建立走时残差目标函数后，采用步骤7中网格搜索法对微地震震源进行定位，得到基于走时目标函数方法的定位结果(图4a)。依据公式(1)至(6)建立本发明的目标函数，采用步骤7中网格搜索法对微地震震源进行定位，得到本发明方法的定位结果(图4b)。
[0059] 从图4a与图4b对比可以看出，在一定初至误差下，图4b的定位结果比图4a更加收敛，即本发明定位误差要比基于走时目标函数的定位方法误差小的多。
[0060] 综上所述，本发明基于初至时差和波形叠加相结合的思路，对微地震震源进行定位，可以提高定位方法的抗噪能力及收敛性，可改善微地震震源定位的精度。
[0061] 本说明书未作详细描述的内容属于本领域专业技术人员公知的现有技术。
