基于地球物理与水文数据融合的地下水水质评估方法与流程

1.本发明属于区域地下水资源调查技术领域，具体涉及一种地下水水质评估方法。


背景技术：

2.地下水咸化是沿海平原区和半干旱地区普遍存在的地质环境问题，过量开采导致的现代海水入侵和古海水越流混合咸化是加重地下水咸化的主要因素。查明地下水水质空间分布现状，包括咸淡水分界面和矿化度，对合理规划不同水质地下水用途和防止地下水进一步咸化意义重大。
3.传统方法多通过钻井的矿化度分析数据和电测井来获取点状矿化度和咸淡水分界数据，插值分别获取各含水层地下水水质空间分布；但钻井受限于高昂的费用和周期，数据量非常有限，不能全空间精细刻画地下水水质，并且如钻井工艺不达标可能会连通地下咸水和淡水，导致地下水局部咸化。
4.地球物理技术有成本低、点位密度高以及无损的优点，广泛应用于区域地下水资源调查。地面地球物理技术按物理性质分为：重力勘探、磁法勘探、电法勘探和地震勘探，淡水和咸水以及不同矿化度的地下水存在明显的导电性差异，使电法勘探成为划分咸淡水界面，评估地下水矿化度方面首选的、应用最为广泛的地球物理技术。常用的电法勘探分为电磁法和直流电法，电磁法又细分为频率域电磁法和时间域电磁法，直流电法又细分为电阻率剖面法和电阻率测深法。电阻率测深法可以同时获取水平方向和垂直方向电阻率特征，适用于地下水水质全空间调查。但现有的地球物理技术既没有通过测井数据来约束电阻率测深法，也没有按多层含水层的系统来分层评估，仅提供了一些笼统的地下水质分布结果，无法得到高精度的量化结果。


技术实现要素：

5.本发明所要解决的技术问题在于克服现有技术不足，提供一种基于地球物理与水文数据融合的地下水水质评估方法，将低精度、高密度的地球物理数据和高精度、低密度的水文数据融合，使多解的、笼统的三维地球物理参数对应至明确的、具体的三维水文地质模型参数。
6.本发明具体采用以下技术方案解决上述技术问题：
7.基于地球物理与水文数据融合的地下水水质评估方法，包括以下步骤：
8.步骤a、用电阻率测深法获取研究区域的换算电阻率数据，并获取研究区域中一组钻井的水文数据，包括水文地质分层数据、电测井数据和矿化度分析数据；
9.步骤b、利用所述电测井数据对换算电阻率数据进行校正，得到研究区域的校正换算电阻率数据；
10.步骤c、基于所述水文地质分层数据构建研究区域的多层含水层结构三维模型，并利用校正换算电阻率数据对所述多层含水层结构三维模型的空间分辨率进行优化；
11.步骤d、根据校正换算电阻率数据确定研究区域的咸淡水分界面；
12.步骤e、根据钻井的矿化度分析数据和钻井旁的平均校正换算电阻率，确定各含水层的矿化度
‑
电阻率相关数学模型，并根据所述数学模型，由各含水层的平均校正换算电阻率反演矿化度，建立各含水层矿化度空间分布。
13.优选地，根据换算电阻率数据对研究区域进行电性分区，并按照电性分区布置所述钻井，以保证每个电性分区都有至少一个钻井。
14.优选地，按照以下方法对换算电阻率数据进行校正：分别找出每个钻井的电测井数据中电阻率
‑
深度曲线以及钻井旁测点的换算电阻率
‑
深度曲线中电阻率变化最大的前n个拐点，按照深度从小到大的排序方式，用电测井数据中电阻率
‑
深度曲线中第i个拐点的深度与钻井旁测点的换算电阻率
‑
深度曲线中第i个拐点的深度之间的比值作为钻井旁测点的换算电阻率
‑
深度曲线中第i
‑
1个拐点至第i个拐点间深度数据所对应的深度校正系数，n为大于等于1的整数，i＝1,2,
…
,n；使用插值或拟合的方法获得研究区域中其它测点换算电阻率
‑
深度曲线的深度校正系数；最后用各电阻率测深法测点的深度校正系数对各测点的换算电阻率
‑
深度曲线中的相应深度数据进行校正，得到研究区域的校正换算电阻率数据。
15.进一步优选地，n为1或2。
16.优选地，具体使用以下方法构建研究区域的多层含水层结构三维模型：利用所述水文地质分层数据，通过插值或拟合计算出各含水层、隔水层以及基岩界面，从而构建多层含水层结构三维模型。
17.优选地，所述利用校正换算电阻率数据对所述多层含水层结构三维模型的空间分辨率进行优化，具体如下：对于所述多层含水层结构三维模型中校正换算电阻率数据存在明显电性差异的相邻含水层，结合校正换算电阻率数据构建空间分辨率更高的所述相邻含水层之间的含水层界面。
18.优选地，所述换算电阻率通过以下公式计算：
[0019][0020][0021][0022][0023]
式中，ρ(i)为第i层换算电阻率，ρ
s
(i)为第i个极距实测视电阻率，h(i)为第i层厚度，r(i)为第i个ab/2极距，k(i)为第i层反射系数，是双对数坐标下电测深曲线在第i
‑
1到第i个极距的斜率，k
′
(i)为斜率校正后的反射系数。
[0024]
相比现有技术，本发明技术方案具有以下有益效果：
[0025]
(1)本发明构建地球物理和水文数据融合流程，将高精度、低密度的水文数据(水
文地质分层、电测井、矿化度分析)和低精度、高密度的地球物理数据(电阻率测深法)融合，以较低成本获取高分辨率地下水水质空间分布，并形成与水文地质模型契合的定量结果。
[0026]
(2)本发明将直接反演方法用于咸淡水共存区的电阻率法反演，与迭代反演法相比，避免了等值现象引起的反演病态，使全区数据处理具有稳定性和统一性。
[0027]
(3)本发明利用电测井数据进行全区深度校正，使电阻率法在区内钻井的约束下，与水文地质结构相匹配。
[0028]
(4)本发明与水文系统中按含水层水质化验和矿化度分析一致，按含水层建立多层矿化度
‑
电阻率数学模型，分别获取每个含水层的矿化度分布。
附图说明
[0029]
图1是本发明方法的原理示意图；
[0030]
图2是具体实施例中地球物理和水文数据测量点位分布图，包含面积性布置的大量电阻率测深法测点和少量水文钻井；
[0031]
图3是具体实施例中由视电阻率计算换算电阻率的单点结果示例；
[0032]
图4是具体实施例中利用电测井数据和井旁换算电阻率求取校正因子并得到井旁校正换算电阻率的部分数据；
[0033]
图5是具体实施例中全区校正因子等值线图，用于非井旁校正换算电阻率；
[0034]
图6是具体实施例中建立的多层含水层结构三维模型；
[0035]
图7是具体实施例中以剖面方式批量解译咸淡水分界线的示例；
[0036]
图8是具体实施例中咸水上界面埋深图；
[0037]
图9是具体实施例中第ⅰ承压含水层平均电阻率；
[0038]
图10是具体实施例中建立的第ⅰ承压和第ⅱ承压含水层矿化度
‑
电阻率相关数学模型；
[0039]
图11是具体实施例中第ⅰ承压矿化度平面分布。
具体实施方式
[0040]
针对现有技术不足，本发明的解决思路是将低精度、高密度的大量地球物理数据和高精度、低密度的有限水文数据相融合，从而以较低的成本获得高分辨率地下水水质空间分布。
[0041]
具体而言，本发明所提出的基于地球物理与水文数据融合的地下水水质评估方法，其基本原理如图1所示，具体包括以下步骤：
[0042]
步骤a、用电阻率测深法获取研究区域的换算电阻率数据，并获取研究区域中一组钻井的水文数据，包括水文地质分层数据、电测井数据和矿化度分析数据：
[0043]
在研究区域面积性布置大量电阻率测深法测点，测点布置方式既可以采用平均分布的规则网，也可以采用灵活的不规则网；测点密度可根据地下水文地质条件复杂程度和比例尺要求确定，所需空间密度要求越高，地下水质空间分布越复杂，测点越密集。
[0044]
电阻率测深法采集的数据为视电阻率数据，需要将其反演为换算电阻率数据，可以采用现有的各种技术。考虑到传统的迭代反演法容易产生等值现象引起的反演病态，为了防止这一现象，从而使全区数据处理具有稳定性和统一性，本发明进一步提出了一种改
进的直接反演方法，基于采集到的视电阻率数据，通过以下公式直接反演出换算电阻率数据：
[0045][0046][0047][0048][0049]
式中，ρ(i)为第i层换算电阻率，ρ
s
(i)为第i个极距实测视电阻率，h(i)为第i层厚度，r(i)为第i个ab/2极距，k(i)为第i层反射系数，是双对数坐标下电测深曲线在第i
‑
1到第i个极距的斜率，k
′
(i)为斜率校正后的反射系数。
[0050]
在研究区域设置少量通过水文钻井得到水文数据，钻井的布置方式既可以与电阻率测深法同步进行，均匀分布于研究区，也可以在电阻率测深法的基础上，对钻井的布置进行优化：根据换算电阻率数据对研究区域进行电性分区，并按照电性分区布置所述钻井，以保证每个电性分区都有至少一个钻井。通过这些钻井得到的数据为水文地质分层数据、电测井数据和矿化度分析数据。通过水文钻井方式获取以上数据为现有成熟技术，为节省篇幅起见，此处不再赘述，更详细内容可参见《dd2019
‑
03水文地质调查技术要求(1：50000)》
[0051]
步骤b、利用所述电测井数据对换算电阻率数据进行校正，得到研究区域的校正换算电阻率数据：
[0052]
地球物理数据所具有的模糊性特点为深度、指向以及参数的模糊性。深度的模糊性表现在地球物理反演结果深度与地下物性分布深度通常存在整体偏离，即换算电阻率数据和地下实际电性分布存在深度的差异，这会给后续定量解释带来困难。为解决这一问题，本发明以测井数据作为联系电阻率测深法数据和钻井水文数据的纽带，以电测井数据作为基准对电阻率测深法所得到的换算电阻率数据进行校正，这是本发明实现水文和地球物理数据融合的第一步，也是关键一步。
[0053]
本发明具体按照以下方法对换算电阻率数据进行校正：分别找出每个钻井的电测井数据中电阻率
‑
深度曲线以及井旁电阻率测深法测点的换算电阻率
‑
深度曲线中电阻率变化最大的前n个拐点，按照深度从小到大的排序方式，用电测井数据中电阻率
‑
深度曲线中第i个拐点的深度与钻井旁测点的换算电阻率
‑
深度曲线中第i个拐点的深度之间的比值作为钻井旁测点的换算电阻率
‑
深度曲线中第i
‑
1个拐点至第i个拐点间深度数据所对应的深度校正系数，n为大于等于1的整数，i＝1,2,
…
,n；使用插值或拟合的方法获得研究区域中其它测点换算电阻率
‑
深度曲线的深度校正系数；最后用各测点的深度校正系数对各测点的换算电阻率
‑
深度曲线中的相应深度数据进行校正，得到研究区域的校正换算电阻率数据。
[0054]
通过大量实际测量研究发现，多数地质环境下所需测井深度拐点的数量n取1或2即可，当然也可以根据实际需要分段确定多个深度校正系数。
[0055]
通过以上的换算电阻率数据校正，可以将地球物理方法所获得的全部换算电阻率数据进行深度校正，使电阻率测深法数据在区内钻井的约束下，与水文地质结构相匹配，从而使整个研究区域的解释深度在钻井框架约束下。
[0056]
步骤c、基于所述水文地质分层数据构建研究区域的多层含水层结构三维模型，并利用校正换算电阻率数据对所述多层含水层结构三维模型的空间分辨率进行优化：
[0057]
地球物理数据的指向模糊性表现在地球物理数据是物理性质的分布，本身不具备地质含义，即换算电阻率不能明确对应根据时代、岩性所确定的水文地质分层。建立优化的多层含水层结构三维模型是本发明实现水文和地球物理数据融合的第二步，可实现校正换算电阻率与水文结构相对应。
[0058]
具体地，首先利用所述水文地质分层数据，通过差值或拟合计算出各含水层、隔水层以及基岩界面，从而构建多层含水层结构三维模型，此为现有技术，此处不再赘述；
[0059]
进一步地，对于所构建的多层含水层结构三维模型，如果其中某个含水层与相邻含水层的校正换算电阻率普遍存在明显的电性差异时，可以结合校正换算电阻率数据构建这两个含水层之间空间分辨率更高的含水层界面；对于不具备明显电性差异的相邻含水层，则仅使用水文地质分层数据。这样就实现了多层含水层结构三维模型的空间分辨率优化。
[0060]
步骤d、根据校正换算电阻率数据确定研究区域的咸淡水分界面：
[0061]
通过校正换算电阻率低阻层的界面即可确定咸淡水分界面，划分的依据是单点电阻率由低变高和由高变低的起跳点(即特征点)深度，通常以剖面的形式分批解译。汇总所有电阻率测深法的测点结果，可得到研究区域的全区咸淡水分界面。从而可在多层含水层结构三维模型的基础上，分析各含水层咸淡水分布。
[0062]
步骤e、根据钻井的矿化度分析数据和钻井旁的平均校正换算电阻率，确定各含水层的矿化度
‑
电阻率相关数学模型，并根据所述数学模型，由各含水层的平均校正换算电阻率反演矿化度，建立各含水层矿化度空间分布：
[0063]
参数的模糊性表现在地球物理参数为单维的，而水文地质参数是多维的，因此应采用主成分地质统计反演方法。水文地质参数反演是本发明水文和地球物理数据融合的第三步，在校正换算电阻率、多层含水层结构三维模型以及钻井矿化度分析数据的基础上，计算水文地质参数主成分—矿化度。
[0064]
根据钻井的矿化度分析数据和钻井旁的平均校正换算电阻率，可确定各含水层的矿化度
‑
电阻率相关数学模型；在多层含水层结构三维模型框架约束下，由校正换算电阻率可计算出各含水层平均电阻率；根据上述数学模型，将各含水层平均电阻率代入，即可由含水层平均电阻率反演矿化度，进而建立各含水层矿化度空间分布。
[0065]
为了验证本发明技术方案的效果，采用本发明方法对江苏沿海地区地下水水质进行评价，评价过程及结果如下：
[0066]
江苏省南通市地处长江三角洲沿海冲积平原，在长江和海洋共同作用下沉积了厚层的第四系和新近系松散岩，松散岩类孔隙含水层具有显著的分层特征，按照时代划分为五大含水岩组，分别为潜水含水层、ⅰ承压含水层、ⅱ承压含水层、ⅲ承压含水层、ⅳ承压含
水层。由于发生多次海水进退，导致区内大范围咸水分布，形成咸淡水共存的水质分布特点。
[0067]
如图2所示，研究区布置电测深测点共2000余个，西部网度1
×
0.5km，东部网度2
×
1km，ab/2供电极距最小5m，最大1000m，mn/2测量极距为ab/2的1/10。电极距系列为5/7.5/10/15/22.5/35/50/75/100/120/150/225/280/350/500/750/1000m；同步布置钻井26个，根据含水砂层的成因时代、埋藏分布、水力联系划分水文地质分层，并进行了测井，按含水层取水样进行水质化验。
[0068]
电阻率测深法所得到的视电阻率数据，使用本发明所改进的公式直接反演出换算电阻率数据，如图3所示。
[0069]
利用电测井数据对钻井旁换算电阻率数据求取深度校正因子并得到钻井旁的校正换算电阻率数据，如图4所示；然后使用插值的方法获得研究区域中非钻井旁测点的深度校正系数；最后用各测点的深度校正系数对各测点的换算电阻率
‑
深度曲线中的相应深度数据进行校正，得到研究区域的校正换算电阻率数据，如图5所示。
[0070]
基于所述水文地质分层数据构建研究区域的多层含水层结构三维模型，并利用校正换算电阻率数据对所述多层含水层结构三维模型的空间分辨率进行优化，优化后的多层含水层结构三维模型如图6所示。
[0071]
根据校正换算电阻率数据，以剖面批量解译咸淡水分界线，所得到的单条剖面如图7所示；汇总所有电阻率测深法的测点结果，得到如图8所示的全空间咸淡水分界面。
[0072]
在多层含水层模型顶底界面约束下，由校正换算电阻率计算各含水层的平均校正换算电阻率，如图9所示；根据钻井的矿化度分析数据和钻井旁的平均校正换算电阻率，确定各含水层的矿化度
‑
电阻率相关数学模型，如图10所示；并根据所述数学模型，由各含水层的平均校正换算电阻率反演矿化度，建立图11所示的各含水层矿化度空间分布。