（地球探测与信息技术专业论文）局部波数在磁异常解释应用中的方法技术讨论.pdf

摘要 二维局部波数法可以用来反演磁性体的特征点位置，顶界埋深和中心埋深，并已得 到了一定的应用。该方法的优点是不依赖于磁化强度和磁倾角的大小，缺点是只适用于 特定的场源模型，即倾斜的断层，倾斜的薄板和水平圆柱体。数值计算技术将影响反演 结果的准确性。 本文在前人研究基础上，对局部波数反演计算过程中存在的一些方法技术进行了讨 论。一是频率域处理中扩边方法的选择，通过对余弦扩边方法和最小曲率扩边方法之间 的对比，认为最小曲率扩边方法优于余弦扩边方法。二是在频率域中计算磁异常高阶导 数时存在高频放大的影响，为了克服该影响，通常的做法是选择一个低通滤波器，而低 通滤波器参数的选择难于掌握，为了解决该问题，提出了向上延拓某一高度来计算高阶 导数或局部波数的方法技术，不但解决了高阶导数计算过程中的高频放大问题，而且也 解决了低通滤波器参数的选择问题。利用f o r t r a n9 0 标准编制了相应程序，试算了同一 模型不同埋深、不同剖面长度、不同点距、不同磁化强度大小和方向对反演结果的影响； 试算了不满足二维局部波数反演条件的三个特殊模型时反演结果的影响。最后，用我国 南海北部的航空磁测数据，反演了异常体的埋深，并对反演结果进行了解释，得出了比 较满意的结果。 关键词：局部波数，反演，磁异常，频率域，最小曲率，向上延拓 a b s t r a c t t h et w o - d i m e n s i o nl o c a lw a v e n u m b e rm e t h o dh a sb e e nd e f i n i t e l yu s e dt oi n v e r s e b o u n d a r yp o s i t i o na n db u r i a ld e p t ho ft o pb o u n d a r yo ft h em a g n e t i cs u b s t a n c e t h ea d v a n t a g e o ft h i sm e t h o di st h a ti td o e s n tg oa l o n gw i t ht h em a g n e t i z a t i o na n dt i p ，a n dt h ed i s a d v a n t a g e i st h a ti ti so n l yf i tf o rp a r t i c u l a rs o u r c em o d e l s ，s u c ha ss l o p p i n gf a u l tm o d e l ，d i p p i n gt h i n s h e e tm o d e la n dh o r i z o n t a lm o d e l s o m ec o m p u t a t i o n a lt e c h n i q u e so fn u m e r i c a lv a l u ew i l l d i r e c t l ya f f e c tt h ea c c u r a c yo ft h ei n v e r s i o nr e s u l t a n td u d n g t h ei n v e r s i o n t h eo r i g i n a li sb a s e do nt h er e s e a r c ho fp r e d e c e s s o r s ，a n dt a l k sa b o u ts o m em e a n sa n d t e c h n i q u e so ft h el o c a lw a v e n u m b e ri n v e r s i o n f i r s t l y , h o wt oc h o o s ee x t e n s i o nm e t h o do f 仃e q u e n c yd o m a i np r o c e s s i n g ：b ym e a n so fc o m p a r i n gc o s i n ef u n c t i o ne x t e n d i n ge d g ea n d m i n i m u mc u r v a t u r ee x t e n d i n ge d g e ，w et h i n km i n i m u mc u r v a t u r ee x t e n d i n ge d g ei ss u p e r i o r t oc o s i n ef u n c t i o ne x t e n d i n ge d g e s e c o n d l y , t h e r ew i l lb ea f f e c t i o no fe n h a n c e dh i g h 1 j r e q u e n c yw h i l ew ec a l c u l a t et h eh i g h e rd e r i v a t i v eo fm a g n e t i ca n o m a l y , f o rt h es a k eo f o v e r c o m i n gt h ea f f e c t i o n ，w ep r o p o s e am e t h o dw h i c hu p c o n t i n u e ss o m ea l t i t u d et oc a l c u l a t e h i g h e rd e r i v a t i v eo rl o c a lw a v e n u m b e r a n dt h i sm e t h o dn o to n l yr e s o l v e s t h ep r o b l e mo f e n h a n c e dh i g h 仃e q u e n c y ，b u ta l s or e s o l v e st h ep r o b l e mo fh o wt oc h o o s et h ep a r a m e t e r so f l o wp a s sf i l t e r w eu s ef o r t r a n9 0t op r o g r a mt h eh o m o l o g o u sp r o g r a m s ，a n dt r yt oc a l c u l a t e t h ea f f e c t i o no fd i f f e r e n tb u r i a ld e p t h ，s e c t i o nl e n g t h ，d o tp i t c ha n dm a g n e t i ci n t e n s i t yt ot h e i n v e r s i o no ft h es a m em o d e l t h e nw et r yt oc a l c u l a t et h ea f f e c t i o no ft h et h r e ep a r t i c u l a r m o d e l sw h i c hd o n tm e e tt h ec o n d i t i o no ft h et w o d i m e n s i o nl o c a lw a v e n u m b e rm e t h o d f i n a l l y , w eu s et h en o r t h e r na i r b o r n em a g n e t i cd a t ao fs o u t hc h i n as e a a n di n v e r s et h eb u r i a l d e p t h o ft h ea n o m a l o u sb o d y ；a l s ow eo b t a i nc o m p a r a t i v e l ys a t i s f a c t o r yr e s u l t k e yw o r d s ：l o c a lw a v e n u m b e r ；i n v e r s i o n ；m a g n e t i c ；f r e q u e n c yd o m a i n ；m i n i m u m c u r v a t u r e ；u p w a r dc o n t i n u a t i o n n 论文独创性声明 本人声明：本人所呈交的学位论文是在导师的指导下，独立进行研究工 作所取得的成果。除论文中已经注明引用的内容外，对论文的研究做出重 要贡献的个人和集体，均已在文中以明确方式标明。本论文中不包含任何 未加明确注明的其他个人或集体已经公开发表的成果。 本声明的法律责任由本人承担。 论文作者签名： 饵彩7 砗易月；日 论文知识产权权属声明 本人在导师指导下所完成的论文及相关的职务作品，知识产权归属学 校。学校享有以任何方式发表、复制、公开阅览、借阅以及申请专利等权 利。本人离校后发表或使用学位论文或与该论文直接相关的学术论文或成 果时，署名单位仍然为长安大学。 ( 保密的论文在解密后应遵守此规定) 论文作者签名：但彩) m 缉矿夕月衫日 导师签名： 锄饰跏狰以月。弓日 长安人学硕士论文 第一章绪论 1 1 本论文研究的背景和意义 已知磁场的空间分布特征来确定地下所对应的场源体特征，如确定磁性体的赋存空 间位置、磁性磁化强度的大小和方向等。这些通常被称为磁异常的反演问题。显然，这 是解决各种地质勘探等实际问题的重要环节。 磁异常反演过程一般分为定性及半定量解释和定量解释两个阶段，这两者是互相关 联，互相辅助的。只有通过正确的定性及半定量解释取得初步对磁性体形状、产状、及 其所引起的地质原因进行判断之后，才能合理地选用定量计算的公式和方法，以便进一 步得到更完善的解释结果。因此，这两个步骤实质上前者为后者的基础，后者又是前者 的继续和深入。 由于磁性体的空间赋存状态和磁化特征( 包括均匀和非均匀磁化) 的不同。因而， 在磁异常解释中存在以下两个较为普遍的问题【4 】。 ( 1 ) 场源体非均匀磁性问题 在自然界中因场源所处的地质环境不同，其非均匀磁性是较为普遍的现象。而反演 解释均依赖于场和场源体之间一定的函数关系式求解，非均匀性无论是对简单规则形体 或复杂形体均会使其函数关系式更加复杂。因此，通常采取两种方法处理：一是对具有 一定埋深的场源在反演解释过程中假设其为均匀磁化；另一是对一些非均匀磁化场源进 行分块切割处理，对每一块视为相对均匀磁化，并且彼此之间不磁化。 ( 2 ) 反演的多解性问题 多解性问题是地球物理勘探反演解释中共同存在的问题。而磁异常反演的多解性要 比重力探勘更为复杂。因为决定磁异常特征的两个主要因素不仅与磁化场大小、方向有 关，还与场源的形态有关。当这些因素不同组合可以获得类同的磁异常分布特征。同样， 利用线性化后的线性反演法求近似解，虽便于计算机自动实现，但其多解性问题仍然不 能克服，而造成多解的因素可能更多。为了避免或者减少磁异常反演的多解性问题，必 须要遵循磁异常解释步骤，并充分利用地质资料和其它地球物理资料进行综合推断解 释。 近年来，我国进行了广泛的区域地面重力、地面磁测及航空磁测，获得了大量高精 度、高分辨率的重、磁数据，面积几乎覆盖了大部分国土，因此充分有效的利用这些资 料，从中提取更多的地质信息，进一步提高解释的准确性和精确度，具有重要的现实意 第一章绪论 义。又由于反演场源体非均匀磁化问题和反演的多解性问题，磁异常的反演方法一直是 地球物理学研究的热点问题。 1 2 国内外研究现状 首先讨论简单规则形体磁异常反演方法。大都以其对应的磁场解析表达式为基础， 分析磁源体的几何形状及磁性体参量与磁异常分布特征之间的关系，导出反演计算公 式。因此，在应用这些计算公式进行反演时必须注意以下几个问题【1 0 】： ( 1 ) 注意用于定量计算剖面的选择，对二度体异常，要选择与异常走向垂直的接 近异常中心的磁异常剖面；对似二度异常应选择与其异常长轴垂直的中心剖面。 ( 2 ) 必须在初步判断磁性体形状之后，才能为合理地选择用相应形体的反演计算 公式提供依据。 ( 3 ) 其它与正演问题的假设条件一致。 目前磁异常的反演方法有以下几类i l o l ： ( 1 ) 特征点法 它是利用磁异常曲线上一些特征点，如极大值、半极值，1 4 极值，拐点，零值点 及极小值等坐标位置和坐标之间的距离，求解磁源体参量的方法称为特征点法。该方法 只能是针对几个规则形体，在自然界很难有这样的规则形体，故实用性不大。 ( 2 ) 切线法 它是利用过异常曲线上的一些特征点，如极值点和拐点的切线之间的交点坐标间的 关系来计算磁性体产状要素的方法。该方法简便，快速、受正常场选择影响小，在航磁 丁异常的定量解释中曾得到广泛应用。但该方法是一种近似的经验公式，另外受使用者 选点位置的影响，一般得出的只是近似解。 ( 3 ) 磁异常梯度的积分法 积分法求解场源反演问题的内容较广泛，通常有两类：一类是对各种简单规则形体， 利用z 。及h 。异常曲线与x 轴之间所围面积或部分面积来实现求解反问题；另一类是对 任意形状有限截面的二度、三度体，利用磁异常的一阶或二阶矩的积分，如卜”z 。d x 型 的积分求解反问题。 ( 4 ) 由磁异常的希尔伯特( h i l b e r t ) 变换求二度磁性体参量：它是利用磁异常2 。b ) 的希尔伯特变换h b ) 曲线与原2 。g ) 曲线，即它们的水平导数异常曲线仨。g ) 、h ，g ) ) 的一些特征点来进行反演解释，可以得到多种参量的反演结果。 2 长安大学硕士论文 复杂磁异常的反演方法，根据其计算原理可分为频率域和空间域两类方法，如空间 域、频率域的最优化选择法及磁性界面的反演。起伏地形下的磁异常及磁异常频谱直接 反演多种参量的解释法。磁异常的线性反演等【4 】。 ( 1 ) 快速自动反演方法 2 0 世纪8 0 年代以来，发展了应用微机对磁异常剖面与面积资料自动反演深度的技 术，可反演出剖面或测区内一系列深度点，再结合地质或其它地球物理资料勾绘地质断 面；也可用此深度资料作初值，提供最优化反演的初始模型之用。磁异常剖面自动反演 方法有w e r n e r 反褶积法、总梯度模法、局部波数法；面积性的磁异常自动反演方法有 欧拉法。 ( 2 ) 最优化方法 这类方法的实质就是选择法，即使模型体引起的场值理论曲线与实测曲线尽可能地 拟合，取得拟合得最好的模型参数作为解释的结果。 ( 3 ) 三维物性反演方法 此方法是指利用观测磁异常，在确定产生这些异常的构造范围内按纵向与横向划分 成许多小的单元组合，选择反演方法，反演每个单元的磁性参数。通过单元组合的磁性 变化，勾画出地质体的轮廓，并给出不同深度的磁性层所反映的地质体深度切片，故具 有层析意义，也有人称为“物性成像”。 1 3 局部波数法的研究现状及论文的研究内容 本论文所研究的反演方法是快速自动化方法中的局部波数法，即二维( 剖面) 磁测 数据确定场源位置的自动化反演方法。上世纪九十年代该方法得到了广泛的关注，近十 年来此方法得到了很大的发展。局部波数法适用于特定的场源模型，即倾斜的断层，倾 斜的薄板和水平圆柱体，利用二维磁测数据的导出参数，可以反演异常体的特征点和埋 深。 这种方法利用了局部波数的振幅通常不依赖于磁化强度和磁倾角的特点，根据 s m i t h ( 1 9 9 8 年) 等人研究结果认为，一阶和二阶局部波数不需先验信息就可求出磁性 体的埋深。因此在确定磁性体位置上比磁场值更直接，并且这种自动化方法，在磁测解 释过程中非常的快捷简单，对于前期的勘探工作具有非常重要的意义。 ( 1 ) 频率域快速f o u r i e r 变换时数据扩边方法的选择； ( 2 ) 频率域各阶导数计算时高频放大问题的解决方法； 3 第一章绪论 ( 3 ) 局部波数计算中的影响因素讨论。 1 4 本论文的研究成果 完成算法程序，并对程序进行了改进。演算了大量的模型，并对结果进行对比，分 析总结了每种模型在算法中的参数设定。根据实际南海资料处理结果，分析了研究区域 异常体的基本信息。 ( 1 ) 对二度体磁异常局部波数反演方法在实际计算中存在的问题进行了讨论，并 给出了改进措施。为了提高计算精度和克服g i b b s 效应，应采用效果好的数据扩边方法， 这里采用了最小曲率扩边方法。 ( 2 ) 在计算过程中为了压制高频干扰对计算结果的影响，我们试用了两种滤波方 式消除高频干扰影响，一种是低通滤波器，另一种是向上延拓。对这两种方法的滤波效 果进行比较后，采用了消除高频干扰效果更好的向上延拓方法。 ( 3 ) 通过三种理论模型的演算，反演得到了相符合实际的埋深和特征点； ( 4 ) 通过大量的单一模型和组合演算证实了局部波数法不依赖于磁化强度和磁倾 角，分析总结了单一模型和组合模型不同深度，不同剖面长度、不同点距的反演规律； ( 5 ) 对南海实际资料进行了试算，试算结果与实际地质情况吻合，解决了南海地 区斜磁化对磁异常解释带来的影响。 4 长安大学硕士论文 第二章二维局部波数法的基本原理 2 1 二维局部波数的定义 2 1 1 一阶局部波数 对二度体磁异常丁g ，z ) ，解析信号的模陋 ，z ) i 定义为1 3 1 ： p ，z ) i = 解析信号的相位口g ，z ) 定义为【1 】： ( 2 1 1 ) 捌一以_ 筹) ( 2 工2 ) 水平局部波数七，g ，z ) 被定义为解析信号相位口g ，z ) 的变化率1 1 1 ，其表达式为： t g ，z ) 掣 d x 1 + a t l o z ) 2 ( 分( 誓) 2 ( 2 1 3 ) 2 1 2 二阶局部波数 日b ，z ) 的定义式为： 蚺伽- 1 糕 亿， 二阶局部波数七。g ，z ) 被定义为口b ，z ) 的变化率【1 7 1 ，其表达式为： 5 第二章二维局部波数法的基本原理 k g ，z ) 。掣 1 + ( 罴) 2 ( 塞) 2 ( 塞) 2 + ( 害) 2 a 2 丁a 3 ra 2 丁a 3 r 0 x o z 以a z 2以20 2 0 x 2 2 + 阿 2 2 三种模型磁异常正演公式 ( 2 1 5 ) 2 2 1 倾斜断层模型磁异常正演计算公式 n a b i g h i a n ( 1 9 7 2 年) 给出了斜断层模型的水平和垂直梯度f 1 4 1 ，表达式为： n 2 脓s 证d 卜c o s 。埘) + s m 。划) l n 昙】 ( 2 - 2 1 ) 其中k 是磁化荤，f 是地磁场，c 。1 一c o s zf s i n ：口，f 是地磁倾角，口是磁偏角， ，。a r c t a n f 堕1 ，d 是断层倾斜面与上界面的夹角，h c 是断层上界面的埋深，测线坐 、c 0 s 口， 标系的原点g = o ) 是斜断层的特征点( 图2 1 ) 。 图2 1 倾斜断层模型示意图 6 长安大学硕士论文 2 2 2 倾斜薄板模型磁异常正演计算公式 r e f o r d ( 1 9 6 4 年) 给出了斜薄板的磁异常的表达式为【1 8 】： 删。2 k f c 地号骞芈型 ( 2 2 2 ) 其中k 是磁化率，f 是地磁场，c = 1 一c o s 2i s i n 2 口，i 是地磁倾角，口是磁偏角， ，。a r c t a n f 堕1 ，d 是薄板倾斜面与上界面的夹角，是斜薄板的厚度， 。是斜薄板 c o s ( y 上界面的埋深，测线坐标系的原点g 0 ) 是斜断层的特征点( 图2 2 ) 。 在磁法勘探中薄板和厚板的概念是相对的，当 一抛出 丁b ) ；f - 1 防o ) 】；渺撕咖 ( 3 1 1 ) 频率域中的求导公式( 这里z 坐标方向向下为正) i s f l 掣l ；捌订口g ，z ) = l c u ) l o r , i 吖掣卜2 1 忡洲= q 挈 = 抚话 掣卜驰砷酬坷“， 吖警卜z 万l u l f 掣卜2 神列出“) 叫等 姊p 掣 甜呻删枣“) 吖掣 _ - 砒2 可掣 - _ 酬帅酬办“) 吖掣卜叫掣卜舻桃力m f 等卜叫掣卜叫忡删础比) ( 3 化) 1 0 长安大学硕士论文 3 1 2 频翠域求导步骤 第1 步：对网格化的数据t ( x ，z ) 进行扩边处理； 第2 步：对扩边后的数据进行快速f o u r i e r 变换，给出z 0 ，z ) 的频谱f p g ，z ) 】； 第3 步：计算各阶导数的频率响应d ( u ，i ，乏，t ) ； 第4 步：计算f 口 ，y ，z ) 。 ，亏，艺，乏) _ f 【器】： 第5 步：进行u d e r 反变换得到罴o to t o t 2 f。1【f【器ot o t o t 】； 12 。 ii 12 。i i 第6 步：去掉扩边部分，输出! 盟即可。 o t o t 2 o t n 3 2 扩边问题讨论 以2 为基的快速f o u r i e r 变换【8 1 要求磁异常数据z k ) 的个数m 2 七，足为整数，如果 m 不足2 个，则需要将数据丁似) 做扩边处理。如果数据丁b ) 的个数，l 一2 ，m 满f f z 2 七 个，但为了消除f o u r i e r 变换中的g i b b s 现象【2 4 1 ，我们也需要将数据扩边。目前常用的 扩边方法是余弦衰减扩边方法，但该方法在接边处数据不光滑( 导数不连续) ，导数计 算结果误差很大。 3 2 1 扩边方法及计算公式 3 2 1 1 余弦扩边方法及计算公式 对磁异常数据t ( m ) 进行扩边，将m 个数据扩为咒个数据( 聆为2 的幂次方) ， z g 。：厅：一1 ) 为左边要扩边的数据，z ( 咒：刀，) 为实际数据，z g ，一l ：n 。) 为右边要扩边的数 据， 咒一以4 一n l + 1 。 若，l 。；o ，咒。2 面l o g m + l ，厅：；鱼专型，l ，- n 2 + m - l , n 4 ；万一1 。余弦扩边公式为： 丁o 。) ；生峰型 r g 。) 一丁g 。) 丁锄1 ) - c o s 蹦和哪：) r n ) f f it 一 蹦和4 ) _ 丁忆) 1 1 咒1 ，l ? 1 2 咒3 n 厅4( 3 2 2 ) 第三章计算方法及程序设计 3 2 1 2 最小曲率扩边方法及计算公式 利用余弦扩边后，计算出r 在x 和z 方向各阶偏导数图形边部出现了高频干扰，而 且随着偏导数的阶次越高，高频干扰越严重。为了解决这个问题，我们采取了最小曲率 扩边方法【1 0 1 ，最小曲率方法的使用，可以减少丁在x 和z 方向的各阶偏导数图形边部出 现的高频干扰。 最小曲率一维微分方程 尹v 2 t ( x ) ) = 妄( 訾) = 訾一圭f t a ( x - ，2 习 边界条件 塑掣：0 a 3 t _ ( x ) ；0( 3 2 4 ) 3 x 。缸 、7 其差分迭代格式 z o ) 。吾仁o + 2 ) + t ( i 一2 ) 一4 p o + 1 ) + r ( i 一1 ) b ( i - 1 , 2 ，m ) ( 3 2 5 ) 其差分边界条件 r ( i 一1 ) 一2 t ( i ) - t ( i + 1 ) ( f 一1 ) r ( + 1 ) 一2 t ( i ) - t ( i - 1 ) ( i m ) 1 f 三g 一苎? 5t ( ，i + 2 ，) - 2 r t ( ，i + 1 ) 一t ( 一1 ) ( = 1 ( 3 2 6 ) i 丁o + 2 ) = 丁( f 一2 ) + 2 f o + 1 ) 一丁( f 一1 ) jo m ) r 。一7 3 2 2 余弦和最小曲率扩边后的z 的导数图形对比 扩边问题解决后我们就可以根据公式 ( 3 1 2 ) 编写程序对磁异常r 在x ，z 方向进行 求导。如图3 1 所示为一个倾斜薄板的磁异 常剖面图，倾斜薄板模型的特征点位于剖面 第6 4 个测点，埋深在地下2 0 0 m 处。图中x 轴为测点号。分别使用余弦扩边和最小曲率 扩边进行扩边后计算各阶导数。 1 2 图3 1 磁异常rg z ) 剖面图 长安大学硕士论文 ( a ) 余弦扩边计算结果 ( b ) 最小曲率扩边计算结果 图3 2 不同扩边方法计算出的坚剖面图 a x ( a ) 余弦扩边计算结果 ( b ) 最小曲率扩边计算结果 图3 3 不同扩边方法计算出的望剖面图 a z 图3 2 和图3 3 分别为一阶水平和垂直导数的计算结果，磁异常丁在x ，z 方向的一阶 导数图形曲线平滑，余弦扩边和最小曲率扩边的效果基本相同。 第三章计算方法及程序设计 ( a ) 余弦扩边计算结果 ( b ) 最小曲率扩边计算结果 图3 4 不同扩边方法计算出的0 2 了t 剖面图 o x 。 ( a ) 余弦扩边计算结果 ( b ) 最小曲率扩边计算结果 图3 5 不同扩边方法计算出的兰剖面图 1 z o x o z 1 4 长安大学硕士论文 ( a ) 余弦扩边计算结果 ( b ) 最小曲率扩边计算结果 图3 6 不同扩边方法计算出的0 2 t 剖面图 a z 图3 4 ，图3 5 和图3 6 分别是磁异常z 在x 方向的二阶导数，z 方向的二阶导数以 及x ，z 方向的二阶混合导数。对比三幅图可以看出，余弦扩边后的图形曲线边部出现了 不光滑现象，有一些高频干扰，这些会对计算结果产生明显的影响，而最小曲率扩边的 图形曲线边部稍有不光滑，相对余弦扩边的计算结果要好。 ( a ) 余弦扩边计算结果 最小曲率扩边计算结果 图3 7 扩边后的三 a x a z 第三章计算方法及程序设计 ( a ) 余弦扩边计算结果( b ) 最小曲率扩边计算结果 图3 8 不同扩边方法计算出的三剖面图 a x 讹 ( a ) 余弦扩边计算结果( b ) 最小曲率扩边计算结果 图3 9 不同扩边方法计算出的0 3 7 t 剖面图 a z 图3 7 ，图3 8 和图3 9 分别是磁异常丁在z ，x ，z 方向的三阶混合导数，毛z ，z 方向的 三阶混合导数和z 方向的三阶导数。对比三幅图可以看出，余弦扩边后的图形曲线边部 出现了严重的不光滑现象，高频干扰严重，最小曲率扩边的图形曲线边部也出现了不光 滑现象，但比余弦扩边的结果要好。 通过以上试验，扩边方法的选择只是提高计算精度的一个方面。仍需要考虑频率放 大问题。 根据公式( 2 1 3 ) 和( 2 1 5 ) 分别用余弦扩边和最小曲率扩边方法的计算结果，计 算一阶局部波数k 。和二阶局部波数k 。，并对图形进行比较。 1 6 长安大学硕士论文 ( a ) 余弦扩边计算结果( b ) 最小曲率扩边计算结果 图3 1 0 一阶局部波数k ，的计算结果剖面图 图3 1 0 是余弦扩边和最小曲率扩边后计算的一阶局部波数七，对比两图可以看出， 由于余弦扩边后，磁异常丁在x ，z 方向二阶导数的高频干扰，导致了k ，图形曲线边部严 重的不光滑现象。最小曲率扩边后计算的尼，图形边部，由于磁异常丁在x ，z 方向二阶导 数也存在高频干扰，所以使得图形曲线边部也出现了不光滑现象，但没有余弦扩边的严 重。 ( a ) 余弦扩边计算结果 u 鹏 0 0 3 o 0 2 一 i 、 o o l m6k4。r刈。 1 f u 0 2 0 4 0 6 08 01 0 0 1 2 0 1 4 0 ( b ) 最小曲率扩边计算结果 图3 n 二阶局部波数k 。的计算结果剖面图 图3 1 1 是余弦扩边和最小曲率扩边后计算的二阶局部波数k ，对比两图可以看出， 由于余弦扩边后，磁异常r 在x ，z 方向二阶、三阶导数的高频干扰，特别是三阶导数严 重的高频干扰，导致了尼。图形曲线边部严重的不光滑。最小曲率扩边后计算的七。图形 边部，由于磁异常z 在x ，z 方向二阶、三阶导数也存在高频干扰，特别是三阶导数出现 了比较明显的高频干扰，所以使得的图形曲线边部也出现了不光滑现象，但没有余弦扩 1 7 第三章计算方法及程序设计 边的严重。 比较余弦扩边和最小曲率扩边这两种方法，能够得到最小曲率扩边方法明显优于余 弦扩边方法，所以我们选择最小曲率方法对磁异常丁的数据进行扩边处理。分析最小曲 率扩边后丁在x ，z 方向上的各阶导数，一阶导数曲线比较光滑没有明显的高频干扰，二 阶导数曲线边部出现了稍许的高频干扰，三阶导数曲线边部出现比较严重的高频干扰。 正是这些高频干扰造成了一阶和二阶局部波数计算结果的严重变形，其整条曲线出现不 光滑现象高频干扰，这使得无法计算k 。一k ，的极大值，即无法反演正确的埋深吮， 无法准确得n - 度体的特征点和构造指数n 。 由于在实际计算过程中出现的高频问题，使得我们无法实现局部波数法的反演结 果，所以下面我们要分析讨论消除高频干扰的办法。 3 3 高频干扰解决办法 虽然利用最小曲率方法代替余弦扩边，可以减少偏导数边部的高频干扰，但是边部 的高频干扰依然比较严重，特别是z 的三阶导数，由于这个原因使得异常体的埋深计算 结果精度很低，为了解决这个问题我们采用了加低通滤波器的方法。 3 3 1 低通滤波器 加低通滤波器就是给多阶导数乘以滤波因子，根据不同的导数阶次加不同的滤波因 子。低通滤波器的使用可以使得导数边部高频干扰降低，提高计算精度。 3 3 1 1 低通滤波因子 6 】 一阶水平和垂直导数滤波因子为： 办0 ) 。o x p - 锄2 a b2 ) j ( 3 3 1 ) 二阶导数滤波因子为： 驴： ) = 阮0 ) 】2( 3 3 2 ) 三阶导数滤波因子为： 丸0 ) 一队0 ) r ( 3 3 3 ) 公式中a 的取值越大，噍0 ) 随着z l 增大而下降的速度越快。快速下降可以使高频成 分较彻底的压制，但同时低频有益成分也受到一定程度的压制而畸变，即低频成分畸变 也很厉害。因此，原则上是a 值越大，压制高频成分同时使低频成分不受畸变，这样的 1 8 长安大学硕上论文 滤波因子效果最好【翻。 3 3 1 2 最小曲率扩边和加低通滤波因子的丁的导数图形对比 利用图3 1 给出的磁异常数据和最小曲率扩边方法，下面给出不加低通滤波器和加 低通滤波器的对比结果。 ( a ) 未加低通滤波因子计算结果 ( b ) 加低通滤波因子计算结果 图3 1 2t 的一阶导数坚计算结果对比 缸 、 厂 。， 7 | v ( a ) 未加低通滤波因子计算结果 ( b ) 加低通滤波因子计算结果 图3 1 3t 的一阶导数竺计算结果对比 彪 图3 1 2 和图3 1 3 为一阶导数计算结果，通过对比可以看出，加了低通滤波因子的 图形和没有加低通滤波因子的图形差别不明显。 1 9 5 2 5 l 5 o 5 l 5 2 5 互 l 仉 吨 一叱 一吃 第三章计算方法及程序设计 ( a ) 未加低通滤波因子计算结果 ( b ) 加低通滤波因子计算结果 图3 1 4z 的二阶导数! 萋计算结果对比 缸 ( a ) 未加低通滤波因子计算结果( b ) 加低通滤波因子计算结果 图3 1 5r 的二阶导数喜! 计算结果对比 d x z ( a ) 未加低通滤波因子计算结果( b ) 加低通滤波因子计算结果 图3 1 6z 的二阶导数重；计算结果对比 a z 长安大学硕士论文 由图3 1 4 ，图3 1 5 和图3 1 6 可以看出，加了低通滤波因子的磁异常z 的二阶导数， 图形边部比较光滑，减少了原有的高频干扰，比不加低通滤波因子要好。 c a ) 未加低通滤波因子计算结果 ( b ) 加低通滤波因子计算结果 图3 1 7t 的三阶导数! 三计算结果对比 a x a z ( a ) 未加低通滤波因子计算结果 ( b ) 加低通滤波因子计算结果 图3 1 8z 的三阶导数三计算结果对比 o z a x 2 1 第三章计算方法及程序设计 ( a ) 未加低通滤波因子计算结果 ( b ) 加低通滤波因子计算结果 图3 1 9z 的三阶导数o a _ t 计算结果对比 a z 由图3 1 7 ，图3 1 8 和图3 1 9 可以看出，加了低通滤波因子的磁异常r 的三阶导数， 图形边部的高频干扰明显改善，减少了原有的高频干扰，比不加低通滤波因子要好。 ( a ) 未加低通滤波因子计算结果 ( b ) 加低通滤波因子计算结果 图3 2 0 一阶局部波数七。计算结果对比 由图3 2 0 可以看出，通过加低通滤波因子后磁异常z 在x ，z 方向各阶导数高频有了 明显的消减，所以一阶局部波数k ，图形边部高频也有了明显的消减，但是仍然还有一些 高频干扰。 长安大学硕士论文 u u i 0 0 3 o 0 2 、 0 o l 4t l 。ik r 。w u i _ o u 02 0 4 06 08 0 1 0 0 1 2 01 4 0 u i q 0 1 2 o 1 0 0 8 0 0 6 0 0 4 0 0 2 l 。 1 、 ii i i i | | l f ，、。，删n u 02 04 06 08 01 0 01 2 01 4 0 第三章计算方法及程序设计 i u 9 8 7 6 5 4 3 2 0 l l v v 、叫yn u 02 0 4 06 0 8 0 1 0 01 2 01 4 0 ( a ) 未加低通滤波因子计算结果 ( b ) 加低通滤波因子计算结果 图3 2 3 n 七计算结果对比 由图3 2 3 可以看出，未加低通滤波因子计算出构造指数咒。图形曲线严重不光滑， 说明高频干扰非常严重，加低通滤波因子图形的高频干扰明显减小，但因为曲线的不光 滑，所以无法判断极大值，即无法计算出构造指数厅。 根据以上结果可以得出结论，加低通滤波因子对于消除高频干扰的确有明显的效 果，但是到了磁异常丁的三阶导数，其图形曲线仍然会出现一些高频现象，这导致了无 法反演结果。分析原因是因为低通滤波因子中a 取值非常难以确定，使得消除高频干扰 的效果难以把握。所以下面我们考虑采用另一种选取参数较为方便的低通滤波方式，即 向上延拓。 3 3 2向上延拓 向上延拓可以压制高频成分，突出低频成分，这也属于低通滤波。向上延拓的作用 可以消除随机干扰的影响，将局部异常展平，使曲线相对光滑：向上延拓将使局部异常 的幅值降低，因而达到了曲线相对光滑去除高频干扰的目的【2 3 1 。 3 3 2 1向上延拓因子 向上延拓因子表达式为【3 l ： h 0 ) ：e - 2 z r i “陋 ( 3 3 4 ) 其中u 为频率因子，h 为向上延拓高度。 3 3 2 2 向上延拓高度的确定问题 为了选取更合理的延拓高度，我们用不同模型不同埋深进行试算，下面模型的点距 都为1 0 m ，剖面长度为1 0 0 0 m 。 长安大学硬论文 ( i ) 倾斜薄板模型进行试算，模型特征点位置在第5 0 号测点，埋深为2 0 0 m ，延拓 不同高度结果如下： 0i 02 03 04 0 5 0 6 0 7 08 09 0i 0 0 图3 2 4 向t 延拓不同高度的i ( k 。- k ，) 图3 2 s 延拓埋深和反演埋深图3 2 6 构造指数 1 ) 由图32 4 町知，延拓高度从o m 到i o o mm 隔为2 0 m 。延拓不同高度1 ( k 。一 。) 的极小值都在5 0 号测点，则异常体的特征点就在此处。延拓高度为o m 时特征点的高频 干扰严重。延拓高度为2 0 m ，4 0 m 时，不能有效地消除高频干扰。 2 ) f l i 罔32 5 可知，反演埋深的值都趋于一个定值时候就是磁性体的真实埋深，所 以磁性体的埋深麻该约为2 0 0 m ，这与真实值2 0 0 m 相符。2 0 m l o o m 的反演埋深平均值的 精度为0 0 4 6 ，符合精度5 要求。 3 ) 由图3 2 6 可知，构造指数在极大值处等于1 0 0 2 ，则可以推断是倾斜的薄板。 度一 0 o o o 叭 雾。到 喜!咖瑚舌|螂蓦|耋言啪。 幽盈 二 0 第三市计算方挂厦程序设计 ( 2 ) 水平圆柱体模型半径为5 0 m ，测点到中心距离为l o o m ，水平圆柱体中心位置 在第5 0 号测点，延拓不同高度结果如下； 1 0 0 0 9 0 0 9 删 8 0 0 瓢u厂研 延拓高度 7 0 0 l n 、 【o 2 0 6 0 0 弋： 4 0 i 6 。 5 0 0 赢j 卜 l 一8 0 4 0 0 3 0 0 2 0 0 1 0 0 1 0 0 j l o2 04 06 08 01 0 0 图3 2 7 向上延拓不同高度的1 耻。一t ，) 图3 2 8 延拓埋深和反演埋深图3 2 9 构造指数 1 ) 由图3 2 7 可知，延拓高度从o m 到l o o m 间隔为2 0 m 。延拓不同高度1 似，一 ，) 的极小值都在5 0 号测点则异常体的特征点就在此处。延拓高度为o m 时特征点的高频 干扰严重。 2 ) 由图3 2 8 可知，反演埋深的值都趋于一个定值时就是磁性体的真实埋深，所咀 磁性体的埋深应该约为l o o m ，这与真实值l o o m 相符。2 0 m l o o m 的反演埋深平均值的精 度为0 0 7 8 ，符合精度5 要求。 3 ) 由图3 2 9 可知，构造指数在极大值处等于20 0 3 ，则可以推断是水平圆柱体。 ( 3 ) 水平圆柱体模型半径为l o o m ，测点到中心距离为2 0 0 m ，水平圆柱体中心位 蚕 长安大学顾士论空 置在第5 0 号测点，延拓不同高度结果如下 02 04 06 08 0l o o1 2 0 图3 3 0 向上延拓不同高度的1 ( k 。一k ，) 0 o 1 0 。 图3 3 1 延拓埋深和反演埋深 延拓高度 | 二刻 j 二刮 图3 3 2 构造指数 i ) 由图33 0 可知，延拓高度从0 m 到l o o m 间隔为2 0 m 。延拓不同高度1 ( k 。一k ，) 的极小值都在5 0 号测点，则异常体的特征点就在此处。延拓高度为o m 时特征点的高频 干扰严重。 2 ) 由圈3 3 1 可知，反演埋深的值都趋于一个定值时就是磁性体的真实埋深，所以 磁性体的埋深应该约为2 0 0 m ，这与真实值2 0 0 m 相符。2 0 m l o o m 的反演埋深平均值的精 度为03 3 ，符合精度5 要求。 3 ) 由图33 2 可知，构造指数在撮大值处等于2 0 1 5 ，则可以推断是水平圆柱体。 ( 4 ) 水平圆柱体模型半径为l o o m ，测点到中心距离为3 0 0 r a ，水平圆柱体中心位 置在第5 0 号测点，延拓不同高度结果如下： 咖湖瑚垂|枷啪卿m o o 第j 章计算方法及程序设计 圈3 3 3 向上延拓不同高度的1 肚。- k ，) 图3 3 4 延拓埋深和反演埋深图3 舶构造指数( 廷拓l o o m ) 1 ) 由图3 3 3 可知，延拓高度从o m 到l o o m 间隔为2 0 m 。延拓不同高度1 ，( 七。一k 。) 的极小值都在5 0 号测点，则异常体的特征点就在此处。延拓高度为0 i r l 时特征点的高频 干扰j 咂重。延拓高度为2 0 m 时，不能有效地消除高频干扰。 2 ) 由图33 4 可知，反演埋深的值都趋于一个定值时就是磁性体的真实埋深，所以 磁性体的埋深应该约为3 0 0 m ，这与真实值3 0 0 m 相符。2 0 ml o o m 的反演埋深平均值的精 度为09 7 ，符合精度5 要求。 3 ) 由图3 3 5 可知，延拓1 0 0 m 的构造指数在极大值处等于20 4 3 ，则可以推断是水 平圆柱体。 ( 5 ) 水平圆柱体模型半径为l o o m ，测点到中心距离为4 0 0 m ，水平圆柱体中心位 置在第5 0 号测点，延拓不同高度结果如下： 二 一 一。l 长安 学硕士论文 1 0 0 0 嚣 舻 6 7 0 0 0 0n 黧 2 1 。0 。0 * a a i a6 i 02 04 06 08 01 0 0 图3 3 6 向上廷拓不同高度的1 肚。一t ，) 延拓高度 0 一：l 6 0 i 一8 0 i 0 0 图3 3 7 延拓埋深和反演埋深图3 3 8 构造指数 1 ) 由图3 , 3 6 可知，延拓高度从o m 到1 0 0 m 间隔为2 0 m 。延拓不同高度1 ，罅。一女。) 的极小值都在5 0 号测点，则异常体的特征点就在此处。延拓高度为o r a 时特征点的高频 干扰严重。延拓高度为2 0 m 时，不能有效地消除高频干扰。 2 ) 由罔3 3 7 可知。反演埋深的值都趋于一个定值时就是磁性体的真实埋深，所以 磁性体的埋深应该约为4 0 0 m ，这与真实值4 0 0 m 相符。2 0 m - l o o m 的反演埋深平均值的精 度为23 ，符合精度5 要求。 3 ) 由图3 3 8u r 知，延拓l o o m 的构造指数在极大值处等于2 0 9 9 ，则可以推断是水 平圆柱体。 ( 6 ) 水平圆柱体模型半径为l o o m ，测点到中心距离为5 0 0 r a ，水平圆柱体中心位 置在第5 0 号测点，延拓不同高度结果如下： 第! 章计算方法厦程序设* 图3 3 9 向上延拓不同高度的1 h 。一k ，) 延拓高度j 0 2 0 4 0 6 0 二到 ，、 、 图3 4 0 延拓埋深和反请埋深图3 a 1 构造指数 1 ) 山图3 3 9 可知，延拓高度从0 m 到1 0 0 m 间隔为2 0 m 。延拓不同高度1 (