R.Mereu B.Roy S.Winardhi
（Dept.of Earth Sciences,University of Western Ontario,London,Ontario,Canada,N6A 5B7）
摘要 过去几年里，在加拿大东部的加拿大地盾开展了一系列长距离地震折射/宽角反射以及近垂直重合反射法的实验研究工作[2]。地震近垂直反射法是一种很有潜力的方法，它使我们能够进一步了解地壳内的详细结构。上述实验的不足之处仍是信噪比较差的问题。这通常是由于地壳内局部侧向不均匀性引起的散射和射线路径延长效应所造成的。在本项研究中，采用了更有效的图像识别法代替传统的CDP映像法。该方法产生了意想不到的效果，大大提高了对地下各主要反射层的分辨能力。地震折射/宽角反射技术利用了组合爆炸，并采用数百个接收器接收。这一方法的实施使我们得以从一个不同的视野去观察和研究地壳。为了优化地下结构的映像，在对数据进行层析时，又发展了适用于侧向非均匀模型参数化的新方法，即三角单元法与延时最小平方反演技术相结合。为了揭示异常的速度结构，从传统的速度解中减去了区域垂直速度梯度，从而获得了具有明显改善的地下结构映像。
关键词 图像识别 层析 地壳 莫霍面 反射 折射 地震映像
1 引言
过去10年里，在加拿大东部穿越加拿大地盾开展了一系列长距离地震折射/宽角反射以及近垂直重合反射法的实验研究工作。主要的研究区域为中央沉积变质带，GrenvilleFront构造带、Sudbury盆地、Midcontinent断裂系等。在1986年进行的Glimpce实验中，对横跨Great lakes的地震测线采用了空气枪震源技术，源距为60m。在最近一次，即1992年的Lithoprobe Abitibi-Grenville的陆地地震折射实验中，采用了44个炮点爆破，平均炮间距为30km。每次爆破由一组415仪器记录，记录点之间的距离为1.0～1.5km。剖面长度为180～640Km。上述实验构成了一个密集的射线覆盖区，因而使我们能够应用地震层析成像技术更详细地研究我们所感兴趣的构造区的地震波速变化。大多数折射测线都由近垂直反射测线做补充。这些垂直反射测线在陆地上使用了连续震动源，而在湖中使用了空气枪震源技术。上述实验的详细情况及分析方法见Green等[1]、Epili和Mereu[2]、Hamilton和Mereu[3]以及Winardhi和Mereu[4]的文章。
2 分析工具
为了有效地分析处理在上述实验中获得的大量数据，特别是地震折射/反射数据，我们设计并研制了一系列的工具。它们是：①能够在喷墨打印机上绘制大型地震剖面的打印机－绘图仪程序；②一套数字滤波器，包括用于三分量数据的极化滤波器和用于CDP数据的图像识别滤波器；③采用自动及人工方法在屏幕上显示并且识别初至波和续至波的计算机交互程序；④在侧向非均匀模型中发射任何类型的射线（P射线、S射线或转换射线）以及多次波射线的交互型射线追踪程序。
3 利用图像识别技术提高信噪比
在大部分开展地壳近垂直地震反射实验的地区，过去都曾经历过某种构造形变，因而地层不再满足均匀层状的条件，即沿侧向和垂向都是非均匀性的。①地震能量在这类介质内传播经过很长射线路径之后，以稍有不同之时间到达各接收站，因而其走时曲线很不光滑。采用某种叠加技术以求提高信噪比的常规处理方法，在许多情形下因非同相叠加，实际上会破坏信号，以致地下界面成像不能令人非常满意。这种现象在一些炮点道集及CDP道集有明显反射体的情况下也会发生。②本项研究着重于以图像识别方法为基础的信号增强技术，该项技术应用于叠前处理流程，用以改善信噪比。Roy和Mereu[5]，根据信号的侧向连续性、振幅、波形、信号频率等特征来鉴别信号和噪声，利用该技术消除CDP道集的噪声，大大降低了叠加剖面上的噪声干扰水平。在信号叠加的方法上也与过去略有不同，即通过在一个小的移动时窗上计算出具有相同极性的信号能量，这样，就连那些非同相性的信号也都考虑在内了。图1a、图1b给出了一条反射剖面的实例。其中图1a是利用传统处理方法得到的映像，而图1b是利用图像识别法得到的改进了的映像。该图所反映的是位于加拿大地盾东南部中央沉积变质带地震剖面中倾斜剪切带的一部分。
4 层析成像分析
利用层析成像技术对大量的地震折射/宽角反射数据进行了分析，从而得到地壳的速度结构映像。具体步骤如下：①建立一个合理的初始模型；②根据该模型用求正问题的方法构制一组理论走时曲线；③用最优化方法对模型进行细调，直至观测走时曲线与理论走时曲线之间的误差达到最小。初始模型的建立是通过对完整的初至Pg波组进行最小二乘法分析后得到的，于是便产生了一个平均的“标准”地壳模型。在我们的实例中，地表波速为6.13km/s，随深度的增加波速呈线性增大，至地下40Km处波速为7.10km/s。研究结果表明，在走时曲线上出现的很多局部分散点都是由近地表岩性沿侧向及垂向的变化所致。我们利用延时法将走时数据转换为近地表的速度结构，然后再利用平滑的走时曲线对较深层的结构做层析成像。
图1 横跨加拿大地盾中央沉积变质带的剖面实例
a—传统处理方法的CDP叠加剖面；b—增强地震信号的图像识别技术
有很多模型参数化的方法，例如Cerveny和Psencik[6],Zelt和Smith[7]。每种方法都有其各自的优缺点。在本研究中，我们选择了最新版的“三角单元法”（Mereu,[8]），该方法是将模型剖分为速度梯度为常数的一系列三角单元。之所以选择三角形做模型的基本单元，是因为它对模型结构的侧向变化及倾斜断层非常有效。沿侧向和垂向的速度变化可以表示为
V（x,z）=ax+bz+c
GradV=ai+bj
其中：V（x,z）是位于模型中的任意一点的地震波速度，GradV是每个三角单元内为常数的速度梯度。常数a、b、c是由每个三角形各节点上的速度值确定的。对每个三角形而言，这三个常数，其速度梯度值是不同的。速度场的梯度值为常数，这样就确保了在每个三角单元内追踪的所有的射线路径均为圆弧，因而对整个模型能够很快地追踪射线。射线追踪是通过求解具有三角形边界线性方程的射线弧方程，并在每个交点上应用斯奈尔定律来实现的。三角单元剖分技术的优点是使我们处理的数据量相对减少。当地下结构比较简单时，三角单元可以剖分得大一些；而当地下结构比较复杂时，则需要将三角单元剖分得小一点。因此，模型剖分的方式并不是唯一的。多次的试验结果表明，一般而言，走时数据对模型剖分的不唯一性并不很敏感。其它试验则表明，如果根据炮间距来确定三角单元的数目，如图2所示，则可以获得最佳的单元数。根据时间延迟、自由表面效应、几何扩展效应以及射线在传播路径中的复传递及复反射系数计算出传输函数，就可以构制出合成地震图。计算机自动搜寻模型中所有主要的走时曲线，并对其做数值编码，确定其终点，进而在模型中以适当的角度进行射线追踪，以便计算出波至时间和振幅。地壳走时数据的反演是个非线性反演问题。对于任何自动反演方法，建立模型时都容易出现不稳定解或不规则路径的情况。为了克服这个问题，我们将原始程序修改成为交互式程序。
图2 采用三角单元法对速度模型参数化的实例
该模型中含有一个断裂带，如图3a所示
新程序的特点是用户通过在计算机键盘上简单地敲几个键就可以修改模型的三角形单元及其速度值。对反演的监控是通过即时显示射线路径以及理论走时和观测走时曲线来实现的。在反演过程中，采用了分层剥离法或扩展单元法[9]（White）对数据进行模拟，反演由表层开始，先构制出浅部地壳结构映像，然后逐渐深入到深部地壳，直至莫霍面。
对莫霍面深度上的模型构制，主要利用了观测的PmP和Pn波走时数据，在那些观测不到PmP波或PmP走时曲线的支很短的地区，莫霍面表现为一个厚的过渡带。为了改善地下结构的映像，将速度剖面减去一个平均的标准地壳速度梯度，即可得到异常速度剖面。此外，我们还发展了一种反映模型不确定性的技术，它以射线覆盖的相对密度为基础，并与模拟退火分析相结合。
图3a显示了一个速度沿侧向变化的模型实例，一条大断层将模型中的高速区与低速区隔开。图2表示了如何利用三角单元法对这样一个模型参数化。在该例中利用了10个炮点的数据。图3b为速度异常映像图，该图明显地改进了地下结构映像，突出地表现了断层的存在。图3c为图3a模型的速度不确定剖面。由图可见，最精确的速度值是在射线覆盖最密集处获得的，而在地壳深部和模型两侧边角的射线稀疏处，估算出的波速值准确性较差。
图3 层析成像图
a—速度模型：一条断层将模型划分为一个高速区和一个低速区；b—速度异常映像图：从（a）减去一个标准地壳，该地壳在地表速度为6.13km/s，莫霍面附近速度为7.10km/s;c—速度不确定剖面，该不确定值取决于射线覆盖的密度
最近，我们将层析成像法应用于1992年穿越加拿大地盾的Lithoprobe Abitibi-Grenville 4条地震折射剖面分析，取得了较好效果（Winardhi和Mereu,1997）。研究结果中最令人感兴趣的是获得了Grenville Front构造带（Superior和Grenville之间地质单元边界）的构造映像图。该图揭示出一个异常的速度梯度带向东南方向倾伏，穿透地壳，直达莫霍面。
5 结论
现代地震勘探采用了近垂直反射法和宽角反射/折射法。本文介绍了两个非常有效的技术——图像识别法和速度异常层析技术。这些技术的应用改善了地壳的速度结构映像。
致谢 感谢LSPF的K.Vasudevan和R.Maier，以及西安大略大学的B.Dunn和J.Brunet在计算机硬件方面所给予的技术上的帮助，并提供了地震资料处理软件包。本项研究得到Lithoprobe和NSERC的共同资助，资助号为A1793。
（江涛译，许云校）
参考文献
[1] A.G.Green,B.Milkereit,A.Davidson,C.Spencer,D.R.Hutchinson,D.R.Cannon,W.F.Lee,M.W.Agena,J.C.Behrendt,W.J.Hinze.Crustal structure of the Grenville Front and adjacent terranes.Geology,1988,16:688～792.
[2] D.Epili and R.F.Mereu.The Grenville front tectonic zone:Results from the 1986 Great Lakes onshore seismic wide-angle reflection experiment.Joural of Geophysical research,1991,96:16335～16348.
[3] D.Hamilton and R.F.Mereu.2-D tomographic imaging across the North American mid-continent rift system.Geophysical Joural international,1993,112:344～358.
[4] S.Winardhi and R.F.Mereu.Crustal velocity structure of the Superior and Grenville Provinces of the Southeastern Canadian Shield.Canadian Journal of Earth Sciences,In press,1997.
[5] B.Roy and R.F.Mereu.Signal enhancement using pattern recognition techniques with application to near vertical crustal seismic reflection experiments.Geophysical research letters,1996,23,1849～1852.
[6] V.Cerveny and I.Psencik.Gaussian beams in two dimensional laterally varying layered structures.Royal Astronomical Society,Geophysical Journal,1984,78:65～91.
[7] C.A.Zelt and R.B.Smith.Seismic travel-time inversion for 2-D crustal velocity structure.Geophysical Joural international,1992,108,16～34.
[8] R.F.Mereu.An interpretation of the seismic-refraction data recorded along profile SJ-6:Morro Bay-Sierra Nevada,California.In:interpretation of the SJ-6 seismic reflection/refraction profile,south-central California,USA.In:Proceedings of the 1985 CCSS Workshop on interpretation of seismic Wave Propagation in laterally Heterogeneous Terranes,USGS Open File Report 87-73,A.W.Walter and W.D.Mooney（Eds.）,1987,20～37.
[9] D.White.Two dimensional seismic refraction tomography.Geophysical Joural international,1989,97,223～245．