速度层析成像正反演技术

/ )
l] -
ldld
。
理论上可以证明: 只要投影采用间隔无限小, 并
Hale Waihona Puke Baidu能得到所有扫描积分测量的投影值时, 就能精确地 得到物性分布 f ( r, )。
由波动理论可知, 各种波动在时间和空间上是
连续的, 不同的物质对某一频率的波的反射、折射或
吸收都是不一样的, 波穿过不同密度的物质时, 波的 传播速度和方向会发生变化, 能量被吸收或衰减的
2期
杨晓弘等:速度层析成 像正反演技术研究
∀ 219∀
图 3 模型 一及反演结果
3. 2 模型二 模型二为一个长 100 m, 深 100 m 的正方形地
质区域, 它的背景速度为 2 000 m / s, 在它的中间部 分有一个长 20 m # 20 m 的 % L&形高速异常体, 波速 为 3 000 m / s( 图 4a) , 在模型区域两边放置发射和 接收装置, 仍按上述方法剖分区域。
法进行模型的反演计算, 计算结 果表明, 反演结果与正演模型十分逼近, 计算精度和速度都能达到满意的效果。
关键词: 速度层析成像; 正反演; 地质模型 ; 慢度
中图分类号: P631. 4
文献标识码 : A
文章编号: 1000- 8918( 2009) 02- 0217- 03
层析成 像技术 ( com puter ized tom ography) 是一 种用数学方法把许多射线路径得到的信息组合成射 线在其中传播的介质图像的技术, 医学上利用层析 成像技术可清 晰地重 现身体 内部 的构造 [ 1- 2] 。目 前, 国内外关于层析成像方法的研究已经取得了丰 硕的成果, 地球物理中的层析成像技术则是利用对 象的各种物理性质和物性参数来重建地质体的内部 图像, 为其他的地球物理资料处理方法提供精确的 速度模型 [ 3] 。自 20世纪 70年代以来, 速度层析成 像技术在油气田勘探开发、地质灾害预报、无损检测 等领域中的应用日益广泛, 已从直射线层析成像发 展到弯曲射线层析成像, 目前存在的主要问题仍然 是计算效率和精度相对较低 [ 4- 6] 。根据 B ackus G il bert理论 [ 7 - 8 ] , 分辨率和精度之间存在折中关系, 在 增加射线数目的情况下, 正演和反演的关键就变成 了求解一个大型稀疏线性方程组的问题。笔者用变 带宽存储方法对正演模型计算, 通过反演, 准确得到 了模型图像, 兼顾了正反演速度和精度的要求。
面, 具有节省计算机内存和运算速度快等显著优点。
2 模型的建立和计算
2. 1 正演模型建立的原理和方法 速度层析成像首先确定走时 T 与慢度 D 及射
线长度矩阵 A 的关 系, 如 果介质 的波速 变化 不大 ( < 10% ), A 可以认为从波发射点到接收点是沿直
线传播, 此时 T = D ds是一个线性系统。
通过编写的程序可以计算出模型中每条射线从 发射点到接收点所需的时间。
由于在实际测量中传播媒质对波能量的吸收, 射线不可能按直线传播, 所以模拟时采用线性回归 的原则以便和实际情况符合。在模型成像工作前, 再进行反演工作, 程序迭代次数取为 1 000, 精度 E 取到 0. 000 01。
利用 SURFER8. 0成图 软件, 对所 得的慢度值 进行成 像 处理 后, 就 得到 模 型一 的 反演 结果 ( 图 3b) , 为方便计算将单位取为 m s。从上图能清楚地 看到速度异常体的位置和大小, 并判断出它为低速 异常体。
程度不一样, 因此只要测得携带地质体内部信息的
波的有关参数, 应用适当反演方法, 就能够重建粗物 体内部的图像。
设区域中的波速分布用 V = V ( x, y )来表示, 其
中 x 为沿两孔方向的水平坐标, y 为孔深方向的坐
标, 当各测线的发射源和接收点坐标一定时, 该测线 上的传播时间 T k 可由下式确定
[ 3] 任云生, 孙德有, 吴国学, 等. 黑河 市八 车力金 钼矿床 地质特 征
及找矿标志 [ J] . 矿物学报, 2007, 27 ( 增 ). [ 4 ] 吕军, 王建民, 王洪波, 等. 土壤 地球化 学测量 在三 道湾子 金矿
床的应用 [ J] . 物探与化探, 2005, 29 ( 6) . [ 5 ] 魏永强. 使用 Surfer软件绘制 地质图件 和处理地 质数据的 方法
D (i)
=
D ( i- 1)
-
D(
i-
1) w i wi wi
P i wi,
这里, w i ( i = 1, ∃, P ), P 为第 i 个方程的投影系数
向量。当投影进行到最后一个方程所表示的超平面
并得到 D ( Q ) 时, 为一 次 完整 迭 代。 依 次迭 代, 有 lim D kQ = D 。当 Q < P 时, 上述迭代过程收敛到解
参考文献:
[ 1] 裴正林. 井间地 震层析成 像的现状 与进展 [ J ]. 地球物理 学进 展, 2001, 16 ( 3) : 91.
[ 2 ] 周灿灿, 王昌学. 水平井测井解释技术综述 [ J] . 地球物 理学进 展, 2006, 21 ( 1) : 152.
[ 3] 毛立峰, 王绪 本, 高永才. 大 地电磁概 率成像的 效果评价 [ J] . 地球物理学报, 2005, 48 ( 2) : 429.
由于受射线条数和区域网格数量的限制, 所得 到的图像在成像区域的上下方有部分失真, 表现为 在异常体的上下部位有速度异常值区域存在, 这种 失真可以通过增加射线条数和细分剖分区域的办法
图 4 模型二及反演结果
来克服。速度层析成像技术在无损检测和研究各种 地下特征、确定地下资源的分布位置等方面有着广 泛的应用。
第 33卷第 2期 2009年 4月
物探与化探
G EO PHY SICAL & GEO CHEM ICA L EXPLORAT ION
V o.l 33, N o. 2 A pr. , 2009
速度层析成像正反演技术
杨晓弘, 何继善, 谢冬琪
(中南大学 信息物理工程学院, 湖南 长沙 410083)
摘 要: 利用层析成像的基本原 理对地质异常体进行速度层析成 像正、反 演研究, 计算结 果可以 确定出 地质模 型中
∀ 218∀
物探与化探
33卷
图 1 区域离散化示意
式中, aij为第 i 条射线通过第 j个网格化区域单元时 的距离, D ij为单元中的慢度值, P、Q 分别为区域中 的网格总数和射线总数。式 ( 2)写成矩阵形式为
T = A∀ D , 其中, T 为 Q 维投影数据向量, D 为 P 维向量, A 为
[ 4 ] 赵爱华, 丁志峰. 宽角反射地震波走时模拟的双重网 格法 [ J] . 地球物理学报, 2005, 48 ( 5) : 1141.
[ 5] 龙海丽, 郝锦绮. 自电位层 析成像 的理论 与实验 研究 [ J] . 地 球物理学报, 2005, 48 ( 6) : 1343.
[ 6 ] 蒋鸿亮, 陈湛文, 陈小宏. 高分辨率 AVO 反演技术研 究 [ J] . 地 球物理学进展, 2006, 21 ( 2) : 478.
速度异常体的大小、位置和慢度 值等物性参数。首先通过建 立正演 模型, 将含有 异常体 的模型 区域剖 分成规 则的
矩形网格, 在区域两边的钻孔中分别安置 发射和接受装置, 采用多点发射、多点接收的方 法获得已 知地质条件 下的
每条速度射线的走时值, 而后根据线性回 归的原则对所得的时间数据加入 20% 的噪声, 采用医学 中成熟的 ART 算
将正演时间数据加噪后进行反演, 可得到模型 二的反演结果 ( 图 4b) , 由于是在钻孔左右两侧发射 和接收, 所以模型上下部位明显出现低速带, 这是由 于射线在速度异常体上下部位密度稀少的缘故。
4 结论
利用速度层析成像技术对模型正演的时间数据 加噪后, 利用 ART 反演得到了地质区域模型中正方 形、% L &形的速度异常体成像结果, 图像结果很好地 反映出不同速度异常体的大小、形状、位置。与所建 模型十分逼近, 取得了理想的效果。
空间中的一个点。当式 ( 1) 存在唯一解时, 其解必 为这 Q 个超平面的焦点, 迭代过程是从初值 D ( 0) 开
始的, 将 D ( 0) 投影到式 ( 1)中第一个方程所表示的
超平面上, 得到 D ( 1) , 再将 D ( 1 ) 投影到第二个方程所
表示的超平面, 得到 D ( 2) , 此过程用数学公式表示为
k
D ( s) , 且有 m in |D ( 0) - D ( s) |。
ART 算法可以 避免在求解大型 稀疏方程组时
遇到的内存不够和耗时过长的问题。
3 模型计算
3. 1 模型一 模型一为一个长 100 m, 深 100 m 的矩形地质
区域, 它的背景波速为 3 000 m / s, 在它的中间部分 有 1个长 20 m # 20 m 的正方形速度异常体, 波速为 2 000 m / s( 图 3a), 将模型区域按上述方法剖分, 左 右钻孔中各放置 110个发射和接收探头。
i
。
2. 2 反演模型的计算
在得到时间数据之后就可以利用 ART 算法对
此模型进行 反演工作, ART ( algebraic reconstract ion
techno logy) 是代数重建方法的 总称, 其基本思想是
K acam a提出的 % 投影方法 &。在 P 维空间中, 式 ( 1)
的每个方程代表了一个超平面, 而图像 D 则是 P 维
P #Q 阶投影矩阵, 写成矩阵的形式为
a11 a12 ∃ a1P
x1
t1
a21 a22 ∃ a2P ∀ x2 = t2 。
! !!! !
!
aQ 1 aQ 2 ∃ aQP
xP
tQ
通过求解上式, 可得到每个网格单元内波速慢 度值。笔者采用代数重建法进行反演, 这种算法在
用于解大型稀疏线性方程 组及生成线性方 程组方
[ 7] 朱留方, 沈建国. 从阵列声 波测井波 形处 理地层 纵、横波 时差 的新方法 [ J] . 地球物理学进展, 2006, 21( 2 ): 483.
[ 8] 赫尔曼 G T. 由投影 重建图像 ∋ ∋ ∋ CT 的理论基 础 [ M ] . 北京: 科学出版社, 1985.
[ 9] 王振东. 浅层地震勘探应用技术 [ M ] . 北京: 地质出版社, 1988.
下转 223页
2期
蔡朝阳等: M apG IS与 Surfer相结合在黑河市化探工作中的应用
∀ 223∀
参考文献:
[ 1] 丛丽娟, 廖蕾. 几种物化探制图软件与 M apG IS之间的数 据转换 [ J] . 内蒙古地质, 2002, ( 1 ).
[ 2] 何明华. M apG IS 数字高 程模拟 系统在 化探 数据处 理中 的应 用 [ J] . 矿山测量, 2005, ( 2 ).
Tk
=
Rk
ds V ( x, y )
=
D ( x,
Rt
y ) ds
( 1)
其中, Rk 为 k 条测线波的传播路径; D 为速度的倒
数, 即慢度值。区域离散化如图 1 所示。式 ( 1) 还
可表示为
QP
! ! T k =
aijD ( i, j)
( 2)
i= 1 j= 1
收稿日期: 2007- 12 - 03
1 速度层析成像的基本原理和方法
层析技术的图像重建算法数学 基础首先是由 R adon于 1917年完成的 [ 9] , 在二维空间中用极坐标 ( r, )表示某被测物理量 f ( r, )沿直线 L 的积分为
+
p ( l, ) = f l2 + s2, + arctan( s /l ) ds ,
-
记为 [ Rf ] ( l, ), 它是 f ( r, ) 在 方向上的投影, 称 为函数 f ( r, )的 Radon变换, 其逆变换为
[ R- 1p ] ( r, 即
+
)=
1
2
2 0-
[ p( rcos(
l, -
)
/ )
l] -
ldld
,
+
f ( r,
)
=
1
2
2 0-
[ p ( l, ) r cos( -
建立一个射线跨孔 (井 )层析勘 探的观测系统 (图 2), 在包含异常体在内的区域两边钻孔中分别 放置声波发射和接收装置, 采用多点发射、多点接收
图 2 区域间发射接收示意
的方式采集速度走时数据。每条射线的走时为
P
! T i = aijD ( i, j ) , j= 1
求和得到总的时间
Q
T=
!T
i= 1