白 天, 张敬丽, 易 娇, 吴雅文
(云南农业大学 园林园艺学院, 云南 昆明 650201)
城镇化改变了城市下垫面地表景观特征,地表景观类型、格局、分布等改变,显著降低了城市下垫面雨水截留和入渗能力。城市地表景观改造彻底改变了原有自然水文过程,加速了城市地表径流局部汇集,使城市区域易受到历时短、强降雨的侵袭而引起城市内涝,给城市带来灾难性的后果。城市地形改造和用地类型是引起城市局部及周围环境雨洪灾害发生的重要因素。随着系统计算理论与技术的发展,基于淹没的理论模型,包括SWMM,MIKE URBAN,HSPF,SWAT,InfoWorks CS等模型广泛运用于中国城市的淹没分析与内涝管理,在城市排水、防洪、环境治理方面取得一定成果[1-2]。但是城市雨洪格局的分析方法仍然有很大研究潜力,以满足城市地表景观、城市地形和土地利用等规划、设计与管理的需求。因此,针对中国城市发展需求的城市地表径流空间分布研究成为了海绵城市建设的重点。汤鹏等[3]对扬州江都区进行产流研究以分析地表径流的分布的特征;姚磊等[4]通过分析北京市产流,得到了地表径流的空间分布规律 ;李孝永等[5]评估了北京市土地利用景观对地表产流与雨洪调节服务的影响。云南省昆明市自2016年已经开始海绵城市建设,但城区夏季雨洪灾害仍然显著。目前关于昆明市地表径流网络与空间分析研究尚未见报道。地表径流网络特征研究目前还处于初级阶段,本次研究将MCR模型与空间自相关结合,构建地表径流网络进行空间分析,以期为城市土地利用规划和海绵城市建设提供科学依据。
昆明市地处云贵高原滇池流域坝区,102°10′—103°40′E,24°23′—26°22′N,受印度洋西南暖湿气流影响雨水充沛,全年降水区间377.1~1 455.6 mm,主要集中在夏秋季5—10月,据1951至2020年统计,5—10月平均降雨量876 mm,约占全年88%。结合2000—2016年昆明国家气候基准站的气象资料,仅2016年昆明市发生城市内涝44例[6-7]。虽然城市积水中心由2016年的40个下降至2020年31个,但夏季暴雨条件下,既使城市37.92%的绿地率与41.62%绿化覆盖率渗透环境,仍未能有效改善城市地表积水的现状。2020年雨季昆明城市暴雨造成大量积水,影响道路交通安全和居民正常活动,地表径流涌入地下车库导致大量车辆受损,漂浮污染物扩散,引起地表水资源的再次污染。
采用ArcGIS 10.5将Spot-5高清卫星影像(2020年2月13日,云量0%),进行融合拼接、地理配准、几何纠正处理,从而获得分辨率为0.5 m的昆明市中心城区54.88 km2正射影像,影像清晰无植被遮挡,可有效分辨城市各类用地边界;将WorldView 2高程精度为1 m数字表面模型 (DSM)转换为城市竖向二维数据,构建GIS数据库[8-9]。
2.2.1 目视解译与汇水区划分 采用目视解译将研究区土地覆盖/利用 (LUCC)类型划分为2大类,7小类,即透水地表P:城市绿地P1(市政绿化用地,公园绿地),水体P2(城市水系,湿地,洼地,开沟),草地P3(公园草坪、运动场草坪),闲置用地P4(裸露土地,城市施工场地);不透水地表I:道路I1(城市交通车型行道路、人行道路),建筑I2(居住、工业、商业区建筑屋顶),综合用地I3(铺装广场、停车场)(图1a)。
采用ArcGIS 10.5将WorldView 2高程精度为1 m数字表面模型 (DSM)转换为城市竖向二维数据(图1b)。子汇水区以城市地表径流单元的地形、河道、街道、市政管网、建筑位置等物理参数为基础进行划分。基于DSM地形构建泰森多边形,根据节点划分汇水区,再根据实际道路、市政网络、建筑位置等信息进行局部细化,明确子汇水区域 ,共划分为51个子汇水区(图1c)[10-11]。
2.2.2 AHP指标与评价 采用层次分析法 (AHP)利用AHP Khaskia软件,对各景观影响因素进行权重分析,分为3层[12-13],即目标层 (地表景观类型),标准层(地形因素与景观因素),因子层 (LUCC、坡度、坡向、起伏度、粗糙度等11项因子)。Towsif等[9],Silva等[14]和Chigbu等[15]学者的研究为子汇水区间阻力影响要素强度的确定提供了参考依据。对各层指标两两比较,得到判断矩阵,并进行一致性检验 (CR<0.1)。
2.2.3 最小阻力累计耗费距离模型 (MCR) 潜在径流从“源”到“汇”的运动过程中存在景观阻力差异,采用最小阻力累计耗费距离模型 (minimum cumulative resistance, MCR)计算空间内流体质点塑性流动方向与路径,并获得从源头到末端最小成本路径[16]。此方法广泛运用于地表径流特征,流域水土资源和生态河网保护等流动扩散研究[17-22]。首先,通过以上AHP法进行权重分析,得到景观阻力消费面 (cost surface);其次,统计各子汇水区景观阻力面的最低值点作为潜在径流的“源”,而从一个子汇水区通过阻力影响汇入下一个子汇水区的过程为潜在径流的“汇”过程,最后,统计潜在径流的流动路径,该路径符合液体流动规律是径流扩散的最易路径[23-24]。公式为:
(1)
式中:Dij为源单元i到汇单元j的空间距离;Ri为源单元i到汇单元j过渡过程中存在的阻力系数;Ri由i的值决定,从i到j单元的路径产生不同的电阻值,当i确定后,计算MCR需要选择空间距离中电阻值最小的路径。
注:P1城市绿地; P2水体; P3草地; P4闲置用地; I1道路; I2建筑; I3公共用地。
2.2.4 潜在径流空间分布的卡方检验 卡方检验是用途广泛的假设检验方法,用于统计观测样本与模拟值之间的拟合程度[25]。在潜在径流分布图中随机选取36个汇流点,统计流经汇流点潜在径流数量(N),与降雨1 h对应汇流点的径流深度(D)进行卡方检验,试验分别于2020年7月1日,7月23日和8月17日进行3次重复测量,验证模拟值与测量值之间的拟合程度,并计算其显著性(p)。公式为:
(2)
2.2.5 重力模型 (gravity model) 建立各子汇水区间作用力矩阵,计算子汇水区之间的相互作用强度[26]。根据矩阵结果,筛选不同强度下的径流路径单元密度,从而得到研究区地表径流分布路径、网络结构、时间顺序[27]。公式为:
式中:Gab是子汇水区间雨水径流路径汇水区a和汇水区b之间的相互作用;Na和Nb分别是a,b汇水区之间的权重值;Dab是a,b子汇水区间潜在径流路径阻力的标准化值;Pa,Pb分别是子汇水区a,b的阻力值;Sa,Sb分别是子汇水区a,b的面积;Lab是子汇水区a到b雨水径流路径累计阻力值;Lmax是研究区潜在雨水径流路径最大阻力值[28]。
2.2.6 地表径流空间自相关分析 空间自相关分析是对地理空间变量分布中相邻位置间的相关性检验的统计方法,以揭示区域化变量取值的空间分布特征,Moran’sI为相关系数,取值范围 [-1, 1],显著性水平下,Moran’sI为正时,表示观测值之间存在显著正相关,具有聚集性,当 Moran’sI为负时,表示观测值之间存在显著负相关,具有离散性[29-30],公式为:
(4)
(5)
地表径流形成的网络化路径是由“源”和“目标”的质量、阻力值和距离共同决定,城市下垫面地形因素和城市景观对于雨水径流的流动过程起着重要作用。试验对不同下垫面产生沿程阻力和局部阻力的因子进行赋值 (表1),通过AHP法进行权重分析(表2),得到景观阻力消费面(图2)。
图2 昆明市中心城区单因素和综合因素景观阻力消费面特征
表2 AHP昆明市中心城区阻力因子权重分析
在ArcGIS中grid模块下,MCR模型能够确定“源”和“目标”之间的最小消耗路径,符合地表径流的流动规律,即地表径流在城市中扩散的最易路径,构建出地表径流路径与网络结构。
由图3a可知,潜在径流数一共有1 274条,呈网络状分布,连通子汇水区上游—下游且具有明确方向。51个子汇水区中潜在地表径流数量n(条)和长度L(km)进行Pearson相关性分析,发现潜在径流数量与长度呈显著正相关(r= 0.911**,p=0.00<0.01,N=51),其回归方程为:y=1.657+0.004x(R2=0.84)。式中:y为潜在径流长(km);x为潜在径流数量(条)。回归直线对观测值拟合程度较好。核密度分析可直观展示潜在径流空间分布的特征,结果表明潜在径流核心位于研究区域中心位置,主要沿道路成带状分布 (图3b)。
图3 昆明中心城区潜在地表径流网络分布特征
研究发现流经汇流点的潜在径流数量(N)与实际测量值径流深度(p)的Pearson卡方值和对应显著性值(p)分别为:①119.883,0.205,>0.05;②86.779,0.063,>0.05;③86.779,0.113,>0.05(表3),即接受原假设,说明原假设模拟值与测量值的总体有效率相等,即两者之间的有效率差异不具有统计学意义,模拟值与观测值之间并无明显差异,存在相关性的假设成立[31]。
表3 潜在径流空间分布卡方检验
模拟条件下潜在径流在相邻和不相邻子汇水区域构成网络。通过重力模型 (gravity model)计算子汇水区间作用力矩阵 (表4)。
表4 基于重力模型的子汇水区域间相互作用矩阵
结果显示,子汇水区之间的相互作用强度存在差异,成本比指数反映网络的复杂程度,成本比指数越接近1网络结构越复杂,其公式为:cost ratio=1-L/d(L为径流数量;d为径流长度),研究区域成本比指数0.83成本较高、结构复杂。51个子汇水区域中的潜在径流数量随着子汇水区间重力强度的增加逐渐减少,但趋势逐渐下降。当重力系数大于19.16时,子汇水区间潜在径流出现断裂,因此,重力系数为19.16是子汇水区潜在径流连通的最大阈值。试验分别选取强度为0,19.16,30,50,100,500这6个重力系数,计算满足该重力条件的子汇水区特征,即重力系数越大,子汇水区间连接越强,越易发生汇流。结果表明,研究区域径流路径单元密度随重力系数的增加而逐渐降低,不同子汇水区径流形成存在显著差异,子汇水区之间相互作用越强、距离越近,越容易为潜在径流的产生和联通提供可能,因此可以直观反映实际不同降雨条件下子汇水区地表径流汇流的先后顺序。
空间自相关分析的6个场景中,Moran’sI都大于0,且Z>2.58,p<0.01,说明潜在径流路径单元密度分布与汇水区地理分布显著正相关呈聚集性,并由中心区域向边缘呈从高到低递减趋势;随着重力系数增加,不满足重力条件的子汇水区间的潜在径流消失,子汇水区潜在径流路径单元密度Moran’sI系数呈现下降趋势 (表5)。
重力系数增加会影响径流地理空间分布,潜在径流逐渐密度减少 (表5):①潜在径流路径单元密度的分布呈显著多核心分布、多级序列级分布规律,并随着子汇水区间的重力增加,呈现出由多核心向单一核心变化趋势;②潜在径流随着子汇水区域之间关系减弱,呈现出由外向内收缩趋势,即在径流形成初期,首先满足重力系数较高的区域,同时低密度径流区域聚集也呈下降趋势(图4)。
表5 潜在径流分布特征分析
注:a为全部潜在径流路径单元密度分布; b为重力系数为19.16时在径流路径单元密度分布; c为重力系数为30时潜在径流路径单元密度分布; d为重力系数为50时潜在径流路径单元密度分布; e为重力系数为100时潜在径流路径单元密度分布; f为重力系数为500时潜在径流路径单元密度分布。
昆明市中心城区土地利用呈现层次性递变特征,城市不透水表面比例由城市中心向外逐步递减,呈同心圆式发展趋势,道路为地表径流提供了扩散通道,试验结果表明昆明中心城区是道路汇流型区域。其中,昆明中心城区地形与城市景观不同程度的影响了潜在地表径流的汇—散关系,昆明中心城区的地形差异较大,五华山、圆通山、虹山等分布,导致局部山地区域汇流减少,而低洼区域汇流增加,引起潜在地表径流的分配不均;而昆明中心城区的城市公园,翠湖公园、圆通山公园等公共绿地空间,以及城市河网水系,盘龙江、大观河等为城市区域提供了地表阻力以减少地表径流的汇流,有效地减缓了局部汇流的现象。
基于高清遥感影像通过地形要素、高程要素、景观要素等多因素分析,并采取AHP法划分阻力值等级,构建景观阻力消费面,形成中心城区二维阻力空间,为模拟试验提供了径流空间分布模拟的基础条件。通过卡方检验发现,潜在径流与实际观测径流深度无明显差异,存在相关性。潜在径流在阻力空间中呈线状分布,连通上游—下游子汇水区且具有明确方向。研究区内域潜在径流呈显著聚集性,主要扩散路径为道路,表现为高密度到低密度由中心向四周逐渐递减的趋势,伴随子汇水区间的重力逐渐降低,能够满足汇流条件的径流逐渐减少,径流路径单元密度由外向内收缩,由多核心分布向单一核心变化。其特征表现为3方面:①潜在径流的空间分布与城市化发展趋势一致,由城市中心向四周逐渐降低,其中径流分布与道路、建筑等不透水景观密度变化趋势一致,而与公园绿地、草地、水体则相反;②子汇水区之间的重力强度决定了潜在径流汇流程度,可推测现实过程中,城市暴雨不同阶段,随着径流量的增加,径流克服子汇水区之间重力的能力增强,反映出地表径流形成过程的先后顺序、聚集程度、有效路径和分布规律;③通过潜在径流的线性结构与方向性能有效的捕捉出研究区域潜在径流的入水口、出水口,明确区域内不同城市用地的进水口、排水口的确切位置,能够较直观反映出径流形成与城市综合环境的相互关系和交流通道,为有效城市引流提供相对科学的依据。
高清遥感技术的发展在大比例尺(10~100 km2)测量任务中表现快速高效、精细准确、成本低等优势,并在地图测绘、应急救灾、国土监测等方面得到广泛应用,本次研究选用Spot-5高清卫星正射影像和WorldView 2 DSM二维栅格影像为参考依据得到可靠的数据基础;城市建设用地面积和结构有严格的建设标准,较自然地貌相对单一,结合MCR模型与空间自相关分析方法,有利于分析城市潜在径流的地理空间分布特征和趋势,能够较全面、直观地展示昆明中心城区地表径流分布特征,并直观、迅速、准确地确定潜在灾害发生的区域与汇流入口节点与出口节点的位置,为有效城市导流、疏浚、提供科学参考,从而有的放矢的进行“生态修复、城市修补”和雨洪资源的循环利用,减少城市防涝的成本。
本研究还处于城市地表景观与城市径流时空分配的持续探究阶段,对系统模型研究还需要进一步的深入探索。研究团队将不断改进研究方法,构建系统化耦合模型和应用系统并不断提高计算精确度,将研究成果运用到城市雨水资源实际管理中。