保角变换在区域面波群速度反演中的应用

[文章编号]　0001-5733(2001)01-0064-08

[中图分类号]　P351保角变换在区域面波群速度反演中的应用

朱良保　许　庆　陈晓非

(北京大学地球物理系,北京　100871)

[摘　要]　通过保角变换,把一个球面局部区域扩展到球面上更大的球台域,再由球台域解析延拓到整个球面域.作为约束条件,变换过程中面波速度保持不变.在变换后的球面域上用球谐函数来拟合速度函数,达到降低球谐系数阶数的目的.面波群速度的反演变成了球谐系数的线性化反演,由球谐系数计算反演分辨核以及方差.该方法的优点是计算速度快,等值线光滑,构造界限清晰,该方法同样适用于各种局部球面区域的分析.

[关键词]　球谐函数,保角变换,群速度反演.

1　引 言面波群速度反演有网格法与无网格法两种.网格法是把研究区域网格化,每一网格中的速度是均匀的,由实际频散数据来反演每一网格的速度.无网格法是把研究区域的速度函数用一组函数展开或用一组函数来拟合.无网格方法又分为两种,一种是选择一组先验函数,一般选择一组正交函数;另一种是根据路径的分布,由具体路径来确定一组函数,主要有Tarantola 非线性反演方法[1,2]和Ditmar 2Yanovskaya 方法[3,4].文献[3,4]的方法是Backus 2G ilbert [5,6]方法在面波反演中的推广.本文讨论的是选择先验函数的方法.

一个函数的球谐展开要求该函数定义在整个球面上.局部区域的反演,速度函数只定义在球面的局部区域上.如果把局部球面区域直接延拓到整个球面区域,原区域保持不变,理论上说,在延拓后的球面上作球谐展开拟合原区域上的速度函数是可能的.但是这种做法可能需要展开到很高的球谐阶数.由于实际资料的数目有限,在反演过程中一般是不定解问题,因此,用反演方法做这种拟合实际上是行不通的.为了在反演方法中应用球谐函数来拟合球面局部区域的速度函数,可以把球面局部区域通过保角变换扩展到球面上更大的区域上.在扩展后的区域上作球谐函数的拟合可以极大地降低球谐阶数,使之成为超定的病态线性反演问题.作者在文献[7]中对球面局部区域的群速度反演的球谐函数法作了一些讨论,并作了理论上的数值对比.

[收稿日期]　2000-01-10收到,2000-07-28收到修定稿.

[基金项目]　国家自然科学基金项目(49625406),国家重点基础研究发展规划项目(G 1998040702).

[作者简介]　朱良保,男,生于1956年,博士,1982年毕业于中国科学技术大学地球和空间科学系.现为北京大学地

球物理系副教授,主要从事横向非均匀介质中高频渐进地震波理论、面波层析成像、地壳上地幔结构的

研究.E 2mail :zhulb @https://www.doczj.com/doc/d311710358.html,

第44卷第1期

　2001年1月 地球物理学报CHIN ESE JOURNAL OF GEOPHYSICS

Vol.44,No.1　Jan.,2001

2　球面局部区域到球台域的保角变换

2.1　球极投影 设球面局部区域定义为φ1,φ2;;1,;2(球面极坐标),其中φ2待定.群速度函数为v (φ,;).作球极投影,使之映射到复平面的一个区域φ1,φ2;r 1,r 2(平面极坐标),群速度函数为V (r ,φ),这种投影满足保角变换.作为约束条件,我们假设在以下的任何一步保角变换中,速度始终保持不变.图1(a )是这种投影的示意图

.由球极投影得

r R 0sin (π-;)=2R 0R 0+R 0cos (

π-;),r =2R 0sin ;1-cos ;,<=φ,V (r ,<)=v (φ,;).

(1)

根据球极投影变换的性质:“平面上任一圆周都变成球面上的圆周,反之亦然”,可以推知,球极投影后,该区域在复平面上的形状是部分圆环,如图1(b )所示.

图1

(a )球极投影;(b )球极投影后,球面局部区域变为复平面内的部分圆环;(c )复平面内的旋转变换;

(d )对数变换后,部分圆环变为矩形;(e )矩形的伸缩;(f )指数变换后,矩形变成有割线的圆环.

Fig.1

(a )Stereographic projection ;(b )After stereographic projection a regional area of a globe becomes a fraction

of an annulus in a complex plane ;(c )Rotation transformation ;(d )After logarithm transformation a fraction

of an annulus becomes a square ;(e )Extension transformation to a square ;(f )After exponent transformation

a square becomes an annulus.

2.2　旋转变换

α1=r e i φ?e -i φ1=r e i (φ-φ1

),旋转变换后图1(b )变为图1(c ).561期 朱良保等:保角变换在区域面波群速度反演中的应用

2.3　对数变换

α

2

=ln r e i(φ-φ1)=ln r+i(φ-φ1),对数变换后图1(c)变为图1(d).

2.4　伸缩变换

α

3=c1(ln r+i(φ-φ1))=ln r c

1+i c1(φ-φ1),

c1=

2π

φ

2

-φ1

.

(2)

伸缩变换后图1(d)变为图1(e),由于φ2为待定量,所以c1是待定常量.

2.5　指数变换

α

4

=e(ln r c1+i c1(φ-φ1))=r c1e i c1(φ-φ1).

指数变换后图1(e)变为图1(f),它是有一割线的圆环.

2.6　再伸缩变换

α

5

=R?e iχ,

其中R=c2r c1,χ=c1(φ-φ1),(3) c2为待定的伸缩系数.

2.7　球极投影

把平面圆环变回到球面上.设变换后的球面坐标为<,θ,

令<=χ,(4)

根据球极投影公式(1)得R=2R0sinθ

1-cosθ,

(5)

由(1)-(5)式得2R0sinθ

1-cosθ

=c2

2R0sin;

1-cos;

c

1

,(6) <=c1(-φ1),

由三角公式可把(6)式化简为　tan θ

2

=c3(2R0)1-c1tan

;

2

c

1

,(7)

其中c3=c-12.

c3由下列条件决定,

当;=;1时θ=θ0,(8)

;=;2时θ=π-θ0,(9)亦即使变换后的球台域关于两极对称.

由式(7)-(9)式得c3=(2R0)c1-1tan ;1

2

tan

;2

2

-

c

1

2

,(10)

把(10)式代入(7)式后得　tan θ

2

=tan

;1

2

tan

;2

2

-

c

1

2

tan

;

2

c

1

.(11)

c1值的选择比较自由,只要它满足下列条件

1≤c1≤2π/(φ32-φ31).

一种较好的选择是c1=π/(;2-;1),

其中φ31,φ32,;1,;2是覆盖所有路径的区域边界.

66地　球　物　理　学　报 44卷

一旦确定了c 1,由式(2)得φ2-φ1=2πc 1

,(12)因此φ2也被确定.如果φ2-φ1没有覆盖所有路径,可加大;2-;1然后重新计算φ2-φ1直到能达到要求为止.

最后的变换可总结如下

tan θ2=tan ;12tan ;22

-c 12tan ;2c 1,(13)<=c 1(φ-φ1),

(14)V (<,θ

)=v (φ,;),(15)c 1=π/(;2-;1),φ2-φ1=2

πc 1.

(13)与(14)式为球面局部区域到球台域的坐标变换,(15)式表示在坐标变换过程中速度

保持不变.当φ2-φ1=2π时,c 1=1.令;1=ε,;2=π-ε,ε→0,则(13)与(14)式退化为

球面自身的变换.

3　面波群速度反演

变换后的球台域并非整个球域,但可以解析延拓到整个球域.如此就可在延拓后的全球域作速度函数的球谐展开.由于(15)式保持速度不变,面波路径的积分可在变换前的原局部区域上进行.

群速度走时为

T i =∫L i d s v (;,φ)　.由于变换后的速度不变,可得T i =∫L i

d s V (θ(;),<(φ)),(16)其中L i 代表第i 条大圆弧路径.应该注意到,积分仍然在原局部区域的路径上进行,将速度函数在变换及延拓后的整个球面作球谐函数展开V (θ,<)=a 0+

∑N

n =1a n P n (cos θ)+∑N n =1∑n m =1(a m n cos m <+b m n sin m <)P m n (cos

θ).(17)

问题变为反演速度函数V 的球谐系数.由(16)式得

δT i =-∫L i δV V 2d s.(18)

将(17)式代入(18)式得到关于球谐系数增量的线性方程,具有如下形式

δT i =

∑n =1f in δc n ,(19)

其中c n 为第n 个球谐系数,f in 为(18)式中第n 个球谐系数增量项前的函数在第i 条路径上的积分.应用线性反演的方法可求出(19)式中球谐系数增量的解.

7

61期 朱良保等:保角变换在区域面波群速度反演中的应用

86地　球　物　理　学　报 44卷4　分辨与误差

以下的讨论中采用下脚标求和规则.

4.1　分辨核

设群速度的球谐展开具有如下形式

U g(r)=f i(r)c i,(20) c i为(17)式中第i个球谐系数,f i(r)为(17)式中第i个归一化的球谐函数,即

δ

=∫f i(r)f j(r)d r.

ij

设反演得到的球谐系数增量为δc′i,真实的球谐系数增量为δc i,反演的分辨矩阵为R ij,则δc′i=R ijδc j,

由上式得f i(r)δc′i=f i(r)R ijδc j,

即δU′g(r)=f i(r)R ijδjkδc k=f i(r)R ij∫f j(r′)f k(r′)d r′δc k,

=∫f i(r)R ij f j(r′)δU g(r′)d r′,

所以分辨核为f i(r)R ij f j(r′).

4.2　方差

由(20)式得ΔU g(r)=f i(r)Δc i,

ΔU

(r)2=Δc iΔc j f i(r)f j(r).

g

群速度反演的方差为〈(ΔU g(r))2〉=f i(r)〈Δc iΔc j〉f j(r),

其中〈Δc iΔc j〉为线性反演中球谐系数增量的协方差矩阵.

5　数值对比

图2(a)是用本文方法计算的63.5s周期的实际资料反演的群速度结果.它清楚的反映出青藏高原块体的轮廓以及物质向东南方向流动的图像.图3(a)是同样分辨率时用Tarantola的非线性方法计算的结果.图2(b),3(b)为相应的分辨核.图4为Ditmar2 Yanovskaya方法计算的结果.

应该指出,用本文计算的结果其分辨率是非均匀的.区域中心的分辨率高于区域的边界,这是该方法本征具有的特点,不同点的分辨率是可以计算的.一般来说,实际资料往往是区域中心的路径密度大于边界的密度,本文介绍的方法正好能反映这一特征,从一定意义上说,这种特征是本方法的优点.另一方面,通过坐标的旋转,可以把研究区域的中心移动到赤道的附近,从而改善分辨的均匀性.通过比较,不同方法计算的结果的差别主要反映在路径较稀的区域,这些区域的约束较差,如印度板块.但总的趋势以及主要结果是一致的.然而,用Tarantola方法或Ditmar2yanovskaya方法计算图3(a)和图4在作者的个人电脑上需要大约9h,用本文的方法计算图7(a)只需大约5min.

96 1期 朱良保等:保角变换在区域面波群速度反演中的应用

图4　Ditmar 2Y anovskaya 方法计算的周期为63.5s 的Rayleigh 波群速度分布图

Fig.4　The map of Rayleigh wave group velocity at the period of 63.5seconds com puted

using Ditmar 2Y anovskaya ’s method

6　结 论

利用保角变换把球面上的局部区域扩展到球面上更大的球台域,在解析延拓后的球面上作速度函数的球谐展开,在相同分辨的条件下极大地降低了球谐展开的阶数,使全球域的球谐展开方法应用于球面局部区域的面波群速度反演中成为可能.数值计算表明,该方法反演速度快,等值线光滑,构造界限清晰.本方法虽然是在面波群速度反演的问题中提出的,但最后的坐标变换式(13),(14)式却与具体的物理问题无关,它代表的是球面局部区域到球台区域的保角变换.因此,变换(13)与(14)式同样适用于任何球面局部区域场的分析.

参　考　文　献

[1]　Tarantola A ,Valette B.G eneralized nonlinear inverse problems solved using the least square criterion.Rev.Geophys.

S pace Phys.,1982,20:219—232.

[2]　Tarantola A ,Nersessian A.Three 2dimensional inversion without blocks.Geophys.J.R.ast r.Soc.,1984,76:299—

306.

[3]　Ditmar P G ,Y anovskaya T B.A generalization of the Backus 2G ilbert method for estimation of lateral variations of sur 2

face wave velocity.Iz v.A kad.N auk SSS R ,Fiz.Zemli ,1987,6:30—60.

[4]　Y anovskaya T B ,Ditmar P G.Smoothness criteria in surface wave tomography.Geophys.J Int.,1990,102:63—72.

[5]　Backus G ,G ilbert F.The resolving power of gross earth data.Geophys.J.R.ast r.Soc.,1968,16:169—205.07地　球　物　理　学　报 44卷

[6]　Backus G ,G ilbert F.Uniqueness in the inversion of inaccurate gross earth data.Philos.Trans.R.Soc.,1970,226:

123—192.

[7]　朱良保.区域面波群速度反演的球谐函数法.地球物理学报,1997,40(4):503—511.

ZHU Liang 2Bao.Regional surface 2wave group velocity inversion by spherical harmonic function method.Chi nese J.Geophys.(in Chinese ),1997,40(4):503—511.

APPL ICATION OF CONFORMAL MAPPING IN REGIONAL

SURFACE WAVE GROUP VE LOCIT Y INVERSION

ZHU L IAN G 2B AO XU Q IN G CHEN X IAO 2F EI

(Depart ment of Geophysics ,Peking U niversity ,Beijing 100871,China )

[Abstract]　Through conformal mapping ,in which the group velocities keep unchanged ,a re 2gional spherical area can be extended to a much lar ger spherical segment area ,and then to the whole sphere by analytical continuation.Representing velocities by spherical harmonic expan 2sion on the sphere after transformation can greatly reduce the number of spherical harmonic coefficients needed.Thus the inversion of surface 2wave velocities turns into the linear inversion of spherical harmonic coefficients.The resolution kernel and the root square deviation of the inversion solution can be obtained by analyzing the spherical harmonic coefficients.The algo 2rithm has advantages in a few aspects :fast computing speed ,smooth contours ,and clear tec 2tonic boundaries.Theoretically ,this method is applicable not only to the surface 2wave group velocity inversion ,but also to various analyses of regional fields.

[K ey w ords]　Spherical harmonic function ,Conformal mapping ,Group velocity inversion ,

Resolution.171期 朱良保等:保角变换在区域面波群速度反演中的应用