地球物理学报  2016, Vol. 59 Issue (5): 1861-1874   PDF    
引用本文
王显祥, 底青云, 王妙月等. 2016. 基于m伪随机序列的电磁法抗噪能力分析. 地球物理学报, 59(5): 1861-1874, doi: 10.6038/cjg20160529   
Wang X X, Di Q Y, Wang M Y, et al. 2016. A study on the noise immunity of electromagnetic methods based on m pseudo-random sequence. Chinese J. Geophysics, 59(5): 1861-1874, doi: 10.6038/cjg20160529   
基于m伪随机序列的电磁法抗噪能力分析
王显祥1,2, 底青云2, 王妙月2, 邓居智1    
1. 东华理工大学放射性地质与勘探技术国防重点学科实验室, 南昌 330013;
2. 中国科学院地质与地球物理研究所 中国科学院页岩气与地质工程重点实验室, 北京 100029
收稿日期 2015-05-16, 2016-01-18 收修定稿
基金项目: 国家重大科研装备研制项目"深部资源探测核心装备研发"(ZDYZ2012-1)-05子项目"多通道大功率电法勘探仪",国防科工委科研项目(科工技[2013]969号)联合资助.
作者简介: 王显祥,男,山东巨野人,1986年生,现东华理工大学讲师,主要从事电磁法的正反演工作.E-mail:wangxianxiang09@163.com
摘要: 电磁法在金属矿、油气资源等矿产的勘探中发挥重要作用,同时电磁干扰也越来越强,在强干扰区利用传统的抗噪方法很难再获取高信噪比信号.鉴于伪随机序列的抗噪能力,近年来国内正在掀起伪随机编码仪研制的高潮.伪随机序列的抗干扰能力与从伪随机响应中提取方法息息相关,本文介绍了两种最常用的提取大地脉冲信号的方法:第一种通过将时间域信号转到频率域,对电流归一化处理后,再转回时间域求取大地脉冲响应;第二种直接在时间域解Wiener-Hopf方程求取.研究表明,在无干扰情况下两种提取方法都可获得高精度大地脉冲信号.随后在具有噪声干扰的情况下,对两种方法提取结果的精度进行了对比分析,表明方法2的抗噪声能力明显优于方法1的,对造成这种状况的原因进行了解析.最后研究了提高伪随机序列抗噪能力的方法,主要包括提高伪随机序列阶数、针对不同噪音使用不同的编码方式、改变循环次数的方法等,以达到抑制噪声的目的.
关键词: 伪随机序列     电磁法     抗噪分析    
A study on the noise immunity of electromagnetic methods based on m pseudo-random sequence
WANG Xian-Xiang1,2, DI Qing-Yun2, WANG Miao-Yue2, DENG Ju-Zhi1     
1. Fundamental Science on Radioactive Geology and Exploration Technology Laboratory, East China Institute of Technology, Nanchang 330013, China;
2. Key Laboratory of Shale Gas and Geoengineering, Institute of Geology and Geophysics, Chinese Academy of Sciences, Beijing 100029, China
Abstract: Electromagnetic methods play an important role in the exploration of metal ore mineral, oil and gas resources. It is very difficult to obtain a high signal to noise ratio with increasingly strong electromagnetic noise using conventional anti-noise methods. In recent years, the pseudo random sequence attracts a lot of attention because of its anti-noise performance. This paper describes two common methods used to extract the impulse response from the pseudo-random sequence response. Studies have shown that the two methods are capable to get high-precision impulse response without noise. The extraction accuracies of the results using the two methods were compared with data contained in the noise, that shows that the results of method 2 is better than 1, then the reasons for this situation was analuzed. Finally, a study is made to improve the ability of the anti- pseudo-random sequence of the methods, including to improve the order of the pseudo-random sequence, use a different encoding, change the number of cycles in order to achieve the purpose of noise suppression.
Key words: Pseudo random sequence     Electromagnetic method     Anti-noise analysis    
1 引言

电磁法勘探具有探测深度深、分辨率高、方便快捷、成本低的特点，已成为资源探测领域不可或缺的一种地球物理方法.近年来该方法在石油、金属矿、天然气、地热、地下水的勘探中发挥着越来越重要的作用(底青云等，2002，2006，2008；何继善，1990；陆其鹄和孙进忠，2004；石昆法，1999；汤井田和何继善，2005).

然而，在电磁勘探中电磁信号强度随着收发距的增大迅速衰减，在接收端电场信号强度一般只能达到mV/km这个级别.另一方面，随着经济社会和人文活动发展，电磁干扰信号越来越强，特别是在危机矿山，工业电、民用电、机械干扰、随机干扰、尖峰干扰等各种干扰混杂在一起.为了克服干扰信号的影响，传统的做法主要通过发射端增加发射功率，接受端仪器陷波、数字滤波和多次叠等方法实现对干扰信号的压制.当干扰信号不是很强时可以取得不错的效果，但在强干扰区利用现有的仪器和方法很难再获得高信噪比的电磁信号.针对此种情况，电法工作者将伪随机序列引入电法勘探中，利用发射端信号与接收端信号相关提取大地电磁信号，在抑制干扰信号方面取得了不错效果.

加拿大多伦多大学的Duncan等人(1980)将m伪随机序列应用到电法勘探中，并通过选用合适的频带宽带和收发距，分别对浅层500 m和深层40 km的目标体进行了探测.2002年前后，Ziolkowski A、Wright D和Bruce A H等人将m伪随机序列应用到MTEM法中，在油气资源的勘探中取得了不错的效果(Wright et al.，2001，2002，2006；Wright，2003; Ziolkowski and Wright，2007; Ziolkowski，2007，2010).赵璧如等人(2002)成功研制了伪随机样机PS100，并开展了相关电法工作.2008年汤井田、罗维斌完成了基于逆重复m序列伪随机电磁法的相关研究，讨论了编码电磁测深原理及参数提取方法，提出了在频域辨识大地系统传输函数的方法(汤井田和罗维斌，2008；罗维斌等，2012).罗先中等(2014)开始了伪随机码KGR-1仪的研制工作，并开展了相关的水槽实验和野外勘探工作，取得了不错的效果.Velikin(2009)针对Duncan等(1980)方案中脉冲发生器引进的误差带来的分辨率降低的问题，提出了改进方案，提高了探测精度.Ziolkowski(2013)利用Wiener 滤波估计大地脉冲响应，并应用到地震和电磁勘探当中，取得了不错效果，但Wiener滤波法的实质仍然是收发互相关的方法.Ilyichev和Bobrovsky(2015)从数值模拟和野外实验的角度分别研究了m伪随机编码的抗噪能力，研究表明伪随机编码的抗噪能力较传统方法可以使信噪比提升100倍左右.武欣等(2015)提出一种通过收发互相关高精度提取大地冲激响应的数学方法，同时对以m序列为发射波形的勘探系统相关参数选择进行了分析.

在前人研究的基础上，本文将较系统地开展提高m伪随机序列抗干扰能力的研究，包括系统地总结了从伪随机响应中提取电磁响应的方法，并对比了不同反褶积算法的抗干扰能力；系统地研究了提高伪随机序列抗噪能力的方法，包括提高伪随机序列阶数、针对不同噪音使用不同的编码方式、改变循环次数的方法等.

2 求取大地脉冲响应主要方法简介

当利用伪随机序列采集电磁场数据时，需要将采集到的伪随机信号转换到时间域或频率域，一方面是由于伪随机序列激发的场并不能直观地反映地电结构的变化，给解释工作带来了诸多的不便；另一方面用伪随机序列激发的场直接反演虽然在理论是可行的，但无论数据叠加还是反演过程中都会带来很多不便，并且会把噪声带入反演过程中.通常的做法是将伪随机序列激发的场转到时间域或者频率域中进行处理，不管求取时间域还是频率域信号思路都是一致的.

伪随机序列激发场在转换到频率或时间域的过程中，也是消除干扰信号的重要过程.为此，本节将首先介绍两种比较流行的求取大地脉冲响应的方法，第一种方法将发射端信号与接收端时域信号分别转到频率域，通过在频率域上的一系列处理取得单位电流的电磁响应，再将电磁响应转换到时间域获得大地脉冲响应(Ziolkowski et al.，2007)，推导过程详见附录A；第二种方法直接在时间域进行处理，主要是通过解Wiener-Hopf方程获取大地脉冲响应(Ziolkowski，2013)，推导过程详见附录B.

为了检验两种方法提取大地脉冲响应的效果，以均匀半空间模型为例进行说明，模型电阻率为100 Ωm，收发距为1000 m，m序列阶数为8阶，具体参数如表 1所示.伪随机序列电流如图 1a所示，伪随机序列激发的场可以通过下式求取(殷长春等，2013)：

式中，*表示褶积，I(t)为伪随机发射电流波形，g(t)为大地脉冲响应，Es(t)为大地阶跃场响应.均匀半空间下，上阶跃场响应Es(t)表达式为(方文藻等，1993)
式中，，I为电流强度，dl为电偶极源长度，ρ为电阻率，r为收发距，erf为误差函数，μ0为磁导率，t为时间.利用(2)式求取上阶跃场，本模型中求取的上阶跃场响应如图 1b所示，通过数值拟合求电流导数，经过(1)式褶积求取伪随机序列激发的场，如图 1c所示.利用上文中介绍的方法1与方法2分别从伪随机序列激发的场中提取大地脉冲响应，提取结果如图 1d所示，通过分析图 1d可知：从早期到晚期，方法1与方法2获得的大地脉冲响应与解析解得到的理论脉冲响应三者吻合度都较高，表明在不存在噪音干扰时方法1与方法2都是正确、可行的.
表 1 均匀半空间模型伪随机模拟参数 Table 1 The pseudo-random simulation parameters for uniform half-space model

图 1 伪随机序列响应模拟及两种提取方法结果对比
(a) 伪随机序列电流； (b) 均匀半空间阶跃场响应； (c) 伪随机序列激发的场； (d) 通过两种方法提取的大地脉冲响应. Fig. 1 Simulation of pseudo-random sequence response and results comparison of two methods
(a) Pseudo-random sequence current； (b) The step response； (c) Pseudo-random sequence response； (d) Results comparison of two methods.
3 两种方法抗噪能力分析

为了比较两种方法在有干扰情况下提取大地脉冲响应的精度，仍沿用表 1中的参数模拟伪随机响应，但循环次数由2次变为5次，利用(1)、(2)式求取的伪随机响应场如图 2a所示.在伪随机序列激发的场中加入振幅为10-9的50 Hz干扰，即10-9cos(2π×50t)，其中t为时间.50 Hz干扰信号如图 2b所示.通过方法1和方法2分别求取的大地脉冲响应曲线如图 2c所示，虚线为解析解，标有下三角的曲线为通过方法1得到的大地脉冲场，实线为由方法2得到的大地脉冲响应，比较方法1与方法2提取的大地脉冲响应可知，由方法2提取的大地脉冲响应更加接近理论脉冲响应曲线.为了更好地量化提取结果与理论脉冲响应之间的差别，定义误差率如下：

式中，E2(t)为t时刻通过反褶积得到的大地脉冲响应，E1(t)为t时刻理论大地脉冲响应.由两种方法提取得到的大地脉冲响应误差率曲线如图 2d所示，由图 2d可知由方法2得到的大地脉冲响应误差率控制在10%以内，而由方法1得到的大地脉冲响应误差率一般都超过10%，在晚期段甚至达到40%.通过以上分析可知，在抗干扰能力方面，方法2要明显优于方法1.下面将分析造成这种结果的原因，在分析的过程中也是对伪随机序列抗干扰原因的解析.
图 2 加入50 Hz干扰后伪随机响应及两种提取方法对比
(a) 伪随机序列激发的场； (b) 50 Hz干扰信号； (c) 通过两种方法获得大地脉冲响应； (d) 误差比. Fig. 2 Simulation of pseudo-random sequence response with 50 Hz noise and results comparison of two methods
(a) Pseudo-random sequence response； (b) The 50 Hz noise； (c) Results comparison of two methods； (d) Error rate.
3.1 方法1不抗噪原因分析

方法1在处理过程中，会把时间域信号通过下式转到频率域：

式中，N1为接收端记录数，N2为发射端记录数.当发射端信号不存在干扰噪声时，Y(ω)与Vs(ω)可以通过(4)、(5)式取得精确值.当接收端存在干扰时(由于发射端测量电压时非常靠近发射源，信号非常强，发射端噪声干扰可以忽略)，利用(4)式求取的频率域信号就会受到干扰，而无法获得十分精确的值.一旦得不到精确的频率域信号，噪声会通过下式进行放大：
当码元宽度为1/6000 s时的8阶m伪随机序列功率谱密度如图 3所示，由图可知在低频段时(伪随机序列带宽范围内)功率谱密度较大，一旦频率超过带宽功率谱会迅速衰减.通过(6)式将反滤波器B(ω)作用于观测信号Y(ω)时，Y(ω)由信号G(ω)和噪音N(ω)组成.当频率在带宽以内时，N(ω)被大大压制，而当频率在带宽边缘或以外时，N(ω)不仅压制不了噪音，噪音反而会被放大.
图 3 当码元宽度为1/6000 s时8阶m码功率谱 Fig. 3 Power spectrum when the code width of 1/6000 s and 8 order
3.2 方法2抗噪原因解析

在方法2中要获取高精度的大地脉冲响应，就要使(B3)式中发射端电流vs(t)与噪声n(t)的互相关函数Rnv(t)要尽可能地小.而m伪随机序列之所以叫伪随机序列，就是因为其性质接近随机信号，而随机信号与其他任何信号的相关性都很差，因此发射电流与干扰信号的互相关函数会很小.另外(B3)式中的有效信号部分包含发射端测量电流的自相关函数Rvv(t)，使有效信号得以增强.以上两点是方法2能够获得高精度大地脉冲响应的原因，实际上这也是m伪随机序列能够抗噪的根本原因之一.

通过上文分析可知，m伪随机序列的随机性是其能够抗噪的重要原因，下面本文将分析什么样的m伪随机序列的随机性才是比较好的.一个序列的随机特性包括以下三个方面的内容(林可祥和汪一飞，1978)：(1)序列中两种不同元素出现的次数大致相等；(2)若把n个同种元素连续性出现叫做一个长度为n的元素游程，则序列中长度为n的元素比长度为n+1的元素游程多一倍；(3)序列具有白噪声的自相关函数.

m伪随机序列的特性为(林可祥和汪一飞，1978)：(1)在每一周期N=2n-1内，0出现2n-1-1次，而1出现2n-1次，0比1少出现一次；(2)m序列中，在每个周期内共有2n-1个元素的游程，其中一半的长度为1，四分之一的长度为2，八分之一的长度为3，即长度为k(1≤k≤n-2)的游程出现的比例为2-k；(3)m序列自相关函数是类似于δ函数的三角波，m序列自相关函数表达式为

m序列自相关函数曲线如图 4所示，当N越大越接近于白噪声.通过对比序列的随机特性的定义和m伪随机序列的特性，m序列的周期越长，即N越大，m序列中两种元素出现的次数越大致相等，并且m序列自相关函数越接近于白噪声.故m序列的阶数越高，其随机性越好，其与干扰信号的相关性也越差，从而能够保证提取更高精度的大地脉冲响应信号.
图 4 m码波形自相关曲线
图中Ap为幅值，Rmm表示自相关函数. Fig. 4 The autocorrelation of pseudo-random sequence
Ap represent amplitude，Rmm represent autocorrelation function.
4 提高m伪随机抗噪方法的数值研究

本节将研究以下两个方面内容：首先对于上节中得出结论，以具体伪随机序列进行检验；随后对确定的伪随机序列，探究提高抗噪的方法.

模型参数仍沿用表 1中的模型参数，但伪随机序列的阶数分别取7、8、9阶.干扰信号也沿用上文中的50 Hz干扰，但强度提升10倍，振幅为10-8，干扰信号表达式为10-8cos(2π×50t).使用不同阶数的编码源电流和50 Hz干扰信号的互相关函数的最大幅值的大小来衡量伪随机编码电流对抗50 Hz规则干扰的能力，幅值越小，抗50 Hz干扰的能力越强.图 5b、5c阶数相同，但由不同本原多项式求得.图 5a为7阶伪随机电流与干扰信号的互相关函数，分析图 5a可知，互相关函数忽高忽低，有很强的震荡性，但对信号的影响取决于它的幅值.因此，为了表征干扰信号对求取结果精度的影响程度，本文选用所有振幅数据中最大值，当m为7阶时振幅最大值为2.448×10-6；对于同阶数的伪随机序列来说，本原多项式可以存在多个，即可以求取多个形式的伪随机序列，图 5b、5c就是当阶数为8阶时，分别利用不同本原多项式求取的伪随机序列电流，再与50 Hz干扰信号求取的互相关函数.通过对比可知，与图 5a一样，两个互相关函数也是存在震荡和周期的特性，但是极大值是不同的，图 5b中的极大值为3.852×10-6，图 5c中的极大值为1.420×10-6，表明即使使用同一阶数的伪随机序列，而利用不同本原多项式确定的序列对干扰信号的抗噪的能力也是不同的；图 5d为当伪随机序列阶数为9时，伪随机序列电流与干扰信号的互相关函数，极大值为1.365×10-6，比7阶、8阶的极大值都要小.通过以上分析可知，随着伪随机序列阶数的增大，m伪随机序列抗干扰能力有增强的趋势，但并不是绝对的；对于同阶的m伪随机序列，利用不同本原多项式求取的序列，对干扰信号的抗干扰能力是不同的.

图 5 不同阶数电流与50 Hz干扰互相关函数
(a) 当阶数为7阶时电流与50 Hz干扰互相关函数； (b) 当阶数为8阶时伪随机电流序列一与50 Hz干扰互相关函数； (c) 当阶数为8阶时伪随机电流序列二与50 Hz干扰互相关函数； (d) 当阶数为9阶时电流与50 Hz干扰互相关函数. Fig. 5 Correlation function of different orders of pseudo-random sequence current and 50 Hz nois
(a) Correlation function of 7 orders； (b) Correlation function of 8 orders； (c) Correlation function of 8 orders； (d) Correlation function of 9 orders.

进一步研究编码电流中编码的循环次数对抗50 Hz周期干扰情况，模拟参数仍沿用表 1中的参数，只改变循环次数.当循环次数变化时，求取的发射电流与干扰信号的互相关函数变化情况如图 6所示，图 6(a—d)是循环次数分别为2、4、8、20次时求取的互相关函数，互相关函数具有一定的震荡性；利用(B4)式求取大地脉冲响应时，利用的互相关数据是有限的，式中的N2是一定的，当N2取500~1000时已经满足求解精度要求，因此只要关心发射电流与接收端数据互相关函数中的一小部分数据的污染情况就可以了.因为中部的值最大，求解精度相对也较高，一般只取这部分的数据，故以后的分析中只关心这部分区域数据受干扰情况即可(数据的中部用虚线标出(图 6))；通过分析图 6可知，循环次数的不同，中间区域的数据受噪声的影响程度不同，这样就可以选择合适的循环次数以达到求取高信噪比大地脉冲响应的目的.

图 6 不同循环次数时电流与50Hz干扰互相关函数
(a) 循环次数为2； (b) 循环次数为4； (c) 循环次数为8； (d) 循环次数为20次. Fig. 6 Correlation function of different cycles of pseudo-random sequence current and 50 Hz noise
(a) 2 cycles； (b) 4 cycles； (c) 8 cycles； (d) 20 cycles.

通过8次循环，伪随机电流与加入50 Hz噪声后伪随机序列激发场的互相关函数如图 7a所示.伪随机电流与干扰信号的互相关函数如图 7b所示，分析此图可知我们关心区域内的互相关函数较小，表明受50 Hz干扰信号的影响较小.通过8次循环得到的大地脉冲响应曲线如图 8所示，图 8a为提取到的大地脉冲响应与解析解对比，图 8b为提取到的大地脉冲响应与解析解的误差率曲线，分析图 8可知，通过8次循环求取的大地脉冲响应曲线误差率控制在6%以内，与解析解吻合度较高，基本不受50 Hz干扰信号的影响.当把循环次数增加一次即增大到9次时，伪随机电流与加入50 Hz噪声后的伪随机序列激发场的互相关函数、伪随机电流与50 Hz干扰信号的互相关函数分别如图 9a、9b所示，通过分析此图可知，在求解区域内电流与干扰信号的互相关函数较强.通过9次循环获取的大地脉冲响应曲线如图 10所示，不难发现：虽然增加了循环次数，但求取的大地脉冲响应与解析解的误差率基本在40%以上，个别点甚至超过120%，其精度远远低于8次循环.在选择循环次数时，并不是越多越好，而是选择合适的循环次数，以达到求取区域内的数据受信号影响较小的目的，从而可获取高精度的数据.

图 7 互相关函数
(a) 伪随机序列电流与伪随机响应互相关函数； (b) 伪随机序列电流与干扰信号互相关函数. Fig. 7 Correlation function
(a) Correlation function of pseudo-random sequence current and pseudo-random sequence response； (b) Correlation function of pseudo-random sequence current and noise.

图 8 提取到的大地脉冲响应与解析解对比图
(a) 通过方法2获得大地脉冲响应； (b) 提取大地脉冲响应误差率. Fig. 8 The comparison of the impulse response extracted from
pseudo-random sequence response with the analytical solutions (a) Results comparison with the analytical solutions； (b) Error rate.

图 9 互相关函数
(a) 伪随机序列电流与伪随机响应互相关函数； (b) 伪随机序列电流与干扰信号互相关函数. Fig. 9 Correlation function
(a) Correlation function of pseudo-random sequence current and pseudo-random sequence response； (b) Correlation function of pseudo-random sequence current and noise.

图 10 计算的大地脉冲响应与解析解对比图
(a) 通过方法2获得大地脉冲响应； (b) 提取大地脉冲响应误差率. Fig. 10 The comparison of the impulse response extracted from
pseudo-random sequence response with the analytical solutions (a) Results comparison with the analytical solutions； (b) Error rate.

通过以上分析可知，伪随机序列激发的场受信号干扰时，伪随机序列的阶数n、利用不同本原多项式求取的伪随机序列、循环次数都会影响大地脉冲响应的求取精度.除随着伪随机阶数的增加，伪随机序列抗干扰能力有明显增加趋势外，其他参数的影响非常复杂.如何选取合适的伪随机序列和循环次数，以达到最佳的求解精度，是接下来要解决的问题.作者通过大量的数值模拟实验，得出了以下规律：

(1)为了表述方便，本文把能够获得最佳大地脉冲响应精度的循环次数，简称为最佳循环次数，当取得最佳循环次数时，提取大地脉冲响应的精度与干扰信号周期的初始相位无关.

(2)同阶伪随机序列有多个本原多项式，即有多个形式的伪随机序列(或多种编码形式)，这些伪随机序列对同一噪声的抗噪能力是不同的，如果某一序列对某一噪声的抗噪能力是最好的，对另一种噪声其抗噪能力不一定是最好的，应该对不同的噪声进行具体分析，但它们的最佳循环次数是相同的.

(3)伪随机序列的周期不能为干扰信号周期的倍数，但干扰信号可以是伪随机序列的倍数.

(4)对于干扰信号频率为f0，长度为N，码元宽度为Δt的伪随机序列，其最佳循环次数n应满足以下条件：

式中，m为正整数或非常接近正整数的值.

为了检验以上结论的正确性，以均匀半空间模型为例进行说明，模型参数仍沿用表 1中的参数，伪随机参数如表 2所示.通过(1)、(2)式计算伪随机序列激发的场，并在场中加入振幅为2×10-7的25 Hz干扰信号，干扰信号表达式为

表 2 均匀半空间模型伪随机模拟参数 Table 2 The pseudo-random simulation parameters for uniform half-space model

伪随机序列激发的场、加入干扰信号后的伪随机序列响应分别如图 11a、11b所示.通过分析以上两图可知，伪随机序列响应最大值约为2×10-8，而噪声的振幅为2×10-7，信噪比约为-10 dB.伪随机电流与伪随机序列响应、25 Hz干扰信号的互相关函数分别如图 12a、12b所示，以上条件都符合上文中总结的规律，可以看到关心区域内的数据受噪声影响较小.通过以上数据获取的大地脉冲响应如图 11c所示，与解析解的误差率如图 11d所示，误差率最大不超过8%，与解析解吻合度较高.

图 11 加入50 Hz干扰后伪随机响应及提取结果对比
(a) 伪随机序列激发的场； (b) 加入干扰信号后伪随机序列激发的场； (c) 通过方法2获得大地脉冲响应； (d) 提取大地脉冲响应误差率. Fig. 11 Simulation of pseudo-random sequence response with 50 Hz noise and results comparison of two methods
(a) Pseudo-random sequence response； (b) Pseudo-random sequence response with 50 Hz noise； (c) Results comparison with the analytical solutions； (D) Error rate.

图 12 互相关函数
(a) 伪随机序列电流与伪随机响应互相关函数； (b) 伪随机序列电流与干扰信号互相关函数. Fig. 12 Correlation function
(a) Correlation function of pseudo-random sequence current and pseudo-random sequence response； (b) Correlation function of pseudo-random sequence current and noise.
5 结论

m伪随机序列的抗噪能力与提取方法息息相关，本文介绍了从m序列响应中提取大地脉冲信号的两种方法.第一种方法将时域信号转到频率域，通过在频率域上的一系列处理取得单位电流电磁响应，再将单位电流电磁响应转换到时间域获得大地脉冲响应.第二种方法直接在时间域进行处理，通过解Wiener-Hopf方程获得.在没有噪声干扰的情况下，两种方法从伪随机响应中获得大地脉冲响应信号的精度都较高.

加入噪声后，第一种方法提取大地脉冲信号的精度要低于第二种方法.造成第一种方法不抗噪的原因是，在时间域转到频率域的过程中，会因噪声的加入造成频率信号不准，而电流归一化的过程中，又因带宽范围以外功率谱密度较低，会把噪声对带宽范围以外频率域信号的影响放大，从而造成信号的畸变.第二种方法能够抗噪的原因是，在互相关的过程中可以利用伪随机序列的随机性，压制干扰信号的影响.

本文还针对m伪随机序列发射源分析了提高抗噪能力的方法，研究表明可通过提高伪随机序列阶数及改变编码方式和循环次数等方式实现对噪音的抑制.对于周期为f0的噪声信号，使循环次数n0、伪随机序列长度N、码元宽度Δt与周期f0的乘积为整数，可以大大降低噪声的影响，至少可以保证在-10 dB噪声条件下工作.

致谢 两位匿名审稿人所提出的修改意见非常细致、中肯，对本文的完善起到了重要作用，在此表示诚挚的感谢！ 附录A 反褶积方法一

采集到得信号需要将发射端仪器响应、接收端仪器响应、大地脉冲响应、噪声干扰等考虑进来，在接收端采集到的信号可表示为

式中，*表示褶积，i(t)为电流变化，是已知的，hs(t)为发射端系统响应，这部分是由发射机内电路、发射端A、B供电极和接地长导线等的等效电路产生的脉冲响应，由于等效电路的电阻抗(包括电阻、电感、电容)是未知的，所以hs(t)是未知的，hr(t)为接收端仪器电路与接收电极MN组成的等效电路产生的脉冲响应，也是未知的，g(t)为未知的大地脉冲响应，n(t)为噪声.记a(t)=i(t)*hs(t)*hr(t)，则(A1)式右端成为a(t)和g(t)的褶积再加上n(t).褶积也称滤波，a(t)称滤波函数.反褶积的目的是从(A1)式提取大地脉冲响应g(t).反褶积也称反滤波，记b(t)为反滤波函数，则
由(A1)和(A2)式的频率表示可知，
(A3)式表示滤波函数的谱A(ω)和反滤波函数的谱B(ω)互为倒数，因此B(ω)很容易由A(ω)求出.知道B(ω)后(A2)式中的反滤波函数b(t)也可由傅里叶变换得到.(A2)式表明为了由观测资料y(t)得到大地脉冲响应，需要知道b(t)或者B(ω)，(A3)式表明B(ω)可由A(ω)的倒数计算，但在a(t)或A(ω)中，hs(t)，Hs(ω)，hr(t)，Hr(ω)是未知的，需要在非常靠近发射端的地方测量电压变化vs(t)(距离单位为cm级.也可在发射端测量电流变化，数据处理方式与发射端测量电压变化并没有本质的区别，在本文中不再赘述.)而得到，发射端测量电压的仪器响应须与接收端测量电压的仪器响应一致，由于非常靠近发射端，所以只有接收端回路产生的脉冲响应.由于非常靠近发射端，电磁信号较强，干扰噪声完全可以忽略：
将(A4)式转换到频率域，并假设发射端仪器响应与接收端仪器响应一致，则有
求上式倒数即可得B(ω).为了避免Vs(ω)为零，并增加稳定性，有
式中，*代表取共轭复数，ε为特别小的常数.考虑伪随机序列带宽，当信号在低频时信噪比较高，在高频时信噪比较差，特设
式中，，Δt为码元宽度.将(A7)式转到时间域可表示为
d(t)可近似为δ(t)：

将(A6)、(A7)式与Y(ω)相乘可得

将其转换到时间域可得
这样就取得了大地脉冲响应信号g(t).附录B 反褶积方法二

方法2主要是借助Wiener-Hopf方程，也是从伪随机序列响应中提取大地脉冲响应的一种常用方法(Ziolkowski，2013).与方法1相同，在接收端采集到的信号可表示为

式中，i(t)为电流变化，g(t)为大地脉冲响应，hs(t)为发射端仪器响应，hr(t)为接收端仪器响应，n(t)为噪声.为了得到发射端电流变化和仪器响应，一般用和接收端系统响应一致的仪器在发射端测量电压变化：
式中，hs(t)为发射端仪器响应.基于Wiener-Hopf方程，将发射端测量的电压与接收端测量到信号相关可得
式中，Ryv(t)、Rvv(t)、Rnv(t)分别为发射端电压与接收端信号的互相关函数、发射端电压的自相关函数、发射端电压与接收端噪声的互相关函数.(B3)式为连续信号公式，为了便于计算需将其离散化.设采样间隔为Δt(采样间隔Δt越小，互相关函数的计算精度越高，故采样间隔应满足相关项的采样定律要求)，采样时间内采样总数为N1，大地脉冲响应点数为N2，则(B3)式离散化后可得
式中，
忽略Rnv项，并将(B4)式转化为矩阵，得
式中，
分别表示观测资料列向量、矩阵、大地脉冲函数列向量，通过解上式即可获得大地脉冲响应.需要说明的是，(B8)式通常是不适定的，需要正则化处理：
通过解代数方程组(B9)式就可以得到G，从而可以得到大地脉冲响应g(m).
参考文献
[1] Di Q Y, Wang M Y, Shi K F, et al. 2002. An applied study on prevention of water bursting disaster in mines with the high resolution V6 system. Chinese J. Geophys. (in Chinese), 45(5): 744-748.
[2] Di Q Y, Wang G J, An Z G, et al. 2006. Geophysical exploration of a long deep tunnel on the west route of South-to-North water diversion project. Chinese J. Geophys.(in Chinese), 49(6): 1836-1842.
[3] Di Q Y, Wang R, et al. 2008. Controlled Source Audio-Frequency Magneto Tellurics (in Chinese). Beijing: Science Press.
[4] Duncan P M, Hwang A, Edwards R N, et al. 1980. The development and applications of a wide band electromagnetic sounding system using a pseudo-noise source. Geophysics, 45(8): 1276-1296.
[5] Fang W Z, Li Y G, Li X. 1993. Theory of Transient Electromagnetic Sounding (in Chinese). Xi′an: North-West Industry University Press.
[6] He J S. 1990. Controlled Source Audio Frequency Magnetotellurics Method (in Chinese). Changsha: Central South University of Technology Press.
[7] Ilyichev P V, Bobrovsky V V. 2015. Application of pseudonoise signals in systems of active geoelectric exploration (Results of mathematical simulation and field experiments). Seismic Instruments, 51(1): 53-64.
[8] Lin K X, Wang Y F. 1978. Principles and Applications of Pseudo-Random (in Chinese). Beijing: People Post Press.
[9] Lu Q H, Sun J Z. 2004. Developing geophysical technology and strengthening its basement-Review on the symposiums on the internationalization of geophysical instruments and the application of geophysical technology in engineering. Progress in Geophysics (in Chinese), 19(4): 726-729.
[10] Luo W B, Li Q C, Tang J T. 2012. Coded source electromagnetic sounding method. Chinese Journal of Geophysics (in Chinese), 55(1): 341-349, doi: 10.6038/j.issn.0001-5733.2012.01.035 .
[11] Luo X Z, Li D W, Peng F P, et al. 2014. Implementation and applications of an coded electrical instrument with anti-interference ability. Progress in Geophysics (in Chinese), 29(2): 944-951, doi: 10.6038/pg20140263.
[12] Shi K F. 1999. The Theory and Applications of Control Source Audio Frequency Magnetotellurics Method (in Chinese). Beijing: Science Press.
[13] Tang J T, He J S. 2005. Controlled Source Electromagnetic Method and Its Applications (in Chinese). Changsha: Central South University Press.
[14] Tang J T, Luo W B. 2008. Pseudo-random electromagnetic exploration based on Invert-Repeated m-Sequence correlation identification. Chinese Journal of Geophysics (in Chinese), 51(4): 1226-1233.
[15] Velikin A B. 2009. Patent No. 2354999.
[16] Wright D, Ziolkowski A, Hobbs B. 2002. Hydrocarbon detection and monitoring with a multicomponent transient electromagnetic (MTEM) survey. The Leading Edge, 21(9): 852-864.
[17] Wright D, Ziolkowski A, Hall G. 2006. Improving signal-to-noise ratio using pseudo random binary sequences in multi transient electromagnetic (MTEM) data.//68th EAGE Conference & Exhibition. SPE, EAGE.
[18] Wright D A, Ziolkowski A, Hobbs B A. 2001. Hydrocarbon detection with a multi-channel transient electromagnetic survey.//71st Annual International Meeting, SEG, Expanded Abstracts, 1435-1438.
[19] Wright D A. 2003. Detection of hydrocarbons and their movement in a reservoir using time-lapse multi-transient electromagnetic MTEM data[Ph. D. thesis]. United Kingdom: University of Edinburgh.
[20] Wu X, Xue G Q, Di Q Y, et al. 2015. Accurate identification for the electromagnetic impulse response of the earth with pseudo random coded waveforms. Chinese Journal of Geophysics (in Chinese), 58(8): 2792-2802, doi: 10.6038/cjg20150814.
[21] Yin C C, Huang W, Ben F. 2013. The full-time electromagnetic modeling for time-domain airborne electromagnetic systems. Chinese Journal of Geophysics (in Chinese), 56(9): 3153-3162, doi: 10.6038/cjg20130928.
[22] Zhao B R, Zhao J, Zhang H K, et al. 2006. The PS100 high precision earth-eletrictity instrument system (IP to IP) with controllable source—application of CDMA technology to the measurement of earth-resistivity for the first time. Progress in Geophysics (in Chinese), 21(2): 675-682.
[23] Ziolkowski A, Wright D. 2007. Removal of the airwave in shallow-marine transient EM data.//77th Annual International Meeting, SEG Expanded Abstracts, 26: 534-538.
[24] Ziolkowski A, Hobbs B A, Wright D. 2007. Multitransient electromagnetic demonstration survey in France. Geophysics, 72(4): F197-F209.
[25] Ziolkowski A, Parr R, Wright D, et al. 2010. Multi-transient electromagnetic repeatability experiment over the North Sea Harding field. Geophysical Prospecting, 58(6): 1159-1176.
[26] Ziolkowski A. 2013. Wiener estimation of the Green′s function. Geophysics, 78(5): W31-W44.
[27] 底青云, 王妙月, 石昆法等. 2002. 高分辨V6系统在矿山顶板涌水隐患中的应用研究. 地球物理学报, 45(5): 744-748.
[28] 底青云, 王光杰, 安志国等. 2006. 南水北调西线千米深长隧洞围岩构造地球物理勘探. 地球物理学报, 49(6): 1836-1842.
[29] 底青云, 王若等. 2008. 可控源音频大地电磁数据正反演及方法应用. 北京: 科学出版社.
[30] 方文藻, 李予国, 李貅. 1993. 瞬变电磁测深法原理. 西安: 西北工业大学出版社.
[31] 何继善. 1990. 可控源音频大地电磁法. 长沙: 中南工业大学出版社.
[32] 林可祥, 汪一飞. 1978. 伪随机码的原理与应用. 北京: 人民邮电出版社.
[33] 陆其鹄, 孙进忠. 2004. 发展技术, 增强基础—地球物理仪器国际化与地球物理技术在工程上的应用研讨会综述. 地球物理学进展, 19(4): 726-729.
[34] 罗维斌, 李庆春, 汤井田. 2012. 编码电磁测深. 地球物理学报, 55(1): 341-349, doi: 10.6038/j.issn.0001-5733.2012.01.035.
[35] 罗先中, 李达为, 彭芳苹等. 2014. 抗干扰编码电法仪的实现及应用. 地球物理学进展, 29(2): 944-951, doi: 10.6038/pg20140263.
[36] 石昆法. 1999. 可控源音频大地电磁法理论与应用. 北京: 科学出版社.
[37] 汤井田, 何继善. 2005. 可控源音频大地电磁法及其应用. 长沙: 中南大学出版社.
[38] 汤井田, 罗维斌. 2008. 基于相关辨识的逆重复m序列伪随机电磁法. 地球物理学报, 51(4): 1226-1233.
[39] 武欣, 薛国强, 底青云等. 2015. 伪随机编码源电磁响应的精细辨识. 地球物理学报, 58(8): 2792-2802, doi: 10.6038/cjg20150814.
[40] 殷长春, 黄威, 贲放. 2013. 时间域航空电磁系统瞬变全时响应正演模拟. 地球物理学报, 56(9): 3153-3162, doi: 10.6038/cjg20130928.
[41] 赵璧如, 赵健, 张洪魁等. 2006. PS100型IP到端可控源高精度大地电测仪系统——CDMA技术首次在地电阻率测量中的应用. 地球物理学进展, 21(2): 675-682.