第十四讲地震波层析成像

1、第六章 地震层析成像 Seismic wave is currently the only effective tool that can penetrate the entire earth Structural information of the EarthFrom IRISSeismic waves1939: Jeffreys & Bullen First travel-time tables：Jeffreys-Bullen Seismological Tables 1D Earth model地球内部结构PREM一维全球速度模型1. PREM参考地球模型： Prelimina

2、ry Reference Earth Model(Dziewonski & Anderson, 1981)2. IASP91速度模型(Kennett and Engdahl, 1991)3. AK135 速度模型Kennett et al. (1995)；Montagner and Kennett (1996).另外另外 还有：还有： MC35， STW105， TNA/SNA模型等模型等 莫霍诺维齐（前南斯拉夫）莫霍诺维齐（前南斯拉夫）莫霍面（壳幔边界）莫霍面（壳幔边界） 地壳与地幔分界面的发现者 古登堡（美国）古登堡（美国）核幔边界核幔边界 地幔与地核分界面的发现者 莱曼（丹麦，女

3、）莱曼（丹麦，女） （内外核边界）（内外核边界） 外核与内核分界面的发现者Topography地球一地球一维维结构结构是是远远不够的！远远不够的！mantle convectionPlate tectonicsShearer, 2009Travel time table from ak135 modelTravel time picks3-D variations of Earths Structure from Seismic TomographySeismic waves in the Earth3-D wave speedsTraveltime/waveformInverse probl

4、emResearchers at MIT and Harvard, led by Keiti Aki and Adam Dziewonski in late 1970s and 1980s, pioneered the technique of seismic tomography.PREM模型给出了地球的一维结构，而地球内部三维结构需要更精细的刻画。地震层析成像方法是给出地球内部三维结构的最重要的方法。某种意义上说，地震是照亮地球内部的明灯。地震层析成像方法可以给出： 全球地球结构的横向不均匀性； 典型地球动力学过程的三维结构：俯冲带、地幔柱、大洋中脊等； 小尺度的构造（断层等）； 地震分布

5、特征。全球地震层析成像全球地震层析成像地幔柱地幔柱区域地震层析成像区域地震层析成像区域地震层析成像区域地震层析成像大洋中脊大洋中脊区域地震层析成像区域地震层析成像俯冲带俯冲带局部地震层析成像局部地震层析成像Where it all began: Radon transform: (Johan Radon, 1917): integral of function over a straight line segment. p(s,) f (x,y)(xcos ysin s)dxdywhere p is the radon transform of f(x, y), and is a Dirac

6、Delta Function (an infinite spike at 0 with an integral area of 1)p is also called sinogram, and it is a sine wave when f(x, y) is a point value. Radon transformTomo Greek for “tomos” (body), graphy - study or subject地震层析成像的基础地震层析成像的基础Radon 变换变换f (x,y) p(xcos ysin),)d0Back projection of the function

7、 is a way to solve f() from p() (“Inversion”):Shepp-Logan Phantom(human cerebral)InputRadon ProjectedRecovered(output)Different “generations” of X-Ray Computed Tomography (angled beams are used to increase resolution). Moral: good coverage & cross-crossingCunningham & Jurdy, 2000A few of the

8、 early medical tomo setupsParallel beamFan beam, Multi-receiver,Moves in big stepsBroader fan beam, Coupled, moving source receivers, fast movingBroader fan beam, Moving source, fixed receivers, fast moving(1976)18Brain Scanning Cool Fact:According to an earlier report, the best valentines gift to y

9、our love ones is a freshly taken brainogram. The spots of red shows your love, not your words!Present Generation of models: Dense receiver sets, all rotating, great coverage and crossing rays. 1920 世纪60 年代初期，美国科学家Cormack 从数学和实验结果证实了根据X射线的投影可以唯一地确定人体内部结构,从而奠定了医学诊断上图像重建的理论基础,即X 射线CT(X Ray Computer Tom

10、ography)。 地震波层析成像的发展历史在数学方法上出现了本质上与奥地利数学家1917年提出的Rndon逆变换方法相同的褶积投影方法，Chapman首先从理论上证明了地震学中的-P 变换即是Radon 变换(Chapman , 1981)。地震波层析成像首先由Aki等提出，并给出了小尺度（Aki and Lee, 1976）和区域尺度（Aki et al., 1977）远震体波层析成像（Teleseismic body-wave tomography）。Dziewonski等在1977年给出了全球尺度的体波层析成像成果（Dziewonski et al., 1977）。体波层析成像：主流，

11、对于远震距离，其可沿着多种穿过地球核%幔等路径传播面波层析成像：仅能对上地幔及其以上部分提供约束，适合于台站稀疏和缺乏地震的区域，按成像的范围来分主要有三种：全球层析成像，区域层析成像和局部层析成像。全球层析成像全球间断面，全球波速场分布等；区域层析成像研究地幔柱，洋中脊，板块俯冲带等局部层析成像石油、天然气、煤矿、矿产等勘探，地下水流向等1. 体波层析成像2. 面波层析成像3. 噪音层析成像4. 衰减层析成像5. 有限频层析成像按成像的方法可以分为：体波速度：P，SV，SH面波速度：Rayleigh，Love类似于面波成像品质因子，衰减系数非射线理论，考虑地震波频率1. 体波走时层析成像体波

12、走时层析成像2.有限频带走时层析成像有限频带走时层析成像Fat-ray tomography (Husen and Kissling, 2001)Banana-doughnut tomography (Dahlen et al., 2000)这两种方法主要是利用了地震波传播时，影响地震波传播的不是简单的射线，而是形状与香蕉类似（即具有一定直径的弯曲传播路径）的区域，而震源和记录台站是该香蕉的两个端点。这两种方法优点是在依据地震波传播理论的基础上，考虑到了地震波频率与模型介质尺度间的影响，同时也用到震相走时信息。3. 地震波衰减层析成像方法地震波衰减层析成像方法 (Attenuation tom

13、ography, Rietbrock, 2001)该方法是利用了地震波在传播过程中，地震波的固有频率会随着周期变化成幂指数衰减的关系（即传播介质的品质因子）及其与波速的关系，对周围介质进行成像，其优点是能够给出近地表介质的更为精细结构特征。体波具有相对小的振幅和尖脉冲特征（Shearer，1999）finite frequency kernels for travel time perturbations对地震波的传统认识认为，震发生后所记录的有效地震事件波形为有效信号。除此外，其余部分被称为噪音。这样，形成噪音的声源广泛地分布在地球表面的每个角落，如海浪、风暴和公路上的车流等，统被称为随机分

14、布的波场( Yang and Ritzwoller, 2008)在一个特定环境中，只要是记录这些随机分布的波场足够长，这些噪音也是有规律的，就可以很好地认识区域噪音规律，对区域进行结构成像，跟面波成像方法相似。4. 噪音层析成像噪音层析成像（Ambient seismic noise tomography, Shapiro et al., 2005; Sabra et al., 2005; Pollitz and Flecher, 2005）Ambient noise is enriched at short periods.Better constraints on crustal and

15、uppermost mantle structure than information from earthquakes.Particularly useful in aseismic areas; e.g., continental interiors.For temporary deployments - do not have to wait for earthquakes to occur.Measurements are repeatable: rigorous uncertainty estimates.R. Weaver,Science, 2005Processing Steps

16、: Remove instrument response, de-mean, de-trend, bandpass filter, time-domain normalization, spectral whitening Cross-correlation: 1 day at a time. Stack over many days. Waveform selection (SNR) for tomographytime (s)16.3 Month StackStation Y12CStation 109CMike Ritzwoller et al., 2008体波走时层析成像体波走时层析成

17、像1. 块体模型方法（Aki and Lee, 1976）该方法是将整个反演模型用多个均匀的六面体表示，而每个六面体中心的速度来表示该六面体的整体速度。块体=BlockRoecher（1982）提出了模型空间由多个尺度大小不同的六面体来描述，即可变块体方法。该方法在我国首先由刘福田等（1989）运用。后来者在块体模型中引入了不连续界面（Wintlinger et al., 1998,2004;胥颐等， 2000）2. 格点模型方法Thurber（1983）提出了六面体内的速度变化用其八个顶点的速度变化来表示的技术。之后有人提出利用四个顶点组成的四面体代替六面体（Lin and Roecker,

18、 1997）在六面体格点模型的基础上，通过引入Snell定律和射线伪弯曲法，Zhao (1992,1994) 利用相关区域的先验信息资料（如莫霍面、康拉德面以及板块边界等）构建较为接近真实的模型结构。Travel time (or slowness) inversions:01 tdsv20011= vtvdsdsvvvvv慢度的扰动，作为反演的基本变量30Menke, 1989Least-Squares SolutionsSuppose we have a simple set of linear equationsTT() ()E AXdAXd (AX d)2A X = dWe can d

19、efine a simple scalar quantity EMean square error (or total error)Error functionWe want to minimize the total error, to do so, find first derivative of function E and set to 0.So, doEX 0, we should haveT2()0AAXdATAX ATdThis is known as the system of normal equations.X (ATA )1ATdSo this involves the

20、inversion of the term ATA, This matrix is often called the inner-product matrix, or Toeplitz matrix. The solution is called the least-squares solution, while X=A-1 d is not a least squares solution.最小二乘法适合于观测数据个数多于未知数个数。Left multiply by( ATAI ) 1Pre-conditioning for ill-conditioned inverse problem (

21、damping, smoothing, regularization. Purpose: Stabilize, enhance smoothness/simplicity)Lets use the same definitionETAX d2Define: An objective function J whereJ EX2where m is the damping or regularization parameter.JXEX(X2)XT2()2AA XdI XI is identity matrixMinimize the above by JX 0(ATA I)X ATdX (ATA

22、 I)1ATd(Damped Least Squares solution)阻尼最小二乘法，适合于观测数据个数少于未知数个数。阻尼最小二乘法，适合于观测数据个数少于未知数个数。Solve for the inverse problem(1) Standard Least Squares Solution GTG may be singular or ill-conditioned singular value decomposition (SVD)(2) Damped Least Squares Solution minimizeRay-based traveltime tomographyT

23、21()TmG GIG dminimizeL: Laplacian operator Smooth modelSolution: m = (GTG+2LTL)-1GTd3. Smooth modelCombined norm and Laplacian regularization misfit functionCrust correctionData misfitModel roughnessModel norm在Zhao 等( 1992, 1994) 方法基础上, 双差走时定位( double difference location; Waldhauser and Ellsw orth,

24、2000; Wolfe, 2002) 方法被引入到地震波层析成像反演中, 即地震波双差走时层析成像方法( Zhang and Thurber, 2003, 2005)基本概念1. 震相对：如果两个时间到同一个站台有相同震相的走时数据，那么对于这两个事件来说这两个震相就是一个震相对。2. 事件对：两个事件有一定数目的震相对，并且震中之间的距离在一定的范围内。Double difference methodcaljkikobsjkikijkttttdr)()(Difference in observed arrival time for stations i and jDifference in

25、calculated arrival time for stations i and jDouble difference for event k aim to minimize this residual汪锐等，2013参考模型（Reference model）：计算理论走时。地震波走时自动拾取：获得观测走时.Crust correction: using 3-D Crust 2.0 as the reference crust modelCrust 2.0 Input model1-D crust reference model3-D crust reference model检测版试验是为了检查模型的分辨率。Checkerboard resolution testsCheckerboard and synthetic res