孙黎明,蔡 红,严 俊,李维朝,肖建章
(中国水利水电科学研究院 岩土工程研究所,北京 100048)
尾矿库通常由筑坝拦截谷口或围地形成,用以堆存尾矿或其他工业废渣。中国有数量巨大的尾矿库,早年低标准修筑或者疏于管理的一些尾矿库存在溃坝风险,一旦失事,将对人民生命和财产安全产生重大威胁,造成重大事故[1-4]。巴西Brumadinho尾矿库和中国襄汾尾矿库失事都造成了巨大的人员财产损失,因此需要研究尾矿库溃坝相关抢险和救灾技术,特别是在溃决后,需快速获取溃决前后的地形数据、尾矿坝、尾矿库等建立现状地形模型,用于快速评估灾情,为抢险救灾提供信息化支撑[5]。
在尾矿库溃决后,首先需要获取溃决前后的尾矿库内地形数据,其中以尾矿坝溃前和溃后的三维几何模型及溃后下游影响区地形为重点,为尾矿坝的数值模拟分析提供几何模型基础。结合GIS(地理信息系统)空间分析技术,可以快速计算溃决的体积、区域面积、损失程度和生境等信息。开展尾矿库溃坝前后模型重建相关的专业研究[6-8],有助于在灾情发生的第一时间快速恢复相关三维几何模型,为虚拟仿真和数值计算提供计算模型基础。过去针对尾矿库相关的分析计算研究,主要是利用卫星、无人机、GNSS(全球定位系统)测量等手段进行数据感知,可以在尾矿坝溃后第一时间获取DEM(数字地形模型)、DOM(数字正射影像)等数据,构建地形模型和尾矿坝体三维几何模型,并基于GIS平台,实现尾矿库溃决过程的精准化虚拟仿真,然而这些只是针对常规矿区溃后的地形模型和尾矿库地形构建方案。此外,在一些偏远地区或者闭库的资料匮乏区域,溃决前的地形数据和坝体详细资料在短时间内较难获取。目前,针对多源数据融合的快速建模或针对少资料地区的尾矿库模型恢复重建的相关研究还相对较少,无法为实际抢险救灾提供支撑,因此,针对少资料地区尾矿库溃决后的地形模型和尾矿库模型恢复方法,是当前亟需研究的重点课题之一[9]。
目前,针对尾矿库的地形获取方法主要集中于利用无人机或者卫星遥感手段直接获取DEM,后续再基于此建立尾矿库影响区的三维模型,进行尾矿库安全的预测预警和分析评估[10-11],以及通过已有的基础图纸,结合GIS、CAD、CAE类软件构建坝区的三维模型,将三维模型用于风险评估、溃决模拟和损失情况等数值模拟分析、空间计算分析和三维虚拟仿真[12-13]。然而,目前针对少资料情况下尾矿库的多源异构数据感知与融合,以及溃坝前后三维模型构建方法的相关研究还较少[14]。对此,本文提出一种适用于少资料地区的基于多源空间数据的尾矿库三维模型重建方法,以光学遥感卫星、雷达卫星、无人机等方式获取到的多种数据为基础数据源,结合尾矿坝和尾矿库的数据结构与几何分布特征,利用克里金插值方法或样条插值方法结合模型转换融合方法,实现栅格数据的高精度融合,完成尾矿库溃坝前后的三维模型和下游地形模型重建,为后续的计算分析评估提供坝体三维几何模型、网格模型和地形模型支持。本文以2008年山西襄汾陶寺乡980沟尾矿库溃坝事件为研究案例,采用此方法构建了襄汾尾矿库溃坝前后的高精度三维模型[15-16],可以为后续少资料地区尾矿库三维模型的快速重建提供参考。
尾矿库溃决后需要尽快进行灾情评估和抢险,而全面的数据是进行评估的基础[17]。研究尾矿库数据获取的方法和标准有助于在尾矿库溃决后第一时间获取数据,实现数据高效融合利用[18-19]。目前,尾矿库溃决后的信息感知有以下几个难题:① 抢险救灾所需数据感知方式分散;② 针对尾矿库失事所需的感知信息时效性不高;③ 多源异构数据不能有效融合;④ 偏远少资料地区的尾矿库在溃后较难重建出满足需求的三维模型。
近些年,关于尾矿库的信息快速感知方面的研究重点集中在采用无人机、卫星和GNSS等手段,获取某一时刻尾矿库所在区域的地形模型,但对多种方法融合应用的研究还较少;且这些模型多用于三维虚拟浏览,而针对少资料地区特别是缺少溃决前高精度地形数据的地区,直接在溃后获取的地形数据不能满足坝体溃决和泥石流损失评估的需求(需要将溃坝前的坝体、库区和下游地形与溃决后的情况进行对比,得到2个时间段的同一区域的地形模型数据)。此外,还需要构建坝体的内部结构三维模型,为坝体的数值模拟分析提供网格结构模型。以下为几种常见的尾矿库空间信息感知方法。
(1) 遥感卫星。光学遥感卫星可以利用卫星对地固定区域的周期性访问,获取地面的地理光学影像。商业光学卫星的分辨率通常可至20 cm,而倾斜摄影技术和雷达卫星可以得到全球大区域范围内的地形数据,这对于少资料地区或者偏远区域的尾矿库溃决研究有较大意义[20]。此外,卫星感知是全天候、全球化、受地面条件影响少的一种方式,特别在突然发生事故的情况下,卫星数据可第一时间得到全面的历史资料,获取溃决前一段时间的地形和影像数据,这对抢险的信息化支持至关重要。
(2) InSAR(合成孔径雷达)卫星。InSAR卫星技术可分析大区域的高精度形变,利用InSAR卫星的高时间分辨率和干涉特征,可以结合地面的常规DEM数据,分析地表形变和尾矿坝的形变特征,进而关注溃决成因和发生时间等关键参数,精度可达厘米和毫米级。InSAR技术还可以普查大区域范围内的滑坡形变及应用于坝体形变分析。
(3) 无人机。无人机主要用于快速、高精度获取小范围内的高精度影像和DEM。随着旋翼无人机的技术发展和成本大幅降低,其成为高精度地形和影像数据获取的最常用方法。在发生灾害第一时间,通常采用无人机数据,其中固定翼无人机体积更大、飞行高度更高,数据获取效率更高,但是在高山峡谷区飞行困难,比较适合大区域范围内的测绘作业;旋翼无人机更小、灵活性更高,比较适合小范围更低高度的数据获取,这样得到的数据精度和时效性都更高[21]。但是,需要对无人机获取的影像进行后期的软件空间解算和处理才能得到地理数据,这个过程较为耗时,无法实时得到数据。
(4) LiDAR。LiDAR主要借助于搭载在无人机上的小型三维激光扫描设备进行地形或其他结构的外部点云获取。点云数据的特点是垂直精度比水平精度高,这是其相比于光学无人机的主要优点。同时,点云不需要复杂的后期计算,只需要与基站进行精度矫正,因而可以实时得到高精度地形数据,在精细化的坝体模型构建和计算部分可以快速应用。
(5) 监测数据。GNSS主要是监测布置在坝体及其附近的监测仪器数据。GNSS基准点和形变监测数据可以为尾矿库的三维模型分析和计算提供高精度位置参考数据。
(6) 图纸资料。尾矿库的地形图、设计图纸和计算资料等可以用于后期的三维建模,特别是在坝体内部结构建模时,用于坝体模型与库区模型的离散构建。但是,通常设计图纸的坐标体系和地理坐标系不同,因而需要考虑图纸和地理数据的融合。
通过感知得到的多源异构数据,其格式和精度通常不同,并且不同时间、不同设备、不同方法,其数据范围也不完全重叠。为了充分利用这些数据,必须进行数据融合,将现有多分辨率、多结构、不同范围的数据融合成为一个范围更大、结构更统一、精度更高的数据,为尾矿库的三维模型构建提供可靠的数据来源[22]。根据数据需求分以下3种情况。
DOM主要指通过卫星或无人机等得到的带有地理投影的航空像片或遥感影像[23]。利用卫星的时序性特征,可以收集低分辨率的尾矿库溃决发生前的DOM数据,其范围通常相对较大,包含了溃决前的房屋、道路等生境特征;溃决后得到的高分辨率无人机影像可以与溃前的卫星数据进行融合。在两者进行投影一致化的基础上,影像融合需要结合特征级融合和决策级融合两种方法,特征级融合的重点是寻找尾矿坝的边界轮廓、房屋的房角、道路等清晰边界进行两张影像的校正和匹配,最终两张图像进行影像镶嵌,完成影像的融合,得到最大范围且多个时间点的DOM影像。这个过程可以借助于QGIS、GRASS、ILWIS等相关GIS和遥感软件中的数据融合功能完成。
地形数据是计算和可视化中最重要的基础数据,通常也是由不同数据结构存储的,包括地形图(等高线CAD数据)、DEM、测点和点云的离散点数据(由LiDAR或者GNSS测量得到),因此需要研究可以融合地形数据的方法。融合方法主要分两步:① 首先进行投影坐标的统一和位置逐点对照,可以选择高精度的GNSS点作为控制点坐标,以高精度的点云数据为基础数据,然后确定坝体形变测量基准点、道路交叉点等进行数据空间校正;② 按照最终需要的数据分辨率进行上采样或者下采样,通过一些常见的算法实现低分辨率采样或者高分辨率插值,对于点云稀疏的地区或者遮挡区域,也需要用插值方法进行加密,常用的插值方法主要是地学领域的空间插值方法,包括克里金插值、反距离加权插值、样条插值等[24],在插值时需要考虑数据分布的各向异性,进行针对性的加权以提升精度,图1(a),(b)分别为低分辨率卫星地形和高分辨率无人机DEM;③ 将不同格式的数据进行结构转换以完成整个地形数据融合过程,这个过程可以借助于相关GIS和遥感软件[25]。
图1 低分辨率和高分辨率DEMFig.1 Low and high resolution DEM
根据三维无人机倾斜摄影模型或BIM等三维模型,可以修正地形模型和影像数据。因为三维模型通常有更详细的信息,而且分辨率更高,因此可利用三维模型中的桥梁、坝体、洪痕等高精度区域数据更新地形模型,融合过程采用关键点匹配更新的方式。数据融合的具体过程如下。
(1) 统一参考体系。将规则格网数据和模型的坐标系调整到一个坐标系,将地理坐标系转换为投影坐标系。
(2) 特征点匹配。选定2个数据的关键特征点,特征点应精度高且容易识别,通常可以参考无人机飞行的控制点,或者道路、房屋的角点等。
(3) 空间校准。利用GIS软件将不同数据进行空间几何校正,相同特征点进行一一对照。
(4) 空间插值。按照最高精度的数据,对低分辨率的数据进行空间插值计算,得到统一高精度的数据。
(5) 模型融合。将规则格网数据转换为模型的三角网格式,然后将已有坝体模型与地形模型进行无缝融合,利用CAD、CAE等3D建模类软件完成。
针对缺乏溃前高精度地形和设计资料的情况,需要构建溃坝前三维模型,用于和溃后模型进行对比分析计算,为数值模拟和三维虚拟仿真提供模型。为了构建溃坝前后尾矿库的三维模型,需要研究基于低分辨率地形数据构建三维坝体模型的方法。本文提出了一种基于样条曲面方法的三维模型恢复法,可根据低精度的地形数据构建尾矿库和影响区域的三维模型。首先从低分辨率的数据上提取溃坝前地形的几何特征,然后用于溃后的无人机高精度数据上,构建坝体和库区的三维模型,最后根据两次数据恢复溃前的地形数据。
在少资料地区的尾矿库溃坝模型构建过程中,前期通常只有卫星的地形数据,分辨率较低(几米),无法满足分析计算误差要求。在少资料的情况下进行模型构建,需要尽可能充分利用、融合多个数据,因而需要从溃坝前的各种资料提取所需要的重点几何特征,如表1所示。
表1 尾矿库的几何特征分类Tab.1 Classification of geometric characteristics of tailings reservoir
(1) 坝体模型构建。构建尾矿库坝体的三维模型通常需要融合溃坝前的设计资料和地形。在缺乏设计资料的前提下,溃坝前的卫星数据分辨率很低,溃决后使用的无人机影像分辨率高,却没有坝体的完整模型,只有溃决后的地形。构建坝体模型需要先根据坝体的高程、宽度、马道和轮廓构建一个基础断面,将断面与地形进行融合,然后对断面沿着坝轴线进行三维拉伸进而完成模型构建,将尾矿库恢复后得到的模型与高精度地形模型进行一致化融合。
(2) 库内构建。通常需要根据低分辨率的地形DEM数据,利用样条曲面函数进行空间插值以构建库区内模型。样条曲面插值适合于三维模型的构建,其可以从已有数据得到所在范围内未知位置点的值,获取高分辨率的平滑数据[26-27]。以坝体内侧边界到山体边界为库内尾矿砂的模型表面。样条曲面函数的计算方法如式(1)~(2)所示:
(1)
(2)
根据上文通过数据提取的几何特征,为得出下游影响区在溃坝前的地形,需要在溃决后的地形上进行恢复。恢复需要通过对洪痕高程插值进行模型构建,另外在河流和道路的边界区域,通过测量溃坝洪痕点所在的地形数据,结合道路、河流等约束条件构建DEM地形曲面;对于房屋模型,需要在DEM地形曲面的基础上进行构建,将房屋模型嵌入到三维地形模型之上;三维房屋桥梁等的建筑模型的构建采用freecad,sketch up,3d max和blender等建模软件完成。
根据本文提出的方法,按照图2所示的技术流程,以2008年山西省襄汾980沟尾矿库溃坝事件为研究对象,构建了三维模型,以用于后续的风险评估、数值分析计算虚拟仿真性。
图2 总体技术流程Fig.2 Overall technical flow
由于尾矿库溃决时已经是闭库状态,且管理不规范,缺乏溃坝前的详细设计资料,因而构建溃坝前后的三维模型存在较大困难。溃坝后过泥面积为0.302 km2,最大纵深2.5 km,最大横宽350 m,卸库体积26.8万m3,造成了重大的人员伤亡。其中,利用旋翼低空无人机进行了倾斜摄影测量获取DEM和DOM,襄汾尾矿库的DOM见图3,无人机飞行面积1.84 km2,航线长度24 km,覆盖范围的长度3 100 m,宽度600 m。利用无人机倾斜摄影技术得到了高精度5 cm分辨率的DOM和10 cm分辨率DSM。此外,收集到在尾矿库发生前,2007年相关区域20 m的卫星遥感DEM数据,通过Google earth获取到了溃坝前多个时间点的DOM数据,但是DEM的分辨率也都较低,无法满足溃坝计算和虚拟仿真的需要。
图3 尾矿库遥感影像DOMFig.3 Remote sensing DOM of tailings reservoir
利用洪痕采集记录溃坝后尾矿砂覆盖的范围和典型断面的堆积厚度,用于更新溃坝前后的地形和模型,为后期数值计算结果准确性提供校核依据。现场共采集洪痕点107处(图4),针对现场测量的房屋和道路的洪痕(图5),采用GPS静态测量完成,精度在毫米级。基于本文提出的方法,首先根据卫星、无人机、GNSS得到的地形和影像进行了多源异构数据的融合,调查了溃决后的房屋、桥梁和河流的损失情况,融合了洪痕点数据,然后在数据基础上建立了尾矿坝、库区和下游影响区的三维模型(图6),模型的精度满足溃坝计算的要求。根据该方法得到溃坝后的DEM与地形图(图7)及溃坝前的DEM与不规则三角网模型(图8)。因此本文提出的方法可为后续的虚拟仿真和数值模拟分析提供了三维模型基础,包括下游冲击损失评估和溃坝模拟分析等。将构建好的模型转换为三角形网格单元,利用光滑粒子流体动力学(Smoothed Particle Hydrodynamics,SPH)流体模拟方法,三维模型作为地形约束曲面,输入条件后模拟溃坝过程,可以详细计算尾矿砂堆积过程、范围和厚度,以此详细进行损失评估和尾矿库溃坝的风险和模式研究,同时验证了方法的可行性。
图4 洪痕点分布Fig.4 Flood marks distribution
图5 洪痕点高程测量Fig.5 Flood marks elevation measurement
图6 基础DEM地形数据和与坝体融合后的三维模型Fig.6 Basic DEM terrain data and the 3D model integrated with dam data
图7 融合后得到的襄汾尾矿库的溃后DEM和地形图Fig.7 DEM data and topographic map of Xiangfen tailings reservoir after dam failure
图8 溃坝前坝体与地形融合模型Fig.8 3D modeling of dam and the terrain integration before dam failure
尾矿库溃坝后的数值分析和三维虚拟仿真计算为抢险救灾提供了重要支撑,而少资料地区通常需要在多源多分辨率异构数据的基础上快速构建尾矿库溃坝前后的三维模型。本文提出了一种基于多源异构数据融合的三维模型高精度构建方法,实现了在卫星、无人机、GNSS等多源多分辨率异构空间数据融合的基础上,提取尾矿库三维模型的几何特征,结合房屋、洪痕、损失调查等,利用样条曲面插值方法,提出了溃坝前后尾矿坝、库区和下游影响区溃坝前后三维模型的恢复方法,实现了最大化的数据利用,为少资料地区尾矿库提供了溃坝前后三维模型快速构建支持,未来可以继续在数据融合与几何特征提取方面进行更多研究,提升整体模型构建的精确度。