CN103592682A - 一种地震波场能量均衡方法及处理装置 - Google Patents

Abstract

本发明提供一种地震波场能量均衡方法，包括：S1：处理装置读入经过预处理的地震数据；S2：处理装置对输入的地震数据，用异常振幅消除技术逐炮进行单炮能量均衡；S3：处理装置对单炮能量均衡处理后的地震数据，进行多炮集记录的炮间能量均衡处理。上述能量均衡方法，从能量角度来消除这两个环节的差异，没有涉及到地震波场的相对关系，在进行地震数据能量均衡后，得到的结果是浅层能量强于深层能量、近道能量强于远道能量，保持了波场的相对关系。此外本发明还提供一种地震波场能量均衡处理装置实施例，可以实现保持地震波场特征的能量均衡处理。

Description

一种地震波场能量均衡方法及处理装置

技术领域

本发明涉及一种地震数据处理方法，特别涉及一种地震波场能量均衡方法及处理装置。

背景技术

地震勘探是一种利用地下介质的弹性和密度的差异，通过观测和分析大地对人工激发地震波的响应，推断地下构造形态和岩层性质的方法。地震勘探是钻探前勘测油气资源的重要手段之一，通过地震勘探，可以了解地下的地质状况，寻找和查明石油与天然气资源(即油气资源)。同时，在煤田和工程地质勘察、区域地质研究和地壳研究等方面，地震勘探也得到广泛应用。地震勘探通常分为折射波地震勘探和反射波地震勘探，目前在石油和天然气资源勘查领域，主要采用反射波地震勘探。

地震勘探的过程主要分为三个阶段，分别是：野外地震数据采集、地震数据处理和地震资料解释。第一个阶段为野外地震数据采集，这一阶段的主要工作内容是：在勘探区域，布置二维或三维测线；使用炸药震源或可控震源激发地震波，炸药震源或可控震源所在点被称为炮点；沿着测线等间距布置多个检波器来接收地震波信号，检波器或检波器组的数量可以根据需要进行设定，各检波器组之间排列方式可以不同，如中间放炮排列、端点放炮排列等，也可以根据需要进行设定；检波器接收到地震波信号后，以等时间间隔离散采样地震数据，并以数字形式记录在磁带上。第二个阶段为地震数据处理，这一阶段的主要工作内容是：以地震波传播理论为基础，利用计算机和地震资料处理软件，加工处理第一阶段中野外采集的地震数据，将地震数据变成能够反映地下构造的地震剖面图和能够反映地下岩性变化的地震波振幅、频率和传播速度等信息剖面图。第三个阶段为地震资料解释，这一阶段的主要工作内容是：是分析解释第二阶段中地震数据处理所得到的地震剖面图等资料，根据石油地质学原理和地震波传播理论，确定地下岩层的构造，从而找出含油气的有利区带并提出钻井位置。

现阶段对野外采集的地震数据的处理主要是以叠前波动方程偏移为核心，该理论主要利用整个地震波场信息来实现复杂构造成像，可以获得空间分辨率较高的地震资料。地震波在地下介质传播的过程中，由于受到扩散和大地吸收的影响，随着传播路程的增加，反射能量逐渐变弱；另外，受激发条件和接收条件等因素的影响，记录到的各炮的能量也存在一定差别。这些变化与地下地质信息无关，容易使第三阶段的地震资料解释陷入误区。为了使地震资料的能量变化真实地反映出地下储层的构造变化，有必要对地震资料进行能量均衡处理。

目前，主要是用球面扩散补偿和自动增益控制来进行能量均衡处理。具体的，球面扩散补偿是补偿地震波向下传播过程中由于球面扩散而造成的能量衰减，使浅层、中层和深层能量得到均衡。地震波传播是球面扩散波。离震源越远，能量就越弱，这样浅层反射回来的地震波，就强于深层反射回来的地震波。最基本的球面扩散补偿方式就是以距离为加权值进行加权，而地震记录是旅行时-振幅值，所以用旅行时乘以地震波速度作为加权值，乘以该旅行时振幅，这样得到的振幅是浅层加权值小、深层加权值大，这样就补偿了深层的能量。这是一种单道处理方式。自动增益控制(AGC)是一种算法，常用于地震记录的显示，这里用来实现能量均衡，也是一种单道处理方式。对于每一道，计算一段时窗内振幅能量，然后用这个能量的倒数作为加权值，乘以该时窗内每个样点的振幅值。由于深层的能量弱，那么能量的倒数就大，于是深层的加权值就大，于是补偿了深层的能量。

经过球面扩散补偿，主要表现是恢复了单道的中深层振幅，但是没有考虑横向空间的球面扩散影响。实际上球面扩散是在空间上，那么单道处理的球面扩散补偿从一维角度考虑，不涉及到道与道之间的均衡，就破坏了波场空间的相对关系，尤其是破坏了波场在空间方向上相对强弱关系。同样自动增益控制也是单道处理，对于一个炮集，不区分远道和近道，同样进行补偿，也是从一维角度考虑，而没有从空间上考虑波场的相对关系，所以就破坏了波场的相对关系，不但破坏了波场的空间强弱关系，最明显的是破坏了波场在时间方向上的强弱关系。

由于以叠前波动方程偏移为核心的数据处理主要是利用整个地震波场信息来实现复杂构造成像，因此在进行地震数据处理过程中都需要尽量保持波场本身在时间方向和空间方向上的相对强弱关系，上述的波场本身在时间方向和空间方向上的相对强弱关系称为波场相对关系。实际地震记录的资料本身是具有波场相对关系的，利用上述的球面扩散补偿和自动增益控制进行能量均衡时，球面扩散补偿破坏了波场在空间方向上相对强弱关系，自动增益控制破坏了波场在时间方向上相对强弱关系，不能实现在进行能量均衡过程中保持地震波场的相对关系。

发明内容

本发明的目的是提供一种能量均衡方法及处理装置，以实现在进行能量均衡处理后保持地震波场特征，为所有基于波场的处理技术提供必备的地震数据。

本发明提供一种地震波场能量均衡方法，包括：

S1：处理装置读入经过预处理的地震数据；

S2：处理装置对读入的地震数据，用异常振幅消除技术逐炮进行单炮能量均衡处理；

S3：处理装置对单炮能量均衡处理后的地震数据，进行多炮集记录的炮间能量均衡处理；

其中：

所述S2包括：

S21：输入炮集中一组地震道x(i，j)，计算复数地震道的虚部y(i，j)；

S22：根据地震道x(i，j)和复数地震道虚部y(i，j)，计算地震道的包络；

S23：计算地震道包络A(i，j)的中值M(i)；

S24：根据中值M(i)，设计检测异常振幅的准则；

S25：对经过S24检测的地震道数据应用权值，衰减异常振幅的能量；

S26：依次移动道窗，并重复执行S23至S25，直至处理完炮集中的所有地震道数据；

所述S3包括：

S31：计算每个炮集的总能量；

S32：根据S31计算出的炮集总能量，对每个炮集的能量进行归一化。

所述S23中的中值M（i)采用alpha-trim算法计算均值替代，具体计算方法如下：将样点值从小到大进行排序，然后放弃前p％个样点和后p％个样点，将剩下的样点值求和取平均，作为要求的M(i)；其中，p是可选参数，p的取值范围为0-30。

所述S25中衰减异常振幅的能量，具体通过下式进行：

$\hat{x}(i,j)=x(i,j)/l(i,j),(i=\mathrm{0,1,2}...,\mathrm{nt}-1,j=\mathrm{tracew}/2)$

其中，x(i，j)是地震道；i是时间采样点；nt是道长，nt为定值；j是空间采样点；tracew是窗口内道数，tracew取奇数；

是异常振幅消除后的地震道；l(i，j)是权值，权值l(i，j)计算方法如下：

$l(i,j)=\left\{\begin{array}{cc}\frac{A(i,j)}{M\left(i\right)\·\mathrm{\α}}& A(i,j)>\mathrm{threshold}\·M\left(i\right)\\ 1& A(i,j)\≤\mathrm{threshold}\·M\left(i\right)\end{array}\right.$

其中，threshold是门槛值，threshold>1；α是衰减系数，0<α≤1。

所述S31中计算每个炮集的总能量是计算每个炮集的均方根总能量，或者是通过计算炮集的包络计算炮集总能量，或者是用绝对值方式来计算总能量，或者是用P模方式来计算总能量。

所述S21中计算复数地震道的虚部是根据希尔伯特变换进行计算的。

本发明一种地震波场能量均衡处理装置，包括：地震数据读入单元、单炮能量均衡处理单元、炮间能量均衡处理单元。其中：

地震数据读入单元，用于读入经过预处理的地震数据；

单炮能量均衡处理单元，用于对读入的地震数据，逐炮进行异常振幅消除，实现单炮能量均衡处理；

炮间能量均衡处理单元，用于对单炮能量均衡后的地震数据进行炮间能量均衡；

所述单炮能量均衡处理单元，包括：地震道虚部计算单元、地震道包络计算单元、包络中值计算单元、异常振幅检测单元、异常振幅能量衰减单元；其中：地震道虚部计算单元，用于根据读入的地震数据，计算复数地震道的虚部；地震道包络计算单元，用于根据读入的地震数据和计算得到的复数地震道虚部，计算地震道的包络；包络中值计算单元，用于根据计算所得的地震道包络计算地震道包络的中值；异常振幅检测单元，用于根据地震道包络的中值判别地震数据是否属于异常振幅；异常振幅能量衰减单元，用于对判别出的属于异常振幅的地震数据进行能量衰减；

所述炮间能量均衡处理单元可以包括：能量计算单元和能量归一化单元；其中：所述能量计算单元，用于计算每个炮集的总能量；所述能量归一化单元，用于根据能量计算单元计算出的炮集总能量对每一炮的能量进行归一化。

所述包络中值计算单元替换为alpha-trim算法均值计算单元；所述alpha-trim算法均值计算单元，用于计算包络中的均值，并将该均值作为M(i)的值。

所述能量计算单元，具体是用于计算每个炮集的均方根总能量，或者是用于通过计算炮集的包络计算炮集总能量，或者是用于通过绝对值方式来计算总能量，或者是用于通过P模方式来计算总能量。

本发明提供的一种地震波场能量均衡方法，首先用异常振幅消除技术进行单炮集中能量均衡，从而消除了接收因素引起的振幅差异；然后再对多炮集记录采用炮间能量均衡技术，消除了由激发因素引起的能量不均衡。采用的上述两种技术都是从能量角度来消除这两个环节的差异，没有涉及到地震波场的相对关系，在进行地震数据能量均衡后，得到的结果是浅层能量强于深层能量、近道能量强于远道能量，保持了波场的相对关系。本发明还提供一种地震波场能量均衡处理装置，所述的地震波场能量均衡处理装置，可以实现保持地震波场特征的能量均衡处理。

附图说明

图1是本发明地震波场能量均衡方法的流程图；

图2是本发明地震波场能量均衡处理装置实施例；

图3是本发明能量均衡处理装置实施例中的单炮能量均衡处理单元；

图4是本发明能量均衡处理装置实施例中的炮间能量均衡处理单元；

图5是米桩号的示意图；

图6是经过预处理后的砾石区共炮集记录；

图7是经过预处理后的山区共炮集记录；

图8是图6经过S2处理后的地震数据；

图9是图7经过S2处理后的地震数据；

图10是各单炮经过能量均衡后的多炮集记录；

图11是图10经过炮间能量均衡处理后的多炮集记录。

具体实施方式

为了使本技术领域的人员更好地理解本申请中的技术方案，下面将结合本申请实施例中的附图，对本申请实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本申请一部分实施例，而不是全部的实施例。基于本申请中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都应当属于本发明保护的范围。

在理想情况下的地震记录，比如用计算机模拟或用物理模拟方式得到的地震记录，都是能量均衡的，因为模拟时认为接收因素相同、激发因素相同。而实际地震资料能量不均衡就是来源于接收和激发这两个环节的差异。只要消除了这两个环节产生的差异，就能实现能量均衡。

图1所示为本发明地震波场能量均衡方法的流程图。下面根据图1所示的流程图来具体说明本发明的第一实施例。如图1所示，本发明地震波场能量均衡方法包括：

S1：处理装置读入经过预处理的地震数据；

S2：处理装置对输入的地震数据，用异常振幅消除技术逐炮进行单炮能量均衡；

S3：处理装置对单炮能量均衡处理后的地震数据，进行多炮记录的炮间能量均衡处理。

其中，

所述S2，具体包括：

S21：输入炮集中一组地震道x(i，j)，计算复数地震道的虚部y(i，j)；

S22：根据地震道x(i，j)和复数地震道虚部y(i，j)，计算地震道的包络；

S23：计算地震道包络A(i，j)的中值M(i)；

S24：根据中值M（i)，设计检测异常振幅的准则；

S25：对经过S24检测的地震道数据应用权值，衰减异常振幅的能量；

S26：依次移动道窗，并重复执行S23至S25，直至处理完炮集中的所有地震道数据。

所述S3，具体包括：

S31：计算每个炮集的总能量；

S32：根据S31计算出的炮集总能量，对每个炮集的能量进行归一化。

图2所示为本发明地震波场能量均衡处理装置实施例。如图2所示，所述能量均衡处理装置，包括：地震数据读入单元1、单炮能量均衡处理单元2、炮间能量均衡处理单元3。其中：

地震数据读入单元1，用于读入经过预处理的地震数据；

单炮能量均衡处理单元2，用于对读入的地震数据，逐炮进行异常振幅消除，实现单炮能量均衡；

炮间能量均衡处理单元3，用于对单炮能量均衡后的地震数据进行炮间能量均衡。

图3所示为本发明地震波场能量均衡处理装置实施例中的单炮能量均衡处理单元2，如图3所示，所述单炮能量均衡处理单元2可以包括：地震道虚部计算单元21、地震道包络计算单元22、包络中值计算单元23、异常振幅检测单元24、异常振幅能量衰减单元25。所述地震道虚部计算单元21用于根据读入的地震数据，计算复数地震道的虚部；所述地震道包络计算单元22用于根据读入的地震数据和计算得到的复数地震道虚部，计算地震道的包络；所述包络中值计算单元23用于根据计算所得的地震道包络计算地震道包络的中值；所述异常振幅检测单元24用于根据地震道包络的中值判别地震数据是否属于异常振幅；所述异常振幅能量衰减单元25用于对判别出的属于异常振幅的地震数据进行能量衰减。

图4所示为本发明地震波场能量均衡处理装置实施例中的炮间能量均衡处理单元3，如图4所示，所述炮间能量均衡处理单元3可以包括：能量计算单元31和能量归一化单元32。所述能量计算单元31可以用于计算每个炮集的总能量；所述能量归一化单元32用于根据计算出的炮集总能量对每一炮的能量进行归一化。

下面利用图2所示的本发明地震波场能量均衡处理装置实施例，按照图1所示的流程对实际施工中采集到的地震数据进行处理。

S1：处理装置读入经过预处理的地震数据。

在地质工作和其他物探工作初步确定的油气勘探有利区带上，按照二维或三维观测系统布置测线和炮点位置，使用炸药震源或可控震源激发地震波，并通过检波器和地震仪以时间离散采样方式把地震波场记录下来。用上述手段得到的就是野外采集的地震数据。在上述野外采集的地震数据中，每个检波点上记录的地震波数据称为地震道。一个检波点上的地震道为一个单道。地震记录时每激发一次地震(也称为单炮)，一般有成千上万个检波器在接收信号，接收的信号的集合叫做单炮集记录。检波器记录的多次地震的地震道的集合称为多炮集记录。多炮集记录是由多个单炮集记录在同一张地震数据上并列显示。在野外采集时，需要记录炮点和接收点的大地坐标和高程信息，并且对每炮和每道都进行编号，炮点的编号叫做炮号，接收点编号叫做道号。

上述野外采集的地震数据一般为时序地震记录。上述地震记录格式由地震道道头字和地震记录数据组成，地震道道头字是指每道地震记录中前240个字符，地震道道头字中记录了该地震道的所有信息，包括炮号、道号、大地坐标、高程等信息。同一炮所有接收道所记录的数据集合形成一个炮集。将每个炮集中按照道号由小到大排序重新进行排列，将时序地震记录重排形成道序地震记录。由于施工时用大地坐标表示地震测线，实际的地震测线往往不是正南正北，为了处理方便，通常将大地坐标换算为以米为单位的桩号(即米桩号)，这样可以把测线旋转为正南或正东方向，便于数据处理。图5所示为米桩号的示意图，如图5所示，横坐标表示米桩号，纵坐标表示高程。

在实际施工过程中，有时会因为激发问题，导致整炮记录中不包含地下信息，这种记录就是废炮记录；有时因为接收问题，导致该道记录无效，这道记录就是废道记录。处理前必须去掉这些无效的记录。

另外，还需要用陷波技术消除50HZ的单频工业干扰，并切除初至波之前的背景波场。

处理装置中的地震数据读入单元1读入经过上述预处理的地震数据。

图6所示为经过预处理后的砾石区单炮集记录，横坐标表示道号，纵坐标表示时间，灰度表示振幅。由图6可以看出，该炮集记录存在单道异常振幅和倾斜分布的异常振幅。图7所示为经过预处理后的山区单炮集记录，横坐标表示道号，纵坐标表示时间，灰度表示振幅。图7中，该炮集记录存在各道之间能量不均匀，有些地震道能量强、有些地震道能量弱的现象。

S2：处理装置对输入的地震数据，用异常振幅消除技术逐炮进行单炮能量均衡。

处理装置中的单炮能量均衡处理单元2，对读入的地震数据，逐炮进行异常振幅消除，实现单炮能量均衡。

由于不均匀地表的影响和检波器接收因素的差异，实际地震炮集记录中存在单道异常振幅及各道之间能量不均匀的情况，出现有些地震道能量强、有些地震道能量弱的现象，这些现象与地下构造信息无关，仅仅是由于近地表因素造成的，所以必须消除这些异常振幅。

复数地震道的包络为x(i，j)+j□y(i，j)，复数的模这是一个实数。其中，x(i，j)代表质点位移的动量，y(i，j)代表质点位移的势能，则复数的模代表质点位移的能量，所以这里可以用包络来衡量该质点的能量大小。如果包络的值太大，就表示该点能量是异常能量，相应的振幅为异常振幅。对于异常振幅，可以利用异常振幅消除技术进行单炮的能量进行均衡处理。具体包括：

S21：输入炮集中一组地震道x(i，j)，计算其复数地震道的虚部y(i，j)。

地震道虚部计算单元21根据读入的地震数据，计算复数地震道的虚部。具体地，先输入炮集中一组地震道x(i，j)，其中，i为时间采样点(i＝0，1，2...，nt-1，nt为道长)，道长是指每道记录的时间长度，常用采样点数来表示，采样点数乘以采样间隔，就等于记录的时间长度；j为空间采样点(j=1，2，..，tracew，tracew为窗口内道数)。再根据输入的地震道x(i，j)，利用希尔伯特(Hilbert)变换计算其复数地震道虚部y(t，j)。

如果把炮集中一组地震道x(i，j)的第j道看成实值函数f(t)的离散值，那么这个实值函数f(t)的希尔伯特变换表示为：

$\hat{f}\left(t\right)={\&Integral;}_{-\&infin;}^{+\&infin;}\frac{f\left(\mathrm{\&tau;}\right)}{\mathrm{\&pi;}(t-\mathrm{\&tau;})}\mathrm{d\&tau;}$

这是函数f(t)与1／tπ的褶积，可以表示为

$\hat{f}\left(t\right)=f\left(t\right)*\frac{1}{\mathrm{t\&pi;}}$

对1／tπ进行傅里叶变换，得

$F[1/\mathrm{t\&pi;}]=-\mathrm{jsgn}\left(\mathrm{\&omega;}\right)=e^{-j\frac{\mathrm{\&pi;}}{2}\mathrm{sgn}\left(\mathrm{\&omega;}\right)}$

其中，ω表示圆频率。

所以希尔伯特变换等效于±π／2相移，对正圆频率产生-π／2相移，对负圆频率产生π／2相移。或者说，对时域信号中每一频率成分移位1／4波长。因此希尔伯特变换又称为90度移相器。所以计算复数地震道虚部y(i，j)有两种方法。

一是用时间域褶积方法，即 $y(i,j)=x(i,j)*\frac{1}{\mathrm{\&pi;i}}=\frac{1}{\mathrm{\&pi;}}\underset{k=1}{\overset{n}{\mathrm{\&Sigma;}}}\frac{x(k,j)}{i-k},$ 其中，n为褶积算子长度，n一般取奇数，如59个点；

另一种方法是在频域内实施，即将每道记录进行傅里叶变换，将变换到频域的数据乘以-jsgn(ω)，其中

即对正圆频率数据的实虚部进行交换，然后将所得的频域数据乘以(-1)，再反傅里叶变换到时间域，即可得到复数地震道虚部y(i，j)。

S22：根据地震道x(i，j)和复数地震道虚部y(i，j)，计算地震道的包络；。

地震道包络计算单元22根据读入的地震数据和计算得到的复数地震道虚部，计算地震道的包络。具体地，在S21中，用地震数据作为实部，代替地震响应中质点位移的动量，那么经过希尔伯特变换后，进行了90度相移，即得到了正交道，该道就相当于质点位移的势能。这样就构成了复数地震道x(i，j)+j□y(i，j)，由此就可以计算地震包络，计算公式如下：

$A(i,j)=\sqrt{x^2(i,j)+y^2(i,j)}$

S23：计算地震道包络A(i，j)的中值M(i)。

包络中值计算单元23根据计算所得的地震道包络计算地震道包络的中值。在消除异常振幅之前需要先找到地震数据中的异常振幅。中值是指在一组数据中居于中间的数(这组数是经过升序排列的)，即这组数据中有一半数据比它大，有一半数据比它小。所以，一般中值能够用来衡量一组数中大部分元素的平均值。

实际地震数据处理工作中，在道数为tracew的道窗中沿道计算A(i，j)的中值M(i)。即将S22求得的地震数据的包络作为一组样点值，对这一组样点值从小到大进行排序，当样点个数为奇数时，取居于中间位置上的样点值作为中值；当样点个数为偶数时，取中间位置上两个样点值的均值作为中值。

S24：根据中值M(i)，设计检测异常振幅的准则。

异常振幅检测单元24根据地震道包络的中值判别地震数据是否属于异常振幅。

地震记录是一种波形，其振幅由正数、零和负数组成。其中，振幅特别大和振幅特别小的值都属于异常振幅。由于S23中将振幅值变成了包络，而包络只有正数和零，那么负数也变成了正数。因此特别大的包络对应于特别大的正数振幅或特别小的负数振幅。可以将特别大的包络值作为异常振幅。

在实际地震数据处理工作中，设计门槛值threshold，若A(i，j)>threshold·M(i)，则认为该包络对应的振幅为异常振幅；否则为正常振幅。其中，一般threshold>1，threshold的值具体可以根据实际资料试验来确定，例如取1.3；M(i)为S23计算所得。

S25：对经过S24检测的地震道数据应用权值，衰减异常振幅的能量。

异常振幅能量衰减单元25对判别出的属于异常振幅的地震数据进行能量衰减。具体地，根据S24中对异常振幅的判别，通过下式对窗内的中间道应用权值l(i，j)，衰减强能量：

$\hat{x}(i,j)=x(i,j)/l(i,j),(i=\mathrm{0,1,2}...,\mathrm{nt}-1,j=\mathrm{tracew}/2)$

其中，tracew是窗口内道数，一般取奇数，例如11；

是异常振幅消除后的地震道；l(i，j)是权值，权值l(i，j)计算方法如下：

$l(i,j)=\left\{\begin{array}{cc}\frac{A(i,j)}{M\left(i\right)\&CenterDot;\mathrm{\&alpha;}}& A(i,j)>\mathrm{threshold}\&CenterDot;M\left(i\right)\\ 1& A(i,j)\&le;\mathrm{threshold}\&CenterDot;M\left(i\right)\end{array}\right.$

其中，threshold为门槛值；α为衰减系数，0<α≤1。α越大，获得的权值就越小，反之权值就越大。α取1时，得到的衰减后的振幅值等于原来的振幅值

倍。对于正常振幅值，权值为1，衰减后振幅不变。α可以根据实际试验来确定，例如为0.8。

S26：依次移动道窗，并重复执行S23至S25，直至处理完共炮集中的所有地震道数据。

依次移动道窗，每一次移动一道，并重复执行S23至S25。对于第一个窗口中间道左边的道和最后一个窗口中间道右边的道，由于无法应用权值，就用当前窗口中各点的中值来代替包络值。

图8所示为图6经过S2处理后的地震数据。图9所示为图7经过S2处理后的地震数据。由图8和图9可以看出，经过S2处理后的地震数据，异常振幅得到了消除，能量更加均匀，同时保持了波场的特征信息。

S3：处理装置对单炮能量均衡处理后的地震数据，进行多炮记录的炮间能量均衡处理。

处理装置中的炮间能量均衡处理单元3，对单炮能量均衡后的地震数据进行炮间能量均衡。

由于激发环境和激发能量的差异，导致不同的炮有不同的能量。因此，在S2对多炮记录中的每一个单炮记录进行能量均衡处理后，还必须进行炮间能量均衡，从而消除激发能量和激发环境差异导致的各炮集中地震波场能量之间的差异。炮间能量均衡技术具体包括：

S31：计算每个炮集的总能量。

能量计算单元31计算每个炮集的总能量。具体地，计算每个炮集的均方根总能量A_rms(k)，计算公式如下：

$A_{\mathrm{rms}}\left(k\right)=\frac{1}{\mathrm{nt}\&times;\mathrm{nr}\left(k\right)}\sqrt{\underset{j=1}{\overset{\mathrm{nr}\left(k\right)}{\mathrm{\&Sigma;}}}\underset{i=1}{\overset{\mathrm{nt}}{\mathrm{\&Sigma;}}}{A}^2(i,j,k)},k=\mathrm{1,2},...,\mathrm{ns}$

其中A(i，j，k)是第k炮、第j道、第i个采样点的振幅值，ns为地震数据中总炮数，nr(k)是第k炮炮集的道数，nt是每道时间采样点个数，上述各参量均为已知量，且这里假设nt保持不变。

S32：根据S31计算出的炮集总能量，对每个炮集的能量进行归一化。

能量归一化单元32根据计算出的炮集总能量对每一炮的能量进行归一化。用归一化方式，就是将各炮的均方根总能量都化为1，这样各炮的总能量就被调整到同一级别，则可以实现炮间能量均衡。

能量归一化计算可以根据下式实现：

其中，A(i，j，k)是第k炮、第j道、第i个采样点的振幅值，为已知量；A_rms(k)是炮集的均方根总能量，由S31计算得到；

就是能量均衡后的多炮集记录。

图10所示为各单炮经过能量均衡后的多炮集记录，横坐标表示道号，纵坐标表示时间，灰度表示振幅。由图10可以看出，在2881～3361道位置上的这一炮记录能量明显弱于其他炮能量，说明不同炮的能量有强有弱。图11所示为图10经过炮间能量均衡处理后的多炮集记录，由图11可以看出，进行炮间能量均衡处理后，各炮的能量趋于一致。

上述本发明第一实施例中，采用异常振幅消除技术消除了接收因素引起的能量不均衡，采用炮间能量均衡技术消除了激发因素引起的能量不均衡。所述两种技术都是从能量角度来消除差异，没有涉及到地震波场的相对关系，得到的结果是浅层能量就应该强于深层能量、近道能量就该强于远道能量。因此本发明的技术方案在进行能量均衡的同时，能够保持波场的相对关系。

需要说明的是，上述实施例S31中，计算炮集总能量也可以是通过计算炮集的包络计算炮集总能量，或用绝对值方式来计算总能量，或用P模方式来计算总能量。可以根据计算要求选取，本发明对此并不作出限定。

以下介绍本发明的第二实施例。

本发明的第二实施例与第一实施例的区别在于，第一实施例S23中的M(i)可以采用alpha-trim算法计算均值来获得。具体实施时，包括：将S22所得的包络作为样点值，对这一组样点值从小到大进行排序，然后放弃前p％个样点和后p％个样点，然后将剩下的样点值求和取平均，作为要求的M(i)。其中，p为可选参数，一般p的取值范围为0-30，可以根据实际资料试验确定，例如取15。

本发明的第二实施例与第一实施例的区别还在于，本发明保持地震波场特征的能量均衡处理装置实施例的单炮能量均衡处理单元2中，包络中值计算单元23可以替换为alpha-trim算法均值计算单元。所述alpha-trim算法均值计算单元可以用于计算包络中的均值，并将该均值作为M(i)。

本发明第二实施例的其他部分与第一实施例完全相同。

本实施例中求得的均值为去除了前p％个样点和后p％个样点后的剩余样点值的均值。

本发明提供的能量均衡方法，首先用异常振幅消除技术进行单炮集中能量均衡，从而消除了接收因素引起的振幅差异；然后再对多炮集记录采用炮间能量均衡技术，消除了由激发因素引起的能量不均衡。两种技术都是从能量角度来消除这两个环节的差异，没有涉及到地震波场的相对关系，在进行地震数据能量均衡后，得到的结果是浅层能量强于深层能量、近道能量强于远道能量，保持了波场的相对关系。

本发明提供的地震波场能量均衡处理装置，可以实现保持地震波场特征的能量均衡处理。

需要说明的是计算M(i)的两种方法，即中值计算方法和alpha-trim算法均值计算方法，可以根据需要进行选取，本发明对此并不作出限定。

Claims (8)

1.一种地震波场能量均衡方法，其特征在于，包括：

S1：处理装置读入经过预处理的地震数据；

S2：处理装置对读入的地震数据，用异常振幅消除技术逐炮进行单炮能量均衡处理；

S3：处理装置对单炮能量均衡处理后的地震数据，进行多炮集记录的炮间能量均衡处理；

其中，

所述S2包括：

S21：输入炮集中一组地震道x(i，j)，计算复数地震道的虚部y(i，j)；

S22：根据地震道x(i，j)和复数地震道虚部y(i，j)，计算地震道的包络；

S23：计算地震道包络A(i，j)的中值M(i)；

S24：根据中值M(i)，设计检测异常振幅的准则；

S25：对经过S24检测的地震道数据应用权值，衰减异常振幅的能量；

S26：依次移动道窗，并重复执行S23至S25，直至处理完炮集中的所有地震道数据；

所述S3包括：

S31：计算每个炮集的总能量；

S32：根据S31计算出的炮集总能量，对每个炮集的能量进行归一化。

2.如权利要求1所述的地震波场能量均衡方法，其特征在于，所述S23中的中值M(i)采用alpha-trim算法计算均值替代，具体计算方法如下：

将样点值从小到大进行排序，然后放弃前p％个样点和后p％个样点，将剩下的样点值求和取平均，作为要求的M(i)；其中，p是可选参数，p的取值范围为0-30。

3.如权利要求1所述的地震波场能量均衡方法，其特征在于，所述S25中衰减异常振幅的能量，具体通过下式进行：

$\hat{x}(i,j)=x(i,j)/l(i,j),(i=\mathrm{0,1,2}...,\mathrm{nt}-1,j=\mathrm{tracew}/2)$

其中，x(i，j)是地震道；i是时间采样点；nt是道长，nt为定值；j是空间采样点；tracew是窗口内道数，tracew取奇数；

是异常振幅消除后的地震道；l(i，j)是权值，权值l(i，j)计算方法如下：

$l(i,j)=\left\{\begin{array}{cc}\frac{A(i,j)}{M\left(i\right)\&CenterDot;\mathrm{\&alpha;}}& A(i,j)>\mathrm{threshold}\&CenterDot;M\left(i\right)\\ 1& A(i,j)\&le;\mathrm{threshold}\&CenterDot;M\left(i\right)\end{array}\right.$

其中，threshold是门槛值，threshold>1；α是衰减系数，0<α≤1。

4.如权利要求1所述的地震波场能量均衡方法，其特征在于，所述S31中计算每个炮集的总能量是计算每个炮集的均方根总能量，或者是通过计算炮集的包络计算炮集总能量，或者是用绝对值方式来计算总能量，或者是用P模方式来计算总能量。

5.如权利要求1所述的地震波场能量均衡方法，其特征在于，所述S21中计算复数地震道的虚部是根据希尔伯特变换进行计算的。

6.一种地震波场能量均衡处理装置，其特征在于，包括：地震数据读入单元、单炮能量均衡处理单元、炮间能量均衡处理单元；其中：

地震数据读入单元，用于读入经过预处理的地震数据；

单炮能量均衡处理单元，用于对读入的地震数据，逐炮进行异常振幅消除，实现单炮能量均衡处理；

炮间能量均衡处理单元，用于对单炮能量均衡后的地震数据进行炮间能量均衡；

所述单炮能量均衡处理单元，包括：地震道虚部计算单元、地震道包络计算单元、包络中值计算单元、异常振幅检测单元、异常振幅能量衰减单元；其中，

地震道虚部计算单元，用于根据读入的地震数据，计算复数地震道的虚部；

地震道包络计算单元，用于根据读入的地震数据和计算得到的复数地震道虚部，计算地震道的包络；

包络中值计算单元，用于根据计算所得的地震道包络计算地震道包络的中值；

异常振幅检测单元，用于根据地震道包络的中值判别地震数据是否属于异常振幅；

异常振幅能量衰减单元，用于对判别出的属于异常振幅的地震数据进行能量衰减；

所述炮间能量均衡处理单元可以包括：能量计算单元和能量归一化单元；其中，

所述能量计算单元，用于计算每个炮集的总能量；

所述能量归一化单元，用于根据能量计算单元计算出的炮集总能量对每一炮的能量进行归一化。

7.如权利要求6所述的地震波场能量均衡处理装置，其特征在于，所述包络中值计算单元替换为alpha-trim算法均值计算单元；所述alpha-trim算法均值计算单元，用于计算包络中的均值，并将该均值作为M(i)的值。

8.如权利要求6所述的地震波场能量均衡处理装置，其特征在于，所述能量计算单元，具体是用于计算每个炮集的均方根总能量，或者是用于通过计算炮集的包络计算炮集总能量，或者是用于通过绝对值方式来计算总能量，或者是用于通过P模方式来计算总能量。

Images