基于自适应镶嵌的地震剖面数据可视化方法①

陈 雷 石家琦 石太昆 张 婷

(*北京信息科技大学计算机学院 北京 100101)

(**国勘数字地球(北京)科技有限公司 北京 100096)

(***北京工业大学信息学部 北京 100124)

0 引言

地震数据剖面可视化技术对于地震勘探中的解释环节有着关键性的作用[1-2]。地震数据直观、正确、特征层次分明的可视化可以极大地方便地质科研人员对相关数据做出准确的分析、得出正确的结论,对油气勘探、地质勘察等具有重要意义。

目前,地震数据可视化技术主要有2 类:一类是二维可视化[3-6],另一类是三维可视化[1-2,7-16]。地震数据二维可视化是将地震剖面数据转换成一幅二维图像进行显示的技术,又可分为波形显示、波形+变面积显示、灰度显示和彩色变密度显示等4 个类别。这种传统可视化方法一定程度上展示了地震数据的空间变化特征,其不足是动态范围小、细节丢失严重、不能直观显示微观地质结构、不便于地震解释和地质分析等。近年来,随着计算机图形学的发展和显卡中图形处理单元(graphics processing unit,GPU)计算能力的提高,地震数据三维可视化获得了广泛研究,取得了较大进展。它是将地震数据在三维空间中进行可视化的技术,是为解决二维可视化的不足而提出的。但现有地震数据三维可视化基本都是对三维地震数据体进行可视化,再按需抽取某一地震剖面的灰度或彩色变密度显示[1-2,7-16],这种可视化方法对地震微观结构的揭示能力并没有实质提高。文献[15]提出了将二维地震剖面进行三维可视化的方法,提高了对微观结构的显示能力,但在顶点镶嵌方面存在镶嵌方式不能根据地震样点倾斜方式自适应选择的问题。

针对文献[15]中可视化方法的网格顶点镶嵌问题,本文提出了一种新的低倾斜相关自适应镶嵌方法,以尽可能使网格顶点根据相应地震样点的倾斜方式来自适应镶嵌。

1 算法分析

本节先对文献[15]中的地震剖面三维可视化方法进行分析,然后再对所提镶嵌方法进行论证。

1.1 地震剖面数据三维可视化方法分析

文献[15]提出的二维地震剖面三维可视化方法具有以下特征。

(1) 对二维地震剖面数据S(x,y),x表示道号,y表示时间深度值,以(x,y) 处的地震样点幅度S(x,y) 作为三维空间中的z坐标,即地震样点对应的三维空间中的点为P(x,y,S(x,y))。

(2) 除要显示的最后一道地震数据外,将每一道上相邻2 个地震样点和其下一道上时间深度相同的相邻2 个样点连接为四边形,由同一道上3 个地震样点和下一道上1 个样点连接成的2 个三角形单元组成,以此将地震剖面数据镶嵌为网格曲面。

(3) 创建一幅由一种颜色渐变到另一种颜色的调色板纹理,将经过光栅化插值之后的地震幅度值S(x,y) 通过线性函数映射为纹理坐标,在GPU 像素着色器中根据此纹理坐标采样调色板纹理为网格顶点和三角形单元内部的像素点着色。

(4) 将相邻道上的4 个地震样点镶嵌为四边形时,镶嵌有2 种方式,即上倾斜镶嵌和下倾斜镶嵌;不正确的镶嵌方式将产生错误的结果(脊或者槽)。对此问题,采用了低倾斜相关自适应镶嵌,在GPU像素着色器中计算以当前样点为中心的某一长度的样点序列与下一道上一定偏移范围内的同样长度序列的相关系数,将最大相关系数所对应的偏移值作为当前样点的倾斜方式,通过渲染到纹理将数值化的倾斜方式存储到纹理中,用来定义网格顶点的索引值。

(5) 将坐标变换、光照计算、纹理坐标映射等内存消耗较大的计算都转移到GPU 中进行,将地震剖面数据填充到一幅浮点纹理中,在GPU 中通过采样浮点纹理实现对地震剖面数据的访问。

1.1.1 网格顶点着色及颜色插值

在固定渲染管线和可编程渲染管线中,顶点数据经过初始装配之后的处理如图1 所示。

图1 渲染管线中经过初始装配之后的顶点处理过程

初始装配之后的顶点数据需进行光栅化与插值。固定渲染管线中,光栅化与插值输出的像素经过光栅器后直接送往帧缓冲等待显示;可编程渲染管线中,光栅化与插值输出的像素先经过可编程像素着色器处理,再送往光栅器直至显示。在像素着色器中,可根据顶点着色器输出的信息对顶点颜色、纹理坐标、光照等进行处理,然后渲染输出。

地震剖面三维可视化中,对不同幅度的样点所对应的网格顶点及其之间的像素用不同的颜色进行着色,可有效提高对地震事件的分辨能力。着色有2 种方式:一种是将颜色作为顶点属性分量为顶点及其之间的像素着色,可称为顶点属性着色方式;另一种是将光栅化与插值输出的像素传送到GPU 可编程像素着色器中,在像素着色器中根据顶点着色器输出的地震幅度对调色板纹理进行采样,在其中,先将线性插值后的地震样点幅度线性映射为纹理坐标,然后再采样调色板纹理,以实现网格顶点及顶点之间三角形单元内部像素的正确着色,这种方式可称为线性渐进调色板着色方式。文献[15]采用第2种着色方式有效避免了地震数据可视化中网格顶点之间像素的颜色交叉问题。

1.1.2 网格顶点低倾斜相关自适应镶嵌

网格顶点镶嵌是地震剖面数据三维可视化中着色之外的另一个关键问题,即如何在三维空间中将地震样点对应的顶点连接成三角形单元,以形成网格曲面。文献[15]所用方法是除要显示的最后一道地震数据外,将每一道上相邻2 个地震样点和其下一道上时间深度相同的相邻2 个样点连接为四边形(此四边形由同一道上2 个地震样点和下一道上1个样点连接成的2 个三角形单元组成),以将地震剖面镶嵌为网格曲面。对以上4 个顶点,不同镶嵌方式会产生不同的效果,如图2 所示。

图2 相同4 个顶点的2 种镶嵌方式

采用图2(a)所示的上倾斜镶嵌方式会产生“脊”结构如图2(c)所示;采用图2(b)所示的下倾斜镶嵌方式会产生“槽”结构如图2(d)所示;若统一采用上倾斜镶嵌,则对下倾斜地震事件不合适,会在事件顶部产生“锯齿”现象,反之亦然。

为解决倾斜事件顶部的“锯齿”问题,文献[15]提出了一种低倾斜相关自适应镶嵌方式,在GPU 像素着色器中计算以当前地震样点为中心的某一长度的样点序列与下一道上同样长度序列在某一滑动窗口(从上向下滑动)内的相关系数,并找出最大相关系数对应的滑动值。如果此滑动值大于等于0,则采用下倾斜镶嵌,用-1 表示,即以当前顶点为控制点(四边形左上角的顶点)的四边形采用下倾斜镶嵌;反之,如果此滑动值小于0,则以当前顶点为控制点的四边形采用上倾斜方式镶嵌,用+1 表示。计算出镶嵌方式之后,再把数值化的镶嵌方式通过渲染到纹理存储到纹理元素中,用以填充网格顶点索引。低倾斜是指在求2 个样点序列的互相关函数值时,只允许滑动窗口有上下一个样点间隔的滑动范围,这样作为四边形控制点的当前道上的地震样点只能与下一道上上倾斜一个采样间隔或下倾斜一个采样间隔的样点相连。没有采用高倾斜相关自适应镶嵌方式的原因是高倾斜情况下,镶嵌方式过于复杂,很难处理自适应镶嵌带来的网格漏洞和三角形单元重叠问题。

文献[15]所提低倾斜相关自适应镶嵌方式的不足是:只有当四边形控制点的倾斜方式为0(此时采用下倾斜镶嵌)或-1(下倾斜)时,组成此四边形的2 个三角形单元的连接方式才与控制点的倾斜方式一致;如果控制点的倾斜方式为1,那么此控制点所控制的4 个顶点都要按照同一地震道上控制点之下一个样点的上倾斜方式来镶嵌,无论下一个样点的倾斜方式是上倾斜、不倾斜还是下倾斜。这样很大一部分网格顶点的镶嵌方式不能与地震样点的倾斜方式相一致。如图3 所示。

图3 低倾斜相关自适应镶嵌方式的不足说明

文献[15]中的方法以四边形为基本单元镶嵌地震数据网格,如图3 所示,以四边形ABCD为例,A点为四边形的控制点,四边形的镶嵌方式由A点的倾斜方式决定。当A为下倾斜(或不倾斜)时采用图3(a)所示的下倾斜镶嵌方式镶嵌四边形ABCD,此时镶嵌方式与A点的倾斜方式一致;当A为上倾斜时采用图3(b)所示的上倾斜方式镶嵌四边形ABCD,而此时按相关系数最大化原则,A点与N点相关度最大,A点应与N点相连,生成四边形MACN。而文献[15]中方法生成了四边形ABCD,这等于忽略了A点的上倾斜方式,同时将B点的倾斜方式强制设为上倾斜,当B点的倾斜方式不是上倾斜时,相当于A、B2 个点都没有根据其倾斜方式进行镶嵌,故文献[15]所提低倾斜相关自适应镶嵌存在镶嵌错误问题,需进一步改进。

1.2 本文算法

为解决1.1.2 小节中低倾斜相关自适应镶嵌方式存在的问题,尽可能使得每一个顶点都按照相应样点的倾斜方式来镶嵌(对于边界顶点和产生漏洞的区域,为了保证网格的完整性与连续性,顶点镶嵌方式可能与地震样点的倾斜方式不一致),避免文献[15]中低倾斜相关自适应镶嵌方法的不足。本文提出了新的低倾斜相关自适应镶嵌方法,具体如下。

首先在GPU 像素着色器中计算相邻地震道上2个样点序列的互相关函数值,以求取地震剖面中各地震样点的倾斜方式及相应的镶嵌方式。GPU 像素着色器中2 个样点序列的互相关函数值的计算原理,如图4 所示。

图4 计算互相关函数值原理图

互相关函数和相关系数的计算公式分别为

其中,式(1)是计算序列X和Y的互相关函数值RXY(τ) 的公式,下标XY表示2 个序列,X、Y的顺序表示序列X未移动,序列Y移动τ个单位;式(2)中表示序列X和Y的相关系数大小。τ是移位参数,T为序列周期,N表示序列X和Y中元素编号的最大值,k表示元素序号,k∈[0,N]。计算出相关系数后,找出最大值,根据其对应的τ值来判断当前地震样点是上倾斜还是下倾斜或者不倾斜,上倾斜用1 表示,下倾斜用-1 表示,倾斜幅度较小时视为不倾斜,用0 表示。也就是说将最大相关系数对应的τ值取反就可得地震样点的倾斜方式。然后,使用渲染到纹理这一操作,把代表倾斜方式的数值存储到纹理中,因纹理元素有R、G、B、A等多个分量,这里只用R分量来存储倾斜方式数值,其他分量均为0。

获得每个地震样点的倾斜方式之后,根据倾斜方式的不同对地震样点进行自适应镶嵌,镶嵌时分为3 种情况,具体如下。

(1) 顶部边缘样点均采用下倾斜镶嵌,如图5所示。图5(a)为上边缘点镶嵌示意图,图5(b)为下边缘点镶嵌示意图,T和B为当前要镶嵌的地震样点,它们分别处于地震剖面的上边缘(时间深度值最小)和下边缘(时间深度值最大)。对边缘点采用统一的镶嵌方式是为了避免当同一地震道上顶部2 个或者底部2 个地震样点的倾斜方式相反时,由于其之间的矩形区域内没有镶嵌三角形而形成漏洞,接下来就可以对其余非上下边缘点进行自适应镶嵌。

图5 边缘点镶嵌示意图

(2) 当非上下边缘点的倾斜方式为上倾斜时,为避免同一矩形内的4 个样点被重复镶嵌而产生三角形交叉重叠现象,先查看其上面相邻地震样点的倾斜方式是否为下倾斜,如果不是,则进行上倾斜镶嵌,如果是,则不进行镶嵌,如图6 所示。

图6 非上下边缘点上倾斜镶嵌示意图

图6 中,B点为当前样点,图6(a)为当前一个样点A的倾斜方式为下倾斜时产生的交叉现象。此时不对B点进行镶嵌,即不给B、D和C分配顶点索引以形成三角形BDC;图6(b)为A点为上倾斜时的情况。

(3) 当非上下边缘点的倾斜方式为下倾斜时,先进行下倾斜镶嵌。然后,为了避免产生漏洞现象,查看其上面相邻样点的倾斜方式是否为上倾斜,如果是,则对上面地震样点和当前样点之间的矩形区域进行镶嵌,以对漏洞进行缝补;如果上面样点的倾斜方式不是上倾斜,则不用进行缝补,如图7 所示。

图7 非上下边缘点下倾斜镶嵌示意图

图7 中,E为当前样点,图7(a)为当B点为上倾斜时产生的漏洞现象,图7(b)和(c)是可以采用的2 种漏洞缝补方法(此2 种方法任选一种即可)。

本文提出的镶嵌方法,与文献[15]每个四边形单元中的4 个样点的镶嵌方式都由其左上角样点的倾斜方式决定的方法不同,可最大程度地根据每一个地震样点的倾斜方式对其进行真正的自适应镶嵌。

2 实验结果及分析

为验证本文算法的有效性,在个人计算机上对新疆某工区的一条二维线中的标准segy 数据进行了实验。运行环境为Windows 7 旗舰版,处理器是i3-3240@3.40 GHz,显卡为ATI Radeon HD 6450A/7450A(1 GB/联想),内存为4 GB(海力士DDR3L 1600 MHz),实验平台是Microsoft Visual Studio 2010,图形库是Dirtect3D 10,着色器shader 4.0。

本文基于以上环境开发的地震剖面三维可视化的系统流程图,如图8 所示。

图8 本文地震剖面三维显示系统流程图

图9～13 所示为不同镶嵌方式对从标准segy剖面数据中选取的一块含有下倾斜和上倾斜结构的区域进行可视化的结果,大小为200×80 样点。

图10 固定模式下倾斜镶嵌

对比图9～13 可知,固定模式镶嵌没有考虑相邻道样点之间的相关性,采用固定上倾斜或下倾斜镶嵌方式镶嵌顶点时,分别会在下倾斜事件和上倾斜事件的顶部产生锯齿(槽或脊),如图9 和10 所示。文献[15]采用低倾斜相关自适应镶嵌方法在一定程度上解决了此问题,如图11 所示。虽然相比固定镶嵌模式,在消除锯齿效果上有了一定的改进,但由于其会出现镶嵌方式与样点倾斜方式不一致的情况,存在镶嵌错误问题,所以锯齿现象仍然很明显。如果根据其方法对每个样点都根据其倾斜方式进行自适应镶嵌,则会产生交叠和漏洞现象,如图12所示。

图11 文献[15]低倾斜相关自适应镶嵌

图12 对每个样点都用文献[15]方法进行自适应镶嵌

图13 是用本文提出的新的低倾斜相关自适应镶嵌方法得到的结果,该方法最大程度地保证了每个顶点的镶嵌方式都与相应样点的倾斜方式一致,只有在产生三角形单元交叠或网格漏洞的情况下,少数顶点的镶嵌方式才可能会与相应样点的倾斜方式不一致,而这里的不一致并不会造成倾斜事件顶部的锯齿现象发生,因为其并未处于倾斜事件上。

图13 本文新的低倾斜相关自适应镶嵌

综合以上对比分析可知,本文所提出的新的低倾斜相关自适应镶嵌方法效果较好,其有效性和可行性得到了验证。

3 结论

本文针对二维地震剖面数据的三维可视化问题进行研究,对于地震数据网格的顶点镶嵌问题,提出了一种新的低倾斜相关自适应镶嵌方法,并通过实验验证了算法的有效性和实用价值。进一步研究的方向是对网格的细分方法进行研究,以求将地震数据可视化的同时,尽可能减少由镶嵌带来的锯齿等现象。