工程科学与技术   2018, Vol. 50 Issue (3): 158-168
汶川县震后潜在泥石流危险性评价研究
韩用顺, 王晶, 吴淼, 李乐林, 张东水, 张金平     
湖南科技大学 资源环境与安全工程学院,湖南 湘潭 411201
基金项目: 交通运输部科技计划资助项目(20153161T906);岩土力学与工程安全湖南省重点实验室开放基金资助项目(16GES06)
摘要: 对震后潜在泥石流灾害进行危险性评价是灾害预测预报与防治规划的重要内容和技术基础,也是目前泥石流领域研究的薄弱环节。以汶川县为典型研究区,通过资料收集、野外考察和遥感调查,遴选出震后泥石流频次、流量和地貌信息熵3个关键指标,改进震后泥石流流量计算方法,建立DNWT模型动态监测并快速提取震后泥石流物源信息;基于此,构建评价指标体系,以小流域为评价单元,采用GIS技术和灰色系统理论与方法,建立评价模型,进行研究区震后潜在泥石流危险性评价与区划。结果表明:1)提出DNWT震后泥石流物源提取模型,能实现震后泥石流判识和崩塌滑坡及其松散堆积体等泥石流物源信息动态监测与自动提取;2)结合野外调查和泥石流物质量计算数据,采用配方法,修正了震后泥石流流量计算公式,与震前相比,堵塞系数增加1.72~3.46倍,流量呈指数增大;3)基于对研究进行小流域划分,引入地貌信息熵,简化震后泥石流危险性评价参数,实现了对子流域地貌发育阶段和沟谷危险程度的判断;4)研究区大部分泥石流沟处于中度及以上危险区,具有极大的潜在危险性,其中高、极高危险性泥石流沟89条,占泥石流总数的39.2%,主要沿岷江干流和渔子溪河谷两岸分布,为防灾减灾的重点区,且资源环境承载力低,不宜进行规模以上开发建设,需加强灾害监测预警与综合防治。研究结果可为震后泥石流遥感动态监测、定量风险评估和防灾减灾管理提供依据与参考。
关键词: 震后泥石流    危险性评价    物源信息提取    流量计算    汶川县    
Hazard Assessment on Potential Post-earthquake Debris Flows in Wenchuan County
HAN Yongshun, WANG Jing, WU Miao, LI Lelin, ZHANG Dongshui, ZHANG Jinping     
School of Resource Environment and Safety Eng., Hunan Univ. of Sci. and Technol., Xiangtan 411201, China
Abstract: The hazard assessment on potential post-earthquake debris flows is an important task and technical basis for disaster prediction and prevention,while related research is limited.Wenchuan was taken as a typical study area.Through data collection,field investigation and remote sensing monitoring,the frequency,the flow quantity and the geomorphic information entropy of post-earthquake debris flow were selected as three key indicators for hazard assessment.Correspondingly,the discharge calculation method of post-earthquake debris flow was improved,the DNWT model was established for dynamic monitoring and quick extraction of the sources of post-earthquake debris flow,and the evaluation index systemfor debris flow hazard was set up.By means of the GIS and the grey system method,the model for evaluating the hazard of post-earthquake debris flows was established,sub-catchments were extracted and selected as hazard assessment units,and the hazard assessment on potential post-earthquake debris flows in the investigated area was conducted.The results show that:1) the proposed DNWT model for extracting debris flow sources can realize the remote sensing identification,dynamic monitoring and automatic extraction of post-earthquake debris flow sources including rockfalls,landslides and their loose accumulations;2) by combining the field investigation with debris-flow discharge calculation data,the discharge calculation formula of debris flow is corrected.The blocking factor of post-earthquake debris flow can increase 1.72-3.46 times,compare with the case before the earthquake,and its discharge increases exponentially;3) based on sub-catchments,geomorphic information entropy is introduced,and meanwhile,it simplifies parameters and procedures for hazard assessment on post-earthquake debris flows and can evaluate the development stage of landforms and the hazard degree of each catchment;4) most of the debris flows in the studied area are in the moderate,high and very high hazard zones and have great potential risks.Among them,there are 89 debris flows in high and very high hazard zones,which are distributed along the main stream of Minjiang River and the valley of Yuzi River,and these areas are crucial for debris flow prevention and mitigation.The research results are consistent with the actual situation and can provide a basis and reference for remote sensing dynamic monitoring,quantitative risk assessment and prevention and mitigation of potential post-earthquake debris flows.
Key words: post-earthquake debris flow    hazard assessment    source information extraction model    corrected discharge calculation method    Wenchuan county    

汶川地震引发的大量崩塌滑坡等次生山地灾害为泥石流形成提供了丰富物源,降低了启动条件[1],改变了震后泥石流灾害的形成、发育、活动与致灾方式,导致其规模、频率、危害等显著增加。位于极震区的汶川县震后连年遭受泥石流灾害:2010年8月13日,映秀镇红椿沟泥石流冲出物质约7 $\times $ 105 m3,堵断岷江并造成72人死亡、8 000余人被迫转移[2];2013年7月10日,都汶公路全线爆发大规模泥石流灾害,堵塞岷江并形成24座堰塞湖,严重损毁道路、桥梁和隧道;G213线震后平均每年发生2~3次大规模泥石流灾害,阻断交通5~10 d,4次被迫重建;震后重建的草坡乡被迫整体异地搬迁并再次重建,等等。震后泥石流规模、数量与活动方式已远超出以往,在形成机制、暴发规模和危险性等方面与震前和非震区泥石流都存在显著差别。由于对震后泥石流特征和规模等缺少科学认识和深入研究,以往的泥石流规模计算和灾害评估方法不能满足震后防灾减灾要求。

国际上,19世纪后半期,俄国开始对泥石流沟的危险性问题进行研究。1930年,瑞典学者Eldeen[3]通过危险性区划图开启泥石流危险性的定量化研究;1977年,日本学者Adachi等[4]首次提出泥石流危险度的判定;1981年,美国Hollingsworth等[5]采用打分和因子叠加方法实现了泥石流危险度评价。1990年以后,随着计算机技术的高速发展,3S技术逐步被引入泥石流危险性评价研究[67],泥石流的空间分析和数值模拟有了很大进展,灰色关联度、模糊数学、神经网络、粗糙集和支持向量机等非线性数学方法也被引入泥石流危险性评价[810],并向基于物理过程模型的定量数值计算与模拟分析趋势发展[1112]。国内,泥石流危险性评价于1980年起全面起步,王礼先、刘希林等[1314]提出泥石流危险度的定量判定方法。1990年以来,随着3S技术、非线性数学模型、人工智能和数值模拟方法的引入,泥石流危险性评价由定性化向定量化发展[1517]。2000年以后,国内学者陆续开展地震泥石流危险性评价研究[18],并在震后泥石流形成机理、运动过程和危险性评价技术与方法方面开展工作,追踪国际前沿,取得了一定成果[1921]。总体而言,有关泥石流危险性评价的原理和方法已趋成熟,并在实践应用中逐步完善和改进,精度和实用性也越来越高。但目前对震后潜在泥石流危险性评价研究仍相对较少,且已有的评价方法需要较多的输入参数,没有考虑震后泥石流流量放大、频率增加等活动特征变化,缺乏对震后泥石流演化过程的动态监测。因此,需要寻找新的科学有效的技术方法进行震后潜在泥石流危险性评价。

以汶川县为典型研究区,通过资料收集、野外考察和遥感调查,提出DNWT震后泥石流物源提取模型,分析震后泥石流发育特征、发生频率和物源规模,改进震后泥石流流量计算方法,引入地貌信息熵,遴选震后泥石流危险性评价关键因子,构建评价指标体系;并以小流域为评价单元,采用GIS技术和灰色系统模型,建立震后潜在泥石流危险性评价方法,进行研究区震后泥石流危险性评价与区划。研究成果可为震后泥石流风险管理与防治规划提供科学依据与技术支持。

1 研究区与数据源 1.1 研究区概况

汶川县位于四川省阿坝藏族羌族自治州东南部,E102°51′~103°44′和N30°15′~31°43′之间,总面积4 084 km2,境内有G213、G317、S303和都汶高速公路,是川西重要交通要道和进藏国防通道。该县位于邛崃山系和龙门山系之间,属典型高山峡谷地貌,地势西北高、东南低,海拔780~6 250 m,沟谷纵横、切割强烈,地形坡度和相对高差大。区内新构造运动和地震活动强烈,地质构造复杂,地层发育较完整,工程地质岩组特征差异大,对震后次生山地灾害形成与发育具有不同的影响作用。研究区属暖温带大陆性半干旱季风气候区,南湿北旱,垂直分带明显;银杏乡苏坡店以南为川西多雨中心区,年均降雨量1 285.1 mm,以北属于岷江上游半干旱河谷区,年均降雨量为526.3 mm。该区属岷江水系,水资源与水能资源丰富,植被和土壤垂直地带性明显,以林地和草地为主。

汶川县位于汶川地震极震区,独特的自然地理环境和复杂的气候条件为滑坡泥石流等山地灾害的形成发育提供了良好条件,震前即为山地灾害易发区。强震严重损伤了岩土体结构,诱发大量崩塌、滑坡、泥石流和不稳定斜坡,其中都汶公路映秀–草坡段,崩塌3 863处、滑坡18处、泥石流沟74条,固体松散物质量1亿余方[22]。汶川县独特的自然条件、强震动力作用和震后强降雨导致震后泥石流频发且危害严重,并具有明显的潜在性、链生性和隐蔽性,是进行震后泥石流危险性评价的理想区域。

1.2 数据源及其处理 1.2.1 基础数据

震后泥石流是地形地貌、地质地震、气象水文、土地利用和人类活动等多种因素共同作用的结果,采集并处理相关数据资料,研究并建立震后泥石流综合数据库是灾害评估的技术基础。地形地貌既是泥石流活动场所,也为泥石流形成发育提供能量与动力条件,其数据来源1∶10 000DLG和16 mDEM,采用ArcGIS10.2提取相对高差、沟床比降等地形因子及居民地、道路等基础地理信息。地质环境是泥石流形成发育的内营力条件,控制着区域地貌形态和地表物质演化与迁移,为泥石流活动提供物质基础和内动力,其数据来源于1∶50 000数字地质图,用于提取地层岩性、构造等因子。降雨是震后泥石流主要激发因素,其数据来源于汶川县震后布设的5个自动气象站、10个雨量计和中国气象科学数据库,通过地形与下垫面耦合,建立高精度降雨插值模型,生成汶川县各类降雨指标及其分布图。水系、径流量、地下水等水文要素是泥石流形成和地貌塑造的重要动力条件,其数据来源于1∶10 000DLG和1∶200 000水文地质图,通过水文模型、格式转换和数字制图生成相关图件。土地利用数据主要通过2008—2016年的TM/ETM+/OLI和ADS40影像遥感解译获取,并结合野外调查和2010年的1∶100 000土地利用图融合与校正。

1.2.2 遥感影像

遥感影像是泥石流灾害研究的重要信息源和数据更新手段,采用2008—2013年的航空ADS40影像及2008—2016年的Quickbird、Spot和TM/ETM+/OLI等多源遥感影像,进行崩塌、滑坡、泥石流等次生山地灾害和土地利用信息提取与动态监测;航空ADS40影像来源于中国科学院,TM/ETM+/OLI影像免费下载,Quickbird和Spot影像则购置。利用ENVI5.3和ArcGIS10.2软件,经影像校正、镶嵌、裁剪、增强等处理,消除影像坐标误差,量化部分增益和偏置的参数,增强目标地物信息及其影像辨识度;再结合1∶10 000DLG和野外调查数据,进行震前震后山地灾害、土地利用及相关专题信息的遥感判识与信息提取,生成并制作相应的灾害专题数据和基础地理信息。

1.2.3 综合数据库

专题数据主要为崩塌、滑坡、泥石流等山地灾害震前震后的类型、分布、规模、频率等灾害发育与致灾特征,是震后泥石流危险性评价的重要因子。其数据来源于研究区2008—2017年30次野外考察、连续9年的多源遥感调查、相关研究报告与文献资料。集成获取的基础数据、多源遥感数据、灾害专题数据和孕灾环境数据,采用ArcSDE和SQL Server2016数据库管理系统,建立研究区震后泥石流综合数据库和基础信息平台。

2 关键指标及其处理 2.1 震后泥石流频次

泥石流频次表征泥石流活跃程度,是震后泥石流危险性评价重要指标之一。作者建立DNWT震后泥石流物源提取模型,对震后泥石流活动进行遥感动态监测,获取并计算震后泥石流频次参数。

2.1.1 随机森林分类算法

灾害遥感判识及其影像分类是震后泥石流遥感信息提取的前提与基础。与常用的最小距离、最大似然、神经网络、决策树和支持向量机等方法相比,随机森林算法是一种更有效的分类器集成算法,具有更高的分类精度和运算速度,无需先验知识,且对噪声数据具有很好的鲁棒性[23]。其基本流程如下(图1):

1)通过Bootstrap技术随机抽取训练样本集D $K$ 个训练子集,未被抽取的数据组成OOB(out-of-bag)数据。

2) $K$ 个训练子集分别构建决策树,每棵决策树节点从 $T$ 个特征中随机抽取 $M$ 个特征( $M$ $T$ ),按照节点不纯洁度最小原则分支生长,节点分裂标准采用基尼指数(式(1))和信息增益两种方法(式(2)):

$Gini(D) = 1 - \sum\limits_{i = 1}^m {\mathop p\nolimits_i^{\rm{2}} } $ (1)
$Entropy(D) = - \sum\limits_{{i} = 1}^m {\mathop p\nolimits_{i}^{} }{\rm lb}\left( {\mathop p\nolimits_i } \right)$ (2)

式中, $m$ 表示训练样本集 $D$ 包含的类别个数, $ p_i$ 表示训练数据属于类别 $C_i$ 的概率。

3)让每棵树最大限度生长,直至所有树生成随机森林, $K$ 个训练子集得到 $K$ 种分类结果,并通过投票表决确定最终分类结果。

图1 随机森林算法流程 Fig. 1 Flow chart of random forest algorithm

2.1.2 震后泥石流物源信息提取

不同地物的波谱曲线形态各异,研究地物光谱特征可进行震后泥石流物源判识与提取。以2009年6月3日的TM数据为例,影像中地物主要为植被、崩塌滑坡泥石流等次生地质灾害、居民地、水系等,另有少量的云、积雪和阴影。根据2008—2017年都汶公路沿线野外考察和高分辨率ADS40遥感解译结果,针对不同地物分别选择一定数量的样本数据,统计TM影像7个波段、缨帽变换3个波段及主成分分析的第1成分(PC1)均值,制作各地物光谱响应曲线(图2)。

图2 汶川县2009年6月3日TM影像典型地物光谱响应曲线 Fig. 2 TM image typical object spectral response curves of Wenchuan county on June 3,2009

1)震后泥石流物源指数模型(DFMSI)

图2可知,TM5和TM6波段上,震后崩塌滑坡及其松散堆积物光谱曲线值的斜率与其他地物差异大,即震后泥石流物源光谱曲线呈下降趋势,其他地物光谱曲线呈上升趋势。据此,建立震后泥石流物源指数模型(式(3)):

$DFMSI = \left( {TM5 - TM6} \right)/\left( {TM5 + TM6} \right)$ (3)

设置Cart树的数量为100,研究区共有各类目标地物样本点6 115个。由于水体面积较小,仅选取715个样本点,其他6类地物各选900个样本;选择4 281个样本点(占样本总量的70%)作为训练样本,剩余1 834个样本点为测试样本。计算可得2009年6月3日震后泥石流物源信息提取结果(图3)。

图3 2009年6月3日TM影像震后泥石流物源信息提取 Fig. 3 Source information extraction of post-earthquake debris flows through TM on June 3,2009

2)震后泥石流物源提取模型(DNWT)

震后泥石流物源区植被损毁严重且松散物质含水量高,可选取归一化植被指数(NDVI)和土壤湿度指数(WET)进行物源区定量遥感分析。纹理特征是区分震后泥石流物源区不同地物的另一个重要特征,考虑变量本质特点、图像灰度指标差和地物纹理特征,选取相关性、熵、对比度、2阶矩等4个纹理特征,通过layer stacking算法,叠加震后泥石流物源关键特征,结合随机森林算法,建立DNWT震后泥石流物源信息提取模型与算法流程(图4),自动提取研究区2008年12月1日、2010年3月8日、2011年8月28日、2013年8月17日、2014年6月17日、2015年2月12日和2016年3月18日的震后泥石流物源信息,并分析其类型、分布、储量和演化规律,获取2008—2016年震后泥石流频次指标数据。

图4 震后泥石流物源信息提取模型与算法流程 Fig. 4 Source information extraction model of post-earthquake debris flows and its algorithm process

2.1.3 精度分析

利用30%的测试样本对DNWT震后泥石流物源提取模型的提取结果进行精度评价,计算各期影像提取结果的生产精度、用户精度和Kappa系数(表1)。

表1中,生产精度和用户精度分别超过88%和90%,且Kappa系数大于0.88。这表明建立的DNWT模型具有较高精度和可操作性,可获取震后泥石流频次等关键指标数据,并满足区域次生山地灾害遥感快速调查和动态监测需要。

表1 基于DNWT模型的震后泥石流物源信息提取结果精度分析 Tab. 1 Accuracy analysis of DNWT-based source information extraction results on post-earthquake debris flows in study area

2.2 震后泥石流流量

震后泥石流规模预测既无成熟的理论,也无可借鉴的技术方法,但可根据已经发生的泥石流状况结合经验公式进行预测。一次泥石流最大冲出量不仅反映泥石流的规模、强度与性质,还决定着防治工程的类型与结构,是判定泥石流规模和危险性最直接的指标,故用该指标代替泥石流规模。研究并通过实地调查和文献资料,获取震后典型泥石流事件(表2)中的泥石流流量(实际值 $Q_{\text F}$ )(表3);结合震区降雨数据和《四川省水文手册》,采用配方法计算上述泥石流事件中的泥石流流量(计算值 $Q_{\text C}$ [24];对比泥石流流量实际值与计算值,可得震后泥石流流量的放大系数( $C_{\text M}$ ),并据此修正配方法计算公式(4);采用修正后的震后泥石流流量模型计算各泥石流沟的一次最大冲出量,并将其作为震后泥石流危险性评价的主要因子。

$ \left\{ \begin{aligned}& \mathop Q\nolimits_{\rm{C}} = \mathop Q\nolimits_{\rm{B}} \times \left( 1 \right. + \left. \phi \right) \times q\text{,} \\& \mathop Q\nolimits_{\rm{B}} = 0.278\left( {\frac{{\mathop S\nolimits_{p} }}{{\mathop \tau \nolimits^n }} - \mu } \right) \times F\text{,} \\& \phi = \frac{{\mathop \gamma \nolimits_{\rm{C}} - 1}}{{\mathop \gamma \nolimits_{\rm{H}} - \mathop \gamma \nolimits_{\rm{C}} }}\text{,} \\& \mathop S\nolimits_{p} = \mathop H\nolimits_{p} \times \mathop t\nolimits^{n - 1} \text{,} \\& \tau = \mathop \tau \nolimits_0 \times \mathop \varphi \nolimits^{\frac{1}{{4 - n}}} \text{,} \\& \mu = 6\mathop F\nolimits^{ - 0.19} \text{,} \\& \mathop \tau \nolimits_0 = {\left( {0.383\frac{\theta }{m}\mathop { \times S}\nolimits_{p}^{ - 0.25} } \right)^{\frac{4}{{4 - n}}}}\text{,} \\& \theta = \frac{L}{{\mathop J\nolimits^{\frac{1}{3}} - \mathop F\nolimits^{\frac{1}{4}} }}\text{,} \\& m = 0.318\mathop \theta \nolimits^{0.204} \end{aligned} \right.$ (4)

式中: $Q_{\text C}$ 为泥石流洪峰值流量,m3/s; $Q_{\text B}$ 为设计频率下洪水洪峰流量,m3/s; $\phi $ 为泥沙修正系数,即单位泥石流体积中固体物质体积与水体体积的比值; $\gamma_{C}$ 为泥石流容重; $\gamma_{\text H}$ 为砂石等固体物质容重; $q$ 为泥石流沟道堵塞系数,其值可查经验表; $S_{P}$ 为雨力,即设计频率 $P$ 下的最大 $t$ 小时暴雨量,mm/h; $\tau $ 为汇流时间,h; $\varphi $ 为洪水径流系数,查水文手册可得; $n$ 为暴雨参数,查水文手册可得; $\mu $ 为产流参数,mm/h; $F$ 为流域面积,km2 ${\theta}$ 为流域特征系数; $J$ 为沟床比降, $L$ 为主沟长度,km; $m$ 为汇流参数。

表2 震后典型泥石流事件的历时与流量调查统计 Tab. 2 Statistical of duration and discharge investigations on typical post-earthquake debris flows

表3 震后典型泥石流事件的流量计算值及其放大系数 Tab. 3 Discharge calculation and amplification coefficient of typical post-earthquake debris flows

利用1∶10 000DEM和多源遥感影像,提取表3中典型泥石流沟相关参数;结合调查数据和《四川省水文手册》,反演震后泥石流沟道堵塞系数的放大值( $C_{\text M}$ )(表3图5)。

图5 震后泥石流流量实测值与计算值关系 Fig. 5 Relationship between measured and calculated values of post-earthquake debris flows

表3图5可知:强震极大改变了震后泥石流形成条件,使得震后泥石流物源增加,流量增大;与震前相比,震后泥石流沟道堵塞系数增加1.72~3.46倍,流量呈指数增大,满足 $ Q_{\rm F}86.02{\rm e}^{0.007Q_{\rm C}}$ 。因此,修正后的震后泥石流流量计算公式:

$\mathop Q\nolimits_{\rm{F}} = 86.02{\rm{\cdot}}\mathop {\rm{e}}\nolimits^{0.007\mathop Q_{\rm{C}} } $ (5)

$\mathop Q\nolimits_{\rm{F}} = \mathop Q\nolimits_{\rm{B}} \times \left( 1 \right. + \left. \phi \right) \times q \times \mathop {C}\nolimits_{\rm{M}} $ (6)

采用式(5)计算研究区各沟道震后泥石流流量,并采用式(6)进行检核。

2.3 地貌信息熵

泥石流强烈的侵蚀切割作用加速了流域地貌形态的演化,是其区别于其他地质灾害的显著特征;同时,沟谷的地貌条件制约着泥石流的形成与运动,影响着泥石流的发生和规模。因此,通过对流域地貌侵蚀形态的定量计算,可在一定程度上反映泥石流沟谷发育阶段,判别其危险性[25]。地貌信息熵作为判定流域地貌侵蚀发育程度及其演化阶段的重要量化指标之一,具体计算过程如下:

1)拟合面积–高差函数

量算每条等高线以上的面积( $a$ )和每条等高线与流域最低点的高差( $h$ ),分别与流域的总面积( $A$ )和流域的高差( $\Delta H$ )建立函数关系式:

$y = f(x)$ (7)

式中, $x$ = $a$ / $A$ $y$ = $h$ / $H$

2)计算Strahler面积–高程积分:

${S} = \int_0^1 {f\left( x \right)} {\rm d}x$ (8)

3)计算地貌信息熵:

$H = S - {\rm ln} S - 1$ (9)

式中, $H$ 为地貌信息熵值, $S$ 为Strahler面积–高程积分值。

利用GIS空间分析功能,对1∶10 000DLG按10 m间隔重分类,再与子流域图层进行叠置分析,得到各子流域内一系列( $x_i$ $y_i$ )值;通过Matlab13.0对面积–高程函数 $f(x)$ 进行拟合;根据拟合函数和式(8)~(9)计算各子流域的 $S$ $H$ 值(图6)。

图6 研究区各子流域地貌信息熵 Fig. 6 Geomorphic information entropy of each sub-catchment in study area

3 研究方法 3.1 评价方法

灾害评价过程中,不可能也没必要采集到所有的指标数据,而已采集的数据所包含的内容与信息时常不完备。在信息不完备或不确定的情况下,震后泥石流灾害评价是一个“灰色”黑箱过程,故可采用灰色系统理论与方法进行震后泥石流危险性评价。灰色系统是基于信息非完备性原则建立的非线性数学模型,能够充分利用系统已有信息,考虑了评价标准和指标的模糊性、灰色性、不确定性和观测资料的随机性,只需对已掌握的部分信息进行合理加工处理,即可较为系统动态地对结果做出科学描述和正确评价,从而在很大程度上解决危险性评价指标权重确定中人为因素造成的偏差问题。因此,本文采用灰色系统原理与方法进行各评价因子的权重计算分析,其主要过程如下:

1)归一化处理评价数据:

$\begin{aligned}[b]{Y}_j^\prime(i) = \frac{{{{Y}_j}(i) - \min {{Y}_j}(i)}}{{\max {{Y}_j}(i) - \min {{Y}_j}(i)}};\\{k = 1,2,\cdots,n;j = 1,2,\cdots,m}\end{aligned}$ (10)

式中, $Y_j(i)$ 为第 $i$ 个评价单元中第 $j$ 个评价指标的真实值; $ {Y}_j^\prime(i)$ ${Y}_j(i)$ 归一化值。

2)计算参考序列与比较序列的绝对差:

${\varDelta _j}(i) = \left| {{Y}_j^\prime(i) - {Y}_j^\prime({\rm{0}})} \right|$ (11)

3)计算绝对差的最大值和最小值:

$\left\{ \begin{aligned}& \varDelta_{{\rm{min}}} = \min (\min {\varDelta _j}(i))\text{,}\\& \varDelta_{\max} = \max (\max {\varDelta _j}(i))\end{aligned} \right.$ (12)

4)计算关联系数:

${\phi _{j}}(i) = \frac{{\varDelta_{{\rm{min}}} + \partial \cdot \varDelta_{{\rm{max}}}}}{{{\varDelta _j}{(i)} + \partial \cdot \varDelta_{{\rm{max}}}}}$ (13)

式中, $ \partial $ 为比例常数。

5)计算关联度:

${{R}_{j}}({Y_{\rm{0}}},{Y_{j}}) = \frac{1}{n}\sum\limits_{i = 1}^n {{\phi _{j}}(i)} $ (14)

6)计算各评价因子权重:

${{w}_{j}} = {M_j}/\sum\limits_{i = 1}^n {{M_j}} $ (15)

式中, $M_{j}$ $R_{j}(Y_{0},Y_{ j})$ 的权数值。

7)计算危险性指数:

为便于确定各评价单元的震后泥石流危险程度及其数值大小,引入危险性指数 $H$ [26],可表达为各评价单元危险性评价因子值与其权重值的函数:

${{H}_i} = \sum\limits_{j = 1}^m {{w_j} \cdot {{{X}}_{{i}j}};\begin{array}{*{20}{c}}{\begin{array}{*{20}{c}}{i = 1,2,\cdots,n;j = 1,2,\cdots,m}\end{array}}\end{array}} $ (16)

式中, $X_{ij}$ 为第 $j$ 个评价因子在第 $i$ 个评价单元的取值。

3.2 评价单元

泥石流的形成、演化与治理通常以小流域为单元,故本文以子流域作为危险性评价单元。利用ArcGIS10.2水文分析功能和1∶10 000DEM,通过洼地填充、汇流方向及其累积量计算、河网提取与分级等处理,自动提取汶川县各水系的子流域;将提取结果与研究区遥感影像和1∶10 000DLG中的水系进行对比分析、拓扑检查和人工编辑,最终将研究区划分为227个子流域。

4 结果与分析 4.1 评价指标体系

评价指标的选取和量化是危险性评价的基础。以震后潜在泥石流形成条件为主导,结合其分布规律与发育特征,从物源、地形地貌、地质地震、气象水文和人类活动等5个方面,选取震后泥石流流量 $X_{1}$ 、震后泥石流频度 $X_{2}$ 、地貌信息熵 $X_{3}$ 、平均坡度 $X_{4}$ 、至断裂带距离 $X_{5}$ 、地震烈度 $X_{6}$ 、1 h雨强 $X_{7}$ 、暴雨日数 $X_{8}$ 、切割密度 $X_{9}$ 和植被覆盖度 $X_{10}$ 等10个因子。其中: $X_{1}$ $X_{2}$ 为泥石流活动特征指标,反映了震后泥石流的规模、频率与强度,为主要评价因子; $X_{3}$ $X_{4}$ 为地形地貌指标,反映了泥石流发育条件; $X_{5}$ $X_{6}$ 为地质指标,反映了强震作用与影响,表征震后泥石流物源条件; $X_{7}$ $X_{8}$ 为气象指标,为震后泥石流激发因子; $X_{9}$ $X_{10}$ 分别反映震后泥石流的水动力条件和诱发因素。

为便于各评价指标的比较与计算,根据各评价因子形成因素贡献率和震后泥石流发育特征,将其统一划分为5个不同的等级区间,分别赋值为5、4、3、2、1,对应震后泥石流的极高、高、中等、低和极低等5个危险等级。基于此,建立研究区震后泥石流危险性评价指标体系见表4

表4 研究区震后泥石流危险性评价指标及其分级赋值 Tab. 4 Hazard assessment indicators and their hierarchical assignment of post-earthquake debris flows in study area

4.2 评价结果及分析

将评价因子图层分别与评价单元进行空间叠置分析,得到各评价单元的每个评价因子属性值,并按表4的分级赋值规则,计算各评价单元中每个评价因子的转换值,如表5所示。

表5 研究区部分评价单元的评价因子计算结果 Tab. 5 Conversion values of each assessment factor in parts of assessment units

采用双系列关联度分析法,选取震后泥石流流量 $X_{1}$ 和频次 $X_{2}$ 两个主要因子为参考序列,将其他评价因子作为比较序列,分别与 $X_{1}$ $X_{2}$ 进行灰色关联度分析。

各评价因子权重值见表6,各评价单元的危险性指数见表7

表6 研究区震后泥石流危险性评价因子权重值 Tab. 6 Weights of each hazard assessment factor on post-earthquake debris flows in study area

表7 研究区震后泥石流危险性指数计算结果 Tab. 7 Calculation results on hazard index to post-earthquake debris flows in study area

采用聚类方法,按照区内程度相似、区间差异最大、地貌单元完整性等原则,进行研究区震后潜在泥石流危险性分级区划与制图。经过试验和实地调查,将研究区震后潜在泥石流分为极高危险区、高危险区、中度危险区、低危险区和极低危险区5个等级(表8图7)。

表8 研究区震后潜在泥石流危险性评估结果分级 Tab. 8 Classification of hazard assessment values of each unit on post-earthquake debris flows in study area

图7 汶川县震后潜在泥石流危险性评价结果 Fig. 7 Hazard assessment on post-Earth-quake potential debris flows in study area

结果表明:汶川县大部分泥石流沟震后处于中度及以上危险区,具有极大的潜在危险性。中度、高、极高危险性泥石流沟共137条,占研究区小流域总数的60.4%;高、极高危险性泥石流沟主要沿岷江干流和渔子溪河谷两岸分布;极低、低、中度、高、极高危险区由外向内大致呈环状分布,受强震、地形和降雨控制作用明显。

1)极高危险区泥石流沟45条,占研究区泥石流总数的19.8%,主要分布在银杏乡、映秀镇和耿达乡。该区位于川西暴雨中心,靠近震中和映秀–北川主断裂带,岩土体严重破碎,大型群发性崩塌、滑坡、泥石流等次生山地灾害广泛发育,灾后恢复重建等人类工程活动强烈,潜在泥石流极易发生且严重威胁灾区工程与公共安全,是研究区人口与社会经济聚集区、恢复重建的高危区和防灾减灾的重点区。该区潜在泥石流危险性极高,破坏力巨大,通常不宜进行规模以上开发建设,人口密度和社会经济规模需控制在资源环境的合理承载能力之内。

2)高危险区共有泥石流沟44条,占研究区泥石流总数的19.4%,涉及三江、草坡、映秀、耿达、绵虒、银杏、水磨、漩口8个乡镇,主要分布在三江、耿达、绵虒和草坡4个乡镇。该区围绕极高危险区,大部分仍位于川西暴雨中心和岷江右岸,区内岩体较破碎,切割强烈,崩塌、滑坡等次生山地灾害点多面广、规模大,震后泥石流具有很高的危险性和潜在性,也是研究区恢复重建和防灾减灾的重点区域。

3)中度危险区共有泥石流沟48条,占研究区泥石流总数的21.2%,主要分布在三江、卧龙、草坡、耿达4个乡镇。该区毗邻极高、高危险区,大部分沿渔子溪和草坡河两岸分布,区内植被覆盖较好,地震影响程度中等,岩土体相对完整,崩塌、滑坡、泥石流灾害以中型为主,总体不甚密集,但峡谷和高陡边坡地带也常有大型崩塌、滑坡等次生山地灾害发生。该区潜在泥石流危险性及其破坏能力中等,耕地面积较大且具有一定肥力,资源环境承载力较好,经勘察、论证或通过山地灾害防治能保障其安全,可开发为人口密度较高的农林牧业区,且能建设中等规模及以下的工矿企业和村镇驻地等。

4)低危险区共有泥石流沟61条,占研究区泥石流总数的26.8%,主要分布在卧龙、草坡、龙溪和耿达4个乡镇。该区毗邻中度危险区,大部分位于高海拔地带,距离震中、主断裂带、干线交通相对较远,各沟道内岩土体较为完整,植被覆盖率高,崩塌、滑坡等不良地质现象较少,且主要受微地貌和小构造控制,以中小型为主。该区潜在泥石流危险性低,破坏能力较小,资源环境承载力较高,不仅可以进行现代化农业建设,也可开展较大规模资源开发和经济建设。

5)极低危险区共有沟道29条,占研究区泥石流总数的12.8%,主要分布在卧龙、龙溪和耿达3个乡镇。该区毗邻低、中度危险区,主要为高山或平缓地带,不具备崩塌、滑坡、泥石流等次生山地灾害形成条件,但受人类工程建设、经济活动或流水侵蚀作用影响,少数流域偶有崩塌、滑坡发生。区内资源开发和经济建设在采取合理地震设防和灾害防治的条件下,可按照安全区进行布局建设。

5 结 论

1)多源遥感技术是进行震后泥石流灾害动态监测与关键评价指标参数获取的重要技术方法。提出的震后泥石流物源信息自动提取模型能实现震后泥石流灾害判识和物源信息自动提取,并分析其类型、分布、储量、频率和演化规律,具有较高精度和可操作性,可满足区域次生山地灾害遥感快速调查和动态监测需要。

2)强震极大改变了震后泥石流形成条件,使得震后泥石流物源增加、临界雨量降低、沟道堵塞严重流量增大。结合野外调查和泥石流物质量计算数据,采用配方法修正了震后泥石流流量计算公式;与震前相比,堵塞系数增加1.72~3.46倍,流量呈指数增大。

3)地貌信息熵蕴含了流域高程、面积、相对高差、沟床比降等较多的重要地形地貌因子,可作为泥石流危险性评价的复合地貌指标。在对研究区进行小流域划分的基础上,引入地貌信息熵,简化了震后泥石流危险性评价参数,实现了对子流域地貌发育阶段和沟谷危险程度的判断。研究区域大部分流域处于地貌侵蚀发育的壮年偏幼期,极少老年期和幼年期,具备发生泥石流的良好地貌条件。

4)根据震后泥石流形成条件、分布规律和发育特征,选取了震后泥石流流量、频度、地貌信息熵等10个重要因子建立评价指标体系;根据形成因素贡献率,对各指标因子进行量化处理和分级赋值,引入灰色系统理论和危险性指数,进行评价因子关联度分析和权重计算,从而构建更为简洁的震后潜在泥石流危险性评价模型。

5)研究区大部分泥石流沟处于中度及以上危险区,具有极大的潜在危险性。高、极高危险性泥石流沟89条,占沟道总数的39.2%,主要沿岷江干流和渔子溪河谷两岸分布,是研究区人口与社会经济聚集区、恢复重建的高危区和防灾减灾的重点区,资源环境承载力低,不宜进行规模以上开发建设,需加强灾害监测预警与综合防治。中度危险区毗邻高、极高危险区,共有泥石流沟48条,占沟道总数的21.2%,大部分沿渔子溪和草坡河两岸分布,次生山地灾害以中型为主且不甚密集,但峡谷和高陡边坡地带,也常有大型崩塌、滑坡等次生山地灾害发生,需要加以处置和防范。低、较低危险区毗邻中度危险区,大部分位于高海拔地带,受地震和人类活动影响较小,次生山地灾害以中小型为主,且受微地貌和小构造控制,泥石流少见,资源环境承载力较高,在合理的地震设防和灾害防治下,可进行较大规模的资源开发与社会经济建设,但人口与资源、环境仍需协调发展。

参考文献
[1]
Cui Peng,Wei Fangqiang,He Siming,et al. Mountain disasters induced by the earthquake of May 12 in Wenchuan and the disasters mitigation[J]. Journal of Mountain Science, 2008, 26(3): 280-282. [崔鹏,韦方强,何思明,等. 5·12汶川地震诱发的山地灾害及减灾措施[J]. 山地学报, 2008, 26(3): 280-282.]
[2]
Xu Qiang. The 13 August 2010 catastrophic deberis flows in Sichuan province:Characteristics,geneticmechanism and suggestions[J]. Journal of Engineering Geology, 2010, 18(5): 596-608. [许强. 四川省8·13特大泥石流灾害特点、成因与启示[J]. 工程地质学报, 2010, 18(5): 596-608.]
[3]
Eldeen M T. Predisaster physical planning:Integration of disaster risk analysis into physical planning—A case study in Tunisia[J]. Disasters, 1980, 4(2): 211. DOI:10.1111/disa.1980.4.issue-2
[4]
Adachi K,Tokuyama K,Nakasuji A,et al. Study on judgment of outbreakability of debris flow[J]. Sabo Gakkaishi, 1977, 30(3): 7-16.
[5]
Hollingsworth R,Kovacs G S. Soil slumps and debris flows:Prediction and protection[J]. Bulletin of the Association of Engineering Geologists, 1981, 38(1): 17-28.
[6]
Carrara A,Cardinali M,Detti R,et al. GIS techniques and statistical models in evaluating landslide hazard[J]. Earth Surface Processes & Landforms, 1991, 16(5): 427-445.
[7]
Chau K T,Lo K H. Hazard assessment of debris flows for Leung King Estate of Hong Kong by incorporating GIS with numericalsimulations[J]. EGS-AGU-EUG Joint Assembly, 2004, 4(1): 103-116.
[8]
Hungr O. A model for the runout analysis of rapid flow slides,debris flows,and avalanches[J]. Canadian Geotechnical Journal, 1995, 32(4): 610-623. DOI:10.1139/t95-063
[9]
Pasuto A,Soldati M. An integrated approach for hazard assessment and mitigation of debris flows in the Italian Dolomites[J]. Geomorphology, 2004, 61(1/2): 59-70.
[10]
Jakob M,Friele P. Frequency and magnitude of debris flows on Cheekye River,British Columbia[J]. Geomorphology, 2010, 114(3): 382-395. DOI:10.1016/j.geomorph.2009.08.013
[11]
Jakob M,Holm K,Weatherly H,et al. Debris flood risk assessment for Mosquito Creek,British Columbia,Canada[J]. Natural Hazards, 2013, 65(3): 1653-1681. DOI:10.1007/s11069-012-0436-6
[12]
Castelli F,Freni G,Lentini V,et al. Modelling of a debris flow event in the Enna area for hazard assessment[J]. Procedia Engineering, 2017, 175: 287-292. DOI:10.1016/j.proeng.2017.01.026
[13]
Wang Lixian. On the torrent classification[J]. Journal of Beijing Forestry College, 1982(3): 94-107. [王礼先. 关于荒溪分类[J]. 北京林学院学报, 1982(3): 94-107.]
[14]
Liu Xilin. Study on judgement on debris flow risk degree[J]. Journal of Catastrophology, 1988(3): 10-15. [刘希林. 泥石流危险度判定的研究[J]. 灾害学, 1988(3): 10-15.]
[15]
Tang Chuan,Liu Xilin,Zhu Jing. The evaluation and application of risk degree for debris flow inundation on alluvial fans[J]. Journal of Natural Disasters, 1993, 2(4): 79-84. [唐川,刘希林,朱静. 泥石流堆积泛滥区危险度的评价与应用[J]. 自然灾害学报, 1993, 2(4): 79-84.]
[16]
Han G,Wang D. Numerical modeling of Anhui debris flow[J]. Journal of Hydraulic Engineering, 1996, 122(5): 262-265. DOI:10.1061/(ASCE)0733-9429(1996)122:5(262)
[17]
Wang Mingwu. Classification of dangerous degree of debris flow based on neural network[J]. Hydrogeology & Engineering Geology, 2000, 27(2): 18-19. [汪明武. 基于神经网络的泥石流危险度区划[J]. 水文地质工程地质, 2000, 27(2): 18-19.]
[18]
Lin C W,Liu S H,Lee S Y,et al. Impacts of the Chi-Chi earthquake on subsequent rainfall-induced landslides in central Taiwan[J]. Engineering Geology, 2006, 86(2): 87-101.
[19]
Tie Yongbo,Tang Chuan,Su Xiaoqin. Shitougou debris flow dynamic characters and its hazard risk assessment in Jiulong of Sichuan provience[J]. Research of Soil and Water Conservation, 2008, 15(5): 168-170. [铁永波,唐川,苏小琴. 四川省九龙县石头沟泥石流动力学特征及其危险性评价研究[J]. 水土保持研究, 2008, 15(5): 168-170.]
[20]
Liu Hongjiang,Lan Hengxing,Cheng Weiming. Numeric simulation for group-occurring debris flows and its risk analysis in Jiegu township after the earthquake of Yushu,Qinghai province[J]. Journal of Mountain Science, 2010, 28(4): 444-452. [刘洪江,兰恒星,程维明. 玉树地震后结古镇群发式泥石流灾害数值模拟及危险性分析[J]. 山地学报, 2010, 28(4): 444-452.]
[21]
Shu Heping,Qi Shi,Ning Na,et al. Risk assessment of debris flow disaster:A case study of wududistrict in the south of Gansu province,China[J]. Journal of Natural Disasters, 2016, 25(6): 34-41. [舒和平,齐识,宁娜,等. 甘肃省南部武都区泥石流灾害风险评价研究[J]. 自然灾害学报, 2016, 25(6): 34-41.]
[22]
Han Yongshun,Cui Peng,Zhu Yingyan,et al. Remote sensing monitoring and assessment of traffical damage by Wenchuan earthquake—A case study in Du-Wen highway[J]. Journal of Sichuan University(Engineering Science Edition), 2009, 41(3): 273-283. [韩用顺,崔鹏,朱颖彦,等. 汶川地震危害道路交通及其遥感监测评估——以都汶公路为例[J]. 四川大学学报(工程科学版), 2009, 41(3): 273-283.]
[23]
Gu Haiyan,Yan Li,Li Haitao,et al. An object-based automatic interpretation method for geographic features based on random forest machine learning[J]. Geomatics and Information Science of Wuhan University, 2016, 41(2): 228-234. [顾海燕,闫利,李海涛,等. 基于随机森林的地理要素面向对象自动解译方法[J]. 武汉大学学报(信息科学版), 2016, 41(2): 228-234.]
[24]
Kang Zhicheng,Zhang Jun. Study and calculation of discharge of debris flow[J]. Soil and Water Conservation in China, 1991(2): 15-19. [康志成,张军. 泥石流洪峰流量的研究与计算[J]. 中国水土保持, 1991(2): 15-19.]
[25]
Ai Nanshan,Yue Tianxiang. Second discussion of the comentropy of drainage system[J]. Journal of Soil and Water Conservation, 1988, 2(4): 1-7. [艾南山,岳天祥. 再论流域系统的信息熵[J]. 水土保持学报, 1988, 2(4): 1-7.]
[26]
Han Y S,Dong S K,Chen Z C,et al. Assessment of secondary mountain hazards along a section of the Dujiangyan-Wenchuan highway[J]. Journal of Mountain Science, 2014, 11(1): 51-65. DOI:10.1007/s11629-012-2516-1