多分量地震数据矢量去噪方法及多分量地震数据矢量去噪
装置
技术领域
背景技术
[0002] 多分量地震数据同时记录了质点在垂直和
水平方向的运动情况，是
地震波矢量场的完整记录。矢量数据不仅同时提供了纵波和横波的走时(速度)、振幅和
频率信息，还提供了同一波在各分量间的相对振幅关系，以及不同波的振幅差异和走时差异等信息。这些矢量波场特有的信息能更加详细、精确地刻画地下介质的构造、岩性、
流体饱和度、孔隙压
力、裂缝等特征。然而由于采集环境、勘探目标和地表条件的复杂性，采集的多分量地震数据
信噪比通常都很低。充分利用矢量波场信息的前提是必须对多分量地震数据进行保幅、保矢量特征的去噪处理。
[0003]
现有技术中，f-x(y)域预测滤波、f-x(y)域Cadzow滤波、中值滤波、SVD奇异值分解、以及基于数学变换的F-K去噪、S变换去噪、
小波变换去噪等方法，这些滤波方法在常规纵波资料处理中已取得很好的去噪效果。但对多分量地震数据进行处理时，现有技术主要是将每个分量当作一个标量场分开来处理，容易破坏多分量地震波场的矢量特征。而且由于多分量地震噪声较为复杂，一些对纵
波数据去噪效果好的方法用于横波数据时去噪效果并不理想。
[0004] 极化滤波是最早的多分量地震数据去噪方法，该方法能很好地识别并分离不相互重叠的不同偏振的波。Naghizadeh和Sacchi提出了一种基于3分量向量自回归(VAR)模型的f-x去噪方法，但他们假设地震反射同相轴是线性的，使其应用范围受到限制，Rodriguez等提出的基于组稀疏约束时频变换的随机噪声衰减方法能对多个分量数据联合进行处理，可较好地保留分量间的振幅相对关系和弱有效
信号，但求稀疏解与时频变换均非常耗时，难以对大批量地震数据进行处理。Wang和寻超等将矢量统计排序滤波扩展应用于多分量地震数据，首先根据相邻道矢量信号间的相关性搜索最优局部同相轴方向，然后沿同相轴选取矢量点进行矢量统计排序滤波。该方法要求地震信号横向连续性要好，对于横向相关性差或只有几道的地震数据无法有效去除噪声。
发明内容
[0005] 本发明的目的在于，针对上述现有技术中的不足，提供一种多分量地震数据去噪方法及多分量地震数据去噪装置，以解决传统滤波方法在去除多分量地震数据噪声的同时，会损坏多分量地震数据的矢量特征，不能对地震数据信号进行有效保护的问题。
[0006] 为实现上述目的，本发明
实施例采用的技术方案如下：
[0007] 第一方面，本发明实施例提供了一种多分量地震数据矢量去噪方法，包括：
[0008] 对多分量
地震道集数据进行分解，得到多个小
块多分量地震数据；对每个所述小块多分量地震数据进行四元傅里叶变换，得到四元频域地震数据；在四元频域中，
抽取频率切片，利用卡佐cadzow滤波方法对所述频率切片进行滤波，得到滤波后的四元频域地震数据；根据滤波后的所述四元频域地震数据，计算反四元傅里叶变换，获取滤波后的各分量地震数据。
[0009] 进一步地，所述对所述多分量地震道集数据进行分解，得到多个小块多分量地震数据，包括：
[0010] 将所述多分量地震道集数据划分为时间长度为T、水平方向道数为N的多个小块多分量地震数据，其中，多个所述小块多分量地震数据中任意相邻的两个小块多分量地震数据之间存在数据重叠，其中，N为大于0、小于地震道集数据总道数的整数，T为大于0、小于地震道集数据时间
采样点数的整数。
[0011] 进一步地，所述对多分量地震道集数据进行分解，得到多个小块多分量地震数据之后，还包括：
[0012] 将所述小块多分量地震数据进行四元数转换，得到转换后的小块多分量地震数据为：q(t,x)＝a(t,x)+b(t,x)i+c(t,x)j+d(t,x)k，其中，x为地震道序数，且x＝1，2，…，N，t为时间采样点数，且t＝1，2，…，T，i、j、k为四元数虚数单位，a(t,x)为压力分量地震数据，b(t,x)为x分量地震数据，c(t,x)为y分量地震数据，d(t,x)为z分量地震数据，当缺少其中某个分量地震数据时用零替代；
[0013] 相应地，所述对多个所述小块多分量地震数据进行四元傅里叶变换，得到四元频域地震数据，包括：
[0014] 对所述转换后的小块多分量地震数据进行四元傅里叶变换，得到四元频域地震数据。
[0015] 进一步地，所述在四元频域中，抽取频率切片，利用卡佐cadzow滤波方法对所述频率切片进行滤波，得到滤波后的四元频域地震数据，包括：
[0016] 在四元频域中，对每个
频率抽取频率切片Q1、Q2、…、QN，其中Qi＝1，2，…N表示四元频域中给定频率下第i道的值，利用频率切片构建四元数汉克尔Hankel矩阵：
[0017]
[0018] 其中M为小于N的正整数；
[0019] 对所述四元数汉克尔Hankel矩阵利用奇异值分解法进行分解；
[0020] 保留k个最大的奇异值，其它奇异值置零，重构后得到降秩矩阵Fk(A)；
[0021] 对所述降秩矩阵Fk(A)沿反对
角线求平均，得到滤波后的频率切片。
[0022] 进一步地，所述根据滤波后的所述四元频域地震数据，计算获取滤波后的各分量地震数据，包括：
[0023] 对所有频率切片均完成滤波后，得到滤波后的四元频域多分量地震数据，对滤波后的所述四元频域多分量地震数据进行反四元傅里叶变换，得到滤波后的各分量地震数据。
[0024] 第二方面，本发明实施例还提供了一种多分量地震数据矢量去噪装置，包括：
[0025] 计算模块，对多分量地震道集数据进行分解，得到多个小块多分量地震数据；对每个所述小块多分量地震数据进行四元傅里叶变换，得到四元频域地震数据；
[0026] 滤波模块，用于在四元频域中，抽取频率切片，利用卡佐cadzow滤波方法对所述频率切片进行滤波，得到滤波后的四元频域地震数据；
[0027] 确定模块，用于根据滤波后的所述四元频域地震数据，计算反四元傅里叶变换，获取滤波后的各分量地震数据。
[0028] 进一步地，所述计算模块，具体用于将所述多分量地震道集数据划分为时间长度为T、水平方向道数为N的多个小块多分量地震数据，其中，多个所述小块多分量地震数据中任意相邻的两个小块多分量地震数据之间存在数据重叠，其中，N为大于0、小于地震道集数据总道数的整数，T为大于0、小于地震道集数据时间采样点数的整数。
[0029] 进一步地，所述计算模块，还用于将所述小块多分量地震数据进行四元数转换，得到转换后的小块多分量地震数据为：q(t,x)＝a(t,x)+b(t,x)i+c(t,x)j+d(t,x)k，其中，x为地震道序数，且x＝1，2，…，N，t为时间采样点序数，且t＝1，2，…，T，i、j、k为四元数虚数单位，a(t,x)为压力分量地震数据，b(t,x)为x分量地震数据，c(t,x)为y分量地震数据，d(t,x)为z分量地震数据，当缺少其中某个分量地震数据时用零替代；
[0030] 相应地，所述计算模块，用于对多个所述小块多分量地震数据进行四元傅里叶变换，得到四元频域地震数据，具体为：对所述转换后的小块多分量地震数据进行四元傅里叶变换，得到四元频域地震数据。
[0031] 进一步地，所述滤波模块，具体用于在四元频域中，对每个频率抽取频率切片Q1、Q2、…、QN，其中Qi＝1，2，…N表示四元f-x域中给定频率下第i道的值，利用频率切片构建四元数汉克尔Hankel矩阵：
[0032]
[0033] 其中M为小于N的正整数；
[0034] 对所述四元数汉克尔Hankel矩阵利用奇异值分解法进行分解；
[0035] 保留k个最大的奇异值，其它奇异值置零，重构后得到降秩矩阵Fk(A)；
[0036] 对所述降秩矩阵Fk(A)沿反对角线求平均，得到滤波后的频率切片。
[0037] 进一步地，所述确定模块，具体用于对所有频率切片均完成滤波后，得到滤波后的四元频域多分量地震数据，对滤波后的所述四元频域多分量地震数据进行反四元傅里叶变换，得到滤波后的各分量地震数据。
[0038] 本发明的有益效果是：通过将四元傅里叶变换和卡佐cadzow滤波方法相结合，在矢量空间中对多分量地震数据进行滤波，在有效压制噪声的同时，较好地保留和恢复原始多分量地震信号的矢量特征。
附图说明
[0039] 为了更清楚地说明本发明实施例的技术方案，下面将对实施例中所需要使用的附图作简单地介绍，应当理解，以下附图仅示出了本发明的某些实施例，因此不应被看作是对范围的限定，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他相关的附图。
[0040] 图1为本
申请实施例提供的多分量地震数据矢量去噪方法流程示意图一；
[0041] 图2为本申请实施例提供的多分量地震数据矢量去噪方法流程示意图二；
[0042] 图3为本申请实施例提供的多分量地震数据矢量去噪方法流程示意图三；
[0043] 图4为本申请一实施例提供的一种合成地震数据加上随机噪声后的x分量、y分量和z分量效果图；
[0044] 图5为本申请提供的利用多分量地震数据矢量去噪方法去噪效果图；
[0045] 图6为本申请提供的利用常规卡佐Cadzow滤波方法进行标量场滤波效果图；
[0046] 图7为本申请一实施例提供的多分量地震数据矢量去噪装置结构示意图；
[0047] 图8为本申请另一实施例提供的多分量地震矢量数据去噪装置结构示意图。
具体实施方式
[0048] 为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。
[0049] 图1为本申请实施例提供的多分量地震数据矢量去噪方法流程示意图一；如图1所示，该多分量地震数据矢量去噪方法包括：
[0050] S101、对多分量地震道集数据进行分解，得到多个小块多分量地震数据。
[0051] 具体地，对多分量地震道集数据进行分解之前，需先进行多分量地震道集数据的采集，此处，读取二维多分量地震数据道集，所谓道集，即就是单个地震
检波器的接收记录，称之为一个地震道，多个地震道的集合，简称道集。
[0052] S102、对每个所述小块多分量地震数据进行四元傅里叶变换，得到四元频域地震数据。
[0053] 需要说明的是，对所述小块多分量地震数据沿时间方向进行四元傅里叶变换，将时域信号转换到频域信号，从而更有利于信号分析，方便对地震数据的分析处理。
[0054] S103、在四元频域中，抽取频率切片，利用卡佐cadzow滤波方法对所述频率切片进行滤波，得到滤波后的四元频域地震数据。
[0055] 具体地，利用卡佐cadzow滤波方法对各个小块多分量地震数据所对应的频率切片进行滤波，获得完整的滤波后的四元频域地震数据。对所述多分量地震数据进行滤波，有效的去除了所述多分量地震数据中的噪声信息，同时很好的保留了矢量特征效果。
[0056] 需要说明的是，多分量地震数据同时记录了质点在垂直和水平方向的运动情况，是地震波矢量场的完整记录。矢量数据不仅同时提供了纵波和横波的走时(速度)、振幅和频率信息，还提供了同一波在各分量间的相对振幅关系，以及不同波的振幅差异和走时差异等信息。这些矢量波场特有的信息能更加详细、精确地刻画地下介质的构造、岩性、流体饱和度、孔隙压力、裂缝等特征。
[0057] S104、根据滤波后的所述四元频域地震数据，计算获取滤波后的各分量地震数据。
[0058] 具体地，对所述滤波后的四元频域地震数据做反四元傅里叶变换，以得到滤波后的各分量地震数据。
[0059] 图2为本申请实施例提供的多分量地震数据矢量去噪方法流程示意图二，如图2所示，进一步地，所述对所述多分量地震道集数据进行分解，得到多个小块多分量地震数据，包括：
[0060] S201、将所述多分量地震道集数据划分为时间长度为T、水平方向道数为N的多个小块多分量地震数据，其中，多个所述小块多分量地震数据中任意相邻的两个小块多分量地震数据之间存在数据重叠，其中，N为大于0、小于道集地震数据总道数的整数，T为大于0、小于道集地震数据时间采样点数的整数。
[0061] 需要说明的是，多个所述小块多分量地震数据中任意相邻的两个小块多分量地震数据之间存在数据重叠，以避免产生边界效应。
[0062] 具体地，一种可实现的方式中，假设多分量地震数据为(1、2、3、4、5、6)，将其分解为多个小块多分量地震数据，多个小块分量地震数据数据之间重叠方式可为：(1、2、3)，(2、3、4)，(3、4、5)，(4、5、6)。
[0063] 进一步地，所述对多分量地震道集数据进行分解，得到多个小块多分量地震数据之后，还包括：
[0064] S202、将所述小块多分量地震数据进行四元数转换，得到转换后的小块多分量地震数据为：q(t,x)＝a(t,x)+b(t,x)i+c(t,x)j+d(t,x)k，其中，x为地震道序数，且x＝1，2，…，N，t为时间采样点数，且t＝1，2，…，T，i、j、k为四元数虚数单位，a(t,x)为压力分量地震数据，b(t,x)为x分量地震数据，c(t,x)为y分量地震数据，d(t,x)为z分量地震数据，当缺少其中某个分量地震数据时用零替代。
[0065] 其中，地震道序数x和时间采样点序数t均为大于零的整数。
[0066] 需要说明的是，四元数理论为矢量
信号处理提供了有力的工具，一个四元数包含一个
实部和三个
虚部，因而非常适合表达三分量或四分量地震数据。通过对所述小块多分量地震数据进行四元数转换，有效的将原始标量数据用矢量进行表达。
[0067] 具体地，i、j、k为四元数虚数单位，三个虚数单元满足的乘法律如下：
[0068] ii＝jj＝kk＝-1
[0069] ij＝-ji＝k
[0070] jk＝-kj＝i
[0071] ki＝-ik＝j
[0072] 因而四元数的乘法不满足交换律。
[0073] 相应地，所述对每个所述小块多分量地震数据进行四元傅里叶变换，得到四元频域地震数据，包括：
[0074] S203、对所述转换后的小块多分量地震数据进行四元傅里叶变换，得到四元频域地震数据。
[0075] 需要说明的是，此处为，对用四元数表示的多个所述小块多分量地震数据沿时间方向进行四元傅里叶变换，以将所述多分量地震数据从时域转换到频域进行分析处理。
[0076] 具体地，四元傅里叶变换有多种定义形式：由四元数的乘积非交换特性就可定义出左型、右型和双边型三种四元数傅里叶变换形式；由于四元数有三个虚数单元，四元傅里叶变换的变换轴可以在三维空间中沿任意方向选取；还可定义出单变换轴和双变换轴的四元傅里叶变换。其中，广义傅里叶变换公式用途较为广泛，其定义为：
[0077]
[0078] 其中f表示频率，μ为变换轴，且可以是任意的单位纯虚四元数。常用的变换轴μ选取方式为
[0079] 图3为本申请实施例提供的多分量地震数据矢量去噪方法流程示意图三，如图3所示，进一步地，所述在四元频域中，抽取频率切片，利用卡佐cadzow滤波方法对所述频率切片进行滤波，得到滤波后的四元频域地震数据，包括：
[0080] S301、在四元频域中，对每个频率抽取频率切片Q1、Q2、…、QN，其中Qi＝1，2，…N表示四元频域中给定频率下第i道的值，利用频率切片构建四元数汉克尔Hankel矩阵：
[0081]
[0082] 其中M为小于N的正整数。
[0083] 需要说明的是，上述是在频域中，对采用四元傅里叶变换后的每一个小块多分量地震数据，其所对应的频率中抽取频率切片。
[0084] 具体地，上述所述汉克尔Hankel矩阵中，M为正整数，为了使矩阵尽可能接近方阵，通常取 符号 表示向上取整。
[0085] S302、对所述四元数汉克尔Hankel矩阵利用奇异值分解法进行分解。
[0086] 具体地，对所述汉克尔Hankel矩阵A采用SVD奇异值分解后，得到A＝USVH，其中，上标H表示共轭转置，U和V的列分别为左奇异向量和右奇异向量，S为实对角线矩阵，对角线元素为奇异值。
[0087] S303、保留k个最大的奇异值，其它奇异值置零，重构后得到降秩矩阵Fk(A)。
[0088] 需要说明的是，保留k个最大的奇异值，其它奇异值置零后，将得到新的对角线矩阵S′，对其进行重构，得到降秩矩阵Fk(A)＝US′VH，k常取为小块地震数据内所含同相轴的个数。
[0089] S304、对所述降秩矩阵Fk(A)沿反对角线求平均，得到滤波后的频率切片。
[0090] 进一步地，根据滤波后的所述四元频域地震数据，计算获取滤波后的各分量地震数据，包括：
[0091] 对所有频率切片均完成滤波后，得到滤波后的四元频域地震数据，对滤波后的所述四元频域地震数据进行反四元傅里叶变换，得到滤波后的各分量地震数据。
[0092] 具体地，所有频率切片均完成滤波后，得到滤波后的四元频域地震数据Q′(f,x)，对其做反四元傅里叶变换即得到滤波后的各分量地震数据q′(t,x)，具体根据反四元傅里叶变换公式进行计算，公式如下所示：
[0093]
[0094] 下面结合具体实验来验证本申请所述多分量地震数据矢量去噪方法的有效性。
[0095] 图4为本申请一实施例提供的一种合成地震数据加上随机噪声后的x分量、y分量和z分量效果图，图5为本申请提供的利用多分量地震数据矢量去噪方法去噪效果图，图6为本申请提供的利用常规卡佐Cadzow滤波方法进行标量场滤波效果图，如图4、图5和图6所示，所述随机噪声较大，使用本发明提供的所述多分量地震数据矢量去噪方法进行矢量滤波，滤波结果如图5所示。为了对比，对各个分量分别采用常规卡佐Cadzow滤波方法进行标量场滤波，结果见图6。
[0096] 由图可知，本发明提供的所述多分量地震数据矢量去噪方法比常规标量去噪方法能更有效地去除随机噪声并更好地保留有效地震信号。特别是信噪比很低的x分量，矢量去噪方法取得更好的效果，使多分量地震数据的矢量特征得到更好地保留。
[0097] 第二实施例
[0098] 图7为本申请一实施例提供的多分量地震数据矢量去噪装置结构示意图，如图7所示，该多分量地震数据矢量去噪装置，包括：计算模块701、滤波模块702以及确定模块703，其中：
[0099] 计算模块701，用于对多分量地震道集数据进行分解，得到多个小块多分量地震数据；对每个所述小块多分量地震数据进行四元傅里叶变换，得到四元频域地震数据。
[0100] 滤波模块702，用于在四元频域中，抽取频率切片，利用卡佐cadzow滤波方法对所述频率切片进行滤波，得到滤波后的多四元频域地震数据。
[0101] 确定模块703，用于根据滤波后的所述四元频域地震数据，计算获取滤波后的各分量地震数据。
[0102] 进一步地，所述计算模块701，具体用于将所述多分量地震道集数据划分为时间长度为T、水平方向道数为N的多个小块多分量地震数据，其中，多个所述小块多分量地震数据中任意相邻的两个小块多分量地震数据之间存在数据重叠，其中，N为大于0、小于道集地震数据总道数的整数，T为大于0、小于道集地震数据时间采样点数的整数。
[0103] 进一步地，所述计算模块701，还用于将所述小块多分量地震数据进行四元数转换，得到转换后的四元频域地震数据为：q(t,x)＝a(t,x)+b(t,x)i+c(t,x)j+d(t,x)k，其中，x为地震道序数，且x＝1，2，…，N，t为时间采样点数，且t＝1，2，…，T，i、j、k为四元数虚数单位，a(t,x)为压力分量地震数据，b(t,x)为x分量地震数据，c(t,x)为y分量地震数据，d(t,x)为z分量地震数据，当缺少其中某个分量地震数据时用零替代。
[0104] 相应地，所述计算模块701，用于对每个所述小块多分量地震数据进行四元傅里叶变换，得到四元频域地震数据，具体为：对所述转换后的小块多分量地震数据进行四元傅里叶变换，得到四元频域地震数据。
[0105] 进一步地，所述滤波模块702，具体用于在四元频域中，对每个频率抽取频率切片Q1、Q2、…、QN，其中Qi＝1，2，…N表示四元f-x域中给定频率下第i道的值，利用频率切片构建四元数汉克尔Hankel矩阵：
[0106]
[0107] 其中M为小于N的正整数。
[0108] 对所述四元数汉克尔Hankel矩阵利用奇异值分解法进行分解。
[0109] 保留k个最大的奇异值，其它奇异值置零，重构后得到降秩矩阵Fk(A)。
[0110] 对所述降秩矩阵Fk(A)沿反对角线求平均，得到滤波后的频率切片。
[0111] 进一步地，所述确定模块703，具体用于对所有频率切片均完成滤波后，得到滤波后的四元频域多分量地震数据，对滤波后的所述四元频域多分量地震数据进行反四元傅里叶变换，得到滤波后的各分量地震数据。
[0112] 上述装置可用于执行上述方法实施例提供的方法，具体实现方式和技术效果类似，这里不再赘述。
[0113] 以上这些模块可以是被配置成实施以上方法的一个或多个集成
电路，例如：一个或多个特定集成电路(Application Specific Integrated Circuit，简称ASIC)，或，一个或多个
微处理器(digital singnal processor，简称DSP)，或，一个或者多个现场可编程
门阵列(Field Programmable Gate Array，简称FPGA)等。再如，当以上某个模块通过处理元件调度程序代码的形式实现时，该处理元件可以是通用处理器，例如
中央处理器(Central Processing Unit，简称CPU)或其它可以调用程序代码的处理器。再如，这些模块可以集成在一起，以片上系统(system-on-a-chip，简称SOC)的形式实现。
[0114] 图8为本申请一实施例提供的多分量地震数据矢量去噪装置结构示意图，如图8所示，该装置包括：处理器801和
存储器802，其中：
[0115] 存储器802用于存储程序，处理器801调用存储器802存储的程序，以执行上述方法实施例。具体实现方式和技术效果类似，这里不再赘述。
[0116] 可选地，本发明还提供一种程序产品，例如计算机可读存储介质，包括程序，该程序在被处理器执行时用于执行上述方法实施例。
[0117] 在本发明所提供的几个实施例中，应该理解到，所揭露的装置和方法，可以通过其它的方式实现。例如，以上所描述的装置实施例仅仅是示意性的，例如，所述单元的划分，仅仅为一种逻辑功能划分，实际实现时可以有另外的划分方式，例如多个单元或组件可以结合或者可以集成到另一个系统，或一些特征可以忽略，或不执行。另一点，所显示或讨论的相互之间的耦合或直接耦合或通信连接可以是通过一些
接口，装置或单元的间接耦合或通信连接，可以是电性，机械或其它的形式。
[0118] 所述作为分离部件说明的单元可以是或者也可以不是物理上分开的，作为单元显示的部件可以是或者也可以不是物理单元，即可以位于一个地方，或者也可以分布到多个网络单元上。可以根据实际的需要选择其中的部分或者全部单元来实现本实施例方案的目的。
[0119] 另外，在本发明各个实施例中的各功能单元可以集成在一个处理单元中，也可以是各个单元单独物理存在，也可以两个或两个以上单元集成在一个单元中。上述集成的单元既可以采用
硬件的形式实现，也可以采用硬件加
软件功能单元的形式实现。
[0120] 上述以软件功能单元的形式实现的集成的单元，可以存储在一个计算机可读取存储介质中。上述软件功能单元存储在一个存储介质中，包括若干指令用以使得一台计算机设备(可以是个人计算机，
服务器，或者网络设备等)或处理器(英文：processor)执行本发明各个实施例所述方法的部分步骤。而前述的存储介质包括：U盘、移动
硬盘、
只读存储器(英文：Read-Only Memory，简称：ROM)、
随机存取存储器(英文：Random Access Memory，简称：RAM)、磁碟或者光盘等各种可以存储程序代码的介质。
