一种基于广义染色算法的粘滞介质靶区成像方法

1.本发明属于数据处理技术领域，具体涉及一种基于广义染色算法的粘滞介质靶区成像方法。


背景技术：

2.从1929年以来美国湾岸地区采用地震勘探方法成功探测油气田以来，地震勘探发展成为油气勘探中应用最广泛和必不可少的方法。经过数十年的勘探实践认识，人们逐渐发现并不是所有油气藏都与背斜构造有关，隐蔽油气藏理论的提出(levorsen，1964)，为油气勘探开发提供了更广阔的前景。勘探重心也由浅层转向深层、由常规地带转向极端地带、由构造油气藏转向隐蔽的岩性油气藏。其中最引人瞩目的是易形成巨型油气藏的盐下区域，然而有限的观测系统难以接收足够多的有效反射波，复杂的上覆结构使深层能量变得复杂且微弱，导致了盐下区域在成像时照明不足，信噪比低，因此发展针对盐下区域的精确成像技术已经成为当今工业界的共识。
3.针对这一问题，国内外学者提出的解决办法大致分为两类：一类以改变观测系统为主，多方位，多尺度接收地震信号，如三维地震测井技术可以极大的丰富地下结构的照明，是一种独特的用于盐下区域成像的方法(chen，1992)。宽方位角三维采集系统可以采集丰富的方位角信息,有利于对深层、高陡构造、各向异性岩体进行成像(唐云，2003)。另一类以对盐下区域信号进行照明补偿为主，如使用对弱照明区域进行加权，同时保留强照明区域的方法，改善了盐下区域的成像效果(gherasim，2014)。多次反射波的反射路径与一次反射波的反射路径不同，所以多次波也被用来增强弱照明的盐下区域的能量信号(李鹏，2006)。多次散射波也被用来对缺少一次反射波照明区域成像(guitton，2002)。多次波成像的结果也被用来修正传统逆时偏移的结果(liu，2011)。值得注意的是一种以目标区域为导向的波动方程偏移速度分析方法,该方法通过使用广义波恩波场建模来合成特定区域的数据集，这个数据集包含速度建模所必须的速度信息，使得这个数据集适合于面向目标的、快速的速度模型构建(tang，2011)。以此类推，在地震波的传播过程中通过标记追踪通过目标区域的波场并将之应用于偏移成像，也可以得到高质量的局部成像结果，chen等(2014)引入源于生物学中“染色”的思想，提出弹性介质复数域染色算法。在此算法中，当地震波接触到染色区域，经过染色区域的地震波将被标记并且在接下来的传播过程中被追踪。应用此算法可以得到一个与目标区域有关的染色波场，并且这个波场与真实波场在传播过程中是同步的。然而，复数域染色算法所得到的染色波场受复速度的影响存在振幅远小于真实波场振幅且波形被扭曲的缺点。li等(2017)提出弹性介质广义染色算法的控制方程，其构建的染色波场与真实波场的匹配度远高于复数域染色算法，可以实现高分辨率的目标区域成像。
4.然而地下介质通常具有粘滞性，地震波在地层传播过程中会产生振幅衰减，相位畸变和频带缩窄等现象，如果忽视地层粘滞性，会导致盐下区域地震波能量更加衰弱，波形扭曲，严重影响下覆构造的成像精度。
5.目前，广泛采用的粘滞声波方程有两类：一是基于标准线性体模型的粘滞波动方程，其最早由carcione等(1988)将广义标准线性体模型(gsls)引入地震波传播理论,采用记忆变量及粘弹性参数表征波场特征；二是基于常q模型的粘滞波动方程，源于kjartansson(1979)提出的常q理论，根据该模型的定义，在地震勘探频段内(《150hz)，品质因子q几乎不随频率变化，更符合物理规律，除此之外相比于广义标准线性体模型常q模型表征的参数更少，计算也更加简洁。结合上述优点采用常q解耦粘声波方程更满足地震勘探的需求。
6.zhu(2014)基于常q模型提出一种分数阶拉普拉斯算子的粘声方程。该方程含有两个拉普拉斯算子分别对应相位频散项和振幅衰减项，可通过快速傅里叶变换进行计算，简化了数值计算，是实现粘声逆时偏移的良好方法。但是，该方法在实施过程中仍然存在诸多问题。第一，由于介质的粘滞效应的不断积累，导致深层同相轴的成像结果振幅能量极弱，以至于无法准确识别地下深层构造。为解决这一问题，zhu(2014)直接改变方程中振幅衰减项的符号进行能量补偿，但是如果使用这种补偿方式会在地震波的传播过程中产生大量高频异常值，降低了成像剖面的信噪比，以至于无法识别有效地质构造。为了解决补偿中出现高频异常值的问题，需要引入低通滤波器。然而在对盐下区域成像时，依旧存在成像时照明不足，信噪比低等问题。


技术实现要素：

7.为了解决背景技术中所述存在的问题，本发明提供了一种基于广义染色算法的粘滞介质靶区成像方法，更能提高盐下复杂油气储层的识别精度，为油气资源的快速、高效开发提供核心技术支撑。
8.本发明采用的技术方案为：一种基于广义染色算法的粘滞介质靶区成像方法，基于广义染色算法的粘滞介质靶区成像方法包括以下步骤：
9.步骤一、设置染色位置，选择与复数域染色算法标记染色位置相同的方法，
10.复数域染色算法方程：
[0011][0012][0013]
方程中：所代表的原波场，所代表的染色波场，为真实速度场，为染色区域速度场，当且仅当时该式成立，当标记的散射点与目标位置重合时，染色波场的波前将与原波场同步传播；标记区域的值只取1或0，1代表标记位置，0代表未标记位置；
[0014]
步骤二、设置震源、检波点、地震记录信息；
[0015]
步骤三、利用分数阶拉普拉斯算子常q解耦粘声波方程和边界条件求取时间步为t的震源地震波场值，并保存波场值；
[0016]
分数阶拉普拉斯算子常q解耦粘声波方程：
[0017][0018][0019]
方程中:c为速度，t为时间，p为粘声波波场；q为品质因子；c0为参考频率ω0处定义的相速度；s为震源；
[0020]
步骤四、将时间步长为t由分数阶拉普拉斯算子常q解耦粘声波方程所计算出的波场值应用于粘滞介质的广义染色算法方程，来计算时间步长为t的染色波场值，保存波场值；
[0021]
粘滞介质的广义染色算法方程：
[0022][0023][0024]
方程中
[0025]
步骤五、同时利用方程(2)和方程(3)计算时间步长为t+1的真实波场与染色波场；
[0026]
步骤六、判断时间循环是否到达最大时间步长t，否则重复步骤三至步骤五，是则结束时间循环；
[0027]
步骤七、加载时间步长为t-t的检波点波场数据；
[0028]
步骤八、利用方程(2)的波动方程和边界条件求解时间步为t-t的检波点地震波场值；
[0029]
步骤九、计算时间步长为t-t-1的检波点地震波场值，判断是否到达时间步长0；否则重复步骤七至步骤九，是则结束循环；
[0030]
步骤十、利用成像条件处理所有炮的染色波场、检波点反传波场信息、去噪后得到地下构造成像结果；
[0031]
成像条件：
[0032]
[0033]
其中r为检波点反传波场，为染色波场。
[0034]
进一步的，采用有限差分法求解等号左边的时间导数，利用伪谱法求取等号右端的空间算子；因此，时间精度为二阶精度，空间精度为谱精度，其计算方程为：
[0035][0036][0037]
其中f和f-1
分别为一维正反傅里叶变化，k为离散波数。
[0038]
进一步的，分数阶粘声方程实现了相位频散项与振幅衰减项的解耦，在波场传播过程中通过改变振幅衰减项的符号来实现补偿；然而补偿过程中补偿算子随波数增加呈指数递增，高频成分被严重放大，导致数值不稳定，为压制数值不稳定现象，施加一个低通滤波器实现稳定传播。
[0039]
本发明的有益效果：提供了一种基于广义染色算法的粘滞介质靶区成像方法，更能提高盐下复杂油气储层的识别精度，为油气资源的快速、高效开发提供核心技术支撑。其主要优点如下：
[0040]
(1)、常规复数域染色法由于极小速度的存在，染色波场振幅值远小于真实波场振幅值，本方法所使用的广义染色算法可以实现保幅的染色波场传播；
[0041]
(2)、常规的弹性波广义染色算法中忽略了地层的吸收衰减作用，不能准确的描述地震波在地下中的传播规律，本方法使用分数阶拉普拉斯算子常q解耦粘声波方程将染色算法推广至粘滞介质，能较为准确的刻画地震波在地下传播过程中的衰减机制，而利用该方程的解耦特性，振幅和相位的解耦，有利于衰减补偿逆时偏移；
[0042]
(3)、在波动方程中引入低通滤波器，缓解了由于直接补偿而引起的高频异常值问题；
[0043]
(4)、基于染色算法的粘声介质靶区逆时偏移成像流程，提高了复杂构造靶区成像剖面的信噪比和分辨率，有利于在地震资料解释中准确识别地下构造及预测储层区间。
附图说明
[0044]
图1是基于广义染色算法的粘声介质靶区成像流程图；
[0045]
图2是水平层速度模型图；
[0046]
图3是水平层q值模型图；
[0047]
图4是复数域染色算法标记位置图；
[0048]
图5是广义染色算法标记位置图；
[0049]
图6是原波场图；
[0050]
图7是粘声复数域染色算法所求染色波场图；
[0051]
图8是原波场和粘声复数域染色算法所求染色波场对比图；
[0052]
图9是粘声广义染色算法所求与染色波场图；
[0053]
图10是原波场图和粘声广义染色算法所求与染色波场图；
[0054]
图11是盐丘速度模型图；
[0055]
图12是盐丘q值模型图；
[0056]
图13是衰减成像图；
[0057]
图14是声波成像图；
[0058]
图15是补偿成像图；
[0059]
图16是衰减复数域染色算法成像图；
[0060]
图17是声波复数域染色算法成像图；
[0061]
图18是补偿复数域染色算法成像图；
[0062]
图19是衰减广义染色算法成像图；
[0063]
图20是声波广义染色算法成像图；
[0064]
图21是补偿广义染色算法成像图；
[0065]
图22是水平位置2.4km处的单道信息图。
具体实施方式
[0066]
实施例一
[0067]
参照各图，
[0068]
一种基于广义染色算法的粘滞介质靶区成像方法，其特征在于：基于广义染色算法的粘滞介质靶区成像方法包括以下步骤：
[0069]
步骤一、设置染色位置，选择与复数域染色算法标记染色位置相同的方法，复数域染色算法方程：
[0070][0071][0072]
方程中：所代表的原波场，所代表的染色波场，为真实速度场，为染色区域速度场，当且仅当时该式成立，当标记的散射点与目标位置重合时，染色波场的波前将与原波场同步传播；标记区域的值只取1或0，1代表标记位置，0代表未标记位置；
[0073]
步骤二、设置震源、检波点、地震记录信息；
[0074]
步骤三、利用分数阶拉普拉斯算子常q解耦粘声波方程和边界条件求取时间步为t的震源地震波场值，并保存波场值；
[0075]
分数阶拉普拉斯算子常q解耦粘声波方程：
[0076]
[0077][0078]
方程中:c为速度，t为时间，p为粘声波波场；q为品质因子；c0为参考频率ω0处定义的相速度；s为震源；
[0079]
步骤四、将时间步长为t由分数阶拉普拉斯算子常q解耦粘声波方程所计算出的波场值应用于粘滞介质的广义染色算法方程，来计算时间步长为t的染色波场值，保存波场值；
[0080]
粘滞介质的广义染色算法方程：
[0081][0082][0083]
方程中
[0084]
步骤五、同时利用方程(2)和方程(3)计算时间步长为t+1的真实波场与染色波场；
[0085]
步骤六、判断时间循环是否到达最大时间步长t，否则重复步骤三至步骤五，是则结束时间循环；
[0086]
步骤七、加载时间步长为t-t的检波点波场数据；
[0087]
步骤八、利用方程(2)的波动方程和边界条件求解时间步为t-t的检波点地震波场值；
[0088]
步骤九、计算时间步长为t-t-1的检波点地震波场值，判断是否到达时间步长0；否则重复步骤七至步骤九，是则结束循环；
[0089]
步骤十、利用成像条件处理所有炮的染色波场、检波点反传波场信息、去噪后得到地下构造成像结果；
[0090]
成像条件：
[0091][0092]
其中r为检波点反传波场，为染色波场。
[0093]
采用有限差分法求解等号左边的时间导数，利用伪谱法求取等号右端的空间算子；因此，时间精度为二阶精度，空间精度为谱精度，其计算方程为：
[0094]
[0095][0096]
其中f和f-1分别为一维正反傅里叶变化，k为离散波数。
[0097]
分数阶粘声方程实现了相位频散项与振幅衰减项的解耦，在波场传播过程中通过改变振幅衰减项的符号来实现补偿；然而补偿过程中补偿算子随波数增加呈指数递增，高频成分被严重放大，导致数值不稳定，为压制数值不稳定现象，施加一个低通滤波器实现稳定传播
[0098]
由于粘声方程的计算量巨大，为了降低算法时间复杂度，满足工业需求，引入多gpu并行加速波场模拟的计算，同时使用cuda命令中的cufft加速求取空间导数的傅里叶变换。
[0099]
从图10中可以得出由粘滞介质广义染色算法方程所计算出的染色波场值与原波场值振幅，相位相同，可以的到保幅的染色波场值；在实际中如果不进行补偿，就会如图13所示，令盐下深层成像效果十分微弱。以声波数据成像作为参考解，图15虽然经过补偿处理，深层能量得到一定程度的补偿，但结构依旧不清晰。应用染色算法后，图18染色波场由于的存在成像振幅、分辨率与图21差距较大。图21相比于图18靶区成像分辨率得到明显提升，靶区构造清晰，断层干脆，断点清晰，能更好的实现对盐下深层构造的准确成像。为了进一步说明染色算法的优越性，抽取了图14、15、18、21水平方向2.4km处单道信息，其中蓝色虚线处为目标区域第一层位。图22也可以验证本文提出的新染色算法与广义染色算法都可得到成像振幅高、分辨率高，构造归位正确的成像结果。