号。将式（5)变换到频率域得：
[0067]
( 6 )
[0068] 其中：D为频率域的炮记录，Hm是m次谐波的相位谱，R为地震震动向大地施加的 外力，am为傅立叶变换后振幅谱的权系数，m为谐波的阶次。将式（6)变形得
[0069]
(7.)
[0070] 其中：D为频率域的炮记录，R为地震震动向大地施加的外力，Hm是m次谐波的相 位谱，am为傅立叶变换后振幅谱的权系数，m为谐波的阶次，Pm是m次谐波振幅谱与基波信 号振幅谱的比值。本发明采用相移法分离地面力信号求取Pm，令
[0071]
(.8)
[007引其中P为预测滤波算子，则谐波干扰的预测和压制由式（9) (10)表示：
[0075] 在步骤102,用伪多道匹配方法对预测的谐波进行自适应修正，结合传统的多道匹 配和约束道均衡方法，得到扩展的伪多道匹配滤波方法。假设初始地震记录为P，有效波记 录为P。，真实的谐波干扰为m，则
[0076] p〇=p-m (11)
[0077] 假设预测的谐波分量为m。，用a表示滤波因子，则真实的谐波干扰可由下式求得
[0078] m=a*m〇 (12)
[0079] 则有效波记录可表示为
[0080] p〇=p-a*m。 (13)
[0081] 通过采用最小二乘的思想最小化目标函数计算滤波因子a:
[0082] e(a) =IIp-MaI12 (14)
[008引其中：e(a)是求解a的最优化问题的函数，其中Μ为预测谐波噪音的矩阵表示形 式。
[0084] 在最小二乘意义下通过计算基于L2模的谐波自适应相减的滤波因子a来最小化 预测的谐波信号与实际谐波信号的差。运需要有效波能量与谐波能量满足正交性，如果不 满足正交性，滤波因子的求取就会出现问题，导致谐波压制不完全或者有效能量的损失。
[0085] 在通常情况下，我们在一个可滑动的窗口中匹配预测的谐波干扰。使用多道方程 代替式（13)。
[0086]
…）
[0087] 其中，p(t)是原始地震记录，nil(t)是预测的第i道谐波干扰，曰1似是多道匹配因 子，N是匹配的道数，Pe(t)为预测滤波后的地震记录。如果我们仅仅使用最小二乘原理而 不加任何约束，匹配的结果将会不准确。
[0088] 我们结合传统的多道匹配和约束道均衡方法，得到扩展的伪多道匹配滤波方法：
[0089]
。色）
[0090] 其中：Pe(t)为预测滤波后的地震记录，w(〇和"'(〇为111似和lAd导数道， fn'ia)为伪多道道频率表达式V(t)表示谐波的希尔伯特变换道。
[0091] 在步骤103,从原始记录中减去预测出的谐波，得到谐波压制后的记录。我们对如 图3所示的国内速度模型进行试算。该速度模型存在较多的起伏地层，同时在中深层存在 一个楠圆形高速异常体。模型大小为5000mX2500m，采用的网格间距为5mX5m。对可控震 源的滑动扫描方式进行模拟，正演模拟的参数如下：采用线性升频扫描方式，起始频率η =化Ζ，终止频率f2 = 8化Ζ，扫描周期Τ= 10s。滑动时间和记录时间都是3s，时间采样间 隔为0.5ms。四个震源车分别分别在1km，2km，3km，4km处激发。正演模拟得到的可控震源 滑动扫描正演模拟多炮记录如图4所示，图4为相关前的原始记录，图中5为图4所示的原 始记录与扫描信号互相关得到的相关记录。从图5中，我们可W清楚的看到，谐波干扰产生 在负时间轴上，与前面炮记录的有效信号混叠在一起，对前面炮记录的有效信号造成干扰， 如图中黑色箭头所示。
[0092] 为了还原较为准确的野外采集情况，我们对真实地面力信号的振幅、相位和频率 随机加一定范围内的误差，并分别采用本文提出的自适应匹配预测滤波方法和传统的预测 滤波算法进行试算，得到的第一炮的谐波干扰的压制结果如图6-9所示。其中图6为谐波 压制后的炮记录，图7为压制的谐波。图6中结果中，谐波干扰被压制的非常干净，而从图7 中也可W看出压制的谐波干扰中不含有有效波，换言之，通过本发明的谐波压制方法，基本 没有造成有效能量的损失，图8为传统预测滤波算法压制谐波的结果，谐波干扰没有压制 干净，而且造成了有效信号的损失。而图9所示的传统预测滤波算法压制谐波的结果中，谐 波干扰没有压制干净，而且造成了有效信号的损失。图10-12为采用本文提出的自适应匹 配预测滤波方法得到的其他炮的谐波压制结果，谐波干扰都得到了明显的消除，并没有造 成有效能量的损失。
[0093] 上述技术方案只是本发明的一种实施方式，对于本领域内的技术人员而言，在本 发明公开了应用方法和原理的基础上，很容易做出各种类型的改进或变形，而不仅限于本 发明上述【具体实施方式】所描述的方法，应此前面描述的方式只是优选的，而并不具有限制 性的意义。
【主权项】
1. 一种基于自适应匹配滤波算子的可控震源谐波压制方法，其特征在于，该基于自适 应匹配滤波算子的可控震源谐波压制方法包括： 步骤1，用传统预测滤波法对谐波干扰进行预测； 步骤2,用伪多道匹配方法对步骤1中预测的谐波进行自适应修正，结合传统的多道匹 配和约束道均衡方法，得到扩展的伪多道匹配滤波方法； 步骤3,从原始记录中减去步骤2中经过修正后的预测出的谐波，得到谐波压制后的记 录。2. 根据权利要求1所述的基于自适应匹配滤波算子的可控震源谐波压制方法，其特征 在于，步骤1包括： A) 采用相移法分离地面力信号求取，令其中P为预测滤波算子，Μ为谐波的阶次，Η"为m次谐波的相位谱，H i为基波的相位谱， P ">为m次谐波振幅与基波振幅的比值； B) 对谐波干扰进行预测和压制：其中，Η是谐波干扰，D为频率域的炮记录，Η"是m次谐波的相位谱，I为基波的相位 谱，R为地震震动向大地施加的外力，aiS基波的振幅。3. 根据权利要求1所述的基于自适应匹配滤波算子的可控震源谐波压制方法，其特征 在于，步骤2包括以下步骤： a) 假设初始地震记录为p，有效波记录为p。，真实的谐波干扰为m，则 P〇= p-m 假设预测的谐波分量为m。，用a表示滤波因子，则真实的谐波干扰可由下式求得： m = a氺m0 则有效波记录可表示为： p0= p-a*m〇 计算滤波因子a: e(a) = I |p-Ma| |2 其中M为预测谐波噪音的矩阵表示形式，e (a)为a的最小化目标函数； b) 在最小二乘意义下通过计算基于L2模的谐波自适应相减的滤波因子a来最小化预 测的谐波信号与实际谐波信号的差，在一个可滑动的窗口中匹配预测的谐波干扰，使用多 道方程式，其中，p(t)是原始地震记录，mi(t)是预测的第i道谐波干扰，ai(t)是多道匹配因子， N是匹配的道数，P(](t)为预测滤波后的地震记录； 根据约束道均衡的方法，将单道扩展为四道，其表达式为：其中：Pc(t)为预测滤波后的地震记录，mH(t)表示谐波的希尔伯特变换道，辦0和 为m(t)和mH(t)导数道，Wl表示每一道的滤波算子； 采用褶积算子来调整上式：其中：Pc(t)为预测滤波后的地震记录，和#?为m(t)和mH(t)导数道，fn(t)为各 道频率表达式，mH(t)表示谐波的希尔伯特变换道； 结合传统的多道匹配和约束道均衡方法，得到扩展的伪多道匹配滤波方法：其中：PD(t)为预测滤波后的地震记录，⑴和V(>)为m(t)和mH(t)导数道，f nil(t)为 伪多道道频率表达式，mH(t)表示谐波的希尔伯特变换道。
【专利摘要】本发明提供一种基于自适应匹配滤波算子的可控震源谐波压制方法，该基于自适应匹配滤波算子的可控震源谐波压制方法包括：步骤1，用传统预测滤波法对谐波干扰进行预测；步骤2，用伪多道匹配方法对步骤1中预测的谐波进行自适应修正，结合传统的多道匹配和约束道均衡方法，得到扩展的伪多道匹配滤波方法；步骤3，从原始记录中减去步骤2中经过修正后的预测出的谐波，得到谐波压制后的记录。该基于自适应匹配滤波算子的可控震源谐波压制方法效果好，精度高，能够很好地达到可控震源数据处理的要求的，并且为后续地震资料成像和解释提供可靠保证。
【IPC分类】G01V1/36
【公开号】CN105277986
【申请号】CN201510278108
【发明人】于富文, 曲英铭, 宋俊海, 胡立新, 莫延钢, 徐雷良, 卢湘鹏, 崔汝国, 时兴文
【申请人】中石化石油工程地球物理有限公司胜利分公司
【公开日】2016年1月27日
【申请日】2015年5月27日