高精度LiDAR 技术在活动断层断错地貌研究中的应用
——以沂沭断裂带莒县至郯城段为例
2024-01-09王华林吴洪斌葛孚刚苏思丽
杨 彬,王华林,吴洪斌,葛孚刚,邹 昊,苏思丽
(山东省地震局,山东 济南 250014)
激光雷达(Light Detection and Ranging,LiDAR)测量技术是近20 年来快速发展的一种新型对地观测手段,是对地貌高程的最直接且精度最高的测量手段,可精确、快速、可靠地获取地面三维数据,在地球科学中的应用日益广泛[1-4]。其最具吸引力的优势在于激光束能够通过植被缝隙到达地表,并能够准确测量地表形态。这些数据可以在适当的尺度和精度上展现断层活动在地表留下的迹象,现已广泛应用于断错地貌定量分析,可以准确地提取和分析活动断裂的构造地貌、几何形态与分布特征[5]。
早在20 世纪90 年代LiDAR 测量技术已开始应用于地形测量,随后该技术得到了空前的发展。美国开展了一系列LiDAR 测量项目,系统采集了包括圣安德烈斯断层南段和圣哈辛托断层的三维地形数据,为南加州地区活动断层研究提供了高精度的基础资料[6-7]。欧洲、日本等也均已经开展沿主要活动断裂带的大规模机载LiDAR 测量。
在国内,虽然LiDAR 技术在其他行业的应用相对较早,而在地形地貌、活动断层的微地貌测量方面起步较晚。2008 年汶川地震时曾开展了唐家山堰塞湖的应急机载激光雷达测量探索研究,但部分技术手段和测量信息并没有完全发挥作用[8-9]。之后LiDAR 技术开始应用于断错微地貌的测量[10-13]。2011 年中国地震局地质研究所沿海原断裂带1920 年海原81/2 级地震120 km 长的地震地表破裂带利用机载激光雷达技术开展了条带状测量,对海原活动断裂带进行了详细的地质地貌研究[14-15],是国内机载激光雷达技术在活动构造研究中的首次应用。
基于LiDAR 数据的量化分析是未来活动构造研究的趋势,与传统方法相比,LiDAR 技术在森林覆盖区和城区的活动断裂填图中具有巨大的优势,在沿活动断裂位移量测量上也更精准,更有效。而震前与震后LiDAR 数据对比也是研究同震变形特征、探索断裂发震模式的重要手段[16]。LiDAR 技术是获取活动断层位移量等定量参数的最有效的新技术手段,在活动构造研究中的应用前景广阔;沿主要活动断裂带的LiDAR 测量将成为未来国内活动断裂研究基础数据获取的重要手段。
郯庐断裂带是我国东部地区最为重要的大型活动构造带。沂沭断裂带位于郯庐断裂带中南段,是断裂和地震活动最为活跃的段落。沂沭断裂带内安丘−莒县断裂是全新世活动断层,沿断裂发生过1668 年郯城8½级地震和公元前安丘7 级地震。该断裂的莒县至郯城段断错地貌特征最清晰,断裂沿线的断层陡坎、冲沟水系错断等现象丰富,是利用LiDAR 技术开展断层测量、获取断层活动定量数据的理想地段。考虑研究地区植被覆盖,采用较为前沿的地基激光雷达与机载激光雷达对沂沭断裂带莒县至郯城段的左山段、钟华山段、岌山段和马陵山段开展了活动断层断错地貌测量,通过LiDAR 点云数据的获取和分析处理,得到了高精度数字高程模型 (Digital Elevation Model,DEM)和定量数据,实现了断裂微地貌形态的高清晰度三维再现,以此更科学、精准地研究了沂沭断裂带莒县至郯城段典型断错地貌和活动特征。
1 构造背景与研究方法
1.1 构造背景
沂沭断裂带规模大,活动性强,具有分段活动特征。自北向南,将沂沭断裂带分为Ⅰ安丘段;Ⅱ莒县−郯城段;Ⅲ新沂−泗洪段(图1),莒县−郯城段是沂沭断裂带中最为活动的段落[17]。断裂带由若干条次级断裂组成,形成“两堑夹一垒”的构造地貌,略呈平行排列,自西向东依次为鄌郚−葛沟断裂、沂水−汤头断裂、白芬子−浮来山断裂、安丘−莒县断裂、昌邑−大店断裂;沂沭断裂带莒县至郯城段长约100 km,表现为4 条右阶斜列的次级断层,分别为左山次级断层、钟华山次级断层、岌山次级断层和马陵山次级断层,并且不连续出露。野外实地考察多处见东侧的白垩系王氏组紫红色砂岩逆冲在西侧的晚第四纪黄土层之上,并错段冲沟水系。1668 年郯城8½级地震就发生在莒县至郯城段[17],遥感影像上极为醒目,控制着沭河水系的发育,形成刀砍般的线性陡坎地貌,反映了该断裂晚第四纪以来的走滑活动性质。
图1 沂沭断裂带莒县至郯城段遥感综合解译与激光雷达(LiDAR)测量分布Fig.1 Maps showing the comprehensive interpretation results of remote sensing images and LiDAR measurement distribution for the Juxian-Tancheng segment of the Yishu fault zone
1.2 构造地貌研究方法
利用单次地震位移与多次地震累积位移的关系可以研究强震的复发习性,所以分析断层的破裂历史是预测地震的有效方法之一[18]。遥感与LiDAR 技术的应用,使得高精度地形地貌测量更加快捷,如高分辨率影像以及DEM 在活断层定量研究中的应用[18-19];LiDAR技术在构造地貌测量中的应用[20-21]。利用高分辨率遥感与高精度地形资料,约束构造地貌位错,重构断层破裂历史,识别古地震事件,明确地震位移分布特征与断层破裂模式,从而估算断层强震复发周期。通过定量分析相关断错地貌参数来开展活动断层构造地貌研究。
本次野外测绘工具分别为奥地利RIEGL 公司研发的VZ-1000 地基激光雷达扫描仪和国产纵横机载JoLiDAR1500 激光雷达;其中VZ-1000 地基激光雷达定位误差仅为0.005 m,纵横机载JoLiDAR1500 激光雷达飞行高度500 m 范围内,定位误差为平面0.01 m,高程0.03 m。经过数据处理后可获得高分辨率的DEM 数据。
冲沟位错可利用基于Matlab 平台开发的LaDiCaoz计算分析软件[22]。该软件是一款专门针对河流及冲沟位错分析测量的专业分析软件,用户可将LiDAR 等设备获取的高精度DEM(图2c)数据导入软件,通过定义断裂迹线位置、地形横断面位置和地貌标志的方向(冲沟、河道或陡坎方向,图2a)获取断裂位移数据。具体操作方法如下:
图2 基于高分辨率DEM 数据的位错量测量方法Fig.2 Method for measuring offset amount based on high-resolution DEM data
(1) 定义断层的位置以及跨冲沟上游和冲沟下游的地形横断面位置(图2a)。
(2) 获取上、下游冲沟横断面的高程剖面信息。
(3) 将红色的剖面形态进行水平移动、垂直移动和拉伸调整(图2d、图2e),最终求得这2 条剖面的高程差绝对值之和的最小值,其对应的水平移动距离即为最佳位移值(图2f)。
(4) 输入得到的最佳位错值,软件即会给出该处直观的原始地貌(图2g)。
经与野外实际测量值比对可知,软件得出的结果在误差范围内与实际相符。
2 地形地貌数据采集与处理
2.1 地基激光雷达数据采集
本文通过实际操作和后期数据的处理及分析,总结出地基激光雷达操作步骤主要分为三个部分(图3)。一是前期准备阶段,其中最主要的任务就是结合野外实际地形,合理制定扫描方案和站点的布设;二是数据采集的扫描阶段;三是数据处理阶段。
图3 地基LiDAR 工作流程Fig.3 Workflow of terrestrial LiDAR
本文采用的仪器RIEGL VZ-1000 激光扫描仪,主要用于采集沿断层典型地貌位错处的三维点云数据,其系统是由激光发射器、纹理采集相机及GPS 天线组成。其主要参数见表1。
表1 RIEGL VZ-1000 主要技术参数Table 1 Primary technical parameters of a RIEGL VZ-1000 terrestrial laser scanner
扫描时主要将测量站点布设在断层破裂带两侧位置上,确保完整采集断层地貌细节。在进行野外勘查时必须确定站点的整体布设情况,需保证相邻站点之间有较高的重叠区以利于后续点云拼接处理,扫描范围要基本涵盖区域内的详尽的地形特征。其中每站扫描激光频率为150 Hz,最远距离为1 200 m,每站扫描的时间5~10 min。图4 为沂沭断裂带莒县至郯城段地貌特征及扫描站点位置影像图。
图4 沂沭断裂带莒县至郯城段地貌特征及扫描站点位置影像图Fig.4 Images showing the geomorphic features and measurement station locations of the Juxian-Tancheng segment of the Yishu fault zone
2.2 机载激光雷达数据采集
无人机机载激光雷达是由纵横CW-25E 固定翼飞机搭载PPK/RTK 差分定位系统(动态后处理差分技术/载波相位差分技术)以及JoLiDAR1500 激光雷达,作业半径达50 km,抗干扰能力强,实现飞机厘米级定位以及双差分定向。
航线首尾会通过绕“8”字的形式,对激光雷达进行雷达惯导标定。无人机飞行时,实时采集惯性传感器IMU 和激光雷达的数据,利用在线标定算法,实时估计IMU 和激光雷达之间的变换关系,并将激光雷达的点云数据和IMU 的惯性测量数据进行融合,以提高数据的定位精度(表2)。
表2 JoLiDAR1500 主要技术参数Table 2 Primary technical parameters of the JoLiDAR1500 measurement system
由于现场山体高耸、落差大,在保证飞行安全的前提下,定高飞行方式获取的数据只能保证最高处的扫描精度,山谷盆地等低处扫描精度会降低。因此,本次无人机激光雷达数据获取采用仿地变高飞行方式,兼顾飞行安全和数据精度,仿地飞行高度500 m。
2.3 数据处理
与光学测量以及其他传统测量方法相比,LiDAR技术最大的优势就是可以去除植被影响。LiDAR 数据包含多重回波信息,可以将末次回波信息保留(一般为最低高程数据),从而更真实地呈现原始地貌信息(图5a 未滤除植被DEM、图5b 滤除植被DEM)。
图5 滤除植被前后的DEM 对比Fig.5 Comparison of the DEM before and after vegetation filtering
LiDAR 数据处理包括三维坐标解算、点云滤波和DEM 生成三个过程。使用Inertial Explorer 软件和点云处理软件进行点云三维坐标解算。导入LiDAR原始点云数据和轨迹数据,计算点云的三维坐标。最后解算得到输出为las 格式点云数据(图6)。
图6 岌山段激光雷达点云解算结果Fig.6 LiDAR point cloud calculation results of the Jishan segment of the Yishu fault zone
由于LiDAR 采集到的点云存在离散和非均匀特性,测区内会存在一些“空洞”,常用的DEM 一般为规格格网,因此,需要通过插值生成DEM。空间插值的常见算法为自然邻域插值、反距离加权法(IDW)、克里金插值法、不规则三角网法等,本研究使用不规则三角网法,最终得到的高精度的DEM 数据。
2.4 地基与机载激光雷达的差异以及数据使用情况
从工作原理上,地基激光雷达与机载激光雷达是一致的,但二者在工作方式上存在差异,工作效率和采集数据的精度也不尽相同:
(1) 地基激光雷达采用的是固定站点式数据采集,需要人工大范围布设站点,工作效率较低,但采集数据的精度高,例如本文VZ-1000 地基激光雷达定位误差为0.005 m。
(2) 机载激光雷达采用的是移动式数据采集,由无人机携带激光雷达在空中飞行,通过持续不断的空中惯导标定,来采集高精度的定位数据,工作效率高,但定位误差要大于地基激光雷达。例如本文纵横CW-25E 搭载的JoLiDAR1500 激光雷达,航高500 m 范围内定位误差为平面0.01 m,高程0.03 m。
概括起来,地基激光雷达和机载激光雷达虽然采集方式、工作效率以及获取的数据精度不同,但数据处理的方式包括数据匹配、点云抽稀、点云滤波、三维建模等过程是一致的,同时都可获得满足工作要求的高分辨率DEM,本次工作左山(一步涧)段、钟华山段以及蒋家岭等地采用地基激光雷达完成野外数据采集,岌山段和马陵山段采用机载激光雷达完成野外数据采集,通过数据处理均可得到的精度为0.1 m 高分辨率的DEM 数据。
3 断错地貌测量分析与讨论
沂沭断裂带莒县至郯城段断错地貌非常发育,由于断裂活动的非均匀性,断错地貌在不同断裂段上的类型及强度都有较大差异。沂沭断裂带莒县至郯城段具有右旋走滑兼逆断层性质,沿断裂带断层陡坎及冲沟位错发育明显。本文通过地基激光雷达(LiDAR)和机载激光雷达(LiDAR)分别获取了莒县至郯城段左山(一步涧)、钟华山、岌山和马陵山等断层段以及蒋家岭等5 处断错地貌典型地段的高精度DEM,对每个地段均进行了右旋冲沟解译,根据高分辨率DEM 解译得到断层断错冲沟的情况分析,断层断错冲沟水平位移的类型(图7)大致可分为以下几种。
图7 地震断层的冲沟位错类型示意Fig.7 Schematic diagram of offset types of earthquake fault-controlled gullies
a 类型,2 次或多次断头冲沟和扭曲冲沟,冲沟的扭曲量代表最新一次水平位错量;b 类型,断层两侧A−A冲沟间的距离代表断层的水平位移量;c 类型,断头冲沟与取直冲沟间的距离为断层的水平位移量,若存在多条断头冲沟说明断层具有多次错动;d 类型,断层一盘冲沟平直另一盘发生扭曲代表AA'间的距离为断层水平错动量;e 类型,断层一盘冲沟平直另一盘冲沟发生多次扭曲,每次扭曲代表一次错动的水平位移,断层总位移为多级扭曲量的总和。
利用冲沟位错分析测量软件LaDiCaoz 获得了71 条冲沟的右旋位错量(表3),对左山(一步涧)、钟华山、岌山和马陵山等断层段水平位移特征以及钟华山和蒋家岭2 段的垂直位移特征进行了分析讨论。
表3 冲沟位错统计Table 3 Statistics of gully offsets
3.1 冲沟的右旋水平位错
3.1.1 左山(一步涧)段水平位错
左山(一步涧)段的断裂走向总体为NNE,该段断层破碎带较宽,冲沟右旋位错最典型,该处共量测获得16 处冲沟位错量。左山北段(图8a)位错量分为4 级:7.1、16.1、29.1 和39.2 m,左山南段(图8b)位错量分为3 级:3.1~3.2、6.2~8.8 和13.2~15.9 m。其中在左山(一步涧)北同一冲沟显示4 级错动,显示该段活动断层是一条多次破裂的地震断层,全新世时期除1668 年8½级地震外,可能还有3 次古地震事件发生,这一结论在开挖的左山(一步涧)探槽中得到了很好的验证。
图8 左山(一步涧)段实测地形及冲沟位错解译Fig.8 Measured topography and interpreted gully offsets in the Zuoshan (Yibujian) segment
3.1.2 钟华山段水平位错
钟华山段紧邻沭河,整体走向NNE,呈现东高西低的格局,断层线平直,线性特征明显。该段共获得14 条冲沟的右旋水平位错量,位错量分为3 级:3.2~5.2、9.5~12.5、16.8~18.6 m(图9a、图9b、图9c),同样说明断层具有多期活动,发生多次古地震事件。
图9 钟华山段实测地形及冲沟位错解译Fig.9 Measured topography and interpreted gully offsets in the Zhonghuashan segment of the Yishu fault zone
3.1.3 岌山、马陵山段水平位错
岌山段和马陵山段是控制岌山和马陵山西边界的断层,主要表现为逆冲性质,倾向南东,断层线总体呈弧形。该段共获得冲沟位错量41 条,其中岌山段20 条,马陵山段21 条。岌山段显示冲沟3.6~5.8、9.8~13.8、18.1~23.8、32.2 和44.1 m 等5 个不同量级的右旋位错,马陵山段也显示了冲沟4.7~8.2、12.5~18.0、25.0~33.0 和51.0~55.0 m 等4 个不同量级的右旋位错。其中岌山段最南侧一条冲沟(图10a)和马陵山段的2 条冲沟(图10b) 均显示有3 级错动,显示了断裂的多期活动。
图10 岌山和马陵山段实测地形及冲沟位错解译Fig.10 Measured topography and interpreted gully offsets in the Jishan and Maolingshan segments of the Yishu fault zone
3.1.4 分析与讨论
地震破裂产生的右旋走滑位错可以被冲沟等地貌保留,而多次地震产生的累计位移也可被冲沟等地貌记录下来,因此,通过冲沟水平位错量分析,可以初步确定断裂的活动期次和活动强度。通过沿断裂带左山(一步涧)段、钟华山段、岌山段和马陵山段的冲沟位错数据统计(表3)分析,初步研究成果表明:①每一段的冲沟水平位错量均存在大致相当的分级特征,说明沂沭断裂带莒县至郯城段晚第四纪以来,可能发生过多期(3~5 次)活动和多次古地震事件;②不同分段间冲沟的水平位错量级存在差异,说明不同段落位错量的差异性,可能与断裂的运动学性质差异性、断层面结构特征、冲沟所处的地质构造部位和地形地貌位置以及冲沟形成过程的等诸多因素相关。因此,冲沟水平位错量的差异性可作为沂沭断裂带郯城至莒县段断层分段的依据。③本次测量所得的最小一级的冲沟右旋水平位错量从一定意义上代表了郯城1668 年8½级地震同震位错量。左山(一步涧)段、钟华山段、岌山段、马陵山段最小一级位错量的最大值分别为7.1、5.2、5.8 和8.2 m。左山北段冲沟最小一级位错量为7.1 m,同时在左山探槽揭露到1668 年郯城8½级地震的最新活动[23-24]。说明冲沟7.1 m 右旋水平位错量作为1668 年郯城8½级地震的位错量是可信的。
综合分析认为,4 个段落最小一级位错量为1668年郯城8½级地震的同震位错量,位错量在5.2~8.2 m之间。这一认识与前人研究结果相吻合[23]。
3.2 垂直位错特征
断层垂直位移的标志主要有冲沟、陡坎和地层,沂沭断裂带莒县至郯城段垂向断错地貌特征主要体现在钟华山、蒋家岭的线性多级坡度陡坎以及左山冲沟的垂直位错。钟华山段总体走向NNE,断裂线性特征明显,发育不同级别的坡度陡坎,剖面1(图11)处累积陡坎高度7.0 m,剖面中可见3 级坡度陡坎,每级坡度陡坎1.5~3.0 m,剖面2(图11)累积陡坎高度5.0 m,剖面中可见2 级坡度陡坎,为多次活动的复合陡坎。蒋家岭段总体走向NNE,2 个剖面发育3 级坡度陡坎,剖面1(图12) 处累积陡坎高度7.5 m,每级坡度陡坎1.5~2.5 m;剖面2(图12)处累积陡坎高度6.8 m,每级坡度陡坎1.5~2.5 m。左山(一步涧)段冲沟右旋位错明显,并伴有4.0 m 的垂直位错(图2g)。累积陡坎高度反映断裂垂直位移情况,多级陡坎反映了断裂的多期活动和多次古地震事件。
图11 钟华山南段断层陡坎DEM 及其剖面图Fig.11 DEM and profiles of fault scarps in the southern Zhonghuashan segment of the Yishu fault zone
表4 给出了断层陡坎剖面反映的垂直位错情况。断层陡坎的多级坡度中断代表断层的多次垂直位移,根据地基激光雷达测得的钟华山和蒋家岭断层陡坎坡度分级位错,同样概括起来可以推测出3 次古地震事件,每次古地震事件垂直位错量为1.5~2.5 m,与1668 年郯城8½级地震的垂直位错量相当。
表4 断层陡坎剖面反映的垂直位错统计Table 4 Statistics of vertical offsets on fault scarp profiles
4 地质构造检验
本次LiDAR 测量成果在山东省防震减灾“十三五”项目临沂市国际生态城地震断层探测与地震危险性评价项目(山东省震灾风险防治中心,2023)中得到了检验[24]。依据本次测量结果对左山(一步涧)段和钟华山段的断层陡坎和冲沟位错进行了地质填图和探槽开挖,取得良好效果。
4.1 左山(一步涧)段地质构造检验
在左山(一步涧)段的断层陡坎处开展了探槽开挖工作(图13)。探槽位置处于带状陡坎附近,探槽以北为人工挖除的第四系最新沉积,探槽以南可见冲沟右旋错动和沿断层方向的陡坎等。通过与文献[25](林伟凡等,1987)在左山探槽揭露的地层对比,认为②褐黄灰黑混杂的含砾亚砂土14C 测定年龄值为230 士85 年,③灰黑色坡积砂砾石层14C 测定年龄值为3 960士395 年,④砂质亚黏土14C 测定年龄值为(1.2~3)万年[24]。在左山(一步涧)探槽(T2) 中发育4 组断面、3个构造楔和崩积楔。根据探槽剖面地层时代、断层错断的地层层位、构造楔和崩积楔等分析,本探槽内揭露出4 次地震事件(包括1668 年郯城8½级),显示了断层的多期活动,断层性质为右旋正断走滑断层。
图13 左山(一步涧)探槽(T2)剖面素描Fig.13 Sketch showing the profile of trench T2 in the Zuoshan (Yibujian) segment of the Yishu fault zone
通过地质填图获得了左山(一步涧)段冲沟位移数据,显示3~7、13~16、25~28 和35~39 m 的分级位移特征,说明断层具有多期活动,发生多次古地震事件。该结论与通过LiDAR 获取冲沟位错特征基本一致。
4.2 钟华山段地质构造检验
在钟华山段的断层陡坎开展探槽开挖(图14)。探槽所处位置地形东高西低,探槽南侧见冲沟右旋水平错动,位错量约4 m。探槽揭露,断层东盘为白垩系灰绿-黄绿的泥岩、砂岩,西盘为河流Ⅰ级阶地。断层东盘向上逆冲,错断了阶地底部的⑦棕褐色粉土和⑥棕褐色亚黏土−黏土,根据与区域地层对比,⑥棕褐色亚黏土−黏土属全新统。分析认为,探槽揭露的断层为全新世活动断层,性质为右旋逆断走滑断层。
图14 钟华山探槽(T1)剖面图Fig.14 Profile of trench T1 in the Zhonghuashan segment of the Yishu fault zone
对通过地质填图获得的钟华山段冲沟位移数据,显示了4、8~10 和16~19 m 的分级位移特征,同样说明断层具有多期活动,发生多次古地震事件。
5 结论
a.利用地基激光雷达和机载激光雷达技术获取了沂沭断裂带莒县至郯城段左山(一步涧)、钟华山、岌山和马陵山等多处典型断错地貌的高精度DEM 数据,实现了对断层陡坎和冲沟位错的详细定量分析研究,真实地再现了沂沭断裂带莒县至郯城段的典型断错地貌特征,获得了断裂水平和垂直位移数据,为断裂运动学和几何学研究提供了重要依据。
b.沂沭断裂带在左山(一步涧)段、钟华山段、岌山段和马陵山段的水平断错地貌具有分级特征。左山段位错量分为4 级:7.1、16.1、29.1 和39.2 m。钟华山段位错量分为3 级:3.2~5.2、9.5~12.5、16.8~18.6 m。岌山段位错量分为5 级:3.6~5.8、9.8~13.8、18.1~23.8、32.2 和44.1 m。马陵山段位错量分为4 级:4.7~8.2、12.5~18.0、25.0~33.0 和51.0~55.0 m。初步分析认为,左山(一步涧)段、钟华山段、岌山段和马陵山段晚第四纪以来可能发生过多期(3~5 次)活动和多次古地震事件,四个段落最小一级位错量为1668 年郯城8½级地震的同震位错量,位错量在5.2~8.2 m。
c.根据地基激光雷达测得的钟华山和蒋家岭断层陡坎坡度分级位错,获得断裂累计垂直位错量为5~7.5 m,一次断裂活动位错量为1.5~2.5 m。初步分析出3 次古地震事件,每次古地震事件位错量1.5~2.5 m,与1668 年郯城8½级地震的位错量相当。
d.通过地质填图和探槽开挖的地质构造检验,认为地基激光雷达和机载激光雷达技术在活动构造研究中有着广泛的应用前景。