左 欣 江 涛 张岱伟

(1.山东科技大学测绘科学与工程学院，山东 青岛 266590；2.中国科学院遥感与数字地球研究所，北京 100101； 3.海岛(礁)测绘技术国家测绘局重点实验室，山东 青岛 266590；4.中国测绘科学研究院，北京 100830)

基于农业生产胁迫下的地下水脆弱性评价与分析*

左 欣1，2江 涛1，3#张岱伟1，4

(1.山东科技大学测绘科学与工程学院，山东 青岛 266590；2.中国科学院遥感与数字地球研究所，北京 100101； 3.海岛(礁)测绘技术国家测绘局重点实验室，山东 青岛 266590；4.中国测绘科学研究院，北京 100830)

以山东省农业大市济宁市为研究区，在充分考虑地下水埋深、含水层介质等7个典型影响因子的同时，将含水层厚度、污染物毒性、污染物迁移性、污染物降解性纳入评价体系，对国际主流DRASTIC模型进行了因子扩充。综合运用地理信息系统(GIS)空间分析、熵值权重等技术，建立了基于农业生产胁迫下的地下水脆弱性评价模型，并进行了定性的评价分析。结果表明：济宁市地下水脆弱性总体呈现由西北向东南递减的空间格局，统计得到脆弱性很高区域占总面积的8.16%；梁山县的农业生产胁迫严重，地下水脆弱性很高；微山县的地下水脆弱性较低，受污染风险较小。

农业生产 地下水污染 脆弱性评价 地理信息系统空间分析DRASTIC

随着人口的急剧增加，社会经济迅速发展，需水量不断上升，人们开始进行大规模的地下水开采。然而，正是农业现代化和工业化的迅速发展，地下水资源面临着越来越大的污染威胁，水质恶化速度也随之加快[1]，越来越多的国家和地区面临着水危机。地下水环境污染已成为世界各国可持续发展的一大障碍。其中，农业生产对地下水环境造成的污染面积较大、治理较难。研究表明，由于治理经费高昂，含水层污染源复杂，所以地下水资源一旦遭到污染，直接修复几乎是不可能的。

法国MGARAT首次提出“地下水脆弱性”这一术语，他希望以此让人类了解地下水污染的危险性，并唤醒人类的环保意识[2]131-132。在后来的几十年，有关“地下水脆弱性”的概念及评价方法层出不穷，不再仅局限于水文地质内部的要素定义，而是慢慢趋向多角度、多方面来定义与评价。目前来说，迭置和指数法、统计分析方法、模糊数学法、过程数学模拟法是主要的地下水水质脆弱性评价方法[3]，国际上以迭置和指数法为主[4-5]，此方法的典型是美国环境保护署AILER于1987年提出的DRASTIC模型。主要对地下水埋深、含水层补给量、含水层介质、土壤类型、地形坡度、包气带介质以及含水层水力传导系数这7个参数进行研究，通过这些参数与地下水脆弱性之间的关系进行层次分析法赋权，并进行分级评价。陈守煜等[6]将含水层固有脆弱性与外界胁迫脆弱性相互作用的复杂关系简化成单一的叠加关系，进行了一系列的理论探索。另外，基于农业生产胁迫利用DRASTIC模型进行地下水脆弱性的分析评价案例也有很多。但是，分析发现，现有的研究方法均未将含水层厚度、污染物毒性、污染物迁移性、污染物降解性纳入评价范畴。

本研究以山东省济宁市为研究区，选用Landsat8-OLI影像数据，在充足的水文地质及气象数据基础上，针对DRASTIC模型进行改进，在研究传统的7个参数的同时，将含水层厚度、污染物毒性、污染物迁移性、污染物降解性纳入评价体系。综合运用地理信息系统(GIS)空间分析法、熵值权重法、统计分析法等多种技术，建立了基于农业生产胁迫下的地下水脆弱性评价模型，并分析地下水污染的高风险区及主要影响因素，为合理的地下水资源保护措施的提出以及农业生产规范化的推进提供可靠的参考依据。

1 研究基础

1.1 研究区概况

济宁市位于山东省的西南部，地处东经115°52′～117°36′、北纬34°26′～35°57′，辖区内共11个县区，属鲁南泰沂低山丘陵与鲁西南黄淮海平原交接地带，地质构造上属华北地区鲁西南断块凹陷区[7]。地形以平原洼地为主，地势东高西低，地貌较复杂；东部山峦绵亘，丘陵起伏。内有中国北方最大的淡水湖——南四湖。此外，济宁市地处暖温带大陆性季风气候区，年平均降水量680.5 mm。全市天然水资源总量约55亿m3，其中地下水天然补给量约21亿m3(占天然水资源总量的38.2%)；可利用水资源总量30.37亿m3，其中地下水12.93亿m3(占可利用水资源总量的42.6%)。并且，地下水年开采量达到94%以上，其中农业生产主要开采60 m以上的浅层地下水。济宁市是山东省农业大市，总土地面积约11 187 km2，其中耕地面积6 113.2 km2，占总面积的54.6%；园地面积97.5 km2，占总面积的0.9%。

1.2 数据获取

本研究使用了Landsat8影像、济宁市数字高程模型(DEM)等遥感数据源以及济宁市地下水埋深、土壤类型等非遥感数据源。所有数据经处理后均采用30 m分辨率、WGS84空间坐标系、UTM投影，具体数据如表1所示。

为助于研究过程的顺利进行，本研究对这些基础数据进行了辐射校正、坏值剔除、局部插值、图幅拼接等预处理操作，并根据需要对遥感影像提取归一化植被指数，对DEM提取地形坡度，对土壤数据进行了重分类得到土壤类型。

2 研究方法

2.1 技术流程

本研究主要考虑地下水固有脆弱性与地下水胁迫脆弱性两个方面，其中固有脆弱性受地下水埋深、含水层介质、土壤类型等8个因子影响；胁迫脆弱性受污染物毒性、迁移性、降解性综合作用的影响，其技术路线图如图1所示。

2.2 影响因子确定与模型建立

地下水脆弱性的影响因子众多，其中一些因子的获得相当困难，这些因子往往具有空间的不连续性，因此很难将影响地下水脆弱性的所有参数量化表示。而且，影响地下水脆弱性的各因子间的关系也错综复杂，很难找出一个统一的模型将所有因子全部纳入其中[8]。因此，现如今的研究多采用参数法，选择地下水脆弱性的主要且较易获得的影响因子对地下水的脆弱性进行评价。DRASTIC模型是一种参数法地下水固有脆弱性评价模型[9]，具有评价因子易于获取、评价方法简单易于掌握等优点，现已被多个国家的学者证实具有较高的精度与可行性[10]。

表1 试验数据列表

图1 技术路线示意图Fig.1 Technology roadmap

但DRASTIC模型也具有一定的局限性，表现在：(1)只针对地下水的固有脆弱性进行评价，没有考虑人类的生产生活对地下水脆弱性的胁迫影响，评价体系不全面；(2)评价因子间具有很高的相关性，各评价因子相互之间不独立，如含水层补给量除了与输入系统的水量有关外还与地形坡度、土壤类型以及包气带介质有关[11]；(3)DRASTIC模型的7个参数中有4个是可以定量获取的，其余3个则需要人为确定，因此评价过程具有较高的人为主观性；(4)各评价因子的得分为离散的值，将本来连续变化的地理事物离散化，使评价结果不够客观。

由于DRASTIC模型具有很多的局限性，因此在使用时不应一味生搬硬套。应当结合研究区的情况以及研究的目的对该模型做一些适当的修改，使模型的评价结果能更好贴近研究区的实际情况。本研究以DRASTIC模型为基础，做了一定的改进。首先，增加了含水层厚度影响因子，表征含水层对污染的稀释作用。其次，地下水的脆弱性是地下水的固有脆弱性与地下水外界胁迫脆弱性共同作用的结果，两者都是影响地下水脆弱性的重要因素。仅靠DRASTIC模型只能评价地下水的固有脆弱性，对地下水脆弱性的评价非常片面。因为在某些人迹罕至的区域，由于没有人类污染的胁迫，即使该地区地下水的固有脆弱性很高，该地区地下水受到污染的可能性也不大[2]134-136。本研究引入污染物的毒性、迁移性、降解性量化指标计算污染源对地下水环境的胁迫指数作为对DRASTIC模型的补充。

表2 可量化因子评分

2.3 地下水固有脆弱性影响因子计算

地下水固有脆弱性8个影响因子中，地下水埋深、含水层补给量、地形坡度、含水层水力传导系数以及含水层厚度5个影响因子可以量化得到，评分如表2所示。

2.3.1 地下水埋深

地下水埋藏得越深污染物就越难到达，地下水就越难受到污染。DRASTIC模型主要用于评价浅层地下水的固有脆弱性，地下水埋深即浅水面的深度。本研究以统计数据进行计算，地下水埋深评分见图2(a)。

2.3.2 含水层补给量

含水层补给量指各种形式的水入渗补给的总量，由于污染必须溶解在水中才能随水入渗污染地下水，因此含水层补给量越大污染物随补给水入渗污染地下水的机会越多。含水层补给量的计算方法如式(1)所示，得到济宁市含水层补给量评分图(见图2(b))。

Rin=P-E-ΔRf-Rs

(1)

式中：Rin为含水层补给量，mm；P为输入水总量，mm；E为蒸发蒸腾量，mm；Rf为田间持水量的变化量，mm；Rs为地表径流量，mm。

由于农业区土壤的田间持水量随时间的波动较小，因此在计算时可忽略田间持水量的变化量。

蒸发蒸腾量的计算分为两部分：(1)考虑植被的地表蒸发蒸腾量计算，采用综合法中的改进彭曼-蒙蒂斯公式(见式(2))[12]；(2)水面蒸发量(E1，mm/d)的计算(见式(3))[13]。另外，入渗率a与地形和土壤的质地有关，土壤越疏松地形越平缓入渗率越高。

ETOP-M=

(2)

(3)

式中：ETOP-M为参考作物蒸发蒸腾量，mm/d；ε为饱和水汽压对温度的导数，kPa/℃；Rn为净辐射量，MJ/(m2·d)；G为土壤热通量，MJ/(m2·d)；γ为湿度表常数，kPa/℃；T为平均气温，℃；U2为2 m高处风速，m/s；ea为饱和水汽压，kPa；ed为实际水汽压，kPa；Δe为自然大水体表面水温对应的饱和水汽压与水面以上150 cm高处空气水汽压的差值，kPa；ΔT为自然大水体表面水温与水面以上150 cm高处气温的差值，℃；r1.5为自然大水体水面以上150 cm高处的相对湿度，%；U1.5为自然大水体水面以上150 cm高处的风速，m/s。

地表径流量计算见式(4)至式(6)。

Rs=ΨgP-ΨgPT

(4)

Ψg=0.52e-860BlD0.036ln(arctanSl)

(5)

(6)

式中：Ψg为产流系数，%；PT为产流阈值，mm；Bl为基岩水力传导系数，mm/h；D为土壤厚度，m；Sl为地形坡度，%；θΔ为雨后含水率，%；θ0为先期含水率，%。

2.3.3 地形坡度

地形坡度表示地表单元陡缓的程度，是坡面垂直距离与水平距离的比值，在栅格计算中使用式(7)计算每个点的地形坡度，可得到研究区地形坡度评分图(见图2(c))。

Sl=ΔHmax/L×100%

(7)

2.3.4 含水层水力传导系数

注：由于东南部是微山湖，故此区域地形坡度、包气带介质、土壤类型、含水层介质定义为0分。图2 影响因子评分Fig.2 Rating for factors

水力传导系数又称渗透系数，是表征流体通过空隙骨架难易程度的参数，水力传导系数越大空隙介质的透水性越强。含水层水力传导系数(K，m/d)的计算方法如式(8)所示。

K=kρg/η

(8)

式中：k为孔隙介质渗透率，m2；ρ为流体密度，kg/m3；g为重力加速度，m/d2；η为动力黏滞性系数，kg/(m·d)。

对于某一种流体而言，重力加速度与流体密度都为常数，其水力传导系数取决于孔隙介质渗透率与动力黏滞性系数。当流体一定时，孔隙介质渗透率与动力黏滞性系数主要与固体骨架的性质有关。在式(8)理论支持下，得到了研究区的含水层水力传导系数评分图(见图2(d))。

2.3.5 含水层厚度

含水层厚度是影响含水层对污染稀释能力的一个重要因素，含水层越厚含水层中赋存的地下水的量越大，对进入含水层的污染物稀释能力也就越强[14]。本研究将含水层厚度作为一个附加因子添加在DRASTIC模型中，作为对DRASTIC模型的扩展与完善，并进行了含水层厚度评分，结果见图2(e)。

2.3.6 含水层介质、土壤类型与包气带介质

含水层介质、土壤类型与包气带介质无法定量化获得，本研究根据各影响因子不同介质类型来确定评分(见表3)，得到各影响因子评分图(见图2(f)、图2(g)和图2(h))。

表3 介质类型评分

2.4 农业胁迫指数计算

胁迫脆弱性是地下水脆弱性的重要组成部分。人类的生产生活中会产生大量的污染物，各种污染物的性质千差万别，对地下水环境的影响也不同。本研究主要针对农业区的地下水进行脆弱性评价，主要考虑农业生产中产生的污染物对地下水环境的影响。农业生产中产生的污染物主要有化肥、农药以及一些常规的综合性指标，如COD、总溶解固体(TDS)等[15]。农业生产中使用的化肥与农药的种类很多，且不同地区不同作物施用的农药与化肥的种类也不同。本研究从我国“水中优先控制污染物”名单中的58种有机毒物及无机毒物筛选出主要的污染物[16]并进行评价。为便于量化的表示污染物对地下水环境的影响，选取污染物的毒性、迁移性与降解性作为评价污染物对环境影响的因子[17-19]。其中，毒性表示一定量的污染物对环境的破坏能力，污染物的毒性越强，污染物对环境的破坏力越大[20]；迁移性表示一定量的污染物在地下水环境中的污染范围，污染物的迁移性越强，污染物在地下水中的污染范围越大；降解性表示污染在地下水中的持续时间，污染物的降解性越弱，污染在地下水中持续的时间越长[21]。针对污染物的3种属性采用参数叠加与层次分析方法，同时考虑污染物排放量，使用式(9)与式(10)计算得到农业生产对济宁市地下水环境的胁迫指数分布。

S=∑CiQi

(9)

Ci=wtDti+wmDmi+wdDdi

(10)

式中：S为农业胁迫指数，即农业生产对地下水环境的危害性；Ci为第i种污染物对地下水环境的危害性；Qi为第i种污染物排放量系数；wt、wm、wd分别为第i种污染物的毒性、迁移性、降解性量化指标的权重；Dti、Dmi、Ddi分别为第i种污染物的毒性、迁移性、降解性量化指标。

为更明显表示出农业生产的胁迫性，研究将污染物毒性、降解性、迁移性归于农业胁迫因子，并得到了研究区农业胁迫指数图(见图2(i))。

3 地下水脆弱性评价与分析

本研究采用熵值权重法结合DRASTIC模型的经验权重，利用3步公式确定各指标的综合权重，从而对各因子间的相关性进行一定的抑制，增强了评价结果的客观性。

(1) 计算第j个指标的信息熵(Cj)：

Cj=-b∑pllnpl

(11)

式中：pl为指标为第l种状态的概率；b为常数。

(2) 计算第j个指标的熵值权重(βj)：

βj=(1-Cj)/∑(1-Cj)

(12)

(3) 确定第j个指标的综合权重(wj)：

wj=ajβj/∑ajβj

(13)

式中：aj为第j个指标的经验权重。

根据式(11)至式(13)计算得出综合权重，结果见表4。使用参数叠加法将影响地下水脆弱性的各因子及其相应综合权重代入式(14)，得到济宁市地下水脆弱性(F)评价图(见图3)。

F=∑wjDj+S

(14)

式中：Dj为地下水固有脆弱性的第j个指标的评分。

表4 地下水脆弱性影响因子综合权重

图3 济宁市地下水脆弱性评价Fig.3 Vulnerability assessment of groundwater in Jining City

为使结果更清晰且易于分析，本研究根据几何间隔进行了脆弱性分级，将研究区划分为脆弱性很高、脆弱性较高、脆弱性中等、脆弱性较低、脆弱性很低5个等级，并得到地下水脆弱性分级图(见图4)。

根据图4，济宁市地下水脆弱性总体呈现由西北向东南递减的空间格局，统计得到脆弱性很高区域占总面积的8.16%。在济宁市西北部的梁山县地下水脆弱性很高，而东南部的微山县地下水脆弱性较低，这种分布情况与当地的水文地质条件相吻合，微山县水域广阔，农业胁迫指数小，并且土壤的碎石体积与含沙量较少，而梁山县农业胁迫指数大，土壤的碎石体积与含沙量较多。由于地下水埋深和土壤类型在脆弱性评价中所占权重较大，故造成此差异。

图4 济宁市地下水脆弱性分级Fig.4 Hierarchical groundwater vulnerability in Jining City

4 结 语

(1) 在现有DRASTIC模型的基础上，对DRASTIC模型做了一定的改进与完善，综合考虑地下水的固有脆弱性与外界胁迫脆弱性，得到较全面准确的基于农业生产胁迫下地下水脆弱性评价结果。

(2) 人类农业生产过程中产生的污染胁迫是地下水脆弱性评价中的重要影响因素。现有的大部分评价方法都只对地下水的固有脆弱性进行评价是非常片面的，这样的评价方法不能客观全面地反映地下水受污染的风险性。

(3) 梁山县在农业生产过程中施用了大量的农药与化肥，并且当地的地下水脆弱性极高，当地的浅层地下水极易受到农业生产的污染。而且，此地的含水层较厚、地下水储量大，发生污染后治理困难。当地政府应采取措施限制梁山县的农药化肥使用量，防止地下水进一步受到污染。

(4) 嘉祥县、泗水县以及济宁市市辖区，由于自然与人为因素的关系，也存在一些脆弱性较高的区域，当地政府也应采取相应的措施限制农业生产中农药化肥的使用量。

(5) 汶上县的地下水脆弱性虽然很低，但此地的地下水固有脆弱性较高，容易受到人类生产生活的污染。而且，此地浅层地下水储量大，发生污染后治理困难。因此，当地政府应注意污染物的排放量，防止该地区地下水受到污染。

(6) 本模型是基于污染物随水入渗理论建立的，不适用于河流、湖区的地下水脆弱性评价。但是，微山湖区由于地表水与地下水相连通，水中的溶解物会相互交换，易造成双向污染，因此要注意湖区水体的污染防治工作。

[1] 郭德伟，陆继东，王卫东，等.济宁市地下水资源面临的问题及保护措施[J].治淮，2014(1)：12-13.

[2] 张树军，张丽君，王学凤，等.基于综合方法的地下水污染脆弱性评价——以山东济宁市浅层地下水为例[J].地质学报，2009，83(1).

[3] 李吉学，孔舒，宋秀真，等.济宁地下水水质脆弱性研究[J].治淮，2011(12)：3-5.

[4] SHUTTLEWORTH W J.Evaporation models in the global water budget[M].Berlin：Springer Netherlands，1983.

[5] DOLMAN A J，GREGORY D.The parametrization of rainfall interception in GCMs[J].Quarterly Journal of the Royal Meteorological Society，2006，118(505)：455-467.

[6] 陈守煜，伏广涛.含水层脆弱性模糊分析评价模型与方法[J].水利学报，2002，33(7)：23-30.

[7] 王营.济宁市城区土地储备制度研究[J].科技经济市场，2013(2)：34-36.

[8] 孙才支，潘俊.地下水脆弱性的概念、评价方法与研究前景[J].水科学进展，1999，10(4)：445-446.

[9] ALLER L，BENNETT T，LEHR J H，et al.DRASTIC：a standardized system for evaluating groundwater pollution potential using hydrogeogic settings[R].San Francisco：U.S.EPA，1987.

[10] 杨贵芳，姜月华，李云.基于DRACSTIC模型的城市地下水脆弱性评价综述[J].地下水，2012，34(1)：5-7.

[11] 钟佐燊.地下水防污性能评价方法探讨[J].地学前缘，2005，12(增刊1)：3-13.

[12] 彭世彰，徐俊增.参考作物蒸发蒸腾量计算方法的应用比较[J].灌溉排水学报，2004，23(6)：5-6.

[13] 闵骞.水面蒸发量计算公式的研究[J].水资源研究，2004，25(1)：38-40.

[14] 鄂建，孙爱荣，钟新永.DRASTIC模型的缺陷与改进方法探讨[J].水文地质工程地质，2010，37(1)：102-107.

[15] ZAPOROZEC A，CONRAD J E，JOHANSSON P O，et al.Groundwater contamination inventory：a methodological guide [R].Paris：United Nations Educational，Scientific and Cultural Organization，2002.

[16] 王昭，石建省，张兆吉，等.我国“水中优先控制有机物”对地下水污染的预警性研究[J].水资源保护，2009，25(1)：90-94.

[17] 陆燕，何江涛，王俊杰，等.北京平原区地下水污染源识别与危害性分级[J].环境科学，2012，33(5)：1526-1528.

[18] 王俊杰，何江涛，陆燕，等.地下水污染风险评价中特征污染物量化方法探讨[J].环境科学，2012，33(3)：771-775.

[19] 罗文彧，金鸣岐.水迁移系数与岩溶发育程度之间关系的初步探讨[J].中国岩溶，1985(增刊1)：20-27.

[20] ARONSON D，BOETHLING R，HOWARD P，et al.Estimating biodegradation half-lives for use in chemical screening[J].Chemosphere，2006，63(11)：1953-1960.

[21] 孟宪萌，束龙仓，卢耀如.基于熵权的改进DRASTIC模型在地下水脆弱性评价中的应用[J].水利学报，2007，38(1)：94-98.

Vulnerabilityassessmentofgroundwaterunderthestressofagriculturalproduction

ZUOXin1，2，JIANGTao1，3，ZHANGDaiwei1，4.

(1.GeomaticsCollege，ShandongUniversityofScienceandTechnology，QingdaoShandong266590；2.InstituteofRemoteSensingandDigitalEarth，ChineseAcademyofSciences，Beijing100101;3.KeyLaboratoryofSurveyingandMappingTechnologyonIslandandReef，StateBureauofSurveyingandMapping，QingdaoShandong266590; 4.ChineseAcademyofSurveyingandMapping，Beijing100830)

Based on the geology，meteorology，hydrology and the remote sensing data，this study regarded the Jining City as the study site which was a city based on agriculture in Shandong Province. Meanwhile seven typical indexes were taken into consideration including the thickness of the aquifer，contaminant toxicity，migratory aptitude and degradation. It was an important supplement of the method of DRASTIC，which was the most popular and advanced method to evaluate the groundwater vulnerability. Moreover，the groundwater vulnerability evaluation model was built based on GIS spatial analysis and the entropy-based weight. The results showed that the groundwater vulnerability gradually present a spatial pattern from northwest to southeast in Jining City，with the highest vulnerability accounted for 8.16% of the total area. Liangshan had the most severe agriculture production stress，and the highest groundwater vulnerability. Moreover，Weishan had lower groundwater vulnerability and lower pollution risk.

agricultural production; groundwater pollution; vulnerability assessment; GIS spatial analysis; DRASTIC

2016-03-14)

左 欣，男，1993年生，本科，研究方向为遥感科学与技术。#

。

*国家自然科学基金资助项目(No.41471331)；国家科技部科技支撑计划课题(No.2012BAH27B04)；国家级大学生创新创业项目(No.201410424010)；山东省高等学校科技计划项目(No.J14LH05)。

10.15985/j.cnki.1001-3865.2017.05.003