一种大地电磁去噪方法及系统与流程

[0001]
本发明涉及数据处理领域，特别是涉及一种大地电磁去噪方法及系统。


背景技术：

[0002]
大地电磁测深法(magnetotelluric sounding，简称mt)在地表观测天然电磁场信号经地下介质产生的分量扰动值，以此来获取地下介质的电性结构。其场源为天然电磁场，有效信号起源于闪电放电、太阳风与地球磁层的相互作用，相较人工场源的电磁测深法，mt具有频带信息丰富、探测深度大等优势，因此被广泛应用于油气田普查勘探、地震预报、矿区找矿等多个领域中。但天然电磁场信号较弱，在测量期间易受到各种电磁噪声的干扰，在矿集区等强干扰环境中受到的影响尤甚。而有效信号与噪声干扰具有明显不一致的性质，有效信号具有随机、非平稳、非高斯等特性，干扰信号往往有着一定的极化方向、信号能量较强，且当强干扰源距测点较近时，干扰信号会完全覆盖有效信号。从而导致数据质量极差，在估算阻抗时发生严重偏差，引起近源干扰现象。因此，实现对噪声的有效压制是mt法能够获得准确地电结构的前提。
[0003]
目前针对大地电磁噪声的去除方法可主要分为以下三类：
①
时域处理方法：该类方法主要针对时间序列中明显存在的噪声形态进行识别并去除，其典型代表有形态滤波法、子空间增强算法、匹配追踪和压缩感知等方法，该类算法对于特定形态的噪声具有较好的处理效果，但对于复杂形态噪声的压制效果并不理想。比如采用矩形结构元素的形态滤波法仅对方波的处理效果最佳，而不能有效去除充放电三角波等其他形态的大尺度噪声与突变噪声。因此，目前常联合多种时域方法进行去噪处理，但相应的处理时间也随之增加；
②
时频转换方法：将时间序列转换到其他域(比如频域、小波域等)，通过筛选信号的谱信息来去除噪声，典型代表有小波分析、经验模态分解(empirical mode decomposition)、变分模态分解(variational mode decomposition)等方法，但具有母小波函数与分解层数选取、模态混叠、预设参数选取等解决难题。
③
频域处理方法：该类方法有最小二乘法、robust估计、远参考法等，且主要与阻抗估计相联系，比如远参考法利用参考点与测点之间的互功率谱数据代替测点本身的自功率谱数据进行阻抗计算，能够去除测点数据中与参考点不相关的噪声，但该方法对参考点数据的信噪比要求较高，且难以确定参考点的位置，当位置选择不恰当时，去除噪声的效果较差。因此当采集数据信噪比较高或存在噪声的形态较规整时，利用上述现有的大地电磁去噪技术能够取得较好的处理效果，但在实际观测时，噪声形态往往复杂多变，在矿集区尤甚，难以满足高信噪比的要求。
[0004]
根据大地电磁有效信号和干扰信号性质的差异，提出利用svd分解(singular value decomposition，奇异值分解)去除噪声。在进行svd分解前，需将一维时序向量转变为多维矩阵，svd分解则是通过多维矩阵的整体特征和每行向量的局部特征提取出不同的奇异值(从大到小排列，在一定程度上可理解为场源值)，通过奇异值计算出时序向量的各分解分量，每个分量为由不同场源激发的电磁波信号，通过判断分量是否为有效信号，来达到去噪的目的。因此该方法在去除噪声时主要与分解的分量个数有关(当分量个数较多时，
信号分解越细致，能将有效信号与噪声完全分离，但处理时间相应增长；当分量个数较少时，处理速度较快，但可能会损失有效信号)，对信噪比无要求，能够应用于矿集区的数据去噪处理。且该方法不依据噪声形态来进行去噪处理，而是通过识别并去除干扰(奇异值)场源。干扰源激发的电磁波形态复杂多变，因此该方法能够同时处理多种类型的噪声，不受噪声形态的限制。但当数据段为纯净信号时，仍进行筛选分量的去噪处理，便会出现有效信号的“过处理”问题。


技术实现要素：

[0005]
本发明的目的是提供一种大地电磁去噪方法及系统，以减少对数据的过处理，提高数据去噪的准确性。
[0006]
为实现上述目的，本发明提供了如下方案：
[0007]
一种大地电磁去噪方法，所述方法包括：
[0008]
利用大地电磁测深法测量大地不同深度的大地电磁脉动信号；
[0009]
对所述大地电磁脉动信号进行重叠分段，获得多个测深数据段；
[0010]
基于近似熵理论，确定每个测深数据段的近似熵；
[0011]
将近似熵小于或等于近似熵阈值的测深数据段确定为待去噪测深数据段，并将近似熵大于近似熵阈值的测深数据段确定为有效测深数据段；
[0012]
根据所述待去噪测深数据段，构造待去噪矩阵；
[0013]
对所述待去噪矩阵进行奇异值分解，获得所述待去噪矩阵的奇异值矩阵、左奇异矩阵和右奇异矩阵；
[0014]
将所述奇异值矩阵的每个对角线元素、所述左奇异矩阵对应的列向量和所述右奇异值矩阵对应的行向量相乘，获得所述待去噪测深数据段的多个场源分量矩阵；
[0015]
将每个场源分量矩阵转变为分量数据段；
[0016]
基于近似熵理论，确定每个分量数据段的近似熵；
[0017]
筛选分量数据段的近似熵大于所述待去噪测深数据段的近似熵的分量数据段，并将筛选出来的分量数据段求和，获得去噪后的测深数据段；
[0018]
将所述去噪后的测深数据段和所述有效测深数据段进行重构，获得重构去噪大地电磁脉动信号。
[0019]
可选的，所述对所述大地电磁脉动信号进行重叠分段，获得多个测深数据段，具体包括：
[0020]
对所述大地电磁脉动信号x1,x2,
…
x
w
进行重叠分段处理，获得每个测深数据段x
v
＝{x
(v-1)*n+1
,x
(v-1)*n+2
,
…
,x
(v-1)*n+n
}；
[0021]
其中，x1,x2,
…
x
w
分别为第1个、第2个、第w个大地电磁脉动信号，w为大地电磁脉动信号的长度，x
(v-1)*n+1
、x
(v-1)*n+2
、x
(v-1)*n+n
分别为第v个测深数据段的第1个、第2个和第n个数据，n为测深数据段长度，v＝1,2，
…
,l，l为测深数据段的个数，l＝(w-αn)/n(1-α)，α为数据重叠率。
[0022]
可选的，所述基于近似熵理论，确定每个测深数据段的近似熵，具体包括：
[0023]
初始化汉克尔矩阵的嵌入维数为h；
[0024]
对每个测深数据段分别构造所述嵌入维数的汉克尔矩阵；所述汉克尔矩阵中的每
一列为一个模式向量；
[0025]
利用公式确定所述汉克尔矩阵中第p个模式向量与所述汉克尔矩阵的第q个模式向量的最大距离；
[0026]
统计所述汉克尔矩阵中第p个模式向量分别与所述汉克尔矩阵的所有模式向量的最大距离小于或等于模式容错阈值的个数；
[0027]
基于近似概率理论，根据所述个数，利用公式c
ph
(r)＝n
ph
(r)/(n-h+1)，确定第p个模式向量在所述汉克尔矩阵的近似概率值；
[0028]
根据所述近似概率值，利用公式确定每个测深数据段在h嵌入维数下的近似概率平均值；
[0029]
令嵌入维数h增加1，返回步骤“对每个测深数据段分别构造嵌入维数为h的汉克尔矩阵”，确定每个测深数据段在h+1嵌入维数下的近似概率平均值；
[0030]
根据每个测深数据段在h嵌入维数下的近似概率平均值和每个测深数据段在h+1嵌入维数下的近似概率平均值，利用公式确定每个测深数据段的近似熵；
[0031]
其中，d
max
[y(p),y(q)]为汉克尔矩阵中第p个模式向量与汉克尔矩阵的第q个模式向量的最大距离，y(p)为汉克尔矩阵中第p个模式向量，y(q)为汉克尔矩阵的第q个模式向量，x(p+δ)为汉克尔矩阵中第p个模式向量的第δ个元素，x(q+δ)为汉克尔矩阵中第q个模式向量的第δ个元素，δ＝1,2,
…
,h，c
ph
(r)为第p个模式向量在汉克尔矩阵的近似概率值，n
ph
(r)为最大距离小于或等于模式容错阈值的个数，r为模式容错阈值，n为测深数据段长度，为测深数据段在h嵌入维数下的近似概率平均值，为测深数据段在h+1嵌入维数下的近似概率平均值，ae
i
为第i个测深数据段的近似熵。
[0032]
可选的，所述根据所述待去噪测深数据段，构造待去噪矩阵，具体包括：
[0033]
根据预设矩阵维数，将第a个待去噪测深数据段x
a
＝{x
(a-1)*n+1
,x
(a-1)*n+2
,
…
,x
(a-1)*n+n
,}等分为多段；
[0034]
根据等分后的每一段，构造待去噪矩阵
[0035]
其中，x
a
为第a个待去噪测深数据段，x
(a-1)*n+1
、x
(a-1)*n+2
、x
(a-1)*n+n
、x
(a-1)*n+1+n
、x
(a-1)*n+2n
、x
a*n+1-n
和x
(a-1)*n+n
分别表示第a个待去噪测深数据段第1个、第2个、第n个、第n+1个、第2n个、第n-n+1个和第n个数据，n为测深数据段长度，n为预设矩阵长度，n＝n/m，m为预设矩阵维数。
[0036]
可选的，所述对所述待去噪矩阵进行奇异值分解，获得所述待去噪矩阵的奇异值矩阵、左奇异矩阵和右奇异矩阵，之后还包括：
[0037]
确定所述奇异值矩阵中每个奇异值与所有奇异值总和的比值；
[0038]
以所述奇异值矩阵中奇异值的个数为横坐标，以所述比值为纵坐标，绘制奇异值
比重折线图；
[0039]
判断所述奇异值比重折线图的尾支对应的纵坐标是否均小于比值阈值，获得判断结果；
[0040]
若所述判断结果表示是，则输出所述待去噪矩阵的奇异值矩阵、左奇异矩阵和右奇异矩阵；
[0041]
若所述判断结果表示否，则更新所述预设矩阵维数，返回步骤“根据预设矩阵维数，将第a个待去噪测深数据段x
a
＝{x
(a-1)*n+1
,x
(a-1)*n+2
,
…
,x
(a-1)*n+n
,}等分为多段”。
[0042]
可选的，所述将所述去噪后的测深数据段和所述有效测深数据段进行重构，获得重构去噪大地电磁脉动信号，具体包括：
[0043]
利用公式w＝ln-(l-1)
·
αn，将所述去噪后的测深数据段和所述有效测深数据段进行重构，获得重构去噪大地电磁脉动信号；
[0044]
其中，w为大地电磁脉动信号的长度，l为测深数据段的个数，α为数据重叠率，n为测深数据段长度。
[0045]
一种大地电磁去噪系统，所述系统包括：
[0046]
大地电磁脉动信号测量模块，用于利用大地电磁测深法测量大地不同深度的大地电磁脉动信号；
[0047]
测深数据段获得模块，用于对所述大地电磁脉动信号进行重叠分段，获得多个测深数据段；
[0048]
测深数据段的近似熵确定模块，用于基于近似熵理论，确定每个测深数据段的近似熵；
[0049]
待去噪测深数据段确定模块，用于将近似熵小于或等于近似熵阈值的测深数据段确定为待去噪测深数据段，并将近似熵大于近似熵阈值的测深数据段确定为有效测深数据段；
[0050]
待去噪矩阵构造模块，用于根据所述待去噪测深数据段，构造待去噪矩阵；
[0051]
奇异值分解模块，用于对所述待去噪矩阵进行奇异值分解，获得所述待去噪矩阵的奇异值矩阵、左奇异矩阵和右奇异矩阵；
[0052]
场源分量矩阵获得模块，用于将所述奇异值矩阵的每个对角线元素、所述左奇异矩阵对应的列向量和所述右奇异值矩阵对应的行向量相乘，获得所述待去噪测深数据段的多个场源分量矩阵；
[0053]
分量数据段获得模块，用于将每个场源分量矩阵转变为分量数据段；
[0054]
分量数据段的近似熵确定模块，用于基于近似熵理论，确定每个分量数据段的近似熵；
[0055]
去噪后的测深数据段获得模块，用于筛选分量数据段的近似熵大于所述待去噪测深数据段的近似熵的分量数据段，并将筛选出来的分量数据段求和，获得去噪后的测深数据段；
[0056]
重构去噪大地电磁脉动信号获得模块，用于将所述去噪后的测深数据段和所述有效测深数据段进行重构，获得重构去噪大地电磁脉动信号。
[0057]
可选的，所述测深数据段获得模块，具体包括：
[0058]
测深数据段获得子模块，用于对所述大地电磁脉动信号x1,x2,
…
x
w
进行重叠分段
处理，获得每个测深数据段x
v
＝{x
(v-1)*n+1
,x
(v-1)*n+2
,
…
,x
(v-1)*n+n
}；
[0059]
其中，x1,x2,
…
x
w
分别为第1个、第2个、第w个大地电磁脉动信号，w为大地电磁脉动信号的长度，x
(v-1)*n+1
、x
(v-1)*n+2
、x
(v-1)*n+n
分别为第v个测深数据段的第1个、第2个和第n个数据，n为测深数据段长度，v＝1,2，
…
,l，l为测深数据段的个数，l＝(w-αn)/n(1-α)，α为数据重叠率。
[0060]
可选的，所述测深数据段的近似熵确定模块，具体包括：
[0061]
嵌入维数初始化子模块，用于初始化汉克尔矩阵的嵌入维数为h；
[0062]
汉克尔矩阵构造子模块，用于对每个测深数据段分别构造嵌入维数为h的汉克尔矩阵；所述汉克尔矩阵中的每一列为一个模式向量；
[0063]
最大距离确定子模块，用于利用公式确定所述汉克尔矩阵中第p个模式向量与所述汉克尔矩阵的第q个模式向量的最大距离；
[0064]
个数统计子模块，用于统计所述汉克尔矩阵中第p个模式向量分别与所述汉克尔矩阵的所有模式向量的最大距离小于或等于模式容错阈值的个数；
[0065]
近似概率值确定子模块，用于基于近似概率理论，根据所述个数，利用公式c
ph
(r)＝n
ph
(r)/(n-h+1)，确定第p个模式向量在所述汉克尔矩阵的近似概率值；
[0066]
近似概率平均值确定子模块，用于根据所述近似概率值，利用公式确定每个测深数据段在h嵌入维数下的近似概率平均值；
[0067]
循环子模块，用于令嵌入维数h增加1，返回步骤“对每个测深数据段分别构造嵌入维数为h的汉克尔矩阵”，确定每个测深数据段在h+1嵌入维数下的近似概率平均值；
[0068]
测深数据段的近似熵确定子模块，用于根据每个测深数据段在h嵌入维数下的近似概率平均值和每个测深数据段在h+1嵌入维数下的近似概率平均值，利用公式确定每个测深数据段的近似熵；
[0069]
其中，d
max
[y(p),y(q)]为汉克尔矩阵中第p个模式向量与汉克尔矩阵的第q个模式向量的最大距离，y(p)为汉克尔矩阵中第p个模式向量，y(q)为汉克尔矩阵的第q个模式向量，x(p+δ)为汉克尔矩阵中第p个模式向量的第δ个元素，x(q+δ)为汉克尔矩阵中第q个模式向量的第δ个元素，δ＝1,2,
…
,h，c
ph
(r)为第p个模式向量在汉克尔矩阵的近似概率值，n
ph
(r)为最大距离小于或等于模式容错阈值的个数，r为模式容错阈值，n为测深数据段长度，为测深数据段在h嵌入维数下的近似概率平均值，为测深数据段在h+1嵌入维数下的近似概率平均值，ae
i
为第i个测深数据段的近似熵。
[0070]
可选的，所述待去噪矩阵构造模块，具体包括：
[0071]
等分子模块，用于根据预设矩阵维数，将第a个待去噪测深数据段x
a
＝{x
(a-1)*n+1
,x
(a-1)*n+2
,
…
,x
(a-1)*n+n
,}等分为多段；
[0072]
待去噪矩阵构造子模块，用于根据等分后的每一段，构造待去噪矩阵
[0073]
其中，x
a
为第a个待去噪测深数据段，x
(a-1)*n+1
、x
(a-1)*n+2
、x
(a-1)*n+n
、x
(a-1)*n+1+n
、x
(a-1)*n+2n
、x
a*n+1-n
和x
(a-1)*n+n
分别表示第a个待去噪测深数据段第1个、第2个、第n个、第n+1个、第2n个、第n-n+1个和第n个数据，n为测深数据段长度，n为预设矩阵长度，n＝n/m，m为预设矩阵维数。
[0074]
根据本发明提供的具体实施例，本发明公开了以下技术效果：
[0075]
本发明公开了一种大地电磁去噪方法及系统，对大地电磁脉动信号进行重叠分段，获得多个测深数据段，确定每个测深数据段的近似熵，将近似熵小于或等于近似熵阈值的测深数据段确定为待去噪测深数据段，并将近似熵大于近似熵阈值的测深数据段确定为有效测深数据段，在采取svd分解进行信噪分离前，选用近似熵对数据段进行筛选，防止了对有效测深数据段的“过处理”；对待去噪测深数据段进行奇异值分解，将奇异值分解后的奇异值矩阵、左奇异矩阵和右奇异矩阵相乘，获得待去噪测深数据段的多个分量矩阵，每个分量矩阵转变为分量数据段，确定每个分量数据段的近似熵，筛选分量数据段的近似熵大于待去噪测深数据段的近似熵的分量数据段，并将筛选出来的分量数据段求和，获得去噪后的测深数据段，对分解的分量采取近似熵筛选，实现了去噪的目的；最后将去噪后的测深数据段和有效测深数据段进行重构，获得重构去噪大地电磁脉动信号，提高了数据去噪的准确性。
附图说明
[0076]
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动性的前提下，还可以根据这些附图获得其他的附图。
[0077]
图1为本发明提供的提供一种大地电磁去噪方法的流程图；
[0078]
图2为本发明提供的提供一种大地电磁去噪方法的原理图。
具体实施方式
[0079]
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
[0080]
本发明的目的是提供一种大地电磁去噪方法及系统，以减少对数据的过处理，提高数据去噪的准确性。
[0081]
为使本发明的上述目的、特征和优点能够更加明显易懂，下面结合附图和具体实施方式对本发明作进一步详细的说明。
[0082]
本发明为了减少有效信号的损失，联合近似熵筛选与svd分解去噪处理，提供了一
种大地电磁去噪方法，如图1和图2所示，方法包括以下步骤：
[0083]
s101，利用大地电磁测深法测量大地不同深度的大地电磁脉动信号。
[0084]
s102，对大地电磁脉动信号进行重叠分段，获得多个测深数据段。
[0085]
s103，基于近似熵理论，确定每个测深数据段的近似熵。
[0086]
s104，将近似熵小于或等于近似熵阈值的测深数据段确定为待去噪测深数据段，并将近似熵大于近似熵阈值的测深数据段确定为有效测深数据段。优选地，近似熵阈值为ae0，ae0一般在0.5-1之间。
[0087]
s105，根据待去噪测深数据段，构造待去噪矩阵。
[0088]
s106，对待去噪矩阵进行奇异值分解，获得待去噪矩阵的奇异值矩阵、左奇异矩阵和右奇异矩阵。
[0089]
s107，将奇异值矩阵的每个对角线元素、左奇异矩阵对应的列向量和右奇异值矩阵对应的行向量相乘，获得待去噪测深数据段的多个场源分量矩阵。
[0090]
s108，将每个场源分量矩阵转变为分量数据段。
[0091]
s109，基于近似熵理论，确定每个分量数据段的近似熵。
[0092]
s110，筛选分量数据段的近似熵大于待去噪测深数据段的近似熵的分量数据段，并将筛选出来的分量数据段求和，获得去噪后的测深数据段。
[0093]
s111，将去噪后的测深数据段和有效测深数据段进行重构，获得重构去噪大地电磁脉动信号。
[0094]
具体过程如下：
[0095]
步骤s102，对大地电磁脉动信号进行重叠分段，获得多个测深数据段，具体包括：
[0096]
对大地电磁脉动信号x1,x2,
…
x
w
进行重叠分段处理，获得每个测深数据段x
v
＝{x
(v-1)*n+1
,x
(v-1)*n+2
,
…
,x
(v-1)*n+n
}。
[0097]
其中，x1,x2,
…
x
w
分别为第1个、第2个、第w个大地电磁脉动信号，w为大地电磁脉动信号的长度，x
(v-1)*n+1
、x
(v-1)*n+2
、x
(v-1)*n+n
分别为第v个测深数据段的第1个、第2个和第n个数据，n为测深数据段长度，v＝1,2，
…
,l，l为测深数据段的个数，l＝(w-αn)/n(1-α)，α为数据重叠率。
[0098]
由于svd分解的处理对象多维矩阵，采用计算机软件对其进行处理时对运行内存需求较大，而此步骤首先将一维向量划分为多个数据段，每个数据段再构造成多维矩阵，缩短了处理对象的长度，因此此步骤能够保证运行的可行性、数据的连续性以及筛选数据段的准确性。
[0099]
s103步骤，基于近似熵理论，确定每个测深数据段的近似熵，具体包括：
[0100]
初始化汉克尔矩阵(hankel矩阵)的嵌入维数为h；
[0101]
对每个测深数据段分别构造嵌入维数的汉克尔矩阵；汉克尔矩阵中的每一列为一个模式向量；
[0102]
利用公式确定汉克尔矩阵中第p个模式向量与汉克尔矩阵的第q个模式向量的最大距离；
[0103]
统计汉克尔矩阵中第p个模式向量分别与汉克尔矩阵的所有模式向量的最大距离小于或等于模式容错阈值的个数；
[0104]
基于近似概率理论，根据个数，利用公式c
ph
(r)＝n
ph
(r)/(n-h+1)，确定第p个模式向量在汉克尔矩阵的近似概率值；
[0105]
根据近似概率值，利用公式确定每个测深数据段在h嵌入维数下的近似概率平均值；
[0106]
令嵌入维数h增加1，返回步骤“对每个测深数据段分别构造嵌入维数为h的汉克尔矩阵”，确定每个测深数据段在h+1嵌入维数下的近似概率平均值；
[0107]
根据每个测深数据段在h嵌入维数下的近似概率平均值和每个测深数据段在h+1嵌入维数下的近似概率平均值，利用公式确定每个测深数据段的近似熵；
[0108]
其中，d
max
[y(p),y(q)]为汉克尔矩阵中第p个模式向量与汉克尔矩阵的第q个模式向量的最大距离，y(p)为汉克尔矩阵中第p个模式向量，y(q)为汉克尔矩阵的第q个模式向量，x(p+δ)为汉克尔矩阵中第p个模式向量的第δ个元素，x(q+δ)为汉克尔矩阵中第q个模式向量的第δ个元素，δ＝1,2,
…
,h，c
ph
(r)为第p个模式向量在汉克尔矩阵的近似概率值，n
ph
(r)为最大距离小于或等于模式容错阈值的个数，r为模式容错阈值，n为测深数据段长度，为测深数据段在h嵌入维数下的近似概率平均值，为测深数据段在h+1嵌入维数下的近似概率平均值，ae
i
为第i个测深数据段的近似熵。
[0109]
采用近似熵在进行信噪分离前对各段数据进行预筛选，噪声对比有效信号具有较强的方向性、分散程度较集中，因此通过该差异借助近似熵预筛选能够排除处理质量较好的数据，防止当数据段为有效纯净信号时，仍进行筛选分量的去噪处理，而导致的有效信号的“过处理”。
[0110]
步骤s105，根据待去噪测深数据段，构造待去噪矩阵，具体包括：
[0111]
根据预设矩阵维数，将第a个待去噪测深数据段x
a
＝{x
(a-1)*n+1
,x
(a-1)*n+2
,
…
,x
(a-1)*n+n
,}等分为多段；
[0112]
根据等分后的每一段，构造待去噪矩阵
[0113]
其中，x
a
为第a个待去噪测深数据段，x
(a-1)*n+1
、x
(a-1)*n+2
、x
(a-1)*n+n
、x
(a-1)*n+1+n
、x
(a-1)*n+2n
、x
a*n+1-n
和x
(a-1)*n+n
分别表示第a个待去噪测深数据段第1个、第2个、第n个、第n+1个、第2n个、第n-n+1个和第n个数据，n为测深数据段长度，n为预设矩阵长度，n＝n/m，m为预设矩阵维数。
[0114]
步骤s106，svd分解：x
a
＝usv
t
，其中u＝[u1…
u
m
]为m
×
m正交的左奇异矩阵，为奇异值矩阵，主对角线为矩阵的奇异值，为x
at
x
a
的n
×
n右奇异矩阵。u1、u
m
分别为左奇异矩阵中的第1个和第m个元素，s1、s
m
分别表示奇异值矩阵中的第1个和第m个元素，v
1t
、v
1t
分别表示右奇异矩阵中的第1个和第n个元素。
[0115]
svd分解时矩阵构造的维数与长度选取十分关键，通过反复试验得到当矩阵长度小于等于噪声长度时，对噪声的压制效果较好。
[0116]
步骤s106，之后还包括：
[0117]
确定奇异值矩阵中每个奇异值与所有奇异值总和的比值；图2中δ
z
为奇异值矩阵中第z个奇异值，δ
sum
为奇异值矩阵中奇异值总和；
[0118]
以奇异值矩阵中奇异值的个数为横坐标，以比值为纵坐标，绘制奇异值比重折线图；
[0119]
判断奇异值比重折线图的尾支对应的纵坐标是否均小于比值阈值，获得判断结果；优选地，比值阈值为0.001。尾支是指奇异值比重折线图中折线的最后一部分。
[0120]
若判断结果表示是，则输出待去噪测深数据段的奇异值矩阵、左奇异矩阵和右奇异矩阵；
[0121]
若判断结果表示否，则更新预设矩阵维数，返回步骤“根据预设矩阵维数，将第a个待去噪测深数据段x
a
＝{x
(a-1)*n+1
,x
(a-1)*n+2
,
…
,x
(a-1)*n+n
,}等分为多段”。
[0122]
根据奇异值比重折线图，通过观察折线图尾支是否小于0.001，来判断矩阵维数的选取是否正确，若折线图尾支数值大于0.001，则重新赋值矩阵维数，直到比重数值小于0.001。
[0123]
步骤s107，根据各个奇异值与左、右特征向量计算出所有的分量矩阵，将各分量矩阵首尾相连，得到不同场源信息发射的信号：
[0124][0125]
其中，x
a1
、x
am
分别为第a个待去噪测深数据段的第1个和第m个分量矩阵。
[0126]
步骤s108，将场源矩阵中的每个场源分量矩阵转变为分量数据段的步骤为步骤s105“根据待去噪测深数据段，构造待去噪矩阵”的相反的步骤。
[0127]
步骤s109，基于近似熵理论，确定每个分量数据段的近似熵的方法与步骤s103的相同。
[0128]
s110，分量数据段的近似熵大于待去噪测深数据段的近似熵的分量数据段被认为有效数据进行保留，小于该值为噪声数据。将保留的有效分量相加求和，实现了待去噪测深数据段的信噪分离。图2中x
b
表示去噪后的测深数据段，x
bk
表示筛选出来的第k个分量矩阵，k＝1,2，
…
,m。
[0129]
步骤s111，将去噪后的测深数据段和有效测深数据段进行重构，获得重构去噪大地电磁脉动信号，具体包括：
[0130]
利用公式w＝ln-(l-1)
·
αn，将去噪后的测深数据段和有效测深数据段进行重构，获得重构去噪大地电磁脉动信号；
[0131]
其中，w为大地电磁脉动信号的长度，l为测深数据段的个数，α为数据重叠率，n为测深数据段长度。
[0132]
本发明的基本原理为：利用近似熵选定待处理的数据段，并基于天然场源与干扰源的性质差异，利用svd分解将时序信号分解为奇异值矩阵与左右特征矩阵向量，再通过svd逆变换构成来自不同场源信息的多个分量信号，对这些分量信号进行近似熵的计算，有
效信号的近似熵值较大(代表着信号分散程度较大)；而噪声的近似熵值较小(具有一定的极化方向)，以此来筛选分量。
[0133]
因此将近似熵与svd分解相结合既能够减少对数据的过处理，又可以准确地选出有效信号所在的分量来重构去噪数据段。在数据质量较低、实测数据存在形态复杂噪声的情况下，基于近似熵和svd分解的大地电磁去噪技术在能够恢复视电阻率-相位曲线的变化趋势，得到连续平滑的曲线形态。
[0134]
本发明采用近似熵在进行信噪分离前对各段数据进行预筛选，噪声对比有效信号具有较强的方向性、分散程度较集中，因此通过该差异借助预筛选能够防止进一步处理质量较好的数据；并且通过信号幅值之间的差异提取出多个奇异值，每一个奇异值能够计算出包含不同信息的信号分量，通过筛选分量达到去噪的目的，因此该技术与数据的信噪比高低无关，与提取的奇异值个数有关(奇异值越多，分量越多，处理速度越慢；奇异值越少，分量越少，可能会损失有效信号)；同理，本发明不依据形态进行去噪，因此可同时对大尺度噪声和突变噪声进行处理，不受噪声形态的限制。
[0135]
本发明还提供了一种大地电磁去噪系统，系统包括：大地电磁脉动信号测量模块、测深数据段获得模块、测深数据段的近似熵确定模块、待去噪测深数据段确定模块、待去噪矩阵构造模块、奇异值分解模块、场源分量矩阵获得模块、分量数据段获得模块、分量数据段的近似熵确定模块、去噪后的测深数据段获得模块和重构去噪大地电磁脉动信号获得模块。
[0136]
大地电磁脉动信号测量模块，用于利用大地电磁测深法测量大地不同深度的大地电磁脉动信号。
[0137]
测深数据段获得模块，用于对大地电磁脉动信号进行重叠分段，获得多个测深数据段。
[0138]
测深数据段的近似熵确定模块，用于基于近似熵理论，确定每个测深数据段的近似熵。
[0139]
待去噪测深数据段确定模块，用于将近似熵小于或等于近似熵阈值的测深数据段确定为待去噪测深数据段，并将近似熵大于近似熵阈值的测深数据段确定为有效测深数据段。
[0140]
待去噪矩阵构造模块，用于根据待去噪测深数据段，构造待去噪矩阵。
[0141]
奇异值分解模块，用于对待去噪矩阵进行奇异值分解，获得待去噪矩阵的奇异值矩阵、左奇异矩阵和右奇异矩阵。
[0142]
场源分量矩阵获得模块，用于将奇异值矩阵的每个对角线元素、左奇异矩阵对应的列向量和右奇异值矩阵对应的行向量相乘，获得待去噪测深数据段的多个场源分量矩阵。
[0143]
分量数据段获得模块，用于将每个场源分量矩阵转变为分量数据段。
[0144]
分量数据段的近似熵确定模块，用于基于近似熵理论，确定每个分量数据段的近似熵。
[0145]
去噪后的测深数据段获得模块，用于筛选分量数据段的近似熵大于待去噪测深数据段的近似熵的分量数据段，并将筛选出来的分量数据段求和，获得去噪后的测深数据段。
[0146]
重构去噪大地电磁脉动信号获得模块，用于将去噪后的测深数据段和有效测深数
据段进行重构，获得重构去噪大地电磁脉动信号。
[0147]
测深数据段获得模块，具体包括：
[0148]
测深数据段获得子模块，用于对大地电磁脉动信号x1,x2,
…
x
w
进行重叠分段处理，获得每个测深数据段x
v
＝{x
(v-1)*n+1
,x
(v-1)*n+2
,
…
,x
(v-1)*n+n
}。
[0149]
其中，x1,x2,
…
x
w
分别为第1个、第2个、第w个大地电磁脉动信号，w为大地电磁脉动信号的长度，x
(v-1)*n+1
、x
(v-1)*n+2
、x
(v-1)*n+n
分别为第v个测深数据段的第1个、第2个和第n个数据，n为测深数据段长度，v＝1,2，
…
,l，l为测深数据段的个数，l＝(w-αn)/n(1-α)，α为数据重叠率。
[0150]
测深数据段的近似熵确定模块，具体包括：
[0151]
嵌入维数初始化子模块，用于初始化汉克尔矩阵的嵌入维数为h。
[0152]
汉克尔矩阵构造子模块，用于对每个测深数据段分别构造嵌入维数为h的汉克尔矩阵。汉克尔矩阵中的每一列为一个模式向量。
[0153]
最大距离确定子模块，用于利用公式确定汉克尔矩阵中第p个模式向量与汉克尔矩阵的第q个模式向量的最大距离。
[0154]
个数统计子模块，用于统计汉克尔矩阵中第p个模式向量分别与汉克尔矩阵的所有模式向量的最大距离小于或等于模式容错阈值的个数。
[0155]
近似概率值确定子模块，用于基于近似概率理论，根据个数，利用公式c
ph
(r)＝n
ph
(r)/(n-h+1)，确定第p个模式向量在汉克尔矩阵的近似概率值。
[0156]
近似概率平均值确定子模块，用于根据近似概率值，利用公式确定每个测深数据段在h嵌入维数下的近似概率平均值。
[0157]
循环子模块，用于令嵌入维数h增加1，返回步骤“对每个测深数据段分别构造嵌入维数为h的汉克尔矩阵”，确定每个测深数据段在h+1嵌入维数下的近似概率平均值。
[0158]
测深数据段的近似熵确定子模块，用于根据每个测深数据段在h嵌入维数下的近似概率平均值和每个测深数据段在h+1嵌入维数下的近似概率平均值，利用公式确定每个测深数据段的近似熵。
[0159]
其中，d
max
[y(p),y(q)]为汉克尔矩阵中第p个模式向量与汉克尔矩阵的第q个模式向量的最大距离，y(p)为汉克尔矩阵中第p个模式向量，y(q)为汉克尔矩阵的第q个模式向量，x(p+δ)为汉克尔矩阵中第p个模式向量的第δ个元素，x(q+δ)为汉克尔矩阵中第q个模式向量的第δ个元素，δ＝1,2,
…
,h，c
ph
(r)为第p个模式向量在汉克尔矩阵的近似概率值，n
ph
(r)为最大距离小于或等于模式容错阈值的个数，r为模式容错阈值，n为测深数据段长度，为测深数据段在h嵌入维数下的近似概率平均值，为测深数据段在h+1嵌入维数下的近似概率平均值，ae
i
为第i个测深数据段的近似熵。
[0160]
待去噪矩阵构造模块，具体包括：
[0161]
等分子模块，用于根据预设矩阵维数，将第a个待去噪测深数据段x
a
＝{x
(a-1)*n+1
,x
(a-1)*n+2
,
…
,x
(a-1)*n+n
,}等分为多段。
[0162]
待去噪矩阵构造子模块，用于根据等分后的每一段，构造待去噪矩阵
[0163][0164]
其中，x
a
为第a个待去噪测深数据段，x
(a-1)*n+1
、x
(a-1)*n+2
、x
(a-1)*n+n
、x
(a-1)*n+1+n
、x
(a-1)*n+2n
、x
a*n+1-n
和x
(a-1)*n+n
分别表示第a个待去噪测深数据段第1个、第2个、第n个、第n+1个、第2n个、第n-n+1个和第n个数据，n为测深数据段长度，n为预设矩阵长度，n＝n/m，m为预设矩阵维数。
[0165]
本发明基于近似熵和svd分解的大地电磁去噪技术在一定情况下能够减少“过处理”数据的现象；对数据信噪比要求较低；能够压制复杂形态的噪声。
[0166]
本说明书中各个实施例采用递进的方式描述，每个实施例重点说明的都是与其他实施例的不同之处，各个实施例之间相同相似部分互相参见即可。对于实施例公开的系统而言，由于其与实施例公开的方法相对应，所以描述的比较简单，相关之处参见方法部分说明即可。
[0167]
本文中应用了具体个例对本发明的原理及实施方式进行了阐述，以上实施例的说明只是用于帮助理解本发明的方法及其核心思想；同时，对于本领域的一般技术人员，依据本发明的思想，在具体实施方式及应用范围上均会有改变之处。综上所述，本说明书内容不应理解为对本发明的限制。