有限频层析成像方法.docx
《有限频层析成像方法.docx》由会员分享，可在线阅读，更多相关《有限频层析成像方法.docx（10页珍藏版）》请在冰豆网上搜索。
有限频层析成像方法
有限频体波层析成像简介
徐小兵赵亮
传统射线走时层析成像基于高频射线理论。
在地球内部的速度异
常变化不大，且异常体的尺寸要远大于地震波的特征波长的情况下，
可以根据费马原理确定地震波沿走时最短的路径，此时可以认为到时
差主要取决于射线路径上的速度结构。
然而，实际的地震波具有有限
频率，当地震波经过尺寸较小或相近（相对地震波波长）的速度异常
体时会发生绕射现象。
随着传播距离的增加，地震波波前面由于非均
匀介质导致的超前或滞后特征能逐渐愈合甚至消失（NoletandDahlen,
2000；Hungetal.,2001）。
地震波的这种有限频率特性表明地震到时
差并不主要取决于射线路径上的速度结构。
针对上述问题，Dahlen
和Hung等人推导了地震到时在有限频地震波传播波场下的灵敏场
函数（或根据其几何形状被简称为“香蕉—甜面圈”理论）（Dahlen
etal.,2000；Hungetal.,2000），进而发展了有限频层析成像方法。
相对传统射线走时层析成像，该方法由于更真实地反映地震波的传播特征，因此能够更好地重建地球深部结构速度图像。
一、有限频层析成像方法原理
地震波的有限频率特性决定了不同频率的地震波的走时对射线
周围三维速度结构具有不同的敏感程度。
Dahlen等（Dahlenetal.,2000）结合体波传播理论和线性Born单次散射近似来模拟通过地球
内部任意一点的散射波到达台站的波形扰动，利用不同频率的体波与
未经散射的直达体波之间的相互干涉结果，得出有限频率走时变化与
中心几何射线附近第一菲涅尔带内的三维速度扰动之间的关系为：
T
Kx,wcx/cxd
3x
（1）
式中，T是经过波形互相关处理得到的散射波与直达波的走时差；
cx/cx为空间点x处的体波速度扰动量；Kx,w是该点三维灵
敏度核，代表了某个散射点x处的速度异常对走时变化的贡献。
从上
式可以看出，敏感场K不仅是空间位置x的函数，同时还是频率w
的函数，即不同频率的地震波对不同空间范围的速度结构进行了采
样。
在射线理论下地震波走时表示为慢度沿射线路径的线积分，而在有限频率灵敏核理论中地震波走时可以更精确地表示为灵敏度核函
数在射线周围一定空间内的体积分，这是有限频率理论与射线理论最本质的差别。
因此，确定在一定速度模型下的有限频率灵敏核函数是进行有限频层析成像的关键所在。
二、有限频率地震体波灵敏度核函数及其性质
获得有限频率地震体波灵敏度核函数是开展有限频层析成像的基础。
图1.有限频体波走时在球对称介质下的几何路径示意图（Hungetal2004）
2
1
3m
sinT'T''T
M'M''M2d
Kc,
1NN'N''
0
c,
2
2
crays'rays''
cr'''
2m
d
0
（2）
公式
（2）表明有限频敏感度核函数是对时间窗内所有单点散射波进
行积分求和的形式。
由于对所有满足条件的单点散射波进行射线追踪
是一个非常耗时的过程，因此在实际成像过程很难实现。
针对此情况，
Dahlen等（2000）在计算敏感场函数时引入旁轴近似，即认为走时残
差只对第一菲涅尔带内的速度扰动比较敏感，从而来有效的减少射线
追踪的运算量。
采用的旁轴近似的公式为：
3Ssyn
2
1
0
sin
Td
K
2
..........（3）
2ccr'''
2
Ssyn
d
0
图2.公式2得出的直达P波灵敏度核函数（Hungetal2000）
图3.根据旁轴近似（公式3）给出的直达P波灵敏度核函数（Hungetal2000）
从图2和图3的对比中我们可以看出通过旁轴近似由公式3得出
的灵敏场可以很好地对由公式2得出的理论敏感场进行近似。
图5.灵敏度核函数与频率的关系（Hungetal.,2004）
由高、中、低三个不同频带内的灵敏度核函数（图5）的形态
可以看出，地震波的频段越高其灵敏度核越窄且强度越大，最接近于
传统的射线；而频带范围越低则灵敏度核越宽，这时的灵敏度核强度
也会降低。
理论分析表明，对宽频带波形分频带处理，不同频带内的
分析结果可以分别约束中心射线路径周围不同范围处的速度结构。
三、有限频理论与射线理论的比较
图6.冰岛下方射线理论与有限频方法的灵敏场函数分布对比。
（a）射线方法的
灵敏场；（b）有限频方法的灵敏场。
图件来自Hungetal.,2004。
图6展示了射线方法和有限频率方法的灵敏场分布：
相对射线方法，有限频率方法对空间的采样更加密集、平滑，且空间范围更大，因此能够充分利用宽频带资料中的信息，减少数据不均匀覆盖的影
响，进而能更加充分对速度异常体进行重建。
如同很多对比实验所示，有限频层析成像能够提高反演结果分辨率。
华北地区的实际观测结果也表明，有限频方法获得的异常体分布尽管大体与射线方法的结果一
致，但重建得到的异常体扰动幅度更大、分辨率更高（图7）。
基于层析成像反演理论，射线走时层析成像和有限频走时层析成
像最终都是求解模型参数与观测数据所构成的方程组，两者的不同之
处在于正演计算模型参数和观测数据的（敏感场泛函）之间的不同。
（a）
（b）
图7.采用射线方法和有限频率方法得到的层析成像结果比较。
（a）采用射线方法得到的150公里深度切片，左图为P波结果，右图为S波结果；（b）有限频方法结果。
四、反演过程
应用网格节点对模型空间参数化后，可将有限频率地震层析成像
的反演问题写成下面的矩阵方程：
dGm
式中，d为观测走时向量；G为灵敏度核函数矩阵，该矩阵第i行第j列的元素Gij是第i个观测走时的灵敏度核函数在三维网格节点
j处一定空间范围内的体积分；m为要求求解的速度扰动向量。
反演问题转化为求解线性方程组Axb的问题，即是求解满足最
小二乘条件minAxb2的解。
A为mn阶矩阵；x为n维矢量，称为模型空间矢量；b为m维矢量，称为数据空间矢量。
一般数据空间矢量的维度m大于模型空间矢量的维度n。
A为大型稀疏矩阵，如果采用传统的求解广义逆的方法来求解，计算机的内存需求很大，计算效率很低。
所以反演过程采用LSQR方法。
LSQR算法运用了三个主要的技术：
（1）Lanczos迭代：
这个迭代方法产生模型空间的正交基矢量v1,k以及数据空间的正交基矢量u1,k。
令数据空间初始矢量为1u1b，则相应的模型空间初始向量为1v1ATu1，利用Lanczos迭代构建正交基，把模型按照
正交基展开，得到k1k阶的下三角矩形阵Bk；如果令模型空间初
始向量为1V1ATb，则相应的数据空间初始向量为1p1
Av1，同样利
用Lanczos迭代构建正交基，把模型按照正交基展开，得到
kk阶上
双对角方阵Rk。
由于所选的正交基是相同的，所以这两个双对角方
阵一定存在特定的关系。
（2）Givens变换
运用一系列非常有规律的Givens平面旋转正交变换实现正交变
换，建立这两个双对角矩阵之间的联系。
由BkRk的过程即是通过
Givens变换对Bk的QR分解。
这个过程实现了最小二乘和QR分解
的完美结合（这也是被称为LSQR的原因）。
（3）模型迭代
通过迭代运算直接计算模型矢量~，逐步逼近最小二乘解。
xk
LSQR算法程序设计的框架流程如下：
1）初始化
1u1
d1
d,1v1ATu11ATu1,
w1
~
0,11,1
1.
v1,x0
2）循环计算fork=1,2,3,...重复（3）-（6）.
3）Lanczos算法，双对角化
k1uk1
Avk
kuk
k1vk1
ATuk1
k1vk
4）正交变换
2
2
k
k
k
1
ck
k
k
sk
k
1
k
k1
sk
k1
k
ck
k
k
1
ckk1
k
1
skk
k1与k1是下一次迭代的起始值。
5）更新~和w
x
~
~
k
kwk
xk
xk1
wk1
vk1
k1
kwk
6）检验是否收敛
如果达到收敛条件，停止迭代。
LSQR算法的优点是，计算机的内存要求低，易于实现并行算法，迭
代收敛快，求解精度高。
缺点是给不出显示的广义逆，也就给不出显
示的模型分辨矩阵和协方差矩阵。
如下是层析成像反演流程图：
开始
选取研究区域
拾取相对走时残差
建立初始速度模型
选用LSQR反演方法调整上一次模
型拟合观测和理论相对走时残差
选相对走时残差是否
小于设定阈值
N
Y
得到最终的速度扰模型
结束
参考文献
1.Dahlen,F.A.,S.H.,Hung,G.Nolet,Frechetkernelsforfinite-frequencytraveltimes-I.Theory[J].GeophysJInt,141:
157-174,2000
2.Hung,S.H.,F.A.,Dahlen,G.Nolet,Frchetékernelsforfinite-frequencytraveltimes-Ⅱ.Example,GeophysJInt,141:
175-203,2000.
3.Hung,S.H.,Y.Shen,L.Y.,Chiao,ImagingseismicvelocitystructurebeneaththeIcelandhotspot:
afinitefrequencyapproach.JGeophysRes,109:
B08305,doi:
10.1029/2003JB002889,2004
4.Zhao,L.,R.M.Allen,T.Zheng,R.Zhu.,High-resolutionbody-wavetomographymodelsoftheuppermantlebeneatheasternChinaandtheadjacentareas,Geochem.Geophys.Geosyst.,13,Q06007,doi:
10.1029/2012GC004119,2012.
5.Zhao,L.,R.M.Allen,T.Y.Zheng,S.Hung（2009）.ReactivationofanArcheancraton:
ConstraintsfromP-andS-wavetomographyinNorthChina,Geophy.Res.Lett.,36,L17306,doi:
10.1029/2009GL039781.
升级会员 