地震正演技术在深反射地震剖面探测中的应用

王光文,王海燕，李洪强，李文辉，庞永香

(1.中国地质科学院 地质研究所岩石圈中心，北京 100037;2.自然资源部 深地动力学重点实验室，北京 100037;3.中国地质科学院，北京 100037;4.中国石油集团测井有限公司 天津分公司，天津 300280)

0 引言

反射地震探测技术是当前地球物理勘探中较为重要的、有效的方法之一。在反射地震勘探中通常需要进行地震波正演模拟来指导野外数据采集、约束数据处理、指导和验证解释，尤其是存在某些特殊地质结构及地面干扰的情况下，地震资料的信噪比和分辨率相对较低，约束资料相对较少，使得地震剖面解释和反演多解性问题突出[1-4]。开展地震波正演模拟研究可以减少多解性，获得合理的地质参数，其中以岩石物性参数为基础的正演模拟所合成的地震记录是研究地下地质结构的桥梁[5-8]。基于石油反射地震技术发展起来的深反射地震剖面技术已被国际地学界公认为是探测地壳上地幔结构的高精技术[9-10]。近几十年来，随着国内外深部探测计划的实施，深反射地震剖面探测技术已经非常成熟，但很少利用正演模拟来约束和指导深地震反射数据采集和处理解释，而地震正演模拟技术在地震勘探中发挥着不可替代的重要作用[11-12]。通过射线追踪及波动方程数值模拟，可以追踪地震波在地下的传播路径，并且通过模型合成的记录可以有效地将地震响应和地质构造联系起来，对地震数据处理与解释起到至关重要的指导作用，使反射地震特征具有明确的地球物理意义及地质意义，并且可以增强人们对地下地质结构的认识。

本文以跨越四川盆地的深反射地震剖面为例，开展深反射地震剖面正演模拟研究。首先根据深反射地震剖面的解释结果，参考区域地球物理资料并结合研究区速度结构设置初始速度和密度，建立初始模型，再根据初始模型和深反射地震剖面的采集参数开展射线追踪正演模拟与波动方程正演模拟，分别将模拟的单炮记录与野外单炮记录、自激自收剖面与偏移叠加剖面进行对比分析，逐步迭代调整初始模型参数，使模型更加合理地反映地下构造情况，最终确定合适的地质模型，为验证和指导地质解释提供依据。同时，深度域地质模型可以为深反射地震剖面构造解释提供深度信息，获得相对准确的地层厚度、地壳厚度和莫霍面埋深等，改善以往仅利用地壳平均速度6 km/s来估算地壳厚度和深度不准的问题。

1 研究区与深地震反射剖面概况

四川盆地位于青藏高原东部、是由挤压作用所形成的沉积盆地之一，由于处在三大洋的结合部位(古亚洲洋、特提斯洋、太平洋)，周边的构造比较复杂[10,13-14]。四川盆地卷入了多期次和多边界的构造变形，是开展盆山耦合作用和深部构造关系研究的理想实验室[10,15,16]，前人在此处做了很多研究工作。为揭示四川盆地的深部结构，中国地质科学院地质研究所在Sionprobe-02项目资助下，完成了跨越四川盆地330 km的深反射地震剖面的数据采集与处理工作。王海燕等利用深地震反射剖面技术得到了四川盆地上地幔精细结构，发现了四川盆地下地壳向上地幔ES方向延伸的倾斜断层，并结合地球化学同位素年龄，认为是扬子克拉通地台西北缘新元古代俯冲的遗迹[10]。康义昌结合钻井资料以及地质资料，认为四川盆地最古老的地层为中震旦系，并且在震旦系沉积地层之下有可能存在更加古老的地层[17]。张岳桥等认为四川盆地北缘和中缘发育的褶皱型构造带为中晚侏罗世，是燕山运动时期陆内造山作用所形成的,其动力学机制与东亚大陆周边板块汇聚碰撞密切相关，板块之间相互碰撞产生的远程效应引起了扬子板块周边造山带的再生复活, 使得四川盆地周缘(华蓥山附近)形成弧形褶皱带，从而奠定了四川盆地的基本构造形态[13]。这些研究让我们对四川盆地的基本构造形态有了一定的认识。

为揭示四川盆地深部地壳结构，研究盆山耦合作用和深部构造关系，在Sinoprobe-02-01项目的资助下，完成了跨越四川盆地总长约330 km的深地震反射剖面(图1)的采集与处理解释工作。该段剖面在2008～2010年期间分三段完成了野外数据采集工作，为获得不同深度地质资料，震源采用不同药量的炸药进行激发，具体采集参数如表1。

图1 深地震反射剖面位置Fig.1 Location of deep seismic reflection profile

表1 跨越四川盆地深反射地震剖面数据采集参数

2 深反射地震正演

2.1 初始模型建立

2.1.1 层位模型建立

四川盆地上地壳覆盖层大致可以划分为中生界地层、古生界地层、前寒武系地层[17-19]，康义昌结合以往这个区域研究的地质资料及沿剖面附近的钻井，并利用本区域一条剖面的深反射信息，研究得出四川盆地震旦系以上的沉积盖层主要可分为四层：侏罗—白垩系陆相红层、三叠系海陆相过渡层、二叠系—三叠系海相上构造层、志留系—中震旦统海相下构造层，并且盆地中部沉积地层存在平行不整合，有部分地层缺失[17]。本文主要依据王海燕等根据深反射剖面建立的解释模型[9-10]及康义昌对此区域的构造解释结果建立初始层位模型，将四川盆地上地壳分为3个层位，中下地壳为一个层位，上地壳的3个沉积层自WN向ES逐渐减薄；在四川盆地中部存在俯冲的莫霍面，测线ES方向存在褶皱带。通过前人对四川盆地结构解释的研究，建立初始地质模型(图2)。

图2 初始地质模型Fig.2 Initial geological model

2.1.2 速度模型构建

目前，研究四川盆地速度结构的文章相当多，研究程度也比较高，如：王椿镛等利用宽角反射/折射地震剖面数据，得到四川盆地地壳平均速度为6.45 km/s，Pn波的速度在8.1～8.2 km/s之间，发现四川盆地地壳平均速度和地幔顶部的Pn波速度偏高[20]；嘉世旭等对反映不同构造单元的震相记录，特别是强震区复杂震相信息进行详细分析和模拟追踪计算，得到龙门山中段褶皱造山带及两侧的横向不均匀地壳速度结构(表2)，显示了地壳速度随深度稳定增加的分层结构特征[21]；李英康等利用深地震反射、宽频地震以及大地电磁的探测结果，将四川盆地Moho面以上的地壳划分为5个层位[22]，并经过震相走时反演和SEIS88程序正演计算，得到了四川盆地二维地壳速度模型(表3)。我们根据以往的钻井、区域地质资料及地球物理特征，参考前人的研究成果，划分研究区Moho面以上层位，并根据测线经过的不同构造位置，对每一层设定不同的速度值，建立初始速度模型(表4)。

表2 四川盆地地壳厚度—速度[21]

表3 四川盆地地壳速度结构[22]

表4 初始速度模型层位速度

2.2 正演模拟

2.2.1 选取地震子波

正演模拟合成的地震记录是地震资料和地质解释之间的桥梁，它对地质层位准确标定起了重要的作用，而子波的选取是合成地震记录的关键[23-24]。挑取不同构造带的野外单炮记录进行频谱分析(图3)，发现单炮频率主要在5～30 Hz之间，主频在10 Hz左右。依据原始数据的频谱分析结果，在进行正演模拟时，将地震子波的频带范围以及主频大小设置为与分析结果一致。对于一般脉冲地震震源(如炸药震源)来说，其原始脉冲基本上是接近最小相位的，因此在模型中选取最小相位地震子波——雷克子波(图4)，然后根据建立的初始地质模型，按照模拟的地震道位置，通过模型中设定的速度和密度值，计算反射系数，最后利用雷克子波与反射系数褶积，合成最终的地震记录[25-26]。

图3 原始单炮频谱分析Fig.3 Spectrum analysis of original single-shot

2.2.2 射线追踪

射线追踪法主要是依据地震波在高频近似条件下其能量沿着射线轨迹传播的原理，进行模拟数据的合成，属于几何地震学方法[27-29]。本文依据建立的初始模型，使用实际地震剖面的采集参数(表1)，进行自激自收和共炮点道集射线追踪正演模拟。依据模型，选定所有的反射界面，模拟射线轨迹及计算射线走时，从射线追踪模拟中(图5)可以直观地看出基于几何地震学原理的射线轨迹，了解地震波在地下构造传播的特点及路径。

a—真实时间;b—振幅谱;c—正交时间;d—相位谱a—real time;b—amplitude spectrum;c—quadrature time;d—phase spectrum图4 子波频谱Fig.4 Wavelet spectrum

a—自激自收记录;b—共炮点道集记录a—zero offset seismic record;b—common shot gather record图5 射线追踪路径Fig.5 Ray tracing routing

2.2.3 共炮点道集记录

深反射地震数据在处理过程中(静校正、偏移)会人为设定一些速度值，如果先将自激自收合成的记录与叠加剖面拟合对比，那么拟合出来的结果可能会存在一定误差，而直接对比接收到的原始单炮记录相对准确。我们先拟合原始单炮记录确定层位的速度值以及模型中层位的构造关系。依据自激自收迭代修正建立的地质模型及野外采集参数，分别在不同构造位置(测线100 km、200 km和280 km处)依据射线追踪合成模拟共炮点道集记录，并将其与原始的单炮记录进行反射波到时对比分析(图6)，依据对比结果，不断迭代修正速度，使单炮中的反射波到时拟合较好为止。

在测线100 km处合成的共炮点道集记录中(图6a1、a2)，发现模拟单炮数据反射波(黑线)的到时偏大，调整速度大小以后，前3层(0～5.5 s)的反射波与实际单炮记录拟合较好。此位置沉积层厚度较深，平均在14 km左右，较厚的沉积层使得反射波到时较大，Moho面反射出现在14 s左右。由于野外震源激发地震波向下传播时吸收衰减严重，Moho面的反射较弱，在原始单炮记录上很难发现，但通过数据模拟，可以得到Moho面在时间剖面上的反射同相轴。

在测线200 km处合成的共炮点道集记录中，模拟数据第一层的反射波与原始数据的直达波叠加在一起，第二、三层的反射波到时比实际单炮数据偏小。调整模型中的速度结构，不断进行迭代拟合，最终使得反射波到时一致(图6b1、b2)。此位置还存在一些弱反射层，并且反射层的到时基本在5 s之内，明显比测线100 km处要小，这是由于此位置沉积层相对较薄所致。莫霍面反射在接近15 s的位置，比测线100 km处的到时偏大，可能是由于地壳到地幔存在倾斜反射层所致。

在测线280 km位置，由于位于华蓥山褶皱带附近，地表起伏较大，地形相对复杂，射线追踪无法到达某些层位，使得合成的单炮记录中，反射波同相轴不能被连续追踪到(图6c1、c2)，并且Moho面位置反射不明显。由于射线追踪方法存在自身的弊端，在华蓥山褶皱带，利用波动方程方法进行正演模拟，此方法基本不会受到地形形态的影响，合成的记录中可以追踪到不同层位深度的反射层[30-33]。从合成的单炮记录中可以发现(图6d1、d2)，此位置沉积层较薄，并且由于层间距较小，使得单炮记录中的浅层反射波到时相差无几，层位拟合不够好，因此需要调试速度大小，使层位拟合效果最佳。通过射线追踪以及波动方程单炮正演模拟，确定测线不同构造段位置的速度值，使拟合的速度接近地下真实速度值。

2.2.4 自激自收地震记录

通过射线追踪及波动方程合成的炮集记录与原始单炮对比分析，对每层的速度进行了初步迭代修正，得到较为合适的速度模型，但是由于原始单炮记录中很难找到Moho面反射层，使得Moho面上下层位的速度值存在不准确性。因此，我们依据共炮点集修正以后的模型，再次利用野外采集参数(表1)，采取逐层迭代修正的方法，进行自激自收合成地震记录与偏移叠加剖面到时信息对比分析，并反复修正模型中各层位的深度、速度等参数(表5)，直至合成的地震记录中层位的到时与叠加剖面中的到时可以基本吻合为止。从偏移叠加剖面(图7)可以看出，浅层反射波同相轴主要分为3层，并且在SE方向浅层褶皱带同相轴明显；Moho面反射清晰，在140～220 km处存在俯冲Moho，俯冲同相轴可以连续追踪。逐层迭代拟合时，主要对比拟合了3个浅层反射同相轴、褶皱带的反射同相轴、Moho面反射同相轴、俯冲Moho反射同相轴位置。

底图—原始单炮记录；覆盖图—模拟单炮记录；a1、a2—100km处射线追踪模拟单炮记录修正前、后对比；b1、b2—200km处射线追踪模拟单炮记录修正前、后对比；c1、c2—280km处射线追踪模拟单炮记录修正前、后对比；d1、d2—280km处波动方程正演模拟单炮记录修正前、后对比 base map—original single shot record；overlay map—simulate single shot record；a1、a2—100km comparison of ray tracing forward simulation before and after single shot correction；b1、b2—200km comparison of ray tracing forward simulation before and after single shot correction；c1、c2—280km comparison of ray tracing forward simulation before and after single shot correction；d1、d2—280km comparison of wave equation forward modeling before and after single shot correction图6 共炮点道集记录对比Fig.6 Comparison of common shot point gather records

图7 偏移叠加剖面(据文献[10]修编)Fig.7 Migration stack profile(revised according to reference [10])

表5 正演模拟层位迭代速度

通过拟合对比，第一层(图8a1、a2)自激自收合成的地震记录与叠加剖面在测线NW方向华蓥山褶皱带位置有一定的差别，模拟数据第一层的到时(双程走时)比叠加剖面偏大，并且测线中部的位置，模拟数据反射波到时偏小。反复修正速度值及反射层位深度关系，使第一层模拟数据与深反射叠加剖面拟合较好为止。依据第一层修正的方法，依次将模型第二、第三层、第四层进行迭代修正(图8b1、b2、c1、c2、d1、d2)，直至合成的地震记录和叠加剖面层位的到时拟合好为止。Moho面(图8d1、d2)经第一次速度迭代，发现到时在俯冲带的位置拟合得不好，俯冲带的同相轴与其东侧Moho面的夹角偏小，并且东部Moho面的到时偏小。调整模型中俯冲角度，减小模型中东部上地幔的速度，经过多次迭代修正，使模拟数据和叠加剖面在俯冲的位置迭代拟合好。通过合成好的自激自收叠加剖面可以看出，俯冲带大约在测线180 km处向下俯冲，到19 s左右时消失，俯冲的到时大约在14～20 s之间。经过逐层速度迭代，得到最终的地质模型(图9)，模型中的速度值也更加接近地下真实速度。

图9 四川盆地最终地质模型Fig.9 Final geological model of Sichuan Basin

3 讨论

本文以深地震反射剖面构造解释结果为基础，参考区域地球物理资料并结合研究区的速度结构设置初始速度和密度，建立初始模型；然后根据初始模型开展射线追踪正演模拟及波动方程模拟，通过对模拟的单炮记录与野外单炮记录、正演的自激自收剖面与叠加剖面进行对比分析，逐步迭代反复调整物性参数，使得在叠加剖面和单炮记录上的拟合结果均达到最佳，从而确定了最终的、相对合适的地质模型。根据该地质模型，验证了深反射地震剖面的构造解释的合理性，并为地质解释提供了约束。

深度域地质模型为深反射地震剖面构造解释提供了深度信息：在龙门山前，沉积层最厚，为14 km，向东南逐渐减薄，在华蓥山前达到最薄约6 km，华蓥山褶皱带沉积层的厚度为5～6 km；莫霍面埋深约为42～44 km，盆地中部(120～200 km范围内)存在ES倾向的反射波组，倾角约30°，认为是扬子板块西北缘新元古代古俯冲的遗迹，俯冲深度达到约60 km。

经过反复修正得到的地质模型和地震剖面构造解释的结果存在差异，地质模型显示古俯冲的角度比实际地震剖面上显示的角度偏大，这可能是由于数据处理过程中偏移处理时所用的偏移速度偏小，导致成像深度变浅，使得俯冲的角度比实际的偏小。正演模型的提出对数据处理过程中参数的设定以及解释模型的合理性具有指导意义。

4 结论

1)射线追踪和波动方程理论建立的正演模型能很好地将地震模型与地质模型联系起来，使地震响应既具有明确的地质意义，又具有地球物理意义，可为深反射地震剖面结果的解释提供一定的指导及约束。

2)根据正演模拟结果，揭示出扬子板块西北缘新元古代俯冲的角度在30°左右，比深反射地震数据的解释模型俯冲角度要大，这可能是在深反射处理过程中，做偏移处理时，给定的偏移速度过小导致同相轴上弯，使得偏移叠加剖面倾斜界面反射波同相轴比实际界面角度偏小。通过正演模拟可以对数据处理过程中参数设定的合理性进行约束。

3)深度域地质模型为深地震反射剖面解释提供了相对准确的深度信息，比利用地壳平均速度 6 km/s 估算的厚度和埋深更合理。

致谢 :中国地质科学院地质所李秋生研究员、卢占武研究员在文章写作构架中给予了诸多启发，中国地质科学院地质所岩石圈中心寒如冰、吴庆宇博士在文章的写作过程中给予了很多专业技术上的支持，并且在模型建立中给出了宝贵建议，在此表示由衷的感谢。