川滇地区位于青藏高原东南缘，在印度-亚欧大陆碰撞构造运动的影响下，形成了剧烈的地形隆升，是理解青藏高原隆起过程中物质运移的关键区域。使用地震学方法，从地震波速的角度对该区域的地壳结构进行描述，是确定青藏高原物质运移过程的重要手段。剧烈的构造运动使得青藏高原的边缘区域具有复杂的地壳结构，需要高质量的地震资料和高精度的成像方法才能细致地确定相应的波速结构。川滇地区具有很高的地震活动性，提供了大量可供使用的天然地震资料。同时，双差地震层析成像方法^[1-2]能够使用地震走时资料，精细反演震源密集区域的速度结构和震源的相对位置。运用双差地震层析成像方法，基于该地区大量的地震资料可以精细地反演该区域的地壳波速结构。另一方面，使用该方法可以得到高精度的地震震中分布，能够更好地刻画断层的几何形态，有助于进一步理解构造运动的状态。

近年来，深度学习方法在地震检测和深度学习领域大量应用，使得自动拾取的精度达到甚至超过人工拾取^[3-6]。同时，在人工拾取中，对于同一段波形记录，不同的拾取者因为拾取经验和主观因素的影响拾取的结果存在差异，而深度学习方法的捡拾结果则有较好的一致性。APP++^[7]是在APP^[6]基础上发展起来的一种台阵地震检测和震相拾取方法，具有很好的P波、S波拾取精度和跨区域的泛化能力。本文利用APP++方法扫描川滇地区2014—2019年间的连续波形所得到的大量的地震目录和到时信息，使用区域双差地震层析成像方法^[1-2]（Regional Double-Difference Seismic Tomography, tomoFDD）反演了该区域的三维P波、S波速度结构和震源位置。

采用的数据来自川滇地区的固定台网（116个）和小江断裂带附件的高密度流动台阵（51个）2014—2019年间6年的连续波形数据（图1）。利用APP++方法扫描该数据，共检测到73 291个地震、537 554个P波记录、471 459个S波记录。为了保证走时数据的质量、适应区域双差层析成像的计算代价并尽量均衡地覆盖整个研究区域，对拾取到的地震进行筛选。首先，选取5个以上台站记录到，与利用taup^[8]方法iasp91模型^[9]计算的理论走时残差在1.5 s以内的地震；其次，将研究区域按经、纬度等分为$ 90\times 90 $的小区域，每个区域内挑选不超过30个被更多台站记录到的地震；最后，挑选出21 160个地震，包含P波到时183 064个、S波到时157 620个（图2）。

图  1  台站分布图

图  2  用于双差层析成像的地震分布图

统计台站的震源距和其波形记录上拾取到的P波和S波走时，两者的频次分布如图3。可以看到，大部分拾取到的P波和S波走时结果与震源距呈现很好的线性关系，弥散也很小，反映出APP++方法拾取的震相走时是准确的。

图  3  震源距-走时关系的频次分布

采用互相关的方法提取相邻（相距0.3°以内，被至少3个相同的台站记录到）2个地震到同一台站的同一震相（P或S）的走时差。将涉及的2个波形记录分别以APP++拾取的到时为原点，截取为−2～3 s和−3～4 s的两段三分量波形，同一分量间做互相关，取最大值作为互相关结果，互相关值最大处的时刻为2个记录的走时差。根据APP++方法拾取的绝对到时数据，也可以得到2个记录间的走时差。通过统计互相关走时差和APP++走时差之间的差别与互相关系数的关系（图4）发现，当互相关系数增加时，两种到时差的差别变小。这可以反映APP++方法的稳定性：当波形相似时，拾取结果接近。最终选取互相关值大于0.6的互相关走时差作为区域双差层析成像的393 558条走时差数据。

图  4  拾取震相到时差-互相关值关系频度分布图

使用Zhang和Thurber提出的区域双差地震层析成像方法^[1-2]，对该区域的速度结构和地震位置进行反演。区域双差地震层析成像方法可以同时使用绝对走时数据和相对走时数据对震源位置和速度结构进行反演，相对走时数据的加入，使得方法能够更加精细地刻画震源密集区内震源的相对位置。为了更好地对区域尺度和存在不连续面的地下介质进行建模和走时正演，该方法将球形的地球结构嵌入到直角坐标系中，并使用有限差分方法进行走时计算和射线追踪。

反演中，将图5展示的一维结构作为区域双差层析成像的初始模型（这里的深度以海平面为原点）。模型垂直方向上设置的节点有−150 km、−2 km、0 km、2.5 km、5 km、7.5 km、10 km、12.5 km、15 km、20 km、30 km、40 km、50 km、60 km、80 km、500 km；水平方向上节点的设置如图6的棋盘格测试所示。反演用到的地震资料为前文提到的APP++拾取的地震目录和走时数据，以及从波形数据中使用互相关方法提取的走时差信息作为该方法的输入。

图  5  一维初始模型

图  6  棋盘格测试模型

区域双差地震层析成像方法的程序中包含了棋盘格测试的部分，具体步骤为：①在给定的参考模型的节点上，周期性地给予+5%或−5%的扰动作为真实模型；②在真实模型的基础上，进行走时正演，生成和所使用地震资料分布相同的模拟数据；③以参考模型作为初始模型，模拟数据作为地震走时资料，用区域双差地震层析成像方法反演区域的速度结构；④反演得到的速度结构与真实模型的异同，可以反映在相同参数设置和约束分布的情况下，反演方法理论上能够有较好分辨率的区域。

图6是用作棋盘格测试的真实模型，其中P波和S波的扰动幅度一致；图7～8展示的是复原得到的P波和S波速度结构。测试结果显示，在研究区域中部存在一个EN向展布的分辨率较好的区域。

图  7  棋盘格测试的P波速度复原结果

图  8  棋盘格测试的S波速度复原结果

使用区域双差地震层析成像反演后，共有17 240个地震被重定位，其震中分布如图9。可以看出，不同地震丛集内的地震更加集中，与断层分布吻合更好，能够更加细致地刻画区域内的断层展布。图10是区域双差地震层析成像反演给出的地震位置在3个方向上的不确定度估计：多数地震的定位误差都在25 m以内，垂直方向的误差比水平方向的误差更大。这说明本文所给的走时数据一致性较好，相对于该区域内基于人工走时数据的双差层析成像^[10]残差更小。

图  9  双差层析成像方法重定位后的地震空间分布

图  10  双差层析成像地震位置误差分布频次图

图11和图12是分层的P波和S波速度反演结果，图13则展示了图2中6条线段所对应垂直剖面的速度结构，结果中同时展示了相应位置附近的地震震中分布。图13中A-A' 剖面中距离A点100 km和200 km的地方，地震呈现很明显倾斜取向，反映出对应位置的断层形态，与龙门山断裂带的断层分布及性质很好对应。100 km的地方，是2个近直立的断层；200 km处是倾向NW的倾斜断层且倾角随深度有变化，反映了龙门山主断裂处由青藏高原向四川盆地仰冲的构造特点。在200 km附近还能看到不太明显的其他取向，可能对应于其他地震活动性稍弱的断层^[11]。图13中的B-B' 剖面中，190 km附近也能观测到龙门山主断裂一个仰冲的断层形态，但是受其他地震影响，取向不够清晰；370 km附近，存在大量的地震丛集，为长山镇断裂上发生的地震。图13中C-C' 剖面是沿龙门山主断裂的一个垂直剖面，地震呈现面状分布的特点。这3条剖面反映的地震分布特点与前人^[11–13]的研究结果一致。D-D’剖面呈现类似现象，100 km处呈现倾斜分布地震条带，200 km处呈现近垂直地震分布条带，分别对应于倾滑断层和走滑断层。

图  11  P波速度反演结果

图  12  S波速度反演结果

图  13  剖面上的P波速度（左侧）和S波速度（右侧）结构及剖面上的断层分布。

E-E' 和F-F' 剖面分别沿纬线和经线方向划分小江断裂带附近的区域。剖面内地震主要集中于25 km以内，与李乐等^[14]、邓山泉等^[10]在该区域的研究结果相似。小江断裂带区域内，断层多、地震密，在E-E' 和F-F' 剖面上体现出一定的丛集，所反映的空间位置取向较为复杂。剖面上，除了地震丛集的地方外，还存在部分弥散分布的地震，说明该地区可能存在许多发育地震的断层。今后研究这些弥散地震的分布，可能有助于理解当地的地质构造及其与地震活动性的关系。

浅层的速度模型表现出与断层很高的一致性。在深度20 km以内的水平结果中，小金河断裂、小江断裂、楚雄断裂、红河断裂、金沙江断裂附近都存在较明显的低速异常，而这5个断层所围绕的区域，即川滇菱形块体的东南部分则速度较高，正好对应于峨眉山大火成岩省。D-D' 剖面经过峨眉山大火成岩省，其相对速度剖面（图13g1～g2）更加显著地刻画了大火成岩对速度结构的影响：在距离D点200 km左右，50 km深度以内呈现明显的高速异常，50 km开始变为低速异常。这部分与Yang et al.^[15]、张智奇等^[16]和郭希等^[17]利用面波频散得到的剪切波速度结构一致。在距离D点100 km附近（图13g1～g2）存在一个垂直的高、低速块体分界面，并分布有与之取向一致的地震丛集，这里可能是大火成岩省的西北边界。

在20 km深度内，龙门山断裂带以东的四川盆地呈现出低速异常，西部的青藏高原为高速异常；当深度大于20 km时，速度结构逐渐向东部高速和西部低速转变。Yang et al.^[15]等人是在20 km左右观测到剪切波高低速异常的转变，雷建设等^[12]则是在45 km观察到P波速度的明显改变，本研究观测到转变的深度则处于两者中间。反演使用数据的差异、拾取精度的差别、建模粒度的不同都有可能影响结果的细节，造成不同研究之间模型的差异。但是总地来说，大尺度的浅部和深部区域的速度异常变化在不同研究中是一致的。由于P波的走时数据较多，故在50 km、60 km深度能够得到较为准确的P波速度结构，而S波的可分辨区域则相对较小。图13中A-A' 的P波速度剖面直观地表现了龙门山断裂带随深度的速度结构变化。青藏高原的地壳比四川盆地的更厚，在50 km左右的深度上，四川盆地已经进入了波速更快的地幔，而青藏高原则还处在波速较之更慢的地壳中，由此导致了这里观察到的速度变化模式。

小江断裂带区域内的E-E' 和F-F' 剖面在近地表附近的P波、S波速度结构（图13e1～f2）观察到了明显的低速，这是由于小江断裂带浅层破裂较多，整体波速偏低。这也与地震弥散分布可能反映的小断层多一致。在20 km深度附近的P波速度结构中观测到明显的低速，存在于小江断裂带的西北部。张智奇等^[16]认为这里低速异常可能来源于青藏高原物质向东南缘的挤出。

由于双差地震层析成像结果严重依赖于地震的空间分布，当某一区域在特定深度上地震分布较少时，这部分速度结构往往不能被很好约束，从而影响整体的速度结构变化模式。在对速度反演结果解释时需要考虑地震分布对分辨率的影响。在E-E' 剖面的速度结构中，观察到东部的低速区域比西部分布更深，但东部区域地震太少，对这部分区域的分辨率有限，因此不予讨论；深度为40 km时，龙门山断裂带西北面的青藏高原部分地区地震分布比较稀疏、分辨率下降，导致部分区域没有呈现出地壳变厚引起的低速异常（图11d），而且穿过这部分区域的B-B' 剖面也相应地没有明显反映出地壳厚度的变化（图13b1～b2）。

运用APP++扫描川滇地区6年的连续数据得到的地震目录和走时数据，使用区域双差地震成像方法获得了川滇地区高精度的地震空间分布和地壳的三维P波、S波速度结构。地震的位置分布和区域内的断层具有很好的一致性，很好地刻画了主要断层的形态。研究所得速度结构与相同区域内的相关研究也基本一致^[18]，部分地反映了研究区域的地壳厚度变化，同时也反映了峨眉山大火成岩省对区域速度结构的影响。

当然，研究结果也显示速度结构的分辨率和台站及地震分布有显著的相关性，再次显示双差层析成像方法需要密集的台站及长期的观测。为了进一步提高分辨率，必须结合尽可能多的高质量地球物理观测资料，进行联合反演。

致谢　中国地震局地震预测研究所为本研究提供地震波形数据；本研究工作得到北京大学高性能计算校级公共平台支持。