各向異性TI介質(zhì)qP反射波走時層析成像

黃光南， ZHOU Bing， 鄧居智， 李紅星， 李澤林， 張華

1 東華理工大學放射性地質(zhì)與勘探技術國防重點學科實驗室， 南昌 330013 2 中國地質(zhì)大學(武漢)地球內(nèi)部多尺度成像湖北省重點實驗室， 武漢 430074 3 中國石油大學(北京)油氣資源與探測國家重點實驗室， 北京 102249 4 Petroleum Geosciences, the Petroleum Institute, P.O. Box 2533, Abu Dhabi, UAE

?

各向異性TI介質(zhì)qP反射波走時層析成像

黃光南1,2,3， ZHOU Bing4， 鄧居智1， 李紅星1， 李澤林1， 張華1

1 東華理工大學放射性地質(zhì)與勘探技術國防重點學科實驗室， 南昌 330013 2 中國地質(zhì)大學(武漢)地球內(nèi)部多尺度成像湖北省重點實驗室， 武漢 430074 3 中國石油大學(北京)油氣資源與探測國家重點實驗室， 北京 102249 4 Petroleum Geosciences, the Petroleum Institute, P.O. Box 2533, Abu Dhabi, UAE

地震走時層析成像是反演地層各向異性參數(shù)分布的有效方法，但是關于地震各向異性介質(zhì)走時層析成像的研究并不多，其技術遠遠沒有達到成熟的階段.在野外數(shù)據(jù)采集時，地表反射波觀測方式相對井間和垂直地震剖面觀測方式的成本更低，利用qP反射波走時反演各向異性參數(shù)具有更加廣泛的實用價值.本文實現(xiàn)的TI介質(zhì)地震走時層析成像方法結合了TI介質(zhì)反射波射線追蹤算法、走時擾動方程和非線性共軛梯度算法，它可以對任意強度的TI介質(zhì)模型進行反演，文中嘗試利用qP反射波走時重建TI介質(zhì)模型的參數(shù)圖像.利用qP反射波對層狀介質(zhì)模型和塊狀異常體模型進行走時反演，由于qP波相速度對彈性模量參數(shù)和Thomsen參數(shù)的偏微分不同，所以可以分別反演彈性模量參數(shù)和Thomsen參數(shù).數(shù)值模擬結果表明：利用qP反射波可以反演出TI介質(zhì)模型的彈性模量參數(shù)與Thomsen參數(shù)，不同模型的走時迭代反演達到了較好的收斂效果，與各向同性介質(zhì)走時反演結果相比較，各向異性介質(zhì)走時反演結果具有較好的識別能力.

各向異性介質(zhì)；各向異性參數(shù)；彈性模量參數(shù)；走時層析成像

1 引言

室內(nèi)研究與野外觀測表明大多數(shù)巖石和地層存在明顯的各向異性.引起地層各向異性的主要因素有：巖石的結構各向異性，地層方向應力引起的各向異性，巖性各向異性，巖石晶體定向排列引起的各向異性，以及巖石裂隙裂縫引起的各向異性(Crampin, 1984; Helbig, 1981; 李芳等, 2012).地震波在各向同性介質(zhì)與各向異性介質(zhì)中的傳播特征不同，在各向異性介質(zhì)有三種體波：qP波，qSV波和qSH波，這三種波各自以自身的傳播速度與相互垂直的極化方向傳播(Zhou and Greenhalgh, 2006).因此，如果地層具有各向異性，那么地震走時層析成像方法必須適用于這種介質(zhì)模型，否則，地震走時反演結果可能會產(chǎn)生錯誤的地質(zhì)解釋.

經(jīng)過40多年的發(fā)展歷程，地震各向同性介質(zhì)走時層析成像方法已經(jīng)比較成熟，很多學者利用各種射線追蹤算法、數(shù)學反演算法和走時擾動方程得到了不同的地震走時層析成像方法.這里列舉一些在勘探地球物理學方面有代表性的地震走時層析成像方法.Zelt等(2006)的有限差分初至波走時層析成像方法得到了廣泛的應用，利用這種初至波走時層析成像方法在地下水污染的圈定和調(diào)查方面取得了良好的應用效果.Tryggvason等(2009)結合初至波走時層析和靜校正方法，將其運用于核廢料處置的勘探選址研究，所得地震剖面的可靠性與分辨率均得到了較大提高.Zhou(2006)提出了地震多尺度可變形地層走時層析成像方法，數(shù)值測試結果表明該方法比基于晶格或者網(wǎng)格節(jié)點的走時層析成像方法具有較高的成像精度.劉玉柱和楊積忠(2014)利用先驗信息在初至波走時層析成像過程中加入正則化約束條件，反演結果使原來隱藏的特征變得更加清晰.

關于地震各向異性介質(zhì)走時層析成像的研究也在不斷升溫，周輝和何樵登(1995)利用層狀TI介質(zhì)的P波和SH波走時，結合雙曲線擬合法和遺傳算法反演得到了五個彈性模量參數(shù)值.張文生和何樵登(1999)利用跨孔qP波和qSV波進行了統(tǒng)計法走時層析成像，提高了計算結果的穩(wěn)定性和精度.盧明輝等(2005)利用正交各向異性介質(zhì)P波走時反演得到了單層介質(zhì)的Thomsen參數(shù).劉玉柱等(2014)利用多參數(shù)聯(lián)合走時層析成像方法對VTI介質(zhì)進行走時反演取得了較好的數(shù)值測試結果.國外對地震各向異性介質(zhì)走時層析成像的研究相對較多，Cerveny(2001)是較早從事各向異性走時層析成像的學者，他運用各向異性介質(zhì)“程函方程”的哈密爾頓形式，推導出了各向異性介質(zhì)走時擾動方程.后續(xù)學者利用該走時擾動方程和不同的觀測系統(tǒng)對各種各向異性介質(zhì)模型進行反演取得了不同的效果.Jech(1988)對該走時擾動方程進行了修改，并且實現(xiàn)了三維各向異性介質(zhì)走時層析成像方法.Chapman和Pratt(1992)、Pratt和Chapman(1992)提出了各向異性介質(zhì)走時反演的線性化方法，利用井間觀測方式進行了走時反演.Wang 和Tsvankin(2013)提出了網(wǎng)格化的各向異性介質(zhì)走時層析成像方法，該方法極大地提高了成像的精度.Wang(2014)對各向異性介質(zhì)彎曲射線法產(chǎn)生的高度非線性問題，提出了修改的牛頓迭代法用于穩(wěn)定非線性問題的解，并利用井間觀測方式反演了各向異性介質(zhì)模型參數(shù).縱觀大多數(shù)各向異性介質(zhì)走時層析成像方法，它們都是采用了弱各向異性介質(zhì)假設條件.實際上，巖石和地層的各向異性強度是不確定的，將它們假設為弱各向異性介質(zhì)具有不合理性.Zhou等(2008)結合各向異性介質(zhì)初至波射線追蹤算法、走時擾動方程和非線性共軛梯度，實現(xiàn)了地震各向異性初至波走時層析成像方法.這種地震走時層析成像方法的優(yōu)點是能夠適用于任意強度的各向異性介質(zhì)模型，它的各向異性介質(zhì)走時擾動方程避免了qSV波和qSH波的奇異值問題，非線性共軛梯度算法可以同時反演出幾個彈性模量參數(shù)(或者Thomsen參數(shù)).

國內(nèi)外關于地震各向異性介質(zhì)走時層析成像的文章不是很多，其技術遠遠沒有達到成熟的階段，因此仍然有很多的理論和數(shù)值模擬研究值得探索.本文以Zhou等(2008)的各向異性初至波走時層析成像方法為基礎，利用反射波射線追蹤算法實現(xiàn)了TI介質(zhì)qP反射波非線性走時層析成像.由于地表反射波觀測方式相比井間和垂直地震剖面觀測方式的成本更低，因此TI介質(zhì)qP反射波非線性走時層析成像方法具有更加廣泛的實用價值.在數(shù)值模擬部分，利用qP反射波走時重建層狀TI介質(zhì)模型和塊狀異常體TI介質(zhì)模型的圖像，用于測試各向異性介質(zhì)參數(shù)重建和成像的可能性.通過分析各種模型的彈性模量參數(shù)與Thomsen參數(shù)反演結果，驗證TI介質(zhì)qP反射波走時層析成像方法的正確性.

2 方法原理

2.1 各向異性介質(zhì)走時擾動方程

Cerveny推導出了各向異性介質(zhì)走時擾動方程的表達式.很多學者對該公式做了不同的修改，Zhou推導出了各向異性介質(zhì)走時擾動方程的三種表達形式，這三種表達式被稱為一階走時擾動方程的波速偏微分形式，其中最為簡單的是走時擾動方程的相速度偏微分形式(Zhou et al., 2008)

(1)

(2)

各向異性模型參數(shù)化后，可以得到相速度vp(m,n)與群速度vg(m,n)模型.

(3)

2.2 雅可比矩陣表達式

當震源與檢波器的數(shù)量比較大時，它們之間存在大量的射線路徑，公式(3)可以寫成非線性反演的矩陣形式(Zhou et al., 2008)

(l=1,2,…)

(4)

式中上標-g表示廣義逆矩陣，Wd和Wm分別代表數(shù)據(jù)與模型參數(shù)的權重矩陣，λ是平衡數(shù)據(jù)擬合度與模型粗糙度之間的變量.(4)式可以利用非線性共軛梯度算法進行求解(Zhouetal., 1992)，其他變量的表達式如下:

(5)

(6)

3 數(shù)值模擬

下面利用各向異性介質(zhì)qP反射波走時層析成像進行數(shù)值模擬.數(shù)值模擬包括兩部分：第一部分是利用qP反射波走時反演TI介質(zhì)模型的彈性模量參數(shù)；第二部分是利用qP反射波走時反演TI介質(zhì)模型的Thomsen參數(shù).測試模型包括層狀介質(zhì)模型和塊狀異常體模型.Zhou對各向異性介質(zhì)初至波走時層析成像做過相應的研究工作，并且指出當?shù)貙訉ΨQ軸傾角不固定時，算法很難反演得到合理的傾角參數(shù)分布.原因是相速度和正演走時均對地層對稱軸傾角非常敏感，因此本文假設模型的地層對稱軸傾角為已知，只反演模型的彈性模量參數(shù)或者Thomsen參數(shù)分布.數(shù)值模擬實驗的觀測走時是射線追蹤得到的正演走時，并且添加了一定百分比的高斯噪聲.在非線性走時反演時，Wd和Wm是根據(jù)數(shù)據(jù)和模型先驗信息確定的權重矩陣，用以減小反演的多解性.文中假設先驗信息的權重相同，即數(shù)據(jù)的權重矩陣和模型的權重矩陣均為單位矩陣(Wd=Id和Wm=Im).阻尼因子λ是通過多次數(shù)值測試，比較反演的收斂情況而選取的數(shù)值大小.在反演迭代過程中，利用了阿爾法中值濾波器剔除每次迭代結果的高波數(shù)噪聲.地震走時層析成像的反演算法是非線性共軛梯度，Zhou等(1992)將其與不同反演算法做過數(shù)值和實際數(shù)據(jù)試驗，證明了它是一種有效而快速的反演算法.

3.1 qP反射波走時反演彈性模量參數(shù)模型

下面分別對層狀介質(zhì)模型和塊狀異常體模型進行TI介質(zhì)反射波非線性走時反演，求取模型的彈性模量參數(shù).

(1)層狀介質(zhì)模型的彈性模量參數(shù)反演.層狀介質(zhì)模型的長度為800 m，深度為600 m，第一層介質(zhì)的彈性模量參數(shù)為c11=8.81 GPa,c13=2.78 GPa,c33=7.34 GPa,c44=2.28 GPa,c66=3.88 GPa，對稱軸傾角為θ0=0°.第二層介質(zhì)的彈性模量參數(shù)為c11=21.49 GPa,c13=10.85 GPa,c33=19.54 GPa,c44=8.71 GPa,c66=8.45 GPa，對稱軸傾角為θ0=90°.第三層介質(zhì)的彈性模量參數(shù)為c11=13.86 GPa,c13=6.23 GPa,c33=12.36 GPa,c44=4.35 GPa,c66=5.09 GPa，對稱軸傾角為θ0=0°.采用地表激發(fā)和接收的觀測方式，炮點與檢波點位置相同，道間距為25 m，炮點與檢波點均為33個，模型的網(wǎng)格間距為50 m.將第一層介質(zhì)參數(shù)作為初始模型用以反演第二層和第三層介質(zhì)的彈性模量參數(shù).為了改善反演成像的效果，我們對模型設置了三個反射界面，充分利用三個反射界面的反射波走時信息.圖1a為層狀彈性模量參數(shù)模型的反射波射線路徑分布圖，射線覆蓋了整個模型范圍，圖1b是利用相同的各向異性介質(zhì)正演走時結合縱波速度走時層析成像所得速度反演結果，可以看出速度反演結果效果不理想，第一層速度值變化較大，第二層下界面位置不正確.圖2為層狀介質(zhì)模型的彈性模量參數(shù)(a—c—e—g)與反演結果(b—d—f—h)的對比圖，左側為真實彈性模量參數(shù)，右側為反演結果，從整體上分析，c11,c13,c33和c44均得到了較好的反演效果.由迭代反演收斂曲線(圖9)可知反演得到了較好的收斂水平，表明反演結果與真實模型相接近.

(2)塊狀異常體模型的彈性模量參數(shù)反演.塊狀異常體模型的長度為800 m，深度為500 m，背景介質(zhì)的彈性模量參數(shù)為c11=12.0 GPa,c13=5.0 GPa,c33=10.0 GPa,c44=4.0 GPa,c66=6.0 GPa，對稱軸傾角為θ0=0°.塊狀異常體的彈性模量參數(shù)為c11=20.0 GPa,c13=10.0 GPa,c33=15.0 GPa,c44=8.0 GPa,c66=12.0 GPa，對稱軸傾角為θ0=45°.采用地表激發(fā)和接收的觀測方式，炮點與檢波點位置重合，水平道間隔為16 m，炮點與檢波點均為51個，模型的網(wǎng)格間距為50 m.將背景介質(zhì)參數(shù)作為初始模型用以反演塊狀異常體的彈性模量參數(shù).為了改善反演成像的效果，我們對模型設置了兩個反射界面，一個是塊狀異常體的上邊界，另一個是模型的底界面.圖3a為塊狀彈性模量參數(shù)模型的反射波射線路徑分布圖，射線覆蓋了塊狀異常體區(qū)域，圖3b是利用相同的各向異性介質(zhì)正演走時結合縱波速度走時層析成像所得速度反演結果，反演所得異常體的形狀能夠得到識別，但是異常體的數(shù)值變化非常大，整體效果不太理想.圖4為塊狀異常體模型的彈性模量參數(shù)(a—c—e—g)與反演結果(b—d—f—h)的對比圖，左側為真實彈性模量參數(shù)，右側為反演結果，c13,c33,c44均得到了較好的反演效果，由于c11的敏感性較弱，它的反演效果相對較差.如果某個參數(shù)的速度偏導數(shù)相對于其他參數(shù)的速度偏導數(shù)不靈敏，那么利用觀測走時數(shù)據(jù)來重建該參數(shù)剖面是不太現(xiàn)實的.國內(nèi)外文獻均是針對敏感性較強的各向異性參數(shù)進行反演，而很難把靈敏度弱的各向異性參數(shù)反演出來.由迭代反演收斂曲線(圖9)可知反演的走時數(shù)據(jù)擬合較好，說明反演結果與真實模型相接近.

3.2 qP反射波走時反演Thomsen參數(shù)模型

下面分別對層狀介質(zhì)模型和塊狀異常體模型進行TI介質(zhì)反射波非線性走時反演，求取模型的Thomsen參數(shù).

(1)層狀介質(zhì)模型的Thomsen參數(shù)反演.層狀介質(zhì)模型的長度為800 m，深度為600 m，第一層介質(zhì)的Thomsen參數(shù)為α0=2.71 km·s-1,β0=1.51 km·s-1,ε=0.10,δ*=-0.13,γ=0.35，對稱軸傾角為θ0=0°；第二層介質(zhì)的Thomsen參數(shù)為α0=4.72 km·s-1,β0=2.89 km·s-1,ε=0.26,δ*=0.17,γ=0.17，對稱軸傾角為θ0=90°；第三層介質(zhì)的Thomsen參數(shù)為α0=3.37 km·s-1,β0=1.83 km·s-1,ε=0.18,δ*=0.02,γ=0.26，對稱軸傾角為θ0=0°.采用地表激發(fā)和接收的觀測方式，炮點與檢波點位置相同，水平道間距為25 m，炮點與檢波點均為33個，模型的網(wǎng)格間距為50 m.將第一層參數(shù)作為初始模型用以反演第二層和第三層的Thomsen參數(shù).為了改善反演成像的效果，我們對模型設置了三個反射界面，充分利用了三個反射界面的反射波走時信息.圖5a為層狀Thomsen參數(shù)模型的反射波射線路徑分布圖，射線覆蓋了整個模型范圍，圖5b是利用相同的各向異性介質(zhì)正演走時結合縱波速度走時層析成像所得速度反演結果，反演所得結果的中間層速度值變化比較大，并且中間層上下兩個邊界位置不易識別.圖6為層狀介質(zhì)模型的Thomsen參數(shù)(a—c—e—g)與反演結果(b—d—f—h)的對比圖，左側為真實Thomsen參數(shù)，右側為反演結果，α0,ε和δ*均得到了較好的反演效果，因為β0的敏感性較弱，它的反演效果稍差，其界面位置很難分辨清楚.由迭代反演收斂曲線(圖9)可知觀測走時與正演走時之間的殘差較小，反演得到的結果與真實模型相接近.

圖1 層狀介質(zhì)彈性參數(shù)模型的射線分布與速度反演結果Fig.1 The raypath distribution and inverted velocity profile of the elastic moduli layered model

圖2 層狀介質(zhì)模型的彈性模量參數(shù)與反演結果Fig.2 The true elastic moduli parameters and inverted results for the elastic moduli layered model

圖3 塊狀異常體彈性參數(shù)模型的射線分布與速度反演結果Fig.3 The raypath distribution and inverted velocity profile of the elastic moduli blocky abnormal body model

圖4 塊狀異常體模型的彈性模量參數(shù)與反演結果Fig.4 The true elastic moduli parameters and inverted results for the elastic moduli blocky abnormal body model

圖5 層狀介質(zhì)Thomsen參數(shù)模型的射線分布與速度反演結果Fig.5 The raypath distribution and inverted velocity profile of the Thomsen layered model

圖6 層狀介質(zhì)模型的Thomsen參數(shù)與反演結果Fig.6 The true elastic moduli parameters and inverted results for the Thomsen layered model

(2)塊狀異常體模型的Thomsen參數(shù)反演.塊狀異常體模型的長度為800m，深度為500m，背景介質(zhì)的Thomsen參數(shù)為α0=3.16 km·s-1,β0=2.0 km·s-1,ε=0.25,δ*=0.3,γ=0.25，對稱軸傾角為θ0=0°；塊狀異常體的Thomsen參數(shù)為α0=3.87 km·s-1,β0=2.83 km·s-1,ε=0.18,δ*=1.14,γ=0.25，對稱軸傾角為θ0=45°.采用地表激發(fā)和接收的觀測方式，炮點與檢波點位置重合，水平道間隔為16 m，炮點與檢波點均為51個，模型的網(wǎng)格間距為50 m.將背景介質(zhì)作為初始模型用以反演塊狀異常體的Thomsen參數(shù).為了改善反演成像的效果，我們對模型設置了兩個反射界面，一個是塊狀異常體的上邊界，另一個是模型的底界面.圖7a為塊狀異常體Thomsen參數(shù)模型的反射波射線路徑分布圖，射線覆蓋了塊狀異常體區(qū)域，圖7b是利用相同的各向異性介質(zhì)正演走時結合縱波速度走時層析成像所得速度反演結果，可以看出速度反演結果效果非常不理想，完全無法識別異常體的形態(tài)和位置.圖8為塊狀異常體模型的Thomsen參數(shù)(a—c—e)與反演結果(b—d—f)的對比圖，左側為真實Thomsen參數(shù)，右側為反演結果，參數(shù)α0,β0和δ*均得到了較好的反演效果，由于塊狀異常體的ε值比背景介質(zhì)的ε值更低，反演算法比較難反演出ε參數(shù).由圖9的收斂曲線可知該模型的走時反演得到了較好的收斂效果，表明反演模型與真實模型的相似性較好.

圖7 塊狀異常體Thomsen參數(shù)模型的射線分布與速度反演結果Fig.7 The raypath distribution and inverted velocity profile of the Thomsen blocky abnormal body model

圖8 塊狀異常體模型的Thomsen參數(shù)與反演結果Fig.8 The true elastic moduli parameters and inverted results for the Thomsen blocky abnormal body model

圖9 四種模型的迭代反演收斂曲線Fig.9 The iterative inversion convergence curve for these four models

4 結論

各向異性現(xiàn)象在地球介質(zhì)中普遍存在，地震走時層析成像需要考慮地層的各向異性.本文的各向異性TI介質(zhì)qP反射波走時層析成像方法適用于任意強度的各向異性介質(zhì)，該方法相比井間和垂直地震剖面走時層析成像方法具有更加廣泛的應用價值.由于qP波相速度對彈性模量參數(shù)和Thomsen參數(shù)的速度偏微分不同，所以利用qP反射波可以分別反演出彈性模量參數(shù)和Thomsen參數(shù).在數(shù)值模擬部分，利用qP反射波分別反演層狀介質(zhì)模型和塊狀異常體模型，均得到了較好的反演結果.這些模型的迭代反演收斂曲線均得到了較好的收斂效果，表明反演結果與真實模型之間相接近.與各向同性介質(zhì)走時反演結果相比較，各向異性介質(zhì)走時反演結果具有較好的識別能力.數(shù)值模擬結果證明了各向異性TI介質(zhì)qP反射波走時層析成像算法的正確性.由于受地表觀測系統(tǒng)的限制，反射波走時層析成像相比井間直達波走時層析成像的反演效果較差，在固定檢波器數(shù)量的情況下，應最大限度地增加qP反射波射線覆蓋范圍以減小走時反演的多解性和非線性.根據(jù)反演理論，已知的信息量越多反演的不確定性越弱，如果能夠綜合利用qP波、qSV波和qSH波進行多波走時聯(lián)合反演，將有助于改善各向異性參數(shù)的反演效果，因此各向異性TI介質(zhì)多波走時聯(lián)合反演值得更深入的研究.

Cerveny V. 2001. Seismic Ray Theory. Cambridge: Cambridge University Press.

Chapman C H, Pratt R G. 1992. Traveltime tomography in anisotropic media—I Theory.GeophysicalJournalInternational, 109(1): 1-19.Crampin S. 1984. Effective anisotropic elastic constants for wave-propagation through cracked solids.GeophysicalJournaloftheRoyalAstronomicalSociety, 76(1): 135-145.

Helbig K. 1981. Systematic classification of layered-induced transverse isotropy.GeophysicalProspecting, 29(4): 550-577.Jech J. 1988. Three-dimensional inversion problem for inhomogeneous transversely isotropic media.StudiaGeophysicaetGeodetica, 32(2): 136-143.Kreyszig E. 1993. Advanced Engineering Mathematics. New York: John Wiley & Sons Inc, 790.

Li F, Cao S Y, Yao J. 2012. Calculation of phase and group velocities in an arbitrary anisotropic medium.ChineseJournalofGeophysics(in Chinese), 55(10): 3420-3426.

Liu Y Z, Wang G Y, Dong L G, et al. 2014. Joint inversion of VTI parameters using nonlinear traveltime tomography.ChineseJournalofGeophysics(in Chinese), 57(10): 3402-3410.

Liu Y Z, Yang J Z. 2014. Offset-weighted seismic tomography aimed at using the first arrival efficiently.GeophysicalProspectingforPetroleum(in Chinese), 53(1): 99-115.Lu M H, Tang J H, Yang H Z, et al. 2005. P-wave traveltime analysis and Thomsen parameters inversion in orthorhombic media.ChineseJournalofGeophysics(in Chinese), 48(5): 1167-1171.Pratt R G, Chapman C H. 1992. Traveltime tomography in anisotropic media—II Application.GeophysicalJournalInternational, 109(1): 20-37.Tryggvason A, Schmelzbach C, Juhlin C. 2009. Traveltime tomographic inversion with simultaneous static corrections—well worth the effort.Geophysics, 74(6): WCB25-WCB33.

Wang X X, Tsvankin I. 2013. Ray-based gridded tomography for tilted transversely isotropic media.Geophysics, 78(1): C11-C23.

Wang Y H. 2014. Seismic ray tracing in anisotropic media: A modified Newton algorithm for solving highly nonlinear systems.Geophysics, 79(1): T1-T7.Zelt C A, Azaria A, Levander A. 2006. 3D seismic refraction traveltime tomography at a groundwater contamination site.Geophysics, 71(5): H67-H78.

Zhang W S, He Q D. 1999. Statistical tomography in anisotropic media by using crosshole traveltimes of qP wave and qSV wave.ChineseJournalofGeophysics(in Chinese), 42(S1): 195-203.

Zhou B, Greenhalgh S A, Sinadinovski C. 1992. Iterative algorithms for the damped minimum norm, least squares and constraining problem in seismic tomography.ExplorationGeophysics, 23: 497-505.

Zhou B, Greenhalgh S A. 2005. Analytic expressions for the velocity sensitivity to the elastic moduli for the most general anisotropic media.GeophysicalProspecting, 53(4): 619-641.Zhou B, Greenhalgh S A. 2006. Raypath and traveltime computations for 2D transversely isotropic media with dipping symmetry axes.ExplorationGeophysics, 37(2): 150-159.Zhou B, Greenhalgh S A, Green A. 2008. Nonlinear traveltime inversion scheme for crosshole seismic tomography in tilted transversely isotropic media.Geophysics, 73(4): D17-D33.Zhou H, He Q D. 1995. Seismic traveltime inversion for transversely

isotropic media.GeophysicalProspectingforPetroleum(in Chinese), 34(3): 63-68.Zhou H W. 2006. Multiscale deformable-layer tomography.Geophysics, 71(3): R11-R19.

附中文參考文獻

李芳, 曹思遠, 姚健. 2012. 任意各向異性介質(zhì)相(群)速度的計算. 地球物理學報, 55(10): 3420-3426, doi: 10.6038/j.issn.0001-5733.2012.10.025.

劉玉柱, 王光銀, 董良國等. 2014. VTI介質(zhì)多參數(shù)聯(lián)合走時層析成像方法. 地球物理學報, 57(10): 3402-3410.

劉玉柱, 楊積忠. 2014. 有效利用初至信息的偏移距加權地震層析成像方法. 石油物探, 53(1): 99-115.

盧明輝, 唐建侯, 楊慧珠等. 2005. 正交各向異性介質(zhì)P波走時分析及Thomsen參數(shù)反演. 地球物理學報, 48(5): 1167-1171, doi: 10.3321/j.issn:0001-5733.2005.05.026.

張文生, 何樵登. 1999. 各向異性介質(zhì)中利用跨孔qP波和qSV波走時的統(tǒng)計法成像. 地球物理學報, 42(S1): 195-203.

周輝, 何樵登. 1995. 橫向各向同性介質(zhì)走時反演. 石油物探, 34(3): 63-68.

(本文編輯 何燕)

Traveltime tomography of qP reflection waves in anisotropic TI media

HUANG Guang-Nan1,2,3， ZHOU Bing4， DENG Ju-Zhi1， LI Hong-Xing1， LI Ze-Lin1， ZHANG Hua1

1FundamentalScienceonRadioactiveGeologyandExplorationTechnologyLaboratory,EastChinaInstituteofTechnology,Nanchang330013,China2HubeiSubsurfaceMulti-scaleImagingKeyLaboratory(SMIL),ChinaUniversityofGeosciences,Wuhan430074,China3StateKeyLaboratoryofPetroleumResourcesandProspecting,ChinaUniversityofPetroleum,Beijing102249,China4PetroleumGeosciences,thePetroleumInstitute,P.O.Box2533,AbuDhabi,UAE

Seismic traveltime tomography is an effective method to invert underground anisotropic parameters. Research about the anisotropic traveltime tomographic method remains insufficient. Most studies on this issue are based on weak anisotropy assumption which was proposed by Leon Thomsen according to previous research. But the degree of anisotropic property of strata and rocks is uncertain in real applications. Therefore, a kind of traveltime tomography by use of qP reflection waves is proposed to invert arbitrary anisotropic TI media.First arrival traveltime tomography is often used for near-surface and crosswell anisotropic parameters inversion. In fact, traveltime tomography of qP reflection waves has a much more extensive practical value, because the surface observation system is more economic and general than the systems of crossholes and vertical seismic profiles. Here, we present a new nonlinear traveltime inversion method for the surface observation system which combines several important features: (1) A robust reflected wave ray tracing method is used for arbitrary anisotropic TI media; (2) The first-order traveltime perturbation equation is not the eigenvector form of Cerveny′s linearized formula, which suffers from a singularity problem for the two quasi-shear waves; (3) An effective computation of the Jacobian matrix is employed for an arbitrary anisotropic TI media; and (4) It adopts a fast, local minimization search style of nonlinear inversion. These features make the traveltime tomography of qP reflection waves can be used to invert anisotropic TI media with arbitrary anisotropic degree.In numerical simulation section, this seismic reflection tomography method was used to invert for a layered model and a blocky abnormal body model, respectively. The elastic moduli parameter and Thomsen parameter of these anisotropic models were obtained, respectively. Firstly, we used the traveltimes of qP reflection waves to invert the elastic moduli parameters of the layered model and blocky abnormal body model. (1) According to the tomograms of the layered model, the velocity image is not reconstructed very well although the ray paths can cover the whole model perfectly. The velocity values of the first and third layers differ from the true velocity value greatly. And the bottom interface of the second layer is not recovered very well. Fortunately, the elastic moduli parameters are all reconstructed very well. We can see these images are all very close to sections of the true elastic moduli parameters. (2) After inversion of the blocky abnormal body model, the velocity image is not reconstructed very well. Because the velocity values of the blocky abnormal body differ from the true velocity values greatly. On the contrary, the elastic moduli parameters are very close to the profiles of true elastic moduli parameters. But thec11parameter is not as good as the other three parameters. Secondly, we used the traveltimes of qP reflection waves to invert the Thomsen parameters of the layered model and blocky abnormal body model. (1) According to the tomograms of the layered model, the velocity image is still not reconstructed very well, because the velocity values of the second layer differ from the true velocity values greatly. The Thomsen parameters are reconstructed very well, except the image ofβ0is slightly not so good. We can see the inverted and true images ofα0,εandδ*are very close to each other. (2) After inversion of the blocky abnormal body model, the velocity image is not reconstructed correctly, because its numerical value and shape differ from the true velocity model greatly. However, the profiles ofα0,β0andδ*are reconstructed successfully, except theεparameter. The images of these three inverted parameters are very close to that of the true parameters. Besides, these four numerical experiments all reach satisfactory convergence levels for different anisotropic models. Hence, the simulation results show that traveltime of qP reflection waves can invert anisotropic parameters correctly.The traveltimes of qP reflection waves have been used to invert anisotropic parameters for the layered model and blocky abnormal body model successfully. The images of the reconstructed results are very close to the images of the true anisotropic parameters. The tomographic method proposed here has superior capability in recovering theoretical models compared with the isotropic traveltime tomographic method. The more information is used, the less uncertainty of the inversion results according to the inverse theory. In other words, the inversion results can be improved greatly if the traveltimes of qP, qSV and qSH waves are used in the traveltime inversion. Therefore, joint traveltime inversion of qP, qSV and qSH waves should worth further research.

Anisotropic media; Anisotropic parameters; Elastic moduli parameters; Traveltime tomography

10.6038/cjg20150618.

東華理工大學博士科研啟動基金(DHBK2013212)，中國地質(zhì)大學(武漢)地球內(nèi)部多尺度成像湖北省重點實驗室基金(SMIL-2015-10)，國家科技重大專項(2011ZX05024-001-02)，國家自然科學基金(41104073, 41004048, 41364004, 41104074, 41304097)和國家科技支撐計劃(2011BAB04B03)聯(lián)合資助.

黃光南，男，1983年生，講師，主要從事地震走時層析成像方法研究.E-mail:bobking2@126.com

10.6038/cjg20150618

P631

2014-06-25，2015-03-16收修定稿

黃光南，ZHOU Bing，鄧居智等.2015. 各向異性TI介質(zhì)qP反射波走時層析成像.地球物理學報，58(6):2035-2045,

Huang G N, Zhou B, Deng J Z，et al. 2015. Traveltime tomography of qP reflection waves in anisotropic TI media.ChineseJ.Geophys. (in Chinese),58(6):2035-2045,doi:10.6038/cjg20150618.