福建地区上部地壳结构深地震测深探测研究

李 培，任丛荣，李海艳,2，王善雄

(1.福建省地震局，福建 福州 350001；2.中国地震局厦门海洋地震研究所，福建 厦门 361116)

福建地处中国大陆东南缘，位于东南沿海地震活动带的北段，因毗邻台湾现代板块活动边界构造带，是新构造活动较为强烈的构造变形和强震活动带[1]。受欧亚板块、太平洋板块、菲律宾海板块相互作用影响，基于断裂交汇及深部热物质的交换，矿产、地热资源丰富、温泉众多[2]，中–新生代岩浆岩广泛出露，为全球构造岩浆活动最活跃的地区之一[3]。福建省处于华南褶皱系东部，地层发育，自上元古界至第四系均有出露。岩石类型繁多，沉积地层和变质地层的总和与燕山期火山岩地层约各占全省陆地总面积的三分之一。由于构造运动的多期性和构造演化的多阶段性，造成了各断代地层及变质岩、火山岩、侵入岩、构造变形带在时空分布上具有较直观的区带特征。福建地质构造活动划分为扬子和加里东、华力西和印支、燕山、喜马拉雅等构造活动期。受华南大陆的多期次构造运动的影响，福建地区断裂构造复杂，主要发育NW 向和NE向2 组断裂，构造地貌形成了“东西分带，南北分块”的格局[4]。根据构造形变、变质作用、沉积建造、岩浆岩活动等，从西到东划分为闽西北隆起带、闽西南拗陷带及闽东火山断拗带(图1)。研究认为，NW 向和NE 向两组断裂控制着福建省新构造运动和第四纪地貌的发育，具有多期次和不同活动性质等特点。因此，开展深地震测深探测研究工作，获得该地区地壳结构，揭示不同构造单元内部包括沉积盖层、结晶基底、地壳的深浅结构与断裂构造特征是阐明福建地区构造形态及其构造演化的关键和基础[5]。

近些年来为研究福建地区的深部结构以及地震孕育的深部构造背景，开展了主动源和被动源的地震探测研究工作并且取得了一定的成果，包括利用远震接收函数研究地壳厚度和泊松比[6-7]，利用主动源地震探测资料建立福建地区速度结构模型[8-12]，利用被动源地震资料建立福建地区三维地壳结构模型[13-18]。这些研究主要针对地壳深部速度结构，较少涉及上地壳上部，尤其是结晶基底与沉积盖层的浅部速度结构与构造特征。2010−2012 年，福建省地震局联合中国地震局地球物理勘探中心在福建地区开展了4 条NW 向和4 条NE 向交叉测线组成的人工爆破深地震测深地震观测，李培等[12]给出了其中4 条NW 向测线的基底结构，分析了测线穿过的NE 向断层的构造特征。但是，横穿NW 向断裂的4 条NE 向测线的上部地壳速度结构、构造特征以及断裂的深部形态等都还不清楚，速度结构与地质构造单元、构造运动、岩性、地热、矿产等的相互关系还需要进一步分析。因此，开展4 条NE 向测线的上部地壳精细结构成像研究，对进一步研究东南缘深部构造关系、断裂带的浅部结构和形成机制等提供地震学证据。

笔者利用福建地区4 条NE 向深地震测深测线的基底折射波(Pg 波)数据，使用有限差分层析成像方法获得了4 条剖面的上部地壳速度精细结构，并通过走时残差、射线数分布、交叉点速度控制等分析了结果的可靠程度。在分析剖面结构对应的构造和上部地壳变形特征的基础上，揭示了测线穿过的NW 向断裂的深部形态和研究区“东西分带，南北分块”的构造格局差异，为进一步研究东南缘结晶基底的构造属性、构造活动的继承性和差异性提供科学依据。

1 数据采集与初至Pg 波波组特征

1.1 数据采集

福建地区4 条NE 向测线和4 条NW 向测线交叉布设，其中4 条NE 向深地震测深测线(L5、L6、L7、L8)沿SW–NE 方向展布，与NW 走向的闽江、沙县−南日岛、永安–晋江、九龙江、上杭–云霄断裂近垂直(图1)。其中，武平–宁化–邵武(L5)剖面位于闽赣交界，南起武平，经长汀、宁化至邵武，全长约340 km，沿剖面布设4 个药量2.5～3.0 t 的爆破点和130 台三分量数字地震仪，观测点距2.5～2.8 km。上杭–永安–建瓯(L6)剖面位于研究区西侧，全长约320 km，沿剖面布设3 个药量2.3～2.6 t 的爆破点和100 台三分量数字地震仪，观测点距3.0～3.4 km。永定–大田–古田(L7)剖面位于研究区中部，经漳平、大田至古田，全长约305 km，沿剖面布设4 个药量2.0～2.8 t 的爆破点和97 台三分量数字地震仪，观测点距3.0～3.2 km。漳浦–南安–罗源(L8)剖面位于闽东沿海，南起诏安，经漳浦、长泰、南安至罗源，全长约360 km，沿剖面布设4 个药量1.8～2.6 t 的爆破点和118 台三分量数字地震仪，观测点距3.0～3.2 km。4 条剖面经过不同的构造带，其中，L5、L6、L7 剖面南段经过闽西南拗陷带，L5、L6 剖面北段经过闽西北隆起带，L7 剖面北段、L8 剖面经过闽东火山断拗带。4 条NE 向剖面的参数见表1，观测系统采用图2 所示的炮点与测点组成的全线追逐相遇观测，通过对获得的4 条NE 向剖面的15个炮点的P 波记录截面进行3～8 Hz 滤波和Pg 震相走时拾取，4 条剖面共拾取704 个Pg 波到时数据(表2)。

表1 4 条深地震测深剖面参数Table 1 Parameters of four deep seismic sounding profiles

表2 4 条深地震测深剖面Pg 波震相统计Table 2 Pg phase statistics of four NE-trending deep seismic sounding survey lines

图2 4 条北东向剖面的观测系统Fig.2 Observation systems of four NE-trending deep seismic sounding profiles

1.2 初至Pg 波波组特征

Pg 波作为结晶基底的折射波或地壳上部的回折波，在较近的炮检距以内是最先到达接收点的初至震相，到时可靠且信噪比较高(图3)。Pg 波走时曲线的变化反映了上部地壳的速度结构分布特征，表明剖面所处的基底结构以及地质构造的差异。一般来说构造隆起(褶皱)区走时曲线超前，表现为高速结构；断陷盆地或构造拗陷内走时曲线滞后，为相对低速区带。本文从4 条剖面中选取穿过不同构造带的4 炮地震记录截面，对Pg 波震相从运动学和动力学特征进行分析(图3)。

图3 部分炮点记录截面Fig.3 Sections of seismic records of some shot points

图3a 为武平–宁化–邵武剖面(L5)SP54 炮的地震折合记录截面图，图3b 为上杭–永安–建瓯剖面(L6)SP63 炮的地震折合记录截面图。SP54 炮整段和SP63炮的北段记录截面Pg 波组反映了南平以北的建瓯、泰宁、邵武一带的闽西北隆起带上部地壳结构。SP54炮和SP63 炮的近炮点20 km 以内，折合到时几乎接近零线，视速度接近5 980 m/s，表明建瓯、邵武附近下方速度存在高值，基底出露。SP54 炮左支记录截面Pg 波组可追踪118 km，其Pg 波震相折合走时较小，其最大折合走时也仅为0.34 s，表明盖层较薄、基底隆起、地表介质速度较高。SP54 炮左支记录截面距炮点100 km 以后，折合走时减小，振幅峰值衰减，显示处于构造带分界附近的Pg 波能量存在衰减作用。

图3c 为永定–大田–古田剖面(L7)SP73 炮的地震折合记录截面图。SP63 炮和SP73 炮的南段记录截面Pg 波组反映了沙县–永安、大田–漳平–永定一带的闽西南拗陷带上部地壳结构。SP73 炮的Pg 波到时最远可追踪到130 km，在距炮点68 km 后折合走时随距离增加而减小，表明龙岩地区的连城−新泉盆地沉积较厚，基底速度低。SP63 炮和SP73 炮地震记录截面经过闽西南拗陷带的折合走时与图3a、图3b 中的闽西北隆起带相比明显滞后，其折合走时为 0.37～0.67 s。SP73 炮左支记录截面Pg 波组的长追踪距离及稍大的折合走时反映了上地壳速度整体偏低，走时曲线幅度较大的锯齿状起伏可能与构造带内的高低速相间结构有关。

图3d 为漳浦–南安–罗源剖面(L8)SP84 炮的地震折合记录截面图。SP73 炮的北段和SP84 炮的整段记录截面Pg 波组反映了尤溪−古田、闽东南沿海的闽东火山断拗带上部地壳结构。Pg 波的最大折合走时为0.5 s 左右，与图3a、图3b 中的闽西北隆起带相比超前，与图3c 中的闽西南拗陷带相比滞后。SP84 炮的左右支的走时曲线变化不同，左支随炮间距增加走时相比右支延迟大，表明基底深度随炮点向南西方向逐渐加深。SP84 炮左支在炮间距60 km 左右走时曲线出现负斜率段，可能与断裂构造有关。SP84 炮的记录截面Pg 波尾波发育，在不同炮间距范围内振幅幅值呈现明显的衰减特性，反映了闽东沿海断块的上部地壳结构的横向非均匀特征。

2 初至Pg 波走时有限差分成像

2.1 初始模型的建立

福建地区4 条NE 向深地震测深剖面经过不同的构造带，穿过了多组NW 向断裂带，所处的地质构造环境和深部构造背景各不相同，因此，建立精准初始速度模型十分重要。结合福建地区已有的地壳一维速度模型[8-11,19]，本次采用分层细化法对初始模型与每一条测线的初至Pg 波走时进行拟合计算，得到与拾取的初至走时拟合较好的一维初始模型(图4)。4 条NE 向剖面的初始模型由不同分层的速度梯度组成，各剖面速度分层存在明显的速度变化，不同的构造单元和断裂剖面的速度模型存在差异性。

图4 4 条深地震测深剖面的初始一维P 波速度模型Fig.4 Initial 1-D P-wave velocity models of four NE-trending deep seismic sounding profiles

2.2 正反演计算

本文使用有限差分层析成像算法对4 条NE 向剖面的初至Pg 波走时数据进行计算，该方法适用于上部地壳速度结构横向变化剧烈的地下介质[20-25]。4 条NE 向剖面观测点距介于2.5～3.5 km，正演计算水平和垂直方向采用1.0 km×1.0 km 的网格单元划分模型空间，通过程函方程的有限差分算子计算初至波走时场，根据走时梯度确定射线路径。反演计算网格尺度水平和垂直向大小分别为4.0 km×1.0 km，使用反投影迭代计算所有射线在网格单元中引起的慢度扰动以逼近真实速度。

2.3 结果可靠性分析

1) 走时残差

有限差分层析成像反演通过迭代过程来实现计算精度，迭代一定次数后走时残差(ERMS)趋于稳定。4 条剖面经过50次迭代的走时残差的均方值如图5a1−图5d1所示，可以看出，30 次迭代收敛良好，均方根走时残差ERMS从初始最大的139 ms 下降到22～31 ms。图5a2−图5d2 是4 条剖面30 次迭代后的拾取走时与理论走时的拟合分布图，图5a3−图5d3 是4 条剖面30 次迭代后的走时残差分布直方图，可以看出，反演后理论走时与观测走时拟合较好，4 条剖面走时残差近似符合正态分布。

图5 4 条剖面反演后的均方根值、走时拟合和走时残差直方图Fig.5 Post-inversion RMS values,travel time fitting,and histograms of travel time residuals of four NE-trending deep seismic sounding profiles

2) 射线分布

图6 给出了4 条NE 向剖面30 次迭代后的射线分布图，沿剖面的绝大部分测段横向上都有射线经过，且射线数至少5 条以上。4 条剖面除南、北两端附近测段外，射线分布较为集中的结晶基底以上的浅部区域密集交叉，射线密集覆盖基底以下的上部地壳部分位置也有射线穿过(图6)。4 条剖面所处的地质情况和构造环境不同，由西到东剖面的射线覆盖深度增加，表明福建地区由内陆到沿海基底以下的上部地壳速度呈增加的趋势。L5−L6 剖面由南到北射线穿透的回折面变浅(图6a−图6b)，与剖面经过的闽西北隆起带与闽西南拗陷带息息相关，表明闽西北基底埋深较浅且速度较高。L8 剖面经过的闽东火山断拗带的射线覆盖的回折面存在明显的横向差异变化(图6d)，L6−L7剖面经过的闽西南拗陷带射线整体穿透深度较深(图6b−图6c)。

图6 4 条NE 向剖面Pg 震相反演的射线覆盖Fig.6 Ray coverage for the inversion of Pg phase in four NE-trending deep seismic sounding profiles

3)交叉验证

福建地区深地震测深的8 条测线呈网格交叉布设，4 条NW 向测线(L1−L4)与4 条NE 向测线(L5−L8)交叉点的爆破点资料可以共用，为获得可靠的剖面基底速度结构奠定了基础。以漳浦–南安–罗源剖面(L8)的4 个炮点为例，SP81 炮点与L4 剖面、SP82 炮点与L3 剖面、SP83 炮点与L2 剖面、SP84 炮点与L1 剖面相互交叉，L8 剖面的4 个炮点与相交剖面交叉点附近的一维速度结构拟合分布如图7 所示，交叉点的速度结构拟合得较好，随深度增加速度拟合误差在±300 m/s 以内。

图7 L8 剖面4 个炮点和L1−L4 剖面交叉点一维速度结构Fig.7 1-D velocity structures at locations where four shot points on L8 profile intersecting the L1-L4 profiles

3 上部地壳结构特征

根据4 条NE 向剖面射线覆盖有效范围(图6)，本文给出了射线经过的有效范围内上部地壳精细速度结构和解释结果(图8)。由图8 可知，4 条剖面速度的横向和纵向变化明显，速度值介于4 700～6 100 m/s，速度等值线随深度的增加由密到疏、形态由平坦到复杂。4 条剖面能够从速度结构上对福建地区的盖层和基底、盆地、平原、岩浆带以及不同构造格架区分开来。参考华南地区浅变质花岗岩系地震P 波速度为5 900～6 000 m/s[26]，通过4 条剖面速度等值线疏密变化、射线覆盖密集情况，本文以速度值为5 900 m/s 的等值线所对应的深度为结晶基底，研究区结晶基底界面埋深介于1.0～6.0 km，4 条剖面基底界面由北往南有所加深，基底面的起伏与构造凸起和凹陷对应。

图8 4 条北东向剖面上地壳速度结构及构造特征分布Fig.8 Velocity structure and structural features of the upper crust along four NE-trending deep seismic sounding profiles

4 条剖面从西到东上部地壳速度总体增加，位于(闽东)东南沿海中生代岩浆带的L7 和L8 剖面的上部地壳速度整体大于(闽西)华夏地块的L5 和L6 剖面。其中L5 剖面(武平–宁化–邵武)位于闽赣交界处，沿剖面结晶基底以上的盖层速度等值线较平缓，基底以下上部地壳总体速度呈高速−低速−高速的变化特征(图8a)，剖面坐标20～65 km、90～275 km 从沉积盖层到结晶基底速度整体偏低。剖面基底埋深介于1.0～6.0 km，其中武平凹陷和泰宁盆地低速层较厚且基底埋深接近5～6 km。剖面坐标65～90 km 的武平以西的桃溪附近地下1 km 处有一高速异常上涌带，结晶基底埋深较浅为1 km，与桃溪隆起中的变质核杂岩金属矿有关[27]。剖面坐标200～230 km 的宁化结晶基底3 km以下的速度呈NW向变化，剖面坐标285～325 km 的邵武基底速度等值线较平缓且基底1 km 以下的上地壳速度明显增加。

L6 剖面(上杭–永安–建瓯)位于(闽西)华夏地块中部，剖面速度从表层的4 700 m/s 增加到6 km 左右的6 000 m/s，剖面沉积盖层速度等值线呈宽缓相间或狭窄弯曲出现，剖面经过的上杭、连城–新泉、沙县等盆地速度整体偏低，剖面坐标140～200 km 速度呈高−低−高明显变化，剖面坐标270～305 km 处速度变化整体偏高。剖面基底埋深介于1.0～5.0 km，剖面中段的永安–南平地区基底面形态起伏较大，剖面坐标150、185、245 km 处的基底埋深较浅为1 km 且基底以下的速度明显增大呈NW 向狭长条带变化(图8b)，与闽西南的推覆构造和宁化–南平岩浆带相关，可以为该区的控矿构造提供深部地球物理学信息[28-29]。

L7 剖面(永定–大田–古田)位于福建中部的政和–大埔断裂带附近，剖面坐标150 km 之前的南段速度变化明显，从表层的4 700 m/s 增加到8 km 左右的6 100 m/s，剖面中南段的永定–尤溪地区盖层速度等值线较平缓，剖面坐标260 km 以后的北段上部地壳速度变化明显，剖面坐标75、175、275 km 附近的局部地区呈高速异常上涌。剖面基底埋深介于1.5～6.0 km，结晶基底面起伏凹凸变化明显，剖面坐标125 km 附近基底界面有明显的下凹特征且深约6.0 km，剖面坐标240 km 后基底界面开始逐渐抬升(图8c)。基底速度结构变化跟所处的政和−大埔分界断裂带内的地幔玄武岩浆上升，喷出地表或者底侵到上地壳浅部有关[30]。

L8 剖面(漳浦–南安–罗源)位于福建东部沿海中生代岩浆带内，右侧是平潭–东山变质带，剖面沉积盖层速度等值线较平缓，剖面坐标130 km 以后的中部和北段速度变化明显，地下8 km 附近存在速度6 100 m/s的高速带。相比L5−L7 剖面的基底速度结果，剖面上部地壳速度整体偏大且纵向变化显著，剖面经过的龙海平原和福州盆地表现为低速下凹。剖面基底埋深介于1.0～3.0 km(图8d)，基底以下速度变化明显呈NW 方向，与岩浆带形成的区域变质和混合岩化作用有关[31]。

4 条剖面的基底精细速度结构与研究区已有的P 波成像结果[6-14,9-11,32-35]相比，提高了福建地区结晶基底与沉积盖层的浅部速度结构的成像精度，给出了NW向断裂带的浅部构造特征，4 条NE 向测线与4 条NW向测线交叉分布，成像结果互为补充，丰富了福建地区上部地壳的精细速度结构。本文结果针对研究区的“东西分带、南北分块”的构造格架从速度结构、基底埋深、构造特征等进行分析，获得NW 向断层的基底构造特征，结合地质构造、地层岩性、地震、地热等资料，为了解该区的结晶基底构造属性以及古太平洋板块向欧亚板块的俯冲作用的东南缘构造特征提供深部地球物理学证据。

4 结果与讨论

L5−L8 4 条剖面自西向东经过(闽西)华夏地块和(闽东)东南沿海中生代岩浆带，从基底速度结构上来看，地表速度整体偏高，结晶基底埋深整体较浅，基底及以下的上部地壳速度偏大，这与4 条剖面所在的东南缘第四纪沉积地层较薄，白垩纪火山岩和花岗岩分布广泛、岩浆活动地表热流值较高有关[36]。

以政和–大埔断裂为界，4 条剖面从西到东经过了闽西北隆起带、闽西南拗陷带及闽东火山断拗带，4 条剖面基底速度结构结果显示，在不同的地质构造单元之间的接触带两侧速度有明显的变化。NE 走向的L5−L7 剖面南段位于闽西南拗陷带，基底速度结构整体表现为低速，且速度横向较为均匀；而L5−L6 剖面北段古老的闽西北隆起带内，上部地壳结构速度较高，且表现为较强的横向非均匀，显示了福建闽西南的凹陷沉积建造与闽西北的抬升构造特征。L8 剖面位于闽东火山断拗带内，发育一系列断陷盆地，基底速度结构呈现明显的高低速差异变化。从4 条剖面的基底埋深来看，位于闽西北隆起带内的L5 和L6 剖面北段基底形态稳定且埋深较浅，表明基底结构形成较早，与该地区的地表出露的前寒武纪变质岩有关，基底以下NE 向展布的高速异常可能与基性玄武岩组成的火山岩带相对应，与晚侏罗世晚期至早白垩纪早期发生的大规模陆相火山岩喷发的地质特征有关[37]。闽西南拗陷带内的L5、L6 和L7 剖面南段所在的盖层埋深变化平缓，推测与闽西南分布的震旦纪–奥陶纪地层变质程度较低有关[38]，基底面起伏明显，最深达6 km 左右且基底深度呈深浅间隔变化，与闽西南拗陷带内的凹陷与隆起多个次级构造单元间隔分布相关[39]。永安、大田附近基底以下存在的高速异常与该地区存在的推覆构造有关，且速度异常的变化方向由NW 向到SE 向，证实上部地壳内确实存在受古太平洋板块的俯冲引起的大型推覆构造[40]。闽东火山断拗带内的L8 剖面基底以下速度高、变化大，与该地区的白垩纪花岗岩和火山岩出露相符合，与该区中生代以来受古太平洋板块俯冲的强烈影响，上地幔的超基性岩受到造山挤压作用，形成的铁镁质岩浆沿断裂构造底侵，深部物质上涌、古老的基底重熔的构造特征有关。

L5−L8 4 条剖面自西向东基底速度逐渐增强，以政和–大埔断裂为界，L7 和L8 剖面所在的(闽东)东南沿海中生代岩浆带基底以下8 km 处存在6 100 m/s 的高速带(图8)，表明华南大陆东部所受的构造活动强于西部，受古太平洋板块向欧亚板块的俯冲所致[16-18]。本次4 条NE 向剖面与研究区已有的4 条NW 向剖面的基底速度结构变化特征来看，本文4 条剖面速度横向变化及界面起伏明显较4 条NW 向剖面的变化小，且4 条NE 向剖面基底速度变化整体呈NW 方向展布，成像结果与研究区东西分带的构造特征和华南地区构造变形的主方向呈NW-SE 向一致[41-42]，本文为该区中生代以来古太平洋板块向欧亚大陆俯冲的方向研究提供深部地球物理学证据。

L5−L8 4 条剖面横穿了研究区的NW 向断裂带,根据剖面断裂两侧的基底速度横向变化和基底埋深起伏情况，结合研究区的地震地质资料，确定的主要基底断裂带有闽江断裂、沙县−南日岛断裂、永安–晋江断裂、九龙江断裂及上杭−云霄断裂，从基底速度结构解释图上来看(图8)，NW 向断裂带所在的位置表现为倾斜舌状低速异常或者过渡带的特征且基底埋深的起伏较大，本次4 条剖面上部地壳结构对NW 向断裂的分布起到约束作用，并且提供了很好的地震学证据。其中闽江断裂控制福州断陷盆地，断层两侧为高速异常、断层近直立的结果与张路等[43]认为的高倾角正断层倾向盆地中心比较符合。沙县–南日岛断裂作为福建地质构造南北分块的分界带，断裂两侧速度结构有明显差异；断裂以北4 条剖面速度偏大，纵向变化明显，基底埋深较浅，基底以下的高速异常与变质岩呈弧状向东突出比较一致。断裂以南基底埋深较大且起伏变化明显，基底以下速度异常与印支期花岗岩呈向西突出相关，可能是受古太平洋板块的俯冲推挤，南北部向西推移的速度差异所致。永安–晋江断裂在剖面所在的位置表现为左右高低速异常明显，断层较陡，永定到大田一带断裂左右两边速度变化剧烈，与该区地震孕育环境密切相关。九江断裂在L7−L8 剖面南段表现为高速特征，与用气枪地震资料反演得到的漳州–高雄断裂上地壳速度结构变化特征比较一致[34]，推断漳州−高雄断裂即是九龙江断裂向东延入台湾海峡的断裂。从L8 剖面速度结果来看，断陷盆地受NW 向和NE向两组断裂控制，第四纪以来NW 向断裂活动强烈且以引张活动为主，因此，在速度结构上经过的新生盆地呈现NW 向展布，向南西掀斜的特征，表明NW 向断裂在东南缘现今板块动力学体系中扮演着十分重要角色。

5 结论

a.福建地区4 条NE 向剖面地表速度整体偏高，结晶基底埋深整体较浅，基底及以下的上部地壳速度偏大，在不同的地质构造单元两侧，上部地壳速度结构横向变化明显。

b.成像结果反映了研究区东西分带的构造特征，表明华南大陆东部所受的构造活动强于西部，与该区的地质构造和构造变形的主方向一致。

c.根据4 条剖面速度横向变化和基底埋深起伏情况，给出的NW 向断裂浅部构造特征研究成果与前人结果具有一定对应关系，表明NW 向断裂构造在东南缘现今板块动力学体系中扮演着十分重要角色。

d.今后的研究需要联合Pg、Sg 数据开展反演工作，或直接进行三维速度结构反演，相互约束提高模型的精度。

致谢：中国地震局地球物理勘探中心诸多同事参与了野外探测工作，感谢参与此项工作的全体人员对文章的帮助！感谢Hole J.A.提供的程序。