地球物理学报  2020, Vol. 63 Issue (11): 3996-4011   PDF    
印度板块下地壳北向俯冲与榴辉岩化的地震学证据:接收函数成像结果
武振波1, 唐国彬2, 徐涛2,3, 梁春涛1,4, 翁雪飞1     
1. 成都理工大学地球物理学院, 成都 610059;
2. 中国科学院地质与地球物理研究所, 岩石圈演化国家重点实验室, 北京 100029;
3. 中国科学院青藏高原地球科学卓越创新中心, 北京 100101;
4. 成都理工大学, 地球勘探与信息技术教育部重点实验室, 成都 610059
摘要:印度地壳与岩石圈地幔的俯冲前缘和俯冲形态,对认识高原构造变形、隆升机制有重要意义.本文基于青藏高原西缘分布的流动宽频带地震台站(TW-80测线和Y2台网)记录的远震波形数据,通过接收函数H-κ网格搜索与CCP叠加方法,对研究区地壳结构进行成像.结果显示:(1)研究区西侧北西—南东向剖面(剖面1,2),狮泉河逆冲断裂带以南,深度67~80 km范围内均观测到连续的Moho界面;40~55 km范围内存在另一组横向上可连续追踪的界面,其形态与之下Moho面横向变化趋势近乎平行;(2)研究区东侧剖面3下方,Moho面从南端喀喇昆仑断裂带下方向北逐渐加深,在雅鲁藏布江缝合带附近增至大约67 km,进入拉萨块体至台站WT20和WT03下方至最深75~80 km,然后向北有所抬升.基于成像结果和岩石学研究成果推测藏南块体下方,自西向东均存在俯冲印度板块下地壳的榴辉岩化现象,可以用来指示印度板块地壳尺度的俯冲前缘,其在青藏高原西部(约80°E)位于班公湖—怒江缝合带附近,向东逐步递减至拉萨块体中部.
关键词: 青藏高原      接收函数      印度板块俯冲      榴辉岩化     
Eclogitization of Indian lower crust during its northward subduction revealed by receiver function depth migration images
WU ZhenBo1, TANG GuoBin2, XU Tao2,3, LIANG ChunTao1,4, WENG XueFei1     
1. College of Geophysics, Chengdu University of Technology, Chengdu 610059, China;
2. State key Laboratory of Lithospheric Evolution, Institute of Geology and Geophysics, Chinese Academy of Sciences, Beijing 100029, China;
3. CAS Center for Excellence in Tibetan Plateau Earth Sciences, Beijing 100101, China;
4. Chengdu University of Technology, Key Laboratory of Earth Exploration and Information Techniques of the China Ministry of Education, Chengdu 610059, China
Abstract: The receiver functions calculated at the broadband seismic stations of TW-80 and Y2 array are used to reveal the main subduction features of the Indian plate beneath western Tibet. The H-κ grid-searching results and CCP depth migration image show that the observed Moho structure at depth 67~80 km is laterally traceable from south to north on the west profiles in the study area, while there is another discontinuity interface at depth 40~55 km almost parallel to the Moho up to the station SQAT nearby the Shiquanhe thrust fault. This feature also can be seen on the CCP migration image on the east profile of the study area. Based on the evidence derived from this and previous petrology studies, we infer that the Indian lower crust front beneath western Tibet is located at Bangonghu-Nujiang suture, while it gradually reduces southward to the internal Lhasa block beneath central and east Tibet.
Keywords: Tibetan plateau    Receiver function    Northward subducting Indian plate    Eclogitization    
0 引言

始新世以来,印度—欧亚大陆的碰撞汇聚造就了世界上规模宏大的青藏高原(Yin and Harrison, 2000).GPS观测结果显示,现今两板块之间的平均汇聚速率约为40 mm·a-1(Gan et al, 2007),而印度与雅鲁藏布江缝合带之间的喜马拉雅区域,其相对汇聚速率约为20 mm·a-1 (Bilham et al., 1997).根据已有研究,持续的碰撞/汇聚已造成印度/欧亚大陆自初始碰撞以来,累积缩短量在高原西部为1800~2100 km,中部为2475 km,东部达到2750~2800 km(Johnson, 2002).高原隆升及不同的东西向缩短量,应该由其内部及周边构造活动、壳幔形变及印度板块自西向东的差异性俯冲协同调节.

关于印度岩石圈板片北向俯冲距离及机制,不同研究方法获得的研究结果存在明显差异.比如,Chen等(2017a)利用体波伴随成像方法认为其向北俯冲至金沙江缝合带附近.Tilmann等(2013)基于层析成像结果推测印度岩石圈板片向北俯冲至班公湖—怒江缝合带,然后倾角变陡近乎垂直地向下进入深部地幔.Kosarev等(1999)基于接收函数成像认为其向北俯冲至羌塘块体下方.而另有部分研究则认为印度岩石圈板片向北只俯冲至藏南块体下方(刘葵等,2006吕庆田等,1998).此外,体波走时层析成像(Li et al., 2008; Bijwaard and Spakman, 2000)、接收函数成像(Zhao et al., 2010)、横波分裂研究(Chen et al., 2015)、Pn波走时层析成像(Li et al., 2008)等研究结果已初步揭示出印度岩石圈板片在青藏高原下方俯冲的前沿位置(图 1)、几何形态和运动方向等存在较为显著的东西向差异,这种差异性可能反映了自西向东俯冲角度的变化,俯冲角度差异加剧了运动学差异,可能进一步造成印度岩石圈板片发生撕裂和断离(Chen et al., 2015; Li and Song, 2018; Liang et al., 2016; Wu et al., 2019a, b; Si et al., 2019).但不同学者提出的撕裂模型也存在差异,比如Chen等(2015)认为印度板块的俯冲角度从东部的陡倾角向西逐渐变缓,而Li和Song(2018)获得的撕裂模型则显示印度板块在青藏高原西部和东部表现为平缓俯冲,而高原中部为陡倾角俯冲.

伴随印度岩石圈板片的向北俯冲,其浅部地壳又有何响应特征?研究表明,其中上地壳通过逆冲断层向南向上折返(Gao et al., 2016)、或以中地壳通道流(Beaumont et al., 2001)的方式转换为喜马拉雅造山带的一部分,而其下地壳则随岩石圈地幔一起向北俯冲.Xie等(2017)根据数条大地电磁剖面观测数据反演的电阻率剖面成像结果,推测印度下地壳俯冲前缘在高原西部约80°E处到达约33.5°N,中部约85°E处到达约31°N,东部约87°E和约92°E处递减至约30.5°N.由分布于青藏高原中部和东部的几条宽频带测线获得的接收函数成像结果显示的“双莫霍”特征,推测印度板块下地壳在向北俯冲过程中可能发生了不同程度的榴辉岩化(Schulte-Pelkum et al., 2005; Xu et al., 2015; Nábělek et al., 2009; Kind et al., 2002; Shi et al., 2015).

在青藏高原西缘,基于TW-80宽频带测线的接收函数成像结果也推测印度大陆板块向北俯冲过程中伴随下地壳发生榴辉岩化(Zhang et al., 2014);此外,结果还显示出剖面下方莫霍面具有分区特征,并与地表主要断裂带的分段性具有较好对应,推测该区域地表的一些主要断裂(如喀喇昆仑断裂带、龙木措—郭扎断裂带)可能一直向下延伸至Moho面附近.然而,该区域Y2台网和羚羊计划西测线的部分宽频带台站获得的接收函数成像图上则显示Moho面在雅鲁藏布江缝合带两侧没有错断,呈现南北向连续性特征,指示喀喇昆仑断裂带并未切穿地壳(Xu et al., 2017).而东侧约81.5°E处的南北向深反射地震剖面进一步认为该断裂带仅向下延伸至15km(Gao et al., 2016).此外,Xu等(2017)Gilligan等(2015)还观测到特提斯—喜马拉雅和拉萨块体的中地壳存在低速异常区,其在Zhang等(2014)的研究结果中未提及.为进一步探讨在青藏高原西缘,喀喇昆仑断裂带向下延伸深度、藏南块体地壳结构以及榴辉岩化印度下地壳的前缘位置,本文统一处理了TW-80剖面和Y2台网记录的远震数据,提取出更多有效的高质量P波接收函数,进行共转换点偏移成像处理,同时结合藏西缘近年来最新的体波层析成像结果和数值模拟研究深入探讨印度大陆板块在此处可能的俯冲模式.

1 构造背景与地震台站分布 1.1 构造背景

从地貌上看(图 1),研究区内高原南北向宽度窄于600 km,印度和欧亚大陆板块在此处相互作用强烈,地表构造十分复杂.从图 2构造上看,该区域由南至北包含的构造块体依次为:特提斯—喜马拉雅块体、拉萨块体、羌塘块体和甜水海块体.由南至北的缝合带/断裂带依次为:喀喇昆仑走滑断裂带、雅鲁藏布江缝合带、狮泉河逆冲断裂、曼东—措北逆冲断裂、班公湖—怒江缝合带、多玛—乌江逆冲断裂、龙木措—郭扎走滑断裂,构造线走向均呈北西-南东向,且都可以和高原东部的主要构造带相接.此外,雅鲁藏布江、班公湖—怒江和狮泉河缝合带/断裂带内均发育有蛇绿岩,代表特提斯大洋岩石圈残片(Lacassin et al., 2004).

图 1 地震台站与印度大陆板块北向俯冲前缘位置 剖面A为Li和Song(2018)基于Pn波成像推测的印度岩石圈板片前缘位置;B、C分别为Bijwaard和Spakman(2000)Li等(2008)基于体波成像推测的200 km深度处印度岩石圈板片前缘位置;D、E、F分别为正文中讨论的Chen等(2017b)Huangfu等(2018)Zhou和Murphy(2005)剖面位置;黑色粗实线为根据接收函数推测的印度下地壳发生榴辉岩化的位置;绿色实心圆点为Xie等(2017)根据大地电磁数据反演结果推测的印度下地壳前缘位置;淡绿色箭头为Zhao等(2010)推测的印度、欧亚板块岩石圈前缘位置.红色三角形:TW-80测线宽频带台站;蓝色菱形:中美合作布设的Y2台阵宽频带台站;灰色菱形:部分已发表成果用到的宽频带探测剖面;红色虚线框:图 2所示研究区域;Tarim basin:塔里木盆地;Qaidm basin:柴达木盆地;SG:松潘—甘孜块体;QB:羌塘块体;LB:拉萨块体;HB:喜马拉雅块体;India plate:印度板块;ATF:阿尔金断裂带;KF:喀喇昆仑断裂带; LGF:龙木措—郭扎断裂带;BNS:班公湖—怒江缝合带;JS:金沙江缝合带;AKMS:阿尼玛卿缝合带; IYS:雅鲁藏布江缝合带. Fig. 1 Map showing broadband seismic stations and northern front of the subducting Indian plate The northern front of subducting Indian plate is shown by profiles A, B, and C, which are respectively derived by Li and Song (2018), Bijwaard and Spakman (2000) and Li et al. (2008). Line D, E, and F denote the seismic profiles later discussed in this study, which respectively referred to the study of Chen et al. (2017), Huangfu et al. (2018) and Zhou and Murphy (2005). The bold black lines represent the horizontally projected location of the partially eclogitized lower crust of the Indian plate inferred from the receiver functions study. Green dots are the varying Indian crustal front revealed by magnetotelluric study of Xie et al. (2017). The light green arrows respectively are the Indian and Asian lithospheric front revealed by receiver functions study of Zhao et al. (2010). Red triangles: TW-80 seismic stations. Blue diamonds: seismic stations of Y2 array. Grey diamonds: other deployed broadband seismic stations. Red dashed box shows the study area in Fig. 2. SG: Songpan-Garzê block; QB: Qiangtang block; LB: Lhasa block; HB: Himalayan block; ATF: Altyn Tagh fault; KF: Karakorum fault; LGF: Longmucuo-Guozha fault; BNS: Bangonghu-Nujiang suture; JS: Jinsha river suture; AKMS: Anyimaqin Kunlun Mustagh Fault. IYS: Indus-Yarlung suture.
图 2 研究区宽频带地震台站与远震事件分布 剖面1, 2, 3为共转换点偏移叠加剖面位置.绿色十字为接收函数在深度80 km处转换点对应的地表投影.ATF:阿尔金断裂带;LGF:龙木措—郭扎断裂带;DWF:多玛—乌江逆冲断裂;MCF:曼冬—措北逆冲断裂;SQHF:狮泉河逆冲断裂;KF:喀喇昆仑断裂带;MCT:主中央逆冲断裂带;MBT:主边界逆冲断裂带;IYS:雅鲁藏布江缝合带;BNS:班公湖—怒江缝合带;JS:金沙江缝合带;TSH:甜水海块体;QB:羌塘块体;LB:拉萨块体;HB:特提斯—喜马拉雅块体. Fig. 2 Broadband seismic stations in the study area and the distribution of teleseismic events Line 1, 2 and 3 respectively show the reference profiles used by receiver function depth migration. The green crosses are the projected positions of piercing points of Ps converted rays at the Moho based on reference depth 80 km. ATF: Altyn Tagh fault; LGF: Longmucuo-Guozha fault; DWF: Domar-Wujiang thrust fault; MCF: Mandong-Cuobei thrust fault; SQHF: Shiquanhe thrust fault; KF: Karakorum fault; MCT: Main Central Thrust; MBT: Main Boundary Thrust; IYS: Indus-Yarlung suture; JS: Jinsha river suture; BNS:Bangonghu-Nujiang suture; TSH: Tianshuihai block; QB: Qiangtang Blcok; LB: Lhasa block; HB: Himalayan block.

从深部结构特征来看,研究区西北部存在塔里木盆地向青藏高原下方俯冲的证据(李秋生等,2000Wittlinger et al., 2004Levin et al., 2013).在西昆仑山脉下方,岩石圈地幔同时存在向南和向北倾的结构(Kao et al., 2001Gao et al., 2000);这种“双向俯冲”特征可能指示青藏高原西北部与塔里木盆地岩石圈存在类似藏南的陆陆碰撞作用.

喀喇昆仑断裂带是分布于青藏高原西缘的一条大型右旋走滑断裂.Tapponnier等(1982)基于物理模型试验提出青藏高原内部不同拼贴块体在印度板块持续性俯冲作用下发生“构造逃逸”,认为北以阿尔金断裂带为界,南以喀喇昆仑—嘉黎断裂带(位置如图 1所示)为界,青藏高原内部刚性块体沿边界走滑断裂带向东挤出.假设该动力学模型成立,相应的观测现象主要应该有两个:(1)高原周边及内部大型走滑断裂带的滑动速率、两侧块体长期以来积累的相对滑动量应该很大;(2)高原内部不同块体之间的边界走滑断裂带均为延伸至岩石圈地幔的深大断裂.那么,青藏高原西缘喀喇昆仑断裂带的滑动速率是多少?累积滑动量是多少?向下延伸有多深?这些问题均关系到上述“构造逃逸”模式的现实合理性.Chevalier(2019)综述了有关喀喇昆仑断裂带的活动构造研究进展,认为从现今的大地测量学尺度到几个百万年的地质学时间尺度,其走滑速率变化范围为3~10 mm·a-1;但也有研究表明其滑动速率高达32 mm·a-1(Valli et al., 2008).沿该断裂带估测的累积滑动量也争议颇大,从几乎没有(Jain and Singh, 2008)到少量的25 km(Searle et al., 2010)、66~150 km(Murphy et al., 2000)、40~150 km(Phillips et al., 2004)再到149~167 km(Robinson, 2009),而Replumaz和Tapponnier(2003)估计约为555 km.关于其深部延伸范围也尚未达成统一,部分研究结果表明其可能为切穿整个地壳的深大断裂(Zhang et al., 2014; Rolland et al., 2009; Klemperer et al., 2013),但也有学者认为其仅局限于中上地壳尺度(Xu et al., 2017; Gilligan et al., 2015; Gao et al., 2016).因此,关于喀喇昆仑断裂带在青藏高原动力学隆升扩展过程中的作用,仍需更多研究和探讨.

1.2 宽频带台站分布

图 2显示本研究用到的宽频带地震台站位置.其中,红色三角形为TW-80剖面的台站,共19台.地震仪型号为CMG3-ESP,数据采集器型号为Reftek-72A,仪器响应的平坦周期范围50 Hz~30/60 s;平均台间距15~20 km,剖面全长400多km,记录时间从2011—2013年,采样频率为40 Hz.中美合作布设的Y2宽频带台阵包含30个地震仪,地震仪类型为STS-2,采集器类型为Q330,仪器响应的平坦周期范围为10 Hz~120 s;平均台间距为40~100 km,采样率为100 Hz;分两期布设,其中10个台站的记录时间为2007—2011年,另20个台站的记录时间为2009—2011年上半年(Gilligan et al., 2015).由于有些台站计算的有效接收函数不多,加之研究区域地下地质结构复杂,为了尽可能避免人为误差,舍弃了难以正确区分和连续追踪有效转换波震相的台站数据.图 2中展示了本文用到的Y2台阵中21个有效地震台站(蓝色菱形).此外,右上角还显示了用到的远震事件分布,其大多为环太平洋地震,故而射线路径也大多从台站东侧入射.基于IASP91参考模型计算的射线在深度80 km处转换点对应的地表投影位置如图 2中绿色十字所示,也大多分布于台站东侧.

2 研究方法 2.1 接收函数计算

接收函数是研究地壳上地幔结构一种典型地球物理探测手段(Kind et al., 2002; Chen and Ai, 2009; 吴庆举和曾融生,1998刘启元等,1996胡家富等,2005李永华等,2006王兴臣等,2015杨传茂等,2019陈俊磊等,2019).针对图 2所示宽频带地震台站记录的远震波形数据,依据USGS-NEIC发布的地震目录挑选MS≥5.0、震中距30°~90°范围内的远震事件.对提取的原始波形数据进行去均值、去线性趋势、剔除三分量记录不全或者波形明显异常的事件,进行0.01~2 Hz频率范围内的巴特沃斯带通滤波.然后,基于IASP91参考模型计算理论P波到时,截取其前20 s和之后85 s的波形数据,采用时间域迭代反褶积算法(Ligorría and Ammon, 1999)计算P波接收函数.其中,高斯系数设为2.5,对应的高频截止频率为1.2 Hz(给定频率域高斯函数,规定G(f)=0.1时,对应频率f即为截止频率).对单台所有事件计算的接收函数进行手动挑选,摒弃转换波震相到时与相邻观测波形明显不同、波形明显异常、直达P波不明显的记录,单台剩余的有效接收函数数目分布于51~377条之间,图 2中所示台站总共获得7953条接收函数,为后续的H-κ网格搜索和共转换点偏移叠加处理提供了良好数据基础.

2.2 动校正

在地震波信号处理当中,叠加技术是一种简单高效增强弱信号的常规处理手段.对于单个台站获得的所有有效接收函数而言,由于远震事件方位角和震中距分布的不均匀性,其相应的射线传播路径、射线参数也各不相同,造成来自台站下方同一深度处速度间断面产生的转换波震相Ps与直达P波的到时差也存在相应变化,因而不能直接叠加.Yuan等(1997)将人工源地震勘探中的动校正技术(Moveout)引入接收函数研究.做动校正之前需要给定参考速度模型和参考慢度值.针对不同射线路径在深度z处产生的Ps转换波相对直达P波的走时,与参考路径(射线参数假定为p0)在深度z处产生的转换波相对于直达P波走时的时间差表示为

(1)

根据此时间差可对不同震中距计算的接收函数做动校正,然后进行线性叠加.同理,也可针对两个多次波PpPs、PsPs+PpSs进行类似动校正.为了进一步阐明该方法,设定一个只包含单层地壳的速度模型(地壳厚度50 km,VP=6.5 km·s-1VS=3.75 km·s-1;地幔纵横波速度分别为8.04 km·s-1、4.47 km·s-1),利用PROGRAMS 330软件包计算理论地震记录,然后同样采用时间域迭代反褶积算法(Ligorría and Ammon, 1999)计算理论P波接收函数,结果如图 3最左侧波形.基于上述方法,分别给定两个参考射线参数:p0=0.00 s·km-1,即假定射线从台站下方垂直入射;p0=0.07 s·km-1;然后对每一个参考值,分别针对莫霍面产生的转换震相Ps及其两个多次波震相做动校正,校正后的接收函数波形如图 3ab中右侧三张图所示;可以看出:(1)针对某一特定震相做动校正后,只有该震相会实现同相对齐,其他震相则由于整个波形沿着时间轴拉伸或压缩发生畸变而层次不齐;(2)在给定速度模型的情况下,针对不同参考射线参数或参考慢度值,动校正后的波形在时间轴上的对齐时间也不同.也就是说,动校正后的接收函数显示的转换波相对于直达P波的走时只对应参考慢度有效.

图 3 合成接收函数针对不同参考慢度和不同震相动校正后的波形 (a)左图为理论合成接收函数波形,最右侧数字为相应的射线参数值;给定参考射线参数p0=0.0 s·km-1,分别对深度50 km处产生的转换波Ps及其两个多次波震相(PpPs、PpSs+PsPs)进行动校正,校正后的波形依次为后面三个图.红色、蓝色及绿色虚线分别为根据参考p0计算的Ps震相及其多次波震相的理论到时;(b)除p0=0.07 s·km-1外,其余与(a)类似. Fig. 3 Waveforms of synthetic receiver functions before and after move-out correction (a) The left panel is synthetic receiver functions with different ray parameters shown on the right of each waveform. The Moho depth is set 50 km. Given the reference ray parameter p0=0.0 s·km-1, which means the vertical incidence of P wave on the earth surface, the moveout corrected waveforms respectively for the phase Ps, PpPs and PpSs+PsPs are sequentially shown on the right panels. The red dashed line indicates theoretical arrival time of Ps, blue for PpPs arrival time, and green for PpSs+PsPs arrival time; (b) All symbols and indications are similar to interpretations of (a), except for p0=0.07 s·km-1.

图 4展示了台站AL06获得的所有有效接收函数,从下往上按照震中距由小往大排列,左图对应原始接收函数,右图为给定p0=0.06 s·km-1、IASP91参考速度模型时,针对Ps震相动校正后的波形.假设台站下方的平均地壳厚度为70 km、VP为6.3 km·s-1、纵横波速度比为1.73,关于三个震相(Ps、PpPs、PsPs+PpSs)计算的理论到时分别如图中红色、蓝色和绿色虚线所示.对比左右图中顶部的叠加波形,动校正后的叠加波形显示Pms震相(位于9~10 s之间)稍有增强,可能是由于该台站下方壳内结构相对简单,计算的原始接收函数波形受噪声干扰小,其上即可连续追踪清晰的Pms震相.

图 4 AL06台站计算的原始接收函数(a)及其动校正后的波形(b),设置p0为0.06 s·km-1 Fig. 4 Receiver functions at station AL06 before(a) and after moveout correction(b)with given p0 0.06 s·km-1
2.3 H-κ网格搜索

Zhu和Kanamori(2000)提出H-κ网格搜索方法来估计接收台站下方的平均地壳厚度和纵横波速度比,并指出估测的地壳厚度对VP变化不敏感,而对VP/VS的取值更敏感.H-κ网格搜索方法的基本原理为设定H和κ的搜索范围和步长、地壳平均P波速度,那么台站记录的每条接收函数对于给定的每组(Hκ),就可以提取三个震相理论到时处的振幅值,记为s1(H, κ)、s2(H, κ)和s3(H, κ).由此获得函数:

(2)

其中wi为各个震相的加权值,且∑wi=1,N为台站记录的接收函数数目.本文设定H搜索范围为30~80 km,间隔0.05 km;κ搜索范围为1.5~2.0,间隔0.002;平均P波速度为6.3 km·s-1;震相Ps、PpPs、PsPs+PpSs的权重分别为0.7、0.2和0.1.图 5显示了台站WT20和AL04的H-κ网格搜索结果.台站WT20位于拉萨块体,其记录的接收函数中Pms及多次波PpPs震相相对清楚,而震相PsPs+PpSs的振幅较弱,扫描获得的平均地壳厚度为75.25 km,κ为1.78.台站AL04位于喀喇昆仑断裂带附近,其下方平均地壳厚度为71.4 km,κ为1.72.可见,H-κ网格搜索方法要求至少有一个较为清楚的多次波震相.实际情况中,由于存在沉积层、倾斜的壳内速度间断面或各向异性介质,对多次波震相干扰较大,使得H-κ网格搜索方法的结果不稳定或无效.本研究中,只获得了32个台站的扫描结果,如图 7a8a中白色圆点所示.

图 5 台站WT20,AL04的H-κ网格搜索结果 第一行为台站WT20的结果;针对不同的(H, κ),对台站记录的每条接收函数可以计算其Ps转换波相对于直达P波的理论走时,提取该走时处的振幅值,然后求和,即可形成图(a);同理,也可以提取两个多次波震相针对不同(H, κ)的理论走时和振幅,并求和,形成图(b)、(c);然后,对三个震相的求和结果进行加权叠加,获得图(d);白色五角星为基于最后的加权叠加结果s(H, κ)扫描获得的振幅最大值点,其对应的H即为平均地壳厚度,κ为平均地壳纵横波速度比.第二行为台站AL04的结果;(e)、(f)、(g)、(h)的含义分别于(a)、(b)、(c)、(d)对应相同. Fig. 5 H-κ grid-searching results at stations WT20 and AL04 The figures on the top are the H-κ grid-searching results at station WT20. Given a pair of (H, κ), the amplitude of Ps phase in each receiver function at one station can be calculated theoretically, then all of which can be stacked to plot the image (a). Similarly images for two multiple phases are generated respectively shown in (b) and (c). Finally, the data used in (a), (b) and (c) can be weighted and stacked to form image (d). The similar results at station AL04 are shown in the bottom.
图 6 台站AL04和WT20计算的有效接收函数 红线分别表示由图 5所示结果计算的每条接收函数中Pms、PpPs及PpSs+PsPs震相的理论到时. Fig. 6 Effective receiver functions obtained at stations AL04 and WT20 The red lines indicate the theoretical arrival time calculated from the results shown in Fig. 5, respectively for Pms, PpPs and PpSs+PsPs phases.
图 7 图 1中剖面1和2的共转换点偏移叠加结果,采用参考速度模型IASP91 (a)剖面1和剖面2的偏移成像图;白色圆点为H-κ网格搜索获得的平均地壳厚度;(b)对应剖面1所用台站的平均叠加接收函数,从左往右按照台站纬度由南至北排列;(c)对应剖面2所用台站的平均叠加接收函数,并按照台站纬度由南至北排列.白色虚线分别示意榴辉岩化印度下地壳的顶界面和Moho界面.KF:喀喇昆仑断裂带;SQHF:狮泉河逆冲断裂带;BNS:班公湖—怒江缝合带;LGF:龙木措—郭扎断裂带. Fig. 7 Images of receiver function depth migration on profiles 1 and 2 based on the reference model IASP91 (a) Depth migration images on profiles 1 and 2. The white circles indicate the average crustal thickness given by H-κ grid-searching beneath each seismic station; (b) The average stacked receiver functions at all stations nearby profile 1 are aligned from south to north based on seismic station latitudes; (c) The average stacked receiver functions at all stations nearby profile 2 are aligned from south to north based on seismic station latitudes. Two dashed white lines schematically indicate the top interface of Indian eclogitized lower crust and Moho structure. LGF: Longmucuo-Guozha fault; SQHF: Shiquanhe thrust fault; KF: Karakorum fault; BNS:Bangonghu-Nujiang suture.
图 8 图 1中剖面3共转换点偏移叠加结果,采用参考速度模型IASP91 (a)剖面3偏移成像图;白色圆点为H-κ网格搜索获得的平均地壳厚度;(b)对应剖面3所用台站的接收函数经过动校正后的叠加波形,且按照台站纬度由南至北排列.白色虚线分别示意榴辉岩化印度下地壳的顶界面和Moho界面.IYS:雅鲁藏布江缝合带;BNS:班公湖—怒江缝合带. Fig. 8 Images of receiver function depth migration on the profile 3 based on the reference model IASP91 (a) Depth migration images on profile 3. The white circles indicate the average crustal thickness given by H-κ grid-searching beneath each seismic station; (b) The average stacked receiver functions at all stations nearby profile 3 are aligned from south to north based on seismic station latitudes. BNS: Bangonghu-Nujiang suture; IYS: Indus-Yarlung suture.
3 接收函数偏移叠加成像

本文主要采用了Zhu等(2000)提出的共转换点偏移叠加成像方法.其基本原理是:首先,提供一个研究区的背景速度模型,可以用全球标准平均速度模型作为参考,也可以采用其他研究方法获得的局部区域性速度模型.其次,通过入射角效应(对振幅而言)的校正,将接收函数时间序列上的每一个离散点都认为对应一个Ps转换震相,基于其与直达P波的走时差,将其向下追踪至相应转换点深度并依次排列在追踪的射线路径上.然后,在三维空间设定叠加区域大小,将每个叠加单元内包含的振幅值进行叠加,即可获得三维叠加数据体.此次研究中,设定叠加单元沿剖面方向和垂直方向均为1 km,垂直剖面方向为50 km.图 2中剖面1和2的成像结果如图 7所示,剖面3的成像结果如图 8所示.

剖面1所涉及台站共26个,每个台站记录的有效接收函数经动校正后的叠加波形由南向北排列,如图 7b所示.整体上看,来自Moho面的转换波相对于直达P波的走时由南至北逐渐增大,AL03台站下方增至最大值10 s左右,然后向北在台站AL04、WT13、AL05、AL06、AL67、WT14下方约为9 s,而台站WT15在9~10 s范围有两个波峰,其北侧SQAH的波形则更为复杂,根据图 2可知此两个台站位于狮泉河逆冲断裂带附近,可能其下方存在倾斜界面等复杂地质结构.狮泉河逆冲断裂带以北,WT16、AL09、AL10的Pms震相走时则逐渐减小至7~8 s,AL11则在7 s和9 s附近各有一个波峰,RUTK位于班公湖—怒江缝合带附近,其在9 s附近观测到清晰的Pms震相.在图 2中,AL12和RUTK台站的地表位置非常接近,但可能由于前者记录时间较短,叠加波形干扰比较大,Pms震相难以准确辨认.

剖面2与剖面1的走向相同,其所涉及台站包括AL15、NOCO、AL16、AL21、AL22和AL24.如图 7c所示,台站AL15的叠加接收函数波形在1~2 s范围内显示很强的浅层转换波震相,可能表明台基下方存在一定厚度的沉积层,该震相的多次波信号对后续来自Moho界面的转换波震相造成严重干扰;剩余其他台站叠加波形的形状也很复杂,可能与此处龙木措—郭扎断裂带的复杂构造运动有关.

剖面3位于研究区域东部,所用台站由南至北依次为PURG、WT04、WT18、WT20、WT03、WT05、WT06、WT09和NOMA,其经动校正后的叠加波形如图 8b所示.Pms震相的走时变化趋势与剖面1结果类似,从南侧的7.8 s左右向北逐渐加深,至台站WT20、WT03下方增至最大值10 s左右,然后向北又逐渐减小为9 s左右.综合这些剖面的接收函数波形,不难看出:在青藏高原西缘,Moho界面的整体变化趋势为从南部的特提斯—喜马拉雅块体下方向北逐渐加深,至拉萨块体中部达到最深,然后向北又逐渐变浅.这一点,与之前观测结果存在一致性(Zhang et al., 2014).

4 讨论 4.1 印度下地壳北向俯冲与榴辉岩化

图 7a来看,综合剖面1和2的叠加成像结果显示:台站SQAH以南,深度67~80 km范围内均观测到连续的Moho界面;深度40~55 km范围内存在另一组横向上可连续追踪的界面,形态与之下的Moho面横向变化趋势近乎平行,类似特征在Zhang等(2014)的接收函数成像图上也有所体现.图 8a中尽管剖面3的成像结果所用台站较少,对地壳内部速度间断面的成像效果变差,导致分辨率降低,但也能够看到Moho界面之上存在另一组近平行的界面.此外,Xu等(2017)利用羚羊计划中西部测线的部分台站对Y2台阵进行了补充,他们的成像结果中接近本文剖面3的位置,显示出Moho之上20 km左右存在一组清楚的壳内速度间断面.那么如何解释这个界面?它有何构造意义?探讨如下.

Searle等(2011)指出印度板块下地壳向北俯冲至60~80 km深度范围内时,由于压力增大,其前寒武纪麻粒岩会向高压相-麻粒岩或榴辉岩相变.而高压相麻粒岩捕虏体在藏南和羌塘块体均有发现(Chan et al., 2009; Hacker et al., 2000);榴辉岩相变则主要来源于远震P波接收函数的间接观测证据.如图 1中黑色粗实线所示,自西向东分布的南北向宽频带台站探测剖面,在其接收函数成像结果中均发现拉萨块体下方Moho界面之上存在一组清晰的壳内间断面(Zhang et al., 2014; Xu et al., 2015, 2017; Nábělek et al., 2009; Schulte-Pelkum et al., 2005; Kind et al., 2002; Shi et al., 2015),结合各向异性分析、横波速度结构、矿物高压实验等多种技术手段,此壳内间断面被解释为向下俯冲的印度下地壳发生榴辉岩化产生的相变顶界面,其向南可进一步与Caldwell等(2013)接收函数成像结果中的MHT连接起来(Zhang et al., 2014).

另外,Xie等(2017)基于在青藏高原实施的数条大地电磁探测剖面所获得的电阻率成像结果,显示此处中下地壳存在不连续分布的低电阻率异常,可能表明相应区域内存在含水流体,而Hetényi等(2007)的研究表明流体含量增加对榴辉岩化有一定的促进作用.武振波等(2016)获得的TW-80测线下方的横波速度结构显示下地壳介质地震波速度偏高,这也与榴辉岩化后带来的岩石密度和地震波传播速度增大相一致.

综合这些观测证据,包括:研究区域内存在的加厚地壳特征(最深至~80 km)会造成下地壳压力增大,低电阻率特征可能表明的含水流体有利于促进榴辉岩相变以及下地壳的高波速特征,我们推测图 7a图 8a中观测到的另一组与下伏Moho面近似平行的界面为向北俯冲的印度下地壳榴辉岩相变顶界面.从图 1可以看出,自西向东分布的宽频带地震剖面均观测到拉萨块体下方印度板块下地壳的榴辉岩化现象,其北端前沿位置与根据大地电磁反演结果推测的印度下地壳俯冲前沿吻合良好,因此可以用来指示印度地壳的北向俯冲前缘,其在高原西部到达班公湖—怒江缝合带附近,向东则位于拉萨块体内部.

4.2 喀喇昆仑断裂带深部延伸尺度

Zhang等(2014)研究结果认为喀喇昆仑断裂带、龙木措—郭扎断裂带均向下延伸到下地壳尺度或Moho界面.然而,Xu等(2017)研究认为喀喇昆仑断裂带仅延伸至中上地壳尺度,Gao等(2016)深反射地震剖面的研究结果进一步认为该断裂带仅向下延伸≤15 km.那么,导致观测结果不同的原因可能有:一般情况下接收函数共转换点偏移成像结果,在处理时存在叠加区域的横向平滑.横向平滑距离越大,成像结果越光滑,但会降低剖面横向分辨率或者形成假象.另一方面,偏移成像结果存在一定程度上的“沿垂直剖面方向上的横向压缩”,也就是说展示的剖面成像结果实际上是剖面两侧一定横向距离范围内的转换点振幅叠加.此外,由于宽频带台站分布并不是严格垂直于构造线方向(如图 2所示),当偏移成像处理时定义的剖面走向不同,成像图也会发生一定变化.Zhang等(2014)定义的偏移剖面近南北走向,而Xu等(2017)定义的偏移剖面走向为北西-南东方向.

为了探讨上述关于喀喇昆仑断裂带深部延伸尺度不一致性,本文中的西侧偏移剖面(图 2中剖面1、2)也取为北西-南东方向,并将沿剖面方向的叠加网格设为1 km,垂直剖面方向设为50 km(小于Xu等设置的200 km).从图 7a成像结果来看,喀喇昆仑断裂带附近台站MONS、AL00和WT08下方的Moho面横向变化连续,但WT08台站下方缺失明显的“印度下地壳顶界面”,可能是由于喀喇昆仑断裂带的影响,不过也不能排除台站本身记录事件波形、噪声干扰、复杂地质结构、偏移叠加成像方法处理等其他因素.如果能够沿着严格垂直于喀喇昆仑断裂带构造走向的方向布设一条密集台站观测剖面,提高台站覆盖率和横向分辨率,或许能够更好地探测该断裂带在地下的延伸产状及深度.

另外,Xu等(2017)Gilligan等(2015)发现在青藏高原西缘的特提斯—喜马拉雅块体和拉萨块体下方存在中地壳低速异常区.本文中剖面1的成像结果上虽然也显示有横向上可连续追踪的负极性界面,但如何排除其是否为滤波带来的旁瓣效应还需要进一步探讨.

4.3 印度板块在青藏高原西缘的俯冲特征

正如引言中的讨论,已有研究结果认为印度板块在青藏高原下方存在东西向差异性俯冲.针对印度板块在青藏高原西缘的俯冲,图 9展示了部分最新的地球物理与数值模拟研究结果,其中图 9a(Chen et al., 2017a)、图 9b(Zhou and Murphy, 2005)、图 9d(Razi et al., 2016)为地震体波成像结果,图 9cHuangfu等(2018)数值模拟研究中用到的综合性模型(来源于Searle et al., 2011; Yin and Harrison, 2000; Zhao et al., 2010; Li et al., 2008),图 9eXu等(2017)综合远震P波和S波接收函数成像结果推测的俯冲模型.

图 9 不同研究成果推测的印度板块在青藏高原西缘可能的俯冲模式 (a) Chen等(2017b)根据体波伴随成像获得的图 1中剖面D下方的速度结构异常;(b) Zhou和Murphy(2005)利用体波走时层析成像结果获得的图 1中剖面F下方的印度板块俯冲模式;(c) Huangfu等(2018)数值模拟研究中用到的青藏高原西缘印度板块俯冲模型(图 1中剖面E);(d) Razi等(2016)基于Y2台阵和体波走时层析成像结果推测的印度板块在青藏高原西缘的演化俯冲模式;其中,eclogite表示印度下地壳发生榴辉岩化,underthrust Indian mantle表示前期俯冲的印度岩石圈板片;delaminated表示前期俯冲板片发生拆沉,Indian mantle表示现今俯冲的印度岩石圈板片;(e) Xu等(2017)基于接收函数研究推测的俯冲模式.其中红色粗线为获得的主要界面,从上往下依次为低速层顶界面(虚线)、印度下地壳榴辉岩化产生的相变顶界面、Moho界面及LAB界面. Fig. 9 Several possible dynamic models of Indian plate northward subducting beneath western Tibet proposed by different studies (a) The abnormal velocity structure beneath profile D in Fig. 1, which is revealed by adjoint tomography by Chen et al. (2017); (b) The dynamic model beneath profile F in Fig. 1 suggested by Zhou and Murphy (2005), which is based on the seismic body wave tomography; (c) The dynamic model beneath profile E in Fig. 1, which is employed in the numerical modelling study of Huangfu et al. (2018); (d) The suggested dynamic model of subducting Indian plate by Razi et al. (2016). The symbol 'eclogite' represents the eclogitized Indian lower crust, 'delaminated' represents that the previous underthrust Indian lithosphere has been delaminated and probably sink into the upper mantle, and 'Indian mantle' represents the current subducting Indian lithospheric mantle. (e) The dynamic model suggested by the receiver function study of Xu et al. (2017). The red bold lines indicate several main velocity discontinuities, respectively representing the top boundary of low velocity zone in the crust, the top interface of Indian eclogitized lower crust, Moho and LAB structures.

由于体波成像方法对结构的分辨率受控于射线路径的交叉覆盖范围,而射线路径越往上,交叉性越小,因而此方法主要用来探测岩石圈地幔结构特征,对地壳分辨率不高.根据对整个青藏高原岩石圈地幔所做的体波成像研究,图 9ab均显示在青藏高原西缘印度岩石圈板片向北俯冲至金沙江缝合带附近,且西北部存在欧亚岩石圈板片向南的俯冲.如果印度大陆板块能够长距离俯冲至整个青藏高原西缘下方,那么此处青藏高原块体下方原来的岩石圈地幔要么呈现弱力学性质或者被逐步剥离,从而不足以抵挡印度岩石圈板片的向北俯冲.另外的问题是,这些拼贴块体的弱岩石圈或被剥离是继承了碰撞前的性质还是与印度板块的俯冲有关?Chen等(2017b)指出沿着印度—欧亚板块碰撞边界,地壳成分以及岩石圈强度存在重大差异,其数值模拟结果显示弱的亚洲地壳有利于应变向北传递,最终在碰撞带后方形成一个宽阔的造山高原.Huangfu等(2018)将青藏高原划分为岩石圈强弱不同的多个块体进行数值模拟研究,结果显示在青藏高原西缘,从印度大陆与欧亚大陆板块碰撞开始,两者之间的汇聚首先引发高原内部弱强度块体的变形,之后其地壳被逐渐剥离拼贴到临近强度大的块体之上,而其下覆岩石圈地幔则相应地剥离沉入软流圈地幔;然后,强度大的块体发生缩短变形,其下覆岩石圈地幔因重力增加发生剥离或拆沉;最后造成高原内部强弱不同的块体均发生岩石圈拆沉,使得印度板块能够长距离向北俯冲至整个青藏高原西缘之下(图 9c).

与之不同的是,图 9dRazi等(2016)基于Y2观测台阵(图 12中蓝色菱形所代表的宽频带地震台站)和体波走时层析成像方法获得的速度结构,推测的印度板块向北俯冲动力学演化模式.他们的速度结构显示两个高速异常区:其中一个异常带近似于垂直,且止于雅鲁藏布江缝合带以南;另一个高速异常呈长条状,在拉萨块体下方向北东方向倾斜,向下延伸至300 km.如图 9d所示,他们认为在印度与欧亚大陆板块初始碰撞以来,随着印度板块的向北俯冲,其下地壳发生榴辉岩化,使得介质密度增大,增加的重力作用诱发岩石圈拆沉,拆离下来的前期俯冲的印度岩石圈板片则对应雅鲁藏布江和班公湖—怒江缝合带之间的高速异常条带.而现今俯冲的印度岩石圈板片,仅位于特提斯—喜马拉雅块体下方,表现为雅鲁藏布江缝合带以南近乎垂直的高速异常区.然而,图 9eXu等(2017)基于Y2台阵和羚羊计划西线部分台站获得的接收函数结果推测的印度板块俯冲模式,该俯冲模式中最底端红色粗线为S波接收函数探测的LAB界面,其在雅鲁藏布江缝合带两侧存在加深但没有断裂且向北延伸至拉萨块体下方,与d中展现的俯冲模式显著不同.

如果存在类似于Huangfu等(2018)Razi等(2016)描述的拆离印度岩石圈板片,其下沉深度如何?是否到达地幔过渡带附近?如果答案是肯定的,则会对地幔过渡带的温度、物质组分及状态产生扰动,可以被地震学研究手段所探测.下一步计划开展有关该区域下方地幔过渡带的结构研究,或许能够为图 9中不同俯冲模式的甄别提供更多参考信息.

5 结论

在青藏高原西缘,本文研究结果及其他部分已发表的接收函数偏移成像结果,均支持印度大陆板块下地壳在俯冲过程中发生了榴辉岩化.根据大地电磁数据反演的电阻率成像结果,推测的印度下地壳俯冲前沿与接收函数方法探测的下地壳榴辉岩化部分的前沿位置具有很好吻合性.本文综合分析认为,印度大陆下地壳的俯冲前缘在高原西部到达班公湖—怒江缝合带附近,在高原中部和东部递减至拉萨块体中南部.印度岩石圈板片向前俯冲的前沿位置及俯冲特征在青藏高原西缘目前尚不统一,有待进一步研究.

致谢  感谢中国科学院地质与地球物理研究所地震实验室及参加野外地震数据采集工作的所有人员.感谢审稿人和中国科学院地质与地球物理研究所苗来成研究员的宝贵建议.感谢IRIS(http://ds.iris.edu/ds/)提供的Y2台阵观测数据.感谢PROGRAMS 330软件包(http://www.eas.slu.edu/eqc/eqc_cps/getzip.html).感谢朱露培老师提供的接收函数叠加成像代码(http://www.eas.slu.edu/People/LZhu/home.html).
References
Beaumont C, Jamieson R A, Nguyen M H, et al. 2001. Himalayan tectonics explained by extrusion of a low-viscosity crustal channel coupled to focused surface denudation. Nature, 414(6865): 738-742. DOI:10.1038/414738a
Bijwaard H, Spakman W. 2000. Non-linear global P-wave tomography by iterated linearized inversion. Geophysical Journal International, 141(1): 71-82. DOI:10.1046/j.1365-246X.2000.00053.x
Bilham R, Larson K, Freymueller J. 1997. GPS measurements of present-day convergence across the Nepal Himalaya. Nature, 386(6620): 61-64. DOI:10.1038/386061a0
Caldwell W B, Klemperer S L, Lawrence J F, et al. 2013. Characterizing the Main Himalayan Thrust in the Garhwal Himalaya, India with receiver function CCP stacking. Earth and Planetary Science Letters, 367: 15-27. DOI:10.1016/j.epsl.2013.02.009
Chan G H N, Waters D J, Searle M P, et al. 2009. Probing the basement of southern Tibet:evidence from crustal xenoliths entrained in a Miocene ultrapotassic dyke. Journal of the Geological Society, 166(1): 45-52. DOI:10.1144/0016-76492007-145
Chevalier M L. 2019. Active Tectonics along the Karakorum Fault, Western Tibetan Plateau:A Review. Acta Geoscientica Sinica, 40(1): 37-54.
Chen J L, Zheng D C, Long F M. 2019. Joint inversion of receiver function and ambient noise:research and application. Progress in Geophysics (in Chinese), 34(3): 862-869. DOI:10.6038/pg2019CC0218
Chen L, Ai Y S. 2009. Discontinuity structure of the mantle transition zone beneath the North China Craton from receiver function migration. Journal of Geophysical Research:Solid Earth, 114(B6): B06307. DOI:10.1029/2008JB006221
Chen L, Capitanio F A, Liu L J, et al. 2017a. Crustal rheology controls on the Tibetan plateau formation during India-Asia convergence. Nature Communications, 8: 15992. DOI:10.1038/ncomms15992
Chen M, Niu F L, Tromp J, et al. 2017b. Lithospheric foundering and underthrusting imaged beneath Tibet. Nature Communications, 8: 15659. DOI:10.1038/ncomms15659
Chen Y, Li W, Yuan X H, et al. 2015. Tearing of the Indian lithospheric slab beneath southern Tibet revealed by SKS-wave splitting measurements. Earth and Planetary Science Letters, 413: 13-24. DOI:10.1016/j.epsl.2014.12.041
Gan W J, Zhang P Z, Shen Z K, et al. 2007. Present-day crustal motion within the Tibetan Plateau inferred from GPS measurements. Journal of Geophysical Research:Solid Earth, 112(B8): B08416. DOI:10.1029/2005JB004120
Gao R, Huang D D, Lu D Y, et al. 2000. Deep seismic reflection profile across the juncture zone between the Tarim Basin and the West Kunlun Mountains. Chinese Science Bulletin, 45(24): 2281-2286. DOI:10.1007/BF02886369
Gao R, Lu Z W, Klemperer L, et al. 2016. Crustal-scale duplexing beneath the Yarlung Zangbo suture in the western Himalaya. Nature Geoscience, 9(7): 555-560. DOI:10.1038/ngeo2730
Gilligan A, Priestley K F, Roecker S W, et al. 2015. The crustal structure of the western Himalayas and Tibet. Journal of Geophysical Research:Solid Earth, 120(5): 3946-3964. DOI:10.1002/2015JB011891
Hacker B R, Gnos E, Ratschbacher L, et al. 2000. Hot and dry deep crustal xenoliths from Tibet. Science, 287(5462): 2463-2466. DOI:10.1126/science.287.5462.2463
Hetényi G, Cattin R, Brunet F, et al. 2007. Density distribution of the India plate beneath the Tibetan plateau:Geophysical and petrological constraints on the kinetics of lower-crustal eclogitization. Earth and Planetary Science Letters, 264(1-2): 226-244. DOI:10.1016/j.epsl.2007.09.036
Hu J F, Zhu X G, Xia J Y, et al. 2005. Using surface wave and receiver function to jointly inverse the crust-mantle velocity structure in the West Yunnan area. Chinese Journal of Geophysics (in Chinese), 48(5): 1069-1076. DOI:10.1002/cjg2.750
Huangfu P P, Li Z H, Gerya T, et al. 2018. Multi-terrane structure controls the contrasting lithospheric evolution beneath the western and central-eastern Tibetan plateau. Nature Communication, 9: 3780. DOI:10.1038/s41467-018-06233-x
Jain A K, Singh S. 2008. Tectonics of the southern Asian Plate margin along the Karakoram Shear Zone:Constraints from field observations and U-Pb SHRIMP ages. Tectonophysics, 451(1-4): 186-205. DOI:10.1016/j.tecto.2007.11.048
Johnson M R. 2000. Shortening budgets and the role of continental subduction during the India-Asia collision. Earth-Science Reviews, 59(1-4): 101-123.
Kao H, Gao R, Rau R J, et al. 2001. Seismic image of the Tarim basin and its collision with Tibet. Geology, 29(7): 575-578. DOI:10.1130/0091-7613(2001)029<0575:SIOTTB>2.0.CO;2
Kind R, Yuan X H, Saul J, et al. 2002. Seismic images of crust and upper mantle beneath tibet:evidence for Eurasian plate subduction. Science, 298(5586): 1219-1221.
Klemperer S L, Kennedy B M, Sastry S R, et al. 2013. Mantle fluids in the Karakoram fault:Helium isotope evidence. Earth and Planetary Science Letters, 366: 59-70. DOI:10.1016/j.epsl.2013.01.013
Kosarev G, Kind R, Sobolev S V, et al. 1999. Seismic evidence for a detached Indian lithospheric mantle beneath Tibet. Science, 283(5406): 1306-1309. DOI:10.1126/science.283.5406.1306
Lacassin R, Valli F, Arnaud N, et al. 2004. Large-scale geometry, offset and kinematic evolution of the Karakorum fault, Tibet. Earth and Planetary Science Letters, 219(3-4): 255-269. DOI:10.1016/S0012-821X(04)00006-8
Levin V, Huang G C D, Roecker S. 2013. Crust-mantle coupling at the northern edge of the Tibetan plateau:Evidence from focal mechanisms and observations of seismic anisotropy. Tectonophysics, 584: 221-229. DOI:10.1016/j.tecto.2012.05.013
Li C, Van Der Hilst R D, Meltzer A S, et al. 2008. Subduction of the Indian lithosphere beneath the Tibetan Plateau and Burma. Earth and Planetary Science Letters, 274(1-2): 157-168. DOI:10.1016/j.epsl.2008.07.016
Li J T, Song X D. 2018. Tearing of Indian mantle lithosphere from highresolution seismic images and its implications for lithosphere coupling in southern Tibet. Proceedings of the National Academy of Sciences of the United States of America, 115(33): 8296-8300. DOI:10.1073/pnas.1717258115
Li Q S, Gao R, Lu D Y, et al. 2000. An explosive seismic sounding profile across the transition zone between west Kunlun Mts. and Tarim Bain. Science China (in Chinese), 44(7): 666-672.
Li Y H, Wu Q J, Zhang A H, et al. 2006. The Poisson ratio and crustal structure across the NE Tibetan Plateau determined from receiver functions. Chinese Journal of Geophysics (in Chinese), 49(5): 1359-1368.
Liang X F, Chen Y, Tian X B, et al. 2016. 3D imaging of subducting and fragmenting Indian continental lithosphere beneath southern and central Tibet using body-wave finite-frequency tomography. Earth and Planetary Science Letters, 443: 162-175. DOI:10.1016/j.epsl.2016.03.029
Ligorría J P, Ammon C J. 1999. Iterative deconvolution and receiver-function estimation. Bulletin of the Seismological Society of America, 89(5): 1395-1400.
Liu Q Y, Kind R, Li S C. 1996. Maximal likelihood estimation and nonlinear inversion of the complex receiver function spectrum ratio. Acta Geophysica Sinica (in Chinese), 39(4): 500-511.
Murphy M A, Yin A, Kapp P, et al. 2000. Southward propagation of the Karakoram fault system, southwest Tibet:timing and magnitude of slip. Geology, 28(5): 451-454. DOI:10.1130/0091-7613(2000)28<451:SPOTKF>2.0.CO;2
Nábělek J, Hetényi G, Vergne J, et al. 2009. Underplating in the Himalaya-Tibet Collision Zone Revealed by the Hi-CLIMB Experiment. Science, 325(5946): 1371-1374. DOI:10.1126/science.1167719
Phillips R J, Parrish R R, Searle M P. 2004. Age constraints on ductile deformation and long-term slip rates along the Karakoram fault zone, Ladakh. Earth and Planetary Science Letters, 226(3-4): 305-319. DOI:10.1016/j.epsl.2004.07.037
Razi A S, Roecker S W, Levin V. 2016. The fate of the Indian lithosphere beneath western Tibet:Upper mantle elastic wave speed structure from a joint teleseismic and regional body wave tomographic study. Physics of the Earth and Planetary Interiors, 251: 11-23. DOI:10.1016/j.pepi.2015.12.001
Replumaz A, Tapponnier P. 2003. Reconstruction of the deformed collision zone between India and Asia by backward motion of lithospheric blocks. Journal of Geophysical Research:Solid Earth, 198(B6): 2285. DOI:10.1029/2001JB000661
Robinson A C. 2009. Evidence against Quaternary slip on the northern Karakorum Fault suggests kinematic reorganization at the western end of the Himalayan-Tibetan orogeny. Earth and Planetary Science Letters, 286(1-2): 158-170. DOI:10.1016/j.epsl.2009.06.025
Rolland Y, Mahéo G, Pêcher A, et al. 2009. Syn-kinematic emplacement of the Pangong metamorphic and magmatic complex along the Karakorum fault (N Ladakh). Journal of Asian Earth Sciences, 34(1): 10-25. DOI:10.1016/j.jseaes.2008.03.009
Schulte-Pelkum V, Monsalve G, Sheehan A, et al. 2005. Imaging the Indian subcontinent beneath the Himalaya. Nature, 435(7046): 1222-1225. DOI:10.1038/nature03678
Searle M P, Elliott J R, Phillips R J, et al. 2011. Crustal-lithospheric structure and continental extrusion of Tibet. Journal of the Geological Society, 168(3): 633-672. DOI:10.1144/0016-76492010-139
Searle M P, Parrish R R, Thow A V, et al. 2010. Anatomy, age and evolution of a collisional mountain belt:the Baltoro granite batholith and Karakoram Metamorphic complex, Pakistani Karakoram. Journal of the Geological Society, 167(1): 183-202. DOI:10.1144/0016-76492009-043
Shi D N, Wu Z H, Klemperer S L, et al. 2015. Receiver function imaging of crustal suture, steep subduction, and mantle wedge in the eastern India-Tibet continental collision zone. Earth and Planetary Science Letters, 414: 6-15. DOI:10.1016/j.epsl.2014.12.055
Si S K, Gao R, Tian X B. 2019. East-west differential underthrusting of the Indian lithospheric plate beneath central Tibet revealed by imaging VP/VS. Journal of Geophysical Research:Solid Earth, 124(9): 9714-9730. DOI:10.1029/2018JB017259
Tapponnier P, Peltzer G, Le Dain A Y, et al. 1982. Propagating extrusion tectonics in Asia:New insights from simple experiments with plasticine. Geology, 10(12): 611-161. DOI:10.1130/0091-7613(1982)10<611:PETIAN>2.0.CO;2
Tilmann F, Ni J, INDEPTH Ⅲ Seismic Team. 2003. Seismic imaging of the downwelling indian lithosphere beneath central Tibet. Science, 300(5624): 1424-1427. DOI:10.1126/science.1082777
Valli F, Leloup P H, Paquette J L, et al. 2008. New U-Th/Pb constraints on timing of shearing and long-term slip-rate on the Karakorum fault. Tectonics, 27(5): TC5007. DOI:10.1029/2007TC002184
Wang X C, Ding Z F, Zhu L P. 2015. Receiver function analysis of the common conversion point stacking method. Progress in Geophysics (in Chinese), 30(6): 2814-2819. DOI:10.6038/pg20150647
Wittlinger G, Vergne J, Tapponnier P, et al. 2004. Teleseismic imaging of subducting lithosphere and Moho offsets beneath western Tibet. Earth and Planetary Science Letters, 221(1-4): 117-130. DOI:10.1016/S0012-821X(03)00723-4
Wu C L, Tian X B, Xu T, et al. 2019a. Deformation of crust and upper mantle in central Tibet caused by the northward subduction and slab tearing of the Indian lithosphere:New evidence based on shear wave splitting measurements. Earth and Planetary Science Letters, 514: 75-83. DOI:10.1016/j.epsl.2019.02.037
Wu C L, Tian X B, Xu T, et al. 2019b. Upper-crustal anisotropy of the conjugate strike-slip fault zone in central Tibet analyzed using local earthquakes and shear-wave splitting. Bulletin of the Seismological Society of America, 109(5): 1968-1984. DOI:10.1785/0120180333
Wu Q J, Zeng R S. 1998. The crustal structure of Qinghai-Xizang plateau inferred from broadband teleseismic waveform. Acta Geophysica Sinica (in Chinese), 41(5): 669-679.
Wu Z B, Xu T, Wu C L, et al. 2016. Crustal shear-wave velocity structure beneath the western Tibetan plateau revealed by receiver function inversions. Chinese Journal of Geophysics (in Chinese), 59(2): 516-527. DOI:10.6038/cjg20160211
Xie C L, Jin S, Wei W B, et al. 2017. Varying Indian crustal front in the southern Tibetan Plateau as revealed by magnetotelluric data. Earth, Planets and Space, 69: 147. DOI:10.1186/s40623-017-0734-z
Xu Q, Zhao J M, Yuan X H, et al. 2015. Mapping crustal structure beneath southern Tibet:Seismic evidence for continental crustal underthrusting. Gondwana Research, 27(4): 1487-1493. DOI:10.1016/j.gr.2014.01.006
Xu Q, Zhao J M, Yuan X H, et al. 2017. Detailed configuration of the underthrusting Indian lithosphere beneath western Tibet revealed by receiver function images. Journal of Geophysical Research:Solid Earth, 122(10): 8257-8269. DOI:10.1002/2017JB014490
Yang C M, Qian Y Y, Wei Z G. 2019. Study on the effects of core-mantle boundary reflected P-wave on receiver functions. Progress in Geophysics (in Chinese), 34(3): 974-985. DOI:10.6038/pg2019CC0102
Yin A, Harrison T M. 2000. Geologic evolution of the Himalayan-Tibetan Orogen. Annual Review of Earth and Planetary Sciences, 28: 211-280. DOI:10.1146/annurev.earth.28.1.211
Yuan X H, Ni J, Kind R, et al. 1997. Lithospheric and upper mantle structure of southern Tibet from a seismological passive source experiment. Journal of Geophysical Research:Solid Earth, 102(B12): 27491-27500. DOI:10.1029/97JB02379
Zhang Z J, Wang Y H, Houseman G A, et al. 2014. The Moho beneath western Tibet:Shear zones and eclogitization in the lower crust. Earth and Planetary Science Letters, 408: 370-377. DOI:10.1016/j.epsl.2014.10.022
Zhao J M, Yuan X H, Liu H B, et al. 2010. The boundary between the Indian and Asian tectonic plates below Tibet. Proceedings of the National Academy of Sciences of the United States of America, 107(25): 11229-11233. DOI:10.1073/pnas.1001921107
Zhou H W, Murphy M A. 2005. Tomographic evidence for wholesale underthrusting of India beneath the entire Tibetan plateau. Journal of Asian Earth Sciences, 25(3): 445-457. DOI:10.1016/j.jseaes.2004.04.007
Zhu L P. 2000. Crustal structure across the San Andreas Fault, southern California from teleseismic converted waves. Earth and Planetary Science Letters, 179(1): 183-190.
Zhu L P, Kanamori H. 2000. Moho depth variation in southern California from teleseismic receiver functions. Journal of Geophysical Research:Solid Earth, 105(B2): 2969-2980. DOI:10.1029/1999JB900322
Chevalier M L. 2019. 青藏高原西部喀喇昆仑断裂活动构造研究进展综述(英文). 地球学报, 40(1): 37-54.
陈俊磊, 郑定昌, 龙飞明. 2019. 基于接收函数与背景噪声联合反演的研究与应用. 地球物理学进展, 34(3): 862-869. DOI:10.6038/pg2019CC0218
胡家富, 朱雄关, 夏静瑜, 等. 2005. 利用面波和接收函数联合反演滇西地区壳幔速度结构. 地球物理学报, 48(5): 1069-1076.
李秋生, 卢德源, 高锐, 等. 2000. 横跨西昆仑-塔里木接触带的爆炸地震探测. 中国科学(D辑), 30(S1): 16-21.
李永华, 吴庆举, 安张辉, 等. 2006. 青藏高原东北缘地壳S波速度结构与泊松比及其意义. 地球物理学报, 49(5): 1359-1368.
刘葵, 赵文津, 江万, 等. 2006. 印度板块的北缘在哪里?. 地质通报, 25(1-2): 43-47.
刘启元, Kind R, 李顺成. 1996. 接收函数复谱比的最大或然性估计及非线性反演. 地球物理学报, 39(4): 500-511.
吕庆田, 姜枚, 许志琴, 等. 1998. 印度板块俯冲仅到特提斯喜马拉雅之下的地震层析证据. 科学通报, 43(12): 1308-1311.
王兴臣, 丁志峰, 朱露培. 2015. 浅析接收函数共转换点叠加成像方法. 地球物理学进展, 30(6): 2814-2819. DOI:10.6038/pg20150647
吴庆举, 曾融生. 1998. 用宽频带远震接收函数研究青藏高原的地壳结构. 地球物理学报, 41(5): 669-679.
武振波, 徐涛, 武澄泷, 等. 2016. 利用接收函数反演青藏高原西部地壳S波速度结构. 地球物理学报, 59(2): 516-527. DOI:10.6038/cjg20160211
杨传茂, 钱韵衣, 危自根. 2019. 核幔边界反射P波对接收函数影响的研究. 地球物理学进展, 34(3): 974-985. DOI:10.6038/pg2019CC0102