CN106772639A - 地下铁质管线埋深磁偶极子构造法优化反演 - Google Patents

Abstract

本发明涉及地下铁质管线埋深磁偶极子构造法优化反演方法，具体涉及一种基于磁偶极子构造法正演的地下铁质管线埋深优化反演计算方法。建立磁异常探测数据预处理、地磁场定位查询及反演搜索模型参数替代，将十二个正演模型参数简化为埋深和磁化率两个未知模型参数。埋深和管道磁化率两个未知模型参数作为全局优化反演求解参数，使用初次反演得到埋深和管道磁化率作为第二次反演参数搜索范围参考点，在新搜索范围内再次进行多次反演搜索，取埋深反演平均值作为地下铁质管线埋深优化反演结果。本发明地下铁质管线埋深磁偶极子构造法优化反演兼顾正演计算速率和精度，反演求解快速且准确度高。

Description

地下铁质管线埋深磁偶极子构造法优化反演

技术领域

本发明涉及一种使用探测磁异常信号反演求解埋地铁质管线埋深的方法，具体涉及地下铁质管线埋深磁偶极子构造法优化反演。

背景技术

地下埋设铁质管线主要有污水、雨水、给水、电缆、燃气、电信、电力和石油管线，为国家工业生产和城镇居民生活中提供基础保障作用。

社会经济的发展，地下铁质管线数量的增多，其建设、维护和管理的成本越来越高，因此地下管线的探测变得越来越重要。随着技术的不断发展和完善，对各种管线都有自身适应的探测方法，其勘测施工难度、探测深度范围、经济投入及准确性各不相同。磁异常铁质管线探测作为一种被动源管线探测方法，使用地球自身的磁场作为激励磁源，不需要人工提供磁源发生装置，使得该方法相对其他管线探测方法而言具有硬件设备简单，操作容易，投入成本低等，探测位置准确优点。

目前，使用磁异常探测地下铁质管线的商业仪器是基于探测总磁异常分布来判断地下管线的有无及大致走向，而不能计算地下铁质管线的埋深。专门针对地下铁质管线磁异常探测，利用数据处理方式来计算管道埋深的方法报道较少，已有技术报道使用磁异常化极估算的方法计算管线埋深精度较低，且容易受噪声影响。另有论文报道使用差分能量谱估计法估算出埋深在印度某小镇附近的一段供水管的深度为3.88±0.10m，而使用磁感应强度大小与磁异常梯度张量的解析信号的比值的方法来估计这段管道的埋深时，其估计埋深为4.09m，且两种方法都未给出准确性分析。另有文献报道使用磁异常梯度的方法探测估计某水管的埋深时，深度估计误差能够小于5％。其它使用磁数据处理估算磁源模型深度的方法主要有Werner反褶积法，解析信号，Euler反褶积法，局部波数法，及源参数成像法等，这些方法尚未在地下铁质管线探测埋深计算方面应用，无法知道其准确性如何。

基于磁异常反演求解地下铁质管线埋深的方法较少，本发明前期研究报道使用遗传算法基于四个反演参数进行过埋深反演，但该方法对于埋深较浅的管线需要进行分块单元划分，反演花费时间较长。其他相关磁异常反演方法尚未应用于地下铁质管线探测埋深反演，无法判断其有效性。

发明内容

本发明的目的是针对地下铁质管线磁异常埋深求解现有的问题，提供了一种地下铁质管线埋深磁偶极子构造法优化反演求解方法。

本发明是通过以下技术方案实现地下铁质管线埋深磁偶极子构造法优化反演，具体步骤如下：

(1)建立正演模型：以地理北为X轴，地理东为Y轴，垂直向下为Z轴，探测地平面为XY平面，建立计算坐标系；坐标系中，设所致探测磁异常有效管道长度为L_p，管道外径为φ，管道厚度为δ，管道磁化率为χ_m；设管道方位角为θ(0≤θ≤π)，管道轴线埋深为h₁；地磁背景场由总强度F，磁倾角I_c及磁偏角D三个独立的地磁要素表征；磁异常正演时，管道轴线位于坐标系原点正下方，管道路径垂直于管道轴线，管道路径的中点在Z轴上，高度为h₂，长度为L_d，管道路径方位角为η；将管道路径上总磁异常大小或总磁异常分量B_d表示为十二个参数表征的函数：

B_d＝f(L_p,φ,δ,θ,h₁,F,I_c,D,χ_m,h₂,L_d,η) ①

(2)测点磁异常计算：将管道分割成有限数量的单元块，将管道分割所得单元视作坐标为单元几何中心的磁偶极子，使用磁偶极子磁场计算公式计算每一个单元在测点处的磁感应强度，使用矢量和计算所有单元在测点处的叠加磁感应强度，该叠加磁感应强度加上该点的地磁场强度即为测点磁异常；设单元块E_i几何中心坐标为(E_ix,E_iy,E_iz)，单元体积为V_i，测点P_j的坐标为(P_jx,P_jy,P_jz)；则有磁偶极子磁场计算公式可得单元块E_i在测点P_j处的磁感应强度B_ij三分量为：

上式②中μ₀为真空磁导率；当V_i为分块单元时当V_i为分块单元时V_i＝πφδ(φ-δ)；上式②中，H_x＝Fcos(I_c)cos(D)，H_y＝Fcos(I_c)sin(D)，H_z＝Fsin(I_c)；则根据磁场叠加原理，地下铁质管线在探测面P_j点的磁感应强度三分量为：

上式③中I为构造地下铁质管道的单元总数；当管道模型参数和单元划分方式确定后，单元块E_i几何中心坐标即随之确定；当探测路径及测点间隔确定后，测点P_j坐标也为随之确定；当使用探测总磁异常大小作为观测值时，则使用模拟总磁异常作为匹配值；当使用探测总磁异常轴向分量作为观测值时，则使用③式中对应的轴向分量作为匹配值；测线上所有测点的总磁异常的大小或总磁异常分量组成的向量为B^sim；

(3)探测计算管道方位角：初步选择几条探测直线路径探测地下铁质管线所致总磁异常分布，确定管道大致走向后，在垂直管道方向上选取两条相隔一定距离的探测路径测量总磁异常，使用曲线拟合的方式找到两条探测曲线的峰值点，使用两个峰值点对应的水平坐标(B_#1x,B_#1y)、(B_#2x,B_#2y)计算管道方位角为：

当θ确认后，测线垂直于管道轴线方位角η等于90°；使用两个峰值点连线初步确定管道的水平位置；

(4)磁异常测量：垂直管道轴线方向布设探测路径，探测路径长度L_d范围在16m到30m之间；探测路径中点设于初步确定的管道的正上方；每隔10cm左右测量记录一次总磁异常或总磁异常分离，得到探测磁异常向量B^obs；

(5)建立反演模型：使用GPS定位得到探测区域坐标，查询国际地磁参考场得到地磁总强度F，磁倾角I_c及磁偏角D；反演时使用有限长管道作为反演管道长度，有效管道长度L_p在18m至30m范围选取；测线高度h₂在探测磁异常时确定；反演时，使用替代管道模型厚度δ与管道外径φ作为确定的参数，管道厚度δ选取范围为5mm至15mm，管道外径φ选取范围为20cm至80cm；此时建立参数反演模型：

(h₁,χ_m)＝f^-1(B^obs) ⑤

(6)反演数据预处理及反演目标函数的确定：使用插值方法找到模拟计算磁异常向量B^sim的最大值或者最小值，使用拟合方式找到探测磁异常向量B^obs拟合曲线上的最大值或者最小值，以上述两个最大值点或两个最小值点作为反演数据定位点，按一定测点间隔选取测线长度为16m至20m的探测磁异常向量B^obs的拟合计算得到目标磁异常向量并选取对应坐标点由计算磁异常B^sim插值计算得到匹配磁异常向量更新⑤式反演模拟为建立反演目标函数：

(7)管道埋深参数优化反演：使用全局优化算法搜索找到式⑥中使得目标函数尽量小的埋深h₁和磁化率χ_m，优化反演得到的埋深h₁即为地下铁质管线的轴线埋深。

对上述内容的进一步补充，使用探测总磁异常的数值大小作为目标求解磁异常向量时，计算总磁异常的数值大小作为匹配磁异常使用探测磁场某一坐标轴向分量，如探测磁异常垂直分量作为目标求解磁异常向量时，计算磁异常某一坐标轴向分量，如计算磁异常垂直分量作为匹配磁异常向量

对上述内容的进一步补充，地下铁质管线埋深磁偶极子构造法优化反演，所述单元划分，将探测路径中点正下方对应的管道轴线上的点作为单元划分的开始端点，沿管道轴线方向向管道两端分割单元；首先将管道从起始点开始分割成长度为管道直径的节单元，当节单元的几何中心到测线中点的距离大于6.5倍的管道直径时，保留该节单元，继续分割节单元，直至长度为探测路径长度L_d的管道单元划分完毕；当所分节单元的几何中心到测线中点的距离小于6.5倍的管道直径时，取消该节单元划分，沿着轴向分割为环，环宽度为δ，将每一个环沿圆周分割为块单元，每一个管道环分割数为π(φ-δ)/δ的向上取整数，每分割一环块单元后再次分割节单元并判断是否满足分割条件，若满足则继续进行节单元分割，如不满足则进行块单元分割，直至长度为探测路径长度L_d的长管道单元划分完毕。

对上述内容的进一步补充，为使管道埋深反演达到较好的准确效果，磁异常探测传感器精度不低于1nT。

对上述内容的进一步补充，为使管道埋深反演准确、快速，反演时确定模型管道长度为20m，测线长度为20m，测点间隔为10cm，管道厚度为10mm，管道外径为50cm。

对上述内容的进一步补充，为使管道埋深反演快速且准确，使用粒子群算法进行埋深参数反演求解；以(h₁,χ_m)参数作为二维搜索粒子，粒子群进化代数p＝40，每一代的粒子数s＝30，参数反演包括两次遗传算法反演搜索，第一次反演管道埋深搜索范围为[0.2,6](m)，磁化率搜索范围为[0.001,100](SI)；第二次反演搜索区间为第一次反演结果的中心扩展区间，以反演结果作为每个参数的搜索中心点向上和向下扩展20％作为新的搜索空间；在新搜索范围内再次进行多次反演搜索，取埋深反演平均值作为地下铁质管线埋深优化反演结果。

本发明相比现有技术具有以下优点：

(1)本发明结合地下铁质管线自身特有的结构及分布特征，使用参数反演的方法求解地下铁质管线的中间轴线埋深，能够直接得出管线轴线深度，准确度高，稳定性好；

(2)本发明基于磁偶极子构造法，使用模拟预处理与模型参数替代方法将十二个磁异常正演参数减少为两个可变模型参数，有效降低了磁异常正演计算量，节约计算时间；

(3)本发明可基于现有地下铁质管线磁异常探测仪器测量数据，应用本发明提出的方法求解地下管线的埋深，对于已有该硬件设备的使用者可以减小硬件投入；

(4)本发明埋深反演具有一定抗噪声能力，相比现有管线埋深反演计算而言误差小。

附图说明

图1铁质管线埋深磁偶极子构造法正演模型。

图2管道分块单元划分示意图。

图3粒子群算法垂直磁场埋深优化反演流程图。

具体实施方式

实施例

磁偶极子构造法经验证为有效的地下铁质管线磁异常正演计算方法，本实施例以磁偶极子构造法基于理论模型计算的垂直磁场数据加上5nT随机噪声作为磁场观测值。地下管线磁异常探测垂直磁场观测值原模型如下：

表1

在被探测区域确定正北方向，先初步选择几条探测直线路径探测地下铁质管线所致总磁异常分布，确定管道大致走向后，在管道正上方选择一点作为原点，以正北方向作为X轴，正东方向作为Y轴；在所建立的坐标系中，在垂直管道方向上选取两条相隔10m的探测路径，并在路径上每隔20cm确定一个探测点并记录探测点坐标，测量每个坐标点的总磁异常，使用10次方曲线拟合的方式找到两条探测曲线的峰值点，使用两个峰值点对应的水平坐标(B_#1x,B_#1y)、(B_#2x,B_#2y)计算管道方位角为

在探测地面连接两条探测曲线的峰值点，垂直于峰值点连线在以上两条探测路径中间重新划定一条长24m的探测线，测线中点位于峰值点连线上；在测线上每隔10cm测量一个垂直磁场分量，所有的垂直磁场分量按坐标顺序构成探测磁异常向量B^obs；基于GPS定位信息通过最新国际地磁参考场查找被探测地点地磁场作为已知信息。

如图1所示，在正演坐标系中，令替代模型中点位于坐标系原点下方，参数管道长度为20m，探测路径长度为20m，测点间隔为10cm，管道厚度为10mm，管道外径为50cm；以X＝(h₁,χ_m)参数作为二维搜索粒子，使用粒子群算法进行粒子搜索，埋深参数初始值为1m，磁化率参数初始值为10SI。

如图1和图2所示，将管道使用式②和式③计算模拟磁异常垂直分量，得到B^sim；使用插值方法找到模拟计算磁异常B^sim的最大值(或者最小值)，使用拟合方式找到探测磁异常Bo^bs拟合曲线上的最大值(或者最小值)，以此两点作为反演数据定位点，按一定测点间隔选取测线长度为18m的探测总磁异常Bo^bs的拟合计算所得向量并选取对应坐标点由B^sim插值计算得到向量

图3为粒子群埋深反演流程图，如图3所示，将管道埋深和磁化率两参数构成的二维向量视为搜索空间中的粒子，每一个粒子对应一个式⑥计算的适应度，使用粒子群速度更新公式和粒子位置更新公式计算粒子速度和粒子在搜索空间中的位置，同时为了防止算法限于局部最优，引入变异方法对第p代第s个粒子的第d维变量进行突变操作；搜索粒子群进化代数p＝40，每一代的粒子数s＝30，第一次反演管道埋深搜索范围为[0.2,6](m)，磁化率搜索范围为[0.001,100](SI)；基本粒子群速度更新与粒子位置计算式如下：

第p代第s个粒子第d维变量变异方法计算式如下：

式⑦⑧⑨中相关参数说明如下：

第p代第s个粒子的搜索速度，p＝1,2，…,P，s＝1,2，…，S；

第p+1代第s个粒子的搜索速度；

c₁、c₂：学习因子，c₁＝c₂＝2.05；

c₃：约束因子，c₃＝0.8；

c₄：约束因子，c₄＝0.2；

w：惯性权重，w＝1.02；

rand(1)：[0,1]范围内的随机值，每次运行该函数时结果都不同；

第p代第s个粒子；

第s个粒子p代及以前的个体最优值；

第p代及以前的全局最优值；

第p代第s个粒子的第d维变量；

X_dmax、X_dmin：第d维变量变化在其范围内的最大及最小值。

反演求解包括两个过程，第二次反演搜索以第一次反演求解一次得到的埋深和磁化率作为中心点向上和向下扩展20％作为新的搜索空间；所述计算个体最优值与全局最优适应值的终止条件为：判断全局最优粒子目标函数值是否小于0.1或判断粒子群进化代数是否达到最大值40，如果任一条件满足，则终止粒子群优化搜索并输出搜索结果；第二次反演对管道埋深进行10次求解，去掉一个最大埋深和一个最小埋深，剩下8次反演埋深如下表2所示：

表2

表2中粒子群管道埋深反演误差均小于5％，如果取8次反演结果平均值作为反演埋深，则埋深反演误差小于1％，此结果说明了地下铁质管线埋深磁偶极子构造法优化反演的有效性。

Claims (6)

1.地下铁质管线埋深磁偶极子构造法优化反演，其特征在于，具体步骤如下：

(1)建立正演模型：以地理北为X轴，地理东为Y轴，垂直向下为Z轴，探测地平面为XY平面，建立计算坐标系；坐标系中，设所致探测磁异常有效管道长度为L_p，管道外径为φ，管道厚度为δ，管道磁化率为χ_m；设管道方位角为θ(0≤θ≤π)，管道轴线埋深为h₁；地磁背景场由总强度F、磁倾角I_c及磁偏角D三个独立的地磁要素表征；磁异常正演时，管道轴线位于坐标系原点正下方，探测路径垂直于管道轴线，探测路径的中点在Z轴上，高度为h₂，长度为L_d，探测路径方位角为η；将探测路径上总磁异常大小或总磁异常分量B_d表示为十二个参数表征的函数：

B_d＝f(L_p,φ,δ,θ,h₁,F,I_c,D,χ_m,h₂,L_d,η)；

(2)测点磁异常计算：将管道分割成有限数量的单元块，将所得单元视作坐标为单元几何中心的磁偶极子，使用磁偶极子磁场计算公式计算每一个单元所在测点处的磁感应强度，使用矢量和计算所有单元在测点处的叠加磁感应强度，该叠加磁感应强度加上该点的地磁场强度即为测点磁异常；设测线上所有测点的总磁异常的大小或总磁异常分量组成的向量为计算磁异常向量B^sim；

(3)探测计算管道方位角：初步选择几条探测直线路径探测地下铁质管线所致总磁异常分布，确定管道大致走向后，在垂直管道方向上选取两条相隔一定距离的探测路径测量总磁异常，使用曲线拟合的方式找到两条探测曲线的峰值点，使用两个峰值点坐标计算管道方位角θ，测线垂直于管道轴线方位角η等于90°,初步确定管道的水平位置；

(4)磁异常测量：垂直管道轴线方向布设探测路径，探测路径长度L_d范围在16m到30m之间，探测路径中点设于初步确定的管道的正上方，每隔10cm左右测量记录一次总磁异常或总磁异常分离，得到探测磁异常向量B^obs；

(5)建立反演模型：使用GPS定位得到探测区域坐标，查询国际地磁参考场得到地磁总强度F，磁倾角I_c及磁偏角D；反演时使用有限长管道作为反演管道长度，有效管道长度L_p在18m至30m范围选取；测线高度h₂在探测磁异常时确定；反演时，使用替代管道厚度δ与管道外径φ作为确定的参数，管道厚度δ选取范围为5mm至15mm，管道外径φ选取范围为20cm至80cm；此时建立参数反演模型：

(h₁,χ_m)＝f^-1(B^obs)；

(6)反演数据预处理及反演目标函数的确定：使用插值方法找到模拟计算磁异常向量B^sim的最大值或者最小值，使用拟合方式找到探测磁异常向量B^obs拟合曲线上的最大值或者最小值，以上述两个最大值点或两个最小值点作为反演数据定位点，按一定测点间隔选取测线长度为16m至20m的探测磁异常向量B^obs的拟合计算得到目标磁异常向量并选取对应坐标点由计算磁异常B^sim插值计算得到匹配磁异常向量更新反演模拟为建立反演目标函数：

(7)管道埋深参数优化反演：使用全局优化算法找到使得目标函数尽量小的埋深h₁和磁化率χ_m，优化反演得到的埋深h₁即为地下铁质管线的轴线埋深。

2.如权利要求1所述地下铁质管线埋深磁偶极子构造法优化反演，其特征在于，使用探测总磁异常的数值大小作为目标求解磁异常向量时，计算总磁异常的数值大小作为匹配磁异常使用探测磁场某一坐标轴向分量，如探测磁异常垂直分量作为目标求解磁异常向量时，计算磁异常某一坐标轴向分量，如计算磁异常垂直分量作为匹配磁异常向量

3.如权利要求1所述地下铁质管线埋深磁偶极子构造法优化反演，其单元划分特征在于，将探测路径中点正下方对应的管道轴线上的点作为单元划分的开始端点，沿管道轴线方向向管道两端分割单元，首先将管道从起始点开始分割成长度为管道直径的节单元，当节单元的几何中心到测线中点的距离大于6.5倍的管道直径时，保留该节单元，继续分割节单元，直至长度为探测路径长度L_d的管道单元划分完毕；当所分节单元的几何中心到测线中点的距离小于6.5倍的管道直径时，取消该节单元划分，沿着轴向分割为环，环宽度为δ，将每一个环沿圆周分割为块单元，每一个管道环分割数为π(φ-δ)/δ的向上取整数，每分割一环块单元后再次分割节单元，并判断是否满足分割条件，若满足则继续进行节单元分割，如不满足则进行块单元分割，直至长度为探测路径长度L_d的管道单元划分完毕。

4.如权利要求2所述地下铁质管线埋深磁偶极子构造法优化反演，其特征在于，磁异常探测传感器精度不低于1nT。

5.如权利要求2所述地下铁质管线埋深磁偶极子构造法优化反演，其特征在于，反演时确定模型管道长度为20m，测线长度为20m，测点间隔为10cm，管道厚度为10mm，管道外径为50cm。

6.如权利要求2所述地下铁质管线埋深磁偶极子构造法优化反演，其特征在于，使用粒子群算法进行埋深参数反演求解；以(h₁,χ_m)参数作为二维搜索粒子，粒子群进化代数p＝40，每一代的粒子数s＝30；参数反演包括两次遗传算法反演搜索，第一次反演管道埋深搜索范围为[0.2,6](m)，磁化率搜索范围为[0.001,100](SI)；第二次反演搜索区间为第一次反演结果的中心扩展区间，以反演结果作为每个参数的搜索中心点向上和向下扩展20％作为新的搜索空间；在新搜索范围内再次进行多次反演搜索，取埋深反演平均值作为地下铁质管线埋深优化反演结果。