一种大距离位场向上延拓计算方法技术

针对传统位场基于快速傅里叶变换算法的波数域向上延拓方法因边界效应而导致向上延拓计算精度低的问题，本发明专利技术公开了一种大距离位场向上延拓计算方法。本发明专利技术通过设定应用场景、加权系数计算，利用加权系数和一种二维离散卷积快速算法，实现位场向上延拓快速、高精度计算，实现了大规模位场数据大距离向上延拓在计算效率和计算精度上的统一。本发明专利技术解决了传统FFT法因边界效应而导致大距离向上延拓计算精度低的问题，为满足大规模位场数据处理、位场分离、大深度位场向下延拓等提供了方法支撑。

【技术实现步骤摘要】
一种大距离位场向上延拓计算方法
本专利技术涉及应用于重磁数据处理领域的位场向上延拓计算方法，特别是大距离位场向上延拓快速、高精度数值计算方法。
技术介绍
重力场和磁场统称位场。位场延拓是重磁数据处理中的一种重要方法，包括向上延拓、向下延拓、曲化平、平化曲和平化平。向上延拓在位场延拓问题中最简单，也是最重要的一种，它是实现其它几类位场延拓问题的基础。传统向上延拓方法是基于快速傅里叶变换算法(FFT)的波数域向上延拓方法，通常称为向上延拓FFT法。FFT法借助快速傅里叶变换算法的高效性，能够快速实现大规模位场数据的向上延拓。但是，由于存在离散傅里叶变换边界效应问题，FFT法大距离向上延拓计算精度低，通常采用某种扩边方法来减小边界效应，但扩边也会引入数据误差，降低延拓结果精度。文献(汪炳柱.用样条函数法求重力异常二阶垂向导数和向上延拓计算.石油地球物理勘探,1996,31(3):415-422.)公开了一种重力异常向上延拓的样条函数法，该方法延拓精度高但是计算效率低。文献(陈龙伟,胡小平,吴美平,马涛.空间域位场延拓新方法研究.地球物理学进展,2012,27(4):1509-1518.)和文献(YileZhang,YauShuWong,YuanfangLin.BTTB-RRCGmethodfordownwardcontinuationofpotentialfielddata.JournalofAppliedGeophysics,2016,126:74-86.)分别公开了一种相同的空间域向上延拓方法，根据系数矩阵具有BTTB矩阵结构的特点，采用BTTB矩阵与向量相乘的快速算法，实现了向上延拓快速计算，并且有效克服了FFT法存在的边界效应问题。该方法存在的问题是向上延拓距离小时，系数矩阵接近奇异，延拓结果不稳定。目前已有的向上延拓方法存在的最大问题是不能同时保证计算效率和计算精度，无法满足快速、高精度、大距离向上延拓的要求。因此，寻找一种计算效率高、同时能保证计算精度的位场延拓方法具有重要的现实意义。
技术实现思路
本专利技术针对传统位场基于快速傅里叶变换算法(FFT)的波数域向上延拓方法因边界效应而导致向上延拓计算精度低的问题，提出了一种大距离位场向上延拓计算方法，其是一种新的空间域大距离位场向上延拓快速、高精度数值计算方法。为解决上述技术问题，本专利技术采用以下技术方案：一种大距离位场向上延拓计算方法，其步骤为：步骤一：设定应用场景建立直角坐标系，z轴垂直向下为正。给定高度为z0的规则网格化位场数据图，规则网格化位场数据图其x方向上等间距分布有M个网格节点，其y方向上等间距分布有N个网格节点，这样规则网格化位场数据图上有M×N个网格节点，M×N个网格节点对应M×N个规则网格化位场数据。给定向上延拓高度Δz，Δz为正数；称高度为z0的水平面为观测面，称高度为z0-Δz的水平面为延拓面；高度为z0的规则网格化位场数据图上的M×N个规则网格化位场数据记为u(ξi,ηj,z0)，i＝1,2,…,M，j＝1,2,…,N；(ξi,ηj,z0)表示观测面网格节点坐标，u(ξi,ηj,z0)表示网格节点坐标(ξi,ηj,z0)处观测面位场值；M表示规则网格化位场数据图其x方向上的网格节点数，N表示规则网格化位场数据图其y方向上的网格节点数。高度为z0-Δz的延拓面上的规则网格化位场数据图上的M×N个规则网格化位场数据记为u(xm,yn,z0-Δz)，(xm,yn,z0-Δz)表示延拓面网格节点坐标，u(xm,yn,z0-Δz)表示网格节点坐标(xm,yn,z0-Δz)处延拓面位场值，m＝1,2,…,M，n＝1,2,…,N。步骤二：根据公式(1)，计算加权系数式(1)中，ω(xm,yn；ξi,ηj)表示加权系数；(xm,yn)表示延拓面网格节点坐标；(ξi,ηj)表示观测面网格节点坐标；arctan()表示反正切函数；Δz表示向上延拓高度；μpq，Xp，Yq，Rpq(p＝1,2；q＝1,2)表示辅助计算变量，其计算式为X1＝xm-ξi+0.5Δx，X2＝xm-ξi-0.5ΔxY1＝ym-ηi+0.5Δy，Y2＝ym-ηi-0.5Δy式中，Δx表示x方向上相邻网格节点间的间距，Δy表示y方向上相邻网格节点间的间距；步骤三：位场向上延拓：根据公式(1)给出的加权系数和一种二维离散卷积快速算法，实现公式(2)给出的位场向上延拓快速、高精度计算式中，(xm,yn,z0-Δz)表示延拓面网格坐标，m＝1,2,…,M，n＝1,2,…,N；(ξi,ηj,z0)表示观测面网格坐标，i＝1,2,…,M，j＝1,2,…,N；z0-Δz表示延拓面高度；z0表示观测面高度；u(xm,yn,z0-Δz)表示网格坐标(xm,yn,z0-Δz)处延拓面位场值；u(ξi,ηj,z0)表示网格坐标(ξi,ηj,z0)处观测面位场值。与现有位场延拓方法相比，采用公式(1)给出的方法计算加权系数，能够提高向上延拓的计算精度，并且保证小距离向上延拓时数值计算的稳定性。步骤三中，根据公式(1)给出的加权系数和一种二维离散卷积快速算法，实现公式(2)给出的位场向上延拓快速、高精度计算，其步骤为：(1)将加权系数ω(x1,y1；ξi,ηj)排列成矩阵t，记为式(3)中，矩阵t中的各矩阵元素用ti,j表示；矩阵t中的矩阵元素ti,j与加权系数ω(x1,y1；ξi,ηj)存在关系ti,j＝ω(x1,y1；ξi,ηj)(4)(2)将矩阵t补零扩展成矩阵将矩阵分成四个块矩阵，记为将块矩阵互换位置，得到矩阵cext(3)将观测面位场数据u(ξi,ηj,z0)(i＝1,2,…,M，j＝1,2,…,N)排列成矩阵g，矩阵元素gi,j与位场数据存在关系gi,j＝u(ξi,ηj,z0)(8)将矩阵g补零扩展成矩阵gext式中，表示Nx×Ny零矩阵；(4)计算式中，fft2()表示二维快速傅里叶变换；(5)计算式中，“.*”表示对应元素相乘运算；(6)计算式中，ifft2()表示二维快速傅里叶反变换；(7)提取矩阵fext的前M行前N列，构成矩阵f，即为二维离散卷积计算结果，也就是位场向上延拓结果。与现有位场向上延拓方法相比，采用上述二维离散卷积快速算法是本专利专利技术点之一，它保证了在大规模位场数据情况下向上延拓的高效性。本专利技术是一个有机整体，根据步骤二提供的一种特殊的加权系数计算公式，采用步骤三运用的二维离散卷积快速算法，实现了大距离位场向上延拓在效率和精度上的统一。本专利技术解决了传统FFT法因边界效应而导致大距离向上延拓计算精度低的问题，为满足大规模位场数据处理、位场分离、大深度位场向下延拓等提供了方法支撑。与现有位场延拓方法相比，本专利技术具有以下优点：(1)能够实现大规模位场数据情况下大距离向上延拓快速、高精度数值计算；(2)能够保证小距离位场向上延拓情况下数值计算的稳定性；(3)大规模位场数据向上延拓时，算法不但计算效率和计算精度高，并且所需计算机内存小。附图说明：图1为大距离位场向上延拓计算流程图；图2为观测面和延拓面网格坐标示意图；图3为棱柱体磁异常组合模型；图4为z＝0m高度面磁异常理论值；图5为z＝20000m高度面磁异常理论值；图6为z＝20000m高度面FFT法向上延拓结果；图本文档来自技高网...

【技术保护点】
一种大距离位场向上延拓计算方法，其特征在于，包括以下步骤：步骤一：设定应用场景建立直角坐标系，z轴垂直向下为正；给定高度为z0的规则网格化位场数据图，规则网格化位场数据图其x方向上等间距分布有M个网格节点，其y方向上等间距分布有N个网格节点，这样规则网格化位场数据图上有M×N个网格节点，M×N个网格节点对应M×N个规则网格化位场数据；给定向上延拓高度Δz，Δz为正数；称高度为z0的水平面为观测面，称高度为z0‑Δz的水平面为延拓面；高度为z0的规则网格化位场数据图上的M×N个规则网格化位场数据记为u(ξi,ηj,z0)，i＝1,2,…,M，j＝1,2,…,N；(ξi,ηj,z0)表示观测面网格节点坐标，u(ξi,ηj,z0)表示网格节点坐标(ξi,ηj,z0)处观测面位场值；M表示规则网格化位场数据图其x方向上的网格节点数，N表示规则网格化位场数据图其y方向上的网格节点数；高度为z0‑Δz的延拓面上的规则网格化位场数据图上的M×N个规则网格化位场数据记为u(xm,yn,z0‑Δz)，(xm,yn,z0‑Δz)表示延拓面网格节点坐标，u(xm,yn,z0‑Δz)表示网格节点坐标(xm,yn,z0‑Δz)处延拓面位场值，m＝1,2,…,M，n＝1,2,…,N；步骤二：根据公式(1)，计算加权系数...

【技术特征摘要】
1.一种大距离位场向上延拓计算方法，其特征在于，包括以下步骤：步骤一：设定应用场景建立直角坐标系，z轴垂直向下为正；给定高度为z0的规则网格化位场数据图，规则网格化位场数据图其x方向上等间距分布有M个网格节点，其y方向上等间距分布有N个网格节点，这样规则网格化位场数据图上有M×N个网格节点，M×N个网格节点对应M×N个规则网格化位场数据；给定向上延拓高度Δz，Δz为正数；称高度为z0的水平面为观测面，称高度为z0-Δz的水平面为延拓面；高度为z0的规则网格化位场数据图上的M×N个规则网格化位场数据记为u(ξi,ηj,z0)，i＝1,2,…,M，j＝1,2,…,N；(ξi,ηj,z0)表示观测面网格节点坐标，u(ξi,ηj,z0)表示网格节点坐标(ξi,ηj,z0)处观测面位场值；M表示规则网格化位场数据图其x方向上的网格节点数，N表示规则网格化位场数据图其y方向上的网格节点数；高度为z0-Δz的延拓面上的规则网格化位场数据图上的M×N个规则网格化位场数据记为u(xm,yn,z0-Δz)，(xm,yn,z0-Δz)表示延拓面网格节点坐标，u(xm,yn,z0-Δz)表示网格节点坐标(xm,yn,z0-Δz)处延拓面位场值，m＝1,2,…,M，n＝1,2,…,N；步骤二：根据公式(1)，计算加权系数式(1)中，ω(xm,yn；ξi,ηj)表示加权系数；(xm,yn)表示延拓面网格节点坐标；(ξi,ηj)表示观测面网格节点坐标；arctan()表示反正切函数；Δz表示向上延拓高度；μpq，Xp，Yq，Rpq表示辅助计算变量，其计算式为：X1＝xm-ξi+0.5Δx，X2＝xm-ξi-0.5ΔxY1＝ym-ηi+0.5Δy，Y2＝ym-ηi-0.5Δyμpq＝(-1)p(-1)q其中：p＝1,2；q＝1,2；Δx表示x方向上相邻网格节点间的间距，Δy表示y方向上相邻网格节点间的间距；步骤三：...

【专利技术属性】
技术研发人员：陈龙伟，杨长福，王宇航，彭健，
申请(专利权)人：桂林理工大学，
类型：发明
国别省市：广西,45