印度板块和欧亚板块于50 Ma前开始碰撞，持续的汇聚导致了青藏高原的隆升和增厚的中下地壳(Molnar and Tapponnier, 1975; Houseman and England, 1993; Royden et al., 2008).龙门山所在的青藏高原东缘具有高原周缘最大的地形梯度，其在几十公里的水平距离内，相对四川盆地抬升了近4 km.GPS观测(Chen et al., 2000; Meade, 2007)和SKS各向异性测量(Wang et al., 2008)显示，地壳尺度(或者是岩石圈和软流圈尺度)物质向东运移，但地壳变形速率很小( < 3 mm·a^-1).2008年5月12日的汶川7.9级地震却恰恰发生在这种似乎稳定的构造背景下(Wang and Meng, 2009)，引发了大量关于龙门山抬升机制的研究和争论.

目前，关于龙门山抬升机制存在许多争议，主要概括为以下两种模型(Hubbard and Shaw, 2009)：脆性地壳增厚模型(Tapponnier et al., 2001)和下地壳流模型(Royden et al., 1997; Clark and Royden, 2000; Zhang et al., 2009; Bai et al., 2010; 王椿镛等, 2006).前者认为岩石圈尺度内伴随大规模走滑的逆冲断裂导致了龙门山的抬升，后者认为来自青藏高原中部的低黏度下地壳物质向外流动，受到坚硬四川盆地的阻挡，推升其上的地壳进而形成陡峭的龙门山.

关于扬子块体和青藏高原的接触关系也存在争议.大量研究认为，龙门山断裂带对应着地形、重力、地壳厚度和速度结构等的陡变带，因此将其定义为高原边界(Clark and Royden, 2000；Cook and Royden, 2008). 崔作舟等(1996)认为四川盆地的坚硬地壳楔形插入松潘—甘孜块体下地壳中；王绪本等(2009)通过大地电磁测深获得的电性结构模型推断可能存在一个西延至若尔盖地块的泛扬子陆块；Guo等(2013)通过深反射地震剖面，指出龙日坝断裂是扬子块体和松甘块体地壳物质的边界.

目前各种观点莫衷一是，存在争议的根源之一在于对藏东缘龙门山断裂带域缺乏精细可靠的地壳结构约束，各种构造和物性因素缺乏统一分析.Jia等(2014)利用近垂直穿越龙门山中段的宽角反射/折射地震数据，获得了该区二维地壳速度结构，揭示了不同构造单元间的地壳结构差异及介质由上到下的岩性变化.本文拟参考该区已有的深地震反射剖面结果(Guo et al., 2013)，利用新发展的起伏地形下高精度走时反演方法，对该剖面数据重新处理解释，以期得到更准确可靠的地壳速度结构，深化认识该区的构造演化动力学过程和强震机制.

龙门山断裂带位于青藏高原东缘，处于青藏高原和扬子块体对接的关键部位，西倚松潘—甘孜地体，东连四川盆地(图 1).在该区域，自西北向东南依次分布有龙日坝断裂、龙门山断裂和龙泉山断裂等断裂构造.

图 1

Fig. 1

图 1 “阿坝—龙门山—遂宁”人工源深地震测深剖面位置(菱形代表炮点，黑点代表接收器)、研究区地震分布 (圆形，另五角星为汶川地震位置)(数据来自国家地震数据中心)及已有深地震测深剖面位置(黑线) ① 龙门山三角剖面(陈学波等，1986)；② 花石峡—简阳剖面(崔作舟等, 1996)；③ 花石峡—简阳剖面南段(Li et al., 2009)；④ 黑水—邵阳北段(Li et al., 2009)；⑤ 唐克—奔子栏剖面(Wang et al., 2007)；⑥ 巴塘—资中剖面(Wang et al., 2007)；⑦ 深地震反射剖面(Guo et al., 2013). Fig. 1 "Aba-Longmenshan-Suning" profile (diamonds are shots and black points are receivers), earthquake distribution of the study area (circles, and the star represents the Wenchuan earthquake) (source: China Earthquake Data Center) and the existing profiles in this area (black lines) ① Longmenshan triangle profile (Chen et al., 1986); ② Huashixia-Shaoyang profile (Cui et al., 1996); ③ Southeastern part of Huashixia-Shaoyang profile (Li et al., 2009); ④ Western part of Heishui-Shaoyang profile; ⑤ Tangke-Benzilan profile (Wang et al., 2007); ⑥ Batang-Zizhong profile (Wang et al., 2007); ⑦ Deep seismic reflection profile (Guo et al., 2003).

青藏高原东缘在中生代和晚新生代经历了强烈的构造变形，急剧抬升.由于青藏高原的隆升及其地壳物质向东运移过程中受到坚硬的四川盆地阻挡，导致松潘—甘孜块体及其上覆复理石沉积层被挤压，形成了一系列叠瓦状构造的逆冲推覆构造带(Chen and Wilson, 1996).扬子块体自晚古生代以来沉积环境比较稳定，具有较厚的未变质沉积盖层.在始新世和渐新世期间，才出现了不同程度的褶皱运动(任纪舜等，1999).四川盆地出现龙门山推覆构造带的最新前陆盆地，盆地中部的龙泉山构造带则是这一前陆盆地的东缘边界构造带(邓起东等，1994).龙门山断裂带地表地质特征主要是广泛分布的前寒武变质岩带，东北向西南依次分布有彭灌杂岩、宝兴杂岩及康定杂岩三个主要杂岩带(Xu et al., 2008).

深部地球物理探测是获得研究区域深部构造环境和孕震环境的根本手段.龙门山断裂带处于青藏高原东缘，是研究青藏高原隆升和扩展动力学过程的重要窗口.尤其汶川地震后，国内外众多学者在该区开展了大量的地震勘探工作，其中Liu等(2014)利用接收函数和面波联合反演的方法获得了该区S波速度结构，给出了松甘块体和四川盆地的碰撞变形关系和速度结构差异.Zhang等(2009)利用接收函数法获得了跨越龙门山的地壳厚度变化和壳内间断面的横向变化特征.Bai等(2011)、Lei和Zhao(2009)和Xu等(2010b)利用天然地震层析成像方法获得了该区P波和S波速度结构.

图 1黑线为在该区所开展的人工源地震探测的测线分布.各观测结果均认为，以龙门山断裂带为界，高原和盆地的地壳厚度、速度结构均存在极大差异.主要表现为高原平均速度低、地壳较厚，盆地平均速度高、地壳较薄.松潘甘孜块体内上地壳是否存在低速层，以及低速层的深度和厚度存在争议(崔作舟等, 1996; 王有学等，2005；Li et al., 2009; Jia et al., 2014).龙门山断裂带下方Moho面被抬升，但Moho面的陡变形态也存在争议.松潘—甘孜块体地壳在南北向也存在差异，以鲜水河断裂为界，以南地壳厚度超过60 km，以北为大约53 km(崔作舟等, 1996; 王椿镛, 2003; Wang et al, 2007).

2010年中国地震局地球物理勘探中心开展了“阿坝—龙门山—遂宁”人工源宽角反射/折射地震探测试验.剖面长度为500 km，沿剖面实施了12次爆破，单炮药量为2吨左右，炮间距7.2~129 km不等.沿剖面分布有450台DAS-1(2)、PDS-1(2) 型便携式三分量数字地震仪器固定观测，观测点距一般0.4~2.5 km.在高分辨探测区段(龙门山断裂带附近)的炮点和观测点间距较小，最小分别可达7.2 km和0.4 km.数据处理过程中，我们选择记录良好的10炮地震数据，沿测线自西向东分别标记为Sp1, Sp2, …, Sp10(图 1).

震相拾取过程中，参照同区域的已有深地震反射剖面(Guo et al., 2013)结果，拾取能量较强震相而忽略部分有歧义震相.所拾取震相符合走时互换和震相分析等测试.

与Jia等(2014)相同，在各炮地震记录图上，Psed或Pg是初至震相，均可清晰辨认.测线西部松潘甘孜块体内，Psed震相仅在0~20 km内出现，视速度和折合到时分别为4.0~5.3 km·s^-1和0~0.8 s.Pg波可在20~150 km范围内接收到，视速度和折合走时分别为5.8~6.3 km·s^-1和0.8~1.0 s.远炮检距清晰接收到的连续Pg震相说明在近地表(深度大约10 km范围内)无低速层存在(Wang et al., 2007).在四川盆地内，Psed出现距离、视速度和折合到时分别为0~50 km、3.5~5.2 km·s^-1和0~1.8 s.Pg波紧随Psed波出现，追踪区间、视速度和折合走时分别是50~100 km，5.8~6.2 km·s^-1和1.5~1.8 s.龙门山断裂带上方炮点附近几乎无Psed震相，Pg震相视速度较高、折合走时较低，分别为5.8~6.3 km·s^-1和0~0.5 s.

在松潘—甘孜块体内有4个壳内反射震相(分别记为P2, P3, P4, P5)，在四川盆地只有3个(无P4)，其中反射震相P4与深地震反射剖面(Guo et al., 2013)中松潘甘孜块体内大约8 s的强反射震相一致.松潘甘孜块体内，中上地壳(P2, P3和P4) 反射能量较强，下地壳(P5) 和Moho面反射能量较弱.各壳内反射震相(P2, P3, P4和P5) 平均速度分别为5.91 km·s^-1、6.00 km·s^-1、6.05 km·s^-1和6.10 km·s^-1.Pm波出现区间和平均速度分别为150~250 km和6.26 km·s^-1.各震相平均速度依次增大，说明无壳内低速层存在.其中，P3和P4震相平均速度相差0.05 km·s^-1, 在地震记录上，Sp1、Sp2和Sp3炮上的右支在近炮点附近(50~100 km)近似平行，但远炮点范围内梯度增大.推测C3和C4界面之间速度变化较小，但并非显著的壳内低速层.四川盆地壳内反射震相(P2, P3和P5) 较为明显，追踪区间为10~150 km，平均速度分别为5.50 km·s^-1、5.80 km·s^-1和6.15 km·s^-1.盆地内Moho面反射震相(Pm)能量较强可准确追踪，出现区间和平均速度分别为100~200 km和6.39 km·s^-1(图 2).

图 2

Fig. 2

图 2 松潘甘孜块体和四川盆地平均速度和平均深度 Fig. 2 Average velocities and depths of each layer in the Songpan-Garzê block and Sichuan basin

地壳浅层初至波(Pg)震相清晰，可准确反演地表浅层速度结构.沿着测线存在较大的海拔高程变化(约4.3 km)，因此浅层结构反演中起伏地表的处理直接影响反演结果的精度及可靠性.本文上地壳结构反演时，采用新发展的起伏地形下高精度有限差分走时反演方法(Hole, 1992；Ma and Zhang, 2014; 张新彦等, 2017)，模型参数化采用与地表起伏变化相一致的贴体网格划分模型.该方法严格考虑了起伏地形对射线分布及反演结果的影响，因此反演结果具有更高的精度.正演过程采用0.5 km×0.5 km的网格单元剖分模型.反演过程中，重新抽样(网格化)因子为(8, 4)，滑动平均滤波器为(7, 3)，以对反演模型进行平滑.

Pg波走时反演基础之上，参考前人成果和一维速度分析结果，构建二维初始速度模型，初始速度模型各界面均为水平界面.用射线走时反演技术Jive3D (Hobro et al., 2003)，线性迭代反演全地壳二维速度结构.该算法模型参数化时采用离散网格点描述速度和界面，网格单元内的速度值和界面深度值通过三次样条插值获得.该模型参数化方法可准确描述介质速度和界面变化，也可直接利用初至波走时反演的结果，满足对地壳结构高精度成像的要求.正演过程用射线法(Xu et al., 2006, 2010a, 2014)获得各震相计算走时，反演过程采用阻尼最小二乘法，通过拟合计算走时和拾取走时，调整界面形态和各层速度分布，获得最终的二维地壳速度结构.

在10炮地震记录上收集了540个初至走时用于有限差分走时反演.初始模型从地表到深度15 km处，速度由4.0 km·s^-1线性增加到6.0 km·s^-1.经8次迭代后，走时残差由1.1 s收敛到0.09 s后保持稳定.图 3为反演获得的上地壳速度结构及相应射线分布图.测线上除龙门山断裂带下方射线穿透深度较浅(大约5 km)外，初至波穿透深度大部分区域能达到10 km(单元内射线覆盖数>10)，因此反演结果的可信度较高.与前人反演结果相比，除构造单元间的速度差异外，各构造单元内也呈现出更为详细的速度变化形态.

图 3

Fig. 3

图 3 初至波有限差分走时反演获得的地壳上部速度结构(a)及射线路径(b) Fig. 3 Velocity structure of upper crust (a) and ray paths (b) obtained by the finite-difference traveltime inversion using first arrivals

二维全地壳模型计算时，通过多次试验，最终正演和反演过程均采用1 km×1 km的网格单元描述介质.初始模型通过震相分析得到(图 2)，其中深度7 km以上部分采用Pg波走时反演结果.采用“剥皮法”，由上到下拟合壳内各折射和反射震相.图 4为所获得的最终二维地壳速度结构.第i层的下界面由P_i震相确定(称为C_i界面，i=2, 3, 4, 5).为节省篇幅，本文仅给出Sp1, Sp2, Sp4, Sp5, Sp8, Sp9炮的地震记录、最终反演模型的走时拟合和射线分布，分别见图 4—图 10.各炮各震相的平均走时残差均控制在0.2 s以内.最终反演模型上的射线分布(图 11)显示该剖面的射线覆盖足够，走时拟合理想，反演的地壳二维速度结构是可靠的.

图 4

Fig. 4

图 4 全地壳二维速度结构 (a)为高程图及各断裂分布，其中蓝线为检波器高程红色五角星为炮点位置. Fig. 4 2D whole crustal velocity model (a) shows elevation and faults along the profile, with the blue line indicating elevation of receivers and red stars as shots.

图 5

Fig. 5

图 5 “阿坝—遂宁”剖面Sp1炮 (a)地震记录和走时拟合；(b)最终反演模型上射线分布. Fig. 5 Shot Sp1 of the "Aba-Suining" profile (a) Seismic records and fitting between the picked and the calculated traveltimes; (b) Ray paths on the final inversion model.

图 6

Fig. 6

图 6 “阿坝—遂宁”剖面Sp2炮 (a)地震记录和走时拟合；(b)最终反演模型上射线分布. Fig. 6 Shot Sp2 of the "Aba-Suining" profile (a) Seismic records and fitting between the picked and the calculated traveltimes; (b) Ray paths on the final inversion model.

图 7

Fig. 7

图 7 “阿坝—遂宁”剖面Sp4炮 (a)地震记录和走时拟合；(b)最终反演模型上射线分布. Fig. 7 Shot Sp4 of the "Aba-Suining" profile (a) Seismic records and the fitting between the picked and the calculated traveltimes; (b) Ray paths on the final inversion model.

图 8

Fig. 8

图 8 “阿坝—遂宁”剖面Sp5炮 (a)地震记录和走时拟合；(b)最终反演模型上射线分布. Fig. 8 Shot Sp5 of the "Aba-Suining" profile (a) Seismic records and the fitting between the picked and the calculated traveltimes; (b) Ray paths on the final inversion model.

图 9

Fig. 9

图 9 “阿坝—遂宁”剖面Sp8炮 (a)地震记录和走时拟合；(b)最终反演模型上射线分布. Fig. 9 Shot Sp8 of the "Aba-Suining" profile (a) Seismic records and the fitting between the picked and the calculated traveltimes; (b) Ray paths on the final version model.

图 10

Fig. 10

图 10 “阿坝—遂宁”剖面Sp9炮 (a)地震记录和走时拟合；(b)最终反演模型射线分布. Fig. 10 Shot Sp9 of the "Aba-Suining" profile (a) Seismic records and the fitting between the picked and the calculated traveltimes; (b) Ray paths on the final inversion model.

图 11

Fig. 11

图 11 最终反演模型上的地震射线路径覆盖 Fig. 11 Ray paths of all shots and phases on the final inversion model

本文识别震相时参考了同区深地震反射剖面结果(Guo et al., 2013)，通过多次叠加的深地震反射资料可给出高精度的壳内界面形态，有助于减少深地震测深资料中震相识别的不确定性.本文精确考虑了起伏地形对成像结果的影响，将浅层速度结构作为全地壳结构反演的初始模型.所得二维地壳速度结构明显表现出与地表断裂相一致的“横向分块”性，沿剖面可分为三个主要构造单元.下面将分别对各构造单元内部结构进行详细描述和分析.

从速度结构上来看，松潘甘孜块体上地壳有低速的三叠纪复理石沉积，厚度可达10 km (Bruguier et al., 1997).在5~10 km的深度范围内，横向分布的高速异常体与地表断裂呈很好的对应关系(图 12)，这些高速异常体可能是中生代的花岗岩侵入，但未出露地表.低速底界面起伏不平则反映了地层的褶皱变形作用.中下地壳速度较低，反映了青藏高原内部地壳较软弱的特点.地壳整体厚度约为57 km，表现出明显的高原地壳性.

图 12

Fig. 12

图 12 (a)由二维地壳结构推测的动力学模型(黑线为构造断裂，红圈为1978年以来该地区发生的MS>4.0地震，黄色五角星为汶川地震)；(b)断裂走向(深度和距离为同一比例显示)，虚线是推测的断裂；(c)深地震反射剖面构造解释图(蓝线标注构造断裂)(Guo et al., 2013) Fig. 12 (a) Dynamic model inferred based on the 2D velocity structure (black lines are faults; red circles are earthquakes occurred along the profile since 1978 with MS>4.0; yellow star is the Wenchuan earthquake); (b) Fault geometry, dashed lines are inferred faults; (c) Tectonic interpretation of deep seismic reflection profile (blue lines mark faults)

相对四川盆地和松甘地体，龙门山断裂带在地表即呈现明显高速，地表速度达5.5 km·s^-1，与地表分布的彭灌杂岩带相一致.P2、P3界面在龙门山断裂带下方被抬升，与龙门山的抬升过程相一致.中下地壳横向速度变化较大，由西部低速过渡为东部高速.Moho面在龙门山断裂带下方不到100 km的水平距离内有近13 km的抬升，但抬升坡度较缓，未出现“陡坎”状的突然变化.

四川中—新生代前陆盆地属于扬子块体的一部分，广泛分布的沉积岩层呈近水平分布.盆地内的龙泉山断裂是龙门山前陆盆地的东缘边界构造带(邓起东等，1994).四川盆地地壳总厚度约为40 km，从浅部到深部表现为“低速—高速”的层状分布特点.若以5.8 km·s^-1为界，浅部的低速沉积层厚度约为11 km，并以龙泉山断裂为界，表现为西厚东薄的特点.另外，盆地西侧表现为清晰的压扭形态(图 3a)，这与青藏高原东缘的逆冲推覆作用和该区所经受的剥蚀作用有关.青藏高原东部降雨充沛，高陡地形条件下剥蚀作用加剧，加之被逆冲抬升的龙门山的重力作用，断裂带与四川盆地的接触带呈现清晰压扭形态.中下地壳平均速度较高，呈现西薄东厚的特点.龙泉山断裂表现为南东倾向，推测是逆冲推覆构造的速度体现.

综上，松潘甘孜块体和四川盆地地壳结构分别属于较软弱的造山带型和较坚硬的地台型(王椿镛等，2003).两者以龙门山断裂带为界，在地壳厚度、平均速度等方面均存在较大差异.

最终模型整体特征与Jia等(2014)相似，但也存在部分差异，主要表现为：(1) 未观察到松潘甘孜块体5~10 km深度范围内的低速层；(2) 中下地壳速度横向上连续变化，壳内高低速异常体较少；(3) 壳内分层较少.出现上述差异的原因，分别因为：(1) 松潘甘孜块体内的地震记录，在远炮检距依然可以接收到连续的初至震相Pg，说明在浅层初至波可以达到的深度不存在低速层；另因本文初至波反演时严格考虑了起伏地形对射线分布及反演结果的影响，因此上地壳反演精度更高；(2) 所拾取的深层震相(P4, P5, Pm)在各个构造单元内，是平稳渐变的；方法方面，本文所采用的反演方法尽量减少人为参与，在给定初始模型和反演参数后，自动反演，以尽量减少主观因素对模型的影响；(3) 所拾取震相与深地震反射剖面对比，尽量拾取较强能量的反射震相，并通过走时互换对比后使得各炮上均有对应震相，忽略在个别单炮上出现而在其他炮集上不存在的壳内小反射界面，所以壳内分层较少.

从整体来看，Moho面为韧性挠曲的构造形态，中下地壳各层速度横向连续变化，不存在纯构造块体间的突变.结合深地震反射剖面结果(Guo et al., 2013)，认为龙日坝断裂以东的中下地壳为俯冲下插到松潘甘孜软弱地壳内的古扬子块体.印支运动之前松潘—甘孜地区为一古洋盆，与东侧的扬子块体为洋壳-陆壳接触关系，即松潘甘孜块体的东缘作为扬子块体西缘，是一个被动大陆边缘(Harrowfield and Wilson, 2005).受印度板块和欧亚板块的碰撞，青藏高原物质向东运移，造成高原东缘的挤压环境，作为被动边界的松潘甘孜地区较扬子块体更为薄弱，更易受到挤压变形，使得地壳相对增厚，地壳速度减小，并在地表形成一组逆冲推覆断裂带.浅部的速度变化剧烈，速度分块较为破碎，反映了脆性上地壳的逆冲推覆过程.而速度横向连续变化的中下地壳以及韧性挠曲的Moho面标志深部以韧性变形为主的构造特征(图 12a).该观点与Guo等(2013)深地震反射结果相一致(图 12c)，而C4界面可能是上地壳脆性变形和中下地壳韧性变形的分界面.

松潘甘孜块体与邻近的四川盆地地壳厚度差大约为17 km、地貌高度差大约为4 km，从而在横向上形成巨大压力差.图 12a显示，龙门山断裂带域绝大多数的震源位置在20 km以上，属于脆性的中上地壳范围.青藏高原地壳物质向东运移过程中，松潘甘孜块体韧性变形的的中下地壳从底部拖曳着被断裂切割的上部脆性地壳，受到坚硬的四川盆地阻挡，导致应变在不同断裂上积累和释放，诱发大量地震.

本文对深地震测深剖面重新处理，以重建该剖面全地壳速度结构.在识别震相时参照同区域深地震反射剖面结果(Guo et al., 2013)，可合理利用多次叠加后的高精度界面信息减少震相识别的不确定性；浅层速度结构反演时用贴体网格描述介质，精确考虑了起伏地形的影响；在全地壳二维反演中，用精细网格点描述介质，使得全地壳模型得以很好继承浅层初至波反演结果.

本文二维全地壳速度结构表明沿剖面具有“横向分块”的特点，可分为三个主要构造单元：松潘甘孜块体、龙门山断裂带和四川盆地.松潘甘孜块体地表被低速的复理石沉积层覆盖，厚度可达大约10 km；5~10 km的深度范围内，与地表断裂相对应的高速异常体可能是中生代的花岗岩侵入但未出露地表；低速层底界面起伏不平则反映了地层的褶皱变形作用；中下地壳速度横向上连续变化，平均速度较低.四川盆地以龙泉山断裂为界，沉积层西厚东薄，并在西侧出现与挤压和剥蚀作用相关的压扭形态；中下地壳西薄东厚，平均速度较高.龙门山断裂带浅表高速且中上地壳内界面(C2和C3) 被抬升，与龙门山的抬升过程相一致；Moho面在龙门山断裂带下方几十公里的范围内自西向东抬升了大约13 km，速度结构和Moho面深度在龙门山断裂带处均发生较大变化，推测龙门山断裂是一条贯穿全地壳的深大断裂.

沿着整个剖面，Moho面呈韧性挠曲形态，中下地壳速度横向上连续变化.结合深地震反射剖面结果(Guo et al., 2013)，认为龙日坝断裂以东的中下地壳为俯冲到松潘甘孜软弱地壳内的古扬子块体.

在青藏高原向东的挤压环境下，松潘甘孜块体下方韧性变形的中下地壳在底部拖曳着被断裂切割的脆性上地壳向东运移，受坚硬四川盆地的阻挡，沿着龙门山断裂向上挤出，应力在不同断裂上积累和释放，诱发大量地震.

非常感谢主持并实施野外地震数据采集的中国地震局地球物理勘探中心共享地震数据；感谢中国地震局地球物理所王椿镛研究员和中国地质科学院地质所李文辉助研究员的有益建议和讨论.