0 引言

大气降水入渗补给是地下水主要的补给来源之一^[1]。大气降水入渗补给量即大气降水通过包气带补给地下水的水量。影响大气降水入渗补给的因素很多，如土地利用类型、降水因素(包括降水强度、降水历时和降水量)、包气带厚度及岩性、地形等^[2]。随着社会经济的发展和城市化进程的推进，以及城市规模和土地利用类型的变化，大气降水入渗补给量也随之改变。因此，如何较为简洁并准确地计算城市化影响下的大气降水入渗补给量具有重要意义。

国内对大气降水入渗补给的研究较多。李亚峰等^[3]利用冉庄8 m定埋深地中蒸渗仪的观测资料，研究了大气降水入渗补给量与地下水水位埋深变化的关系，然而该方法只适用于室内实验研究，且需要较长时间的观测资料，在实际生产中应用难度较大；肖起模等^[4]通过对大气降水入渗补给系数与地层的相关分析，推求大气降水入渗补给量，该方法适用于地下水开发程度不高、植被覆盖度较高的山丘区，对于城市化发展程度较高的平原区并不适用。此外，国内很多研究和地下水调查评价在对区域大气降水入渗补给量进行计算时，未考虑土地利用类型对大气降水入渗补给量的影响^[5-7]。随着现代遥感技术的发展，遥感技术越来越多地被应用到水文水资源的研究领域^[8-11]。遥感数据具有周期短、信息量大和成本低的特点，可为水文水资源研究提供丰富的数据源^[12]。其应用大致可分为两个方面：一是利用遥感资料推求土地利用类型、地形等的变化；二是利用遥感资料进行水文过程中参数的推求^[13]。

本文利用遥感(RS)技术、遥感影像解译(ENVI)和地理信息系统(GIS)技术，在前人研究成果的基础上，将传统大气降水入渗补给量的计算方法与现代遥感技术相结合，通过对卫星遥感影像的解译，将研究区土地利用类型划分为硬化地面、水系及可入渗面三类；综合考虑土地利用类型、降水量、土壤岩性等因素，对大气降水入渗补给的有效区域进行面积划分和参数赋值，构建了一种基于土地利用类型的大气降水入渗补给量计算方法(precipitation infiltration recharge calculation based on land-use type，PIRCL)；以广州市广花盆地为例进行实例计算，并将计算结果与水量平衡法计算结果进行比较，以验证计算结果的准确性。将该研究成果应用到区域地下水资源调查评价中，能快捷准确地计算区域大气降水入渗补给量，在保障计算精度的同时，提升实际应用性；以期为指导地下水资源评价、地下水资源合理开发利用与保护、土地利用类型的合理配置提供科学依据。

1 研究区概况

广花盆地地处广州市西北部(图 1)，地理位置为112°95′15″E-113°43′55″E，23°14′51″N-23°50′14″N，总面积1 092 km²。受副热带高压控制，广花盆地属热带季风海洋气候，夏季盛行南风和东南风，冬季盛行北风和东北风。区内多年平均降水量为1 761.7 mm (图 2)。地势东北高、西南低。境内河流纵横交错，流溪河和白泥水及其支流分别自北、西北向东南贯穿全境。广花盆地地层由老至新有蓟县-青白口系、泥盆系、石炭系、二叠系、三叠系、侏罗系、白垩系、古近系和第四系。盆地地质构造单元属华南褶皱系的粤中凹陷区，盆地东部和北部与丘陵分界，以第四纪沉积地层为主，沉积地层下大面积分布覆盖型碳酸盐岩。广花盆地主要岩性为粉质黏土，局部为砂质黏土和细砂，给水度为0.02~0.06，大气降水入渗补给系数为0.033~0.201。广花盆地土地利用类型丰富，主要为现代化城市建筑区、水库、河流及农田耕地等。

2 大气降水入渗补给量计算 2.1 大气降水入渗补给量计算方法

首先确定研究区域适合的栅格大小，并对该区域进行栅格划分。分析土地利用类型、土壤类型、大气降水入渗补给系数、降水量等控制因素对大气降水入渗补给量的影响^[14]，将这些因素以栅格的形式输入ArcGIS 10.0。计算各栅格大气降水入渗补给量，然后将每个栅格计算得到的大气降水入渗补给量相加，得到整个研究区的大气降水入渗补给量，其计算公式为

其中，

式中：α为大气降水入渗补给系数；P为降水量(m)；F为可入渗面面积(m²)；μ为水位变动带含水层给水度；ΔH_次为降水造成的水位增幅(m)；P_次为ΔH_次对应时段的次降水量(m)。

2.2 基于VBA技术计算区域大气降水入渗补给量

上述方法可用基于VBA(visual basic for applications)技术的Arcobjects进行编程^[15]。过程如下：1) 数据初始化。将遥感解译得到的土地利用类型图(水域、硬化地面、可入渗地面三部分)和计算得到的大气降水入渗补给系数、降水量数据加载到ArcMap中。2) 将所得土地利用类型图输入ArcGIS 10.0系统，采用系统的“重采样”功能将栅格大小设置为500 m×500 m^[16]，同时也将研究区栅格大小确定为500 m×500 m。3) 利用克里金插值将大气降水入渗补给系数插值到整个研究区，栅格转点。栅格转点的目的是将研究区划分为多个栅格，提取每个栅格的中心点，以每个栅格中心点的大气降水入渗补给系数代表每个栅格的平均大气降水入渗补给系数，并以计算得到的每个栅格中心点的大气降水入渗补给量代表每个栅格的平均大气降水入渗补给量，因此栅格转点，以点代面。4) 将降水量栅格图栅格转点，转换到栅格中心点上。5) 对每个栅格的土地利用情况进行判定。如果是水域，则α=0；如果是可入渗地面即植被、裸土部分，则α即为插值得到的该位置点的值。同时考虑到城市中面积较小的花坛、草坪、树木等也接受大气降水入渗补给，但因解译误差均被直接划分为硬化地面，因此硬化地面部分降水入渗补给系数值取0.01^[17]。6) 创建新字段，命名为Q_降，含义为平均大气降水入渗补给量，按照公式(1) 进行计算。

3 广花盆地数据处理和输入 3.1 遥感影像分类及解译结果

在对遥感影像进行解译之前，第1步先对遥感影像进行预处理。遥感影像预处理是遥感影像处理过程中非常重要的环节，主要目的是消除图像中的无关信息，增强有关信息的可检测性，从而改进特征识别和提取的可靠性^[18]。

第2步进行遥感影像几何校正。遥感影像几何校正是利用地面控制点和几何校正数学模型来矫正非系统因素产生的误差，使其在空间上达到匹配，同时也是将图像投影到平面上的过程。达到高精度几何校正的关键是合理选取地面控制点(GCP)，地面控制点应该尽量选地物明显的点(比如道路交叉点、河流弯曲或者分叉处、飞机场等)，并尽可能满幅均匀选取，图像的边缘部分也要选取控制点以免图像发生外推现象。本次研究选取成像质量较好的2013年10月5日广州市10 m×10 m的SPOT-5高分辨率多光谱数据，影像共5个波段，使用321波段进行组合。影像校正方法选用Image to Image几何校正法，采用已经几何校正的本区1:5万地形图作为基准图，通过从地形图上选取同名点(或控制点)来配准高分遥感影像图，使同样的地物出现在校正以后图像的同一个位置。采用二次多项式校正法，选取18个控制点(图 3中标记)，采用双线性内插法进行重采样，完成几何校正，校正总体误差在1个像元以内。

第3步进行遥感数据切割。即通过研究区范围与遥感图像的叠置，切割出所需区域的遥感数据。接着利用ENVI软件，对经过预处理后的遥感影像进行监督分类，划分土地利用类型。本次研究将研究区土地利用类型分为3类，即硬化地面(城市中现代化建筑、公路等不透水面)、水系(河流、湖泊、水塘等)、可入渗面(植被、裸土、田地等)。通过建立各土地利用类型的遥感解译标志，在经过预处理的影像图上进行监督分类。

第4步检查样本代表性。为了评价样本的好坏，需要通过样本的基本统计值(如均值、最值、协方差矩阵等)检查其代表性，若分离度均大于1.8，则样本选取结果较好^[19]。本次研究的分离性检查结果表明，硬化地面与可入渗面、可入渗面与水系、硬化地面与水系可分离度分别为1.846、1.882、1.934，其可分离度均大于1.800，可见样本选取得较好。

第5步为结果分析。利用最大似然法进行监督分类，最终解译结果见图 4。利用混淆矩阵对解译结果进行评价，得到本次解译的总体分类精度达98.18%(表 1)；再对比遥感影像原图(图 3)，可见本次解译结果基本符合实际情况，解译结果合理可信。

3.2 大气降水入渗补给系数

大气降水入渗补给系数的计算：首先将研究区全年降水量与每一个自动监测井年内对应时段的地下水水位绘制成表，观察年内几场典型降雨所引起的水位增幅；再选取次降雨强度较大的若干场次降雨，划定对应次降雨时段的地下水水位增幅；然后利用公式(2) 进行计算，以多次降雨求得的降雨入渗补给系数平均值作为该自动监测井所在位置的降雨入渗补给系数。计算结果如表 2所示。

图 5为自动监测井A1-7自2014-01-22至2014-12-31的水位过程曲线及对应时段的降雨量，可见年内有3次较强降雨(即圈出的3部分)，并伴随着较大的水位波动；再结合钻井柱状图，可发现该自动监测井年内水位均在粉砂层波动，因此给水度取粉砂给水度；然后分别计算各次降雨时段内的降雨入渗补给系数；最后取自动监测井A1-7所在位置的大气降雨入渗补给系数为3次计算结果的平均值。

利用克里金插值法得到广花盆地大气降水入渗补给系数分布图，见图 6。可见广花盆地东南部大气降水入渗补给系数较大，西部中间地区及中部大气降水入渗补给较小。

4 降水入渗补给量计算结果分析 4.1 大气降水入渗补给量计算结果

本文降水量输入值为广州市2005-2014年10年平均降水量，根据上述方法最终计算求得广花盆地2005-2014年多年平均大气降水入渗补给量为14 369万m³。

4.2 计算结果与水量平衡法比较

根据水量平衡法^[20-21]，默认该地区多年地下水均衡差近似为0，即总补给Q_总补与总排泄Q_总排近似相等。广花盆地总排泄项主要由潜水蒸发量Q_潜、地下水开采量Q_采、地下水向河道排泄量Q_河1、地下水向湖库渠排泄量Q_湖库渠1及侧向排泄量Q_侧排 5项构成：

总补给项主要由大气降水入渗补给量Q_降、灌溉入渗补给量Q_灌、河流侧漏补给量Q_河、湖库渠侧漏补给量Q_湖库渠及山前侧向补给量Q_山前 5项构成：

通过计算其他补给项和排泄项反推该地区的多年平均大气降水入渗补给量：

计算结果见表 3。

根据表 2及公式(3) (5)，由水量平衡法计算可得广花盆地多年大气降水入渗补给量为14 036万m³，与PIRCL法的计算结果相比，相对误差约为2.4%，结果相近。2005-2014年逐年大气降水入渗补给量计算结果见表 4。

4.3 广花盆地大气降水入渗补给量空间分布

根据计算结果，绘制广花盆地的大气降水入渗补给量空间分布图以反映广花盆地内任意位置的大气降水入渗补给量(图 7)。因空间分布图反映的是点坐标年大气降水入渗补给量，与上文整个研究区面单位上的大气降水入渗补给量单位万m³有所区别，故单位为mm/a。对比参照图 4和图 6可以看出：淡蓝色区域解译得到了河流及水库，大气降水入渗补给量为0；黑色区域解译结果为建筑区、居民区、机场等，因大气降水入渗系数均取为0.01，故其单位面积的大气降水入渗补给量也相同且较小，主要分布在广花盆地中东部和东南部，该地区地下水主要接受侧向补给；绿色区域解译得到结果为可入渗面，主要分布在广花盆地的北部、西部及东南部分区域。实际结果表明：广花盆地西部和北部大多为植被覆盖下的山区或者耕地林地，有利于大气降水补给地下水；广花盆地东南部可入渗区，大气降水入渗补给系数最大，故单位面积的大气降水入渗补给量也最大。由此可见本次计算结果得到的大气降水入渗补给量空间分布基本合理。利用PIRCL法能够方便直接地观察出大气降水入渗补给量的空间分布情况。

5 结论与展望

1) 本次研究将传统的大气降水入渗补给量的计算与现代遥感技术相结合，提出了一种基于土地利用类型的大气降水入渗补给量的计算方法(PIRCL法)。首先对遥感影像进行解译，得出区域土地利用类型，进而合理分析研究区的有效降水入渗区域，结合大气降水入渗补给系数的空间分布、降水量等因素，对大气降水入渗补给量进行计算。相较于水均衡法，PIRCL法只需区域遥感数据、区域降水量、对应时段地下水位动态变化及区域水文地质参数，而不需要计算其他地下水的补给量和排泄量，简化了工作量。

2) 本次研究通过对卫星遥感影像的解译，将广州市广花盆地土地利用类型划分为硬化地面、水系及可入渗面3类，遥感影像样本选取合理，解译结果精确性较高，为98.18%。综合考虑土地利用类型、降水量、土壤岩性等因素，对大气降水入渗补给的有效区域进行面积划分和参数赋值，最终计算得到广州市广花盆地多年平均大气降水入渗补给量为14 369万m³。与水均衡法计算结果相比，相对误差约为2.4%。本次计算结果的空间分布基本合理，广花盆地西部、中部和东北部植被、耕地面积大，可直接接受大气降水补给面积较大，因此大气降水入渗补给量较大；北部、东中部及东南部城市建筑林立，建筑物占主要面积，降水入渗补给系数较小，因此大气降水入渗补给量较小。

3) 本次研究尚存在很多不足，比如对土地利用类型的划分还不是特别详细，仅分为3类，其原因之一也是由于区域各土地利用类型的降水入渗补给系数限制，以后可在此基础上加以改进和优化。