华北克拉通是中朝地台在中国的部分，由中部造山带分为东西部两个块体^[1].与其他克拉通相同，华北克拉通从其在1.8G 年前形成到侏罗纪时期与南部的扬子板块相互碰撞这一段地质历史时期内，一直处于一个相对稳定的状态^[2].之后在扬子板块､太平洋板块､印度洋板块的作用下，华北克拉通经历了重新改造，岩石圈遭到了破坏，厚度发生了很大的变化，这也使华北克拉通成为地质学家和地球物理学家关注的重点区域^[3-6].华北克拉通的西部陆块较为稳定.中部造山带由山西断陷带和太行山隆起组成，存在大量北东和北北东向断裂^[7].东部陆块被改造的最为强烈，经过中生代和新生代的改造，尤其是进入新生代之后，在太平洋板块向欧亚板块俯冲以及印度板块与欧亚板块碰撞所产生的右旋拉张应力场的作用下^[8]，东部陆块渤海湾地区形成了大规模的裂谷带并最终形成了现今覆盖大量沉积层的渤海湾盆地.

接收函数由于其出色的纵向分辨率被广泛地应用于地壳和上地幔的结构研究之中^[9-11]，前人也有很多利用接收函数研究华北克拉通地下结构的成果^[3,4,12,13].由于华北东部陆块中的渤海湾盆地被大量的沉积层所覆盖，Ps转换波震相被沉积层基底反射产生的大振幅的反射波所覆盖，在计算盆地地区的地壳厚度时显得并不理想.Zheng等^[3]利用接收函数计算了华北克拉通中部和东部的地壳结构的变化，通过接收函数波形拟合的方法来反演的地壳结构.考虑到这种算法在拟合过程中不稳定性较强，本文提出利用扩大共转换点单元的共转换点(CCP)叠加方法来减弱浅层反射信号的干扰，得到盆地地区的地壳厚度.这种方法直接､简便，是在Zhu 等(2006)^[14]中所提出的方法基础上进行的改进.2006年10月-2009年9 月3 年间，中国地震局地球物理研究所在华北地区布设了200个流动地震台站(其中甚宽带台站10个，地震计为CMT-3T，频带范围:50 Hz-120s； 宽频带台站190 个，地震计为CMG-3ESP，频带范围:50 Hz-60s)，覆盖了华北克拉通中部陆块和东西陆块的部分地区.本文选取其中的一条北西向测线(L 测线)作为研究对象.该测线横跨了华北克拉通的三个陆块，台站分布较为密集，台站间距约为15km, 包含宽频带台站35个，有利于通过叠加方法增强地下速度间断面的反射信号.本文利用该测线所收集的地震数据，通过接收函数的H-κ 叠加和CCP叠加方法，得到该测线下方的地壳结构，并与已有的地质和地球物理的结果相对比，探讨其深部构造含义.

使用华北地震台阵中的L 测线宽频带台站记录的远震数据，选取从2006年10月至2009年9月的3年间，震中距在30°-90°､震级大于5.5 级､波形记录性噪比高､方位分布均匀的地震事件.台站分布如图 1 所示，该测线为北西-南东走向，其西北端位于黄河东岸的清水河，向东南穿过太行山进入华北平原的渤海湾盆地，其东南端在德州附近.测线上的台站间距为15km 左右，共35个台站，测线长度约为550km.选用的远震震中分布如图 2所示，共有远震事件895个.

图 1

Fig. 1

图 1 台站分布图(根据[1,3,15,16]修改) JZ:冀中坳陷，CX:沧县隆起，HH:黄骅坳陷，CN:埕宁隆起，JY:济阳坳陷，WB:西部陆块，CB:中部陆块，EB:东部陆块. Fig. 1 The distribution of the station s(modified from [1,3,15,16]) JZ:Jizhong Depression, CX:Cangxian Uplift, HH:Huanghua Depression, CN:Chengning Uplift, JY:Jiyang Depression, WB:Western Block, CB:Central Block； EB:Eastern Block.

图 2

Fig. 2

图 2 震中分布图 Fig. 2 Distribution of the epicenters

将每个台站的地震记录截取直达P波的前50s和后150s, 对所截取的部分进行接收函数的计算，并将计算得到的接收函数以30km深度处的入射点间距不小于0.2km 的原则进行组合叠加.

远震P波穿过地壳到达台站时，会有一部分转换为SV 波在地壳中传播.将台站接收到的地震波转换为水平径向､切向和垂直向后，将水平径向和切向的记录分别与垂直向的地震信号进行反褶积，得到的时间序列即为接收函数.由于所选择的30°-90°震中距范围的远震几乎都是近垂直入射的，所以径向方向的转换波振幅很大，利用径向和垂直向的地震波的频率域响应得到的接收函数可以表示为:

(1)

其中的是高斯低通滤波，c为水准量，用来稳定反褶积运算，σ₀² 为Z(ω)的自相关，来规范化c，从而对于不同震级的地震，可以使c缩小其取值范围；(1+c)用来恢复引用水准量带来的振幅减弱^[17].

反褶积运算也可以通过计算下列公式在时间域中进行^[18]:

即计算

(2)

其中σ_i和c_i分别为Z(t)的自相关值和Z(t)与R(t)的互相关值.(2)式中最左侧的矩阵为托普利兹矩阵.与反褶积引入的水准量功能相同，在此，在(2)式中引入了惩罚因子λ 使反演更为稳定.此时，(2)式变为

(3)

本文采用时间域的计算接收函数的方法，Zhu等^[18]认为此方法的分辨率更高一些.

CCP 叠加方法是由Zhu 根据共反射点叠加(CMP)方法提出来的^[19].首先在背景速度下对接收函数的射线路径进行计算，并将时间域的接收函数的每一个振幅都看作某个深度的界面产生的Ps转换波，在对接收函数进行时深转换和入射角矫正之后，这些振幅就会在其发生的深度上排列.对台站下方进行深度划分，本文是以0.5km 对深度进行划分的，之后在每层内设定共转换点单元和像素点的大小.在对某层某一转换点进行共转换点叠加时，将该转换点周围共转换点单元内的所有转换点对应的振幅进行叠加，之后赋值给该转换点周围像素点大小的范围.通过这样的叠加方法即可对台站下方的结构进行成像，生成的图像即为共转换点叠加的结果.该方法可以很直观地得到一条测线下方的地壳结构的变化.本文中所设定共转换点单元是由菲涅耳带半径加平滑长度所决定的.菲涅耳带半径在深度为z的地方大小为$\sqrt {\lambda z/2} $，λ为波长^[14].本文在计算CCP叠加的时候，针对盆地地区沉积层信号对Moho面转换波信号干扰较强的情况，通过调整平滑长度，增加共转换点单元中参与叠加的射线数目，来增强盆地地区Moho面转换波的信号强度.

本文研究区域横向构造较为复杂，既包括山西断陷盆地和渤海湾盆地，又包括太行山隆起和吕梁山地区.针对这种情况，在确定CCP 叠加的背景模型时，对于架设在基岩上的台站，选取了利用H-κ方法计算所得的地壳结构，而对于架设在沉积层上的台站，则选取了利用相邻算法和波形模拟得到的沉积层的结构，沉积层下方的地壳结构选用的是平均的速度结构模型.

由台站接收到的地震记录计算出的接收函数经过坐标变换，转换为水平径向､切向和垂直向三个方向.接收函数的水平径向主要包括直达P 波，Ps波以及其经过多次反射的PpPs, PsPs和PpSs等震相.利用这些震相和直达P 波之间的到时差，可以通过公式(4)，得到地壳的厚度^[15].(4)式中，t_Ps､t_PpPs和t_PpSs+PsPs分别代表Ps, PpPs和PpSs以及PsPs和直达P 波的到时差，p为射线参数，V_P，V_S 表示P波和S 波的速度.本文设定V_P 的值为6.3km/s,V_S=V_P/κ.可以看出，每一个式子中的一组H和κ都对应一个到时，通过这些到时，可以求出其在接收函数的水平径向分量中对应的振幅.将这些振幅按照(5)式进行叠加，当对一定范围内的H和κ 进行扫描时，得到的最大叠加振幅的地方，即为所求得的H和κ 值，这种方法即为H-κ 叠加.叠加的权重w_i之和为1.在进行叠加的时候，每个震相的权重是不同的.本文对三个震相权重的取值分别为0.6，0.3和0.1.

(4)

(5)

在s(H，κ)取最大值时，对其进行泰勒展开，忽略高阶项，就得到H和κ 的方差^[15]，见公式(6).σs是在叠加过程中估算出的s(H，κ)的方差.

(6)

沉积层结构对于架设在沉积层较厚地区的台站是十分重要的.若不考虑沉积层结构，则利用CCP叠加成像得到的地壳厚度大于实际的地壳厚度.作者通过接收函数波形模拟的方法获得了沉积层的厚度及平均速度(具体方法另文给出).由于CCP 叠加中需要的沉积层结构横向变化不能太剧烈，因此对得到的结果做光滑处理，最后得到了测线台站下方的沉积层结构，如表 1所示.

表 1(Table 1)

表 1 盆地中各台站下方沉积层结构 Table 1 The sedimentary structure under the stations in the basin

表 1 盆地中各台站下方沉积层结构 Table 1 The sedimentary structure under the stations in the basin

利用时间域接收函数的计算方法，对测线上的35个台站所接收到的895个地震进行了计算.并对在地壳30km 处的转换点间距小于0.2km 范围内的接收函数进行了叠加，按方位角分布，得到的结果如图 3a所示.将接收函数按照射线参数排列，并对相同的射线参数的接收函数进行了叠加，得到如图3b所示的结果.位于盆地中的台站，由于有较厚的沉积层的存在，沉积层基底的多次反射波震相能量较强，对莫霍面的转换震相产生干扰，造成Moho面的Ps转换波震相无法识别，如图 4所示的情形.对比计算所得的接收函数(图 3 和图 4)，根据接收函数波形的不同，选取出基岩台站的接收函数和架设在沉积层上台站的接收函数，然后对基岩台站的接收函数进行H-κ 叠加的计算，H的取值范围为20-60km,κ 的取值范围为1.5-2.0.进而得到每个台站下的地壳深度和波速比.同样以L205 台站为例，如图 5所示.对于非盆地的台站，通过这样的计算方式，即可得到台站下方的地壳厚度和波速比.由于有些相邻台站之间的波速比存在较大的差异，不利于CCP叠加的成像，也不符合实际情况，因此将结果进行修正，在κ 的方差范围内，对其进行调整.再利用t_Ps计算得到该台站的地壳厚度.对于L211台，由于其t_Ps明显与周围台站不同，故将其厚度和κ 取周围两台的平均值.通过这样的方法，得到用于CCP叠加的背景结构模型，见表 2.

图 3

Fig. 3

图 3 L205台站计算得的接收函数 (a)按方位角分布；(b)按射线参数分布. Fig. 3 The receiver function of L205 (a) Sorted by azimuth； (b)Sorted by ray-parameter.

图 4

Fig. 4

图 4 位于沉积层上的台站的接收函数(说明同图 3) Fig. 4 The receiver function of the station located on the sedimentary

图 5

Fig. 5

图 5 H-κ 叠加结果 Fig. 5 H-κstacking of the station

表 2(Table 2)

表 2 H-κ 叠加方法得到的基岩台站的地壳结构 Table 2 The crustal structure under the stations located on the rock from theH-κ stacking

表 2 H-κ 叠加方法得到的基岩台站的地壳结构 Table 2 The crustal structure under the stations located on the rock from theH-κ stacking

L211台的接收函数结果如图 6所示.从图中可以看到L211台存在两个可能为Ps转换波的震相，经过H-κ 叠加的结果见图 7，这两个震相对应的t_Ps值分别为4.0s和7.2s.若将其κ 值取为周围台站的平均值1.73，对应的地壳厚度分别为32km和异常两种情况.而对比前人的结果^[20-22]，地壳内的波速比低于1.6基本是不可能，无法得到合理的解释.而对于波速比的高异常，则可能对应于地幔物质的上涌^[23]，故还是高异常的结果更为合理；(2)Moho面不连续，Ps震相不清楚，即图 6 中的两种可能的Ps震相都不是由Moho面转换造成的.对于这种情况，我们对L211 台站的厚度和波速比取了周围台站的平均值作为CCP 叠加的背景模型.然而Moho面的不连续，也可能是地幔物质上涌的通道.所以通过这两种情况的分析，作者认为L211 台站下方地壳可能存在地幔物质上涌的情况.

图 6

Fig. 6

图 6 L211台的接收函数结果 (a)tps 为4.0s时接收函数的三组震相；(b)tps 为7.2s时接收函数的三组震相. Fig. 6 The receiver function results of the station L211 (a)The three phases of the receiver function when thetps is4.0s； (b)The three phases of the receiver function when thetps is 7.2s.

图 7

Fig. 7

图 7 L211接收函数的H-κ 叠加结果 (a)tps 为4.0s时的叠加结果；(b)tps 为7.2s时的叠加结果. Fig. 7 TheH-κ stacking results of the station L211 (a)The stacking results whentps is 4.0s； (b)The stacking results whentps is 7.2s.

利用表 1和表 2的沉积层结构和地壳结构作为CCP叠加的背景速度，对测线所有台站的接收函数进行CCP叠加，得到的结果如图 8所示.从图 8中可以看出，L222 台站以西的地壳结构成像清晰，L222以东进入渤海湾盆地地区，沉积层较厚，而且叠加范围较小，对于标准的1 Hz的接收函数来说，Moho面深度处的菲涅尔带的半径约为5km, 而台站间距为15km, 所以无法使多台站接收函数进行叠加，而同一台站的接收函数叠加无法消除沉积层基底反射产生的较强的信号，无法识别出Moho面的位置.因此对L222 以东盆地内的台站加大了平滑长度，将平滑长度由3km 增加至30km 来增大共反射点单元.通过这种方法，可以使相邻多个台站接收到的接收函数进行叠加，以削弱沉积层基底反射的信号强度，并加强了Ps转换波的信号(见图 9).

图 8

Fig. 8

图 8 以3km 为平滑间距的CCP叠加结果(黑色箭头为重力梯度带所在位置) Fig. 8 CCP resultusing 3km as the smooth length (the black arrow is Gravity Line ament)

从图 9中可以看出，该测线西侧地壳厚度约为45km, 逐渐减薄至东部的32km 左右.L204 台站下方的Moho界面存在一个小的错断.L207-L210之间的盆地对应的地壳厚度略微变浅，L211-L212台站下方Moho面出现不连续.L212-L219之间的地壳厚度较为平缓，约为40-42km.L220-L223之间，地壳厚度迅速减薄至33km 左右.L224 以东进入渤海湾盆地，L227 和L228 台站下方地壳出现上隆，最薄可至29km.与冀中坳陷较深的结晶基底相呼应.L220处为测线与南北向的重力强梯度带的交会处附近，从结果中可以看出，地壳厚度的强烈变化区与重力梯度带的位置吻合得很好.

图 9

Fig. 9

图 9 CCP叠加结果 L222台站以西为以3km 作为平滑长度进行的叠加，以东以30km 作为平滑长度进行的叠加，黑色箭头为南北重力梯度带. Fig. 9 CCP stacking result The west of L222 use 3km as the smooth length while the east of it use 30km； the arrow is the North-South Gravity Lineament.

本文利用华北地震台阵中的一条高密度测线记录的远震事件所获得的接收函数，通过H-κ 叠加方法获得了基岩台阵下方的地壳结构，并结合利用波形模拟算法得到的沉积层结构作为CCP 叠加方法的背景结构模型，通过CCP叠加成像获得了测线下方的地壳结构.成像结果显示华北克拉通西部陆块较厚，平均为42km左右，中部陆块西侧地壳较为平坦，阜平至太行山山前断裂带地壳迅速减薄，前人推测可能与太平洋板块俯冲引起的地幔对流的加强有关^[24].东部陆块地壳厚度大约为33km.这些结果均与已有的研究结果相吻合^[25-29].

由图 9 的CCP 叠加结果可以看出，L201-L203的Moho面深度由44km 逐渐加深到46km左右，L204- L205 之间地壳略有抬升，抬升至43km.Zhao等^[1]根据岩石学､地球化学､结构学和P-T轨迹约束等方法划分了华北克拉通的构造分区.他在文中还推测了早中元古代时期东部和西部陆块整合为超级大陆的过程中，中部造山带由东西部陆块碰撞形成的，并将多伦-大同-华山-离石断层作为西部和中部陆块的分界线，将涿县-石家庄和邢台-安阳断层作为中部和东部陆块的分界线(见图 1).我们推测在L204 和L205 之间可能存在中部陆块和西部陆块的深部分界带，这个与Zhao等^[1]的结果略有差别，而Zhao等^[1]的结果是近地表采样得到的结果，本文得到的结果是深部结构的结果，故不冲突.但是我们的数据量也不够充分，至于这里究竟是不是深部的分界线，还需进一步研究分析.

由于L211台站位于大同盆地的边缘，而根据以往关于山西断陷带的研究，山西断陷带为一系列的断陷和断陷间的断隆组成，且断陷的部分为负均衡异常，隆起的部分为正的均衡异常^[30].层析成像的结果^[31]也显示，L211 台站下方的地壳上地幔处为低速异常区，暗示现在的断陷盆地仍然有上地幔物质上涌.所以L211下方Moho面的间断，可能是地幔热物质上涌过程中造成的地壳底界的不连续，成为热物质上涌的通道.

由表 1 沉积层的厚度结果可知，L227 和L228两个台附近是沉积层厚度最大的地区，与其下方Moho面的隆起相对应.这一现象可能是由渤海湾盆地形成受到的拉张作用力而产生的^[8].L209 和L210台站下方也有此现象，说明渤海湾盆地和山西断陷盆地这两个新生代断陷盆地的形成机制相似.总之，根据这两种接收函数的叠加方法，可以看出，华北克拉通的地壳呈现一个西厚东薄的趋势.Moho面较为平坦，盆地下方的莫霍面有所起伏，与上覆沉积层厚度相呼应.山西盆地下方存在由于热物质上涌而造成的Moho面不连续的现象.西部陆块下方存在Moho面的突变，或许为中晚元古代华北克拉通形成时的西部和中部陆块的边界.