基于小波变换的SEG-Y格式地震数据压缩

汤君健，张正炳，付青青 (长江大学电子信息学院，湖北 荆州 434023)

地震勘探数据压缩方法按照有无信息丢失可分为无损压缩和有损压缩2大类：无损压缩能保证解压数据与原始数据完全一样，特别适合于对解压后还需要做进一步处理的地震数据进行压缩，但其压缩效率低，一般压缩倍数在2倍左右；有损压缩在允许有少量信息损失的情况下可以获得高倍压缩，压缩倍数较高，但存在信息丢失等问题[1]。

小波变换由于其在时域、频域内具有良好的局部分析特性，在地震数据压缩领域一直是学者们关注的热点[2～4]：在地震勘探过程中，获取的地震信号是频率依赖的，受大地滤波的影响，地震波类似于小波的结构，正好与小波的特征相一致；小波变换属于变换域数据压缩方法的一种,它具有很好的局部化特性,地震信号经过小波变换后的能量集中在少数变换系数上,用它代替离散余弦变换并合理地利用其变换系数的特点,可获得较好的压缩效果,同时可克服JPEG方法产生的方块效应；同时地震波信号有高度的非平稳性,所以有良好局部分析能力的小波变换就成为地震数据压缩的主流方法。王喜珍等[5]在研究地震数据的有损压缩过程中，通过试验发现使用bior(3,7)小波函数，分解层数为5时信噪比较高，均方差较小，验证了地震数据压缩的效果与小波变换时滤波器系统长度和分解层数有关；武文波等[6]选用基于小波变换的地震信号压缩方法，当压缩比达到50∶1以上时,经解压恢复后的信息没有明显损失；徐峰涛等[7]提出了一种基于小波变换的零树编码与算术编码相结合的地震数据压缩算法，其算法对叠前地震数据压缩4倍时，信噪比可达50dB以上，选用适用小波基对叠后地震数据压缩16倍时，信噪比仍可达30dB以上；张正炳等[8]针对地震勘探数据压缩改进JPEG 2000算法，对叠前数据压缩4倍以下时，信噪比高于42dB，对叠后数据压缩15倍以下时，信噪比高于30dB。于文茂等[9]提出基于SPIHT改进算法的地震数据压缩方法，对叠后地震数据压缩16倍时，信噪比可达28.78dB。

徐峰涛等[7，8]基于小波变换的压缩算法所得地震数据压缩倍数较高，数据解压效果较好。但其只是对SEG-Y文件中的地震采样数据集进行压缩，忽略了SEG-Y文件的格式信息，这就导致其解压后的数据缺乏格式信息，他人无法读取解析；其研究的是一个具有固定道数和采样点数的地震采样数据集(800道,每道800个采样点)，但是实际地震勘探中获取的地震数据SEG-Y文件的道数与采样点数是任意的，绝大数不满足进行多次小波变换的条件；另外，压缩实际地震数据时，由于数据量极大，如果将按照其压缩方案全部数据作为一个整体进行压缩编码，则需要占用很大的内存，且计算速度慢，十分不利于实际应用。为此，笔者对其压缩方案进行如下改进，将其算法研究向实际应用方向推进一步：读取标准SEG-Y文件，将其分割为卷头、道头数据和地震采样数据集两部分，对卷头、道头数据保留原始信息不压缩；对地震采样数据集进行分块压缩处理，先按采样点数进行数据扩展填充，使其满足多次小波变换条件，采样点数扩展值为划分块的宽度值，这样每块的道数和采样点数都满足多次小波变换条件；比较3种填充方法，探索最佳扩展填充方式。

1 SEG-Y格式地震数据

标准SEG-Y格式，是勘探地球物理学会(SEG)制定的、在地震勘探中最常用的地震数据格式[10]。标准SEG-Y格式文件的构成如图1所示，由以下2部分组成：

1)卷头数据。共3600字节，由一个3200字节的EBCDIC卡和一个400字节的二进制编码头组成，其中(3201-3260)中定义了对整个SEG-Y格式文件有效的信息,包括工作识别号、测线号、卷号、地震道数、辅助道数、采样率、每道样点数、样点数据类型、覆盖次数等27个数据项。

2)地震道数据。由m个地震道数据块组成，每一个地震道数据由一个240字节的道头数据和多个采样数据构成。地震道的每个采样数据值以4个连续字节记录，采用IBM格式，4个字节组成一个32位浮点数，由一个符号位qs、7位阶数qc和24位尾数qf组成。其中，qs表示指定数值是正还是负；qc表示16的幂；qf为一个6个字节的二进制数字，其基数点在有效位的左边[11]。

二维离散地震采样数据表示为:

f(x,y)x=0,1,…,n-1；y=0,1,…,m-1

(1)

式中：n为道数；m为每道采样点数。

笔者采用3种扩展填充方法，填充“0”、按边界对称复制原始数据填充和边界值复制填充。

2 数据压缩/解压方案

基于小波变换的SEG-Y格式地震数据压缩/解压方案如图2所示。

1)SEG-Y格式文件的分割。读取SEG-Y格式文件，将其分割为卷头、道头数据和地震采样数据。

2)地震采样数据格式转化。SEG-Y格式文件中地震采样数据记录格式为工作站的32位IBM浮点数。由于微机和工作站的 CPU 架构不同，微机不能对工作站的IBM浮点数直接进行读取和利用，因此，工作站的IBM浮点型数据与微机的IEEE浮点型数据之间要相互转换[12～13]。

4)对各个地震采样数据块分别进行处理。先进行小波正变换，采用Mallat算法，将其分解为不同分辨率的子带，EZW编码则利用小波变换后不同子带上，相同位置的小波系数之间的空间相似性进行编码[15]。算术编码将EZW编码形成的二进制文件进一步进行无失真压缩，提高压缩效率，形成最终压缩数据文件，压缩过程结束。

5)解码时先分块进行算术解码，后进行EZW解码，再进行小波逆变换，依次得到各块地震采样数据。根据原始地震采样数据尺寸，对解压过程中小波逆变换得到的地震采样数据块进行裁剪，去除掉数据扩展填充部分。

6)将各块地震采样数据块按原始道号顺序拼接，重组地震采样数据部分；之后再将地震采样数据格式由IEEE浮点数转化为IBM浮点数；最后将原始卷头、道头数据和IBM格式地震采样数据按SEG-Y格式文件结构重组，形成解压重构后的SEG-Y格式文件，解压过程结束。

3 试验结果

试验数据为原始叠前地震勘探数据A，大小为4.28MB，其地震采样数据集为856道，每道1250个采样点(见图3)。

由于地震道数小于填充后的样点数，所以可看作按笔者方法分块的最后一块的情况。试验将根据原始地震采样数据(不包含填充部分)与解压重建数据的差值计算得到的信噪比(SNR)作为技术指标，研究填充值对压缩质量的影响。其中，全局信噪比为原始地震856道数据与其解压重建数据的信噪比，边界信噪比为最后一道原始数据与其解压重建数据的信噪比。采用3种扩展填充方式：方法1为“0”值填充；方法2为按边界对称复制原始数据填充；方法3为边界值复制填充。基于已有研究[7]，选用Antonini小波基，对叠前地震勘探数据A做4层小波变换，根据笔者算法将其扩展填充为1264道，每道1264个采样点。

根据3种扩展填充方法压缩倍数与全局信噪比和边界信噪比的关系曲线如图4所示。在相同压缩倍数下，“0”值填充无论是全局还是边界，重构质量都相对最好，对压缩质量影响相对最小，那是因为“0”值填充只扩充了极少量信息；按边界对称复制填充对压缩质量影响相对最大，那是因为其扩充了大量信息，使编码器编码了大量携带扩充区域信息的小波变换重要系数，浪费了编码资源，从而大大降低了压缩质量。

选用Antonini小波基，对叠前地震勘探数据A，扩展时进行“0”值填充，做4层小波变换，如表1所示，在压缩小于4倍时，信噪比大于48dB，如图5(a)所示。图5与图3相比较，与原始数据无明显差别，重构误差相对较小，可用于实际地震数据压缩。在压缩小于16倍时，信噪比大于13dB，如图5(b)所示，与图3相比较，重构误差较大，仅仅可用于质量监督。

4 结语

笔者采用地震数据存储标准SEG-Y格式文件作为压缩原始文件，可以实现IBM与IEEE浮点数据之间的转化，对地震采样数据集进行分块处理，每块数据进行扩展填充操作，使地震采样数据集适合进行多次小波变换。试验结果表明，对SEG-Y格式文件中地震采样数据分块处理，数据扩展填充时使用“0”值填充，对压缩效果影响相对较小，对叠前地震数据仍然保持着较高的压缩效率，叠前地震数据压缩4倍以下，信噪比高于48dB，对后续处理影响不大，有利于后续地震数据进一步有效利用。