第４５卷　第４期２０２３年７月
物探化探计算技术
ＣＯＭＰＵＴＩＮＧＴＥＣＨＮＩＱＵＥＳＦＯＲＧＥＯＰＨＹＳＩＣＡＬＡＮＤＧＥＯＣＨＥＭＩＣＡＬＥＸＰＬＯＲＡＴＩＯＮ
Ｖｏｌ．４５　Ｎｏ．４
Ｊｕｌ．２０２３
收稿日期：２０２２ ０４ ０２
基金项目：中国电建集团贵阳勘测设计有限公司重大专项（ＹＪＺＤ２０２０－０２）第一作者：杜兴忠（１９７３－），男，高级工程师，主要从事地球物理勘探及地震监测相关技术研究工作，Ｅ ｍａｉｌ：ｄｕｘｉｎｚｈｏｎ＿ｇｙｙ
＠ｐ
ｏｗｅｒｃｈｉｎａ．ｃｎ。文章编号：１００１ １７４９（２０２３）０４ ０５０５ ０９半航空瞬变电磁法正反演算法及
岩溶洼地实测数据验证
杜兴忠１，李永铭１，李万荣１，郭　明２，高　嵩２，李登柯２
（１．中国电建集团　贵阳勘测设计研究院有限公司，贵阳　５５００８１；２．成都理工大学　地球勘探与信息技术教育部重点实验室，成都　６１００５９
）摘　要：半航空瞬变电磁法（ｓｅｍｉ－ａｉｒｂｏｒｎｅｔｒａｎｓｉｅｎｔｅｌｅｃｔｒｏｍａｇｎｅｔｉｃ，ＳＡＴＥＭ）是一种新型电磁勘探技术，相比较于传统地面瞬变电磁法，勘探工作不受地形条件制约。特别是岩溶地貌地区，由于地形条件复杂，地面勘探工作很难开展。因此这里将ＳＡＴＥＭ应用于贵州某岩溶地貌地区的储能水电站工程探测。由于ＳＡＴＥＭ系统独特的工作方式，将会导致信号中噪声成分复杂，单点反演结果可靠性差，而过度的滤波会可能导致关键信息的丢失。据此首先根据岩溶地貌的两个经典构造模型设计了正演模拟，并在正演模拟的响应信号中加入噪声以模拟野外实测信号再反演计算，分析了使用单点正则化阻尼最小二乘反演和横向约束的阻尼最小二乘反演法刻画各模型的能力差异，理论上证明横向约束反演方法的优越性。最后将两种反演方法用于解释某储能水电站勘探的实测数据，反演结果与实际地质资料具有较好地对应，表明了ＳＡＴＥＭ用于岩溶地貌勘探的是可行的，且效果明显。
关键词：半航空瞬变电磁法（ＳＡＴＥＭ）
；岩溶地貌；横向约束中图分类号：Ｐ６３１．３  文献标志码：Ａ  犇犗犐：１０．３９６９／ｊ
．ｉｓｓｎ．１００１ １７４９．２０２３．０４．１１０　引言
岩溶地貌广泛分布在我国西南地区，多为碳酸盐岩发育地区，在溶蚀风化作用、地下暗河岩溶塌陷作用下，形成四周多低山和峰林，底部多为含水量较
高的软粘土层的地貌特征［
１－２
］。其特殊的水文地质条件是水利水电工程建设中常遇的难题，往往由于对岩溶区域的勘测计算不够精确，会造成防护工程
的巨大浪费，或是引起巨大的经济损失［３－４
］，所以在
岩溶地貌工程施工时，能准确快速地探测出浅层岩溶管道、深层岩溶管道及深大暗河等不良地质构造
有重要的意义［５－６］。半航空瞬变电磁法是一种地面
激发、
空中采集的工作方式，具有方便、高效、勘探范围大、空间分辨率高等优势，适合在山地、起伏地形、
沼泽地等地面物探难以开展的地区［７－８］
。
ＳＡＴＥＭ作为一种新兴的地球物理勘探方法，近几年在我国已经被有效地应用于矿产资源勘探、
地下水资源勘探以及隧道工程勘探等领域［９－１１］
。
然而，ＳＡＴＥＭ探测系统一直在飞行运动中连续采集数据，这使得该系统不能像地面瞬变电磁法那样经过长时间的数据叠加来提高信噪比，而且还会在
数据中引入其他类型的噪声（如低频运动噪声）［１２］。
因此该系统的较弱的信号和复杂的噪声类型导致数据处理流程复杂，受信号处理人员的主观影响较大，且处理后的数据仍可能残存噪声干扰，常规单点反
演无法保证结果的可靠性，但将横向约束方法应用于瞬变电磁数据反演中却能有效改善反演效
果［１３－１４］，故选择带横向约束的阻尼最小二乘法分别
反演ＳＡＴＥＭ的模拟数据和实测数据。
模拟数据选择两种经典的岩溶地质模型［１５
］，通
过正演模拟得到对应的半航空瞬变电磁响应特征，并且为尽可能模拟野外实测数据，考虑在模拟数据中加入常见的两种噪声（
白噪声和仪器背景噪声）。通过反演添加噪声的模拟数据，分析横向约束反演对岩溶地质模型的刻画能力，
证明了横向约束反演的有效性。实测数据选择贵州岩溶地貌地区的某水库半航空瞬变电磁法勘探试验数据，
反演结果与实际地质资料具有较好的对应，表明了ＳＡＴＥＭ用于岩溶地貌勘探效果明显。
１　正演验证
笔者所用半航空瞬变电磁探测以接地长导线源作为激发源，导线中心位于坐标原点上，沿狓轴向两侧延伸，狕轴以向下为正，半航空瞬变电磁法工作示意图见图１。
在层状模型中，线圈接收回路中心的频域的垂
直分量磁感应强度表示为［１６
］：
犎狕＝犐４π∫犔－
犔狔犚∫
∞
０（１＋狉ＴＥ）犲狌０狕λ２
狌０犑１（
λ犚）ｄλｄ狓＇（１
）式中：犐和犔分别为电流大小和导线的半长；狕表示
飞行高度的坐标值；犚＝（（狓－狓＇）２＋狔２）１／２
；犑１为一阶贝塞尔函数；狉ＴＥ为反射系数。Ｇ
ｕｐｔａｓａｒｍａ提出的Ｈａｎｋｅｌ变换用于计算频域响应［１
７］，王华军［１８］
提出的正弦变换用于将响应虚部变换到时域，垂直磁
场的时间导数的时域表达式为式（２）。犱犎狕（狋）犱狋＝－１
２
２槡
π∑∞
狀＝－∞
犐犕［犎狕（犲狀Δ正则化的约束条件
狋）］·犠（狀Δ）（２
）式中：Δ＝
ｌｎ（１０）２０
为采样间隔；犠为正弦变换滤波系数；狀＝－１４９～１００，即２５０点正弦滤波系数。
为验证本文正演算法的准确性，选择将正演模拟结果和对应解析解对比（图２）。正演参数如下：线源长度为１０００ｍ；发射电流１Ａ；测点狓、狔、狕坐标为（０，２５０，－３０）；均匀半空间电阻率为１００Ω·ｍ；在０．１ｍｓ至１００ｍｓ区间范围内等对数间隔选取３
０个时间点。
图１　半航空瞬变电磁法工作示意图Ｆｉｇ．１　Ｓｃｈｅｍａｔｉｃｄｉａｇｒａｍｏｆｓｅｍｉ－ａｉｒｂｏｒｎｅ  　ｔｒａｎｓｉｅｎｔｅｌｅｃｔｒｏｍａｇ
ｎｅｔｉｃｍｅｔ
ｈｏｄ图２　正演与解析解对比图
Ｆｉｇ．２　Ｃｏｍｐａｒｉｓｏｎｏｆｆｏｒｗａｒｄｍｏｄｅｌｉｎｇａ
ｎｄ  ａｎａｌｙ
ｔｉｃａｌｓｏｌｕｔｉｏｎｓ由图２可以看出，误差曲线显示，在前１０ｓ内正演模拟结果误差极小，仅在１００ｓ处出现了明显地震荡。
２　模拟噪声
为使模拟数据更接近真实数据，还需要在模拟数据中加入噪声干扰，笔者加入白噪声和仪器背景噪声两种。一般半航空瞬变电磁法中，用于反演的数据都是感应电动势经过等对数间隔时窗抽道后的数据，各时窗的宽度和此时窗中心点的延时成正比，且各时窗的瞬态数据都是整个时窗内所有感应电动势的平均值。当环境噪声为白噪声时，
各个时窗的数据受不同时窗宽度的影响会使白噪声与狋－１／２
成正比衰减［１９］。仪器背景噪声用统一标准偏差表示，通常设置为２％，加噪声后模拟响应公式如下［２０］：
６
０５    物探化探计算技术４５卷
图３　单点模拟响应加噪效果图
Ｆｉｇ．３　Ｔｈｅｅｆｆｅｃｔｏｆａｄｄｉｎｇｎｏｉｓｅｔｏｔｈｅｓｉｎｇｌｅ－  ｐｏｉｎｔｓｉｍｕｌａｔｉｏｎｒｅｓｐｏｎｓｅ
犞＇＝犞＋犌（０，１）·犖狅犻２犢＋
犖狅犻犅（）犞
［］２１／２·犞
（３）
犖狅犻犅＝犫·
狋
１·１０－
（）３１／２（４）
式（３）中：犞＇为加噪声后响应；犞为正演模拟响应；犌（０，１）是均值为“０”且标准差为“１”的高斯分布；犖狅犻犅为白噪声；犖狅犻犢为仪器背景噪声。
在单点模拟响应依照上述方式加１０ｎＴ／ｓ白噪声和标准偏差为２％仪器背景噪声后的含噪声响应如图３所示，其中单点正演参数如下：线源长度为１０００ｍ；发射电流１０Ａ；测点狓、狔、狕坐标为（０，５００，－３０）；双层模型阻率分别为３０Ω·ｍ和１００Ω·ｍ，其中第一层厚为６０ｍ；在０．１ｍｓ至１０ｍｓ区间范围内等对数间隔选取５０个时间点。
３　约束反演理论
半航空瞬变电磁法的数据处理和反演解释是一个线性流程，只有经过数据处理后才能反演解释，但是在噪声成分复杂的条件下，数据处理很难完全剔除噪声干扰。对含有噪声的数据进行单点反演将很难得到一个可靠的结果，所以引入了带横向约束的阻尼最小二乘反演方法约束最终的反演结果，减小噪声对反演结果的影响。阻尼最小二乘的目标函数为：
Φ＝‖犠犱（犱犻－犉犻（犿））‖２＋α‖犠１Δ犿‖２
（５）式中：Φ为总目标函数；犱为观测数据向量，其数据个数犻与抽道数量有关；犿为模型参数向量；
α为阻尼因子；Δ犿为模型修正量；犠
犱
为观测数据加权矩
阵，该矩阵为一个对角矩阵；犠
１
为模型修正量数据加权矩阵。
设犉（犿）为模型犿计算出的理论电磁响应，为
使计算简便，选择在初始模型犕
０
处对犉（犿）做一阶泰勒级数展开并忽略高阶项得式（６）。
Δ犱＝犱－犉（犿０）＝犌Δ犿＋犲ｏｂｓ（６）式中：Δ犱为实测数据与在模型参数犿条件下计算
出理论响应的差值；Δ犿＝犿－犿
０
是模型修正量。
横向约束的中心思想是减小相邻测点反演结果中模型参数的差异，使得反演结果中电阻率和层厚度缓慢变化，根据［２１］，可以设：
犚狆犿－犲狉狆＝０（７）
式中：犲
ｒｐ
是相邻点之间电阻率和层厚度的差异，犚
狆是由１和－１组成的横向约束矩阵，该矩阵为一个稀疏矩阵，具体可以表示为：
犚狆＝
１０…０－１０…０００
０１０…０－１０…００
０００…０１０…０－
烄
烆
烌
烎
１犕×犖
（８）
其中：犕＝（犖
狆
－１） （２狀－１）；犖＝犖狆 （２狀－１）。
对式（７）左右两边减去犚
狆
犿０可得：
犚狆Δ犿＝－犚狆犿０＋犲ｒｐ（９）
在实际反演过程中可以根据已知的地质信息，对不同的地层采用不同的横向约束强度，此时可通过引入加权因子来调整各个模型参数光滑度，将式（９）左右两边乘以加权系数矩阵得：
犠狆犚狆Δ犿＝犠狆犚狆犿０＋犠狆犲狉狆（１０）
我们将带有阻尼因子约束的数据拟合方程和模型参数横向约束方程集合起来，可以得到带横向约束的阻尼最小二乘反演的总体方程：
犌
犠狆犚
［］
狆
Δ犿＝
犱－犉（犿０）
－犠狆犚狆犿
［］
０
＋
犲ｏｂｓ
犠狆犲ｒ
［］
ｐ
（１１）引入反演迭代结束的条件为［２２］：①达到最大迭代次数；②拟合误差小于设定拟合误差。拟合差为：
犚犿狊＝
∑
狀
犼＝１
［（犱
犼－犉犼
（犿））／犱
犼
］２
槡犖（１２）
式中：犱
犼
观测数据向量元素；犉
犼
（犿）为理论正演响应值；犖为采样时间个数。
７
０
５
４期杜兴忠，等：半航空瞬变电磁法正反演算法及岩溶洼地实测数据验证   
４　岩溶地形模型建立
溶蚀地貌主要是由可溶岩的溶解、腐蚀，并伴随局部的基岩塌陷形成的。溶蚀洼地被看作是区域内地表水与地下岩溶网络的主要联通方式，其发育程度与该区域内大气降水与地表径流补给密切相关［１］。为尽量模拟实测地质构造的电磁响应特征，笔者对测线上各测点经正演模拟，然后按照前面加噪声的方式，对各测点添加白噪声和仪器背景噪声，最后形成含噪声的多测道电磁剖面。其中，剖面由各测点时间为０．０１ｍｓ～１０ｍｓ内５０个等对数间隔时窗构成。正演模拟参数见表１。
４．１　石林型溶蚀构造模式
溶蚀洼地近地表岩石裂隙在地表水和粘土不断冲刷与腐蚀作用不断向下扩张，由于可溶岩矿物组成的不同以及不同位置的裂隙发育速度不同，在多个位置形成地表水、裂隙水的小型排水中心。这些小型排水中心裂隙进一步向下发展，会将基岩腐蚀为石林型的构造模式。石林型溶蚀洼地构造可主要分为两层：不溶碎屑和细颗粒风化沉积物组成的覆盖层以及碳酸盐岩组成的基岩。对应电性特征：粘土覆盖层呈现相对低阻，基岩呈现相对高阻。本文正演模拟中覆盖层呈“Ｖ”型，在部分区域向下腐蚀，其中覆盖层电阻率取３０Ω·ｍ，基岩电阻率取１００Ω·ｍ，模型如图４所示。往正演模拟所得响应加入１０ｎＴ／ｓ白噪声和２％的仪器背景噪声，最后生成含噪声的多测道电磁响应剖面见图５。
表１　正演模拟参数表
Ｔａｂ．１　Ｆｏｒｗａｒｄｓｉｍｕｌａｔｉｏｎｐａｒａｍｅｔｅｒｔａｂｌｅ
线源长度飞行高度偏移距发射电流波形电流大小１０００ｍ３０ｍ５００ｍ双极方波１
０Ａ
图４　石林型溶蚀模型
Ｆｉｇ．４　Ｓｔｏｎｅｆｏｒｅｓｔｔｙｐｅｄｉｓｓｏｌｕｔｉｏｎｍｏｄｅ
ｌ
图５　石林型溶蚀多测道电磁响应剖面
Ｆｉｇ．５　Ｍｕｌｔｉ－ｍｅａｓｕｒｅｍｅｎｔｅｌｅｃｔｒｏｍａｇｎｅｔｉｃｒｅｓｐｏｎｓｅｐｒｏｆｉｌｅｏｆｓｔｏｎｅｆｏｒｅｓｔｔｙｐｅｄｉｓｓｏｌｕｔｉｏｎ
反演方法选择正则化反演和带横向约束的最小二乘
反演方法，反演初始模型为均匀半空间模型，模型第
一层厚度为２ｍ，层厚按１．０４指数增长，共３０层，
总深度为１０５ｍ。反演迭代结束的条件为：①最大
迭代次数５次；②迭代拟合差１％。最终的反演结
果如图６所示。
石林型溶蚀模型单点正则化反演结果中向基岩
腐蚀的低阻柱体大小和深度反映准确，但覆盖层和
基岩的分界面不清晰，部分反演结果受噪声影响产
生横向突变。横向约束反演结果中覆盖层和基岩的
分界面清晰；向基岩腐蚀的低阻柱体大小和深度反
映准确，但柱体下方的高阻基岩的典型特征反映一
般。
４．２　含隔水层型溶蚀构造模式
当碳酸岩中存在溶解难或者不溶的岩质隔水层
时，在隔水层未被侵蚀破坏前，上方通常存在包气带８
０
５    物探化探计算技术４５卷
图６　石林型溶蚀模型反演结果
Ｆｉｇ．６　Ｉｎｖｅｒｓｉｏｎｒｅｓｕｌｔｓｏｆｔｈｅｓｔｏｎｅｆｏｒｅｓｔｔｙｐ
ｅｄｉｓｓｏｌｕｔｉｏｎｍｏｄｅ
ｌ图７　含隔水层型溶蚀模型
Ｆｉｇ
．７　Ｄｉｓｓｏｌｕｔｉｏｎｍｏｄｅｌｗｉｔｈｗａｔｅｒｂａｒｒｉｅ
ｒ图８　含隔水层型溶蚀多测道电磁响应剖面
Ｆｉｇ．８　Ｍｕｌｔｉ－ｃｈａｎｎｅｌｅｌｅｃｔｒｏｍａｇｎｅｔｉｃｒｅｓｐｏｎｓｅｐｒｏｆｉｌｅｏｆａｑｕｉｆｅｒ－ｂｅａｒｉｎｇｄ
ｉｓｓｏｌｕｔｉｏｎ水，将导致隔水层上方的碳酸岩剧烈溶解，且岩溶地区地形起伏不定，地表水会透过覆盖层向基岩低洼处堆积。当隔水层被侵蚀破坏后，地表水将与地下碳酸岩接触，极有可能形成岩溶排水通道。含隔水层岩溶洼地构造可主要分为四部分：①不溶碎屑和细颗粒风化沉积物组成的覆盖层；②碳酸盐岩组成的基岩；③隔水层；④岩溶排水通道。对应电性特征：覆盖层呈现相对低阻，基岩呈现相对高阻。其中覆盖层电阻率取３０Ω·ｍ，基岩电阻率取１００Ω·ｍ，
隔水层电阻率为５００Ω·ｍ，岩溶排水通道电阻率为１０Ω·ｍ，模型如图７所示。按设计参数正演模拟所得响应再加入１０ｎＴ／ｓ白噪声和２％的仪器背景噪声，
最后生成含噪多测道电磁响应剖面见图８。
反演初始模型为均匀半空间模型，模型第一层厚度为２ｍ，层厚按１．０４指数增长，共３０层，
总深度为１０５ｍ。反演迭代结束的条件为［２２］
：①最大迭
代次数５次；②迭代拟合差１％。最终的反演结果如图９所示。
含隔水层溶蚀模型单点正则化反演结果中岩溶通道有一定地反映，但覆盖层和基岩的分界面不清晰，测线两端反演结果受噪声影响产生横向突变。横向约束反演结果中覆盖层和基岩的分界面较清晰；两条岩溶通道反映明显，但二者连通处岩溶通道反映不明显，只能看出部分倾斜趋势；高阻隔水层无明显反映。
９０５４期
杜兴忠，等：半航空瞬变电磁法正反演算法及岩溶洼地实测数据验证

版权声明：本站内容均来自互联网，仅供演示用，请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198，我们将在24小时内删除。