提高信噪比的数字滤波技术及其vc++实现.doc

提高信噪比的数字滤波技术及其VC++实现 论文提要 滤波是地震资料数据处理中最常用的一种方法。本文首先介绍了地震数据处理中的数字滤波技术的基本原理，然后利用VC++程序设计工具进行滤波程序的编译和调试，在此基础上，选择不同的滤波参数，对滤波结果进行了对比和分析。在软件设计的整个过程中，融入了软件工程的基本思想。所编译的程序能够在一定程度上对地震数据进行简单快速的滤波处理。并具有一定的滤波效果，但由于一维滤波本身的局限性，滤波虽然有一定效果，但存在不足，需要发展更加优越的滤波技术，如二维视速度滤波、非线性滤波、最小平方滤波等。 正文 一、 前言 （一）C程序上的实现来展开讨论。 二、 一维时间域滤波的设计与实现 根据软件工程的理论，像这样一个软件(程序)的生命周期包括以下几个步骤：问题定义、可行性研究、需求分析、总体设计、详细设计、编码和单元测试、综合测试、软件维护这几个阶段。前面已经对方法基本原理、思想进行了介绍，本章将主要按上面几个部分来进行介绍。 这次要实现的是时间域上的一维滤波。之所以选择它来实现，其一，是因为它可以在很大程度上体现出数字滤波基本原理与思路。其二，无论频率滤波还是时间域滤波都是对数据进行频率上的处理，本质思想相同。其三，和其他滤波方法相比，这个更容易在VC++上实现。 (一) 程序的总体设计 在结构上，程序如果要顺利运行出结果，就应该包括以下几个组成部分，即：SEGY读取部分，傅立叶变换部分，滤波部分，输出部分。 下面逐个部分的给出程序组成部分 第一部分：读取SEGY数据 SEGY数据是地震数据的通常标准形式。不同的处理软件所接受的数据格式并不相同，所以一般都先存为SEGY文件，然后再转为可处理的文件格式。这里就是把SEGY转为C语言可以接受的数据文件的步骤。 使用的时候，要注明输入文件与输出文件的名称。SEGY文件中包含有许多信息，包括3200字节的卷头信息，400字节道头信息，最后才是数据。前两项可以根据需要进行取舍，数据才是要读取的核心内容，所以在读取的时候要根据需要跳过一些字节。最后处理完输出到指定的文件之内。 第二部分：傅立叶正变换 将SEGY文件处理完之后，还需要对其进行傅立叶变换，变换为频谱，才能进行频率滤波。不过单纯的傅立叶变换对计算机运行速度非常浪费，所以这里采用的是一种快速傅立叶算法。 第三部分：滤波部分 如图2-1为未经滤波处理原始炮集剖面，从上图可以清晰看到有效信号中有噪音和面波的干扰。而通过滤波所要去掉的就是图中的干扰。 图2-1 原始地震数据图 由于屏幕限制，上面只截取了1到111道，每道中200到900多个样点的图形。 增益后的信号简单的说就是将地震波在地下传播中损失的能量进行补充处理后的信号。具体过程比较复杂，这里不再过多解释，但是从图形可以更明显的看出噪音干扰和面波干扰。 利用OMEGA软件可以先划定出需要处理的区域，然后在此区域进行去面波处理。可以得到精度更高的结果。 至于频率带宽的选取则不尽相同。应当根据所收集数据的地区和以往经验来确定。而且往往需要不断的调整以便得到最佳效果。 第四部分：傅立叶逆变换 由于进行的是时间域的滤波变换，所以在滤波完成后可以不用再进行傅立叶逆变换。但是如果是进行的频率域滤波，则需要再做傅立叶逆变换，然后才能读取出来。 第五部分：数据转换为SEGY文件 处理完的数据一般最后都转变为SEGY的通用格式。 整个程序的结构如图2-2： 图2-2 程序结构图 (二) 程序编译 编译所使用的工具为创天中文VC++.操作系统及WINDOWSXP2。 程序大体结构是按照上面所说的结构来安排的。在主程序中调用读取、傅立叶变换、确定主频、滤波等的子程序。 下面给出所使用参量的对应关系。 表2-1 程序内参量对应关系 Dt tl L Nj LN flag T NR X maxf 采样间隔 采样点数 接收道数 2的N次方 滤波因子长度 数据内包含的道数 正/负傅立叶变换 2的阶数 存放读出来的数据 所选道的主频 此外在各个模块中，zN代表噪音函数。fft子函数用于进行傅立叶正/逆变换。filter子函数用于进行滤波。 具体源码可以参见附属文件，具体在编译中遇见的问题也将在相关的地方有具体说明。 (三) 程序测试 测试内容包括以下几个方面：一个是在编译过程中利用编写工具查看是否存在编写错误等问题，再一个是要用实际数据来验证此程序能否按照编写的意愿进行工作，此外，还测试了本程序适用的范围。 在东方地球物理公司的实习中，我有幸得到了2个单炮数据，其中一个是另一个数据经增益处理后的结果。这两个文件就成为测试程序的基础。具体使用中由于增益文件的效果更明显，所以图片都是在增益文件的基础上做成的。 首先是