（通信与信息系统专业论文）瞬变电磁法矩形线框一维正演及电磁信号去噪的研究.pdf

摘要 瞬变电磁法是探测地下浅层结构的一种有效方法，其基本理论已经比较成熟，特 别是在均匀水平层状介质的情况下，对于偶极子方式和环形发射框的中心方式问题， 人们已经很好的利用了数学方法加以解决，得到了解析的表达式。所以从上个世纪5 0 年代开始，瞬变电磁测深的方法得到了广泛的应用，与之相关的文章也不断涌现。 然而，这些文章大都是以偶极子方式或环形发射框的中心方式问题作为论述的重 点，而有关以固定的矩形线框作为发射源的论述却很少。但解决这样的问题却十分必 要。其原因有二：( 1 ) 由于发射机功率的限制，瞬变电磁法中的接收线圈要放置在靠 近而不是远离发射线圈的位置，这时发射源就不可以再简单地被看作是偶极子了；( 2 ) 对于大尺寸发射线圈的野外铺设来说，矩形要比圆形容易得多，所以野外测量大多数 情况下是以固定的矩形线圈作为发射源的。如果在这两种情况下仍然用传统的偶极子 方式或环形发射框中心方式的方法来分析问题，其结果必然会有所偏差。所以，有必 要对瞬变电磁法中的矩形线圈方式进行研究。 本文系统地讨论了时间域瞬变电磁法( t i m e d o m a i ne l e c t r 。m a g n e t i cm e t h o d ， 简称t e m ) 一维正演的基本理论，着重研究了层状地层上以固定矩形线圈作为发射源 的瞬变响应问题，取得了良好的结果。文中还探时了利用小波多分辨率分析的方法对 电磁信号去噪的问题。 关键词：瞬变电磁场；时间域响应；一维正演；汉克变换；g u p t a s a r m a 算法；小波分 析 a b s t r a c t t h ei d e ao fu s i n gt r a n s i e n te l e c t r o m a g n e t i cf i e l dt od e l i n e a t el a y e r e d e a r t hs t r u c t u r ei s n o tn e w i nf a c t ，m e t h e o r yo ft e m ( t r a n s i e n te l e c t r o m a g n e t i cm e t h o d ) i sv e r yw e l l d e v e l o p e dm a i n l yb e c a u s e a l l e l e c t r o m a g n e t i cp r o b l e m si n s t r a t i f i e dm e d i ai nw h i c h s t r a t i f i c a t i o ni sh o r i z o n t a la n d h o m o g e n e o u s a l ea m e n a b l et o e l e g a n t m a t h e m a t i c a l s o l u t i o n s f r o m1 9 5 0 s ， m a n yp a p e r so nt e m h a v e b e e n p u b l i s h e dc o n t i n u o u s l y i nf e wo ft h ep u b l i c a t i o n s h o w e v e r ，a r et h es o l u t i o n sf o re l e c t r o m a g n e t i cf i e l dc l o s e t oar e c t a n g l el 0 0 0s o u r c eo fc u r r e n to ff i n i t es i z ep r e s e n t e d t h ed e v e l o p m e n to fs u c h s o l u t i o ni s n e c e s s a r yb e c a u s e ( 1 ) t h et r a n s m i t t e rp o w e rr e q u i r e m e n t i sm u c hl e s sf o r s o u n d i n g n e a ral o o pt h a na tad i s t a n c e ( w h e nt h el o o pc a l lb er e g a r d e da sad i p o l e ) a n d ( 2 ) i ti se a s i e rt ol a ya l o o po fl a r g es i z ei nas q u a r e o r r e c t a n g l et h a n ac i r c l e i nt h i st h e s i st h ee s s e n t i a lt h e o r yo f1 一df o r w a r df o rt r a n s i e n te l e c t r o m a g n e t i cf i e l di s s y s t e m a t i c a l l ye x p l a i n e d t h ep r o b l e m o fa r e c t a n g u l a rl o o p s o u r c eo fc u r r e n to n m u l t i l a y e r e de a r t h h a sb e e ns t u d i e da n ds o l v e d g o o dr e s u l t sh a v eb e e no b t a i n e da n d d e m o n s t r a t e d m o r e o v e r ，t h em e t h o dt or e d u c et h en o i s ei ne l e c t r o m a g n e t i cs i g n a lb y w a v e l e tm u l t i s c a l ed e c o m p o s i t i o ni sa l s oe x p l o r e d k e yw o r d s ：t r a n s i e n te l e c t r o m a g n e t i cf i e l d ；t i m e d o m a i nr e s p o n s e ；1 - df o r w a r d ； h a n k e lt r a n s f o r m ；g u p t a s a r m aa l g o r i t h m ；w a v e l e ta n a l y s i s 瞬变电磁法中矩形线框一维正演算法的研究 0 前言 l 、浅部电磁方法的比较 电磁方法是浅部环境与工程探测的重要手段，就目前来说，所有浅部探测方法均 是由相应的深部探测方法发展而来的。电磁方法分两大类：频率域电磁方法( 如m t ) 和时间域电磁方法( 如t e m ) 。 1 ) 频率域电磁方法的特点：对于频率域来说，地下介质一定时，探测深度是由 频率决定的，要达到探测浅层目标体的目的，就必须提高观测信号的频率。例如目前 的新型仪器可以记录l o k h z 的信号，若地下介质的电阻率为5 0 n m ，根据概念性的 穿透深度 8f p 盯f 斗 “2 ( 0 1 ) 估算，其最小探测深度约为3 5 米。实际上真正的有效探测范围在一百多米以下，因 为频率在i k h z 以上的信号需要人工源发射，对于强噪音地区，该频段的效果不甚理 想。我们知道，电磁波的传播特性与频率有关，当频率非常高时，其主要表现为波的 直线传播，这就是探地雷达( g p r ) 的设计原理。探地雷达的频率范围一般为几十千 赫兹到数兆赫兹，由上式估算其探测范围可从地下几厘米到几十米，而且分辨率非常 高，超浅层分辨率可达厘米的量级。实际上探地雷达的最佳探测范围只有不到十米， 特别在具有低阻盖层的地质条件下探测深度更小( 一般为2 5 米) 。原因有二：首先， 低阻盖层的屏蔽作用，高频成分迅速衰减，只有低频成分才有可能到达下层；其次， 穿过低阻层到达下层的低频成分波的绕射和感应特性增强，使探测能力降低。由此可 见，地下2 10 0 米左右的深度范围是频率域方法的一个弱视区( 或半盲区) 。 2 ) 时间域电磁方法的特点：对于时间域来说，地下介质一定时，探测深度是由 观测时间的早晚决定的，要探测浅层目标体，唯一的办法就是缩短延迟的时间( 即电 流关断后开始记录的时间) 。目前国际上最好的浅层t e m 仪器e m 一4 7 的最早记录时间 为6 8 肛s ，对一个5 0 n m 的半空间模型，根据概念扩散深度公式 8 。= 2pt 肛 “2( 0 2 ) 估算，最小探测深度约1 8 米。所以，对于目前市场上的时间域t e m 仪器，从地表到 地下20 米左右是一个半盲区。通过多年的实践，人们已经认识到时间域比频率域具 有多方面的优势。特别是1 千米以上范围是时间域探测的强项。所以在浅部探测中人 们主要研究并采用的是时间域的t e m 方法。因此，解决从地表到20 多米( 实际中要 考虑到5 0 米) 范围的探测理论和技术就成了浅层电磁方法的关键问题。 一 壁窭皇壁鳖堑垄垡焦= 堡至塑苎鲨竺堑壅 2 、时间域瞬变电磁法( t e m ) 的基本原理 时间域瞬变电磁法( 简称t e m 法) 是利用不接地回线( 或接地线源) 向地下发送一次 脉冲磁场( 或电场) ，在一次脉冲磁场( 或电场) 的间歇期间，利用线圈( 或接地电极) 观测 地层响应的二次涡流场的方法。 t x ( r x ) lh 2 忍船 i ( 发射电流) 图0 i 瞬变电磁法原理示意图 f i g 0 1g r a p h i ci l l u s t r a t i o no f t e m 1 0 瞬变场的时间范围 ， 0 1 ( 图0 2 瞬变电磁法中发射线圈j 目 的阶跃电流 f i g 0 2 s t e p c u r r t m tc a n i e d b y t r a n s l i f i t t e r l o o p i n t e m t 时j 可) 如图o 1 所示，在t e m 法中，首先向放置在层状地层表面的发射线圈中加载一个 如图0 2 所示的在t 。时刻截止的阶跃电流。当一个稳定的电流i 。在t ，时刻突然断开 时，由于电磁感应，将会在空间产生一个瞬变的电磁场，这个瞬变的一次场将在地下 感应出一个瞬变的二次场并不断向下扩散，当电磁波遇到不同介质的表面时就会有回 波反射回来，因此瞬变的二次场携带了地下介质结构的信息。如果记录到t - 时刻以后 的二次场，就可得到地下不同深度介质的信息。因为发射电流位于地面，二次场首先 从地面开始向下扩散，理论上讲，随着记录时间的不同，二次场在地下扩散的深度也 瞬变电磁法中矩形线框一维正演算法的研究 不同。也就是说，不同时间记录的电磁信号含有不同深度的地下结构的信息，即早期 信号代表了浅部结构的信息，晚期信号代表了深部结构的信息。所以，瞬变电磁法也 被称作时间域( t i m e d o m a i n ) 的测深方法。由此可以看出如果想得到更高分辨率的 浅层地质结构的信息，就需要早期信号记录得更加准确。 竖奎皇壁堡堕垄垡焦= 丝至遗苎堡塑堕 1 问题的提出 在瞬变电磁法中，我们为什么要研究正演算法? 正演算法在瞬变电磁测深中具有 什么样的作用和意义呢? 这是我们在做进一步研究之前，首先要回答的问题。本章将 从浅层探测分辨率的提高以及对记录到的电磁信号进行反演，从丽得到地质分析结果 两方面，对这一问题进行分析和阐述。 1 1 浅层高分辨率测深中的正演问题 t e m 的原理如图0 2 所示，当一个稳定电流i 。在t 。时刻突然断开后，将在地下产 生一个瞬变的二次场并不断向下扩散，如果记录到t ，时刻以后的二次场，就可得到地 下一定深度范围内介质的信息。因为发射电流位于地面，二次场首先从地面开始向下 扩散，理论上讲，我们可以得到从地面到信号到达的任意深度。然而，实际情况并非 如此。实际上发射电流断开后并不立即下降为o ，丽是要经过一定的时间吼后才下降 为0 。如图1 1 所示。 i ( 发射电流) | | 此段记录受到关 t 0 号。断电流的影口商 i o 一常规仪器的观测范围 0t lt 2t r i 图1 i 实际测量中早期信号受到关断电流的影响 f i g 1 1e a r l ys i g n a li sa f f e c t e db y t u r n o f f c u r r e n t 图1 1 中t 。一t 。这段时间间隔叫关断时间( t u r n o f f t i m e ) ，该时间段内的电流 叫关断电流( t u r n o f fc u r r e n t ) 。实际观测表明，从t ：开始的相当一段时间内( 如 到时刻t ) ，瞬变场严重受到关断电流的一次场的影响而表现出非常复杂的性质，如果 用于解释将产生错误的结果。因此，在实际记录中人们采取了一种忍痛割爱的措施， 即舍去t 。t 这段时间的早期信号，只记录t 时刻以后的晚期信号。这就意味着我们 舍去了浅层信息。 图1 2 ( b ) 是记录的一条试验曲线，曲线所对应的发射电流示于图1 2 ( a o 4 瞬变电磁法中矩形线框维正演算法的研究 兮 e 、- 一 3 2 0 j r j 7 b ： i k 。 一2 0 0口2 0 04 0 0即0 8 0 01 0 0 0 t ( u s ) 图1 2 全程总瞬变场的变化特征 ( a ) 发射电流( 一次场) 的变化 ( b ) 总瞬变场曲线，该曲线的早期信号既包含二次场也包含一次场 f 追1 2c h a n g e s i na l l t i m eg e n e r a lt r a n s i e n te l e c t r o m a g n e t i cf i e l d 图1 2 ( b ) 的曲线具有重要意义，它详实地揭示了瞬变过程与一次场之间的关系， 为我们认识瞬变的内在机理以及采取相应对策提供了第一手依据。它启示我们，如果 能够从总场中剔除关断电流所产生的一次场，那么就可以得到从电流关断时刻以后整 个过程的二次瞬变场，把t e m 的记录和解释扩展到全时域，如果正确的记录到早期信 号，就可以提高浅层测深的分辨率。 那么如何得到关断电流所产生的一次场呢? 实际上，关断电流可以被分成若干小 电流元，每一个小电流元可以近似看作阶跃函数( 图1 3 o 如果我们能够求得t 时刻 阶跃电流源( 阶跃电流为i ) 在层状介质表面所产生的响应b z ( 。) 口) = f ( 1 ，r ) ，则每一 瞬变电磁法中矩形线框一维正演算法的研究 个小电流源在t 时刻所产生的响应为舾：( 。) ( 丁) = ，( ( 由，6 f ) ，r ) ，通过积分运算 曰。( 玎) 汀) = l 。f ( ( 6 1 ，鼙) ，r ) 扫，我们就可以得到关断电流在t 时刻所产生的一次场 b 。( o f f ) ( t ) 。为了得到去除关断电流影响的准确的二次场早期信号，我们可以从记录到 的总场中减掉由关断电流所产生的场，即如陬) ( r ) = b z ( 一j ( 丁) 一易( 甜】( r ) 。由此可 见去除关断电流对早期信号影响的关键在于得到阶跃电流源在层状介质表面所产生 的电磁场的响应，而这正是瞬变电磁法正演算法所要解决的问题。所以提高浅层探测 能力的关键在于正演问题的解决。 图1 3 关断电流的分割( 观测时间t 位于t o 以后) f i g 1 3d i v i s i o n so f t u r n o f f c u r r e n t 1 2 反演算法中的正演问题 1 反演的目的和任务 宇宙万物都在矛盾中存在着、运动着。任何事物都是一分为二的，对立统一的， 既有正面，也有它的反面。在初等数学中有指数和对数，三角函数和反三角函数，双 曲函数和反双曲函数等；在高等数学中有傅氏变换和傅氏反变换，拉普拉斯变换和拉 普拉斯反变换，褶积和反褶积等；在物理学中有正电子和负电子，质子和反质子，物 质和反物质等；同样，在地球物理学中有正( 演) 问题和反( 演) 问题。 如果我们把用已知自变量x 去求函数y = f ( x ) 的值，已知一个函数( x ) 去求其变 换f - - l ( f ( x ) ) ，已知一个地球物理模型去求其响应函数这类问题叫正演问题的话，那 么由f ( x ) 求x ，由l ( f ( x ) ) 求f ( x ) ，由地球物理观测资料去反推地球物理模型，这类 问题就统称为反问题，或者求反演问题。 6 瞬变电磁法中矩形线框一维正演算法的研究 由此可见，地球物理学中反演的目的是根据观测数据求取相应的地球物理模型。 因此，首先必须确定观测数据和地球模型参数之间的函数关系，使地球物理工作者既 可根据给定的模型参数计算相应的观测数据( 即实现正演计算) ，也可以根据观测数 据求取地球物理模型的参数，实现反演映射。显然，正演是反演的前提和条件，只有 解决了正演计算，不管是靠勰析的方法还是数值的方法，才有可能实现反演映射。 2 反演的计算机实现 在上个世纪5 0 年代以前，由于计算工具的限制，用于解释推断反演结果的方法主 要是将测量数据与事先计算好的大量正演结果进行比对的选择法和量板法。1 9 4 7 年第 一台电子数字计算机的诞生标志着科学技术的一场革命，对地球物理数据分析产生了 巨大的推动。5 0 年代地球物理中的正演问题已广泛的使用了电子计算机。到了6 0 年 代，人们开始利用计算机高速运算的特性，试验用计算机代替解释者自动修改地球模 型的参数，实现所谓的计算机自动解释。 反演之前首先要解决的是与之相关的正演问题，即找到一种函数关系，这个函数 关系表示了，在给定一组地球模型参数和试验方案的条件下，应得到的精确测量数据 值。在瞬变电磁测深法中，就是在给定地下各层的厚度和电阻率的条件下，求通有阶 跃电流的线圈在地面所产生的响应。下面我们用简化的形式来说明如何用计算机来实 现自动反演。 在给定模型参数之后，我们可以把精确数据d ，和模型参数m ，的关系表示成以下关 系式( 正演) ： d 。= a i ( m l ，m ) 其中i = l ，n ；j ：1 ，j ，用向量形式记为 埘7 = ( m 1 ，m j ) d 7 = ( 盔，d ) 我们可以根据先验知识给出模型参数的初始猜测 m o t = ( m ? ，m ? ) 然后按照下面的计算流程图( 图1 4 ) 来计算反演问题。 7 瞬变电磁法中矩形线框一维正演算法的研究 图1 ，4 反演的一种迭代流程 f i g 1 4 a f l o wc h a n o f i n v 耐o n 由流程图可以看出，反演实际上就是不断对模型参数修正，进行正演迭代的过程。图 1 5 给出了一个5 层模型进行反演迭代的结果。 由图1 5 可以看出经过不断的修改地球模型( 各层厚度h ，各层的电阻率p ) 正 演的结果( 嗵线) 逐渐逼近实际的测量数据( 圆点) 。当由模型所得到的正演结果与 测量数据之间达到一种很好的匹配时，我们就可以把此时的模型作为地下地质结构的 一种可能分布，从而实现了反演的自动化。 由以上两节的讨论可以看出，无论是想提高瞬变电磁法浅层测深的分辨率，还是 想实现反演的自动化，都是以正演算法作为基础和前提的。在接下来的几章中我们将 会对如何实现正演作专门的讨论。 瞬变电磁法中矩形线框一维正演算法的研究 图15 一个5 层模型的反演迭代结果 圆点表示实际测量的数据， 折线表示修正的地球模型， 曲线表示正演的结果。 f i g 1 5 r e s u l t o f i n v e r s i o n o na f i v e - l a y e r e m h 9 ll#萨柚。l 瞬变电磁法中矩形线框一维正演算法的研究 2 频率域矩形线框的一维正演算法 在对通有阶跃电流的矩形线框在水平均匀层状介质表面上的响应( 正演) 进行讨 论之前，我们先从一种简单的情况人手，研究一下通有简谐电流( e “) 的矩形线框 在水平均匀层状介质表面上的响应。这样做的原因有二，一是由简谐源产生的电磁场 方程可以化为简单的复矢量形式，有利于进一步的数学推导；二是由高等数学中的傅 立叶变换可知，个符合狄氏条件的函数可以表示为由无数简谐函数的线性叠加，所 以，我们可以容易的由简谐源的响应得到阶跃源的响应。这个过程实质上就是一个傅 氏反变换的过程。所以，人们有时也会把简谐源的响应称作频率域的响应，而把阶跃 源的响应称作时问域的响应。 在对通电矩形线框的频率域响应求解的过程中，我们的思路是由简单到复杂，首 先导出通有简谐电流的磁偶极子在水平层状介质表面上的电磁场响应的表达式。然 后，根据电磁场的互惠原理，把磁偶极子的响应的表达式推广到矩形线框响应的表达 式。但是，由于表达式中存在着振荡的b e s s e l 函数项，所以给数值计算带来了麻烦， 为了解决这个问题，我们引入了h a n k e l 积分( 一种计算含b e s s e l 函数积分的快速算 法) ，从而实现了对矩形线框响应计算的计算机化。本章的最后，我们给出了计算结 果，通过和前人所作工作的比较，可以证明本文算法的正确性和准确性。 2 1水平层状介质表面上垂直磁偶极子电磁场表达式 本节将导出水平层状介质表面上电磁场各个分量的表达式。因为一次电场只有切 向分量e ? ，而该分量又和不同电导率地层的分界面不相交，所以不产生电荷积累。 根据这和原因，感应电流都位于水平面内，并根据电磁场的柱对称性，所以，仅有三 个不等于零的分量： e = o ，匕，o ) ，h = 【芹，0 ，h z ) ( 2 1 ) 假定通有简谐电流l e ，磁偶极矩为m 的磁偶极源位于地球表面( 图2 1 ) 柱坐标系 的原点处。根据电磁场的理论，我们用磁矢量位a + 定义电场和磁场为： e = f 甜必v a 和h = 七a + v v 。a 其中，m 为简谐电流角频率，斗为磁导率，k 为波数。 由此可继续推出电磁场个分量的表达式： 驴e ：- 0 ，岛一f 掣警 和 耻篆，驴眦+ 等一- o 1 0 ( 2 3 ) 瞬变电磁法中矩形线框一维正演算法的研究 z 图2 1 层状介质表面上垂直轴磁偶极源模型 f i g 2 1m o d e l o fav e r t i c a lm a g n e t i cd i p o l ea b o v e t h es u r f a c eo fas 仃a f i f i 。d “e d i 8 r 这种矢量位在每一层中都应满足h e l m h o t z 方程式，并用i 表示层位，有 v 2 a ( + 七j a 州= 0 ( 2 4 ) 式中，a ( j 为第i 层的矢量位，由电磁测深中的似稳场条件可得简化的波数 ? ：i t 7 ，芦街，o 。为第i 层的电导率。 由于任一层中的电场分量只有水平分量，并考虑到e = i w v a ，所以我们仅用。 的垂直分量求解；也就是说，利用 a + = ( 0 , 0 ，) 囊熬雾繁錾霖装篇箬翟嬲搿麓墨蓉糯嚣嬲崔 雩数要雪蝥：。耍坐 n 标x e 中霍蒿! 焉雾昙萼兰罘裟曩；齐翁翕毒姜箍要姜 (1 = n h 2 ，lin he 】= m x e 2 】) ，在每一个分界回上，甄1 i j 以雅哥田出了里 位a ：表示的边界条件： 瞬变电磁法中矩形线框一维正演算法的研究 a ；”= a ；什” 娑：半，艨。地 ( 2 6 ) 式中乜为从地面到第i 层和第i + l 层间分界面的深度。这里，我们首先假设该层状介 质为三层，那么，由电磁场的波动方程我们可以解出个层中电磁场的矢量位表示： 如= 鬈r ( p 叫纠+ d 0 一以m r ) 咖，z 0 a 】z = 4 粕0 ( c 1 e - m * z + d l e m t z ) ，。( m r ) 出拜，o z l ( 2 7 ) 如= 篆f ( q e - r z + d 2 氓( m 胁， z 轨 a ；z = 4 m 刀j 0 “c 3 e - n o z ，。( 删咖，z 魄 式中，m = 占可，m 。= 巧蕊，= 0 孺。k ，k 。，b 分别为第一层、第 二层、第三层的波数。应当指出，因为z 在任一单层中都保持有界，所以a 乞和a 乞的 表达式中都包含有正负指数；只有当z 无界时，才能应用指数条件。 由( 2 6 ) 式给定的边界条件，就能得到一个代数方程组，以确定未知系数d o ， c ：，n 和c ，。利用汉克尔积分( 含b e s s e l 函数的积分) 的正交性，得： 1 + d o = c i + d t 一? 托+ 州d o = 一m i c l + m i d l c l e 一”1 + d t e l h = c 2 已一恍 + d 2 e 啦 一，n lc l e 一“ + ，2 ld l 已“l = 一，n 2 c 2 p 一4 2 喝+ ， 2 d 2 p “1 c 2 p 一啦k + d 2 e 啦 2 = c 3 e 一码k m 2 c 2 已一啦 2 + m 2 d 2 e 帆k = m 3 c 3 e 一均k 取以上方程式组中最后两个方程之比： c z e 一“k + d 2 9 ”一m 2 c 2 e “”一d 2 e “ 2m 3 或 l + ( | d 2 c z ) e 2 ”一i 1 1 2 卜( d 2 c 2 ) f 2 “”m 3 1 2 ( 2 8 ) 幽堕鳖主星垄丝焦二竺垩堡茎整! 燮 由于 所以有 d 2 c 2 = e 1 ( 0 2 1 c 2 ) 1 + ( d ，c ，) e 2 愧屯1 + p 2 啦 2 + 如7 岛) l 。-_ 1 一( j d 2 c 2 沁。吒1 一p 2 j ，1 2 2 + h i ( 岛7 c 2 一c 。血m 2 h + ：l n ( 4 慨) - 叫 于是 m 2 h 2 + 1 n ( d 2 c 2 ) = 一c o t h 。1 ( m 2 i m 3 )( 2 9 ) 这里c 。t h ( x ) = e p x ，+ 一p e 。- x = 1 l + 二e 二- 五2 x 是双曲余切函数。 利用方程式组( 2 8 ) 在z = h ，界面匕的边界条件，有： 或 c 劬h + ：l m 刚g ， = 蚴c 。m m ：啊+ 扣峨 由( 2 9 ) 式移项变换得： 1 n ( d 2 c 2 ) = 一m 2 h 2 - c o t h 一渤2 m ；) 把陔式代入( 2 1 0 ) ，有： c o t n h + ：1 m 刚g ) = ( m l m 2 ) c o t l l b 2 h i 所2 h 2 一c o t h 一1 ( 川2 m 3 ) j 或 c 础m 。h 1 + 1 2 ( d r c j ) = 一( m l m 2 ) c o t h m 2 h 2 + c o t h 一1 ( m 2 m 3 ) j 其中，h j ：h 1 ，h ：= h ：一h 。，代表了各层的厚度。 最后，由( 2 8 ) 式中的头两个方程式，得： 1 + d o c l + d tm 1 一d oc l d lm l 或 ( 2 1 0 ) ( 2 1 1 ) 吣丽 尹一， d d+ 一一 竺坝 一 _ l 撑抟 c c 舰m 一一 h m p p d d + 一 h m m p一已 c c 鳖竺皇堡鳖! 堑堑堡堡二堡至塑苎鳖盟堕塞 等署c 。出巴- 删q ) 由( 2 1 1 ) 式求出c o m 丢l n ( d j q ) i 代人( 2 1 2 ) 得： 嚣一卜，一一曙c o t h ( m 2 h z + c o t h - t 。m ，2 ， 于是 器；署尺冲- + = 式中 尺i ：c o t h m l h l + c o t h 一1 ( 优i r n 2 ) c o t h h ，h ，+ c o t h - i f m ，m ，) i 这样，将( 2 1 3 ) 式代回( 2 7 ) 式，则得在由三层水平介质组成的地层的情况下 地面上( z = 0 ) 的矢量位表达式为： 疋：兰r 卫坐塑由 2 万如，”+ m ，r ： 当介质为两层( 譬如，m 2 = m 3 h ：一* ) 时，我们可导出此时： 群= c o t hl h l + c o t h 一1 ( m l m 2 ) j( 2 1 4 ) 当介质为一层时( 譬如，m l = m 2 = m 3 ，h ，_ + ) 时，由于c o t h ( ) ：1 ，c o t h1 ( 1 ) ：，我 们也可导出此时： r i2 i( 2 1 5 ) 对于三层模型，我们可以化为递推的形式， 研一 m 。h t + 训c 马m 2 弛k m ) l j 式中 尺；( 七2 ，k 3 ，m ) = c o t h l n 2 日2 + c o t h 一1 ( m 2 i m 3 ) j( 2 1 7 ) 式( 2 1 7 ) 右边的c o t h k ：日：+ c o t h 一1 ( m 2 m 3 ) 】相当于由参数k ：，k ，和h 2 表征的二层介 质的函数r ，如式( 2 1 4 ) 所示。 根据数学归纳法，n 层介质( 如图2 2 ) 的表达式可写为 r ：= 1 1 4 、，、， 地 2 2 l ( 瞬变电磁法中矩形线框一维正演算法的研究 r ：一l = c o t h m h 日+ c o t h 一1 k m hi m 。) 】 尺? ：c o m b 。h 。+ c o t h 一1 im f m 。) ，聪j | ( 2 1 8 ) r j ：c o t h i n l h l + c o t h 一1 ( m l m 2 ) r il l ：靠，i t ，。l7 妒l ”l _ 一r ? ： 只，h 。 l 螅l “ l 图2 2n 层介质模型 f i g 2 2m o d e l o fa n l a y e re a r t h 将由( 2 1 8 ) 式递推得到的研代回到( 2 7 ) 式中，我们就可以得到在n 层水平层状 介质表面上( z = 0 ) 的矢量位表达式为： 铲篆f 糟 。1 9 ) 由( 2 3 ) 式和( 2 1 9 ) 式所得结果可知，磁偶极源电流和水平层状介质中感应电流 一 一 一 一一一 一 一一 一一一一 一 十，+ 一一一一 瞬变电磁法中矩形线框一维正演算法的研究 在地球表面上严生的总电磁场为： 和罐篆f 器 日，= 篓r 暑石+ ：，砰州m r ) 咖 ( 。z 。) 日z = 篆f m + 东磊_ 以m ，) 咖 至此，我们求得了垂直磁偶极子在水平层状介质表面的响应，在下一节中我们将 会把这一结果推广到矩形线框在水平层状介质表面的响应。本节所导出的磁偶极子响 应公式是参考了文献 i o i 的结果，与之稍有不同的是关于月? 的定义。关于电磁场的理 论和更为详尽的推导可以参考文献i 0 1 。 2 2 水平层状介质表面上矩形线框响应的表达式 本节我们将在上一节所作讨论的基础上，继续对通有简谐电流的矩形线圈在均匀 水平层状介质表面上的响应作进一步的推导。 在推导之前，让我们首先了解一下电磁场中一种有趣的互易( r e c i p r o c i t y ) 现 象。如图2 3 所示， c 图2 3 电磁场中有趣的互易现象 f i g 2 - 3t h ep h e n o m e n o no f r e c i p r o c i t y c 。，cz 是两个线圈，如果我们在c 。中通一简谐电流i e ，由于电磁感应，此时， c ：将感应出电压v ；如果我们反过来做相同的实验，在c ：中通入同样的简谐电流l e ， 那么这时c 。中也将感应出同样的电压v 。这就是电磁场中有趣的互易( r e c i p r o c i t y ) 现象。 下面我们给出对这一互易( r e c i p r o c i t y ) 现象的数学证明， 设线圈c 。，c 2 的面积分别是s ：和s 。，c ，和c 。之间的互感足够小，以至将任一个线 1 6 瞬变电磁法中矩形线框一维正演算法的研究 剧放八另外一个线圈所产生的电磁场中时，对电磁场的场形没有任何影响。 首先在c ，中通人电流i e “，则在c 。处产生的磁场为： h 。= 石2 m 2 f & r ( 1 一妇) 其中，心= s 。i e “为线圈c ：的磁偶极矩，r 为两线圈之间的距离，k 为波数。此时， c 。中感应出的电压为： k = 一等= 一嘎o t 点= 一警= 嘤掣2 - 训 一器( 即f s ：以1 - i k r ) 同理，如c ，中通电流l e “时，则c 。产生的电压为： k 一警( s 2 ，f 即e i k r ( 1 一蚴 因此，v = k 所以，在c ：中通以电流i e “，在c ，处产生电压v ，反过来在c ，中通以同样的电 流在c 。处产生的电压也等于v 。 接下来，我们将运用这一有趣的互易来推导矩形线框在水平层状介质表面处的 响应。如图2 4 所示，我们建立了一个直角坐标系，坐标原点位于矩形线框的中心， z 轴垂直向下指向地球深处。矩形框的长，宽分别为2 a ，2 b 。我们所要求的就是通有 简谐电流的矩形线框在z - - - 0 的地球表面处的任一点( x ，y ) 处，所产生的响应，即磁场 的垂直分量b ：。 为了求取( x ，y ) 处的响应，我们可以在( x t y ) 处放一个小的接收线圈r 。( 如图 x z 2 0p l a n e k u 一uj f u 一 一一2 b f 一 心a i 一 时 o ( 0 0 ) f，h 、 一 、 i l u 一u y 图z 4 放置在地球表面( z - - - o ) 的矩形框 f i g 2 4g e , o m c t t yo ft h er e c m n g u l a rl o o ps e t - u pa b v et h es u r f a c e 1 7 瞬变电磁法中矩形线框一维正演算法的研究 2 4 ) ，我们在其中通入简谐电流如“，由于线圈足够小，可以看作是一个磁偶极子， 由上一节的推导可知，在水平层状介质表面，磁偶极子产生的电场只有e + 分量。如果 我们将e + 平行于矩形四边的分量，即e + s i n0 ( 图2 4 ) 沿着矩形线框的四周进行积分， 就可以的到矩形线框的感应电压v 。由互惠定理可知，如果在矩形线框中通入电流 l e ，则在( x ，y ) 处感生出的电压为v 。由y = 一警s 就可以得到吃= 一去善a 经过如上所述的计算，我们最终得到了通有简谐电流的矩形线框在均匀水平层状 介质表面的响应( b z 分量) 的表达式： b ，= 上竺( + 曰+ c + d ) ( 2 2 1 ) “ 2 万 这里，a ，b ，c ，d 分别沿是矩形四边的积分式： 一( 屺孚f 糟 r = ( x - x ) 2 + ( 6 一y ) 2 口= 一( 口一z ) 了d y f r 黯 r = 、瓜i 万丽 c 叫帅，孚f 器 r 一瓜i 万而 一曲f 6 譬r 糟 ，= ( 口+ 工) 2 + ( y l _ y ) 2 2 3h a n k e i 积分一含b e s s e i 函数积分的陕速算法 式( 2 2 1 ) 中的积分a ，b ，c 和d 的内积分都是含有特殊函数b e s s e l 函数的积分。 这类积分在许多科学研究( 如物理、光学、地球物理、电磁学、电予工程) 中都经常 出现，所以，人们给这类积分起了一个好听名字h a n k e l 积分( 也称傅立叶贝塞尔 积分) 。在利用计算机对其计算时，由于b e s s e l 函数的震荡性，一般的数值积分方法 的收敛速度非常慢，无法使用。所以我们需要找到一种计算h a n k e l 积分的新方法。 下式( 2 2 2 ) ，是一个典型的n 阶h a n k e l 积分等式： 瞬变电磁法中矩形线框一维正演算法的研究 h 。( r ) = r ) j 。( m r ) d i n 我们用x = i n ( r ) ，y = i n ( 1 m ) ，对上式中的变量i l l ，r 做变量的替换， h n ( e x ) = 二 ( e - y ) 陋一，。e - y ) 坳 在等式( 2 2 3 ) 两边同时乘以r ，即e 。，得， 掰。( r ) = e x h 。 0 ，w 。是计算出来的数字滤波器的系数，a ：是与滤波器横坐标相关的对数 等间隔抽样点，其范围为3 1 a 。 5 0 ，i = 1 ，2 ，8 0 1 ，截断范围为l 锄， n z - 8 0 1 ， n 。，n ：的取值由计算结果所要达到的精度决定。 下图是利用求得的数字线性滤波器来计算h a n k e l 积分的个实例。我们采用一个 已知的h a n k e l 积分式( 2 2 5 ) ( 令a = 1 ) 来对滤波器进行验证。这里，我们计算的是 一阶h a n k e l 积分，即n = 1 ，输人函数为 ( m ) = ”1 2e x p ( - m 2 ) ，实际积分的真实结果可 由式”x p ( - h ，算出。图2 5 的上半部，显示了由数字线性滤波器计算的h a n k e l 积 分的结果；下半部则在对数坐标系中显示了计算结果与真实值之间的相对误差。 瞬变电髓法中矩形线框一维正演算法的研究 可 苦 击 重 罡 图2 5 利用数字线性滤波器计算h a n k e l 积分的结果 f i g 2 5e v a l u a t i o n o fh a n k e lt r a n s f o r mu s i n g l i n e a rd i g i t a lf i l t e r 由图我们可以看出计算结果的相对误差很小，可以满足我们实际应用的要求。而 且由于计算中避免了b e s s e l 函数的计算，只有简单的乘积和求和运算，所以运算的 速度是一般数值积分方法无法比拟的。 瞬变电磁珐中矩形线框一维正演算法的研究 2 4 频率域计算结果成像与比较 至此，我们已经推导出了均匀水平层状介质表面矩形线框的频域响应的公式 ( 2 2 5 ) ，式中包含了a ，b ，c 和d 四个积分式。其内层的h a n k e l 积分我们可以由上 一节介绍的数字线性滤波器的方法来求解，而外层的定积分我们则可以采用一般的辛 普森数值积分方法加以解决。为了验证我们推导的公式是否正确。我们利用公式 ( 2 2 5 ) 来计算矩形线框在如下几个给定参数的地球模型表面的频域响应。 三层模型参数设定( 发射框的放置参照图2 4 ) ： 矩形线框的半边长：a = b = 2 5 0m 测量点的x 坐标：x = 0m 第一层电导率： 盯1 = 0 0 1s i e m e n s m 第二层电导率：盯2 ：o 0 3s i e m e n s l m 第三层电导率： 盯 - o 0 0 1s i e m e n s m 第一层的厚度：h f 3m 第二层的厚度：h f 5 ，1 0 ，2 0 ，3 0 ，5 0 ，i 0 0m 线圈电流频率： f = 1 3 4 4h z 图2 6 显示了由上面设定的模型所计算的结果： 4 5 一p)。d三盖山v 瞬变电磁法中矩形线框一维正演算法的研究 0 0 25 50 7 5 1 0 0 y ( m ) 图2 6 通电矩形线框在3 层地球模型表面频域响应曲线。 ( a ) 幅度响应，( b ) 相位响应。 f i g 26 g e o m e t r i cs o u n d i n gc u f v e sf o ri n t e r i o ro f as q u a r el o o p p l a c e do n a i lt h r e e l a y e re a r t h y ( m ) o 加 锄 亭 舶 娜 卿 郴 枷 舶 ph#e(de$一 言p)d暑=duj 壁銮皇堡整! 笙星垡堡= 堡垩塑蔓鎏墼堡塞 05 0 1 0 0 1 5 0 y f m ) 图2 7 m p o d d a r l 9 8 3 年计算得到的矩形线框在3 层模型表面的频域响应 ( a )为幅度响应，( b ) 为相位响应。 f i g 2 7g e o m e t r i cs o u n d i n g c u r v e sf o rt h ei n t e r i o ro fa s q u a r el o o pp l a c e d o n a n t h r e e l a y e re a r t hm a d eb y m 。p o d d a ri n1 9 8 3 图2 7 是由m p o d d a r 在1 9 8 3 年计算出来的在均匀水平层状介质表面的通电矩形 线框响应曲线，他所用来计算的地球模型，与我们在本节开头所给出的模型一样。图 2 6 ( a ) 和图2 7 ( a ) 所示的幅度曲线，都是用相同设置的线框在全空问( 无层状介质 存在的自由空间) 的响应归一化的结果。从图2 6 和2 7 的比较中我们可以看出，无 论是响应曲线的形状还是数值