本公开涉及一种伪随机编码磁性源瞬变电磁系统获取大地电导率等信息底层深度、金属物体的位置的数据预处理方法，尤其是涉及一种基于伪随机编码阈值优化的相关辨识数据预处理方法，属于时间域磁性源伪随机编码瞬变电磁资料处理技术领域。
背景技术：
伪随机编码磁性源瞬变电磁系统的发射、接收装置均位于大地表面，基于瞬变电磁法原理对采集的数据进行预处理，从而可获得大地电导率等信息；可广泛应用于海洋、极地、高峰等含有的金属矿产和油气资源的复杂地形，具有节约成本，保护生态环境等优点。
通过对电流采集设备和响应接收设备采集的响应波形进行预处理，可以得到大地脉冲响应估计值，可较好的反映观测区域的地电信息。但由于电流采集设备带宽有限，只能采集一定频率范围内的信号导致出现频率截断，频率截断引发时域信号产生吉布斯效应，从而无法得到较好的预处理结果。由于吉布斯效应的本质是频率缺失引起的，无法通过滤波法对电流进行很好的优化。
因此，目前亟待解决的技术问题是，如何消除电流中的吉布斯震荡。
技术实现要素：
为了消除时间域磁性源伪随机编码瞬变电磁系统中电流中的吉布斯震荡，本公开提供一种基于伪随机编码阈值优化的相关辨识数据预处理方法，可以有效的消除电流中的吉布斯震荡，具体方案如下。
一种基于伪随机编码阈值优化的相关辨识数据预处理方法，其特征在于，包括：
s1，获取伪随机编码发射电流的电流数据以及所述伪随机编码发射电流激发地下良导体产生的二次磁场响应数据；
s2，选取阈值，对所述电流数据进行阈值法优化；
s3，判断阈值法优化后的伪随机编码发射电流的电流数据中是否存在振荡，若是，则返回步骤s2，重新选择新的阈值，对所述电流数据进行阈值法优化，若否，则进入下一步骤；
s4，对所述二次磁场响应的数据进行降采样操作，所述二次磁场响应数据的采样频率等于阈值法优化后的伪随机编码发射电流的电流数据的采样频率；
s5，对阈值法优化后的伪随机编码发射电流数据对时间求导，得到电流导数，对所述电流导数和所述降采样后的二次磁场响应数据进行相关辨识法运算，得到大地脉冲响应估计值。
根据本公开的一些实施例，所述s2中，包括：
s21，将所述伪随机编码发射电流的电流数据进行分组，每n个点一组，选取每组中绝对值最大的元素组成新的电流序列，其中，n为大于1的整数；
s22，提取所述新的电流序列中的最大值和最小值，选取阈值，将所述新的电流序列中的元素逐一与所述阈值比对，若该元素大于等于所述阈值，则调整该元素为所述最大值，若该元素小于所述阈值，则调整该元素为所述最小值，得到阈值法优化后的伪随机编码发射电流的电流数据。
根据本公开的一些实施例，根据所述新的电流序列的波形图选取阈值，其中，经过所述阈值优化后的所述伪随机编码发射电流的电流数据的码元宽度与所述新的电流序列的码元宽度误差小于20％。
根据本公开的一些实施例，所述s3中，根据阈值法优化后的伪随机编码发射电流的电流数据的波形图判断是否存在振荡，或根据阈值法优化后的伪随机编码发射电流的电流数据的自相关函数旁瓣判断是否存在振荡。
根据本公开的一些实施例，所述s5中，所述相关辨识法运算还包括卷积运算，所述电流导数与降采样后的二次磁场响应数据的卷积关系包括：
其中，u(t)为降采样后的二次磁场响应数据，i(t)为电流波形，g(t)为大地系统磁场脉冲响应，n(t)为非相干性噪声，
根据本公开的一些实施例，所述相关辨识法运算包括，将公式(1)等式两侧对电流对时间的导数
rui(t)＝rii(t)*g(t) rni(t)(2)
其中，rui(t)表示u(t)与
根据本公开的一些实施例，对公式(3)进行优化得到：
rui(t)＝rii(t)*g(t)(3)
其中，噪声n(t)与
根据本公开的一些实施例，还包括对公式(3)进行离散化计算，包括：
其中，tg表示系统响应g(t)持续的时间。
根据本公开的一些实施例，tg表示大地系统磁场脉冲响应g(t)持续的时间，设置采样频率为fs，则采样时间为dτ＝1/fs，g(t)在tg时间段内的点数为
rui≈riig(5)
对公式(6)进行伪逆矩阵求解得到大地脉冲响应估计值，包括：
通过上述技术方案，本公开通过对伪随机编码发射电流的电流数据进行阈值法优化，配合对优化后的电流数据及处理后的二次磁场响应数据进行相关辨识法运算，可以有效地消除电流中的吉布斯震荡，得到大地脉冲响应估计值。
附图说明
图1示意性示出了本公开实施例的基于伪随机编码阈值优化的相关辨识数据预处理方法的流程图；
图2示意性示出了本公开实施例的伪随机编码发射电流的电流中的吉布斯震荡；
图3示意性示出了本公开实施例的基于伪随机编码阈值优化的相关辨识数据预处理方法的码元频率为256hz的伪随机编码电流数据图；
图4示意性示出了本公开实施例的基于伪随机编码阈值优化的相关辨识数据预处理方法的分组降采样后的电流波形图；
图5示意性示出了本公开实施例的基于伪随机编码阈值优化的相关辨识数据预处理方法的阈值法优化后的电流波形图；
图6示意性示出了本公开实施例的基于伪随机编码阈值优化的相关辨识数据预处理方法的二次磁场响应的波形图；
图7示意性示出了本公开实施例未优化电流数据的大地电磁响应估计值；
图8示意性示出了本公开实施例的基于伪随机编码阈值优化的相关辨识数据预处理方法的优化后大地电磁响应估计值。
具体实施方式
为使本发明的目的、技术方案和优点更加清楚明白，以下结合具体实施例，并参照附图，对本发明作进一步的详细说明。
但是应该理解，这些描述只是示例性的，而并非要限制本发明的范围。在下面的详细描述中，为便于解释，阐述了许多具体的细节以提供对本发明实施例的全面理解。然而，明显地，一个或多个实施例在没有这些具体细节的情况下也可以被实施。此外，在以下说明中，省略了对公知技术的描述，以避免不必要地混淆本发明的概念。
在此使用的术语仅仅是为了描述具体实施例，而并非意在限制本发明。在此使用的术语“包括”表明了特征、步骤、操作的存在，但是并不排除存在或添加一个或多个其他特征。
在此使用的所有术语(包括技术和科学术语)具有本领域技术人员通常所理解的含义，除非另外定义。应注意，这里使用的术语应解释为具有与本说明书的上下文相一致的含义，而不应以理想化或过于刻板的方式来解释。例如，瞬变电磁法，是利用不接地回线或接地线源向地下发射一次脉冲磁场，在一次脉冲磁场间歇期间利用线圈或接地电极观测地下介质中引起的二次感应涡流场，从而探测介质电阻率的一种方法；例如，相关辨识法，是相关分析方法的一种，对总体中确实具有联系的标志进行分析，其主体是对总体中具有因果关系标志的分析，具体的包括自相关和互相关。
一种基于伪随机编码阈值优化的相关辨识数据预处理方法，图1示意性示出了本公开实施例的基于伪随机编码阈值优化的相关辨识数据预处理方法的流程图。
如图1所示，一种基于伪随机编码阈值优化的相关辨识数据预处理方法，包括如下步骤。
s1，获取伪随机编码发射电流的电流数据以及伪随机编码发射电流激发地下良导体产生的二次磁场响应数据。可选的，设备于地面或空中通以一伪随机编码电流的发射线圈，并在地下导电岩矿体中产生的二次磁场响应，记录这个伪随机编码电流数据，以及下导电岩矿体中产生的二次磁场响应数据。
s2，选取阈值，对电流数据进行阈值法优化。其目的在于通过阈值法优化筛除电流数据中的吉布斯震荡。
s3，判断阈值法优化后的伪随机编码发射电流的电流数据中是否存在振荡，若是，则返回步骤s2，重新选择新的阈值，对电流数据进行阈值法优化，若否，则进入下一步骤。其目的在于，判断上述优化是否消除了吉布斯震荡，若消除了，则进行下一步，若未消除，则返回上一步重新选取阈值，进行阈值法优化。
s4，对二次磁场响应的数据进行降采样操作，降采样后的二次磁场响应数据的采样频率等于阈值法优化后的伪随机编码发射电流的电流数据的采样频率。其目的在于，保持二次磁场响应数据的采样频率与阈值法优化后的伪随机编码发射电流的电流数据的采样频率一致，便于后面的相关辨识法运算。
s5，对阈值法优化后的伪随机编码发射电流数据对时间求导，得到电流导数，对电流导数和降采样后的二次磁场响应数据进行相关辨识法运算，得到大地脉冲响应估计值。
根据本公开的一些实施例，s2中，包括如下子步骤。
s21，将伪随机编码发射电流的电流数据进行分组，每n个点一组，选取每组中绝对值最大的元素组成新的电流序列，其中，n为大于1的整数。
s22，提取新的电流序列中的最大值和最小值，选取阈值，将新的电流序列中的元素逐一与阈值比对，若该元素大于等于阈值，则调整该元素为最大值，若该元素小于阈值，则调整该元素为最小值，得到阈值法优化后的伪随机编码发射电流的电流数据。
根据本公开的一些实施例，根据新的电流序列的波形图选取阈值，其中，经过阈值优化后的伪随机编码发射电流的电流数据的码元宽度与新的电流序列的码元宽度误差小于20％。如4或图5所示，两条纵线之间的宽度即为码元宽度，其中，选取阈值后码元宽度越小越好。
根据本公开的一些实施例，s3中，根据阈值法优化后的伪随机编码发射电流的电流数据的波形图判断是否存在振荡，具体的，如图3中右上角所示电流波形，即为吉布斯震荡波形。
根据本公开的一些实施例，s5中，相关辨识法运算还包括卷积运算，电流导数与降采样后的二次磁场响应数据的卷积关系包括：
其中，u(t)为降采样后的二次磁场响应数据，i(t)为电流波形，g(t)为大地系统磁场脉冲响应，n(t)为非相干性噪声，
根据本公开的一些实施例，相关辨识法运算包括，将公式(1)等式两侧对电流对时间的导数
rui(t)＝rii(t)*g(t) rni(t)(2)
其中，rui(t)表示u(t)与
根据本公开的一些实施例，对公式(2)进行优化得到：
rui(t)＝rii(t)*g(t)(3)
其中，噪声n(t)与
根据本公开的一些实施例，还包括对公式(3)进行离散化计算，包括：
其中，tg表示系统响应g(t)持续的时间。
根据本公开的一些实施例，tg表示大地系统磁场脉冲响应g(t)持续的时间，设置采样频率为fs，则采样时间为dτ＝1/fs，g(t)在tg时间段内的点数为
rui≈riig(5)
对公式(5)进行伪逆矩阵求解得到大地脉冲响应估计值，包括：
下面结合具体实施例对本公开的技术方案进行阐述，应当理解的是，这些具体实施例仅仅是为了更好地说明本公开的技术方案，而不应该理解为对本公开技术方案保护范围的限制，具体如下。
根据本公开的一些实施例，如图2和图3所示，电流数据采集到的伪随机编码发射电流数据图中，可以明显的观测到多个吉布斯震荡(图3右上角的波形即为吉布斯震荡)。其中，图3中伪随机电流的码元频率为256hz。
根据本公开的一些实施例，电流数据采样频率为192khz。将电流数据每10个点分为一组，筛取每组中绝对值最大的点的数据构成新的电流序列，此时采样频率降低为19.2khz，新的电流序列的电流波形图如图4所示，通过观测(可辅助于网格线)可以得出，阈值选取值在图4中两条横向虚线中间位置的时候，可使得优化后的码元宽度与原电流码元宽度大致相同，试选取阈值为128，进行阈值法优化。
根据本公开的一些实施例，如图4所示，新的电流序列中，最小值为-239，最大值为258，将新的电流序列中的元素逐一与阈值128比对，将大于或等于128的元素修正为258，小于128的元素修正为-239，即可得到图5中的阈值法优化后的电流波形图。
根据本公开的一些实施例，将二次磁场响应数据采样频率fs调整为19.2khz，进行分组降采，用于与优化后电流数据采样频率保持一致，即可得到如图6所示的二次磁场响应的波形图。
根据本公开的一些实施例，对阈值法优化后的伪随机编码发射电流数据对时间求导，得到电流导数，对电流导数和降采样后的二次磁场响应数据进行相关辨识法运算，得到大地脉冲响应估计值，其波形图如图8所示。此外，图7为未优化电流数据的大地电磁响应估计值的波形图，通过图7和图8的对比可以得出，本公开提供的基于伪随机编码阈值优化的相关辨识数据预处理方法可以消除电流中的吉布斯震荡，得到更为准确的大地脉冲响应估计值。
通过上述技术方案，本公开通过对伪随机编码发射电流的电流数据进行阈值法优化，配合对优化后的电流数据及处理后的二次磁场响应数据进行相关辨识法运算，可以有效地消除电流中的吉布斯震荡，得到大地脉冲响应估计值。
至此，已经结合附图对本公开实施例进行了详细描述。需要说明的是，在附图或说明书正文中，未绘示或描述的实现方式，均为所属技术领域中普通技术人员所知的形式，并未进行详细说明。此外，上述对各零部件的定义并不仅限于实施例中提到的各种具体结构、形状或方式，本领域普通技术人员可对其进行简单地更改或替换。
还需要说明的是，在本公开的具体实施例中，除非有所知名为相反之意，本说明书及所附权利要求中的数值参数是近似值，能够根据通过本公开的内容所得的所需特性改变。具体而言，所有使用于说明书及权利要求中表示组成的尺寸、范围条件等等的数字，应理解为在所有情况中是受到“约”的用语所修饰。一般情况下，其表达的含义是指包含由特定数量在一些实施例中±10％的变化、在一些实施例中±5％的变化、在一些实施例中±1％的变化、在一些实施例中±0.5％的变化。
本领域技术人员可以理解，本发明的各个实施例和/或权利要求中记载的特征可以进行多种组合或/或结合，即使这样的组合或结合没有明确记载于本发明中。特别地，在不脱离本发明精神和教导的情况下，本发明的各个实施例和/或权利要求中记载的特征可以进行多种组合和/或结合。所有这些组合和/或结合均落入本发明的范围。
以上所述的具体实施例，对本发明的目的、技术方案和有益效果进行了进一步详细说明，应理解的是，以上所述仅为本发明的具体实施例而已，并不用于限制本发明，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。
1.一种基于伪随机编码阈值优化的相关辨识数据预处理方法，其特征在于，包括：
s1，获取伪随机编码发射电流的电流数据以及所述伪随机编码发射电流激发地下良导体产生的二次磁场响应数据；
s2，选取阈值，对所述电流数据进行阈值法优化；
s3，判断阈值法优化后的伪随机编码发射电流的电流数据中是否存在振荡，若是，则返回步骤s2，重新选择新的阈值，对所述电流数据进行阈值法优化，若否，则进入下一步骤；
s4，对所述二次磁场响应的数据进行降采样操作，所述二次磁场响应数据的采样频率等于阈值法优化后的伪随机编码发射电流的电流数据的采样频率；
s5，对阈值法优化后的伪随机编码发射电流数据对时间求导，得到电流导数，对所述电流导数和所述降采样后的二次磁场响应数据进行相关辨识法运算，得到大地脉冲响应估计值。
2.根据权利要求1所述的基于伪随机编码阈值优化的相关辨识数据预处理方法，其特征在于，所述s2中，包括：
s21，将所述伪随机编码发射电流的电流数据进行分组，每n个点一组，选取每组中绝对值最大的元素组成新的电流序列，其中，n为大于1的整数；
s22，提取所述新的电流序列中的最大值和最小值，选取阈值，将所述新的电流序列中的元素逐一与所述阈值比对，若该元素大于等于所述阈值，则调整该元素为所述最大值，若该元素小于所述阈值，则调整该元素为所述最小值，得到阈值法优化后的伪随机编码发射电流的电流数据。
3.根据权利要求2所述的基于伪随机编码阈值优化的相关辨识数据预处理方法，其特征在于，根据所述新的电流序列的波形图选取阈值，其中，经过所述阈值优化后的所述伪随机编码发射电流的电流数据的码元宽度与所述新的电流序列的码元宽度误差小于20％。
4.根据权利要求1所述的基于伪随机编码阈值优化的相关辨识数据预处理方法，其特征在于，所述s3中，根据阈值法优化后的伪随机编码发射电流的电流数据的波形图判断是否存在振荡，或根据阈值法优化后的伪随机编码发射电流的电流数据的自相关函数旁瓣判断是否存在振荡。
5.根据权利要求1至4任一所述的基于伪随机编码阈值优化的相关辨识数据预处理方法，其特征在于，所述s5中，所述相关辨识法运算还包括卷积运算，所述电流导数与降采样后的二次磁场响应数据的卷积关系包括：
其中，u(t)为降采样后的二次磁场响应数据，i(t)为电流波形，g(t)为大地系统磁场脉冲响应，n(t)为非相干性噪声，
6.根据权利要求5所述的基于伪随机编码阈值优化的相关辨识数据预处理方法，其特征在于，所述相关辨识法运算包括，将公式(1)等式两侧对电流对时间的导数进行相关运算得到：
rui(t)＝rii(t)*g(t) rni(t)(2)
其中，rui(t)表示u(t)与
7.根据权利要求6所述的基于伪随机编码阈值优化的相关辨识数据预处理方法，其特征在于，对公式(2)进行优化得到：
rui(t)＝rii(t)*g(t)(3)
其中，噪声n(t)与
8.根据权利要求7所述的基于伪随机编码阈值优化的相关辨识数据预处理方法，其特征在于，还包括对公式(3)进行离散化计算，包括：
其中，tg表示系统响应g(t)持续的时间。
9.根据权利要求8所述的基于伪随机编码阈值优化的相关辨识数据预处理方法，其特征在于，tg表示大地系统磁场脉冲响应g(t)持续的时间，设置采样频率为fs，则采样时间为dτ＝1/fs，g(t)在tg时间段内的点数为
rui≈riig(5)
对公式(6)进行伪逆矩阵求解得到大地脉冲响应估计值，包括：
