基于主成分分析的地震资料工频干扰自动识别与压制方法与流程
本发明涉及一种地震数据处理去噪的方法，尤其是工频干扰严重的地震记录去噪方法，利用地震记录各道的平均能量值识别干扰道以及利用主成分分析技术压制工频干扰的方法。
背景技术：
：
地震勘探资料采集过程中往往受到工频干扰的影响，采集到的地震记录中常混有工频噪声，当工频干扰很强，且贯穿于一道地震数据的始终时，该地震道常表现为坏道，引起地震资料品质下降。通过对大量坏道数据的分析，很多地震道并不是源于地面附近的强干扰源，往往来自地震检波器本身的问题，与地震记录的时间深度无关，当前对坏道的处理是剔除坏道，这种处理降低了地震数据质量，而本专利提供的技术，能够有效识别由于工频干扰产生的坏道，这为进一步压制工频干扰提供了可能。针对工频干扰当前还没有非常有效的解决方法。目前主要的压制方法有陷波法、工频回归相减抑制法，自适应滤波法三类。文章“地震资料中工频干扰的自动识别与压制”提出陷波法滤除工业干扰，但该方法同时也会滤除50Hz频率的地震信号；文章“强工频干扰波的提取与消除方法”提出时间域内运用共轭梯度算法提取单频干扰，但该方忽略了实际地震资料中工频频率的可变性，研究表明在不同区域，不同时间，工频干扰频率是变化的；文章“压制地震资料中工业电干扰的余弦逼近法的改进及应用”提出基于余弦逼近法滤除工业干扰，但该方法无法实现自动识别干扰道；文章“在检波点域分离50Hz工业干扰”提出炮集转换滤除工业干扰，但该方法效率低并无法自动识别干扰道；文章“基于维纳滤波的50Hz工业干扰去噪方法及应用”提出维纳滤波器来调整参考道的振幅和相位的匹配来预测单频噪声，但该方法不适应于平稳性差的工频干扰；文章“基于独立分量分析的工频干扰消除技术”提出从混合信号中分离出相互独立的各个信号分量达到消除工业噪声的目的，但该方法依靠经验性，缺乏严格的理论依据；
CN104570118A公开的“一种基于双因素的自动识别与去除工业干扰的方法”提出正余弦加权逼近法处理工业干扰，但该方法忽略了工频的频率，相位和幅值的不固定性；
CN103630935A公开的“去除地震数据中交流电干扰信号的方法”提出频率域提取工频频率的方法，但该方法没有考虑工频谐波的存在；
CN101907726A公开的“一种自动识别和消除地震勘探工业电干扰的方法”提出基于自相关理论确定工频干扰，但该方法自动识别工频噪声时依赖于初至时间之前地震资料质量。
技术实现要素：
：
本发明的目的就在于针对上述现有技术的不足，提供一种基于主成分分析的地震资料工频干扰自动识别与压制方法。
本发明的主要思想是：地震采集过程中地震资料往往会受到工频干扰的影响，从而造成地震记录中的大量坏道或强干扰道，严重影响地震勘探数据的质量，本发明是通过自动识别受工频干扰的地震道，然后通过主成分分析技术分离出干扰，再对重构后的数据进行保幅，实现了地震资料工频干扰的自动识别与压制。
本发明是通过以下技术方案实现的：
基于主成分分析的地震资料工频干扰自动识别与压制方法，包括以下步骤：
a、针对待处理的原始地震记录，选取记录时间位于中下段的时窗TW，计算该地震记录在TW内的各道平均能量值，如公式
其中Ei为时窗内地震记录第i道平均能量值，xi(tj)为时窗内地震记录第i道的第j个采样点值，i＝1,2…n，n为地震记录总道数，j＝1,2…k，k为时窗内地震记录的采样点数，时窗TW的长度建议为总记录长度的到时间内；
b、定义第l道特征因子根据λl的取值识别干扰道，El为地震记录第l道在TW内的平均能量值，q为TW内具有最小平均能量值的地震道号，其中l＝1,2,…n，识别干扰道具体步骤如下：
b1、若对所有道均存在λl＜5，表明该地震记录无显著工频干扰，无需后续去噪过程，否则执行步骤b2；
b2、若λl≥5，表明地震记录有显著工频干扰，记第l道为干扰道，记为Ul；
c、提取与Ul相邻且包括Ul的连续5道记录组成子地震干扰道集，记为X，X中每一道为一列，求出X每列的均值，并用X中的每个元素减去该列的均值，得到X处理后的子地震道集，记为X*；
d、计算X*的协方差矩阵Γ，如公式
其中b为X*的列数，X*T为X*的转置矩阵，“·”表示矩阵乘法；
e、利用奇异值分解法计算出协方差矩阵Γ的特征值矩阵Λ和特征向量矩阵R，则存在公式
Γ＝R·Λ·RT (3)
其中Λ为特征值从大到小排列的对角矩阵，R每一列为对应特征值的特征向量，RT为转置矩阵，满足RT·R＝R·RT＝E,E为单位矩阵；
f、X经RT线性映射，得到主成分矩阵，如公式
Φ＝RT·X (4)
其中Φ为主成分矩阵；
g、将Φ第一行置0，得到Φ'，令
X1*＝R·Φ' (5)
其中X1*为重构矩阵；
h、提取X1*中对应干扰道的数据，记为Ul'，对Ul'进行保幅处理，如公式
其中Ul*为保幅处理后干扰道信号，A1为去噪后干扰道信号振幅绝对值的平均值，A2为与Ul相邻两道数据的振幅绝对值平均值之和的平均。
有益效果：经试验，本发明公开的一种基于主成分分析的地震资料工频干扰的自动识别与压制方法能够实现在地震勘探中压制工频干扰源带来的噪声,并且能够自动识别干扰道，该算法能够修复工频干扰引起的坏道，有效提高了地震资料质量，降低了采集成本。
附图说明：
图1实际地震记录
图2压制工频干扰后的地震记录
具体实施方式：
下面结合附图和实施例对本发明做进一步的详细说明：
基于主成分分析的地震资料工频干扰自动识别与压制方法，包括以下步骤：
a、采用实际单炮记录，本例TW＝700ms～800ms，计算TW内的Ei，如公式
其中Ei为地震记录第i道平均能量值，xi(tj)为地震记录第i道的第j个采样点元素，本例中i＝1,2…n，n＝48，j＝1,2…k，k＝100；
b、定义第l道特征因子根据λl的取值识别干扰道，El为地震记录第l道在TW内的平均能量值，本例中q为TW内的第一道，E1＝0.0042，其中l＝1,2,…n，识别干扰道具体步骤如下：
b1、若对所有道均存在λl＜5，表明该地震记录无显著工频干扰，无需后续去噪过程，本例中满足λl＜5的道为第1～7道，第9～42道，第44～48道，否则执行步骤b2；
b2、若λl≥5，表明地震记录有显著工频干扰，记第l道为干扰道，记为Ul，本例中满足λl≥5的道为第8和第43道，记为U8和U43，E8＝0.078，E43＝0.23，λ8＝18.57，λ43＝54.76；
c、以U8为例，提取第6，7，8，9，10道组成子地震干扰道集，记为X，X中每一道为一列，求出X每列的均值，并用X中的每个元素减去该列的均值，得到X处理后的子地震道集，记为X*；
d、计算X*的协方差矩阵Γ，如公式
其中b为X*的列数，X*T为X*的转置矩阵，“·”表示矩阵乘法，式中b＝5；
e、利用奇异值分解法计算出协方差矩阵Γ的特征值矩阵Λ和特征向量矩阵R，则存在公式
Γ＝R·Λ·RT (3)
其中Λ为特征值从大到小排列的对角矩阵，R每一列为对应特征值的特征向量，RT为转置矩阵，满足RT·R＝R·RT＝E,E为单位矩阵；
f、X经RT线性映射，得到主成分矩阵，如公式
Φ＝RT·X (4)
其中Φ为主成分矩阵；
g、将Φ第一行置0，得到Φ'，令
X1*＝R·Φ' (5)
其中X1*为重构矩阵；
h、提取X1*中对应干扰道的数据，记为Ul'，对Ul'进行保幅处理，如公式
其中Ul*为保幅处理后干扰道信号，A1为去噪后干扰道信号振幅绝对值的平均值，A2为与Ul相邻两道数据的振幅绝对值平均值之和的平均,本例中A2＝0.0689,A1＝0.0057。
- 还没有人留言评论。精彩留言会获得点赞！