基于LiDAR数据开展活动断层填图的实验研究
——以新疆独山子背斜-逆冲断裂带为例
2014-07-02魏占玉何宏林徐锡伟甘卫军卫蕾华
魏占玉 何宏林 高 伟 徐锡伟 甘卫军 卫蕾华
(中国地震局地质研究所,活动构造与火山重点实验室,北京 100029)
基于LiDAR数据开展活动断层填图的实验研究
——以新疆独山子背斜-逆冲断裂带为例
魏占玉 何宏林 高 伟 徐锡伟 甘卫军 卫蕾华
(中国地震局地质研究所,活动构造与火山重点实验室,北京 100029)
机载LiDAR技术为描绘活动构造相关构造地貌和最新的地表形变提供更精确的基础数据。如何将LiDAR新技术、新数据应用于活动构造填图和活动断层地震危险性评价等方面,是今后活动构造研究领域的一个重要的发展方向。文中以新疆天山北麓的独山子背斜-逆冲断裂带为试验区,开展了基于LiDAR数据的活动构造填图实验研究。首先,采用机载LiDAR技术进行数据采集,获得点云密度为6.6个/m2、平均点间距为0.39m的LiDAR原始数据;其次,利用试验区内12个测量精度可达mm级的GPS静态测量点评估LiDAR的相对垂直精度为0.12m、均方差值为0.078m;最后,对密度为6.4个/m2的地面点云数据进行DEM最佳分辨率评估,利用反距离权重算法获得0.5m分辨率的数据高程模型(DEM)。该分辨率的DEM数据足以完成独山子背斜-逆冲断裂带的精细构造地貌特征的确定以及高精度的空间解译。文中仅使用DEM可视化工具从不同虚拟的视角、不同色度或其他处理方式来识别微构造地貌、划分地貌面和确定断层位置等,宏观上获得与前人通过航片解译和野外调查一致的断裂分布特征,微观上较前者具有更高的精细程度。此外,数据采集、数据质量检验、数据处理及数据应用等技术和方法适用于其他能够获得LiDAR地形数据的活动断裂研究工作。
LiDAR 独山子逆冲-背斜带 活动断层填图
0 引言
在过去十几年的发展中,机载LiDAR(Light Detection And Ranging)技术已成为可精确、快速、可靠地获取地面3维数据的重要手段,在地球科学中的应用日益广泛。机载LiDAR最吸引人的优势在于激光束能够通过植被缝隙到达地表,并能够准确测量地表形态。这些数据可以在适当的尺度(m)和精度(dm)上展现断层活动在地表上所留下的现象。在本世纪初美国已经开始“B4 LiDAR”项目,该项目利用LiDAR技术(也被称为Airborne Laser Swath Mapping,ALSM)对南加州圣安德烈斯和圣哈辛托断层进行地质调查和地貌填图(Arrowsim th et al.,2008,2009),推动LiDAR技术在构造地貌测量及地貌演化研究中的应用(Crosby et al.,2007)。近年来,日本及欧洲的研究者利用LiDAR数据获得高精度、高分辨率的DEM(Digital Elevation Model)成功地识别地震地表破裂带及精细地貌划分(Cunningham et al.,2006;Kondo et al.,2008;Lin et al.,2013)。
在活动构造地质填图工作中,通常都会使用各种遥感信息进行构造地貌解译,特别是具有m级分辨率的航空照片(何宏林,2011),但是航空摄影测量等技术仅用于2维层面上的断层破裂定位,以及定性的3维空间辅助解译工作(陈桂华等,2006)。与传统的摄影测量方法相比,LiDAR测量是对地表3维坐标的直接测量,避免通过摄影像对解算或雷达干涉测量(In-SAR)间接获得地表3维坐标,从而大大提高测量精度。另外,通过对LiDAR原始数据后期处理,能够获得高精度、高分辨率的DEM数据。这些DEM数据能够精细表达地表地貌结构,为研究和理解活动构造相关构造地质情况、地貌和最新的形变历史提供更精确的基础数据。如何将这些LiDAR新技术、新数据最大效率地用于活动构造的填图工作和活动断层地震危险性评价等方面,无疑是今后活动构造研究领域的一个重要的发展方向。
新疆天山北麓的独山子背斜-逆冲断裂带发育在准噶尔盆地南缘,为天山向准噶尔盆地挤压逆冲在沉积盖层中形成的新生代褶皱带,晚第四纪以来依然在持续活动,而横切独山子背斜-逆冲断裂带的奎屯河发育7级基座阶地,它们的形成与晚第四纪时期独山子背斜、玛纳斯背斜的褶皱隆升相关(杨晓平等,2012)。选取新疆天山北麓的独山子逆冲-背斜带作为试验区,进行机载LiDAR数据采集并开展基于LiDAR数据的活动构造研究,其目的是试验和建立基于LiDAR数据进行活动构造填图的技术流程。LiDAR应用研究是非常广泛的,从数据采集、数据质量检验、数据处理及数据应用各方面都是值得研究的内容。如何检查原始LiDAR覆盖区域的正确性,评估原始LiDAR数据质量以及制作DEM数据是利用LiDAR数据开展构造地貌研究的基础。但是,LiDAR数据使用者往往只关心数据的应用,对于LiDAR原始数据质量检验和数据处理方法重视不足,这将影响使用者对LiDAR数据的理解及有效使用。文中将以新疆独山子背斜-逆冲断裂带的原始LiDAR数据为基础,介绍如何进行LiDAR原始数据质量检验和数据处理,如何以LiDAR数据为基础制作高精度的DEM数据,如何基于高精度的DEM数据进行断层解译、识别和制图。
1 数据采集工作原理
机载LiDAR是一种集激光扫描测距系统(Laser Scanner)、全球定位系统(GPS)和惯性导航系统(INS)于一体的新型航空遥感系统,能够快速、精确地获得地面3维信息。在飞行测量过程中,LiDAR系统将发射高频激光脉冲依次扫描被测目标,利用每个激光脉冲从发射源到被测物体表面再返回系统的时间来计算距离,航摄时配合动态GPS及惯性导航系统,即可获取脉冲发射瞬间摄点的空间坐标及其空间姿态,从而计算出被测物体的3维坐标(图1a)。为提高测量过程中各时刻的3维坐标数据的精度,还需要在地面沿航线布设一定数量的GPS基站。
测量原始数据是具有3维坐标信息的点云(Point Cloud)数据,这些点云数据可以构成被测目标虚拟影像(图1b)。但是这些原始数据只是没有属性的离散点阵数据,不能区分地面点和非地面点(如植被、建筑物);而且,所采集的点云数据包含大量的冗余数据及系统误差,不能直接用于构造地貌识别和测量。鉴于上述原因,在使用LiDAR数据之前必须进行数据质量检验、数据分类和DEM数据提取等。
图1 LiDAR工作原理示意图(a)及LiDAR点云渲染图(b)Fig.1 Diagram of LiDAR operating principle(a)and renderings of LiDAR point-cloud(b).
2 LiDAR测量概况
所选取的LiDAR测量区为位于天山北麓近EW向展布的独山子逆断裂-背斜带,与其东侧的哈拉安德和安集海等逆断裂构成天山北麓第3排逆断裂-褶皱带。独山子逆断裂-背斜是该带南倾主逆断裂和背斜构造发育最好、在地表出露最为完好的段落。该逆断裂东西长约18km,由发育在背斜带北翼山前地带的主逆断裂、后冲断裂和NW向剪切断裂构成(邓起东等,1991)。奎屯河近SN向横切该逆断裂-背斜带。由于背斜的强烈抬升和河流的下切侵蚀,在河谷两侧形成多级阶地,阶地受到褶皱弯曲和断裂切割(邓起东等,1991;Molnar et al.,1994)。由于地表破裂带和河流阶地是重复地震活动和构造运动的记录仪,能够反映该区最新逆断裂-褶皱变形的几何学和运动学特征,因此本次LiDAR测量选取的测量区呈“T”型,分别覆盖逆断裂带和逆冲上盘的河流阶地(图2)。
LiDAR航测选在10月底进行,这样可尽量减少植被对地表的遮挡,以提高激光脉冲到达地面的几率。测量目标区分为两部分,目标区Ⅰ呈EW向覆盖山前地表破裂带,扫描带长约22km,宽约3km;目标区Ⅱ呈SN向覆盖奎屯河阶地,扫描带长约25km,宽4km。LiDAR点云数据存储格式采用美国摄影测量与遥感协会(ASPRS)于2003年发布的面向机载LiDAR数据的标准格式LAS1.0,该点云格式可在不同软、硬件之间实现交互和存储。在航测过程中,同时获取对应地面的彩色数码影像,用于后期正射影像或其他专题图。
图2 独山子背斜-逆冲断裂地质简图及LiDAR测量区Fig.2 Geological sketch of Dushanzi anticline-reverse fault zone and the LiDAR measurement area.红色线代表山前逆断裂,虚线框带表LiDAR测量目标区
3 数据处理流程
原始LiDAR测量数据包括原始LAS点云数据和原始影像两部分,一般是由硬件设备随带的软件进行解算处理的。在航测采集数据之后,首先需要检查LiDAR测量区域是否与目标区一致、是否存在测量空区、航带重叠区与原始点密度是否符合要求等等。另外,原始LiDAR点云数据仅仅是空间离散点坐标数据,需要进行点云数据滤波与分类、去除异常点数据、测量精度评估等。最后,利用地面点云数据生成数据地面模型(DSM)和数字高程模型(DEM)得到地形与地物等信息,才可供后续的应用。
3.1 LiDAR数据检测
在LiDAR航测中,单次飞行只能测量有限宽度的地面范围,而且单次扫描的点云密度难以满足测量要求。由于LiDAR测量的最终成果是基于点云数据构建地物3维模型,点云密度直接影响到DSM和DEM的分辨率和精度(Langridge et al.,2013)。为了提高点云密度达到一定要求,常用的方法是多次重复扫描目标区,增大航测条带间的重叠度。如美国南圣安德烈斯断裂(San Andreas Fault,SAF)LiDAR测量“B4”项目中,在覆盖该断裂带的1.5km的条带上飞行5次获得地面点云密度约为3.6个/m2(Arrowsmith et al.,2009);在日本Neodani断裂上森林覆盖地区,对1.25km宽的条带进行13次重复LiDAR航测,平均原始点云密度可达到约12.7个/m2(Lin等,2013)。
在LiDAR原始数据处理之前,首先利用航带拼接图检验独山子LiDAR航测范围(图3)。航带拼接图显示LiDAR航测范围覆盖目标区(红色外框),不存在扫描空区。目标区Ⅰ宽度约为3km,进行了8次飞行扫描,表明目标区内任何位置被扫描不少于2~3次。目标区Ⅱ宽度约为4km,进行了16次飞行扫锚,目标区内任何位置至少被扫描2次。点云密度是扫描区内单位面积上测量点的数量。利用lastools软件的lasinfo函数统计点云信息,整个扫描区原始点云密度为6.6个/m2,平均点间距为0.39m;目标区Ⅰ原始点云密度为6.9个/m2,平均点间距为0.37m;目标区Ⅱ原始点云密度为6.5个/m2,平均点间距为0.39m。本次LiDAR测量点云密度高于美国“B4”项目的点云密度(Arrowsmith et al.,2008)。
图3 独山子LiDAR测量航带拼接图Fig.3 Mosaic of airborne LiDAR measurements in Dushanzi area.颜色代表LiDAR测量次数;黑十字点为野外GPS地面控制点位置;红线框是目标区范围
LiDAR数据精度决定其测量成果能否作为基础测绘资料直接应用,也是评价LiDAR数据质量好坏的标准之一。因此,LiDAR数据精度评估是开展LiDAR数据应用前必须回答的问题(王晓辉等,2007)。数据误差主要是由观测条件、激光信号发射率和回波信号强度等造成的,包含绝对误差和相对误差2个方面(赵礼剑等,2009)。为了测试独山子LiDAR数据精度,利用测量精度可达mm级的高精度GPS在目标区内静态测量12个观测点。每个观测点的持续测量时间>100min,有效观测卫星总数>5个。这些观测点主要分布在呈三角形的3个地点上(图3中黑点)。为了避免地形起伏引起的2次高程误差的影响,选取的GPS观测位置为视野开阔,地势平坦,周边地表无茂密植被(图4)。绝对垂直误差主要通过对机载LiDAR点与水平地面点的实际高程比较来评定。通过选取每个GPS观测点周边6m×6m范围内所有点云数据,利用不规则三角网计算GPS位置高程值(表1:lidar_Z),然后将该高程值与该GPS实测的高程值(表1:GPS_P)进行比较来评估LiDAR数据的绝对垂直误差。计算结果显示,独山子LiDAR数据的绝对垂直误差平均值为0.812m,均方差值为0.010m(表1:高程差)。
除了绝对垂直误差外,相对垂直误差也是评价数据精度(质量)的重要标准。利用GPS观测点评估独山子LiDAR数据的垂直相对误差,野外获取GPS观测点的垂直相对误差小于Li-DAR数据的相对误差。独山子LiDAR数据的垂直相对误差(Errrel)可以通过GPS点高程差矩阵(DiffGPS)与LiDAR点高程差矩阵(Difflidar)的差值获得:
图4 GPS测量野外现场照片Fig.4 Field views of GPSmeasurement.
表1 GPS观测点的测量数据及其LiDAR测量数据Table 1 Measured data of GPS and LiDAR
公式中:
式中,eij=Zi-Zj。
式中,eij=Pi-Pj。
计算结果的统计分析显示,LiDAR数据的垂直相对误差主要分布在0.1~0.2m之间(图5)。通过概率分布拟合,误差分布服从正态分布,均值为0.12m,均方差值为0.078m。该垂直误差可应用于后期利用LiDAR数据进行地貌测量及量化分析等方面。通过增加测量区点云密度,可以有效减小LiDAR数据的相对误差,如日本Nedani断裂上平均原始点云密度可达到约12.7个/m2,其相对垂直误差可降低到0.085m(Lin等,2013)。
3.2 LiDAR数据分类
由于地表要素复杂多样,所以LiDAR点云数据中包含了地面点和非地面点(如建筑物、树木和植被)。只有准确而可靠地将LiDAR点云数据进行分类处理后,才能有效地利用Li-DAR数据。点云数据分类是将LiDAR系统获取的分布在地表不同地物上的不规则离散点进行处理和识别,最基础的工作是将地面点与非地面点进行区分(隋立春等,2010)。针对机载LiDAR数据分类,国内外学者进行了深入研究,并提出多种分类算法,如不规则三角网渐进加密(Axelsson,1999)、分层线性内插(Kraus,1998)、渐进窗口数学形态学滤波(Zhang et al.,2003)和分层迭代选权滤波方法(刘经南等,2008)等。
图5 LiDAR原始数据相对误差概率分布Fig.5 Probability distribution of the relative error of raw LiDAR data.
由于独山子LiDAR测量目标区内建筑物较少,地表植被稀少,地形坡度变化小,文中使用Lastools软件提供的经典不规则三角网渐进加密方法对LiDAR点云数据进行分类处理。该方法简单且处理速度快,但是需要大量的经验参数,如最大的建筑物尺寸、三角型所在平面的倾角、确定平面点阈值等多个参数,这些参数的最优阈值需随不同场地复杂度进行实时调整。首先,选取一处具有建筑物和树木的区域作为点云数据分类的实验数据(图6)。自然地形是连续、平缓变化的,而地表树木、建筑物等点云数据通常比其相邻点的高程值大1~2m,因此分类方法可以识别并归类。图6显示大部分建筑物(褐色)和树木(绿色)能够被识别出来。从分类图的剖面图可以看出,建筑物的房顶和树木的形态符合自然特征,分类效果较好。对于一些低矮灌木和小型房屋等点云,还需人工去除归类(王明华等,2010)。
图6 LiDAR点云数据分类结果及分类点云横剖面图Fig.6 Classification of LiDAR point-cloud data and the transverse section.灰色代表地面点,绿色代表树木植被,棕色代表建筑物,黑色线AA′代表剖面位置
数字高程模型(DEM)把这些不规则离散的地面点云数据进行插值生成具有一定分辨率的2维网格,因此需要评估地面点密度,从而确定DEM的最佳分辨率。去除地表植被和建筑物等非地面点,分类后的LiDAR点云代表裸露地表。由于扫描区内建筑物和植被较少,分类后地面点云密度与原始点云密度相差不大。分类后地面点云平均密度为6.4个/m2,平均点间距为0.39m。
3.3 生成DEM
为了便于表现地形和量化分析地貌特征,通过插值将不规则离散的地面点转换成能够最佳表现地形表面特征的网格化数据模型(DEM),这一过程称为机载LiDAR点云的网格化(靳克强等,2011)。DEM的网格单元代表DEM的分辨率,越小的网格单元代表越高的分辨率和更多的细节,但是点云密度往往约束着DEM的分辨率(Langridge et al.,2013)。因此根据点云密度确定DEM最佳网格尺寸是利用插值方法生成DEM的必要步骤,然后选取恰当的插值方法形成所需要的DEM。
3.3.1 确定DEM分辨率
在利用点云内插生成DEM过程中,希望获得高分辨率(小尺寸网格)的DEM。但是如果设定网格过小,内插方法将引入一些人为特征,DEM不能准确表现地貌特征;如果设定网格过大,DEM不能表现细小的地貌特征,失去LiDAR数据的优势。目前主要根据经验和反复对比来确定最佳的网格尺寸(Arrowsmith et al.,2009)。地形数据的形态分析表明,地形起伏具有自相似性或自仿射性,并可以被分解成一组正弦组分(Perron et al.,2008),而最小正弦波代表着该地形数据的最高分辨率。根据该理论,文中介绍在一定点云密度下如何确定DEM最佳分辨率。
首先生成一个DEM模拟波状起伏的地形,该模拟DEM分辨率为5mm。该DEM在X方向存在一组波长和振幅相等的正弦波,分别为0.25m、0.5m、1m和2m,而在Y方向不变(图7)。在模拟的DEM中按某一密度(平均密度)随机选取离散点模拟LiDAR的点云分布,该离散点不偏向任何区域,利用这些离散点重构DEM。然后,利用傅里叶分析方法确定X方向上的正弦波的波长和振幅。根据Nyquist定理(采样定理),在进行模拟信号转换过程中,当采样频率大于信号中最高频率的2倍时,即采样波长小于最小波长的1/2时,采样信号可以完整保留原始信号的信息。在文中确定的DEM分辨率中,采样波长就是DEM的网格尺寸,最小波长是X方向上能够准确识别最小正弦波的波长。因此,重构后的DEM在X方向上能够准确识别最小正弦波长的一半就是在该离散点密度下最佳DEM分辨率。
图7 模拟波状起伏地形的DEM模型(a)以及离散点重构DEM模型中识别正弦波分布(b)Fig.7 DEM of simulated undulating terrain(a)and the sine distribution identified from the DEM reconstructed using discrete points(b).a黑色点表示密度为6个/m2的离散点分布;b每个明显峰值代表可以识别的正弦波
图7展示了点云密度为4个/m2、6个/m2和15个/m2的分析结果。结果显示,密度为4个/m2的点云构建的DEM可以较好确定波长为2m的地形起伏,因此该点云密度下DEM最佳网格尺寸为1m;当点云密度为6个/m2时,构建的DEM除了能够识别2m的地形起伏外,还能识别波长为1m的地形起伏,因此该密度下点云数据生成DEM最佳网格尺寸为0.5m;独山子LiDAR点云密度约为6个/m2,因此可以构建网格为0.5m的DEM。通过模拟发现,如果想获得更高分辨率的DEM,需要大大提高点云密度,如15个/m2点云密度可以获得0.25m网格的DEM。Lin等(2013)根据经验利用点密度为6.2个/m2和12.0个/m2的地面点可分别生成网格为0.5m和0.25m的DEM,与本文获得结果一致。
一些学者根据长度-面积关系得到点云密度与DEM网格尺寸的经验关系(Hu,2003),
公式中s为DEM网格大小(单位:m),d为点云密度(单位:个/m2)。Langridge等(2013)利用该公式确定点云密度低于2个/m2时DEM网格尺寸(图8),发现该公式适合于地势平坦的区域,对于植被复杂、地形多变的区域,该公式往往高估DEM分辨率。通过对比,在低密度条件(<3个/m2)下本文方法与经验公式(1)存在较大差别,但是在高密度(>5个/m2)下2种方法确定的DEM网格尺寸相似。
3.3.2 创建DEM
本研究的目的是基于LiDAR数据,通过准确识别构造地貌特征来确定活动断层位置及其活动特征,以此完成活动断裂填图。因此需要所选用的插值方法能够保留和凸显地形上的微小变化,这些微小变化可能指示着断层的最新活动;同时,要抑制插值过程所产生的线性特征,避免将这些人为线性特征解译成构造地貌特征。为了平衡微小地貌和人为特征,应根据点云数据选取恰当的插值方法以及参数,甚至反复实验、比较不同方法和参数(Arrowsm ith et al.,2009)。
由于地形特征依赖于近距离的测量点,而与远距离的测量点相关性低。El-Sheimy等(2005)提出了局部面元方法,该方法利用一定范围内的测量点来估计网格节点的高程值,既可提高计算效率,又可避免远距离测量点对网格节点的影响(图9d)。在这些面元上,DEM节点上的高程值是特定搜索半径(r)内测量点的函数。提高查找半径(r)可以确保面元内存在测量点(否则就会给网格节点分配一个空值null),但是该点的值实际上代表了一个更大范围的平均值。反距离权重(IDW)方法通过赋予近点比远点更高的权重而在一定程度上解决了这种问题(El-Sheimy等,2005)。Kim等(2006)给出了反距离权重(IDW)的最佳搜索半径(r)与DEM网格尺寸(s)的经验关系:
图8 地面点云密度与DEM网格尺寸关系图Fig.8 Relationship between the density of point-cloud and DEM grid size.方块及标注代表数据来源:1 Langridge等(2013);2 Adams等(2011):3 Arrowsmith等(2009);4 Lin等(2013):5 Begg等(2009);6本文结果
图9 使用反距离权重局部面元方法获得不同网格尺寸和搜索半径的DEM的阴影图Fig.9 DEM shadow maps with differentmesh sizes and different radius of search window using the method of IDW local surface element.a局部面元方法生成DEM示意图;b、c和d为网格尺寸分别为0.25m、0.5m和1m的DEM阴影图,依据公式(2)设置搜索半径
Arrowsmith等(2009)使用反距离权重局部面元方法获得0.25~0.5m分辨率、0.8m或者更大搜索半径的数字高程模型(DEM),并讨论了搜索半径对DEM分辨率的影响。Langridge等(2013)使用反距离权重方法,利用多植被、低密度地面点LiDAR数据仍然获得较好分辨率的DEM,能够识别m级的地表破裂特征。本文也采用反距离权重方法生成独山子LiDAR数据的DEM模型。DEM的网格尺寸分别为0.25m、0.5m和1m,依据公式(2)设置搜索半径。图9显示不同网格尺寸和搜索半径的DEM的阴影图。结果显示:1)当网格尺寸为0.25m时,许多网格因在查找半径内找不到测量点而返回空值(图9a中白色);2)网格尺寸为0.5m时,DEM模型中空值网格数迅速降低,并能清晰显示地表断层陡坎及冲沟特征(图9c);3)网格尺寸为1m时,DEM模型中不存在空值网格,并显示较为平滑的陡坎和冲沟特征(图9d)。但是在考虑分辨率的情况下,0.5m的DEM模型是最理想的DEM。
4 基于DEM的解译
4.1 DEM可视化方法
LiDAR数据处理的最终成果之一是数字高程模型(DEM),为利用3维数据可视化方法解译地貌单元和确定断层位置提供基础数据。目前,3维数据的可视化方法很多,需要筛选出适合于构造地貌解译的可视化方法。通过比较,晕渲图(Hill-Shading Map)(图10a)和坡度图(Slope Map)(图10b)适合于表现陡坎地貌。晕渲图通过模拟太阳光对地面照射所产生的明暗程度表现地形起伏,符合人眼观察实际地形的习惯。由于地面的阴影取决于太阳光的入射角和高度角,因此需要不断地调整入射光方向和高度角才能有效地识别微小地貌现象。坡度图尽管在表现坎状地貌(断层坎、阶地坎)方面效果较好,但是缺乏立体感。
图10 DEM数据的不同可视化方法效果图Fig.10 Different visualization renderings of DEM derived from the LiDAR data. a晕渲图;b坡度图;c彩色晕渲图;d赤色图
为了增强3维数据可视化的直观性和立体感,在地貌解译过程中将不同地形参数叠加在一起增强图形的可读性。因此,文中也使用了彩色晕渲图(Colored Hill-shading Map)(图10c)和赤色图(Red Relief Image Map)(图10d)。彩色晕渲图是晕渲图叠加高程设色图而成的,该图除了具有晕渲图的优势外,还能表现3维数据的高程等信息。赤色图是地形坡度、正openness、负openness 3个地形参数叠加融合的平面图,利用彩色和灰度的渐变来表现立体感(Chiba et al.,2008;张玲等,2014)。正、负openness通过计算一定范围内的最大天顶角(zen-ith)和天底角(nadir)来表现不规则表面的凹凸程度,即突显正负地形。正openness代表凸状表面,值越大,凸度越大;而负openness代表凹状表面,值越大,凹度越大(Chiba et al.,2008;张玲等,2014)。赤色图最大的优势在于可以突显地貌上的微小构造特征,如小冲沟、地表破裂带,同时也能准确地勾勒出面状构造的边界,如洪积扇、阶地面等。Lin等(2013)在对日本中部植被覆盖严重的区域进行微小构造地貌识别中,利用赤色图解译出了以前未知的断层陡坎和其他隐伏的地貌现象。
在地貌解译过程中,需综合使用各种DEM可视化工具从不同虚拟的视角、不同色度或其他处理方式来识别微构造地貌、地貌面划分和确定断层位置等(Arrowsm ith et al.,2008)。图10展示了上述4种可视化方法的结果,晕渲图和坡度图清楚显示断层陡坎位置、陡坎坡度和倾向;彩色晕渲图显示地形总体向NE倾,断层坎为逆断层坎;赤色图突出显示各类线性边界和小型阶地面,有助于识别和研究最新的地表变形。
4.2 断层解译结果
基于LiDAR数据生成的高精度、高分辨率的DEM及其阴影图,尝试进行独山子背斜-逆冲断裂地质填图。这是一种完全数字化、以GIS为基础、无野外验证的方法,通过寻找断层坎、错断冲沟,以及其他现象来确定断层位置。独山子逆断裂发育了非常清晰的构造地貌现象,包括清晰的断层坎、冲沟位错和废弃沟等,利用DEM晕渲图可以识别大部分断层坎和冲沟位错。图11对比了基于传统填图技术的独山子背斜-逆冲断裂地质图(邓起东等,2013)和基于LiDAR-DEM解译的独山子背斜-逆冲断裂的分布及主要地貌面的划分。尽管传统方法和基于LiDAR数据解译的结果都显示,独山子背斜-逆冲断裂由西向东逐渐分为2条分支断层,可以分为独山子背斜北翼山前逆断裂带(F1)、背斜东侧剪切断裂带(F2)和后冲逆断裂带(F3)三部分,但是对断层识别的精细程度上后者要远高于前者。
4.2.1 背斜北翼山前逆断裂带(F1)
背斜北翼山前逆断裂带总体走向近EW,空间上构成独山子背斜的北边界。该断裂带在地貌上主要表现为北倾断层崖,并被冲沟切割或被冲积扇覆盖而不连续。在奎屯河东、西两岸Ⅱ级阶地,该断裂均形成断层坎(图12)。在西岸Ⅱ级低阶地上,断层由2条分支断层构成,向西延伸至西岸Ⅱ级阶地上。由于山前冲洪积扇覆盖和人类活动的影响,西岸Ⅱ级阶地上断层陡坎向西逐渐消失。在奎屯河东岸的Ⅱ级阶地面上也发育2条近EW向的北倾断层坎,但延伸至石河子市南侧山前便消失,推测该断裂被冲洪积层覆盖。
4.2.2 背斜东侧剪切断裂带(F2)
背斜东侧剪切断裂带从独山子背斜东端向SE延伸,总体走向NW,地貌表现为多条不连续、NE倾断层崖。根据解译结果,该剪切断裂带在几何结构上由1条较长断层崖和3条较短断层崖呈左阶排列而成,其中每条较短断层崖又由若干次级破裂构成,形成帚状断层崖带(图13)。该剪切带上发育一系列大小不等的冲沟横跨断层崖,一些较大冲沟边界和阶地存在右旋位错(图10),一些小型冲沟因断层活动而废弃或改道。
4.2.3 后冲逆断裂带(F3)
图12 奎屯河东、西两岸Ⅱ级阶地上山前逆断裂带(a)及主要地貌面(b)的解译Fig.12 Interpretation of the piedmont reverse fault zone(a)and main geomorphological surface(b)on the terraceⅡalong the western and eastern banks of the Kuitun River.
图13 背斜东侧剪切断层及主要地貌面的解译结果Fig.13 Interpretation of strike-slip faults and main geomorphological surfaces on the eastern wing of the Dushanzi anticline.
后冲逆断裂带从独山子背斜东端哈里排特他乌沟口附近开始发育,然后向东延伸至依什克他乌沟。相对于其他2条断裂带,该断裂带构造几何形态简单,由14条EW向、相互平行的南倾断层崖组成,在平面上构成一向南突出弧形断层崖带。该构造在微地貌上形成一种形似搓板的棱状地形,朱海之等(1990)称之为搓板地形。该断裂带内各分支断层长短不一,并被近垂直的冲沟切割,但冲沟没有明显水平位错现象(图14)。在哈里排特他乌沟口附近,后冲逆断裂带与背斜山前逆断裂带存在重叠段落,后冲逆断裂发育在山前逆断裂上盘而形成南倾的断层崖地貌。
图14 后冲逆断裂带的断层及主要地貌面的解译结果Fig.14 Interpretation of back-thrust reverse fault zone and main geomorphological surfaces.
文中仅使用DEM可视化工具从不同虚拟的视角、不同色度来识别微构造地貌、划分地貌面和确定断层位置等。与前人通过航片解译和野外考察获得的断层填图成果(邓起东等,1991;Deng et al.,2000)比较,本文的LiDAR-DEM解译结果宏观上获得基本相同的断裂分布特征,但微观上较前者具有更高的精细程度,而且在地貌面划分上更加丰富。更为重要的是,基于LiDAR生成的DEM提供了一个具有准确地理参考系的数据库,因此任何以LiDAR数据为底图的解译成果都是具有可靠的坐标的数字化成果,便于数据共享。
基于LiDAR的高精度DEM的断层解译,是获得详细独山子逆冲断裂的几何展布和地形地貌分类的必要方法。尽管本文的结果是通过目视解译完成的,但是该结果是对构造地貌特征客观分析和评价获得的,结合DEM的可靠地理坐标系和高分辨率,这些构造地貌识别和断层解译结果是可靠的,可以指导野外调查工作。同时,DEM还可以提供构造地貌的垂直量信息,如陡坎高度、阶地面高差等,这些信息是无法从2维的航片上获取的。因此,利用高精度DEM进行活动断层识别和地貌面划分是地质填图中不可或缺的重要技术手段。除了基于DEM的活断层解译,作为数字化的基础3维数据,DEM在地形地貌量化分析方面也存在优势。今后工作应该利用GIS空间分析方法完成地貌特征分析获得断层几何分布,还可以进行地貌形态量化分析获得断层活动历史和特征。
5 结束语
利用机载LiDAR技术获得高精度、高分辨率的DEM数据,为研究和理解活动构造相关构造地质情况、地貌和最新的形变历史提供更精确的基础数据。如何将LiDAR新技术、新数据应用于活动构造填图工作和活动断层地震危险性评价等方面,无疑是今后活动构造研究领域的一个重要的发展方向。
在新疆天山北麓的独山子背斜-逆冲断裂带进行机载LiDAR数据采集,并从数据采集、数据质量检验、数据处理及数据应用等方面介绍基于LiDAR数据进行活动构造填图的方法和流程。通过数据质量检验,LiDAR原始数据的点云密度为6.6个/m2,平均点间距为0.39m,测量精度可达0.12m。通过对地面点云数据进行最佳分辨率评估,利用反距离权重算法可以获得0.5m分辨率的数据高程模型(DEM)。该分辨率的DEM数据足够精细到确定独山子逆冲断裂带的构造地貌特征并对它们进行解析。仅使用DEM可视化工具从不同虚拟的视角、不同色度或其他处理方式来识别微构造地貌、划分地貌面和确定断层位置等。由于LiDAR DEM具有可靠的地理坐标系和高分辨率,获得与前人通过航片解译和野外调查一致的断裂分布特征,但微观上较前者具有更高的精细程度,而且在地貌面划分上更加丰富。作为数字化的基础3维数据,今后工作还应该在地形地貌量化分析方面发挥DEM的优势。
陈桂华,徐锡伟,闻学泽,等.2006.数字航空摄影测量学方法在活动构造中的应用[J].地球科学-中国地质大学学报,31(2):405—410.
CHEN Gui-hua,XU Xi-wei,WEN Xue-ze,et al.2006.Application of digital aerophotogrammetry in active tectonics[J].Earth Science-Journal of China University of Geoscience,31(2):405—410(in Chinese).
邓起东,冯先岳,尤惠川,等.1991.新疆独山子-安集海活动逆断裂-褶皱带的变形特征及其形成机制[A].见:国家地震局地质研究所编.活动断裂研究(1).17—36.
DENG Qi-dong,FENG Xian-yue,YOU Hui-chuan,et al.1991.Characteristics and mechanism of deformation along the Dushanzi-Anjihai active reverse fault and fold zone,Xinjiang[A].In:Institute of Geology,SSB(ed). Research of Active Fault(1).17—36(in Chinese).
邓起东,冯先岳,张培震,等.2013.北天山山前和吐鲁番盆地中央隆起带活动逆断裂-褶皱带地质图(1∶50,000)[Z].北京:地震出版社.
DENG Qi-dong,FENG Xian-yue,ZHANG Pei-zhen,et al.2013.Geologic Map of Active Reverse Fault-Fold Zones in the Northern Piedmont of Tianshan and the Central Uplift Zone in the Turpan Basin[Z].Seismological Press,Beijing(in Chinese).
HE Hong-lin.2011.Some problems of aerial photo interpretation in active faultmapping[J].Seismology and Geology,33(4):938—950(in Chinese).
靳克强,龚志辉,汤志强,等.2011.机载LiDAR技术原理及其几点应用分析[J].测绘与空间地理信息,34(1):144—147.
JIN Ke-qiang,GONG Zhi-hui,TANG Zhi-qiang,et al.2011.Key principles and its applications of airborne LiDAR[J].Geomatics and Spatial Information Technology,34(1):144—147(in Chinese).
刘经南,许晓东,张小红,等.2008.机载激光扫描测高数据分层迭代选权滤波方法及其质量评价[J].武汉大学学报-信息科学版,33(6):551—555.
LIU Jing-nan,XU Xiao-dong,ZHANG Xiao-hong,et al.2008.Adaptive hierarchical and weighted iterative filtering of airborne LiDAR data and its quality assessment[J].Geomatics and Information Science ofWuhan University,33(6):551—555(in Chinese).
隋立春,张熠斌,柳艳,等.2010.基于改进的数学形态学算法的LiDAR点云数据滤波[J].测绘学报,39(4):390—396.
SUI Li-chun,ZHANG Yi-bin,LIU Yan,et al.2010.Filtering of airborne LiDAR point cloud data based on the adaptive mathematical morphology[J].Acta Geodaetica et Cartographica Sinica,39(4):390—396(in Chinese).
王明华,张小红,曾涛,等.2010.机载LiDAR数据滤波预处理方法研究[J].武汉大学学报,35(2):224—227.
WANG M ing-hua,ZHANG Xiao-hong,ZENG Tao,et al.2010.Preprocessing algorithms for filtering airborne LiDAR data[J].Geomatics and Information Science ofWuhan University,35(2):224—227(in Chinese).
王晓辉,胡伍生,刘行波.2007.LiDAR系统测量成果精度检测[J].测绘工程,16(3):67—70.
WANG Xiao-hui,HUWu-sheng,LIU Xing-bo.2007.Precision tests of LiDAR products[J].Engineering of Surveying and Mapping,16(3):67—70(in Chinese).
张玲,杨晓平,魏占玉,等.2014.三维数据的二维可视化方法综述[J].地震地质,36(1):275—284.doi:10.3969/j.issn.0253-4967.2014.02.023.
ZHANG Ling,YANG Xiao-ping,WEI Zhan-yu,et al.2014.Overview of visualization methods of three dimensional topographic data[J].Seismology and Geology,36(1):275—284(in Chinese).
赵礼剑,程新文,李英成,等.2009.记载LIDAR点云高程数据精度检核及误差来源分析[J].地理空间信息,7(1):58—60.
ZHAO Li-jian,CHENG Xin-wen,LIYing-cheng,et al.2009.Accuracy check and analysis of LiDAR elevation data[J].Geospatial Information,7(1):58—60(in Chinese).
Adams T,Brack C,Farrier T,et al.2011.So you want to use LiDAR?An guide on how to use LiDAR in forestry[J].New Zeal J For,55:19—23.
Arrowsmith JR,Zielke O.2008.Slip along the San Andreas Fault associated with the great 1857 earthquake derived from“B4”LiDAR high resolution topographic data[J].Eos Transactions AGU 89(53):Fall Meet Suppl,Abstract T44B-03.
Arrowsmith J R,Zielke O.2009.Tectonic geomorphology of the San Andreas Fault zone,from high resolution topography:An example from the Cholame segment[J].Geomorphology,113:70—81.http:∥dx.doi.org/10.101 6/j.geomorph.2009.01.002.
Axelsson P.1999.DEM generation from laser scanner data using adaptive TIN model[J].International Archives of the Photoframmetry,Remote Sensing and Spatial Information Sciences,33(B4/1):110—117.
Begg J G,Mouslopoulou V.2009.Analysis of late Holocene faulting within an active rift using LiDAR,Taupo Rift[J].New Zeal J Volcanol Geotherm Res,190:152—167.http:∥dx.doi.org/10.101 6/j.jvolgeores. 2009.06.001.
Chiba T,Kaneda S,Suzuki Y.2008.Red Relief Image Map:New visualization method for three dimeasional data[J]. The International Archives of the Photogrammetry,Remote Sensing and Spatial Information Science,2008,Beijing:ISPRS Vol.XXXVII.Part B2,1071—1076.
Crosby C J,Arrowsmith JR.2007.App lication of LiDAR data to constraining a late Pleistocene slip rate and vertical deformation of the Northern San Andreas Fault,Fort Ross to Mendocino,California[D].US,Arizona State University.
Cunningham D,Grebby S,Tansey K,et al.2006.Application of airborne LiDAR to mapping seismogenic faults in forested mountainous terrain,southeastern Alps,Slovenia[J].Geophysical Research Letter,33:L20308.http:∥dx.doi.org/10.1029/2006GL027014.
Hu Y.2003.Automated extraction of digital terrain models,roads and buildings using airborne LiDAR Data[D]. Department of Geomatics Engineering,The University of Calgary,Calgary,Alberta,Canada.
Kim H,Arrowsmith JR,Crosby C J,et al.2006.An efficient implementation of a local binning algorithm for digital elevation model generation of LiDAR/ALSM dataset[J].Eos Trans,AGU,87(52):Fall Meet Suppl,Abstract G53C-0921.
Kondo H,Toda S,Okamura K,et al.2008.A fault scarp in an urban area identified by LiDAR survey:A case study on the Itoigawa-Shizuoka Tectonic Line,central Japan[J].Geomorphology,101:731—739.
Kraus K.1998.Determination of terrain models in wooded areaswith airborne laser scanner data[J].ISPRS Journal of Photogrammetry and Remote Sensing,53:193—203.
Langridge R M,Ries W F,Farrier T,et al.2013.Developing sub 5-m LiDAR DEMs for forested sections of the Alpine and Hope Faults,South Island,New Zealand:Implications for structural interpretations[J].Journal of Structural Geology,1—14.
Lin Z,Kaneda H,Mukoyama S,et al.2013.Detection of subtle tectonic geomorphic features in densely forested mountains by very high-resolution airborne LiDAR survey[J].Geomorphology,2013:104—115.
Molnar P,Brown E T,Buchfie B C,et al.1994.Quaternary climate change and the formation of river terraces across growing anticlines on the north flank of the Tian Shan,China[J].The JGeology,102:583—602.
Perron J T,Kirchner JW,Dietrich W E.2008.Spectral signatures of characteristic spatial scales and nonfractal structure in landscapes[J].Journal of Geophysical Research:Earth Surface,113(F4).
Zhang K,Chen S,Whitman D,et al.2003.A progressive morphological filter for removing non-ground measurements from airborne LiDAR data[J].IEEE Transactions on Geoscience and Remote Sensing,41(4):872—882.
EXPERIMENTAL STUDY ON GEOLOGIC MAPPING OF ACTIVE TECTONICS BASED ON LIDAR DATA—A CASE OF DUSHANZI
ANTICLINE-REVERSE FAULT ZONE IN XINJIANG
WEIZhan-yu HE Hong-lin GAO Wei XU Xi-wei GAN Wei-jun WEI Lei-hua
(Key Laboratory of Active Tectonics and Volcano,Institute of Geology,China Earthquake Adm inistration,Beijing 100029,China)
Airborne LiDAR(Light Detection And Ranging)provides a more advanced technique and more accurate basic data to describe geomorphological features and the latest surface deformation associated with active tectonics.How to app ly this new technique and dataset to mapping of active fault and seismic hazard assessment is an important trend in the field of active tectonics.Taking the Dushanzi anticline-reverse fault zone in Xinjiang as test area,we made an experimental study on geologic mapping of active tectonics based on the LiDAR data.Firstly,we collected raw data using the airborne LiDAR technique,and obtained a raw point-cloud with a point density of 6.6 points/m2and an average space of 0.39m between any two points.Secondly,using twelve ground control points(GCP)which is acquired by static GPSmeasurementwith accuracy up tom illimeter,we evaluated the vertical error of the ground point-cloud data with density of 6.4 points/m2,and the result shows a vertical error of 0.12m,mean square value 0.078m.Finally,using the inverse distance weighting algorithm,we obtained the digital elevation model(DEM)of 0.5m-resolution.The resolution of the DEM is high enough to describe and analyze spatially the fine feature of tectonic landform of the Dushanzi anticline-reverse fault zone.In this paper,we identify the fine tectonic landforms using merely the DEM visualization tools based on different virtual perspectives,different shades or different treatmentmethods.The active tectonics and their distribution identified based on the high resolution DEM derived from LiDAR are not only consistent with previous results identified from airinterpretation and field investigation,but also finer and more precise than the latter.In addition,these methods of data acquisition,quality inspection and data processing introduced in this paper are also applied to other active fault researches in which LiDAR data have been acquired.
LiDAR,Dushanzi anticline-reverse fault zone,mapping of active fault
P315.2
A
0253-4967(2014)03-0794-20
林.2011.活动断层填图中的航片解译问题[J].地震地质,33(4):938—950.
10.3969/j.issn.0253-4967.2011.04.017.
魏占玉,男,1981年生,2010年在中国地震局地质研究所获得构造地质学博士学位,现主要研究方向为活动构造与构造地貌,电话:010-62009031,E-mail:weizhanyu@gmail.com。
doi:10.3969/j.issn.0253-4967.2014.03.019
2014-03-10收稿,2014-07-17改回。
我国地震重点监视防御区活动断层地震危险性评价项目(201308001)资助。