磁法勘探是地质构造研究和矿产资源勘探的主要地球物理勘探方法之一.利用磁异常资料,进行地质体的识别推断,如岩体范围、断裂的位置等场源边界的识别,是进一步深化研究的基础.但由于受到斜磁化、场源埋深、异常叠加以及噪声等多种因素的影响,基于磁异常资料直接进行场源的边界识别往往存在很大的误差和不确定性,所以,人们研究了各种场源边界特征增强的方法技术,以便于进行场源边界识别工作的开展.这里需要明确一下,边界识别与边界特征增强两者在概念上还是有明显差异的.前者常指利用一系列的点或线段明确标定如断裂构造线或地质体的水平投影线位置的表达方式,而后者常指使用某种数值转换方法技术对实测磁异常数据进行处理转换,得到其转换结果数据,目的是希望在转换数据中能突显场源边界附近的数值特征,便于进一步用于场源边界的识别.边界特征增强的方法已经有很多,其理论基础是在地质体边界位置的两侧,由于岩石磁性差异明显,通常其磁异常的变化率也较大,而在地质体边界内部或没有地质体的围岩区域,岩性不变,其对应的磁异常也显得平缓.故而,目前几乎所有的边界特征增强类方法都是利用这一特点进行方法技术设计的,目标是使经过其处理后的磁异常结果,在诸如量值等属性在边界位置出现特征值(例如:极大值或零值等),从而将其作为判别边界位置的依据,进行后续的场源边界识别.因此,场源边界的磁异常特征增强既是磁异常处理的研究热点之一,也是磁异常场源边界识别研究工作的基础.
基于重磁异常及其变换量进行场源边界特征增强的思想出现较早,方法种类也很多,大致可以归为两大类:数理统计类与导数分析类.数理统计类方法主要有归一化标准偏差与小子域滤波两种.归一化标准偏差最初由Cooper和Cowan(2008)提出,该方法计算一个滑动窗口内,垂向一阶导数的标准偏差与三个坐标方向一阶导数标准偏差之和的比值,将该比值记为滑动窗口中心点的转换量值,并用极大值位置标定场源边界.与归一化标准偏差有相似之处,小子域滤波同样利用标准偏差的概念,它最初由杨高印(1995)提出,基本原理是将一个滑动窗口划分成八个不同的小子域,分别计算每个小子域内重磁异常的标准偏差与平均值,然后将标准偏差最小的小子域异常平均值记为滑动窗口中心点的量值.经过该方法处理的异常梯级带更加突出,从而有利于场源边界的标定.后来,小子域滤波经过马涛等(2007)、张凤旭等(2007)和段晓旭(2014)等人的改进,应用效果得到了进一步的改善,但是处理结果存在折线走样失真等问题依然未能得到很好的解决.数理统计类方法的优点主要是可以通过改变窗口的大小来调节对噪声的压制作用,大窗口对噪声压制作用强,但对边界的分辨力较低,小窗口对噪声的压制作用弱,但分辨率较高.
相对于主要运用统计学规律的数理统计类方法,导数分析类方法则主要着眼于异常数据导数的变化特点,而且一直是边界特征增强方法的主要研究方向,方法种类也更多.例如,1936年Evjen(1936)就提出了垂向导数可以突显重力异常的局部特征.Peters(1949)提出了对磁异常向上或向下延拓、求取垂向二阶导数等方法,分离出叠加在低频特征基本场上的异常场,可用于浅部地质构造单元的划分.1965年,Hood和McClure(1965)提出利用磁异常垂向分量的垂向一阶导数与垂向二阶导数的零值位置标定铅垂台阶的水平投影边界,这也就是利用转换异常的特征值(此处为零值)直接识别出场源边界了;同年,Bhattacharyya(1965)提出了频率域计算磁异常垂向导数、延拓、化极以及化极磁异常垂向二阶导数的方法,并用化极磁异常的垂向二阶导数零值位置确定棱柱体的边界,体现了化极的作用.另外还讨论了地磁场方向和磁化方向不一致对化极结果的影响.1979年,Cordell(1979)提出运用重力异常的总水平导数极大值位置确定地质体的边界位置.Cordell和Grauch(1985)将总水平导数法应用到磁异常数据的解释中,提出先将磁异常转换为磁源重力异常,后对其求取总水平导数,并利用极大值位置确定具有垂直侧面的地质体边界.Blakely与Simpson(1986)的实验表明,若将总水平导数方法与最大值自动提取技术结合,可更好地显示场源边界.余钦范和楼海(1994)在Cordell等人的研究基础上,探讨了总水平导数极大值确定地质体边界时存在的问题,认为对单个垂直边界的场源体,总水平导数可以较准确地确定场源体边界,但当场源埋深较大、侧面倾斜或场源横向较薄时,确定的边界就会存在较大偏差.针对总水平导数法分辨率不足的问题,Fedi和Florio(2001)提出了增强总水平导数法,即计算不同阶次垂向导数加权结果的总水平导数,认为该方法通过改变垂向导数的阶次及其相应的权重,可以提高横向分辨率;还有,潘玉(2005)提出了线性构造信号法;Wang等(2009)提出了归一化总水平导数垂向导数法,等等.这些方法在一定程度都提高了横向分辨率.但是,事情往往都有两面性,有些方法也带来了新的副作用,例如明显放大了噪声干扰,这是我们需要注意的问题.
垂向导数法与总水平导数法是两种得到广泛应用的传统的边界特征增强方法,此外,Nabighian(1972,1984)与Roest等(1992)等人提出的解析信号模法也是一种基本方法,可利用其极大值作为识别场源边界的标志.二维情形下解析信号模不受磁化方向的影响,三维情形下解析信号模受磁化方向的影响,只是受影响的程度较低,这是该方法最主要的特征也是优点(管志宁和姚长利,1997;王万银等,2010).但是,解析信号模的横向分辨率尚显不足,为此人们对解析信号模做了多种增强处理.如Hsu等(1996)提出了增强解析信号法,即利用磁异常垂向n阶导数的解析信号模识别场源边界,并认为该方法比传统的解析信号模法分辨率更高;Debeglia和Corpel(1997)提出了利用磁异常垂向n阶导数或水平n阶导数代替原来的磁异常,从而构建一种新的解析信号模法;Bournas和Baker(2001)提出了解析信号模总水平导数的概念,并用于场源边界的识别,在与解析信号模法、导数解析信号模法等方法的模型实验对比后,认为该方法的分辨率更高,识别的边界位置更准确.针对解析信号模法进行场源边界识别时存在的问题,王万银(2012)研究了解析信号模极值位置空间变化规律,认为解析信号模可用来识别浅部地质体的边界,但不能用来识别深部地质体的边界,当场源埋深较大时,解析信号模的极大值对应形体的“中心位置”.
当存在多个埋深不同的场源体时,垂向导数法、总水平导数法与解析信号模法等方法很难探测出深部场源体的边界,因为深部场源的磁异常往往被浅部场源的磁异常掩盖.为此,Miller和Singh(1994)首次提出了倾斜角法的概念,认为该方法的处理结果对于场源的埋深不敏感,因此能同时识别不同埋深场源体的边界,但是受磁化方向影响较大.在倾斜角法的基础上,Verduzco等(2004)提出了倾斜角总水平导数法,认为该方法不受地质体边界倾角的影响,且受场源磁化方向影响较小.
在上述的垂向导数法、总水平导数法、解析信号模法与倾斜角法的基础上,先后衍生出了很多新的方法.如秦葆瑚(1998)与Cooper(2014)提出的解析信号模倾斜角法,该方法用解析信号模代替磁异常进行倾斜角值计算;Wijns等(2005)提出的θ图法;王万银等(2010)提出的解析信号模垂向导数法;Cooper(2009)提出的平衡解析信号模法与剖面曲率法;Alamdar等(2012)提出的倾斜角解析信号模法;Ma(2013)(Ma and Li,2012)提出的归一化总水平导数法与改进的局部相位法.这些方法均促进了边界特征增强方法的发展,对某些情形也有不错的应用效果,但是要作为适用性强、应用效果稳定的方法,还有一些不足之处.
通过上述对前人研究工作的分析介绍,我们知道边界特征增强方法众多,并且各种方法的应用条件与处理效果也差异明显.实际上,一些方法的应用效果还有待深化检验,不仅如此,不同方法之间也缺少应用条件和效果的对比,主要存在验证方法效果的理论模型过于简单,未全面细致地分析出它们的优缺点等问题.为此,为了深入了解这些方法技术的优缺点,便于实际资料处理解释时的方法选取具有针对性,本文选取了目前15种具有代表性的方法,通过针对性设计的多种类型的模型实验,研究分析它们的增强识别效果,并尝试归纳总结出它们的应用条件和优缺点,目的是为基于磁异常的边界识别方法的实际应用,选择合理的方法技术打下基础.
为了便于进行模型计算的比较,我们先对这些方法的基本原理、主要特点进行必要的分析简介.
2 部分边界特征增强方法简介 2.1 垂向导数法(VDR)、总水平导数法(THD)、解析信号模方法(AS)垂向导数法是一种传统的基本方法,它求取磁异常的垂向导数,并利用其零值位置作为确定场源边界的标志.该方法简单直接,最初的处理对象是磁异常垂直分量Za(Hood and McClure,1965).由于磁异常通常包含斜磁化的影响,用其定位的边界位置往往与真实位置有不同程度的偏差.另外,随着场源埋深增大,磁异常变得宽缓,其垂向导数零值点也逐渐外扩,偏离真实边界.
除此之外,总水平导数法也是一种很常用的基本方法,它是对水平导数的一种扩充,其表达式为
(1) |
其中THD表示计算点的磁异常总水平导数,ΔT表示磁异常.场源边界位置的磁异常变化剧烈,其水平导数数值较大,总水平导数法便充分利用这一特点,先分别计算磁异常沿x方向和y方向的导数,后求取水平方向导数向量的模,并用极大值位置标定场源边界,从而实现对边界的显著增强.
解析信号模法(AS)求取磁异常沿三个方向的总导数模,亦称总导数法(Nabighian,1984;Roest et al.,1992),它是垂直导数和总水平导数的综合,其具体表达式为
(2) |
通常利用AS的极大值位置确定场源边界.主要特点是在二维情形下边界增强效果不受磁化方向的影响,在三维情形下受影响相对较小(管志宁和姚长利,1997).其最大的特点就是在斜磁化情形下也能较准确定位并增强场源边界,但其也有缺点,主要表现为横向分辨率不足.为此,该方法有多种方式的扩展,如解析信号模垂向导数法(ASz)等.
上述方法都是最基本的导数分析类方法,是磁异常处理中常用的突出异常变化特征的一类方法.
2.2 归一化总水平导数法(NTHD)、倾斜角法(Tilt)上述导数类方法具有目标明确,计算简单的特点,但存在的主要问题是对弱异常不敏感,对深部地质体的异常信息缺少足够的反映.为此,人们提出一些针对性的改进,具有代表性的是归一化总水平导数法和斜倾角法等.
归一化总水平导数法(Ma and Li,2012)来源于总水平导数法,其表达式为
(3) |
其中THD表示计算点的磁异常总水平导数,THDWmax为包含计算点的移动窗口内各点总水平导数中的最大值.该方法针对总水平导数法对大埋深场源的边界增强不足的特点,提出利用移动窗口内总水平导数的最大值归一化计算点的总水平导数,并利用极大值标定场源边界,以期增强弱异常.我们的实验结果表明该方法往往对边界的增强作用不足.
倾斜角法(Miller and Singh,1994)是近年来应用比较广泛的信息增强方法,其表达式为
(4) |
其中ΔT表示磁异常.该方法设计思路巧妙,磁异常的垂直导数与水平导数本身就对异常的变化敏感,但其随深度增加而衰减的速度比异常本身更快,因而深部信号显得比较弱.为了克服导数随深度加大而迅速衰减的问题,Miller和Singh(1994)提出将垂直导数与水平导数的比值作为新的转换量,以弥补这个衰减问题,并认为该方法受场源埋深影响较小,可实现在场源内部为正值,边界附近为零值,场源外部为负值,从而识别深部和浅部场源体的边界.通过实验,我们验证了该方法具有受场源埋深影响较小,可同时增强异常幅值不同的场源边界的特点,另外发现它受磁化方向影响严重,定位的边界位置误差较大.
2.3 解析信号模与倾斜角法的组合类方法为了保持优点并克服各自的缺点,人们尝试将倾斜角法应用于解析信号模,这里列举主要几种组合方法.
解析信号模倾斜角法(Cooper,2014)是直接将倾斜角法中的磁异常ΔT用解析信号模AS替换,试图吸收两者的优点,进一步改善效果.其表达式为
(5) |
其中AS表示磁异常的解析信号模.该方法先用(2)式计算磁异常的解析信号模AS,后求取解析信号模AS的倾斜角值,通常利用它的极大值位置标定场源体边界.经过模型实验,我们认为该方法是优点较突出、应用效果较好的方法.特别是在斜磁化情况下,AS比ΔT受磁化方向影响更小的优点得到体现.
当然,Cooper(2014)认为解析信号模倾斜角法隐含二阶导数,对噪声较为敏感,为此,进一步提出了只包含一阶导数的“零阶解析信号模”倾斜角法TAS0(Cooper,2014).该方法使用基于希尔伯特变换定义的解析信号AS0替代原来的解析信号AS,即
(6) |
使得TAS0只包含一阶导数,从而降低噪声干扰.我们经过模型实验,验证了该方法受噪声干扰较小的特点,但同时发现它的处理结果存在较多旁瓣现象,已严重影响真实边界的辨别,没有达到想象中的理想效果.
倾斜角解析信号模法(Alamdar et al.,2012)则是直接计算倾斜角法结果值的一阶导数模量的方法,其表达式为
(7) |
表达式中的T表示磁异常的倾斜角值.该方法意图整合倾斜角与解析信号模的优点,并直接利用极大值位置标定场源边界.经过模型对比实验,我们认为该方法的效果并不理想,通常对边界的增强作用不足且对噪声极其敏感,实用效果有限.
倾斜角总水平导数法(Verduzco et al.,2004)则与倾斜角解析信号模法非常相似,其表达式为
(8) |
其中T表示磁异常的倾斜角值.显然,该方法是联合倾斜角法和总水平导数法,即先计算磁异常的倾斜角值,再求取倾斜角值的总水平导数.Verduzco等(2004)认为该方法可以更准确定位浅部地质体的边界,且受地质体倾斜侧面的影响相对较小,但我们认为其效果将与倾斜角解析信号模法类似,实用效果有限.经过模型实验,确实表明该方法通常对边界的增强作用不明显、数值计算不稳定、对噪声极其敏感.
2.4 θ图法(θ)与二阶θ图法(STM)θ图法(Wijns et al.,2005)也是近年来应用广泛的边界增强方法,其设计思路与倾斜角法有相似之处,该方法的具体计算表达式为
(9) |
其中ΔT是磁异常.该方法形式上利用解析信号模归一化总水平导数,本质上是求取解析信号向量与水平面的夹角余弦,其目的与倾斜角法相同,即增强弱异常、降低场源埋深的影响.θ图的数值在0到1之间,它利用极大值标定场源体的边界位置.我们通过模型实验,验证了该方法具有受场源埋深影响较小的特点,但同时发现它受磁化方向影响严重、定位的边界位置不准确、假边界较多、对噪声敏感.
针对θ图法的分辨率稍显不足的问题,Ma和Li(2012)认为其原因是θ图公式中包含的导数阶次较低,为此提出了二阶θ图法STM,即使用磁异常垂向导数代替原来的磁异常,计算θ图值,并利用极大值位置确定边界的位置.我们通过模型实验,验证了二阶θ图法具有分辨率较高,受场源埋深影响较小的特点,但不出意外,同样发现它受磁化方向影响严重、定位的边界位置往往不准确、对噪声极其敏感等问题.
2.5 归一化标准偏差法(NSTD)归一化标准偏差法(Cooper and Cowan,2008)是为数不多的数理统计类方法,受到广泛关注,其表达式为
(10) |
其中ΔT表示磁异常,σ表示一个大小为m×n的移动窗口内,相关量的标准偏差.该方法计算移动窗口内磁异常沿三个方向导数标准偏差的比率,并作为窗口中心点的数值,其本质是利用场源边界位置磁异常变化剧烈,标准偏差较大的特征来定位并增强场源边界.我们通过模型实验,验证了前人提出的该方法具有受场源埋深和磁化方向的影响较小,可增强弱异常的特点,另外发现它也存在不少缺点,例如存在双峰现象与假边界较多,这将严重干扰真实边界的辨别.
2.6 改进的局部相位法(ILP1与ILP2)改进的局部相位法(Ma,2013)意图联合水平一阶导数与水平二阶导数,不使用垂向导数,以提高方法的计算稳定性并同时增强不同幅值大小的异常.我们经过模型实验,认为此两种方法的明显缺点表现为它们对边界的增强不足,且受磁化方向影响严重,无法通过处理结果辨别边界位置.
虽然还有其他方法,但基本上可以认为它们与上述15种方法中的某种技术类似,效果接近,所以,限于篇幅就不再赘述了.
针对上述方法技术的效果对比,仅仅基于理论分析还难以直观感受各自的效果,为此,我们重点通过理论模型的计算,来对其不同效果加以分析评价.
3 模型实验与分析针对磁异常梯级带增强方法应用效果的检验,模型实验是最有效的手段,其具有操作简单、易于观察、可比性强,可单独测试某项影响因素的特点,因而广泛应用于地球物理勘探方法技术的研究中.为了了解各种方法的边界特征增强效果及其影响因素,下面针对性地运用多种模型,测试各种方法在不同磁化方向、不同场源埋深、不同场源形态、异常叠加及噪声干扰等因素影响下的应用效果,以更好地辨明各种方法的应用特点及优缺点.
3.1 单个垂直磁化立方体模型的实验垂直磁化情况下磁性体产生的磁异常通常具有最简单的特征,与场源的对应关系也最直接,各种边界特征增强方法只有对垂直磁化磁异常表现出优异的增强效果,才有可能在斜磁化磁异常上有良好的表现.为此设计下面的模型A,用以测试各种方法所能发挥的最佳边界特征增强能力.模型A为直立立方体,其东西宽100 m,南北长140 m,上顶埋深30 m,下底埋深60 m.磁化强度M=6 A·m-1,网格间距2 m.地磁场方向和磁化方向均垂直向下.
图 1为模型A的磁异常与边界位置叠合图.可以看出,由于模型单一,没有干扰,垂直磁化的磁异常与场源范围及边界已经具有相当好的对应关系了.这里只是基于此简单模型例子,来测试各增强方法在这样“有利”的条件下,是否能表现出良好的增强效果.如果方法的效果不好,就更难以通过复杂模型的测试,也肯定不能应用于更复杂的实际资料的处理了.
图 2为我们选择的十五种边界特征增强方法对模型A磁异常的处理结果.直观来看,解析信号模(AS)、解析信号模垂向导数(ASz)、解析信号模倾斜角(TAS)、垂向导数(VDR)、总水平导数(THD)和倾斜角(Tilt)的增强效果较好,边界的菱角清晰、位置大体正确、幅值均匀,干扰真实边界识别的假边界没有或极少.详细来讲,AS的增强作用明显,无假边界,但分辨率稍显不足;ASz的效果与AS相似,但分辨率提高更明显;TAS处理结果的正常场区域低值与边界高值过渡平缓,稍显边界位置模糊;零阶解析信号模倾斜角(TAS0)的处理结果出现沿着真实边界向外顺势延伸的假边界,即旁瓣较多;改进的局部相位(ILP1和ILP2)、归一化标准偏差(NSTD)、θ图、二阶θ图(STM)与倾斜角总水平导数(THDR)的效果相近,突显了边界位置的数值特征,但存在双峰现象,表现为在边界外围出现高值圈闭;θ图、STM与THDR处理结果的边界峰值较宽,且有外扩迹象,难以准确定位边界位置;倾斜角解析信号模(ASTA)与归一化总水平导数(NTHD)的效果相似,增强作用不足,无法分辨边界位置.总之,对垂直磁化模型,AS、ASz、TAS、VDR、THD与Tilt的图像相似,效果也较好,而其他的方法大多都存在不同程度的双峰现象,影响真实边界的辨别.
由于垂直磁化情况只在地磁极附近才会出现,而绝大多数情况下磁性岩矿体的磁化方向都是倾斜的,即斜磁化才具有普遍性,因此,各种边界特征增强方法只有在斜磁化模型上表现出良好的增强效果,才可能具备实用价值.为此我们这里设计的模型B,几何参数与模型A相同,但磁化方向不同.具体为:地磁场倾角I0=37°、偏角D0=5°,磁化方向倾角I=30°、偏角D=10°,模型磁化强度M=6 A·m-1,网格间距2 m.
图 3为模型B的磁异常与边界位置.与垂直磁化磁异常比较,磁异常与场源边界位置不再具有简单的对应关系了,可以看出磁化方向对磁异常的巨大影响.
图 4为所选十五种边界特征增强方法对模型B磁异常的处理结果.解析信号模类方法受磁化方向影响较小,故而AS与ASz可较准确定位边界位置,且增强效果明显,但增强的边界不完整,例如此处只有东西方向的边界得到增强;ASTA也是解析信号模类方法,但是对边界的增强不足、数值计算不稳定;NSTD受斜磁化影响较小,可相对完整地增强所有边界,但边界位置稍有偏差,且在正常场区域存在很多假边界;TAS的增强作用明显,斜磁化情形下也可识别出较完整的边界,且边界位置准确、假边界极少,虽然横向分辨率稍低;TAS0处理结果的旁瓣较多,严重干扰真实边界的辨别;THD的增强作用明显,但斜磁化情形下存在双峰现象;VDR受斜磁化影响严重,指示边界的零值位置与真实边界位置偏差非常大;Tilt、θ图和STM的效果与VDR相似,受斜磁化影响严重,定位的边界位置误差很大;THDR受磁化方向影响较小,可相对完整地增强所有边界,但存在双峰现象;ILP1、ILP2与NTHD三者都没有明显增强边界位置的数值特征,故而无法通过其处理结果辨别边界位置.总之,斜磁化通常引发各种方法对边界位置的定位不准确,并导致产生假边界.除了AS、ASz与TAS对边界位置的定位较准确之外,其他方法均存在不同程度的缺陷,如定位的边界位置有偏差、假边界较多等.
经过上面的垂直磁化模型与斜磁化模型的实验,我们认为对解析信号模(AS)、解析信号模垂向导数(ASz)、归一化标准偏差(NSTD)、解析信号模倾斜角(TAS)、θ图、垂向导数(VDR)、总水平导数(THD)、倾斜角(Tilt)与倾斜角总水平导数(THDR)九种方法有继续探讨其应用效果的必要,而对于倾斜角解析信号模(ASTA)、改进的局部相位(ILP1与ILP2)、归一化总水平导数(NTHD)、二阶θ图(STM)与零阶解析信号模倾斜角(TAS0)六种方法,由于它们的缺陷过于直接与明显,我们不建议在实际应用中使用,因此,在接下来的模型实验中不再对它们进行成图与分析.
3.3 不同埋深场源体的组合模型实验由于随着磁性体埋深的增加,其磁异常迅速衰减,这就极大地增加了深部场源边界识别的困难,也是对边界特征增强方法的挑战.由于这些增强方法对埋深不同的场源体会有不同的响应,有些方法受场源埋深影响较大,识别结果往往受到浅源大幅值异常的控制,而有些方法则受场源埋深影响相对较小,比较而言,可同时增强深部和浅部场源体的边界特征.为此我们设计组合模型C,以测试各种方法对不同埋深场源体边界的增强效果.模型C包括三个大小相同,但上顶面埋深不同的长方体.长方体东西宽60 m,南北长140 m,上下表面间距40 m.地磁场方向倾角I0=35°、偏角D0=5°,网格间距2 m,磁化强度M=6 A·m-1.其他参数见表 1.
图 5为模型C的磁异常与边界位置.可以看出,由于埋深不同,三个长方体场源的磁异常幅值各异,左边的幅值最小,右边的幅值最大.另外,由于是倾斜磁化,图中异常范围与场源边界位置没有直接的对应关系,这些因素都给边界特征增强加大了难度.
剔除六种有明显不足的方法后,剩余九种边界特征增强方法对模型C磁异常的处理结果见图 6.观察图 6可以发现:AS、ASz和THD三者效果相似,受场源埋深影响较大,对浅部场源边界的增强明显,而对于深部场源边界的增强不足,表现为异常幅值较低;θ图与Tilt受场源埋深影响较小,可同时增强幅值大小不同的异常,但识别的边界位置误差较大,假边界信息很多;TAS可同时增强深部和浅部场源边界,且突显的边界位置准确、完整度高,假边界很少;NSTD也可同时增强埋深不同的场源体边界,且边界较完整、位置相对准确,但在正常场区域存在很多“假极大值”,即假边界;VDR随场源埋深增大,其幅值下降速度比THD等方法稍慢,但零值位置与边界位置偏差很大;THDR受场源埋深影响较小,但结果量值较小、数值计算不稳定.综上分析,通常公式中使用“归一化”的方法(如NSTD、TAS、θ图、Tilt与THDR),其结果受场源埋深影响较小,可同时增强埋深不同的场源边界,而没有使用“归一化”的方法(如AS、ASz、VDR与THD)则受场源埋深影响较大,表现为其结果量值随埋深增大而迅速下降.
对于空间位置上下叠加的组合场源体,它们的磁异常也相互叠加,表现为浅源异常幅值较大且与场源位置对应较好,其边界易于增强,而深源异常受制于浅源异常的掩盖,其边界往往难以定位及增强.为此设计模型D,以测试各种方法对上下叠加场源体边界的增强效果.模型D包含两个长方体,上部长方体尺寸小于下部长方体,另外它们的磁化强度大小和方向均不同.地磁场方向倾角I0=30°、偏角D0=10°,网格间距2 m.其他参数见表 2.
图 7为模型D的侧视图,图 8为模型D的磁异常与边界位置,图 9为九种边界特征增强方法对模型D磁异常的处理结果.观察图 9可以发现,AS、ASz与THD对叠加异常的分辨能力相近,对上部场源边界的增强明显,幅值大、位置准确,而对下部场源边界的增强不足,幅值较低,难以分清边界;θ图与Tilt的增强效果相似,分辨力相当,对上部与下部场源边界都可同时增强,但定位的边界位置误差很大、假边界很多;TAS对叠加异常的分辨能力较强,可同时增强上部和下部场源边界,且突显的边界清晰、完整、位置正确,假边界较少;NSTD受叠加异常影响较小,可同时增强上部和下部场源边界,但存在双峰现象、定位的边界位置有偏差;VDR对叠加异常的分辨能力较弱,外加斜磁化情形下定位的边界位置误差很大,总体效果不好;THDR的数值计算不稳定,无法分辨边界位置.综上分析,TAS对叠加异常的分辨力较强,在增强上部场源边界的同时也可较好地增强下部场源边界,而其他方法的分辨力相对较弱,表现为对下部场源边界的增强不足,如AS、ASz与THD,或定位的边界位置不准确,如θ图、Tilt、NSTD与VDR.
实际地质状况复杂多变,地质体可能是上下叠加、埋深不同、间距变化、形态各异、磁化强度大小与方向均不同的.为了模拟实际地质状况,设计组合模型E,以测试各种方法的综合应用效果.模型E包含三个几何形态不同的磁性体:三棱柱体、圆柱体、四棱柱体.三个形体的顶底面间距相同,为100 m,但上顶面埋深、水平间距、磁化强度大小和方向均不同,且有横向有纵向排列.地磁场方向倾角I0=40°、偏角D0=10°,网格间距2 m.考虑了不同磁化即剩磁明显存在的影响,具体参数见表 3.
图 10为模型E的磁异常与边界位置.可以看出,由于受到磁化强度和场源埋深的影响,磁异常的幅值各异,表现为三角棱柱体投影位置的异常主要为正值,幅值最大,而圆柱体投影位置的异常主要为负值,还有,异常范围与场源边界位置没有直接的对应关系.另外,由于圆柱体与四棱柱体的水平间距较小,两者产生的磁异常存在明显的叠加现象.
图 11为九种边界特征增强方法对模型E磁异常的处理结果.从图 11可以看出,AS、ASz与THD的图像效果相似,对形态、水平间距不同的场源边界均可准确定位并增强,当然,受埋深增大的影响,四棱柱体的边界较为模糊;Tilt的效果比θ图稍好,受斜磁化影响,总体上两者的结果图像都较为混乱,表现为边界不清晰、假边界很多;TAS受场源形态、埋深和水平间距的干扰较小,显示的边界清晰、完整、位置准确,假边界很少,是一种较好的方法;NSTD可增强形态不同的场源边界,且受场源埋深影响较小、显示的边界较完整,但在正常场区域存在较多干扰真实边界辨别的假边界;VDR受斜磁化影响严重,定位的边界位置误差很大,效果不好;THDR受斜磁化影响较小,但对边界的增强不足、结果量值很小、数值计算不稳定.综上分析,通常场源形态、间距对增强效果的影响相对较小,场源埋深、磁化方向才是主要的影响因素.对无噪声数据,TAS对多种影响因素的抗干扰能力较强,适应性较好,而其他方法均不同程度受到场源埋深与磁化方向的影响,定位的边界位置不准确或假边界比较多.具体来讲,AS、ASz与THD受场源埋深影响较大;θ图、VDR与Tilt受斜磁化干扰较强;NSTD在正常场区域出现很多假边界;THDR数值计算不稳定.
实测磁异常总是包含各种干扰的,虽然通常通过滤波处理,如补偿圆滑滤波(侯重初,1979,1981)等方法可以消除部分干扰,但是干扰是不可能完全消除的,它的存在将会干扰场源边界的增强处理.此处使用模型A的磁异常数据(图 1),对其添加5%的随机噪声后进行“补偿圆滑滤波”处理,得到包含有少量噪声的磁异常,并用以测试各种边界特征增强方法对噪声的抗干扰能力.
图 12为对模型A磁异常(图 1)添加5%高斯噪声后的数据,图 13为对图 12磁异常进行补偿圆滑滤波后的结果.观察图 12和图 13可发现,经过补偿圆滑滤波处理,图 13磁异常包含的噪声已明显少于图 12,这说明该消噪处理的作用显著,将有益于后续的边界增强处理.
图 14为九种边界特征增强方法对图 13磁异常的处理结果.分析图 14可见,AS、VDR与THD只包含一阶导数,受噪声干扰较小,正常场区域高频噪声幅值较低,与边界的高幅值区别明显;ASz与TAS包含二阶导数,对噪声较敏感,表现为正常场区域存在较多高频干扰,边界位置出现不规则弯曲与偏移,但总体上噪声的影响可控,正常场区域与边界的区别依然明显;NSTD只包含一阶导数,但由于该方法计算窗口标准偏差,对噪声更为敏感,正常场区域存在较多高频干扰;θ图与Tilt对噪声非常敏感,正常场区域存在较多高频干扰,严重影响真实边界的辨别;THDR受噪声干扰也非常严重,边界没有得到突显,无法辨别.
影响方法对噪声敏感性的原因很多,包含的导数阶次与方法公式中基本项的数学关系(如使用“归一化”)是主要来源.通常高阶导数可以提高方法的横向分辨率,但也会放大噪声的干扰,如ASz、TAS与THDR包含二阶导数,横向分辨率较高对噪声也更敏感,而AS、VDR与THD只包含一阶导数,横向分辨率相对较低对噪声也不敏感.但是,观察ASz与TAS的图像可发现:虽然两者对噪声较敏感,但是由于它们对边界的增强作用充足,边界依然是清晰、易于辨别的,这相当于弥补了噪声干扰的影响.
4 结论通过对主要的边界增强类方法的理论分析,并针对性地设计了从简单到复杂、分门别类的多种模型实验及后续的特征对比分析,我们对这些边界特征增强方法有以下主要认识和结论:
(1) 隐含高阶导数的方法,其横向分辨率较高,如包含二阶导数的二阶θ图(STM)、解析信号模倾斜角(TAS)与包含一阶导数的θ图、“零阶解析信号模”倾斜角(TAS0)比较,前者的横向分辨率要更高.但是,因为高阶导数的作用,往往噪声干扰也会放大.所以,实际应用时,需要事先最大程度可靠地滤除磁异常中的高频噪声.
(2) 所有边界特征增强方法都受到斜磁化的影响,其中有些方法受影响相对较小,比如AS、ASz和TAS等方法,而θ图、Tilt、垂向导数(VDR)等方法受影响较大.解析信号模类方法受磁化方向的影响均较小,如AS和ASz在磁化倾角很小的情况下也能准确识别出边界,但是也有缺点,主要体现在识别的边界往往不够完整.
(3)使用“归一化”方式的方法(如倾斜角法Tilt、θ图法与归一化标准偏差法)受场源埋深影响较小,可同时增强大幅值的强异常和小幅值的弱异常,但是在异常相对平缓的正常场区域会出现很多“假极大值”,即假边界,这将严重干扰真实边界的辨别.相反地,像AS、ASz和总水平导数(THD)等方法没有使用“归一化”的方式,其增强结果存在的假边界较少,但是对大埋深的场源边界往往增强不足.
(4)综合对比所选择的15种方法可以确定:解析信号模垂向导数ASz受磁化方向影响较小,处理的结果边界位置准确、假边界较少、分辨率较高,其不足之处是受场源埋深影响较大;总水平导数THD对边界的增强作用比ASz明显,受噪声干扰较小,其欠缺处是稍有双峰现象;解析信号模倾斜角TAS受磁化方向和场源埋深的影响均较小,突显的边界也较完整、清晰,其不足的方面是受噪声干扰稍大.即便如此,综合评价而言,我们依然认为这三种方法(ASz、THD、TAS)是现今基于磁异常的场源边界特征增强效果相对更好的方法.
Alamdar K, Rouhani A K, Ansari H. 2012. A new edge detection method based on the analytic signal of tilt angle (ASTA) for magnetic anomalies.//Istanbul International Geophysical Conference and Oil & Gas Exhibition.Istanbul, Turkey. | |
Bhattacharyya B K. 1965. Two-dimensional harmonic analysis as a tool for magnetic interpretation. Geophysics , 30 (5) : 829-857. DOI:10.1190/1.1439658 | |
Blakely R J, Simpson R W. 1986. Approximating edges of source bodies from magnetic or gravity anomalies. Geophysics , 51 (7) : 1494-1498. DOI:10.1190/1.1442197 | |
Bournas N, Baker H A. 2001. Interpretation of magnetic anomalies using the horizontal gradient analytic signal. Annali di Geofisica , 44 (3) : 505-526. | |
Cooper G R J, Cowan D R. 2008. Edge enhancement of potential-field data using normalized statistics. Geophysics , 73 (3) : H1-H4. DOI:10.1190/1.2837309 | |
Cooper G R J. 2009. Balancing images of potential-field data. Geophysics , 74 (3) : L17-L20. DOI:10.1190/1.3096615 | |
Cooper G R J. 2014. Reducing the dependence of the analytic signal amplitude of aeromagnetic data on the source vector direction. Geophysics , 79 (4) : J55-J60. DOI:10.1190/geo2013-0319.1 | |
Cordell L. 1979. Gravimetric expression of graben faulting in Santa Fe Country and the Espanola Basin, New Mexico.//New Mexico Geological Society Guidebook, 30th Annual Fall Field Conference.New Mexico:New Mexico Geological. Society,59-64. | |
Cordell L, Grauch V J S. 1985. Mapping basement magnetization zones from aeromagnetic data in the San Juan Basin, New Mexico.//Hinze W J ed. The Utility of Regional Gravity and Magnetic Anomaly Maps. Tulsa, OK:Society of Exploration Geophysics, 181-197. | |
Debeglia N, Corpel J. 1997. Automatic 3-D interpretation of potential field data using analytic signal derivatives. Geophysics , 62 (1) : 87-96. DOI:10.1190/1.1444149 | |
Duan X X. 2014. Research and improvement of small subdomain filtering of gravity and magnetic anomaly[Master thesis](in Chinese).Beijing:China University of Geoscience. | |
Evjen H M. 1936. The place of the vertical gradient in gravitational interpretations. Geophysics , 1 (1) : 127-136. DOI:10.1190/1.1437067 | |
Fedi M, Florio G. 2001. Detection of potential fields source boundaries by enhanced horizontal derivative method. Geophysical Prospecting , 49 (1) : 40-58. DOI:10.1046/j.1365-2478.2001.00235.x | |
Guan Z N, Yao C L. 1997. Inversion of the total gradient modulus of magnetic anomaly due to dipping dike. Earth Science-Journal of China University of Geosciences (in Chinese) , 22 (1) : 81-85. | |
Hood P, McClure D J. 1965. Gradient measurements in ground magnetic prospecting. Geophysics , 30 (3) : 403-410. DOI:10.1190/1.1439592 | |
Hou Z C. 1979. A filtering method in frequency domain for blanking jamming. Geophysical and Geochemical Exploration (in Chinese) , 3 (5) : 50-54. | |
Hou Z C. 1981. Filtering of smooth compensation. Geophysical Prospecting for Petroleum (2) : 22-29. | |
Hsu S K, Sibuet J C, Shyu C T. 1996. High-resolution detection of geologic boundaries from potential-field anomalies:An enhanced analytic signal technique. Geophysics , 61 (2) : 373-386. DOI:10.1190/1.1443966 | |
Ma G Q, Li L L. 2012. Edge detection in potential fields with the normalized total horizontal derivative. Computers & Geosciences , 41 : 83-87. | |
Ma G Q. 2013. Edge detection of potential field data using improved local phase filter. Exploration Geophysics , 44 (1) : 36-41. DOI:10.1071/EG12022 | |
Ma T, Wang T C, Wang Y. 2007. An improved gradient-preserved grid data filtering approach. Oil Geophysical Prospecting (in Chinese) , 42 (2) : 198-201. | |
Miller H G, Singh V. 1994. Potential field tilt-a new concept for location of potential field sources. Journal of Applied Geophysics , 32 (2-3) : 213-217. DOI:10.1016/0926-9851(94)90022-1 | |
Nabighian M N. 1972. The analytic signal of two-dimensional magnetic bodies with polygonal cross-section:Its properties and use for automated anomaly interpretation. Geophysics , 37 (3) : 507-517. DOI:10.1190/1.1440276 | |
Nabighian M N. 1984. Toward a three-dimensional automatic interpretation of potential field data via generalized Hilbert transforms:Fundamental relations. Geophysics , 49 (6) : 780-786. DOI:10.1190/1.1441706 | |
Pan Y. 2005. The theory and application of Linear Structure Signal of gravity anomaly in detecting the fault position[Master thesis] (in Chinese).Xi'an:Chang'an University. | |
Peters L J. 1949. The direct approach to magnetic interpretation and its practical application. Geophysics , 14 (3) : 290-320. DOI:10.1190/1.1437537 | |
Qin B H. 1998. A design of application system of geophysical visualization in scientific computing. Computing Techniques for Geophysical and Geochemical Exploration (in Chinese) , 20 (3) : 261-266. | |
Roest W R, Verhoef J, Pilkington M. 1992. Magnetic interpretation using the 3-D analytic signal. Geophysics , 57 (1) : 116-125. DOI:10.1190/1.1443174 | |
Verduzco B, Fairhead J D, Green C M, et al. 2004. New insights into magnetic derivatives for structural mapping. The Leading Edge , 23 (2) : 116-119. DOI:10.1190/1.1651454 | |
Wang W Y, Pan Y, Qiu Z Y. 2009. A new edge recognition technology based on the normalized vertical derivative of the total horizontal derivative for potential field data. Applied Geophysics , 6 (3) : 226-233. DOI:10.1007/s11770-009-0026-x | |
Wang W Y. 2012. Spatial variation law of the extreme value position of analytical signal amplitude for potential field data. Chinese Journal of Geophysics (in Chinese) , 55 (4) : 1288-1299. DOI:10.6038/j.issn.0001-5733.2012.04.024 | |
Wijns C, Perez C, Kowalczyk P. 2005. Theta map:Edge detection in magnetic data. Geophysics , 70 (4) : L39-L43. DOI:10.1190/1.1988184 | |
Wang W Y, Qiu Z Y, Yang Y, et al. 2010. Some advances in the edge recognition of the potential field. Progress in Geophysics (in Chinese) , 25 (1) : 196-210. DOI:10.3969/j.issn.1004-2903.2010.01.027 | |
Yang G Y. 1995. A new technique for potential-field data processing:Small subdomain filtering. Oil Geophysical Prospecting (in Chinese) , 30 (2) : 240-244. | |
Yu Q F, Lou H. 1994. Locating the boundaries of magnetic or gravity sources using horizontal gravity sources using horizontal gradient anomalies. Computing Techniques for Geophysical and Geochemical Exploration (in Chinese) , 16 (4) : 363-367. | |
Zhang F X, Zhang F Q, Liu C, et al. 2007. A technique for elaborate explanation of faulted structures:Three-directional small subdomain filtering. Chinese Journal of Geophysics (in Chinese) , 50 (5) : 1543-1550. DOI:10.3321/j.issn:0001-5733.2007.05.031 | |
段晓旭. 2014. 重磁异常小子域滤波算法的分析研究与改进[硕士论文]. 北京:中国地质大学(北京). http://cdmd.cnki.com.cn/article/cdmd-11415-1014233762.htm | |
管志宁, 姚长利. 1997. 倾斜板体磁异常总梯度模反演方法. 地球科学-中国地质大学学报 , 22 (1) : 81–85. | |
侯重初. 1979. 一种压制干扰的频率滤波方法. 物探与化探 , 3 (5) : 50–54. | |
侯重初. 1981. 补偿圆滑滤波方法. 石油物探 (2) : 22–29. | |
马涛, 王铁成, 王雨. 2007. 一种改进的网格数据保持梯度滤波方法. 石油地球物理勘探 , 42 (2) : 198–201. | |
潘玉. 2005. 重力异常线性构造信号在识别断裂上的理论及应用[硕士论文]. 西安:长安大学. http://xuewen.cnki.net/CMFD-2006030647.nh.html | |
秦葆瑚. 1998. 用解析信号计算倾角异常. 物探化探计算技术 , 20 (3) : 261–266. | |
王万银, 邱之云, 杨永, 等. 2010. 位场边缘识别方法研究进展. 地球物理学进展 , 25 (1) : 196–210. DOI:10.3969/j.issn.1004-2903.2010.01.027 | |
王万银. 2012. 位场解析信号振幅极值位置空间变化规律研究. 地球物理学报 , 55 (4) : 1288–1299. DOI:10.6038/j.issn.0001-5733.2012.04.024 | |
杨高印. 1995. 位场数据处理的一项新技术——小子域滤波法. 石油地球物理勘探 , 30 (2) : 240–244. | |
余钦范, 楼海. 1994. 水平梯度法提取重磁源边界位置. 物探化探计算技术 , 16 (4) : 363–367. | |
张凤旭, 张凤琴, 刘财, 等. 2007. 断裂构造精细解释技术——三方向小子域滤波. 地球物理学报 , 50 (5) : 1543–1550. DOI:10.3321/j.issn:0001-5733.2007.05.031 | |