一种基于小生境差分进化算法的中值滤波方法

摘要

本发明提供一种基于小生境差分进化算法的中值滤波方法，其特征在于，包括步骤：将待处理数据划分到一系列重叠的时空窗中；获取单个所述时空窗中的同相轴能量曲线，根据所述同相轴能量曲线的局部极大值获得同相轴位置；在所述同相轴位置处，通过所述小生境差分进化算法获得同相轴倾角及其相干值；根据所述同相轴倾角及其相干值筛选所述同相轴倾角；沿着筛选出的倾角依次进行中值滤波。本发明能同时获得所有交叉同相轴倾角的能力和真三维的特性使其能够获得更好的滤波效果。

权利要求

1.一种基于小生境差分进化算法的中值滤波方法,其特征在于,包括步骤:S1将待处理数据划分到一系列重叠的时空窗中；

S2获取单个所述时空窗中的同相轴能量曲线,根据所述同相轴能量曲线的局部极大值获得同相轴位置；

S3在所述同相轴位置处,通过所述小生境差分进化算法获得同相轴倾角及其相干值；

S4根据所述同相轴倾角及其相干值筛选所述同相轴倾角；

S5沿着筛选出的倾角依次进行中值滤波；

其中,所述小生境差分进化算法包括:基于邻域拥挤的差分进化算法、hill‑valley探测函数和Nelder‑Mead单纯形算法；

所述基于邻域拥挤的差分进化算法用于驱使局部解朝着各自的局部最优移动；

所述hill‑valley探测函数用于从所述局部解中获取局部优秀解；

所述Nelder‑Mead单纯形算法用于提高所述局部优秀解的精度；

其中,所述基于邻域拥挤的差分进化算法包括步骤:S31随机生成N个初始个体；

S32遍历每一个所述初始个体；

寻找参数上最相似的m个个体形成一个亚种群；

S33在所述亚种群中使用所述差分进化算法的变异和交叉操作生,成一个试验子个体；

S34在整个种群中计算所述试验子个体与其他所述个体的欧氏距离,S35比较所述试验子个体与同所述试验子个体欧式距离最近个体的相干值,如果所述试验子个体的相干值更小,则用所述试验子个体替换掉所述欧式距离最近的个体；

S36重复M次步骤S32‑S35后算法终止,所述M为最大迭代次数。

2.根据权利要求1所述的一种基于小生境差分进化算法的中值滤波方法,其特征在于,还包括步骤:S6储存所述时空窗中的中值滤波数据。

3.根据权利要求2所述的一种基于小生境差分进化算法的中值滤波方法,其特征在于,还包括步骤:S7确定各个所述时空窗中的待分析数据是否已被处理；

若否,则返回执行步骤S2及步骤S2之后的步骤；

若是,则累加所有时空窗中的所述中值滤波数据。

4.根据权利要求1所述的一种基于小生境差分进化算法的中值滤波方法,其特征在于,所述hill‑valley探测函数用来确定两个所述局部解是否属于同一个小生境,包括步骤:S37所述hill‑valley探测函数在两个所述局部解之间采样,获取所述局部解之间的内部倾角样本的相干值；

S38若任何所述内部倾角样本的相干值大于两个所述局部解的相干值,则确定两个所述局部解属于不同的小生境；

否则,确定两个所述局部解属于同一个小生境；

S39从所述局部解中获取局部优秀解。

说明书

一种基于小生境差分进化算法的中值滤波方法

技术领域

[0001] 本发明涉及地震数据处理领域，尤其涉及一种基于小生境差分进 化算法的中值滤波方法。

背景技术

[0002] 在地震勘探过程中，采集到的地震数据通常含有规则噪声和不规 则噪声，其中不规则噪声即随机噪声。地震数据中常含有许多随机噪 声，如微震、背景干扰等，这些噪声分布很广，会严重影响了地震数 据的信噪比，目前，中值滤波是一种普遍认为较快速高效的随机噪声 压制方法。

[0003] 对于地震数据来说，在应用中值滤波之前需要把同相轴拉平，现 有技术中，拉平操作通过在局部时空窗中估计同相轴的倾角来实现， 即倾角导向的中值滤波方法。在简单情况下，一个数据位置只对应一 个倾角值。但是在一些复杂的地质环境下，如不整合面或断层，两个 或者更多的同相轴会相交，在一个数据位置处会产生多个倾角值。

[0004] 现有应用中用于求取同相轴倾角的方法包括局部倾斜叠加法，平 面波摧毁滤波器(PWD)方法和迭代f‑k方法。局部倾斜叠加法是在 数据点处按照相干性计算公式计算出沿不同方向的相干值，然后挑选 出最大相干值所对应的倾角进行倾角导向中值滤波。PWD方法首先 将隐式有限差分格式应用到局部平面波微分方程，然后用无限脉冲响 应滤波器来近似相移算子将局部平面波微分方程变成倾角的非线性 方程，进而通过迭代反演算法来计算同相轴的倾角。迭代f‑k方法在 一次迭代中首先将单个时空窗中的地震数据变换到径向傅里叶域，然 后沿着径向方向叠加获得峰值图，最终根据径向傅里叶域中选出的最 大峰值获得时空域中的倾角。

[0005] 目前，在简单情况下求取同相轴倾角的方法已经很常见，对于复 杂情况下交叉同相轴倾角的计算，还没有很好的方法。

发明内容

[0006] 本发明提供一种基于小生境差分进化算法的中值滤波方法，用以 解决上述问题，具体的本发明提供一种基于小生境差分进化算法的中 值滤波方法，其特征在于，包括步骤：S1将待处理数据划分到一系 列重叠的时空窗中；S2获取单个所述时空窗中的同相轴能量曲线， 根据所述同相轴能量曲线的局部极大值获得同相轴位置；S3在所述 同相轴位置处，通过所述小生境差分进化算法获得同相轴倾角及其相 干值；S4根据所述同相轴倾角及其相干值筛选所述同相轴倾角；S5 沿着筛选出的倾角依次进行中值滤波。

[0007] 优选的，还包括步骤：S6储存所述时空窗中的中值滤波数据。

[0008] 优选的，还包括步骤：S7确定所述各个时空窗中的待分析数据 是否已被处理；若是，则返回执行步骤S2及步骤S2之后的步骤；若 否，则累加所有时空窗中的所述中值滤波数据。

[0009] 优选的，步骤S2中，获取单个所述时空窗中的同相轴能量曲线， 包括：应用公式确定所述同相轴能量曲线的局部极大值， 根据所述同相轴能量曲线的局部极大值来识别同相轴的位置；其中，若 否则，  E
为所述时空窗中的同相轴能量，t表示计算时间位置，H表示所述 时空窗的空间宽度，d表示待分析数据，T为能量计算范围的半宽， φ为大于等于1的常数，x0的值为H/4、H/2或3H/4，y0的值为H/4、H/2 或3H/4；q为大于等于1的常数。

[0010] 优选的，所述小生境差分进化算法包括：基于邻域拥挤的差分进 化算法、hill‑valley探测函数和Nelder‑Mead单纯形算法；所述基于 邻域拥挤的差分进化算法用于驱使局部解朝着各自的局部最优移动； 所述hill‑valley探测函数用于从所述局部解中获取局部优秀解；所述 Nelder‑Mead单纯形算法用于提高所述局部优秀解的精度；其中，所 述基于邻域拥挤的差分进化算法包括步骤：S31：随机生成N个初始 个体；S32：遍历每一个所述初始个体；寻找参数上最相似的m个个 体形成一个亚种群；S33在所述亚种群中使用所述差分进化算法的变 异和交叉操作生成一个试验子个体；S34在整个种群中计算所述试验 子个体与其他所述个体的欧氏距离；S35比较所述试验子个体与同所 述试验子个体欧式距离最近个体的相干值，如果所述试验子个体的相 干值更小，则用所述试验子个体替换掉所述欧式距离最近的个体；S36 重复M次步骤S32‑S35后算法终止，所述M为最大迭代次数。

[0011] 优选的，所述hill‑valley探测函数用来确定两个所述局部解是否 属于同一个小生境，包括步骤：S37所述hill‑valley探测函数在两个 所述局部解之间采样，获取所述局部解之间的内部倾角样本的相干值； S38若任何所述内部倾角样本的相干值大于两个所述局部解的相干值， 则确定两个所述局部解属于不同的小生境；否则，确定两个所述局部 解属于同一个小生境；S39从所述局部解中获取局部优秀解。

[0012] 优选的，所述Nelder‑Mead单纯形算法用于提高所述局部优秀解 的精度，包括步骤：。S310围绕单个局部优秀解构建具有n+1个顶点 的单纯形，n为解的维数，排序并对n+1个顶点重新标签为x1、x2、…、 xn+1，使f(x1)≤f(x2)≤…≤f(xn+1)。由于目标是最小化f，所以在本发明 的实例中称x1为最好的顶点，xn+1为最差的顶点。S311根据 xr＝xo+ρ(xo‑xn+1)计算反射点xr。式中，xo表示n个最好点的质心，即  计算f(xr)。如果f(x1)≤f(xr)否则将xn+1替换为xr并转到步骤S316。S313如果 f(xn)≤f(xr)

[0013] xic＝xo+γ(xn+1‑xo)计算外部收缩点xic。计算f(xic)。如果f(xic)≤f(xn+1)， 则将xn+1替换为xic并转到步骤S316；否则，转到步骤S315；S315在 n个新的顶点x′i＝x1+γ(xi‑x1),i＝2,…,n+1处计算f(x′i)。将顶点x2、…、 xn+1替换为x′2、…、x′n+1。S316排序并对新的单纯形的n+1个顶点重新 标签为x1、x2、…、xn+1，使f(x1)≤f(x2)≤…≤f(xn+1)。

[0014] 如果 (ε>0是一个预先给定的阈值)或者达 到最大迭代次数，则停止；否则，转到步骤S311。

[0015] 优选的，步骤S4中，根据所述同相轴倾角及其相干值筛选所述 同相轴倾角，包括：S41设定同相轴倾角阈值δ；S42对于获得的k 个所述同相轴倾角，p1、p2、…、pk，以及所述同相轴倾角相应的相 干值，f1、f2、…、fk，对于任意两个同相轴倾角po和pq，计算 |po‑pq|；S43若|po‑pq|<δ，则移除具有较大相干值的同相轴倾角； 其中，o和q为大于0的整数，且o≠q,o≤k,q≤k。

[0016] 本发明提供的基于小生境差分进化算法的中值滤波方法，不需要 沿着三维数据的两个空间方向分别计算视倾角，而是直接获得真正的 三维倾角。由于该算法随机性和智能性，相比于局部倾斜叠加算法(三 维情况下需要沿着两个方向进行均匀采样)，在获得更高计算效率的 同时能够同时获得所有交叉同相轴的倾角。该算法同时获得所有交叉 同相轴倾角的能力和真三维的特性使其能够获得更好的滤波效果。

附图说明

[0017] 为了更清楚地说明本发明或现有技术中的技术方案，下面将对实 施例或现有技术描述中所需要使用的附图作一简单地介绍，显而易见 地，下面描述中的附图是本发明的一些实施例，对于本领域普通技术 人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得 其他的附图。

[0018] 图1是本发明提供的基于小生境差分进化算法的中值滤波方法 的流程示意图之一；

[0019] 图2是本发明提供的基于小生境差分进化算法的中值滤波方法 的时空窗中同相轴位置的识别示意图；

[0020] 图3是本发明提供的基于小生境差分进化算法的中值滤波方法 的hill‑valley探测函数的方案示意图

[0021] 图4是本发明提供的信噪比为20dB的二维含高斯白噪地震数据；

[0022] 图5a是本发明提供的基于迭代f‑k方法的倾角导向中值滤波结 果；图5b是基于小生境差分进化算法的倾角导向中值滤波结果；图 5c图5a同含噪数据的残差；图5d图5b同含噪数据的残差；

[0023] 图6是本发明提供的信噪比为20dB的三维含高斯白噪地震数据；

[0024] 图7是本发明提供的在inline号为130位置处的数据切片：图 7a是基于迭代f‑k方法的倾角导向中值滤波结果；图7b是基于小生 境差分进化算法的倾角导向中值滤波结果；图7c图7a同含噪数据的 残差；图7d图7b同含噪数据的残差；

[0025] 图8是本发明提供的在crossline号为160位置处的数据切片： 图8a是基于迭代f‑k方法的倾角导向中值滤波结果；图8b是基于小 生境差分进化算法的倾角导向中值滤波结果；图8c图8a同含噪数据 的残差；图8d图8b同含噪数据的残差；

[0026] 图9是本发明提供的在时间为0.2s位置处的数据切片：图9a基 于迭代f‑k方法的倾角导向中值滤波结果；图9b基于小生境差分进化 算法的倾角导向中值滤波结果；图9c图9a同含噪数据的残差；图9d 图9b同含噪数据的残差。

具体实施方式

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

[0028] 本发明发明人在对局部倾斜叠加法、PWD方法和迭代f‑k方法研 究的过程中发现，局部倾斜叠加法在处理交叉同相轴时只能计算出最 大相干值所对应的倾角方向。此外，三维情况下该方法需要沿着两个 方向进行均匀采样，如沿每个方向的采样点个数设为100，4
则需要进 行10次的相干计算，会严重影响计算效率。而PWD方法在求解交 叉同相轴的倾角时方程是欠定的，最终的解严重依赖初始条件及正则 化条件，该方法对强噪音敏感，而强噪音在实际数据中广泛存在，导 致用二维PWD方法沿着三维数据的两个空间方向进行计算难以准确 获得三维数据的实际倾角。对于迭代f‑k方法，虽然交叉同相轴的个 数不受限制，但需预先设定同相轴交叉的个数(即迭代次数)，且由 于采取迭代的方式，势必影响次序靠后的同相轴倾角的计算精度。对 于三维数据，受f‑k变换的限制，一般采用沿着三维数据的两个空间 方向来分别计算视倾角。而这同样难以准确获得三维数据的实际倾角， 从而影响倾角导向中值滤波的滤波效果。

[0029] 在研究过程中，本发明的发明人通过观察交叉同相轴沿不同倾角 方向的相干曲线(二维)或曲面(三维)，发现它是多峰的局部极值 准确对应着相应同相轴的倾角。为此，基于沿不同倾角方向交叉同相 轴的相干性所表现出的多峰特性，本发明创造性的将交叉同相轴倾角 的求取转化成智能优化领域中的多峰优化问题，进而采用小生境差分 进化算法来计算交叉同相轴的倾角，将其应用于倾角导向的中值滤波。 小生境差分进化算法可以分布式将种群的子集驱使到搜索空间的不 同区域，进而在一次运行中同时分布式收敛到多峰空间中的多个吸引 域，从而同时获得各个交叉同相轴所对应的倾角。

[0030] 此外，对于三维地震数据，小生境差分进化算法不需要沿着三维 数据的两个空间方向分别计算视倾角，而是直接获得真正的三维倾角。 由于该算法随机性和智能性，相比于局部倾斜叠加算法(三维情况下 需要沿着两个方向进行均匀采样)，在获得更高计算效率的同时能够 同时获得所有交叉同相轴的倾角。该算法同时获得所有交叉同相轴倾 角的能力和真三维的特性使其能够获得更好的滤波效果。

[0031] 下面结合图1‑图3描述本发明的一种基于小生境差分进化算法 的中值滤波方法，其特征在于，包括步骤S1‑S5，具体的：

[0032] S1将待处理数据划分到一系列重叠的时空窗中。

[0033] 在本发明的一个实施例中，所述待处理数据是地震数据，进一步 的，所述地震数据是二维或三维地震数据。所述时空窗是沿着时间空 间方向的取样窗口。

[0034] 在本发明的一个实施例中，将待处理数据划分到重叠的时空窗中 时，可以选择时间和空间都重叠50％，在其它实施例中，也可以选择 其它重叠方式。

[0035] S2获取单个所述时空窗中的同相轴能量曲线，根据所述同相轴 能量曲线的局部极大值获得同相轴位置。

[0036] 在本发明的一个实施例中，在二维情况下，计算同相轴能量的公 式为：

[0037]

[0038] 式中，t表示计算时间位置，H表示空间窗的空间宽度，x0的值分别 为H/4、H/2和3H/4，E表示计算位置的能量，d表示待处理数据，T为 能量计算范围的半宽，|*|表示对数据取绝对值，φ为大于等于1的常 数，在本发明的一个实施例中取4，q为大于等于1的常数，在本发 明的一个实施例中取10。

[0039] 根据同相轴能量曲线的局部极大值来识别同相轴的位置。如图2 所示，是二维情况下x0取H/2时，某一时空窗中同相轴位置的识别示 意图，其中地震数据由五个同相轴构成，同相轴有交叉，右侧蓝色曲 线为在H/2处根据计算同相轴能量的公式获得的能量曲线，红色五角 星为识别出的能量曲线的局部极大值，即同相轴位置。

[0040] 在本发明的一个实施例中，在三维情况下，计算同相轴能量的公式为：

[0041] 其中，若 否则， 

[0042] E为所述时空窗中的同相轴能量，t表示计算时间位置，H表示 所述时空窗的空间宽度，d表示待处理数据，T为能量计算范围的半 宽，φ为大于等于1的常数，在本发明的一个实施例中取4，x0的值 为H/4、H/2或3H/4，y0的值为H/4、H/2或3H/4，x0,y0表示计算的位 置，即二维情况下四等分点的中间三个点，三维情况下就是平面上的 九个点；q为大于等于1的常数，在本发明的一个实施例中取10。

[0043] S3在所述同相轴位置处，通过所述小生境差分进化算法获得同 相轴倾角及其相干值。

[0044] 其中所述小生境差分进化算法包括：基于邻域拥挤的差分进化算 法、hill‑valley探测函数和Nelder‑Mead单纯形算法；

[0045] 所述基于邻域拥挤的差分进化算法用于驱使局部解朝着各自的 局部最优移动；所述hill‑valley探测函数用于从所述局部解中获取局 部优秀解；所述Nelder‑Mead单纯形算法用于提高所述局部优秀解的 精度。

[0046] 其中，所述基于邻域拥挤的差分进化算法包括步骤：

[0047] S31：随机生成N个初始个体。

[0048] 所述个体是某一个倾角种群，表示所有倾角的集合，其中倾角是 优化算法的解。

[0049] S32：遍历每一个所述初始个体。

[0050] 寻找参数上最相似的m个个体形成一个亚种群，在本发明的一 个实施例中，所述最相似的个体是两个倾角间的距离(即欧氏距离) 最近的个体。

[0051] S33在所述亚种群中使用所述差分进化算法的变异和交叉操作生 成一个试验子个体。

[0052] 在本发明的一个实施例中，变异公式为

[0053]

[0054] 式中，vi＝[vi,1,vi,2,……,vi,D]是变异向量(D表示解的维数)，G表示迭代 次数，r1、r2和r3是[1，m]之间三个互斥的随机索引，它们也应同当前 目标向量 的索引不同，F∈[0,2]是缩放因子。

[0055] 交叉公式为

[0056]

[0057] 式中，CR∈[0,1]是控制参数，randj是在[0,1]区间均匀分布的随机数， k是[1,D](G+1)之间任意随机选择的自然数，用它来确保mi 至少从 获 得一个分量。

[0058] S34在整个种群中计算所述试验子个体与其他所述个体的欧氏距 离。

[0059] S35比较所述试验子个体与同所述试验子个体欧式距离最近个 体的相干值，如果所述试验子个体的相干值更小，则用所述试验子个 体替换掉所述欧式距离最近的个体。

[0060] S36重复M次步骤S32‑S35后算法终止，所述M为最大迭代次 数。

[0061] 在本发明的一个实施例中，所述亚种群的个体数m为25，所述 最大迭代次数N为25，针对交叉同相轴倾角的计算问题，基于邻域 拥挤的差分进化算法能够保证每个最优解附近出现一个解。

[0062] 其中，如图3所示，所述hill‑valley探测函数用来确定所述时空 窗中的两个所述局部解是否属于同一个小生境，包括步骤：

[0063] S37所述hill‑valley探测函数在两个所述局部解之间采样，获取 所述局部解之间的内部倾角样本的相干值。

[0064] S38若任何所述内部倾角样本的相干值大于两个所述局部解的相 干值，则确定两个所述局部解属于不同的小生境；

[0065] 否则，确定两个所述局部解属于同一个小生境。

[0066] S39从所述局部解中获取局部优秀解。

[0067] 在本发明的实施例中，对步骤S获得局部优秀解用单纯形算法继 续优化，所述单纯形算法是一种最流行的无导数非线性优化方法。所 谓单纯形是由n维空间中的n+1个顶点(x1、x2、…、xn+1)和它们的 函数值构成的一个几何形状。通过反射、扩张、收缩和缩边操作，选 出的优秀解可以确定性收敛到相应的局部最优。Nelder‑Mead单纯形 算法在每一次迭代过程中保持一个非退化的单纯形。Nelder‑Mead单 纯形算法需要设定四个常量参数，即反射系数ρ(ρ>0)、扩张系数χ (χ>1)、收缩系数γ(0<γ<1)及缩边系数σ(0<σ<1)。具体的，所 述单纯形算法用于提高所述局部优秀解的精度，包括步骤：

[0068] S310围绕单个局部优秀解构建具有n+1个顶点的单纯形，n为解 的维数，并排序并对n+1个顶点重新标签为x1、x2、…、xn+1，使 f(x1)≤f(x2)≤…≤f(xn+1)。由于目标是最小化f，所以在本发明的实例中 称x1为最好的顶点，xn+1为最差的顶点。

[0069] S311根据

[0070] xr＝xo+ρ(xo‑xn+1)

[0071] 计算反射点xr。式中，xo表示n个最好点的质心，即 计算f(xr)。如果f(x1)≤f(xr)

[0072] S312如果f(xr)

[0073] xe＝xo+χ(xr‑xo)

[0074] 计算扩张点xe。计算f(xe)。如果f(xe)

[0075] S313如果f(xn)≤f(xr)

[0076] xoc＝xo+γ(xr‑xo)

[0077] 计算外部收缩点xoc。计算f(xoc)。如果f(xoc)≤f(xr)，则将xn+1替换 为xoc并转到步骤S315；

[0078] S314如果f(xr)≥f(xn+1)，则进行内部收缩。根据

[0079] xic＝xo+γ(xn+1‑xo)

[0080] 计算外部收缩点xic。计算f(xic)。如果f(xic)≤f(xn+1)，则将xn+1替换 为xic并转到步骤S316；否则，转到步骤S315；

[0081] S315在n个新的顶点

[0082] x′i＝x1+γ(xi‑x1),i＝2,…,n+1

[0083] 处计算f(x′i)。将顶点x2、…、xn+1替换为x′2、…、x′n+1。

[0084] S316排序并对新的单纯形的n+1个顶点重新标签为x1、x2、…、 xn+1，使f(x1)≤f(x2)≤…≤f(xn+1)。

[0085] 如果 (ε>0是一个预先给定的阈值)或者达 到最大迭代次数，则停止；否则，转到步骤S311。

[0086] 工程应用需要综合考虑算法的稳定性、效率和精度。毫无疑问， 在地震数据处理的应用中需要对交叉同相轴的倾角进行稳定高精度 的计算。但是，稳定高精度的计算需要精细的搜索，这意味着更多的 函数调用，从而极大地影响计算效率。当完全只使用一种小生境方法 时，由于其内在的随机性，这个问题就会变的更加突出。本发明将 hill‑valley探测函数和Nelder‑Mead单纯形算法同基于邻域拥挤的差 分进化算法级联，避免了只采用基于邻域拥挤的差分进化算法而由于 其内在的随机性导致的过度搜索问题。

[0087] S4根据所述同相轴倾角及其相干值筛选所述同相轴倾角。

[0088] 根据所述同相轴倾角及其相干值筛选所述同相轴倾角，包括：

[0089] S41设定同相轴倾角阈值δ；

[0090] S42对于获得的k个所述同相轴倾角，p1、p2、…、pk，以及所 述同相轴倾角相应的相干值，f1、f2、…、fk，对于任意两个同相轴 倾角po和pq，计算|po‑pq|；

[0091] S43若|po‑pq|<δ，则移除具有较大相干值的同相轴倾角；

[0092] 其中，o和q为大于0的整数，且o≠q,o≤k,q≤k。

[0093] S5沿着筛选出的倾角依次进行中值滤波。

[0094] 优选的，在本发明的一个实施例中，所述的一种基于小生境差分 进化算法的中值滤波方法还包括步骤：

[0095] S6储存所述时空窗中的中值滤波数据。

[0096] S7确定所述各个时空窗中的待处理数据是否已被处理；

[0097] 若是，则返回执行步骤S2及步骤S2之后的步骤；

[0098] 若否，则累加所有时空窗中的所述中值滤波数据。

[0099] 如图4所示是信噪比为15dB的二维含高斯白噪地震数据，该数 据波场复杂，同相轴交叉现象严重。本发明分别使用基于迭代f‑k方 法和基于小生境差分进化算法的倾角导向中值滤波对随机噪声进行 压制，并对二者进行了对比。图5a是基于迭代f‑k方法的倾角导向中 值滤波的滤波结果。图5b是基于小生境差分进化算法的倾角导向中 值滤波的滤波结果。图5c是图5a同含噪数据的残差。图5d是图5b 同含噪数据的残差。

[0100] 从滤波结果和残差数据中可以看出基于小生境差分进化算法的 倾角导向中值滤波可以很好地保护同相轴，残差结果几乎是不相干的。 而基于迭代f‑k方法的倾角导向中值滤波对于交叉轴的滤波效果不好， 弱能量的同相轴会被遗漏掉。而且由于迭代f‑k方法对于局部倾角不 精确的估计，还会引入一些异常振幅。

[0101] 图6是信噪比为20dB的三维含高斯白噪地震数据。该地震数据 由于同相轴交叉现象的存在，波场很复杂。分别使用基于迭代f‑k方 法和基于小生境差分进化算法的倾角导向中值滤波对随机噪声进行 压制。图7a和7b分别是基于迭代f‑k方法和基于小生境差分进化算 法的倾角导向中值滤波在inline号为130位置处的滤波数据切片。图 7c和7d分别是它们同含噪数据的残差。图8a和8b分别是这两种滤 波方法在crossline号为160位置处的滤波数据切片。图8c和8d分别 是它们同含噪数据的残差。图9a和9b分别是这两种滤波方法在时间 为2.0s处的滤波数据切片。图9c和9d分别是它们同含噪数据的残差。

[0102] 可以看到对于基于小生境差分进化算法的倾角导向中值滤波的 滤波结果，同相轴得到了很好的保护，在残差数据中基本没有相干信 号。但是对于基于迭代f‑k方法的倾角导向中值滤波的滤波结果，一 些相干同相轴却被遗漏了。对于基于迭代f‑k方法的倾角导向中值滤 波，需要预先设定迭代滤波次数。当时空窗中同相轴的倾角个数少于 滤波次数时，会引入额外的噪音(如图7a)。而对于基于小生境差分 进化算法的倾角导向中值滤波，滤波次数是自适应的。基于小生境差 分进化算法的倾角导向中值滤波同时获得所有交叉同相轴倾角的能 力和真三维的特性使得它相比于基于迭代f‑k方法的倾角导向中值滤 波能够获得更好的滤波效果。

[0103] 本发明提供的基于小生境差分进化算法的中值滤波方法，不需要 沿着三维数据的两个空间方向分别计算视倾角，而是直接获得真正的 三维倾角。由于该算法随机性和智能性，相比于局部倾斜叠加算法(三 维情况下需要沿着两个方向进行均匀采样)，在获得更高计算效率的 同时能够同时获得所有交叉同相轴的倾角。该算法同时获得所有交叉 同相轴倾角的能力和真三维的特性使其能够获得更好的滤波效果。

[0104] 最后应说明的是：以上实施例仅用以说明本发明的技术方案，而 非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领 域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技 术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修 改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方 案的精神和范围。